CN106023146A - 用于摄影测量中的场相关单边自标定光束平差方法 - Google Patents

Abstract

本发明提供了一种用于摄影测量中的场相关单边自标定光束平差方法，包括步骤：a)建立线性的针孔成像模型，并附加场相关的非线性畸变模型；b)确定场相关的单边自标定光束平差理论模型；c)建立场相关单边自标定成像模型的误差方程；d)求取误差方程中场相关成像模型对于外方位参数的偏导数；e)求取误差方程中场相关模型对于空间坐标的偏导数；f)求取误差方程中场相关模型对于单边内参数的偏导数；g)对误差方程雅可比矩阵中的各项进行自适应的比例调整；h)通过分块方式快速计算法方程中的各项；i)求解单边自标定中的各项参数；j)参数比例逆调整，消除比例变化；k)计算中误差，并进行单边自标定各项参数的不确定度估计。

Description

用于摄影测量中的场相关单边自标定光束平差方法

技术领域

本发明涉及一种场相关单边自标定光束平差方法，尤其涉及大尺寸多站位的静态摄影测量系的场相关单边自标定光束平差方法。

背景技术

摄影测量使用测量网络中多站位、多角度拍摄的被测物图片，利用光线束交会定位空间点坐标，大量冗余信息的使用降低了测量的系统误差。相比于测量范围，摄影测量的普遍坐标测量相对误差在[1/100000，1/200000]，长度测量相对误差在[1/50000，1/100000]^[1]。摄影测量除了具有非接触、成本低、高效率的优点外，在精度方面又不亚于其他测量设备。

摄影测量的高精度得益于两个关键因素：回光反射目标点的使用以及像面目标点的高精度定位；具有自标定能力的光束平差技术。自标定光束平差是指在测量过程中完成相机内参数标定的光束平差技术，自标定的意义在于：一方面降低了实验室标定的要求，另一方面使得相机内参数的标定结果与拍摄环境(温度、压力、光照、材质等)相适应，能够在复杂的测量环境中得到最精确的测量结果。

当人们最开始研究成像系统时，已经注意到象差同放大率之间的关系，但是很长时间并未在摄影测量相机模型中应用。直至上世纪70年代，Duane Brown提出了简便的数学描述和标定方法，将相关模型引入摄影测量，Fryer和Fraser也认为物距对于畸变参数有不可忽略的影响。但是这种模型在自标定光束平差中的应用受到限制。

国内对于摄影测量的研究和开发没有深入探讨目标点空间距离对于其成像带来的影响，光束平差模型和算法也不考虑空间坐标、外方位参数与畸变参数之间的相关性。虽然平差结果给出较为满意的空间坐标误差估计，但事实上，相关研究指出，空间误差远大于平差结论。

发明内容

为了解决上述技术问题，本发明提供一种摄影测量相机的场相关单边自标定光束平差方法，包括步骤：a)建立线性的针孔成像模型，并附加场相关的非线性畸变模型，用于描述物空间点、外方位参数、内方位参数之间的数学关系；b)确定场相关单边自标定光束平差理论模型；c)建立场相关单边自标定光束平差理论模型的误差方程；d)求取误差方程中场相关成像模型对于外方位参数的偏导数；e)求取误差方程中场相关成像模型对于空间坐标的偏导数；f)求取误差方程中场相关成像模型对于单边内参数的偏导数；g)对误差方程雅可比矩阵中的各项进行自适应的比例调整，以降低法方程的病态；h)通过分块方式快速计算法方程中的各项；i)通过最小二乘技术迭代求解单边自标定光束平差理论模型中的各项参数；j)参数比例逆调整，消除步骤g)中的人为定义的比例变化；k)计算中误差，并进行单边自标定光束平差理论模型各项参数的不确定度估计。

优选地，所述步骤a)中建立线性成像模型的步骤为：

a1)确定物空间坐标系至相机空间坐标系之间的刚体旋转矩阵：

$\left[\begin{array}{ccc}\cos \left(Ro\right)\cos \left(Az\right)+\sin \left(Ro\right)\sin \left(El\right)\sin \left(Az\right)& \cos \left(Ro\right)\sin \left(Az\right)-\sin \left(Ro\right)\sin \left(El\right)\cos \left(Az\right)& \sin \left(Ro\right)\cos \left(El\right)\\ -\sin \left(Ro\right)\cos \left(Az\right)+\cos \left(Ro\right)\sin \left(El\right)\sin \left(Az\right)& -\sin \left(Ro\right)\sin \left(Az\right)-\cos \left(Ro\right)\sin \left(El\right)\cos \left(Az\right)& \cos \left(Ro\right)\cos \left(El\right)\\ \cos \left(El\right)\sin \left(Az\right)& -\cos \left(El\right)\cos \left(Az\right)& -\sin \left(El\right)\end{array}\right]$

其中，Az,El,Ro为空间坐标系的三个旋转角；

a2)通过旋转矩阵和平移向量描述物空间至相机空间之间的刚体变换：

$\left[\begin{array}{c}\stackrel{\‾}{X}\\ \stackrel{\‾}{Y}\\ \stackrel{\‾}{Z}\end{array}\right]=R(\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]-\left[\begin{array}{c}{X}_0\\ Y_0\\ Z_0\end{array}\right])=\left[\begin{array}{ccc}{a}_1& b_1& c_1\\ a_2& b_2& c_2\\ a_3& b_3& c_3\end{array}\right](\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]-\left[\begin{array}{c}{X}_0\\ Y_0\\ Z_0\end{array}\right])$

a3)在相机坐标系下，经小孔成像，将相机空间中的点投影至像面，得到对应的像点坐标：

$x_l=-f\frac{\stackrel{\‾}{X}}{\stackrel{\‾}{Z}}=-f\frac{a_1(X-X_0)+b_1(Y-Y_0)+c_1(Z-Z_0)}{a_3(X-X_0)+b_3(Y-Y_0)+c_3(Z-Z_0)}$

$y_l=-f\frac{\stackrel{\‾}{Y}}{\stackrel{\‾}{Z}}=-f\frac{a_2(X-X_0)+b_2(Y-Y_0)+c_2(Z-Z_0)}{a_3(X-X_0)+b_3(Y-Y_0)+c_3(Z-Z_0)}$

优选地，所述步骤a)中建立非线性畸变模型的步骤为：

a4)计算位于某一坐标(x，y)处的像点畸变量：

$\mathrm{\Δ}x=\stackrel{\‾}{x}(K_{1{\mathrm{ss}}^{\′}}{r}^2+K_{2{\mathrm{ss}}^{\′}}{r}^4+K_{3{\mathrm{ss}}^{\′}}{r}^6)+P_1(2{\stackrel{\‾}{x}}^2+r^2)+2P_2\stackrel{\‾}{xy}+B_1\stackrel{\‾}{x}+B_2\stackrel{\‾}{y}$

$\mathrm{\Δ}y=\stackrel{\‾}{y}(K_{1{\mathrm{ss}}^{\′}}{r}^2+K_{2{\mathrm{ss}}^{\′}}{r}^4+K_{3{\mathrm{ss}}^{\′}}{r}^6)+P_2(2{\stackrel{\‾}{y}}^2+r^2)+2P_1\stackrel{\‾}{xy}$

$\begin{array}{ccc}\stackrel{\‾}{x}=x-x_p& \stackrel{\‾}{y}=y-y_p& r=\sqrt{{\stackrel{\‾}{x}}^2+{\stackrel{\‾}{y}}^2}\end{array}$

其中，相机的主点偏移为x_p和y_p，对应此空间点的畸变参数为K_1ss’,K_2ss’,K_3ss’,P₁,P₂,B₁,B₂。

a5)获取空间点在该站位相机像面上最终的像点坐标：

x＝x_l+x_p-Δx

优选地，所述步骤b)中确定场相关单边自标定光束平差理论模型的步骤为：

b1)标定两个距离s₁和s₂上径向畸变参数；

b2)推导任意其他物距s′上的畸变参数:

$k_{1{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^2}{{C_{s_1}}^2}{k}_{1{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^2}{{C_{s_2}}^2}{k}_{1{\mathrm{ss}}_2}$

$k_{2{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^4}{{C_{s_1}}^4}{k}_{2{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^4}{{C_{s_2}}^4}{R}_{2{\mathrm{ss}}_2}$

$k_{3{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^6}{{C_{s_1}}^6}{k}_{3{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^6}{{C_{s_2}}^6}{k}_{3{\mathrm{ss}}_2}$

${\mathrm{\α}}_{s^{\′}}=\frac{S_2-S^{\′}}{S_2-S_1}\·\frac{S_1-f}{S^{\′}-f}$

其中，s是成像系统对焦距离，以及是两个距离上径向畸变参数的标定结果，和分别是两个距离对应高斯成像公式的像距；s′为物空间点到相机xoy平面的距离。

b3)将两个距离s₁和s₂中任一距离的径向畸变参数作为已知量，对另一端的径向畸变参数以及其他参数进行平差。

优选地，所述步骤c)中建立场相关单边自标定成像模型的误差方程为：

$\begin{array}{c}{\hat{l}}_{ij}=l_{ij}+v_{ij}=K_{ij}{\hat{x}}_{ij}\\ \⇒\left[\begin{array}{c}{x}_{ij}-\left({X_{ij}}^0\right)_x\\ y_{ij}-\left({X_{ij}}^0\right)_y\end{array}\right]+\left[\begin{array}{c}{v}_{xij}\\ v_{yij}\end{array}\right]=A_{ij}\left[\begin{array}{c}{\mathrm{\ΔX}}_{0i}\\ {\mathrm{\ΔY}}_{0i}\\ {\mathrm{\ΔZ}}_{0i}\\ {\mathrm{\ΔAz}}_i\\ {\mathrm{\ΔEl}}_i\\ {\mathrm{\ΔRo}}_i\end{array}\right]+B_{ij}\left[\begin{array}{c}{\mathrm{\ΔX}}_j\\ {\mathrm{\ΔY}}_j\\ {\mathrm{\ΔZ}}_j\end{array}\right]+C_{ij}\left[\begin{array}{c}{\mathrm{\Δx}}_0\\ {\mathrm{\Δy}}_0\\ \mathrm{\Δ}f\\ {\mathrm{\ΔK}}_{1ss2}\\ {\mathrm{\ΔK}}_{2ss2}\\ {\mathrm{\ΔK}}_{3ss2}\\ {\mathrm{\ΔP}}_1\\ {\mathrm{\ΔP}}_2\\ {\mathrm{\ΔB}}_1\\ {\mathrm{\ΔB}}_2\end{array}\right]=A_{ij}{\mathrm{\δ}}_i+B_{ij}{\stackrel{\·}{\mathrm{\δ}}}_j+C_{ij}\stackrel{\·\·}{\mathrm{\δ}}\end{array}$

其中，用(x_ij,y_ij)表示第i张相片对于空间第j点的成像坐标；X_ij ⁰是所有与成像坐标(x_ij,y_ij)相关的参数，(v_xij,v_yij)是残余误差，A_ij，B_ij和C_ij分别是步骤a5中获取的空间点在该站位相机像面上最终的像点坐标对于第i图片的外方位参数、第j空间点坐标以及相机成像参数的雅可比矩阵。

优选地，所述步骤d)至步骤f)中求取外方位参数、空间坐标以及内参数的偏导数的方法为：利用描述的场相关参数与畸变参数之间的直接关系，通过复合函数求导方法，利用中间量减少偏导数求解的复杂性和计算量。

优选地，所述步骤g中对误差方程雅可比矩阵中的各项进行自适应的比例调整的步骤为：

g1)在每次平差迭代过程开始，首先统计误差方程中对应单边内方位参数雅可比矩阵各列的数量级M_j及调整系数k_j，计算方法如下：

$M_j=average\left(round\right(-\ln \left(A\left(:,j\right)\right)\left)\right)\⇒k_j={10}^{M_j}$

g2)将对应单边内方位参数的雅可比矩阵的每一列计算结果乘以相应列的比例系数，得到调整后的雅可比矩阵为：

