時序Sentinel-1TOPS模式SAR數據精配準

馬張烽，蔣 彌，丁琪瑄

1. 河海大學地球科學與工程學院，江蘇 南京 211100； 2. 中山大學測繪科學與技術學院，廣東 珠海 519082

Sentinel-1合成孔徑雷達(SAR)衛星具有重訪周期短(A、B組合為6 d)和覆蓋范圍廣等優點，近年來已被廣泛應用于成像大地測量和地球物理學研究。其中，干涉寬幅模式采用逐行掃描地形觀測(terrain observation with progressive scanning,TOPS)[1-2]作為標準SAR成像模式，可以同時獲得3個子條帶，幅寬達250 km。然而在TOPS模式下，相鄰burst影像(圖1)高多普勒質心頻率差(約5000 Hz)對方位向配準誤差十分敏感。理論上，TOPS模式數據方位向配準精度需小于0.001像素，才能避免干涉圖相位不連續或相位跳變現象。現有研究表明，采用精密軌道的幾何粗配準[3]和增強譜分集(enhanced spectral diversity,ESD)技術[4]是實現高精度TOPS模式數據配準，消除相位跳變的有效手段。

圖1 TOPS成像Fig.1 Description of TOPS imaging geometry

受地形起伏影響，傳統低階多項式擬合影像配準方法無法滿足配準偏移量為恒量的假設前提，導致非平穩信號的出現[5]。隨著GNSS定軌技術的成熟以及SAR影像聚焦理論的發展，聯合SAR精密軌道信息和外部數字高程模型(digital elevation model,DEM)進行粗配準是相對有效的方法[6-9]。假設忽略電離層效應和方位向地表形變的相位貢獻，粗配準后的方位向配準誤差主要與軌道誤差、SAR傳感器鐘差誤差相關[10]。這些誤差均產生系統性的方位向偏移量，需要進一步精配準操作消除殘余的方位向偏移量。

ESD技術是消除系統性方位向偏移量的有效方法。從已有研究來看，文獻[11]將ESD技術與譜分集[4]技術進行了比較，利用TerraSAR-X TOPS模式數據證明了ESD的優越性能。文獻[10，12]對Sentinel-1TOPS模式數據配準進行了詳盡分析，論證ESD技術具有普適性。在國內，文獻[6—7,13]詳細闡述了ESD技術的誤差敏感度，并定量證明了精密軌道對ESD相位估計的重要性。與分頻干涉法[14]類似，ESD利用上下頻帶(每個burst影像長度的8%)之間的重疊區干涉相位差來估計方位向殘余偏移量。因此，ESD的估計精度十分依賴干涉相干性[15]。就時序SAR配準而言，由于干涉對遭受的時空去相關作用不同，ESD的估計精度也不同。在低相干性場景，ESD技術可能完全歪曲估計結果并形成粗差，從而在數據堆棧中造成誤差傳播，降低整個數據集的配準精度。

最大化時間網絡相干性是改善ESD估計精度的有效手段。類似于InSAR時序分析，小基線技術仍然是降低時空去相關的主流方法。例如，文獻[10]提出了一種序貫時間網絡ESD方法(network ESD，NESD)，即通過時間去相關指數衰減模型將最鄰近時刻的影像相連，以便減少時間去相關對時序配準誤差的影響，最后對網絡進行整體平差。然而，這種方法無法適用于因環境變化導致的相干性突變(由于氣候條件的影響，如雪和降水)。文獻[15]提出的MST-ESD(minimum spanning tree based ESD)方法通過構造最小生成樹來克服這種突發的相干損失，但由于缺少冗余觀測，平差過程中解的穩定性可能依賴于研究區的整體相干性和其他估計手段。

綜上所述，目前的時序TOPS模式數據配準方法無法兼顧相干性和多余觀測。雖然小基線方法可以提供多余觀測，但是網絡中邊(干涉對)的數量是根據閾值試探性確定的，過大的閾值易增加低相干干涉對，而過小的閾值易形成隔離子集。為解決這一問題，本文在圖論框架下使用Bellman-Ford單源最短路徑算法最大化時間網絡的相干性并確保多余觀測。考慮到經驗性閾值選取的局限性，本文提出在先驗序貫網絡的基礎上構筑任意兩個TOPS模式影像的最短路徑，在包含N幅影像的堆棧中，從N(N-1)/2個干涉集合中挑選出相干性較高的差分干涉子集，并確保網絡的全連通性。通過模擬和實測數據證明，本文方法可顯著提高ESD估計精度，進而改善整個時序TOPS模式數據集配準精度。

1 增強譜分集技術

