[发明专利]一种基于定焦数码相机的旋转全景摄影测量方法和装置

权利要求书

1.一种基于定焦数码相机的旋转全景相机的摄影测量方法，其特征在于，包括以下步骤：

步骤1：在全景三维控制场中选择M个摄站，在每一摄站用定焦数码相机拍摄含控制点的影像，共M张影像；

步骤2：计算定焦数码相机的内方位元素(f,x₀,y₀)、畸变参数(K₁,K₂,P₁,P₂)；

提取M张影像上控制点的像点坐标(x,y)，利用已知的控制点像点坐标(x,y)、控制点在全景三维控制场的物方空间坐标(X,Y,Z)和共线方程式，见公式(1)，采用最小二乘的平差方法，解求定焦数码相机的内方位元素(f,x₀,y₀)、外方位元素(X_s,Y_s,Z_s,a₁,a₂,a₃,b₁,b₂,b₃,c₁,c₂,c₃)、畸变参数(K₁,K₂,P₁,P₂)，f,x₀,y₀和K₁,K₂,P₁,P₂的结果作为真值参与后续计算；

公式(1)如下所示：

$\left\{\begin{array}{c}x-x_0-\mathrm{\Δ}x=-f\frac{a_1(X-X_s)+b_1(Y-Y_s)+c_1(Z-Z_s)}{a_3(X-X_s)+b_3(Y-Y_s)+c_3(Z-Z_s)}\\ y-y_0-\mathrm{\Δ}y=-f\frac{a_2(X-X_s)+b_2(Y-Y_s)+c_2(Z-Z_s)}{a_3(X-X_s)+b_3(Y-Y_s)+c_3(Z-Z_s)}\end{array}\right.---\left(1\right)$

式中，(x₀,y₀)为像主点坐标，Δx,Δy为像点畸变差，

$\left\{\begin{array}{c}\mathrm{\Δ}x=(x-x_0)(K_1{r}^2+K_2{r}^4)+P_1(r^2+2{\left(x-x_0\right)}^2)+2P_2(x-x_0)\·(y-y_0)\\ \mathrm{\Δ}y=(y-y_0)(K_1{r}^2+K_2{r}^4)+P_2(r^2+2{\left(y-y_0\right)}^2)+2P_1(x-x_0)\·(y-y_0)\end{array}\right.,$

式中，K₁,K₂为径向畸变参数，P₁,P₂为偏心畸变参数，r为M张影像上控制点的像点到像主点的距离，r＝(x-x₀)²+(y-y₀)²，f为定焦数码相机的焦距，X_s,Y_s,Z_s为投影中心S在全景三维控制场物方空间坐标系(G-X,Y,Z)中的坐标，a₁,a₂,a₃,b₁,b₂,b₃,c₁,c₂,c₃为从像空间坐标系到全景三维控制场物方空间坐标系的旋转矩阵中的元素；

步骤3：在全景三维控制场中选择N个摄站，在每一摄站上用旋转全景相机拍摄360°水平角的全景影像，共N张全景影像；

步骤4：计算旋转全景相机的旋转参数(U_S,V_S,W_S)和R_S；

提取N张全景影像上控制点的像点坐标(x₁,y₁)；根据步骤2中的(f,x₀,y₀)、(K₁,K₂,P₁,P₂)、全景三维控制场中控制点的像点坐标(x₁,y₁)、全景三维控制场中控制点的物方空间坐标(X,Y,Z)和公式(2)，采用最小二乘的平差方法解求旋转中心在全景三维控制场物方空间坐标系中的坐标(X_O,Y_O,Z_O)、旋转坐标系在全景三维控制场物方空间坐标系中的初始角度α₀和旋转全景相机的旋转参数；旋转全景相机的旋转参数包括定焦数码相机投影中心与旋转中心的偏移量(U_S,V_S,W_S)和相对姿态R_S；

公式(2)如下所示：

$\left\{\begin{array}{c}{x}_1-x_0-\mathrm{\Δ}x=-f\frac{\begin{array}{c}(a_1^{\′}{\mathrm{cos\α}}_0-c_1^{\′}{\mathrm{sin\α}}_0)(X-X_O)+b_1^{\′}(Y-Y_O)\\ +(a_1^{\′}{\mathrm{sin\α}}_0+c_1^{\′}{\mathrm{cos\α}}_0)(Z-Z_O)-a_1^{\′}{U}_S-b_1^{\′}{V}_S-c_1^{\′}{W}_S\end{array}}{\begin{array}{c}(a_3^{\′}{\mathrm{cos\α}}_0-c_3^{\′}{\mathrm{sin\α}}_0)(X-X_O)+b_3^{\′}(Y-Y_O)\\ +(a_3^{\′}{\mathrm{sin\α}}_0+c_3^{\′}{\mathrm{cos\α}}_0)(Z-Z_O)-a_3^{\′}{U}_S-b_3^{\′}{V}_S-c_3^{\′}{W}_S\end{array}}\\ y_1-y_0-\mathrm{\Δ}y=-f\frac{\begin{array}{c}(a_2^{\′}{\mathrm{cos\α}}_0-c_2^{\′}{\mathrm{sin\α}}_0)(X-X_O)+b_2^{\′}(Y-Y_O)\\ +(a_2^{\′}{\mathrm{sin\α}}_0+c_2^{\′}{\mathrm{cos\α}}_0)(Z-Z_O)-a_2^{\′}{U}_S-b_2^{\′}V-c_2^{\′}{W}_S\end{array}}{\begin{array}{c}(a_3^{\′}{\mathrm{cos\α}}_0-c_3^{\′}{\mathrm{sin\α}}_0)(X-X_O)+b_3^{\′}(Y-Y_O)\\ +(a_3^{\′}{\mathrm{sin\α}}_0+c_3^{\′}{\mathrm{cos\α}}_0)(Z-Z_O)-a_3^{\′}{U}_S-b_3^{\′}{V}_S-c_3^{\′}{W}_S\end{array}}\end{array}\right.---\left(2\right)$

式中，(x₀,y₀)为像主点坐标，Δx,Δy为像点畸变差，

$\left\{\begin{array}{c}\mathrm{\Δ}x=(x_1-x_0)(K_1{r}_1^2+K_2{r}_1^4)+P_1(r_1^2+2{\left(x_1-x_0\right)}^2)+2P_2(x_1-x_0)\·(y_1-y_0)\\ \mathrm{\Δ}y=(y_1-y_0)(K_1{r}_1^2+K_2{r}_1^4)+P_2(r_1^2+2{\left(y_1-y_0\right)}^2)+2P_1(x_1-x_0)\·(y_1-y_0)\end{array}\right.,$

式中，K₁,K₂为径向畸变参数，P₁,P₂为偏心畸变参数，r₁为N张影像上控制点的像点到像主点的距离，r₁＝(x₁-x₀)²+(y₁-y₀)²，f为定焦数码相机的焦距，X′,Y′,Z′为全景三维控制场内控制点的物方空间坐标，(U_S,V_S,W_S)为投影中心S在旋转坐标系(O-U,V,W)中的坐标，a₁′,a₂′,a₃′,b₁′,b₂′,b₃′,c₁′,c₂′,c₃′为从像空间坐标系到旋转坐标系的旋转矩阵中的元素，(X_O,Y_O,Z_O)为旋转坐标系与全景三维控制场物方空间坐标系的相对位置，α₀为旋转坐标系在全景三维控制场物方空间坐标系中的初始角度；

步骤5：在测量现场布设3个以上的控制点Q₁,Q₂,Q₃……；

步骤6：在测量现场选择T个摄站，在每一摄站上拍摄含Q₁,Q₂,Q₃……和待测点的测量现场的全景影像，共T张全景影像；

步骤7：计算测量现场中旋转全景相机的旋转中心在测量现场物方空间坐标系中的坐标(X′_O,Y′_O,Z′_O)和旋转坐标系在测量现场物方空间坐标系中的初始角度α′₀；

提取测量现场的T张全景影像中控制点Q₁,Q₂,Q₃……的像点坐标(x₂,y₂)；根据步骤2中的(f,x₀,y₀)、(K₁,K₂,P₁,P₂)、步骤4计算得到的全景相机的旋转参数U_S,V_S,W_S、R_S、控制点Q₁,Q₂,Q₃……的像点坐标(x₂,y₂)、控制点在测量现场物方空间坐标系中的坐标(X_Q,Y_Q,Z_Q)和公式(3)，采用最小二乘的平差方法解求旋转全景相机在测量现场的旋转中心坐标(X′_O,Y′_O,Z′_O)、旋转坐标系在测量现场物方空间坐标系中的初始角度α′₀；

公式(3)如下所示：

$\left\{\begin{array}{c}{x}_2-x_0-\mathrm{\Δ}x=-f\frac{\begin{array}{c}(a_1^{\′}{\mathrm{cos\α}}_0^{\′}-c_1^{\′}{\mathrm{sin\α}}_0^{\′})(X_P-X_O^{\′})+b_1^{\′}(Y_P-Y_O^{\′})\\ +(a_1^{\′}{\mathrm{sin\α}}_0^{\′}+c_1^{\′}{\mathrm{cos\α}}_0^{\′})(Z_P-Z_O^{\′})-a_1^{\′}{U}_S-b_1^{\′}{V}_S-c_1^{\′}{W}_S\end{array}}{\begin{array}{c}(a_3^{\′}{\mathrm{cos\α}}_0^{\′}-c_3^{\′}{\mathrm{sin\α}}_0^{\′})(X_P-X_O^{\′})+b_3^{\′}(Y_P-Y_O^{\′})\\ +(a_3^{\′}{\mathrm{sin\α}}_0^{\′}+c_3^{\′}{\mathrm{cos\α}}_0^{\′})(Z_P-Z_O^{\′})-a_3^{\′}{U}_S-b_3^{\′}{V}_S-c_3^{\′}{W}_S\end{array}}\\ y_2-y_0-\mathrm{\Δ}y=-f\frac{\begin{array}{c}(a_2^{\′}{\mathrm{cos\α}}_0^{\′}-c_2^{\′}{\mathrm{sin\α}}_0^{\′})(X_P-X_O^{\′})+b_2^{\′}(Y_P-Y_O^{\′})\\ +(a_2^{\′}{\mathrm{sin\α}}_0^{\′}+c_2^{\′}{\mathrm{cos\α}}_0^{\′})(Z_P-Z_O^{\′})-a_2^{\′}{U}_S-b_2^{\′}V-c_2^{\′}{W}_S\end{array}}{\begin{array}{c}(a_3^{\′}{\mathrm{cos\α}}_0^{\′}-c_3^{\′}{\mathrm{sin\α}}_0^{\′})(X_P-X_O^{\′})+b_3^{\′}(Y_P-Y_O^{\′})\\ +(a_3^{\′}{\mathrm{sin\α}}_0^{\′}+c_3^{\′}{\mathrm{cos\α}}_0^{\′})(Z_P-Z_O^{\′})-a_3^{\′}{U}_S-b_3^{\′}{V}_S-c_3^{\′}{W}_S\end{array}}\end{array}\right.---\left(3\right)$

式中，(x₂,y₂)为像点坐标，(x₀,y₀)为像主点坐标，Δx,Δy为像点畸变差，K₁,K₂为径向畸变参数，P₁,P₂为偏心畸变参数，r₂为T张全景影像中控制点Q₁,Q₂,Q₃,……的像点到像主点的距离r₂＝(x₂-x₀)²+(y₂-y₀)²，f为定焦数码相机的焦距，X_P,Y_P,Z_P为测量现场内控制点的物方空间坐标，(U_S,V_S,W_S)为投影中心S在旋转坐标系(O-U,V,W)中的坐标，a′₁,a′₂,a′₃,b′₁,b′₂,b′₃,c′₁,c′₂,c′₃为从像空间坐标系到旋转坐标系的旋转矩阵中的元素，(X′_O,Y′_O,Z′_O)为旋转坐标系与测量现场物方空间坐标系的相对位置，α′₀为旋转坐标系在测量现场物方空间坐标系的初始角度；

步骤8：由α′₀和旋转角度β＝(360/C)°，根据公式(4)计算每张像片的旋转角度α；

公式(4)如下所示：

α＝α′₀+(i-1)β (i＝1,2,…,C) (4)

其中，C为转动定焦数码相机的次数；

步骤9：计算待测点在测量现场物方空间坐标系中的坐标(X_A,Y_A,Z_A)；

提取测量现场的T张全景影像上待测点的像点坐标(x₃,y₃)；利用已知的f,x₀,y₀、K₁,K₂,P₁,P₂、(U_S,V_S,W_S)、R_S、(x₃,y₃)、(X′_O,Y′_O,Z′_O)、α和公式(5)，采用最小二乘的平差方法解求待测点在测量现场物方空间坐标系中的坐标(X_A,Y_A,Z_A)；

公式(5)如下所示：

$\left\{\begin{array}{c}{x}_3-x_0-\mathrm{\Δ}x=-f\frac{\begin{array}{c}(a_1^{\′}{\mathrm{cos\α}}_0^{\′}-c_1^{\′}{\mathrm{sin\α}}_0^{\′})(X_A-X_O^{\′})+b_1^{\′}(Y_A-Y_O^{\′})\\ +(a_1^{\′}{\mathrm{sin\α}}_0^{\′}+c_1^{\′}{\mathrm{cos\α}}_0^{\′})(Z_A-Z_O^{\′})-a_1^{\′}{U}_S-b_1^{\′}{V}_S-c_1^{\′}{W}_S\end{array}}{\begin{array}{c}(a_3^{\′}{\mathrm{cos\α}}_0^{\′}-c_3^{\′}{\mathrm{sin\α}}_0^{\′})(X_A-X_O^{\′})+b_3^{\′}(Y_P-Y_O^{\′})\\ +(a_3^{\′}{\mathrm{sin\α}}_0^{\′}+c_3^{\′}{\mathrm{cos\α}}_0^{\′})(Z_A-Z_O^{\′})-a_3^{\′}{U}_S-b_3^{\′}{V}_S-c_3^{\′}{W}_S\end{array}}\\ y_3-y_0-\mathrm{\Δ}y=-f\frac{\begin{array}{c}(a_2^{\′}{\mathrm{cos\α}}_0^{\′}-c_2^{\′}{\mathrm{sin\α}}_0^{\′})(X_A-X_O^{\′})+b_2^{\′}(Y_A-Y_O^{\′})\\ +(a_2^{\′}{\mathrm{sin\α}}_0^{\′}+c_2^{\′}{\mathrm{cos\α}}_0^{\′})(Z_A-Z_O^{\′})-a_2^{\′}{U}_S-b_2^{\′}V-c_2^{\′}{W}_S\end{array}}{\begin{array}{c}(a_3^{\′}{\mathrm{cos\α}}_0^{\′}-c_3^{\′}{\mathrm{sin\α}}_0^{\′})(X_A-X_O^{\′})+b_3^{\′}(Y_A-Y_O^{\′})\\ +(a_3^{\′}{\mathrm{sin\α}}_0^{\′}+c_3^{\′}{\mathrm{cos\α}}_0^{\′})(Z_A-Z_O^{\′})-a_3^{\′}{U}_S-b_3^{\′}{V}_S-c_3^{\′}{W}_S\end{array}}\end{array}\right.---\left(5\right)$

式中，(x₃,y₃)为像点坐标，(x₀,y₀)为像主点坐标，Δx,Δy为像点畸变差，K₁,K₂为径向畸变参数，P₁,P₂为偏心畸变参数，r₃为T张全景影像上待测点的像点到像主点的距离，r₃＝(x₃-x₀)²+(y₃-y₀)²，f为定焦数码相机的焦距，X_A,Y_A,Z_A为待测点在测量现场物方空间坐标系中的坐标，(U_S,V_S,W_S)为投影中心S在旋转坐标系(O-U,V,W)中的坐标，a′₁,a′₂,a′₃,b′₁,b′₂,b′₃,c′₁,c′₂,c′₃为从像空间坐标系到旋转坐标系的旋转矩阵中的元素，(X′_O,Y′_O,Z′_O)为旋转坐标系与测量现场物方空间坐标系的相对位置，α′₀为旋转坐标系在测量现场物方空间坐标系的初始角度；

步骤10：结束。