$A^{\′}=A\left[\begin{array}{cccc}{k}_1& 0& \mathrm{...}& 0\\ 0& k_2& \mathrm{...}& 0\\ .& .& .& .\\ .& .& .& .\\ .& .& .& .\\ 0& 0& \mathrm{...}& k_u\end{array}\right]$

优选地，所述步骤h)通过分块方式快速计算法方程中的各项的步骤为：

h1)计算A_ij，B_ij和C_ij，分别是步骤c中获取的成像模型对于第i图片的外方位参数、第j空间点坐标以及相机成像单边内参数的雅可比矩阵。

h2)将求解得到的A_ij，B_ij和C_ij，按照如下公式对法方程中的各项进行有规则的分块描述：从而在对图像、目标点的索引过程中逐步得到法方程中的矩阵：

$A^{T}l=\left[\begin{array}{c}\underset{j=1}{\overset{n}{\Σ}}{A}_{1j}^T{l}_{1j}\\ \underset{j=1}{\overset{n}{\Σ}}{A}_{2j}^T{l}_{2j}\\ .\\ .\\ .\\ \underset{j=1}{\overset{n}{\Σ}}{A}_{mj}^T{l}_{mj}\\ \underset{i=1}{\overset{m}{\Σ}}{A}_{i1}^T{l}_{i1}\\ \underset{i=1}{\overset{m}{\Σ}}{B}_{i2}^T{l}_{i2}\\ .\\ .\\ .\\ \underset{i=1}{\overset{m}{\Σ}}{B}_{in}^T{l}_{in}\\ \underset{i=1}{\overset{m}{\Σ}}\underset{j=1}{\overset{n}{\Σ}}{C}_{ij}^T{l}_{ij}\end{array}\right]=\left[\begin{array}{c}\underset{\left(6m\×1\right)}{w_1}\\ \underset{\left(3n\×1\right)}{w_2}\\ \underset{\left(10\×1\right)}{w_3}\end{array}\right]$

h3)获取最终的法方程的分块模型：

$\left[\begin{array}{ccc}\underset{\left(6m\×6m\right)}{N_{11}}& \underset{\left(6m\×3n\right)}{N_{12}}& \underset{\left(6m\×10\right)}{N_{13}}\\ \underset{\left(3n\×6m\right)}{N_{21}}& \underset{\left(3n\×3n\right)}{N_{22}}& \underset{\left(3n\×10\right)}{N_{23}}\\ \underset{\left(10\×6m\right)}{N_{31}}& \underset{\left(10\×3n\right)}{N_{32}}& \underset{\left(10\×10\right)}{N_{33}}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}\underset{\left(6m\×1\right)}{w_1}\\ \underset{\left(3n\×1\right)}{w_2}\\ \underset{\left(10\×1\right)}{w_3}\end{array}\right]\⇒N\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=W$

优选地，所述步骤i)通过最小二乘技术迭代求解单边自标定中的各项参数的步骤为：

i1)利用N₂₂项的分块对角性，将法方程进一步改写为：

$\left[\begin{array}{ccc}{N}_{11}-N_{12}{N}_{22}^{-1}{N}_{21}& 0& N_{13}-N_{12}{N}_{22}^{-1}{N}_{23}\\ N_{21}& N_{22}& N_{23}\\ N_{31}-N_{32}{N}_{22}^{-1}{N}_{21}& 0& N_{33}-N_{32}{N}_{22}^{-1}{N}_{23}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}{w}_1-N_{12}{N}_{22}^{-1}{w}_2\\ w_2\\ w_3-N_{32}{N}_{22}^{-1}{w}_2\end{array}\right]$

i2)利用上式第一、三等式：

$\left[\begin{array}{cc}{M}_{11}& M_{12}\\ M_{21}& M_{22}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}{T}_1\\ T_2\end{array}\right]$

求解外方位参数和单边内方位参数的增量：

$\mathrm{\δ}={\left(M_{11}-M_{12}{M}_{22}^{-1}{M}_{21}\right)}^{-1}\left(T_1-M_{12}{M}_{22}^{-1}{T}_2\right)$

$\stackrel{\·\·}{\mathrm{\δ}}=M_{22}^{-1}\left(T_2-M_{21}\mathrm{\δ}\right)$

i3)将i2)的计算结果带入i1)中求解空间坐标增量其中利用N₂₂的分块对角性进行简化：

$\stackrel{\·}{\mathrm{\δ}}={\left(N_{22}\right)}^{-1}\left(w_2-N_{21}\mathrm{\δ}-N_{23}\stackrel{\·\·}{\mathrm{\δ}}\right)$

优选地，所述步骤j)参数比例逆调整，消除步骤g)中的人为定义的比例变化的方法为：确定未知量δ′与需要求解的δ之间的关系为：

$\mathrm{\δ}=\left[\begin{array}{cccc}{k}_1& 0& \mathrm{...}& 0\\ 0& k_2& \mathrm{...}& 0\\ .& & & .\\ .& & & .\\ .& & & .\\ 0& 0& \mathrm{...}& k_u\end{array}\right]{\mathrm{\δ}}^{\′}$

优选地，所述步骤k)计算中误差，并进行单边自标定各项参数的不确定度估计的步骤为：

k1)计算单位权中误差：

${\hat{s}}_0=\sqrt{\frac{V^{T}PV}{n-u}}$

其中，V是所有观测坐标的残余误差；n是观测坐标数目；u是未知变量数目。

k2)通过协方差矩阵的误差传播定律，计算外方位参数、内参数以及空间坐标的协方差矩阵：

$\underset{\left(u\×u\right)}{C}={\hat{s}}_0^2{\underset{\left(u\×u\right)}{N}}^{-1}$

其中，N为系数矩阵。

总结上述描述，本发明的用于摄影测量中的场相关单边自标定光束平差方法适用于摄影测量的大型稀疏矩阵的分块运算技术已经非常成熟，本发明的相关研究也展示出了场相关的畸变模型对于测量精度的有效提高。如何运用摄影测量传统光束平差中的分块运算技术实现场相关模型的扩展并取得更高精度的测量结果，在国内是一个空白。

本发明以准确的数学描述为基础，从算法上实现场相关的自标定光束平差，解决了以下技术问题：

1、分析场相关的相机成像模型，确定同畸变量相关的相机参数；

2、建立场相关的单边自标定误差方程，对误差方程、法方程中相关项进行准确描述；

3、分块光束平差算法、不确定度估计及自适应的数值比例调节。

附图说明

为了更清楚的说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单的介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他附图。

图1为本发明的用于摄影测量中的场相关单边自标定光束平差方法的流程图；

图2a为本发明的用于摄影测量中的场相关单边自标定光束平差方法的实验中测量场；

图2b为实验中测量场中的定向靶。

具体实施方式

通过参考示范性实施例，本发明的目的和功能以及用于实现这些目的和功能的方法将得以阐明。然而，本发明并不受限于以下所公开的示范性实施例，可以通过不同形式来对其加以实现。说明书的实质仅仅是帮助相关领域技术人员综合理解本发明的具体细节。

如图1所示，为本发明的用于摄影测量中的场相关单边自标定光束平差方法的流程图，步骤为：

101)建立线性的针孔成像模型，并附加场相关的非线性畸变模型，来描述物空间点、外方位参数、内方位参数之间的数学关系；

空间点通过坐标系转换、线性成像以及非线性畸变，投影成为一个像面点，其中，建立线性成像模型的步骤为：

a1)在世界坐标系下，令空间点坐标为[X Y Z]^T，相机位置和指向通过[X₀ Y₀ Z₀AzEl Ro]^T描述，其中前三项表示相机位置坐标，后三项表示相机的方位、俯仰及旋转角，可以通过下式计算旋转矩阵：

$\left[\begin{array}{ccc}\cos \left(Ro\right)\cos \left(Az\right)+\sin \left(Ro\right)\sin \left(El\right)\sin \left(Az\right)& \cos \left(Ro\right)\sin \left(Az\right)-\sin \left(Ro\right)\sin \left(El\right)\cos \left(Az\right)& \sin \left(Ro\right)\cos \left(El\right)\\ -\sin \left(Ro\right)\cos \left(Az\right)+\cos \left(Ro\right)\sin \left(El\right)\sin \left(Az\right)& -\sin \left(Ro\right)\sin \left(Az\right)-\cos \left(Ro\right)\sin \left(El\right)\cos \left(Az\right)& \cos \left(Ro\right)\cos \left(El\right)\\ \cos \left(El\right)\sin \left(Az\right)& -\cos \left(El\right)\cos \left(Az\right)& -\sin \left(El\right)\end{array}\right]$

