時(shí)丕旭，楊志強(qiáng)，計(jì)國(guó)鋒

(長(zhǎng)安大學(xué)地質(zhì)工程與測(cè)繪學(xué)院，陜西 西安710054)

一、引 言

GAT磁懸浮陀螺全站儀是一種采用高速旋轉(zhuǎn)的磁懸浮陀螺敏感地球自轉(zhuǎn)角動(dòng)量，并利用力矩平衡反饋法對(duì)陀螺角動(dòng)量軸線與子午線夾角進(jìn)行的精確測(cè)定以確定真北方向的一種精密定向儀器［1］，其觀測(cè)數(shù)據(jù)為兩組電流值，分別為定子電流數(shù)據(jù)和轉(zhuǎn)子電流數(shù)據(jù)。通常應(yīng)用于地下工程，如隧道、礦山、地鐵等的貫通測(cè)量工作，為地下工程提供精確的方位角檢核條件［2］。然而，由于施工現(xiàn)場(chǎng)環(huán)境的復(fù)雜性，在實(shí)際觀測(cè)過(guò)程中，通過(guò)傳感器采集的信號(hào)，受到諸如溫度、濕度、振動(dòng)等對(duì)傳感器的干擾，使得信號(hào)受到噪聲的污染，因此要得到高精度的陀螺數(shù)據(jù)通常需要對(duì)信號(hào)進(jìn)行去噪處理，進(jìn)而提高磁懸浮陀螺全站儀觀測(cè)結(jié)果的精度。

本文針對(duì)磁懸浮陀螺數(shù)據(jù)中的隨機(jī)噪聲，提出了一種新的時(shí)頻分析方法——希爾伯特-黃變換，引入了固有模態(tài)函數(shù)(IMF)的概念，在經(jīng)驗(yàn)?zāi)Ｊ椒纸?EMD)的基礎(chǔ)上，對(duì)每個(gè)IMF進(jìn)行Hilbert變換得到瞬時(shí)頻率，從而將信號(hào)精確表示為頻率-時(shí)間-能量(或者幅度)的分布，稱(chēng)為Hilbert譜。HHT是基于信號(hào)局部特征的，能對(duì)信號(hào)進(jìn)行自適應(yīng)、高效分解，特別適用于分析非線性、非平穩(wěn)信號(hào)，能最大程度去除信號(hào)的干擾成分，還原信號(hào)的有用信息。

二、希爾伯特-黃變換原理

希爾伯特-黃變換方法(HHT)處理數(shù)據(jù)主要包含兩個(gè)步驟:第1步是對(duì)原始數(shù)據(jù)進(jìn)行經(jīng)驗(yàn)?zāi)Ｊ椒纸?EMD)，分解完成后會(huì)得到有限階的固有模態(tài)函數(shù)(IMF);第2步是對(duì)分解后的每一階IMF函數(shù)進(jìn)行Hilbert變換，得到每階IMF在時(shí)頻平面上的能量分布圖，即Hilbert譜。

1．黃變換

黃變換是基于這樣的假設(shè):任何復(fù)雜信號(hào)都是由簡(jiǎn)單的固有模態(tài)函數(shù)組成的，每個(gè)模態(tài)可以是線性和平穩(wěn)的，也可以是非線性和非平穩(wěn)的，但這些模態(tài)都有一個(gè)共同的特點(diǎn)，那就是在整個(gè)信號(hào)長(zhǎng)度內(nèi)，每個(gè)模態(tài)具有相同數(shù)目的極值點(diǎn)和過(guò)零點(diǎn)，并且滿(mǎn)足“零均值”條件。這種方法的本質(zhì)是通過(guò)信號(hào)特征尺度來(lái)獲取固有模態(tài)函數(shù)。特征尺度定義為相鄰極值點(diǎn)的時(shí)間跨度。根據(jù)特征尺度的概念，黃變換采用“篩”的過(guò)程，將原始信號(hào)中各固有模態(tài)信號(hào)層層分解出來(lái)。

固有模態(tài)函數(shù)(IMF)必須滿(mǎn)足兩個(gè)條件:①在整個(gè)數(shù)據(jù)集中，極值數(shù)和過(guò)零點(diǎn)數(shù)相等或最多差1;②在任一點(diǎn)，有局部極大值構(gòu)成的包絡(luò)和有局部極小值構(gòu)成的包絡(luò)均值為零。經(jīng)驗(yàn)?zāi)Ｊ椒纸?EMD)是將信號(hào)分解成IMF的組合，經(jīng)驗(yàn)?zāi)Ｊ椒纸馐墙柚粩嘀貜?fù)的篩選來(lái)逐步找出IMF［3］，基本過(guò)程如下。

