Activity 5: Area Estimation for Images with Defined Edges

Finding an area of a 2-dimensional image has many practical applications. (cancer research, remote sensing, automated product inspection [1] ) The problem arises when the image is non-uniform and asymmetric. Think of a blob or a poorly drawn circle. One cannot measure the area so easily. The best way to measure the area of such images is with the use of Green's Theorem.

Green's Theorem relates a double integral to a line integral. Let Region R with contour taken in the counterclockwise direction with F1 and F2 be continuous everywhere and with continuous partial derivatives:

                      

Which can be reduced  by letting  pairs F1 = 0 and F2 = x and second pair  F1 = -y and F2 = 0 via averaging to:

The Green's Theorem now becomes summation of edge point products of xi, yi+1 and xi+1, yi.

In summary, the area of any closed figure regardless of shape can be computed via Green's Theorem.


To obtain the edges of the images Scilab has a library of functions which finds edges of single channel image:

edge(im,method)


The methods are: 

1. sobel- discrete differentiation operator, computing an approximation of the gradient of the image intensity function; the result of the Sobel operator is either the corresponding gradient vector or the norm of the vector. [2]
2. prewitt - similar to sobel operator
3. log - 
4. fftderiv - detects zero-crossing of the second-order directional derivative in the gradient direction [3]
5. canny - function described by the sum of four exponential terms, but can also be approximated by the first derivative of a Gaussian.[4]

The following methods will be examined. The edge-finding algorithm which computes closest to the theoretical image area will also be determined.

So here I created a 2-dimensional image of a circle: (from Activity 1 actually)

Figure 1

This image is 300 x 300 with a circle radius of 105 pixels.
The edge detection algorithms are applied on the image above and the following results are observed:

Figure 2

In terms of pixel count, the theoretical area of the circle is 3.4636e+4 pixels.
Now I will discuss the very interesting method on how I was able to compute the area of the circle via Green's Method. First I will show the algorithm/code snippet of the code:

(1)Edge1 = edge(Image, 'canny'); //sobel, prewitt,log,fftderiv,canny
(2)imshow(Edge1);
(3)imwrite(Edge1,'canny.bmp');
(4)[i,j] = find(Edge1); //locates coordinates of the edges
(5)xc = sum(i)/(size(i,2));  // locates center coodinates
(6)yc = sum(j)/(size(j,2));  // locates center coodinates
(7)x =j-xc;
(8)y =yc-i; //subtracts center coordinates to all edge coordinates 
(9)r = sqrt(x^2 + y^2) //computes radius
(10)theta = atan(y,x)*(180/3.14159);  //computes theta
(11)z = [theta; x; y; r]';// creates matrix containing theta, x and y edge coordinates, and radius and transposes them
(12)[C, I] = gsort(z(:,1),'r','i'); //sorts Theta 
(13)sorted = z(I,:) // creates Index and sorts Theta with x and y following it
(14)xiy = sum(sorted(1:size(sorted,1)-1,2).*sorted(2:size(sorted,1),3)) //Applies the Green Function
(15)xyi = sum(sorted(2:size(sorted,1),2).*sorted(1:size(sorted,1)-1,3))//Applies the Green Function
(16)Area = 0.5*(xiy - xyi)//Applies the Green Function Area Formula
(17)Area2 = 3.14158*(nx*0.35)^2 //Area of the Circle
(18)Error = 100*sqrt((Area-Area2)**2)/Area2

The code basically goes like this:
I locate the edge points (white pixels) and their coordinates. The coordinates are not x-y cartesian and uses row x column. Lines 7 and 8 compensate for this. I then compute the corresponding theta of the x and y coordinates. The list of angles are then matched with their respective x and y coordinates as in line 11, and the matrix is transposed for convenicence. Then here lines 12 and 13 are used. I found this code snippet online http://www.equalis.com/forums/posts.asp?group=&topic=302065&DGPCrPg=1&hhSearchTerms=&#Post302949, which allows me to sort rows by order of their column. Specifically I need the x and y coordinates sorted by increasing angle.
Finally, lines 14 and 15 computes the Green's function using equation 3. The sorted x and y coordinates are computed. I had a hard time implementing the Green's function. But in my opinion, this is the most elegant and simple method to do so. All thanks to the operator ".*" This operator multiplies the matrix element by element. Without the dot, an error occurs. 

After which the two areas are compared and here are the results:

Theoretical Area:  3.4636e+4
Area - % Error:

sobel-  33830.75 - 2.32%

prewitt - 34334.25 - 0.87%

log - 34386.25 -0.72%

fftderiv - 34334.25 -  0.87%

canny - 34328.28 - 0.89%

The results show that the "log" method seem to be the best method in determining the area of a given image. "sobel" method seem to be the worst. As one can see, the circle image formed in the "log" method is well-defined and most importantly continuous as compared to the "sobel" method which has some gaps in the circumference of the circle. For all intents and purposes, I will use "log" and "prewitt" method.

Here is my attempt to draw a PERFECT circle =P

Figure 3

Here is the edge finding method using "log:"

Figure 4

Here is the edge finding method using "prewitt:"

Figure 5

Area for the circle is 

43601.23 pixels -"log" method

43562.70 pixels - "prewitt" method

Conclusion:

I need to practice drawing my circles. XD

Now I want to measure the land area of our house in Cavite. So I take a google map image of our lot, which looks like this:

Figure 6

Erase the land area off the image via Paint:

Figure 7

Crop, then perform edge finding methods (prewitt and log) on the cropped image (142 x 104):

Figure 8

The new image is now this, which area can be easily computed:

Note that 76 pixels is equivalent to 10 meters (from the original image). Which means there is 0.1315789 meters per pixel  or  0.0173130 sq. meters per  square pixels.  The computed pixel area is 

10800.88 sq. pixels *0.0173130 = 186.99  sq. meters.

I cleaned off the area and recomputed:

10390.70 sq. pixels * 0.0173130 = 179.89419 sq. meters.
Our land area according to the blueprints from our Architect says its 180 sq. meter!!!!
Yay!!!
One of my favorite activities, I must say. I enjoyed coding this part...

References: 

1. Jing Soriano A5- Area estimation of images with defined edges 2013
2-5. wikipeida.org