Introduction to Computer and Human Vision Exercise Set 2

Document Sample
Introduction to Computer and Human Vision Exercise Set 2 Powered By Docstoc
					            Introduction to Computer and Human Vision
                            Exercise Set 2
                Deblurring using Natural Image Priors

                          Submission: individual or in pairs
                       Due date: Thursday, December 31st , 2009

    In this assignment we will investigate the solution of inverse problems using natural
image priors. In particular we will focus on the deblurring problem and explore solutions
using Gaussian and sparse priors.
Brief overview: In class we discussed the problem of deblurring an m × n image y. We
assume that the blurred image y is generated by
                                       y =x k+n
where k is a blur kernel and n ∼ N (0, η 2 ) is zero mean Gaussian noise.
   Introducing a prior on natural images
                                         1                           y
                                                       |gij (x)|α +|gij (x)|α
                              p(x) ∝ e− σ2       i,j

where σ 2 is a parameter of the prior and gij (x) and gij (x) are the responses at coordinate
(i, j) of the image to the filters (−1, 1) and (−1, 1)T respectively. We then formulated the
problem of deblurring using a prior on the distribution of the edges in the image as finding
a minimizer of
                         1                 1                              y
                L(x) = 2 y − Ax 2 + 22
                                                           |gij (x)|α + |gij (x)|α
                        2η                2σ

Where A is the mn × mn blur matrix corresponding to k, y is the column-stacked blurred
image and η 2 is the variance of the noise in the blurred image. Setting α = 2 we get the
Gaussian prior and using α = 0.8 we get the sparse prior.
   The “iteratively reweighted least squares” algorithm is an iterative procedure to find a
minimizer of a function of the form
                                  E(x) =           ρi (|bi x − zi |)

Given an estimate of a minimizer xt define an iteration by
                              xt+1 = arg min               q it (|bi x − zi |)

Where q it (c) is a quadratic (scaler) function, which coincides with ρi at ct = bi xt − zi and
which is strictly greater than ρi at all other points. From this we can derive three conditions
on q it
  1. q it (ct ) = ρi (ct )
  2. q it (ct ) = ρi (ct )
  3. q it is symmetric around the y-axis (thus q it (c) = ac2 + d for some constants a, d).
These conditions determine q it . Since all we care about is finding the minimizer of q it we
can neglect the additive term d and use
                                                           ρ (ct ) 2
                                           q it (c) =             c
This is a Matlab programming assignment.
   • It is very important that you thoroughly document your code. In particular each
     function should begin with a short description of what it does and what are the input
     and output.
   • Remember that Matlab is efficient when the code is vectorized, so try to avoid using
     loops. Specifically linear systems should be formulated in matrix notation and solved
     without using loops (looping over IRLS iterations is fine).
   • You will most probably want to use the following Matlab functions: help, imread,
     imwrite, imshow, figure, conv2, reshape, spdiags...
   • In this exercise you will have to use Matlab’s sparse matrix representation. If you are
     not familiar with sparse matrices please consult the Matlab premier (you can find a
     link on the “links” page of the course website). The important thing to note here is
     that if you don’t keep all your matrices sparse you will run out of memory.
   • You can use the function getConvMat.m which is posted on the website to generate
     the convolution matrices.
   • Save each matlab function in a separate file.
1. Write a function
function x_star = solveGaussianPrior(A,y,eta,sigma)
which computes x_star the minimizer of
                              1            2        1                                  y
                   L(x) =         y − Ax   2   +                          x
                                                                        |gij (x)|2 + |gij (x)|2
                             2η 2                  2σ 2

Where A is an mn × mn matrix (which could correspond to deblurring super resolution or
any linear operator), y is an m × n image and eta and sigma are scalars as defined above.
The output is of size m × n. The optimization should be done in the spatial domain, solving
a linear system of equations using Matlab’s backslash.
2. Write a function

function x_hat = solveSparsePrior(A,y,eta,sigma,x0,n_iteration)
which computes x_hat an (approximate) minimizer of
                         1             2        1                                    y
               L(x) =        y − Ax    2   +                          x
                                                                    |gij (x)|0.8 + |gij (x)|0.8
                        2η 2                   2σ 2

by running n_iteration iterations of “iteratively reweighted least squares” with an initial
guess of x0 (also of size m × n). Thus the iteration should be
                                    1 1t                    1                            y
    xt+1 = arg min                    2
                                        q (|yij − aij x|) + 2 q 2t (|gij (x)|) + q 2t (|gij (x)|)
                                   2η                      2σ

where q 1t (c) = c2 , q 2t (c) = |ct |−1.2 c2 and aij is the ij’th row of A. Or equivalently
                           1               2        1                                     y y
         xt+1 = arg min        y − Ax      2   +                         x x
                                                                        wij |gij (x)|2 + wij |gij (x)|2
                          2η 2                     2σ 2

                                 y      y
        x     x
