CN106885585B - 一种基于光束法平差的星载摄影测量系统一体化检校方法 - Google Patents

Abstract

本发明公开了一种基于光束法平差的星载摄影测量系统一体化检校方法，该方法利用卫星摄影测量系统自身较稳定的内部结构几何关系，构建系统所包含的星载对地相机、卫星定位观测设备、卫星姿态观测设备之间的几何关联，得出卫星多种原始观测值的误差方程，从而建立卫星摄影测量系统一体化检校平差模型，实现系统结构几何参数的高精度估求，从而有效补偿系统参数误差在对地目标定位中的影响。

Description

一种基于光束法平差的星载摄影测量系统一体化检校方法

技术领域

本发明涉及光学测绘卫星数据几何处理技术领域，特别涉及一种基于光束法平差的星载摄影测量系统一体化检校方法，应用于卫星摄影测量系统的在轨几何检校及对地目标定位。

背景技术

近年来，随着我国光学遥感卫星的快速发展，国产卫星影像数据量大幅增长，卫星影像几何处理成为卫星数据深度挖掘、提取空间地理信息的重要基础，而卫星摄影测量技术是卫星影像高精度几何处理的关键支撑技术。

卫星摄影测量技术的基本任务之一，是通过星载摄影测量系统的精确几何检校，消除在轨系统误差影响，使卫星对地目标定位精度达到卫星系统的应用潜力。由于星载系统对定位精度影响最大的是角元素系统误差，目前常用的卫星摄影测量检校方法是在建立星载系统成像几何模型的基础上，利用适量控制点对角元素系统误差进行检校，如姿态角常差检校法、相机偏置角检校法、CCD侧视角检校法、探元指向角法等。

尽管以上检校方法在实际应用中取得了较好的效果，但仍存在不足：(1)现有模型主要针对系统的角元素误差进行校正，往往忽略其他几何误差源，而这些误差源的影响在高精度应用中不可忽视；(2)现有模型依赖像点坐标观测值，检校参数与相机内外方位元素具有强相关性，导致求解参数的法方程病态，参数解算的精度较差；(3)星载系统的其他观测值如卫星定位观测值和卫星姿态观测值在检校模型中不是作为直接观测值使用，而是作为模型参数初值，既造成了不必要的精度损失，又不能与像点观测值形成有效互补，数据价值未得到充分利用。

目前在实际工程应用中，为取得精确的星载摄影测量系统检校结果，要求在数据处理时提供大量高精度控制点，检校作业过程繁复，因此亟待寻求一种既准确又简单可靠的系统检校方法。

发明内容

本发明提出一种基于光束法平差的星载摄影测量系统一体化检校方法，可充分地消除由于现有技术的限制和缺陷导致的一个或多个问题。

本发明另外的优点、目的和特性，一部分将在下面的说明书中得到阐明，而另一部分对于本领域的普通技术人员通过对下面的说明的考察将是明显的或从本发明的实施中学到。通过在文字的说明书和权利要求书及附图中特别地指出的结构可实现和获得本发明目的和优点。

本发明提供了一种基于光束法平差的星载摄影测量系统一体化检校方法，其特征在于，所述方法包括以下步骤：

步骤1，获取基于惯性系的以归一化四元数形式表示的卫星姿态数据；

步骤2，将通过步骤1获得的所述卫星姿态数据的表达形式由归一化四元数的形式转换为欧拉角形式；

步骤3，将通过步骤2获得的欧拉角形式的卫星姿态数据从惯性坐标系转换至地球固定直角坐标系；

步骤4，获取基于地球固定直角坐标系的卫星定位观测值；

步骤5，构建卫星摄影测量一体化检校平差模型；

步骤6，基于步骤5建立的星载摄影测量系统一体化检校平差模型，按最小二乘原理求检校参数估计值。

本发明针对目前卫星摄影测量检校方法存在的不足，通过引入卫星系统结构几何参数，建立星载摄影测量系统一体化检校光束法平差模型。结合后续给出的具体实施例可以看出，本方法以单个影像的光线束为处理对象，将原始像点坐标、卫星定位观测值、卫星姿态观测值归入一个统一的平差模型进行一体化平差，同时实现对系统结构几何参数及内方位元素的精确求解。其主要优点在于：

(1)基于严密的光束法平差理论，直接对卫星摄影测量系统的各类观测值进行严格平差，避免中间环节的误差因素影响，能够充分地利用卫星原始观测成果。

(2)通过引入卫星定位观测误差方程和姿态观测方程，可使卫星摄影测量的平差几何条件显著增强，检校结果精度更高。

(3)一体化平差可放宽对地面控制点数量和分布的要求，有利于提高卫星摄影测量系统检校的作业效率。

附图说明

图1为根据本发明实施例的、基于光束法平差的星载摄影测量系统一体化检校方法的流程图。

具体实施方式

下面对本发明进行更全面的描述，其中说明本发明的示例性实施例。

如图1所示，本发明所提出的基于光束法平差的星载摄影测量系统一体化检校方法，所述方法具体包括以下步骤：

