閆 滿 劉國(guó)林 郝華東 陶秋香

(1)山東科技大學(xué)測(cè)繪科學(xué)與工程學(xué)院,青島 266510

2)黑龍江科技學(xué)院計(jì)算機(jī)與信息工程學(xué)院,哈爾濱 150027)

二維元胞自動(dòng)機(jī)相位解纏算法的改進(jìn)＊

閆 滿1,2)劉國(guó)林1)郝華東1)陶秋香1)

(1)山東科技大學(xué)測(cè)繪科學(xué)與工程學(xué)院,青島 266510

2)黑龍江科技學(xué)院計(jì)算機(jī)與信息工程學(xué)院,哈爾濱 150027)

提出一種改進(jìn)的二維元胞自動(dòng)機(jī)相位解纏算法。通過將原有元胞自動(dòng)機(jī)相位解纏算法中的 4鄰域迭代規(guī)則改進(jìn)為最近鄰域規(guī)則,以干涉圖的邊緣作為始點(diǎn)進(jìn)行迭代,逐一減少干涉條紋的數(shù)目,達(dá)到相位解纏的目的。用模擬和真實(shí)的干涉圖進(jìn)行實(shí)驗(yàn),結(jié)果表明,該算法減少了原有元胞自動(dòng)機(jī)相位解纏算法的迭代次數(shù),降低了計(jì)算機(jī)內(nèi)存的占用率,提高了相位解纏的精度。

元胞自動(dòng)機(jī);相位解纏;最近鄰域規(guī)則;4鄰域規(guī)則;干涉圖

1 引言

合成孔徑雷達(dá)干涉測(cè)量(InSAR)利用雷達(dá)干涉數(shù)據(jù)的相位信息獲取高精度數(shù)字地形信息。其基本原理是通過兩部天線或兩次重復(fù)軌道對(duì)地面的同一地區(qū)進(jìn)行觀測(cè),獲取兩幅復(fù)雷達(dá)圖像,利用兩幅復(fù)圖像相干運(yùn)算得到的相位差反演地物的相對(duì)高程。但是實(shí)際上這樣得到的干涉相位差是真實(shí)相位差的卷疊,即 InSAR干涉圖中與地面位置直接相關(guān)的相位是以 2為模的,所以為了計(jì)算每一點(diǎn)的高程必須給每一個(gè)相位測(cè)量值加上整數(shù)倍的相位周期,這種求解 2模糊性問題的技術(shù)稱為相位解纏[1]。

傳統(tǒng)的二維相位解纏算法有 Goldstein枝切法[2]、區(qū)域增長(zhǎng)法[3]等,它們是一種局部算法,其優(yōu)點(diǎn)是可以隔絕相位不連續(xù)點(diǎn),阻止局部相位誤差在整個(gè)積分區(qū)域的傳播,計(jì)算速度較快,在相干性較好的區(qū)域可以獲得精確的解纏相位,但是在強(qiáng)噪聲條件下,很難獲得最佳積分路徑,容易造成誤差傳遞或無法解纏的孤立區(qū)域。除此之外還有最小LP范數(shù)法,其中應(yīng)用較為廣泛的是無權(quán)重的最小二乘法[4,5]、加權(quán)的最小二乘法[6,7],它們都是一種全局算法,其優(yōu)點(diǎn)是運(yùn)算穩(wěn)定性好,不需要識(shí)別殘差點(diǎn)。元胞自動(dòng)機(jī)算法也屬于全局算法,1987年 Ghigial等人提出了將元胞自動(dòng)機(jī)算法用于相位解纏。其基本思想是通過檢驗(yàn)內(nèi)部 4鄰域相位值的不一致性產(chǎn)生的規(guī)則來進(jìn)行迭代,在每一次迭代后相位的變化范圍增大,向相位真值的方向趨近,經(jīng)過多次迭代,最后得到連續(xù)的相位。但是該算法也存在如下問題：在總體迭代中,需要存儲(chǔ)振蕩前后的狀態(tài)值,這使得計(jì)算機(jī)內(nèi)存的占用率高;且每一次總體迭代中又包含局部迭代,隨著干涉條紋的增加,局部迭代次數(shù)嚴(yán)重影響解纏的速度。為了解決上述問題,本文將原有的 4鄰域規(guī)則改進(jìn)為最近鄰域規(guī)則,同時(shí)將每一點(diǎn)設(shè)置為零權(quán)值進(jìn)行相位解纏,最后分別采用模擬和真實(shí)的干涉圖進(jìn)行實(shí)驗(yàn),來驗(yàn)證改進(jìn)算法的有效性。

