Tuesday, June 25, 2013

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

    Monday, June 17, 2013

    Activity 4: Image Types and Formats


    In this activity we will be dealing with different Image Types and Formats. Image Types come in four different types: 
    a. Binary Images 
    b. Grayscale Images: 
    c. Truecolor Images
    d. Indexed Images



    Details:
    Item Type: JPEG Image
    Size: 1.10 MB
    Dimensions: 1920 pixels x 1200 pixels
    Horizontal Resolution: 96 dpi
    Vertical resolution: 96 dpi
    Bit Depth: 24

    Anyway here is a camera image of my family's dog, Tami. The details are as follows.


    Details:
    Item Type: JPEG Image
    Size: 3.00 MB
    Dimensions: 4608 pixels x 3456 pixels
    Horizontal Resolution: 180 dpi
    Vertical resolution: 180 dpi
    Bit Depth: 24
    Resolution unit: 2
    Color representation: sRGB
    Compressed bits/pixel: 3

    Camera maker: Canon
    Camera model: Canon PowerShot A4000 IS
    F-stop: f/3
    Exposure time: 1/8 sec.
    ISO speed: ISO-400
    Exposure bias: 0 step
    Focal length: 5 mm
    Max. aperture: 3.15625
    Metering mode: Pattern
    Flash mode: No flash, compulsary


    In this activity, we will be experimenting with different file formats. There are many file formats. File Formats are categorized as either Raster or Vector. Raster types describes the characteristics of each individual pixels. Vector formats, on the other hand contain a geometric description that can be rendered smoothly at any display size. [2]  Some popular Raster formats include: JPEG, TIFF, RAW, GIF, BMP, PNG, PAM, WEBP. Vector Format images need to be "rasterized" be displayed on digital monitors. Cathode Ray tube technology (such as in radars, video games, medical monitors) make heavy use of Vector formats. Examples of Vector formats include CGM (computer graphics metafile), SVG, and Gerber File formats.
    The scope of this activity will unfortunately fall only on Raster file formats, specifically, PNG, JPEG, TIFF, and BMP.

    Here is an image of Juggernaut Dota 2 hero originally in PNG file file, courtesy of http://kkcdn-static.kaskus.co.id/images/2012/07/30/1101061_20120730125032.png






    Here is the same file converted to
    a. JPEG
        Dimensions: 481 x 432
        Size: 39.3 KB
        Bit Depth: 24
    b. GIF
        Dimensions: 481 x 432
        Size: 42.7 KB
        Bit Depth: 8
    c. Bitmap Image - 24-bit
        Dimensions: 481 x 432
        Size: 205 KB
        Bit Depth: 8
    d. Bitmap Image -  16 Color-bit
        Dimensions: 481 x 432
        Size: 103 KB
        Bit Depth: 4
    e. Bitmap Image - Monochrome bitmap
        Dimensions: 481 x 432
        Size: 27 KB
        Bit Depth: 1

    PNG or Portable Network Graphics supports lossless data compression. It was created to improve and replace GIF (Graphics Interchange Format). PNG is used for the 2D interlacing method, the cross-platform gamma correction, and the ability to make anti-aliased balls, buttons, text and other graphic elements. It is currently the most lossless data compression image file format on the web. PNG started in 1977 and 1978 when two Israeli researchers, Jacob Ziv and Abraham Lempel, first published a pair of papers on a new class of lossless data-compression algorithms, now collectively referred to as ``LZ77'' and ``LZ78.'[3]
    In 1983, Terry Welch of Sperry (which later merged with Burroughs to form Unisys) developed a very fast variant of LZ78 called LZW.[3]
    JPEG or Joint Photographic Experts Group images allow tradeoff between quality image and file size. JPEG allows 10:1 compression ratio. JPEG utilizes Discrete Cosine Transform which mathematically converts each frame from spatial (2D) domain to frequency domain.
    GIF (Graphics Interchange Format) is a bitmap image format created in 1987. It was originally called 87a. Later, 89a added support for animation delays, transparent background color, and storage of application-specific metadata. [4] As of now all relevant patents have expired. GIF's are now widely used for its simple animation found in websites such as 9gag.

    Sources:
    1.http://web.archive.org/web/20080610170124/http://www.codersource.net/csharp_color_image_to_binary.aspx
    2. A4 - Image Types and Formats 2013.pdf, Soriano, Jing
    3. History of the Portable Network Graphics, Greg Roelofs, http://linuxgazette.net/issue13/png.html
    4. Graphics Interchange Format, Version 89a, 1 July 1990 Retrieved 6 March 2009

    Thursday, June 13, 2013

    Activity 3 Scilab Basics

    So our activity for today was to get familiar with the programming language that is Scilab. The objective was to create synthetic images out of test algorithms. 

    Centered Square Aperture
    Code Snippet in Scienotes:

    nx = 100; ny = 100; //defines the the number of elements along x and y
    x = linspace(-1,1,nx); //defines the range
    y = linspace(-1,1,ny);
    [X,Y] = ndgrid(x,y); //creates two 2D arrays with x and y coordinatesnx = 100;
    l = abs(X); //length
    w = abs(Y); //width
    A = zeros (nx,ny);
    A (find(l<0.5 & w<0.5) ) = 1;
    f = scf();
    grayplot(x,y,A);
    f.color_map = graycolormap(32);

    As one can see I could have made a rectangle just by simply editing the length and width. It can be edited on lines 5 and 6.


    Sinusoid along the x-direction (corrugated roof)
    nx = 100; ny = 100; //defines the the number of elements along x and y
    x = linspace(-1,1,nx); //defines the range
    y = linspace(-1,1,ny);
    [X,Y] = ndgrid(x,y); //creates two 2D arrays with x and y coordinatesnx = 100;
    A = zeros(nx,ny);
    A = 10*sin(20*X); //defines the sinusoid with frequency and amplitude
    f = scf();
    grayplot(x,y,A);
    f.color_map = graycolormap(32);

    This was achieved by simply letting the matrix (A) have the assigned values of sine with varying frequency on line 6 of the code snippet. The frequency is 20 while the amplitude is 10.


    Grating along the x-direction

    nx = 100; ny = 100//defines the the number of elements along x and y
    = linspace(-1,1,nx)//defines the range
    = linspace(-1,1,ny);
    [X,Y] = ndgrid(x,y)//creates two 2D arrays with x and y coordinatesnx = 100;
    = zeros (nx,ny);
    = 10*sin(20*X);
    (find(g<0)) = 1;
    = scf();
    grayplot(x,y,A);

    f.color_map = graycolormap(32);
      


    Annulus
    nx = 100; ny = 100//defines the number ofelements along x and y
    x = linspace(-1,1,nx); //defines the range
    y = linspace(-1,1,ny);
    [X,Y] = ndgrid(x,y); //creates two 2-D arrays of xand y coordinates
    r= sqrt(X.^2 + Y.^2); //note element-per-element squaring of X and Y
    A = zeros (nx,ny);
    A (find(r<0.7 & r>0.3) )= 1;
    f = scf();
    grayplot(x,y,A);

    f.color_map = graycolormap(32);

    Circular Aperture with graded transparency (Gaussian transparency)

    nx = 100; ny = 100//defines the number ofelements along x and y
    = linspace(-1,1,nx)//defines the range
    = linspace(-1,1,ny);
    [X,Y] = ndgrid(x,y)//creates two 2-D arrays of xand y coordinates
    = exp(-X.^2 -Y.^2);
    = scf();
    grayplot(x,y,Z);
    f.color_map = graycolormap(32);




     Here is a product of the Sine function and Gaussian function.

    Code Snippet:
    nx = 100; ny = 100; //defines the the number of elements along x and y
    x = linspace(-1,1,nx); //defines the range
    y = linspace(-1,1,ny);
    [X,Y] = ndgrid(x,y); //creates two 2D arrays with x and y coordinatesnx = 100;
    G = sin(10*X)*exp(-X.^2 - Y.^2);
    f = scf();
    plot3d(x,y,G);
    f.color_map = graycolormap(32);

    Pretty cool huh?