CN101609551A - 基于线阵推扫式异步采样卫星影像几何模型的正射校正方法 - Google Patents

Abstract

本发明公开了一种基于线阵推扫式异步采样卫星影像几何模型的正射校正方法，采用双中心线投影的投影模式来描述线阵推扫式异步采样卫星影像的几何变形规律，通过恢复双中心线投影的各项参数，建立卫星影像像素坐标与对应的地面坐标之间的关系，从而实现对线阵推扫式异步采样卫星影像正射校正，包括以下步骤：(1)确定双中心线投影模型参数初始值；(2)误差方程系数的计算与模型参数改正数的解算；(3)计算输出正射影像的范围；(4)正射影像格网点灰度值的内插计算。

Description

基于线阵推扫式异步采样卫星影像几何模型的正射校正方法

技术领域

本发明涉及一种正射校正方法，特别是一种基于双中心线投影的线阵推扫式异步采样卫星影像正射校正方法，属于卫星影像处理领域。

背景技术

从上世纪九十年代末至今，由于高分辨率商业遥感卫星影像的大量获取和分发，人们开始将这些卫星影像用于地理空间数据的采集和建库。卫星遥感影像不可避免地存在各种变形，地理空间数据库要求消除变形，得到精确的二维或者三维地理空间信息。高分辨率光学遥感卫星大都采用线阵推扫式传感器，获取的每一行影像有独立的外方位元素，不能采用传统的中心投影共线方程进行整体几何处理。因此，对线阵推扫式卫星影像的几何建模一直是学术界研究的一个重要课题。

线阵推扫式遥感卫星的成像模式主要有两种：同步采样成像模式和异步采样成像模式。在同步采样模式模式下，传感器的地面采样速度与卫星的地速相等，在地图投影坐标系中传感器的外方位元素尽量保持不变。在异步采样模式下，传感器的地面采样速度与卫星的地速不相等。

针对线阵推扫式遥感卫星不同的成像模式，学术界主要提出了严密几何模型、简化几何模型和主观模型。

严密几何模型根据成像光束在传感器线阵方向符合中心投影原理，对每一行影像使用中心投影共线方程：

$\left.\begin{array}{c}0=a_{i1}(X-X_{\mathrm{si}})+b_{i1}(Y-Y_{\mathrm{si}})+c_{i1}(Z-Z_{\mathrm{si}})\\ y=-f\frac{a_{i2}(X-X_{\mathrm{si}})+b_{i2}(Y-Y_{\mathrm{si}})+c_{i2}(Z-Z_{\mathrm{si}})}{a_{i3}(X-X_{\mathrm{si}})+b_{i3}(Y-Y_{\mathrm{si}})+c_{i3}(Z-Z_{\mathrm{si}})}\end{array}\right\}---\left(1\right)$

其中

在严密模型中不同行之间的外方位元素的变化采用一般多项式进行描述。严密模型需要传感器的物理和几何参数以及卫星星历数据，并且存在参数之间的高度相关的现象，导致参数解不稳定。

在严密模型的基础上假设传感器在成像过程中姿态不变，卫星沿直线轨道运行，或者舍弃一些在具体条件下不重要的参数，得到简化的几何模型，其中包括平行推扫式成像几何模型和平行投影模型。简化模型只能描述卫星影像的主要变形因素，在地面起伏较大以及卫星姿态不稳定的条件下，简化么模型得到的几何精度较低。

主观模型不考虑成像的几何原理和传感器的物理参数，直接采用一般多项式或者有理多项式来建立影像坐标与地面坐标的关系，其参数与实际的成像几何没有直接的对应关系。采用一般多项式作为模型时，与确定的地面控制点拟合很好，但在其它点的内插值则可能有明显的偏离，即可能在某些点处产生振荡，导致用多项式近似计算的中误差明显超出平均误差。有理多项式模型也是一种主观模型，从理论上讲有理多项式模型所描述的数学关系只在控制点处成立，而在其它地方都是近似的。因此，有理多项式模型的精度与控制点的精度、分布、数量以及纠正范围密切相关。

