文章信息
- 吴文豪, 李陶, 龙四春, 周志伟, 祝传广
- WU Wenhao, LI Tao, LONG Sichun, ZHOU Zhiwei, ZHU Chuanguang
- RADARSAT-2卫星TOPS模式影像干涉处理
- Interferometric Processing of RADARSAT-2 TOPS Data
- 武汉大学学报·信息科学版, 2020, 45(1): 134-140
- Geomatics and Information Science of Wuhan University, 2020, 45(1): 134-140
- http://dx.doi.org/10.13203/j.whugis20180212
-
文章历史
收稿日期: 2018-12-09
2. 武汉大学卫星导航定位技术研究中心, 湖北 武汉, 430079;
3. 中国科学院测量与地球物理研究所, 湖北 武汉, 430077
2. GNSS Research Center, Wuhan University, Wuhan 430079, China;
3. Institute of Geodesy and Geophysics, Chinese Academy of Sciences, Wuhan 430077, China
加拿大合成孔径雷达遥感卫星(RADARSAT-2)搭载的C波段传感器具有多种极化方式和成像模式,2013年欧空局联合加拿大航天局基于RADARSAT-2卫星宽幅扫描模式(ScanSAR)成功进行了宽幅干涉模式(terrain observation by progressive scans,TOPS)成像,为即将发射的哨兵合成孔径雷达遥感卫星(Sentinel-1)进行了技术验证。RADARSAT-2卫星TOPS模式采用3个子测绘带,侧视角范围为[27.5°, 42.2°],斜视角最大为±0.6°,发射脉冲带宽为30 MHz[1]。尽管其聚焦影像标称质量与Sentinel-1卫星基本相同,但由于RADARSAT-2卫星轨道精度较低,影像定位精度低于Sentinel-1卫星,导致其干涉处理时难以达到0.001个像元量级的配准精度[2]。
卫星轨道精度不但决定了几何配准精度,也影响着TOPS模式影像的配准方案。Sentinel-1卫星精密轨道精度约为0.05 m,直接进行几何配准即可达到较高的配准精度,利用增强谱分集可进一步提高方位向配准精度[3-8]。即使是精度较低的实时轨道,只需采用相关配准方法对几何配准结果进行优化即可满足增强谱分集对初配准的精度要求[9-11]。然而RADARSAT-2卫星轨道精度为0.65 m(置信度99.99%),而且轨道夹角较大,几何配准精度远低于Sentinel-1卫星,如果只采用相关配准方法改正几何配准结果,将导致增强谱分集相位缠绕。因此需要在Sentinel-1卫星干涉处理方法基础上,针对RADARSAT-2卫星影像进一步完善其干涉处理流程,为低精度轨道卫星TOPS影像干涉处理积累经验[12-13]。
1 TOPS成像TOPS模式成像时天线波束由后向前扫描,确保了每个目标的波束方向图、多普勒模糊的一致性,延长了连续脉冲序列(Burst)成像时间,降低了回波信号的边缘效应[14-15]。但天线波束的旋转效应导致每个Burst雷达回波相位历程在慢时间维存在明显差异,如图 1所示。图 1中,横轴为方位向时间η,纵轴表示多普勒频率fη。绿色阴影区为雷达回波有效信号,红色阴影区为Burst聚焦后的方位向频谱,目标P1、P2分别为第一个和最后一个有效Burst回波信号。雷达回波进行零多普勒聚焦时,每个目标对应的回波多普勒中心与聚焦影像多普勒中心并未发生变化,而是将回波时间转换至零多普勒时间,这导致了聚焦影像持续时间Tfocused大于回波持续时间Tb,聚焦Burst影像之间存在重叠,如目标P1原始回波起始时间为ηRaw,对应的零多普勒起始时间为ηZD。RADARSAT-2卫星Burst间方位向重叠度设计为10%,不仅能确保地表覆盖的连续性,也为增强谱分集配准提供了足够的观测单元。
多普勒中心是雷达成像与干涉的重要参数,理论上聚焦影像多普勒中心分布可以根据卫星星历计算得出。鉴于RADARSAT-2卫星轨道精度较差,其轨道有效速度应优先基于方位调频率ka进行计算。因此雷达回波多普勒中心调频率ks(η)可由波束旋转角速度kψ、雷达波长λ、零多普勒斜距时间τ0以及电磁波速度c等参数得到:
${k_s}\left( \eta \right) = \sqrt {\left| {\frac{{{k_a}c{\tau _0}}}{\lambda }} \right|} \cdot \frac{{\rm{ \mathit{ π} }}}{{180}} \cdot {k_\psi }$ | (1) |
聚焦影像方位向多普勒中心调频率kt可进一步表示为:
${k_t} = \frac{{{k_a} \cdot {k_s}}}{{{k_a} - {k_s}}}$ | (2) |
RADARSAT-2卫星影像提供了每组Burst参考时间点的方位调频率和多普勒中心参数,根据式(2)可得到聚焦Burst影像多普勒中心分布,如图 2所示,其聚焦影像多普勒在方位向发生较大变化,此时Burst聚焦方位线具有时变信号特性。如果多普勒中心参数存在误差,干涉处理时将直接影响到图像重采样、谱分集和增强谱分集估计结果[16-17]。当对雷达影像进行插值或重采样时需要根据多普勒中心估计结果对插值核进行调制或者对插值信号进行去调频处理,使插值信号与内插核频谱相匹配,否则将导致相位信息泄露,干涉相位误差增大。谱分集对影像进行方位向频带分解时,需要根据影像多普勒中心参数将方位线解调至基带信号;增强谱分集需要依据重叠区域的多普勒中心参数影像的残余偏移量。准确获取聚焦影像多普勒中心是TOPS模式进行高精度干涉处理的前提条件。
2 多普勒中心估计尽管RADARSAT-2卫星在椭圆轨道上运行时每个历元姿态和速率均存在差异,且基于星历参数的多普勒中心估计精度有限,但每组Burst已根据多普勒中心进行了初步去调频处理,影像方位线已接近于基带信号,此时进行残余多普勒中心估计无需考虑多普勒中心模糊问题[18]。因此本文采用相位增量法进一步估计残余多普勒中心参数。相位增量法通过测量雷达回波或聚焦影像复数据的相位变化来估计基带信号的多普勒中心,具体方法如下:假设s(η)为雷达影像方位线上的某一采样单元,
$C\left( \eta \right) = \mathop \sum \nolimits^ {s^{\rm{*}}}\left( \eta \right)s\left( {\eta + {\rm{\Delta }}\eta } \right)$ | (3) |
平均相位增量类似于求取延迟为1方位线的自相关函数,计算效率较高。为保证残余多普勒中心估计结果的精度和可靠性,运算时需对雷达回波数据实部和虚部进行符号化处理,以消除雷达数据中噪声以及较强后向散射目标对平均相位增量估计结果的干扰。在频谱对称前提下,残余多普勒中心
$f{'_{\eta c}} = \frac{{{F_a}}}{{2{\rm{ \mathit{ π} }}}}{\rm{arg}}\left( {c\left( \eta \right)} \right)$ | (4) |
式中,arg表示复数的辐角;Fa为RADARSAT-2卫星聚焦影像方位向采样频率,一般为400 Hz。鉴于雷达波束姿态变化是连续的,多普勒中心在影像内一般不会出现跳变,因此本文将Burst影像划分为多个格网单元,如图 3所示。
基于式(4)获取影像中每个格网的残余多普勒中心,然后进行解缠处理,获得残余多普勒中心在Burst中的分布和残余多普勒中心调频率Δk t。从图 3中可以看出残余多普勒中心最大值为60 Hz,也就意味着基于星历参数的多普勒中心估计误差是不能忽略的,因此影像最终的多普勒分布应为残余多普勒中心的估计结果与基于星历参数的多普勒中心估计结果之和,如图 4所示。本文也将多普勒分布称之为联合多普勒估计结果。
相位增量法估计多普勒中心克拉美-罗(Cramer-Rao)下界精度为:
$\sigma = 0.340\;7 \times \frac{{{F_a}}}{{\sqrt N }}$ | (5) |
式中,N为分辨率单元数量。理论上分辨率单元数量越多,影像残余多普勒中心估计精度越高。但格网并非越大越好,还需顾及到格网多普勒中心解缠时的梯度。本文对Burst影像处理时依据过采样率以及方位线距离线间的比例,将每个格网设置为200×1 000像素,根据式(5)可知其残余多普勒中心估计精度约为0.413 6 Hz,对应格网梯度只有5 Hz,远低于理论梯度极限值200 Hz要求,因此残余多普勒中心估计结果是可靠的。多普勒中心应为两次多普勒中心估计之和,如图 4所示。相对应方位向多普勒中心斜率
$k_t^{{\rm{EST}}}\left( \eta \right) = {k_t}\left( \eta \right) + {\rm{\Delta }}{k_t}\left( \eta \right)$ | (6) |
因此对Burst影像进行去调频处理时,需重新按照
$\Delta \omega = \frac{{\left| {k_t^{{\rm{EST}}}} \right|}}{{{B_{az}}}} \approx 6{\rm{Hz}}$ | (7) |
RADARSAT-2卫星TOPS模式影像Δω约为6 Hz,相对于229 Hz的多普勒带宽,该值可以忽略不计,所以去调频时无需对影像进行过采样处理。本文进一步根据影像方位向功率谱验证多普勒中心估计结果。理论上,均匀场景下方位向功率谱和天线波束方向图是一致的,方位功率谱的峰值位置对应于多普勒中心。
图 5为RADARSAT-2卫星TOPS模式瞬时功率谱与多普勒中心对比示意图,其横坐标为聚焦影像所对应的零多普勒时间,纵坐标为功率谱频率值,背景颜色表示影像功率谱,黄色对应功率谱峰值。TOPS模式影像聚焦处理后通常对Burst影像进行降采样处理,以确保3个子测绘带方位向采样频率一致,使方位向采样频率仅是目标多普勒带宽的1.1~1.3倍,远低于方位向信号带宽,因此功率谱将出现多次周期性的折叠。图 5(a)红色折线为基于卫星星历和成像参数估计的多普勒中心结果,图 5(b)红色折线为多普勒中心联合估计结果。实验结果证明了联合估计多普勒中心参数更加准确。
3 TOPS模式影像干涉处理RADARSAT-2卫星Burst方位向重叠区域多普勒中心基线约为5 000 Hz,理论上初始配准精度需达到0.05个像素以内,才能确保增强谱分集不会出现相位缠绕问题。但RADARSAT-2卫星轨道精度差,几何配准通常无法满足增强谱分集对初始配准精度的要求。因此该卫星TOPS模式影像干涉处理的关键是利用相关配准或谱分集改善几何配准精度。
由于每个Burst持续时间只有3.6 s,其几何配准误差可认为是一恒定值,残余偏移量估计窗口为整个Burst影像,即15 185×1 433像素。依据相关配准和谱分集配准精度Cramer-Rao公式即可得到每个Burst影像极限配准精度,如图 6所示。从图 6中可以看出谱分集配准精度优于相关配准,而且谱分集在相干性很低条件下其精度也优于0.05个像素。但考虑到相关配准对初始配准精度没有要求,本文首先在几何配准的基础上,采用相关配准对配准结果进行一次方位向偏移量校正,消除大尺度方位向残余偏移量,然后采用谱分集进行二次校正,以保证增强谱分集估计的正确性,具体操作流程如图 7所示。基于优化后的配准结果完成影像重采样处理后,即可按照常规的干涉方法及流程进行去平地效应、干涉图滤波等处理。
尽管相关配准理论上可以达到较高精度,但是实际操作中其配准结果易受到影像地形、多项式阶数等因素的影响。本文分别通过提高单个窗口的配准精度和增加观测量提高影像整体配准精度。相关配准时选择以点目标作为配准窗口位置中心,使窗口的特征信息也参与影像匹配,提高配准的可靠性和精度。尽管TOPS模式影像每个Burst均独立成像,但在同一轨道面内均进行零多普勒聚焦,因此相关配准结果可按照配准窗口的零多普勒时间进行拼接处理,以增加多项式拟合的观测量。图 8、图 9分别为同一子测绘带内14组Burst相关配准获取的偏移量距离向、方位向拼接结果,图中颜色代表了不同的Burst。14组Burst影像偏移量的正确拼接相当于将子测绘带配准虚拟化为条带模式影像配准,可以避免Burst间的几何配准校正值出现跳变问题,从而确保增强谱分集配准上下频带干涉图初始配准精度保持一致。
RADARSAT-2卫星TOPS影像即使通过几何配准和相关配准处理,已消除大尺度的方位向残余偏移量,但仍无法确保其满足增强谱分集的要求。图 10为南美洲墨西哥高原地区TOPS影像第3子测绘带偏移量估计结果。几何配准误差在短时间内可视为系统误差,如果增强谱分集估计正确,其结果应该是连续的。如果不采用谱分集校正相关配准结果,增强谱分集会出现严重的相位缠绕,其估计出的残余偏移量将出现明显跳变,如图 10(a)所示。针对这一问题,本文增加谱分集方法对几何配准结果做进一步校正。为确保几何配准校正值连续,谱分集估计时每个Burst依然需要进行方位向拼接,并采用低阶多项式对估计结果进行拟合处理。谱分集配准完成后,紧接着进行增强谱分集配准,最终完成影像配准和干涉处理。图 10(b)为几何配准后依次使用相关配准、谱分集、增强谱分集偏移量改正结果,可以看出相关配准虽然消除了大尺度方位向残余偏移量,但影像残余偏移量依然达到0.2个像素,无法满足增强谱分集的要求,这也解释了图 10(a)中增强谱分集出现缠绕的原因。因此对于RADARSAT-2卫星利用谱分集优化几何配准结果是必要的,谱分集配准可进一步消除残余偏移量,确保增强谱分集估计结果正确。
图 11为增加谱分集配准处理前后的干涉图对比结果。图 11(a)为未采用谱分集配准所得到的干涉图,可以发现Burst干涉图在方位向存在相位跳变,这是增强谱分集出现严重的相位缠绕,影像间的残余偏移量估计错误造成的。然而如果按照本文所提出的方法增加谱分集配准,如图 11(b)所示,其最终Burst干涉图间相位连续,与之相对应的图 10(b)中增强谱分集估计结果也是连续的,符合几何配准误差分布规律。上述实验证明了本文所提方法的可行性。
4 结语本文首先探讨了RADARSAT-2卫星TOPS模式影像干涉处理时的多普勒中心估计问题。由于仅依据星历参数难以准确估计多普勒中心,因此本文在初步去调频处理的基础上采用相位增量进行残余多普勒中心估计;然后针对RADARSAT-2卫星轨道精度过低的特点,在几何配准的基础上,采用相关配准和谱分集依次进行方位向配准偏移量的校正;最后采用增强谱分集估计残余偏移量并实施配准误差的校正,从而实现RADARSAT-2卫星TOPS模式影像高精度配准和干涉处理,为低精度轨道卫星TOPS影像干涉处理积累经验。
致谢: 感谢中南大学地球科学与信息物理学院的段梦、浙江环球星云遥感科技有限公司的殷幼松提出的宝贵建议。[1] |
Attema E, Davidson M, Floury N, et al. Sentinel-1 ESA's New European Radar Observatory[C]. Synthetic Aperture Radar (EUSAR), Friedrichshafen, Germany, 2008 https://www.researchgate.net/publication/224232803_Sentinel-1_ESAAs_New_European_Radar_Observatory
|
[2] |
Prats-Iraola P, Rodriguez-Cassola M, Zan F D, et al. Role of the Orbital Tube in Interferometric Spaceborne SAR Missions[J]. IEEE Geoscience and Remote Sensing Letters, 2015, 12(7): 1 486-1 490. DOI:10.1109/LGRS.2015.2409885 |
[3] |
Sansosti E, Berardino P, Manunta M, et al. Geometrical SAR Image Registration[J]. IEEE Transactions on Geoscience and Remote Sensing, 2006, 44(10): 2 681-2 870. |
[4] |
Wu Wenhao, Zhou Zhiwei, Li Tao, et al. A Study of Sentinel-1 TOPS Mode Co-registratgion[J]. Acta Geodaetica et Cartographica Sinica, 2017, 46(9): 1 156-1 164. (吴文豪, 周志伟, 李陶, 等. 精密轨道支持下的哨兵卫星TOPS模式干涉处理[J]. 测绘学报, 2017, 46(9): 1 156-1 164. ) |
[5] |
Wu Wenhao. TOPS Interferometry with Sentinel-1[D]. Wuhan: Wuhan University, 2006 (吴文豪.哨兵雷达卫星TOPS模式干涉处理研究[D].武汉: 武汉大学, 2016) http://cdmd.cnki.com.cn/Article/CDMD-10486-1016124261.htm
|
[6] |
Prats P, Marotti L, Wollstadt S, et al. TOPS Interferometry with TerraSAR-X[J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(8): 3 179-3 188. DOI:10.1109/TGRS.2011.2178247 |
[7] |
Martinez N Y, Prats P, Gonzalez F R, et al. Interferometric Processing of Sentinel-1 TOPS Data[J]. IEEE Transactions on Geoscience and Remote Sensing, 2016, 54(4): 2 220-2 234. DOI:10.1109/TGRS.2015.2497902 |
[8] |
Fattahi H, Agram P, Simons M. A Network-Based Enhanced Spectral Diversity Approach for TOPS Time-Series Analysis[J]. IEEE Transactions on Geoscience and Remote Sensing, 2016, 55(2): 777-786. |
[9] |
Meta A, Prats P, Steinbrecher U, et al. First TOPSAR Image and Interferometry Results with TerraSAR-X[J]. American Heart Journal, 2008, 36(6): 943-949. |
[10] |
Urs W, Charles W, Tazio S, et al. Sentinel-1 Support in the GAMMA Software[J]. Procedia Computer Science, 2016, 80(3): 100-105. |
[11] |
Xu B, Li Z, Feng G, et al. Continent-Wide 2-D Co-Seismic Deformation of the 2015 Mw 8.3 Illapel, Chile Earthquake Derived from Sentinel-1A Data: Correction of Azimuth Co‑Registration Error[J]. Remote Sensing, 2016, 8(5): 1-12. |
[12] |
Scheiber R, Jager M, Prats P, et al. Speckle Tracking and Interferometric Processing of TerraSAR-X TOPS Data doe Mapping Nonstationary Scenarios[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2014, 8(4): 1 709-1 720. |
[13] |
Cassola M R, Prats P, Zan F D, et al. Doppler-related Distortions in TOPS SAR Images[J]. IEEE Transactions on Geoscience and Remote Sensing, 2015, 53(1): 25-35. DOI:10.1109/TGRS.2014.2313068 |
[14] |
Wu Wenhao, Wang Mingzhou, Li Sha, et al. Analysis of Sliding Spotlight SAR Interferometry[J]. Geomatics and Information Science of Wuhan University, 2015, 40(12): 1 588-1 593. (吴文豪, 王明洲, 李沙, 等. 滑动聚束模式SAR影像干涉处理方法[J]. 武汉大学学报·信息科学版, 2015, 40(12): 1 588-1 593. ) |
[15] |
Bai Xia, Sun Jinping, Hong Wen, et al. On the TOPS Mode Spaceborne SAR[J]. Scientia Sinica Informationis, 2010, 40(3): 468-480. (白霞, 孙进平, 洪文, 等. 星载合成孔径雷达TOPS模式研究[J]. 中国科学:信息科学, 2010, 40(3): 468-480. ) |
[16] |
Huang P, Li S, Xu W. Investigation on Full Aperture Multichannel Azimuth Data Processing in TOPS[J]. IEEE Geoscience and Remote Sensing Letters, 2014, 11(4): 728-732. DOI:10.1109/LGRS.2013.2278183 |
[17] |
Prats P, Scheiber R, Mittermayer J, et al. Processing of Sliding Spotlight and TOPS SAR Data Using Baseband Azimuth Scaling[J]. IEEE Transactions on Geoscience and Remote Sensing, 2010, 48(2): 770-780. DOI:10.1109/TGRS.2009.2027701 |
[18] |
Xu W, Huang P P, Wang R, et al. Processing of Multichannel Sliding Spotlight and TOPS Synthetic Aperture Radar Data[J]. IEEE Transactions on Geoscience and Remote Sensing, 2013, 51(8): 4 417-4 429. DOI:10.1109/TGRS.2013.2265306 |