内视场光学分割型大面阵CCD影像几何拼接方法

摘要

本发明涉及一种内视场光学分割型大面阵CCD航摄相机子影像几何拼接合成方法，利用该方法生成的大幅面数字影像可满足航空摄影测量对高摄影效率和高航测成图精度要求。基本思路是：首先高精度量取各子影像重叠区域同名点影像坐标，并获取由相机厂商提供的子影像基准像点在实际焦面上的位置控制坐标；其次，基于拼接关系模型建立观测方程组；然后，采用最小二乘法平差技术求解各子影像到最终合成影像的拼接参数；最后利用拼接参数对子影像进行变换和重采样，从而完成各子影像到最终合成影像的拼接。实验表明本发明具有很好的可靠性和良好的拼接精度，对于内视场光学分割型数字相机，本发明能够无缝拼接生成高精度无几何错位的大幅面航摄影像。

说明书

技术领域

本发明涉及一种内视场光学分割型大面阵CCD航摄相机的子影像几何拼接合 成方法，属于摄影测量与遥感技术领域。

背景技术

到目前为止，随着数字航摄相机的兴起，以前使用的胶片型模拟航摄相机已被完 全取代。但由于受到CCD制造工艺的限制，单张数字航摄影像的像幅规模始终无法达 到原有模拟相机的像面尺寸，这很大程度影响了摄影效率和测图精度等数字航摄相机 的具体应用。

但是在旺盛的应用需求推动下，国内外出现了多种复合型大面阵航摄相机，以解 决数字航摄影像像幅规模受限于单个CCD面阵大小的问题，这些复合型相机可分为外 视场分割型、内视场分割型两大类。其中外视场分割型相机是利用多个单台中小面阵 CCD相机捆绑而成，典型代表有美国Intergraph公司生产的DMC型相机、国家测绘遥感 院制造的SWDC相机等；内视场分割型相机包括Microsofot/Vexcel生产的UltraCam相 机、南京尖兵遥感信息技术有限公司制造的DMZ相机等，这两种内视场分割型相机也 存在构造上的差异，前者是利用多镜头分时摄影模式实现的内视场光学分割，而后者 是直接在单个镜头的后节点视场内利用光学半反半透棱镜一次性进行整体分割。内、 外视场两大类相机各有优缺点，但此处需要注意的是，内视场光学分割型相机的明显 优势是具有良好的拼接精度以及拼接后影像具有良好的量测性能，这对满足后期测图 精度要求至关重要。

本申请正是针对内视场光学分割型这一类相机，提出一种内视场光学分割大面阵 CCD影像拼接方法，以解决高精度子影像拼接、高量测性能的大幅面影像复合生成问 题，国内外与之相同或相似的大面阵影像合成方法未见出现。

发明内容

本发明需要解决的技术问题是提供一种能够用于内视场分割型航摄相机大幅面 影像高精度复合生成方法，利用该方法生成后的大幅面数字影像可满足航空摄影测量 对摄影效率和航测成图精度的应用需求。

本发明提出的内视场分割型航摄相机大幅面影像高精度复合生成方法的基本思 路是：首先高精度量取各子影像重叠区域的同名点影像坐标，并获取由相机厂商提供 的某个/些子影像基准像点在实际焦面上的位置控制坐标；其次，基于一定的拼接关 系模型——如二维等形变换或仿射变换，建立观测方程组；然后，采用最小二乘法平 差技术求解各个子影像到最终合成影像的拼接参数；最后利用拼接参数对子影像进行 变换和重采样，从而完成各个子影像到最终合成影像的拼接。

本发明的技术解决方案是：

内视场光学分割型大面阵CCD影像几何拼接方法，其特征在于该方法包括以下 几个步骤：

（1）量取各相邻个子影像重叠区域内同名像点影像坐标，获取相机厂商提供 的子影像基准点的焦面位置坐标；

（2）基于拼接关系模型建立观测方程组

若确立拼接关系模型为二维等形变换，公式为

X=ax-by+T_X   （1）

Y=ay+bx+T_Y

