CN103697864B - 一种基于大虚拟相机的窄视场双相机影像拼接方法 - Google Patents

Abstract

一种基于大虚拟相机的窄视场双相机影像拼接方法，根据两个单相机建立大虚拟相机，根据两个单相机的几何成像参数，建立大虚拟相机的几何成像参数；根据两个单相机和大虚拟相机的几何成像参数，建立各自的几何成像模型；解算并输出大虚拟相机对应的有理多项式模型系数；根据基于几何成像模型的坐标正算过程和坐标反算过程对两个单相机的图像分别进行间接法几何纠正，得到大虚拟相机图像坐标系下的两个图像，得到拼接后的大虚拟相机图像。本发明巧妙地借用了大虚拟相机的概念，在实现窄视场双相机影像的高精度拼接的同时，还提供与之相对应的有理多项式模型；且处理过程是全自动的、无需人工干预，适用于地面预处理过程。

Description

一种基于大虚拟相机的窄视场双相机影像拼接方法

技术领域

本发明属于航天与航空摄影测量领域，涉及两台窄视场线阵相机同时推扫成像的情况下，一种基于大虚拟相机的窄视场双相机影像拼接方法。

背景技术

线阵推扫成像方式是目前获取高分辨率光学卫星影像的主要传感器。为了提高光学影像的空间分辨率，常采用长焦镜头；而长焦镜头导致观测视场变窄；为增加观测视场角，采用多片CCD（电荷耦合元件图像传感器）拼接或多台相机同时观测的方式。在多台相机同时观测的情况下，每台相机有一套独立的光学系统，遵循各自的几何成像模型，这给后续几何处理带来额外的工作。

常规的双相机图像处理方法有两种。一种是针对单相机图像进行高精度几何处理，即：对两个相机的图像独立进行几何纠正，做成正射影像之后再进行拼接和匀光处理；这种方法几何精度较高，但纠正选点和拼接处理的工作量都比较大。另一种是不考虑原始单相机成像的几何条件，仅根据双相机重叠区域内的图像，基于同名点匹配进行影像拼接处理；这种方法得到的拼接图像缺少严格的物理成像模型，几何精度比较差，难以满足测绘等领域的需求。由此可见，一种既能保证几何精度又不增加常规后续处理工作量的双相机图像拼接方法是迫在眉睫的。

发明内容

本发明所要解决的问题是：针对同时成像的两台线阵推扫的窄视场相机图像，通过基于大虚拟相机图像的几何纠正，实现全自动拼接处理，同时输出大虚拟相机图像对应的高精度有理多项式模型系数。

本发明的技术方案为一种基于大虚拟相机的窄视场双相机影像拼接方法，根据两个单相机建立大虚拟相机，所述大虚拟相机的焦距介于两个单相机的焦距之间，视场为两个单相机视场之和，主光轴位于两个单相机主光轴的中间；基于大虚拟相机进行以下步骤，

步骤1，根据两个单相机的几何成像参数，建立大虚拟相机的几何成像参数，所述几何成像参数包括相机参数和辅助数据，大虚拟相机的辅助数据和两个单相机的辅助数据相同，对大虚拟相机的相机参数建立包括求取本体坐标系下的视线方向如下，

设两个单相机分别记为相机A、B，相机A、B分别具有N₁、N₂个探元，相机A、B在本体坐标系上的重叠探元数为N₀，则大虚拟相机的总探元数=N₁+N₂-N₀，

设相机A的CCD的两端点A0,A1在本体坐标系中的平面投影坐标为（x_A0,y_A0）和（x_A1,y_A1），相机B的CCD的两端点B0,B1在本体坐标系中的平面投影坐标为（x_B0,y_B0）和（x_B1,y_B1），大虚拟相机的CCD两端点在本体坐标系下的平面投影坐标（x_C0,y_C0）（x_C1,y_C1）通过下式计算，

$x_{C0}=\frac{x_{A0}+x_{B0}}{2},$ y_C0＝y_A0, $x_{C1}=\frac{x_{A0}+x_{B0}}{2},$ y_C1＝y_B1

大虚拟相机的CCD的探元编号为s=0、1…N₁+N₂-N₀-1，经线性内插，根据下式得到大虚拟相机的CCD每个探元s在本体坐标系中的视线方向[x_C(s) y_C(s) 1]^T，

x_C(s)＝x_C0, $y_C\left(s\right)=y_{C0}+\frac{y_{C1}-y_{C0}}{N_1+N_2-N_0-1}s$

步骤2，根据两个单相机和大虚拟相机的几何成像参数，建立各自的几何成像模型；解算并输出大虚拟相机图像对应的有理多项式模型系数；根据基于几何成像模型的坐标正算过程和坐标反算过程对两个单相机的图像分别进行间接法几何纠正，得到大虚拟相机图像坐标系下的两个图像；

步骤3，对步骤2所得大虚拟相机图像坐标系下的两个图像进行基于坐标的拼接，得到拼接后的大虚拟相机图像。

而且，步骤1中，设卫星的本体坐标系o-X_bY_bZ_b中，原点o位于卫星质心，X_bY_bZ_b轴分别为卫星的滚动轴、俯仰轴和偏航轴，按以下方式求取相机A的CCD的两端点A0,A1在本体坐标系中的平面投影坐标（x_A0,y_A0）和（x_A1,y_A1），以及相机B的CCD的两端点B0,B1在本体坐标系中的平面投影坐标（x_B0,y_B0）和（x_B1,y_B1），

设相机A和B的内方位元素分别为a_0A,a_1A,a_2A,a_3A，b_0A,b_1A,b_2A,b_3A，a_0B,a_1B,a_2B,a_3B，b_0B,b_1B,b_2B,b_3B，则像点(s,l)的光线在相机A、B的相机坐标系中的矢量分别为[x_A(s) y_A(s) 1]^T、[x_B(s) y_B(s) 1]^T，

x_A(s)＝a_0A+a_1As+a_2As²+a_3As³   x_B(s)＝a_0B+a_1Bs+a_2Bs²+a_3Bs³

