鄭偉員，梁子安，周 俊，張穎華，金占領(lǐng)

(1 南京農(nóng)業(yè)大學(xué) 工學(xué)院,江蘇 南京 210031；2 連云港雙亞機(jī)械有限公司,江蘇 連云港 222000)

水田整地是水稻種植前的重要環(huán)節(jié)，精細(xì)的攪漿平地技術(shù)可以使稻田水層深淺均勻，節(jié)約農(nóng)業(yè)灌溉用水，滿足水稻種植“寸水不過田，寸水不漏泥”的要求；可以提高化肥、農(nóng)藥和除草劑的利用率，減少農(nóng)業(yè)生產(chǎn)成本，促進(jìn)水稻生產(chǎn)可持續(xù)發(fā)展；可以提高撒播、拋秧和插秧的作業(yè)質(zhì)量，提高勞動(dòng)生產(chǎn)率，是保證水稻高產(chǎn)、穩(wěn)產(chǎn)的重要措施[1-4]。

全球衛(wèi)星導(dǎo)航系統(tǒng)(Global navigation satellite system,GNSS)控制的水田攪漿平地技術(shù)是目前先進(jìn)的土地平整技術(shù)，相較于激光平地技術(shù)，GNSS 技術(shù)精細(xì)平整，受風(fēng)雪、雨霧等外界環(huán)境影響小，易于進(jìn)行大面積水田作業(yè)，能夠?qū)崿F(xiàn)地形地勢(shì)的測(cè)量[5]。但是受多路徑效應(yīng)誤差、拖拉機(jī)和攪漿機(jī)振動(dòng)誤差的影響，GNSS 接收到的高程數(shù)據(jù)存在不可避免的偏差，導(dǎo)致水田地形的測(cè)量結(jié)果與真實(shí)情況不一致，影響基準(zhǔn)面的準(zhǔn)確度，進(jìn)而降低水田攪漿平地機(jī)的作業(yè)質(zhì)量。

為提高GNSS 高程數(shù)據(jù)的精度與可信度，有必要對(duì)地形測(cè)量過程中采集到的數(shù)據(jù)進(jìn)行降噪處理。夏友祥等[6]分析了旱地農(nóng)田平整系統(tǒng)GNSS 定位數(shù)據(jù)誤差源，分別利用卡爾曼和小波變換濾波對(duì)多路徑效應(yīng)和振動(dòng)誤差進(jìn)行處理，提高了農(nóng)田平整GNSS定位精度；梁友斌等[7]采用小波硬閾值和卡爾曼融合濾波算法對(duì)水田旋耕平地機(jī)GNSS 傾角信號(hào)進(jìn)行去噪處理，減少了傾角信號(hào)中的噪聲成分。雖然小波變換降噪和卡爾曼濾波在GNSS 數(shù)據(jù)處理方面有一定的效果，但是小波基的選擇不同會(huì)導(dǎo)致降噪效果不一樣，小波基的選擇原則需要進(jìn)一步研究，且采用卡爾曼濾波算法需要建立精確的噪聲模型才能取得較優(yōu)的濾波效果，而對(duì)于水田平地機(jī)GNSS 高程數(shù)據(jù)來說很難建立準(zhǔn)確的噪聲模型。胡煉等[8]通過限幅加權(quán)遞推平均濾波算法減少了GPS高程數(shù)據(jù)誤差，濾波后的數(shù)據(jù)可以滿足平地要求，但用于濾波的GPS 高程定位數(shù)據(jù)的個(gè)數(shù)和權(quán)值的選擇缺少理論支撐。由于經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical mode decomposition，EMD)算法在處理非線性、非平穩(wěn)信號(hào)中具有優(yōu)勢(shì)，近年來其用于降噪的研究逐漸增多。戴吾蛟等[9]根據(jù)高斯白噪聲信號(hào)的本征模態(tài)函數(shù)(Intrinsic mode functions，IMFs)分量的能量密度與其平均周期的乘積為一常量這一特性來確定分解級(jí)數(shù)，并將其用于GPS 多路徑效應(yīng)研究，有效削弱了隨機(jī)誤差，得到了更準(zhǔn)確的多路徑效應(yīng)改正模型。羅飛雪等[10]分別利用EMD 低通濾波去噪法和小波分解對(duì)動(dòng)態(tài)變形監(jiān)測(cè)的GPS 數(shù)據(jù)進(jìn)行濾波處理，驗(yàn)證EMD 算法在GPS 動(dòng)態(tài)變形監(jiān)測(cè)數(shù)據(jù)去噪中的可行性與有效性。蔣卓蕓[11]利用EMD 去噪法對(duì)高光譜數(shù)據(jù)去噪，效果良好。以上EMD 降噪方法將信號(hào)分解為若干個(gè)頻率由高到低的IMFs，其中，高頻分量一般包含的是信號(hào)的噪聲成分，將其剔除并將剩余IMF 分量重構(gòu)即完成對(duì)信號(hào)的濾波。然而中頻分量中仍可能含有噪聲成分，因此傳統(tǒng)的EMD 降噪方法在一定程度上會(huì)造成濾波不徹底，信號(hào)仍然含有少量噪聲[12]。

