find_extrema.m
上传用户:trade789
上传日期:2018-05-10
资源大小:603k
文件大小:3k
- %/////////////////////////////////////////////////////////////////////////////////////////////
- %
- % find_extrema - finds local maxima within a grayscale image. Each point is
- % checked against all of the pixels within a given radius to be a local max/min.
- % The magnitude of pixel values must be above the given threshold to be picked
- % as a valid maxima or minima.
- %
- % Usage: m = find_extrema(img,thresh,radius,min_separation)
- %
- % Parameters:
- % img : image matrix
- % thresh : threshold value
- % radius : pixel radius
- %
- % Returns:
- %
- % m - an nX2 matrix of row,column coordinates of selected points
- %
- % Author:
- % Scott Ettinger
- % scott.m.ettinger@intel.com
- %
- % May 2002
- %/////////////////////////////////////////////////////////////////////////////////////////////
- function [mx] = find_extrema(img,thresh,radius)
- %img = abs(img);
- [h,w] = size(img);
- % get interior image subtracting radius pixels from border
- p = img(radius+1:h-radius, radius+1:w-radius);
- %get pixels above threshold
- [yp,xp] = find(p>thresh);
- yp=yp+radius; xp=xp+radius;
- pts = yp+(xp-1)*h;
- %create offset list for immediate neighborhood
- z=ones(3,3);
- z(2,2)=0;
- [y,x]=find(z);
- y=y-2; x=x-2;
- if size(pts,2)>size(pts,1)
- pts = pts';
- end
- %test for max within immediate neighborhood
- if size(pts,1)>0
- maxima=ones(length(pts),1);
- for i=1:length(x)
- pts2 = pts + y(i) + x(i)*h;
- maxima = maxima & img(pts)>img(pts2);
- end
-
- xp = xp(find(maxima)); %save maxima
- yp = yp(find(maxima));
- pts = yp+(xp-1)*h; %create new index list of good points
- end
-
- %create offset list for radius of pixels
- [y,x]=meshgrid(-20:20,-20:20);
- z = (x.^2+y.^2)<=radius^2 & (x.^2+y.^2)>(1.5)^2; %include points within radius without immediate neighborhood
- [y,x]=find(z);
- x=x-21; y=y-21;
- %create offset list for ring of pixels
- [y2,x2]=meshgrid(-20:20,-20:20);
- z = (x2.^2+y2.^2)<=(radius)^2 & (x2.^2+y2.^2)>(radius-1)^2;
- [y2,x2]=find(z);
- x2=x2-21; y2=y2-21;
- maxima = ones(length(pts),1);
- %test within radius of pixels (done after first test for slight speed increase)
- if size(pts,1)>0
- for i = 1:length(x)
- pts2 = pts + y(i) + x(i)*h;
- maxima = maxima & img(pts)>img(pts2); %test points
- end
- xp = xp(find(maxima)); %save maxima from immediate neighborhood
- yp = yp(find(maxima));
- pts = yp+(xp-1)*h; %create new index list
- mx = [yp xp];
-
- else
- mx = [];
- end
-