董小桐，韓春明，岳昔娟，趙迎輝

(1 中國科學(xué)院大學(xué)， 北京 100049； 2 中國科學(xué)院遙感與數(shù)字地球研究所， 北京 100094)(2018年5月4日收稿； 2018年6月1日收修改稿)

SAR影像上地物的位置信息是地形測(cè)繪、災(zāi)害監(jiān)測(cè)、資源普查、變化檢測(cè)等空間對(duì)地觀測(cè)技術(shù)應(yīng)用的基礎(chǔ)信息[1]。在攝影測(cè)量遙感領(lǐng)域，高精度的三維定位技術(shù)一直是研究的重點(diǎn)和難點(diǎn)[2]。機(jī)載重軌干涉合成孔徑雷達(dá)(repeat pass interferometric synthetic aperture radar, RP-InSAR)[3-6]，突破了雙天線系統(tǒng)受載機(jī)平臺(tái)尺寸的限制，實(shí)現(xiàn)不同時(shí)間重復(fù)觀測(cè)相同區(qū)域，在局部地區(qū)的地形測(cè)繪和形變監(jiān)測(cè)中具有不可替代的作用。機(jī)載重軌測(cè)量中，滿足多波段測(cè)繪要求，基線長(zhǎng)度可增至幾百米。然而基線長(zhǎng)度越大，獲取圖像像元高精度三維信息的難度越大。

機(jī)載重軌InSAR要求載機(jī)飛行軌跡高度重合，保持在與預(yù)設(shè)理想航跡相差±5 m的范圍內(nèi)，天線波束指向在±20°內(nèi)變化，精度維持在1°以內(nèi)[7]。干涉圖像對(duì)的高相干性保證DEM的高精度。重軌InSAR的飛行軌跡必須保持在空間基線附近的一個(gè)很小的圓柱范圍內(nèi)。不同于機(jī)載雙天線InSAR的剛性基線結(jié)構(gòu)，在現(xiàn)有的GPS系統(tǒng)和慣性導(dǎo)航系統(tǒng)條件下，機(jī)載重軌InSAR測(cè)量?jī)纱蔚娘w行軌跡，難以保證重合性。

機(jī)載重軌InSAR通?；€較長(zhǎng)，飛機(jī)飛行航跡嚴(yán)重偏離理想航跡，機(jī)載重軌InSAR圖像幾何畸變明顯，干涉圖像對(duì)之間相對(duì)形變明顯，同名像元對(duì)應(yīng)不同地面大小，造成干涉處理困難。基線越長(zhǎng)，干涉圖像對(duì)之間的相對(duì)形變?cè)矫黠@，尤其在遠(yuǎn)斜距端，失配越發(fā)嚴(yán)重，干涉圖像對(duì)的相干性越低，無法獲取高精度DEM，亦無法實(shí)現(xiàn)高精度三維定位。

多項(xiàng)式模型法[8-9]、共線方程法[10]和距離多普勒模型法(range-Doppler model，RD)[11]是主要的傳統(tǒng)機(jī)載SAR圖像定位方法。

嚴(yán)格意義上，多項(xiàng)式模型和共線方程模型，不能實(shí)現(xiàn)圖像的自動(dòng)實(shí)時(shí)定位，不存在明確的物理意義，并不符合SAR側(cè)視成像原理。多項(xiàng)式模型和共線方程模型在反演模型參數(shù)，實(shí)現(xiàn)定位校正時(shí)，嚴(yán)重依賴控制點(diǎn)信息[12]。

距離多普勒模型符合SAR側(cè)視成像原理，具有明確的幾何關(guān)系和物理意義[13]，因而目前在SAR影像定位中使用較多。

