matlab实现TSAI标定及其误差分析
1.结果截图
2.结果讨论
用FASI标定的结果得到的结果图如图三所示,用HALCON得到的标定结果如图四所示。由于是用手机进行拍照,得到的相机内部参数不全,于是本次实验结果仅讨论焦距。手机拍照时的真实焦距为25mm,比较两次实验结果的相对误差可得,FASI与真实焦距的相对误差为7.228%,HALCON与真实焦距的相对误差为11.384%。就焦距这一项数据而言,FASI的准确率更高。
3.食用方法
将提取到的特征点写入ImageCoordinate.txt文件,将特征点所对应的真实世界坐标写入WordCoordinate.txt,运行testData.m即可。
4.代码
专业课老师给的代码
testData.m
clc
clear all
x=load('WordCoordinate.txt');
X=load('ImageCoordinate.txt');
Xf=X(:,1); %图像横坐标,表示像素
Yf=X(:,2);%图像纵坐标,表示像素
xw=x(:,1);
yw=x(:,2);
[M,N]=size(x);
zw=zeros(M,1);
Ncx=1;%可理解为平面像素的采样频率,一般设置为1
Nfx=1;
Cx=640; %Cx,Cy是像平面中心像素位置,即U0,V0
Cy=512;
dx=0.0052;%每个像素的物理尺寸
dy=0.0052;
sx=1;%相当于误差项,也设置为1
[R, T, f, k1] = Tsai(Xf, Yf, xw, yw, zw, Ncx, Nfx, dx, dy, Cx, Cy, sx);
a=atan(-R(2,3)/R(3,3));
a1=a*180/pi;
c=atan(-R(1,2)/R(1,1));
c1=c*180/pi;
q=cos(c);
b=atan(R(1,3)*q/R(1,1));
b1=b*180/pi;
Tsai.m
% [R, T, f, k1] = Tsai (Xf, Yf, xw, yw, zw, Ncx, Nfx, dx, dy, Cx, Cy, sx)
%
% **********************************************************************************************
% ******* Calibrating a Camera Using a Monoview Coplanar Set of Points *******
% **********************************************************************************************
% 6/2004 Simon Wan
% simonwan@hit.edu.cn
%
% Note: Xf, Yf, xw, yw, zw are all column vectors
%
% (xw, yw, zw) is the 3D coordinate of the object point P in the 3D world coordinate system
% (x, y, z) is ths 3D coordinate of the object point P in the 3D camera coordinate system
% (X, Y) is the image coordinate system centered at Oi where is the intersection of the optical center axis z and the front plane
% (Xu, Yu) is the image coordinate of (x, y, z) if a perfect pinhole
% camera model is used(利用相机坐标表示的图像坐标)
% Xu = f * x / z (4a)
% Yu = f * y / z (4b)
% (Xd, Yd) is the actual image coordinate which differs from (Xu, Yu)
% due to lens distortion实际的图像坐标(加了畸变)
% (Xf, Yf) is the coordinate used in the computer, is the number of
% pixels for the discrete image in the frame memory(像素坐标)
% R is the 3*3 rotation matrix
% = [r1, r2, r3; r4, r5, r6; r7, r8, r9]; (2)
% [x, y, z]' = R * [xw, yw, zw]' + T (1)
% T is the translation vector
% = [Tx, Ty, Tz]' (3)
% f is the effective focal length 是有效的焦距
% Dx = Xd*( k1*r^2 + k2*r^4 + ... ) P327
% Xd+Dx=Xu (5a)
% Dy = Yd*( k1*r^2 + k2*r^4 + ... ) P327
% Yd+Dy=Yu (5b)
% r = (Xd^2 + Yd^2)^(0.5) P327
% k1 is the distortion coeffient
% Xf = sx * dxp^(-1) * Xd + Cx (6a)
% Yf = dy^(-1) * Yd + Cy (6b)
% dxp = dx * Ncx / Nfx (6d)
% dx is the center to center distance between adjacent sensor elements in X (scan line) diretion是X(扫描线)方向上相邻传感器元件之间的中心距离
% dy is the center to center distance between adjacent CCD sensor in the Y direction是相邻CCD传感器在Y方向上的中心距离
% Ncx is the number of sensor elements in the X direction是X方向上的传感器元件的数量,平面像素的采样频率
% Nfx is the number of pixels in a line as sampled by the computer是计算机采样的行中像素数,表明平面像素的采样频率
% sx is the uncertainty image scale factor是不确定性图像比例因子,相当于误差项,也设置为1
% X = (Xd * Nfx) / (dx * Ncx) P328
% X = Xf - Cx P328
% Y = Yf - Cy P328
% sx^(-1)*dxp*X + sz^(-1)*dxp*X*k1*r^2 = f*x/z (7a)
% dxp*Y + dy*Y*k1*r^2 = f*y/z (7b)
% r = ( ( sx^(-1)*dxp*X )^2 + (dx*Y)^2 )^(0.5)
% sx^(-1)*dxp*X + sx^(-1)*dxp*X*k1*r^2 = f*(r1*xw + r2*yw + r3*zw +
% Tx) / (r7*xw + r8*yw + r9*zw +Tz) (8a)
% dy*Y + dy*Y*k1*r^2 = f*(r1*xw + r2*yw + r3*zw +
% Tx) / (r7*xw + r8*yw + r9*zw +Tz) (8b)
% Since the calibration points are on a common plane, the (xw, yw, zw) coordinate system can be chosen such that zw=0 and the
% corigin is not lose to the center of the view or y axis of the camera coordinate system. Since the (xw, yw, zw) is user-defined
% and the origin is arbitrary, it is no problem setting the origin of (xw, yw, zw) to be out of the field of view and not close
% to the y axis. the purpose for the latter is to make sure that Ty is not exactly zero.
%
% REF: "A versatile camera calibration technique for high-accuracy 3D machine
% vision metrology using off-the-shelf TV cameras and lens"
% R.Y. Tsai, IEEE Trans R&A RA-3, No.4, Aug 1987, pp 323-344.
%
function [R, T, f, k1] = Tsai(Xf, Yf, xw, yw, zw, Ncx, Nfx, dx, dy, Cx, Cy, sx)
% Stage 1 --- Compute 3D Orientation, Position (x and y):
% a) Compute the distored image coordinates (Xd, Yd) Procedure:
dxp = dx * Ncx / Nfx;
X = Xf - Cx;
Y = Yf - Cy;
Xd=sx^(-1)*dxp*(Xf-Cx);
Yd=dy*(Yf-Cy);
% b) Compute the five unknowns Ty^(-1)*r1, Ty^(-1)*r2, Ty^(-1)*Tx, Ty^(-1)*r4, Ty^(-1)*r5
% r1p=Ty^(-1)*r1;
% r2p=Ty^(-1)*r2;
% Txp=Ty^(-1)*Tx;
% r4p=Ty^(-1)*r4;
% r5p=Ty^(-1)*r5;
A=[Yd.*xw Yd.*yw Yd -Xd.*xw -Xd.*yw];
B=Xd;
C=A\B;
r1p=C(1);
r2p=C(2);
Txp=C(3);
r4p=C(4);
r5p=C(5);
clear A B C;
% c) Compute (r1,...,r9,Tx,Ty) from (Ty^(-1)*r1, Ty^(-1)*r2, Ty^(-1)*Tx, Ty^(-1)*r4, Ty^(-1)*r5):
% 1) Compute |Ty| from (Ty^(-1)*r1, Ty^(-1)*r2, Ty^(-1)*Tx, Ty^(-1)*r4, Ty^(-1)*r5):
C=[r1p, r2p; r4p, r5p];
Sr=r1p^2 + r2p^2 + r4p^2 + r5p^2;
if rank(C)==2
Ty2=( Sr - (Sr^2-4*(r1p*r5p-r4p*r2p)^2)^(0.5) )/(2*(r1p*r5p-r4p*r2p)^2);
else
z = C(abs(C) > 0);
Ty2 = 1.0 / (z(1)^2 + z(2)^2);
end
Ty = sqrt(Ty2);
clear C Sr Ty2 z
% 2) Determine the sign of Ty:
[ymax i] = max(Xd.^2 + Yd.^2);
r1 = r1p*Ty;
r2 = r2p*Ty;
r4 = r4p*Ty;
r5 = r5p*Ty;
Tx = Txp*Ty;
x = r1*xw(i) + r2*yw(i) + Tx;
y = r4*xw(i) + r5*yw(i) + Ty;
% if (sign(x) == sign(Xf(i))) & (sign(y) == sign(Yf(i))),
if (sign(x) == sign(X(i))) && (sign(y) == sign(Y(i))),
Ty = Ty;
else
Ty = -Ty;
end
clear ymax i x y
% 3) Compute the 3D rotation matrix R, or r1, r2,...,r9
r1 = r1p*Ty;
r2 = r2p*Ty;
r4 = r4p*Ty;
r5 = r5p*Ty;
Tx = Txp*Ty;
s = -sign(r1*r4 + r2*r5);
R=[r1, r2, (1-r1^2-r2^2)^(0.5); r4, r5, s*(1-r4^2-r5^2)^(0.5)];
R = [R(1:2,:); cross(R(1,:), R(2,:))];
r7 = R(3,1);
r8 = R(3,2);
r9 = R(3,3);
y = r4*xw+r5*yw+Ty;
w = r7*xw+r8*yw;
z = [y -dy*Y] \ [dy*(w.*Y)];
f = z(1);
if f < 0
R(1,3) = -R(1,3);
R(2,3) = -R(2,3);
R(3,1) = -R(3,1);
R(3,2) = -R(3,2);
end
r3 = R(1,3);
r6 = R(2,3);
r7 = R(3,1);
r8 = R(3,2);
clear s y w z
% 2) Stage 2 --- Compute Effective Focal Length, Distortion Coefficients, and z Position:
% d) Compute an approximation of f and Tz by ignoring lens distortion:
y = r4*xw+r5*yw+Ty;
w = r7*xw+r8*yw;
z = [y -dy*Y] \ [dy*(w.*Y)];
f = z(1);
Tz = z(2);
% Compute the exactly solution for f, Tz, k1:
params_const = [r4 r5 r6 r7 r8 r9 dx dy sx Ty];
params = [f, Tz, 0]; % add initial guess for k1
[x,fval,exitflag,output] = fminsearch( @Tsai_8b, params, [], params_const, xw, yw, zw, X, Y);
f = x(1);
Tz = x(2);
k1 = x(3);
T=[Tx, Ty, Tz]';
% fval the value of the objective function fun at the solution x.
fval
% exitflag that describes the exit condition of fminsearch
% >0 Indicates that the function converged to a solution x.
% 0 Indicates that the maximum number of function evaluations was exceeded.
% <0 Indicates that the function did not converge to a solution.
exitflag
% output that contains information about the optimization
% output.algorithmThe algorithm used
% output.funcCountThe number of function evaluations
% output.iterationsThe number of iterations taken
output
Tsai_8b.m
% f = Tsai_8b(params, params_const, xw, yw, zw, X, Y)
%
% **********************************************************************************************
% ******* Calibrating a Camera Using a Monoview Coplanar Set of Points *******
% **********************************************************************************************
% 6/2004 Simon Wan
% simonwan@hit.edu.cn
%
% Note: This is not called directly but as a function handle from the "fminsearch "
%
function f = Tsai_8b(params, params_const, xw, yw, zw, X, Y)
% unpack the params
f = params(1);
Tz = params(2);
k1 = params(3);
% unpack the params_const
r4 = params_const(1);
r5 = params_const(2);
r6 = params_const(3);
r7 = params_const(4);
r8 = params_const(5);
r9 = params_const(6);
dx = params_const(7);
dy = params_const(8);
sx = params_const(9);
Ty = params_const(10);
rsq = (dx*X).^2 + (dy*Y).^2;
res = (dy*Y).*(1+k1*rsq).*(r7*xw+r8*yw+r9*zw+Tz) - f*(r4*xw+r5*yw+r6*zw+Ty);
f = norm(res, 2);
5.结果:
运行代码后,结果及其误差分析会在当前目录下生成txt文件。
matlab实现TSAI标定及其误差分析相关推荐
- 基于matlab测量物体直径,基于MATLAB的零件尺寸检测误差分析的软件设计
ISSN1672-4305 CN12-1352/N 实 验 室 科 学 LABORATORY SCIENCE 第 14 卷 第 4 期 2011 年 8 月 Vol. 14 No. 4 Aug. 20 ...
- matlab双目相机标定校正_基于双目视觉的无人机避障算法(一)
讲述在10月到12月所做的所有工作 对于一个无人机自主避障来说,存在着以下流程: 感知:障碍物检测.行人检测.目标检测 SLAM:为无人机提供位置估计,构建稀疏环境地图 路径规划:规划一条从当前位置到 ...
- 利用Matlab进行相机标定并使用openCV进行简单三维重建
注:本文主要针对Matlab和OpenCV跨平台进行相机标定.单相机三维重建工作的实现,因为我发现网上竟然没有一篇博客径直指出这两者在进行图像处理时的巨大差异(坐标系完全不同),不然我也不会走了很多弯 ...
- matlab 光平面标定代码,线结构光视觉传感器/线激光深度传感器标定工具(matlab)...
线结构光视觉传感器/线激光深度传感器标定工具(matlab) 发布时间:2018-03-13 17:32, 浏览次数:655 , 标签: matlab 线结构光视觉系统有着结构简单.使用灵活.受周围光 ...
- 三维扫描仪[4]——标定·Matlab单目标定
前半部分由孤旅人翻译,非常感谢. 原文链接相机校准:Toolbox for Matlab 这部分带你们完成一个基于20张平面棋盘图的完整的标定例程. 这个例程将让你们学会怎么使用关于这个工具箱所有的功 ...
- matlab双目相机标定校正_双目相机的标定过程详解!-----MATLAB
基于双目视觉的测距.三维重建等过程中的第一步就是要进行标定.双目相机的标定过程在网上有很多资料,但是基本都没有matlab官方网址讲的好.所以请参考MATLAB官方文档:https://ww2.mat ...
- MATLAB 绘制圆点标定板
前言 本文主要利用Matlab中的rectangle函数,绘制圆点标定板.主要描述如何理解和设计圆点的大小.数量,如何建立标定板的原点坐标系.使得标定板正中央在原点上. 正文 建立圆点参数 如何使用M ...
- 相机标定的理解及采用opencv和matlab工具箱的标定方法
一.相机标定的目的 确定空间物体表面某点的三维几何位置与其在图像中对应点之间的相互关系,建立摄像机成像的几何模型,这些几何模型参数就是摄像机参数. 二.通用摄像机模型 世界坐标系.摄像机坐标系和像平面 ...
- matlab相机标定工具箱进行相机标定,相机标定的理解及采用opencv和matlab工具箱的标定方法...
一.相机标定的目的 确定空间物体表面某点的三维几何位置与其在图像中对应点之间的相互关系,建立摄像机成像的几何模型,这些几何模型参数就是摄像机参数. 二.通用摄像机模型 世界坐标系.摄像机坐标系和像平面 ...
最新文章
- 【MySQL】记一次线上重大事故:二狗子竟然把线上数据库删了!!
- (无奈之举)2011年01月26号日志
- 客户端JavaScript加密数据,服务端Java解密数据
- python epoll多路复用技术_python:多路复用+零拷贝
- python降温了吗_为“Python将纳入高考”降温
- VsCode常用设置,新手必备!
- 荆门市掇刀石中学2021高考成绩查询,2021年荆门中考所有高中学校排名 荆门重点高中分数线...
- iphone ios编译ffmpeg
- asp.net中获得客户端的IP地址
- 《常微分方程教程》习题2.3.6
- Red Hat 5.6-64位使用yum升级PHP
- 在OpenSSL中添加自定义加密算法
- python怎么看自己安装的第三方包_安装第三方包查看python版本/第三方包版本
- C# 读取Word中的书签
- 济南:建设中国软件名城
- 美国纽约摄影学院摄影教材 学习笔记1
- Druid (大数据实时统计分析数据存储)
- FFMPEG AvFilter使用实例(实现视频缩放,裁剪,水印等)
- png格式的icon使用方法
- 论文笔记 | code pretraining(代码预训练系列)
热门文章
- win7蓝屏0x000000f4修复_关于近期财务电脑蓝屏的处理和预防方法,必看!
- 30本世界名著浓缩成的经典话语
- 获取浏览器原生滚动条宽度的方法
- Nginx - 正反向代理,nginx.conf配置
- 25幅精美绝伦的光涂鸦摄影作品欣赏
- PHP爬取微信公众号文章(可做为扩展类直接使用)
- ssh连接linux网段不同_ssh 多台服务器之间连接(linux)
- 基于共享存储的Harbor高可用-Docker部署方案
- 操作系统之共享存储区详解(Linux shm)
- python导入excel到access_如何使用Python访问(读取和写入)Google表格电子表格?