摘 要：本實驗主要使用維納濾波法（又名為最小均方誤差濾波）實作圖像複原與重建。首先我們通過對一幅圖像加入運動污損濾波和高斯噪聲，然後從噪聲中提取出原始圖像信号。在各種估計方法中，維納濾波是一種最基本的方法，相比起逆濾波，它更适用于從噪聲中分離出有用信号的圖像。
一、技術論證
1 圖像複原與重建
圖像複原技術的主要目的是以預先确定的目标來改善圖像。圖像複原試圖利用退化現象的某種先驗知識來複原被退化的圖像。因而複原技術是面向退化模型的，并且采用相反的過程進行處理，以便恢複出原圖像。
圖1給出了圖像退化/複原過程的模型。圖像複原處理是建立在圖像退化的數學模型基礎上的，這個退化數學模型能夠反映圖像退化的原因。退化過程中可以被模型化為一個退化函數和一個加性噪聲項，處理一幅輸入圖像
産生一幅退化圖像。給定
和關于退化函數
的一些知識以及外加噪聲項
，圖像複原的目的是獲得關于原始圖像的近似估計
。通常我們希望這一估計盡可能接近原始輸入圖像，并且
和
的資訊知道得越多，所得到的
就會越接近
。
如果系統H是一個線性、位置不變性的過程，那麼在空間域中給出的退化圖像可由下式給出：
其中，
是退化函數的空間描述；“*”表示空間卷積。由于空間域上的卷積等同于頻域上的乘積，是以，可以把退化圖像的模型函數寫成等價的頻域下的描述：
公式中各項是式1中相應項的傅裡葉變換。本實驗使用的退化函數為：
該函數表示圖像在X軸和Y軸分别以給定的速度
做勻速直線運動。實驗中設定：
圖1圖像退化/複原過程的模型
2 高斯噪聲
由于高斯噪聲在空間和頻域中數學上的易處理性，這種噪聲(也稱為正态噪聲)模型經常被用于實踐中。事實上，這種易處理性非常友善，使高斯模型經常用于臨界情況下。高斯随機變量z的PDF由下式給出：
圖2高斯函數的波形
3 維納濾波
維納濾波綜合了退化函數和噪聲統計特性兩個方面進行複原處理，該方法是建立在認為圖像和噪聲是随機過程的基礎上，目标是找一個未污染圖像
的估計值
，使它們之間的均方誤差最小。誤差度量由下式給出：
其中
是參數的期望值。這裡假設噪聲和圖像不相關，其中一個或另一個有零均值，且估計的灰階級是退化圖像灰階級的線性函數。在這些條件下，
中誤差函數的最小值在頻域用下清單達式計算：
這裡，我們應用了這樣一個事實：一個複數量與它的共轭的乘積等于該複數量幅度的平方。這個結果就是維納濾波，是N.Wien
er[1942]首次提出的概念。由方括号裡邊的項所組成的濾波器通常還叫做最小均方誤差濾波器，或最小二乘方誤差濾波器。式中各項的意義如下：
二、實驗結果讨論
本實驗使用岡薩雷斯主編數字圖像處理第三版中的素材Fig5.07(a).jpg進行+45度方向，T=1的污損濾波，并對污損圖像進行高斯噪聲加噪處理。其次，驗證維納濾波法的有效性。從實驗結果可以看出維納濾波對于噪聲的複原效果是很好的。
編寫函數gaussian_noise(ima,a,b)産生高斯噪聲，其中ima為輸入圖像，a，b分别是高斯噪聲的均值和方差，函數中主要調用了函數randn()，該函數的效果是産生正态分布的随機數或矩陣的函數。要産生一個随機分布的指定均值和方差的矩陣，隻需執行以下公式：
式中(m,n)為圖像的尺寸。為了驗證函數效果，圖3給出均值為0，方差分别為20，200，1000的加性噪聲圖像。本實驗過程中要求生成均值為0、方差為10的高斯噪聲。
圖3均值為0，方差分别為20，200，1000的加性高斯噪聲圖像
編寫函數motion(T,a,b)實作圖像模糊化，其中使用的核心公式為：
使用函數進行+45度方向且T=1 的污損濾波，對比MATLAB自帶函數的效果，如圖4所示，有巨大差異，但圖像确實在x和y軸都有同幅度的位移，事實證明隻有核心公式還無法做到示例的輸出效果。
圖4兩種函數的運動模糊輸出圖像
編寫函數實作對運動模糊且存在噪聲的圖像的重建和複原，其中使用的退化函數為：
圖5和圖6分别給出了維納濾波恢複結果，結果顯示維納濾波能有效重建院圖像。而對于圖6這種噪聲污染嚴重的圖像，使用該濾波方法也能得到有效改善。
圖5 圖像經過運動模糊、加噪處理（0均值，10方差）和維納濾波後的效果對比
圖6圖像經過運動模糊、加噪處理（0均值，100方差）和維納濾波後的效果對比
以下給出主要部分的MATLAB代碼：
% 編寫函數給圖像加入高斯噪聲
% ima為輸入圖像，a和b分别為高斯函數的均值和方差
function imt=gaussian_noise(ima,a,b)
% a = 0;
% b = 1000;
ima = imread('Fig5.26(a).jpg');
[m,n] = size(ima);
x = uint8(a + sqrt(b) * randn(m,n)); % 産生一個m×n的高斯噪聲矩陣
imt = ima + x; % 給原圖疊加高斯噪聲
% imshow(imt);
return
% 編寫函數給圖像加入模糊污損噪聲
% ima為輸入圖像，a和b分别為x軸和y軸的運動速度參數，T為圖像位移時間
function imt=motion(ima,T,a,b)
% ima = imread('Fig5.26(a).jpg');
%
% a = 0.1;
% b = 0.1;
% T = 1;
[m,n] = size(ima);
imt = fftshift(fft2((ima)));
for u=1:m
for v=1:n
Y(u,v) = u*a + v*b;
H(u,v) = T * sin(pi*Y(u,v)) * exp(-1j*pi*Y(u,v)) / (pi*Y(u,v));
imk(u,v) = imt(u,v)* H(u,v);
end
end
imt = ifft2(fftshift(imk));
imt = uint8(abs(imt))
return
%維納濾波函數
F0=fftshift(fft2(imt)); % 對運動模糊圖像（未加高斯噪聲）作傅立葉變換
LA = fftshift(fft2(ima)); % 對原始圖像作傅立葉變換
K=0.1; % 加在退化函數的所有項上的特定常數
for u=1:m
for v=1:n
H(u,v) = F0(u,v)/LA(u,v); % 加噪圖像作傅立葉變換除以原始圖像作傅立葉變換，得到退化函數
% 計算維納濾波函數
H0(u,v)=(abs(H(u,v)))^2;
H1(u,v)=H0(u,v)/(H(u,v)*(H0(u,v)+K));
F2(u,v)=H1(u,v)*F0(u,v); % 恢複、計算原始圖像的傅立葉變換估計
end
end
J2=ifft2(fftshift(F2));
J2=uint8(abs(J2));
subplot(2,2,4),imshow(J2,[]),title('維納濾波複原圖');