Lab (5) - Image Processing
Image Restoration

Objectives:
In this lab you will learn about:
Creating noisy image
Estimation of the noise

Please DO NOT Print this lab as it will not appear correctly

Preparations:
Copy all the files shown below in your work directory.  If you need more files, save them from the link on the notes web page.  You will prepare 2 MS Word files one for the inlab part and the other for the postlab.  E-mail the inlab file at the end of  the lab and the postlab by the due date.  You need some files that are not included in the Toolbox.  Save these files in your work directory.
eb01.bmp
histroi.p
statmoments.p
imnoise3.p
imnoise2.p
noisy_org.bmp

A Model of the Image Degradation/Restoration Process
The degradation process is modeled as a degradation function and an additive noise operating on an input image f(x,y) to produce a degraded image g(x,y):

Where H is the degradation function and the last term is the additive noise. Giving g(x,y) and some knowledge about H and the additive noise, we can come up with an estimate for f for the original image.  If H is a linear, spatially invariant process, the degraded image can be obtain as:

where h(x,y) is the spatial representation of the degradation function and "*" denotes the convolution.  In frequency domain, this can be written as:

When everything are shown in Fourier Transform.  We first look into modeling the noise and assume H is the unity matrix. Basically, for now assume it has no effect.

In MATLAB, you can add noise to an image using the imnoise function.  The syntax for this function is:
g = imnoise(f, type, parameters)
Where f denotes the image, type is the type of noise we wish to add, and parameters are the parameters required for the noise type we are using.  Note: The imnoise converts the image in the range [0 1].

One of the noise types we can add is salt & pepper.  The syntax for this type of noise is:
g = imnoise(f, 'salt & pepper', d)

Where d is the noise density which is the percentage of the image containing noise values.  Thus, approximately d*numel(f) pixels are affected by this type of noise.  The default is 0.05.  See page 143 for more information and details.

Figure (1)
You have stored the image shown in Figure (1) in bmp format at the beginning.  If you haven't done so, right click here (eb01) to save the image in your work directory.

Lab Activity (1)

Corrupt this image by the percentages given below and save the resulting image in your MS Word file.  Mark the images with the corresponding percentages.

Act1A) 5% Salt and Pepper
Act1B) 25% Salt and Pepper
Act1C) 50% Salt and Pepper
Act1D) 75% Salt and Pepper
Act1E) 100% Salt and Pepper

Let me ask you a question. You know what you have added in each part, can you take the additive noise back to obtain the original image back?

Similarly, you can corrupt the image using Gaussian noise.  The syntax for this type of noise is:
g = imnoise(f, 'gaussian', m, var)

Where m is the mean and var is the variance.  The default value of m is 0 and that of var is 0.01.
Note: Due to the fact that imnoise converts the image to the range [0 1], you need to convert m and var based on the class type of the image.  For example, for an image in uint8 class type, a Gaussian noise with mean 64 and variance of 400 will have m = 64/255 and var = 400/(2552).

Lab Activity (2)
In this lab activity, you are required to corrupt the original image with 4 combinations of m and var (4 cases ONLY).  To help you out with the type of value you may want to select, I will explain how you can find the mean and variance of the original image.  This may give you an idea of what kind of values you can choose for the reasonable mean and variance.

>> mean(mean(f))   % will compute the mean of the image
>> sqrt(std(f(:)))   % will compute the variance of the image which is the sqrt of standard deviation

Use imshow to display the 4 resulting noisy images and include them in your inlab file.  Mark the images clearly.

Once again, don't forget imnoise works in [0 1] range, so you need to convert m and var.

Generating Spatial Random Noise with a Specified Distribution
In MATLAB, function imnoise2 can be used to generate random noise based on a given Cumulative Distribution Function (CDF).  The syntax for this function is:
R = imnoise2(type, M, N, a, b)
This generates an array R of size M-by-N whose elements are random numbers of the specified type with parameters a and b (See Page 148 for more details).

Example:  The following command line:
r = imnoise2('gaussian', 100000, 1, 0, 1);
produces 100000 random values (M is 100000 and N is 1) with a Gaussian distribution with mean of 0 and standard deviation of 1.

You can display the histogram of this distribution using:
p = hist(r, bins);  % where you set the number of bins, let's say 50
bar(p);

With 50 bins, I have produced:

Lab Activity (3)
Create two more histogram using two other types as defined on Page 148.  I let you choose two different types as you wish. You can use the default values for the parameters that you need in each case.

Here are the types:
1) uniform - uniform random numbers in the interval (a, b), default (0, 1).
2) salt & pepper - number of 0 with probability a, and number of 1 with probability b, default for both is 0.05.
3) lognormal - lognormal numbers with offset a and shape parameter b. Default is a = 1 and b = 0.25.
4) rayleigh - Rayleigh noise with parameter a and b.  Defaults are a = 0 and b = 1.
5) exponential - Exponential random numbers with parameters a and b.  Defaults are a = 1, b does not matter.
6) erlang - Erlang (gamma) random numbers with parameters a and b.  B must be positive integer.  Defaults are a = 2 and b = 5.

Estimating Noise Parameters
In MATLAB, there is a function called  roipoly that allows you to select a Region Of Interest  (ROI) from an image to use as the noise matrix.  You can use the mouse's left button to select the area of image that you wish to use as the noise. Remember that the less feature you have the closer the selected area is to the actual noise.

Example:
>> clear  % this clears the MATLAB memory
% now use the left mouse click to select a polygon that best describe the noise, once done with the selection of last corners, press Enter
>> [B, c, r] = roipoly(f);

After you press Enter, you will get:
B which is the same size image as the original image but in binary, r are pairs holding the coordinates of the corners of the polygon.  To view the mask, you can try:
>> imshow(B);

This mask marks the polygon that you have selected to represent the noise.  The white area contains 1's and the black area contains 0's.  In theory, if you multiply this mask by the original image, it will produce the noise matrix of the same size as the original image with non-zero data appearing in the polygon area.

To display the histogram of the selected area, we can use:
[ p  npix] = histroi(f, c, r);  % find out what p and npix are
figure, bar(p, 1);

The next thing is to use the standard noise functions on Page 148 and play with parameters until you find a histogram that is a match to this one or is very close to it. If we are lucky enough, then we have found a model for our noise.  Remember, the parameters are very important in such cases.

By the way, the above histogram may not reveal much on the noise because it is taken from the background area of the image which sometimes is filled after the image is taken to make it to a suitable size.  It was best to choose the mask somewhere inside the image where there is not much feature.

Lab Activity (4)
Try to select a polygon inside the image itself where you think is a good area with minimum feature. Then show your mask and histogram in your inlab file.

Postlab
1. Repeat the same procedure in the last part of the lab to find out the histogram for the noise in the following image and use the standard noise functions given in the textbook and to come up with a noise model that matches the noise in this image.  Please use this link to store the image you haven't done that already.

You need to submit the noise mask, the histogram of the noise, and the histogram of the standard noise functions you have used with the parameters to come up with the matching noise.  Note, you may need to spend sometimes to find the histogram that matches. Find the best match and write down the parameters.  See Page 156-157 for more information.

2. Use the procedure explained in section 5.3 to create a noiseless image (as noiseless as possible) of the above image.