CN114114181A

CN114114181A - 基于轨道误差相位基的星载sar干涉基线矫正方法

Info

Publication number: CN114114181A
Application number: CN202210106363.8A
Authority: CN
Inventors: 杨晴月; 王樱洁; 王宇; 郑明洁
Original assignee: Aerospace Information Research Institute of CAS
Current assignee: Aerospace Information Research Institute of CAS
Priority date: 2022-01-28
Filing date: 2022-01-28
Publication date: 2022-03-01
Anticipated expiration: 2042-01-28
Also published as: CN114114181B

Abstract

本发明公开一种基于轨道误差相位基的星载SAR干涉基线矫正方法，涉及导航和雷达测量领域。该方法用景中心基线长度和基线变化速率平行、垂直雷达视线方向的分量共四个参数来描述干涉基线，并根据建立的轨道误差相位观测方程，将由残差相位近似的轨道误差相位投影到四个能够反映基线误差参数信息的轨道误差相位基以及经地理编码到雷达坐标系的外部数字高程模型上，再根据投影系数矫正各个基线参数。方法通过中值池化、最小二乘解算和迭代来保证算法的精度和鲁棒性。该方法普适性强，获得的基线精度较高，收敛速度快，且在实现上更简单、直观。应用本方法可以在根源上缓解基线误差对干涉处理带来的影响，提高干涉、差分干涉、时序干涉处理的精度。

Description

基于轨道误差相位基的星载SAR干涉基线矫正方法

技术领域

本发明属于导航和雷达测量领域，是一种星载合成孔径雷达（SyntheticAperture Radar，SAR）干涉基线矫正方法，可应用于干涉、差分干涉、时序干涉处理中。

背景技术

合成孔径雷达干涉（Interferometric SAR，InSAR）和差分干涉（DifferentialInSAR，D-InSAR）由于具有全天时、全天候、分辨率高、覆盖范围广的工作能力已成为对地观测的重要工具。InSAR和D-InSAR利用SAR图像的相位信息获取地表高程、地表在雷达视线方向的形变信息，广泛应用于地形测绘、地表沉降监测以及自然灾害监测。

以重轨干涉为例，干涉基线是SAR照射地表同一目标时主卫星天线指向辅卫星天线的空间矢量，是干涉处理的重要参数。若基线估计误差较大则会严重影响干涉处理的精度。具体表现在两个方面：基线误差造成平地相位计算误差，即轨道误差相位，所产生的影响通常称为“相位倾斜”，基线误差还导致地形相位与高程转换系数的误差，从而影响InSAR生成数字高程模型（Digital Elevation Model，DEM）的准确性和D-InSAR处理中的地形相位去除过程。

已经有许多研究致力于缓解基线误差造成的影响，这些方法可以分为两大类：第一类方法专注于缓解“相位倾斜”的影响，以某种模型拟合轨道误差相位，再从差分相位中减去，以重建真实的差分相位。另一类方法通过精确地估计基线来从根本上解决基线误差对干涉处理的影响。

在第一类方法中，常用的模型有定义在雷达结构即斜距和方位时间坐标下的线性模型、二次模型。线性模型通常用于较小范围的近似，然而在整幅干涉图中，线性模型就不足以拟合轨道误差相位，这时常用到二次模型。这一类方法操作简单，可以有效缓解相位倾斜问题，已经被广泛应用于差分干涉和时序干涉处理中。但轨道误差相位并不是严格多项式形式的，应用这种方法会带来一个小的偏差。另外，即便相位倾斜被很好地补偿了，地形相位与高程的转换系数仍不准确，这会使得差分干涉去除地形相位不准确，尤其在山地地区会带来较大影响。

第二类方法根据所依赖的信息可以分为基于图像自身信息的基线估计方法和基于外部信息的基线估计方法。基于图像自身信息的基线估计方法通常用到基线与配准的关系、基线与干涉图频率条纹的关系来解算基线。这类方法由于不依赖其它信息，因此对干涉数据本身的质量、特点会提出要求，而所提出的要求在大部分情况下不一定能满足。另外，处理过程中一些参数的选取也会给基线估计带来一定的不确定性。基于外部信息的基线估计方法常用到两种外部信息。第一种是地面控制点（Ground Control Points，GCP）信息，通过建立干涉观测方程，结合若干GCP信息解算观测方程参数，进而得到基线信息。这一方法理论精度较高，然而在很多情况下，我们无法获得感兴趣区域的地面控制点，而且在控制点信息不够精确的情况下，这一方法的效果将大打折扣。第二种是外部DEM信息，给定基线初值后，这类方法要么通过不断调节基线使得借助DEM模拟的干涉数据的某一特性逼近原始干涉数据的某一特性，要么借助DEM来获取残差相位，再根据基线与之关系构建观测方程以矫正基线。这类方法所得到的基线参数一般比较精确，而且外部DEM很容易获取，但现有方法迭代收敛速度慢、所构建的方程复杂度较高，使得该类方法计算量较大、计算效率较低，尤其用于大数据量的时序干涉处理中时，这些缺点会更加凸显。

