实习报告

学生:邬代杰

院系:遥感信息工程学院

学号:2019302130329

班级:19D9

文章目录

实习报告
- 一. 基本原理
- - 1. 寻找控制点坐标
  - 2. 对控制点进行编号
  - 3. 计算相机外部参数
  - 4. 计算控制点三维坐标
- 二. 实现过程
- - 1. 寻找控制点
  - 2. 寻找圆盘中心点
  - 3. 对控制点分组
  - 4. 对控制点编号
  - 5. 优化外部参数
  - 6. 前方交会反算控制点世界坐标以检验外方位元素计算准确性
- 三. 操作说明
- - 1. 运行文件位置说明
  - 2. 文件路径替换
  - 3. MATLAB工作区变量说明
- 四. 结论分析
- 五. 实习成果

一. 基本原理

1. 寻找控制点坐标

控制点就是圆盘上黑色区域的中心坐标，体现在图像上就是黑色像素集合区域的中心，那我们首先要做的就是提取出黑色像素。由于原图像为RGB彩色图像，佛像与背景、控制点黑色区域与圆盘其他部分对比度较高，故考虑将图像二值化后再进行操作，具体操作应用MATLAB自带函数就可以解决，不再赘述。

2. 对控制点进行编号

提取出控制点坐标后，接下来的任务就是利用已知的控制点编号规则在图像上对控制点进行编号。

注意到控制点是在一个圆盘之上的，即是说：
任意过同一编号组控制点的直线一定都相交于圆盘圆心且过圆心的直线至多与4个控制点中心相交

又有：投影（相机摄像）不会改变物体的相对位置，并且图片已经消除了畸变，不会出现同一组控制点不在一条直线上的情况。

基于以上事实，我们可以得到将控制点编号的方法：

选取两组控制点中各两点，列直线方程；
计算两直线相交点作为圆盘中心点；
逐点进行遍历并分组：首先任意选取一点作为起始点，然后应用直线拟合对其他点进行判定（即是否为同一组控制点）并记录同一组的控制点的像素坐标，然后不断选取未分组的点进行分组操作，直到全部点分组完毕；
将分组后的控制点从外到内进行排序：在分好组的控制点间计算其与中心点之间的棋盘距离，采用冒泡排序以便于同时更新控制点坐标；
编号：根据控制点区域大小赋给控制点不同的值（0或1），4个点得到一个二进制数（不满4个点的控制点组舍弃），依此序列计算组号。以0100为例，组号为：0×23+1×22+0×21+0×20=40\times 2^3+1\times2^2+0\times2^1+0\times2^0=40×23+1×22+0×21+0×20=4. 若组号为nnn，从外到内的编号依次为：
N=4(n−1)+i−1(i是从外到内点的位置)N=4(n-1)+i-1\quad(i是从外到内点的位置) N=4(n−1)+i−1(i是从外到内点的位置)

至此，编号结束。

3. 计算相机外部参数

采取后方交会计算左右相片各自的外部参数。

后方交会原理:

共线方程
x−x0=−fa1(X−XS)+b1(Y−YS)+c1(Z−ZS)a3(X−XS)+b3(Y−YS)+c3(Z−ZS)=−fX‾Z‾y−y0=−fa2(X−XS)+b2(Y−YS)+c2(Z−ZS)a3(X−XS)+b3(Y−YS)+c3(Z−ZS)=−fY‾Z‾\begin{aligned} x-x_0&=-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)}=-f\frac{\overline{X}}{\overline{Z}}\\ y-y_0&=-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)}=-f\frac{\overline{Y}}{\overline{Z}} \end{aligned} x−x0y−y0=−fa3(X−XS)+b3(Y−YS)+c3(Z−ZS)a1(X−XS)+b1(Y−YS)+c1(Z−ZS)=−fZX=−fa3(X−XS)+b3(Y−YS)+c3(Z−ZS)a2(X−XS)+b2(Y−YS)+c2(Z−ZS)=−fZY
其中
R=[a1a2a3b1b2b3c1c2c3]=[cos⁡φcos⁡κ−sin⁡φsin⁡ωsin⁡κ−cos⁡φsin⁡κ−sin⁡φcos⁡ωcos⁡ωsin⁡κcos⁡ωcos⁡κ−sin⁡ωsin⁡ωcos⁡κ+cos⁡φsin⁡ωsin⁡κ−sin⁡φsin⁡κ+cos⁡φsin⁡ωcos⁡κcos⁡φcos⁡ω]\begin{aligned} R&=\begin{bmatrix} a_1 & a_2& a_3\\ b_1 & b_2& b_3\\ c_1 & c_2& c_3 \end{bmatrix}\\ &=\begin{bmatrix} \cos\varphi\cos\kappa-\sin\varphi\sin\omega\sin\kappa & -\cos\varphi\sin\kappa & -\sin\varphi\cos\omega\\ \cos\omega\sin\kappa & \cos\omega\cos\kappa & -\sin\omega\\ \sin\omega\cos\kappa+\cos\varphi\sin\omega\sin\kappa & -\sin\varphi\sin\kappa+\cos\varphi\sin\omega\cos\kappa & \cos\varphi\cos\omega \end{bmatrix} \end{aligned} R=⎣⎡a1b1c1a2b2c2a3b3c3⎦⎤=⎣⎡cosφcosκ−sinφsinωsinκcosωsinκsinωcosκ+cosφsinωsinκ−cosφsinκcosωcosκ−sinφsinκ+cosφsinωcosκ−sinφcosω−sinωcosφcosω⎦⎤
(X,Y,Z)(X,Y,Z)(X,Y,Z)为像点世界坐标, (XS,YS,ZS)(X_S,Y_S,Z_S)(XS,YS,ZS)为相片外部参数, (x,y)(x,y)(x,y)为像点坐标, (x0,y0)(x_0,y_0)(x0,y0)为像点中心坐标。

