張桂芳 郭利民 單新建 屈春燕

(中國地震局地質(zhì)研究所，地震動力學(xué)國家重點實驗室，北京 100029)

0 引言

雷達遙感中的角反射器是指能夠?qū)走_波束產(chǎn)生強烈反射，并且使其能夠逆向反射的地物。雷達角反射器可分為兩大類:一是天然角反射器，即永久散射體(PS點)，一般指靠向雷達波束入射方向的裸露山峰、山脊、巨石等天然地物，還包括尺寸不統(tǒng)一和形狀不規(guī)則的人工地物，如靠近水面的金屬塔、路面、堤岸等;二是人工角反射器(CR點)，特指那些利用導(dǎo)電性能和導(dǎo)磁性能良好，電容率大的金屬材料做成的一種點狀人工地物目標(biāo)。角反射器最初被用于雷達系統(tǒng)的外部定標(biāo)(Ulander，1991;Mohr et al.，2001)。

CR-InSAR技術(shù)是對危險地段特定結(jié)構(gòu)(活動構(gòu)造、大壩、橋梁和核電站等)形變監(jiān)測的有效手段，在進行構(gòu)造微形變監(jiān)測時，安裝在活動斷層兩側(cè)的角反射器可補充構(gòu)成高精度控制點網(wǎng)絡(luò)(張景發(fā)等，2006)。利用CR-InSAR技術(shù)對形變進行監(jiān)測的文獻很多:Xie等(2002，2004，2008)最早利用人工角反射器方法對滑坡進行監(jiān)測;程滔等(2007)提出了CR和PS兩種高相干散射點聯(lián)合解算算法;Renaud等(2010)用2臺CR監(jiān)測冰川移動速率;邢學(xué)敏等(2011)利用CR-InSAR與PS-InSAR聯(lián)合算法探測區(qū)域線性沉降;許才軍等(2012)利用CR-InSAR技術(shù)進行鮮水河斷層的形變速率監(jiān)測。盡管國內(nèi)外學(xué)者開展了CR技術(shù)用于形變監(jiān)測的研究，但都側(cè)重于論述CR點形變獲取方法，而對CR點像素級精確定位提及的甚少。在能查閱到的幾篇相關(guān)文獻中(Sakurai-Amano et al.，1999;范景輝，2008;Renaud et al.，2010)，均利用搜索最大幅度值定位角反射器像素。但由于角反射器散射幅度受背景影響大，加之聚焦性差，會使復(fù)雜背景條件下的角反射器出現(xiàn)定位偏差。

角反射器微形變研究需要解決的第1個難題是角反射器的精確定位。如果角反射器定位錯誤，那么基于角反射器的形變研究除了得到錯誤結(jié)論外，就沒有任何意義了。本文將通過分析角反射器影像特征，確定角反射器精確定位的最佳參數(shù)。

1 數(shù)據(jù)介紹及其預(yù)處理

研究區(qū)位于北京西北部的延懷盆地及其周邊地區(qū)，圖1中三角形給出的是研究區(qū)內(nèi)架設(shè)的三角形三面角角反射器，角反射器(特指三角形三面角角反射器，下同)實物照片見圖2。文中用到2種SAR數(shù)據(jù):一種數(shù)據(jù)是c波段的ENVISAT ASAR的升軌數(shù)據(jù)，波長為5.6cm，獲取時間是2008年12月12日，其覆蓋范圍見圖1中黑色實線，有18臺角反射器被其覆蓋;另一種數(shù)據(jù)是x波段的TerraSAR升軌數(shù)據(jù)，波長是3.1cm，獲取時間是2012年12月6日，其覆蓋范圍見圖1中黑色虛線，有12臺角反射器被其覆蓋;有7臺角反射器同時被上述2種數(shù)據(jù)覆蓋(表1)。

圖1 Envisat數(shù)據(jù)(黑色實線方框)和TerraSAR數(shù)據(jù)(黑色虛線方框)覆蓋范圍以及角反射器分布圖(白色三角)Fig.1 The coverage area of Envisat and TerraSAR data，and CRs distribution map.