综合以上分析可以看出，第二类方法相比第一类方法明显复杂，但由于矫正了基线误差，从而在根源上缓解了基线误差带来的影响，可以大大提高干涉处理的精度。但为了更好地满足干涉处理的需求，仍需要普适性强、精度高、复杂度低、收敛速度快的基线估计算法。

发明内容

有鉴于此，本发明的主要目的在于提出基于轨道误差相位基的星载SAR干涉基线矫正方法，以从根源上缓解基线误差造成影响。本发明的方法仅需要外部DEM数据，对干涉数据的质量和特点无特殊要求，且能够在保证所需基线精度（本发明以有效去除基线误差引起的轨道误差相位体现）和算法鲁棒性的前提下，收敛速度快，操作简单。

为达到上述目的，本发明的技术方案是这样实现的：一种基于轨道误差相位基的星载SAR干涉基线矫正方法，包括如下步骤：

步骤1、计算初始基线；

步骤2、提取轨道误差相位基；

步骤3、差分干涉处理得到残差相位；

步骤4、对残差相位、轨道误差相位基、DEM数据进行中值池化；

步骤5、残差相位向轨道误差相位基和DEM投影并矫正基线；

步骤6、在基线参数矫正量均小于阈值的情况下输出矫正后的基线，否则以当前矫正后的结果为基线初值迭代进行步骤3-5。

有益效果：

本发明的方法所提算法对干涉数据的质量和特点无特殊要求，且能够在保证所需基线精度和算法鲁棒性的前提下，收敛速度快，操作简单。

附图说明

图1为基于轨道误差相位基的星载SAR干涉基线矫正方法流程图；

图2为二维平面内基线的不同表示方式示意图；（a）基线表示为基线长度与基线倾角，（b）基线表示为平行基线和垂直基线；

图3 为轨道误差相位基；（a）db _||=160m（b）

（c）db _⊥=2m（d）

；

图4为残差相位；（a）基线矫正前的残差相位，（b）基线矫正后的残差相位；

图5为高程与高程误差；（a）基线矫正前的高程，（b）基线矫正前的高程误差，（c）基线矫正后的高程，（d）基线矫正后的高程误差。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整的描述，显然，所描述的实施例仅为本发明的一部分实施例，而不是全部的实施例，基于本发明中的实施例，本领域的普通技术人员在不付出创造性劳动的前提下所获得的所有其他实施例，都属于本发明的保护范围。

结合本发明的一个实施例，提出一种基于轨道误差相位基的星载SAR干涉基线矫正方法，如图1所示，包括如下步骤：

步骤101：确定基线参数、给出轨道误差相位基的定义并建立轨道误差相位观测方程

以重轨干涉为例，如图2所示，在二维平面中，基线可以以基线长度B(t)和基线倾角α(t)表示（图2的（a）），还可以以平行基线B _||(t,θ)和垂直基线B _⊥(t,θ)表示（图2的（b）），二者的相互转换关系为：

(1)

其中，t是方位时间变量，θ是主卫星的雷达视角。

根据干涉原理，干涉图上任一方位时间t和斜距R处的干涉相位可以近似表示为：

(2)

其中，λ表示雷达波长，R表示主卫星与地面目标的距离（R是θ的函数，后文将直接以变量θ取代变量R）。

如果方位时间为t，雷达视角为θ时基线存在误差dB _||(t,θ)，相应地会产生干涉相位误差：

(3)

对于平地，式（3）符合轨道误差相位的定义，因此将dφ(t,θ)记为φ _orbit(t,θ)。φ _orbit(t,θ)随t和θ变化，将其在整幅干涉图的时间范围和雷达视角范围内关于景中心坐标(t _c,θ _c)进行泰勒展开并忽略高阶项可得：

(4)

