视觉测量实验报告

实验报告
实验名称典型零件尺寸的图像处理与视觉测量专业班级测控
学号2009013003
姓名侯爵
指导教师刘力双
实验目的：
通过本实验，给学生提供边学习边实践的机会，应用平时所学的知识，通过一些编程软件处理一些图像，观察处理之后的效果，可以巩固所学知识，提高学生的动手能力和编程水平。

实验要求：
已修过《数字图像处理》和《视觉测量技术》课程。

实验报告要求：
每人提交一份电子版实验报告，报告内容包括：封皮，实验要求，实验思路和过程，实验结果，实验遇到的问题及解决方法，附上程序。

实验内容：
采用MATLAB或VC++编制程序，实现图像中尺寸参数的测量。

必做：
（1）求图1中直线的参数方程，按照y=kx+b的形式给出。

（难度：**）（思路：先提取出边缘点，然后用边缘点拟合直线）
图1 求出图中直线的参数
clc;
clear;
close all;
I=imread('straight_line.jpg'); %读取图像
I1=im2double(I); %将彩图序列变成双精度
I2=rgb2gray(I1); %将彩色图变成灰色图
[thr, sorh, keepapp]=ddencmp('den','wv',I2);
I3=wdencmp('gbl',I2,'sym4',2,thr,sorh,keepapp); %小波除噪
%I4=medfilt2(I3,[9 9]); %中值滤波
%I5=imresize(I4,0.3,'bicubic'); %图像大小
%I6=imresize(I2,0.3,'bicubic'); %图像大小
BW1=edge(I3,'sobel'); %sobel图像边缘提取
[y,x]=find(BW1);
m=polyfit(x,y,1);
yy=polyval(m,x);
%k=(yy(1)-m)/x(1);
func=['y=',num2str(m(2)),'+',num2str(m(1)),'*x'];
subplot(1,2,1); %图划分为一行三幅图，第一幅图
imshow(I2); %绘图
title('原始图像 '); %标注
subplot(1,2,2); %第二幅图
imshow(I3);title('消噪后图像 ');
%subplot(1,3,3); %第三幅图
%imshow(I4);title('中值滤波图像 ');
figure;生成窗口
subplot(1,2,1);
imshow(BW1);
title('Sobel 算子 ');
subplot(1,2,2);
imshow(I3);
hold on;
plot(x,yy,'r-','linewidth',2);
xlabel('x');
ylabel('y');
text(x(2),yy(1)+5,func);
title('拟合直线 '); %标注
（2）按照拟合法求图2中圆的参数，即给出圆心位置和半径。

（难度：**）（思路：先提取出边缘点，然后用边缘点拟合圆）
图2 求出图中圆的圆心位置和半径
clc;
clear;
close all;
I=imread('circle.jpg'); %读取图像
I1=im2double(I); %将彩图序列变成双精度
I2=rgb2gray(I1); %将彩色图变成灰色图
[thr, sorh, keepapp]=ddencmp('den','wv',I2);
I3=wdencmp('gbl',I2,'sym4',2,thr,sorh,keepapp); %小波除噪BW1=edge(I3,'sobel'); %sobel图像边缘提取
[y,x]=find(BW1);
subplot(1,2,1); %图划分为一行三幅图，第一幅图
imshow(I2); %绘图
title('原始图像 '); %标注
subplot(1,2,2); %第二幅图
imshow(I3);title('消噪后图像 ');
figure;生成窗口
subplot(1,2,1);
imshow(BW1);
title('Sobel 算子 ');
N=length(x);
x1 = 0;
x2 = 0;
x3 = 0;
y1 = 0;
y2 = 0;
y3 = 0;
x1y1 = 0;
x1y2 = 0;
x2y1 = 0;
for i = 1 : N
x1 = x1 + x(i);
x2 = x2 + x(i)*x(i);
x3 = x3 + x(i)*x(i)*x(i);
y1 = y1 + y(i);
y2 = y2 + y(i)*y(i);
y3 = y3 + y(i)*y(i)*y(i);
x1y1 = x1y1 + x(i)*y(i);
x1y2 = x1y2 + x(i)*y(i)*y(i);
x2y1 = x2y1 + x(i)*x(i)*y(i);
end
C = N * x2 - x1 * x1;
D = N * x1y1 - x1 * y1;
E = N * x3 + N * x1y2 - (x2 + y2) * x1;
G = N * y2 - y1 * y1;
H = N * x2y1 + N * y3 - (x2 + y2) * y1;
a = (H * D - E * G)/(C * G - D * D);
b = (H * C - E * D)/(D * D - G * C);
c = -(a * x1 + b * y1 + x2 + y2)/N;
A = a/(-2); %x坐标
B = b/(-2); %y坐标
R = sqrt(a * a + b * b - 4 * c)/2; subplot(1,2,2);
imshow(I3);
hold on;
X=(A-R):0.1:(A+R);
for i=1:(2*10*R+1);
Ya(i)=sqrt(R^2-(X(i)-A)^2)+B ;
Yb(i)=-sqrt(R^2-(X(i)-A)^2)+B ;
end
a=length(X);
for i=1:a
xx(i)=X(a+1-i);
yy(i)=Yb(a+1-i);
end
X=[X xx];
Y=[Ya yy];
plot(X,Y,'r-','linewidth',2);
func=[num2str(R),'^2=(x-',num2str(A),')^2+(y-',num2str(B),')^2'];
text(5,140,func);
text(5,155,['R=',num2str(R)]);
text(5,170,['圆心坐标（',num2str(A),',',num2str(B),')']);
xlabel('x');
ylabel('y');
title('拟合圆');
（3）按照拟合法求图3中圆弧的参数，即给出圆心位置和半径。