由于线阵推扫式异步采样遥感卫星出现得比较晚，其成像几何尚未引起足够的重视，1B级异步采样影像的几何模型还未建立。

发明内容

本发明目的就是针对线阵推扫式异步采样卫星的1B级影像的几何建模。

实现本发明目的采用的技术方案是：利用双中心线投影模型描述线阵推扫式异步采样卫星1B级影像的影像坐标与对应点地面坐标之间的数学关系。线阵推扫式异步采样卫星影像正射纠正方法包括以下步骤：

(1)双中心投影模型参数初始值的确定；

(2)误差方程系数的计算与模型参数改正数的解算：将误差方程系数矩阵标准化，将常数项中心化，用标准化的误差方程解算标准化的模型参数改正值，再采用标准化的逆过程求出模型参数改正值；

(3)输出正射影像范围的计算；

(4)正射影像格网点灰度值的内插计算。

本发明通过下式建立1B级影像的影像坐标与对应的地面坐标之间的关系：

$\left.\begin{array}{c}F=\frac{a_3(X-X_o)+b_3(Y-Y_o)+c_3(Z-Z_o)+f^{\′}}{a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}}-\frac{a_1(X-X_o)+b_1(Y-Y_o)+c_1(Z-Z_o)}{a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)}\\ G=\frac{a_3(X-X_o)+b_3(Y-Y_o)+c_3(Z-Z_o)+f^{\′}}{a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}}-\\ \frac{[a_2(X-X_o)+b_2(Y-Y_o)+c_2(Z-Z_o)][a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)-f]}{[a_2^{\′}(x-x_o)+b_2^{\′}(y-y_o)].[a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}]-[a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)](f+f^{\′})\mathrm{tg\α}}\\ \frac{[a_2^{\′}(x-x_o)+b_2^{\′}(y-y_o)](f+f^{\′})-[a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)](f+f^{\′})\mathrm{tg\α}}{[a_2^{\′}(x-x_o)+b_2^{\′}(y-y_o)].[a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}]-[a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)](f+f^{\′})\mathrm{tg\α}}=0\end{array}\right\}$

上式1B级影像双中心线投影共线方程，该方程为隐函数方程组，除了像点观测值x和y以外，还有x_o、y_o、X_o、Y_o、Z_o、

ω、κ、

ω′、κ′、f、f′、α共14个参数，只要有7个以上的像点坐标和对应点的地面坐标，就可以用最小二乘法唯一确定上述参数的解。

本发明具有以下优点：

双中心线投影模型不但适用于处理经过几何改正的线阵推扫式异步采样卫星影像，而且可以通过用来对中心投影影像、线阵推扫式同步采样的原始卫星影像进行几何处理。

附图说明

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

图2为异步采样对称位置成像光束交会示意图。

图3为线阵方向与运行轨道垂直时成像光束的交会示意图。

具体实施方式

中心投影是广为熟悉的一种投影模式，实现中心投影必须具备2个要素：投影中心、投影面。从物点发出的光线通过投影中心，在投影面上形成影像，从而实现中心投影。在中心投影中，投影中心是几何上的一个点，光学上以一个小孔替代小孔。设想该几何点为两条直线的交点，光学上可以用两个相交的狭缝来实现，投影面不变，中心投影仍然成立。将上述假设中的两个相交的狭缝分离，使它们既不平行，也不相交，同时保持投影面不变，从物点发出的光线依次通过两个狭缝，在投影面上仍然可以成像。这里把这种由物点发出的光线依次通过两个既不平行、也不相交的狭缝，然后在投影面上形成影像的过程称作双中心线投影。几何上可以把这两个既不平行、也不相交的狭缝看作两条既不平行、也不相交的直线，这里把它们称做投影中心线。由此可见，双中心线投影包含3个要素：2条投影中心线和1个投影面。

线阵推扫式传感器异步采样的双中心线投影原理

