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.

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...

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.