自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法
【專利摘要】本發(fā)明公開(kāi)了一種自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，該方法包括將干涉相位變換到復(fù)數(shù)域并對(duì)復(fù)數(shù)域的實(shí)部和虛部做非局部濾波，對(duì)得到的復(fù)干涉相位的實(shí)部和虛部計(jì)算，得到干涉相位，選出殘差點(diǎn)最少的濾波后的干涉相位并記錄迭代次數(shù)參數(shù)m＝1；增長(zhǎng)搜索窗口、鄰域窗口尺寸，m＝m+1，將殘差點(diǎn)最少的第m-1次濾波后干涉相位作為待濾波圖像進(jìn)行非局部濾波，選出殘差點(diǎn)最少的第m次濾波后干涉相位，獲取最佳鄰域窗口值；判斷是否繼續(xù)迭代，若繼續(xù)迭代，則增長(zhǎng)搜索窗口、鄰域窗口尺寸；若停止迭代，則輸出濾波后的干涉相位。
【專利說(shuō)明】自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法
【技術(shù)領(lǐng)域】
[0001]本發(fā)明涉及干涉合成孔徑雷達(dá)信號(hào)處理領(lǐng)域，尤其涉及一種自適應(yīng)迭代的非局部干涉相位濾波方法。
【背景技術(shù)】
[0002]干涉合成孔徑雷達(dá)(InterferometricSynthetic Aperture Radar, InSAR)利用兩個(gè)通道雷達(dá)回波的干涉相位信息提取地形高程、運(yùn)動(dòng)目標(biāo)速度、地表形變等各種信息，將合成孔徑雷達(dá)(SAR)的測(cè)量拓展到三維空間。其具有全天時(shí)、全天候、高精度的特點(diǎn)，因此在地形測(cè)繪、冰川研究、海洋測(cè)繪以及地面沉降監(jiān)測(cè)等多個(gè)領(lǐng)域都有廣泛的應(yīng)用。
[0003]干涉SAR在獲取相位信息的過(guò)程中，干涉相位會(huì)受到熱噪聲去相干、基線或幾何去相干、時(shí)間去相干等隨機(jī)誤差的影響，具體表現(xiàn)為相位噪聲。低質(zhì)量的干涉條紋將會(huì)影響后續(xù)的干涉相位展開(kāi)及DEM反演的準(zhǔn)確性，
[0004]最終給高程測(cè)量帶來(lái)誤差，影響測(cè)繪精度。因此，為獲取高質(zhì)量的干涉相位，必須對(duì)干涉條紋進(jìn)行濾波處理。
[0005]針對(duì)該問(wèn)題，前人提出了多種相位濾波方法，這些算法可以分為兩類:空域?yàn)V波算法，如圓周期均值濾波(Eichel P H, Ghiglia D C, Jakowatz C VJr, etal..SpotlightSAR interferometry for terrain elevation mapping and interferometric changedetection.Sandia National Labs Tech Report, Washingdon DC,1993, US-D0E.)、Lee濾波(Lee J S, Papathanassiou P, Ainsworth T L, etal..A new technique for noisefiltering of SAR interferometric phase images.1EEE Trans.0n Geoscience andRemote Sensing, 1998, 36 (5):1456-1465.)等，以及變換域?yàn)V波算法，如 Goldstein 濾波(Goldstein R M and Werner C L.Radar interferogram filtering for geophysicalapplications.Geophysical Research Letters, 1998, 25 (21):4035-4038.)、小波變換濾波(Martinez C L and Fabregas X.Modeling and reduction ofSAR mterferometricphase noise in the wavelet domain.1EEE Trans, on Geoscience and Remote Sensing,2002,40(12):2553-2566.)等。其中，變換域?yàn)V波算法在實(shí)際中應(yīng)用更為廣泛，其基于信號(hào)與噪聲在變換域中所處的位置不同這一假設(shè)來(lái)進(jìn)行噪聲抑制，在很多情況下可以得到滿意結(jié)果。然而，在相干性較差或地形變化劇烈的情況下，由于細(xì)節(jié)信息與噪聲信息在變換域中所處位置無(wú)法完全區(qū)分，因此上述算法不能同時(shí)滿足去除噪聲和保持相位細(xì)節(jié)兩方面的要求。為得到高質(zhì)量的干涉相位，有必要研究一種適用于低相干及地形變化劇烈區(qū)域的相位濾波方法。

【發(fā)明內(nèi)容】