其中，x，y为子影像坐标系坐标，X,Y为拼接后合成影像坐标系坐标。a,b,T_X,T_Y为4个变换参数，以此对子影像基准像点和相邻子影像同名像点进行拼接关系设计

·第一种拼接关系设计——关于子影像基准像点

对于VI号子影像上的基准像点C为有

$\left(\begin{array}{c}{X}_C=a_{\mathrm{VI}}{x}_{C_{\mathrm{VI}}}-b_{\mathrm{VI}}{y}_{C_{\mathrm{VI}}}+T_{X_{\mathrm{VI}}}\\ Y_C=a_{\mathrm{VI}}{y}_{C_{\mathrm{VI}}}+b_{\mathrm{VI}}{x}_{C_{\mathrm{VI}}}+T_{Y_{\mathrm{VI}}}\end{array}\right)---\left(2\right)$

其中，为子影像VI上基准点C在其子影像坐标系下的像点坐标， X_C,Y_C则为基准点C在焦面上控制点位置的合成后影像坐标，为子影像VI的拼接参数，对于其它I、IV、VII、IX、X、XII子影像上的A、B、D、 E、F、G基准点同样可按上式写出一组对应的拼接关系方程组。

·第二种拼接关系设计——关于相邻子影像同名像点

利用相邻子影像上同名点所对应的合成后影像像点坐标相等这一条件列出拼接 关系方程，对于相邻子影像I和子影像II上的同名点1为有：

$\left(\begin{array}{c}{a}_I{x}_{1_I}-b_I{y}_{1_I}+T_{X_I}-(a_{\mathrm{II}}-x_{1_{\mathrm{II}}}-b_{\mathrm{II}}{y}_{1_{\mathrm{II}}}+T_{X_{\mathrm{II}}})=0\\ a_I{y}_{1_I}+b_I{x}_{1_I}+T_{Y_I}-(a_{\mathrm{II}}{y}_{1_{\mathrm{II}}}+b_{\mathrm{II}}{x}_{1_{\mathrm{II}}}+T_{Y_{\mathrm{II}}})=0\end{array}\right)---\left(3\right)$

对于其它每对相邻子影像，重叠区中出现的每个同名点可按上式写出一组对应的 拼接关系方程组。

这样，将所有基准点及同名像点所列出对应拼接关系观测方程组进行组合，最后 得到的观测方程组矩阵形式为

_mAⁿ_nX¹=_mL¹+_mV¹   （4）

其中，系数矩阵A按照两种拼接关系组建，X为所有子影像的拼接参数，L为由 元素0和基准点焦面位置坐标组成的常量矩阵，V为残差向量，m大小为同名点和基 准点个数的两倍，n大小为子影像数量的四倍；

（3）利用最小二乘法平差技术求解各个子影像到最终合成影像的拼接参数， 并给出解算精度

所有观测方程构建完毕后，每张子影像对应的拼接参数（a,b,T_X,T_Y）可通过最 小二乘平差解求，这些参数确定了所有子影像到合成后影像的几何拼接关系

_mAⁿ_nX¹=_mL¹+_mV¹，X=(A^TA)^-1(A^TL)

对拼接参数X的内部精度按如下公式计算：

V=AX-L（5）

单位权标准差为：

$S_0=\sqrt{\frac{V^{T}V}{r}},r=m-n$

为(A^TA)^-1的第i行、第i列元素，

$S_{x_i}=S_0\sqrt{Q_{x_i{x}_i}}$

即为X中第i个拼接参数的内部精度；

（4）最后利用拼接参数对子影像进行变换和重采样，从而完成各个子影像到 最终合成影像的拼接

实际操作中应将这些拼接参数“反向”应用，公式（1）的矩阵形式如下：

$\left(\begin{array}{c}X\\ Y\end{array}\right)=\left(\begin{array}{cc}a& -b\\ b& a\end{array}\right)\left(\begin{array}{c}x\\ y\end{array}\right)+\left(\begin{array}{c}{T}_X\\ T_Y\end{array}\right)---\left(6\right)$

反向求解，得