為解決水田平地機(jī)測(cè)量地形時(shí)GNSS 高程數(shù)據(jù)存在噪聲進(jìn)而影響基準(zhǔn)面建立準(zhǔn)確性的問題，本文分析了基于GNSS 控制的水田攪漿平地機(jī)測(cè)量高程數(shù)據(jù)時(shí)的誤差來源，在傳統(tǒng)EMD 濾波的基礎(chǔ)上提出了一種EMD 與S-G 聯(lián)合降噪的方法，S-G濾波器同樣適用于分析非平穩(wěn)信號(hào)，且去除噪聲的同時(shí)可以最大限度保持原始信號(hào)波形峰值的高度和寬度。聯(lián)合降噪方法首先采用EMD 將原始高程信號(hào)分解成若干IMFs，利用歸一化自相關(guān)函數(shù)和相關(guān)系數(shù)將其細(xì)分為噪聲IMFs、混合IMFs 和有效IMFs，然后利用S-G 算法對(duì)混合IMFs 進(jìn)行濾波，最后將S-G 濾波后的IMFs 與有效IMFs 進(jìn)行重構(gòu)，得到最終降噪后的數(shù)據(jù)。通過靜態(tài)試驗(yàn)與拉繩位移傳感器測(cè)量數(shù)據(jù)進(jìn)行對(duì)比，驗(yàn)證算法的可行性，并進(jìn)行了田間試驗(yàn)進(jìn)一步說明聯(lián)合濾波算法的有效性，以期減小水田平地機(jī)在測(cè)量地勢(shì)時(shí)存在的多路徑效應(yīng)和機(jī)器振動(dòng)影響，從而提高水田平地作業(yè)質(zhì)量。

1 系統(tǒng)平臺(tái)

基于GNSS 的水田攪漿平地機(jī)的系統(tǒng)硬件組成如圖1 所示，系統(tǒng)主要由懸掛架、攪漿機(jī)、平地機(jī)構(gòu)、液壓系統(tǒng)、自動(dòng)控制系統(tǒng)5 部分組成。懸掛架和拖拉機(jī)三點(diǎn)懸掛裝置掛接，平地機(jī)構(gòu)設(shè)置于攪漿機(jī)后方，平地機(jī)構(gòu)為平行四邊形結(jié)構(gòu)，液壓系統(tǒng)固定于懸掛架上方，自動(dòng)控制系統(tǒng)主要包括RTK GNSS 系統(tǒng)、集成控制柜和觸摸顯示屏，其中，RTK GNSS 系統(tǒng)包括基站、移動(dòng)站接收機(jī)、主副衛(wèi)星接收天線和差分信號(hào)天線。

圖1 系統(tǒng)硬件組成Fig.1 System hardware components