y_A(s)＝b_0A+b_1As+b_2As²+b_3As³   y_B(s)＝b_0B+b_1Bs+b_2Bs²+b_3Bs³

其中，变量s为探元的编号，变量l为图像的行号；

又设相机A和B相对于本体坐标系的安装矩阵分别为：R_BSa、R_BSb，则相机A和B的CCD端点A0,A1,B0,B1在本体坐标系中的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A0}}& \stackrel{\‾}{y_{A0}}& \stackrel{\‾}{z_{A0}}\end{array}\right]}^T,$ ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A1}}& \stackrel{\‾}{y_{A1}}& \stackrel{\‾}{z_{A1}}\end{array}\right]}^T,$ ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B0}}& \stackrel{\‾}{y_{B0}}& \stackrel{\‾}{z_{B0}}\end{array}\right]}^T,$ ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B1}}& \stackrel{\‾}{y_{B1}}& \stackrel{\‾}{z_{B1}}\end{array}\right]}^T$ 如下，

$\left[\begin{array}{c}\stackrel{\‾}{x_{A0}}\\ \stackrel{\‾}{y_{A0}}\\ \stackrel{\‾}{z_{A0}}\end{array}\right]=R_{\mathrm{BSa}}\left[\begin{array}{c}{x}_A\left(0\right)\\ y_A\left(0\right)\\ 1\end{array}\right],$ $\left[\begin{array}{c}\stackrel{\‾}{x_{A1}}\\ \stackrel{\‾}{y_{A1}}\\ \stackrel{\‾}{z_{A1}}\end{array}\right]=R_{\mathrm{BSa}}\left[\begin{array}{c}{x}_A(N_1-1)\\ y_A(N_1-1)\\ 1\end{array}\right],$ $\left[\begin{array}{c}\stackrel{\‾}{x_{B0}}\\ \stackrel{\‾}{y_{B0}}\\ \stackrel{\‾}{z_{B0}}\end{array}\right]=R_{\mathrm{BSb}}\left[\begin{array}{c}{x}_B\left(0\right)\\ y_B\left(0\right)\\ 1\end{array}\right],$ $\left[\begin{array}{c}\stackrel{\‾}{x_{B1}}\\ \stackrel{\‾}{y_{B1}}\\ \stackrel{\‾}{z_{B1}}\end{array}\right]=R_{\mathrm{BSb}}\left[\begin{array}{c}{x}_B(N_2-1)\\ y_B(N_2-1)\\ 1\end{array}\right]$

将本体坐标系下的光线矢量进行Z_b轴归一化，得到各端点的平面投影坐标（x_A0,y_A0）（x_A1,y_A1）（x_B0,y_B0）（x_B1,y_B1）如下，

$x_{A0}=\frac{\stackrel{\‾}{x_{A0}}}{\stackrel{\‾}{z_{A0}}},$ $y_{A0}=\frac{\stackrel{\‾}{y_{A0}}}{\stackrel{\‾}{z_{A0}}},$ $x_{A1}=\frac{\stackrel{\‾}{x_{A1}}}{\stackrel{\‾}{z_{A1}}},$ $y_{A1}=\frac{\stackrel{\‾}{y_{A1}}}{\stackrel{\‾}{z_{A1}}},$ $x_{B0}=\frac{\stackrel{\‾}{x_{B0}}}{\stackrel{\‾}{z_{B0}}},$ $y_{B0}=\frac{\stackrel{\‾}{y_{B0}}}{\stackrel{\‾}{z_{B0}}},$ $x_{B1}=\frac{\stackrel{\‾}{x_{B1}}}{\stackrel{\‾}{z_{B1}}},$ $y_{B1}=\frac{\stackrel{\‾}{y_{B1}}}{\stackrel{\‾}{z_{B1}}}$

其中，x_A0、y_A0表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A0}}& \stackrel{\‾}{y_{A0}}& \stackrel{\‾}{z_{A0}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量；x_A1、y_A1表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A1}}& \stackrel{\‾}{y_{A1}}& \stackrel{\‾}{z_{A1}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量；x_B0、y_B0表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B0}}& \stackrel{\‾}{y_{B0}}& \stackrel{\‾}{z_{B0}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量；x_B1、y_B1表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B1}}& \stackrel{\‾}{y_{B1}}& \stackrel{\‾}{z_{B1}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量。

而且，步骤2中，基于几何成像模型的坐标正算过程实现如下，

设R_t、R_GF、R_FB、R_BSa、R_BSb分别为t时刻下从J2000惯性坐标系到地固地心直角坐标系的旋转矩阵、从轨道坐标系到J2000惯性坐标系的旋转矩阵、从本体坐标系到轨道坐标系的旋转矩阵、从相机坐标系A到本体坐标系的旋转矩阵、从相机坐标系B到本体坐标系的旋转矩阵,且[X_t Y_t Z_t]^T为t时刻卫星质心在地固地心直角坐标系下的坐标矢量，

相机A的坐标正算过程为，对于相机A图像上的某像点(s,l)，相机坐标系下的光线矢量[x_A(s) y_A(s) 1]^T，则物方坐标[X Y Z]^T由下列共线方程给出，

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+m{R}_t{R}_{\mathrm{GF}}{R}_{\mathrm{FB}}{R}_{\mathrm{BSa}}\left[\begin{array}{c}{x}_A\left(s\right)\\ y_A\left(s\right)\\ 1\end{array}\right]$

相机B的坐标正算过程为，对于相机B图像上的某像点(s,l)，相机坐标系下的光线矢量[x_B(s) y_B(s) 1]^T，则物方坐标[X Y Z]^T由下列共线方程给出，

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+m{R}_t{R}_{\mathrm{GF}}{R}_{\mathrm{FB}}{R}_{\mathrm{BSb}}\left[\begin{array}{c}{x}_B\left(s\right)\\ y_B\left(s\right)\\ 1\end{array}\right]$

大虚拟相机坐标正算过程为，对于大虚拟相机图像上的某像点(s,l)，首先求出本体坐标系下的光线矢量[x_C(s) y_C(s) 1]^T，则物方坐标[X Y Z]^T由下列共线方程给出：

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+m{R}_t{R}_{\mathrm{GF}}{R}_{\mathrm{FB}}\left[\begin{array}{c}{x}_C\left(s\right)\\ y_C\left(s\right)\\ 1\end{array}\right]$

其中，m为比例系数。

而且，步骤2中，基于几何成像模型的坐标反算过程根据有理多项式模型实现，相机A、B和大虚拟相机坐标的有理多项式模型的系数分别通过以下步骤得到，

首先，对相机所得图像分别划分规则格网，对物方高程划分多个高程面，利用相机的几何成像模型的坐标正算过程，计算所有虚拟立体格网点的物方坐标（X,Y,Z），并转换为WGS84地理经纬度坐标（B,L,H）；

然后，将上述各虚拟立体格网点作为控制点，对有理多项式模型的系数列出误差方程式，基于最小二乘准则迭代求解，解算出有理多项式模型的系数。

而且，步骤2中，两个单相机的图像分别进行间接法几何纠正的实现方式为，以大虚拟相机图像坐标系为坐标系②，以相机A和B的图像坐标系分别作为坐标系①执行以下处理，

1)将坐标系①的原始图像四个角点坐标，基于相应单相机的几何成像模型的坐标正算过程和大虚拟相机的几何成像模型的坐标反算过程，得到输出图像在坐标系②中的范围；