1)找出原始信號(hào)s(t)中的所有局部極大值以及局部極小值，接著利用三次樣條，分別將局部極大值串連成上包絡(luò)線與局部極小值串聯(lián)成下包絡(luò)線。

2)求出上下包絡(luò)線的平均值，得到均值包絡(luò)線m1(t)。

3)原始信號(hào)s(t)與均值包絡(luò)線相減，得到第一個(gè)分量h1(t)

4)檢查h1(t)是否符合IMF的條件。如果不符合，則返回步驟1)并且將h1(t)當(dāng)作原始信號(hào)，進(jìn)行第二次篩選。即

重復(fù)篩選k次

直到hk(t)符合IMF的條件，即得到第一個(gè)IMF分量c1(t)，即

5)原始信號(hào)s(t)減去c1(t)可得剩余量r1(t)，即

6)將r1(t)當(dāng)作新的資料，重新執(zhí)行步驟1)至步驟5)，得到新的剩余量r2(t)。如此重復(fù)n次

當(dāng)?shù)趎個(gè)剩余量rn(t)已成為單調(diào)函數(shù)，將無(wú)法再分解IMF時(shí)，整個(gè)EMD的分解過(guò)程完成。原始信號(hào)s(t)可以表示成n個(gè)IMF分量與一個(gè)平均趨勢(shì)分量rn(t)的組合，即

這樣，原始信號(hào)便分解成n個(gè)IMF和一個(gè)趨勢(shì)函數(shù)，即可將IMF做希爾伯特變換進(jìn)行瞬時(shí)頻率的分析。

2．希爾伯特變換

希爾伯特變換是信號(hào)分析中的重要工具，當(dāng)給定一個(gè)連續(xù)的時(shí)間信號(hào)x(t)，其Hilbert變換定義為

對(duì)每個(gè)IMF進(jìn)行Hilbert變換后，可得

這里省略了殘余函數(shù)r，Re表示取實(shí)部。展開(kāi)式為Hilbert幅值譜，簡(jiǎn)稱(chēng)Hilbert譜

HHT中的邊際譜可利用Hilbert譜進(jìn)行定義，即

而表示信號(hào)s(t)的第k個(gè)分量的邊際譜。

三、陀螺數(shù)據(jù)處理分析

陀螺儀觀測(cè)數(shù)據(jù)的核心是轉(zhuǎn)子電流與定子電流采樣數(shù)據(jù)。根據(jù)目前磁懸浮陀螺儀內(nèi)置程序設(shè)計(jì)，在一個(gè)尋北測(cè)回中，共采集40 000組轉(zhuǎn)子電流與定子電流數(shù)據(jù)，其中回轉(zhuǎn)前后的兩個(gè)精尋位置各采集20 000組數(shù)據(jù)，真北方向計(jì)算的關(guān)鍵數(shù)據(jù)就是這40 000組數(shù)據(jù)，由轉(zhuǎn)子電流數(shù)據(jù)和定子電流數(shù)據(jù)共同確定指向力矩來(lái)最后計(jì)算陀螺旋轉(zhuǎn)軸的北向偏角。為了驗(yàn)證希爾伯特-黃變換方法的可靠性，筆者在不同的觀測(cè)環(huán)境下采集磁懸浮陀螺數(shù)據(jù)，分析不同分布形態(tài)的轉(zhuǎn)子電流的濾波效果。

1．定子電流處理與分析

在觀測(cè)過(guò)程中，定子電流數(shù)據(jù)與轉(zhuǎn)子電流數(shù)據(jù)同步采樣。與轉(zhuǎn)子電流數(shù)據(jù)相比較，定子電流數(shù)據(jù)具有顯著的平穩(wěn)性，大量的數(shù)據(jù)分析顯示，定子電流數(shù)據(jù)對(duì)測(cè)量環(huán)境變化不敏感，從數(shù)據(jù)分布形態(tài)到數(shù)量級(jí)都具有明顯的一致性，表現(xiàn)出近似正態(tài)分布或近似白噪聲特征［4］。定子電流的這種特征與采集定子電流的定子固定于殼體上不能轉(zhuǎn)動(dòng)的物理特性是相符的。