使用120 馬力的拖拉機(jī)作為牽引動(dòng)力，采用連云港雙亞機(jī)械有限公司生產(chǎn)的攪漿機(jī)(作業(yè)幅寬為3 m)，試驗(yàn)中使用的導(dǎo)航系統(tǒng)是上海華測(cè)導(dǎo)航公司生產(chǎn)的P3-DT 北斗高精度定位測(cè)向接收機(jī)，水平和垂直測(cè)量精度分別為(0.80±0.01)和(1.50±0.01) cm，接收數(shù)據(jù)頻率為10 Hz。選用深圳市米蘭特科技有限公司的MPS 拉繩位移傳感器，量程為1 000 mm，輸出信號(hào)為RS485 數(shù)字信號(hào)，線性精度為 ± 0.1%，重復(fù)性精度為 ± 0.01%，在程序中設(shè)置GNSS 和位移傳感器同步接收數(shù)據(jù)。

2 誤差分析與算法描述

2.1 誤差來源

影響全球衛(wèi)星導(dǎo)航系統(tǒng)定位精度的原因有多種，如衛(wèi)星星歷誤差、衛(wèi)星鐘差、電離層與對(duì)流層折射誤差、多路徑效應(yīng)等。其中，前4 種相關(guān)性較強(qiáng)的誤差可以通過差分技術(shù)和誤差糾正模型消除，但相關(guān)性較弱的多路徑效應(yīng)難以通過上述方法消除[13-14]。因此，多路徑效應(yīng)是目前影響GNSS 高定位精度的主要誤差源[15]。

基于GNSS 自動(dòng)控制的水田攪漿平地機(jī)，移動(dòng)站衛(wèi)星接收天線通過強(qiáng)磁鐵吸附于平地機(jī)構(gòu)兩端，拖拉機(jī)和攪漿機(jī)啟動(dòng)后產(chǎn)生的無規(guī)則振動(dòng)噪聲會(huì)疊加到接收機(jī)采集的高程數(shù)據(jù)中，引起測(cè)量誤差[16]。為了提高系統(tǒng)的精度，需要對(duì)多路徑效應(yīng)引起的數(shù)據(jù)漂移和機(jī)械振動(dòng)產(chǎn)生的隨機(jī)噪聲進(jìn)行濾波處理。

2.2 算法及原理

2.2.1 EMD 算法 經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical mode decomposition，EMD)是一種自適應(yīng)的時(shí)頻分析方法，該方法根據(jù)信號(hào)自身的特點(diǎn)將信號(hào)分解成若干個(gè)頻率由高到低的本征模態(tài)函數(shù)(Intrinsic mode functions，IMFs) 和1 個(gè)殘余分量，其中，高頻IMFs 分量一般為信號(hào)的尖銳部分，即噪聲；低頻IMFs 分量一般為信號(hào)的有效成分。EMD 算法適用于分析非線性、非平穩(wěn)GNSS 高程信號(hào)[17-19]。

提取本征模態(tài)函數(shù)步驟如下[20]：找出原始信號(hào)x(t)的所有極大值點(diǎn)和極小值點(diǎn)，利用3 次樣條函數(shù)擬合極大值點(diǎn)和極小值點(diǎn)分別得到上包絡(luò)線emax和下包絡(luò)線emin，計(jì)算均值包絡(luò)m1(t)、原信號(hào)x(t)與均值包絡(luò)m1(t)的差值h1(t)。

判斷h1(t)是否滿足以下2 個(gè)條件：一是其極值點(diǎn)個(gè)數(shù)和過零點(diǎn)個(gè)數(shù)相差不大于1；二是其任意點(diǎn)處的上下包絡(luò)均值為0。若不滿足，則重復(fù)上述步驟，直至得到滿足IMFs 條件的模態(tài)分量c1(t)。

用原信號(hào)x(t)減去c1(t)，得到1 個(gè)去除高頻成分的新信號(hào)r1(t)：

