How can local maxima be pinpointed on a surf plot?

Hello everyone. I posted a question a few days ago and I repost it now since I got no answers (maybe because Ι wasn't clear about what I needed). So I have a surf plot (I include an image of it) and I need to pinpoint the local maxima on it. The range of the local maxima that will be pinpointed I want to be controlled by me. I also include my data if that helps. Thanks in advance!

[Paschalis Garouniatis]

In reference to the range of the local maxima. What I mean is that I need to pinpoint local maxima points from a certain value and higher. Hope that clarifies the meaning of range.

Answer 1

[KSSV]

You can get your maximum value using max along with indices. To get your range values use find.

Answer 2

[Walter Roberson]

Calculate a mask image, Data>threshold. regionprops asking for PixelIdxList. For each region returned use the indices to index the data and take max of that subset of the data with the two output max to find the index of the max relative to to that subset; use that to index the the PixelIdxList to find the linear index of the local max. You can then convert linear to 2d coordinate. Caution: the local max might be duplicated which your problem description did not account for.

[Paschalis Garouniatis]

How can I calculate a mask image?

[Walter Roberson]

Mask=TheArray>TheThreshold;

You can then regionprops the results

[Paschalis Garouniatis]

Thank you for your answer. I will give it a shot and see if I can make it.

[Paschalis Garouniatis]

Using the following code I got a struct(292x1). How do I proceed from here? I read something about a step function but got me nowhere.

    Mask=Zq>0.2
    H=regionprops(Mask,'pixelidxlist')

[Walter Roberson]

for K = 1 : length(H)
  this_region_idx = H(K).PixelIDxList;
  this_region_vals = Zq(this_region_idx);
  [~, relative_max_idx] = max(this_region_vals);
  this_region_idx_of_max = this_region_idx(relative_max_idx);
  [r, c] = ind2sub(Zq, this_region_idx_of_max);
  text(c, r, '+', 'color', 'red');   %notice r is y and c is x
end

[Image Analyst]

If he needs to find "local maxima points from a certain value and higher." he can simply threshold to find the values above a threshold, then AND that with imregionalmax() to find local max that are within the above-threshold regions only:

localMaxAboveThreshold = (grayImage > threshold) & imregionalmax(grayImage);

[Walter Roberson]

A couple of questions have been left open so far:

• what do you want to do if there are multiple locations in the region that share the same maximum?
• what do you want to do if there are multiple peaks within the same region, such as if you had

     A
    * *   B
    * *  * *
    * **** ***
   *          *
   ------------  threshold

where A and B are separated and not necessarily equal, but are still both relatively high peaks? Especially when signals are noisy this is quite common to encounter. It also can happen with narrow gaussians with high phases: depending exactly where the gaussian is sampled relative to its phased peak, the data can peak, drop a lot, then peak again.

[Paschalis Garouniatis]

For the first question the answer is that if there are multiple locations with the same maximum then I won't ignore any of them. For the second question the answer is the same. In your previous comment the command

 text(c, r, '+', 'color', 'red');

does not have impact on the surf plot. Should I find the corresponding Xq and Yq of the variable this_region_vals (if I am not mistaken it contains the local maxima over the given threshold)and then use scatter3 over the surf plot?

[Image Analyst]

Did you see my comment? Just use find() on it to get the rows and columns

localMaxAboveThreshold = (grayImage > threshold) & imregionalmax(grayImage);
[r, c] = find(localMaxAboveThreshold);

[Walter Roberson]

Hmmm, yes, you would need to convert row and column to x y z coordinates.

[Paschalis Garouniatis]

Below I am posting the code I wrote. I tried to run it and seems to need some improvement on the plotting section. What am I missing? My aim is to place on the surf plot markers for the local maxima.

 % Create data for surf plot
X=-30.0:0.06006:30.0
Y=20:50:420
Z=[perc1{1,1}';perc1{2,1}';perc1{3,1}';perc1{4,1}';perc1{5,1}';perc1{6,1}';perc1{7,1}';perc1{8,1}';perc1{9,1}']
Xq=-30.0:0.006006:30.0
Yq=20:5:420
[Xq,Yq]=meshgrid(Xq,Yq)
Zq=interp2(X,Y,Z,Xq,Yq)

 % Find indices of local maxima
data_plottlr=[{Xq},{Yq},{Zq}]
grayImage=data_plottlr{1,3}
Threshold=0.2;
localMaxAboveThreshold=(grayImage>Threshold)
idx={}
for i=1:size(localMaxAboveThreshold,1)
    count=1;
    for j=1:size(localMaxAboveThreshold,2)
        if true(localMaxAboveThreshold(i,j))
            idx{i,1}(count,:)=[data_plottlr{1,1}(i,j), data_plottlr{1,2}(i,j), data_plottlr{1,3}(i,j)];
            count=count+1;
        else
            continue;
        end
    end
end

 % Surf plot with demonstration of local maxima
figure;
hold on
surf(Xq,Yq,Zq)
for k=1:size(idx,1)
    for l=1:size(idx{k,1},1)
        scatter3(idx{k,1}(l,1),idx{k,1}(l,2),idx{k,1}(l,3))
    end
end
shading interp
hold off

[Walter Roberson]

Move

      count=1;

to before the "for" loop.

Replace

if true(localMaxAboveThreshold(i,j))

with

if localMaxAboveThreshold(i,j)

However, that entire section can be replaced:

   V = @(M) M(:);    %helper

idx = [V(data_plottlr{1,1}(localMaxAboveThreshold)), V(data_plottlr{1,2}(localMaxAboveThreshold)), V(data_plottlr{1,3}(localMaxAboveThreshold))];    %no loops

However your code gets all of the points above the threshold, not just the maxima.

[Image Analyst]

He must not have the Image Processing Toolbox and so can't use imregionalmax(), which would make it trivial.

[Walter Roberson]

"For the first question the answer is that if there are multiple locations with the same maximum then I won't ignore any of them."

How would you like to choose, then? Do you want to select any point which is greater than all of its immediate neighbours? If so then Image Analyst's suggestion of imregionalmax would probably help.

The difference between what I was suggesting earlier and imregionalmax: if you slightly modify the diagram I posted earlier, to make

    A
    * *   B
    * *  * * C
    * **** ***
   *          *         D
   -----------------------  threshold

then imregionalmax would include C and D because they are greater than all of their surrounding neighbours, even though they are obviously not the peak of the entire region, and might well be just a small bit of noise above their immediate neighbours. (If you choose continuous points randomly independently from each other, in a 2D array, then on average 1/9 of them would be considered regional maxima.)

[Paschalis Garouniatis]

Sorry for my absence but I had a lot of work to pull off. I tried what ImageAnalyst suggested (imregionalmax) and it worked. Didn't answer earlier because I was trying some codes. Walter I made a mistake about what I need to do about the maxima.

            D                              L
       C*********E           H        K********M
     B****************F   G******I  J************N
 A**************************************************O
 -------------------------------------------------------- threshold

This schematic is an example of my data. With imregionalmax or your suggestion all peaks are included. I need to keep only D,H,L. How can this be done? Again sorry for my mistake and the absence.

[Image Analyst]

What suggestion of mine did you say worked? Using imregionalmax() or using an arrow to point to some points?

[Paschalis Garouniatis]

I was refering to Imregionalmax().

[Walter Roberson]

Sorry, I missed seeing the update here. Could you include a bit of data with this kind of structure that imregionalmax is not working for? Looking at the diagram it looks to me as if imregionalmax would work.

Answer 3

[Image Analyst]

I don't know what "pinpoint" means. Perhaps you want a 3-D arrow with the tip of it pointing to the max. If so, there are 3-D arrow functions in the File Exchange. Please attach a screenshot that's an image file, like PNG, NOT a .fig file which requires me to right click, save it to disk, switch to MATLAB, browse to the file I saved, and open it. I did not do any of that. With a PNG file you can jsut insert it here and save me a lot of work.

By the way, local max can be found with the function imregionalmax().

[Paschalis Garouniatis]

Thank you for your response Image Analyst. I want to show on the surf plot the local maxima with a symbol like "+" or something like that and their values. That's what I meant with the word "pinpoint". Hope this helps.

[Walter Roberson]

The code I showed will do that '+' .

[Image Analyst]

I'd probably use plot3() to put a marker at the top peak location. Let's say you somehow found the x,y,z location(s) of the peaks, then say

plot3(x, y, z, 'r.', 'MarkerSize', 30);

Of course the arrow in 3D is also a good alternative way.