岳昔娟等[14]基于差分GPS和慣性導(dǎo)航系統(tǒng)數(shù)據(jù)，在無控制點(diǎn)的情況下，根據(jù)RD模型和DEM數(shù)據(jù)，導(dǎo)出一種進(jìn)行機(jī)載SAR影像主動(dòng)定位的數(shù)學(xué)模型。并試驗(yàn)驗(yàn)證此數(shù)學(xué)模型的正確性,分析主要系統(tǒng)誤差源及系統(tǒng)誤差的改正方法。吳穎丹[15]針對(duì)星載合成孔徑雷達(dá)對(duì)地目標(biāo)的定位問題,根據(jù)距離-多普勒模型和地球模型進(jìn)行對(duì)地定位的理論和不同計(jì)算方法，系統(tǒng)地分析RD模型現(xiàn)有的解算方法。劉佳音等[16]提出一種新的SAR圖像斜距多普勒定位模型的直接解法，推導(dǎo)出斜距多普勒定位模型的明確數(shù)學(xué)解析解，并與常用的數(shù)值迭代求解方法進(jìn)行比較分析，仿真結(jié)果驗(yàn)證了推導(dǎo)的正確性。

張紅敏等[17]針對(duì)無地面控制情況下的單幅SAR圖像定位問題,設(shè)計(jì)一種基于DEM和圖像仿真的單幅SAR圖像無控制定位方案。在外部DEM數(shù)據(jù)的支持下，將實(shí)際SAR圖像與仿真SAR圖像進(jìn)行匹配來提取控制點(diǎn)以完成定向參數(shù)解算。張紅敏等[18]還針對(duì)稀少控制下斜側(cè)視SAR圖像高精度定位難題，設(shè)計(jì)利用單個(gè)地面控制點(diǎn)的SAR圖像立體定位方案。邱春平等[19]使用斜距與航向之間的夾角表示RD模型中的多普勒方程，推導(dǎo)SAR距離-斜視角模型的線性化形式,設(shè)計(jì)相應(yīng)的定向參數(shù)解算和立體定位方案。

為了解決機(jī)載重軌InSAR三維定位問題，本文提出一種機(jī)載重軌InSAR高精度三維定位方法。研究和推導(dǎo)機(jī)載重軌InSAR定位模型，利用POS數(shù)據(jù)和少量地面控制點(diǎn)信息，從精準(zhǔn)SAR成像幾何關(guān)系出發(fā)，實(shí)現(xiàn)干涉圖像對(duì)精準(zhǔn)的相對(duì)定位，生成高精度DEM；并利用高精度DEM，絕對(duì)定位圖像像元，最終實(shí)現(xiàn)三維定位。

1 機(jī)載重軌InSAR形變分析

機(jī)載單景SAR圖像定位存在的主要問題是，載機(jī)飛行高度和姿態(tài)不穩(wěn)定產(chǎn)生的圖像畸變和雷達(dá)圖像近距壓縮的幾何特點(diǎn)。一般來說，為了保證機(jī)載SAR圖像同時(shí)實(shí)現(xiàn)方位向和距離向的高分辨率，載機(jī)飛行過程中應(yīng)盡可能保持勻速直線運(yùn)動(dòng)。而實(shí)際情況是，受氣流、GPS系統(tǒng)和慣性導(dǎo)航系統(tǒng)測(cè)量精度等的影響，載機(jī)在作非勻速飛行的同時(shí)會(huì)產(chǎn)生上下左右的偏轉(zhuǎn)和抖動(dòng)，極大地降低了成像質(zhì)量和定位精度。非勻速運(yùn)動(dòng)的航跡會(huì)造成方位向采樣不均勻；非直線運(yùn)動(dòng)的航跡會(huì)造成天線波束中心偏離，產(chǎn)生不規(guī)則幾何畸變。

機(jī)載重軌InSAR三維定位，不僅需要解決單幅SAR圖像的定位問題，還要解決干涉圖像對(duì)之間相對(duì)定位不準(zhǔn)確的問題。受POS系統(tǒng)實(shí)時(shí)精度的影響，重軌飛行航跡不完全重合且存在較大差異。造成自共同測(cè)繪帶內(nèi)測(cè)繪獲得的干涉復(fù)圖像對(duì)，距離向像元數(shù)不等。在長(zhǎng)基線條件下，從SAR圖像近距端到遠(yuǎn)距端，圖像對(duì)失配情況俞加明顯，難以保持相干性。

