|
|
EDA365欢迎您登录!
您需要 登录 才可以下载或查看,没有帐号?注册
x
FFT是离散傅立叶变换的快速算法,可以将一个信号变换4 B& H6 o8 l. y
到频域。有些信号在时域上是很难看出什么特征的,但是如
: c7 |( ] s+ g% `2 n. _/ V6 Q果变换到频域之后,就很容易看出特征了。这就是很多信号
* R. ~1 t! P# ]: ?! V分析采用FFT变换的原因。另外,FFT可以将一个信号的频谱' k% m+ i- A7 W
提取出来,这在频谱分析方面也是经常用的。
4 Y$ P; T l/ ]* u. a- q 虽然很多人都知道FFT是什么,可以用来做什么,怎么去
- }* y: a' k/ T: L0 |* w做,但是却不知道FFT之后的结果是什意思、如何决定要使用5 B5 B) D8 x# h+ a
多少点来做FFT。% X/ K% T) d0 F( x [
, g+ v5 s$ @" J8 ?2 C 现在圈圈就根据实际经验来说说FFT结果的具体物理意义。
/ o- v8 e) T0 d" Y0 ~9 n一个模拟信号,经过ADC采样之后,就变成了数字信号。采样
$ v9 m1 E8 {3 i0 W5 t) J( e定理告诉我们,采样频率要大于信号频率的两倍,这些我就( N- {& ^+ J3 h" O: e$ S3 K" |
不在此罗嗦了。
' Y7 c# [4 S7 d/ A% i- e
: j0 F( A7 S. i; S3 {# X1 y0 h4 a 采样得到的数字信号,就可以做FFT变换了。N个采样点,
$ h1 H1 y6 U* Z" W# C经过FFT之后,就可以得到N个点的FFT结果。为了方便进行FFT
+ \: I/ @& l7 N运算,通常N取2的整数次方。
) x' N3 [; s/ I3 |/ N. i0 p
) M! D4 V% p5 a ^$ s 假设采样频率为Fs,信号频率F,采样点数为N。那么FFT
2 z _( R' r3 ?- _; B% D之后结果就是一个为N点的复数。每一个点就对应着一个频率$ q" E4 U4 D9 ]2 z! y- ]$ c
点。这个点的模值,就是该频率值下的幅度特性。具体跟原始+ X- k8 D7 F# ~+ n. i( W5 }) B
信号的幅度有什么关系呢?假设原始信号的峰值为A,那么FFT
( M$ G$ ]& \0 {- Q' c的结果的每个点(除了第一个点直流分量之外)的模值就是A! K, m7 |$ B1 y% ^5 `& Z
的N/2倍。而第一个点就是直流分量,它的模值就是直流分量4 ?( `4 v3 B0 R. |! O* w
的N倍。而每个点的相位呢,就是在该频率下的信号的相位。) |/ G: x7 }& U% c7 l5 q: \
第一个点表示直流分量(即0Hz),而最后一个点N的再下一个
/ O: s2 O: ~* i, d4 w$ d2 p点(实际上这个点是不存在的,这里是假设的第N+1个点,也 o1 A; X% N6 w1 d, O& z3 }" g
可以看做是将第一个点分做两半分,另一半移到最后)则表示8 e% D4 m8 S7 ~* S0 }
采样频率Fs,这中间被N-1个点平均分成N等份,每个点的频率$ H! h0 h) c) t& ~( C- n
依次增加。例如某点n所表示的频率为:Fn=(n-1)*Fs/N。
2 h! W J, m/ Y由上面的公式可以看出,Fn所能分辨到频率为为Fs/N,如果( y/ o8 r. s6 S
采样频率Fs为1024Hz,采样点数为1024点,则可以分辨到1Hz。6 p1 ^. t. V4 `
1024Hz的采样率采样1024点,刚好是1秒,也就是说,采样1秒
; A7 G2 Z. [! [/ W3 a; ]时间的信号并做FFT,则结果可以分析到1Hz,如果采样2秒时& S" ~' o0 q6 c+ [6 a
间的信号并做FFT,则结果可以分析到0.5Hz。如果要提高频率! x: O( D7 }' W. m& L. V
分辨力,则必须增加采样点数,也即采样时间。频率分辨率和* s2 ~8 ]* d0 O9 B) H6 D
采样时间是倒数关系。1 q9 h+ v t# H
假设FFT之后某点n用复数a+bi表示,那么这个复数的模就是( ]9 R- e3 z% W/ W
An=根号a*a+b*b,相位就是Pn=atan2(b,a)。根据以上的结果,4 n5 v- H+ V$ G$ A
就可以计算出n点(n≠1,且n<=N/2)对应的信号的表达式为:# _ `7 f' Q3 {1 } d+ i6 t
An/(N/2)*cos(2*pi*Fn*t+Pn),即2*An/N*cos(2*pi*Fn*t+Pn)。
3 N2 _4 J8 G4 @$ t; {/ Y对于n=1点的信号,是直流分量,幅度即为A1/N。
" K i, k6 N6 s0 L; u 由于FFT结果的对称性,通常我们只使用前半部分的结果,4 G9 A# Q6 a1 }2 P) J- b
即小于采样频率一半的结果。
4 l* s0 L7 W$ u4 c& L; h4 N) _ {/ k2 u8 B- F( b) @+ k
好了,说了半天,看着公式也晕,下面圈圈以一个实际的
& s3 \- F7 O, X, ~7 S信号来做说明。: S, ?( y- S& @! A( D5 s- i% |
, X! j2 a( j2 { ?9 o: \/ Z 假设我们有一个信号,它含有2V的直流分量,频率为50Hz、; ]- s3 D5 ^: J0 Z
相位为-30度、幅度为3V的交流信号,以及一个频率为75Hz、
5 b& U1 ?" c: M& F相位为90度、幅度为1.5V的交流信号。用数学表达式就是如下:
6 a; m q1 N& I V2 U) r1 @8 v6 r" M0 d0 D
S=2+3*cos(2*pi*50*t-pi*30/180)+1.5*cos(2*pi*75*t+pi*90/180)
6 e# G% }. m/ v) X1 v* g4 t
2 O/ R7 I5 c2 M) r* F 式中cos参数为弧度,所以-30度和90度要分别换算成弧度。
- G: C& b: A2 {+ F我们以256Hz的采样率对这个信号进行采样,总共采样256点。
( N6 ^) ?6 u% V4 J( @5 s按照我们上面的分析,Fn=(n-1)*Fs/N,我们可以知道,每两个
$ }" Q! a* u9 l- A6 T$ ^2 S点之间的间距就是1Hz,第n个点的频率就是n-1。我们的信号
" o) P1 N0 E8 d& T有3个频率:0Hz、50Hz、75Hz,应该分别在第1个点、第51个点、: T z' z# g4 D4 W2 z3 ^. a, i
第76个点上出现峰值,其它各点应该接近0。实际情况如何呢?
- k+ i1 H; Y" I我们来看看FFT的结果的模值如图所示。# U. d% }1 L5 l$ W3 f8 b! C8 a
' L" X* d9 A( l& p X+ X' w# t( h" G( L! \2 G2 Q
4 l% J0 {5 Y& `
图1 FFT结果/ [' Q; O8 F) g( V4 J9 U
从图中我们可以看到,在第1点、第51点、和第76点附近有
2 ] X ^- c S6 ?7 b比较大的值。我们分别将这三个点附近的数据拿上来细看:
1 T1 Q f5 h4 s! R, J/ B+ I* |1点: 512+0i. N3 h4 ^( {* ]; I
2点: -2.6195E-14 - 1.4162E-13i: ~( L0 _+ J+ V( b% p" b! }: |
3点: -2.8586E-14 - 1.1898E-13i2 f- b) z3 W. h; _# a
$ ]3 O, O0 d( W. E
50点:-6.2076E-13 - 2.1713E-12i
3 ]6 T$ v1 n* z5 k( a8 i51点:332.55 - 192i
; r; u' V+ b1 C52点:-1.6707E-12 - 1.5241E-12i- G9 R- @4 l) x3 f: Q
$ \* G7 o( J' `
75点:-2.2199E-13 -1.0076E-12i
* R8 Q I8 J8 R4 n76点:3.4315E-12 + 192i
% v1 `. L- C7 V+ X! V7 g77点:-3.0263E-14 +7.5609E-13i# c4 D x3 B& @) ?) l, e
4 p9 B Y/ Z% S# \ 很明显,1点、51点、76点的值都比较大,它附近的点值/ N4 A) s1 i9 M
都很小,可以认为是0,即在那些频率点上的信号幅度为0。8 W1 P4 @3 R6 b2 J# E. @. d. _6 C
接着,我们来计算各点的幅度值。分别计算这三个点的模值,
- A9 U2 p% {, z/ J2 C结果如下:
" B9 i* `0 ^/ ]% m' `( p- f* L1点: 512
+ Q& O; \! N* Q4 [9 d% ^51点:384
+ q6 [( R: X: g76点:192( m3 F/ r* w" V! n0 v
按照公式,可以计算出直流分量为:512/N=512/256=2;
9 x8 t/ K2 s+ ]8 ~50Hz信号的幅度为:384/(N/2)=384/(256/2)=3;75Hz信号的
8 H: ^: U4 ?1 Q( E" q! \幅度为192/(N/2)=192/(256/2)=1.5。可见,从频谱分析出来6 z7 X1 o9 d2 W7 Q) F0 m
的幅度是正确的。
% z. D) K8 q% x9 O' u: E q+ h# [ 然后再来计算相位信息。直流信号没有相位可言,不用管0 t+ c9 a* @ y) \1 {3 E0 E
它。先计算50Hz信号的相位,atan2(-192, 332.55)=-0.5236, |5 ]- F( a5 u
结果是弧度,换算为角度就是180*(-0.5236)/pi=-30.0001。再' A- ?8 d3 a9 h1 e
计算75Hz信号的相位,atan2(192, 3.4315E-12)=1.5708弧度,
; o( [' K# t# V) {3 w换算成角度就是180*1.5708/pi=90.0002。可见,相位也是对的。
( u7 v, J+ L7 L& r/ A7 M根据FFT结果以及上面的分析计算,我们就可以写出信号的表达
Y8 ^8 Q! b$ G9 c# }式了,它就是我们开始提供的信号。
I* z2 b- W8 y0 ?2 G4 f1 ~) @$ E+ |% H, K# {0 R) G4 U+ Y. y
总结:假设采样频率为Fs,采样点数为N,做FFT之后,某
$ O0 p" K$ N$ ~* `4 {% C6 r: I4 S一点n(n从1开始)表示的频率为:Fn=(n-1)*Fs/N;该点的模值
3 ?7 m0 s: p8 T/ |除以N/2就是对应该频率下的信号的幅度(对于直流信号是除以
e8 `7 p7 |8 I6 G- \- {N);该点的相位即是对应该频率下的信号的相位。相位的计算
5 q+ {- _/ ^# c _9 T: m8 Y' y可用函数atan2(b,a)计算。atan2(b,a)是求坐标为(a,b)点的角3 r% u" Q& {4 T0 [5 m/ } v( f1 x3 g
度值,范围从-pi到pi。要精确到xHz,则需要采样长度为1/x秒% h% k) _, M# F: I. `% O3 ^5 s! q
的信号,并做FFT。要提高频率分辨率,就需要增加采样点数,
; _" F* O0 @5 V! n0 h这在一些实际的应用中是不现实的,需要在较短的时间内完成$ o" a. Y, K: Y& f' b9 q% G
分析。解决这个问题的方法有频率细分法,比较简单的方法是0 i. S1 B+ i9 E! Z
采样比较短时间的信号,然后在后面补充一定数量的0,使其长度1 W! M4 P) M) h5 n* N
达到需要的点数,再做FFT,这在一定程度上能够提高频率分辨力。
3 g9 ?0 d; {1 J4 A8 I, P& H+ q具体的频率细分法可参考相关文献。
4 d# J; F; P }0 p, m& }) F
- Y3 t( y1 h" L" s$ I6 F! D[附录:本测试数据使用的matlab程序]
) L* p! F" P; T4 j. u0 xclose all; %先关闭所有图片
( \ K( y. Y5 `- Z, Y6 {+ pAdc=2; %直流分量幅度
0 @% z" u! }( Z0 f* xA1=3; %频率F1信号的幅度
! a% p0 [+ ?1 A% WA2=1.5; %频率F2信号的幅度0 l4 V- ^$ `: i/ x
F1=50; %信号1频率(Hz)
. D% x3 V4 W: b$ ]' ~F2=75; %信号2频率(Hz)
8 V7 |1 _7 v8 ]- J' q( T# o0 E; T, E7 w3 `Fs=256; %采样频率(Hz)
0 x, F* P) L# A/ LP1=-30; %信号1相位(度)
7 }( }" R- G8 G6 L5 D" ~) { U# \P2=90; %信号相位(度)
- s% l! O; j. @- Q( N; f5 PN=256; %采样点数
7 T& z6 e- X7 i. u! R. C: Ot=[0:1/Fs:N/Fs]; %采样时刻. R, t+ m4 ]: j; s$ f" J+ I2 u: g
: O. \, P# _- Q( n+ c! Q%信号
W- E0 f1 z0 n/ d# G+ B7 k) p) @" ?7 ]' HS=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);
) x. P/ Z1 U1 N$ |9 X% P%显示原始信号
2 y) @% f6 z3 ?, Y! Oplot(S);
0 m% T& p3 S# `9 o9 L, {6 ftitle('原始信号');1 u; X" E+ n; q6 H* @8 g! W
: x9 @: X8 e5 ]3 d& k, o% U Q& ^3 G
figure;9 y5 |9 ~- H6 k. [/ \! U
Y = fft(S,N); %做FFT变换$ H+ d3 p0 d r0 X- V& x
Ayy = (abs(Y)); %取模- p% y5 g$ r+ _# N0 m: Q, e
plot(Ayy(1:N)); %显示原始的FFT模值结果$ A9 P% o$ @. B" C
title('FFT 模值');5 n& @) l1 l$ R+ F! L
9 w+ K6 ~2 T7 q7 W$ T/ Rfigure;( n2 {0 u+ P8 t! x
Ayy=Ayy/(N/2); %换算成实际的幅度
, j3 S+ o" x( ^7 z' {' z+ tAyy(1)=Ayy(1)/2;6 G" x" K4 g# J% P( c2 F) J
F=([1:N]-1)*Fs/N; %换算成实际的频率值
" c9 U: U; m6 p/ f( q& [plot(F(1:N/2),Ayy(1:N/2)); %显示换算后的FFT模值结果# H/ d7 p1 _; G( b) T9 z$ t5 l
title('幅度-频率曲线图');
7 [5 g& U% N6 k( V ?, g% ?. C- {; r5 r. |' _1 Z
figure;" I/ W6 i& d, _7 X2 t% \
Pyy=[1:N/2];* H4 `2 F$ A1 p+ B5 T4 R" F+ ~/ Q
for i="1:N/2"- b3 o) |7 u, H7 o
Pyy(i)=phase(Y(i)); %计算相位
3 L. Z$ r' v4 M: C2 M) n. }$ T6 NPyy(i)=Pyy(i)*180/pi; %换算为角度 m3 O+ Q1 g) ^! A- i9 B
end;, X% T$ O& R1 |
plot(F(1:N/2),Pyy(1:N/2)); %显示相位图
- A" P. m M( I1 z, V2 ftitle('相位-频率曲线图');4 h- N! w0 {0 Q' R
8 R5 j1 A6 F+ ~8 E看完这个你就明白谐波分析了。+ p( F- Q# y9 T1 |5 k% m; x k
|
|