%detect pupil with bwtraceboundary
%Imput: the name of the picture
%Output:center of pupil:[xc,yc]; radius
%Jason Wang
%2005.4.7
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [xc,yc,radius]=detect_pupil(imagename);
imagename=('0020_000.bmp');%test
I_o=imread(imagename);
% imshow(I_o);figure;
I=I_o(100:380,100:540);
TWO=twovalue(I,26);%im2bw(I,0.1);
BW = ~im2bw(TWO);%objects of interest must be white
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%median filter
I_m=medfilt2(BW,[7 7]);
% %close operation,remove eyelash
% se = strel('disk',4);
% I_c = imopen(BW,se);
%fill spot on the pupil
% I_f = imfill(I_c,'holes');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%multi-boundaries
figure;imshow(~I_m);title('multi-boundaries');
hold on;
boundaries = bwboundaries(I_m);
maxlen = 0;
maxk = 0;
for k=1:size(boundaries)
   b = boundaries{k};
   len = length(b);
   if len>maxlen
       maxlen = len;
       maxk=k;
   end
   plot(b(:,2),b(:,1),'g','LineWidth',3);
end
%----------------------------------------%
% %尋找single邊緣
% dim = size(I_f);
% row = round(dim(1)/2);
% col = max(find(I_f(row,:)));
% boundary = bwtraceboundary(I_f,[row, col],'N');%找到瞳孔邊緣
% figure;
% imshow(I);
% title('detect_pupil');
% hold on;
% plot(col,row,'xy');%start point
% plot(boundary(:,2),boundary(:,1),'w','LineWidth',2);%畫出實際邊緣,綠色線
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%fit a circle
[xc,yc,radius] = fit_circle(boundaries{maxk});

%view step results
figure;imshow(~BW);title('二值化圖像');
figure;imshow(~I_m);title('中值濾波圖像');
% figure;imshow(~I_c);title('開操作');
% figure;imshow(~I_f);title('填充圖像');