其中，Az,El,Ro为空间坐标系的三个旋转角；

a2)通过旋转矩阵和平移向量描述物空间至相机空间之间的刚体变换：

$\left[\begin{array}{c}\stackrel{\‾}{X}\\ \stackrel{\‾}{Y}\\ \stackrel{\‾}{Z}\end{array}\right]=R(\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]-\left[\begin{array}{c}{X}_0\\ Y_0\\ Z_0\end{array}\right])=\left[\begin{array}{ccc}{a}_1& b_1& c_1\\ a_2& b_2& c_2\\ a_3& b_3& c_3\end{array}\right](\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]-\left[\begin{array}{c}{X}_0\\ Y_0\\ Z_0\end{array}\right])---\left(1\right)$

a3)在相机坐标系下，经过成像系统投影至像面上，得到线性投影坐标：

$x_l=-f\frac{\stackrel{\‾}{X}}{\stackrel{\‾}{Z}}=-f\frac{a_1(X-X_0)+b_1(Y-Y_0)+c_1(Z-Z_0)}{a_3(X-X_0)+b_3(Y-Y_0)+c_3(Z-Z_0)}$

$y_l=-f\frac{\stackrel{\‾}{Y}}{\stackrel{\‾}{Z}}=-f\frac{a_2(X-X_0)+b_2(Y-Y_0)+c_2(Z-Z_0)}{a_3(X-X_0)+b_3(Y-Y_0)+c_3(Z-Z_0)}---\left(2\right)$

其中，建立非线性畸变模型的步骤为：

a4)计算位于某一坐标(x，y)处的像点畸变量：

$\mathrm{\Δ}x=\stackrel{\‾}{x}(K_1{r}^2+K_2{r}^4+K_3{r}^6)+P_1(2{\stackrel{\‾}{x}}^2+r^2)+2P_2\stackrel{\‾}{xy}+B_1\stackrel{\‾}{x}+B_2\stackrel{\‾}{y}$

$\mathrm{\Δ}y=\stackrel{\‾}{y}(K_1{r}^2+K_2{r}^4+K_3{r}^6)+P_2(2{\stackrel{\‾}{y}}^2+r^2)+2P_1\stackrel{\‾}{xy}$

$\begin{array}{ccc}\stackrel{\‾}{x}=x-x_p& \stackrel{\‾}{y}=y-y_p& r=\sqrt{{\stackrel{\‾}{x}}^2+{\stackrel{\‾}{y}}^2}\end{array}---\left(3\right)$

其中，相机的主点偏移为x_p和y_p，对应此空间点的畸变参数为K_1ss’,K_2ss’,K_3ss’,P₁,P₂,B₁,B₂。

a5)于是空间点在此站位相机像面上最终的像点坐标为：

x＝x_l+x_p-Δx

102)确定场相关单边自标定光束平差理论模型；根据本发明的一个实施例，确定场相关单边自标定光束平差理论模型的步骤为：

b1)在回光反射共面直线阵列标定场内标定两个距离s₁和s₂上的径向畸变参数：

b2)推导任意其他物距s′上的畸变参数:

$k_{1{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^2}{{C_{s_1}}^2}{k}_{1{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^2}{{C_{s_2}}^2}{k}_{1{\mathrm{ss}}_2}$

$k_{2{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^4}{{C_{s_1}}^4}{k}_{2{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^4}{{C_{s_2}}^4}{k}_{2{\mathrm{ss}}_2}$

$k_{3{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^6}{{C_{s_1}}^6}{k}_{3{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^6}{{C_{s_2}}^6}{k}_{3{\mathrm{ss}}_2}---\left(5\right)$

${\mathrm{\α}}_{s^{\′}}=\frac{S_2-S^{\′}}{S_2-S_1}\·\frac{S_1-f}{S^{\′}-f}---\left(6\right)$

其中，s是成像系统对焦距离，和分别是两个距离对应高斯成像公式的像距；s′为物空间点到相机xoy平面的距离。

b3)将s₁和s₂两端的径向畸变参数区别对待。本发明定义，固定一端的径向畸变参数而平差另一端参数的方法叫单边自标定光束平差。将s₁端的径向畸变参数作为已知量，对s₂端的径向畸变参数 以及其他参数进行平差。

103)建立场相关单边自标定光束平差理论模型的误差方程；

用(x_ij,y_ij)表示第i张相片对于空间第j点的成像坐标，可以得到如下面公式描述的误差方程。其中，X_ij ⁰是所有与成像坐标(x_ij,y_ij)相关的参数，(v_xij,v_yij)是残余误差，A_ij，B_ij和C_ij分别是成像模型(4)对于第i图片的外方位参数、第j空间点坐标以及相机成像参数的雅可比矩阵。

$\begin{array}{c}{\hat{l}}_{ij}=l_{ij}+v_{ij}=K_{ij}{\hat{x}}_{ij}\\ \⇒\left[\begin{array}{c}{x}_{ij}-\left({X_{ij}}^0\right)_x\\ y_{ij}-\left({X_{ij}}^0\right)_y\end{array}\right]+\left[\begin{array}{c}{v}_{xij}\\ v_{yij}\end{array}\right]=A_{ij}\left[\begin{array}{c}{\mathrm{\ΔX}}_{0i}\\ {\mathrm{\ΔY}}_{0i}\\ {\mathrm{\ΔZ}}_{0i}\\ {\mathrm{\ΔAz}}_i\\ {\mathrm{\ΔEl}}_i\\ {\mathrm{\ΔRo}}_i\end{array}\right]+B_{ij}\left[\begin{array}{c}{\mathrm{\ΔX}}_j\\ {\mathrm{\ΔY}}_j\\ {\mathrm{\ΔZ}}_j\end{array}\right]+C_{ij}\left[\begin{array}{c}{\mathrm{\Δx}}_0\\ {\mathrm{\Δy}}_0\\ \mathrm{\Δ}f\\ {\mathrm{\ΔK}}_{1ss2}\\ {\mathrm{\ΔK}}_{2ss2}\\ {\mathrm{\ΔK}}_{3ss2}\\ {\mathrm{\ΔP}}_1\\ {\mathrm{\ΔP}}_2\\ {\mathrm{\ΔB}}_1\\ {\mathrm{\ΔB}}_2\end{array}\right]=A_{ij}{\mathrm{\δ}}_i+B_{ij}{\stackrel{\·}{\mathrm{\δ}}}_j+C_{ij}\stackrel{\·\·}{\mathrm{\δ}}\end{array}---\left(7\right)$

将公式(5)描述的畸变参数模型应用至自标定光束平差的关键在于将s₁和s₂两端的径向畸变参数区别对待和处理。本发明中定义，固定两端径向畸变参数的方法叫做非自标定光束平差；固定一端的径向畸变参数而平差另一端参数的方法叫单边自标定光束平差；两端的径向畸变参数都参与平差的方法叫做双边自标定光束平差。本发明中只介绍单边自标定光束平差。

由上可知，径向畸变参数同空间坐标、外方位参数紧密相关，下面的分析均将s₁端的径向畸变参数作为已知量，对s₂端的径向畸变参数 以及其他参数进行平差。

104)求取场相关成像模型的误差方程对于外方位参数的偏导数；

详细推导成像模型对于X₀偏导数，其他项偏导数可以按照类似方法求解。各外方位参数及空间坐标通过s′与畸变参数产生联系，先求解畸变参数对于s′的偏导数：

$\frac{\∂C_{s^{\′}}}{\∂s^{\′}}=-\frac{f^2}{{\left(s^{\′}-f\right)}^2}$

$\frac{\∂{\mathrm{\α}}_{s^{\′}}}{\∂s^{\′}}=-\frac{\left(S_1-f\right)\left(S_2-f\right)}{\left(S_2-S_1\right){\left(S^{\′}-f\right)}^2}$

$\frac{\∂k_{1{\mathrm{ss}}^{\′}}}{\∂s^{\′}}=\frac{\∂{\mathrm{\α}}_{s^{\′}}}{\∂s^{\′}}\frac{{C_{s^{\′}}}^2}{{C_{s_1}}^2}{k}_{1{\mathrm{ss}}_1}+2\frac{\∂C_{s^{\′}}}{\∂s^{\′}}\frac{{\mathrm{\α}}_{s^{\′}}{C}_{s^{\′}}}{{C_{s_1}}^2}{k}_{1{\mathrm{ss}}_1}+(-\frac{\∂{\mathrm{\α}}_{s^{\′}}}{\∂s^{\′}})\frac{{C_{s^{\′}}}^2}{{C_{s_2}}^2}{k}_{1{\mathrm{ss}}_2}+2\frac{\∂C_{s^{\′}}}{\∂s^{\′}}\frac{(1-{\mathrm{\α}}_{s^{\′}})C_{s^{\′}}}{{C_{s_2}}^2}{k}_{1{\mathrm{ss}}_2}$

$\frac{\∂k_{2{\mathrm{ss}}^{\′}}}{\∂s^{\′}}=\frac{\∂{\mathrm{\α}}_{s^{\′}}}{\∂s^{\′}}\frac{{C_{s^{\′}}}^4}{{C_{s_1}}^4}{k}_{2{\mathrm{ss}}_1}+4\frac{\∂C_{s^{\′}}}{\∂s^{\′}}\frac{{\mathrm{\α}}_{s^{\′}}{C_{s^{\′}}}^3}{{C_{s_1}}^4}{k}_{2{\mathrm{ss}}_1}+(-\frac{\∂{\mathrm{\α}}_{s^{\′}}}{\∂s^{\′}})\frac{{C_{s^{\′}}}^4}{{C_{s_2}}^4}{k}_{2{\mathrm{ss}}_2}+4\frac{\∂C_{s^{\′}}}{\∂s^{\′}}\frac{(1-{\mathrm{\α}}_{s^{\′}}){C_{s^{\′}}}^3}{{C_{s_2}}^4}{k}_{2{\mathrm{ss}}_2}$

$\frac{\∂k_{3{\mathrm{ss}}^{\′}}}{\∂s^{\′}}=\frac{\∂{\mathrm{\α}}_{s^{\′}}}{\∂s^{\′}}\frac{{C_{s^{\′}}}^6}{{C_{s_1}}^6}{k}_{3{\mathrm{ss}}_1}+6\frac{\∂C_{s^{\′}}}{\∂s^{\′}}\frac{{\mathrm{\α}}_{s^{\′}}{C_{s^{\′}}}^5}{{C_{s_1}}^6}{k}_{3{\mathrm{ss}}_1}+(-\frac{\∂{\mathrm{\α}}_{s^{\′}}}{\∂s^{\′}})\frac{{C_{s^{\′}}}^6}{{C_{s_2}}^6}{k}_{3{\mathrm{ss}}_2}+6\frac{\∂C_{s^{\′}}}{\∂s^{\′}}\frac{(1-{\mathrm{\α}}_{s^{\′}}){C_{s^{\′}}}^5}{{C_{s_2}}^6}{k}_{3{\mathrm{ss}}_2}$

可知：

$\frac{\∂s^{\′}}{\∂X_0}=a_3$

$\frac{\∂\mathrm{\Δ}x}{\∂s^{\′}}=\stackrel{\‾}{x}(\frac{\∂k_{1{\mathrm{ss}}^{\′}}}{\∂s^{\′}}{r}^2+\frac{\∂k_{2{\mathrm{ss}}^{\′}}}{\∂s^{\′}}{r}^4+\frac{\∂k_{3{\mathrm{ss}}^{\′}}}{\∂s^{\′}}{r}^6),\frac{\∂\mathrm{\Δ}y}{\∂s^{\′}}=\stackrel{\‾}{y}(\frac{\∂k_{1{\mathrm{ss}}^{\′}}}{\∂s^{\′}}{r}^2+\frac{\∂k_{2{\mathrm{ss}}^{\′}}}{\∂s^{\′}}{r}^4+\frac{\∂k_{3{\mathrm{ss}}^{\′}}}{\∂s^{\′}}{r}^6)$

$\frac{\∂\mathrm{\Δ}x}{\∂X_0}=\frac{\∂\mathrm{\Δ}x}{\∂s^{\′}}\frac{\∂s^{\′}}{\∂X_0},\frac{\∂\mathrm{\Δ}y}{\∂X_0}=\frac{\∂\mathrm{\Δ}y}{\∂s^{\′}}\frac{\∂s^{\′}}{\∂X_0}$

通过公式(2)的描述，可以求解x_l和y_l对于X₀的偏导数为：

$\begin{array}{cc}\frac{\∂x_l}{\∂X_0}=\frac{f\left(a_1\stackrel{\‾}{Z}-a_3\stackrel{\‾}{X}\right)}{{\stackrel{\‾}{Z}}^2}& \frac{\∂y_l}{\∂X_0}=\frac{f\left(a_2\stackrel{\‾}{Z}-a_3\stackrel{\‾}{Y}\right)}{{\stackrel{\‾}{Z}}^2}\end{array}$

于是，结合公式(4)，观测值对于外方位参数中X₀项的偏导数为：

$\begin{array}{cc}\frac{\∂x}{\∂X_0}=\frac{\∂x_l}{\∂X_0}-\frac{\∂\mathrm{\Δ}x}{\∂X_0}& \frac{\∂y}{\∂X_0}=\frac{\∂y_l}{\∂X_0}-\frac{\∂\mathrm{\Δ}y}{\∂X_0}\end{array}---\left(8\right)$

同理，可以得到观测值对于其他外方位参数的偏导数：

$\begin{array}{cc}\frac{\∂x}{\∂Y_0}=\frac{\∂x_l}{\∂Y_0}-\frac{\∂\mathrm{\Δ}x}{\∂Y_0}& \frac{\∂y}{\∂Y_0}=\frac{\∂y_l}{\∂Y_0}-\frac{\∂\mathrm{\Δ}y}{\∂Y_0}\end{array}---\left(9\right)$

$\begin{array}{cc}\frac{\∂x}{\∂Z_0}=\frac{\∂x_l}{\∂Z_0}-\frac{\∂\mathrm{\Δ}x}{\∂Z_0}& \frac{\∂y}{\∂Z_0}=\frac{\∂y_l}{\∂Z_0}-\frac{\∂\mathrm{\Δ}y}{\∂Z_0}\end{array}---\left(10\right)$

线性项x_l和y_l对于各角度的偏导数描述比较复杂，可以参照相关文献。这里仅就Δx和Δy对于各角度的偏导数进行分析。

$\frac{\∂s^{\′}}{\∂Az}=-\left(\cos \left(El\right)\cos \left(Az\right)\left(X-X_0\right)+\cos \left(El\right)\sin \left(Az\right)\left(Y-Y_0\right)\right)$

$\begin{array}{cc}\frac{\∂\mathrm{\Δ}x}{\∂Az}=\frac{\∂\mathrm{\Δ}x}{\∂s^{\′}}\frac{\∂s^{\′}}{\∂Az}& \frac{\∂\mathrm{\Δ}y}{\∂Az}=\frac{\∂\mathrm{\Δ}y}{\∂s^{\′}}\frac{\∂s^{\′}}{\∂Az}\end{array}$

$\begin{array}{cc}\frac{\∂x}{\∂Az}=\frac{\∂x_l}{\∂Az}-\frac{\∂\mathrm{\Δ}x}{\∂Az}& \frac{\∂y}{\∂Az}=\frac{\∂y_l}{\∂Az}-\frac{\∂\mathrm{\Δ}y}{\∂Az}\end{array}---\left(11\right)$

$\begin{array}{cc}\frac{\∂x}{\∂El}=\frac{\∂x_l}{\∂El}-\frac{\∂\mathrm{\Δ}x}{\∂El}& \frac{\∂y}{\∂El}=\frac{\∂y_l}{\∂El}-\frac{\∂\mathrm{\Δ}y}{\∂El}\end{array}---\left(12\right)$

$\begin{array}{cc}\frac{\∂x}{\∂Ro}=\frac{\∂x_l}{\∂Ro}& \frac{\∂y}{\∂Ro}=\frac{\∂y_l}{\∂Ro}\end{array}---\left(13\right)$

105)求取误差方程中场相关成像模型对于空间坐标的偏导数；

由于：

$\begin{array}{ccc}\frac{\∂s^{\′}}{\∂X}=-\frac{\∂s^{\′}}{\∂X_0}& \frac{\∂s^{\′}}{\∂Y}=-\frac{\∂s^{\′}}{\∂Y_0}& \frac{\∂s^{\′}}{\∂Z}=-\frac{\∂s^{\′}}{\∂Z_0}\end{array}$

$\begin{array}{ccc}\frac{\∂x_l}{\∂X}=-\frac{\∂x_l}{\∂X_0}& \frac{\∂x_l}{\∂Y}=-\frac{\∂x_l}{\∂Y_0}& \frac{\∂x_l}{\∂Z}=-\frac{\∂x_l}{\∂Z_0}\end{array}$

$\begin{array}{ccc}\frac{\∂y_l}{\∂X}=-\frac{\∂y_l}{\∂X_0}& \frac{\∂y_l}{\∂Y}=-\frac{\∂y_l}{\∂Y_0}& \frac{\∂y_l}{\∂Z}=-\frac{\∂y_l}{\∂Z_0}\end{array}$

因此：

$\begin{array}{ccc}\frac{\∂x}{\∂X}=-\frac{\∂x}{\∂X_0}& \frac{\∂x}{\∂Y}=-\frac{\∂x}{\∂Y_0}& \frac{\∂x}{\∂Z}=-\frac{\∂x}{\∂Z_0}\end{array}---\left(14\right)$

$\begin{array}{ccc}\frac{\∂y}{\∂X}=-\frac{\∂y}{\∂X_0}& \frac{\∂y}{\∂Y}=-\frac{\∂y}{\∂Y_0}& \frac{\∂y}{\∂Z}=-\frac{\∂y}{\∂Z_0}\end{array}---\left(15\right)$

106)求取误差方程中场相关成像模型对于单边内参数的偏导数；

由上可知，参与平差的内参数是x_p,y_p,f,P₁,P₂,B₁,B₂以及S₂单边径向畸变参数

$\frac{\∂\mathrm{\Δ}x}{\∂x_p}=-(K_{1{\mathrm{ss}}^{\′}}{r}^2+K_{2{\mathrm{ss}}^{\′}}{r}^4+K_{3{\mathrm{ss}}^{\′}}{r}^6)+\stackrel{\‾}{x}(-2\stackrel{\‾}{x})(K_{1{\mathrm{ss}}^{\′}}+2K_{2{\mathrm{ss}}^{\′}}{r}^2+3K_{3{\mathrm{ss}}^{\′}}{r}^4)+P_1(-6\stackrel{\‾}{x})+2P_2(-\stackrel{\‾}{y})-B_1$

$\frac{\∂\mathrm{\Δ}y}{\∂x_p}=\stackrel{\‾}{y}(-2\stackrel{\‾}{x})(K_{1{\mathrm{ss}}^{\′}}+2K_{2{\mathrm{ss}}^{\′}}{r}^2+3K_{3{\mathrm{ss}}^{\′}}{r}^4)+P_2(-2\stackrel{\‾}{x})+2P_1(-\stackrel{\‾}{y})$

$\frac{\∂\mathrm{\Δ}x}{\∂y_p}=\stackrel{\‾}{x}(-2\stackrel{\‾}{y})(K_{1{\mathrm{ss}}^{\′}}+2K_{2{\mathrm{ss}}^{\′}}{r}^2+3K_{3{\mathrm{ss}}^{\′}}{r}^4)+P_1(-2\stackrel{\‾}{y})+2P_2(-\stackrel{\‾}{x})-B_2$

$\frac{\∂\mathrm{\Δ}y}{\∂y_p}=-(K_{1{\mathrm{ss}}^{\′}}{r}^2+K_{2{\mathrm{ss}}^{\′}}{r}^4+K_{3{\mathrm{ss}}^{\′}}{r}^6)+\stackrel{\‾}{y}(-2\stackrel{\‾}{y})(K_{1{\mathrm{ss}}^{\′}}+2K_{2{\mathrm{ss}}^{\′}}{r}^2+3K_{3{\mathrm{ss}}^{\′}}{r}^4)+P_2(-6\stackrel{\‾}{y})+2P_1(-\stackrel{\‾}{x})$

$\begin{array}{cc}\frac{\∂x}{\∂x_p}=1-\frac{\∂\mathrm{\Δ}x}{\∂x_p}& \frac{\∂y}{\∂x_p}=-\frac{\∂\mathrm{\Δ}y}{\∂x_p}\end{array}---\left(16\right)$

$\begin{array}{cc}\frac{\∂x}{\∂y_p}=-\frac{\∂\mathrm{\Δ}y}{\∂y_p}& \frac{\∂y}{\∂y_p}=1-\frac{\∂\mathrm{\Δ}y}{\∂y_p}\end{array}---\left(17\right)$

$\begin{array}{cc}\frac{\∂x}{\∂f}=-\frac{\stackrel{\‾}{X}}{\stackrel{\‾}{Z}}& \frac{\∂v}{\∂f}=-\frac{\stackrel{\‾}{Y}}{\stackrel{\‾}{Z}}\end{array}---\left(18\right)$

$\begin{array}{ccc}\frac{\∂K_{1{\mathrm{ss}}^{\′}}}{\∂K_{1{\mathrm{ss}}_2}}=(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^2}{{C_{s_2}}^2}& \frac{\∂K_{2{\mathrm{ss}}^{\′}}}{\∂K_{2{\mathrm{ss}}_2}}=(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^4}{{C_{s_2}}^4}& \frac{\∂K_{3{\mathrm{ss}}^{\′}}}{\∂K_{3{\mathrm{ss}}_2}}=(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^6}{{C_{s_2}}^6}\end{array}$

$\begin{array}{ccc}\frac{\∂x}{\∂K_{1{\mathrm{ss}}_2}}=-\stackrel{\‾}{x}{r}^2\frac{\∂K_{1{\mathrm{ss}}^{\′}}}{\∂K_{1{\mathrm{ss}}_2}}& \frac{\∂x}{\∂K_{2{\mathrm{ss}}_2}}=-\stackrel{\‾}{x}{r}^4\frac{\∂K_{2{\mathrm{ss}}^{\′}}}{\∂K_{2{\mathrm{ss}}_2}}& \frac{\∂x}{\∂K_{3{\mathrm{ss}}_2}}=-\stackrel{\‾}{x}{r}^6\frac{\∂K_{3{\mathrm{ss}}^{\′}}}{\∂K_{3{\mathrm{ss}}_2}}\end{array}---\left(19\right)$

$\begin{array}{ccc}\frac{\∂y}{\∂K_{1{\mathrm{ss}}_2}}=-\stackrel{\‾}{y}{r}^2\frac{\∂K_{1{\mathrm{ss}}^{\′}}}{\∂K_{1{\mathrm{ss}}_2}}& \frac{\∂y}{\∂K_{2{\mathrm{ss}}_2}}=-\stackrel{\‾}{y}{r}^4\frac{\∂K_{2{\mathrm{ss}}^{\′}}}{\∂K_{2{\mathrm{ss}}_2}}& \frac{\∂y}{\∂K_{3{\mathrm{ss}}_2}}=-\stackrel{\‾}{y}{r}^6\frac{\∂K_{3{\mathrm{ss}}^{\′}}}{\∂K_{3{\mathrm{ss}}_2}}\end{array}---\left(20\right)$

$\begin{array}{cc}\frac{\∂x}{\∂P_1}=-(2{\stackrel{\‾}{x}}^2+r^2)& \frac{\∂y}{\∂P_1}=-\end{array}2\stackrel{\‾}{xy}---\left(21\right)$

$\begin{array}{cc}\frac{\∂x}{\∂P_2}=-2\stackrel{\‾}{xy}& \frac{\∂y}{\∂P_2}=-(2{\stackrel{\‾}{y}}^2+r^2)\end{array}---\left(22\right)$

$\begin{array}{cc}\frac{\∂x}{\∂B_1}=-\stackrel{\‾}{x}& \frac{\∂y}{\∂B_1}=0\end{array}---\left(23\right)$

$\begin{array}{cc}\frac{\∂x}{\∂B_2}=-\stackrel{\‾}{y}& \frac{\∂y}{\∂B_2}=0\end{array}---\left(24\right)$

107)对误差方程雅可比矩阵中的各项进行自适应的比例调整，以降低法方程的病态；