通過(guò)大量的試驗(yàn)和分析可以得出，定子電流數(shù)據(jù)是一個(gè)均值非常平穩(wěn)、均方差較平穩(wěn)的數(shù)據(jù)序列;定子電流數(shù)據(jù)可以理解為是一個(gè)常值加白噪聲的分布形態(tài)。但定子電流中含有大量的高頻噪聲，對(duì)此本文采用希爾伯特-黃變換的濾波方法對(duì)定子電流數(shù)據(jù)進(jìn)行處理，來(lái)濾除其高頻噪聲。筆者取試驗(yàn)中的定子電流數(shù)據(jù)進(jìn)行EMD分解，得到15個(gè)IMF函數(shù)，從低階到高階依次濾除高頻噪聲，并計(jì)算去噪后的均值和均方差的值，見(jiàn)表1，各階的濾波結(jié)果對(duì)比圖如圖1所示。

圖1 各階濾波結(jié)果對(duì)比圖

表1 濾去各階IMF后的均值和均方差

從表1可以看出，無(wú)論去至哪一階的IMF，電流均值保持不變，這也說(shuō)明了電流均值的平穩(wěn)性。從圖1各階濾波結(jié)果對(duì)比圖和表1的均方差可以看出，依次濾除各階IMF后，電流分布圖逐漸變得平滑，均方差逐漸減小。最后濾波結(jié)果趨近于電流值為－0．063 2 A的直線，這和采集定子電流數(shù)據(jù)的定子運(yùn)動(dòng)規(guī)律和物理狀態(tài)是相符合的。通過(guò)濾波，筆者看到了剝離噪聲后定子電流數(shù)據(jù)的本質(zhì)，這說(shuō)明了EMD濾波［5］方法對(duì)于定子電流去噪的有效性。

因此，在最后計(jì)算方位角時(shí)，對(duì)定子電流首先進(jìn)行數(shù)據(jù)差分，剔除粗差，然后進(jìn)行去噪處理，利用去噪后的數(shù)據(jù)進(jìn)行方位角計(jì)算。

2．轉(zhuǎn)子電流的處理與分析

力矩器轉(zhuǎn)子位于靈敏部的下方，轉(zhuǎn)子電流會(huì)隨著靈敏部懸浮狀態(tài)和陀螺旋轉(zhuǎn)軸的變化而變化。因此，轉(zhuǎn)子電流的變化將反映反向力矩的變化。在不同的環(huán)境下，由于干擾力矩的不同導(dǎo)致獲取的反向力矩不同，從而采集的轉(zhuǎn)子電流數(shù)據(jù)時(shí)間序列以不同的走勢(shì)表現(xiàn)出來(lái)。常見(jiàn)的轉(zhuǎn)子電流的分布有3種情況:具有顯著周期的轉(zhuǎn)子電流數(shù)據(jù);具有不明顯周期和受外界干擾的轉(zhuǎn)子電流數(shù)據(jù);類(lèi)似于定子電流分布的轉(zhuǎn)子電流數(shù)據(jù)。

相對(duì)于定子電流數(shù)據(jù)的平穩(wěn)性及其白噪聲(近似白噪聲)影響的情形，轉(zhuǎn)子電流數(shù)據(jù)序列表現(xiàn)出明顯的系統(tǒng)運(yùn)動(dòng)特征與白噪聲(近似白噪聲)影響特征。更為嚴(yán)重的問(wèn)題是，由于轉(zhuǎn)子運(yùn)動(dòng)系統(tǒng)對(duì)環(huán)境因素變化的高度敏感性，甚至可能發(fā)生在觀測(cè)環(huán)境條件突變時(shí)，產(chǎn)生觀測(cè)數(shù)據(jù)突然跳變，表現(xiàn)出顯著的非線性特征。因此，本文采用針對(duì)非平穩(wěn)非線性數(shù)據(jù)的希爾伯特-黃變換處理方法。以下分3種情況進(jìn)行說(shuō)明。

(1)具有顯著周期的轉(zhuǎn)子電流數(shù)據(jù)

取試驗(yàn)中具有顯著周期的轉(zhuǎn)子電流數(shù)據(jù)進(jìn)行說(shuō)明，其原始數(shù)據(jù)分布如圖2所示。

圖2 具有顯著周期的轉(zhuǎn)子電流數(shù)據(jù)分布圖

從圖2中可以看出轉(zhuǎn)子電流數(shù)據(jù)具有明顯的周期性，說(shuō)明當(dāng)時(shí)的觀測(cè)條件良好，受到的干擾較小。數(shù)據(jù)表現(xiàn)出低頻波幅的特征，同時(shí)含有大量高頻噪聲。首先對(duì)其進(jìn)行EMD分解，濾除高頻噪聲。然后對(duì)分解產(chǎn)生的IMF函數(shù)進(jìn)行希爾伯特變換，建立每一階IMF函數(shù)所對(duì)應(yīng)的邊際譜，分析數(shù)據(jù)的波動(dòng)情況。

