【图像复原方法综述】典型的图像复原方法有哪些

图像复原方法综述

1、摘要

图像是人类视觉的基础，给人具体而直观的作用。图像的数字化包括取样和量化两个步骤。数字图像处理就是将图像信号转换成数字格式，并利用计算机进行加工和处理的过程。

图像复原是图像处理中的一个重要问题，对于改善图像质量具有重要的意义。解决该问题的关键是对图像的退化过程建立相应的数学模型，然后通过求解该逆问题获得图像的复原模型并对原始图像进行合理估计。

本文主要介绍了图像退化的原因、图像复原技术的分类和目前常用的几种图像复原方法，详细的介绍了维纳滤波、正则滤波、LR 算法和盲区卷积，并通过实验证明了该方法的可行性和有效性。

关键词：图像退化、图像复原、维纳滤波、正则滤波、LR 算法、盲区卷积、

2、图像复原概述

在图像的获取、传输以及保存过程中，由于各种因素，如大气的湍流效应、摄像设备中光学系统的衍射、传感器特性的非线性、光学系统的像差、成像设备与物体之间的相对运动、感光胶卷的非线性及胶片颗粒噪声以及电视摄像扫描的非线性等所引起的几何失真，都难免会造成图像的畸变和失真。通常，称由于这些因素引起的质量下降为图像退化。

图像退化的典型表现是图像出现模糊、失真，出现附加噪声等。由于图像的退化，在图像接受端显示的图像已不再是传输的原始图像，图像效果明显变差。为此，必须对退化的图像进行处理，才能恢复出真实的原始图像，这一过程就称为图像复原[1]。

图像复原技术是图像处理领域中一类非常重要的处理技术，与图像增强等其他基本图像处理技术类似，也是以获取视觉质量某种程度的改善为目的，所不同的是图像复原过程实际上是一个估计过程，需要根据某些特定的图像退化模型，对退化图像进行复原。简言之，图像复原的处理过程就是对退化图像品质的提升，并通过图像品质的提升来达到图像在视觉上的改善。

由于引起图像退化的因素众多，且性质各不相同，目前没有统一的复原方法，众多研究人员根据不同的应用物理环境，采用了不同的退化模型、处理技巧和估计准则，从而得到了不同的复原方法。

图像复原算法是整个技术的核心部分。目前，国内在这方面的研究才刚刚起步，而国外

却已经取得了较好的成果。早期的图像复原是利用光学的方法对失真的观测图像进行校正，而数字图像复原技术最早则是从对天文观测图像的后期处理中逐步发展起来的。其中一个成功例子是NASA 的喷气推进实验室在1964年用计算机处理有关月球的照片。照片是在空间飞行器上用电视摄像机拍摄的，图像的复原包括消除干扰和噪声，校正几何失真和对比度损失以及反卷积。另一个典型的例子是对肯尼迪遇刺事件现场照片的处理。由于事发突然，照片是在相机移动过程中拍摄的，图像复原的主要目的就是消除移动造成的失真[2]。

早期的复原方法有：非邻域滤波法，最近邻域滤波法以及效果较好的维纳滤波和最小二乘滤波等。随着数字信号处理和图像处理的发展，新的复原算法不断出现，在应用中可以根据具体情况加以选择。

目前国内外图像复原技术的研究和应用主要集中于诸如空间探索、天文观测、物质研究、遥感遥测、军事科学、生物科学、医学影象、交通监控、刑事侦察等领域。如生物方面，主要是用于生物活体细胞内部组织的三维再现和重构，通过复原荧光显微镜所采集的细胞内部逐层切片图，来重现细胞内部构成；医学方面，如对肿瘤周围组织进行显微观察，以获取肿瘤安全切缘与癌肿原发部位之间关系的定量数据；天文方面，如采用迭代盲反卷积进行气动光学效应图像复原研究等。

3、图像退化模型

