Segmenting a grayscale image

Question

I am having trouble achieving the correct segmentation of a grayscale image:

The ground truth, i.e. what I would like the segmentation to look like, is this:

I am most interested in the three components within the circle. Thus, as you can see, I would like to segment the top image into three components: two semi-circles, and a rectangle between them.

I have tried various combinations of dilation, erosion, and reconstruction, as well as various clustering algorithms, including k-means, isodata, and mixture of gaussians--all with varying degrees of success.

Any suggestions would be appreciated.

Edit: here is the best result I've been able to obtain. This was obtained using an active contour to segment the circular ROI, and then applying isodata clustering:

There are two problems with this:

• The white halo around the bottom-right cluster, belonging to the top-left cluster
• The gray halo around both the top-right and bottom-left cluster, belonging to the center cluster.

Predilection

It looks like you have a tri-modal histogram. Check out my answer on dsp.stackexchange.com -> dsp.stackexchange.com/questions/3643/…. In any case, your question should be moved to that site. – Matthia

I have tried histogram-based segmentation. The problem with this approach is that the values of pixels surrounding the right-most cluster are most similar to the pixels of the left-most cluster, resulting in a "halo" around the right-most cluster. – Predilection

Have you tried bwboundaries, or bwlabel? They may work for you. But if the halo effect is too much, you might not get the results you want. Still, with some manipulation and cleaning before and after, you might get what you want. It might also help to use something like hough circles, or other methods like this one, so that you know the bounds of your circle. – Leninist

I should have mentioned this originally, but I use an active contour to segment the circular ROI. I've updated my question with a picture of the results. I'm not sure how bwboundaries and bwlabel can help here, since they require a binary image. – Predilection

Answer 1

Here's a starter... use circular Hough transform to find the circular part. For that I initially threshold the image locally.

 im=rgb2gray(imread('Ly7C8.png'));
 imbw = thresholdLocally(im,[2 2]); % thresold localy with a 2x2 window
 % preparing to find the circle
 props = regionprops(imbw,'Area','PixelIdxList','MajorAxisLength','MinorAxisLength');
 [~,indexOfMax] = max([props.Area]);
 approximateRadius =  props(indexOfMax).MajorAxisLength/2;
 radius=round(approximateRadius);%-1:approximateRadius+1);
 %find the circle using Hough trans.
 h = circle_hough(edge(imbw), radius,'same');
 [~,maxIndex] = max(h(:));
 [i,j,k] = ind2sub(size(h), maxIndex);
 center.x = j;     center.y = i;

 figure;imagesc(im);imellipse(gca,[center.x-radius  center.y-radius 2*radius 2*radius]);
 title('Finding the circle using Hough Trans.');

select only what's inside the circle:

 [y,x] = meshgrid(1:size(im,2),1:size(im,1));
 z = (x-j).^2+(y-i).^2;
 f = (z<=radius^2);
 im=im.*uint8(f);

EDIT:

look for a place to start threshold the image to segment it by looking at the histogram, finding it's first local maxima, and iterating from there until 2 separate segments are found, using bwlabel:

  p=hist(im(im>0),1:255);
  p=smooth(p,5);
  [pks,locs] = findpeaks(p);

  bw=bwlabel(im>locs(1));
  i=0;
  while numel(unique(bw))<3
     bw=bwlabel(im>locs(1)+i); 
     i=i+1;
  end


 imagesc(bw);

The middle part can now be obtained by taking out the two labeled parts from the circle, and what is left will be the middle part (+some of the halo)

 bw2=(bw<1.*f);

but after some median filtering we get something more reasonble

 bw2= medfilt2(medfilt2(bw2));

and together we get:

 imagesc(bw+3*bw2); 

The last part is a real "quick and dirty", I'm sure that with the tools you already used you'll get better results...

Allargando

Very nice, thanks for taking the time to come up with this. However, the result depends heavily on the threshold value, in this case 186. For example, a value of 183 significantly distorts the result. I realize I didn't specify this in my question, but this is only a typical example of the problem; there are many other images that require similar segmentation. Therefore, I was hoping to come up with an automated, general solution. Perhaps there is some way to "guess" the threshold parameter. – Predilection

There are several ways to do this, see my edit for a quick a dirty way to automatically set a threshold that is not far from what I had before... – Allargando

Answer 2

One can also obtain an approximate result using the watershed transformation. This is the watershed on the inverted image -> watershed(255-I) Here is an example result:

Another Simple method is to perform a morphological closing on the original image with a disc structuring element (one can perform multiscale closing for granulometries) and then obtain the full circle. After this extracting the circle is and components withing is easier.

se = strel('disk',3);
Iclo = imclose(I, se);% This closes open circular cells.
Ithresh = Iclo>170;% one can locate this threshold automatically by histogram modes (if you know apriori your cell structure.)
Icircle = bwareaopen(Ithresh, 50); %to remove small noise components in the bg

Ithresh2 = I>185; % This again needs a simple histogram.

Ethnocentrism