步骤1，获取基于惯性系的以归一化四元数形式表示的卫星姿态数据；

步骤1具体包括：通过抽取测绘卫星下传的对应于某个检校区域的卫星姿态数据，得到一组基于惯性系的离散的以归一化四元数形式表示的卫星姿态数据。

根据本发明的实施例，可以抽取光学测绘卫星对应任意测区的卫星姿态数据；优选的，本发明的实施例抽取的是资源三号02星下传的卫星姿态数据，但并不限于此。随机抽取一轨卫星下传的姿态数据，该轨数据为发射初期的一轨姿态数据，获取时间为2016年5月31日，数据包含1202个离散的以归一化四元数形式表示的姿态数据值。

资源三号02星下传的该轨姿态数据具体如表1所示：

表1

表1中序号列表示姿态数据的次序，其对应于卫星姿态获取的时间，两个序号之间的时间间隔为0.25s，即，以0.25s的时间间隔获取卫星姿态数据。

所述卫星姿态数据优选的用归一化四元数表示为：

其中，q为一个卫星姿态值，

q₁、q₂和q₃表示归一化四元数的虚部，也是表1中对应的q₁、q₂和q₃列数据，其中，一组q₁、q₂和q₃值描述了卫星的一种姿态，q₄表示归一化四元数的实部，可由虚部计算得到，计算公式为

在实际操作中，为了避免缀余，可以不下传实部部分。

步骤2，将通过步骤1获得的所述卫星姿态数据的表达形式由归一化四元数的形式转换为欧拉角形式；

多数遥感卫星采用四元数的形式进行姿态数据的描述，四元数形式描述的姿态在物理意义表达和直观性方面存在缺陷，而欧拉角形式描述的姿态直接反映卫星三轴姿态变化量，物理意义明显，在摄影测量中广泛应用；为了便于检校平差模型的建立，本发明首先将姿态描述形式统一到欧拉角。

根据本发明的优选实施例，根据转换公式(11)将卫星姿态数据的表达形式由归一化四元数的形式转换为欧拉角形式：

其中，q₁、q₂和q₃表示归一化四元数的虚部，也是表1中对应的q₁、q₂和q₃列数据，q₄表示归一化四元数的实部，可由虚部计算得到，计算公式为

Roll、Pitch和Yaw分别为以欧拉角形式表示的横滚角、俯仰角和偏航角，A_eular表示转换后的三轴欧拉角矩阵。

步骤3，将通过步骤2获得的欧拉角形式的卫星姿态数据从惯性坐标系转换至地球固定直角坐标系；

资源三号卫星下传的姿态数据的坐标系为J2000惯性坐标系，为了实施对地目标定位，需要将其转换到地球固定直角坐标系(卫星摄影测量一般采用地球固定直角坐标系作为物方坐标系)。

设某摄影光线在卫星本体坐标系下的方向矢量为

该摄影光线在惯性坐标系下的方向矢量为

式中R_GCRS为卫星在惯性坐标系下的旋转矩阵。设卫星在惯性坐标系中的三轴姿态角为A_{eular_GCRS}＝[RollPitchYaw]^T，则R_GCRS中的元素就是由三个姿态角Roll、Pitch、Yaw计算的方向余弦。

同理，该摄影光线在地球固定直角坐标系下的方向矢量为：

式中R_ITRS为卫星在地球固定直角坐标系下的旋转矩阵。设卫星在地球固定直角坐标系中的三轴姿态角为A_ITRS＝[phiomegakappa]^T，则R_ITRS中的元素就是由三个姿态角phi、omega、kappa计算的方向余弦。

顾及极移、岁差、章动等现象的存在，同一向量在惯性坐标系和地球固定直角坐标系之间的转换关系可表示为：

其中，

和

分别为同一向量在地球固定直角坐标系和惯性坐标系中的坐标。W(t)，R(t)和Q(t)分别对应极移，自转和岁差章动转换矩阵。由此可以得出，

由上式可得到

R_ITRS＝H(t)·R_GCRS (12)

根据上式，将卫星在惯性坐标系中的三轴姿态Roll、Pitch、Yaw转换为在地球固定直角坐标系下的姿态角phi、omega、kappa。

即，步骤3中根据下述公式将卫星姿态数据从惯性坐标系转换至地球固定直角坐标系：

R_ITRS＝H(t)·R_GCRS

其中，R_ITRS为卫星在地球固定直角坐标系下的旋转矩阵，R_ITRS中的元素是由基于地球固定直角坐标系的三个姿态角phi、omega、kappa计算的方向余弦；R_GCRS为卫星在惯性坐标系下的旋转矩阵，R_GCRS中的元素是由基于惯性坐标系的三个姿态角Roll、Pitch、Yaw计算的方向余弦；H(t)＝W(t)·R(t)·Q(t)，W(t)，R(t)和Q(t)分别对应极移，自转和岁差章动转换矩阵。

步骤4，获取基于地球固定直角坐标系的卫星定位观测值；

根据本发明的实施例，卫星定位观测值和卫星姿态观测值是由星上同步下传的数据，通过抽取与卫星姿态数据同步获取的卫星定位观测数据来获得所述卫星定位观测值。