圖2 研究區(qū)內(nèi)角反射器照片F(xiàn)ig.2 Photograph of CR.

SAR系統(tǒng)能以圖像的形式記錄地物后向散射電磁波的幅度和相位信息，根據(jù)雷達波與地物所產(chǎn)生的不同散射機制，可以將地面目標(biāo)分為兩大類:點目標(biāo)和分布目標(biāo)。點目標(biāo)是在像素分辨單元內(nèi)占主導(dǎo)地位的散射體，在進行距離向和方位向多視處理時，其光譜能量基本保持不變?？紤]點目標(biāo)這一特性，基于SLC(Single Look Complex)圖像的頻譜差異(spectral diversity)性，設(shè)置距離向和方位向視數(shù)為4，用3×3的像素窗口統(tǒng)計生成強度圖、光譜相關(guān)圖和光譜離散指數(shù)圖。

表1 角反射器初始定位像素Table 1 The initial pixels of CRs

2 角反射器初始定位

利用GPS定標(biāo)得到的實測角反射器三維信息一般是以WGS84橢球體為基準(zhǔn)，用大地坐標(biāo)系下的緯度、經(jīng)度、高程來表示，即(Φ，Λ，h)，而ENVISAT衛(wèi)星的狀態(tài)矢量是在以WGS84橢球體為基礎(chǔ)的地心旋轉(zhuǎn)(Earth Centered Rotating，ECR)坐標(biāo)系中進行描述的，也稱地心直角坐標(biāo)系(范景輝，2008)。要在影像上對角反射器定位，首先需要進行坐標(biāo)轉(zhuǎn)換，也就是大地坐標(biāo)向ECR坐標(biāo)轉(zhuǎn)換，即P(Φ，Λ，h)→P(X，Y，Z)。在DORIS精確軌道進行校正后的SLC參數(shù)文件中提供了有限的衛(wèi)星位置數(shù)據(jù)，為了獲得SAR影像中任意一行對應(yīng)的衛(wèi)星位置矢量，需對精軌數(shù)據(jù)進行插值，得到衛(wèi)星位置矢量和速度矢量隨成像時間t的函數(shù)。對于已處理至零多普勒成像參數(shù)的SLC數(shù)據(jù)，一行像元對應(yīng)著同一成像時間tl，在與SLC對應(yīng)的參數(shù)文件中，給出了SLC數(shù)據(jù)成像的開始時刻t0，以及每行數(shù)據(jù)對應(yīng)的時間間隔Δtl。

作為一種利用了合成孔徑技術(shù)的雷達系統(tǒng)，SAR不僅可以通過記錄地物反射信號的時間信息求得傳感器與地物目標(biāo)間的距離，還能記錄下地物反射信號的多普勒信息?；赟AR數(shù)據(jù)的這些特點，研究者發(fā)展了距離-多普勒模型(R-D模型)(Curlander，1982)，用以描述衛(wèi)星狀態(tài)矢量、地面目標(biāo)狀態(tài)矢量、成像參數(shù)以及參考橢球體之間的關(guān)系。利用R-D模型，采用牛頓迭代法(李慶揚等，2000;范景輝，2008)可以求解出角反射器瞬時成像時間t。則角反射器所在像素的行l(wèi)和列s可分別用下面的公式來計算:

式(2)中:ceil(l)表示對l向上取整;r_s表示距離向采樣數(shù)。表1中分別列出了研究區(qū)內(nèi)角反射器在Envisat數(shù)據(jù)和TerraSAR數(shù)據(jù)中初始像素位置。

3 角反射器輻射特性研究

3.1 RCS定義

雷達橫截面是描述目標(biāo)回波強度的物理量，簡稱RCS。RCS的正式定義為(Merrill，2008):

式(3)中:R是目標(biāo)到雷達的距離;E0是目標(biāo)接收到的入射波電場強度;Es是目標(biāo)散射波的電場強度。

簡單角反射器RCS的近似值可表示為