2)对坐标系②输出图像范围内的每个像素，通过大虚拟相机的几何成像模型的坐标正算过程和相应单相机的几何成像模型的坐标反算过程，得到在坐标系①的原始图像坐标；

3)根据坐标系①中的原始图像坐标，通过灰度重采样，得到坐标系②相应输出图像上的每个像素的灰度值。

本发明针对线阵推扫的窄视场双相机影像，构造一个线阵推扫的大虚拟相机，旨在实现双相机图像拼接的同时，还能提供拼接图像的高精度几何模型。通过以上步骤就可以实现窄视场双相机图像的高精度拼接。该方法借助大虚拟相机概念，使拼接图像具有与之对应的高精度有理多项式模型，从而简化了后续图像几何处理与应用。该方法适用于同时推扫成像的窄视场双相机图像，前提条件是需要对两个单相机进行严格的几何标定，处理过程是全自动的、无需人工干预。这种双相机拼接工作可以在地面预处理阶段通过全自动批处理完成，因此，可以减少后续处理的工作量。

附图说明

图1为本发明实施例的两个单相机图像与大虚拟相机图像关系示意图；

图2为本发明实施例的间接法几何纠正流程示意图。

图3为本发明实施例的流程图。

具体实施方式

以下结合附图和实施例详细说明本发明技术方案。

参见图3，实施例针对同时推扫成像的窄视场双相机图像，执行步骤如下，可采用计算机软件技术实现自动运行流程：

步骤1，根据两个单相机的几何成像参数，建立大虚拟相机的几何成像参数。

本发明提出大虚拟相机概念，并计算每个探元在本体坐标系的视线方向。

首先介绍卫星本体坐标系与单相机坐标系。

卫星本体坐标系o-X_bY_bZ_b：原点o位于卫星质心，X_bY_bZ_b轴分别为卫星的滚动轴、俯仰轴和偏航轴。

单相机坐标系o-X_cY_cZ_c：原点o位于卫星质心，Z_c轴指向主光轴方向，与卫星本体坐标系的Z_b轴有一个夹角，称为安装角。X_cY_c平面为焦平面，与Z_c轴垂直，X_cY_cZ_c构成右手系。其中：X_c指向卫星飞行方向，Y_c轴指向CCD方向。

两个单相机的关系可以从两个方面来看。1）两个单相机的焦距、CCD探元数、视场角、CCD线阵在焦平面的位置等相机参数都比较接近，因此，用三次多项式表示的单相机内方位元素也都比较接近。2）两个单相机在卫星本体坐标系下的安装角方向相反、大小相近，两个主光轴之间的夹角略小于单相机的视场角，这样得到的两个单相机图像在地面覆盖范围有N₀个像素的重叠。设两个单相机分别记为相机A、B。

对于大虚拟相机，假设它是安装在卫星本体坐标系下的，焦距介于两个单相机的焦距之间（因为一般相机A、B焦距会有差异，可等于相机A、B焦距的均值），视场则扩大到两个单相机视场之和（相当于CCD探元数翻了一倍），主光轴位于两个单相机主光轴的中间（相当于安装角接近于0）。因此，从大虚拟相机坐标系到本体坐标系的安装矩阵可以假设为单位矩阵，这样，大虚拟相机坐标系下的每个像素的视线方向就等于本体坐标系下的视线方向。

几何成像参数包括相机参数和辅助数据。由于两个单相机是同时推扫成像，因此，拥有相同的辅助数据，包括成像时间范围、卫星星历和卫星姿态数据等。对于大虚拟相机，直接使用辅助数据即可。

相机参数包括相机坐标系下视线方向和相机坐标系与卫星本体坐标系的安装矩阵，两者合在一起就是本体坐标系下的视线方向。为了求出大虚拟相机在本体坐标系下的视线方向，需要将两个单相机在相机坐标系下的视线方向变换到本体坐标系下，并投影到一个平面上。

设相机A、B分别具有N₁、N₂个探元，A和B内方位元素分别为：a_0A,a_1A,a_2A,a_3A，b_0A,b_1A,b_2A,b_3A，a_0B,a_1B,a_2B,a_3B，b_0B,b_1B,b_2B,b_3B，则像点(s,l)的光线在A、B相机坐标系中的矢量分别为：[x_A(s) y_A(s) 1]^T、[x_B(s) y_B(s) 1]^T：

x_A(s)＝a_0A+a_1As+a_2As²+a_3As³

y_A(s)＝b_0A+b_1As+b_2As²+b_3As³

（1）

x_B(s)＝a_0B+a_1Bs+a_2Bs²+a_3Bs³

y_B(s)＝b_0B+b_1Bs+b_2Bs²+b_3Bs³

（2）

