# Minimization of The Blind Deconvolution Functional

I want to minimize the functional of the Blind Deconvolution model as given in: Total Variation Blind Deconvolution by Chan and Wong.

Their model is given by:

$$z = h \ast u + \eta$$

Where $$\ast$$ is the convolution operator, $$h$$ is the blurring kernel, $$u$$ is the sharp noiseless image and $$\eta$$ is additive white gaussian noise (AWGN).

The functional to be minimized is given by (Assuming the Blurring Kernel is known):

$$\underset{u}{\min} f \left( u \right ) = \underset{u}{\min} \left \{ \frac{1}{2} {\left \| h \ast u - z \right \|}^{2}_{L_{2}} + \alpha \int \left | \nabla u \right | dx dy \right \}$$

Where $$\alpha$$ is the smoothing term.

Yet, in Blind Deconvolution the Kernel isn't known and the minimization functional is given by:

$$\underset{u}{\min} f \left( u \right ) = \underset{u}{\min} \left \{ \frac{1}{2} {\left \| h \ast u - z \right \|}^{2}_{L_{2}} + {\alpha}_{1} \int \left | \nabla u \right | dx dy + {\alpha}_{2} \int \left | \nabla h \right | dx dy \right \}$$

Now, The Euler Lagrange equations I calculated and given in the article are:

\begin{align} \frac{\delta L}{\delta h} & = \left ( u \ast h - z \right ) \ast u \left( -x, -y \right ) - {\alpha}_{2} \nabla \cdot \left( \frac{\nabla h}{\left | \nabla h \right |} \right) \\ \frac{\delta L}{\delta u} & = \left ( u \ast h - z \right ) \ast h \left( -x, -y \right ) - {\alpha}_{1} \nabla \cdot \left( \frac{\nabla u}{\left | \nabla u \right |} \right) \end{align}

What I'm not sure about is how can I use it to solve the problem.
At the article they suggest the alternating method, namely once solve for $$h$$ and then for $$u$$.

Yet I don't see how to write in in MATLAB code (Or any other pseudo code).
It should be some kind of a Gradient Descent step.

The blind deconvolution should be a minimization on both $u$ and $h$:

$$\underset{u,h}{\min} f \left( u, h \right ) = \underset{u,h}{\min} \left \{ \frac{1}{2} {\left \| h \ast u - z \right \|}^{2}_{L_{2}} + {\alpha}_{1} \int \left | \nabla u \right | dx dy + {\alpha}_{2} \int \left | \nabla h \right | dx dy \right \}$$

The standard way to do this is to alternately hold one variable fixed and minimize on the other, ie.

\begin{cases} u^{k+1} = \underset{u}{\text{argmin}} f(u,h^k) \\ h^{k+1} = \underset{h}{\text{argmin}} f(u^k,h) \\ \end{cases}

Back when Chan and Wong wrote their paper the way to minimize total variation was to write out the Euler-Lagrange equations which leads to those complicated and slow to solve PDEs. Now there are faster and easier ways to solve each of those subproblems (eg Split-Bregman) which will simplify things for you.

A bit of searching led to me to this 2014 paper: http://www.cvg.unibe.ch/dperrone/tvdb/ which has more information as well as publicly-available MATLAB code. I'm not familiar with their work but it's not easy to publish in CVPR so it's probably worth a look...

• Thank You for your answer. 2 Questions: 1. What if I want to implement the straight way, using Gradient Descent, I don't understand how create the corrects dimensions of the signals (At least in the iteration solving for the Blur Kernel). 2. Can I use MATLAB's fminunc function for the minimization? Thank You.
– Royi
Dec 22, 2014 at 22:45
• Since you have a TV-penalty on h gradient descent won't work since the TV functional is not differentiable. There is some research on smoothing out the TV-term using a parameter (usually called epsilon) that goes under a square root somewhere but the newer Bregman stuff is easier to implement. I don't know how fminunc works but I suspect it won't be a good choice for this difficult problem. After a first pass over the linked CVPR paper I suggest taking a look at their research and code. Dec 22, 2014 at 23:43
• Hi, I looked at their paper and they do use Gradient Descent (Am I wrong?). They use it alternatively (Holding the Kernel, calculating the Image and vice versa).
– Royi
Dec 23, 2014 at 6:28
• Huh, yeah I see it now in alg. 1. In the following paragraph they remark that they only do one iteration - I don't know how that changes things. Also, they have alpha_2 = 0 which simplifies things a bit. In any event, computing the div*(nabla/|nabla|) (the curvature of a level curve) isn't easy. The linked paper has a function gradTVcc(f,epsilon) in their code that handles the computation. Dec 23, 2014 at 21:32