Image Restoration using kalman filter

Question

I have been trying to restore an image that was blurred with a known Point Spreading Function and corrupted with noise using a kalman filter. I have looked at theory and have a basic understanding of how kalman filters work. But I can't find any good material on how I can use it to restore an image. I have tried a method where I considered that the image was the state to be predicted with no control signal and a observation model equal to the blur function with a known measurement noise i.e,

$$X(k+1) = X(k)$$ $$Z(k) = H*X(k) + R$$ where,

X(k) is the 2D image
Z(k) is the observed noisy blurred image
H is the PSF
R is the noise

here is the code is used,

I = im2double(imread('E:\Wallpapers\2.jpg'));
I = I(:,:,1);

LEN = 2;
THETA = 5;
PSF = fspecial('motion', LEN, THETA);
blurred = imfilter(I, PSF, 'conv', 'circular');

SNR = 20;
noisy = awgn1(blurred,SNR);

noise = noisy - blurred;
R = cov(noise);
H = fft2(PSF,size(I,1),size(I,2));
K = 1e5*eye(size(I));
P = 1e5*eye(size(I));
X = zeros(size(I));
S = zeros(size(I));

for i = 1 : 1000
    S = R + H*P*H';
    K = P*H'/S;
    X = X + K*(noisy - H*X);
    P = (I - K*H)*P*(I - K*H)' + K*R*K';
end

MSE_noise = sum(sum((noisy-I).^2));
MSE1_restored = sum(sum((X-I).^2));

The image is filled with NaNs which is because of K = P*H'/S being badly scaled. Where am I going wrong? Is there a problem with the code or am I supposed to change the model to restore the image?

Answer 1

Concerning the code of Kalman filter all is OK in your code excepted the update of the estimate covariance noted P where I should be the identity matrix, not the input image.

Answer 2

I suggest you to investigate Image Restoration using deconvolution since it is more appropriate for your application. Kalman filtering is designed to remove or ignore noisy data in a signal but it is not suited for deblurring (as far as I know).

Noting $\epsilon$ the additive noise, $H$ the convolutional mask (your Point Spreading Function), $x$ the original image you want to retrieve and $y$ the blurred image you have:

$y = H*x + \epsilon$

Deconvolution allows you to retrieve $x$ knowing $y$ and $H$.

Using Matlab you can use the following function:

deconvreg(image, PSF)

where PSF is an estimate of the Point Spread filter and image is the blurred and noisy image.

If you want to know more about the theory, this is called an ill-posed problem that you solve by minimizing a Least Square Error through an iterative process:

$$x_{estimated} = \arg\min_{x}||y - Hx||^2$$

Under some approximations this can be computed in a fast way through FFT.

This does not answer your initial question about Kalman Filtering but I hope it will help you! If you are interested I have documents about image restoration using inverse and Wiener filtering.