式(4)中:Aeff是貢獻多重內(nèi)部反射的有效面積，可以通過對每個反射面在視線方向的投影面積求和來確定(Kbittm，1981);λ是入射波波長。

由式(4)可知:角反射器RCS主要與其有效面積和接收的入射波波長有關(guān)。除了用m2為單位反映雷達散射截面外，另一種更通用的方法是用雷達散射截面對數(shù)值的10倍來表示，符號是σdBsm，單位是分貝平方米(dBsm)，即σdBsm=10lgσ。

3.2 角反射器影像特征分析

SAR以圖像的形式記錄角反射器后向散射電磁波的幅度和相位信息，將幅度值定標(biāo)并平方后就可得角反射器的強度圖，在強度圖上可求得角反射器RCS，單位是m2。角反射器對SAR的后向散射信號可看作是脈沖響應(yīng)函數(shù)(Curlander，1991)，則角反射器的回波信號包絡(luò)面可表示為(薛笑榮等，2009):

式(5)中:x為方位角(azimuth coordinate)，t為斜距時間(slant range time)，Ls為合成孔徑長度;τp為脈沖長度;ρr為斜距向分辨率;ρa為方位向分辨率。

由式(5)可知，在理想狀態(tài)下，角反射器強度在斜距向和方位向都呈現(xiàn)為sinc函數(shù)形狀，其模擬圖像是對稱的十字型。

當(dāng)波長、入射角等成像參數(shù)為同一模式的SAR傳感器對角反射器成像時，角反射器表現(xiàn)出強而穩(wěn)定的后向散射特性，由于回波信號較強，角反射器在強度圖上表現(xiàn)為比較亮的色調(diào)。在進行角反射器初定位后，截取大小為25×25像素的影像圖，在ENVISAT數(shù)據(jù)的強度圖中18臺角反射器表現(xiàn)為亮點(圖3)，而在TerraSAR數(shù)據(jù)的強度圖中12臺角反射器均表現(xiàn)出十字型特征(圖4)。角反射器模擬強度圖像是十字型，但在實際的SAR強度影像中，角反射器的影像特征主要受其信號強弱的影響。由式(4)可知，當(dāng)角反射器安裝固定后，其后向散射的強度主要受波長影響，波長越短，角反射器后向散射的強度值越大，所以角反射器在TerraSAR數(shù)據(jù)的強度影像呈十字型。

圖3 角反射器在Envisat數(shù)據(jù)中的強度圖像Fig.3 The intensity image of CRs in Envisat data.

圖4 角反射器在TerraSAR數(shù)據(jù)中的強度圖像Fig.4 The intensity image of CRs in TerraSAR data.

MSR圖像是光譜信息在距離向和方位向上，取不同帶寬進行多視處理后統(tǒng)計的均值和標(biāo)準(zhǔn)差的比率。在ENVISAT數(shù)據(jù)的MSR圖像上，由于角反射器的MSR統(tǒng)計值與背景值差別很小，所以18臺角反射器在MSR圖像中無法分辨，圖像都為黑色(圖5)。在TerraSAR數(shù)據(jù)的MSR圖像上(圖6)，盡管角反射器MSR影像特征明顯，但表現(xiàn)在多個像素上。

相關(guān)圖像記錄的是在進行多視處理時光譜的相關(guān)性信息。相關(guān)性指標(biāo)的優(yōu)勢在于點目標(biāo)的高相關(guān)值通常僅聚焦在一個像素內(nèi)。圖7和圖8分別是ENVISAT數(shù)據(jù)和TerraSAR數(shù)據(jù)5×5像素大小的角反射器相關(guān)圖像，從每個角反射器對應(yīng)的相關(guān)圖像上能夠直接分辨出中間像素的相關(guān)值比較高。

圖5 角反射器在ENVISAT數(shù)據(jù)中的MSR圖像Fig.5 The MSR image of CRs in Envisat data.

圖6 角反射器在TerraSAR數(shù)據(jù)中的MSR圖像Fig.6 The MSR image of CRs in TerraSAR data.

