5

I have a image(png format) in hand. The lines that bound the ellipses (represent the nucleus) are over straight which are impractical. How could i extract the lines from the image and make them bent, and with the precondition that they still enclose the nucleus.

The following is the image:
enter image description here

After bending

enter image description here

EDIT: How can i translate the Dilation And Filter part in answer2 into Matlab language? I can't figure it out.

4
  • Are you talking about the black lines? Nov 1, 2011 at 8:04
  • And what do you mean by "bent", in this case. Can you perhaps choose one line in the image and draw how it should look after "bending" it? Nov 1, 2011 at 8:06
  • 3
    @Ivy Depends on how many of these you need, but honestly the best phantom for this type of thing would seem to be a thresholded/simplified image derived from a real image like your example.
    – John Colby
    Nov 1, 2011 at 15:51
  • @John I think it can be done simulating some kind of random growing from the core ellipses, but I still can not get it wrking. Nov 2, 2011 at 3:46

3 Answers 3

7

Ok, here is a way involving several randomization steps needed to get a "natural" non symmetrical appearance.

I am posting the actual code in Mathematica, just in case someone cares translating it to Matlab.

(* A preparatory step: get your image and clean it*)
i = Import@"https://i.stack.imgur.com/YENhB.png";
i1 = Image@Replace[ImageData[i], {0., 0., 0.} -> {1, 1, 1}, {2}];
i2 = ImageSubtract[i1, i];
i3 = Inpaint[i, i2]

enter image description here

(*Now reduce to a skeleton to get a somewhat random starting point.  
The actual algorithm for this dilation does not matter, as far as we 
get a random area slightly larger than the original elipses *)
id = Dilation[SkeletonTransform[
             Dilation[SkeletonTransform@ColorNegate@Binarize@i3, 3]], 1] 

enter image description here

(*Now the real random dilation loop*)
(*Init vars*)
p = Array[1 &, 70]; j = 1;
(*Store in w an image with a different color for each cluster, so we 
can find edges between them*)
w = (w1 = 
      WatershedComponents[
       GradientFilter[Binarize[id, .1], 1]]) /. {4 -> 0} // Colorize;
(*and loop ...*)
For[i = 1, i < 70, i++,
 (*Select edges in w and dilate them with a random 3x3 kernel*)
 ed = Dilation[EdgeDetect[w, 1], RandomInteger[{0, 1}, {3, 3}]];
 (*The following is the core*)
 p[[j++]] = w =
   ImageFilter[  (* We apply a filter to the edges*)
    (Switch[
          Length[#1],  (*Count the colors in a 3x3 neighborhood of each pixel*)
          0, {{{0, 0, 0}, 0}},          (*If no colors, return bkg*)
          1, #1,                        (*If one color, return it*)
          _, {{{0, 0, 0}, 0}}])[[1, 1]] (*If more than one color, return bkg*)&@
      Cases[Tally[Flatten[#1, 1]], 
       Except[{{0.`, 0.`, 0.`}, _}]] & (*But Don't count bkg pixels*),
    w, 1,
    Masking -> ed,       (*apply only to edges*)
    Interleaving -> True (*apply to all color chanels at once*)]
 ]

The result is:

enter image description here

Edit

For the Mathematica oriented reader, a functional code for the last loop could be easier (and shorter):

NestList[
 ImageFilter[  
   If[Length[#1] ==  1, #1[[1, 1]], {0, 0, 0}] &@
     Cases[Tally[Flatten[#1, 1]], Except[{0.` {1, 1, 1}, _}]] & , #, 1,
   Masking      -> Dilation[EdgeDetect[#, 1], RandomInteger[{0, 1}, {3, 3}]],  
   Interleaving -> True ] &,
 WatershedComponents@GradientFilter[Binarize[id,.1],1]/.{4-> 0}//Colorize, 
5]
8
  • The algorithm you used above, if my understanding is correct, you assign different colors to the watershed regions and once the dilation meet each other,it will return the edge in background color which is in black. It works great in Mathematica. Now i'm trying to apply this method in matlab, hope i can meke it. Thanks.
    – Elsie
    Nov 9, 2011 at 9:46
  • @Ivy Yes. Beware of the randomization part. If you don't, you will end up with a beautiful squared tartan patchwork. Good luck! Nov 9, 2011 at 13:39
  • @Ivy The user Amro is a master in the art of translating Mathematica to Matlab ... just in case Nov 9, 2011 at 14:12
  • thanks for the vote of confidence but that's a difficult one to translate (I am still learning Mathematica myself)... Anyways I will take a look as soon as I can and report back. very nice solution btw +1
    – Amro
    Nov 9, 2011 at 18:28
  • @Amro You largely deserve my praise! Anyway, if you are learning Mma, remember to check out the Mma tag questions and answers in SO. There a lot of knowledgeable users to learn from. Nov 11, 2011 at 4:45
7

What you have as input is the Voronoi diagram. You can recalculate it using another distance function instead of the Euclidean one.

Here is an example in Mathematica using the Manhattan Distance (i3 is your input image without the lines):

ColorCombine[{Image[
   WatershedComponents[
    DistanceTransform[Binarize@i3, 
      DistanceFunction -> ManhattanDistance] ]], i3, i3}]

enter image description here

Edit

I am working with another algorithm (preliminary result). What do you think?

enter image description here

7
  • The OP wants more life-like cells. i.e. less like the Voronoi diagram. Look in the example photo and you'll see some cells have more pressure or whatever and bulge out and squish the others.
    – John Colby
    Nov 1, 2011 at 15:49
  • @John I understand that. I am trying to keep things simple. I'll try another shot when work allows. Nov 1, 2011 at 16:20
  • Nice! I will follow this post. It's a neat problem.
    – John Colby
    Nov 1, 2011 at 16:43
  • 2
    @John What do think of my new drawing? I'm still working out some details! Nov 4, 2011 at 5:12
  • That's really awesome. I love it!
    – John Colby
    Nov 4, 2011 at 15:41
6

Here is what I came up with, it is not a direct translation of @belisarius code, but should be close enough..

%# read image (indexed image)
[I,map] = imread('https://i.stack.imgur.com/YENhB.png');

%# extract the blobs (binary image)
BW = (I==1);

%# skeletonization + dilation
BW = bwmorph(BW, 'skel', Inf);
BW = imdilate(BW, strel('square',2*1+1));

%# connected components
L = bwlabel(BW);
imshow(label2rgb(L))

%# filter 15x15 neighborhood
for i=1:13
    L = nlfilter(L, [15 15], @myFilterFunc);
    imshow( label2rgb(L) )
end

%# result
L(I==1) = 0;                %# put blobs back
L(edge(L,'canny')) = 0;     %# edges
imshow( label2rgb(L,@jet,[0 0 0]) )

myFilterFunc.m

function p = myFilterFunc(x)
    if range(x(:)) == 0
        p = x(1);                %# if one color, return it
    else
        p = mode(x(x~=0));       %# else, return the most frequent color
    end
end

The result:

screenshot

and here is an animation of the process:

animation

2
  • +1 Nice. The masking was used in my answer to provide an extra randomization step. Anyway, I think this is good enough for a start. Congrats! Nov 13, 2011 at 0:18
  • @Amro:Neighboring section cost time...whatever, that's really cool!
    – Elsie
    Nov 16, 2011 at 9:33

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.