探讨异步采样卫星遥感影像成像的几何原理，设想如图3所示，卫星首先从S₁点开始，对某目标区域进行扫描，此时卫星传感器投影中心是S₁，主点为o₁，主光轴S₁o₁向前倾斜一个角度φ，传感器线阵列为a₁b₁，传感器的焦距为f。卫星在向前运行的过程中，对目标区域进行扫描，同时为了提高传感器的感应时间，卫星以一个接近均匀的角速度向后旋转，到达位置S的时候，卫星传感器线阵列为ab，此时主光轴So垂直于S₁S。卫星继续向前运行，到达位置S₂时，传感器线阵列为a₂b₂，主光轴S₂o₂向后倾斜φ。假设卫星平台向后转的速度是是对称的，卫星从S₁运行到S所用的时间和从S运行到S₂所用的时间相等，又因为卫星运行的速度也是一个常数，所以有

S₁S＝S₂S    (3)

这里继续采用一个基本假设，即在一景影像范围内，卫星平台质心的运行轨道近似为一条直线。由于卫星在进行异步采样时，要不停的向后转动，因此，在采样过程中，传感器的投影中心相对于卫星质心来说，会有一定的变化。卫星本身的尺寸较小，与扫描一景影像的过程中卫星的运行距离相比，传感器投影中心相对于卫星质心位置变化是一个很小的量，所以，在一景影像范围内，传感器投影中心的运行轨迹S₁S₂也可以看作近似为一条直线。

卫星平台在采集一景影像的过程中，由于运行时间较短，可以认为卫星平台的侧偏角和自旋角不发生变化，也就是说，卫星传感器在对目标区域的扫描过程中，传感器线阵方向不发生变化，而是保持相互平行的状态。

在上述两个基本前提条件下，以图3中的S为原点，以oS方向为w轴正方向，以垂直于面Sab且与SS₂成锐角的方向为u轴正方向，建立起右手直角坐标系S-uvw作为传感器参考坐标系。设S₁在S-uvw坐标系中的坐标为S₁(-Δu，-Δv，0)，因为S₁、S₂相对于S是对称的，则有S₂在S-uvw坐标系中的坐标为S₂(Δu，Δv，0)，因此，传感器在S₁、S、S₂三点的外方位元素分别为：

将上述参数代入式(1)的第一个等式和(2)分别得到

$\left.\begin{array}{c}\cos \mathrm{\φ}.(u+\mathrm{\Δu})+\sin \mathrm{\φ}.w=0\\ \cos \mathrm{\φ}.(u-\mathrm{\Δu})-\sin \mathrm{\φ}.w=0\\ u=0\end{array}\right\}---\left(7\right)$

将上述第三式与第一式联立求解，得到：

$\left.\begin{array}{c}w=-\mathrm{ctg\φ}.\mathrm{\Δu}\\ u=0\end{array}\right\}---\left(8\right)$

将(7)的第三式与第二式联立求解，也得到(8)。

从式(8)可以看出，它实际上是一个直线方程，该直线在u轴方向和w轴方向的分量为零。由此证明了Sab、S₁a₁b₁、S₂a₂b₂三个面交于同一条直线，该直线为平行于ab，与S的距离取决于SS₁和SS₂在垂直于Sab的直线方向的分量以及φ的大小。

由以上的分析可以看出，卫星在扫描成像的过程中，投影中心处于点S的时候，以及处于与S成对称位置的S₁点和S₂的时候，成像光束交会于一条平行于传感器线阵的直线。进一步可以证明，当卫星运行轨道与传感器线阵相互垂直时，经过主点的成像光束交会于一个点，与主光束对称的成像光束也分别交会于各自的点，这些点构成一条平行于传感器线阵的直线，如图3所示，So、So₁、So₂交会于点S′，Sa₁、Sa₂交会于点S₁′，Sb₁、Sb₂交会于点S₂′。

