cftbyfft.m

来自「本代码是相应的快速傅立叶以及反傅立叶算法的相应实现」· M 代码 · 共 30 行

function [AW,f]=cftbyfft(wt,t,flag)
%cftbyfft.m
%本程序采用FFT计算连续时间Fourier变换。输出幅频谱数据对(f,AW)。
%	输入量(wt,t)为已经窗口化了的时间函数wt(t)，它们分别是长度为N的向量。
%	对于"非平凡"取值时段有限的情况，应使该时段与窗口长度相比足够小。以
%	提高频率分辨率。
%	对于"非平凡"取值时段无限的情况，窗口长度的选取应使窗口外的函数值小
%	到可忽略，以提高近似精度。
%	输入宗量flag控制输出CFT的频率范围。
%		flag取非0时（缺省使用），频率范围在[0,fs）；
%		flag取0时，频率范围在[-fs/2,fs/2)。
if nargin==2;flag=1;end
N=length(t);					%采样点数，应为2的幂次，以求快速。
T=t(length(t))-t(1);			%窗口长度
dt=T/N;						%时间分辨率。
W0=fft(wt);					%施行FFT变换								<16>
W=dt*W0;						%算得[0,fs)上的N点CFT值
df=1/T;						%频率分辨率
n=0:1:(N-1);
%把以上计算结果改写到[-fs/2,fs/2]范围
if flag==0
   n=-N/2:(N/2-1);
   W=fftshift(W);			%产生满足式（4.13.3.1-6）的频谱
end
f=n*df;					%频率分度向量
AW=abs(W);					%福频谱数据向量
if nargout==0
	plot(f,AW);grid,xlabel('频率f');ylabel('|w(f)|')
end