對(duì)r1(t)重復(fù)上述過程，得到第2 個(gè)IMF 分量c2(t)，反復(fù)進(jìn)行，直到最后一階分量(殘余分量) 至多存在1 個(gè)極值點(diǎn)，經(jīng)驗(yàn)?zāi)B(tài)分解結(jié)束。最終原始信號(hào)x(t)被分解為n個(gè)IMFs 分量和1 個(gè)殘余分量：

2.2.2 S-G 算法 S-G 濾波法連續(xù)選取噪聲信號(hào)x[n]中2M+1 個(gè)數(shù)據(jù)(濾波窗口)，利用一個(gè)階次為N的多項(xiàng)式進(jìn)行擬合，多項(xiàng)式的系數(shù)根據(jù)最小二乘法擬合誤差最小原則確定，將多項(xiàng)式曲線在濾波窗口中心點(diǎn)處的值作為降噪后的值，從左到右依次移動(dòng)窗口，直至完成所有數(shù)據(jù)的濾波[21-22]。構(gòu)造擬合多項(xiàng)式：

擬合殘差為：

令ε對(duì)ak的所有偏導(dǎo)數(shù)為0，即：

化簡(jiǎn)得：

引入矩陣A=其中an,i=ni,-M≤n≤M,0 ≤i≤N，再引入輔助矩陣B=ATA，令：

可得：

式中，H第1 行即為所求系數(shù)。

2.2.3 EMD 與S-G 聯(lián)合濾波算法 EMD 與S-G聯(lián)合濾波算法如圖2 所示，主要分為3 個(gè)步驟。

圖2 EMD 與S-G 聯(lián)合濾波算法流程Fig.2 Flow chart of EMD and S-G combined filter algorithm

利用EMD 算法將含有噪聲的高程數(shù)據(jù)分解為若干IMFs 和1 個(gè)殘余分量；通過求解IMFs 的歸一化自相關(guān)函數(shù)以及與原始信號(hào)的相關(guān)系數(shù)，將分解得到的IMFs 分為噪聲IMFs、混合IMFs 和有效IMFs；將噪聲IMFs 舍棄，對(duì)混合IMFs 進(jìn)行S-G 濾波，最后將S-G 濾波后的“干凈”信號(hào)和有效IMFs 進(jìn)行重構(gòu)得到最終降噪高程信號(hào)。

本文利用隨機(jī)噪聲和一般信號(hào)的歸一化自相關(guān)函數(shù)特性的不同來區(qū)分噪聲信號(hào)和其他信號(hào)的IMFs 分量界限[23]。如圖3 所示，隨機(jī)噪聲信號(hào)在零點(diǎn)處歸一化自相關(guān)函數(shù)值最大，在其他點(diǎn)迅速衰減至0 附近，而對(duì)于一般信號(hào)，其最大值也出現(xiàn)在零點(diǎn)，但由于信號(hào)間存在相關(guān)性，在其他點(diǎn)并沒有迅速衰減到很小的值。

圖3 不同信號(hào)歸一化自相關(guān)函數(shù)Fig.3 Normalized autocorrelation function for different signals

本文利用相關(guān)系數(shù)(r) 區(qū)分有效IMFs 與其他IMFs，通過計(jì)算每個(gè)IMF 分量與原始信號(hào)的相關(guān)系數(shù)，把第1 個(gè)r大于0.5 及后面的低頻IMFs 作為有效IMFs[24]。r表達(dá)式為：

式中，cov(,xi)表示與xi的協(xié)方差，和 σx分別表示與xi的標(biāo)準(zhǔn)差。

3 試驗(yàn)與結(jié)果分析

3.1 靜態(tài)驗(yàn)證試驗(yàn)