从图3可以看出，在扫描一景影像的过程中，所有的成像光束必须经过直线S₁S₂，S₁S₂的作用类似于中心投影中的投影中心，这里把它称作第一投影中心线。为了探讨异步采样卫星遥感影像的整体构像方程，这里再进一步假设卫星从S₁开始，扫描到主光轴与卫星运行轨迹相垂直的位置S，然后继续扫描，直到与S₁对称的位置S₂，完成一景影像的扫描，整个过程中的成像光束全都交会于一条位于Sab平面内，且平行于ab的直线，也就是说，整景影像的成像光束除了经过S₁S₂直线外，还同时经过另外一条直线。由于该条线的作用也类似于中心投影中的投影中心，这里把它称作第二投影中心线。从几何原理可知，过一个已知空间点和两条既不平行也不相交的空间直线，可以唯一地确定一条空间直线。因为过地面点且同时经过两条投影中心线的直线与过像点且同时经过两条投影中心线的直线是同一成像光线，它们应该共线，可形成共线方程。由于上述结论是从成像光束的空间关系得出的，因此不会因为影像空间参照系的改变而改变。设想已知一定数量的点的像点坐标和像点对应的地面点的地面坐标，就可以用空间后方交会法，求出这两条直线的方程，从而确立从像点坐标到地面坐标的对应关系，实现异步采样卫星遥感影像的几何处理。

本发明利用双中心线投影模型进行正射纠正的步骤：

(1)双中心线投影模型参数初始值的确定

基于双中心线投影模型的EROS A1卫星的1B级影像的模型参数的解算步骤如下：

根据双中心线模型的定义，参数中f实际上是以影像比例尺表示的传感器焦距，由于EROS A1卫星的1B级影像是根据卫星星历和传感器姿态参数，将原始影像投影到WGS84 UTM坐标系后得到的，在处理EROS A1卫星的1B级影像时，可以用影像中心对应的地面点与卫星之间的距离作为f的初始值。而f′是传感器线阵与传感器投影中心组成的平面在不同位置的交线与虚拟像平面之间的距离，可以根据卫星地面采样弧长、卫星轨道弧长在地面的投影长度以及f等3个参数按下式计算：

$\frac{f^{\′}+f}{D}=\frac{f^{\′}}{d}---\left(9\right)$

式中，D为卫星的飞行弧长，d为地面采样距离，f为第一焦距。上述参数可以从影像的元数据中得到或者根据影像元数据推算出来。

由于EROS A1的1B级影像是经过几何纠正的，并且投影到WGS84 UTM坐标系，而正射纠正也在WGS84 UTM坐标系中进行，所以主点o的坐标可以直接在EROS A1的1B级影像上量取。具体做法是找出影像扫描行方向长度最短的一条线的中点，量出该点的X、Y坐标，作为该景影像的主点坐标X_o和Y_o的初始值。而影像扫描行方向可以由EROS A1的1B级影像的上边缘方向确定。

转角

是卫星处于轨道平面内轨道的垂线方向与WGS84 UTM坐标系的Z轴的夹角，实际上反映了卫星轨道与水平面的夹角，一般情况下近似为零。转角ω是卫星成像时主光轴与本地大地法线之间的夹角，可以用元数据中传感器的侧偏角作为初始值。转角κ是卫星轨道方向在WGS84 UTM坐标系中的方位角，可以由EROS A1的1B级影像的上边缘的向下的垂线方向与WGS84 UTM坐标系的X轴之间的夹角来确定。而转角α是卫星实际飞行轨迹与卫星轨道方向的夹角，是地球自转和传感器线阵方向与卫星运行方向不垂直两个因素共同造成的，其初始值可以根据卫星的轨道参数、地球自转速度以及影像所处的地理位置来计算；也可以在影像上直接量取，其做法是在EROS A1卫星的1B级影像上画出第一行影像和最后一行影像的中点之间的连线，量取该连线与第一行影像的垂线之间的夹角β，然后利用下面的近似公式计算α：

$\mathrm{\α}=\frac{d}{D}\×\mathrm{\β}---\left(10\right)$

(2)误差方程系数的计算与模型参数改正数的解算