$\left(\begin{array}{c}x\\ y\end{array}\right)={\left(\begin{array}{cc}a& -b\\ b& a\end{array}\right)}^{-1}\left(\begin{array}{c}X-T_X\\ Y-T_Y\end{array}\right)=\frac{1}{a^2+b^2}\left(\begin{array}{cc}a& b\\ -b& a\end{array}\right)\left(\begin{array}{c}X-T_X\\ Y-T_Y\end{array}\right)---\left(7\right)$

使用尚未生成的拼接后影像某像素坐标（X，Y）及其对应子影像变换参数，通 过公式（7）可计算出该像素所对应子影像上相应的像素位置，然后可采样得到灰度 信息，对于重叠区域所对应的拼接后影像，可在包含重叠区的任意子影像上采样得 到灰度值或进行灰度值的加权平均即可。

本发明可以求解内视场分割型数字航摄相机各个子影像的拼接参数，从而完成 大幅面影像高精度拼接生成，大幅面合成影像可满足航空摄影测量对高摄影效率和 高影像量测性能的需求。本发明有望增强国内光学内视场分割型数字航摄相机子影 像高精度拼接合成这一技术短板，为该类型航摄相机的后期航测作业提供必要的技 术支撑。

附图说明

图1是1分12内视场分割型DMZ相机I～XII号CCD机背外形安置

图2是1分12内视场分割型DMZ相机内视场分割器及I～XII号CCD分布

图3是图2下对应子影像拓扑关系以及假定的同名像点和基准点分布

图4是系数矩阵A的元素填充效果

图5是拼接后影像上位于不同类型区域的像素点

具体实施方式

以南京尖兵遥感信息技术有限公司的DMZ航摄相机为例，图1是I～XII号共计 12个CCD在相机机背上的安置关系，对应的内视场光学分割器及I～XII号CCD安装 示意如图2，I～XII号CCD所对应的子影像拓扑关系以及假定的同名像点和基准点分 布如图3，子影像标号编码同样为：I～XII，同名点标号编码为：1～33，基准点标号编 码为：A～G，下面具体说明本发明提出的内视场光学分割型大面阵CCD影像几何拼 接方法的实施方式。

（1）量取各个子影像重叠区域同名点影像坐标，获取相机厂商提供的子影像基 准点的焦面位置坐标，定义子影像、同名像点及基准点的编号为：子影像标号编码采 用希腊字母，如I、II、III……，同名点标号编码采用阿拉伯数字，如1、2、3……， 基准点标号编码采用英文大写字母，如A、B、C……

关于同名像点坐标的量测，可利用目前商用全数字摄影测量像点量测工具 Eardas、PCI、ENVI及Ermaper等进行像点的高精度影像匹配量测，依据表1中各个 同名像点及与之对应的相邻子影像关系，量测可得各相邻子影像对应的坐标点串：

$\left(\begin{array}{c}I,\mathrm{II}:(x_{1_I},y_{1_I}),(x_{1_{\mathrm{II}}},y_{1_{\mathrm{II}}});(x_{4_I},y_{4_I}),(x_{4_{\mathrm{II}}},y_{4_{\mathrm{II}}});(x_{9_I},y_{9_I}),(x_{9_{\mathrm{II}}},y_{9_{\mathrm{II}}})\\ ...\\ \mathrm{XI},\mathrm{XII}:(x_{{25}_{\mathrm{XI}}},y_{{25}_{\mathrm{XI}}}),(x_{{25}_{\mathrm{XII}}},y_{{25}_{\mathrm{XII}}});(x_{{30}_{\mathrm{XI}}},y_{{30}_{\mathrm{XI}}}),(x_{{30}_{\mathrm{XII}}},y_{{30}_{\mathrm{XII}}});(x_{{33}_{\mathrm{XI}}},y_{{33}_{\mathrm{XI}}}),(x_{{33}_{\mathrm{XII}}},y_{{33}_{\mathrm{XII}}})\end{array}\right)$

关于基准点焦面位置坐标，其由相机厂商提供，A～G基准点坐标串形式为：

(X_A,Y_A),(X_B,Y_B)…(X_G,Y_G)

表1各个同名像点及与之对应的相邻子影像关系

