chuangxinyanjiu.asv

实现对小孔的高斯滤波

function pic
%r为圆的半径，【rx，ry】为圆心坐标
 f=imread('4ms.bmp');
 subplot(221)
 imshow(f);
 title('the original image')
 [M,N]=size(f);
 F=fft2(f);
 a=20*log(1+abs(fftshift(F)));
 subplot(222)
 imshow(a,[])
 title('the frequency domain of the original image')
 sig=10;
 H=lpfilter('gaussian',M,N,sig);%生成高斯低通滤波器
 G=H.*F;
 g=real(ifft2(G));
 subplot(223)
 imshow(g,[])
 title('the image after gaussian lowpass filter')
 h=edge(g,'sobel');%采用sobel算子进行边缘提取，方向默认为both.
 subplot(224)
 imshow(h,[])
 title('the edge of the image')
 h1=im2uint8(h);
 hadd=f+h1;
 figure
 subplot(221)
 imshow(hadd)%检验所得的圆与原始圆的位置是否发生了变化
 title('the add of the edge and the original image')
 ma=zeros(sum(h(:)),2);%ma用来存放边缘轮廓的点坐标
 k=1;
 for x=1:480
     for y=1:640
         if h(x,y)>0
             ma(k,:)=[x,y];
             k=k+1;
         end
     end
 end
 %横坐标存放在第一列，纵坐标存放在第二列
 
 %用最小二乘法计算
 N=sum(h(:));%所有边缘点的个数
 x1=ma(:,1);%横坐标
 y1=ma(:,2);%纵坐标
 A1=(sum(x1)/N).^2-sum(x1.^2)/N;
 A2=(sum(x1)/N)*(sum(y1)/N)-sum(x1.*y1)/N;
 A3=sum(x1.^2)/N*(sum(x1)/N)+(sum(x1)/N)*(sum(y1.^2)/N)-sum(x1.^3)/N-sum(x1.*(y1.^2))/N;
 B1=(sum(x1)/N)*(sum(y1)/N)-sum(x1.*y1)/N;
 B2=(sum(y1)/N)^2-sum(y1.^2)/N;
 B3=(sum(x1.^2)/N)*(sum(y1)/N)+(sum(y1)/N)*(sum(y1.^2)/N)-(sum(x1.^2)/N)*(sum(y1)/N)-sum(y1.^3)/N;
 a=(A3*B2-B3*A2)/(2*A1*B2-2*A2*B1)
 b=(B3*A1-A3*B1)/(2*A1*B2-2*A2*B1)
 r=sqrt(a^2-2*(sum(x1)/N)*a+b^2-2*(sum(y1)/N)*b+sum(x1.^2)/N+sum(y1.^2)/N)