一種高分辨率InSAR時(shí)序分析反演地物高程與地面沉降量的方法
【專利摘要】一種高分辨率InSAR時(shí)序分析反演地物高程與地面沉降量的方法，它有五大步驟：步驟一：高分辨率SAR數(shù)據(jù)選??；步驟二：初始高程相位估計(jì)；步驟三：選取InSAR時(shí)序分析干涉像對數(shù)據(jù)集并提取相干目標(biāo)候選點(diǎn)；步驟四：高程相位校正和線性形變分量估計(jì)；步驟五：相干目標(biāo)的形變序列估計(jì)。本發(fā)明能夠克服城區(qū)高分辨率InSAR形變監(jiān)測數(shù)據(jù)處理中無高分辨率DEM去除高層建筑物附加相位，有效地將城市高層建筑物在高分辨率InSAR中的附加相位與形變相位進(jìn)行分離，從而大大提高高分辨率InSAR在城區(qū)地面沉降中的監(jiān)測精度。
【專利說明】-種高分辨率InSAR時(shí)序分析反演地物高程與地面沉降量 的方法

【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及一種高分辨率InSAR時(shí)序分析反演地物高程與地面沉降量的方法，屬 于合成孔徑雷達(dá)干涉測量技術(shù)（InSAR)領(lǐng)域。它能夠有效地將城市高層建筑物在高分辨率 InSAR干涉圖中的附加相位與形變相位進(jìn)行分離，可以大大提高高分辨率InSAR在城區(qū)地 面沉降中的監(jiān)測精度。

【背景技術(shù)】
[0002] 利用InSAR求解高程（DEM)或者形變量，相位解纏是必須的步驟之一，這一步驟在 差分干涉處理中稱為解纏，而在時(shí)序分析方法中則稱為參數(shù)估計(jì)，其實(shí)質(zhì)是求解纏繞相位 的整周數(shù)，基本過程是求解相鄰點(diǎn)間的相位差，然后按照特定的路徑或網(wǎng)絡(luò)以一定的約束 條件進(jìn)行積分，求解觀測范圍整體的解纏相位，進(jìn)而反演形變場或高程。相位解纏的前提是 干涉圖連續(xù)分布且變化平緩，滿足相位差小于η的約束條件，整個(gè)相位場為無旋場，解纏 結(jié)果不隨路徑而變。而實(shí)際上干涉圖受噪聲或不連續(xù)相位（如不連續(xù)的陡坡）的影響，難 以滿足解纏條件，因而需要按照給定路徑求解后再進(jìn)行連接。高分辨率條件下城市高層建 筑物所引起的相位變化類似于非連續(xù)的陡坡相位。對于求解高程而言，條紋密度隨基線的 增大而密集，且存在與形變條紋混合的可能。因而，如何同時(shí)解算高程和形變相位，提高形 變量估計(jì)的精度是高分辨率InSAR應(yīng)用面臨的主要難題。
[0003] 地形相位補(bǔ)償是InSAR時(shí)序分析處理過程中差分相位求取的基本步驟，同時(shí)地物 高程也是InSAR時(shí)序分析的主要解算結(jié)果。中等分辨率SAR數(shù)據(jù)中表現(xiàn)為一個(gè)或者幾個(gè)點(diǎn) 目標(biāo)的像元在高分辨率雷達(dá)數(shù)據(jù)中以數(shù)個(gè)點(diǎn)目標(biāo)表現(xiàn)出來，每個(gè)散射體的散射中心不同， 如一棟建筑物，在中分SAR圖像中表現(xiàn)為一個(gè)點(diǎn)目標(biāo)，而在高分SAR圖像中以多個(gè)不同的散 射體表現(xiàn)出來。對高分辨率InSAR地形相位補(bǔ)償?shù)淖钪苯臃椒ㄊ谦@取每個(gè)散射中心高精度 的DEM數(shù)據(jù)，利用Lidar、TanDEM-X等獲取的DEM數(shù)據(jù)實(shí)現(xiàn)近似同等分辨率差分干涉處理的 高程相位補(bǔ)償。然而這種高分辨率的DEM數(shù)據(jù)很多情況下難以獲得。因此，這種通過外部 DEM(如SRTM)補(bǔ)償整個(gè)像元高程的方法在高分辨率數(shù)據(jù)處理中將難以適用，需要研究單個(gè) 目標(biāo)高程精確計(jì)算的方法。
[0004] 在InSAR干涉圖中，干涉條紋的頻率與基線的關(guān)系如圖1所不，同等1?程下，基線 越短，條紋越稀疏，利于相位解纏或參數(shù)估計(jì)。當(dāng)?shù)匦巫兓癁? π時(shí)對應(yīng)的高程差為：
[0005]

