郭 剛,李 釗,李建軍,胡禮勇,董巍巍

(1.第二炮兵青州士官學(xué)校,山東青州262500;2.第二炮兵駐石家莊地區(qū)軍事代表室,河北石家莊050081)

0 引言

目前,以單片機為核心的手持式振動信號采集器配有簡潔的信號調(diào)理和顯示電路,能夠完成低頻振動信號的采集、儲存和傳輸,成本低廉,在機器振動狀況檢測領(lǐng)域應(yīng)用廣泛。但受單片機硬件資源的制約,其一般不具備現(xiàn)場對采樣數(shù)據(jù)進行FFT的能力[1,2]。然而,FFT是振動信號分析最常用的方法,對許多現(xiàn)場采集的信號數(shù)據(jù),只需進行少量采樣點的FFT,即可算出機器振動的主頻,進而快速判斷運轉(zhuǎn)狀態(tài),實現(xiàn)對機器運行狀態(tài)的快速檢測和故障初診,具有很大的實用價值。下面將著重闡述在以ATmega64L單片機為核心的手持式振動信號采集器上實現(xiàn)128點FFT的過程和方法。

1 按時間抽取基2FFT算法[3,4]

離散傅里葉變換(DFT)的本質(zhì)就是建立了以時間為自變量的信號與以頻率為自變量的頻譜函數(shù)之間的變換關(guān)系。設(shè)X(n)為N點有限長序列,其DFT公式為:

式中,

按時間抽取的基2FFT算法,利用旋轉(zhuǎn)因子的對稱性、周期性和可約性,根據(jù)輸入序列時間上的奇偶,將一個N點的DFT分解成2個N/2點的DFT,再把2個N/2點的DFT又可按照奇偶分成2個N/4的DFT,直到最終分解為N/2個2點的DFT運算。本算法的基本運算單元是2點蝶形運算結(jié)構(gòu)。當(dāng)N=2M時(這里M為自然數(shù)),整個運算共有M級,每級又有N/2個蝶形運算組成。蝶形運算計算流程如圖1所示。

圖1 蝶形運算單元

蝶形運算方程式為:

在處理采集到的實際數(shù)據(jù)時,輸入信號序列通常是實數(shù)。理論上使用N點復(fù)數(shù)FFT完成2N點實數(shù)序列FFT運算,可以提高運算速度并且減小存儲空間。但是變換完成后的結(jié)果轉(zhuǎn)化計算量較大,導(dǎo)致單片機整體運算速度不能顯著提高,又鑒于ATmega 64L單片機自帶 4 KB的 RAM,完全滿足128個32位浮點數(shù)進行FFT的存儲空間需求,所以此處直接使用128點復(fù)數(shù)FFT算法完成128個實數(shù)采樣數(shù)據(jù)的FFT運算,在運算時要把原實數(shù)數(shù)據(jù)作為復(fù)數(shù)的實部,而對應(yīng)復(fù)數(shù)虛部需全部置零。

2 ATmega64L單片機采集系統(tǒng)

該系統(tǒng)采用ATmega64L單片機作為控制核心。它是基于RISC結(jié)構(gòu)的8位低功耗CMOS微處理器,大部分指令都可以在一個時鐘周期內(nèi)完成,最大可以支持16 MHz的系統(tǒng)時鐘,片載64 KB的可編程Flash,2B的EEPROM,4KB的SRAM,53個通用I/O口線,32個通用工作寄存器和2個串口、SPI接口和TWI接口等豐富的外圍接口電路。整個振動信號采集系統(tǒng)的硬件框圖如圖2所示。

圖2 采集系統(tǒng)硬件

3 128點采樣數(shù)據(jù)實時FFT的實現(xiàn)

3.1 進行復(fù)數(shù)FFT計算的步驟及流程

使用ATmega64L單片機進行128點采樣數(shù)據(jù)FFT計算時,應(yīng)首先把輸入的128個實信號虛部置零,使其變?yōu)閺?fù)數(shù),然后將數(shù)據(jù)序列按碼位倒置的方式反序排列,再進行按時間抽取基2FFT運算,最后計算轉(zhuǎn)換后各復(fù)數(shù)的模并比較分析振動信號主頻,具體FFT計算流程圖如圖3所示。

圖3 128點 FFT流程

3.2 在ICCAVR中編程實現(xiàn)FFT的要點

3.2.1 反序排列算法的實現(xiàn)

反序排列是由FFT算法的性質(zhì)決定的。首先將128個數(shù)據(jù)點自然排列的十進制數(shù)轉(zhuǎn)換成7位二進制數(shù),再將這些二進制數(shù)的首位至末尾的順序顛倒過來重新?lián)Q算成十進制數(shù),即得到反序排列的結(jié)果。為減小單片機的運算負(fù)擔(dān),同時充分利用ATmega64L程序存儲器大的優(yōu)勢,使用MATLAB軟件算出128點反序排列結(jié)果,并以數(shù)組形式存放在單片機程序存儲器中,FFT運算時,單片機只需進行簡單的查表,即可快速完成數(shù)據(jù)反序操作。單片機存儲128個序列碼為:

const unsigned char inv_tab[128]={0,64,32,96,16,80,……,47,111,31,95,63,127}。

3.2.2 旋轉(zhuǎn)因子的計算和儲存

進行FFT的子函數(shù)運算時,實時計算旋轉(zhuǎn)因子的值,對單片機而言是十分繁重的任務(wù)。根據(jù)旋轉(zhuǎn)因子的復(fù)數(shù)形式,預(yù)先算出 cos(2πkn/N)和sin(2πkn/N)的值(取小數(shù)點后4位精度),并以數(shù)組的形式存儲在單片機的程序存儲器內(nèi),供FFT運算過程中隨時調(diào)用。單片機存儲旋轉(zhuǎn)因子實部和虛部的形式為:

