基于Hough变换的人眼虹膜定位

Github个人博客：https://joeyos.github.io

clear all; 
close all;
i=imread('yanjing.bmp'); 
imshow(i); 
iii=i; 
%把输入图象二值化，用canny算法返回阈值
sigma=3.0;
thresh=[0.03,0.09];
bw_1=i>70;
edgerm=edge(bw_1,'canny',thresh,sigma); 
figure,imshow(edgerm);
t1=280;
s=0;
while t1>10
t2=1;
while t2<310
%查找第一个边缘点
if edgerm(t1,t2)==1 
         u1=t1;
         u2=t2;
         s=1;
end
if s==1
   break;
end
  t2=t2+1;  
end
t1=t1-1;
end
po=1;
sum2=0;
%第一个边缘点
o1=u1; 
o2=u2;
hang=zeros(0,0);
lie=zeros(0,0);
while (po==1)
   while (po==1)
         sum1=0;
         for t3=1:5
            for t4=1:5
               % 第一个边缘点的左上方5个像素内有边缘点
               if edgerm(u1-t3+1,u2+t4-1)==1                  
% 第一个边缘点周围的边缘点个数
sum1=sum1+1; 
                  sum2=sum2+1;
                  % 第sum1个边缘点位置x
hang(sum1,1)=u1-t3+1;
                  % 第sum1个边缘点位置y
hang(sum1,2)=u2+t4-1;
                  lie(sum2,1)=u1-t3+1;
                  lie(sum2,2)=u2+t4-1;
               end
            end
         end
         % 边缘点只有一个
if sum1==1 
            po=0;
         % 没有边缘点
elseif sum1==0 
            po=0;
         else
            % 以最后的边缘点为起点，进行下一轮搜索
u1=hang(sum1,1); 
            u2=hang(sum1,2);
            po=1;
         end
      end
      % 边缘点个数小于30个
if sum2<30 
         u1=o1;
         u2=o2+1;
         po=1;
         sum2=0;
      % 横坐标不变，改变纵坐标值得到边缘点
while (edgerm(u1,u2)~=1)   
         while (edgerm(u1,u2)~=1)&(u2<310)
            % 不是边缘点，纵坐标加1
u2=u2+1; 
         end 
         % 没有得到边缘点
if u2==310 
            u1=u1-1;
            u2=1;
         end
      end
      % x不变，改变y重新得到边缘点
o1=u1; 
      o2=u2;
      else
         break;
      end  
   end 
% 边缘点个数
a1=size(lie); 
w1=lie(a1(1),1);
w2=lie(a1(1),2);
po1=1;
      while (po1==1)
         sum1=0;
         for t1=1:3
            for t2=1:5
               % 边缘点向左方3个像素，上方5个像素
if edgerm(w1-t1+1,w2-t2+1)==1 
                  sum1=sum1+1;
                  sum2=sum2+1;
                  lie(sum2,1)=w1-t1+1;
                  lie(sum2,2)=w2-t2+1;
                  hang(sum1,1)=w1-t1+1;
                  hang(sum1,2)=w2-t2+1;
               end
            end
         end   
         % 边缘点只有一个
if sum1==1 
            po1=0;
         else
            po1=1;
            w1=hang(sum1,1);
            w2=hang(sum1,2);
         end
      end

 po2=1;
 while (po2==1)
         sum1=0;
         for t1=1:7
            for t2=1:15
               if edgerm(w1+t1-1,w2-t2+1)==1 
                  sum1=sum1+1;