其中：变量s为探元的编号，变量l为图像的行号。x_A(s)、_yA(s)分别代表A相机坐标系中沿轨、垂轨方向瞬时视场角的正切值，x_B(s)、_yB(s)分别代表B相机坐标系中沿轨、垂轨方向瞬时视场角的正切值。相机A的CCD的探元编号为s=0、1…N₁-1，即s=0、s=N₁-1对应于相机A的CCD的两端点A0,A1，相机B的CCD的探元编号为s=0、1…N₁-1，即s=0、s=N₂-1对应于相机B的CCD的两端点B0,B1。

又设A和B相机相对于本体坐标系的安装矩阵分别为：R_BSa、R_BSb，则可计算出A和B相机的CCD端点A0,A1,B0,B1在本体坐标系中的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A0}}& \stackrel{\‾}{y_{A0}}& \stackrel{\‾}{z_{A0}}\end{array}\right]}^T,$ ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A1}}& \stackrel{\‾}{y_{A1}}& \stackrel{\‾}{z_{A1}}\end{array}\right]}^T,$

${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B0}}& \stackrel{\‾}{y_{B0}}& \stackrel{\‾}{z_{B0}}\end{array}\right]}^T,$ ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B1}}& \stackrel{\‾}{y_{B1}}& \stackrel{\‾}{z_{B1}}\end{array}\right]}^T:$

$\left[\begin{array}{c}\stackrel{\‾}{x_{A0}}\\ \stackrel{\‾}{y_{A0}}\\ \stackrel{\‾}{z_{A0}}\end{array}\right]=R_{\mathrm{BSa}}\left[\begin{array}{c}{x}_A\left(0\right)\\ y_A\left(0\right)\\ 1\end{array}\right],$ $\left[\begin{array}{c}\stackrel{\‾}{x_{A1}}\\ \stackrel{\‾}{y_{A1}}\\ \stackrel{\‾}{z_{A1}}\end{array}\right]=R_{\mathrm{BSa}}\left[\begin{array}{c}{x}_A(N_1-1)\\ y_A(N_1-1)\\ 1\end{array}\right],$ $\left[\begin{array}{c}\stackrel{\‾}{x_{B0}}\\ \stackrel{\‾}{y_{B0}}\\ \stackrel{\‾}{z_{B0}}\end{array}\right]=R_{\mathrm{BSb}}\left[\begin{array}{c}{x}_B\left(0\right)\\ y_B\left(0\right)\\ 1\end{array}\right],$ $\left[\begin{array}{c}\stackrel{\‾}{x_{B1}}\\ \stackrel{\‾}{y_{B1}}\\ \stackrel{\‾}{z_{B1}}\end{array}\right]=R_{\mathrm{BSb}}\left[\begin{array}{c}{x}_B(N_2-1)\\ y_B(N_2-1)\\ 1\end{array}\right]---\left(3\right)$

再将本体坐标系中的光线矢量投影到平面Z_b＝1上，即将本体坐标系下的光线矢量进行Z_b轴归一化，得到X_bY_b两个方向的分量，即各端点的平面投影坐标（x_A0,y_A0）（x_A1,y_A1）（x_B0,y_B0）（x_B1,y_B1）如下：

$x_{A0}=\frac{\stackrel{\‾}{x_{A0}}}{\stackrel{\‾}{z_{A0}}},$ $y_{A0}=\frac{\stackrel{\‾}{y_{A0}}}{\stackrel{\‾}{z_{A0}}},$ $x_{A1}=\frac{\stackrel{\‾}{x_{A1}}}{\stackrel{\‾}{z_{A1}}},$ $y_{A1}=\frac{\stackrel{\‾}{y_{A1}}}{\stackrel{\‾}{z_{A1}}},$ $x_{B0}=\frac{\stackrel{\‾}{x_{B0}}}{\stackrel{\‾}{z_{B0}}},$ $y_{B0}=\frac{\stackrel{\‾}{y_{B0}}}{\stackrel{\‾}{z_{B0}}},$ $x_{B1}=\frac{\stackrel{\‾}{x_{B1}}}{\stackrel{\‾}{z_{B1}}},$ $y_{B1}=\frac{\stackrel{\‾}{y_{B1}}}{\stackrel{\‾}{z_{B1}}}---\left(4\right)$

