digitalimages.m

Objectives The purpose of this notebook is to give you a brief introduction to the DiscreteWav

% conversion.  The NTSC map is
%
% gray = .299*red + .587*green +  .114*blue
%
% In the cell below, load and plot the legos image that comes with the
% package.   Next create two grayscale images from this color image. Build
% the first by averaging the three channels and construct the second using
% the NTSC map.  Call the first image gray1 and the second gray2.  Plot 
% both images - do you notice a difference? 
%
% Repeat the exercise with some other color images included with the
% package or one from the internet.

%%
% Place your answer here.  To type an answer, start a line with %.  Matlab
% commands are entered as usual.

%%
% *Exercise 4*
%
% Gamma correction is a simple image processing tool that can be used to
% lighten or darken an image.  The process is quite simple.  Suppose A is a
% grayscale image matrix.  We first divide all elements of A by 255 to 
% obtain values in the interval [0, 1] and then raise each value to some 
% positive exponent |gamma|.  We then multiply each element by 255 and 
% round the result to the nearest integer.
%
% In Matlab we can do "incorrect" mathematical operations.  If A is
% our image matrix, then A/255 makes some sense, but (A/255).^gamma 
% actually raises each element of A to the gamma power.  In the cell 
% below, using the |round| function, write code that will perform gamma 
% correction on the image from Exercise 1.  Call the output B.  Try several
% values of |gamma| (I would suggest naming a variable gamma and then 
% changing its values) and plotting the result.  
%
% Can you characterize what different values of gamma do to the image?
%
% Feel free to load other grayscale images in the package and perform gamma
% correction on them.

%%
% Place your answer here.  To type an answer, start a line with %.  Matlab
% commands are entered as usual.

%%
% *Exercise 5*
%
% In an example above we computed the negative of a color image.  In the
% original image, each pixel was represented by a linear combination of 
% vectors representing red, green, and blue.  For the image negative, what
% vectors are used to form the linear combination for each pixel?  What are
% the associated colors?

%%
% Place your answer here.  To type an answer, start a line with %.  Matlab
% commands are entered as usual.


%% Cumulative Energy
%
% One of the important skills students develop in this course is the
% ability to write modules or functions.  The cumulative energy function is 
% an excellent starting point for this development. 
%
% To compute the cumulative energy of vector v, we sort the absolute values
% of the elements largest to smallest and then square the components of the
% resulting vector.  Call this vector y.  We find the kth element of the 
% cumulative energy vector by summing the first k elements of y and 
% dividing the result by the norm of v squared.
%
% The first two steps of this process are easy.  Let's use a small vector
% as an example.

v = [-3 -2 0 0 5 -8 3 1 1 2]
y = fliplr(sort(abs(v)))
y = y.^2

%%
% Next comes the cumulative sums.  We could write a loop for k = 1 through
% 10 and an inner loop j = 1 to k to perform the task, but we instead 
% encourage our students to take advantage of built-in functions that 
% accomplish tasks.  In most cases, these built-in functions are much 
% faster than attempting to extract individual elements from a vector or 
% list.  For cumulative sums, the built-in command |cumsum| is well-suited
% for our needs. 