3.1.1 試驗(yàn)方案 靜態(tài)試驗(yàn)前，將GNSS 基站放置于田間地頭，移動(dòng)站衛(wèi)星接收天線對(duì)稱磁吸于平地機(jī)構(gòu)兩側(cè)，差分信號(hào)天線吸附于液壓閥塊頂部，位移傳感器安裝在調(diào)高油缸上，如圖4 所示。采用靜態(tài)調(diào)高性能測(cè)試試驗(yàn)采集的數(shù)據(jù)進(jìn)行驗(yàn)證，在調(diào)節(jié)過程中，GNSS 和線性位移傳感器同時(shí)記錄數(shù)據(jù)，將線性位移傳感器采集的數(shù)據(jù)經(jīng)過轉(zhuǎn)換后的高度作為真值，若經(jīng)過聯(lián)合濾波后的GNSS 高程曲線和線性位移對(duì)應(yīng)的高程吻合度高，則說明EMD 與S-G聯(lián)合濾波可以有效地剔除GNSS 原始高程數(shù)據(jù)中的多路徑效應(yīng)和拖拉機(jī)振動(dòng)的影響。

圖4 靜態(tài)驗(yàn)證試驗(yàn)平臺(tái)Fig.4 Static verification test platform

3.1.2 數(shù)據(jù)算法處理 根據(jù)“2.2”，首先將主副衛(wèi)星接收天線的原始中間點(diǎn)高程數(shù)據(jù)采用EMD 算法進(jìn)行分解得到7 個(gè)IMFs 分量與1 個(gè)趨勢(shì)項(xiàng)，如圖5所示；然后對(duì)EMD 分解的7 個(gè)IMFs 分量進(jìn)行歸一化自相關(guān)函數(shù)求解，結(jié)果如圖6 所示，從圖6 可以看出IMF1 和IMF2 為噪聲分量，將其剔除；其次分別計(jì)算IMF1～I(xiàn)MF7 分量與原始高程信號(hào)的相關(guān)系數(shù)，結(jié)果依次為0.116 2、0.122 1、0.202 4、0.187 6、0.262 0、0.615 9 和0.657 5，由計(jì)算結(jié)果可知，IMF6 和IMF7 為有效分量，對(duì)IMF3～I(xiàn)MF5 進(jìn)行SG 濾波，最后將濾波后的IMF3～I(xiàn)MF5 分量與有效IMF 分量及殘余分量重構(gòu)得到聯(lián)合濾波后的高程信號(hào)。為更好地說明聯(lián)合濾波的有效性，另外設(shè)置單獨(dú)采用EMD 濾波方案和單獨(dú)采用S-G 濾波作為對(duì)比。3 種方案的降噪結(jié)果如圖7 所示。

圖5 原始高程數(shù)據(jù)進(jìn)行EMD 分解Fig.5 EMD decomposition of raw elevation data

圖6 IMFs 的歸一化自相關(guān)函數(shù)Fig.6 Normalized autocorrelation function of IMFs

圖7 3 種方案降噪結(jié)果對(duì)比圖Fig.7 Denoising results of three schemes

3.1.3 結(jié)果分析 由圖7a、7b 可以看出單獨(dú)使用EMD 和S-G 濾波雖然可以抑制一部分多路徑效應(yīng)誤差和振動(dòng)誤差，但數(shù)據(jù)仍有較多噪聲，信號(hào)依然存在“毛刺”現(xiàn)象，誤差較大。由圖7c 可以看出EMD 與S-G 聯(lián)合濾波降噪較充分，數(shù)據(jù)最為平滑，大大提高了GNSS 高程數(shù)據(jù)的精度。選用均方根誤差(Root mean square error，RMSE) 和信噪比(Signal to noise ratio，SNR)進(jìn)行濾波效果評(píng)價(jià)，結(jié)果如表1 所示。

表1 不同算法的降噪效果統(tǒng)計(jì)Table 1 Statistics of noise reduction effect of different algorithms

RMSE 是指濾波后信號(hào)與真實(shí)信號(hào)方差的平方根，用來衡量濾波后信號(hào)與真實(shí)信號(hào)之間的偏差，RMSE 越小表示降噪效果越好，其表達(dá)式為：

