Tuesday, June 21, 2011

A2 - Scilab Basics

Solving real world problems is a very complex task. The observable and non-observable interactions between components influence the over-all characteristic of a certain real world phenomena. However, the complicated tasks can be partially solved by computer programming. Programming as a scientific tool is very useful in recreating the physical world by implementing rules and variables which can somewhat simulate the present governing rules. In relation to this, our task is to familiarize ourselves with Scilab, powerful scientific programming language which is a free counterpart of Matlab.

The first thing to do is install Scilab and an image processing toolbox called SIP. However, for current versions of Scilab, SIVP (Signal, Image and Video Processing) Toolbox will be more compatible.

For beginners, Scilab can be learned by just reading its documentation and tutorials at Scilab Tutorial and practicing code creation.

In order to test the depth of our understanding in using Scilab, an exercise on creating patterns is followed. These patterns can be used in masking and as apertures in optical systems. 

The first thing to do is to used the code given by Dr. Soriano to create a centered circle which can simulate a circular aperture or pinhole. The resulting image is shown below(figure 1) and the code to produce it follows after.

 
Figure 1.  Centered circle/circular aperture.

Code for a Circular Aperture:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates
r = sqrt(X.^2 + Y.^2);   //note element-per-elementsquaring of X and Y
A = zeros(nx,ny);   //creating a matrix of zeros
A(find(r<0.7)) = 1;  //replacing the zeros by 1 when the condition r<0.7 is satisfied
//imshow(A, []);   //showing the image produced
imwrite(A,'centered_circle_500.png');   //saving the image

The next task is to create codes to produce synthetic images described as a centered square aperture, sinusoid along the x-direction(corrugated roof), grating along the x-direction, annulus and a circular aperture with graded transparency(a gaussian transparency to be exact).

A. Centered Square Aperture 
A centered square aperture can be created as simple as creating a centered circle(figure 1) but instead of implementing a condition to replace the zeros of A by one for a radius r, a condition is changed to -0.7<X<0.7 and -0.7<Y<0.7. The resulting image is in figure 2 and the corresponding code follows.

Figure 2.  Centered square/rectangular aperture.

Code for a Rectangular Aperture:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates

A = zeros(nx,ny);   //creating a matrix of zeros
A(find(X<0.7 & X>-0.7 & Y<0.7 & Y>-0.7)) = 1;  //replacing the zeros by 1 when the condition is satisfied
//imshow(A, []);   //showing the image produced
imwrite(A,'centered_square_aperture_500.png');   //saving the image


B. Sinusoid along the x-direction(corrugated roof) 
In a corrugated roof, the only thing to do is to operate a sine function to the X points. The resulting image is shown in figure 3 while the code follows.

Figure 3.  Sinusoid along the x-direction(corrugated roof).

Code for a Corrugated Roof:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates

A = sin(6 * %pi * X);;  //operating a sine function on X
//imshow(A, []);   //showing the image produced
imwrite(A,'corrugated_roof_500.png');   //saving the image


C. Grating along the x-direction 
To create a grating along the x-direction, it is as if creating a rectangular aperture but varying the locations of the ones. Here, I divided a 500x500 pixel area into ten rectangular strips, 5 of which have a value of 1 and 5 have a value of 0.
Figure 4.  Grating along the x-direction.

Code for a Grating:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates

A = zeros(nx,ny);   //creating a matrix of zeros   
for i = 0:4
  A((1 + 100*i):(nx/10 + 100*i), 1:ny) = 1;   //
replacing the zeros by 1 when the condition is satisfied
end
//imshow(A, []);   //showing the image produced
imwrite(A,'grating_500.png')   //saving the image

  
D. Annulus 
An annulus is the same as a circular aperture, the only difference is that the condition that r<0.7 is changed to 0.4<r<0.7.
Figure 5.  Annulus.


Code for an Annulus:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates
r = sqrt(X.^2 + Y.^2);   //note element-per-elementsquaring of X and Y
A = zeros(nx,ny);   //creating a matrix of zeros
A(find(r<0.7 & r>0.4)) = 1;  //replacing the zeros by 1 when the condition 0.4<r<0.7 is satisfied
//imshow(A, []);   //showing the image produced
imwrite(A,'annulus_500.png');   //saving the image


E. Circular Aperture with Graded Transparency(gaussian transparency) 
In this part, the same technique used to create a circular aperture is used. However, this time, instead of replacing the zeros with ones for certain radius r, a gaussian function of the form exp(-r^2) is used where I consider the standard deviation b, as 1 in the original gaussian distribution of exp(-x^2/2*b^2).
 
Figure 6. Circular aperture with gaussian transparency.

Code for a Circular Aperture with Gaussian Transparency:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates
r = sqrt(X.^2 + Y.^2);   //note element-per-elementsquaring of X and Y
A = exp(-10*r.^2);   //operating a gaussian function
//imshow(A, []);   //showing the image produced
imwrite(A,'gaussian_transparency_500.png');   //saving the image


I then created two more images which are often used in real world optical systems, the double slit(figure 7) and a cross slit(figure8).

Figure 7. Double Slit

Figure 8. Cross Slit.

Code for a Double Slit:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates

A = zeros(nx,ny);   //creating a matrix of zeros
A(find(Y<-0.1 & Y>-0.2 & X<0.8 & X>-0.8)) = 1;   //replacing the zeros by 1 when the condition is satisfied

A(find(Y<0.2 & Y>0.1 & X<0.8 & X>-0.8)) = 1;
//imshow(A, []);   //showing the image produced
imwrite(A,'double_slit_500.png');   //saving the image


Code for a Cross Slit:

nx = 500; ny = 500;   //defines 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 2-D arrays of x and y coordinates

A = zeros(nx,ny);   //creating a matrix of zeros
A(find(Y<0.1 & Y>-0.1 & X<0.8 & X>-0.8)) = 1;   //replacing the zeros by 1 when the condition is satisfied

A(find(X<0.1 & X>-0.1 & Y<0.8 & Y>-0.8)) = 1;
//imshow(A, []);   //showing the image produced
imwrite(A,'cross_500.png');   //saving the image



In summary, the goals of this activity are to:
  1. Familiarize ourselves with Scilab.
  2. Create codes that can simulate apertures in optical systems using Scilab.
Based on the goals shown above, I would give myself a score of 12. 10 for being able to do the activity alone and creating images which can serve as optical tools in optical systems using the power of Scilab; and a bonus of 2 for taking the initiative of creating additional images for other kinds of apertures.



No comments:

Post a Comment