图像复原问题的有效性关键之一取决于描述图像退化过程模型的精确性。要建立图像的退化模型，则首先必须了解、分析图像退化的机理并用数学模型表现出来。在实际的图像处理过程中，图像均需以数字离散函数表示，所以必须将退化模型离散化[3]。

对于退化图像g (x , y ) ：

+∞+∞

g (x , y ) =

-∞-∞⎰⎰f (α, β) h (x -α, y -β) d αd β+n (x , y ) (1)

如果上式中f ，h ，n ，g 按相同间隔采样，产生相应的阵列[f (i , j ) ]AB 、[h (i , j ) ]CD 、[n (i , j ) ]AB 、[g (i , j ) ]AB ，然后将这些阵列补零增广得到大小为M ⨯N 的周期延拓阵列，为了避免重叠误差，这里M ≥A +C -1，N ≥B +D -1。由此，当k=0,1,L,M-1;l=0,1,L,N-1时，即可得到二维离散退化模型形式：

M -1N -1

g e (k , l ) =∑∑f (i , j ) h (k -i , l -j ) +n (k , l ) (2) e e e

i =0j =0

如果用矩阵表示上式，则可写为：

g =Hf +n (3)

其中，f , g , n 为一个行堆叠形成的MN ⨯1列向量，H 为MN ⨯MN 阶的块循环矩阵。 现实中造成图像降质的种类很多，常见的图像退化模型及点扩展函数有如下情景[15]：

(1) 线性移动降质

在拍照时，成像系统与目标之间有相对直线移动会造成图像的降质。水平方向线性移动可以用以下降质函数来描述：

⎧1⎪h (m , n ) =⎨d ⎪⎩0若0≤m ≤d and n =0 (4) 其他

式中，d 是降质函数的长度。在应用中如果线性移动降质函数不在水平方向，则可类似地定义移动降质函数。

(2) 散焦降质

当镜头散焦时，光学系统造成的图像降质相应的点扩展函数是一个均匀分布的圆形光斑。此时，降质函数可表示为：

⎧1⎪h (m , n ) =⎨πR 2

⎪⎩0若m 2+n 2=R 2其他 (5)

式中，R 是散焦半径。

(3) 高斯(Gauss)降质

Gauss 降质函数是许多光学测量系统和成像系统最常见的降质函数。对于这些系统，决定系统点扩展函数的因素比较多。众多因素综合的结果总是使点扩展函数趋于Gauss 型。典型的系统可以举出光学相机和CCD 摄像机、γ相机、CT 相机、成像雷达、显微光学系统等。Gauss 降质函数可以表达为：

⎧K exp[-α(m 2+n 2)]若(m , n ) ∈C (6) h (m , n ) =⎨0其他⎩

式中，K 是归一化常数，α是一个正常数，C 是h (m , n ) 的圆形支持域。

4、几种较经典的复原方法介绍

图像复原算法有线性和非线性两类。线性算法通过对图像进行逆滤波来实现反卷积，这类方法方便快捷，无需循环或迭代，直接可以得到反卷积结果，然而，它有一些局限性，

比如无法保证图像的非负性。而非线性方法通过连续的迭代过程不断提高复原质量，直到满足预先设定的终止条件，结果往往令人满意。但是迭代程序导致计算量很大，图像复原时耗较长，有时甚至需要几个小时。所以实际应用中还需要对两种处理方法综合考虑，进行选择[4]。

(1) 维纳滤波法

维纳滤波法是由Wiener 首先提出的，应用于一维信号处理，取得了很好的效果。之后，维纳滤波法被用于二维信号处理，也取得了不错的效果，尤其在图像复原领域，由于维纳滤波计算量小，复原效果好，从而得到了广泛的应用和发展。

维纳滤波器寻找一个使统计误差函数

e =E {(f -f ) } (7) 最小的估计f 。E 是期望值操作符，f 是未退化的图像。该表达式在频域可表示为 ∧2∧2