const float cos _tab[64]={0,-0.0491,-0.0980,-0.1467,…,-0.9808,-0.9892,-0.9952,-0.9988};

const float sin _tab[64]={1,0.9988,0.9952,…,-0.1467,-0.0980,-0.049}。

3.2.3 FFT各級蝶形運算的實現(xiàn)

根據(jù)FFT算法的特性,使用3級循環(huán)完成各級蝶形運算[5]。

for(L=1;L＜=7;L++){//L控制蝶形運算的級數(shù)

b=1;i=L-1;while(i＞0)

{b=b*2;i--;}//b控制蝶形運算2點之間的距離

for(j=0;j＜=b-1;j++)

{p=1;i=7-L;

while(i＞0)/*p=pow(2,7-L)*j;

{p=p*2;i--;}

p=p*j;//計算旋轉(zhuǎn)因子虛實部在數(shù)組中的對應(yīng)序號p

for(k=j;k＜128;k=k+2*b)//完成本級N/2個蝶形運算

{TR=dataR[k];TI=dataI[k];temp=dataR[k+b];

dataR[k]=dataR[k]+dataR[k+b]*cos_tab[p]+dataI[k+b]*sin_tab[p];

dataI[k]=dataI[k]-dataR[k+b]*sin_tab[p]+dataI[k+b]*cos_tab[p];

dataR[k+b]=TR-dataR[k+b]*cos_tab[p]-dataI[k+b]*sin_tab[p];

dataI[k+b]=TI+temp*sin_ tab[p]-dataI[k+b]*cos_ tab[p];}

}

}

3.2.4 FFT結(jié)果的處理[6]

根據(jù)采樣定理,只要采樣頻率大于信號頻率的2倍,采樣得到的數(shù)字信號就可以做FFT變換了。N個采樣點經(jīng)過FFT之后就可得到N個點的FFT結(jié)果。為了方便進行FFT運算,通常N取2的整數(shù)次方。假設(shè)采樣頻率為fs,信號頻率為f,采樣點數(shù)為N。那么FFT之后結(jié)果就是一個為N點的復(fù)數(shù),每個點就對應(yīng)有一個頻率點。這個點的模值就是該頻率值下的幅度,而每個點的相位就是該頻率下的信號的相位。FFT運算后第一個點表示直流分量(即0 Hz),第n個點所表示的頻率為:fn=(n-1)*fs/N。

4 實驗分析

為了驗證程序的正確性,使用手持式振動信號采集系統(tǒng)對FG1617函數(shù)波形信號發(fā)生器產(chǎn)生的標(biāo)準(zhǔn)9 Hz矩形波進行數(shù)據(jù)采集,采樣頻率設(shè)為128 Hz(該系統(tǒng)可以在10～1 000 Hz之間靈活設(shè)置采樣頻率)。然后截取采集到的128個連續(xù)數(shù)據(jù)點,使用Lab VIEW編制的專用FFT處理程序和ATmega64L內(nèi)嵌的FFT程序分別進行分析,得到的幅頻結(jié)果如表1所示。

表1 FFT結(jié)果比較表

通過比較分析發(fā)現(xiàn),ATmega64L單片機對128點采樣數(shù)據(jù)進行FFT,計算所得到的諧波主頻和幅值與使用專用數(shù)據(jù)處理程序計算得到的主頻和幅值基本吻合,達(dá)到設(shè)計要求精度。

5 結(jié)束語

FFT是信號采集領(lǐng)域的重要數(shù)據(jù)處理工具。上述充分利用高性能低成本的ATmega64L單片機硬件資源和高速運算能力,結(jié)合查表手段,可以在220 ms時間內(nèi)完成128點的32位精度浮點FFT運算,而且經(jīng)過實驗驗證轉(zhuǎn)換結(jié)果正確、精度夠用,證明了該型號單片機實現(xiàn)FFT運算的可行性,為FFT實現(xiàn)方法開辟了新的應(yīng)用領(lǐng)域。在低成本手持振動信號采集儀器上集成FFT功能,通過對現(xiàn)場采集的振動信號進行實時FFT運算,算出被檢測機器的振動主頻和幅值,實現(xiàn)了對采樣數(shù)據(jù)的現(xiàn)場分析,滿足了設(shè)備運行狀態(tài)檢測和故障初判的要求。

[1]肖 飛.手持式振動測試系統(tǒng)的軟件開發(fā)[D].南京:東南大學(xué),2005:1-5.

[2]賈民平,劉玉春,鐘秉林,等.便攜式數(shù)據(jù)采集與工況監(jiān)測分析系統(tǒng)的研制[J].東南大學(xué)學(xué)大報,1997,27(2):99-102.

[3]胡廣書.數(shù)字信號處理—理論、算法與實現(xiàn)[M].北京:清華大學(xué)出版社,2003:177-198.

[4]文其林,白曉東,周 洪,等.2048點FFT在TMS320C240x定點DSP上的實現(xiàn)[J].微計算機信息,2006(5-2):159-160.

[5]肖宛昂.嵌人式系統(tǒng)中FFT算法研究[J].單片機與嵌入式系統(tǒng)應(yīng)用,2003(1):69-71.

[6]李全利,劉長亮.CCS上FFT運算的實現(xiàn)[J].自動化技術(shù)____與應(yīng)用,2009,28(2):61-63.