Difference between fft and manually coding for a fourier transform

9 Ansichten (letzte 30 Tage)
Hello everyone!
I've been trying to introduce the fast fourier transform function fft into my code, to replace my manually coded fourier transform. The results I get are different and I have no idea why. The function I am transforming (called cx) is an approximation of a double peak delta function - the peaks are at different places each time but always symmetrical about the origin. e.g. one peak at x=0.1 and the other at x=-0.1. Here's my code:
% Variables N=512; x=linspace(-0.3,0.3,N); L=0.6; dx=L/N; k=linspace(-N/(2*L),N/(2*L),N);
% Manual Fourier Tranform (approximation of the integral FT) % cx is the function to be transformed - it is the same size as x
for g=1:N F(g) = sum(exp(-2*pi*i*k(g)*x).*cx)*dx; end
% FFT
FF=fft(cx)*dx; FFS=fftshift(FF);
The two delta functions are at different positions each time I run the code (it's modelling the positions of particles), and at some positions the fft and manual FT give exactly the same results, whereas at other positions they give different results. The difference is usually that one of the FTs has a peak in the centre but the other does not.
Thanks a lot!
Paul

Akzeptierte Antwort

Dr. Seis
Dr. Seis am 25 Mai 2012
You are not correctly defining the frequencies "k" that are associated with the amplitudes in the frequency domain, but it starts out with incorrectly defining "dx".
In general (assuming your time/space domain data are stored in "cx" and time/space samples stored in "x"):
N = length(cx);
dx = (x(end)-x(1))/(N-1); % Time increment
Nyq = 1/(2*dx); % Nyquist frequency
df = 1/(N*dx); % Frequency increment
if mod(N,2) == 0 % N is even
k = -Nyq : df : Nyq-df;
else % N is odd
k = [sort(-1*(df:df:Nyq)) (0:df:Nyq)];
end
  2 Kommentare
Paul
Paul am 28 Mai 2012
Thanks for your answer Elige - this works well!
Nihal Pai
Nihal Pai am 5 Mai 2017
Dr.Seis, could you please expand on this specific line? k = [sort(-1*(df:df:Nyq)) (0:df:Nyq)];
Why do you sort it in this way?

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Kategorien

Mehr zu Fourier Analysis and Filtering finden Sie in Help Center und File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by