以该轨姿态数据同步的卫星定位数据为例，该轨数据包含301个离散的卫星定位数据值，如表2所示。

表2

表2中序号列表示定位数据的次序，其对应于卫星空间位置获取的时间，两个序号之间的时间间隔为1.0s，即，以1.0s的时间间隔获取卫星定位数据。Xg，Yg,Zg即卫星在地固直角坐标系中的三维坐标。由于进行一体化检校的空间基准是地固直角坐标系，因此无需再对Xg，Yg,Zg进行转换处理。

步骤5，构建星载摄影测量一体化检校平差模型；

根据本发明的实施例，一体化检校平差模型由像点观测方程、卫星姿态观测方程和卫星定位观测方程三类方程组成，这三类方程分别建立了像点观测值(通过一般遥感影像量测方式得到)、卫星姿态观测值(通过步骤2和3得到)、卫星定位观测值(通过步骤4得到)与待检校参数的函数关系。

步骤5具体包括下述子步骤：

步骤5.1，建立像点观测方程；

建立像点观测方程就是建立像点观测值与待求参数的函数关系。由于星载对地相机一般为线阵推扫式相机，所成影像为动态影像，即每个成像历元对应不同的摄影外方位元素，本发明采用定向片内插模型，通过若干合理分布的定向片外方位元素以内插出任意时刻的相机外方位元素。

像点观测方程即公式(1)给出了像点坐标与成像时刻成像相机内外方位元素的函数关系：

