Does anyone know of a good method of finding "peaks" in 2D images? A peak would be defined as a region of higher gray level than the surrounding area. For example, let's say you had topography of the Rocky Mountains in Colorado. Well there are millions of local peaks. Even little peaks on the sides of much larger peaks. Is there any way to find the, say, three "largest" peaks. What algorithms are there? I'm sure it can get quite complicated as to what exactly defines a peak and what is too small to be a real peak and should just be considered noise on a larger peak. But I'm open to ideas.
I've used imregionalmax() and imextendedmax() in MATLAB and they don't really seem to do what I want, plus I don't really understand the algorithm and they don't explain it. For example, they say: "BW = imextendedmax(I,H) computes the extended-maxima transform, which is the regional maxima of the H-maxima transform. H is a nonnegative scalar. Regional maxima are connected components of pixels with a constant intensity value, and whose external boundary pixels all have a lower value." Yep, clear as mud. Anyway it didn't work because it found small dots in my largest peak and found a larger base on my smallest peaks - very counter-intuitive.
> Does anyone know of a good method of finding "peaks" in 2D images? A > peak would be defined as a region of higher gray level than the > surrounding area. For example, let's say you had topography of the > Rocky Mountains in Colorado. Well there are millions of local peaks. > Even little peaks on the sides of much larger peaks. Is there any way > to find the, say, three "largest" peaks. What algorithms are there? > I'm sure it can get quite complicated as to what exactly defines a > peak and what is too small to be a real peak and should just be > considered noise on a larger peak. But I'm open to ideas.
> I've used imregionalmax() and imextendedmax() in MATLAB and they don't > really seem to do what I want, plus I don't really understand the > algorithm and they don't explain it. For example, they say: > "BW = imextendedmax(I,H) computes the extended-maxima transform, which > is the regional maxima of the H-maxima transform. H is a nonnegative > scalar. > Regional maxima are connected components of pixels with a constant > intensity value, and whose external boundary pixels all have a lower > value." Yep, clear as mud. Anyway it didn't work because it found > small dots in my largest peak and found a larger base on my smallest > peaks - very counter-intuitive.
> Thanks, > ImageAnalyst
Why don't you just mark pixels with luminance >= epsilon, starting with epsilon = 255 and working downwards until you have three segregated regions of whatever size you consider not noise?
> Does anyone know of a good method of finding "peaks" in 2D images? A > peak would be defined as a region of higher gray level than the > surrounding area. For example, let's say you had topography of the > Rocky Mountains in Colorado. Well there are millions of local peaks. > Even little peaks on the sides of much larger peaks. Is there any way > to find the, say, three "largest" peaks. What algorithms are there? > I'm sure it can get quite complicated as to what exactly defines a > peak and what is too small to be a real peak and should just be > considered noise on a larger peak. But I'm open to ideas.
> I've used imregionalmax() and imextendedmax() in MATLAB and they don't > really seem to do what I want, plus I don't really understand the > algorithm and they don't explain it. For example, they say: > "BW = imextendedmax(I,H) computes the extended-maxima transform, which > is the regional maxima of the H-maxima transform. H is a nonnegative > scalar. > Regional maxima are connected components of pixels with a constant > intensity value, and whose external boundary pixels all have a lower > value." Yep, clear as mud. Anyway it didn't work because it found > small dots in my largest peak and found a larger base on my smallest > peaks - very counter-intuitive.
I think that regional maximum just means that it can have several pixels of the same value in the maximum. So largest peaks are probably more sharp and have less pixels having same values. But Matlab is not my tool, so just some simple ideas.
To get just locations of the largest peaks, you should get a list of the locations and heights of the local maxima, sort it according to the height value and select the largest one. Then you define a buffer around the peak, perhaps making the buffer size proportional to the height of the peak and find the next largest maximum that is not within the buffer. And so on.
If you need the mountains as well, then one way is low pass filtering to get rid of the little peaks on the sides, and watershed segmentation of the negative filtered image to get segments around the local minima in the negative image. It may be useful to select the size of the low pass filter according to the height and filter more heavily on the high areas.
function [m,i] = localmax(x, w) % LOCALMAX Find indices and amplitudes of local maxima % [m,i] = localmax(x, w) returns the indices and maxima defined % over a local window of size 2w+1 given by w points on either % side of the point being considered as a local maximum.. % inputs: % x: vector (double) of data to find local max % w: window % outputs: % m: data value at local maxima % i: indices of local maxima % % by P.E.McSharry % These routines are made available under the GNU general public license. % If you have not received a copy of this license, please download from % http://www.gnu.org/ % % Please distribute (and modify) freely, commenting where you have % added modifications. The author would appreciate correspondence % regarding corrections, modifications, improvements etc. % % G. Clifford : g...@ieee.org
N = length(x); w = round(w); % change-- colon operator used below requires integer limits k = 2*w+1; %-- not needed y = zeros(k,1); m = []; i = []; l = 0; for j=w+1:N-w y = x(j-w:j+w); [ymax,imax] = max(y); if imax == w+1 l = l+1; m(l) = ymax; i(l) = j; end end
> On Mar 26, 2:10 pm, ImageAnalyst <imageanal...@mailinator.com> wrote:
> > Does anyone know of a good method of finding "peaks" in 2D images? A > > peak would be defined as a region of higher gray level than the > > surrounding area. For example, let's say you had topography of the > > Rocky Mountains in Colorado. Well there are millions of local peaks. > > Even little peaks on the sides of much larger peaks. Is there any way > > to find the, say, three "largest" peaks. What algorithms are there? > > I'm sure it can get quite complicated as to what exactly defines a > > peak and what is too small to be a real peak and should just be > > considered noise on a larger peak. But I'm open to ideas.
> > I've used imregionalmax() and imextendedmax() in MATLAB and they don't > > really seem to do what I want, plus I don't really understand the > > algorithm and they don't explain it. For example, they say: > > "BW = imextendedmax(I,H) computes the extended-maxima transform, which > > is the regional maxima of the H-maxima transform. H is a nonnegative > > scalar. > > Regional maxima are connected components of pixels with a constant > > intensity value, and whose external boundary pixels all have a lower > > value." Yep, clear as mud. Anyway it didn't work because it found > > small dots in my largest peak and found a larger base on my smallest > > peaks - very counter-intuitive.
> > Thanks, > > ImageAnalyst
> Why don't you just mark pixels with luminance >= epsilon, starting > with epsilon = 255 and working downwards until you have three > segregated regions of whatever size you consider not noise?
> Also, a key word is "hillocks".- Hide quoted text -
> - Show quoted text -
--------------------------------------------------------------------------- -------------------------- aruzinsky: It's more difficult than that. Let's take a 1D example, let's say your signal was this: [ 0 255 0 240 0 255 255 0 100 100 100 100 100 100 0] The method you suggested would find the first 3 peaks first, but the third peak would have less volume (510) than the 4th peak (600). So you can't quit when you have your required number. Now let's take another case: [0 255 30 100 100 220 100 200 100 200 100 0] Going down from 255, you'd first have 1 peak, then you'd have 2, then you'd have 4, and then (at less than 100) you'd be back down to only 2 peaks again, and finally (less than 30) you'd have only 1 peak. So it can get pretty tricky, particularly when you have peaks on top of peaks. Now I know there are lots of 1D peak detector algorithms developed by the mass-spectrometry and chromatography community, but I'm looking for a 2D one that makes sense.
I'm currently looking into morphological methods such as: L. Vincent. (1993). Morphological Grayscale Reconstruction in Image Analysis: Applications and Efficient Algorithms. IEEE Transactions on Image Processing, vol.2 (2) 176-201 and http://www.afscet.asso.fr/resSystemica/Crete02/Halkiotis.pdf It seems like they have potential but are a bit tricky to understand. Doing a watershed and then calculating volume in each zone (region) may also work as long as I can get a decent separation of peaks and not have it over-segmented.
> > Does anyone know of a good method of finding "peaks" in 2D images? A > > peak would be defined as a region of higher gray level than the > > surrounding area. For example, let's say you had topography of the > > Rocky Mountains in Colorado. Well there are millions of local peaks. > > Even little peaks on the sides of much larger peaks. Is there any way > > to find the, say, three "largest" peaks. What algorithms are there? > > I'm sure it can get quite complicated as to what exactly defines a > > peak and what is too small to be a real peak and should just be > > considered noise on a larger peak. But I'm open to ideas.
> > I've used imregionalmax() and imextendedmax() in MATLAB and they don't > > really seem to do what I want, plus I don't really understand the > > algorithm and they don't explain it. For example, they say: > > "BW = imextendedmax(I,H) computes the extended-maxima transform, which > > is the regional maxima of the H-maxima transform. H is a nonnegative > > scalar. > > Regional maxima are connected components of pixels with a constant > > intensity value, and whose external boundary pixels all have a lower > > value." Yep, clear as mud. Anyway it didn't work because it found > > small dots in my largest peak and found a larger base on my smallest > > peaks - very counter-intuitive.
> I think that regional maximum just means that it can have several pixels > of the same value in the maximum. So largest peaks are probably more > sharp and have less pixels having same values. But Matlab is not my > tool, so just some simple ideas.
> To get just locations of the largest peaks, you should get a list of the > locations and heights of the local maxima, sort it according to the > height value and select the largest one. Then you define a buffer around > the peak, perhaps making the buffer size proportional to the height of > the peak and find the next largest maximum that is not within the > buffer. And so on.
> If you need the mountains as well, then one way is low pass filtering to > get rid of the little peaks on the sides, and watershed segmentation of > the negative filtered image to get segments around the local minima in > the negative image. It may be useful to select the size of the low pass > filter according to the height and filter more heavily on the high areas.
> Juho- Hide quoted text -
> - Show quoted text -
------------------------------------------------------------ Juho: Thanks for your reply. I've seen some papers that have had some success via the watershed method and I might try that. ImageAnalyst
> function [m,i] = localmax(x, w) > % LOCALMAX Find indices and amplitudes of local maxima > % [m,i] = localmax(x, w) returns the indices and maxima defined > % over a local window of size 2w+1 given by w points on either > % side of the point being considered as a local maximum.. > % inputs: > % x: vector (double) of data to find local max > % w: window > % outputs: > % m: data value at local maxima > % i: indices of local maxima > % > % by P.E.McSharry > % These routines are made available under the GNU general public > license. > % If you have not received a copy of this license, please download > from > %http://www.gnu.org/ > % > % Please distribute (and modify) freely, commenting where you have > % added modifications. The author would appreciate correspondence > % regarding corrections, modifications, improvements etc. > % > % G. Clifford : g...@ieee.org
> N = length(x); > w = round(w); % change-- colon operator used below requires integer > limits > k = 2*w+1; > %-- not needed y = zeros(k,1); > m = []; > i = []; > l = 0; > for j=w+1:N-w > y = x(j-w:j+w); > [ymax,imax] = max(y); > if imax == w+1 > l = l+1; > m(l) = ymax; > i(l) = j; > end > end
--------------------------------------------------------------------------- ------------------ Bob: Thanks for the link. I did try that code but for noisy data, it did a terrible job. Essentially finding every other point in the whole 1D array - since about half the values in a noisy array are local max and half are local mins. That algorithm didn't do as good a job as the one I'm currently using for 1D peak detection: http://billauer.co.il/peakdet.html There are lots and lots of 1D peak detection algorithm: http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=2631518 reviews several of them. It appears most of them come from the chemistry community where they have a variety of instruments that generate 1D signals (mass sepctrometry, etc.)
2D peak detection is a lot more complicated, particularly if you don't just want to find the location of the peak's very tip top but want to find the region that its base covers. I might take a look deeper into watershed, such as http://cmm.ensmp.fr/~beucher/wtshed.html or else morphological methods like promoted by Dr. Luc Vincent (currently at Google). L. Vincent. (1993). Morphological Grayscale Reconstruction in Image Analysis: Applications and Efficient Algorithms. IEEE Transactions on Image Processing, vol.2 (2) 176-201
> On Mar 26, 5:01 pm, aruzinsky <aruzin...@general-cathexis.com> wrote:
> > On Mar 26, 2:10 pm, ImageAnalyst <imageanal...@mailinator.com> wrote:
> > > Does anyone know of a good method of finding "peaks" in 2D images? A > > > peak would be defined as a region of higher gray level than the > > > surrounding area. For example, let's say you had topography of the > > > Rocky Mountains in Colorado. Well there are millions of local peaks. > > > Even little peaks on the sides of much larger peaks. Is there any way > > > to find the, say, three "largest" peaks. What algorithms are there? > > > I'm sure it can get quite complicated as to what exactly defines a > > > peak and what is too small to be a real peak and should just be > > > considered noise on a larger peak. But I'm open to ideas.
> > > I've used imregionalmax() and imextendedmax() in MATLAB and they don't > > > really seem to do what I want, plus I don't really understand the > > > algorithm and they don't explain it. For example, they say: > > > "BW = imextendedmax(I,H) computes the extended-maxima transform, which > > > is the regional maxima of the H-maxima transform. H is a nonnegative > > > scalar. > > > Regional maxima are connected components of pixels with a constant > > > intensity value, and whose external boundary pixels all have a lower > > > value." Yep, clear as mud. Anyway it didn't work because it found > > > small dots in my largest peak and found a larger base on my smallest > > > peaks - very counter-intuitive.
> > > Thanks, > > > ImageAnalyst
> > Why don't you just mark pixels with luminance >= epsilon, starting > > with epsilon = 255 and working downwards until you have three > > segregated regions of whatever size you consider not noise?
> > Also, a key word is "hillocks".- Hide quoted text -
> > - Show quoted text -
> --------------------------------------------------------------------------- -------------------------- > aruzinsky: > It's more difficult than that. Let's take a 1D example, let's say > your signal was this: [ 0 255 0 240 0 255 255 0 100 100 100 100 100 > 100 0] > The method you suggested would find the first 3 peaks first, but the > third peak would have less volume (510) than the 4th peak (600). So > you can't quit when you have your required number. Now let's take > another case: [0 255 30 100 100 220 100 200 100 200 100 0] > Going down from 255, you'd first have 1 peak, then you'd have 2, then > you'd have 4, and then (at less than 100) you'd be back down to only 2 > peaks again, and finally (less than 30) you'd have only 1 peak. So it > can get pretty tricky, particularly when you have peaks on top of > peaks. Now I know there are lots of 1D peak detector algorithms > developed by the mass-spectrometry and chromatography community, but > I'm looking for a 2D one that makes sense.
> I'm currently looking into morphological methods such as: > L. Vincent. (1993). Morphological Grayscale Reconstruction in Image > Analysis: Applications and Efficient Algorithms. IEEE Transactions on > Image Processing, vol.2 (2) 176-201 > andhttp://www.afscet.asso.fr/resSystemica/Crete02/Halkiotis.pdf > It seems like they have potential but are a bit tricky to understand. > Doing a watershed and then calculating volume in each zone (region) > may also work as long as I can get a decent separation of peaks and > not have it over-segmented.
> Thanks, > ImageAnalyst- Hide quoted text -
> - Show quoted text -
Regarding, "Going down from 255, you'd first have 1 peak, then you'd have 2, then you'd have 4, and then (at less than 100) you'd be back down to only 2 peaks again, and finally (less than 30) you'd have only 1 peak.", you didn't put a minimum size restriction on the regions as I suggested. Putting a size restriction >=2, and working down, the first peak of
[0 255 30 100 100 220 100 200 100 200 100 0]
is
[100 100 220 100 200 100 200 100]
and the second peak is
[255 30].
Also, noise can be reduced by general methods and then region size disregarded. Suppose we apply a median filter (which would be appropriate for shot noise) of length three to your examples, reflecting at borders, and disregard region size.