待求未知量在数量级上相差巨大，尤其是内参数，例如K₁、P₁、P₂、B₁和B₂项一般具有10^-5，K₂项具有10^-7，而K₃具有10^-11，这种数量级上的巨大差异带来了误差方程中矩阵A的病态，由于数值计算中普遍存在的机器精度、大数吃小数以及大数除小数等现象，会严重影响法方程的计算结果，往往会出现矩阵缺秩等问题。为了统一误差方程中各项的数量级，这里设计了自适应的比例调整技术。

在每次平差迭代过程开始，首先统计误差方程中各列的数量级M_j及调整系数k_j，计算方法如下：

$M_j=average\left(round\left(-\ln \left(A\left(:,j\right)\right)\right)\right)\⇒k_j={10}^{M_j}---\left(37\right)$

然后将每一列的计算结果乘以相应列的比例系数，得到调整后的误差方程矩阵：

$A^{\′}=A\left[\begin{array}{cccc}{k}_1& 0& \mathrm{...}& 0\\ 0& k_2& \mathrm{...}& 0\\ .& & & .\\ .& & & .\\ .& & & .\\ 0& 0& \mathrm{...}& k_u\end{array}\right]$

108)通过分块方式快速计算法方程中的各项；

h1)对于第i图片的外方位参数、第j空间点坐标以及相机成像单边内参数分别求解偏导数：

$A_{ij}=\left[\begin{array}{cccccc}\frac{\∂x_{ij}}{\∂X_{0i}}& \frac{\∂x_{ij}}{\∂Y_{0i}}& \frac{\∂x_{ij}}{\∂Z_{0i}}& \frac{\∂x_{ij}}{\∂{\mathrm{Az}}_i}& \frac{\∂x_{ij}}{\∂{\mathrm{El}}_i}& \frac{\∂x_{ij}}{\∂{\mathrm{Ro}}_i}\\ \frac{\∂y_{ij}}{\∂X_{0i}}& \frac{\∂y_{ij}}{\∂Y_{0i}}& \frac{\∂y_{ij}}{\∂Z_{0i}}& \frac{\∂y_{ij}}{\∂{\mathrm{Az}}_i}& \frac{\∂y_{ij}}{\∂{\mathrm{El}}_i}& \frac{\∂y_{ij}}{\∂{\mathrm{Ro}}_i}\end{array}\right]$

$B_{ij}=\left[\begin{array}{ccc}\frac{\∂x_{ij}}{\∂X_j}& \frac{\∂x_{ij}}{\∂Y_j}& \frac{\∂x_{ij}}{\∂Z_j}\\ \frac{\∂y_{ij}}{\∂X_j}& \frac{\∂y_{ij}}{\∂Y_j}& \frac{\∂y_{ij}}{\∂Z_j}\end{array}\right]$

$C_{ij}=\left[\begin{array}{cccccccccc}\frac{\∂x_{ij}}{\∂x_p}& \frac{\∂x_{ij}}{\∂y_p}& \frac{\∂x_{ij}}{\∂f}& \frac{\∂x_{ij}}{\∂K_{1{\mathrm{ss}}_2}}& \frac{\∂x_{ij}}{\∂K_{2{\mathrm{ss}}_2}}& \frac{\∂x_{ij}}{\∂K_{3{\mathrm{ss}}_2}}& \frac{\∂x_{ij}}{\∂P_1}& \frac{\∂x_{ij}}{\∂P_2}& \frac{\∂x_{ij}}{\∂B_1}& \frac{\∂x_{ij}}{\∂B_2}\\ \frac{\∂y_{ij}}{\∂x_p}& \frac{\∂y_{ij}}{\∂y_p}& \frac{\∂y_{ij}}{\∂f}& \frac{\∂y_{ij}}{\∂K_{1{\mathrm{ss}}_2}}& \frac{\∂y_{ij}}{\∂K_{2{\mathrm{ss}}_2}}& \frac{\∂y_{ij}}{\∂K_{3{\mathrm{ss}}_2}}& \frac{\∂y_{ij}}{\∂P_1}& \frac{\∂y_{ij}}{\∂P_2}& \frac{\∂y_{ij}}{\∂B_1}& \frac{\∂y_{ij}}{\∂B_2}\end{array}\right]$

h2)设一个具有m张图片的测量网络对于n个空间点进行摄影测量，那么误差方程组为：

$\begin{array}{c}\left[\begin{array}{c}{v}_{11}\\ v_{12}\\ .\\ .\\ .\\ v_{1n}\\ v_{21}\\ v_{22}\\ .\\ .\\ .\\ v_{2n}\\ .\\ .\\ .\\ v_{m1}\\ v_{m2}\\ .\\ .\\ .\\ v_{mn}\end{array}\right]+\left[\begin{array}{c}{l}_{11}\\ l_{12}\\ .\\ .\\ .\\ l_{1n}\\ l_{21}\\ l_{22}\\ .\\ .\\ .\\ l_{2n}\\ .\\ .\\ .\\ l_{m1}\\ l_{m2}\\ .\\ .\\ .\\ l_{mn}\end{array}\right]=\left[\begin{array}{ccccccccc}{A}_{11}& 0& \mathrm{...}& 0& B_{11}& 0& \mathrm{...}& 0& C_{11}\\ A_{12}& 0& \mathrm{...}& 0& 0& B_{12}& \mathrm{...}& 0& C_{12}\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ A_{1n}& 0& \mathrm{...}& 0& 0& 0& \mathrm{...}& B_{1n}& C_{1n}\\ 0& A_{21}& \mathrm{...}& 0& B_{21}& 0& \mathrm{...}& 0& C_{21}\\ 0& A_{22}& \mathrm{...}& 0& 0& B_{22}& \mathrm{...}& 0& C_{22}\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ 0& A_{2n}& \mathrm{...}& 0& 0& 0& \mathrm{...}& B_{2n}& C_{2n}\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ 0& 0& \mathrm{...}& A_{m1}& B_{m1}& 0& \mathrm{...}& 0& C_{m1}\\ 0& 0& \mathrm{...}& A_{m2}& 0& B_{m2}& \mathrm{...}& 0& C_{m2}\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ 0& 0& \mathrm{...}& A_{mn}& 0& 0& \mathrm{...}& B_{mn}& C_{mn}\end{array}\right]\\ \⇒v+l=A\mathrm{\δ}\end{array}\left[\begin{array}{c}{\mathrm{\δ}}_1\\ {\mathrm{\δ}}_2\\ .\\ .\\ .\\ {\mathrm{\δ}}_m\\ {\stackrel{\·}{\mathrm{\δ}}}_1\\ {\stackrel{\·}{\mathrm{\δ}}}_2\\ .\\ .\\ .\\ {\stackrel{\·}{\mathrm{\δ}}}_n\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]---\left(25\right)$

