www.pudn.com > eyedemo.rar > hysthresh.asv

% HYSTHRESH - 后继阈值处理
% 使用: bw = hysthresh(im, T1, T2)
%
% 参数:
%             im  - 待处理的图像
%             T1  - 高阈值
%             T2  - 低阈值
%
% 输出:
%             bw  - the thresholded image (containing values 0 or 1)
%
% 值高于T1的像素被标记为边缘，之后被连接为边缘。高于T2的也被标记为边缘。
% 假设输入图像非负的


function bw = hysthresh(im, T1, T2)

if (T2 > T1 | T2 < 0 | T1 < 0)  % Check thesholds are sensible
  error('T1 must be >= T2 and both must be >= 0 ');
end

[rows, cols] = size(im);    % 预先计算一些值
rc = rows*cols;
rcmr = rc - rows;
rp1 = rows+1;

bw = im(:);                 % 图片存为列向量
pix = find(bw > T1);        % 找值大于 T1的像素点
npix = size(pix,1);         % 所得像素数量

stack = zeros(rows*cols,1); % 创建堆栈数组，使得运算不会溢出

stack(1:npix) = pix;        % 边缘点放入堆栈
stp = npix;                 % 设置堆栈指针
for k = 1:npix
    bw(pix(k)) = -1;        % 标记点为边缘
end


% 预先计算一个数组，任何的值都与其周围八个点的值相关. 注意到图像已经被转换为向量,所以如果将数组重塑为图像，将会是这个样子:
%              n-rows-1   n-1   n+rows-1
%
%               n-rows     n     n+rows
%                     
%              n-rows+1   n+1   n+rows+1

O = [-1, 1, -rows-1, -rows, -rows+1, rows-1, rows, rows+1];

while stp ~= 0            % 堆栈不为空
    v = stack(stp);         % 索引入栈
    stp = stp - 1;
    
    if v > rp1 & v < rcmr   % Prevent us from generating illegal indices
			    % Now look at surrounding pixels to see if they
                            % should be pushed onto the stack to be
                            % processed as well.
       index = O+v;	    % Calculate indices of points around this pixel.	    
       for l = 1:8
	   ind = index(l);
	   if bw(ind) > T2   % if value > T2,
	       stp = stp+1;  % push index onto the stack.
	       stack(stp) = ind;
	       bw(ind) = -1; % mark this as an edge point
	   end
       end
    end
end



bw = (bw == -1);            % Finally zero out anything that was not an edge 
bw = reshape(bw,rows,cols); % and reshape the image