一种多场景干涉SAR图像的三维定位方法与流程

本发明涉及三维定位方法技术领域，具体涉及一种多场景干涉SAR图像的三维定位方法。

背景技术：

干涉SAR技术是合成孔径雷达技术与干涉技术的结合，干涉SAR通过两副不同位置处的天线来观测同一场景，利用雷达回波获取并进行数据处理得到场景的SAR图像和干涉相位数据，然后对干涉相位数据进行干涉定标处理，利用干涉定标处理后的干涉相位数据来反演SAR像元的三维位置。

三维定位是干涉SAR图像走向应用的必要前提，所谓三维定位是指获取SAR图像像元所对应地物目标的东西向位置(经度)、南北向位置(纬度)和高程。由于参数误差、航迹误差以及系统非理想特性的存在，导致干涉SAR三维重建几何模型中的参数取值并不准确，需要通过干涉定标来校准各参数，以提高几何定位的精度。因此，干涉定标处理是实现SAR图像精确定位的关键。

干涉定标处理的常规思路是利用三维位置已知的地面控制点(GroundControl Point，GCP)，来闭环修正各干涉参数。为满足干涉SAR的测图作业的要求，需要在少量控制点的条件下实现大区域多场景的干涉定标，采用多场景联合定标方法可实现上述要求。联合定标能保证各场景的三维位置精度，又能保证相邻场景之间的三维位置衔接性。但现有的多场景联合定标方法，通常采用同名点(Tie Point，TP，连接点，又称同名点，即同一地物目标在不同SAR场景中所分别对应的像元)位置信息的传递来实现，实现过程较为繁琐，不便于自动化处理；并且仅面向高程定位，较少涉及三维位置同时定位，现有的关于多场景干涉SAR图像的三维定位方法多是面向定位精度的提升，而鲜有涉及处理流程的优化改进，不便于利用计算机进行自动化、批量化处理。

技术实现要素：

有鉴于此，本发明提供了一种多场景干涉SAR图像的三维定位方法，该定位方法的流程简洁、高效，便于自动化批量化的计算机处理，对三维位置进行定位，通过构建多场景三维联合敏感度矩阵，能够避免现有方法中设计同名点传递路径的繁琐，本发明是通过以下技术方案实现的：

步骤1：确定参与联合定标的干涉SAR场景的总数目I，确定每个场景中GCP点的数目M_i，确定所有GCP点的三维地理位置信息和其在SAR图像中的行列位置信息；

步骤2：在所有参与联合定标的干涉SAR场景中，确定在位置上互相邻接并部分重合的场景组数J；确定第j组互相邻接并部分重合的场景所包含的两个场景的编号a_j、b_j，其中a_j、b_j∈{1，2...I}且a_j≠b_j，j＝1，2...J；在各组场景中选取TP点，确定第j组场景中的TP点数目N_j，确定所有TP点在SAR图像中的行列位置；

步骤3：确定第i个干涉SAR场景的三维位置东向x_i、北向y_i和高程h_i反演的几何表达式

其中i＝1，2...I，自变量u_i，1，u_i，2，...u_i，k...，u_i，K为待标定参数组，u_i，k表示第i个干涉SAR场景中的第k个待定标参数u_i，k，k＝1，2...K，K为待标定参数总个数；

步骤4：构建多场景参数矩阵U，矩阵U的维度为I·K行1列，其构成为U＝[U₁ U₂…U_i…U_I]^T，其中，U_i＝[u_i，1 u_i，2…u_i，k…u_i，K]^T

步骤5：构建关于GCP点的多场景联合敏感度矩阵G，矩阵G的维度为行I·K列，其构成为G＝diag(G₁，G₂，…G_I)，其中表示第i个SAR场景中所有GCP点处三维位置关于待定标参数的敏感度，的矩阵结构如下：

其中a＝x，y，h；

步骤6：构建关于TP点的多场景联合敏感度矩阵T：

其中，矩阵T的维度为行I·K列，和分别表示第a_j和b_j个SAR场景中所有TP点处三维位置关于待定标参数的敏感度，a_j和b_j表征位置上相互邻接并且共同拥有TP点的第j组的两个场景的编号，同一编号的或位于同一列；共同拥有TP点的第j组的两个场景在同一行；

步骤7：构建多场景参数误差矩阵ΔU：矩阵ΔU的维度为I·K行1列，其构成为ΔU＝[ΔU₁ ΔU₂ … ΔU_i … ΔU_I]^T，其中，ΔU_i＝[Δu_i，1 Δu_i，2 … Δu_i，k … Δu_i，K]^T，ΔU为待求解参量，表征待定标参数修正量矩阵；