圖1顯示選自2014年10月在海南陵水地區(qū)進(jìn)行的重軌測(cè)繪中的4條載機(jī)飛行航跡，4條航跡均嚴(yán)重偏離預(yù)設(shè)航跡，差異明顯。圖中，航高和距離向坐標(biāo)間隔1 m；方位坐標(biāo)間隔0.16 m。

圖1 4條重軌航跡示意圖Fig.1 Drawing of four repeat paths

圖2中的散點(diǎn)表示統(tǒng)計(jì)得到的，一對(duì)重軌干涉SAR圖像對(duì)上同名像元的距離坐標(biāo)差值及其隨斜距的變化趨勢(shì)。隨著斜距的增長(zhǎng)，同名像元主輔圖像的距離坐標(biāo)差異逐漸增大，在近斜距端變化較快，在遠(yuǎn)斜距端變化漸緩。該現(xiàn)象說明在重軌機(jī)載InSAR中，存在干涉圖像對(duì)隨斜距增加，相干性逐漸降低的現(xiàn)象。

圖2 同名像元坐標(biāo)差值及其變化趨勢(shì)Fig.2 The coordinate differences and its trends of identical points

2 干涉幾何關(guān)系

假設(shè)基線長(zhǎng)度為B，基線傾角為α，r1為主天線斜距，r2為輔天線斜距，θ為斜視角。根據(jù)幾何關(guān)系和余弦定理可得

(1)

即

(2)

(3)

對(duì)于機(jī)載重軌干涉InSAR，Q=2。

因而，可以得出

(4)

P點(diǎn)的高度信息

h=H-r1cosθ,

(5)

水平距離

(6)

3 同名像元位置確定

3.1 RD模型

RD主動(dòng)定位模型包括:距離方程、多普勒方程和橢球方程。

(7)

式中：R=R0+m·j為目標(biāo)點(diǎn)的斜距；R0為初始斜距；m為距離分辨率，j為目標(biāo)點(diǎn)距離向坐標(biāo)；(Xs,Ys,Zs)、(Vx,Vy,Vz)分別表示天線相位中心坐標(biāo)和速度；fdol為多普勒頻率偏移；λ為波長(zhǎng)；(X,Y,Z)是目標(biāo)點(diǎn)地面坐標(biāo)；Re為地球橢球赤道半徑；Rp=(1-f)·(Re+h)，f為地球扁率；h為地面點(diǎn)高程。

3.2 定位參數(shù)改正

將RD主動(dòng)定位模型線性化，構(gòu)建誤差方程，代入地面控制點(diǎn)坐標(biāo)(X,Y,Z)，并根據(jù)最小二乘原理，求解dΧ：

dΧ=(BTB)-1(BTL).

(8)

式中，

(9)

3.3 主輔圖像配準(zhǔn)

將改正值代入新的初始斜距和航高，假設(shè)存在地面點(diǎn)P，則點(diǎn)P的主圖像地距為

(10)

式中：m1為圖上斜距分辨率；j1是點(diǎn)P的圖上斜距坐標(biāo)；H1new=H1+Δh1，H1為主圖像航高，Δh1為改正值；r1new=r1+Δr1，r1為主圖像初始斜距，Δr1為改正值。

航跡2相對(duì)航跡1的地距偏移為ΔY。點(diǎn)P在輔圖像上的斜距坐標(biāo)j2new為

(11)

式中：H2new=H2+Δh2，H2為輔圖像航高，Δh2為改正值；r2new=r2+Δr2，r2為輔圖像初始斜距，Δr2為改正值。

根據(jù)j2new，重采樣輔圖像，獲得像元一一對(duì)應(yīng)的干涉圖像對(duì)，完成機(jī)載重軌InSAR圖像配準(zhǔn)。

