本發(fā)明涉及信息濾波應(yīng)用中干涉相位展開(kāi)算法領(lǐng)域，具體涉及一種基于無(wú)味信息濾波的相位展開(kāi)算法。

背景技術(shù)：

干涉合成孔徑雷達(dá)(InSAR)可以高精度、高可靠性地獲取地表三維信息和高程變化信息，被廣泛應(yīng)用在大地測(cè)繪、海洋監(jiān)測(cè)、火山監(jiān)測(cè)和地震檢測(cè)等領(lǐng)域。干涉相位展開(kāi)是干涉合成孔徑雷達(dá)(lnSAR)技術(shù)應(yīng)用中尤為關(guān)鍵的環(huán)節(jié)，其相位展開(kāi)精度直接影響著干涉合成孔徑雷達(dá)(InSAR)系統(tǒng)高程測(cè)量精度，故干涉圖展開(kāi)問(wèn)題受到越來(lái)越多的關(guān)注。相位展開(kāi)技術(shù)能從干涉圖模2π映射的相位主值區(qū)間恢復(fù)出真實(shí)干涉相位，廣泛應(yīng)用于數(shù)字全息顯微、磁共振成像處理、散斑成像、光學(xué)干涉測(cè)量、自適應(yīng)光學(xué)、合成孔徑雷達(dá)干涉測(cè)量等領(lǐng)域，目前基于數(shù)據(jù)融合的相位展開(kāi)算法主要包括擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)以及無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)，這些算法雖然改善了相位展開(kāi)算法的精度問(wèn)題，但精度還待進(jìn)一步提高。

技術(shù)實(shí)現(xiàn)要素：

針對(duì)上述現(xiàn)有技術(shù)，本發(fā)明要解決的技術(shù)問(wèn)題是如何進(jìn)一步提高相位展開(kāi)算法的精度。

為解決上述技術(shù)問(wèn)題，本發(fā)明提供的技術(shù)方案是一種基于無(wú)味信息濾波的相位展開(kāi)算法，包括以下過(guò)程：

(1)建立無(wú)味信息濾波(UIF)相位展開(kāi)算法遞推狀態(tài)估計(jì)模型，利用基于修正矩陣束模型(AMPM)的相位梯度估計(jì)方法獲取上述遞推狀態(tài)估計(jì)模型所需的相位梯度信息，隨后利用LEVENBERG-MARQUARDT方法優(yōu)化上述遞推狀態(tài)估計(jì)模型，提高算法收斂性；

(2)引入堆排序的快速質(zhì)量圖引導(dǎo)策略，把已展開(kāi)像元的鄰接纏繞像元作為待展開(kāi)像元嵌入堆數(shù)組，根據(jù)待展開(kāi)像元質(zhì)量值調(diào)整堆數(shù)組為最大堆；

(3)在每一展開(kāi)步驟中利用無(wú)味信息濾波相位展開(kāi)算法遞推狀態(tài)估計(jì)模型展開(kāi)堆數(shù)組根結(jié)點(diǎn)處的最佳待展開(kāi)像元，隨后從堆數(shù)組中刪除該像元，并調(diào)整堆數(shù)組為最大堆，直至最終完成所有纏繞像元的相位展開(kāi)；

其具體步驟如下：

步驟1定義干涉像元質(zhì)量圖；

步驟2創(chuàng)建基于完全二叉樹(shù)的堆數(shù)組，存儲(chǔ)待展開(kāi)像元；

步驟3選取相位質(zhì)量最高的非邊界像元作為起始像元，把其纏繞相位作為狀態(tài)估計(jì)值，并在(0，1)范圍內(nèi)設(shè)定其估計(jì)誤差方差，以起始像元為中心，把其鄰接的上下左右四個(gè)纏繞像元標(biāo)記為待展開(kāi)像元，分別把待展開(kāi)像元嵌入堆數(shù)組，并根據(jù)待展開(kāi)像元的相位質(zhì)量值，即把待展開(kāi)像元的相位質(zhì)量值作為堆排序時(shí)的關(guān)鍵字，把無(wú)序堆數(shù)組調(diào)整為最大堆；