對(duì)該數(shù)據(jù)進(jìn)行EMD分解后共產(chǎn)生12個(gè)IMF函數(shù)，并對(duì)每一階的IMF建立邊際譜，如圖3所示。從上到下，頻率依次從高到低，時(shí)間尺度依次增大，邊際譜中的能量也依次減小，最后一個(gè)IMF為趨勢(shì)項(xiàng)，代表了數(shù)據(jù)序列的整體趨勢(shì)。

圖3 EMD分解后的IMF函數(shù)和邊際譜

圖3(續(xù))EMD分解后的IMF函數(shù)和邊際譜

邊際譜中某一頻率的能量存在則反映出在整個(gè)時(shí)間軸上有這樣一個(gè)頻率的振動(dòng)波在局部出現(xiàn)過(guò)，即某一頻率的能量出現(xiàn)時(shí)就表示一定有該頻率的振動(dòng)波出現(xiàn)［6］。從圖3可以看出，從上到下邊際譜中的能量依次減弱，到9階IMF函數(shù)的邊際譜時(shí)，已看不到其頻率的能量出現(xiàn)，即不存在振動(dòng)波。因此可以斷定，到9階IMF函數(shù)時(shí)，其噪聲已不存在，只需濾去前8階的IMF函數(shù)即可。

試驗(yàn)中的轉(zhuǎn)子電流數(shù)據(jù)表現(xiàn)出低頻周期波幅圖像，從上到下依次濾去高頻的IMF函數(shù)，濾波結(jié)果如圖4所示。從圖4中可以看出信號(hào)的噪聲逐漸被濾掉，信號(hào)變得越來(lái)越光滑，越來(lái)越清晰。每濾掉一階IMF函數(shù)，對(duì)剩下的信號(hào)進(jìn)行重構(gòu)［7］，并求取重構(gòu)信號(hào)的均值和均方差。

圖4 逐階濾波各層效果對(duì)比圖

表2 依次濾去各階IMF后重構(gòu)信號(hào)的均值與均方差

從表2中可以看出，去掉前8階IMF后，其均方差達(dá)到了±1．011 9×10－5，并且前8個(gè)均值基本保持一致。當(dāng)繼續(xù)再往下濾波時(shí)其均值發(fā)生了很大的變化，并且Matlab圖像也發(fā)生了很大的變化，與原始圖像的走勢(shì)發(fā)生很大的背離。因此對(duì)于該數(shù)據(jù)，筆者可以認(rèn)為前8階IMF函數(shù)為高頻噪聲數(shù)據(jù)序列，后4階低頻IMF函數(shù)包含了原始數(shù)據(jù)的絕大數(shù)有用信號(hào)。濾波最終結(jié)果如圖4所示。

(2)具有不明顯周期和受外界干擾的轉(zhuǎn)子電流數(shù)據(jù)

本文以試驗(yàn)中受干擾的轉(zhuǎn)子電流數(shù)據(jù)為例，該數(shù)據(jù)受到外界環(huán)境干擾較大，電流數(shù)據(jù)分布不規(guī)則。對(duì)數(shù)據(jù)進(jìn)行EMD分解，逐階對(duì)它進(jìn)行濾波，去掉相關(guān)的高頻噪聲數(shù)據(jù)。并逐階計(jì)算重構(gòu)信號(hào)的均值與均方差值，隨時(shí)觀察濾波后的圖像，并與原始圖像進(jìn)行對(duì)比。轉(zhuǎn)子電流數(shù)據(jù)經(jīng)EMD分解后產(chǎn)生了14個(gè)IMF函數(shù)(這里分解后的IMF函數(shù)和邊際譜不再顯示)。同樣，從低階到高階，時(shí)間尺度依次增大，代表了原始數(shù)據(jù)序列頻率從高頻到低頻的分布。從產(chǎn)生的邊際譜分析可知，到9階IMF函數(shù)的邊際譜時(shí)，振動(dòng)波已不存在，說(shuō)明已沒(méi)有噪聲，剩余的為有用信號(hào)，因此，去掉前8階的IMF函數(shù)即可。

圖5 濾波結(jié)果對(duì)比圖

圖5顯示了原始數(shù)據(jù)，去掉2階IMF函數(shù)，去掉5階IMF函數(shù)以及去掉8階IMF函數(shù)的轉(zhuǎn)子電流分布圖。從圖5中可以看到在數(shù)據(jù)波動(dòng)較大的區(qū)域，濾波結(jié)果的對(duì)比效果比較明顯，這說(shuō)明在數(shù)據(jù)波動(dòng)較大的采樣區(qū)域，噪聲嚴(yán)重干擾了原始信號(hào)，掩蓋了原始信號(hào)所要表達(dá)的信息。經(jīng)過(guò)EMD去噪，去除了高頻噪聲對(duì)原始信號(hào)的掩蓋，信號(hào)的有用成分得以還原［8］。

