ECE563 Image Processing
Chapter 6 Project
Wiener Filter Image Restoration Implemented with FFTs

1. The problem

During image acquisition images are invariable corrupted by blurring processes and
noise. Here we will explore the restoration of an image corrupted by the degradations
modeled in the midterm exam project. In particular, we will consider diffraction limited
blur, detector integration effects, and atmospheric turbulence. Use the exact same
parameters from the midterm project for the degradation here. Assume a single Fried
parameter of 0 0.05r  meters. In addition to that blurring, here you are to also model
additive Gaussian noise with standard deviation of 1 digital unit (DU). Use the same
“ideal” image as you used for the midterm project.

In this project, you are to create a custom user-defined function in MATLAB called
'wiener_filter_2d.m' that will restore an image like our observed image given some
relevant information. This function should work for the current project, but be general
purpose, so it can be used for any similar restoration problem where the blurring PSF is
known. The Wiener filter should be the FFT-based IIR Wiener filter discussed in class
with the simple constant noise-to-signal ratio model. The Wiener filter design and image
processing are to be done using FFTs. Apply the DC gain adjustment as shown in the
examples in class so that the DC gain of the Wiener filter is 1.

Your Wiener filter function should follow the function template provided at the end of
this project statement. You may include additional optional inputs and/or outputs if you
choose. In addition to the Wiener filter function, prepare a script file called
'run_wiener_filter_2d.m'. The script file will load the input image, generate the
artificially degraded image as described above (blur and noise). Next, the script file
should call your new restoration function 'wiener_filter_2d.m' to generate a restored
version of the image. Use 'dif.m' to quantitatively compare the ideal image with the
restored image (or sue your own code). Also, compute the SSIM metric using the built-in
ssim() function. Also, compute the error for the observed image relative to the ideal
image for reference. Make note of the mean absolute error (MAE), mean squared error
(MSE), and structural similarity image metric (SSIM) for all of your results.

A tuning parameter required by the Wiener filter is the noise-to-signal (NSR) ratio. In
practice this will not be known in advance. It is usually selected in a real application
based on subjective evaluation of the restored image using trial and error. In this case,
however, we know the ideal image so we can explore the sensitivity of the Wiener filter
to this parameter by quantitatively by computing the restoration error as a function of the
NSR used in the Wiener filter. Put your Wiener filter in a loop (yes, it is OK to use a
loop here) and evaluate the MSE for a suitable range of NSR values. What is a suitable
range you ask? That is for you to determine. Find the NSR that provides approximately
the lowest MSE and show the restored image for that one NSR. Graph the error on the
vertical axis versus the NSR found using your loop, showing where the minimum lies.
Use semilogy instead of plot to make the range of MSE values easier to display. As
always, be sure to fully title and label all images and plots and include correct units.
Compute and display the MAE, MSE and SSIM for a restored image using this optimum
Wiener NSR. You do not need to optimize the NSR for each metric (use MSE only for
choosing NSR).

2. What to turn in

Once complete, publish your script file run_wiener_filter_2d.m, to a .pdf file. Upload the
output .pdf file and your wiener_filter_2d.m function file. The script file should generate
the outputs enumerated below in order, properly titled, and with all axes fully labeled.
Each plot must be on its own figure so they can all be inspected after execution (no
overwriting figures). Your script file should:

1. Show the ideal input image.

2. Show a mesh plot of the degradation frequency response as 1 2| ( , ) |H f f . Label
axes with units of cycles/sample. Plot this from -0.5 to 0.5 cyc/sample on each
axis. This can be accomplished using the midterm project code (or you can use
different code).

3. Show the artificially degraded image with blur and noise. The degraded image
can be generated using the midterm project code (or you can use different code).

4. Graph the MSE versus the NSR used by your Wiener filter function. Use
semilogy and fully and properly label the graph.

5. Show a mesh plot of the magnitude Wiener filter frequency response 1 2| ( , ) |WH f f
similar to that above in Item 2. Show only the Wiener filter frequency response
corresponding to the filter that gives the lowest MSE found in Part 4. Use abs()
to get the magnitude and be sure to zero-center your frequency response samples.

6. Show the Wiener filtered image (restored). Show only the Wiener filter output
corresponding to the filter that gives the lowest error found in Part 4.

7. Display the MAE, MSE, and SSIM for the Wiener filter with NSR that gives the
lowest error found in Part 4 (you don’t need to optimize for each metric). You
can use fprintf( ) to make a nice output statement for each error metric.

Your plots should be clearly and correctly labeled and presented in order. Make sure all
plots are have properly labeled axes with units and appropriate titles. Your MATLAB
code should be well commented and clearly presented. You must use cells to organize
your script file. A cell is a block of code separated using %%. Experiment for how the
comments come out after publishing and make it look nice. Your function should be
formatted like the template at the end of this document.

Look at the code provided in chp6/run_wiener_filter.m and the slides in
chp6/Wiener_filters.pdf for helpful information.

You do not need to prepare a written report (just submit the published script file and your
function file wiener_filter_2d.m). This project is to be done independently and the
submitted code must be uniquely your own.

function out = wiener_filter_2d( in, psf, NSR, pad )
% out = wiener_filter_2d( in, psf, NSR, pad )
% 2D IIR Wiener filter implemented with FFTs. Uses the constant noise to
% signal PSD ratio approximation.
% Prepared for ECE 563 Image Processing with Dr. Hardie
% in Input image (2D array)
% psf Degradation point spread function
% NSR Noise to signal PSD ratio (constant)
% pad Padding size for the FFT processing to minimize border effects
% out.img Output image (restored), same size as in (and aligned)
% out.HW Wiener filter frequency response
% out.f1 Frequency axis in x for out.HW in cyc/sample
% out.f2 Frequency axis in y for out.HW in cyc/sample
% Author: ----------------
% Date: ----------------
% COPYRIGHT © 2021 ------------. ALL RIGHTS RESERVED.
% This is UNPUBLISHED PROPRIETARY SOURCE CODE of -------------; the
% contents of this file may not be disclosed to third parties, copied or
% duplicated in any form, in whole or in part, without the prior written
% permission of --------------.