Different results when compearing DFT from fft to the "real" fourier transformation

1 次查看（过去 30 天）
Lillebror Sagmen-Andersson
编辑： Paul 2022-5-21
I'm trying to compare the DFT computed from matlabs fft, to the "real" fourier transformation for the signal that can derived from the figure above. However I can't seem to get the the signals to match.
I do not have access to syms for this.
It is known that T = 3, and a = 1. EDIT : T_0 = 3, || T =/= 3
I have tried doing it likes this:
x=@(t) exp(-a.*t).*((t>=0)&(t<=3)); % a = 1
T=T_0/N_0; %T_0 & N_0 may be choosen to achive T = 3.
omega=linspace(-pi/T,pi/T,4097);
X = (1-exp(-(a+1i*omega)*T))./(a+1i*omega);
t=(0:T:T*(N_0-1))';
xf=T*x(t);
xf(1)=T*(x(T_0)+1)/2;
X_r=fft(xf);
r=[-N_0/2:N_0/2-1]';
omega_r=r*2*pi/T_0;
then using
subplot(211);
plot(omega,abs(X),'r',omega_r,fftshift(abs(X_r)),'bo');
xlabel('\omega');ylabel('|X(\omega)|')
%and
legend('True FT',strcat('DFT with T_0 = ',num2str(T_0),' , N_0 = ',num2str(N_0)));
subplot(212);
plot(omega,angle(X),'r',omega_r,fftshift(angle(X_r)),'bo');
xlabel('\omega'); ylabel('\angle X(\omega)')
legend('True FT',strcat('DFT with T_0 = ', num2str(T_0),' , N_0 = ',num2str(N_0)));
However as you can see the plots do not agree, and I don't see whats is going wrong.

请先登录，再回答此问题。

采纳的回答

Lillebror Sagmen-Andersson
I made an error when implementing the manualy calculated X(W)
in the original code it was
X = (1-exp(-(a+1i*omega)*T))./(a+1i*omega);
While it should have been
X = (1-exp(-(a+1i*omega)*T_0))./(a+1i*omega); %T --> T_0
and T =/= 3, but rather T_0 = 3.
Many thanks to Paul for moving my attention to the mixup in T

更多回答（3 个）

Lillebror Sagmen-Andersson
Soultion found in using different T values for the FFT etc, and the "real" transformation. While changing omega for a longer? axis:
a=1;
T_0=768; N_0=256;
T=T_0/N_0;
omega=linspace(-20.*pi/T,20.*pi/T,4097); %Change
Numerator = 1-exp(-(a+1i*omega)*T);
Denumerator = a+1i*omega;
X = (Numerator./Denumerator); % fouriertransform <-> x(t)
T_0=3; N_0=500; T=T_0/N_0;%Change
x=@(t) exp(-a.*t).*((t>=0)&(t<=3));
t=(0:T:T*(N_0-1))';
xf=T*x(t);
xf(1)=T*(x(T_0)+1)/2;
X_r=fft(xf);
r=[-N_0/2:N_0/2-1]';
omega_r=r*2*pi/T_0;
%----------------------PLOT STUFF-----------------------%
figure
grid on
subplot(211);
plot(omega,abs(X),'k',omega_r,fftshift(abs(X_r)),'ko');
xlabel('\omega');ylabel('|X(\omega)|')
axis([-0.1 20 -0.5 2]);
legend('True FT',strcat('DFT with T_0 = ',num2str(T_0),' , N_0 = ',num2str(N_0)));
subplot(212);
plot(omega,angle(X),'k',omega_r,fftshift(angle(X_r)),'ko');
xlabel('\omega'); ylabel('\angle X(\omega)')
axis([-0.1 20 -2 2]);
legend('True FT',strcat('DFT with T_0 = ', num2str(T_0),' , N_0 = ',num2str(N_0)));
now the graphs agree as expected.
Paul
Paul 2022-5-20
编辑：Paul 2022-5-21
Ran in:
Hi Lillebror,
Referring only to the code in the original question, it looks like there is a mix up between the variables used to define the duration of the signal and the sample period to generate the samples of the same, but it's hard to say because the entire code is not posted in the question. I think you're looking for something like this.
Define the signal of interest
T = 3; % T defines the duration of the signal
a = 1;
xfunc = @(t) exp(-a.*t).*((t >= 0) & (t <= T)); % a = 1, T = 3
Pick a sampling period Ts. Though not really necessary, choose Ts such that T/Ts is a nice integer
N = 31;
t = linspace(0,T,N);
Ts = t(2); % Ts is the sampling period
Define the frequency vector for the CTFT of x(t) and then compute the CTFT
omega_c = linspace(-pi/Ts,pi/Ts,4097);
XCTFT = (1 - exp(-(a + 1i*omega_c)*T)) ./ (a + 1i*omega_c);
Samples of x(t)
x = xfunc(t);
Adjust the endpoint to account for the effect of impulse sampling at discontinuities. I think this is similar to your adjustment of xf(1).
x(1) = x(1)/2;
x(end) = x(end)/2;
Compute the DFT
XDFT = fft(x);
Do the fftshift and get the associated frequency vector for N odd
XDFT = fftshift(XDFT);
omega_n = (-(N-1)/2 : (N-1)/2)*2*pi/N/Ts;
Compare, note the scaling by Ts on the DFT
figure
subplot(211);hold on
plot(omega_c,abs(XCTFT));
stem(omega_n,abs(Ts*XDFT));
subplot(212);hold on
plot(omega_c,180/pi*angle(XCTFT));
stem(omega_n,180/pi*angle(T*XDFT));
xlabel('rad/sec')
Michal
Michal 2022-5-18
Can you clarify the purpose of this operation (why multiply by T?):
t=(0:T:T*(N_0-1))';
xf=T*x(t)
I'm not following your operations with the time axis, sorry

类别

在 Help Center 和 File Exchange 中查找有关 Spectral Measurements 的更多信息

标签

产品


版本

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by