对应的法方程为：

(A^TA)δ＝A^Tl (26)

鉴于A的稀疏性，A^TA及A^Tl可以通过对于图片和目标点的索引进行有规则的分块描述：

$\begin{array}{c}{A}^{T}A=\left[\begin{array}{ccccccccc}\underset{j=1}{\overset{n}{\Σ}}{A}_{1j}^T{A}_{1j}& 0& \mathrm{...}& 0& A_{11}^T{B}_{11}& A_{12}^T{B}_{12}& \mathrm{...}& A_{1n}^T{B}_{1n}& \underset{j=1}{\overset{n}{\Σ}}{A}_{1j}^T{C}_{1j}\\ 0& \underset{j=1}{\overset{n}{\Σ}}{A}_{2j}^T{A}_{2j}& \mathrm{...}& 0& A_{21}^T{B}_{21}& A_{22}^T{B}_{22}& \mathrm{...}& A_{2n}^T{B}_{2n}& \underset{j=1}{\overset{n}{\Σ}}{A}_{2j}^T{C}_{2j}\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ 0& 0& \mathrm{...}& \underset{j=1}{\overset{n}{\Σ}}{A}_{mj}^T{A}_{mj}& A_{m1}^T{B}_{m1}^T& A_{m2}^T{B}_{m2}& \mathrm{...}& A_{mn}^T{B}_{mn}& \underset{j=1}{\overset{n}{\Σ}}{A}_{mj}^T{C}_{mj}\\ B_{11}^T{A}_{11}& B_{21}^T{A}_{21}& \mathrm{...}& B_{m1}^T{A}_{m1}& \underset{i=1}{\overset{m}{\Σ}}{B}_{i1}^T{B}_{i1}^T& 0& \mathrm{...}& 0& \underset{i=1}{\overset{m}{\Σ}}{B}_{i1}^T{C}_{i1}\\ B_{12}^T{A}_{12}& B_{22}^T{A}_{22}& \mathrm{...}& B_{m2}^T{A}_{m2}& 0& \underset{i=1}{\overset{m}{\Σ}}{B}_{i2}^T{B}_{i2}^T& \mathrm{...}& 0& \underset{i=1}{\overset{m}{\Σ}}{B}_{i2}^T{C}_{i2}\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ .& .& .& .& .& .& .& .& .\\ B_{1n}^T{A}_{1n}& B_{2n}^T{A}_{2n}& \mathrm{...}& B_{mn}^T{A}_{mn}& 0& 0& \mathrm{...}& \underset{i=1}{\overset{m}{\Σ}}{B}_{in}^T{B}_{in}^T& \underset{i=1}{\overset{m}{\Σ}}{B}_{in}^T{C}_{in}\\ \underset{j=1}{\overset{n}{\Σ}}{C}_{1j}^T{A}_{1j}& \underset{j=1}{\overset{n}{\Σ}}{C}_{2j}^T{A}_{2j}& & \underset{j=1}{\overset{n}{\Σ}}{C}_{mj}^T{A}_{mj}& \underset{i=1}{\overset{m}{\Σ}}{C}_{i1}^T{B}_{i1}& \underset{i=1}{\overset{m}{\Σ}}{C}_{i2}^T{B}_{i2}& \mathrm{...}& \underset{i=1}{\overset{m}{\Σ}}{C}_{in}^T{B}_{in}& \underset{i=1}{\overset{m}{\Σ}}\underset{j=1}{\overset{n}{\Σ}}{C}_{ij}^T{C}_{ij}\end{array}\right]\\ =\left[\begin{array}{ccc}\underset{\left(6m\×6m\right)}{N_{11}}& \underset{\left(6m\×3n\right)}{N_{12}}& \underset{\left(6m\×10\right)}{N_{13}}\\ \underset{\left(3n\×6m\right)}{N_{21}}& \underset{\left(3n\×3n\right)}{N_{22}}& \underset{\left(3n\×10\right)}{N_{23}}\\ \underset{\left(10\×6m\right)}{N_{31}}& \underset{\left(10\×3n\right)}{N_{32}}& \underset{\left(10\×10\right)}{N_{33}}\end{array}\right]\end{array}---\left(27\right)$

$A^{T}l=\left[\begin{array}{c}\underset{j=1}{\overset{n}{\Σ}}{A}_{1j}^T{l}_{1j}\\ \underset{j=1}{\overset{n}{\Σ}}{A}_{2j}^T{l}_{2j}\\ .\\ .\\ .\\ \underset{j=1}{\overset{n}{\Σ}}{A}_{mj}^T{l}_{mj}\\ \underset{i=1}{\overset{m}{\Σ}}{A}_{i1}^T{l}_{i1}\\ \underset{i=1}{\overset{m}{\Σ}}{B}_{i2}^T{l}_{i2}\\ .\\ .\\ .\\ \underset{i=1}{\overset{m}{\Σ}}{B}_{in}^T{l}_{in}\\ \underset{i=1}{\overset{m}{\Σ}}\underset{j=1}{\overset{n}{\Σ}}{C}_{ij}^T{l}_{ij}\end{array}\right]=\left[\begin{array}{c}\underset{\left(6m\×1\right)}{w_1}\\ \underset{\left(3n\×1\right)}{w_2}\\ \underset{\left(10\×1\right)}{w_3}\end{array}\right]---\left(28\right)$

h3)于是，对于(26)的求解演变成对于下式的求解：

$\left[\begin{array}{ccc}\underset{\left(6m\×6m\right)}{N_{11}}& \underset{\left(6m\×3n\right)}{N_{12}}& \underset{\left(6m\×10\right)}{N_{13}}\\ \underset{\left(3n\×6m\right)}{N_{21}}& \underset{\left(3n\×3n\right)}{N_{22}}& \underset{\left(3n\×10\right)}{N_{23}}\\ \underset{\left(10\×6m\right)}{N_{31}}& \underset{\left(10\×3n\right)}{N_{32}}& \underset{\left(10\×10\right)}{N_{33}}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}\underset{\left(6m\×1\right)}{w_1}\\ \underset{\left(3n\×1\right)}{w_2}\\ \underset{\left(10\×1\right)}{w_3}\end{array}\right]\⇒N\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=W---\left(29\right)$

109)通过最小二乘技术迭代求解单边自标定光束平差理论模型中的各项参数；

i1)充分利用N₂₂项的分块对角性，将法方程进一步描述为：

$\left[\begin{array}{ccc}{N}_{11}-N_{12}{N}_{22}^{-1}{N}_{21}& 0& N_{13}-N_{12}{N}_{22}^{-1}{N}_{23}\\ N_{21}& N_{22}& N_{23}\\ N_{31}-N_{32}{N}_{22}^{-1}{N}_{21}& 0& N_{33}-N_{32}{N}_{22}^{-1}{N}_{23}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}{w}_1-N_{12}{N}_{22}^{-1}{w}_2\\ w_2\\ w_3-N_{32}{N}_{22}^{-1}{w}_2\end{array}\right]---\left(30\right)$

i2)利用上式第一、三等式：

$\left[\begin{array}{cc}{M}_{11}& M_{12}\\ M_{21}& M_{22}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}{T}_1\\ T_2\end{array}\right]---\left(31\right)$

求解外方位参数和内参数的增量：

$\mathrm{\δ}={\left(M_{11}-M_{12}{M}_{22}^{-1}{M}_{21}\right)}^{-1}\left(T_1-M_{12}{M}_{22}^{-1}{T}_2\right)---\left(32\right)$

$\stackrel{\·\·}{\mathrm{\δ}}=M_{22}^{-1}\left(T_2-M_{21}\mathrm{\δ}\right)---\left(33\right)$

i3)将i2)的计算结果带入1)中求解空间坐标增量其中利用N₂₂的分块对角性进行简化：

$\stackrel{\·}{\mathrm{\δ}}={\left(N_{22}\right)}^{-1}\left(w_2-N_{21}\mathrm{\δ}-N_{23}\stackrel{\·\·}{\mathrm{\δ}}\right)---\left(34\right)$

110)参数比例逆调整，消除步骤g)中的人为定义的比例变化；

未知量δ′与需要求解的δ之间的关系为：

$\mathrm{\δ}=\left[\begin{array}{cccc}{k}_1& 0& \mathrm{...}& 0\\ 0& k_2& \mathrm{...}& 0\\ .& & & .\\ .& & & .\\ .& & & .\\ 0& 0& \mathrm{...}& k_u\end{array}\right]{\mathrm{\δ}}^{\′}---\left(38\right)$

111)计算中误差，并进行单边自标定光束平差理论模型各项参数的不确定度估计。

公式(32)、(33)和(34)计算的外方位参数、单边内参数以及空间坐标的不确定度可以通过其误差的协方差矩阵进行描述，误差协方差矩阵通过误差传播定律给出。单位权中误差为：

${\hat{s}}_0=\sqrt{\frac{V^{T}PV}{n-u}}---\left(35\right)$

其中，V是所有观测坐标的残余误差；n是观测坐标数目；u是未知变量数目。

各未知量的协方差矩阵为：

