一种无控制点的干涉SAR基线矢量估计方法与流程

本发明涉及雷达信号处理技术领域，具体涉及一种无控制点的干涉sar基线矢量估计方法。

背景技术：

合成孔径雷达干涉测量(interferometricsyntheticapertureradar,insar)技术是随着信息技术、摄影测量技术、数字信号处理技术等相关技术的发展而迅速发展起来的一种高精度对地观测新技术。利用干涉sar技术快速获取高精度数字高程模型(digitalelevationmodel,dem)是目前干涉sar技术的主要应用之一。干涉sar获取dem的基本原理是利用合成孔径雷达(syntheticapertureradar,sar)系统的两副天线(或一副天线重复观测)，来获取同一地区具有一定视角差的两幅具有相干性的单视复数sar图像，并根据其干涉相位信息来提取地表的高程信息和重建dem。

在干涉sar进行dem反演时，基线矢量是影响高程dem重建精度的重要因素。为了降低基线矢量误差的影响，必须对基线矢量进行标定。目前，常常采用的是基于地面控制点的干涉定标方法。但是，这种方法要求在每景干涉像对中实地布放足量且分布合理的地面控制点，工作量大、作业效率低，而且野外布放地面控制点受地形条件限制，存在某些测区如荒山、沼泽等难以实现地面控制点的布放。因此，在无地面控制点时，需要利用信号处理技术实现对干涉基线矢量的精确估计。

技术实现要素：

本发明的目的在于提供一种无控制点的干涉sar基线矢量估计方法，实现了无控制点情况下的干涉sar基线矢量估计，能够有效提高干涉sar大区域测绘的作业效率。

为了达到上述目的，本发明通过以下技术方案实现：

一种无控制点的干涉sar基线矢量估计方法，其特征是，包含以下步骤：

s1、根据干涉sar几何关系，建立基线矢量误差与地物目标高程重建误差的关系；

s2、分别建立基线长度和基线倾角随着方位向时间的数学模型，并代入根据基线矢量误差与地物目标高程重建误差的关系中；

s3、在两次飞行的重叠区域，基于地物目标真实高度不变建立两次干涉测高过程中基线矢量误差之间的关系，结合步骤s2的结果获得关于基线矢量误差的线性方程组；

s4、通过引入权值来区分重叠区域中不同点处的相位质量差异，该权值与相干系数和位置分布相关；

s5、根据步骤s3和步骤s4的结果利用加权最小二乘法求解线性方程组，获得估计的基线矢量误差；

s6、在原始基线矢量基础上叠加估计的基线矢量误差即为最终的基线矢量，利用其进行dem高程的高精度重建。

上述的无控制点的干涉sar基线矢量估计方法，其中，所述的步骤s1具体包含：

根据干涉sar的几何关系，设h为主天线a1的高度，主天线a1和副天线a2之间的基线长度为b，基线倾角为α；点p为干涉测高区域中的一点，其高程为h；主天线a1与点p之间的斜距为r，且其与垂直方向的夹角为θ；

由干涉测高原理得到：

(r+δr)²＝r²+b²+2rbsin(θ-α)(1)

h＝h-rcosθ(3)

式中，λ为波长，δr为主天线a1和副天线a2与点p之间斜距的差值；δφ为主天线a1和副天线a2与点p之间的斜距差引起的相位差；

联合求解上述三个方程得到地物目标的高程；

对式(3)求基线矢量的偏导数，可得：

建立基线矢量误差与高程重建误差之间的关系为：

其中，hobs为干涉sar反演得到的目标高程。

上述的无控制点的干涉sar基线矢量估计方法，其中，所述的步骤s2具体包含：

建立基线长度和基线倾角随方位向实践的数学模型：

δb(t)＝δb0+mt(7)

δα(t)＝δα0+nt(8)

其中，δb0和δα0分别为基线长度常数误差和基线倾角常数误差，t为方位向时间，δb(t)和δα(t)分别为基线长度和基线倾角随着时间变化的误差，m和n分别为基线长度误差和基线倾角误差随着时间变化的系数值；

根据式(7)、(8)将式(6)改写为：

上述的无控制点的干涉sar基线矢量估计方法，其中：

对于某一测绘区域，利用干涉sar对其进行作业两次，根据式(9)分别获得干涉基线矢量误差与高程重建误差之间的关系，即：

将地物目标的高程认为是不变的，由式(11)和(12)得到：

对重叠区域内k个同名点分别按照式(12)建立方程，可以得到如下线性方程组：

δh＝ax(13)

其中，

x＝[δb01m1δα01n1δb02m2δα02n2]^t(15)

上述的无控制点的干涉sar基线矢量估计方法，其中，所述的步骤s4具体包含：

假设权值矩阵为p，则其可表示为

p＝diag[p1,p2,…,pk](17)

鉴于各点处的相位噪声各不相同，各同名点处的式(12)所示高程约束方程的误差各不相同，其误差具体表达如下：

式中，l1和l2分别为同名点在两次干涉测绘时干涉处理的多视视数，φ1和φ2分别为同名点在两次干涉测绘时干涉处理后的绝对干涉相位；

因此，第k个同名点处的权值可设计为的反比例函数，即：

上述的无控制点的干涉sar基线矢量估计方法，其中，所述的步骤s5具体包含：

在获得式(13)所示的线性方程组和式(17)所示的权值矩阵后，利用加权最小二乘法求解，可得：

x＝(a^tpa)^-1a^tpδh(20)；

估计得到基线矢量误差后，根据以下两式分别获得两次干涉测高时的基线长度和基线倾角估计值：