步骤8：构建多场景位置误差矩阵ΔP：矩阵ΔP的维度为行1列，其构成为ΔP＝[ΔP_G，1 ΔP_G，2 … ΔP_G，i … ΔP_G，I；ΔP_T，1 ΔP_T，2 … ΔP_T，j … ΔP_T，J]^T，其中ΔP_G,i＝[Δx_G,i Δy_G,i Δh_G,i]^T，列向量Δx_G,i Δ_yG,i Δh_G,i分别表示第i个SAR场景中所有GCP点处三维位置误差，即在当前的参数数值下计算出的各GCP点三维位置数值减去其真实三维地理位置数值所得差值；ΔP_T,j＝[Δx_T,j Δy_T,j Δh_T,j]^T，列向量Δx_T,j Δy_T,j Δh_T,j分别表示第j组SAR场景中所有TP点处三维位置误差即同一地物点在2个场景中的三维位置坐标差值，即在当前的参数数值下计算出的TP点在第a_j个场景中的三维位置数值减去该点在第b_j个场景中的三维位置数值所得差值；

步骤9：构建多场景联合敏感度方程：

并求解此方程获得参数修正量矩阵ΔU，计算公式为其中[·]⁺表示矩阵[·]的Moor-Penrose逆；

步骤10：更新矩阵G、T、ΔP：利用步骤9求解出的参数修正量矩阵ΔU，U与ΔU相加，获得更新后的待定标参数矩阵U，依据更新后的U更新矩阵G、矩阵T，然后按照步骤3中的干涉SAR场景的三维位置反演几何表达式计算各GCP、TP点处的三维位置，并通过比较GCP点的真实三维、TP点在一组场景中的三维位置差异，获得更新后的位置误差矩阵ΔP；

步骤11：判断位置误差矩阵ΔP是否满足精度要求或满足收敛条件，若满足精度要求或已收敛，则将步骤10计算出的各GCP、TP点处的三维位置作为各干涉SAR场景中各像元处的三维位置，完成多场景干涉SAR图像的三维几何定位，否则，转步骤9，利用更新后的G、T、ΔP求解步骤9中的方程(1)。

其中，所述步骤3中的待定标参数包括交轨基线长度、顺轨基线长度、基线倾角、干涉相位偏置、回波延时、天线相位中心初始位置和多普勒中心频率。

进一步地，所述步骤10中的收敛条件为为最近2次迭代获得的ΔP值的差异小于系统要求的三维位置精度的取值。

有益效果：

本发明提供了一种多场景干涉SAR图像的三维定位方法，通过同时对多个场景的SAR数据进行数学建模表达，实现多场景同时定位，优化了面向多场景干涉SAR图像地三维定位流程，便于自动化与批量化的计算机处理；

本发明的定位方法包含东向位置(经度)、北向位置(纬度)和高程的三维位置同时定位，通过同时对多个场景的SAR数据进行三维位置干涉定标，实现三维定位；

本发明通过构建多场景三维联合敏感度矩阵，能够使多个干涉SAR场景的所有待定标参数同时进行迭代计算，保证各场景的定位的精度以及场景邻接处的位置衔接性，并且避免了现有方法中设计同名点传递路径的繁琐。

附图说明

图1为多场景干涉SAR图像的三维定位方法流程图。

图2为有两个场景，三个邻接位置时，联合定标场景邻接分布图。

具体实施方式

下面结合附图并举实施例，对本发明进行详细描述。

本发明提供了一种多场景干涉SAR图像的三维定位方法，所述方法包括：

步骤1：确定参与联合定标的干涉SAR场景的总数目I，确定第i个场景中GCP点的数目M_i，i＝1，2...I，确定所有GCP点的三维地理位置信息和其在SAR图像中的行列位置信息。

步骤2：在所有参与联合定标的干涉SAR场景中，确定在位置上互相邻接并部分重合的场景组数J，每组包括2个场景；确定第j组场景所包含的2个场景的编号a_j、b_j，其中a_j、b_j∈{1，2...I}，j＝1，2...J；在各组场景中选取TP点，确定第j组场景中的TP点数目N_j，确定所有TP点在SAR图像中的行列位置信息。

步骤3：根据系统的照射关系来确定干涉SAR的三维位置反演的几何表达式。

第i个干涉SAR场景的三维位置东向x_i、北向y_i和高程h_i反演的几何表达式为

其中自变量(u_i，1，u_i，2，…u_i，k…，u_i，K)为待标定参数组，u_i，k表示第i个干涉SAR场景中的第k(k＝1，2…K)个待定标参数u_i，k。从第1个到第K个待定参数的具体物理量依据实际需求具体确定，干涉SAR待定标参数u_i，k的物理量通常包括交轨基线长度、顺轨基线长度、基线倾角、干涉相位偏置、回波延时、天线相位中心初始位置和多普勒中心频率。

