fft原理详解.docx
《fft原理详解.docx》由会员分享,可在线阅读,更多相关《fft原理详解.docx(18页珍藏版)》请在冰点文库上搜索。
fft原理详解
FFT算法
FFT是离散傅立叶变换的快速算法,可以将一个信号变换到频域。
有些信号在时域上是很难看出什么特征的,但是如果变换到频域之后,就很容易看出特征了。
这就是很多信号分析采用FFT变换的原因。
另外,FFT可以将一个信号的频谱提取出来,这在频谱分析方面也是经常用的。
虽然很多人都知道FFT是什么,可以用来做什么,怎么去做,但是却不知道FFT之后的结果是什意思、如何决定要使用多少点来做FFT。
现在圈圈就根据实际经验来说说FFT结果的具体物理意义。
一个模拟信号,经过ADC采样之后,就变成了数字信号。
采样定理告诉我们,采样频率要大于信号频率的两倍,这些我就不在此罗嗦了。
采样得到的数字信号,就可以做FFT变换了。
N个采样点,经过FFT之后,就可以得到N个点的FFT结果。
为了方便进行FFT运算,通常N取2的整数次方。
假设采样频率为Fs,信号频率F,采样点数为N。
那么FFT之后结果就是一个为N点的复数。
每一个点就对应着一个频率点。
这个点的模值,就是该频率值下的幅度特性。
具体跟原始信号的幅度有什么关系呢假设原始信号的峰值为A,那么FFT的结果的每个点(除了第一个点直流分量之外)的模值就是A的N/2倍。
而第一个点就是直流分量,它的模值就是直流分量的N倍。
而每个点的相位呢,就是在该频率下的信号的相位。
第一个点表示直流分量(即0Hz),而最后一个点N的再下一个点(实际上这个点是不存在的,这里是假设的第N+1个点,也可以看做是将第一个点分做两半分,另一半移到最后)则表示采样频率Fs,这中间被N-1个点平均分成N等份,每个点的频率依次增加。
例如某点n所表示的频率为:
Fn=(n-1)*Fs/N。
由上面的公式可以看出,Fn所能分辨到频率为为Fs/N,如果采样频率Fs为1024Hz,采样点数为1024点,则可以分辨到1Hz。
1024Hz的采样率采样1024点,刚好是1秒,也就是说,采样1秒时间的信号并做FFT,则结果可以分析到1Hz,如果采样2秒时间的信号并做FFT,则结果可以分析到。
如果要提高频率分辨力,则必须增加采样点数,也即采样时间。
频率分辨率和采样时间是倒数关系。
假设FFT之后某点n用复数a+bi表示,那么这个复数的模就是An=根号a*a+b*b,相位就是Pn=atan2(b,a)。
根据以上的结果,就可以计算出n点(n≠1,且n<=N/2)对应的信号的表达式为:
An/(N/2)*cos(2*pi*Fn*t+Pn),即2*An/N*cos(2*pi*Fn*t+Pn)。
对于n=1点的信号,是直流分量,幅度即为A1/N。
由于FFT结果的对称性,通常我们只使用前半部分的结果,即小于采样频率一半的结果。
好了,说了半天,看着公式也晕,下面圈圈以一个实际的信号来做说明。
假设我们有一个信号,它含有2V的直流分量,频率为50Hz、相位为-30度、幅度为3V的交流信号,以及一个频率为75Hz、相位为90度、幅度为的交流信号。
用数学表达式就是如下:
S=2+3*cos(2*pi*50*t-pi*30/180)+*cos(2*pi*75*t+pi*90/180)
式中cos参数为弧度,所以-30度和90度要分别换算成弧度。
我们以256Hz的采样率对这个信号进行采样,总共采样256点。
按照我们上面的分析,Fn=(n-1)*Fs/N,我们可以知道,每两个点之间的间距就是1Hz,第n个点的频率就是n-1。
我们的信号有3个频率:
0Hz、50Hz、75Hz,应该分别在第1个点、第51个点、第76个点上出现峰值,其它各点应该接近0。
实际情况如何呢我们来看看FFT的结果的模值如图所示。
图1FFT结果
从图中我们可以看到,在第1点、第51点、和第76点附近有比较大的值。
我们分别将这三个点附近的数据拿上来细看:
1点:
512+0i
2点:
-
3点:
-
50点:
-
51点:
-192i
52点:
-
75点:
76点:
+192i
77点:
+
很明显,1点、51点、76点的值都比较大,它附近的点值都很小,可以认为是0,即在那些频率点上的信号幅度为0。
接着,我们来计算各点的幅度值。
分别计算这三个点的模值,
结果如下:
1点:
512
51点:
384
76点:
192
按照公式,可以计算出直流分量为:
512/N=512/256=2;50Hz信号的幅度为:
384/(N/2)=384/(256/2)=3;75Hz信号的幅度为192/(N/2)=192/(256/2)=。
可见,从频谱分析出来的幅度是正确的。
然后再来计算相位信息。
直流信号没有相位可言,不用管它。
先计算50Hz信号的相位,atan2(-192,=,结果是弧度,换算为角度就是180*/pi=。
再
计算75Hz信号的相位,atan2(192,=弧度,换算成角度就是180*pi=。
可见,相位也是对的。
根据FFT结果以及上面的分析计算,我们就可以写出信号的表达式了,它就是我们开始提供的信号。
总结:
假设采样频率为Fs,采样点数为N,做FFT之后,某一点n(n从1开始)表示的频率为:
Fn=(n-1)*Fs/N;该点的模值除以N/2就是对应该频率下的信号的幅度(对于直流信号是除以N);该点的相位即是对应该频率下的信号的相位。
相位的计算可用函数atan2(b,a)计算。
atan2(b,a)是求坐标为(a,b)点的角度值,范围从-pi到pi。
要精确到xHz,则需要采样长度为1/x秒的信号,并做FFT。
要提高频率分辨率,就需要增加采样点数,这在一些实际的应用中是不现实的,需要在较短的时间内完成分析。
解决这个问题的方法有频率细分法,比较简单的方法是
采样比较短时间的信号,然后在后面补充一定数量的0,使其长度达到需要的点数,再做FFT,这在一定程度上能够提高频率分辨力。
具体的频率细分法可参考相关文献。
[附录:
本测试数据使用的matlab程序]
实例一:
S=2+3cos(2pi*50t-pi/6)+(2pi*75t+pi/2)
closeall;%先关闭所有图片
Adc=2; %直流分量幅度
A1=3; %频率F1信号的幅度
A2=;%频率F2信号的幅度
F1=50; %信号1频率(Hz)
F2=75; %信号2频率(Hz)
Fs=256;%采样频率(Hz)
P1=-30;%信号1相位(度)
P2=90; %信号相位(度)
N=256; %采样点数
t=[0:
1/Fs:
N/Fs];%采样时刻
%信号
S=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);
%显示原始信号
plot(S);
title('原始信号');
figure;
Y=fft(S,N);%做FFT变换
Ayy=(abs(Y));%取模
plot(Ayy(1:
N));%显示原始的FFT模值结果
title('FFT模值');
figure;
Ayy=Ayy/(N/2); %换算成实际的幅度
Ayy
(1)=Ayy
(1)/2;
F=([1:
N]-1)*Fs/N;%换算成实际的频率值
plot(F(1:
N/2),Ayy(1:
N/2)); %显示换算后的FFT模值结果
title('幅度-频率曲线图');
figure;
Pyy=[1:
N/2];
fori=1:
N/2
Pyy(i)=phase(Y(i));%计算相位
Pyy(i)=Pyy(i)*180/pi;%换算为角度
end;
plot(F(1:
N/2),Pyy(1:
N/2)); %显示相位图
title('相位-频率曲线图');
实例一:
S=1+(2pi*20t)+(2pi*60t)
closeall;%先关闭所有图片
Adc=1; %直流分量幅度
A1=; %频率F1信号的幅度
A2=;%频率F2信号的幅度
F1=20; %信号1频率(Hz)
F2=60; %信号2频率(Hz)
Fs=256;%采样频率(Hz)
P1=0;%信号1相位(度)
P2=0; %信号相位(度)
N=256; %采样点数
t=[0:
1/Fs:
N/Fs];%采样时刻
%信号
S=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);
%显示原始信号
plot(S);
title('原始信号');
figure;
Y=fft(S,N);%做FFT变换
Ayy=(abs(Y));%取模
plot(Ayy(1:
N));%显示原始的FFT模值结果
title('FFT模值');
figure;
Ayy=Ayy/(N/2); %换算成实际的幅度
Ayy
(1)=Ayy
(1)/2;
F=([1:
N]-1)*Fs/N;%换算成实际的频率值
plot(F(1:
N/2),Ayy(1:
N/2)); %显示换算后的FFT模值结果
title('幅度-频率曲线图');
figure;
Pyy=[1:
N/2];
fori=1:
N/2
Pyy(i)=phase(Y(i));%计算相位
Pyy(i)=Pyy(i)*180/pi;%换算为角度
end;
plot(F(1:
N/2),Pyy(1:
N/2)); %显示相位图
title('相位-频率曲线图');
实例三
closeall;%先关闭所有图片
Adc=2;%直流分量幅度
A1=3;%频率F1信号的幅度
F1=50;%信号1频率(Hz)
P1=-30;%信号1相位(度)
A2=;%频率F2信号的幅度
F2=75;%信号2频率(Hz)
P2=90;%信号相位(度)
Fs=512;%采样频率(Hz)
N=1024;%采样点数
t=[0:
1/Fs:
(N-1)/Fs];%采样时刻
%信号
S=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);
%%显示原始信号
plot(S);
title('原始信号');
%%FFT变换后
figure;
Y=fft(S,N);%做FFT变换
Ayy=(abs(Y));%取模
plot(Ayy(1:
N));%显示原始的FFT模值结果
title('FFT模值');
%%幅度频率曲线图
figure;
Ayy=Ayy/(N/2);%换算成实际的幅度
Ayy
(1)=Ayy
(1)/2;
F=([1:
N]-1)*Fs/N;%换算成实际的频率值
plot(F(1:
N/2),Ayy(1:
N/2));%显示换算后的FFT模值结果
title('幅度-频率曲线图');
%%相位频率曲线图
figure;
Pyy=[1:
N/2];
fori=1:
N/2
Pyy(i)=phase(Y(i));%计算相位
Pyy(i)=Pyy(i)*180/pi;%换算为角度
end;
plot(F(1:
N/2),Pyy(1:
N/2));%显示相位图
title('相位-频率曲线图');
实例四
关于FFT的相位谱
(2011-07-1311:
41:
56)
转载▼
标签:
相位谱正弦信号延拓进行it
分类:
机械技术
先看一下我收到的程序,作为研究对象的信号是这样产生的:
T=128;
N=128;
dt=T/N;
t=dt*(1:
N);
x=2*cos(2*t-pi/4);
...
(我觉得这个信号存在一点问题,因为t是从1开始的,所以它的初相应该和-pi/4有点差别吧。
)
为什么进行FFT,用angle得到相位-频率特性却不能反映这个信号的初始相位?
胡广书的《数字信号处理-理论、算法与实现(第二版)》第三章第八节《关于正弦信号抽样的讨论》,得出了关于正弦信号抽样的六个结论,最后总结了一个总的原则:
抽样频率应为信号频率的整数倍,抽样点数应包含整周期。
书中的结论五与采样频率和抽样点数有很大的关联。
结论五主要说只有满足了上面的那个总的原则,频谱泄漏才不会发生。
我想不光是幅度谱的频谱泄漏现象,抽样频率和抽样点数同样会对相位谱产生影响。
考虑一个无限长的正弦信号(相当于初相为-90°),如果我们截取它的整数个周期,然后对截短的信号进行周期延拓,则得到的延拓的信号与原无限长正弦没有区别。
现在我们再次对这个无限长的正弦进行截短,长度为个周期,然后对截短信号进行周期延拓,看看我们得到了什么?
下图、无限长正弦:
下图、截短信号
下图、对截短信号周期延拓:
可以看出,此时进行周期延拓得到的信号与原来的正弦信号大相径庭。
新的周期信号是一个周期的偶函数,原无限长正弦是一个周期的奇函数,两者奇偶性都不一样了,因此不能指望利用新的信号的DFT求出原信号的初相。
exp(-jωt)=cos(ωt)-jsin(ωt),进行变换的时候,若f(t)为实偶函数,则f(t)sin(ωt)就是奇函数,对一个奇函数在对称区间内积分只能得到0,因此实偶函数的傅立叶变换肯定是实的,对一个实数用angle求相位,当然相位是0。
而原正弦肯定是初相为-90°。
我想这就是问题所在,DFT就是DFS,只不过DFT先将有限长信号进行周期延拓,然后求DFS,再截取一个周期。
使用DFT,在有限的观测时间内采集信号的信息。
如若观测时间内正好得到了整数个正弦周期,则DFT的周期延拓可以不失真的表示原正弦,可是如果观测时间内得到的信号不是整数个周期,那么问题随之而来,就像上面的例子,观测时间内得到了个周期的正弦,然后进行周期延拓,显然乱了套。
如果满足了胡广书老师所总结的抽样条件,则对正弦的DFT谱无疑可以很好地反映初相,我写了两个例子:
第一个例子,信号只包含一个正弦:
t=linspace(0,,16);
x=cos(2*pi*t+pi/4);
X=fft(x);
stem(abs(X));
figure;
stem(angle(X)/pi*180);
幅度谱:
相位谱:
可以看见DFT相位谱第三个点对应正弦的相位,刚好是45°。
第二个例子信号中包含两个正弦:
t=linspace(0,,16);
x=cos(2*pi*t+pi/4)+2*cos(2*pi**t+pi/8);
X=fft(x);
stem(abs(X));
figure;
stem(angle(X)/pi*180);
幅度谱:
相位谱:
可以看见DFT相位谱第二个和第三个点对应两个正弦的相位,刚好是°和45°。
如果没有满足上面所说的条件,就会得到不准确的结果,有兴趣可试试下面的代码:
t=linspace(0,
x=cos(2*pi*t+pi/4);
X=fft(x);
stem(abs(X));
figure;
stem(angle(X)/pi*180);
如何克服这个问题我觉得这非常困难。
在不能预知信号频率的情况下,无法确定采样频率和观测点数。
也许可以先进行一次观测,通过幅度谱估计出正弦的频率,然后根据频率调整抽样频率,重新对信号进行采样,使采样符合上面所述的条件。
但是这样做有很多的问题,例如硬件可能不好实现。
而且虽然第二次调整了采样频率和抽样点数,可是初始相位已经无法得到了,因为第二次采样不可能再从零时刻开始。
Sandygreta同学说可以这样做,先以较高的抽样频率对信号进行采样,通过FFT幅度谱估算出正弦信号的频率,然后计算出满足抽样条件的最佳的抽样频率和观测时间,使抽样频率为正弦频率的整数倍(大于2倍),且观测时间内能正好得到整数个正弦周期。
然后对刚才采集的信号样本进行插值,接着使用计算出来的采样频率和观测时间对插值的结果重新采样,计算FFT,得到初始相位。
实例五
x=0:
.001:
1;fs=1000;
y=sin(2*pi*50*x);
N=64;
fork=1:
4
N=N*2;
M=fft(y,N);
Py=abs(M)*2/N;
f=fs*(0:
N/2)/N;%fs采样率
subplot(4,1,k);
stem(f,Py(1:
N/2+1));
xlim([0100]);
title(['N='num2str(N)]);
end
实例六
同频率信号滤波的问题
实例⑦
小波变换尺度相关性去噪程序
实例⑧
关于低通滤波器设计的问题