RegionProps Major / Minor axis determination Problem

Hi, I am getting a strange result using regionprops. I have an image that I have the locations of the centroids (x1,y2) and my aim is to determine the elipticity somehow. I am using the major and minor length as my first approach, via regionprops.

So I create a binary image knowing the centroids and dilate

%Create Binary image from the already found centroid locations (x1,y1)
      binary = false(size(IM));
      binary(sub2ind(size(IM), round(y1), round(x1))) = true;

      %Dilate slighlty
      n=8;
      se = strel('disk',n);    %2 = default
      binary=imdilate(binary,se,'same');
      subplot(1,3,2);  imshow(binary); title('Binary Image')
      label = bwlabel(binary); %// label each object
      s = regionprops(label, double(IM),'WeightedCentroid', 'Area','Orientation', 'MajorAxisLength', ...
      'MinorAxisLength', 'Eccentricity','MeanIntensity');
      centroids = cat(1, s.WeightedCentroid);
      x1=centroids(:,1); y1=centroids(:,2);    %Use these weighted centroids now rather the first guess x1,y1.
      hold on;   plot(x1,y1,'r.')

      %Now get major and minor lengths.
      major=cat(1,s.MajorAxisLength)
      minor=cat(1,s.MinorAxisLength)

     for k = 1:length(s)
      a(k,1) = s(k).MajorAxisLength/2;
      b(k,1) = s(k).MinorAxisLength/2;
    end

    [a(:) b(:)]   %Display in column form

But this is where something is going wrong - they are all the same value, here are the values reported:

ans =

      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905
      7.6905    7.6905

Answer 1

[Steven Lord]

To my eyes, your binary image and your raw image don't look the same. The objects in the binary image look circular, while the objects in your raw image do look more elliptical. That's because of the way you created the binary image. You set one pixel in the binary image for each object, then you dilated with a disk.

For each of the points in the binary image, when you dilate you move the center of the disk around the "perimeter" of the point. What you receive from that operation is the circular disk itself. A circle is an ellipse with both focal points in the same place and in that case the semi-major axis and the semi-minor axis are both radii of the circle. Since they're both radii of the circle they must have the same length.

So the results you receive seem correct for the image on which you're computing, but the image on which you're computing is not representative of your real data.

[Jason]

Hi, so I thought the binary image would act as a kind of mask which regionprops would then operate on the real image that isn't masked out.

[Image Analyst]

No.

[Steven Lord]

So in the documentation for regionprops you're using the "stats = regionprops(_,I,properties)" syntax? Look at the matrix label that you pass into regionprops as the first input. I suspect that the labeled areas as exactly those circular regions you created by dilation.

Ah, I think I understand where you may be confused. You were expecting MATLAB to look at regions in the original image that contain the labeled sections plus any elements connected to them, so the circular disks would identify but not be the full extent of the regions on which regionprops would act?

I don't believe that's how it works. I believe if an element is not in any labeled region in that first input, regionprops ignores it even if it's adjacent to an element that is in a labeled region. Close doesn't count.

I think what you want instead is the series of steps in this example from the Image Processing Toolbox documentation.

Answer 2

[Image Analyst]

You already had a binary image with elliptical blobs - that's how you got the centroids, right? So just use that. In fact, you can get the centroids, weighted centroids, and eccentricity all in the same call:

props = regionprops(labeledImage, 'Centroid', 'WeightedCentroid', 'Eccentricity');

No need to do that other stuff you did with dilation. Not only is it wrong, but it's unnecessary.

[Jason]

You already had a binary image with elliptical blobs - that's how you got the centroids, right?

No, I used normakused cross correlation to give me the first pass at the centroids. I then created a binary from this

[Image Analyst]

Why? It looks like a simple global threshold would work. Even with your way, didn't your initial binary image have elliptical blobs?

[Jason]

OK, so I have realised that the binary image itself has to reflect the shape of my real object. I seemed to have this bit worked out I think and have attempted to get the average major and minor axis and plot them. Im not convinced they are correct, they look too long and Im not sure about the angle.

%Create Binary image from these lcoations
    binary = false(size(IM));
    binary(sub2ind(size(IM), round(y1), round(x1))) = true;

      %Dilate slighlty    
      se = strel('disk',20);    %2 = default
      binary=imdilate(binary,se,'same');

      IM(~binary)=0;   %Put all values in IM to zero that don't pass thru binary
      level=mean(IM(:))+1.5*std2(double(IM));  
      binaryImage = IM >level;   
      labeledImage = bwlabel(binaryImage);
      area_measurements = regionprops(labeledImage,'Area');
      allAreas = [area_measurements.Area];
      %ax area, areas below are discarded
      AllowedBlobIndex1 = find(allAreas >= min);
  %    AllowedBlobIndex=AllowedBlobIndex1(allAreas<max);  THIS DOESNT WORK
  %    FOR SOME REASON

      keeperBlobsImage = ismember(labeledImage, AllowedBlobIndex1);
      s = regionprops(keeperBlobsImage,double(IM),'WeightedCentroid', 'Area','MajorAxisLength','MinorAxisLength','Orientation','Eccentricity')
      centroids = cat(1, s.WeightedCentroid);
      x=centroids(:,1); y=centroids(:,2);
      boundaries = bwboundaries(keeperBlobsImage);

    for i=1:length(boundaries)
     blobBoundary = boundaries{i};
     plot(blobBoundary(:,2), blobBoundary(:,1), 'g-', 'LineWidth', 1);
    end  
    hold off; 

    major=cat(1,s.MajorAxisLength)
    minor=cat(1,s.MinorAxisLength)
    orientation=cat(1,s.Orientation)
    Eccen=cat(1,s.Eccentricity)

    dMajor=mean(major)
    dMinor=mean(minor)
    Theta=mean(orientation)
    Ec=mean(Eccen)

    [ysize,xsize]=size(IM);
    xc=round(xsize/2)
    yc=round(ysize/2)

    xMinor=xc + [-1 1]*dMinor*sind(Theta);
    yMinor=yc + [-1 1]*dMinor*cosd(Theta);
    xMajor=xc + [-1 1]*dMajor*-cosd(Theta);
    yMajor=yc + [-1 1]*dMajor*sind(Theta);

    subplot(1,3,3)
    hold on;
    line(xMinor,yMinor,'Color','r','LineWidth',2);
    line(xMajor,yMajor,'Color','r','LineWidth',2);
    hold off

[Image Analyst]

Attach your original image by itself if you want me to try anything.

[Jason]

This is my Image data.

Also to explain a bit more about my method.

My raw image (middle one) isn't great and has some high intensity values off the grid. So using normalised cross correlation (or "template matching" via

c=normxcorr2(template,Image);

And then thresholding this (and accounting for padding) does identify quite well the centroids. I then create a binary image with holes at these locations and dilated so all the off grid pixels will be set to zero. I find this then helps with creating another binary for regionprops via a global threshold. I just plot the major and minor axis on the RHS image at the centre - but it still looks wrong

[Image Analyst]

Well attach whatever image you'd like to start with, like the middle one, or the one you passed in to normxcorr2().

[Jason]

Hi, my images are tiffs.I did attach an ascii version. Is this what you wanted?

Here is the tiff zipped up

[Image Analyst]

I don't want an ascii version. I'd have to write code to read it and turn it into an image. A tiff or PNG file is preferred. Or even jpg if it's not too compressed/corrupted.

[Jason]

OK understood, did the zipped up tiff work?