why fft (exp(-t)) don't fit my analytique expression?
4 次查看(过去 30 天)
显示 更早的评论
Hey guys, I am trying to find the fft of exp (-t).
Also, I know that's the fft of (exp (-t)) is --> 1/(1+i*2*pi *f)
I tried to compare the fft (exp (-t)) with 1/(1+i*2*pi *f) and I didn't get the same result, when I zoom in I have almost the same shape but not the same amplitude.
can someone help me ?
thank you in advance
my code is :
clear all
clc
dt=0.1;
t=[-5:dt:5];
z=length(t);
f= (-z/2 : z/2 -1)*1/dt;
h_t=exp(-t);
h_f1=fftshift(fft(h_t))/z;
h_f2=1./(1+1i*2*pi*f);
subplot(1,2,1)
plot(t,h_t,'r')
title('h_t = e^-^t')
subplot(1,2,2)
plot(f,real(h_f1),'b',f,real(h_f2),'r')
grid on
title('h_f=TF(h_t)= 1/(1+i2\pif)')
legend('h_f_1=TF(h_t)', 'h_f_2=1/(1+2\pi i f)')
0 个评论
采纳的回答
Matt J
2022-2-3
编辑:Matt J
2022-2-3
dt=0.001;
t=-5:dt:5;
N=numel(t);
f= t/N/dt^2;
h_t=exp(-t).*(t>=0); %<---- truncate to zero for t<0
h_f1=fftshift(fft(ifftshift(h_t)))*dt; %<---- scale by dt, not 1/z
h_f2=1./(1+1i*2*pi*f);
subplot(1,2,1)
plot(t,h_t,'r')
title('h_t = e^-^t')
subplot(1,2,2)
plot(f,real(h_f1),'xb',f,real(h_f2),'r')
grid on
title('h_f=TF(h_t)= 1/(1+i2\pif)')
legend('h_f_1=TF(h_t)', 'h_f_2=1/(1+2\pi i f)')
xlim([-2,2])
2 个评论
Matt J
2022-2-4
编辑:Matt J
2022-2-4
Regarding f, when working with FFTs, you have to make sure dt,df, and N satisfy N*dt*df=1. So,
f=df*(t/dt)=t/(N*dt^2)
Regarding h_f1, you are using the FFT to approximate the continuous Fourier transform integral. In the formula for the Fourier transform, we can clearly see that the exponentials need to be weighted by both the signal samples x(t) and dt..
更多回答(1 个)
William Rose
2022-2-3
The difference is due to differences in normalization of the discrete Fourier transform (DFT) of x(t), where x(t) is discrete, and the continuous Fourier transform (cFT) of x(t), where x(t) is continuous.
The code below compares the absolute value of the DFTs of three discrete signals to the Fourier transform of the continuous signal. The three discrete signals are x=e(-at), sampled at different rates and durations.
The deviation between the DFT and cFT at high frequencies (where high means approaching the Nyquisy frequency) is due to the fact that the DFT is the convolution in frequency domain, or multiplication in the time domain, of a boxcar sequence with x(t). Another way of thinking of it is that the DFT must produce a signal that repeats over and over. In this case, that means restarting the decaying exponential after 5 or 10 secnds, depending on which version of x(t) we're looking at. That requires an almost discontinuous jump, which requires power at high frequencies to produce. So the DFT has power at high frequencies while the cFT does not.
The deviaiton between the DFT and cFT at f=0 is due to the fact that the DFT's value at DC, X(0), equals the sum of all the points in the signal. Put another way, it equals the mean value of the signal times N, the number of points. The cFT=1 at f=0. COnsider the three discrete signals X1, X2, X3. The mean value is approximately 1/a (=the area under a decaying exponential with rate constant a) divided by the signal duration in seconds. So, for the three signals, that is
mean =(1/a)/Tmax=(1/a)/(N*dt)=fs/(a*N)
X(f=0)=mean*N=[fs/(a*N)]*N=fs/a
The plots of abs(X(f)) show that this is in fact true. and it epxlains the discrepancy at DC between the DFT and the cFT.
%fourierTransformComparison.m
%Compare FFT of discrete signal to the theoretical Fourier transform.
%W. Rose 2022-02-03
clear;
N1=100;
fs1=10; %sampling rate (Hz)
t1=(0:N1-1)/fs1; %time vector (s)
f1=(0:N1-1)*fs1/N1; %frequency vector (Hz)
a=1.0; %rate constant (1/s)
x1=exp(-a*t1);
X1=fft(x1);
N2=50;
fs2=10; %sampling rate (Hz)
t2=(0:N2-1)/fs2; %time vector (s)
f2=(0:N2-1)*fs2/N2; %frequency vector (Hz)
a=1.0; %rate constant (1/s)
x2=exp(-a*t2);
X2=fft(x2);
N3=100;
fs3=20; %sampling rate (Hz)
t3=(0:N3-1)/fs3; %time vector (s)
f3=(0:N3-1)*fs3/N3; %frequency vector (Hz)
a=1.0; %rate constant (1/s)
x3=exp(-a*t3);
X3=fft(x3);
Xtheo=1./(a+1i*2*pi*f3);
figure;
subplot(121)
plot(t1,x1,'-r.',t2,x2,'-gx',t3,x3,'-bo');
grid on; xlabel('Time (s)'); ylabel('x(t)');
legend('N=100,fs=10','N=50,fs=10','N=100,fs=20');
subplot(122)
plot(f1,abs(X1),'-r.',f2,abs(X2),'-gx',f3,abs(X3),'-bo',f3,abs(Xtheo),'-k^');
grid on; xlabel('Frequency (Hz)'); ylabel('|X(f)|');
legend('N=100,fs=10','N=50,fs=10','N=100,fs=20','Theo');
Experiment with the code above.
2 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Spectral Measurements 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!