由间接平差模型
Av=Bx−lAv=Bx-l Av=Bx−l
对共线方程线性化后有
xT=[ΔXSΔYSΔZSΔφΔωΔκ]B=[a11a12a13a14a15a16a21a22a23a24a25a26]l=[lxly]x^T=\begin{bmatrix}\Delta X_S&\Delta Y_S&\Delta Z_S&\Delta\varphi&\Delta\omega&\Delta\kappa\end{bmatrix}\\ B=\begin{bmatrix} a_{11}&a_{12}&a_{13}&a_{14}&a_{15}&a_{16}\\ a_{21}&a_{22}&a_{23}&a_{24}&a_{25}&a_{26} \end{bmatrix}\\ l=\begin{bmatrix} l_x\\l_y \end{bmatrix} xT=[ΔXSΔYSΔZSΔφΔωΔκ]B=[a11a21a12a22a13a23a14a24a15a25a16a26]l=[lxly]
其中
{a11=1Z‾[a1f+a3(x−x0)]a12=1Z‾[b1f+b3(x−x0)]a13=1Z‾[c1f+c3(x−x0)]a21=1Z‾[a2f+a3(y−y0)]a22=1Z‾[b2f+b3(y−y0)]a23=1Z‾[c2f+c3(y−y0)]\left\{ \begin{aligned} a_{11}&=\frac{1}{\overline{Z}}[a_1f+a_3(x-x_0)]\\ a_{12}&=\frac{1}{\overline{Z}}[b_1f+b_3(x-x_0)]\\ a_{13}&=\frac{1}{\overline{Z}}[c_1f+c_3(x-x_0)]\\ a_{21}&=\frac{1}{\overline{Z}}[a_2f+a_3(y-y_0)]\\ a_{22}&=\frac{1}{\overline{Z}}[b_2f+b_3(y-y_0)]\\ a_{23}&=\frac{1}{\overline{Z}}[c_2f+c_3(y-y_0)] \end{aligned} \right. ⎩⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎧a11a12a13a21a22a23=Z1[a1f+a3(x−x0)]=Z1[b1f+b3(x−x0)]=Z1[c1f+c3(x−x0)]=Z1[a2f+a3(y−y0)]=Z1[b2f+b3(y−y0)]=Z1[c2f+c3(y−y0)]

{a14=(y−y0)sin⁡ω−{(x−x0)f[(x−x0)cos⁡κ−(y−y0)sin⁡κ]+fcos⁡κ}cos⁡ωa15=−fsin⁡κ−x−x0f{(x−x0)sin⁡κ+(y−y0)cos⁡κ}a16=+(y−y0)a24=−(x−x0)sin⁡ω−{(y−y0)f[(x−x0)cos⁡κ−(y−y0)sin⁡κ]−fcos⁡κ}cos⁡ωa25=−fsin⁡κ−y−y0f{(x−x0)sin⁡κ+(y−y0)cos⁡κ}a26=−(x−x0)\left\{ \begin{aligned} a_{14}&=(y-y_0)\sin\omega-\{\frac{(x-x_0)}{f}[(x-x_0)\cos\kappa-(y-y_0)\sin\kappa]+f\cos\kappa\}\cos\omega\\ a_{15}&=-f\sin\kappa-\frac{x-x_0}{f}\{(x-x_0)\sin\kappa+(y-y_0)\cos\kappa\}\\ a_{16}&=+(y-y_0)\\ a_{24}&=-(x-x_0)\sin\omega-\{\frac{(y-y_0)}{f}[(x-x_0)\cos\kappa-(y-y_0)\sin\kappa]-f\cos\kappa\}\cos\omega\\ a_{25}&=-f\sin\kappa-\frac{y-y_0}{f}\{(x-x_0)\sin\kappa+(y-y_0)\cos\kappa\}\\ a_{26}&=-(x-x_0)\\ \end{aligned} \right. ⎩⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎧a14a15a16a24a25a26=(y−y0)sinω−{f(x−x0)[(x−x0)cosκ−(y−y0)sinκ]+fcosκ}cosω=−fsinκ−fx−x0{(x−x0)sinκ+(y−y0)cosκ}=+(y−y0)=−(x−x0)sinω−{f(y−y0)[(x−x0)cosκ−(y−y0)sinκ]−fcosκ}cosω=−fsinκ−fy−y0{(x−x0)sinκ+(y−y0)cosκ}=−(x−x0)
计算xxx的公式为
x=(BTPB)−1BPlx=(B^TPB)^{-1}BPl x=(BTPB)−1BPl
其中PPP矩阵为权阵，在本次实习中，各控制点间相互独立且精度相同，故权阵为单位矩阵，变为求
x=(BTB)−1Blx=(B^TB)^{-1}Bl x=(BTB)−1Bl
在MATLAB中表示为
x=(B′∗B)B∗lx=(B'*B)\text{\\}B*l x=(B′∗B)B∗l
在计算出改正值后，不断迭代更新外方位元素平差值直至改正值在限差之内即可。

4. 计算控制点三维坐标

原理：前方交会

已知：
相机内参：x0,y0,f,width,height相片外参：XS1,YS1,ZS1,φ1,ω1,κ1XS2,YS2,ZS2,φ2,ω2,κ2同名像点坐标：(x1,y1),(x2,y2)\begin{aligned} 相机内参：&x_0,y_0,f,width,height\\ 相片外参：&X_{S1},Y_{S1},Z_{S1},\varphi_1,\omega_1,\kappa_1\\ &X_{S2},Y_{S2},Z_{S2},\varphi_2,\omega_2,\kappa_2\\ 同名像点坐标：&(x_1,y_1),(x_2,y_2) \end{aligned} 相机内参：相片外参：同名像点坐标：x0,y0,f,width,heightXS1,YS1,ZS1,φ1,ω1,κ1XS2,YS2,ZS2,φ2,ω2,κ2(x1,y1),(x2,y2)
待求：
控制点三维坐标：Xp,Yp,Zp控制点三维坐标：X_p,Y_p,Z_p 控制点三维坐标：Xp,Yp,Zp
过程：

计算左右片在地辅坐标系中旋转矩阵的方向余弦：该步骤与后方交会时计算旋转矩阵一致；
R=[a1a2a3b1b2b3c1c2c3]R′=[a1′a2′a3′b1′b2′b3′c1′c2′c3′]\begin{aligned} R&=\begin{bmatrix} a_1 \quad a_2\quad a_3\\ b_1 \quad b_2\quad b_3\\ c_1 \quad c_2\quad c_3 \end{bmatrix}\\ R'&=\begin{bmatrix} a_1' \quad a_2'\quad a_3'\\ b_1' \quad b_2'\quad b_3'\\ c_1' \quad c_2'\quad c_3' \end{bmatrix} \end{aligned} RR′=⎣⎡a1a2a3b1b2b3c1c2c3⎦⎤=⎣⎡a1′a2′a3′b1′b2′b3′c1′c2′c3′⎦⎤
计算基线分量
BX=Xs2−Xs1BY=Ys2−Ys1BZ=Zs2−Zs1B _ { X } = X _ { s2 } - X _ { s1 } \\ B _ { Y } = Y _ { s 2 } - Y _ { s1 } \\ B _ { Z } = Z _ { s 2 } - Z _ { s 1 } BX=Xs2−Xs1BY=Ys2−Ys1BZ=Zs2−Zs1
计算像点的像空间辅助坐标
[XYZ]=R[x1y1−f],[X′Y′Z′]=R[x2y2−f]\begin{bmatrix} X\\Y\\Z \end{bmatrix} =R\begin{bmatrix} x_1\\y_1\\-f \end{bmatrix}, \begin{bmatrix} X'\\Y'\\Z' \end{bmatrix} =R\begin{bmatrix} x_2\\y_2\\-f \end{bmatrix} ⎣⎡XYZ⎦⎤=R⎣⎡x1y1−f⎦⎤,⎣⎡X′Y′Z′⎦⎤=R⎣⎡x2y2−f⎦⎤
计算投影系数
N=BXZ′−BZX′XZ′−ZX′N′=BXZ−BZXXZ′−ZX′N=\frac{B_X Z'-B_Z X'}{XZ'-ZX'}\\ N'=\frac{B_X Z-B_Z X}{XZ'-ZX'} N=XZ′−ZX′BXZ′−BZX′N′=XZ′−ZX′BXZ−BZX
计算地面点的左像辅系坐标
ΔX=NXΔY=(NY+N′Y′+BY)/2ΔZ=NZ\begin{aligned} \Delta X&=NX\\ \Delta Y&=(NY+N'Y'+B_Y)/2\\ \Delta Z&=NZ \end{aligned} ΔXΔYΔZ=NX=(NY+N′Y′+BY)/2=NZ
计算地面点的地面坐标
Xp=Xs1+ΔXYp=Ys1+ΔYZp=Zs1+ΔZ\begin{aligned} X_p&=X_{s1}+\Delta X\\ Y_p&=Y_{s1}+\Delta Y\\ Z_p&=Z_{s1}+\Delta Z\\ \end{aligned} XpYpZp=Xs1+ΔX=Ys1+ΔY=Zs1+ΔZ

二. 实现过程

以下均以左边图像为例

1. 寻找控制点

main.m第52~74行:

RGB转灰度图–>灰度图转二值图–>对二值图取质心及面积–>储存于结构体ldftInfo.coodInfo中

2. 寻找圆盘中心点

main.m第76~85行

选定第1,13,28,40点计算中心点坐标–>计算交点(中心点坐标)–>存于ldftInfo.centreInfo

3. 对控制点分组

main.m第87~119行

flag数组用于标记是否分组–>将选择出来的未分组的控制点与中心点进行一次幂拟合,得到系数–>遍历所有控制点,寻找在一条直线上的点–>对有的不足4个元素的数组进行补充–>重复迭代直至遍历所有控制点

4. 对控制点编号

main.m第122~158行

距离中心点的距离从小到大进行排序–>求二进制编码组–>编号–>存于leftInfo.coodInfo第4列(Inf为图中有但未进行编号的控制点,其原因是一组中控制点由于佛像遮挡有丢失,无法编号)

编号校验无误.

5. 优化外部参数

main.m第355~392行
dealA.m计算B矩阵，dealR.m计算旋转矩阵

6. 前方交会反算控制点世界坐标以检验外方位元素计算准确性

main.m第437~469行

寻找左右图同名点后根据前方交会原理一步一步计算即可，不再赘述。

三. 操作说明

1. 运行文件位置说明

请将main.m、dealR.m以及dealA.m放在同一文件夹下，并将该文件夹更改为工作文件夹，然后再运行main.m文件。

2. 文件路径替换

main.m文件中路径需要根据自己的情况更改

第49行：左边图像路径

第167行：右边图像路径

第295行：控制点坐标文件

第303行：内部参数文件

第317行：外部参数初值文件

3. MATLAB工作区变量说明

innerInfo内部参数
leftBI rightBI二值化图像
leftRGB rightRGB读取的图像
leftInfo rightInfo左右图片点的信息，包括.coodInfo （像点坐标）以及centreInfo（中心点坐标）
lOutInfo rOutInfo外部参数
point控制点信息，包括.data(数据)和.textdata(文件头)
其余多余变量均已清除。

着重介绍一下point.data结构，它各列代表的信息如下表：

列数	1	2	3	4	5
含义	控制点编号	文件读取的X	文件读取的Y	文件读取的Z	控制点类别
6	7	8	9	10	11
前方交会的X’	前方交会的Y’	前方交会的Z’	X-X’	Y-Y’	Z-Z’

其中第6~11列中的控制点，若没有同名点与之对应，则全部赋值为0.

第61行第9~11列数据为累积误差。

四. 结论分析

考虑到求解控制点像点坐标时开运算、闭运算带来的位置偏差和矩阵运算时的精度误差，这个结果是在可以接受的范围内的。在计算过程中，要注意到各个参数的意义，套用公式时要明白其中的原理，特别是各坐标系之间的转换关系，不然很可能会出现某些结果异常的情况，找bug时是十分痛苦的（T x T）

最后附上代码结构

五. 实习成果

实习MATLAB代码：点击下载
博客园链接:点击跳转

WHU2021计算机视觉实习项目相关推荐

机器学习与计算机视觉入门项目——视频投篮检测（二）
机器学习与计算机视觉入门项目--视频投篮检测(二) 一.手工特征与CNN特征在上一次的博客中,介绍了计算机视觉和机器学习的关系.篮球进球检测的基本问题和数据集的制作.这次的我们主要介绍如何从原始图像 ...
计算机视觉一些项目实战技术（续）
计算机视觉一些项目实战技术(续) PROTO-OBJECT BASED SALIENCY 在本项目中,提出一种新的方法来完成显著目标侦测的任务.与以往基于聚光灯注意理论的显著目标检测器相比,遵循基于对 ...
计算机视觉一些项目实战技术
计算机视觉一些项目实战技术 SELECTIVE SEARCH FOR OBJECT LOCALISATION 需要多种策略来查找上述图像中的所有对象.勺子在桌子上的沙拉碗里.因此,图像本质上是层次性的 ...
直方图均衡视觉显著_计算机视觉一些项目实战技术（续）
计算机视觉一些项目实战技术(续) 1. PROTO-OBJECT BASED SALIENCY 在本项目中,提出一种新的方法来完成显著目标侦测的任务.与以往基于聚光灯注意理论的显著目标检测器相比,遵循 ...
软件实习项目2——贪吃喵（猫吃鱼版贪吃蛇）（代码实现）
软件实习项目2--贪吃喵(猫吃鱼版贪吃蛇)(代码实现) 类变量的定义以及类的初始化__init__ 一.游戏的逻辑 1.猫头的生成 2.鱼的生成 3.猫头和鱼骨的移动 4.按下键盘,改变方向二.主窗 ...
软件实习项目2——贪吃喵（猫吃鱼版贪吃蛇）（实验准备与设计）
软件实习项目2--贪吃喵(猫吃鱼版贪吃蛇)(实验准备与设计) 实验内容编程语言以及开发环境的选择实验思路(游戏设计) 一.游戏的逻辑设计 1.猫头的生成 2.鱼的生成 3.猫头和鱼骨的移动 4.按 ...
软件实习项目2——贪吃喵（猫吃鱼版贪吃蛇）（成品展示）
软件实习项目2--贪吃喵(猫吃鱼版贪吃蛇)(成品展示) 成品展示 1.开始游戏界面 2.游戏主界面 3.结束游戏界面视频演示成品展示 1.开始游戏界面速度选择: 猫咪类型选择: 2.游戏主界面 ...
介绍ReLU6！计算机视觉实习面经(京东/商汤/思谋/依图/图森/字节/腾讯)
点击下方卡片,关注"CVer"公众号 AI/CV重磅干货,第一时间送达本文作者:大白鲨 | 来源:知乎(已授权) https://zhuanlan.zhihu.com/p/34 ...
实习证明| 大数据在线实习项目意义
实习,在实践中学习,任何的知识源于实践,归于实践,最终也要通过实践进行检验,它是衔接学生在校学习与步入社会两个阶段的重要桥梁.一个好的实习机会,会帮助学生积累工作经验.验证职业选择.提升综合素质.但是 ...

WHU2021计算机视觉实习项目