圖7 角反射器在ENVISAT數(shù)據(jù)中的光譜相關(guān)圖像(5×5)Fig.7 The spectral correlation image of CRs in Envisat data.

圖8 角反射器在TerraSAR數(shù)據(jù)中的光譜相關(guān)圖像(5×5)Fig.8 The spectral correlation image of CRs in TerraSAR data.

在強度圖像中角反射器影像特征明顯，但角反射器的強度體現(xiàn)在相鄰多個像素中;角反射器影像特征在MSR圖像中表現(xiàn)不明顯;在相關(guān)圖像上角反射器影像特征沒有強度影像特征明顯，與強度圖像相比需要放大25倍才能識別出角反射器相關(guān)信息，但角反射器的相關(guān)信息聚焦在一個像素內(nèi)，而且相關(guān)值明顯高于周邊像素，這使角反射器定位容易而且準(zhǔn)確。

4 角反射器精定位最優(yōu)參數(shù)的選擇

選用反映角反射器輻射特性的3個參數(shù):強度σ、相關(guān)系數(shù)γ、平均值與標(biāo)準(zhǔn)差的比率MSR來研究角反射器精確定位。Envisat數(shù)據(jù)中角反射器初定位結(jié)果在行值上有4行的系統(tǒng)誤差，以系統(tǒng)誤差校正后的行列值為中心，在5×5窗口內(nèi)搜索，最大值σ和γ值對應(yīng)的像素位置定為角反射器行列值。因為角反射器參數(shù)MSR與背景差別極小(圖5)，所以對參數(shù)MSR不做數(shù)值分析。表2中列出了強度σ和相關(guān)系數(shù)γ的最大值，以及分別用這2個參數(shù)確定的角反射行列值。在表2中粗體字給出的是2個參數(shù)定位結(jié)果不同的7臺角反射器，這7臺角反射器在行或列上有1個像素差別。圖9描述了這7臺角反射器不同參數(shù)定位結(jié)果對應(yīng)像素的相位對比結(jié)果。圖中每一組中間的綠色柱是角反射器所在像素的相位，左右2個柱狀圖是與角反射器同行相鄰2個像素的相位;上圖是用參數(shù)σ定位的角反射器及同行相鄰像素的相位值，下圖是用參數(shù)γ定位的角反射器及同行相鄰像素的相位值;上下圖對比可知，下圖每組同行相鄰的3個像素相位值符號相同，而且數(shù)值差別小于上圖每組3個像素之間相位的差別。其他11臺角反射器2個參數(shù)定位結(jié)果相同，角反射器像素和同行相鄰像素的相位差別也很小(圖10)。表2中角反射器CR1和CR2相距約50m，CR3和CR4相距約90m，這2組角反射器因距離近，強度值彼此影響的原因，所以選用強度參數(shù)定位的結(jié)果是同一像素。

表2 Envisat數(shù)據(jù)參數(shù)最大值及其定位行列值Table 2 The max value of parameters and the pixels of CRs located by these parameters in Envisat data

圖9 Envisat數(shù)據(jù)中不同參數(shù)定位角反射器相位對比圖(縱軸為弧度)Fig.9 The phase contrast figure of CRs located with different parameters in Envisat data(Y-axis is radians).

圖10 Envisat數(shù)據(jù)中角反射器及同一行左右相鄰像素相位圖(縱軸為弧度)Fig.10 The phase of CR and its adjacent pixels at same line in Envisat data(Y-axis is radians).

在Envisat數(shù)據(jù)中，通過對角反射器定位后的強度、相關(guān)性和相位特性的分析可知，角反射器在Envisat數(shù)據(jù)中的強度受背景干擾大;角反射器的相關(guān)性聚焦性好，不受背景影響;用相關(guān)性參數(shù)定位后的角反射器對應(yīng)相位穩(wěn)定，表現(xiàn)在角反射器相位值和同行相鄰的2個像素相位值差別小。而在TerraSAR中，角反射器散射幅度強，強度信息受背景干擾性小，所以在用強度參數(shù)和相干性參數(shù)定位的角反射器位置基本一致，只有表3中粗體字所示1個角反射器像素位置在行上差一個像素;而用MSR定位的角反射器像素位置與表3中前2個參數(shù)定位的像素行列值都有差別。圖11是TerraSAR數(shù)據(jù)用參數(shù)γ定位角反射器及其周圍8鄰域像素的相位柱狀圖，可見角反射器及其周圍像素的相位值基本相等，這是角反射器強散射效應(yīng)的結(jié)果，據(jù)此說明了用參數(shù)γ定位角反射器像素位置的準(zhǔn)確性。