其中，上标^·及^··分别代表对时间求一阶导和二阶导。由于单幅SAR图像成像时的视角变化范围很小，成像时间很短，所以该近似在整幅干涉图范围内是合理的。另外，假设基线分量随方位时间线性变化，则有：

(5)

这样，干涉基线误差就可以平行、垂直雷达视线方向坐标系下的景中心基线长度误差

以及景中心基线变化速率误差

四个参数表示。也就是说，景中心基线长度

以及景中心基线变化速率

四个参数足以描述干涉基线。

根据式（5）可得：当仅存在景中心平行基线误差db _||时，产生轨道误差相位分量φ _para（式（6），图3中的（a））；当仅存在景中心平行基线变化速率误差

时，产生轨道误差相位分量φ _parav（式（7），图3中的（b））；当仅存在景中心垂直基线误差db _⊥时，产生轨道误差相位分量φ _perp（式（8），图3中的（c））；当仅存在景中心垂直基线变化速率误差

时，产生轨道误差相位分量φ _perpv（式（9），图3中的（d））。这四个分量正好对应式（5）等式右边的四项，而db _||、

、db _⊥和

也可以看作调节轨道误差相位分量幅度的可变系数。此外，由于式（5）为泰勒展开而得，这四个轨道误差相位分量是相互正交的，具有基函数的性质，因此本发明将其称作“轨道误差相位基”。根据式（6）到（9）可知，轨道误差相位基只与系统参数（t、θ和λ）及基线误差相关而与真实基线无关，这意味着在真实基线未知的情况下就可以获取轨道误差相位基。

(6)

(7)

(8)

(9)

进一步观察式（5）可以发现，轨道误差相位与基线误差有线性关系，改写式（5）形式为：

(10)

以整幅干涉图来看，则有：

(11)

其中：

(12)

符号“^—”代表矩阵转变为列向量的操作，用以满足矩阵相乘时对于矩阵大小的要求。

综合以上分析，基线矫正问题就可以描述为：将轨道误差相位投影到四个轨道误差相位基上，投影系数乘以与轨道误差相位基对应的基线误差参数值就为对应的基线参数的改正数。这样，基线矫正就成为一个线性问题，或者说是一个求解线性系数的问题。由该原理所得的基线矫正方法的流程可以用图1表述。具体描述见以下步骤。

步骤102：计算初始基线（基线一般用基线参数描述，因此“基线”与“基线参数”等价，只有特别需要强调参数的数值时称为“基线参数”）。首先根据已有的基线估计方法估计初始基线

和

。可以选取最简单的根据轨道状态向量信息计算基线的方法，如果没有卫星的精轨数据，则计算所得的基线误差较大，将在以下步骤中用本发明的方法进行矫正。

步骤103：提取轨道误差相位基

轨道误差相位基的提取方法为：分别给基线长度初值

加上平行基线误差db _||和垂直基线误差db _⊥，给基线变化速率初值

加上平行基线变化速率误差

和垂直基线变化速率误差

。再分别用根据加了已知误差的基线计算所得的平地相位减去根据初始基线计算所得的平地相位，即可准确地获取轨道误差相位基φ _para、φ _perp、φ _parav、φ _perpv。

为了保证接下来残差相位投影的可靠性，每一轨道误差相位基的量程（最大值减去最小值）应该相当。因此每一轨道误差相位基对应的基线误差参数的大小建议可以按以下规则设定：

(13)

其中，k为正实数作为放大倍数，Δθ为生成干涉图的雷达视角范围，Δt为生成干涉图的时间范围。以上公式含义均为：当轨道误差相位基的量程为2kπ时对应的基线误差参数的值。根据上式或者直接根据实验观测调整即可对轨道误差相位基对应的基线误差参数值进行设置。

步骤104：差分干涉处理得到残差相位

由于无法准确地提取轨道误差相位，因此只能用残差相位φ _res来近似轨道误差相位。残差相位可由二轨差分干涉处理得到，具体为: 借助经地理编码到雷达坐标系的外部DEM数据（后文简称为DEM数据）以及初始基线模拟干涉相位，从真实的干涉相位中减去模拟干涉相位得到缠绕的残差相位，对其进行解缠后进而得到解缠的残差相位，简称为残差相位。

步骤105：对残差相位、轨道误差相位基、DEM数据进行中值池化

