Matlab algorithm for finding ROI

Would anyone be able to help me code my algorithm:

WHILE entire image has not been examined by 16 16 window
MOVE 16 16 window to next position
RECORD position and grey level value of pixel with
largest grey level in window

IF pixels surrounding the largest pixel are as bright as the
largest pixel grey level value
AND outer pixels are darker than the largest pixel grey
level value
THEN largest pixel position is the center pixel of a microcalcification
area
END IF
END WHILE

I have done the first dour lines of the algorithm as follows but I am confused about the coding for the rest of it.

N = 16;
  info = repmat(struct, ceil(size(Z, 1) / N), ceil(size(Z, 2) / N)); %creates an array of structs consisting of an m-by-n tiling 
  for row = 1:N:size(Z, 1)%loop through each pixel in the image matrix
      for col = 1:N:size(Z, 2)
        r = (row - 1) / N + 1;% store the row index for each window
        c = (col - 1) / N + 1;%store the col index for each window 
          imgWindow = Z(row:min(end,row+N-1), col:min(end,col+N-1));
          largest = max(imgWindow(:));
          [rLarg, cLarg] = find(imgWindow == largest, 1, 'first');
          window(r, c).largest = largest;
          window(r, c).row = rLarg + row - 1;
          window(r, c).col = cLarg + col - 1;
      end
  end

Any help or guidance would really be appreciated! Thanks!

[Image Analyst]

Define "next position" - is that one pixel over, or 16 pixels over? It matters because it says whether you need to use conv2() or blockproc().

Anyway, let's say that the brightest pixels is at row 16, column 4 - at the bottom row of the window. So are the "surrounding pixels" the (up to) 3 by 3 pixels immediately adjacent to the brightest one but still in the window? So in my example, at (15,3), (15,4), (15,5), (16,3) and (16,5)? Those are the 5 pixels surrounding the pixel at (16,4). Now, what are the "outer pixels"? Are those all the pixels at the perimeter of the 16x16 window, EXCEPT for the pixels which may be classified at being neighbors of the brightest pixel?

As it is, this algorithm is not stated accurately enough to determine how to proceed.

[Joseph Cheng]

I was thinking about this a bit more. Wouldn't a centroiding/center of mass calculation on a binary [0 max] image give you center pixel? Would have to do some sub 16x16 ROI if there are more than 1 blob. This would also be pretty useful without having to deal with boundaries conditions.

[Image Analyst]

Yes, but who said it was a binary image ? It's a gray level image, with presumably more than 2 gray levels.

What I find unusual about this algorithm is that it slides a window along like normal, but it doesn't necessarily operate on the center pixel (which it actually doesn't even have because it's not an odd number of pixels wide) but it could be operating on a pixel anywhere inside the window. It's not saying that the center pixel needs to be the brightest. As the original poster states the algorithm, if the pixel in the lower left is the brightest pixel and it's surrounded by darker pixels, then the center pixel is a microcalcification , which just seems wrong to me. What if the center pixel is really dark??? I don't think the algorithm is being relayed accurately.

[Joseph Cheng]

well from the pseudocode he/she is looking for just the max/"brightest" in the 16x16 block window. If "we" take the current 16x16 and am looking for only the max values and have blobs then it should work. However i do see the case where the blob is in the shape of a doughnut.

I haven't dealt with mammogram images for a while but given it's an X ray the center of the calcification deposit probably has the highest density in the center. In the original post by the user you did suggest to dialate and fill right?

[Image Analyst]

The code creates an array of structures. For each pixel, there is a structure that holds the max gray level and the location of that pixel relative to the center. So the code does not do what the algorithm spells out and ends up with an array of structures that I don't know what to do with. So now we have two mysteries (1) was the algorithm relayed correctly, and (2) what does his algorithm do, since it doesn't do the algorithm and ends up with something that is not used and we don't know what to do with.

Answer 1

[Joseph Cheng]

Here is an example of what you can do to determine the surroundings using convolution.

x= randi(5,16,16); %generate random 16x16 image
x(4:8,7:9)=20;     %make portion of it bright
X=[x(:,1) x x(:,end)];   %Pad column.
PaddedX=[X(1,:);X;X(end,:)]; %pad it again for rows

Nmask = [0 1 0;1 1 1;0 1 0]; % Convolution mask for neighbors.  

Check = conv2(PaddedX,Nmask/5,'valid'); %here conv will place the mask over the now 18x18 image and element by element multiply this against the 3x3 it is over and add the result.  in the Check.    

[r c]= find(Check == max(x(:)))  %find the row and column index for where check equals the max value of the 16x16 block.

As you can see the convolution mask is a + sign which when element by element multiplied and summed and divide by 5. you are basically checking to see if the average of the pixels surrounding the center point is. I'll leave it to you to play around with this to see how you can use a similar approach to check the pixel surroundings