馬嘯宇，張金生，*，郝亮亮，李婷，王少博，李琳

（1.火箭軍工程大學(xué) 精確制導(dǎo)與仿真實(shí)驗(yàn)室，西安710025； 2.火箭軍96863部隊(duì)，洛陽471000）

在近年來的導(dǎo)航技術(shù)發(fā)展過程中，衛(wèi)星導(dǎo)航和慣性導(dǎo)航技術(shù)一直是飛行器導(dǎo)航技術(shù)的主流發(fā)展方向。但衛(wèi)星導(dǎo)航的抗干擾能力較差，慣性導(dǎo)航的誤差隨時(shí)間不斷累積。地磁場隨時(shí)間變化緩慢且不受氣候、地域影響，且由于其無需主動(dòng)發(fā)射探測信號、抗干擾能力強(qiáng)，地磁導(dǎo)航技術(shù)作為一種輔助導(dǎo)航方式逐漸受到關(guān)注［1-2］。

地磁匹配導(dǎo)航技術(shù)主要使用實(shí)時(shí)采集的地磁場數(shù)據(jù)與該區(qū)域地磁基準(zhǔn)圖進(jìn)行匹配［3-4］，從而實(shí)現(xiàn)定位。因此地磁基準(zhǔn)圖的構(gòu)建精確度對定位精度的好壞起著決定性作用。目前構(gòu)建地磁基準(zhǔn)圖的方法主要有2種：一是根據(jù)現(xiàn)有的地磁場物理模型進(jìn)行構(gòu)建；二是根據(jù)實(shí)測地磁場數(shù)據(jù)，構(gòu)建網(wǎng)格化地磁基準(zhǔn)圖［5］。現(xiàn)有國際地磁參考磁場（IGRF）和世界磁場模型（WMM）主要描述的是地球主磁場模型［6］。地磁場模型依賴衛(wèi)星監(jiān)測數(shù)據(jù)，地球內(nèi)部異常場隨海拔高度升高而衰減，因此地磁場模型中幾乎不含異常場信息，用于刻畫局部地磁場時(shí)誤差很大，通過實(shí)測數(shù)據(jù)獲得的地磁基準(zhǔn)圖具有更大優(yōu)勢。目前，地磁基準(zhǔn)圖的構(gòu)建主要集中于插值法的研究，常用的方法有：三次樣條插值、Kriging插值［7］、基于粒子群優(yōu)化（PSO）算法的Kriging（PSO-Kriging）插值［8］等。

現(xiàn)有插值方法存在的共同問題是在實(shí)測數(shù)據(jù)量較小時(shí)誤差較大，而壓縮感知理論能夠在采集少量數(shù)據(jù)的同時(shí)幾乎完全重構(gòu)出原始信號，可以彌補(bǔ)在低采樣率下插值方法的缺陷。因此本文以現(xiàn)有的地磁場實(shí)測數(shù)據(jù)為基礎(chǔ)，利用信號處理中的壓縮感知理論進(jìn)行地磁基準(zhǔn)圖的構(gòu)建，研究適用于地磁信息特征的壓縮感知算法，并驗(yàn)證其在實(shí)測數(shù)據(jù)量較小時(shí)相比于傳統(tǒng)插值法的優(yōu)越性。

1 壓縮感知理論

壓縮感知是2006年由Candes［9］和Donoho［10-11］等提出的具有革命性的信號處理理論。壓縮感知理論指出，當(dāng)信號在某個(gè)變換域上是稀疏的，就可以投影為低維觀測向量，以遠(yuǎn)低于奈奎斯特理論的采樣頻率進(jìn)行信號采集，同時(shí)能夠高精度甚至完全重構(gòu)出原始信號［12］。其核心思想是在信號發(fā)送端進(jìn)行壓縮采樣，并在信號接收端進(jìn)行高精度重建的過程，目的是在保證信號精度的同時(shí)大幅度提升信號采集的效率并降低采集成本［13］。壓縮感知理論主要分為稀疏表示和信號重建2個(gè)部分。設(shè)一個(gè)信號x（x∈RN×1），在基底Ψ 下表示為信號α，即x＝Ψα。若α中只有K個(gè)元素非零或其他元素趨近于零，則稱該信號在基底Ψ 下是K稀疏的，Ψ 稱為稀疏字典。

若能找到合適的觀測矩陣Φ∈RM×N（K＜M?N）對信號x進(jìn)行采樣，則可得到M 維觀測信號為

式中：A為傳感矩陣。即已知觀測信號y，求解出稀疏信號α，從而得到初始信號x。