where wij = |gij (xt )|−1.2 and wij = |gij (xt )|−1.2 . To avoid division by zero we modify the
weights to
                                    x             x
                                  wij = max(|gij (xt )|, )−1.2
                                     y          y
                                    wij = max(|gij (xt )|, )−1.2
In this exercise we use = 1e − 5 and for the initial estimate x0 we use the solution of
the Gaussian prior minimization. Technical note: assume that the input to the function
is σ 2 = 2obsvar (which is the correct estimate for the Gaussian case). Use this value
to compute x0 . Then to compensate for the inaccuracy of this estimation for the non-
Gaussian case multiply σ by 5 and proceed with the IRLS iterations. In all the questions
in this exercise we will perform 3 iterations of IRLS for the sparse deblurring.
 3.     Download the file ex2_q3.mat, this file contains three blur kernels ker_1,ker_2
and ker_3. It also contains 4 versions of the same image. im_1 the original image and
y1 y2 y3 image im_1 blurred with the corresponding kernels plus noise. Use the clean
image to estimate σ 2 . A good estimate is 2obsvar where obsvar is the mean of the squared
filter responses (gij (xt ))2 and (gij (xt ))2 . You can also use the clean image to estimate the
noise η 2 .

 (a) Estimate the noise for each one of the noisy images and the prior parameter σ 2 . Include
      the estimated values in your report along with a brief explanation of how you obtained
      the noise estimators.

 (b) Deblur the blurry images using a Gaussian image prior and then a sparse image prior.
     For each blur kernel j = 1, 2, 3 generate two images im_1_deblur_ker_j_Gauss and
     im_1_deblur_ker_j_sparse where j is and index which should be substituted with
     the value it assumes. Also generate a (3 × 3) image array of the following form.
     Row 1 consists of the three blurred images. Row 2 shows the results of the Gaussian
     deblurring and row 3 the results of sparse deblurring.
      Save this image array as a .png image called ex2_q3_b.png.

 (c) Investigate the effect of under and over estimating the noise parameter on debluring.
      Use y2 to generate a (1 × 3) image array of the following form. The left most image
      shows an example of deblurring when the noise is underestimated (multiply η by 0.2).
      The central image is deblurring with the correct noise estimation and the rightmost
      should show a result of deblurring with overestimation (multiply η by 10).
     Save this image array as a .png image called ex2_q3_c.png. Also include a brief
     discussion which explains your results and why you think they turn out as they do.

 (d) Investigate the effect of using different blur kernels on debluring. Generate a (2 × 2)
      image array of the following form. The left most image shows deblurring of y1 with
      ker_1. The next image shows the results of deblurring y1 with ker_3. The next row
      is deblurring of y3 with ker_1 and y3 with ker_3
     Save this image array as a .png image called ex2_q3_d.png. Also include a brief
     discussion which explains your results and why you think they turn out as they do.

 4.    Download the file ex2_q4.mat, this file contains a blur kernel and a blurry-noisy
image called ker and y respectively. Deblur using a Gaussian prior and a sparse prior. Use
η = 0.01 and σ = 0.09. Save the results as im_2_deblur_sparse and im_2_deblur_Gauss
Also save a (1×3) image array (blurry, Gaussian, sparse) as a .png image called ex2_q4.png.
    In this exercise we perform the debluring in the spatial domain, thus the computations
are pretty demanding. In particular you will need more memory than is usually available
on a home PC. Thus in this question you will have to perform the deblurring in overlapping
blocks which should be stitched together to form the complete image. There should be
an overlap so as to avoid boundary problems within the image (of course you will still get
boundary effects on the boundary of the stitched image). Another thing to note is that it
is important to give the function solveSparsePrior an initial guess which comes from the
final stitched image of the Gaussian deblurring part.
Submission instructions: Submit one .zip archive called by sending an email to with the subject “CVfall2009 ex2”. This zip archive
should contain

  1. A discussion of the results in pdf, ps or rtf format. The first page should include the
     names of the submitters and your ID numbers.

  2. A .mat file called ex2.mat which includes the deblurred images from question 3
     (im_1_deblur_ker_j_Gauss and im_1_deblur_ker_j_sparse for j = 1, 2, 3) and
     from question 4 (im_2_deblur_sparse and im_2_deblur_Gauss).

  3. all the .png images namely: ex2_q3_b.png ex2_q3_c.png ex2_q3_d.png ex2_q4.png.

  4. Two .m files implementing the functions from questions 1 and 2. and the scripts /
     functions which you used to answer questions 3 and 4.

  5. Please submit a hard copy of your code and report to the “computer vision” mailbox.