（难度：**）（思路：先提取出边缘点，然后用边缘点拟合圆）
图3 求出图中圆弧的圆心位置和半径
clc;
clear;
close all;
I=imread('half_circle.jpg'); %读取图像
I1=im2double(I); %将彩图序列变成双精度
I2=rgb2gray(I1); %将彩色图变成灰色图
[thr, sorh, keepapp]=ddencmp('den','wv',I2);
I3=wdencmp('gbl',I2,'sym4',2,thr,sorh,keepapp); %小波除噪
BW1=edge(I3,'sobel'); %sobel图像边缘提取
[y,x]=find(BW1);
subplot(1,2,1); %图划分为一行三幅图，第一幅图
imshow(I2); %绘图
title('原始图像 '); %标注
subplot(1,2,2); %第二幅图
imshow(I3);title('消噪后图像 ');
figure;
subplot(1,2,1);
imshow(BW1);
title('Sobel 算子 ');
N=length(x);
x1 = 0;
x2 = 0;
x3 = 0;
y1 = 0;
y2 = 0;
y3 = 0;
x1y1 = 0;
x1y2 = 0;
x2y1 = 0;
for i = 1 : N
x1 = x1 + x(i);
x2 = x2 + x(i)*x(i);
x3 = x3 + x(i)*x(i)*x(i);
y1 = y1 + y(i);
y2 = y2 + y(i)*y(i);
y3 = y3 + y(i)*y(i)*y(i);
x1y1 = x1y1 + x(i)*y(i);
x1y2 = x1y2 + x(i)*y(i)*y(i);
x2y1 = x2y1 + x(i)*x(i)*y(i);
end
C = N * x2 - x1 * x1;
D = N * x1y1 - x1 * y1;
E = N * x3 + N * x1y2 - (x2 + y2) * x1;
G = N * y2 - y1 * y1;
H = N * x2y1 + N * y3 - (x2 + y2) * y1;
a = (H * D - E * G)/(C * G - D * D);
b = (H * C - E * D)/(D * D - G * C);
c = -(a * x1 + b * y1 + x2 + y2)/N;
A = a/(-2); %x坐标
B = b/(-2); %y坐标
R = sqrt(a * a + b * b - 4 * c)/2;
subplot(1,2,2);
imshow(I3);
hold on;
X=(A-R):0.1:(A+R);
for i=1:(2*10*R+1);
Ya(i)=sqrt(R^2-(X(i)-A)^2)+B ;
Yb(i)=-sqrt(R^2-(X(i)-A)^2)+B ;
end
a=length(X);
for i=1:a
xx(i)=X(a+1-i);
yy(i)=Yb(a+1-i);
end
X=[X xx];
Y=[Ya yy];
plot(X,Y,'r-','linewidth',2);
func=[num2str(R),'^2=(x-',num2str(A),')^2+(y-',num2str(B),')^2'];
text(5,205,func);
text(5,225,['R=',num2str(R)]);
text(5,247,['圆心坐标（',num2str(A),',',num2str(B),')']);
xlabel('x');
ylabel('y');
title('拟合圆弧');
（4）求图4中光斑的参数，即给出光斑中心和光斑大小。

（难度：**）（思路：采用重心法计算光斑中心，然后采用二阶矩法计算光斑大小）
图4 求出图中光斑的中心和大小clc;
clear;
close all;
I=imread('light_spot.jpg'); %读取图像
I1=im2double(I); %将彩图序列变成双精度
I2=rgb2gray(I1); %将彩色图变成灰色图
[thr, sorh, keepapp]=ddencmp('den','wv',I2);
I3=wdencmp('gbl',I2,'sym4',2,thr,sorh,keepapp); %小波除噪
subplot(1,2,1); %图划分为一行三幅图，第一幅图
imshow(I2); %绘图
title('原始图像'); %标注
subplot(1,2,2); %第二幅图
imshow(I3);title('消噪后图像');
[a,b]=size(I3);
x=1:b;
for i=1:a
X(i,:)=x;
end
y=1:a;
y=y';
for j=1:b
Y(:,j)=y;
end
A=I3.*X;
B=I3.*Y;
A_SUM=sum(A);
B_SUM=sum(B);
C_SUM=sum(I3);
Xc=A_SUM/C_SUM
Yc=B_SUM/C_SUM
J=I3.*((X-Xc).^2);
K=I3.*((Y-Yc).^2);
J_SUM=sum(J);
K_SUM=sum(K);
ox=sqrt(J_SUM/C_SUM); oy=sqrt(K_SUM/C_SUM); Wx=2*ox;
Wy=2*oy;
R=(Wx+Wy)/2;
Xt=(Xc-R):0.1:(Xc+R); for i=1:(2*10*R+1);
Ya(i)=sqrt(R^2-(Xt(i)-Xc)^2)+Yc ;
Yb(i)=-sqrt(R^2-(Xt(i)-Xc)^2)+Yc ;
end
a=length(Xt);
for i=1:a
xx(i)=Xt(a+1-i);
yy(i)=Yb(a+1-i);
end
XX=[Xt xx];
YY=[Ya yy];
figure;
subplot(1,2,1);
imshow(I3);
hold on;
plot(XX,YY,'r-','linewidth',2);
func=[num2str(R),'^2=(x-',num2str(Xc),')^2+(y-',num2str(Yc),')^2']; text(250,20,func);
text(250,40,['R=',num2str(R)]);
text(250,60,['圆心坐标（',num2str(Xc),',',num2str(Yc),')']);
xlabel('x');
ylabel('y');
title('由所求参数得到的圆');
9。