基于Matlab實(shí)驗(yàn)仿真對(duì)經(jīng)典功率譜估計(jì)性能分析

功率譜估計(jì)在分析平穩(wěn)各態(tài)遍歷隨機(jī)信號(hào)頻率成分領(lǐng)域被廣泛使用，并且已被成功應(yīng)用到雷達(dá)信號(hào)處理、故障診斷等實(shí)際工程中。本文給出了經(jīng)典功率譜估計(jì)的幾類(lèi)方法，并通過(guò)Matlab的實(shí)驗(yàn)仿真對(duì)經(jīng)典功率譜估計(jì)方法性能進(jìn)行了分析，最后說(shuō)明了經(jīng)典功率譜估計(jì)法的局限性和造成這種局限性的原因。

1.引言

給定一個(gè)標(biāo)準(zhǔn)的正弦信號(hào)，我們可以通過(guò)傅里葉變換來(lái)分析它的頻率成分。然而，實(shí)際工程應(yīng)用中，由于存在著各種干擾、噪聲，我們得到的信號(hào)往往不是理想的，如圖1-1這種信號(hào)，具有不確定性，幅度不能預(yù)知，非周期，但往往服從一定的統(tǒng)計(jì)特性，這種信號(hào)叫作隨機(jī)信號(hào)。需要注意的是，本文所說(shuō)的隨機(jī)信號(hào)是指平穩(wěn)各態(tài)遍歷的隨機(jī)信號(hào)，關(guān)于非平穩(wěn)隨機(jī)信號(hào)的分析方法[1]本文不予討論。

圖1-1 一種隨機(jī)信號(hào)時(shí)域形式

對(duì)于圖1-1的隨機(jī)信號(hào)，我們可以通過(guò)功率譜來(lái)分析它的頻率成分，如圖1-2所示為圖1-1隨機(jī)信號(hào)的功率譜。實(shí)際過(guò)程中，我們只能獲得隨機(jī)信號(hào)的一些離散數(shù)據(jù)點(diǎn)（假設(shè)為N個(gè)），本文將討論如何利用這N個(gè)數(shù)據(jù)點(diǎn)，來(lái)得到一個(gè)"非精確"的功率譜來(lái)對(duì)真實(shí)隨機(jī)信號(hào)的功率譜進(jìn)行估計(jì)，并討論如何更好的估計(jì)，即在下一章要講述的幾個(gè)經(jīng)典的功率譜估計(jì)法。

圖1-2 上圖所示的隨機(jī)信號(hào)功率譜

2.經(jīng)典功率譜估計(jì)法

上一章我們已經(jīng)知道功率譜估計(jì)法是通過(guò)利用已經(jīng)獲得的N個(gè)數(shù)據(jù)點(diǎn)，來(lái)得到一個(gè)"非精確"的功率譜對(duì)真實(shí)隨機(jī)信號(hào)的功率譜進(jìn)行估計(jì)，所以在給出具體的方法之前，如何來(lái)評(píng)價(jià)我們得出的這個(gè)"非精確"的功率譜的好壞呢？

評(píng)價(jià)功率譜性能好壞的標(biāo)準(zhǔn)有很多，本文只給出兩個(gè)影響最大的標(biāo)準(zhǔn)：分辨率和方差。分辨率即功率譜上能夠區(qū)分的最小相鄰頻率成分，分辨率越高，我們觀察信號(hào)的頻率成分越清晰；方差大小則反映到功率譜波動(dòng)性的大小，如果方差太大，功率譜波動(dòng)性大，則很容易造成有用的頻率成分被噪聲淹沒(méi)。所以，我們希望得到的這個(gè)"非精確"的功率譜，分辨率越高越好，方差越小越好。

2.1.2周期圖法性能（Matlab仿真）

上一小節(jié)我們已經(jīng)給出了周期圖法的原理。本節(jié)將通過(guò)Matlab仿真給出數(shù)據(jù)點(diǎn)數(shù)N對(duì)功率譜性能好壞的影響，正如上文所述，將通過(guò)對(duì)所得功率譜的分辨率和方差兩方面進(jìn)行分析。