表3 TerraSAR數(shù)據(jù)定位參數(shù)及其定位行列值Table 3 The max value of parameters and the pixels of CRs located by these parameters in TerraSAR data

圖11 TerraSAR數(shù)據(jù)中角反射器及其8鄰域像素相位圖(縱軸為弧度)Fig.11 The phase of CR and its 8 neighborhood pixels in TerraSAR data(Y-axis is radians).

5 結(jié)論與討論

(1)分析總結(jié)了Envisat、TerraSAR兩種數(shù)據(jù)對應(yīng)的角反射器強度、相關(guān)性和MSR影像特征:1)角反射器強度影像特征明顯，在多個像素中呈現(xiàn)，Envisat數(shù)據(jù)中角反射器強度表現(xiàn)為亮點，而在TerraSAR數(shù)據(jù)中則呈現(xiàn)明顯的十字形狀;2)相反角反射器MSR影像特征不明顯，僅在TerraSAR數(shù)據(jù)中有影像特征;3)在2種數(shù)據(jù)中，角反射器的相關(guān)性都聚焦在一個像素內(nèi)，而且相關(guān)值明顯高于周圍像素。

(2)對比分析不同參數(shù)定位出的角反射器像素位置及其對應(yīng)的相位信息可知:1)不同參數(shù)定位角反射器像素位置稍有不同，但對應(yīng)相位差別很大;2)參數(shù)γ定位角反射器像素的相位比較穩(wěn)定，體現(xiàn)在Envisat數(shù)據(jù)中角反射器像素和同行上相鄰像素的相位值差別小，而TerraSAR數(shù)據(jù)中角反射器及其8鄰域像素的相位值基本相同，這是角反射器強散射效應(yīng)抑制噪聲的結(jié)果;3)Envisat數(shù)據(jù)中角反射器像素和同行上相鄰像素的相位值差別小，表明角反射器效應(yīng)在方位向(也就是衛(wèi)星飛行方向)表現(xiàn)明顯，這可能與傳感器按行記錄回波信息有關(guān)。

(3)本文結(jié)果表明參數(shù)γ定位角反射器像素位置較其他參數(shù)準(zhǔn)確，下一步的研究工作將在時序影像集中進一步檢驗這一結(jié)論。

程滔，單新建.2007.CR、PS干涉測量聯(lián)合解算算法研究［J］.地震，27(2):64—71.

CHENG Tao，SHAN Xin-jian.2007.Research on algorithm of CRInSAR and PSInSAR combined calculation ［J］.Earthquake，27(2):64—71(in Chinese).

范景輝.2008.基于相干目標(biāo)的DInSAR技術(shù)地表形變監(jiān)測研究與應(yīng)用［D］:［學(xué)位論文］.北京:中國科學(xué)院遙感應(yīng)用研究所.82—87.

FAN Jing-hui.2008.Research on coherent targets based DInSAR technique and its applications to surface deformation monitoring［D］.PhD Dissertation.Institute of Remote Sensing Applications，Chinese Academy of Sciences，Beijing.82—87(in Chinese).

李慶揚，關(guān)治，白峰杉.2000.數(shù)值計算原理［M］.北京:清華大學(xué)出版社.

LI Qing-yang，GUAN Zhi，BAI Feng-shan.2000.Numerical Calculation Principle［M］.Tsinghua University Press，Beijing(in Chinese).