误差方程的系数是根据上一个步骤计算出的模型参数初始值、量测的像点坐标、对应于像点的地面点坐标等三组数据计算的。由于EROS A1的1B级影像经过了几何改正，并且投影到WGS84 UTM坐标系中，所以像点坐标与地面坐标同属于一个坐标系，像点坐标也用WGS84 UTM坐标表示。将模型参数初始值、像点坐标、对应于像点的地面点坐标代入双中心线投影共线方程，就可以得到误差方程的系数矩阵和常数项。设一景影像中量测的控制点有n个，就可以得到2n个误差方程，当n＞4时，需要用最小二乘法求解模型参数的改正值。

由于模型中各个参数值在数量级上差别相当大，直接将误差方程法化后进行平差计算会遇到法方程病态化问题。为了使法方程各个元素在同一个数量级，需要将误差方程系数矩阵标准化，将常数项中心化，用标准化的误差方程解算标准化的模型参数改正值，再采用标准化的逆过程求出模型参数改正值，具体步骤如下。

将上面求出误差方程系数矩阵记为

$A=\left[\begin{array}{cccccccc}{a}_{\mathrm{1,1}}& a_{1,2}& a_{\mathrm{1,3}}& a_{\mathrm{1,4}}& a_{\mathrm{1,5}}& a_{\mathrm{1,6}}& a_{\mathrm{1,7}}& a_{\mathrm{1,8}}\\ a_{\mathrm{2,1}}& a_{\mathrm{2,2}}& a_{\mathrm{2,3}}& a_{\mathrm{2,4}}& a_{\mathrm{2,5}}& a_{\mathrm{2,6}}& a_{\mathrm{2,7}}& a_{\mathrm{2,8}}\\ a_{\mathrm{3,1}}& a_{\mathrm{3,2}}& a_{\mathrm{3,3}}& a_{\mathrm{3,4}}& a_{\mathrm{3,5}}& a_{\mathrm{3,6}}& a_{\mathrm{3,7}}& a_{\mathrm{3,8}}\\ a_{\mathrm{4,1}}& a_{\mathrm{4,2}}& a_{\mathrm{4,3}}& a_{\mathrm{4,4}}& a_{\mathrm{4,5}}& a_{\mathrm{4,6}}& a_{\mathrm{4,7}}& a_{\mathrm{4,8}}\\ & & & .& & & & \\ & & & .& & & & \\ & & & .& & & & \\ a_{2n-\mathrm{1,1}}& a_{2n-\mathrm{1,2}}& a_{2n-\mathrm{1,3}}& a_{2n-\mathrm{1,4}}& a_{2n-\mathrm{1,5}}& a_{2n-\mathrm{1,6}}& a_{2n-\mathrm{1,7}}& a_{2n-\mathrm{1,8}}\\ a_{2n,1}& a_{2n,2}& a_{2n,3}& a_{2n,4}& a_{2n,5}& a_{2n,6}& a_{2n,7}& a_{2n,8}\end{array}\right]---\left(11\right)$

常数项记为

W＝[w₁ w₂ … w_2n]′           (12)

求出误差方程系数矩阵每一列的平均值

${\stackrel{\‾}{a}}_j=\frac{1}{2n}\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{a}_{i,j}---\left(13\right)$

引进每一列的缩放系数

$R_j=\sqrt{\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{(a_{i,j}-{\stackrel{\‾}{a}}_j)}^2}---\left(14\right)$

然后将矩阵A的每一列元素减去该列元素的平均值，再除以该列的缩放系数R_j得到

$\stackrel{\‾}{A}=\left[\begin{array}{cccc}\frac{(a_{\mathrm{1,1}}-\stackrel{\‾}{a_1})}{R_1}& \frac{(a_{\mathrm{1,2}}-\stackrel{\‾}{a_2})}{R_2}& ...& \frac{(a_{\mathrm{1,8}}-\stackrel{\‾}{a_8})}{R_8}\\ \frac{(a_{\mathrm{2,1}}-\stackrel{\‾}{a_1})}{R_1}& \frac{(a_{\mathrm{2,2}}-\stackrel{\‾}{a_2})}{R_2}& ...& \frac{(a_{\mathrm{2,8}}-\stackrel{\‾}{a_8})}{R_8}\\ .& .& & .\\ .& .& & .\\ .& .& & .\\ \frac{(a_{2n,1}-\stackrel{\‾}{a_1})}{R_1}& \frac{(a_{2n,2}-\stackrel{\‾}{a_2})}{R_2}& ...& \frac{(a_{2n,8}-\stackrel{\‾}{a_8})}{R_8}\end{array}\right]---\left(15\right)$