圖2-1 實(shí)驗(yàn)所用的隨機(jī)信號(hào)

當(dāng)數(shù)據(jù)點(diǎn)數(shù)N分別為128、256、512和1024時(shí)，得到的功率譜分別如圖2-2、圖2-3、圖2-4和圖2-5所示。分辨率能夠直觀的通過(guò)功率譜圖形看出，方差的數(shù)值由表2-1給出。

表2-1 不同N值得到功率譜的方差值

通過(guò)上面實(shí)驗(yàn)結(jié)果的比較，我們很容易發(fā)現(xiàn)，周期圖法得到的功率譜隨著數(shù)據(jù)點(diǎn)數(shù)N的增大，分辨率變大、方差變也大。

2.1.3平均周期圖法

周期圖法得到的功率譜與我們所期望的"分辨率大、方差小"是矛盾的。為了進(jìn)一步降低方差，將N個(gè)觀測(cè)樣本數(shù)據(jù)點(diǎn)分為L(zhǎng)段，每段數(shù)據(jù)長(zhǎng)度為M, 分別對(duì)每段數(shù)據(jù)求周期圖功率譜估計(jì)，然后求平均值，這種方法稱(chēng)平均周期圖法。

2.1.4平均周期圖法性能（Matlab仿真）

當(dāng)數(shù)據(jù)點(diǎn)數(shù)N為1024，分段數(shù)分別為8、4、2時(shí)，平均周期圖法得到的功率譜分別如圖2-6、圖2-7、圖2-8所示。分辨率能夠直觀的通過(guò)功率譜圖形看出，方差的數(shù)值由表2-2給出。

表2-2 不同L值得到功率譜的方差值

L=1時(shí)，平均周期圖法退化為周期圖法。通過(guò)上面實(shí)驗(yàn)結(jié)果的比較，我們很容易發(fā)現(xiàn)，平均周期圖法得到的功率譜隨著分段數(shù)L變大，方差變小，但分辨率變小。

當(dāng)觀測(cè)樣本序列數(shù)據(jù)個(gè)數(shù)N固定時(shí)，要降低方差需要增加分段數(shù)L。當(dāng)N不大時(shí)分段長(zhǎng)度M取值較小，則功率譜分辨率降低到較低的水平。若分段數(shù)L固定時(shí)，增加分辨率需要分段長(zhǎng)度M，則需要采集到更長(zhǎng)的檢測(cè)數(shù)據(jù)序列。實(shí)際中恰恰是檢測(cè)樣本序列長(zhǎng)度不足。

2.1.5修正的平均周期圖法

上一節(jié)已經(jīng)提到實(shí)際中檢測(cè)樣本序列長(zhǎng)度是有限的。對(duì)現(xiàn)有數(shù)據(jù)長(zhǎng)度N，如果能獲得更多的段數(shù)分割，將會(huì)得到更小的方差。允許數(shù)據(jù)段間有重疊部分，來(lái)得到更多的段數(shù)。對(duì)段間重疊長(zhǎng)度的選取，最簡(jiǎn)單是取為段長(zhǎng)度M的一半。由式（2-5）可知更多的段數(shù)可以進(jìn)一步降低方差。

數(shù)據(jù)截?cái)嗟倪^(guò)程中相當(dāng)于數(shù)據(jù)加矩形窗，矩形窗幅度較大的旁瓣會(huì)造成"頻譜泄漏"。我們分段時(shí)采取的窗函數(shù)更為多樣(三角窗,海明窗等), 以減小截?cái)鄶?shù)據(jù)（加矩形窗）窗函數(shù)帶來(lái)的影響[2]

2.1.6修正的平均周期圖法性能（Matlab仿真）

利用修正平均周期圖法，分別使用矩形窗、Blackman窗和Hamming窗得到的功率譜如圖2-9所示。

圖2-9 不同窗函數(shù)的修正平均周期圖法得到的功率譜

可以發(fā)現(xiàn)，矩形窗的分辨率最高，但是方差也最大，這是由于矩形窗頻譜主瓣最窄，分辨率因此最高，旁瓣也高，導(dǎo)致頻譜泄漏最嚴(yán)重，方差最大。