残差相位中不仅包含轨道误差相位还包含大气相位、形变相位、由基线误差造成的地形相位去除残留误差、高程误差相位、解缠误差、噪声相位等。其中，轨道误差相位是空间低频分量，大气湍流分量偶尔也会表现出空间低频特性，除此之外，在不考虑大尺度形变（如地震造成的形变等）的情况下，其余相位都表现出空间高频特性。通过对残差相位池化处理，可以有效突出空间低频分量并削弱空间高频分量。另外，池化后的数据量减小，可以减小后续步骤的计算量。

池化窗口的大小由经验决定，在使得残差相位的高频分量得到有效抑制和尽可能保留轨道误差相位特征之间进行折中。池化后的取值可以取池化窗口内的中值或平均值等，通常中值能比较稳定地反映池化窗口内的数据特性，因此本发明采用中值，这一处理也称为中值池化。为了保持矩阵大小的一致性，对四个轨道误差相位基和DEM也要进行相同尺度的中值池化操作，其中池化后的DEM的作用将在下一步骤中体现。

需要注意的是，可以在干涉图中发生局部较大形变的区域和地形坡度十分陡峭的区域进行掩膜处理，以避免形变和相位解缠误差对解算结果造成影响。

步骤106：残差相位向轨道误差相位基和DEM投影并矫正基线

根据观测方程（11）将轨道误差相位投影到四个轨道误差相位基上。正如前面提到的，轨道误差相位不能直接获取，只能由残差相位近似，因此观测方程（11）应做出相应修正。其次，大气相位的垂直分层分量φ _{atm_stra}近似表现为与高程的线性关系，可以在投影过程中一并考虑。最后，提取的残差相位由于相位解缠的参考点选取不同而相差一个常数。结合以上三点，（11）式可以改写为：

(14)

其中，上标MP代表步骤105中值池化处理后的结果，φ _I是与中值池化后的轨道误差相位基同样大小的元素全为1的矩阵，c ₀为对应于φ _I的系数，φ _noise是φ _res中包含的除-φ _orbit、φ _{atm_stra}和常数相位外的相位成分，均视为噪声相位，h是地形高程。分别用矩阵符号A和C代替式（14）中的设计矩阵和系数向量，则式（14）又可以写为：

(15)

根据最小二乘原理可以得到系数解：

(16)

其中，

是C=[c ₀ c ₁ c ₂ c ₃ c ₄]^T的最小二乘估计结果，此时有：

(17)

V是投影残差相位。

根据式（16）得到系数向量

后，接下来可以对基线参数进行矫正得到矫正后的基线参数

和

：

(18)

步骤107：基线参数矫正量均小于阈值的情况下输出矫正后的基线，否则以当前矫正后的结果为基线初值迭代进行步骤104-106；

由于轨道误差相位是由残差相位近似的，其中包含误差。当基线逼近真实基线的时候，模拟的地形相位越逼近真实地形相位，残差相位也就越逼近轨道误差相位。因此，有必要进行迭代来提高基线矫正结果的精度。具体方法为：首先判断每个基线参数的矫正量（绝对值）是否小于设定的阈值，如果小于则输出矫正后的基线，否则将经步骤106矫正的结果当做基线初值重复进行步骤104到106，直到每个基线参数的矫正量小于阈值时输出当前迭代矫正后的基线作为最终的基线矫正结果。

本发明基于轨道误差相位基的星载SAR干涉基线矫正方法适用性分析，为：

本发明仅需要外部粗精度DEM数据，对干涉数据的质量和特点没有特殊要求。由于以残差相位近似轨道误差相位，因此残差相位中可能包含的大尺度形变（如地震造成的形变）以及大气的低频分量可能会对投影结果造成影响。所以本发明不适用于大尺度形变的场景，而大气的低频分量对本发明的影响在大部分应用中可以接受。如果在大气影响比较严重的情况下对基线精度的需求较高，那么需要在基线估计前进行去大气操作，否则需要对本发明进一步改进以克服大气低频分量的影响。

下面结合具体实施例对本发明的技术方案作进一步详细说明。

实施例1

本发明采用C波段国产卫星高分三号（Gaofen-3，GF-3）在嵩山地区的数据进行实验验证。GF-3是我国自主研发的首颗具备干涉能力的民用全极化星载合成孔径雷达卫星，但GF-3卫星实时定轨精度较差，基线估计问题是应用GF-3数据做干涉处理的一个重要难题。因此，本发明以GF-3数据来说明所提基线矫正算法的有效性。外部DEM采用美国航天飞机雷达地形测绘任务（Shuttle Radar Topography Mission，SRTM）所获取的30m分辨率DEM数据。