2 二維元胞自動(dòng)機(jī)相位解纏算法

二維元胞自動(dòng)機(jī)相位解纏是一種檢驗(yàn) 4鄰域相位值的不一致而逐步解纏和修正的過程,并且要求在無方向性的條件下進(jìn)行并行計(jì)算。具體二維元胞自動(dòng)機(jī)相位解纏算法如下[8]：

1)以相位圖中的一點(diǎn)為基準(zhǔn),分別求出當(dāng)前點(diǎn)和與它正交 4鄰域點(diǎn)的相位差;

2)相鄰兩點(diǎn)的權(quán)值定義為兩點(diǎn)之間展開相位所需 2π的整數(shù)倍,將每一點(diǎn)對(duì)應(yīng)的 4個(gè)權(quán)值相加;

3)如果權(quán)值之和為正,則該中心點(diǎn)的相位值加上 2π,為負(fù)則減去 2π;

4)如果 4個(gè)權(quán)值都為零,則該點(diǎn)和 4鄰域點(diǎn)的相位差都不超過π,該點(diǎn)值不變。如果正的權(quán)值和負(fù)的權(quán)值剛好抵消,則該點(diǎn)加上一個(gè) 2π;

5)步驟 1)到 4)的一次重復(fù)稱為一次局部迭代。經(jīng)過一次局部迭代,干涉圖中每一點(diǎn)的相位值同時(shí)變化到了下一時(shí)刻的狀態(tài)。重復(fù)步驟 1)到4),直到出現(xiàn)周期為 2的振蕩,即相鄰兩個(gè)周期的相位圖完全相同。這時(shí)算法陷入周期性的重復(fù),將相鄰的兩個(gè)相位圖相加求平均。兩個(gè)振蕩狀態(tài)的平均稱為總體迭代。然后判斷相位是否已完全解纏,如果是則結(jié)束,否則轉(zhuǎn)到步驟 1),重新開始局部迭代。

該算法從相位圖的邊緣開始進(jìn)行解纏,為了確定是否達(dá)到周期為 2的振蕩狀態(tài),需要將前后兩個(gè)時(shí)刻的狀態(tài)值都存儲(chǔ)在計(jì)算機(jī)中,這樣大大占用了計(jì)算機(jī)的內(nèi)存。而且在進(jìn)行相位解纏時(shí),每一次總體迭代中包含若干次局部迭代,隨著干涉條紋的增加,局部迭代次數(shù)的增加速度難以控制。為了解決上述問題,本文提出一種改進(jìn)的元胞自動(dòng)機(jī)算法。

3 改進(jìn)的二維元胞自動(dòng)機(jī)相位解纏算法

在改進(jìn)的算法中將每一點(diǎn)設(shè)置為零權(quán)值,并且將判斷 4鄰域點(diǎn)相位差的規(guī)則改進(jìn)為最近鄰域規(guī)則：

其中：a1∈{-2,-1,1,2},a2∈{-2,-1,1,2};b1∈{-2,-1,0,1,2},b2∈{-2,-1,0,1,2},b3∈{-2,-1,0,1,2};c1∈{-2,-1,0,1,2},c2∈{-2,-1,0,1,2},c3∈{-2,-1,0,1,2};d1∈{-2,-1,1,2},d2∈{-2,-1,1,2}。方程 (2)和 (3)被稱為最近鄰域規(guī)則。

對(duì)文獻(xiàn)[8]的元胞自動(dòng)機(jī)相位解纏算法,利用最近鄰域規(guī)則進(jìn)行改進(jìn),得改進(jìn)的算法如下：

1)以相位圖中的每一點(diǎn)為基準(zhǔn),分別求出當(dāng)前點(diǎn)和與它相鄰點(diǎn)的相位差;

2)如果相位差滿足方程 (2)或 (3),即在相鄰相位之間存在相位跳變,則的值為 1,也就是在該點(diǎn)值上加上一個(gè) 2π,否則為 0,該點(diǎn)值不變;

3)按照方程(1)計(jì)算出每一點(diǎn)的相位值;

