FFT and IFFT Problem (numerically and analytically result is not same)

12 views (last 30 days)
ashish gupta on 6 Jul 2016
Commented: Jan Orwat on 7 Jul 2016
Hello every one i stuck in one problem i have two signal first one is h(t)=exp(-a*t) and other one is x(t)=1-exp(-a*t). i convoluted these two function and estimated y(t)= h(t) convolution x(t) or simply multiplication in frequency domain. my numerical estimation of y(t) is not same as i got analytical expression of these two function convolution can any one help me what is problem in numerically estimation in matlab. is it matlab fft and ifft problem? my code is also attached with this
ts=1e-4;
point=2000;
t=ts*(0:point-1);
fs=1/ts;
f1=fs*(0:point-1)/point;
a=4e2;
h_t=a*exp(-a*t);
b=3e2;
x_t=1-exp(-b*t);
fft_x=fft(x_t)*ts;
fft_h=fft(h_t)*ts;
y_t=real(ifft(fft_h.*fft_x)/ts);
figure(1)
subplot(3,1,1)
plot(t,x_t)
subplot(3,1,2)
plot(t,h_t)
subplot(3,1,3)
plot(t,y_t)
and analytical expression of y(t) is 1+1/(-b+a)*(-exp(-b*t)*a+b*exp(-a*t));

Jan Orwat on 6 Jul 2016
This is caused by side effects of digital domain
1. You are calculating circular convolution instead of normal one. To make convolution via fft/ifft you have to use zero padding.
2. After convolution you will get vector longer than those first two vectors. Because your vector x(t)=1-exp(-a*t) should contain values near 1 when we are counting towards infinity and we used zero padding to be able to convolve via fft it is wise to abandon some of this data.
3. No need to multiply and then divide by ts.
ts=1e-4;
point=2000;
t=ts*(0:point-1);
fs=1/ts;
f1=fs*(0:point-1)/point;
a=4e2;
h_t=a*exp(-a*t);
b=3e2;
x_t=1-exp(-b*t);
%%fft convolution
n = 2*point-1;
y_t = real(ifft(fft(x_t, n).*fft(h_t, n)));
y_t = y_t(1:point)*ts;
figure(1)
subplot(3,1,1)
plot(t,x_t)
subplot(3,1,2)
plot(t,h_t)
subplot(3,1,3)
plot(t,y_t)
Jan Orwat on 7 Jul 2016
You are welkome.
Regarding your question. Well, it's yes and no. Your algebraic solution/logic is correct. Your numeric intuition is correct too. As you can see, change to your original code was minimal. Here the "no" part comes.
Note, those algebraic functions and their convolution via Fourier transform are defined in ℝ. But in this MATLAB simulation their representation is changed from continous to digital where we have set of quantised, discrete-time values. This digital representation has some properties we have to deal with. That's the "price" we have to pay for easy/cheap signal processing and storage.
Check out the differences between Fourier transform and its discrete version at Wikipedia. Look into the boxes on the right side. I'm sorry I'm not very familiar with proper english literature to recommend for a good start in this topic.