步驟4在最大堆根結(jié)點(diǎn)處獲取相位質(zhì)量最高的最佳待展開(kāi)像元，標(biāo)記為像元x，并利用基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)此纏繞像元；從堆數(shù)組中去除像元x，并將堆數(shù)據(jù)組中剩余元素調(diào)整為最大堆；以像元x為中心，把其鄰接像元x的上下左右四個(gè)像元中的纏繞像元標(biāo)記為待展開(kāi)像元，分別把上述待展開(kāi)像元嵌入堆數(shù)組，并根據(jù)待展開(kāi)像元的相位質(zhì)量值，把無(wú)序堆數(shù)組調(diào)整為最大堆；

步驟5堆數(shù)組中是否存在元素？否，結(jié)束；否則轉(zhuǎn)步驟4。

過(guò)程(1)中所述建立無(wú)味信息濾波(UIF)相位展開(kāi)算法遞推狀態(tài)估計(jì)模型，利用基于修正矩陣束模型(AMPM)的相位梯度估計(jì)方法獲取上述遞推狀態(tài)估計(jì)模型所需的相位梯度信息，過(guò)程如下：

無(wú)味信息濾波器(UIF)是由無(wú)味卡爾曼濾波器(UKF)發(fā)展到信息濾波領(lǐng)域而來(lái)，同樣采用Sigma點(diǎn)來(lái)捕捉非線性系統(tǒng)狀態(tài)矢量的后驗(yàn)均值和方差，具有比無(wú)味卡爾曼濾波器(UKF)更高的計(jì)算效率。用一維坐標(biāo)k代替二維坐標(biāo)(m,n)，則一維相位展開(kāi)系統(tǒng)模型可表述為如下：

上式中，x(k)表示干涉圖k像元真實(shí)干涉相位，和w(k)分別為干涉圖k像元的相位梯度估計(jì)值及估計(jì)誤差，可利用基于修正矩陣束模型(AMPM)的局部相位梯度估計(jì)算法^[5]來(lái)獲?。粃(k+1)和v(k+1)分別為干涉圖k+1像元觀測(cè)值及其附加噪聲，v₁(k+1)和v₂(k+1)分別為附加在z(k+1)虛部和實(shí)部的噪聲。針對(duì)上述相位展開(kāi)系統(tǒng)模型，一維UIF相位展開(kāi)算法遞推狀態(tài)估計(jì)按如下進(jìn)行：

設(shè)干涉圖中k像元干涉相位的狀態(tài)估計(jì)值及其估計(jì)誤差方差為和其相應(yīng)的Sigma點(diǎn)為：

其中，n_x表示系統(tǒng)狀態(tài)矢量的維數(shù)，本文中n_x＝1；指是均方根矩陣的第j列向量，參數(shù)λ＝α²(n+κ)-n_x，α＝0.01，κ＝0。則干涉圖k+1像元干涉相位可按如下遞推估計(jì)：

(1)傳播Sigma點(diǎn)并計(jì)算狀態(tài)預(yù)測(cè)值和狀態(tài)預(yù)測(cè)方差：

其中，和分別為相應(yīng)的調(diào)節(jié)權(quán)值系數(shù)，其參考取值為現(xiàn)有技術(shù)；Q(k)為干涉圖k像元局部相位梯度估計(jì)誤差方差，其相應(yīng)計(jì)算取值為現(xiàn)有技術(shù)。

(2)計(jì)算信息矩陣Y(k+1)及其相應(yīng)的信息狀態(tài)向量

(3)傳播Sigma點(diǎn)并計(jì)算觀測(cè)預(yù)測(cè)值和互協(xié)方差

(4)計(jì)算信息狀態(tài)貢獻(xiàn)度及其相應(yīng)的信息矩陣：

其中，R(k+1)為干涉圖k+1像元量測(cè)噪聲方差，其相應(yīng)計(jì)算取值為現(xiàn)有技術(shù)。

(5)更新k+1像元的狀態(tài)估計(jì)值及其相應(yīng)的估計(jì)誤差方差：

其中，為干涉圖k+1像元狀態(tài)估計(jì)值，為干涉圖k+1像元估計(jì)誤差方差；

過(guò)程(1)中所述利用LEVENBERG-MARQUARDT方法優(yōu)化遞推狀態(tài)估計(jì)模型,其過(guò)程如下：

為了提高一維無(wú)味信息濾波(UIF)相位展開(kāi)算法的收斂性，用LEVENBERG-MARQUARDT方法對(duì)無(wú)味信息濾波(UIF)狀態(tài)預(yù)測(cè)方差進(jìn)行優(yōu)化，即通過(guò)參數(shù)u調(diào)整一維無(wú)味信息濾波(UIF)狀態(tài)預(yù)測(cè)方差；因此，一維無(wú)味信息濾波(UIF)相位展開(kāi)算法中的公式③需修正為如下：