觀測信號的維度M大于α的稀疏度K，因此可證明存在唯一解，即在已知壓縮矩陣Ψ 和觀測矩陣Φ時(shí)，根據(jù)低維信號y必然可解出原始高維信號x。因此實(shí)際上該問題可轉(zhuǎn)化為在一定稀疏度條件下求解最小L0范數(shù)解的問題，即

式中：b為原始信號；A＝ΦΨ。顯然Ψ 和Φ 共同決定了變換后信號的稀疏性，Ψ 與Φ 的差異越大，得到的解α就越稀疏。因此定義Ψ 與Φ 的互相關(guān)性為

2 稀疏字典

現(xiàn)有常用的稀疏字典有離散傅里葉基、離散余弦基、離散小波變換基、Curvelet基、冗余字典等［15］。

2.1 離散傅里葉變換

離散傅里葉變換將二維信號利用一系列周期函數(shù)進(jìn)行表示，實(shí)際上是從時(shí)域轉(zhuǎn)化到頻域進(jìn)行分析，離散傅里葉基可表示為

2.2 離散余弦變換

2.3 離散小波變換

離散小波變換是對于連續(xù)小波的尺度和平移進(jìn)行離散化，在圖像處理方面有顯著的應(yīng)用效果。本文采用Haar小波變換作為小波基。二維小波變換過程如圖1所示。

圖1 二維小波變換過程Fig.1 Two-dimensional wavelet transform process

3 觀測矩陣

3.1 高斯隨機(jī)矩陣

高斯隨機(jī)矩陣是在壓縮感知中應(yīng)用最廣泛的一種觀測矩陣，由于其具有很強(qiáng)的隨機(jī)性，只需觀測維度M滿足：

式中：c為常數(shù)即可高概率地滿足RIP條件。因此在壓縮感知算法仿真中一般使用高斯隨機(jī)矩陣作為觀測矩陣以達(dá)到較好的實(shí)驗(yàn)效果。但由于構(gòu)造矩陣時(shí)數(shù)據(jù)的連續(xù)性和隨機(jī)性，硬件上實(shí)現(xiàn)困難，因此難以得到實(shí)際應(yīng)用。

3.2 隨機(jī)伯努利矩陣

隨機(jī)伯努利矩陣的性質(zhì)類似于高斯隨機(jī)矩陣，由于其中每一個(gè)元素都服從伯努利分布，M≥cK lg（N／K）時(shí)以接近1的概率滿足RIP條件。相比于高斯隨機(jī)矩陣，由于其元素都為±1，更有利于硬件的實(shí)現(xiàn)。

3.3 稀疏隨機(jī)矩陣

稀疏矩陣構(gòu)造方法為：在一個(gè)M×N的全零陣基礎(chǔ)上，將每一行隨機(jī)d個(gè)位置置1，其中d∈｛4，8，10，16｝，在保證一定的重建精度的條件下易于進(jìn)行構(gòu)造和實(shí)際實(shí)現(xiàn)。

3.4 循環(huán)矩陣

循環(huán)矩陣的構(gòu)造方法為：首先構(gòu)造隨機(jī)向量β，即β＝（β1，β2，…，βN）∈RN，β中的元素取±1，而后經(jīng)過M 次移位操作構(gòu)造出完整的矩陣，同時(shí)矩陣中每一個(gè)元素服從伯努利分布。由于可通過移位寄存器實(shí)現(xiàn)矩陣構(gòu)造，具有較好的應(yīng)用前景。

4 仿真結(jié)果與性能分析

為研究適用于地磁基準(zhǔn)圖構(gòu)建的壓縮感知方法，本文比較了多種稀疏基和觀測矩陣在地磁基準(zhǔn)圖構(gòu)建中的效果。稀疏基選擇離散小波基和離散余弦基，觀測矩陣選擇高斯隨機(jī)矩陣、隨機(jī)伯努利矩陣、稀疏隨機(jī)矩陣、循環(huán)矩陣以及單位矩陣。本文以美國地質(zhì)調(diào)查局（United States Geological Survey，USGS）于2002年發(fā)布的北美地區(qū)地磁場異常數(shù)據(jù)［16］為基礎(chǔ)，結(jié)合IGRF計(jì)算得到的地磁主磁場強(qiáng)度，從而得到地磁總磁場實(shí)驗(yàn)數(shù)據(jù)。實(shí)驗(yàn)樣本為10張128×128網(wǎng)格具有不同起伏程度的地磁基準(zhǔn)圖，選取區(qū)域大小都為約50 km×50 km，以保證基本涵蓋全部類型的地磁特征；使用觀測矩陣維度、圖像處理中的峰值信噪比（Peak Signal to Noise Ratio，PSNR）以及絕對誤差（MAE）和均方根誤差（RMSE）作為評價(jià)標(biāo)準(zhǔn)。由于在信號處理端無需考慮實(shí)時(shí)性問題，信號重構(gòu)算法采用精度較高的壓縮采樣匹配追蹤（Co-SaMP）算法，通過不斷更新迭代原子提高精度。PSNR值計(jì)算公式為

