% Mini-Project #2

% Tuneup #7

% Load the image.  Download the image from 
% https://users.ece.utexas.edu/~bevans/courses/signals/homework/echar512.matLinks to an external site.
% and place it in the current directory or a directory on the Matlab path.

load echar512.mat

% The load command will define a Matlab matrix echart.

% Display the image. 

hAxes = axes(figure); 
hImage = imshow(echart, [], 'Parent', hAxes);
title(hAxes, 'Original echart image');

% 2.5

% Apply a first-order difference filter along the rows.
% Display the resulting image.

bdiffh = [1, -1];
yy1 = conv2(echart, bdiffh);
hAxes1 = axes(figure); 
hImage1 = imshow(yy1, [], 'Parent', hAxes1);
title(hAxes1, 'Image filtered with [1 -1] along rows');

% Apply a first-order difference filter along the columns.
% Display the resulting image.

bdiffh = [1, -1];
yy2 = conv2(echart, bdiffh');
hAxes2 = axes(figure); 
hImage2 = imshow(yy2, [], 'Parent', hAxes2);
title(hAxes2, 'Image filtered with [1 -1] along columns');

% Not asked
% Apply a first-order difference filter along the rows
% and columns.
% Display the resulting image.

bdiffh = [1, -1];
yy3 = conv2(echart, bdiffh);
yy4 = conv2(yy3, bdiffh');
hAxes4 = axes(figure); 
hImage4 = imshow(yy4, [], 'Parent', hAxes4);
title(hAxes4, 'Image filtered with [1 -1] along rows and columns');

% 3.2
load echar512.mat
hAxes = axes(figure); 
hImage = imshow(echart, [], 'Parent', hAxes);
title(hAxes, 'Original echart image');

% 3.2.1
% Apply two 2-D FIR filters
% h1[n] has non-zero coefficients of 1 and -q:
q = 0.9;
h1 = [ 1 -q ];
yy5 = conv2(echart, h1);
ech90 = conv2(yy5, h1');
hAxes6 = axes(figure); 
hImage6 = imshow(ech90, [], 'Parent', hAxes6);
title(hAxes6, 'Image filtered with FIR filter #1 along rows and columns');

% h2[n] = r^n for 0 <= n <= M and zero elsewhere:
r = 0.9;
M = 22;
l = 0 : M;
h2 = r .^ l;
yy7 = conv2(ech90, h2);
yy8 = conv2(yy7, h2');
hAxes8 = axes(figure); 
hImage8 = imshow(yy8, [], 'Parent', hAxes8);
title(hAxes8, 'Image filtered with FIR filter #1 and #2 along rows and columns');

% Worst case error in the restored image vs. the echart image
% on a pixel-by-pixel basis.  We have to remove the black border
% in the restored image so the restored images are same size
% and then normalize the number gray levels.
[numrows numcols] = size(echart);
yy8cropped = yy8(1:numrows, 1:numcols);
yy8min = min(min(yy8cropped));
yy8max = max(max(yy8cropped));
echartmin = min(min(echart));
echartmax = max(max(echart));
yy8normalized = (yy8cropped - yy8min) * (echartmax / (yy8max - yy8min));
worstCaseError = max(max(abs(yy8normalized - echart)));

% Image quality measures
% We have to make sure the two images have the same dimensions
% The order of the arguments do not matter
ech90cropped = ech90(1:512, 1:512);
yy8cropped = yy8(1:512, 1:512);

ech90immse = immse(echart, ech90cropped)
yy8immse = immse(echart, yy8cropped)

ech90psnr = psnr(echart, ech90cropped)
yy8psnr = psnr(echart, yy8cropped)

ech90ssim = ssim(echart, ech90cropped)
yy8pssim = ssim(echart, yy8cropped)