how to find the right frequency axis when we take the Fourier transform of a function?

81 visualizzazioni (ultimi 30 giorni)
I want to know how to find the correct frequency axis when I am taking the fourier transfor. To make it clear, we know that for sinc function, we have the below relationship between the lenght of the sinc function and its fourier (the rect function):
based on this fact, each time I tried different frequency axix and then finally I find the correct one. But, actually I still do not know how to find the correct axis. This was a simple function that I already knew about the lenght of the function and its transforemed but for other functions I do not have any idea about the correct frequency axis.
Any help will be appreciated

Risposta accettata

Paul
Paul il 22 Giu 2022
Modificato: Paul il 23 Giu 2022
Hi Donya,
By default, the frequency variable for fourier in the Symbolic Math Toolbox is in rad/sec.
Also, Matlab definition of sinc in the SMT is sinc(x) = sin(pi*x)/(pi*x).
So with the example in the question we have:
syms t omega f T real
s(t,T) = sinc(2*t/T);
S_omega(omega,T) = simplify(fourier(s(t,T),t,omega));
Now, if we want the CTFT expressed in Hz
S_f(f,T) = simplify(S_omega(2*sym(pi)*f,T))
S_f(f, T) = 
Now make some plots for T = 2
figure
fplot(s(t,2),[-3 3]),grid
xlabel('t (sec)')
figure
fplot(S_f(f,2),[-3 3]),grid
xlabel('f (Hz)')
As expected and consistent with the plot in the Question, the mainlobe of s(t) has a span of T the CTFT spans -1/T to 1/T.
See fourier and sympref for more discusson the how the SMT parameterizes the CTFT and its inverse, and how to control the parameterization.
  13 Commenti
2NOR_Kh
2NOR_Kh il 25 Giu 2022
However, in my work I totally ignore the phase, but it was a great point to have in mind in case I wanted to work on the phase one day.
Thank you again Paul, you helped me a lot.
Paul
Paul il 26 Giu 2022
Modificato: Paul il 26 Giu 2022
for normalizing the DFT coefficients, we should divide them to the length of the vector,
Instead of "normalizing," I thnk a better word is "scaling." Whether or not we "should" scale the DFT coefficients depends on a) what we want the DFT coefficients to represent, and b) how the DFT is defined for the particular implementation being used.
One motivation for dividing the DFT coefficients by the length of the input vector, N, can be illustrated as follows. Consider again the case of the sine wave with period P and we take the DFT of time domain samples from 5 periods.
P = 0.1;
Ts = P/10;
Fs = 1/Ts;
N = 50;
t = (0:(N-1))*Ts;
y = sin(2*pi/P*t);
figure
stem(t,y)
Compute the FFT and apply fftshift()
Y = fft(y);
Y = fftshift(Y);
Compute the associated frequency vector in Hz for the fftshifted samples of Y.
if mod(N,2) == 0 % N is even
k = ( (-N/2) : ((N-2)/2) )/N*Fs;
else % N is odd
k = ( (-(N-1)/2) : ((N-1)/2) )/N*Fs;
end
Now plot
figure
stem(k,abs(Y)),grid
xlim([-15 15])
xlabel('Frequency (Hz)')
We see that, using Matlab's fft() implementation, the DFT amplitudes at +-10 Hz are 25. Suppose now we want to take those DFT coefficients and reconstruct the original 5 periods of the sine wave. We know that the DFT coefficients multiply complex exponentials of the form exp(1j*2*pi*k*n*Ts). Note that n*Ts = t and, in this case, we only need to worry about k = -10 and k = 10. After the fftshift, Y(21) corresponds to k = -10 and Y(31) to k = 10
yrecon = Y(21)*exp(1j*2*pi*(-10)*t) + Y(31)*exp(1j*2*pi*10*t);
figure
stem(t,yrecon)
This plot is the same as that above, except the amplitude is too high by a factor of 50, which not coincidentally is the value of N
N
N = 50
So, if we want to reconstruct the sine wave, we can divide Y by N
Y = Y/N;
figure
stem(k,abs(Y))
yrecon = Y(21)*exp(1j*2*pi*(-10)*t) + Y(31)*exp(1j*2*pi*10*t);
figure
stem(t,yrecon)
After scaling by N, the DFT coeffcients have amplitide of 1/2 that of the sine wave.
An alternative approach to reconstructing the sine wave is to divide by N in the reconstruction equation instead of scaling the DFT coeffeicients
Y = Y*N; % get back the original DFT
yrecon = (1/N) * ( Y(21)*exp(1j*2*pi*(-10)*t) + Y(31)*exp(1j*2*pi*10*t) );
which we won't plot because we know what it will look like.
Now look the documenation for fft, particularly the equation for the ifft. We see that equation includes the 1/N factor in front of the summation of the complex exponentials, exactly as shown in the yrecon equation immediately above.
In summary, using Matlab's implementation of fft/ifft, you can scale the DFT coefficients output from fft by 1/N if you want them to be "consistent" in some sense with the amplitude of the exponentials that would be used to reconstruct the signal. I think some people call this "amplitude scaling." IMO, I don't think you should ever explicitly scale Y by N, i.e., do not do
Y = fft(y)/N
because then you run the risk of forgetting later in the analysis thay Y is scaled and something else won't work, e.g., ifft(Y) yields the wrong answer if you forget to multiply the result by N. Instead, just divide by N on the fly as needed, e.g.,
plot(k.abs(Y)/N)
if you decide you want to see the "amplitude scaled" coefficients.
The CTFT, DTFT, and DFT are all related to each other, so if you change the definition of the DFT by scaling it by 1/N, you may need to reconsider how to maintain its relationshihp to the CTFT and DTFT as well if want all three to be consistent with each other.
Finally, before you do any scaling you really need to know how the DFT is defined for the particular software you're using. As we've seen, Matlab scales the DFT sum by 1 and the iDFT sum by 1/N (I think this convention is very common). However, the only requirement is that the product be 1/N, so other s/w might have the 1/N scale factor on the DFT sum and unity scale factor on the iDFT sum (or 1/sqrt(N) on both!). In this case, obviously, you wouldn't want to divide Y by N.

Accedi per commentare.

Più risposte (0)

Prodotti


Release

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by