[0006](一 )要解決的技術(shù)問(wèn)題
[0007]為解決上述問(wèn)題，本發(fā)明提供了一種自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，以克服傳統(tǒng)的基于像素的處理方法不能兼顧去噪和細(xì)節(jié)保持的缺點(diǎn)，從而提聞干涉相位的精度。
[0008]( 二)技術(shù)方案
[0009]根據(jù)本發(fā)明的一個(gè)方面，提供一種自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，包括如下步驟:
[0010]步驟S1:將干涉相位Φ變換到復(fù)數(shù)域；
[0011]步驟S2:設(shè)定初始搜索窗口尺寸參數(shù)S,初始鄰域窗口參數(shù)N,對(duì)復(fù)數(shù)域的實(shí)部進(jìn)行非局部濾波處理，得到多個(gè)濾波后復(fù)干涉相位的實(shí)部& ;
[0012]步驟S3:對(duì)復(fù)數(shù)域的虛部進(jìn)行非局部濾波處理，得到多個(gè)濾波后的復(fù)干涉相位的虛部Im ；
[0013]步驟S4:對(duì)多個(gè)濾波后的復(fù)干涉相位的實(shí)部和虛部進(jìn)行計(jì)算，得到多個(gè)濾波后的干涉相位為fzarctam.1rn/lic)，從多個(gè)濾波的干涉相位的結(jié)果中選出殘差點(diǎn)最少的濾波后的干涉相位，同時(shí)得到該殘差點(diǎn)最少的第I次去除噪聲的濾波后的干涉相位對(duì)應(yīng)的鄰域窗口尺寸參數(shù)值η,并記錄迭代次數(shù)參數(shù)m = I ;
[0014]步驟S5:設(shè)定搜索窗口尺寸參數(shù)增長(zhǎng)步長(zhǎng)為r，鄰域窗口尺寸增長(zhǎng)步長(zhǎng)為t，令S=S+r,N = (n,n+t,n+2t) ,m = m+1,將得到的殘差點(diǎn)最少的第m-Ι次去除噪聲的濾波后的干涉相位作為待濾波的圖像，計(jì)算待濾波的圖像的殘差點(diǎn)數(shù)目RES1，并對(duì)待濾波的圖像進(jìn)行非局部濾波處理，并選出殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位，確定殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位對(duì)應(yīng)的鄰域窗口參數(shù)值NO ；
[0015]步驟S6:考察鄰域窗口參數(shù)值N = n+k大于殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位所對(duì)應(yīng)的鄰域窗口參數(shù)值NO = n+2t，k > 2t且為偶數(shù)時(shí)的第m次去除噪聲的濾波后的干涉相位，獲取最佳鄰域窗口值；如果殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位所對(duì)應(yīng)的鄰域窗口參數(shù)值NO古n+2t,則記錄該殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位的殘差點(diǎn)數(shù)目RES2，本次迭代結(jié)束，轉(zhuǎn)到步驟S7 ；
[0016]步驟S7:判斷是否繼續(xù)迭代，若繼續(xù)迭代，則轉(zhuǎn)到步驟S5 ;若停止迭代，則輸出濾波后的干涉相位。
[0017](三)有益效果
[0018]本發(fā)明基于非局部去噪的思想，利用圖像的冗余信息進(jìn)行濾波，根據(jù)干涉相位的自身特性，利用干涉相位噪聲標(biāo)準(zhǔn)差對(duì)權(quán)值計(jì)算中的衰減系數(shù)進(jìn)行估計(jì)，實(shí)現(xiàn)了衰減系數(shù)的自適應(yīng)，之后通過(guò)調(diào)節(jié)窗口尺寸及自動(dòng)迭代，在濾除噪聲的同時(shí)能夠更好的利用相位的紋理信息，使得干涉相位細(xì)節(jié)結(jié)構(gòu)不被破壞，從而克服了傳統(tǒng)的基于像素的濾波方法不能兼顧去噪和細(xì)節(jié)保持的缺點(diǎn)，實(shí)現(xiàn)高精度的干涉相位濾波，為高精度的干涉測(cè)量提供了條件，為干涉合成孔徑雷達(dá)干涉相位濾波提供了一種有效的方法。
【專利附圖】