$\underset{\left(u\×u\right)}{C}={\hat{s}}_0^2{\underset{\left(u\×u\right)}{N}}^{-1}---\left(36\right)$

其中，N是h3)描述的系数矩阵。

以下为实验验证本发明的用于摄影测量中的场相关单边自标定光束平差方法的技术效果。

实验为用一台AVT工业相机，具有1600万像素的全画幅传感器，搭配35mm尼康定焦镜头及商用闪光灯。数据通过千兆网络传输，并直接进入编写的单边场相关自标定光束平差软件算法进行处理，图2a为本发明的用于摄影测量中的场相关单边自标定光束平差方法的实验中测量场，如图所示，由编码点201a和普通点201b、定向靶203以及背景202组成测量场；其中如图2b所示，定向靶203的两端布置有编码点或/和普通点。

测量之前通过特殊标定场和特殊标定方法标定相机在4.489米处的各项畸变参数(主要是径向畸变参数)，在单边场相关自标定光束平差的过程中计算自标定另一距离上的各项畸变参数，结果如下表所示。

表1场相关单边自标定光束平差相机内参数

	距离1	距离2
			S(mm)	4489.4	1909.5
xp(mm)	-0.1267	-0.2057
			yp(mm)	0.0278	0.1973
f(mm)		36.0997
			K1	6.6597E-05	6.5463E-05
K2	-5.2602E-08	-3.0186E-08
			K3	-1.1628E-11	-4.6766E-11
P1	8.1816E-06	6.9022E-06
			P2	-1.2852E-05	-3.4025E-06
B1	-7.9364E-03	1.1967E-04
			B2	-1.7230E-04	-2.8629E-05

可以看出，自标定结果体现出了不同距离上径向畸变参数的差别。空间点坐标不确定度估计的统计结果如表2所示。

表2场相关模型光束平差统计结果(单位：mm)

为了实现场相关及场无关自标定模型测量精度的比较，设计了相关实验。在测量场中不同位置放置等长度基准尺，用两种模型分别进行摄影测量，最后统计长度测量的不确定度，统计结果如下表所示。

表3不同光束平差模型长度测量不确定度的统计结果

	场相关模型	场无关模型
			长度不确定度(mm)	0.012	0.020
长度测量相对误差	1/79000	1/47500

本发明详细介绍了场相关的径向畸变模型在摄影测量中的应用方法。推导了场相关的单边自标定光束平差中的误差方程；解决了如何利用误差方程的稀疏性分块求解法方程及未知数改正数，在平差结束给出未知参数的不确定度估计；通过自适应的比例调节消除平差运算中的病态问题；最后通过实验验证了场相关单边自标定技术的可行性，并验证了这种模型对于摄影测量长度测量精度的提高。

以上只是本发明较佳的实例，并非来限制本发明实施范围，故凡依本发明申请专利范围所述的构造、特征及原理所做的等效变化或修饰，均应包括于本发明专利申请范围内。

Claims (10)

1.一种用于摄影测量中的场相关单边自标定光束平差方法，所述步骤包括：

a)建立线性的针孔成像模型，并附加场相关的非线性畸变模型，用于描述物空间点、外方位参数、内方位参数之间的数学关系；

b)确定场相关单边自标定光束平差理论模型；

c)建立场相关单边自标定光束平差理论模型的误差方程；

d)求取误差方程中场相关成像模型对于外方位参数的偏导数；

e)求取误差方程中场相关成像模型对于空间坐标的偏导数；

f)求取误差方程中场相关成像模型对于单边内参数的偏导数；

g)对误差方程雅可比矩阵中的各项进行自适应的比例调整，以降低法方程的病态；

h)通过分块方式快速计算法方程中的各项；

i)通过最小二乘方法迭代求解单边自标定光束平差理论模型中的各项参数；

j)参数比例逆调整，消除步骤g)中的人为定义的比例变化；

k)计算中误差，并进行单边自标定光束平差理论模型各项参数的不确定度估计。

2.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤a)中建立线性成像模型的步骤为：

a1)确定物空间坐标系至相机空间坐标系之间的刚体旋转矩阵：

$\left[\begin{array}{ccc}\cos \left(Ro\right)\cos \left(Az\right)+\sin \left(Ro\right)\sin \left(El\right)\sin \left(Az\right)& \cos \left(Ro\right)\sin \left(Az\right)-\sin \left(Ro\right)\sin \left(El\right)\cos \left(Az\right)& \sin \left(Ro\right)\cos \left(El\right)\\ -\sin \left(Ro\right)\cos \left(Az\right)+\cos \left(Ro\right)\sin \left(El\right)\sin \left(Az\right)& -\sin \left(Ro\right)\sin \left(Az\right)-\cos \left(Ro\right)\sin \left(El\right)\cos \left(Az\right)& \cos \left(Ro\right)\cos \left(El\right)\\ \cos \left(El\right)\sin \left(Az\right)& -\cos \left(El\right)\cos \left(Az\right)& -\sin \left(El\right)\end{array}\right]$

其中，Az,El,Ro为空间坐标系的三个旋转角；

a2)通过旋转矩阵和平移向量描述物空间至相机空间之间的刚体变换：

$\left[\begin{array}{c}\stackrel{\‾}{X}\\ \stackrel{\‾}{Y}\\ \stackrel{\‾}{Z}\end{array}\right]=R(\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]-\left[\begin{array}{c}{X}_0\\ Y_0\\ Z_0\end{array}\right])=\left[\begin{array}{ccc}{a}_1& b_1& c_1\\ a_2& b_2& c_2\\ a_3& b_3& c_3\end{array}\right](\left[\begin{array}{c}X\\ Y\\ Z\end{array}\right]-\left[\begin{array}{c}{X}_0\\ Y_0\\ Z_0\end{array}\right])$

a3)在相机坐标系下，经小孔成像，将相机空间中的点投影至像面，得到对应的像点坐标：

$x_l=-f\frac{\stackrel{\‾}{X}}{\stackrel{\‾}{Z}}=-f\frac{a_1(X-X_0)+b_1(Y-Y_0)+c_1(Z-Z_0)}{a_3(X-X_0)+b_3(Y-Y_0)+c_3(Z-Z_0)}$

$y_l=-f\frac{\stackrel{\‾}{Y}}{\stackrel{\‾}{Z}}=-f\frac{a_2\left(X-X_0\right)+b_2\left(Y-Y_0\right)+c_2\left(Z-Z_0\right)}{a_3\left(X-X_0\right)+b_3\left(Y-Y_0\right)+c_3\left(Z-Z_0\right)}$

3.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤a)中建立非线性畸变模型的步骤为：

a4)计算位于某一坐标(x，y)处的像点畸变量：

$\mathrm{\Δ}x=\stackrel{\‾}{x}(K_{1{\mathrm{ss}}^{\′}}{r}^2+K_{2{\mathrm{ss}}^{\′}}{r}^4+K_{3{\mathrm{ss}}^{\′}}{r}^6)+P_1(2{\stackrel{\‾}{x}}^2+r^2)+2P_2\stackrel{\‾}{xy}+B_1\stackrel{\‾}{x}+B_2\stackrel{\‾}{y}$

$\mathrm{\Δ}y=\stackrel{\‾}{y}(K_{1{\mathrm{ss}}^{\′}}{r}^2+K_{2{\mathrm{ss}}^{\′}}{r}^4+K_{3{\mathrm{ss}}^{\′}}{r}^6)+P_2(2{\stackrel{\‾}{y}}^2+r^2)+2P_1\stackrel{\‾}{xy}$

$\begin{array}{ccc}\stackrel{\‾}{x}=x-x_p& \stackrel{\‾}{y}=y-y_p& r=\sqrt{{\stackrel{\‾}{x}}^2+{\stackrel{\‾}{y}}^2}\end{array}$

其中，相机的主点偏移为x_p和y_p，对应此空间点的畸变参数为K_1ss’,K_2ss’,K_3ss’,P₁,P₂,B₁,B₂。

a5)获取空间点在该站位相机像面上最终的像点坐标：

x＝x_l+x_p-Δx

y＝y_l+y_p-Δy

4.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤b)中确定场相关单边自标定光束平差理论模型的步骤为：

b1)标定两个距离s₁和s₂上径向畸变参数；

b2)推导任意其他物距s′上的畸变参数:

