# Estimate the Image Using Multi Many Realizations of Its Convolution with a Known Filters Using Wiener Filter

Suppose we have a corrupted image $$Y = H*X + \epsilon$$ that is formed by taking an image $$X$$, convolving it with a point-spread function $$H$$, and adding gaussian noise $$\epsilon$$. Then we know that the Wiener Filter can compute the MMSE estimator of $$X$$ given $$H$$ and the signal-to-noise ratio (SNR).

Given a set of $$n$$ images $$Y_i = H_i * X + \epsilon_i$$, is there a generalized Wiener Filter estimate of $$X$$ given the $$H_i$$'s and the SNRs?

• Hi: You may be better off sending that to dsp.stackexchange since the people on that list will be more familiar with the notation, terminology and the weiner filter. Of course, that's not say that somewhere here can't answer. I'm not sure if you're supposed to cross-post ? I don't think so though. Aug 14 '21 at 22:07

This is a nice question.

The math is actually pretty simple once you embrace the method I derived the Wiener Filter in - How Is the Formula for the Wiener Deconvolution Derived?

So, here is the model:

$$\boldsymbol{y}_{i} = \boldsymbol{h}_{i} \ast \boldsymbol{x} + \boldsymbol{w}, \; i = 1, 2, \ldots, n$$

Where $$\boldsymbol{w}$$ is an additive white gaussian noise which is independent of the signal and $$\boldsymbol{x} \sim N \left( 0, {\sigma}_{x}^{2} \right)$$.

Then the optimization model becomes (The MMSE Estimator):

$$\arg \min_{\boldsymbol{x}} \frac{1}{ {\sigma}_{n}^{2} } \sum_{i = 1}^{n} {\left\| H \boldsymbol{x} - \boldsymbol{y} \right\|}_{2}^{2} + \frac{1}{2 {\sigma}_{x}^{2}} {\left\| \boldsymbol{x} \right\|}_{2}^{2}$$

Where $$H = \begin{bmatrix} {H}_{1} \\ {H}_{2} \\ \vdots \\ {H}_{n} \end{bmatrix}$$ is the model matrix where $${H}_{i}$$ is the matrix form of the $$i$$ -th filter and $$\boldsymbol{y} = \begin{bmatrix} \boldsymbol{y}_{1} \\ \boldsymbol{y}_{2} \\ \vdots \\ \boldsymbol{y}_{n} \end{bmatrix}$$ is a concatenation of the output images in a column stack form.

Then the solution is given by:

$$\hat{\boldsymbol{x}} = {\left( {H}^{T} H + \frac{ {\sigma}_{n}^{2} }{ {\sigma}_{x}^{2} } I \right)}^{-1} {H}^{T} \boldsymbol{y}$$