x = cumsum(y)
ce = x/(v*v')

%%
%
% We can plot the cumulative energy using the Matlab command |plot|.
% Note the x-axis is the element number of ce.

plot(ce);

%% The Cumulative Energy Function
%
% In Matlab, it is very easy to write modules or functions.  In the My
% Documents folder, you will find a copy of the DiscreteWavelets Toolbox -
% the folder is called DiscreteWavelets.  In this folder, you should find a
% folder called MyFiles.  In this folder, open the m-file called MyCE.m
%
% Using this file, we will write a function for computing the cumulative
% energy of a vector (or matrix!) together.

%%
% 
% The cell below creates a vector of 20 random integers in the range 
% 0,...,10, and then calls the |MyCE| function and the DiscreteWavelets 
% command |CE| to compute the cumulative energy.  Don't execute the cell 
% until you have completed and executed the cell above.

v=round(rand(1,20)*10)
MyCE(v)
CE(v)

%% Exercises
%
%%
% *Exercise 1*
%
% The |MyCE| function above will not work on matrices.  The function is 
% expecting a vector.  How can you modify the |MyCE| function so that it will
% work on matrices?  (Hint: Check the command |reshape| under Help.)

%%
% Place your answer here.  To type an answer, start a line with %.  Matlab
% commands are entered as usual.

%%
% *Exercise 2*
%
% Load the grayscale image of the chess pieces from the DiscreteWavelets
% package, find its cumulative energy using the |MyCE| function, and then 
% plot the cumulative energy of the image.  How many elements of the image 
% constitute 90% of the energy?  How many zeros are in the image?

%%
% Place your answer here.  To type an answer, start a line with %.  Matlab
% commands are entered as usual.

%%
% *Exercise 3*
%
% Describe the plot of the cumulative energy of a (nonzero) constant
% vector.  (Hint: If you wish, you can generate a constant vector using the
% |ones| command - see Help for more information.)

%%
% Place your answer here.  To type an answer, start a line with %.  Matlab
% commands are entered as usual.

%% Entropy
%
% The DiscreteWavelets command for entropy is |Entropy|.  The function 
% takes either a vector or matrix as input.  Here are some example calls:

v = 1:16
str=sprintf('The entropy of vector v is %f.', Entropy(v));
disp(str);
A = eye(8);
str=sprintf('The entropy of the 8x8 identity matrix is %f.', Entropy(A));
disp(str);

%% Exercises
%

%% 
% *Exercise 1 (Challenge)*
%
% Under the MyFiles folder, open the m-file MyEntropy.m.  In this file,
% create a function that computes the entropy of either a vector or a
% matrix.  Useful Matlab commands are |histc|, |uniq|, |length|, and
% |log2|.
%
% Use the cell below to test your code.

v = round(rand(1,10)*10);
Entropy(v)
MyEntropy(v)

A = round(rand(8,8)*10);
Entropy(A)
MyEntropy(A)

%% Peak Signal-To-Noise Ratio
%
% The DiscreteWavelets command for peak signal-to-noise ratio is |PSNR|.  
% The input values are two matrices of equal size.

img = ImageNames('ImageType','GrayScale','ListThumbnails','True');
A = ImageRead(img{4});
[rows cols] = size(A);

A1 = A + round(rand(rows,cols)*20 - 10);
A2 = A + round(rand(rows,cols)*80 - 40);

ImagePlot(A);
figure;
ImagePlot(A1);
figure;
ImagePlot(A2);

psnr1=PSNR(A,A1);
psnr2=PSNR(A,A2);

str=sprintf('The PSNR of A and A1 is %f and the PSNR of A and A2 is %f.',psnr1,psnr2);
disp(str);

%% Exercises
%

%%
% *Exercise 1*
%
% Write a function to compute the PSNR of two matrices.  Use the m-file
% MyPSNR that can be found in the MyFiles folder.  Useful Matlab commands
% are |sum| and |log10|.
%
% Here is some code to test your module :

mypsnr1 = MyPSNR(A, A1)
mypsnr2 = MyPSNR(A, A2)

%% Huffman Codes
%
% The DiscreteWavelets Toolbox includes two functions that are useful for
% creating and analyzing Huffman codes.
%

%% MakeHuffmanCodes
%
% The DiscreteWavelets function |MakeHuffmanCodes| can be used to generate 
% Huffman codes for a list of integers, a string, or a matrix of integers. 
% Note that integer input must be nonnegative.
%
% The routine returns five pieces of information.  The first is a list of 
% unique integers or characters that appeared in the input, the second is
% list of the relative frequencies of each unique value, and the third is 
% a list of Huffman codes for each unique value.  The fourth output is the
% original bitstream length of the input and the last return value is the
% new bitstream length.
%

[uniq, freq, codes, origlen, newlen] = MakeHuffmanCodes('pitterpatter')

%%
% We can also find the Huffman codes of an image:

img = ImageNames('ImageType', 'GrayScale', 'ListThumbnails', 'True');
A = ImageRead(img{1});
ImagePlot(A);
[uniq, freq, codes, origlen, newlen] = MakeHuffmanCodes(A);
str=sprintf('The original bitstream length is %i and the new bistream length is %i.',origlen,newlen);
disp(str);

%% HuffmanTree
%
% For small input, we can use |HuffmanTree| to visualize the output of
% |MakeHuffmanCodes|.  The input of |HuffmanTree| is the first three 
% outputs of |MakeHuffmanCodes|.  There are several graphics options for 
% |HuffmanTree|.  See Help for more information.

[uniq, freq, codes, origlen, newlen] = MakeHuffmanCodes('tennessee');
figure;
HuffmanTree(uniq,freq,codes);

%% Exercises
%
% *Exercise 1*
%
% Load the surfer image (small version) from the DiscreteWavelets package
% and find the bits per pixel that results from the Huffman coded version
% of the image.

%%
% Place your answer here.  To type an answer, start a line with %.  Matlab
% commands are entered as usual.

%%
%
% *Exercise 2 (Challenge)
%
% Suppose you pass a string |s| to |MakeHuffmanCodes|.  Write a function 
% that will take the codes returned by |MakeHuffmanCodes| and |s| and 
% return the new bitstream.  Call the function |MakeBitstream|.  You will 
% have to add this m-file to the MyFiles folder.
%
% For example:

s = 'boohoo';
figure;
[uniq, freq, codes, origlen, newlen] = MakeHuffmanCodes(s);
HuffmanTree(uniq, freq, codes);

%%
%
% Using the tree above and |s|, we see that the new bitstream for 
% "boohoo" is 01110011.

%%
close all;
displayEndOfDemoMessage(mfilename)