如圖1所示，連續burst之間獨立成像且存在8%重疊區域。因此在干涉處理時，將方位向上兩個相鄰burst的重疊區域分為前后視影像，即將上一個burst影像的重疊區域視為前視影像，將下一個burst影像的重疊區域視為后視影像。通過計算前后視影像干涉相位差，并結合估計的多普勒質心頻率差解算出burst的方位向殘余偏移量，進而實現影像精配準[16]。其中，連續burst之間的多普勒質心頻率差Δfovl≈5000 Hz，這使得式(1)中方位向配準時間誤差Δt的微小變化會引起干涉相位差(ESD相位)φESD的大幅變化

φESD=2πΔfovlΔt

(1)

該現象說明配準誤差對相位變化十分敏感，可以通過反算ESD相位來獲得高精度配準時間誤差。任一像元點P的ESD相位φESD,P可通過兩次差分得到

(2)

式中，mi,p和si,p分別表示第i個主輔burst單視復數影像；mi+1,p和si+1,p分別表示第i+1個主輔burst單視復數影像；上標*表示復共軛算子。基于Sentinel-1高精度軌道姿態控制，ESD相位通常小于一個周期的相位變化，所以像元點P的ESD相位φESD,p無須進行相位解纏，而是通過周期圖進行估算[12]

(3)

(4)

假設兩個burst干涉圖的重疊區域具有相同去相干機制，則Δx的精度可以通過相位噪聲標準差σφ來衡量[10,15,17]

(5)

σΔx和γ的關系如下

(6)

圖2展示了配準精度σΔx與樣本數n和干涉相干性γ的關系。可以看出，當重疊區樣本數n固定時，相干性將明顯影響TOPS模式數據配準精度。

圖2 ESD在40萬、400萬和4000萬個相干像元點下的配準偏移量理論精度Fig.2 Theoretical azimuth co-registration accuracy of ESD under four hundred thousand,four million and forty million coherent pixels respectively

2 時序TOPS模式數據精配準算法

2.1 Bellman-Ford單源最短路徑算法

最短路徑是圖論研究的經典問題，旨在尋找圖中兩節點間的最短路徑并求得該路徑下各邊的最小權值總和。針對這類問題，本文采用Bellman-Ford單源最短路徑算法。其基本思想是運用迭代法對所有邊進行松弛操作，并由此更新起始節點到其他節點的最短路徑值，直到獲得最優解時迭代停止[18-22]。

在最短路徑問題中，常用無向圖G=(V,E,w)表示不含負回路的網絡，其中V={v1,v2,…,vN}為節點集(影像序列)，E為G的邊集(干涉組合)，w為G中的權矩陣(與相干性有關)。Bellman-Ford算法通過松弛操作推測V中所有節點最短路徑的層次，并逐層生成以起始節點v1為根的最短路徑樹。松弛操作過程中，邊集E的每次遍歷是基于前一次遍歷對當前遍歷的部分節點的最短路徑值進行更新。在第k次遍歷中，最短路徑層次為1～k-1的節點均已確定，此時最短路徑層次為k的節點由最短路徑層次為k-1的節點獲得。由于最短路徑不包含回路，(v1,vj)最短路徑層次最高為N-1，因此算法最多循環N-1次。

Bellman-Ford算法的具體步驟如下：

(1) 對起始節點v1與節點vj的距離進行初始化。令起始節點v1到自身的距離d(v1)=0，起始節點v1到其他節點vj的距離d(vj)=+∞。此時(v1,vj)最短路徑值d1(j)恰好等于起始節點v1與節點vj相連的單邊距離

d1(j)=w1j(1≤j≤N)

(7)

(2)進入循環。循環下標為k，節點vj在第k次循環中到起始節點v1的距離為dk(j)。在循環內部，遍歷邊集E的每條邊并進行松弛計算

(8)

(3) 當所有節點均滿足dk(j)=dk+1(j)(k≤N)時，dk(j)即為(v1,vj)的最短路徑值，此時循環結束。由式(8)可以看出，距離的度量實際上通過權w來實現，而權w就是與相關性有關的量。

2.2 InSAR時間網絡

用上述Bellman-Ford方法構建InSAR時間網絡的方法由以下3個步驟組成：①確定起點和終點;②從邊集合中挑選出每個起點和終點的最短路徑;③刪除重復的邊并形成網絡。

在步驟①中，起點和終點可以通過小基線集合事先確定。本文采用頂點度數為4的序貫網絡中所有邊的起點和終點(vi,vj)作為先驗(圖3)。假設幾何配準后的N景SAR影像按照時間升序排列V={v1,v2,…,vN}。步驟②從帶權函數w:E→R的加權無向圖G=(V,E)開始，找到從起點vi到最終vj之間最小權重的路徑。E是所有可能邊的集合,即N(N-1)/2個干涉圖的平均相干性。根據Bellman-Ford算法，可以從E中挑選出經過(vi,vj)相干性最高的干涉組合。在步驟③中，重復步驟②直至所有起點vi和終點vj被遍歷。在刪除重復的干涉對之后，獲得優化的干涉子集。圖3概括了時間網絡優化的具體過程。可以看出，增強后的時序網絡能夠使低質量的影像相互連接，因而比小基線中的經驗閾值更加可靠，避免后者可能產生的孤立子集。與最小生成樹方法相比，使用Bellman-Ford算法可以提供更多的多余觀測或形成閉合環，這保證了在最小二乘平差過程中獲得更高的估計精度。

