Lab: Image Formation

CSC 295 - Computer Vision - Weinman

Summary:

We acquire some of our own images in a laboratory setting and use them to estimate the "noise" in the sensors.

Due:

Friday 2/5

Deliverables

• The script you use to process the images and generate the data in your report as a single .m file
  • Note that you can create a Matlab script by using File -> Blank M-File.
  • You may either add your commands directly to the file as you go, or copy them from your Command History.
• (5 points) One of the original (or stock) images you captured in C.1 (saved in D.3)
• (10 points) The adjusted noise image from F.3
• (10 points) Answer to F.2
• (15 points) The estimated noise plot from F.11
• (20 points) Answers to F.(12,13,14).
• (10 points) Professionalism of write-up

Preparation

Exercises

(A,B,C). Alternate

1. You may use a set of pre-captured images located on the MathLAN in

   /home/weinman/courses/CSC295/images/raw
2. Go to Part D.

A. Getting Started (5 minutes)

1. If you are not already, log in to the MathLAN on one of three machines with the image acquisition software installed (bush, verea, or galler). Ask the instructor to point them out if you cannot find them.
2. Plug the camera into one of the USB ports on the machine.
3. Start the capture software by using the Gnome panel menu and navigating from the main menu (at the lower left) to Applications > Point Grey Research > FlyCapture.
4. The software should bring up a dialog box allowing you to choose a device. There should be only one, so select it and click OK.

B. Acquiring Images (5 minutes)

1. Take the lens cap off of the camera if you haven't already. Move the camera around and notice how the brightness of the image seems to adapt to the amount of light in the scene being captured. We will need to fix this so it does not affect our experiments.
2. Settle the camera somewhere it (or the surface it is resting on) is unlikely to be bumped. Typically this will be a table adjacent to the one your workstation is located on.
3. Point the camera at a scene that will not change. Pay particular attention to shadows and moving wires and/or rollable objects (not to mention others in the lab).
4. When you are sure you have the camera fixed at the unchanging scene you will use for your experiments, click Camera Control. This brings up a control panel with a host of options.
5. In this dialog box, uncheck the "Auto" box for the following settings (leaving their current selections unchanged): Exposure, Pan, Tilt, Shutter, and Gain.
6. Close the window by using the X in the upper right. Now we are ready to acquire some data!

C. Saving Images (10 minutes)

1. When you are satisfied with your image of the scene, select File > Save As or use Ctrl+S to save the image.
2. If the dialog box that comes up is not at your home directory or another place you have write permissions to, click on your home directory in the side of the box under "Places."
3. Important Step!! From the pull down menu in the lower-right that reads "All files (*.*)" select "Raw data (*.raw)" This is important for making sure you store the data in a raw, unprocessed format.
4. Type pic1 in the "Name:" field (omit the extension, .raw is automatically appended for you)
5. Click "Save As"
6. Repeat this process for pic2 through pic10, acquiring a total of 10 images from your scene. Please be very careful to in naming your pictures, saving them in the appropritae place and saving them as raw files. It will be very easy to generate bogus data by moving your camera in the process.
7. When you are through, place the lens cap back on the camera and exit the Flycap program.

D. Loading Images (5 minutes)

1. Start up Matlab. Change Matlab's workspace so that its working directory is the place your images are stored.
2. Whereas before you used imread to load images, those you just stored are nothing but pure bytes. Thus, Matlab does not know how to read them. Instead, you will use a custom procedure called rawread to load them. For instance,

   X = rawread('pic1.raw');

   Load one of your images as above and display it using imshow.
3. Write a generally readable copy (by other programs) of this image in the PNG format using imwrite.
4. Do you see strange plaid checkerboard patterns? This is primarily due to aliasing, which was briefly mentioned in your reading. We will discuss this more shortly. It is also due to another reason. The image you see is the readout (discretized to 0-255) from the sensors. These are color cameras and have an array of sensors that are more sensitive to different wavelengths of light. Zoom in (using the magnifying class) until you can see the individual raw pixels as larger squares. Do they look checkered? Make sure you understand why.
5. You captured 10 images, and we are going to store all of these in a single array, rather than separate arrays, in order to make it easier for us to let Matlab do the work of estimating the image noise. Although you don't need to allocate space for an array beforehand, sometimes this can make your code much quicker. We can force Matlab to pre-allocate 10 slices for our images by assigning one value to the last (lower right) pixel of the last image slice we want:

   X(end,end,10) = 0;

   Do this. Now examine the size of X. Is the result what you expect?
6. Your images can easily be loaded into this array with a for loop and a little creative string processing in Matlab. Actually, there are no strings in Matlab, only arrays of a char type. Thus, strings can be appended in just the same way that arrays can. Assuming you named them as suggested, the loop below will read the kth image into the corresponding slice of X using procedure num2str to convert from a numeric type to a char "string."

   for  k=1:10

     X(:,:,k) = rawread(['pic' num2str(k) '.raw']);

   end;

   Use this to load all of your images, but be sure you understand how it works.

E. Processing Images (10 minutes)