H (u , v ) 1]G (u , v ) (8) F (u , v ) =[H (u , v ) H (u , v ) 2+S η(u , v ) /S η(u , v ) ∧2

其中，

H (u , v ) 表示退化函数

H (u , v ) =H *(u , v ) H (u , v ) 2

H *(u , v ) 表示H (u , v ) 的复共轭

S η(u , v ) =N (u , v ) 表示噪声的功率谱

S f (u , v ) =F (u , v ) 表示未退化图像的功率谱

比率S η(u , v ) /S η(u , v ) 称为信噪功率比。在IPT 中维纳滤波使用函数deconvwnr 来实现的。 模拟实验结果如下：

22

(2) 正则滤波法

另一个容易实现线性复原的方法称为约束的最小二乘方滤波，在IPT 中称为正则滤波，并且通过函数deconvreg 来实现。

在最小二乘复原处理中，常常需要附加某种约束条件。例如令Q 为f 的线性算子，那么最小二乘方复原的问题可以看成使形式为Q f

∧2∧2的函数，服从约束条件

g -H f =n 的最小化问题，这种有附加条件的极值问题可以用拉格朗日乘数法2

来处理。

寻找一个f ，使下述准则函数为最小： ∧

W (f ) =Q f ∧∧2

+λg -H f ∧2

-n (9) 2

式中λ叫拉格朗日系数。通过指定不同的Q ，可以得到不同的复原目标。

模拟实验结果如下：

(3)Lucy-Richardson算法

L-R 算法是一种迭代非线性复原算法，它是从最大似然公式印出来的，图像用泊松分布加以模型化的。当下面这个迭代收敛时模型的最大似然函数就可以得到一个令人满意的方程：

f k +1(x , y ) =f k (x , y )[h (-x , -y ) *

∧∧∧g (x , y ) h (x , y ) *f k (x , y ) ∧] (10)

*代表卷积，f 代表未退化图像的估计，g 和h 和以前定义一样。

在IPT 中，L-R 算法由名为deconvlucy 的函数完成的。

模拟实验结果如下：

(4)盲去卷积

在图像复原过程中，最困难的问题之一是，如何获得PSF 的恰当估计。那些不以PSF 为基础的图像复原方法统称为盲区卷积。

它以MLE 为基础的，即一种用被随机噪声所干扰的量进行估计的最优化策略。工具箱通过函数deconvblind 来执行盲区卷积。

模拟实验结果如下：

图像分别迭代5次，10次，20次的结果。

5、参考文献

[1] 罗军辉冯平等。MA TLAB7.0在图像处理中的应用[M]，北京：机械工业出版社 2006，202~204，257~278。

[2] 刘维一，于德月，王肇圻等。用迭代法消除数字图像放大后的模糊[J]。光电子.

激光，2002，13（4）：398~400。

[3] 陶洪。数字共焦显微技术及其图像复原算法研究四川大学硕士学位论文 CNKI::CDMD:10610.2.2003.6632。

[4] 冈萨雷斯等著。数字图像处理（MATLAB 版）[M]，电子工业出版社2009。

6、程序附录清单

(1) 维纳滤波

clc

clear all

close all

I=imread("d:\1.bmp");

subplot(231);

imshow(I);

title("原图");

I=rgb2gray(I);

subplot(232);

imshow(I);

title("灰度图像");

noise=0.1*randn(size(I));

PSF=fspecial("motion",21,11);

Blurred=imfilter(I,PSF,"circular");

BlurredNoisy=im2uint8(Blurred);

NP=abs(fftn(noise)).^2;

NPOW=sum(NP(:)/prod(size(noise)));

NCORR=fftshift(real(ifftn(NP)));

IP=abs(fftn(I)).^2;

IPOW=sum(IP(:)/prod(size(noise)));

ICORR=fftshift(real(ifftn(IP)));

ICORR1=ICORR(:,ceil(size(I,1)/2));

NSR=NPOW/IPOW;

subplot(233);imshow(BlurredNoisy,[]);

title("A=Blurred and Noisy");

subplot(234);imshow(deconvwnr(BlurredNoisy,PSF),[]);

title("deconbwnr(A,PSF,NSR)");

subplot(235);imshow(deconvwnr(BlurredNoisy,PSF,NCORR,ICORR),[]); title("deconbwnr(A,PSF,NCORR,ICORR)");

subplot(236);imshow(deconvwnr(BlurredNoisy,PSF,NPOW,ICORR1),[]); title("deconbwnr(A,PSF,NPOW,ICORR_1_D)");

(2) 正则滤波

clc

clear all

close all

I=imread("d:\1.bmp");

subplot(231);

imshow(I);

title("原图");

I=rgb2gray(I);

subplot(232);

imshow(I);

title("灰度图像");

PSF=fspecial("gaussian",7,10);

V=0.01;

BlurredNoisy=imnoise(imfilter(I,PSF),"gaussian",0,V);

NOISEPOWER=V*prod(size(I));

[J LAGRA]=deconvreg(BlurredNoisy,PSF,NOISEPOWER); subplot(233);imshow(BlurredNoisy);

title("A=Blurred and Noisy");

subplot(234);imshow(J);

title("[J LAGRA]=deconvreg(A,PSF,NP)");

subplot(235);imshow(deconvreg(BlurredNoisy,PSF,[],LAGRA/10)); title("deconvreg(A,PSF,[],0.1*LAGRA)");

subplot(236);imshow(deconvreg(BlurredNoisy,PSF,[],LAGRA*10)); title("deconvreg(A,PSF,[],10*LAGRA");

(3) L-R算法

clc

clear all

close all

I=imread("C:\Documents and Settings\Administrator\桌面\原始图.jpg"); subplot(231);

imshow(I);

I=rgb2gray(I);

subplot(232);

imshow(I);

PSF=fspecial("gaussian",7,10);

V=.0001;

BlurredNoisy=imnoise(imfilter(I,PSF),"gaussian",0,V);

BlurredNoisy=double(BlurredNoisy);

WT=zeros(size(I));

WT(5:end-4,5:end-4)=1;

J1=deconvlucy(BlurredNoisy,PSF);

J2=deconvlucy(BlurredNoisy,PSF,20,sqrt(V));

J3=deconvlucy(BlurredNoisy,PSF,20,sqrt(V),WT);

subplot(233);

imshow(mat2gray(BlurredNoisy));

title("A=Blurred and Noisy");

subplot(234);

imshow(mat2gray(J1));

title("deconvlucy(A,PSF)");

subplot(235);imshow(mat2gray(J2));

title("deconvlucy(A,PSF,NI,DP)");

subplot(236);imshow(mat2gray(J3));

title("deconvlucy(A,PSF,NI,DP,WT)");

(4) 盲去卷积

clc

clear all

close all

I=imread("D:\1.bmp");

subplot(231);

imshow(I);

I=rgb2gray(I);

subplot(232);

imshow(I);

PSF=fspecial("gaussian",7,10);

V=.0001;

BlurredNoisy=imnoise(imfilter(I,PSF),"gaussian",0,V); BlurredNoisy=double(BlurredNoisy);

WT=zeros(size(I));

WT(5:end-4,5:end-4)=1;

INITPSF=ones(size(PSF));

FUN=inline("PSF+P1","PSF","P1");

[J P]=deconvblind(BlurredNoisy,INITPSF,5,10*sqrt(V), WT,FUN,0);

[K P]=deconvblind(BlurredNoisy,INITPSF,10,10*sqrt(V), WT,FUN,0);

[L P]=deconvblind(BlurredNoisy,INITPSF,20,10*sqrt(V), WT,FUN,0); subplot(233);imshow(mat2gray(BlurredNoisy)); title("A=Blurred and Noisy");

subplot(234);imshow(mat2gray(J));

title("True PSF");

subplot(235);imshow(mat2gray(K));

title("Deblured Image");

subplot(236);imshow(mat2gray(L));

title("Recovered PSF");