circfit.m

基于最小二乘估计的圆心拟合

i=imread('1.bmp');                    %      输入图片
g=rgb2gray(i);                        %      转为灰度图像
bw=edge(g,'sobel');                   %       边缘检测
imshow(bw)
pixval                                %      显示坐标
j=1;                                  %      同一圆周内的搜索边缘坐标
for m=20:1:50
   for n=12:1:42
       k=impixel(bw,m,n);              %     边缘点上的impixel值为1，不在边缘上的impixel值为0
       if k==1                         %      循环检测，直至搜索完整个区域
           x11(j)=m;y11(j)=n;
           j=j+1;
       end
   end
end
x11,y11
function [a,b,r]=circfit(x,y)                            %          a,b为圆心坐标，r为半径
n=length(x);
c1=sum(x);c2=sum(y);c3=sum(x.*x);c4=sum(y.*y);             %最小二乘估计
c5=sum(x.*y);c6=sum(x.*x.*x);c7=sum(y.*y.*y);
c8=sum(x.*y.*y);c9=sum(x.*x.*y);
d1=2*(c3-c1.*c1/n);d2=2*(c5-c1*c2/n);d3=c6+c8-c1*(c3+c4)/n;
d4=2*(c4-c2*c2/n);d5=c9+c7-c2*(c3+c4)/n;
a=(d3*d4-d2*d5)/(d1*d4-d2*d2);
b=(d1*d5-d2*d3)/(d1*d4-d2*d2);
r=sqrt(((c3+c4-2*a*c1-2*b*c2)/n)+a*a+b*b);