4)步驟 1)到 3)的一次重復(fù)稱為一次迭代。重復(fù)步驟 1)到 3),直到所有的干涉條紋均消失為止。

4 實(shí)驗(yàn)結(jié)果與分析

利用改進(jìn)算法進(jìn)行實(shí)驗(yàn),分別給出 3組實(shí)驗(yàn)結(jié)果：第一組以無噪聲的模擬干涉圖作為實(shí)驗(yàn)對(duì)象,第二組以有噪聲的模擬干涉圖作為實(shí)驗(yàn)對(duì)象,第三組以真實(shí)的干涉圖作為實(shí)驗(yàn)對(duì)象。

4.1 無噪聲干涉圖相位解纏

采用 200×200的無噪聲干涉圖 (圖 1)。利用改進(jìn)的元胞自動(dòng)機(jī)算法,將干涉圖中的每一點(diǎn)利用方程(1)進(jìn)行計(jì)算,如果滿足方程(2)或(3),則加上2π,一直到所有的纏繞相位都消失在干涉圖的中心區(qū)域?yàn)橹埂２捎迷械亩S元胞自動(dòng)機(jī)相位解纏算法和改進(jìn)算法得到了同樣的解纏結(jié)果(圖 2)。改進(jìn)算法經(jīng)歷 80次迭代就能得到比較理想的結(jié)果,而原有的元胞自動(dòng)機(jī)方法需經(jīng)歷 400次迭代才能得到相同的結(jié)果。

圖1 無噪聲干涉圖Fig.1 Interferogram without noise

圖2 解纏結(jié)果Fig.2 Result of phase unwrapping

4.2 有噪聲干涉圖相位解纏

在原有的干涉圖中,加上相干系數(shù)為 0.8的噪聲,產(chǎn)生的干涉圖如圖 3所示,選擇干涉圖中左上角的點(diǎn)作為相位解纏的起始點(diǎn),利用式 (1)進(jìn)行計(jì)算,得到的解纏結(jié)果如圖 4所示。由于噪聲的存在,改進(jìn)算法經(jīng)歷了 85次迭代得到了解纏結(jié)果。而原有的二維元胞自動(dòng)機(jī)相位解纏算法需要 420次迭代才能得到同樣的結(jié)果?？梢钥闯?在有噪聲情況下,本文算法能較快地進(jìn)行相位解纏,而原有的元胞自動(dòng)機(jī)算法解纏速度較慢,甚至在噪聲較多時(shí),不能正確進(jìn)行解纏。

圖3 有噪聲干涉圖Fig.3 Interferogram with noise

圖4 解纏結(jié)果Fig.4 Result of phase unwrapping

4.3 真實(shí)干涉圖相位解纏

采用以 2003-12-03和 2004-01-07的兩景伊朗Bam地區(qū) ENV ISAT-ASAR圖像為實(shí)驗(yàn)的主副圖像。

對(duì)完成配準(zhǔn)的主副圖像做干涉處理,經(jīng)過去除平地效應(yīng)和濾波等處理,得到干涉圖像,其部分干涉相位如圖 5所示。選擇從干涉圖的右下方開始進(jìn)行相位解纏,圖 6、7、8、9分圖的右下方逐漸向左上方移動(dòng),并且隨著迭代次數(shù)的增加,干涉條紋隨之減少,直到經(jīng)過 216次迭代后所有的干涉條紋全部消失,即相位解纏過程結(jié)束。而原有的元胞自動(dòng)機(jī)算法,由于需要確定相鄰兩個(gè)振蕩狀態(tài)值是否一致,使計(jì)算機(jī)內(nèi)存占有率增加,從而使解纏速度明顯降低。另外原有的元胞自動(dòng)機(jī)算法在總體迭代中包括多次局部迭代,使得迭代次數(shù)遠(yuǎn)遠(yuǎn)超過改進(jìn)算法,而改進(jìn)算法的迭代次數(shù)不會(huì)超過干涉圖的大小(不會(huì)超過250)。改進(jìn)算法在干涉圖存在不連續(xù)的情況下,也能得到較好的解纏結(jié)果。

圖5 真實(shí)干涉圖Fig.5 Real interfergram

圖6 第8次迭代結(jié)果Fig.6 Result after 8th iteration

圖7 第64次迭代結(jié)果Fig.7 Result after 64th iteration

圖8 第104次迭代結(jié)果Fig.8 Result after 104th iteration