經(jīng)過該處理之后，InSAR圖像對(duì)之間的相對(duì)定位更加準(zhǔn)確。

4 機(jī)載重軌InSAR圖像定位

機(jī)載重軌InSAR圖像定位的流程如圖3所示，具體過程如下：

1)獲取模型參數(shù)，根據(jù)主圖像像元坐標(biāo)(i,j)，讀取其天線相位中心的瞬時(shí)位置向量(Xs,Ys,Zs)、速度向量(Vx,Vy,Vz)和斜距R。

2)計(jì)算地心直角坐標(biāo)系坐標(biāo)原點(diǎn)，從地圖頭文件中讀取SAR圖像中心的大地經(jīng)緯度(BC,LC)，并轉(zhuǎn)換坐標(biāo)系

(12)

式中：(X0,Y0,Z0)為地心直角坐標(biāo)系初始值，Have是圖像測(cè)繪區(qū)域的平均高程，N為地球曲率半徑。

3)主輔圖像配準(zhǔn)

RD模型解算并改正初始斜距和航高，計(jì)算主圖像像元在輔圖像上的距離向坐標(biāo)，使用sinc插值法對(duì)輔圖像進(jìn)行重采樣，獲得像元大小一致且一一對(duì)應(yīng)的干涉圖像對(duì)。

4)干涉處理

干涉處理分為相干計(jì)算、干涉濾波以及相位解纏繞等，最后獲得解纏繞的相位圖。

5)高程計(jì)算

根據(jù)干涉幾何關(guān)系，通過干涉相位計(jì)算主圖像像元對(duì)應(yīng)的地面高程。

6)坐標(biāo)轉(zhuǎn)換

地心直角坐標(biāo)系坐標(biāo)(X,Y,Z)向高斯平面直角坐標(biāo)系坐標(biāo)(x,y)轉(zhuǎn)換。

7)計(jì)算像元灰度和高程

對(duì)應(yīng)像元坐標(biāo)(i,j)的灰度值σ(i,j)和高程值h(i,j)賦值給坐標(biāo)轉(zhuǎn)換后的像元。

5 試驗(yàn)與分析

試驗(yàn)所用機(jī)載InSAR重軌圖像為2014年10獲取的，海南省陵水地區(qū)的測(cè)量結(jié)果。采用正側(cè)式成像，分辨率0.5 m，所用波段為C波段，極化方式為VV極化，由POS AV 610系統(tǒng)獲取POS數(shù)據(jù)。圖像基本參數(shù)如表1所示。選擇距離向和方位向均勻分布的9個(gè)D級(jí)控制點(diǎn)，通過RD模型解算，改正初始斜距和航高后，計(jì)算圖像像元在主輔圖像的距離向坐標(biāo)差值，即Δj=j2new-j1，如圖2中實(shí)線所示。

圖2中散點(diǎn)和實(shí)線在數(shù)值和變化趨勢(shì)上的一致性說明，使用新到的初始斜距和航高，可以準(zhǔn)確地定位主圖像像元在輔圖像上的位置，干涉圖像對(duì)之間的相對(duì)定位更加精準(zhǔn)。

圖3 機(jī)載重軌InSAR圖像定位的流程Fig.3 Flow chart of 3D location of airborne repeat pass InSAR

表1 基本參數(shù)Table 1 Basic parameters

為了確定本文方法對(duì)機(jī)載重軌InSAR圖像的處理效果。選取幾何畸變嚴(yán)重的遠(yuǎn)斜距區(qū)域進(jìn)行相關(guān)處理，結(jié)果如圖4、圖5所示。

圖4 圖像對(duì)相干性Fig.4 Coherence of image pairs

圖5 干涉處理結(jié)果Fig.5 Result of interfereometric processing