邢學(xué)敏，丁曉利，朱建軍，等.2011.CRInSAR與PSInSAR聯(lián)合探測區(qū)域線性沉降研究［J］.地球物理學(xué)報，54(5):1193—1204.doi:10.3969/j.issn.0001-5733.2011.05.008.

XING Xue-min，DING Xiao-li，ZHU Jian-jun，et al.2011.Detecting the regional linear subsidence based on CRInSAR and PSInSAR integration［J］.Chinese J Geophys，54(5):1193—1204(in Chinese).

許才軍，何平，溫?fù)P茂，等.2012.利用CR-InSAR技術(shù)進行鮮水河斷層的形變速率監(jiān)測［J］.武漢大學(xué)學(xué)報· 信息科學(xué)版，37(3):302—305.

XU Cai-jun，HE Ping，WEN Yang-mao，et al.2012.Crustal deformation monitoring of Xianshuihe Fault by CRInSAR［J］.Geomatics and Information Science of Wuhan University，37(3):302—305(in Chinese).

薛笑榮，王愛民，曾琪明.2009.基于圖像模擬的SAR圖像角反射器檢測方法研究［J］.安陽師范學(xué)院學(xué)報，5:69—72.

XUE Xiao-rong，WANG Ai-min，ZENG Qi-ming.2009.Research on SAR image corner reflector detection method［J］.Journal of Anyang Normal University，5:69—72(in Chinese).

張景發(fā)，龔利霞，姜文亮.2006.PS-InSAR技術(shù)在地殼長期緩慢形變監(jiān)測中的應(yīng)用［J］.國際地震動態(tài)，6:1—6.

ZHANG Jing-fa，GONG Li-xian，JIANG Wen-liang.2006.Application of PSInSAR technique to measurement of long-term crustal deformation［J］.Recent Developments in World Seismology，6:1—6(in Chinese).

Curlander J C.1982.Location of space-borne SAR imagery［J］.IEEE Transactions on Geoscience and Remote Sensing，20(3):359—364.doi:10.1109/TGRS.1982.350455.

Curlander J C，McDonough R N.1991.Synthetic Aperture Radar［M］.Wiley Interscience，New York，p396.

Kbittm E F.1981.A tool for predicting the radar cross section of an arbitrary corner reflector［R］.IEEE Southeastcon’81 Conference，17—20.

Mohr J J，Madsen S N.2001.Geometric calibration of ERS satellite SAR images ［J］.IEEE Transactions on Geoscience and Remote Sensing，39(4):842—850.doi:10.1109/36.917909.

Merrill I S.2008.Radar Handbook(Third Edition)［M］.Publishing House of Electronics Industry.545.

Renaud F，Olivier H，Emmanuel T，et al.2010.Monitoring temperate glacier displacement by multi-temporal TerraSAR-X images and continuous GPS measurements［J］.IEEE Journal of Selected Topics in Applied Earth Observation and Remote Sensing，1—26(submitted).doi:10.1109/JSTARS.2010.2096200.

Sakurai-Amano T，Kobayashi S，F(xiàn)ujii N.1999.Detection of singular corner reflectors in residential and mountainous areas from SAR images［A］.In:IEEE Proc，1999 International，2:1454—1456.doi:10.1109/IGARSS.1999.774662.

Ulander L M.1991.Accuracy of using point targets for SAR calibration［J］.IAE，27:139—148.doi:10.1109/7.68156.

Xie Ye，Kaufmann H，Xiaofang G.2002.Differential SAR Interferometry Using Corner Reflectors［A］.In:IEEE Proc，2002 International Geoscience and Remote Sensing Symposium，2:1243—1246.doi:10.1109/IGARSS.2002.1025902.

Xie Ye，Kaufmann H，Xiaofang G.2004.Landslide monitoring in the Three Gorges area using D-InSAR and corner reflectors［J］.Photogrammetric Engineering and Remote Sensing，70(10):1167—1172.

Xie Ye.2008.CR-Based SAR-Interferometry for landslide monitoring［A］.In:IEEE Proc，2008 International Geoscience and Remote Sensing Symposium，2:1239—1242.doi:10.1109/IGARSS.2008.4779226.