2014-10-11 91 views
1

我想將exp(-t^2)的FFT與函數的解析傅里葉變換exp( - (w^2)/ 4)/ sqrt(2 ),在-3到3的頻率範圍內。在Matlab中比較FFT函數和分析FT解決方案

我已經寫了下面的matlab代碼,並且現在迭代了很多次,但沒有成功。

fs = 100; %sampling frequency 
dt = 1/fs; 
t = 0:dt:10-dt; %time vector 
L = length(t); %number of sample points 
%N = 2^nextpow2(L); %necessary? 

y = exp(-(t.^2)); 
Y=dt*ifftshift(abs(fft(y))); 

freq = (-L/2:L/2-1)*fs/L; %freq vector 

F = (exp(-(freq.^2)/4))/sqrt(2); %analytical solution 

%Y_valid_pts = Y(W>=-3 & W<=3); %compare for freq = -3 to 3 
%npts = length(Y_valid_pts); 
% w = linspace(-3,3,npts); 
% Fe = (exp(-(w.^2)/4))/sqrt(2); 

error = norm(Y - F) %L2 Norm for error 

hold on; 
plot(freq,Y,'r'); 
plot(freq,F,'b'); 
xlabel('Frequency, w'); 
legend('numerical','analytic'); 
hold off; 

你現在可以看到,我只是試圖讓這兩個地塊看起來相似。最終,我想找到一種方法來做兩件事:1)找到最小採樣率, 2)找到最小採樣數, 達到一個誤差(定義爲二者之差的L2範數解決方案)的10^-4。

我覺得這很簡單,但我似乎甚至看不到這兩個圖形在視覺上一致。 如果有人能讓我知道我要出錯的地方以及我如何處理上述兩點(最小採樣頻率和最小採樣數),我會非常感激。

由於

回答

2

要注意的第一件事是,傅立葉變換對的功能exp(-t^2)在+/-無窮大的範圍,如可以從tables of Fourier transforms導出實際上是:

Gaussian FT: exp(-t^2) <==> exp(-(w^2)/4)*sqrt(pi)

最後,在生成函數exp(-t^2)時,您將t的範圍限制爲正值(而不是採用整個+/-無窮大範圍)。 對於保持關係,你將因此必須產生exp(-t^2)的東西,例如:

t = 0:dt:10-dt;  %time vector 
t = t - 0.5*max(t); %center around t=0 
y = exp(-(t.^2)); 

然後,可變w表示其中通過相關的歸一化頻率freq弧度角頻率:

w = 2*pi*freq; 

因此,

F = (exp(-((2*pi*freq).^2)/4))*sqrt(pi); %analytical solution