其中，x_A0、y_A0表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A0}}& \stackrel{\‾}{y_{A0}}& \stackrel{\‾}{z_{A0}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量；x_A1、y_A1表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{A1}}& \stackrel{\‾}{y_{A1}}& \stackrel{\‾}{z_{A1}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量；x_B0、y_B0表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B0}}& \stackrel{\‾}{y_{B0}}& \stackrel{\‾}{z_{B0}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量；x_B1、y_B1表示本体坐标系下的光线矢量 ${\left[\begin{array}{ccc}\stackrel{\‾}{x_{B1}}& \stackrel{\‾}{y_{B1}}& \stackrel{\‾}{z_{B1}}\end{array}\right]}^T$ 在X_bY_b坐标轴方向上的分量。

如图1所示，设大虚拟相机的CCD位于A和B相机中间，则其两端点C0,C1在本体坐标系下的平面投影坐标（x_C0,y_C0）（x_C1,y_C1）可通过下式计算：

$x_{C0}=\frac{x_{A0}+x_{B0}}{2},$ y_C0＝y_A0, $x_{C1}=\frac{x_{A0}+x_{B0}}{2},$ y_C1＝y_B1   （5）

由于镜头畸变等原因，A、B两个单相机的CCD不一定严格排列在一条直线上，但这并不影响，大虚拟相机的CCD仍然可以排列在一条直线上。设两个单相机在本体坐标系(平面Z_b＝1)上的重叠探元数为N₀，则大虚拟相机的总探元数=N₁+N₂-N₀，大虚拟相机的CCD的探元编号为s=0、1…N₁+N₂-N₀-1，经线性内插，可得到大虚拟相机的每个探元s在本体坐标系中的视线方向[x_C(s) y_C(s) 1]^T：

x_C(s)＝x_C0, $y_C\left(s\right)=y_{C0}+\frac{y_{C1}-y_{C0}}{N_1+N_2-N_0-1}s---\left(6\right)$

步骤2，对两个单相机图像进行几何纠正，得到大虚拟相机的图像坐标系下的两个图像，同时解算并输出大虚拟相机图像对应的有理多项式模型系数，包括以下子步骤，

步骤2.1，根据两个单相机和大虚拟相机的几何成像参数，建立各自的几何成像模型，包括：坐标正算（从各自的图像坐标系坐标到物方平均高程面坐标）和坐标反算（从物方平均高程面坐标到各自的图像坐标系坐标）两个过程。

基于几何成像模型的坐标正、反算过程均为现有技术。一般坐标正算过程基于共线方程模型实现，而坐标反算过程基于有理多项式模型实现。以下先描述共线方程，再在此基础上计算虚拟立体格网点坐标，用以解算有理多项式模型系数。其中单相机图像的坐标反算过程不限于有理多项式模型、也可以基于共线方程模型实现；而大虚拟相机图像对应的有理多项式模型系数是需要解算并输出的。

假设R_t、R_GF、R_FB、R_BSa(R_BSb)分别为t时刻下从J2000惯性坐标系到地固地心直角坐标系的旋转矩阵、从轨道坐标系到J2000惯性坐标系的旋转矩阵、从本体坐标系到轨道坐标系的旋转矩阵、从A(或B)相机坐标系到本体坐标系的旋转矩阵,且[X_t Y_t Z_t]^T为t时刻卫星质心在地固地心直角坐标系下的坐标矢量。

单相机坐标正算过程如下：对于A(或B)相机图像上的某像点(s,l)，先根据步骤1中式（1）或（2）所得相机坐标系下的光线矢量[x_A(s) y_A(s) 1]^T(或[_xB(s) y_B(s) 1]^T)，则物方坐标[X Y Z]^T由下列共线方程给出：

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+m{R}_t{R}_{\mathrm{GF}}{R}_{\mathrm{FB}}{R}_{\mathrm{BSa}}\left[\begin{array}{c}{x}_A\left(s\right)\\ y_A\left(s\right)\\ 1\end{array}\right]$ 或 $\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+m{R}_t{R}_{\mathrm{GF}}{R}_{\mathrm{FB}}{R}_{\mathrm{BSb}}\left[\begin{array}{c}{x}_B\left(s\right)\\ y_B\left(s\right)\\ 1\end{array}\right]---\left(7\right)$

其中，m为比例系数，可利用现有技术，通过光线与参考椭球的平均高程面H相交，求解一元二次多项式，得到m的值，再代入上式，求得物方坐标[X Y Z]^T。

大虚拟相机坐标正算过程如下：对于大虚拟相机图像上的某像点(s,l)，根据步骤1中式（6）所得本体坐标系下的光线矢量[x_C(s) y_C(s) 1]^T，则物方坐标[X Y Z]^T由下列共线方程给出：

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+m{R}_t{R}_{\mathrm{GF}}{R}_{\mathrm{FB}}\left[\begin{array}{c}{x}_C\left(s\right)\\ y_C\left(s\right)\\ 1\end{array}\right]---\left(8\right)$

同样，m为比例系数，可利用现有技术，通过光线与参考椭球的平均高程面H相交，得到m的值，再代入上式，求得物方坐标[X Y Z]^T。

由物方坐标反算像方坐标的过程可以通过有理多项式模型计算，属于现有技术。其中有理多项式模型的系数可以利用虚拟立体格网点，将其作为控制点，用最小二乘法求解得到。虚拟立体格网点的坐标计算和有理多项式模型的系数解算都是现有技术。对相机A、B、大虚拟相机，可以分别采用以下两步求取系数：。

首先，对相机所得图像分别划分规则格网，对物方高程划分多个高程面，利用相机的几何成像模型，即坐标正算公式（7）或（8），计算所有虚拟立体格网点的物方坐标（X,Y,Z），并转换为WGS84地理经纬度坐标（B,L,H）；

然后，将上述各虚拟立体格网点作为控制点，对有理多项式模型的系数列出误差方程式，基于最小二乘准则迭代求解，可解算出有理多项式模型的系数。解释过程可参见现有最小二乘准则，为便于实施参考起见，提供解算过程如下：

1）设当前的迭代次数为k，令k的初始值为1，

设有理多项式模型为：

$s^{\′}=\frac{{\mathrm{Num}}_S(U,V,W)}{{\mathrm{Dem}}_S(U,V,W)},$ $l^{\′}=\frac{{\mathrm{Num}}_L(U,V,W)}{{\mathrm{Den}}_L(U,V,W)}---\left(9\right)$

其中，

Num_S(U,V,W)＝a₁+a₂V+a₃W+a₄U+a₅VU+a₆VW+a₇UW+a₈V²+

a₉U²+a₁₀W²+a₁₁UVW+a₁₂V³+a₁₃VU²+a₁₄VW²

+a₁₅V²U+a₁₆U³+a₁₇UW²+a₁₈WV²+a₁₉U²W+a₂₀W³

Den_S(U,V,W)＝b₁+b₂V+b₃W+b₄U+b₅VU+b₆VW+b₇UW+b₈V²+

b₉U²+b₁₀W²+b₁₁UVW+b₁₂V³+b₁₃VU²+b₁₄VW²

+b₁₅V²U+b₁₆U³+b₁₇UW²+b₁₈WV²+b₁₉U²W+b₂₀W³

Num_L(U,V,W)＝c₁+c₂V+c₃W+c₄U+c₅VU+c₆VW+c₇UW+c₈V²+

c₉U²+c₁₀W²+c₁₁UVW+c₁₂V³+c₁₃VU²+c₁₄VW²

+c₁₅V²U+c₁₆U³+c₁₇UW²+c₁₈WV²+c₁₉U²W+c₂₀W³

Den_L(U,V,W)＝d₁+d₂V+d₃W+d₄U+d₅VU+d₆VW+d₇UW+d₈V²+

d₉U²+d₁₀W²+d₁₁UVW+d₁₂V³+d₁₃VU²+d₁₄VW²

+d₁₅V²U+d₁₆U³+d₁₇UW²+d₁₈WV²+d₁₉U²W+d₂₀W³

Num_S(U,V,W)表示列方向的分子多项式，Den_S(U,V,W)表示列方向的分母多项式，Num_L(U,V,W)表示行方向的分子多项式，Den_L(U,V,W)表示行方向的分子多项式。

多项式中的a_j,b_j,c_j,d_j,(j＝1,…,20)为有理多项式模型的系数，（s′,l′）为归一化的图像坐标，（U,V,W）为归一化的经纬度坐标。

设共有n个虚拟立体格网控制点，其中第i个控制点的图像坐标(s_i,l_i)，WGS84地理经纬度坐标为（B_i,L_i,H_i），则第i个控制点的归一化的图像坐标（s_i′,l_i′）和归一化的经纬度坐标（U_i,V_i,W_i）可通过下式求出：

$s_i^{\′}=\frac{s_i-\mathrm{SampleOff}}{\mathrm{SampleScale}},$ $l_i^{\′}=\frac{l_i-\mathrm{LineOff}}{\mathrm{LineScale}}$ $U_i=\frac{B_i-\mathrm{LatOff}}{\mathrm{LatScale}},$ $V_i=\frac{L_i-\mathrm{LonOff}}{\mathrm{LonScale}},$ $W_i=\frac{H_i-\mathrm{HeiOff}}{\mathrm{HeiScale}}---\left(10\right)$

i＝1,…,n

其中：平移系数SampleOff,LineOff,LatOff,LonOff,HeiOff为虚拟立体格网控制点的坐标均值；缩放系数SampleScale,LineScale,LatScale,LonScale,HeiScale为虚拟立体格网控制点的最大坐标值与最小坐标值之差。

2）为了求解未知数（即：有理多项式模型的系数）a_j,b_j,c_j,d_j,(j＝1,…,20)，先对公式（9）变换形式，得到有理多项式模型的另一种表达形式（11）：

F_S＝s′×Den_S(U,V,W)-Num_S(U,V,W)=0   （11）

F_L＝l′×Den_L(U,V,W)-Num_L(U,V,W)=0

3）更新有理多项式模型的系数：

当k=1时，未知数a_j,b_j,c_j,d_j,(j＝1,…,20)的初值a_j ⁰,b_j ⁰,c_j ⁰,d_j ⁰根据经验给出，当k＞1时，未知数a_j,b_j,c_j,d_j,(j＝1,…,20)的初值为上一次迭代结果a_j ^k-1,b_j ^k-1,c_j ^k-1,d_j ^k-1。

公式（11）按泰勒级数展开，取一次项得：

其中：由未知数初值代入公式（11）求得（是关于控制点坐标的函数），第一次执行到3）时，k=1，采用a_j ⁰,b_j ⁰,c_j ⁰,d_j ⁰作为第k次迭代的初值代入公式（11），后续执行到3）时采用上一轮迭代（第k-1次迭代）所得结果作为本次迭代的初值代入公式（11）。未知数的初值a_j ⁰,b_j ⁰,c_j ⁰,d_j ⁰可由本领域技术人员预先设置。为公式（11）中的F_S,F_L分别对有理多项式模型系数a_j,b_j,c_j,d_j求偏导数得到的系数（也是关于控制点坐标的函数）。da_j,db_j,dc_j,dd_j表示a_j,b_j,c_j,d_j的误差。

设va_j,vb_j,vc_j,vd_j为未知数a_j,b_j,c_j,d_j的改正数，将每个控制点的归一化坐标代入（12）式，可建立如下误差方程式：

A_2n×80X_80×1=L_2n×1   （13）

其中：n为控制点数，A_2n×80为系数矩阵，X_80×1为有理多项式模型系数的改正数，L_2n×1为常数向量。

X_80×1=[va₁ va₂ … vb₁ vb₂ … vc₁ vc₂ … vd₁ vd₂ …]^T

其中，为第i个控制点的归一化坐标代入得到的值；

（i＝1,…,n）为第i个控制点的归一化坐标代入得到的值。

基于最小二乘准则，可解算出有理多项式模型系数的改正数：

X_80×1＝（A_2n×80 ^TA_2n×80）^-1（A_2n×80 ^TL_2n×1）   （15）

更新有理多项式模型的系数，得到本次迭代的结果a_j ^k,^b _j ^k,c_j ^k,^d _j ^k,(j＝1,…,20)：

a_j ^k＝a_j ^k-1+va_j,b_j ^k＝b_j ^k-1+vb_j,c_j ^k＝c_j ^k-1+vc_j,d_j ^k＝d_j ^k-1+vd_j

4）判断，当前迭代次数k是否达到一预设定值或有理多项式模型系数的改正数向量的每个分量的绝对值是否都小于预设阈值，“是”则结束迭代；“否”则令k=k+1，返回3）进行下一次迭代。

通过上述流程可解算出A、B单相机的有理多项式模型系数和大虚拟相机的有理多项式模型系数，从而实现基于有理多项式模型的坐标反算（已知物方坐标求图像坐标的过程）。

具体实施时，坐标反算也可以采用其他方式，例如单相机可采用共线方程模型。

步骤2.2，对单相机图像进行间接法几何纠正，得到大虚拟相机图像坐标系下的两个图像。

间接法几何纠正属于现有技术，以相机A和B的图像坐标系分别作为坐标系①执行处理，处理流程如图2所示，以下分三步简述从坐标系①纠正到坐标系②的过程（坐标系①、②分别为单相机图像坐标系和大虚拟相机图像坐标系）：

1)将坐标系①的原始图像四个角点坐标，通过步骤2.1的单相机的坐标正算和大虚拟相机的坐标反算，得到输出图像在坐标系②中的范围；