4.1 稀疏矩陣分析

以高斯隨機(jī)矩陣作為觀測矩陣，檢驗(yàn)離散小波變換和離散余弦變換在不同維度觀測矩陣下的重建效果。圖2為離散小波變換和離散余弦變換重建效果隨觀測矩陣維度變化曲線，PSNR值為10張樣本圖分別重構(gòu)后的平均值?？梢钥闯?，當(dāng)觀測矩陣維度較大時(shí)離散小波變換重建效果較高，但隨著觀測矩陣維度的降低，其重建性能迅速降低，M＜50時(shí)幾乎無法重建信號；而對于離散余弦變換，M對重建效果影響較小，即使降到70，PSNR值依然可以保持在80 dB以上。

由于取的信號稀疏度K＝M／4，當(dāng)M 減小時(shí)信號稀疏度也隨之減小。經(jīng)過離散小波變換后的信號難以保證足夠的稀疏度，因此當(dāng)觀測矩陣維度降低時(shí)會(huì)丟失大量的信息，造成重建信號的嚴(yán)重失真；而離散余弦變換可以保證足夠的壓縮率，能夠在保證一定重建精度的條件下降低觀測信號的維度。

圖2 離散小波變換和離散余弦變換重建效果隨觀測矩陣維度的變化Fig.2 Variation of reconstruction effect of discrete wavelet transform and discrete cosine transform with dimension of measurementmatrix

4.2 觀測矩陣分析

為尋找適合地磁基準(zhǔn)圖的觀測矩陣，以離散余弦變換作為稀疏矩陣，分別以高斯隨機(jī)矩陣、隨機(jī)伯努利矩陣、稀疏隨機(jī)矩陣、循環(huán)矩陣、單位矩陣作為觀測矩陣檢驗(yàn)重建效果，實(shí)驗(yàn)效果如圖3所示。

圖3 不同觀測矩陣重建效果隨觀測矩陣維度的變化Fig.3 Variation of reconstruction effect with dimension of matrix for differentmeasurementmatrices

從圖3可以看出，當(dāng)觀測矩陣維度較大時(shí)，循環(huán)矩陣、隨機(jī)伯努利矩陣和高斯隨機(jī)矩陣都具有較好的信號重構(gòu)效果；但隨著維度的降低，高斯隨機(jī)矩陣、隨機(jī)伯努利矩陣、稀疏隨機(jī)矩陣、循環(huán)矩陣的重構(gòu)性能都迅速降低；在M＜55時(shí)，單位矩陣由于性能基本不受觀測矩陣維度影響逐漸體現(xiàn)出其明顯優(yōu)勢。因此在保證一定重構(gòu)精度條件下使用單位矩陣作為觀測矩陣能夠減小觀測信號的維度。

為檢驗(yàn)重建算法的穩(wěn)定性，將采樣率設(shè)為30%，對每種觀測矩陣進(jìn)行100次重建，重建效果如圖4所示。

以PSNR值的標(biāo)準(zhǔn)差（SD）作為判定重建算法穩(wěn)定性的標(biāo)準(zhǔn)，計(jì)算結(jié)果如表1所示，在30%采樣率下，利用單位矩陣作為觀測矩陣，在重建精度和穩(wěn)定性都優(yōu)于其他觀測矩陣。

圖4 不同觀測矩陣重建穩(wěn)定性Fig.4 Reconstruction stability of different measurementmatrices

表1 不同觀測矩陣重建100次PSNR標(biāo)準(zhǔn)差Table 1 PSNR standard deviation of 100 tim es reconstruction for d ifferentm easurem ent m atrices dB

4.3 與其他重建算法比較

為檢驗(yàn)CoSaMP算法的效果，分別使用匹配追蹤（MP）算法、正交匹配追蹤（OMP）算法、分段正交匹配追蹤（StOMP）算法對30%采樣率下的地磁基準(zhǔn)圖進(jìn)行重建，以單位矩陣作為觀測矩陣，以PSNR值作為評價(jià)標(biāo)準(zhǔn)，實(shí)驗(yàn)結(jié)果如表2所示。

4.4 壓縮感知與插值法比較

為驗(yàn)證壓縮感知在地磁基準(zhǔn)圖構(gòu)建中的實(shí)際效果，將離散余弦變換及單位矩陣構(gòu)成的傳感矩陣和三次樣條插值、Kriging插值、PSO-Kriging插值分別在25%、18.75%、12.5%、6.25%采樣率下的基準(zhǔn)圖構(gòu)建效果進(jìn)行對比。

