Monte Carlo Simulation to estimate area from image

Questions
1

Hello all

I have a little problem which I thought I could solve with Mathematica 9.

I am interested in calculating, or shall I say, approximating the area of a non defined geometric shape.

My input is an image of a small “canal” created in the human body. This canal is colored in black due to limitations of the imaging device. I need to calculate the area of this canal, so I can take another image a month later, and see if is still the same size.

My idea was to use Monte Carlo simulation. To be more specific, what I want to do, is to somehow draw a rectangle to contain the shape in it. The rectangle should have known dimensions, so I can know its area. Then, I want by simulation, to take n random points inside this rectangle. For each point, I thought to check if it’s black or not (the rectangle can be white inside). If the point is black, it’s in the shape, if not, it is not. Knowing the area of the rectangle, I thought I could approximate the area of the non defined shape.

In order to to this, I thought to try it first with known shapes, i.e., to create a rectangle and inside it to put a circle with a known radius, and then to run a simulation.

My problem - I do not even know where to start . Can this be implemented with Mathematica ? If yes, how ?

Thanks for any idea !

2

Try this and see if you can adapt it to your problem.

In[1]:= g=Graphics[{Red, Rectangle[{5,10},{20,30}], Black, Polygon[{{8,15},{10,18},{18,22},{17,21}}]}]
Out[1]= <snip>

In[2]:= pixels = ImageData[g];

In[3]:= Dimensions[pixels]
Out[3]= {432, 328, 3}

In[4]:= p = Partition[Flatten[pixels], 3];

In[5]:= Count[RandomSample[p, 10^4], {0., 0., 0.}]/10^4 // N
Out[5]= 0.026

RandomSample chooses distinct points so the 10^4 must be less than the number of pixels in your image.
{0.,0.,0.} is “pure black” and if you are processing a captured image you may need to adjust that threshold to get all you accept as black.
Try this repeatedly on constructed images where you can manually verify that the result is approximately correct before depending on this.

If you look at this, which attempts to have half red and half black pixels and which finds a list of {color,number of pixels} for each color from ImageData,

In[6]:= g=Graphics[{Red,Rectangle[{5,10},{20,30}],Black,Polygon[{{5,10},{20,10},{20,20},{5,20}}]}];
    pixels = ImageData[g];
    p = Partition[Flatten[pixels], 3];
    Count[RandomSample[p, 10^4], {0., 0., 0.}]/10^4 // N
    Map[{#, Count[p, #]} &, Union[p]]

Out[7]= 0.4619

Out[8]= {{{0.,0.,0.}, 65208}, {{0.12549,0.,0.}, 418}, {{0.74902,0.74902,0.74902}, 312}, {{0.780392,0.780392,0.780392}, 2}, {{1.,0.,0.}, 65626}, {{1.,1.,1.}, 10130}}

it finds 65208 black pixels and 65626 red pixels, but also finds a small number of other colors and 10130 white pixels, so ImageData appears to be actually including more than just the contents of that colored rectangle and you should be careful using this.