2)对坐标系②输出图像范围内的每个像素，通过步骤2.1的大虚拟相机的坐标正算和单相机的坐标反算，得到其在坐标系①的原始图像坐标；

3)根据坐标系①中的原始图像坐标，通过灰度重采样，得到坐标系②输出图像上的每个像素的灰度值。

图2中，表示原始单相机图像，表示输出大虚拟相机坐标系下的单相机图像，粗箭头表示从坐标系①到坐标系②的坐标正算过程，细箭头表示从坐标系②到坐标系①的坐标反算过程；由于反算到原始单相机图像上的坐标不是整像素值，图中dx和dy表示坐标的小数部分，根据反算的坐标进行灰度重采样，再将内插得到的灰度值赋给输出图像上的像素。

步骤3，对大虚拟相机图像坐标系下的两个输出图像进行基于坐标的拼接，得到拼接后的大虚拟相机图像。基于坐标的拼接具体实现可采用现有技术。

通过以上步骤就可以实现窄视场双相机图像的拼接并输出高精度的有理多项式模型系数。该方法借助大虚拟相机的概念，使拼接图像具有与之对应的高精度有理多项式模型，不仅简化了后续图像几何处理与应用，而且处理过程是全自动的、无需人工干预。本发明方法中的大虚拟相机参数不需要对用户公开，提供给用户的是有理多项式模型系数，这使得图像产品具备较好的保密性和通用性。本发明方法适用于同时推扫成像的窄视场双相机图像的地面预处理，前提条件是需要对两个单相机进行严格的几何标定，以保证两个单相机图像的相对几何精度。

