基于两步法的数码相机定位-毕业论文

(23)算出的含参数的测量值和实际距离值进行比较，便可计算出值，从而求得矩阵。模型和算法的优缺点分析及改进方向由于数码照相机模型的参数太多，而且是一个复杂的非线性方程，因此一步求解不仅繁冗，而且稳定性差；本文所采用的两步法确定像坐标和两相机相对位值的方法，不仅考虑了径向畸变，还考虑了切向畸变因素，求得较为准确的像坐标和相机相对位置，精度高，稳定性好，解的收敛速度快。由于考虑了相机的径向畸变和切向畸变即非线性畸变，计算程序略显复杂，计算量稍大。由于BP神经网络在逼近能力，收敛速度等方面都优于其它的非线性优化算法，很适宜于相机的标定。将神经网络技术用于相机的标定，不需要预先确定摄像机具体内外参数及有关模型的先验知识，利用BP神经网络对复杂非线性映射关系强大的逼近能力来对视觉系统中的摄像机进行隐式标定，算法标定精度较高，可靠性好，实用性较强，是相机定位的一种趋势。参考文献ADDINNE.Bib[1]刘红宁.基于双目立体视觉拟人机器人定位关键技术的研究[D].河北工业大学,2006.[2]马颂德,张正友.计算机视觉[M].北京:科学技术出版社,1998.[3]毛剑飞,邹细勇,诸静.改进的平面模板两步法标定摄像机[J].中国图象图形学报,2004,9(7):847~852.[4]张磊,刘小力,罗翔,田梦倩,史金飞.一种考虑畸变的摄像机线性标定简化方法[J].工业仪表与自动化装置,2007(1):3~7.[5]陈杰.Matlab宝典[M].北京:电子工业出版社,2007.[6]李海滨,单文军,刘彬.双目立体视觉测距系统误差模型的研究[J].光学技术,2006(01):13~16.[7]张健新,段发阶,叶声华.”两步法“求取双目视觉传感器中摄象机位置关系[J].光电工程,1998,25(5):37~41.[8]李雪鹃,蔡怀宇,黄战华.基于双目视差的测距系统设计[J].电子测量技术,2007(08):31~34.附录基于Canny算子边界检测源程序1%基于Canny算子求边界坐标%参数name为图片名function[x,y]=tp(name)I=imread(name);%读入图片I=double(I);I=I(1:768,1:1024);G=edge(I,'Canny');%基于Canny算子求边界imshow(G);I=uint8(I);G=uint8(G);subplot(1,2,1)imshow(I);%输出图像subplot(1,2,2)imshow(G);%输出图像x=[];y=[];fori=100:600%求边界坐标forj=200:700ifG(i,j)==0）y=[y,i];x=[x,j];endendend基于简化模型求取像坐标初值源程序2%基于简化模型求解像圆心%参数x，y为物的圆心坐标u为物距function[x0,y0]=jiandan(x,y,u)x0=[];y0=[];v=1577;%给像距赋初值h1=24;f=ff(h1,h2,v)%调用求解焦距的函数，n=length(u);forj=1:n;k=f/(u(j)+v)*3.78;%求比例系数a=x(j)*k;%求像坐标并保存b=y(j)*k;x0=[x0,a];y0=[y0,b];end两步法求解像坐标源程序3%两步法求解像坐标%读入数据并整理成可处理的形式%M为靶标数据%m为靶标的像数据%求解焦距functionf=ff(h1,h2,v)n=length(h2);fork=1:nu(k)=h1*v/h2(k);%求每一个像高对应的物距endu=u*3.78;fork=1:nF(k)=1/(1/v+1/u(k));%求每一个物距对应的焦距endf=mean(F);%求所以焦距的算术平平均clearall;closeall;M=imread('tu.bmp');%读入图片m1=imread('tu1.bmp');m2=imread('tu2.bmp');m3=imread('tu3.bmp');m4=imread('tu4.bmp');m5=imread('tu5.bmp');MM=[M(:,1:2);M(:,3:4);M(:,5:6);M(:,7:8)];%提取数据mm1=[m1(:,1:2);m1(:,3:4);m1(:,5:6);m1(:,7:8)];mm2=[m2(:,1:2);m2(:,3:4);m2(:,5:6);m2(:,7:8)];mm3=[m3(:,1:2);m3(:,3:4);m3(:,5:6);m3(:,7:8)];mm4=[m4(:,1:2);m4(:,3:4);m4(:,5:6);m4(:,7:8)];mm5=[m5(:,1:2);m5(:,3:4);m5(:,5:6);m5(:,7:8)];M=[M(:,1:2);M(:,3:4);M(:,5:6);M(:,7:8)];m1=[m1(:,1:2);m1(:,3:4);m1(:,5:6);m1(:,7:8)];m2=[m2(:,1:2);m2(:,3:4);m2(:,5:6);m2(:,7:8)];m3=[m3(:,1:2);m3(:,3:4);m3(:,5:6);m3(:,7:8)];m4=[m4(:,1:2);m4(:,3:4);m4(:,5:6);m4(:,7:8)];m5=[m5(:,1:2);m5(:,3:4);m5(:,5:6);m5(:,7:8)];M=M';m(:,:,1)=m1';m(:,:,2)=m2';m(:,:,3)=m3';m(:,:,4)=m4';m(:,:,5)=m5';%M为靶标数据%m为靶标的像数据functionliangbufa(M,m)%M=[X,Y]’转化为M=[X,Y,1]'%m=[u,v]'转化为m=[u,v,1]'[rows,npts]=size(M);matrixone=ones(1,npts);M=[M;matrixone];num=size(m,3);fori=1:numm(3,:,i)=matrixone;endH=jiandan(x,y,u)%调用简单模型求解初值V=[];forflag=1:num%迭代v12(:,:,flag)=[H(1,1,flag)*H(2,1,flag),H(1,1,flag)*H(2,2,flag)+H(1,2,flag)*H(2,1,flag),H(1,2,flag)*H(2,2,flag),H(1,3,flag)*H(2,1,flag)+H(1,1,flag)*H(2,3,flag),H(1,3,flag)*H(2,2,flag)+H(1,2,flag)*H(2,3,flag),H(1,3,flag)*H(2,3,flag)];v11(:,:,flag)=[H(1,1,flag)*H(1,1,flag),H(1,1,flag)*H(1,2,flag)+H(1,2,flag)*H(1,1,flag),H(1,2,flag)*H(1,2,flag),H(1,3,flag)*H(1,1,flag)+H(1,1,flag)*H(1,3,flag),H(1,3,flag)*H(1,2,flag)+H(1,2,flag)*H(1,3,flag),H(1,3,flag)*H(1,3,flag)];v22(:,:,flag)=[H(2,1,flag)*H(2,1,flag),H(2,1,flag)*H(2,2,flag)+H(2,2,flag)*H(2,1,flag),H(2,2,flag)*H(2,2,flag),H(2,3,flag)*H(2,1,flag)+H(2,1,flag)*H(2,3,flag),H(2,3,flag)*H(2,2,flag)+H(2,2,flag)*H(2,3,flag),H(2,3,flag)*H(2,3,flag)];V=[V;v12(:,:,flag);v11(:,:,flag)-v22(:,:,flag)];endk=V'*V;%将Vb=0转化为kb=V'*V*b=0,这样k就是一个对称正定矩阵，[u,v,d]=svd(k);[e,d2]=eig(k);b=d(:,6);v0=(b(2)*b(4)-b(1)*b(5))/(b(1)*b(3)-b(2)^2);s=b(6)-(b(4)^2+v0*(b(2)*b(4)-b(1)*b(5)))/b(1);alpha_u=sqrt(s/b(1));alpha_v=sqrt(s*b(1)/(b(1)*b(3)-b(2)^2));skewness=-b(2)*alpha_u*alpha_u*alpha_v/s;u0=skewness*v0/alpha_u-b(4)*alpha_u*alpha_u/s;A=[alpha_uskewnessu00alpha_vv0001];D=[];%求解所有外部参数d=[];Rm=[];forflag=1:nums=(1/norm(inv(A)*H(1,:,flag)')+1/norm(inv(A)*H(2,:,flag)'))/2;rl1=s*inv(A)*H(1,:,flag)';rl2=s*inv(A)*H(2,:,flag)';rl3=cross(rl1,rl2);RL=[rl1,rl2,rl3];[U,S,V]=svd(RL);RL=U*V';TL=s*inv(A)*H(3,:,flag)';TL是位移矩阵t(外参)RT=[rl1,rl2,TL];%H=A[r1r2t]XY=RT*M;%M是modelplane点的坐标UV=A*XY;%sm=A[Rt]M,UV是等式的右边UV=[UV(1,:)./UV(3,:);UV(2,:)./UV(3,:);UV(3,:)./UV(3,:)];XY=[XY(1,:)./XY(3,:);XY(2,:)./XY(3,:);XY(3,:)./XY(3,:)];endk=inv(D'*D)*D'*d;para=[Rm,k(1),k(2),alpha_u,skewness,u0,alpha_v,v0];options=optimset('LargeScale','off','LevenbergMarquardt','on');%用最小二乘法求解[x,resnorm,residual,exitflag,output]=lsqnonlin(@simon_HHH,para,[],[],options,m,M);k1=x(num*6+1)%显示结果k2=x(num*6+2)A=[x(num*6+3)x(num*6+4)x(num*6+5);0x(num*6+6)x(num*6+7);0,0,1]基于Canny算子边界检测确定圆心源程序4%基于Canny算子求出的边界坐标来拟合圆心%参数x，y为边界坐标function[x0,y0]=yuanxin(x,y)x1=(min(x)+max(x))/2;%估计圆心的坐标y1=(min(y)+max(y))/2;X=[];Y=[];fora=1:length(x)-2%任意取三个边界点确定一个圆心的坐标forb=a+1:length(x)-1forc=b+1:length(x)x0=((x(b)^2-x(a)^2+y(b)^2-y(a)^2)*(y(b)-y(a))/2-(x(a)^2-x(c)^2+y(a)^2-y(c)^2)*(y(a)-y(b))/2)/((x(a)-x(b))*(y(a)-y(c))-(x(a)-x(c))*(y(a)-y(b)));y0=((x(a)^2-x(c)^2+y(a)^2-y(c)^2)*(x(a)-x(b))/2-(x