(三) 时域信号的频谱分析
首先要将从外设输入或采集的时域波形数据经抽样量化后,通过CFile类的Open(……)、Read(……)等成员函数将其读取到缓存中,并将其转化为复变量存放于复变量数组A中,同时需要验证以下数据量的长度是否为2的整数次幂,如若不是则必须用0来补齐,否则无法用"蝴蝶图"进行分解运算。下面代码用于完成对原始采样时域序列的快速傅立叶变换,A、M分别表示指向原始采样数据数组的指针和序列长度的2的整数次幂:
…… Complex U,W,T; int LE,LE1,I,J,IP; int N=(int)pow(2,M); //在此采用的是时间抽选奇偶分解方式,所以在参加运算前首先要对时间序列进行倒序 ReverseOrder(A,N); int L=1; while(L<=M) { LE=(int)pow(2,L); LE1=LE/2; U.Re=1.0f; U.Im=0.0f; W.Re=(float)cos(PI/(1.0*LE1));//计算W算子的值 W.Im=(float)-1.0*sin(PI/(1.0*LE1)); if(abs(W.Re)<1.0e-12) W.Re=0.0f; if(abs(W.Im)<1.0e-12) W.Im=0.0f; J=1; while(J<=LE1) { I=J; while(I<=N) { IP=I+LE1; T.Re=(float)A[IP-1].Re*U.Re-A[IP-1].Im*U.Im;//计算复数运算A*U T.Im=(float)A[IP-1].Re*U.Im+A[IP-1].Im*U.Re; A[IP-1].Re=(float)A[I-1].Re-T.Re;//计算复数运算A-T A[IP-1].Im=(float)A[I-1].Im-T.Im; A[I-1].Re+=T.Re;//计算复数运算A+T A[I-1].Im+=T.Im; I+=LE; } float temp=U.Re; U.Re=(float)U.Re*W.Re-U.Im*W.Im;//计算复数运算U*W U.Im=(float)temp*W.Im+U.Im*W.Re; J++; } L++; } …… |
上述代码执行完毕时,原先存放着时域数值的复变量数组内存放的就是经过分析后的频域值了,对此数据可以通过绘图将频域波形直观的显示出来,也可以将其存成数据文件,以备进一步使用。
四、 测试及运算结果分析
编译运行程序,打开一三角脉冲的数据文件,并将分析结果保存,该三角脉冲幅度为1,持续时间2毫秒,采样时抽样时间间隔是20微秒,延拓周期(数据记录长度)为10毫秒,采样点数目500点,取2的整数次幂512个样点。下附该三角脉冲频谱的计算结果及误差分析:
频率(Hz) FFT求得 X(f) 误差 0.00 1.00006E-03 1.00000E-03 6.10352E-08 100.00 9.67593E-04 9.67531E-04 6.14332E-08 200.00 8.75203E-04 8.75150E-04 6.25092E-08 300.00 7.36904E-04 7.36849E-04 6.39413E-08 400.00 5.72852E-04 5.72787E-04 6.52926E-08 500.00 4.05351E-04 4.05285E-04 6.61362E-08 600.00 2.54638E-04 2.54572E-04 6.61847E-08 700.00 1.35403E-04 1.35338E-04 6.53870E-08 800.00 5.47602E-05 5.46963E-05 6.39612E-08 900.00 1.20072E-05 1.19448E-05 6.23453E-08 1000.00 6.10719E-08 1.17757E-32 6.53870E-08 1100.00 8.05672E-06 7.99613E-06 6.05985E-08 1200.00 2.43706E-05 2.43095E-05 6.11450E-08 1300.00 3.93026E-05 3.92400E-05 6.25965E-08 1400.00 4.68226E-05 4.67581E-05 6.45128E-08 1500.00 4.50979E-05 4.50316E-05 6.62543E-08 1600.00 3.58664E-05 3.57992E-05 6.71930E-08 1700.00 2.30135E-05 2.29466E-05 6.69399E-08 1800.00 1.08697E-05 1.08042E-05 6.55073E-08 1900.00 2.74348E-06 2.68014E-05 6.33390E-08 2000.00 6.11826E-08 1.17757E-32 6.11826E-08 2100.00 2.25379E-06 2.19395E-06 5.98376E-08 2200.00 7.29243E-06 7.23256E-06 5.98625E-08 2300.00 1.25974E-05 1.25360E-05 6.13467E-08 2400.00 1.59746E-05 1.59107E-05 6.38421E-08 2500.00 1.62779E-05 1.62114E-05 6.64915E-08 2600.00 1.36254E-05 1.35571E-05 6.83226E-08 2700.00 9.16539E-06 9.09679E-06 6.86075E-08 2800.00 4.53216E-06 4.46500E-06 6.71550E-08 2900.00 1.21487E-06 1.15945E-06 6.44190E-08 注:在此,FFT运算结果都倍乘了系数10毫秒(0.01秒)。
在分析结果中产生了误差,是由于待分析的连续时间信号不具备离散性或周期性,也可能有无限长度。为了适应FFT方法的需要,对波形进行了抽样和截断,这样再用程序分析采样数据必然会引入误差,从分析结果可以看出,频率越高,误差波动也越大,此分析结果产生的误差在允许范围之内,是一个可以满意的近似。实践证明,本程序的算法是正确可靠的。
小结:
DFT尤其是FFT的应用已遍及各个科学领域,"DFT的应用"与 "FFT的应用"几乎成为同义语。通过本文介绍和程序示例可以清楚的看到FFT方法在直接处理离散信号数据的作用,而且也可以很好的用于对连续时间信号分析的逼近。本程序在Windows 2000 Professional下、由Microsoft Visual C++ 6.0编译通过。 
说明:本教程来源互联网或网友上传或出版商,仅为学习研究或媒体推广,wanshiok.com不保证资料的完整性。
2/2 首页 上一页 1 2 |