上式中，I表示單位矩陣；用經(jīng)LEVENBERG-MARQUARDT方法優(yōu)化后的一維無(wú)味信息濾波(UIF)相位展開(kāi)算法逐行或逐列地展開(kāi)干涉相位，便可同時(shí)執(zhí)行纏繞相位展開(kāi)和相位噪聲抑制，完成干涉圖展開(kāi)工作。

基于堆排序的快速質(zhì)量圖引導(dǎo)策略結(jié)合，過(guò)程如下：

將經(jīng)LEVENBERG-MARQUARDT方法優(yōu)化后的一維無(wú)味信息濾波(UIF)相位展開(kāi)算法與基于堆排序的快速質(zhì)量圖引導(dǎo)策略結(jié)合起來(lái)實(shí)現(xiàn)快速二維相位展開(kāi)，即利用基于堆排序的快速質(zhì)量圖引導(dǎo)策略指導(dǎo)無(wú)味信息濾波相位展開(kāi)算法(UIFPU)從干涉圖中具有較高質(zhì)量干涉相位的區(qū)域到低質(zhì)量干涉相位的區(qū)域快速完成干涉相位的展開(kāi)工作，盡可能減小相位展開(kāi)過(guò)程中的誤差傳遞效應(yīng)。同時(shí)為了保證待展開(kāi)干涉相位的連續(xù)性，待展開(kāi)像元(m,n)可用其相鄰8像元中已展開(kāi)像元的狀態(tài)估計(jì)值來(lái)預(yù)測(cè)：

其中，像元(a,s)表示待展開(kāi)像元(m,n)相鄰8像元中的已展開(kāi)像元，G表示像元(m,n)的相鄰8像元中已展開(kāi)像元的集合；和分別表示像元(a,s)狀態(tài)估計(jì)值及其估計(jì)誤差方差；SNR(a,s)表示(a,s)像元的信噪比；為干涉圖(m,n)像元與(a,s)像元之間的相位梯度估計(jì)值，可利用基于修正矩陣束模型(AMPM)的局部相位梯度估計(jì)算法來(lái)獲取；κ(a,s)指已展開(kāi)像元(a,s)在待展開(kāi)像元預(yù)測(cè)中的權(quán)重；和分別指像元(m,n)一步預(yù)測(cè)估計(jì)及其估計(jì)誤差方差?？紤]到本文中狀態(tài)變量維數(shù)n_x＝1，在一維相位展開(kāi)算法中的公式②和公式③需做相應(yīng)的修正，其中公式②更改為：

其中公式③修正為：

其中，Q_(m,n)|(a,s)為干涉圖(m,n)像元與(a,s)像元之間的相位梯度估計(jì)誤差方差，其相應(yīng)計(jì)算取值為現(xiàn)有技術(shù)；此外，需要注意在二維相位展開(kāi)中需把一維相位展開(kāi)中的像元序號(hào)k+1變?yōu)槎S序號(hào)(m,n)。利用堆排序的快速質(zhì)量圖引導(dǎo)策略引導(dǎo)UIFPU快速地沿高質(zhì)量區(qū)域到低質(zhì)量區(qū)域的路徑對(duì)纏繞相位進(jìn)行遞推估計(jì)。

采用本發(fā)明的技術(shù)方案算法展開(kāi)精度明顯高于現(xiàn)有技術(shù)方案采用的常用算法如質(zhì)量引導(dǎo)算法、擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)以及無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)。

附圖說(shuō)明

圖1為基于無(wú)味信息濾波的相位展開(kāi)算法步驟框圖；

圖2為預(yù)濾波噪聲干涉圖；

圖3為預(yù)濾波質(zhì)量圖法展開(kāi)結(jié)果；

圖4為預(yù)濾波擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)展開(kāi)結(jié)果；

圖5為預(yù)濾波無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)展開(kāi)結(jié)果；

圖6為預(yù)濾波基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)結(jié)果；

圖7為無(wú)濾波噪聲干涉圖；

圖8為無(wú)濾波質(zhì)量圖法展開(kāi)結(jié)果；

圖9為無(wú)濾波擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)展開(kāi)結(jié)果；

圖10為無(wú)濾波無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)展開(kāi)結(jié)果；