圖3 時間網絡優化流程Fig.3 Flow chart of temporal network optimization

2.3 BF-ESD算法

將ESD技術與上述時間網絡聯合，可得到時序TOPS模式數據集配準算法，主要步驟為：

(1) 根據公共主影像的SAR幾何，采用外部DEM和精密軌道信息將N景SAR數據進行幾何配準。

(2) 用基于第2類統計的相干性估計量估計burst影像重疊區的樣本相干性[21]，獲得N(N-1)/2景干涉圖的平均相干性。

(3) 用Bellman-Ford算法建立兼顧相干性最大化和多余觀測的時間網絡。

(4) 在使用式(2)—式(4)估計時間網絡上每個干涉對的方位向偏移量之后，用加權最小二乘估計器獲得時序數據集關于公共主影像的偏移量。

在后續的試驗章節中本文的方法簡寫為BF-ESD算法。

3 結果與討論

本文通過模擬數據和真實數據對BF-ESD的性能進行了評估和精度驗證，并與現有的NESD和MST-ESD兩種時序配準方法進行對比分析。

3.1 模擬試驗

3.1.1 時序TOPS模式數據模擬

首先采用蒙特卡羅模擬試驗評估時間網絡對方位向精配準的影響。與前述研究類似，在復圓高斯模型的假設下，本文通過復圓高斯矩陣和復協方差矩陣的哈達馬內積來模擬SLC時間序列[12,23-25]。在SLC模擬過程中，采用指數衰減去相關模型模擬復協方差矩陣中的非對角元素

(9)

式中，γi,j表示第i與j個SLC影像之間的相干系數。γ0=0.5和γ∞=0.1分別表示短期和長期相干性。γ0的物理意義是當時間基線ti,j趨于0時相干性的初值，而γ∞表示隨著時間基線趨于無窮大，相干性衰減到的最小值。數值上，γ0=0.5表征了低相干場景，而γ∞=0.1模擬了相干性估計的偏差，即在純噪聲的環境下的相干性估計值并不為0。t=50 d表示與相干性衰減速度有關的量。為了凸顯真實場景中突發的相干性損失，本文對協方差矩陣的非對角元素添加了[0,1]區間的隨機信號。

式(2)需要兩個burst干涉圖的差分，因此在試驗中用相同的協方差矩陣和不同的復高斯隨機數模擬兩組N=30的時間序列，并在差分干涉圖中隨機添加區間在[-0.02,0.02]像素之間的相位貢獻。時序影像的獲取間隔定義為12 d，其中假設第1景影像為主影像。為滿足ESD技術的假設條件，線性形變速率設置為0。在每次模擬試驗中，用于估算協方差的樣本設置為100，從而減少樣本復協方差矩陣估計對網絡解算的干擾，強調時間網絡的重要性。為公平比較，除了時間網絡的差異，3種配準方法的相位和相干性估計量均相同，所使用的估計參數也相同。

3.1.2 分析與討論

圖4是10 000次蒙特卡羅的試驗結果。將各配準方法估計的偏移量的方差除以基于單主影像配準的偏移量的方差之后轉換成分貝(dB)，獲得相對方差。相對方差值越小，說明配準精度越高。可以看出，本文提出的BF-ESD方法由于最大化相干性并包含多余觀測，得到的方差最小。相比之下，MST-ESD中的最小生成樹的估計結果最差。結果是可以預料的，因為在低相干場景中缺少多余觀測，無法獲得穩健的平差結果。NESD中的序貫網絡雖然增加了多余觀測，但包含了相對低質量的干涉對，因而精度低于BF-ESD。BF-ESD在所有時間序列數據中都有最小方差，說明該方法具有普適性，可獲得更高的配準精度。

圖4 10 000次蒙特卡洛模擬試驗的配準偏移量相對方差Fig.4 The relative variance of azimuth offset after 10 000 times Monte Carlo simulation

3.2 真實數據驗證

3.2.1 Sentinel-1 TOPS模式數據處理