（2）基于第一种、第二种拼接关系建立观测方程组

以拼接关系模型公式

X=ax-by+T_X   （1）

Y=ay+bx+T_Y

为基础，设计建立两种拼接关系，对于VI号子影像上的基准像点C为有

$\left(\begin{array}{c}{X}_C=a_{\mathrm{VI}}{x}_{C_{\mathrm{VI}}}-b_{\mathrm{VI}}{y}_{C_{\mathrm{VI}}}+T_{X_{\mathrm{VI}}}\\ Y_C=a_{\mathrm{VI}}{y}_{C_{\mathrm{VI}}}+b_{\mathrm{VI}}{x}_{C_{\mathrm{VI}}}+T_{Y_{\mathrm{VI}}}\end{array}\right)---\left(2\right)$

对于相邻子影像I和子影像II上的同名点1为有

$\left(\begin{array}{c}{a}_I{x}_{1_I}-b_I{y}_{1_I}+T_{X_I}-(a_{\mathrm{II}}-x_{1_{\mathrm{II}}}-b_{\mathrm{II}}{y}_{1_{\mathrm{II}}}+T_{X_{\mathrm{II}}})=0\\ a_I{y}_{1_I}+b_I{x}_{1_I}+T_{Y_I}-(a_{\mathrm{II}}{y}_{1_{\mathrm{II}}}+b_{\mathrm{II}}{x}_{1_{\mathrm{II}}}+T_{Y_{\mathrm{II}}})=0\end{array}\right)---\left(3\right)$

并结合表1建立关于所有关于同名点和基准点的观测方程组，系数矩阵₈₀A⁴⁸， 未知数向量₄₈X¹及常数项₈₀L¹的结构如下：

矩阵A中的元素以I、II号相邻子影像重叠区域的1、4同名点，XI、XII号相邻 子影像重叠区域的30、33同名点以及基准点A、G为代表样例进行填充，其它同名 点和基准点可参考进行，而后剩余的大部分元素为0，该稀疏矩阵的元素填充效果见 图4，每个填充块的大小为2行×4列，代表同名点填充块，色代表基准点 填充块。

${}_{48}X^1=\left(\begin{array}{c}{a}_I\\ b_I\\ T_{X_I}\\ T_{Y_I}\\ .\\ .\\ .\\ a_i\\ b_i\\ T_{X_i}\\ T_{Y_i}\\ .\\ .\\ .\\ a_{\mathrm{XII}}\\ b_{\mathrm{XII}}\\ T_{X_{\mathrm{XII}}}\\ T_{Y_{\mathrm{XII}}}\end{array}\right),{}_{80}L^1=\left(\begin{array}{c}0\\ .\\ .\\ .\\ 0\\ X_A\\ Y_A\\ .\\ .\\ .\\ X_G\\ Y_G\end{array}\right)$

至此，观测方程组的矩阵形式_mAⁿ_nX¹=_mL¹+_mV¹组建完毕，其中，m=（33+7） ×2=80，n=12×4=48。

（3）利用最小二乘法平差技术求解各个子影像到最终合成影像的拼接参数，并 给出解算精度

所有观测方程构建完毕后，每张子影像对应的拼接参数（a,b,T_X,T_Y）可通过最 小二乘平差解求，这些参数确定了所有子影像到合成后影像的几何拼接关系

_mAⁿ_nX¹=_mL¹+_mV¹，X=(A^TA)^-1(A^TL)

对拼接参数X的内部精度按如下公式计算：

V=AX-L   （5）

单位权标准差为：

$S_0=\sqrt{\frac{V^{T}V}{r}},r=m-n$

为(A^TA)^-1的第i行、第i列元素，

$S_{x_i}=S_0\sqrt{Q_{x_i{x}_i}}$

即为X中第i个拼接参数的内部精度；

（4）最后利用拼接参数对子影像进行变换和重采样，从而完成各个子影像到最 终合成影像的拼接

实际操作中应将这些拼接参数“反向”应用，公式（1）的矩阵形式如下：