圖11為無(wú)濾波基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)結(jié)果；

圖12為實(shí)測(cè)數(shù)據(jù)Etna火山地形干涉圖；

圖13為實(shí)測(cè)數(shù)據(jù)基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)結(jié)果。

具體實(shí)施方式

下面結(jié)合附圖對(duì)本發(fā)明作進(jìn)一步的說(shuō)明。

為了驗(yàn)證算法性能，將與常用算法，包括質(zhì)量圖引導(dǎo)算法、擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)以及無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)在模擬和實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)中進(jìn)行比較。

圖1為基于無(wú)味信息濾波的相位展開(kāi)算法步驟框圖。

圖2(a)為模擬真實(shí)干涉相位圖，圖2(b)為含噪聲纏繞干涉圖，其信噪比為4.94dB。在相位展開(kāi)前用均值濾波器(3×3)對(duì)含噪聲干涉圖進(jìn)行濾波，濾波后干涉圖見(jiàn)圖2(c)。

圖3-圖6分別為質(zhì)量圖引導(dǎo)算法、擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)、無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)和基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)圖2(c)的結(jié)果。

從圖3-圖6中，可以看出本發(fā)明和無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)的展開(kāi)結(jié)果較好，其誤差遠(yuǎn)小于質(zhì)量圖引導(dǎo)算法和擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)的誤差。質(zhì)量圖引導(dǎo)算法、無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)、基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)所采用的引導(dǎo)展開(kāi)路徑的相位質(zhì)量圖均為相位微分偏差圖。

表1各算法的均方根誤差

表1列出了各算法展開(kāi)不同信噪比干涉圖的結(jié)果，可以看出基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)的誤差不僅遠(yuǎn)小于質(zhì)量圖引導(dǎo)算法，且同樣小于同屬數(shù)據(jù)融合框架(亦即貝葉斯框架)下的擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)、無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)。

表2各算法的運(yùn)行時(shí)間

表2列出了上述算法在同一MATLAB計(jì)算環(huán)境下展開(kāi)圖1(c)所消耗的時(shí)間，可以看出本專利算法消耗的時(shí)間小于質(zhì)量圖引導(dǎo)算法、擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)以及無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)。故本發(fā)明能在較小的時(shí)間消耗范圍內(nèi)獲得較高精度的展開(kāi)相位圖。

無(wú)預(yù)濾波的干涉圖展開(kāi)實(shí)驗(yàn)：

為了進(jìn)一步檢驗(yàn)各算法性能，利用上述算法直接展開(kāi)噪聲纏繞相位圖；

圖7(a)為真實(shí)干涉相位圖(與圖2(a)相同)，圖7(b)為信噪比為3.01dB的噪聲纏繞相位圖。

圖8-圖11分別為質(zhì)量圖法、擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)、無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)和基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)圖7(b)的結(jié)果。

由于質(zhì)量圖引導(dǎo)算法不能有效抑制纏繞相位中相位噪聲，所以相位展開(kāi)誤差較大，見(jiàn)圖8；擴(kuò)展卡爾曼濾波相位展開(kāi)算法(EKFPU)同樣對(duì)噪聲比較敏感，故其相位展開(kāi)結(jié)果不理想，見(jiàn)圖9；

無(wú)味卡爾曼濾波相位展開(kāi)算法(UKFPU)和基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)從噪聲干涉圖中獲得了較好結(jié)果，且基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)相位展開(kāi)誤差更小，見(jiàn)圖10和圖11。

表3各算法的均方根誤差

表3列出了上述算法直接展開(kāi)不同信噪比的噪聲纏繞相位圖的結(jié)果，可以看出即使在干涉圖的信噪比較低的情況下，基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)仍具有較高的展開(kāi)精度。

實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)：

圖12(a)為部分Etna火山干涉圖。

圖12(b)為干涉圖相位殘差點(diǎn)。

圖13(a)和圖13(b)為基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)結(jié)果及其重纏繞相位圖。

可以看出基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)展開(kāi)相位圖平滑連續(xù)，其重纏繞相位圖干涉條紋與圖12(a)中的干涉條紋保持一致，且干涉條紋清晰，幾乎不存在相位噪聲，這表明基于無(wú)味信息濾波的相位展開(kāi)算法(UIFPU)不僅獲得了較好的展開(kāi)結(jié)果，且可有效地抑制干涉圖中的相位噪聲。