選用2019年1月至2020年1月的31景Sentinel-1 TOPS模式影像(track84,frame110)，并選用覆蓋了關中平原大部分區域的IW2和IW3兩個子帶。覆蓋區屬于溫帶季風氣候區和亞熱帶季風氣候區的過渡地帶，年平均降水量約561 mm，平均濕度為69.6%，植被較茂盛，地形包含平原峽谷，具有部分相干性條件優良的裸地。選取該研究區的主要依據是氣候與地貌條件使得不同干涉對存在干性差異，進而在執行圖論算法時能夠挑選出優質干涉組合、最大化時間網絡的相干性。相反，相干性過高(如沙漠、城市)或完全去相關的區域(如中國西南部地區)均不能體現出算法的優勢。這是因為時間網絡上的干涉組合相干性對比度不強，各路徑的最小權值接近，因此不能從本質上改善網絡質量。

在幾何配準之后，將所有SLC影像進行干涉組合(圖5)，包括NESD中的序貫網絡(n=4，圖5(a))、MST-ESD中的最小生成樹(圖5(b))，以及本文提出的Bellman-Ford網絡(圖5(c))。對圖5中每個干涉對，采用式(2)—式(4)估計ESD相位和方位向偏移量，并用加權最小二乘平差獲得關于公共主影像(第1景)的偏移量。與模擬試驗相同，除時間網絡之外，所有方法均使用了相同的相干性估計器、相位估計器以及估計參數(7×7規則窗)。

圖5 NSED、MST-ESD和BF-ESD的時序網絡Fig.5 Temporal network of NESD,MST-ESD and BF-ESD

為定量評定各算法的精度，首先利用各時間網絡估計30個輔影像的方位向偏移量。然后對其校正并重采樣獲得配準序列。最后將配準序列再次進行ESD估計，獲得方位向偏移量。運用此法，對每一種網絡進行50次重復試驗，每次試驗都基于前一次配準序列。最終獲得的偏移量標準偏差越小，說明算法越穩定，精度越高。

3.2.2 分析與討論

由圖5可以看出，本次試驗研究區的相干性低，平均相干性約為0.2。序貫網絡包含的干涉對最多，導致網狀結構復雜且包含相干性低的邊。最小生成樹挑選了所有干涉對中相干性最高的子集，然而受限于高相干點數量，最小生成樹不能形成閉合的網絡。相比之下，Bellman-Ford網絡不僅剔除了序貫網絡中相干性低的邊，同時保證了網絡的閉合性，這在2019年4月至2019年8月的子集中更為顯著。

圖6展示了BF-ESD方法在方位向偏移量糾正前后的對照。可以看出，圖6(a)中左下角的相位不連續在圖6(b)BF-ESD算法糾正之后被消除，說明BF-ESD算法能夠有效去除方位向偏移量誤差。圖7展示了各網絡下時間序列的方位向偏移量的標準偏差。從中可以看出，多余觀測仍然是保證精度的主要前提，這可以從基于最小生成樹的偏移量標準偏差中得以證明。其次，干涉對的質量仍然會影響偏移量估計的精度。在確保多余觀測并改善觀測質量的情況下，BF-ESD在所有配準的圖像中均能得到最小的不確定性。定量地，基于最小生成樹網絡的偏移量標準偏差為0.001 9像素，NESD為0.001 3像素，BF-ESD為0.000 8像素。相比之下，BF-ESD方法估計精度遠高于NESD方法以及MST-ESD方法，證明了本文方法的改進性能。

圖6 配準誤差糾正前后的干涉圖對比(2019-04-24—2019-07-05，干涉圖經過濾波處理以供視覺判讀)Fig.6 Interferograms of 2019-04-24—2019-07-05 before and after correction of co-registration error using BF-ESD (both interferograms have been filtered for visual interpretation)

圖7 時序數據集配準殘余偏移量標準偏差Fig.7 Standard deviation of residual mis-registration

4 結論與展望

針對目前時序TOPS模式數據配準因去相干問題導致的精度衰減，本文提出了一種改善時間網絡質量的方法。該方法利用圖論中Bellman-Ford算法反復修正先驗網絡中低相干性的干涉組合，進而獲得相干性更高的時間網絡，并保留多余觀測。高質量的時間網絡能夠減少因ESD估計誤差導致的最小二乘誤差傳播，因而獲得更精確的時間配準偏移量。本文結合模擬數據和真實數據對算法進行了可靠性評估并與現有時序配準算法的時間網絡進行了對比分析。結果表明，相比NESD中的序貫網絡以及MST-ESD中的最小生成樹網絡，利用高質量網絡的BF-ESD算法具有更高的估計精度，不確定性相應減少。

雖然本文方法僅從時序TOPS模式數據配準給出解決方案，但是基于Bellman-Ford的時間網絡能夠在不改變算法的條件下延伸至小基線集的InSAR時序分析，從而改善用經驗性閾值確定干涉子集的缺點。

致謝：Sentinel-1數據由ESA/Copernicus提供。