上式即为标准化的误差方程系数矩阵，用标准化的误差方程计算出了法方程系数矩阵实际上就是模型之间的相关系数矩阵。同样将误差方程常数项中心化，记

$\stackrel{\‾}{w}=\frac{1}{2n}\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{w}_i---\left(16\right)$

W＝[w₁-w w₂-w…w_2n-w]′           (17)

标准化的误差方程系数矩阵和中心化的误差方程常数项构成新的误差方程

V＝A·δX+W                       (18)

用最小二乘法解上式，得到经过缩放后的参数改正量。

δX＝[A′·P·A]^-1A′·P·W       (19)

参数的实际改正量为

δx_j＝δx_j×R_j                    (20)

在式(19)中，P为观测值权阵。如果像控点的地面坐标测量误差对空间后方交会平差精度有比较显著的影响，可以对地面点坐标列立误差方程，给定其权矩阵，采用上面一样的方法解算。

解算完毕后，将各项参数改正量与对应的参数初始值相加，得到改正后的参数值。将改正后的参数值、像点坐标以及像点对应的地面点坐标代入双中心线投影共线方程，重新计算误差方程系数项和常数项，按照前面的过程重新平差计算，这样反复迭代，直到各项参数的改正量小于允许值为止。

平差结束后，用各项参数最后的值和像点坐标、地面坐标重新计算出误差方程的常数项，得到每个像点坐标的残差。利用像点坐标的残差统计像点坐标测量精度，作为精度评定的依据。在控制点较多的情况下，可以采用一部分点作为控制点，另外一部分点作为独立检查点，利用检查点的地面坐标和模型参数平差值，求出这些检查点的像点坐标，将实际量测的检查点的像点坐标与计算出的像点坐标进行比较，以评定模型定向的精度。

(3)输出正射影像范围的计算

纠正后图像边界范围的确定过程如下：

3-1求解原始图像的四个角点1、2、3、4对应的地面点在地图投影坐标系中的坐标值，得到8个坐标值：(X₁，Y₁)、(X₂，Y₂)、(X₃，Y₃)、(X₄，Y₄)。

由于双中心线投影模型从影像的二维坐标到地面的三维坐标变换需要用到地面点的高程，可以用EROS A1 1B级影像四个角点在WGS84 UTM坐标系中的坐标值，直接从DEM中内插出对应的地面的高程值。将角点的影像坐标、对应的高程值和解算出的模型参数值代入共线方程，求解出每一个角点对应的地面点的坐标，此时算出的地面点坐标与从影像上直接量取的坐标值不会相同。用解算出的地面点坐标再次从DEM中内插出对应的高程值，然后将角点在影像上的坐标值、模型参数和内插出的地面点高程再次代入共线方程，再次求出角点对应的地面点坐标。采用上述方法反复迭代，直到两次计算出的地面坐标相等，即得出影像四个角点对应的地面点坐标。

3-2利用8个坐标值计算正射影像的范围

$\left.\begin{array}{c}{X}_{\min }=d\×\mathrm{int}\left[\frac{\min (X_1,X_2,X_3,X_4)}{d}\right]-d\\ Y_{\min }=d\×\mathrm{int}\left[\frac{\min (Y_1,Y_2,Y_3,Y_4)}{d}\right]-d\\ X_{\max }=d\×\mathrm{int}\left[\frac{\min (X_1,X_2,X_3,X_4)}{d}\right]+d\\ Y_{\max }=d\×\mathrm{int}\left[\frac{\min (Y_1,Y_2,Y_3,Y_4)}{d}\right]+d\end{array}\right\}---\left(21\right)$