实测数据初始基线（见表1，由GF-3数据产品中的轨道状态向量信息计算而得）对应的残差相位如图4中的（a）所示，由于基线误差较大，轨道误差相位比较明显，其余相位成分不明显。利用本发明所提出的算法进行基线矫正。经过4次迭代，矫正后的基线如表1所示，其对应的残差相位如图4中的（b）所示，此时的相位倾斜已经被很好的补偿掉了，因此大气干扰凸显出来，可以看到湍流分量在大气相位中占主导。

表1. 初始基线与矫正后的基线

用初始基线和矫正后的基线生成的高程及高程误差（以SRTM DEM为标准）如图5中的（a）-（d）所示，在基线存在较大误差时，干涉所生成高程完全失真，使得干涉处理失效，而基线矫正后，干涉生成的高程接近SRTM DEM。用高程误差标准差来衡量的高程精度见表2，矫正基线并去除所估计的大气垂直分量影响后，标准差为16.863m，考虑到SRTM DEM本身的绝对高程精度小于16m以及该区域有比较严重的大气影响，这一结果是比较精确的。

表2. 高程误差标准差

尽管上面对本发明说明性的具体实施方式进行了描述，以便于本技术领域的技术人员理解本发明，且应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。

Claims

1.一种基于轨道误差相位基的星载SAR干涉基线矫正方法，其特征在于，包括如下步骤：

步骤1、计算初始基线；

步骤2、提取轨道误差相位基；

步骤3、差分干涉处理得到残差相位；

步骤5、残差相位向轨道误差相位基和DEM投影并矫正基线；

2.根据权利要求1所述的一种基于轨道误差相位基的星载SAR干涉基线矫正方法，其特征在于，在步骤1之前，确定基线参数、给出轨道误差相位基的定义并建立轨道误差相位观测方程，具体如下：

将干涉基线关于场景中心对方位时间和雷达视角进行泰勒展开得到四个基线参数：景中心平行基线、景中心垂直基线、景中心平行基线变化率、景中心垂直基线变化率；根据干涉原理，轨道误差相位与基线误差具有线性关系，因此轨道误差相位被分成分别与四个基线参数相关的四个正交分量，将其称作轨道误差相位基；有了轨道误差相位基的定义后，根据干涉原理就建立起进行基线矫正所依据的轨道误差相位观测方程，即轨道误差相位与基线误差和轨道误差相位基的关系方程。

3.根据权利要求1所述的一种基于轨道误差相位基的星载SAR干涉基线矫正方法，其特征在于，所述步骤1计算初始基线，具体如下：

根据轨道状态向量信息计算初始基线。

4.根据权利要求1所述的一种基于轨道误差相位基的星载SAR干涉基线矫正方法，其特征在于，所述步骤2提取轨道误差相位基，具体如下：

分别给初始基线加上平行基线误差、垂直基线误差、平行基线变化速率误差和垂直基线变化速率误差，再分别用根据加了已知误差的基线计算所得的平地相位减去根据初始基线计算所得的平地相位，即可准确地获取轨道误差相位基。

5.根据权利要求1所述的一种基于轨道误差相位基的星载SAR干涉基线矫正方法，其特征在于，所述步骤3差分干涉处理得到残差相位，具体如下：

借助经地理编码到雷达坐标系的外部DEM数据以及初始基线模拟干涉相位，从真实的干涉相位中减去模拟干涉相位得到缠绕的残差相位，对其进行解缠后进而得到解缠的残差相位，简称为残差相位。

6.根据权利要求1所述的一种基于轨道误差相位基的星载SAR干涉基线矫正方法，其特征在于，所述步骤4对残差相位、轨道误差相位基、DEM数据进行中值池化，具体如下：

对残差相位、轨道误差相位基、DEM数据划分大小相同的网格，并提取网格内的中值，得到比原始数据尺寸更小的数据以代替原始数据，这一操作即为中值池化。

7.根据权利要求1所述的一种基于轨道误差相位基的星载SAR干涉基线矫正方法，其特征在于，所述步骤5残差相位向轨道误差相位基和DEM投影并矫正基线，具体如下：

将由残差相位近似的轨道误差相位投影到四个轨道误差相位基及DEM上，投影系数乘以与轨道误差相位基对应的基线误差参数值即为对应的基线参数的改正数。