1. Our images are currently stored as bytes. Since we'll be doing math with them, we want to convert them into doubles (usually in the range of 0-1 so we can still display them). Do this using

   X = im2double(X);
2. Remember from your earlier Matlab lab that sum can add things up along a particular dimension. To calculate the standard deviation of each pixel, we first must calculate the mean of each pixel. Use sum (along the third dimension) to add up the values of the pixels in your image array along the slices and divide by the number of slices to give a "mean" image.

   Xmean = sum( ...

   The size of Xmean should be 960×1280. You can verify this using the workspace variables window, the command whos, or the command size.
3. To calculate the standard deviation, we need to take another sum, but we are adding the sum of (the square of) the difference between the mean image and one of the slices. What happens if you try to do this in a vectorized fashion?

   X - Xmean
4. Unfortunately, while we can subtract a scalar from an array, using syntax like X-5, if we are to subtract two arrays, they must be the same size. Since we want to subtract the matrix Xmean from each slice of X, we can use the procedure repmat to REPeat the MATrix across the slices of X. The syntax is

   B = repmat(A,[m n p]);

   where m, n, and p specify how many times to "tile" along the rows, columns, and slices, respectively.
   How many times do we want to tile along rows? One: we need only one copy in the vertical direction. Same with the columns. How many times do we want to tile across the slices? As many slices as there are in X. Use this to create a tiled version of the mean image replicated across the slices so that it is the same size as X.

   Xmeantile = repmat(Xmean, ...
5. Now we are ready to finally calculate the standard deviation. Begin by calculating the difference between your tiled mean and the entire image.

   Xmeandiff = ...
6. Next we need to square the result. We are not taking a matrix power, but calculating the element-wise square of each value in the array. Recall that the .^ operator does this (carefully note the dot before the caret, which makes the operatation element-wise).

   Xdiffsq = ...
7. Now that we have the squared difference, you can calculate the sum using the same technique you used in finding Xmean, only this time note that the denominator of Eq. (2.17) is n-1. This quantity is called the variance.

   Xvar = sum(Xdiffsq, ...
8. Finally, you can take the sqrt of the result to give an "image" containing the standard deviation of each pixel. The result should once again be 960×1280.

   Xstd = ...

F. Analyzing Noise (15 minutes)

1. You can display the deviation image using imshow. However, by default this procedure clips the values with black at 0 and white at 1. In order to show the lowest value in the image as black and the highest value as white, we can specify an empty display range as the second argument:

   imshow(A,[])

   Use this syntax to display your deviation image.
2. Do you notice more noise in certain regions of the image? Describe qualitatively what relationship there is, if any, between the image content and the amount of noise.
3. Save this view adjusted image using imwrite and imadjust:

   imwrite(imadjust(A), 'filename.png');

   where A is the image that you want to have re-scaled to [0-1].
4. We'll try to get a better handle on the noise by looking at a curve along one row of pixels using the plot command, which expects a 1-D array of data to display. Recall that you can extract all the columns for an entire row using indexing like A(row,:). Using this, we can plot the middle row of the of the average image as follows

   figure; % Open a new figure window

   plot(Xmean(end/2,:));

   Is the curve bumpier than you expect? This is because the brightness values plotted next to each other in the curve correspond to different color channels, thus registering different sensitivities to the light hitting the sensor.
5. We can extract just one of the color channels by looking only at the odd columns of one row. Write a command that plots only every other column in the same row of X, starting from 1. Does this look smoother?
6. If you didn't open a new figure, you may have noticed that plot overwrites whatever you had plotted previously. We are going to add two more curves to the plot you just made. In order to do this, we need to tell the figure to keep what was already there:

   hold on;

   Use this command to "hold on" to the current plot.
7. Next, plot the same brightness curve but with the standard deviation added to it. That is, if m is a brightness point of a particular pixel on the plot, you need to plot m+s, where s is the standard deviation of that pixel. To do this, you will extract the same array (or vector) as in the previous command and add to it a similar array extracted from the standard deviation image. (Make sure you are still only taking every other column of the data.)
8. Add to your plot the analogous curve for m-s (brightness minus standard deviation).
9. Use the magnifying glass to zoom into an interesting (there should be some change), but relatively stable (no wild swings) region of your new graph so that you can see all three curves.
10. Use the xlabel and ylabel commands to label the axes with something like 'Column' and 'Brightness'
11. Save your figure as a PNG file using the print command:

    print('-dpng','filename.png');

    Unless you specify a directory, print will deposit the file in the Matlab desktop's current working directory. When you are writing up your lab, you may find that the figure has an excessive amount of whitespace. You can trim this automatically at the (login shell) command line using the rather handy ImageMagick program convert. The shell syntax to auto-crop such figures is:

    $ convert -trim in-file out-file

    The input and output file can have the same name.
12. Roughly how many discrete gray levels does the noise span? (Hint: Multiply the difference between the top and bottom of the curves-twice the standard deviation-by 256.)
13. What is the [estimate of the] average (over all the pixels) of the noise?
14. What is the [estimate of the] worst-case acquisition noise? (Hint: Matlab has a max function.)