【附圖說(shuō)明】
[0019]圖1為本發(fā)明實(shí)施例干涉合成孔徑雷達(dá)自適應(yīng)迭代的非局部相位濾波的流程圖。
[0020]圖2為本發(fā)明干涉合成孔徑雷達(dá)自適應(yīng)迭代的非局部相位濾波的實(shí)施例流程圖。
[0021]圖3為意大利Etna火山的實(shí)測(cè)干涉合成孔徑雷達(dá)干涉相位圖。
[0022]圖4a-圖4d為幾種濾波方法對(duì)圖3的干涉相位進(jìn)行濾波后的結(jié)果?！揪唧w實(shí)施方式】
[0023]為使本發(fā)明的目的、技術(shù)方案和優(yōu)點(diǎn)更加清楚明白，以下結(jié)合具體實(shí)施例，并參照附圖，對(duì)本發(fā)明進(jìn)一步詳細(xì)說(shuō)明。
[0024]需要說(shuō)明的是，在附圖或說(shuō)明書描述中，相似或相同的部分都使用相同的圖號(hào)。且在附圖中，以簡(jiǎn)化或是方便標(biāo)示。再者，附圖中未繪示或描述的實(shí)現(xiàn)方式，為所屬【技術(shù)領(lǐng)域】中普通技術(shù)人員所知的形式。另外，雖然本文可提供包含特定值的參數(shù)的示范，但應(yīng)了解，參數(shù)無(wú)需確切等于相應(yīng)的值，而是可在可接受的誤差容限或設(shè)計(jì)約束內(nèi)近似于相應(yīng)的值。
[0025]如圖1示出本發(fā)明提供的一種自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法包括如下步驟:
[0026]步驟S1:將干涉相位Φ變換到復(fù)數(shù)域；[0027]步驟S2:設(shè)定初始搜索窗口尺寸參數(shù)S,初始鄰域窗口參數(shù)N,對(duì)復(fù)數(shù)域的實(shí)部進(jìn)行非局部濾波處理，得到多個(gè)濾波后復(fù)干涉相位的實(shí)部Kc ;
[0028]步驟S3:對(duì)復(fù)數(shù)域的虛部進(jìn)行非局部濾波處理，得到多個(gè)濾波后的復(fù)干涉相位的虛部
[0029]步驟S4:對(duì)多個(gè)濾波后的復(fù)干涉相位的實(shí)部和虛部進(jìn)行計(jì)算，得到多個(gè)濾波后的干涉相位為^ = arcunK Iiii / Re )，從多個(gè)濾波的干涉相位的結(jié)果中選出殘差點(diǎn)最少的濾波后的干涉相位，同時(shí)得到該殘差點(diǎn)最少的第1次去除噪聲的濾波后的干涉相位對(duì)應(yīng)的鄰域窗口尺寸參數(shù)值η,并記錄迭代次數(shù)參數(shù)m = I ;
[0030]步驟S5:設(shè)定搜索窗口尺寸參數(shù)增長(zhǎng)步長(zhǎng)為r，鄰域窗口尺寸增長(zhǎng)步長(zhǎng)為t，令S=S+r,N = (n,n+t,n+2t) ,m = m+1,將得到的殘差點(diǎn)最少的第m-Ι次去除噪聲的濾波后的干涉相位作為待濾波的圖像，計(jì)算待濾波的圖像的殘差點(diǎn)數(shù)目RES1，并對(duì)待濾波的圖像進(jìn)行非局部濾波處理，并選出殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位，確定殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位對(duì)應(yīng)的鄰域窗口參數(shù)值NO ；
[0031]步驟S6:考察鄰域窗口參數(shù)值N = n+k大于殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位所對(duì)應(yīng)的鄰域窗口參數(shù)值NO = n+2t，k > 2t且為偶數(shù)時(shí)的第m次去除噪聲的濾波后的干涉相位，獲取最佳鄰域窗口值；如果殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位所對(duì)應(yīng)的鄰域窗口參數(shù)值NO古n+2t,則記錄該殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位的殘差點(diǎn)數(shù)目RES2，本次迭代結(jié)束，轉(zhuǎn)到步驟S7 ；
[0032]步驟S7:判斷是否繼續(xù)迭代，若繼續(xù)迭代，則轉(zhuǎn)到步驟S5 ;若停止迭代，則輸出濾波后的干涉相位。
[0033]圖2為本發(fā)明自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法的實(shí)施例流程圖。本發(fā)明中的參數(shù)選取可選擇但不限于該實(shí)施例中的參數(shù)選取情況。如圖2所示，本實(shí)施例包括以下步驟:
[0034]步驟S1:將干涉相位Φ變換到復(fù)數(shù)域= cos<i)+jsin<i)中，分別取復(fù)數(shù)域的實(shí)部Re {e.^} = cos Φ、虛部InUe.^} = sin Φ ； j為虛數(shù)，Φ為濾波前的干涉相位；
[0035]步驟S2:設(shè)定窗口尺寸參數(shù)(S，N) = (1，1)，(1，2)，(1，3)，其中，所述搜索窗口尺寸參數(shù)S對(duì)應(yīng)搜索窗口大小(2S+1) X (2S+1)，且參數(shù)S的初始值為小于4的自然數(shù)；所述鄰域窗口尺寸參數(shù)N對(duì)應(yīng)鄰域窗口大小(2N+1) X (2N+1)，且參數(shù)N取值為多個(gè)逐漸增大的自然數(shù)。[0036]對(duì)復(fù)數(shù)域的實(shí)部進(jìn)行非局部濾波處理，得到三個(gè)濾波后復(fù)干涉相位的實(shí)部Se;
[0037]上述對(duì)復(fù)數(shù)域的實(shí)部進(jìn)行非局部濾波處理，得到濾波后的坐標(biāo)i處的像素值NL ( υi)表示為:
【權(quán)利要求】
1.一種自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，包括如下步驟: 步驟S1:將干涉相位Φ變換到復(fù)數(shù)域； 步驟S2:設(shè)定初始搜索窗口尺寸參數(shù)S，初始鄰域窗口參數(shù)N，對(duì)復(fù)數(shù)域的實(shí)部進(jìn)行非局部濾波處理，得到多個(gè)濾波后復(fù)干涉相位的實(shí)哥Re 步驟S3:對(duì)復(fù)數(shù)域的虛部進(jìn)行非局部濾波處理，得到多個(gè)濾波后的復(fù)干涉相位的虛部Im ； 步驟S4:對(duì)多個(gè)濾波后的復(fù)干涉相位的實(shí)部和虛部進(jìn)行計(jì)算，得到多個(gè)濾波后的干涉相位為^ a_n(lm/liC)，從多個(gè)濾波的干涉相位的結(jié)果中選出殘差點(diǎn)最少的濾波后的干涉相位，同時(shí)得到該殘差點(diǎn)最少的第1次去除噪聲的濾波后的干涉相位對(duì)應(yīng)的鄰域窗口尺寸參數(shù)值η,并記錄迭代次數(shù)參數(shù)m = I ; 步驟S5:設(shè)定搜索窗口尺寸參數(shù)增長(zhǎng)步長(zhǎng)為r，鄰域窗口尺寸增長(zhǎng)步長(zhǎng)為t，令S =S+r, N = (n, n+t, n+2t), m = m+1,將得到的殘差點(diǎn)最少的第m_l次去除噪聲的濾波后的干涉相位作為待濾波的圖像，計(jì)算待濾波的圖像的殘差點(diǎn)數(shù)目RES1，并對(duì)待濾波的圖像進(jìn)行非局部濾波處理，并選出殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位，確定殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位對(duì)應(yīng)的鄰域窗口參數(shù)值NO ； 步驟S6:考察鄰域窗口參數(shù)值N = n+k大于步驟S5得到的殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位所對(duì)應(yīng)的鄰域窗口參數(shù)值NO = n+2t, k > 2t且為偶數(shù)時(shí)的第m次去除噪聲的濾波后的干涉相位，獲取最佳鄰域窗口值；如果步驟S5得到的殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位所對(duì)應(yīng)的鄰域窗口參數(shù)值NO古n+2t，則記錄該殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位的殘差點(diǎn)數(shù)目RES2，本次迭代結(jié)束，轉(zhuǎn)到步驟S7 ； 步驟S7:判斷是否繼續(xù)迭代，若繼續(xù)迭代，則轉(zhuǎn)到步驟S5 ;若停止迭代，則輸出濾波后的干涉相位。
2.根據(jù)權(quán)利要求1所述的自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，其特征在于，所述復(fù)數(shù)域= cos+jsin<i)中，分別取復(fù)數(shù)域的實(shí)部Re {e.^} = cos Φ、虛部InUe.^} = sin(J) ；j為虛數(shù)，Φ為濾波前的干涉相位。
3.根據(jù)權(quán)利要求1所述的自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，其特征在于，所述搜索窗口尺寸參數(shù)S，其對(duì)應(yīng)的搜索窗口大小為(2S+1) X (2S+1),且參數(shù)S的初始值為小于4的自然數(shù)；鄰域窗口尺寸參數(shù)N，其對(duì)應(yīng)鄰域窗口大小(2Ν+1) X (2Ν+1)，且參數(shù)N取值為多個(gè)逐漸增大的自然數(shù)。
4.根據(jù)權(quán)利要求1所述的自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，其特征在于，所述對(duì)復(fù)數(shù)域的實(shí)部進(jìn)行非局部濾波處理，得到濾波后的坐標(biāo)i處的像素值NL ( υ⑴)表示為:
5.根據(jù)權(quán)利要求1所述的自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，其特征在于，所述根據(jù)待濾波的干涉相位的相位噪聲及鄰域窗口大小設(shè)置非局部濾波處理中的衰減系數(shù)h的值的步驟為: 在條紋較稀疏即干涉相位每變化2 π所經(jīng)歷的像素在30個(gè)像素以上的情況下，衰減系數(shù)h值在5到幾百的范圍內(nèi)得到的濾波結(jié)果變化不大，衰減系數(shù)h值取值為經(jīng)典值10 σ n ;而當(dāng)條紋較密集即干涉相位每變化2 π所經(jīng)歷的像素在30個(gè)像素以下時(shí)，衰減系數(shù)h值的最優(yōu)值與噪聲標(biāo)準(zhǔn)差及鄰域窗口大小近似呈線性，且衰減系數(shù)h值在最優(yōu)值附近一個(gè)±1到±5的范圍內(nèi)得到的濾波結(jié)果變化不大，即與理想干涉相位的均方誤差值變化在0.1以內(nèi)，即具有一定的魯棒性，為采用擬合近似提供了條件，通過(guò)多項(xiàng)式擬合，最優(yōu)衰減系數(shù)h值表示為h = 1.3ση+0.5.(2Ν+1)+C3，其中，ση表示噪聲標(biāo)準(zhǔn)差，C3為常數(shù)項(xiàng)，C3取值在-0.5到-1之間。
6.根據(jù)權(quán)利要求1所述的自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，其特征在于，所述當(dāng)NO = n+2t時(shí)，繼續(xù)考察N > NO時(shí)的第m次去除噪聲的濾波后的干涉相位以獲取最佳鄰域窗口值的具體操作為:計(jì)算N = n+k，k > 2t且為偶數(shù)時(shí)的第m次去除噪聲的濾波后的干涉相位，其中，k表示鄰域窗口邊長(zhǎng)增加了 2k個(gè)像素，更新NO值及殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位，直至殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位對(duì)應(yīng)的鄰域窗口大小不是目前的最大值即滿足NO古n+k時(shí)，或者隨鄰域窗口增大，濾波效果基本不變即滿足N = n+k，n+(k-2)，n+(k-4)時(shí)的殘差點(diǎn)數(shù)目不變時(shí)，記錄該NO值對(duì)應(yīng)的第m次去除噪聲的濾波后的干涉相位的殘差點(diǎn)數(shù)目RES2，本次迭代結(jié)束，轉(zhuǎn)到步驟S7。
7.根據(jù)權(quán)利要求1所述的自適應(yīng)迭代的非局部干涉合成孔徑雷達(dá)干涉相位濾波方法，其特征在于，所述判斷是否繼續(xù)迭代的具體操作為:若4 < S < 10且待濾波圖像與殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位的殘差點(diǎn)數(shù)目RES1, RES2滿足:(RES1-RES2) / RESl < C，其中，C為常數(shù)，且O <C< 0.5，停止迭代，輸出殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位作為最終結(jié)果；若4 < S < 10且待濾波圖像與殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位的殘差點(diǎn)數(shù)目RES1，RES2滿足:(RES1-RES2) /RESl > C，則繼續(xù)迭代，η = NO，轉(zhuǎn)到步驟S5 ;若S≤10，則停止迭代，輸出殘差點(diǎn)最少的第m次去除噪聲的濾波后的干涉相位作為最終結(jié)果。
【文檔編號(hào)】G01S13/90GK103823219SQ201410095231
【公開(kāi)日】2014年5月28日 申請(qǐng)日期:2014年3月14日 優(yōu)先權(quán)日:2014年3月14日
【發(fā)明者】林雪, 李曾璽, 李芳芳, 胡東輝, 丁赤飚 申請(qǐng)人:中國(guó)科學(xué)院電子學(xué)研究所