(3)類(lèi)似于定子電流分布的轉(zhuǎn)子電流數(shù)據(jù)

以試驗(yàn)中受干擾的轉(zhuǎn)子電流數(shù)據(jù)為例，該數(shù)據(jù)分布與定子電流類(lèi)似，說(shuō)明觀測(cè)條件良好，觀測(cè)過(guò)程中受到的干擾很小。通過(guò)對(duì)數(shù)據(jù)進(jìn)行EMD分解，分解為14個(gè)IMF函數(shù)，并建立每個(gè)IMF函數(shù)的邊際譜(這里分解后的IMF函數(shù)和邊際譜不再顯示)，從產(chǎn)生的邊際譜分析可知，到9階IMF函數(shù)的邊際譜時(shí)，振動(dòng)波已不存在，說(shuō)明已沒(méi)有噪聲，剩余的為有用信號(hào)，因此，去掉前8階的IMF函數(shù)即可。濾波結(jié)果如圖6所示。

圖6 濾波結(jié)果對(duì)比圖

雖然該轉(zhuǎn)子電流數(shù)據(jù)與定子電流分布圖相似，但卻與定子電流數(shù)據(jù)有本質(zhì)的區(qū)別。從濾波結(jié)果圖可以看出，類(lèi)似定子電流的轉(zhuǎn)子電流數(shù)據(jù)分布曲線以很小的幅度在均值附近做上下擺動(dòng)，而并不是像定子電流那樣以直線的形態(tài)完全集中于均值。這就是它和定子電流的本質(zhì)區(qū)別。

四、結(jié) 論

本文分析了磁懸浮陀螺數(shù)據(jù)的特性，針對(duì)轉(zhuǎn)子電流非線性非平穩(wěn)性的特征，提出了一種希爾伯特-黃變換的方法。該方法對(duì)IMF函數(shù)建立的邊際譜能很直觀地判斷噪聲的存在性，對(duì)于去除高頻噪聲，還原信號(hào)的有用成分是十分有效的。希爾伯特-黃變換非常適合于分析非線性非平穩(wěn)數(shù)據(jù)，完全自適應(yīng)，并且不受Heisenberg測(cè)不準(zhǔn)原理的制約。對(duì)陀螺數(shù)據(jù)的處理結(jié)果可以看出，希爾伯特-黃變換方法能有效去除隨機(jī)噪聲，還原信號(hào)的有用成分，最終提高陀螺的定向精度。

［1］ 石震，楊志強(qiáng)，楊帥，等．基于度盤(pán)配置法的陀螺定向新方法及相關(guān)問(wèn)題研究［J］．測(cè)繪通報(bào)，2011(2):87-89．

［2］ 李慧茹，楊志強(qiáng)，石震，等．Vondrak濾波在磁懸浮陀螺全站儀數(shù)據(jù)處理中的應(yīng)用［J］．地球科學(xué)與環(huán)境學(xué)報(bào)，2012，34(4):107-110．

［3］ 黃光南．希爾伯特-黃變換及其在地震資料分析處理中的應(yīng)用［D］．青島:中國(guó)海洋大學(xué)，2009．

［4］ 楊建華．GAT磁懸浮陀螺全站儀自動(dòng)化與陀螺儀觀測(cè)數(shù)據(jù)分析［D］．西安:長(zhǎng)安大學(xué)，2011．

［5］ 賈嶸，徐其惠，田錄林，等．基于經(jīng)驗(yàn)?zāi)B(tài)分解和固有模態(tài)函數(shù)重構(gòu)的局部放電去噪方法［J］．電工技術(shù)學(xué)報(bào)，2008，23(1):13-18．

［6］ 鐘佑明，秦樹(shù)人，湯寶平．希爾伯特-黃變換中邊際譜的研究［J］．系統(tǒng)工程與電子技術(shù)，2004，26(9):1323-1326．

［7］ 王英，曾光宇．基于小波的EMD去噪法應(yīng)用于心電信號(hào)去噪［J］．?dāng)?shù)值計(jì)算與計(jì)算機(jī)應(yīng)用，2011，32(4):274-282．

［8］ 張愛(ài)能．希爾伯特-黃變換方法在GAT磁懸浮陀螺以數(shù)據(jù)處理中的應(yīng)用［D］．西安:長(zhǎng)安大學(xué)，2014．