分別計(jì)算輔圖像重采樣之前和之后與主圖像的相關(guān)性。設(shè)定相關(guān)窗口4×4；搜索窗口16×16，搜索窗口數(shù)目1 000，均勻分布在整個(gè)圖像。評(píng)價(jià)參數(shù)是，全部搜索窗口內(nèi)相關(guān)系數(shù)極大值的分布和窗口內(nèi)部相關(guān)系數(shù)極大值點(diǎn)即匹配點(diǎn)的斜距坐標(biāo)偏移量。

如圖4所示，圖4(a)、4(b)分別表示原圖像對(duì)和新圖像對(duì)的相關(guān)系數(shù)極大值分布直方圖，4(c)、4(d)分別表示匹配的距離向坐標(biāo)差值。圖4(b)的值分布明顯高于圖4(a)，峰值出現(xiàn)在0.9附近，說明新圖像對(duì)的相干性整體更高；由圖4(c)顯示匹配的絕大多數(shù)分布于窗口邊緣，說明大部分干涉圖像對(duì)的同名點(diǎn)分布于搜索窗口范圍之外，而圖4(d)相對(duì)于圖4(c)變化更加明顯，斜距坐標(biāo)偏移量絕大多數(shù)分布于0值附近，說明干涉圖像對(duì)的同名點(diǎn)基本分布于搜索窗口中心，新圖像對(duì)的匹配性非常高，遠(yuǎn)超過原圖像對(duì)。

圖5為試驗(yàn)數(shù)據(jù)的干涉處理結(jié)果。

試驗(yàn)區(qū)域地表信息保存完整，圖像對(duì)相干性高，解纏繞相位圖連續(xù)性好。高程變化敏感度高，在高度不同的地物類型分布區(qū)域，具有十分明顯的高程變化。

上述試驗(yàn)結(jié)果表明，本文方法在實(shí)現(xiàn)InSAR復(fù)圖像對(duì)精準(zhǔn)相對(duì)定位的同時(shí)可以保持圖像對(duì)之間的高度相干性，特別適合相對(duì)幾何畸變隨斜距增長(zhǎng)而增大的長(zhǎng)基線機(jī)載重軌InSAR圖像對(duì)。使用檢查點(diǎn)評(píng)價(jià)像元的三維定位，實(shí)測(cè)數(shù)據(jù)和計(jì)算數(shù)據(jù)比對(duì)如表2所示。

表2 定位精度分析Table 2 Analysis of location accuracy

對(duì)比檢查點(diǎn)實(shí)測(cè)三維坐標(biāo)和本文方法計(jì)算得到的三維坐標(biāo)顯示，方位向定位誤差控制在0.32 m以內(nèi)；距離向定位誤差控制在0.8 m以內(nèi)；高程誤差控制在0.6 m之內(nèi)，符合高精度三維定位的要求。

6 結(jié) 論

本文提出一種利用地面控制點(diǎn)，標(biāo)定初始斜距和航高，獲得高精度DEM數(shù)據(jù)，最后實(shí)現(xiàn)機(jī)載重軌InSAR圖像三維定位的方法。通過重采樣輔圖像，獲得像元一一對(duì)應(yīng)的InSAR圖像對(duì)，使得圖像對(duì)之間的相對(duì)定位更加精確。相干性更加良好，干涉處理后獲得高精度的DEM。實(shí)驗(yàn)結(jié)果表明本文方法可以實(shí)現(xiàn)機(jī)載重軌InSAR圖像高精度三維定位。但是，本文試驗(yàn)數(shù)據(jù)是正側(cè)視模型下的成像數(shù)據(jù)，還需驗(yàn)證對(duì)斜側(cè)視模型下成像數(shù)據(jù)的處理效果。試驗(yàn)區(qū)域地形平坦，無需外源DEM輔助，在地形復(fù)雜地區(qū)需要引入外源DEM輔助處理[20]。另外，本文未考慮大氣效應(yīng)[21-22]和土壤含水量[23]的影響，且本文方法對(duì)于長(zhǎng)條帶的重軌機(jī)載InSAR圖像的有效性，還需進(jìn)一步試驗(yàn)。