$\left(\begin{array}{c}X\\ Y\end{array}\right)=\left(\begin{array}{cc}a& -b\\ b& a\end{array}\right)\left(\begin{array}{c}x\\ y\end{array}\right)+\left(\begin{array}{c}{T}_X\\ T_Y\end{array}\right)---\left(6\right)$

反向求解，得

$\left(\begin{array}{c}x\\ y\end{array}\right)={\left(\begin{array}{cc}a& -b\\ b& a\end{array}\right)}^{-1}\left(\begin{array}{c}X-T_X\\ Y-T_Y\end{array}\right)=\frac{1}{a^2+b^2}\left(\begin{array}{cc}a& b\\ -b& a\end{array}\right)\left(\begin{array}{c}X-T_X\\ Y-T_Y\end{array}\right)---\left(7\right)$

使用尚未生成的拼接后影像某像素坐标（X，Y）及其对应子影像变换参数，通 过公式（7）可计算出该像素所对应子影像上相应的像素位置，然后可采样得到灰度 信息，当然具体采样方法可采取邻近值或双线性等内插方法，对于重叠区域的拼接 后影像，可在包含重叠区的任意一个子影像上采样得到灰度值或进行多个灰度值的 加权平均即可，图5中最外围矩形即是拼接后影像，其中示意了位于三种不同类型 区域的像素点，像素点L位于“无效区”，像素点J只能在I号子影像上找到对应像 点，像素点K位于I号、2号相邻子影像的重叠区域内。

对于像素点L：使用任何一个子影像拼接参数计算所得像点位置都不会落在相 应子影像的有效像坐标范围内，即是无效的，这时拼接影像像素点L处的灰度值可 用特定代表无效灰度的某一数值填充，比如数值0；

对于像素点J：只能使用I号子影像拼接参数计算而得到位于I号子影像坐标系 下的有效像点，这时拼接影像像素点J处的灰度值可直接用I号子影像坐标 处的采样灰度值进行填充

$\left(\begin{array}{c}{x}_I^J\\ y_I^J\end{array}\right)=\frac{1}{a_I^2+b_I^2}\left(\begin{array}{cc}{a}_I^2& b_I^2\\ -b_I^2& a_I^2\end{array}\right)\left(\begin{array}{c}{X}_J-T_{X_I}\\ Y_J-T_{Y_I}\end{array}\right)$

对于像素点K：可使用I号、2号子影像拼接参数计算而分别得到对应位于I 号、2号子影像坐标系下的两个有效像点，这时拼接影像像素点K处的灰度值可用 I号子影像坐标和II号子影像坐标两处的灰度采样值进行 加权平均处理后填充

$\left(\begin{array}{c}{x}_I^K\\ y_I^K\end{array}\right)=\frac{1}{a_I^2+b_I^2}\left(\begin{array}{cc}{a}_I^2& b_I^2\\ -b_I^2& a_I^2\end{array}\right)\left(\begin{array}{c}{X}_K-T_{X_I}\\ Y_K-T_{Y_I}\end{array}\right),$

$\left(\begin{array}{c}{x}_{\mathrm{II}}^K\\ y_{\mathrm{II}}^K\end{array}\right)=\frac{1}{a_{\mathrm{II}}^2+b_{\mathrm{II}}^2}\left(\begin{array}{cc}{a}_{\mathrm{II}}^2& b_{\mathrm{II}}^2\\ -b_{\mathrm{II}}^2& a_{\mathrm{II}}^2\end{array}\right)\left(\begin{array}{c}{X}_K-T_{X_{\mathrm{II}}}\\ Y_K-T_{Y_{\mathrm{II}}}\end{array}\right)$

至此，对于拼接后影像的每个像素，其必然属于点L、J或K中的某一种像素 点，依据上述对应的方法进行像素灰度值的采样和填充，直至最终拼接影像的生成， 从大量实验结果可以看出，对于内视场光学分割型航摄相机的子影像，运用本发明 提出的内视场光学分割型大面阵CCD影像几何拼接方法，能够无缝拼接生成高精 度无几何错位的大幅面航摄影像，从而满足航测作业对数字航摄相机高摄影效率和 高几何量测性能的实际应用需求。