步骤4：构建多场景参数矩阵U，矩阵U的维度为I·K行1列，其构成为U＝[U₁ U₂ … U_i … U_I]^T，其中，U_i＝[u_i，1u_i，2 … u_i，k … u_i，K]^T

步骤5：构建关于GCP点的多场景联合敏感度矩阵G，矩阵G的维度为行I·K列，其构成为G＝diag(G₁，G₂，…G_i…G_I)。其中表示第i个SAR场景中所有GCP点处三维位置关于待定标参数的敏感度。具体地，

其中分别表示第i(i＝1，2...I)个干涉SAR场景中第m(m＝1，2...M_i)个GCP点位置处三维位置x_i、y_i、h_i关于第k(k＝1，2...K)个待定标参数u_i，k的敏感度，通过步骤3中的反演的几何表达式以及第m个GCP点对应的干涉SAR场景参数的初值获得。例如矩阵元素为将第一个GCP点对应的相应物理参数值代入函数x(u_i，1，u_i，2，…u_i，k…，u_i，K)对u_i，l求偏导后得到的偏导函数得到的数值。

步骤6：构建关于TP点的多场景联合敏感度矩阵T。

矩阵T的维度为行I·K列，其矩阵结构如下

其中，矩阵T的维度为行I·K列，和分别表示第a_j和b_j个SAR场景中所有TP点处三维位置关于待定标参数的敏感度，a_j和b_j表征位置上相互邻接并且共同拥有TP点的第j组的两个场景的编号，同一编号的或位于同一列；共同拥有TP点的第j组的两个场景在同一行；

如图2所示的联合定标场景邻接分布图，图中画出的场景编号为No1、No2、No3、NoI…，第1组的两个场景的编号为1和2，第2组的两个场景的编号为1和3，第3组的两个场景的编号为2和I，按编号确定T₁、T₂、T₃、T_I所属的第1、2、3、I列，该分布图对应的矩阵T为

对于具体地，

其中和分别表示第a_j、b_j个干涉SAR场景中第n个TP点位置处三维位置x_i、y_i、h_i关于第k(k＝1，2...K)个待定标参数和的敏感度，和的表达式通过步骤3中的公式以及第n个GCP点干涉SAR场景参数的初值获得。

步骤7：构建多场景参数误差矩阵ΔU。矩阵ΔU的维度为I·K行1列，其构成为ΔU＝[ΔU₁ ΔU₂ … ΔU_i … ΔU_I]^T。其中，ΔU_i＝[Δu_i，1 Δu_i，2 … Δu_i，k … Δu_i，K]^T。ΔU表征待定标参数修正量矩阵，为待求解参量。

步骤8：构建多场景位置误差矩阵ΔP。矩阵ΔP的维度为行1列，其构成为ΔP＝[ΔP_G，1 ΔP_G，2 … ΔP_G，i … ΔP_G，I；ΔP_T，1 ΔP_T，2 … ΔP_T，j … ΔP_T，J]^T。

其中ΔP_G，i＝[Δx_G，i Δy_G，i Δh_G，i]^T，

列向量Δx_G,iΔy_G,iΔh_G,i分别表示第i个SAR场景中所有GCP点处三维位置误差，即在当前的参数数值下计算出的各GCP点三维位置数值减去其真实三维地理位置数值所得差值，列向量Δx_T,jΔy_T,jΔh_T,j分别表示第j组SAR场景中所有TP点处三维位置误差，即在当前的参数数值下计算出的TP点在第a_j个场景中的三维位置数值减去该点在第b_j个场景中的三维位置数值所得差值。

步骤9：构建多场景联合敏感度方程并求解此方程获得参数修正量ΔU，计算公式为其中[·]⁺表示矩阵[·]的Moor-Penrose逆。

步骤10：更新矩阵G、T、ΔP：利用步骤9求解出的参数修正量矩阵ΔU，U与ΔU相加，获得更新后的待定标参数矩阵U，依据更新后的U更新矩阵G、矩阵T，然后按照步骤3中的干涉SAR场景的三维位置反演几何表达式计算各GCP、TP点处的三维位置，并通过比较GCP点的真实三维、TP点在一组场景中的三维位置差异，获得更新后的位置误差矩阵ΔP；

步骤11：判断位置误差矩阵ΔP是否满足精度要求或满足收敛条件，若满足精度要求或已收敛，则将步骤10计算出的各GCP、TP点处的三维位置作为各干涉SAR场景中各像元处的三维位置，完成多场景干涉SAR图像的三维几何定位，否则，转步骤9，利用更新后的G、T、ΔP求解步骤9中的方程(1)。

所述收敛条件为最近2次迭代获得的ΔP值的差异小于设定要求，系统设定要求通常结合制图的三维位置精度需求进行设定，其取值即为三维位置精度的取值。

综上所述，以上仅为本发明的较佳实施例而已，并非用于限定本发明的保护范围。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。