本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (5)

1.一种基于大虚拟相机的窄视场双相机影像拼接方法，其特征在于：根据两个单相机建立大虚拟相机，所述大虚拟相机的焦距介于两个单相机的焦距之间，视场为两个单相机视场之和，主光轴位于两个单相机主光轴的中间；基于大虚拟相机进行以下步骤，

步骤1，根据两个单相机的几何成像参数，建立大虚拟相机的几何成像参数，所述几何成像参数包括相机参数和辅助数据，大虚拟相机的辅助数据和两个单相机的辅助数据相同，对大虚拟相机的相机参数建立包括求取卫星的本体坐标系下的视线方向如下，

设两个单相机分别记为相机A、B，相机A、B分别具有N₁、N₂个探元，相机A、B在本体坐标系上的重叠探元数为N₀，则大虚拟相机的总探元数＝N₁+N₂-N₀，

设相机A的CCD的两端点A0,A1在本体坐标系中的平面投影坐标为(x_A0,y_A0)和(x_A1,y_A1)，相机B的CCD的两端点B0,B1在本体坐标系中的平面投影坐标为(x_B0,y_B0)和(x_B1,y_B1)，大虚拟相机的CCD两端点在本体坐标系下的平面投影坐标(x_C0,y_C0)(x_C1,y_C1)通过下式计算，

$x_{C0}=\frac{x_{A0}+x_{B0}}{2},y_{C0}=y_{A0},x_{C1}=\frac{x_{A0}+x_{B0}}{2},y_{C1}=y_{B1}$

大虚拟相机的CCD的探元编号为s＝0、1…N₁+N₂-N₀-1，经线性内插，根据下式得到大虚拟相机的CCD每个探元s在本体坐标系中的视线方向[x_C(s) y_C(s) 1]^T，

$x_C\left(s\right)=x_{C0},y_C\left(s\right)=y_{C0}+\frac{y_{C1}-y_{C0}}{N_1+N_2-N_0-1}s$

步骤2，根据两个单相机和大虚拟相机的几何成像参数，建立各自的几何成像模型；解算并输出大虚拟相机图像对应的有理多项式模型系数；根据基于几何成像模型的坐标正算过程和坐标反算过程对两个单相机的图像分别进行间接法几何纠正，得到大虚拟相机图像坐标系下的两个图像；

步骤3，对步骤2所得大虚拟相机图像坐标系下的两个图像进行基于坐标的拼接，得到拼接后的大虚拟相机图像。

2.根据权利要求1所述基于大虚拟相机的窄视场双相机影像拼接方法，其特征在于：步骤1中，设卫星的本体坐标系o-X_bY_bZ_b中，原点o位于卫星质心，X_b、Y_b、Z_b轴分别为卫星的滚动轴、俯仰轴和偏航轴，按以下方式求取相机A的CCD的两端点A0,A1在本体坐标系中的平面投影坐标(x_A0,y_A0)和(x_A1,y_A1)，以及相机B的CCD的两端点B0,B1在本体坐标系中的平面投影坐标(x_B0,y_B0)和(x_B1,y_B1)，

设相机A和B的内方位元素分别为a_0A,a_1A,a_2A,a_3A，b_0A,b_1A,b_2A,b_3A，a_0B,a_1B,a_2B,a_3B，b_0B,b_1B,b_2B,b_3B，则像点(s,l)的光线在相机A、B的相机坐标系中的矢量分别为[x_A(s) y_A(s) 1]^T、[x_B(s) y_B(s) 1]^T，

x_A(s)＝a_0A+a_1As+a_2As²+a_3As³    x_B(s)＝a_0B+a_1Bs+a_2Bs²+a_3Bs³

y_A(s)＝b_0A+b_1As+b_2As²+b_3As³    y_B(s)＝b_0B+b_1Bs+b_2Bs²+b_3Bs³

其中，变量s为探元的编号，变量l为图像的行号；

又设相机A和B相对于本体坐标系的安装矩阵分别为：R_BSa、R_BSb，则相机A和B的CCD端点A0,A1,B0,B1在本体坐标系中的光线矢量 如下，

$\left[\begin{array}{c}\stackrel{\‾}{x_{A0}}\\ \stackrel{\‾}{y_{A0}}\\ \stackrel{\‾}{z_{A0}}\end{array}\right]=R_{BSa}\left[\begin{array}{c}{x}_A\left(0\right)\\ y_A\left(0\right)\\ 1\end{array}\right],\left[\begin{array}{c}\stackrel{\‾}{x_{A1}}\\ \stackrel{\‾}{y_{A1}}\\ \stackrel{\‾}{z_{A1}}\end{array}\right]=R_{BSa}\left[\begin{array}{c}{x}_A(N_1-1)\\ y_A(N_1-1)\\ 1\end{array}\right],\left[\begin{array}{c}\stackrel{\‾}{x_{B0}}\\ \stackrel{\‾}{y_{B0}}\\ \stackrel{\‾}{z_{B0}}\end{array}\right]=R_{BSb}\left[\begin{array}{c}{x}_B\left(0\right)\\ y_B\left(0\right)\\ 1\end{array}\right],\left[\begin{array}{c}\stackrel{\‾}{x_{B1}}\\ \stackrel{\‾}{y_{B1}}\\ \stackrel{\‾}{z_{B1}}\end{array}\right]=R_{BSb}\left[\begin{array}{c}{x}_B(N_2-1)\\ y_B(N_2-1)\\ 1\end{array}\right]$