Kriging插值所采用的變異函數(shù)模型為指數(shù)模型，其具體形式為

式中：T0為基臺值；T為塊金值；h為到插值點(diǎn)的距離；3a為變程。PSO-Kriging插值變異函數(shù)同樣采用指數(shù)模型，PSO算法粒子總個(gè)數(shù)設(shè)置40個(gè)，迭代次數(shù)設(shè)置100次。

具體采樣方法為：設(shè)128×128的原始地磁基準(zhǔn)圖為P＝pij（i，j＝1，2，…，128），在25%采樣率（M／N）時(shí)，壓縮感知傳感矩陣維度選擇32×128，即M＝32在6.25%采樣率時(shí)，傳感矩陣維度選擇8×128，插值所用基準(zhǔn)圖由隨機(jī)下采樣獲得，即保證所有方法采樣數(shù)據(jù)量相同。以PSNR值作為評價(jià)指標(biāo)，結(jié)果如表3所示。當(dāng)采樣率為25%時(shí)，壓縮感知與插值法的PSNR值基本相近，壓縮感知精度略高于插值法；當(dāng)采樣率降到6.25%時(shí)，壓縮感知依然可以保持70 dB以上的PSNR值，明顯優(yōu)于同采樣率下的插值法。

圖5為用壓縮感知和Kriging插值重建后的基準(zhǔn)圖效果，可以看出Kriging插值后失去了很多細(xì)節(jié)信息，經(jīng)過壓縮感知后圖像與原始圖相似度更高。

為充分驗(yàn)證實(shí)驗(yàn)結(jié)果的有效性，利用文獻(xiàn)［5］中的絕對誤差、均方根誤差和標(biāo)準(zhǔn)差3種性能指標(biāo)進(jìn)行分析。

絕對誤差為

表2 不同重建算法結(jié)果Tab le 2 Results of differen t reconstruction algorithm sdB

表3 不同采樣率下各方案構(gòu)建基準(zhǔn)圖的PSNR值Tab le 3 PSNR values of reference m ap reconstructed by different schem es at differen t sam p ling rates

圖5 Kriging插值法和壓縮感知重構(gòu)效果對比Fig.5 Comparison of Kriging interpolation method and compressed sensing reconstruction effect

均方根誤差為

式中：ξi為真值與估計(jì)值之間的誤差。

為驗(yàn)證各種方法構(gòu)建效果隨采樣率變化情況，本文選取了25%、18.75%、12.5%、6.25%這4種采樣率進(jìn)行實(shí)驗(yàn)，計(jì)算結(jié)果分別如表4～表7所示。

從表4～表7可以看出，在各種采樣率下利用壓縮感知理論重構(gòu)地磁基準(zhǔn)圖性能都優(yōu)于插值法，且采樣率越低優(yōu)勢越明顯。

表4 25%采樣率下各性能指標(biāo)數(shù)據(jù)Table 4 Perform ance index data at 25% sam p ling rate nT

表5 18.75%采樣率下各性能指標(biāo)數(shù)據(jù)Table 5 Perform ance index data at 18.75%sam p ling rate nT

表6 12.5%采樣率下各性能指標(biāo)數(shù)據(jù)Table 6 Perform ance index data at 12.5%sam p ling rate nT

表7 6.25%采樣率下各性能指標(biāo)數(shù)據(jù)Tab le 7 Perform ance index data at 6.25%sam p ling rate nT

5 結(jié) 論

1）驗(yàn)證了壓縮感知理論在地磁基準(zhǔn)圖構(gòu)建中的可行性以及穩(wěn)定性。

2）確定了以離散余弦變換作為稀疏基、單位矩陣作為觀測矩陣，以CoSaMP算法作為重構(gòu)算法的基于壓縮感知的高精度地磁基準(zhǔn)圖構(gòu)建方法。

3）仿真結(jié)果表明，在6.25%采樣率下通過壓縮感知理論重構(gòu)的地磁基準(zhǔn)圖PSNR值、平均誤差、均方根誤差指標(biāo)都明顯優(yōu)于插值方法。

4）根據(jù)壓縮感知理論，要想完成對信號的壓縮感知，首先要進(jìn)行壓縮采樣，通過一個(gè)隨機(jī)矩陣控制采樣點(diǎn)選擇是否采集信號。這就需要重新設(shè)計(jì)一種地磁信息采集方案。

5）如何將壓縮感知理論實(shí)際應(yīng)用到地磁傳感器設(shè)計(jì)中并設(shè)計(jì)低成本、高重構(gòu)精度同時(shí)利于硬件實(shí)現(xiàn)的觀測矩陣將是今后的研究方向。