b(t)＝b0+δb0+mt(21)

α(t)＝α0+δα0+nt(22)

其中，b0和α0分别代表基线长度和基线倾角的初始值。

本发明与现有技术相比具有以下优点：

(1)不需要在进行测绘时布放地面控制点，大大降低测绘工作量，提高干涉sar进行测绘的作业效率；

(2)使得在对某些如荒山、沼泽等难以布放地面控制点的测区进行测绘时仍然能够精确估计基线矢量，获得高精度的地物目标高程；

(3)基线矢量沿着方位向时变，能够更加准确的提高dem高程的重建精度。

附图说明

图1为本发明的方法流程图；

图2为本发明中干涉sar的几何关系图。

具体实施方式

以下结合附图，通过详细说明一个较佳的具体实施例，对本发明做进一步阐述。

如图1示出本发明方法的流程图，本发明所提出的一种无控制点的干涉sar基线矢量估计方法，包括以下步骤：

s1、根据干涉sar几何关系，建立基线矢量误差与地物目标高程重建误差的关系；

s2、分别建立基线长度和基线倾角随着方位向时间的数学模型，并代入根据基线矢量误差与地物目标高程重建误差的关系中；

s3、在两次飞行的重叠区域，基于地物目标真实高度不变建立两次干涉测高过程中基线矢量误差之间的关系，结合步骤s2的结果获得关于基线矢量误差的线性方程组；

s4、通过引入权值来区分重叠区域中不同点处的相位质量差异，该权值与相干系数和位置分布相关；

s5、根据步骤s3和步骤s4的结果利用加权最小二乘法求解线性方程组，获得估计的基线矢量误差；

s6、在原始基线矢量基础上叠加估计的基线矢量误差即为最终的基线矢量，利用其进行dem高程的高精度重建。

所述步骤s1具体包含：

如图2所示，为干涉sar的几何关系图，根据干涉sar的几何关系，h为主天线a1的高度；主天线和副天线都是干涉sar上的，主天线a1和副天线a2之间的基线长度为b，基线倾角为α；点p为干涉测高区域中的一点，其高程为h；主天线与点之间的斜距为r，且其与垂直方向的夹角为θ。由干涉测高原理可知：

(r+δr)²＝r²+b²+2rbsin(θ-α)(1)

h＝h-rcosθ(3)

因此，联合求解上述三个方程即可得到地物目标的高程，式中，λ为波长，δr为主天线a1和副天线a2与点p之间斜距的差值；δφ为主天线a1和副天线a2与点p之间的斜距差引起的相位差。

为了建立基线矢量误差与高程重建误差之间的关系，对式(3)求基线矢量的偏导数，可得：

因此，基线矢量误差与高程重建误差之间的关系可以表示为：

其中，hobs为干涉sar反演得到的目标高程。

所述的步骤s2具体包含：

在传统的基于地面控制点的干涉基线矢量估计中，通常假设基线矢量在干涉测高时间内是不变的。而实际上，干涉基线矢量受到天线杆姿态抖动等因素的影响，存在着缓慢的变化。因此，为了更加准确的对基线矢量进行估计，特对其建立以下模型：

δb(t)＝δb0+mt(7)

δα(t)＝δα0+nt(8)

其中，δb0和δα0分别为基线长度常数误差和基线倾角常数误差，t为方位向时间，δb(t)和δα(t)分别为基线长度和基线倾角随着时间变化的误差，m和n分别为基线长度误差和基线倾角误差随着时间变化的系数值。

此时，式(6)可改写为

所述的步骤s3具体包含：

对于某一测绘区域，如果利用干涉sar对其进行作业两次，则根据式(9)分别获得干涉基线矢量误差与高程重建误差之间的关系，即：

其中，下标1和2分别代表第1次干涉测高和第2次干涉测高时的系统参数。

由于两次干涉测高的时间间隔很短，地物目标的高程可以认为是不变的，因此，由式(11)和(12)可得：

因此，对重叠区域内k个同名点分别按照式(12)建立方程，可以得到如下线性方程组：

δh＝ax(13)

其中，

x＝[δb01m1δα01n1δb02m2δα02n2]^t(15)

所述的步骤s4具体包含：

假设权值矩阵为p，则其可表示为：

p＝diag[p1,p2,…,pk](17)

鉴于各点处的相位噪声各不相同，各同名点处的高程约束方程(如式(12)所示)的误差各不相同，其误差可具体表达如下：

式中，l1和l2分别为同名点在两次干涉测绘时干涉处理的多视视数，φ1和φ2分别为同名点在两次干涉测绘时干涉处理后的绝对干涉相位；

因此，第k个同名点处的权值可设计为的反比例函数，即：

所述的步骤s5具体包含：

在获得式(13)所示的线性方程组和式(17)所示的权值矩阵后，利用加权最小二乘法求解，可得：

x＝(a^tpa)^-1a^tpδh(20)

上式(20)估计得到的值再根据式(7)和式(8)计算可以得到基线矢量误差，在估计得到基线矢量误差后，可以根据以下两式分别获得两次干涉测高时的基线长度和基线倾角估计值：

b(t)＝b0+δb0+mt(21)

α(t)＝α0+δα0+nt(22)

其中，b0和α0分别代表基线长度和基线倾角的初始值。

尽管本发明的内容已经通过上述优选实施例作了详细介绍，但应当认识到上述的描述不应被认为是对本发明的限制。在本领域技术人员阅读了上述内容后，对于本发明的多种修改和替代都将是显而易见的。因此，本发明的保护范围应由所附的权利要求来限定。