其中，(x,y)为像点坐标(其中_y坐标为影像列坐标，x坐标为影像行坐标(因每一行对应一个成像时刻，故有x＝0))；(x₀,y₀)为相机的像主点坐标，即摄影中心在像面投影点相对于行中心的偏移量；f为成像相机主距，即摄影中心到像面的距离；由于采用地球固定直角坐标系作为物方坐标系，(X_S',Y_S',Z_S')为成像相机外方位线元素，即，成像相机摄影中心在地固坐标系中的空间位置；a_j,b_j,c_j(j＝1,2,3)为成像相机外方位角元素构成的旋转矩阵方向余弦值；(X_T,Y_T,Z_T)为该像点对应的地面点在地固坐标系下的空间坐标。

公式(13)表示成像时刻成像相机摄影中心坐标X_S',Y_S',Z_S'与基准相机摄影中心坐标X_S,Y_S,Z_S的数学关系：

式中，u',v',w'为成像相机摄影中心相对于基准相机摄影中心的空间偏移向量。R_B为成像时刻的基准相机三轴姿态角

ω,κ构成的旋转矩阵。

公式(14)表示成像时刻成像相机三轴姿态角

ω',κ'与基准相机三轴姿态角

ω,κ的数学关系：

式中，R_A为成像相机外方位角元素

ω',κ'构成的旋转矩阵；R_AB为成像相机与基准相机的相对旋转角

ω_I,κ_I构成的旋转矩阵。

公式(15)给出了成像时刻基准相机外方位元素与定向片时刻基准相机外方位元素的函数关系：

其中，

W_i为定向片i对成像时刻t的拉格朗日内插系数；t_i、t_k分别表示第i张定向片和第k张定向片的时刻；n为大于1的整数，优选的，n＝3；X_S,Y_S,Z_S为成像时刻的基准相机摄影中心坐标；

ω,κ为成像时刻基准相机的三轴姿态角。X_Si,Y_Si,Z_Si为第i张定向片时刻的基准相机摄影中心坐标，

ω_i,κ_i为第i张定向片时刻的基准相机三轴姿态角。

结合式(1)、式(13)、式(14)、式(15)可以得出像点坐标与定向片时刻基准相机外方位元素、待检校参数(包括成像相机内方位元素、成像相机相对于基准相机的摄影中心空间偏移量增量、成像相机相对于基准相机的相对旋转角等)及地面点坐标的函数关系，由于该函数是非线性函数，需进行线性化求解。

将式(1)进行线性化，得到像点坐标误差方程为：

V_X＝A_X1t₁+A_X2t₂+Ii+Jj+Bx+Cc-L_X (2)

其中：V_X为像点坐标误差，t₁为基准相机外方位线元素增量；t₂为基准相机外方位角元素增量；x为地面点物方坐标增量；c为相机内方位元素参数增量；A_X1表示像点坐标对基准相机外方位线元素的偏导系数矩阵；A_X2表示像点坐标对基准相机外方位角元素的偏导系数矩阵；I表示像点坐标对相机间摄影中心空间偏移量的偏导系数矩阵；J表示像点坐标对相机间相对旋转角的偏导系数矩阵；B表示像点坐标对地面点物方坐标的偏导系数矩阵；B表示像点坐标对地面点物方坐标的偏导系数矩阵；C表示像点坐标对相机内方位元素的偏导系数矩阵；L_X为像点坐标观测值的残差向量；对于前(后)视相机所成像点，i为前(后)视相机相对于正视相机的摄影中心空间偏移量增量；j为前(后)视相机本体相对于正视相机本体的相对旋转角增量；对于正视相机所成像点，有i,j＝0。

上述线性化的具体推导过程为：设像点坐标的精确值为X，其近似值为X⁰，参数的精确值为Y，其近似值为Y⁰，像点观测方程可简写为

X＝F₁(Y)

像点坐标精确值X可看作像点坐标观测值X⁰(观测值作为近似值)与像点坐标观测误差V_X之和，参数的精确值Y看作参数近似值Y⁰与参数增量y之和，对函数F₁(Y)按泰勒级数展开至一次项，有

X⁰+V_X＝F₁(Y⁰)+F₁'y

上式可写为：

V_X＝F₁'y-L_X

其中，F'为像点坐标对参数的偏导系数矩阵，L_X＝X⁰-F₁(Y⁰)。对于像点观测方程，参数为成像相机内方位元素、定向片时刻基准相机外方位元素、成像相机相对于基准相机的摄影中心空间偏移量、成像相机相对于基准相机的相对旋转角以及地面点坐标。

y＝[t₁ t₂ i j x c]^T为参数增量向量。其中：t₁为基准相机外方位线元素增量、t₂为基准相机外方位角元素增量、x为地面点物方坐标增量、c为成像相机内方位元素参数增量、i为成像相机相对于基准相机的摄影中心空间偏移量增量、j为成像相机相对于基准相机的相对旋转角增量。

像点坐标(x,y)对基准相机外方位线元素的偏导系数矩阵子块

对基准相机外方位角元素的偏导系数矩阵子块

对摄影中心空间偏移量的偏导系数矩阵子块

对相对旋转角的偏导系数矩阵子块

对地面点坐标的偏导系数矩阵子块

对相机内方位元素偏导系数矩阵子块

步骤5.2，建立卫星定位观测方程：

建立卫星定位观测方程就是建立步骤4获得的卫星定位观测值Xg，Yg,Zg与待求参数的函数关系。设O、O'分别为基准相机摄影中心与卫星定位观测中心，它们之间的空间偏移向量可分解为平行于卫星摄影测量系统坐标系三轴的三个分量[u_g v_g w_g]^T，即卫星定位观测中心在卫星摄影测量系统坐标系中的三维坐标，从而得到卫星定位观测方程如式(16)所示，该式表示了基准相机摄影中心坐标[X_S Y_S Z_S]^T和卫星定位观测值[X_g Y_g Z_g]^T之间的数学关系：

式中，[X_g Y_g Z_g]^T表示卫星定位观测值，[X_S Y_S Z_S]^T表示基准相机摄影中心坐标，[u_g v_g w_g]^T表示卫星定位观测中心在卫星摄影测量系统坐标系中的三维坐标，R_B表示基准相机外方位角元素

构成的旋转矩阵。

在建立卫星定位接收机与基准相机的几何关系后，卫星定位观测方程的基本形式可表示为(即，可将式(16)简化为式(3))：

G＝S+R_B·u (3)

其中，G为卫星定位观测值向量；S为基准相机摄影中心空间坐标向量；R_B表示基准相机外方位角元素

构成的旋转矩阵；u为卫星定位接收机观测中心相对于基准相机摄影中心的空间偏移向量。

将式(3)进行线性化，得到卫星定位观测误差方程：

V_G＝A_G1t₁+A_G2t₂+Dd-L_G (4)

其中，V_G表示卫星定位观测误差，A_G1表示卫星位置坐标对基准相机外方位线元素的偏导系数矩阵子块，A_G2表示卫星位置坐标对基准相机外方位角元素的偏导系数矩阵子块，t₁为基准相机外方位线元素增量，t₂为基准相机外方位角元素增量，D为卫星位置坐标对观测中心空间偏移量的偏导系数矩阵子块，d为卫星定位观测中心空间偏移量增量，L_G为卫星定位观测值的残差向量。

上述线性化的具体推导过程为：设卫星位置坐标的精确值为G，其近似值为G⁰，参数的精确值为Y，其近似值为Y⁰，卫星定位观测方程函数形式为

G＝F₂(Y)

卫星位置坐标精确值G看作卫星定位观测值G⁰(观测值作为近似值)与卫星定位观测误差V_G之和，参数的精确值Y看作参数近似值Y⁰与参数增量y之和，对函数F₂(Y)按泰勒级数展开至一次项，则有

G⁰+V_G＝F₂(G⁰)+F₂'y

上式可写为：

V_G＝F₂'y-L_G

其中，F₂'为卫星位置坐标对参数的偏导系数，L_G＝G⁰-F₂(Y⁰)。对于卫星定位观测方程，涉及的参数包括定向片时刻基准相机外方位元素、卫星定位观测中心相对于基准相机的摄影中心空间偏移量。

y＝[t₁t₂d]^T为参数增量向量。其中：t₁为基准相机外方位线元素增量、t₂为基准相机外方位角元素增量、d为卫星定位观测中心空间偏移量增量。

F₂'＝[A_G1 A_G2 D]。其中，卫星位置坐标对基准相机外方位线元素的偏导系数矩阵子块

卫星位置坐标对基准相机外方位角元素的偏导系数矩阵子块

卫星位置坐标对观测中心空间偏移量的偏导系数矩阵子块

步骤5.3，建立卫星姿态观测方程；

建立卫星姿态观测方程就是建立步骤3获得的卫星姿态观测值phi、omega、kappa与待求参数的函数关系。卫星姿态输出设备本体坐标系与系统坐标系之间的相对旋转角用

ω_A,κ_A表示，载荷本体坐标系可看作由系统坐标系O-UVW分别绕其V,U,W轴依次旋转

ω_A,κ_A后再平移至O'所得到的空间坐标系。卫星姿态输出设备本体空间姿态角构成的旋转矩阵可表示为

上式(式(5))即卫星姿态观测方程的表达式，其中

为相对旋转角构成的旋转矩阵，为绕系统坐标系V轴旋转

对应的旋转矩阵，

为绕系统坐标系U轴旋转ω_A对应的旋转矩阵，

为绕系统坐标系W轴旋转κ_A对应的旋转矩阵；R'为卫星姿态观测值phi、omega、kappa构成的旋转矩阵，R_B为基准相机三轴姿态角即外方位角元素

构成的旋转矩阵。

将式(5)进行线性化，得到卫星姿态观测误差方程形式为：

V_S＝A_St₂+Mm-L_S (6)

其中，V_S表示卫星姿态观测误差，A_S表示卫星姿态对基准相机外方位角元素的偏导系数矩阵子块，t₂为基准相机外方位角元素增量，M为卫星姿态对相对旋转角的偏导系数矩阵子块，m为姿态输出设备本体与基准相机本体的相对旋转角增量，L_S为卫星姿态观测值的残差向量。

所述线性化的具体推导过程为：设卫星三轴的精确值为S，其近似值为S⁰，参数的精确值为Y，其近似值为Y⁰，卫星姿态观测方程的函数形式为

S＝F₃(Y)

卫星姿态精确值S看作卫星姿态观测值S⁰(观测值作为近似值)与卫星姿态观测误差V_S之和，参数的精确值Y看作参数近似值Y⁰与参数增量y之和，对函数F₃(Y)按泰勒级数展开至一次项，则有

S⁰+V_S＝F₃(S⁰)+F₃'y

上式可写为：

V_S＝F₃'y-L_S

其中，F₃'为卫星姿态对参数的偏导系数，L_S＝S⁰-F₂(Y⁰)。对于卫星姿态观测方程，涉及的参数包括定向片时刻基准相机外方位角元素、卫星姿态输出设备本体坐标系与基准相机之间的相对旋转角用

ω_A,κ_A。

y＝[t₂m]^T为参数增量向量。其中：t₂为基准相机外方位角元素增量、m为姿态输出设备本体与基准相机本体的相对旋转角增量。

F₃'＝[A_S M]。其中，卫星姿态对基准相机外方位角元素的偏导系数矩阵子块

卫星姿态对相对旋转角的偏导系数矩阵子块

步骤5.4，联立所述像点坐标误差方程、卫星定位观测误差方程和卫星姿态观测误差方程，建立星载摄影测量系统一体化检校平差模型；

对于式(1)、式(3)、式(5)三类观测方程，利用待求参数初值进行线性化，得到相应的误差方程，即式(2)、式(4)、式(6)。以基准相机外方位线元素、基准相机外方位角元素、地面点物方坐标、成像相机内方位元素、成像相机相对于基准相机的摄影中心空间偏移量、成像相机相对于基准相机的相对旋转角、卫星定位观测中心空间偏移量、姿态输出设备本体与基准相机本体的相对旋转角作为待求参数，其中，成像相机内方位元素、成像相机相对于基准相机的摄影中心空间偏移量、成像相机相对于基准相机的相对旋转角、卫星定位观测中心空间偏移量、姿态输出设备本体与基准相机本体的相对旋转角是待检校参数，将这三类方程同时进行求解，得到星载摄影测量系统一体化检校平差模型如下式所示：

其中：

E为单位矩阵，其对角线元素全为1，非对角线元素全为0；P_G为卫星定位观测值的权矩阵，其对角线元素的数值为

(σ₀为像点坐标观测中误差，σ_G为卫星定位观测中误差)，非对角线元素全为0；P_S为卫星姿态观测值的权矩阵，其对角线元素的数值为

(σ_S为卫星姿态观测中误差)，非对角线元素全为0；

V_X为像点坐标误差，t₁为基准相机外方位线元素增量；t₂为基准相机外方位角元素增量；x为地面点物方坐标增量；c为相机内方位元素参数增量；A_X1表示像点坐标对基准相机外方位线元素的偏导系数矩阵；A_X2表示像点坐标对基准相机外方位角元素的偏导系数矩阵；I表示像点坐标对相机间摄影中心空间偏移量的偏导系数矩阵；J表示像点坐标对相机间相对旋转角的偏导系数矩阵；B表示像点坐标对地面点物方坐标的偏导系数矩阵；B表示像点坐标对地面点物方坐标的偏导系数矩阵；C表示像点坐标对相机内方位元素的偏导系数矩阵；L_X为像点坐标观测值的残差向量；对于前(后)视相机所成像点，i为前(后)视相机相对于正视相机的摄影中心空间偏移量增量；j为前(后)视相机本体相对于正视相机本体的相对旋转角增量；对于正视相机所成像点，有i,j＝0；

V_G表示卫星定位观测误差，A_G1表示卫星位置坐标对基准相机外方位线元素的偏导系数矩阵子块，A_G2表示卫星位置坐标对基准相机外方位角元素的偏导系数矩阵子块，t₁为基准相机外方位线元素增量，t₂为基准相机外方位角元素增量，D为卫星位置坐标对观测中心空间偏移量的偏导系数矩阵子块，d为卫星定位观测中心空间偏移量增量，L_G为卫星定位观测值的残差向量；

V_S表示卫星姿态观测误差，A_S表示卫星姿态对基准相机外方位角元素的偏导系数矩阵子块，t₂为基准相机外方位角元素增量，M为卫星姿态对相对旋转角的偏导系数矩阵子块，m为姿态输出设备本体与基准相机本体的相对旋转角增量，L_S为卫星姿态观测值的残差向量。

根据平差模型，将上述星载摄影测量系统一体化检校平差模型表示为：

V＝Ky-l,P

其中，V＝Ky-l为平差函数模型，P为平差随机模型，即观测值的权矩阵。

具体地，

即观测值改正数向量；

y＝[t₁ t₂ i j x c d m]^T，即待求参数增量向量；其中，t₁为基准相机外方位线元素增量，t₂为基准相机外方位角元素增量，i为成像相机相对于基准相机的摄影中心空间偏移量增量；j为成像相机本体相对于基准相机本体的相对旋转角增量，x为地面点物方坐标增量，c为相机内方位元素参数增量，d为卫星定位观测中心空间偏移量增量，m为姿态输出设备本体与基准相机本体的相对旋转角增量；

即观测值对待求参数的偏导系数矩阵；

即观测值与其近似值之差的向量；

为观测值的权矩阵；其中，E为单位矩阵，其对角线元素全为1，非对角线元素全为0；P_G为卫星定位观测值的权矩阵，其对角线元素的数值为

(σ₀为像点坐标观测中误差，σ_G为卫星定位观测中误差)，非对角线元素全为0；P_S为卫星姿态观测值的权矩阵，其对角线元素的数值为(σ_S为卫星姿态观测中误差)，非对角线元素全为0。

步骤6，基于步骤5建立的星载摄影测量系统一体化检校平差模型，按最小二乘原理求检校参数估计值；

在步骤5中，已经建立了观测值与检校参数的平差数学模型，在此基础上，步骤6的任务是求出检校参数的估计值。

在步骤5得出的卫星摄影测量一体化检校平差模型中，观测值个数即方程个数大于待求参数个数，理论上待求参数存在多个解，采用最小二乘解作为参数的估计值。

按最小二乘原理，构建卫星摄影测量一体化检校平差的法方程，其形式为

其中，K为观测值对待求参数的偏导系数矩阵、P为观测值的权矩阵、l为观测值与其近似值之差的向量；

代表待求参数增量的最小二乘估计值。

由此可得到待求参数增量的最小二乘估计值：

待求参数的估计值为

上式中，Y⁰为待求参数(包括基准相机外方位线元素、基准相机外方位角元素、地面点物方坐标、成像相机内方位元素、成像相机相对于基准相机的摄影中心空间偏移量、成像相机相对于基准相机的相对旋转角、卫星定位观测中心空间偏移量、姿态输出设备本体与基准相机本体的相对旋转角)的初始值。

为取得足够的参数估计值精度，可进行迭代计算，即将上一次平差得到的参数估计值作为下一次平差的参数初值，重新列出误差方程、构造法方程并计算出新的参数估计值，直至参数收敛。

至此，本发明得到了星载摄影测量系统检校参数的高精度估计值，这些系统检校参数包括成像相机内方位元素、成像相机相对于基准相机的摄影中心空间偏移量、成像相机相对于基准相机的相对旋转角、卫星定位观测中心空间偏移量、姿态输出设备本体与基准相机本体的相对旋转角(注意，检校参数是所有待求参数的其中一部分)。

另外，本发明还提供了一种计算机程序，当执行所述计算机程序时，执行上述任意一个实施例所述的方法。

本发明构造了与像点观测值、卫星定位观测值、卫星姿态观测值相应的三类观测方程。一体化检校模型以卫星摄影测量系统的正视相机为基准，确保其他对地相机与正视相机像点的同名光束约束关系，体现了光束法平差的严密性。此外，模型也从数学上反映了卫星姿态输出设备、卫星定位设备与基准相机的几何关联，因此称之为基于光束法平差的星载摄影测量系统一体化检校方法。相比传统方法，一体化检校方法相比传统方法的改进之处在于：

(1)传统检校模型主要针对角元素系统误差进行校正，而本发明提出的一体化检校方法不仅可求出角度系统误差(成像相机相对于基准相机的相对旋转角、姿态输出设备本体与基准相机本体的相对旋转角)，还可以求出内方位元素误差和线偏移量误差，进一步提高了系统检校精度；

(2)传统模型依赖像点坐标观测值，只能利用像点观测方程求解参数，在这种情况下，参数之间往往具有强相关性，参数解算的精度无法保证；而一体化检校方法在像点观测方程基础上，增加了卫星定位观测方程和卫星姿态观测方程。由于三种观测值对待求参数分别具有不同的误差传递特性，因此联立三类方程对检校参数进行一体化解算，可以克服传统方法难以解决的参数强相关，使法方程具有更优良的数学性态，解算稳定性提高，因而可以取得更可靠的检校结果。

(3)一体化检校方法将卫星定位观测值和卫星姿态观测值作为带权观测值，与像点观测值一起纳入平差模型，直接进行参数联合解算，最大限度地发挥了这些观测值的作用，使星载摄影测量系统采集的多源数据价值都得到充分利用。由于一体化检校模型平差几何条件较强，无需按传统方法要求大量控制点参与平差，也大大减少了控制点外业测量的繁重工作，可显著降低检校成本。

以上内容仅为本发明的较佳实施例，对于本领域的普通技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，本说明书内容不应理解为对本发明的限制。

Claims (7)

1.一种基于光束法平差的星载摄影测量系统一体化检校方法，其特征在于，所述方法包括以下步骤：

步骤1，获取基于惯性系的以归一化四元数形式表示的卫星姿态数据；

步骤2，将通过步骤1获得的所述卫星姿态数据的表达形式由归一化四元数的形式转换为欧拉角形式；

步骤3，将通过步骤2获得的欧拉角形式的卫星姿态数据从惯性坐标系转换至地球固定直角坐标系；

步骤4，获取基于地球固定直角坐标系的卫星定位观测值；

步骤5，构建卫星摄影测量一体化检校平差模型；

其中，步骤5具体包括以下子步骤：步骤5.1，建立像点观测方程，所述像点观测方程如下式(1)所示：

其中，(x,y)为像点坐标；(x₀,y₀)为相机的像主点坐标；f为成像相机主距；(X_S',Y_S',Z_S')为成像相机外方位线元素；a_j,b_j,c_j(j＝1,2,3)为成像相机外方位角元素

构成的旋转矩阵方向余弦值；(X_T,Y_T,Z_T)为该像点对应的地面点在地固坐标系下的空间坐标；

将式(1)进行线性化，得到像点坐标误差方程为：

V_X＝A_X1t₁+A_X2t₂+Ii+Jj+Bx+Cc-L_X (2)

其中，V_X为像点坐标误差，t₁为基准相机外方位线元素增量；t₂为基准相机外方位角元素增量；x为地面点物方坐标增量；c为相机内方位元素参数增量；A_X1表示像点坐标对基准相机外方位线元素的偏导系数矩阵；A_X2表示像点坐标对基准相机外方位角元素的偏导系数矩阵；I表示像点坐标对相机间摄影中心空间偏移量的偏导系数矩阵；J表示像点坐标对相机间相对旋转角的偏导系数矩阵；B表示像点坐标对地面点物方坐标的偏导系数矩阵；B表示像点坐标对地面点物方坐标的偏导系数矩阵；C表示像点坐标对相机内方位元素的偏导系数矩阵；L_X为像点坐标观测值的残差向量；对于前后视相机所成像点，i为前后视相机相对于正视相机的摄影中心空间偏移量增量；j为前后视相机本体相对于正视相机本体的相对旋转角增量；对于正视相机所成像点，有i,j＝0；

步骤5.2，建立卫星定位观测方程，所述卫星定位观测方程如下式(16)所示：

式中，[X_g Y_g Z_g]^T表示卫星定位观测值，[X_S Y_S Z_S]^T表示基准相机摄影中心坐标，[u_g v_gw_g]^T表示卫星定位观测中心在卫星摄影测量系统坐标系中的三维坐标，R_B表示基准相机外方位角元素

构成的旋转矩阵；

在建立卫星定位接收机与基准相机的几何关系后，将式(16)简化为式(3)

G＝S+R_B·u (3)

其中，G为卫星定位观测值向量；S为基准相机摄影中心空间坐标向量；R_B表示基准相机外方位角元素

构成的旋转矩阵；u为卫星定位接收机观测中心相对于基准相机摄影中心的空间偏移向量；

将式(3)进行线性化，得到卫星定位观测误差方程：

V_G＝A_G1t₁+A_G2t₂+Dd-L_G (4)

其中，V_G表示卫星定位观测误差，A_G1表示卫星位置坐标对基准相机外方位线元素的偏导系数矩阵子块，A_G2表示卫星位置坐标对基准相机外方位角元素的偏导系数矩阵子块，t₁为基准相机外方位线元素增量，t₂为基准相机外方位角元素增量，D为卫星位置坐标对观测中心空间偏移量的偏导系数矩阵子块，d为卫星定位观测中心空间偏移量增量，L_G为卫星定位观测值的残差向量；

步骤5.3，建立卫星姿态观测方程，所述卫星姿态观测方程如下式(5)所示：

其中，

为相对旋转角构成的旋转矩阵，

为绕系统坐标系V轴旋转

对应的旋转矩阵，

为绕系统坐标系U轴旋转ω_A对应的旋转矩阵，

为绕系统坐标系W轴旋转κ_A对应的旋转矩阵；R'为卫星姿态观测值phi、omega、kappa构成的旋转矩阵，R_B为基准相机三轴姿态角即外方位角元素

构成的旋转矩阵；

将式(5)进行线性化，得到卫星姿态观测误差方程如下式(6)所示：

V_S＝A_St₂+Mm-L_S (6)

其中，V_S表示卫星姿态观测误差，A_S表示卫星姿态对基准相机外方位角元素的偏导系数矩阵子块，t₂为基准相机外方位角元素增量，M为卫星姿态对相对旋转角的偏导系数矩阵子块，m为姿态输出设备本体与基准相机本体的相对旋转角增量，L_S为卫星姿态观测值的残差向量；

步骤5.4，联立所述像点坐标误差方程、卫星定位观测误差方程和卫星姿态观测误差方程，建立星载摄影测量系统一体化检校平差模型，如下式所示：

其中，所述步骤5.4还包括：

根据平差模型，将上述星载摄影测量系统一体化检校平差模型表示为：

V＝Ky-l,P

其中，V＝Ky-l为平差函数模型，P为平差随机模型，即观测值的权矩阵；

具体地，

即观测值改正数向量；

y＝[t₁ t₂ i j x c d m]^T，即待求参数增量向量；

即观测值对待求参数的偏导系数矩阵；

即观测值与其近似值之差的向量；

为观测值的权矩阵；其中，E为单位矩阵，其对角线元素全为1，非对角线元素全为0；P_G为卫星定位观测值的权矩阵，其对角线元素的数值为

σ₀为像点坐标观测中误差，σ_G为卫星定位观测中误差，非对角线元素全为0；P_S为卫星姿态观测值的权矩阵，其对角线元素的数值为σ_S为卫星姿态观测中误差，非对角线元素全为0；

步骤6，基于步骤5建立的星载摄影测量系统一体化检校平差模型，按最小二乘原理求待检校参数估计值；其中，所述步骤6具体包括：

按最小二乘原理，构建卫星摄影测量一体化检校平差的法方程，其形式为

其中，K为观测值对待求参数的偏导系数矩阵，P为观测值的权矩阵，l为观测值与其近似值之差的向量；

表示待求参数增量的最小二乘估计值；

由此可得到参数增量的最小二乘估计值：

则待求参数的估计值为

上式中，Y⁰为待求参数的初始值。

2.根据权利要求1所述的方法，其特征在于，在所述步骤1中，所述卫星姿态数据优选的用归一化四元数表示为：

其中，q为一个卫星姿态值，

q₁、q₂和q₃表示归一化四元数的虚部，q₄表示归一化四元数的实部，其中，

3.根据权利要求1所述的方法，其特征在于，在所述步骤2中，根据下述转换公式(11)将卫星姿态数据的表达形式由归一化四元数的形式转换为欧拉角形式：

其中，Roll、Pitch和Yaw分别为以欧拉角形式表示的横滚角、俯仰角和偏航角，A_eular表示转换后的三轴欧拉角矩阵。

4.根据权利要求1所述的方法，其特征在于，在所述步骤3中，根据下述公式将卫星姿态数据从惯性坐标系转换至地球固定直角坐标系：

R_ITRS＝H(t)·R_GCRS

其中，R_ITRS为卫星在地球固定直角坐标系下的旋转矩阵，R_ITRS中的元素是由基于地球固定直角坐标系的三个姿态角phi、omega、kappa计算的方向余弦；R_GCRS为卫星在惯性坐标系下的旋转矩阵，R_GCRS中的元素是由基于惯性坐标系的三个姿态角Roll、Pitch、Yaw计算的方向余弦；H(t)＝W(t)·R(t)·Q(t)，W(t)，R(t)和Q(t)分别对应极移、自转和岁差章动转换矩阵。

5.根据权利要求1所述的方法，其特征在于，所述待检校参数是所述待求参数的一部分。

6.根据权利要求1所述的方法，其特征在于，其特征在于，所述待求参数包括基准相机外方位线元素、基准相机外方位角元素、地面点物方坐标、成像相机内方位元素、成像相机相对于基准相机的摄影中心空间偏移量、成像相机相对于基准相机的相对旋转角、卫星定位观测中心空间偏移量、姿态输出设备本体与基准相机本体的相对旋转角；所述待检校参数包括成像相机内方位元素、成像相机相对于基准相机的摄影中心空间偏移量、成像相机相对于基准相机的相对旋转角、卫星定位观测中心空间偏移量、姿态输出设备本体与基准相机本体的相对旋转角。

7.一种介质，存储有计算机程序，当执行该计算机程序时，执行上述权利要求1-6中的任意一项所述的方法。

Images