2.1.7總結(jié)

周期圖法獲得的功率譜隨著樣本點(diǎn)數(shù)越多，分辨率越大、方差越大；平均周期圖法以犧牲分辨率來(lái)進(jìn)一步改善方差；修正的平均周期圖法允許段的重疊來(lái)進(jìn)一步增大分段數(shù)、或者分段數(shù)相同，每段樣本點(diǎn)數(shù)變多。無(wú)論是哪種方法都沒(méi)有徹底結(jié)局方差與分辨率之間的矛盾。

2.2相關(guān)功率譜估計(jì)法-BT法

正如我們之前介紹的，要提高功率譜估計(jì)的分辨率，必須增加數(shù)據(jù)序列的長(zhǎng)度N，但是較長(zhǎng)的數(shù)據(jù)序列，由噪聲引起的隨機(jī)性得到更為充分的體現(xiàn)－較大的方差。事實(shí)上，當(dāng)N無(wú)窮大時(shí)，方差為一非零常數(shù)。即周期圖法無(wú)法實(shí)現(xiàn)功率譜的一致估計(jì)。而這節(jié)講述的相關(guān)功率譜估計(jì)法（下文稱(chēng)作BT法），是一致估計(jì)。

2.2.1 BT法的原理

維納辛欽定理指出，隨機(jī)信號(hào)的相關(guān)函數(shù)與它的功率譜是一對(duì)傅里葉變換對(duì)。BT法就是基于這個(gè)原理。先由觀測(cè)數(shù)據(jù)估計(jì)出自相關(guān)函數(shù)，然后求自相關(guān)函數(shù)的傅立葉變換，以此變換作為對(duì)功率譜的估計(jì)，也稱(chēng)為間接法。BT法要求信號(hào)長(zhǎng)度N以外的信號(hào)為零，這也造成BT法的局限性。

2.2.2 BT法的性能（Matlab仿真）

數(shù)據(jù)點(diǎn)數(shù)N分別為128、256、512和1024的BT法，得到的功率譜如圖2-10、圖2-11、圖2-12和圖2-13所示。

圖2-13 N=1024時(shí)，BT法得到的功率譜

由上面實(shí)驗(yàn)可以發(fā)現(xiàn)，M隨著N的增大而增大時(shí)，分辨率提高，方差變大。BT法仍然沒(méi)有解決分辨率與方差之間的矛盾，但是BT法得到的功率譜當(dāng)N為無(wú)窮大時(shí)，方差會(huì)趨向于零，即為一致估計(jì)[2]。

2.2.3 周期圖法與BT法的關(guān)系

結(jié) 論

本文通過(guò)Matlab仿真，以一個(gè)具體的隨機(jī)信號(hào)為例，簡(jiǎn)單介紹了周期圖法、平均周期圖法、修正的平均周期圖法以及BT法的基本原理，并對(duì)這些方法的性能進(jìn)行分析?？梢钥闯?，無(wú)論是周期圖法及其改進(jìn)算法還是BT法都沒(méi)有從根本上解決分辨率與方差的矛盾。經(jīng)典功率譜估計(jì)是利用傅里葉變換估計(jì)功率譜，而我們之前分析隨機(jī)信號(hào)不滿(mǎn)足傅里葉變換的條件，所以經(jīng)典功率譜估計(jì)方法不得不從無(wú)限長(zhǎng)數(shù)據(jù)點(diǎn)截取有限長(zhǎng)數(shù)據(jù)點(diǎn)，加入限制條件（周期圖法實(shí)際上假定N點(diǎn)外數(shù)據(jù)周期重復(fù)、BT法假定N點(diǎn)外數(shù)據(jù)為零）來(lái)"強(qiáng)制"作傅里葉變換，這也是造成它局限性的原因。

參考資料

[1]朱哲,鐘宏偉. 非平穩(wěn)隨機(jī)信號(hào)分析處理方法研究[J] 安徽電子信息技術(shù)學(xué)院學(xué)報(bào)2008.6:28-28

[2]皇甫堪.現(xiàn)代數(shù)字信號(hào)處理[M].電子工業(yè)出版社