【權(quán)利要求】
1. 一種商分辨率InSAR時(shí)序分析反演地物商程與地面沉降量的方法，其特征在于：該 方法具體步驟如下： 步驟一：高分辨率SAR數(shù)據(jù)選取 以TerraSAR - X和COSMO - Skymed為代表的高分辨率雷達(dá)衛(wèi)星系統(tǒng)為開展InSAR 反演城區(qū)地物高程和形變精細(xì)化監(jiān)測提供了數(shù)據(jù)源，高分辨率InSAR相對于中等分辨率 InSAR技術(shù)，其總體優(yōu)勢體現(xiàn)在兩個(gè)方面，S卩（i )高密度相干點(diǎn)目標(biāo)和短周期4-16天； (ii )對地面點(diǎn)目標(biāo)的準(zhǔn)確定位；高分辨率InSAR數(shù)據(jù)處理的數(shù)據(jù)選取要滿足空間上盡可 能覆蓋完整城區(qū)的軌道數(shù)據(jù)，在時(shí)間上數(shù)據(jù)能夠連續(xù)接收； 步驟二：初始高程相位估計(jì) InSAR形變時(shí)序分析所構(gòu)建的二維參數(shù)估計(jì)模型中，考慮到大氣的空間相關(guān)性，對相鄰 兩點(diǎn)求差以削弱大氣的影響；相干目標(biāo)i和j差分干涉相位的互差為：
(3) 上式中，CB為與垂直基線相關(guān)的系數(shù)，T為時(shí)間基線，Λ ε為相對高程誤差，Λν為相 對形變速率，μ m為非線性形變量，α為大氣相位，η為噪聲，k表示干涉圖個(gè)數(shù)，與干涉圖 序列的組合有關(guān)；構(gòu)建目標(biāo)函數(shù)如下：
(,1) 將上式從相位互差式（3)中減去，得到殘余相位為：
(5) 顯然，當(dāng)目標(biāo)函數(shù)的參數(shù)△ ε和Λν在準(zhǔn)確估計(jì)時(shí)，殘余相位將最小化； 采取將差分干涉圖相位解纏的基礎(chǔ)上進(jìn)行的估計(jì)，這時(shí)式（3)轉(zhuǎn)換為二維線性函數(shù)， 通過建立Delanay三角網(wǎng)或利用冗余網(wǎng)構(gòu)建更為復(fù)雜的連接關(guān)系強(qiáng)化對待解算方程組的 約束，利用鄰近法則將所有距離滿足大氣相關(guān)距離的相干目標(biāo)連接起來，在求解完成相鄰 點(diǎn)間的互差后，通過最小二乘或加權(quán)平均的方法求解每個(gè)目標(biāo)相對于參考點(diǎn)的高程值和形 變速率場； 當(dāng)干涉圖時(shí)間間隔和空間基線足夠小時(shí)，變形相位可忽略，則可以將二維模型式（3) 轉(zhuǎn)為一維模型式（6)，進(jìn)而解算得到DEM ;
(6) (1)選取超短時(shí)間和空間基線干涉像對數(shù)據(jù)集 依據(jù)上述理論，根據(jù)建筑物的高度、干涉圖基線及地形相位條紋密度的關(guān)系，選取研究 區(qū)平均建筑物高度作為高程約束，確定計(jì)算原始地物高程時(shí)需要的干涉像對序列；選取時(shí) 間基線和空間基線都較短的干涉像對序列作為初始求解的數(shù)據(jù)集，設(shè)定初始時(shí)間和空間基 線閾值分別為1\和&，對數(shù)據(jù)集中的干涉像對進(jìn)行干涉處理，無需外部DEM進(jìn)行地形相位 校正，綜合利用快速傅里葉變換估計(jì)和多項(xiàng)式擬合估計(jì)的方法去除干涉圖中的軌道誤差和 趨勢性干涉條紋，生成超短時(shí)空基線干涉圖數(shù)據(jù)集、相干圖數(shù)據(jù)集和強(qiáng)度圖數(shù)據(jù)集；對每一 個(gè)干涉圖進(jìn)行相位解纏，按照一維模型多次迭代修正高程誤差相位估計(jì)得到原始DEM ;需 要說明的是，在地理編碼之前所有的操作都是在雷達(dá)坐標(biāo)系下進(jìn)行的； (2) 利用點(diǎn)目標(biāo)識別方法提取相干目標(biāo)候選點(diǎn) 針對上述所有生成的干涉圖序列，綜合采用幅度離散指數(shù)即Amplitude Dispersion Index和相干系數(shù)coherence來篩選相干目標(biāo)候選點(diǎn)； 幅度離散指數(shù)的計(jì)算公式為：
其中，σ A和％分別為像素幅度值的標(biāo)準(zhǔn)差和均值；給定一適當(dāng)閾值
DA低于閾值的 像元確定為相干目標(biāo)候選點(diǎn)； 雷達(dá)干涉相位圖的相干系數(shù)估計(jì)公式為：
(8) 根據(jù)各像元點(diǎn)在相干圖中的相干系數(shù)序列、和給定的相干系數(shù)閾值f，如果 mean
那么則將該像元確定為相干目標(biāo)候選點(diǎn)； (3) 多次迭代修正高程誤差相位估計(jì)地物高程 對于每一相干目標(biāo)點(diǎn)，利用超短時(shí)空基線干涉圖數(shù)據(jù)集并按照式（6)所述的一維模型 多次迭代修正高程誤差相位，將多次估計(jì)的高程誤差修正相位相加，得到相干目標(biāo)點(diǎn)處的 高程，然后，將點(diǎn)矢量轉(zhuǎn)換為柵格數(shù)據(jù)，生成研究區(qū)的原始DEM; 步驟三：選取InSAR時(shí)序分析干涉像對數(shù)據(jù)集并提取相干目標(biāo)候選點(diǎn) 將空間基線和時(shí)間基線不超過某一設(shè)定閾值的所有干涉像對進(jìn)行干涉處理，將步驟二 生成的原始高程模擬地形相位用于所有干涉圖的高程相位補(bǔ)償，生成差分干涉圖序列；針 對單一差分干涉圖中出現(xiàn)的軌道誤差和趨勢性干涉條紋，綜合利用快速傅里葉變換估計(jì)和 多項(xiàng)式擬合估計(jì)的方法予以去除，最終生成用于InSAR時(shí)序分析的差分干涉圖數(shù)據(jù)集和相 應(yīng)的相干圖數(shù)據(jù)集以及強(qiáng)度圖數(shù)據(jù)集；對上述所有生成的差分干涉圖序列，綜合采用幅度 離散指數(shù)和相干系數(shù)來篩選相干目標(biāo)候選點(diǎn)，降低高分辨率SAR條件下相干目標(biāo)數(shù)量冗 余； 步驟四：高程相位校正和線性形變分量估計(jì) 解纏 InSAR時(shí)序分析差分干涉圖序列中的每一幅差分干涉圖；選取時(shí)間和空間基線閾 值分別為T2和B2的干涉像對數(shù)據(jù)集，重新利用一維模型二次估計(jì)高程殘余值dhOTi ;將上述 一維模型所得高程殘余值dhOTi作為形變時(shí)序分析二維參數(shù)估計(jì)模型的初始高程，以Λ T、 ΛΒ和Adh作為時(shí)間、空間基線和高程修正的步長，將Τ2+Λ?^ΡΒ2+ΛΒ內(nèi)的干涉圖參與到 二維參數(shù)估計(jì)的序列中，設(shè)定初始最大高程修正閾值DH max OTi，并初次估計(jì)高程修正值和形 變速率；分別以步長ΛΤ和ΛΒ逐步增加用于迭代估算的干涉圖數(shù)量，并以步長Adh逐步 減小最大高程修正閾值，通過多次迭代處理逐步修正高程估計(jì)值和形變速率；直至滿足時(shí) 間和空間基線以及最大高程修正閾值的所有差分干涉圖都參與計(jì)算，估計(jì)最終的高程改正 和形變速率；由于逐次迭代求解得到的高程改正是對已經(jīng)去除高程值的估計(jì)，因而最終的 (9) 地物高程估計(jì)值為逐次迭代修正值之和，即： 式中Η為迭代終止時(shí)的地物高程估計(jì)值，u為總迭代次數(shù)； 步驟五：相干目標(biāo)的形變序列估計(jì) 對于具有顯著非線性形變過程，仍需對殘余相位進(jìn)行更為復(fù)雜的處理，以提取非線性 形變量；處理的前提仍是基于不同相位的時(shí)空頻率特性，從差分相位中去除高程和線性速 率后的殘余相位為：
(10) 由于已解算出整周相位，因而殘余相位為相位主值，其大小滿足：
(11) 這里，首先要實(shí)現(xiàn)對單個(gè)差分干涉圖中殘余相位的濾波處理；由于大氣表現(xiàn)出空間低 頻特性，而非線性變形的空間變化范圍較小，相對大氣相位而言表現(xiàn)為高通特性，因而，對 殘余相位圖進(jìn)行空間低通濾波進(jìn)一步弱化大氣的影響；在短基線集條件下求解的非線性 形變量，不同于PSInSAR的處理策略，不直接對大氣進(jìn)行估計(jì)，而是通過濾波減弱大氣的影 響；殘余相位經(jīng)過空間高通濾波后，大氣分量已經(jīng)減弱，此時(shí)，進(jìn)一步的處理是求解與雷達(dá) 數(shù)據(jù)對應(yīng)的不同時(shí)刻的非線性變形量，包含噪聲影響；根據(jù)主輔影像的關(guān)系，將相位解纏后 的殘余相位分解為：
(12) 式中，P和q表示生成第k張差分干涉圖的主輔影像的獲取時(shí)間；由于采用短基線原 貝1J，在求解每個(gè)時(shí)刻對應(yīng)的殘余相位時(shí)，會出現(xiàn)秩虧方程組，解決這一問題的辦法就是奇異 值分解法即SVD ;由此，在奇異值分解的基礎(chǔ)上，對Μ張殘余相位進(jìn)行時(shí)域低通濾波處理，以 提取最終的非線性形變量，表示為：
(13) 在求解出非線性形變量后，每個(gè)相干目標(biāo)對應(yīng)的形變序列為：
(14) 其中，vest為解算得到線性形變速率。
【文檔編號】G01S13/90GK104123464SQ201410353107
【公開日】2014年10月29日 申請日期:2014年7月23日 優(yōu)先權(quán)日:2014年7月23日
【發(fā)明者】葛大慶, 王艷, 劉斌, 張玲, 李曼, 郭小方 申請人:中國國土資源航空物探遙感中心