圖9 第160次迭代結(jié)果Fig.9 Result after 160th iteration

圖10 解纏結(jié)果Fig.10 Phase unwrapping result

5 結(jié)論

通過對(duì)原有的元胞自動(dòng)機(jī)相位解纏算法中的 4鄰域迭代規(guī)則改進(jìn)為最近鄰域規(guī)則,分別對(duì)模擬的無噪聲干涉圖,有噪聲干涉圖和真實(shí)的干涉圖進(jìn)行相位解纏,實(shí)驗(yàn)結(jié)果說明改進(jìn)算法能夠得到令人滿意的結(jié)果。由于不必保存振蕩前后的狀態(tài)值,計(jì)算機(jī)的內(nèi)存被極大地降低了,同時(shí),改進(jìn)算法不必計(jì)算相鄰兩個(gè)振蕩狀態(tài)的平均值,在總體迭代中沒有局部迭代,提高了相位解纏的效率。即使在有噪聲的情況下,改進(jìn)算法的迭代次數(shù)也明顯少于原有算法的迭代次數(shù),甚至在干涉圖中有不連續(xù)的相位時(shí),也能得到令人滿意的結(jié)果。

1 王超,張紅,劉智.星載合成孔徑雷達(dá)干涉測(cè)量[M].北京：科學(xué)出版社,2002.(Wang Chao,Zhang Hong and Liu Zhi.Spaceborne synthetic aperture radar interferometry[M]. Beijing：Science Press,2002)

2 Goldstein R M,Zerber H A and Werner C L.Satellite radar interferometry：Two-dimensional phase unwrapping[J].Radio Science,1988,23(4)：713-720.

3 XuW and Cumming I.A region-growing algorithm for inSAR phase unwrapping[J].IEEE Trans Geosci Remote Sens., 1999,37(3)：124-133.

4 PrittM D and Shipman J S.Least-squares two-dimensional phase unwrapping using FFTs[J].IEEE Trans Geosci Remote Sens.,1994,32(3)：706-708.

5 Pritt M D.Phase unwrapping by means of multigrid techniques for interferometric SAR[J].IEEE Trans Geosci Remote Sens.,1996,34(3)：728-738.

6 Flynn T J.Two dimensional phase unwrapping with minimum weighted discontinuity[J].J Opt Soc Am.,1997, (14)：2 692-2 701.

7 Ghiglia D C and Romero L A.Robust two-dimensional weighted and unweighted phase unwrapping that uses fast transforms and iterative methods[J].J Opt Soc Am., 1994,A11(1)：107-117.

8 Ghiglia D C,Mastin GA and Romero L A.Cellular-automata method for phase unwrapping[J].J Opt Soc Am.,1987,4 (1)：267-80.

IM PROVEM ENT OF PHASE UNW RAPPING M ETHOD FOR TWO-D IM ENSI ONAL CELLULAR AUTOMATA

YanMan1,2),Liu Guolin1),Hao Huadong1)and Tao Qiuxiang1)

(1)Geom atics College,Shandong University of Science and Technology,Q ingdao 266510 2)Computer and Infor mation Engineering College,Heilongjiang Institute of Science and Technology,Haerbin 150027)

An improved phase unwrappingmethod for cellular automata is proposed.It revises the four neighborhood of the originalmethod to the nearest neighborhood,begins from the boundary of the interferogram to iterate,eliminates the phase fringes one by one to achieve the purpose of phase unwrapping.The experi ments by using the si mulated and real interferograms aremade.The result shows that thismethod can decrease the original iterative times,reduce the occupancy rate of the memory and improve the accuracy of the phase unwrapping.

cellular automata;phase unwrapping;nearest neighborhood rule;four neighborhood rule;interferogram

1671-5942(2010)04-0142-05

2010-01-23

國(guó)家自然科學(xué)基金(40874001);國(guó)家 863計(jì)劃項(xiàng)目(2009AA12Z147);山東省泰山學(xué)者建設(shè)工程專項(xiàng)(TSXZ-0520);黑龍江科技學(xué)院 2006年引進(jìn)人才科研啟動(dòng)基金(06115)

閆滿,女,1979年生,吉林長(zhǎng)春,博士生,講師,研究方向：微波遙感技術(shù)及應(yīng)用.E-mail：yanmansdust@126.com

P225.1

A