式中X_min，Y_min为数字正射影像的起始点(左下角)坐标，X_max，Y_max为数字正射影像的终止点(右上角)坐标，d为数字正射影像的地面分辨率。

正射影像格网点对应于原影像中像点坐标的计算

求出数字正射影像的范围后，即可计算出正射影像上每一个像元中心对应的地面坐标值。以正射影像的左下角开始，行坐标向上取正值，列坐标向左取正值，对应第i行第j列像元中心的坐标值为

$\left.\begin{array}{c}{X}_{\mathrm{ij}}=X_{\min }+d\×i\\ Y_{\mathrm{ij}}=Y_{\min }+d\×j\end{array}\right\}---\left(22\right)$

利用上式计算出X_ij和Y_ij坐标后，就可根据DEM内插出对应的高程值Z_ij。将数字正射影像每个格网点的三维坐标(X_ij，Y_ij，Z_ij)代入共线方程，即可求出对应的像点坐标。

(4)正射影像格网点灰度值的内插计算

求出数字正射影像格网点中心在原影像中对应的像点坐标后，即可内插出该格网点的影像灰度值。影像灰度值的内插方法可以采用SINC函数法、双三次卷积重采样法和最邻近像元法等。

如果计算出的影像坐标超出原影像的范围，可采用空白值填充。采用上述方法将正射影像上每一个格网点的灰度值计算出来，并赋给该像元，即得到完整的数字正射影像，整个正射纠正过程结束。

Claims (2)

1.一种基于线阵推扫式异步采样卫星影像几何模型的正射校正方法，其特征在于采用双中心线投影的投影模式来描述线阵推扫式异步采样卫星影像的几何变形规律，通过恢复双中心线投影的各项参数，建立卫星影像像素坐标与对应的地面坐标之间的关系，从而实现对线阵推扫式异步采样卫星影像正射校正，包括以下步骤：

(1)确定双中心线投影模型参数初始值；

(2)误差方程系数的计算与模型参数改正数的解算；

(3)计算输出正射影像的范围；

(4)正射影像格网点灰度值的内插计算。

2.根据权利要求1所述基于线阵推扫式异步采样卫星影像几何模型的正射校正方法，其特征在于通过下式建立1B级影像的影像坐标与对应的地面坐标之间的关系：

$\left.\begin{array}{c}F=\frac{a_3(X-X_o)+b_3(Y-Y_o)+c_3(Z-Z_o)+f^{\′}}{a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}}-\frac{a_1(X-X_o)+b_1(Y-Y_o)+c_1(Z-Z_o)}{a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)}=0\\ G=\frac{a_3(X-X_o)+b_3(Y-Y_o)+c_3(Z-Z_o)+f^{\′}}{a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}}\\ \frac{[a_2(X-X_o)+b_2(Y-Y_o)+c_2(Z-Z_o)][a_3^{\′}(x-x_0)+b_3^{\′}(y-y_o)-f]}{[a_2^{\′}(x-x_o)+b_2^{\′}(y-y_o)]\·[a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}]-[a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)](f+f^{\′})\mathrm{tg\α}}\\ \frac{[a_2^{\′}(x-x_o)+b_2^{\′}(y-y_o)](f+f^{\′})-[a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)](f+f^{\′})\mathrm{tg\α}}{[a_2^{\′}(x-x_o)+b_2^{\′}(y-y_o)]\·[a_3^{\′}(x-x_o)+b_3^{\′}(y-y_o)+f^{\′}]-[a_1^{\′}(x-x_o)+b_1^{\′}(y-y_o)](f+f^{\′})\mathrm{tg\α}}=0\end{array}\right\}$

上式1B级影像双中心线投影共线方程，该方程为隐函数方程组，除了像点观测值x和y以外，还有x_o、y_o、X_o、Y_o、Z_o、

ω、κ

ω′、κ′、f、f′、α共14个参数，只要有7个以上的像点坐标和对应点的地面坐标，就可以用最小二乘法唯一确定上述参数的解。

Images