$k_{1{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^2}{{C_{s_1}}^2}{k}_{1{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^2}{{C_{s_2}}^2}{k}_{1{\mathrm{ss}}_2}$

$k_{2{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^4}{{C_{s_1}}^4}{k}_{2{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^4}{{C_{s_2}}^4}{k}_{2{\mathrm{ss}}_2}$

$k_{3{\mathrm{ss}}^{\′}}={\mathrm{\α}}_{s^{\′}}\frac{{C_{s^{\′}}}^6}{{C_{s_1}}^6}{k}_{3{\mathrm{ss}}_1}+(1-{\mathrm{\α}}_{s^{\′}})\frac{{C_{s^{\′}}}^6}{{C_{s_2}}^6}{k}_{3{\mathrm{ss}}_2}$

${\mathrm{\α}}_{s^{\′}}=\frac{S_2-S^{\′}}{S_2-S_1}\·\frac{S_1-f}{S^{\′}-f}$

其中，s是成像系统对焦距离，以及是两个距离上径向畸变参数的标定结果；

和分别是两个距离对应高斯成像公式的像距；

s′为物空间点到相机xoy平面的距离；

α_s′为中间计算量。

b3)将两个距离s₁和s₂中任一距离的径向畸变参数作为已知量，对另一端的径向畸变参数以及其他参数进行平差。

5.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤c)中建立场相关单边自标定成像模型的误差方程为：

$\begin{array}{c}{\hat{l}}_{ij}=l_{ij}+v_{ij}=K_{ij}{\hat{x}}_{ij}\\ \⇒\left[\begin{array}{c}{x}_{ij}-\left({X_{ij}}^0\right)_x\\ y_{ij}-\left({X_{ij}}^0\right)_y\end{array}\right]+\left[\begin{array}{c}{v}_{xij}\\ v_{yij}\end{array}\right]=A_{ij}\left[\begin{array}{c}{\mathrm{\ΔX}}_{0i}\\ {\mathrm{\ΔY}}_{0i}\\ {\mathrm{\ΔZ}}_{0i}\\ {\mathrm{\ΔAz}}_i\\ {\mathrm{\ΔEl}}_i\\ {\mathrm{\ΔRo}}_i\end{array}\right]+B_{ij}\left[\begin{array}{c}{\mathrm{\ΔX}}_j\\ {\mathrm{\ΔY}}_j\\ {\mathrm{\ΔZ}}_j\end{array}\right]+C_{ij}\left[\begin{array}{c}{\mathrm{\Δx}}_0\\ {\mathrm{\Δy}}_0\\ \mathrm{\Δ}f\\ {\mathrm{\ΔK}}_{1ss2}\\ {\mathrm{\ΔK}}_{2ss2}\\ {\mathrm{\ΔK}}_{3ss2}\\ {\mathrm{\ΔP}}_1\\ {\mathrm{\ΔP}}_2\\ {\mathrm{\ΔB}}_1\\ {\mathrm{\ΔB}}_2\end{array}\right]=A_{ij}{\mathrm{\δ}}_i+B_{ij}{\stackrel{\·}{\mathrm{\δ}}}_j+C_{ij}\stackrel{\·\·}{\mathrm{\δ}}\end{array}$

其中，用(x_ij,y_ij)表示第i张相片对于空间第j点的成像坐标；X_ij ⁰是所有与成像坐标(x_ij,y_ij)相关的参数，(v_xij,v_yij)是残余误差，A_ij，B_ij和C_ij分别是步骤a5中获取的空间点在该站位相机像面上最终的像点坐标对于第i图片的外方位参数、第j空间点坐标以及相机成像参数的雅可比矩阵。

6.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤g中对误差方程雅可比矩阵中的各项进行自适应的比例调整的步骤为：

g1)在每次平差迭代过程开始，首先统计误差方程中对应单边内方位参数雅可比矩阵各列的数量级M_j及调整系数k_j，计算方法如下：

g2)将对应单边内方位参数的雅可比矩阵的每一列计算结果乘以相应列的比例系数，得到调整后的雅可比矩阵为：

$A^{\′}=A\left[\begin{array}{cccc}{k}_1& 0& ...& 0\\ 0& k_2& ...& 0\\ \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}\\ 0& 0& ...& k_u\end{array}\right]$

7.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤h)通过分块方式快速计算法方程中的各项的步骤为：

h1)计算A_ij，B_ij和C_ij，分别是步骤c中获取的成像模型对于第i图片的外方位参数、第j空间点坐标以及相机成像单边内参数的雅可比矩阵。

h2)将求解得到的A_ij，B_ij和C_ij，按照如下公式对法方程中的各项进行有规则的分块描述：从而在对图像、目标点的索引过程中逐步得到法方程中的矩阵：

$\begin{array}{c}{A}^{T}A=\left[\begin{array}{ccccccccc}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{1j}^T{A}_{1j}& 0& ...& 0& A_{11}^T{B}_{11}& A_{12}^T{B}_{12}& ...& A_{1n}^T{B}_{1n}& \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{1j}^T{C}_{1j}\\ 0& \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{2j}^T{A}_{2j}& ...& 0& A_{21}^T{B}_{21}& A_{22}^T{B}_{22}& ...& A_{2n}^T{B}_{2n}& \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{2j}^T{C}_{2j}\\ \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}\\ 0& 0& ...& \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{mj}^T{A}_{mj}& A_{m1}^T{B}_{m1}& A_{m2}^T{B}_{m2}& ...& A_{mn}^T{B}_{mn}& \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{mj}^T{C}_{mj}\\ B_{11}^T{A}_{11}& B_{21}^T{A}_{21}& ...& B_{m1}^T{A}_{m1}& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{i1}^T{B}_{i1}^T& 0& ...& 0& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{i1}^T{C}_{i1}\\ B_{12}^T{A}_{12}& B_{22}^T{A}_{22}& ...& B_{m2}^T{A}_{m2}& 0& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{i2}^T{B}_{i2}^T& ...& 0& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{i2}^T{C}_{i2}\\ \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}\\ B_{1n}^T{A}_{1n}& B_{2n}^T{A}_{2n}& ...& B_{mn}^T{A}_{mn}& 0& 0& ...& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{in}^T{B}_{in}^T& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{in}^T{C}_{in}\\ \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{C}_{1j}^T{A}_{1j}& \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{C}_{2j}^T{A}_{2j}& & \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{C}_{mj}^T{A}_{mj}& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{C}_{i1}^T{B}_{i1}& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{C}_{i2}^T{B}_{i2}& ...& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{C}_{in}^T{B}_{in}& \underset{i=1}{\overset{m}{\mathrm{\Σ}}}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{C}_{ij}^T{C}_{ij}\end{array}\right]\\ =\left[\begin{array}{ccc}\underset{\left(6m\×6m\right)}{N_{11}}& \underset{\left(6m\×3n\right)}{N_{12}}& \underset{\left(6m\×10\right)}{N_{13}}\\ \underset{\left(3n\×6m\right)}{N_{21}}& \underset{\left(3n\×3n\right)}{N_{22}}& \underset{\left(3n\×10\right)}{N_{23}}\\ \underset{\left(10\×6m\right)}{N_{31}}& \underset{\left(10\×3n\right)}{N_{32}}& \underset{\left(10\×10\right)}{N_{33}}\end{array}\right]\end{array}$

$A^{T}l=\left[\begin{array}{c}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{1j}^T{l}_{1j}\\ \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{2j}^T{l}_{2j}\\ \begin{array}{c}.\\ .\\ .\end{array}\\ \underset{j=1}{\overset{n}{\mathrm{\Σ}}}{A}_{mj}^T{l}_{mj}\\ \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{i1}^T{l}_{i1}\\ \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{i2}^T{l}_{i2}\\ \begin{array}{c}.\\ .\\ .\end{array}\\ \underset{i=1}{\overset{m}{\mathrm{\Σ}}}{B}_{in}^T{l}_{in}\\ \underset{i=1}{\overset{m}{\mathrm{\Σ}}}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{C}_{ij}^T{l}_{ij}\end{array}\right]=\left[\begin{array}{c}\underset{\left(6m\×1\right)}{w_1}\\ \underset{\left(3n\×1\right)}{w_2}\\ \underset{\left(10\×1\right)}{w_3}\end{array}\right]$

h3)获取最终的法方程的分块模型。

$\left[\begin{array}{ccc}\underset{\left(6m\×6m\right)}{N_{11}}& \underset{\left(6m\×3n\right)}{N_{12}}& \underset{\left(6m\×10\right)}{N_{13}}\\ \underset{\left(3n\×6m\right)}{N_{21}}& \underset{\left(3n\×3n\right)}{N_{22}}& \underset{\left(3n\×10\right)}{N_{23}}\\ \underset{\left(10\×6m\right)}{N_{31}}& \underset{\left(10\×3n\right)}{N_{32}}& \underset{\left(10\×10\right)}{N_{33}}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}\underset{\left(6m\×1\right)}{w_1}\\ \underset{\left(3n\×1\right)}{w_2}\\ \underset{\left(10\×1\right)}{w_3}\end{array}\right]\⇒N\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=W$

8.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤i)通过最小二乘技术迭代求解单边自标定中的各项参数的步骤为：

i1)利用N₂₂项的分块对角性，改写法方程；

$\left[\begin{array}{ccc}{N}_{11}-N_{12}{N}_{22}^{-1}{N}_{21}& 0& N_{13}-N_{12}{N}_{22}^{-1}{N}_{23}\\ N_{21}& N_{22}& N_{23}\\ N_{31}-N_{32}{N}_{22}^{-1}{N}_{21}& 0& N_{33}-N_{32}{N}_{22}^{-1}{N}_{23}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·}{\mathrm{\δ}}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}{w}_1-N_{12}{N}_{22}^{-1}{w}_2\\ w_2\\ w_3-N_{32}{N}_{22}^{-1}{w}_2\end{array}\right]$

表示为：

$\left[\begin{array}{cc}{M}_{11}& M_{12}\\ M_{21}& M_{22}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}\\ \stackrel{\·\·}{\mathrm{\δ}}\end{array}\right]=\left[\begin{array}{c}{T}_1\\ T_2\end{array}\right]$

i2)求解外方位参数和单边内方位参数的增量；

$\mathrm{\δ}={\left(M_{11}-M_{12}{M}_{22}^{-1}{M}_{21}\right)}^{-1}\left(T_1-M_{12}{M}_{22}^{-1}{T}_2\right)$

$\stackrel{\·\·}{\mathrm{\δ}}=M_{22}^{-1}\left(T_2-M_{21}\mathrm{\δ}\right)$

i3)将i2)的计算结果带入i1)中求解空间坐标增量其中利用N₂₂的分块对角性进行简化；

$\stackrel{\·}{\mathrm{\δ}}={\left(N_{22}\right)}^{-1}\left(w_2-N_{21}\mathrm{\δ}-N_{23}\stackrel{\·\·}{\mathrm{\δ}}\right)$

9.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤j)参数比例逆调整，消除步骤g)中的人为定义的比例变化的方法为：确定未知量δ′与需要求解的δ之间的关系为：

$\mathrm{\δ}=\left[\begin{array}{cccc}{k}_1& 0& ...& 0\\ 0& k_2& ...& 0\\ \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}& \begin{array}{c}.\\ .\\ .\end{array}\\ 0& 0& ...& k_u\end{array}\right]{\mathrm{\δ}}^{\′}$

10.根据权利要求1所述的场相关单边自标定光束平差方法，其特征是：所述步骤k)计算中误差，并进行单边自标定各项参数的不确定度估计的步骤为：

k1)计算单位权中误差：

${\hat{s}}_0=\sqrt{\frac{V^{T}PV}{n-u}}$

其中，V是所有观测坐标的残余误差；n是观测坐标数目；u是未知变量数目；

k2)通过协方差矩阵的误差传播定律，计算外方位参数、内参数以及空间坐标的协方差矩阵：

$\underset{\left(u\×u\right)}{C}={\hat{s}}_0^2{\underset{\left(u\×u\right)}{N}}^{-1}$

其中，N为系数矩阵。