# Deconvolution Question on Article "Deriving Intrinsic Images from Image Sequences" by Yair Weiss

there are n derivative filters: $f_i$, and denote $f_i^r$ as $f_i$'s reverse filter such that

$$f_i(x,y)=f_i^r(-x, -y)$$

$r_i, f_i$ given, to find $r$ from the equations: $$f_i * r = r_i, (1 \leq i \leq n)$$

Professor Weiss recovers $r$ using the following method: $$r = g * \left(\sum_{i=1}^n f_i^r*r_i\right)$$ and $g$ satisfies: $$\delta = g * \left(\sum_{i=1}^n f_i^r*f_i\right)$$

1. Why can the method solve the equation? What books should I refer to?
2. I want to know how to find the $g$?

Professor Weiss uses these two derivative filters: $[0\ 1 -1]$ and $[0; 1; -1]$. and in his code the function invDel() returns $g$, which I could not understand. Article link, implementation link.

Professor Weiss's MATLAB code for getting $g$:

function [invK]=invDel2(isize)    % isize is 2 * max(imgWidth, imgHeight)
K=zeros(isize);
K(isize/2,isize/2)=-4;
K(isize/2+1,isize/2)=1;
K(isize/2,isize/2+1)=1;
K(isize/2-1,isize/2)=1;
K(isize/2,isize/2-1)=1;

Khat=fft2(K);
I=find(Khat==0);
Khat(I)=1;
invKhat=1./Khat;
invKhat(I)=0;
invK=ifft2(invKhat);
invK=-real(invK);
invK=conv2(invK,[1 0 0;0 0 0;0 0 0],'same');% shift by one

• Any chance you review my answer?
– Royi
Aug 5 at 14:36

I didn't have time to fully read the article yet I'll try to answer according to your question.

If you take the first equation and convolove it on the left by ${f}_{i}$ and use the second identity you can see this specific $g$ holds the equation.

Regarding how to get this $g$, I think he uses "Deconvolution" in the Fourier domain.

I'm struggling a little with this notation because I don't usually do any significant image processing. Many LTI filters have an inverse and it seems like this is where your discussion is leading.

I'm assuming your image filter with three values could be expressed something like this: $$h[n] = h_0\delta[n+1] + h_1\delta[n] + h_2\delta[n-1]$$

If so, getting the inverse is pretty simple. You start with a Z-transform, take the inverse, and convert it back to an impulse reponse.

\begin{align} H(z) &= h_0z + h_1 + h_2z^{-1}\\ H_i(z) &= \frac{1}{H(z)} = \frac{1}{h_0z + h_1 + h_2z^{-1}}\\ &= \frac{1}{z(h_0 + h_1z^{-1} + h_2z^{-2})}\\ &= \frac{z^{-1}}{(1-\alpha z^{-1})(1-\beta z^{-1})}, \textrm{ROC}^{\dagger}\\ \alpha &= \frac{-h_1+\sqrt{h_1^2-4h_0h_2}}{2h_0}\\ \beta &= \frac{-h_1+\sqrt{h_1^2-4h_0h_2}}{2h_0}\\ \end{align} You can use a Z-transform table to convert back. The Region Of Convergence needs to be such that it follows the available patterns for a possible ROC: $$|z|<\operatorname{min}(|\alpha|, |\beta|)$$ or $$|z|>\operatorname{max}(|\alpha|, |\beta|)$$ or $$|\alpha|<|z|<|\beta|$$ or $$|\beta|<|z|<|\alpha|$$ and is stable ($$|z|=1$$ must be in the ROC). You convert back using one form if $$z$$ is greater than the pole magnitude, the other if it is less. I don't think the 1-sample delay in the numerator is actually needed and can probably be removed.