式中，xi表示真實(shí)信號(hào)，表示濾波后信號(hào)。

SNR 是指真實(shí)信號(hào)能量與噪聲能量的比值，信噪比愈高表示降噪效果愈好，其表達(dá)式為：

式中，Ps、Pn分別表示真實(shí)信號(hào)功率和噪聲信號(hào)功率。

由表1 可知經(jīng)EMD 與S-G 聯(lián)合濾波后的均方根誤差最小，信噪比最大，其次為S-G 濾波方法，EMD 相較于前2 種濾波方法的降噪效果差。表1說明經(jīng)聯(lián)合降噪后的高程信號(hào)與拉繩位移傳感器高程波形保持了較高的相似，聯(lián)合濾波具有較優(yōu)的效果。

3.2 田間試驗(yàn)結(jié)果

在江蘇省灌云縣五圖河農(nóng)場(chǎng)選取一塊面積為20 m × 50 m 的水田，使用GNSS 控制的水田攪漿平地機(jī)進(jìn)行水田高程信息的采集。作業(yè)前，將設(shè)備按照“3.1.1”的方法安裝，測(cè)量水田耕作層高程時(shí)，通過改變液壓系統(tǒng)中電磁閥閥芯的位置將調(diào)高油缸設(shè)置為浮動(dòng)模式，采用蛇形行走的方式進(jìn)行離散點(diǎn)高程數(shù)據(jù)的采集，在地頭拐彎處不采集數(shù)據(jù)。對(duì)高程數(shù)據(jù)濾波分幅寬進(jìn)行，限于篇幅，只選擇第1 幅寬數(shù)據(jù)進(jìn)行分析，結(jié)果如圖8 所示。

圖8 第1 幅寬原始高程數(shù)據(jù)和聯(lián)合濾波處理效果圖Fig.8 The raw elevation data and the effect of combined filtering of the first width

由圖8 分析可知，經(jīng)EMD 與S-G 濾波聯(lián)合降噪后的高程數(shù)據(jù)得到較優(yōu)平滑效果的同時(shí)保持了與原始數(shù)據(jù)幾乎相同的變化趨勢(shì)，說明噪聲去除充分，且數(shù)據(jù)沒有出現(xiàn)失真現(xiàn)象，信號(hào)中的有效信息得到保留。聯(lián)合濾波后最大值與最小值均有所改善，最大值由17.107 m 下降到17.100 m，最小值由17.006 m 上升到17.011 m，高程差由10.1 cm 下降到8.9 cm，數(shù)據(jù)波動(dòng)范圍減小了11.9%。通過對(duì)田間測(cè)量數(shù)據(jù)的分析進(jìn)一步說明了 EMD 和S-G 濾波聯(lián)合對(duì)水田平地機(jī)GNSS 高程信號(hào)降噪的優(yōu)越性。

4 結(jié)論

本文基于水田平地機(jī)測(cè)量的GNSS 實(shí)際高程數(shù)據(jù)，分析了GNSS 高程數(shù)據(jù)誤差來源，提出了一種EMD 與S-G 聯(lián)合降噪的方法削弱多路徑效應(yīng)和振動(dòng)對(duì)定位精度的影響，更加準(zhǔn)確地反映水田信息。

靜態(tài)驗(yàn)證試驗(yàn)聯(lián)合濾波后的RMSE 為3.95 mm，比濾波前降低了36.9%，SNR 為67.465 3 dB，比濾波前提高了6.3%。田間測(cè)量試驗(yàn)濾波后的數(shù)據(jù)波動(dòng)范圍減少了11.9%。試驗(yàn)結(jié)果表明，本文所提的EMD 與S-G 聯(lián)合濾波算法有效削減了多路徑效應(yīng)誤差和振動(dòng)誤差，改善了數(shù)據(jù)的平滑度，對(duì)于水田平地前基準(zhǔn)面的建立及提高水田平地作業(yè)質(zhì)量具有現(xiàn)實(shí)意義。