将本体坐标系下的光线矢量进行Z_b轴归一化，得到各端点的平面投影坐标(x_A0,y_A0)(x_A1,y_A1)(x_B0,y_B0)(x_B1,y_B1)如下，

$x_{A0}=\frac{\stackrel{\‾}{x_{A0}}}{\stackrel{\‾}{z_{A0}}},y_{A0}=\frac{\stackrel{\‾}{y_{A0}}}{\stackrel{\‾}{z_{A0}}},x_{A1}=\frac{\stackrel{\‾}{x_{A1}}}{\stackrel{\‾}{z_{A1}}},y_{A1}=\frac{\stackrel{\‾}{y_{A1}}}{\stackrel{\‾}{z_{A1}}},x_{B0}=\frac{\stackrel{\‾}{x_{B0}}}{\stackrel{\‾}{z_{B0}}},y_{B0}=\frac{\stackrel{\‾}{y_{B0}}}{\stackrel{\‾}{z_{B0}}},x_{B1}=\frac{\stackrel{\‾}{x_{B1}}}{\stackrel{\‾}{z_{B1}}},y_{B1}=\frac{\stackrel{\‾}{y_{B1}}}{\stackrel{\‾}{z_{B1}}}$

其中，x_A0、y_A0表示本体坐标系下的光线矢量在X_b、Y_b坐标轴方向上的分量；x_A1、y_A1表示本体坐标系下的光线矢量在X_b、Y_b坐标轴方向上的分量；x_B0、y_B0表示本体坐标系下的光线矢量在X_b、Y_b坐标轴方向上的分量；x_B1、y_B1表示本体坐标系下的光线矢量在X_b、Y_b坐标轴方向上的分量。

3.根据权利要求2所述基于大虚拟相机的窄视场双相机影像拼接方法，其特征在于：步骤2中，基于几何成像模型的坐标正算过程实现如下，

设R_t、R_GF、R_FB、R_BSa、R_BSb分别为t时刻下从J2000惯性坐标系到地固地心直角坐标系的旋转矩阵、从轨道坐标系到J2000惯性坐标系的旋转矩阵、从本体坐标系到轨道坐标系的旋转矩阵、从相机坐标系A到本体坐标系的旋转矩阵、从相机坐标系B到本体坐标系的旋转矩阵,且[X_t Y_t Z_t]^T为t时刻卫星质心在地固地心直角坐标系下的坐标矢量，

相机A的坐标正算过程为，对于相机A图像上的某像点(s,l)，相机坐标系下的光线矢量[x_A(s) y_A(s) 1]^T，则物方坐标[X Y Z]^T由下列共线方程给出，

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+{\mathrm{mR}}_t{R}_{GF}{R}_{FB}{R}_{BSa}\left[\begin{array}{c}{x}_A\left(s\right)\\ y_A\left(s\right)\\ 1\end{array}\right]$

相机B的坐标正算过程为，对于相机B图像上的某像点(s,l)，相机坐标系下的光线矢量[x_B(s) y_B(s) 1]^T，则物方坐标[X Y Z]^T由下列共线方程给出，

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+{\mathrm{mR}}_t{R}_{GF}{R}_{FB}{R}_{BSb}\left[\begin{array}{c}{x}_B\left(s\right)\\ y_B\left(s\right)\\ 1\end{array}\right]$

大虚拟相机坐标正算过程为，对于大虚拟相机图像上的某像点(s,l)，首先求出本体坐标系下的光线矢量[x_C(s) y_C(s) 1]^T，则物方坐标[X Y Z]^T由下列共线方程给出：

$\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]=\left[\begin{array}{c}{X}_t\\ Y_t\\ Z_t\end{array}\right]+{\mathrm{mR}}_t{R}_{GF}{R}_{FB}\left[\begin{array}{c}{x}_C\left(s\right)\\ y_C\left(s\right)\\ 1\end{array}\right]$

其中，m为比例系数。

4.根据权利要求3所述基于大虚拟相机的窄视场双相机影像拼接方法，其特征在于：步骤2中，基于几何成像模型的坐标反算过程根据有理多项式模型实现，相机A、B和大虚拟相机坐标的有理多项式模型的系数分别通过以下步骤得到，

首先，对相机所得图像分别划分规则格网，对物方高程划分多个高程面，利用相机的几何成像模型的坐标正算过程，计算所有虚拟立体格网点的物方坐标(X,Y,Z)，并转换为WGS84地理经纬度坐标(B,L,H)；

然后，将上述各虚拟立体格网点作为控制点，对有理多项式模型的系数列出误差方程式，基于最小二乘准则迭代求解，解算出有理多项式模型的系数。

5.根据权利要求4所述基于大虚拟相机的窄视场双相机影像拼接方法，其特征在于：步骤2中，两个单相机的图像分别进行间接法几何纠正的实现方式为，以大虚拟相机图像坐标系为坐标系②，以相机A和B的图像坐标系分别作为坐标系①执行以下处理，

1)将坐标系①的原始图像四个角点坐标，基于相应单相机的几何成像模型的坐标正算过程和大虚拟相机的几何成像模型的坐标反算过程，得到输出图像在坐标系②中的范围；

2)对坐标系②输出图像范围内的每个像素，通过大虚拟相机的几何成像模型的坐标正算过程和相应单相机的几何成像模型的坐标反算过程，得到在坐标系①的原始图像坐标；

3)根据坐标系①中的原始图像坐标，通过灰度重采样，得到坐标系②相应输出图像上的每个像素的灰度值。