数字图像处理6PPT课件

第6章 图像恢复 6.1 退化模型及恢复技术基础 6.2 空间域滤波恢复 6.3频率域滤波恢复 6.4 逆滤波 6.5 最小均方误差滤波器-维纳滤波 第1页/共71页前 言 图像恢复和图像增强一样，都是为了改善图像视觉效果，以及便于后续处理。只是图像增强方法更偏向主观判断，而图像恢复则是根据图像畸变或退化原因，进行模型化处理。本章安排如下：首先介绍退化模型及恢复技术基础；随后分析空间域滤波恢复技术；最后介绍频率域图像恢复技术及应用。 第2页/共71页6.1退化模型及恢复技术基础- 退化的原因 成像系统的象差、畸变、带宽有限等造成图像图像失真；由于成象器件拍摄姿态和扫描非线性引起的图像几何失真；运动模糊，成象传感器与被拍摄景物之间的相对运动，引起所成图像的运动模糊；灰度失真，光学系统或成象传感器本身特性不均匀，造成同样亮度景物成象灰度不同；辐射失真，由于场景能量传输通道中的介质特性如大气湍流效应、大气成分变化引起图像失真；图像在成象、数字化、采集和处理过程中引入的噪声等。第3页/共71页6.1退化模型及恢复技术基础-退化模型 实际所得退化图像 可简化为yxnyxfHyxg, yxnyxhyxfyxnddyxhfyxg, 第4页/共71页6.1退化模型及恢复技术基础-恢复技术的概念及分类 定义：图像恢复是根据退化原因，建立相应的数学模型，从被污染或畸变的图像信号中提取所需要的信息，沿着使图像降质的逆过程恢复图像本来面貌。 第5页/共71页6.1退化模型及恢复技术基础-恢复技术的概念及分类 图像恢复技术的分类： (1)在给定退化模型条件下，分为无约束和有约束两大类;(2)根据是否需要外界干预，分为自动和交互两大类;(3)根据处理所在域，分为频域和空间域两大类。 第6页/共71页6.2 噪声模型（Noise Models） 数字图像的噪声的来源 噪声的频率特性 噪声的空间特性第7页/共71页6.2.1 噪声类型-高斯噪声 定义：222/)(21)(zezp第8页/共71页6.2.1 噪声类型-均匀分布噪声 定义：其它01)(bzaabzp第9页/共71页6.2.1 噪声类型-脉冲噪声（椒盐噪声） 定义其它0)(bzPazPzpba第10页/共71页例6.1A=imread(fig606a.jpg); %读取图像figure,imshow(A); %显示图像A=rgb2gray(A);figure,hist(double(A),10); %求出A的直方图并显示B=imnoise(A,gaussian,0.05); %对A附加高斯噪声figure,imshow(B); %显示附加高斯噪声后的图像Bfigure,hist(double(B),10); %求出B的直方图并显示C=imnoise(A,speckle,0.05); %对A附加均匀分布噪声figure,imshow(C); %显示附加均匀噪声后的图像Cfigure,hist(double(C),10); %求出C的直方图并显示D=imnoise(A,salt & pepper,0.05); %对A附加椒盐噪声figure,imshow(D); %显示附加椒盐噪声后的图像Dfigure,hist(double(D),10); %求出D的直方图并显示第11页/共71页6.2.1 噪声类型-示例 (a) 原图 (b) 高斯噪声图 第12页/共71页6.2.1 噪声类型-示例 (c) 均匀分布噪声 (d) 椒盐噪声 第13页/共71页6.3 空间域滤波恢复 定义： 空间域滤波恢复即是在已知噪声模型的基础上，对噪声的空间域滤波 第14页/共71页6.3.1 均值滤波器 采用均值滤波模板对图像噪声进行滤除第15页/共71页6.2.2 均值滤波-类型 算术均值滤波器： 几何均值滤波器 1( , )( , )( , )xymns tSf x yg s t Ststsgmnyxf),(),(1),(第16页/共71页6.2.2 均值滤波-类型 谐波均值滤波器 逆谐波均值滤波器 xyStstsgmnyxf),(),(1),(xyxyStsQStsQtsgtsgyxf),(),(1),(),(),(第17页/共71页6.2.2 均值滤波-示例 (a) 输入图像； (b)高斯噪声污染图像；(c) 用均值滤波结果 第18页/共71页6.2.2 均值滤波-示例 (d) 几何均值滤波（e）Q1.5的逆谐波滤波 (f) Q=1.5滤波的结果第19页/共71页6.3.2 顺序统计滤波 1.中值滤波 其中，其中，g为输入图像，sxy为滤波窗口。修正后的阿尔法均值滤波器 xyStstsgMedyxf),(),(),(Stsrtsgdmnyxf),(),(1),(第20页/共71页img=imread(cameraman.tif); imshow(img); %显示图像 img_noise=double(imnoise(img,gaussian,0.06);%对图像附加高斯噪声figure,imshow(img_noise,); %显示加噪图像img_mean=imfilter(img_noise,fspecial(average,3);%对附加有高斯噪声的图像实行算术均值滤波figure; imshow(img_mean,); %显示算术均值滤波后的图像img_mean=exp(imfilter(log(img_noise),fspecial(average,3);%对附加有高斯噪声的图像实行几何均值滤波figure; imshow(img_mean,); %显示几何均值滤波后的图像Q=-1.5; %对高斯噪声图像实行Q取负数的逆谐波滤波img_mean=imfilter(img_noise.(Q+1),fspecial(average,3)./imfilter(img_noise.Q,fspecial(average,3);figure; imshow(img_mean,); %显示逆谐波滤波后的图像Q=1.5; %对高斯噪声图像实行Q取正数的逆谐波滤波img_mean=imfilter(img_noise.(Q+1),fspecial(average,3)./imfilter(img_noise.Q,fspecial(average,3);figure; imshow(img_mean,); %显示逆谐波滤波后的图像第21页/共71页例6.3img=rgb2gray(imread(football.jpg); %读取图像并转换成灰度图像figure; imshow(img); %显示图像img_noise=double(imnoise(img,salt & pepper,0.06);%加椒盐噪声figure,imshow(img_noise,); %显示附加椒盐噪声后的图像img_mean=imfilter(img_noise,fspecial(average,5);%对附加有椒盐噪声的图像实行算术均值滤波figure; imshow(img_mean,); %显示算术均值滤波后的图像img_median=medfilt2(img_noise); %对附加有椒盐噪声的图像实行中值滤波figure; imshow(img_median,); %显示中值滤波后的图像img_median2=medfilt2(img_median); %对中值滤波处理后的图像再次实行中值滤波figure; imshow(img_median2,); %显示再次中值滤波后的图像第22页/共71页第23页/共71页中值滤波示例 (a)椒盐噪声污染的图像 (b) 均值滤波结果； 第24页/共71页中值滤波示例(续） (c) 中值滤波结果 (d)对c图再次中值滤波 第25页/共71页最大/最小滤波2.最大/最小滤波1）最大值滤波器为： 2）最小值滤波器为： xyStstsgyxf),(),(max),(xyStstsgyxf),(),(min),(第26页/共71页最大/最小滤波示例 (a)噪声图像 (b) 最大滤波结果 (c) 最小滤波结果第27页/共71页6.4频率域滤波恢复 原理： 时域卷积相当于频域乘积。因此可以在频率域中直接设计滤波器，对信号进行恢复处理。 分类：常用的图像恢复方法有带阻滤波器、带通滤波器 、陷波滤波器等第28页/共71页6.4.1 带阻滤波器 1.理想带阻滤波器2.巴特沃思带阻滤波器 00001,( , )2( , )0,( , )221,( , )2WD u vDWWH u vDD u vDWD u vD22201( , )( , )1( , )nH u vD u v WD u vD第29页/共71页6.4.1 带阻滤波器 高斯带阻滤波器 2220( , )12( , )( , )1Du vDD u v WH u ve (a)理想带阻滤波器;(b)巴特沃思带阻滤波;(c)高斯带阻滤波器 第30页/共71页6.4.1 带阻滤波器示例(a)(b)(c)(d) (a) 被正弦噪声污染的图像；(b) 图(a)的频谱； (c) 巴特沃思带阻滤波器；(d) 滤波效果图 第31页/共71页 I=imread(woman1.bmp); %读取图像I=rgb2gray(I); %转换成灰度图像M,N=size(I); %得到图像的高度和宽度P=I;for i=1:M for j=1:N P(i,j)=P(i,j)+20*sin(20*i)+20*sin(20*j); %添加周期噪声 endendfigure,imshow(I); %显示原图像figure,imshow(P); %显示加噪图像IF=fftshift(fft2(I); %对原图像作傅里叶变换，并将原点移至中心IFV=log(1+abs(IF); %原图像的频谱PF=fftshift(fft2(P); %对加噪图像作傅里叶变换，并将原点移至中心PFV=log(1+abs(PF); %加噪图像的频谱第32页/共71页figure,imshow(IFV,); %显示原图像的频谱figure,imshow(PFV,); %显示加噪图像的频谱freq=50; %设置带阻滤波器中心频率width=5; %设置带阻滤波器频带宽度ff = ones(M,N);for i=1:M for j=1:N ff(i,j) = 1-exp(-0.5*(i-M/2)2+(j-N/2)2)-freq2)/(sqrt(i.2+j.2)*width)2); %高斯带阻滤波器 endendfigure,imshow(ff,); %显示高斯带阻滤波器out = PF.* ff; %矩阵点乘实现频域滤波out = ifftshift(out); %原点移回左上角out = ifft2(out); %傅里叶逆变换out = abs(out); %取绝对值out = out/max(out(:); %归一化figure,imshow(out,); %显示滤波结果第33页/共71页第34页/共71页第35页/共71页6.4.2 带通滤波器 带通滤波器执行与带阻滤波器相反的操作 可用全通滤波器减去带阻滤波器来实现带通滤波器),(1),(vuHvuHbsbp第36页/共71页6.4.3 陷波滤波器 陷波滤波器被用于阻止(或通过)事先定义的中心频率领域内的频率 由于傅立叶变换时对称的,因此陷波滤波器必须以关于原点对称的形式出现。10200( , )( , )( , )1D u vDD u vDH u v或其他第37页/共71页6.3.3 陷波滤波器(a)(b)(c)(a)理想陷波滤波器;(b)巴特沃思陷波滤波器;(c)高斯陷波滤波器 第38页/共71页6.5 估计退化函数( , ) ( , )( , )g x yH f x yx y 退化模型退化模型:(1) 如果如果:1212( , )( , )( , )( , )H af x ybfx yaH f x ybH fx y则系统则系统H是一个线性系统是一个线性系统.( , ) ( , )( , ),:g x yH f x yf x y若系统对于任意和有 (,)(,)H f xyg xy则系统则系统H称为位置不变系统称为位置不变系统(或空间不变系统或空间不变系统).(2)如果退化模型为线性和位置不变的如果退化模型为线性和位置不变的,其可表示为其可表示为:( , )( ,) (,)( , )g x yfh xyd dx y 即即:( , )( , )*( , )( , )g x yh x yf x yx y( , )( , ) ( , )( , )G u vH u v F u vN u v( , , ,).h xyH为系统 的冲激响应(点扩散函数)第39页/共71页6.5 估计退化函数估计退化函数许多退化类型可以近似表示为线性的位置不变过程许多退化类型可以近似表示为线性的位置不变过程.非线性的与位置有关的技术难以求解非线性的与位置有关的技术难以求解.由于退化模型为卷积的结果由于退化模型为卷积的结果,且图像复原需要滤波器且图像复原需要滤波器,应此术语应此术语”图像去卷积图像去卷积”常用于表示线性图像复原常用于表示线性图像复原,而用于复原处理的滤波器称为而用于复原处理的滤波器称为”去卷积滤波器去卷积滤波器”.第40页/共71页6.5 估计退化函数退化函数通常未知退化函数通常未知,因此在复原之前需要估计退化函数因此在复原之前需要估计退化函数.估计退化函数的方法估计退化函数的方法:(1)观察法观察法(2)实验法实验法(3)数学建模法数学建模法( , )( , )*( , )( , )g x yh x yf x yx y( , )( , ) ( , )( , )G u vH u v F u vN u v第41页/共71页6.5 估计退化函数(1) 观察法观察法( , )( , )( , )sssG u vH u vF u v收集图像自身的信息来估计退化函数收集图像自身的信息来估计退化函数.例如例如: 对于模糊图像对于模糊图像,选择一小部分图像选择一小部分图像,强信号区强信号区,减少噪声影响减少噪声影响.并构建一个不退化的图像并构建一个不退化的图像( , )sgx y( , )sfx y( , )( , )sH u vH u v然后根据推出第42页/共71页6.5 估计退化函数( , )g x y(2) 试验估计法试验估计法使用与获取退化图像的设备相似的装置使用与获取退化图像的设备相似的装置,得到准确的退化估计得到准确的退化估计.小亮点小亮点成像系统成像系统H( , )( , )G u vH u vA由于冲激的傅立叶变换为常数由于冲激的傅立叶变换为常数A,可得可得:实验估计模型如下实验估计模型如下:第43页/共71页6.5 估计退化函数冲激特性的退化估计冲激特性的退化估计(a)一个亮脉冲一个亮脉冲(b) 图像化的图像化的(退化的退化的)冲激冲激第44页/共71页6.5 估计退化函数(3) 模型估计法模型估计法建立退化模型建立退化模型,模型要把引起退化的环境因素考虑在内模型要把引起退化的环境因素考虑在内.22 5/6()( , )k uvH u ve例如退化模型例如退化模型就是基于大气湍流的物理特性而提出来的就是基于大气湍流的物理特性而提出来的,其中其中k为常数为常数,与湍流特性相关与湍流特性相关.第45页/共71页6.5 估计退化函数大气湍流模型的解释大气湍流模型的解释(a)可忽略的湍流可忽略的湍流(b) 剧烈湍流剧烈湍流,k=0.0025(c)中等湍流中等湍流,k=0.001(d) 轻微湍流轻微湍流,k=0.00025另外也可以从基本原理开始推导出退化模型另外也可以从基本原理开始推导出退化模型.如匀速直线运动造成的模如匀速直线运动造成的模糊就可以运用数学推导出其退化函数糊就可以运用数学推导出其退化函数.第46页/共71页2. 运动模糊模型 当成像传感器与被摄景物之间足够快的相对运动时，所摄取的图像就会出现“运动模糊”，运动模糊是场景能量在传感器拍摄瞬间内在像平面上的非正常积累。0( ),( )tyyt dtT00g(x,y)=fx-x第47页/共71页例6.7I=imread(i_camera.bmp); %读取图像I=rgb2gray(I); %转换为灰度图像figure,imshow(I); %显示图像LEN=25; %设置线性运动位移THETA=11; %设置旋转角度PSF=fspecial(motion,LEN,THETA); %图像线性运动Blurred=imfilter(I,PSF,circular,conv); %图像被线性运动模糊figure,imshow(Blurred); %显示运动模糊后的图像第48页/共71页第49页/共71页6.6 逆滤波 1. 逆滤波原理： 图像退化模型： 傅立叶变换 逆滤波恢复 yxnyxhyxfyxg, vuNvuFvuHvuG,vuHvuNvuHvuGvuF,第50页/共71页逆滤波原理 分析： 1.实验证明，当退化图像的噪声较小，即轻度降质时，采用逆滤波恢复的方法可以获得较好的结果。通常，在离频率平面原点较远的地方数值较小或为零，因此图像恢复在原点周围的有限区域内进行，即将退化图像的傅立叶谱限制在没出现零点而且数值又不是太小的有限范围内。 2. 当噪声作用范围很大时，逆滤波不能从噪声中提取图像。第51页/共71页逆滤波示例 (a)原图 (b)退化图像 (c) 逆滤波结果第52页/共71页6.7 最小均方误差滤波器-维纳滤波 逆滤波恢复方法对噪声极为敏感，要求信噪比较高，通常不满足该条件。 为了解决高噪声情况下的图像恢复问题，可采用最小均方滤波器来解决，其中，用得最多的是维纳滤波器第53页/共71页6.7 最小均方误差滤波器-维纳滤波 目标函数:采用拉格朗日乘数法，在有噪声条件下，从退化图像g(x,y)复原出f(x,y)的估计值，该估计值符合一定的准则。 维纳滤波器结果上式称为维纳滤波，括号中的项组成的滤波器通常称为最小均方误差滤波器,或最小二乘方误差滤波器。 222)(minnfHgQffJ),(),(/ ),(),(),(),(2*vuGvuSvuSvuHvuHvuFffnn第54页/共71页6.7 最小均方误差滤波器-维纳滤波示例 (a) 运动模糊退化图像 (b)7次循环 (c) 15次循环 第55页/共71页例6.8I=imread(i_camera.bmp); %读取图像figure(1);imshow(I,); %显示图像PSF=fspecial(motion,25,11); %运动模糊函数，运动位移是25像素，角度是11Blurred=imfilter(I,PSF,conv,circular); %对图像运动模糊处理Noise=0.1*randn(size(I); %正态分布的随机噪声BlurredNoisy=imadd(Blurred,im2uint8(Noise);%对退化后的图像附加噪声figure(2);imshow(BlurredNoisy,); %显示运动模糊且加噪声后图像WI1=deconvwnr(BlurredNoisy,PSF); %不带参数的维纳滤波（逆滤波）复原figure(3); imshow(WI1,); %显示逆滤波复原结果NSR=sum(Noise(:).2)/sum(im2double(I(:).2); %计算噪信比WI2=deconvwnr(BlurredNoisy,PSF,NSR); %带噪信比参数的维纳滤波复原figure(4);imshow(WI2,); %显示带噪信比参数维纳滤波复原结果NP=abs(fftn(Noise).2;NCORR=real(ifftn(NP); %计算噪声的自相关函数IP=abs(fftn(im2double(I).2;ICORR=real(ifftn(IP); %计算信号的自相关函数WI3=deconvwnr(BlurredNoisy,PSF,NCORR,ICORR); %带自相关函数的维纳滤波复原 figure(5);imshow(WI3,); 第56页/共71页第57页/共71页6.8 几何失真校正 在图像的获取或显示过程中往往会产生几何失真，例如成像系统有一定的几何非线性。这主要是由于视像管摄像机及阴极射线管显示器的扫描偏转系统有一定的非线性，因此会造成如图6.20所示的枕形失真或桶形失真。 第58页/共71页6.8.1 空间变换（Spatial Transformation） 假设一幅图像为 ，经过几何失真变成了 ，这里的 表示失真图像的坐标，它已不是原图像的坐标了。上述变化可表示为：这里 ， 和是空间变换，产生了几何失真图像 。 ( , )f x y( , )g u v( , )u v( , )ur x y( , )vs x y( , )r x y( , )s x y( , )g u v第59页/共71页1已知 和 条件下的几何校正（1）直接法。先由 推出 ，然后依次计算每个像素的校正坐标值， ( , )r x y( , )s x y( , )( , )ur x yvs x y( , )( , )xr u vys u v第60页/共71页（2）间接法。设恢复的图像像素在基准坐标系统为等距网格的交叉点，从网格交叉点的坐标(x,y)出发算出在已知畸变图像上的坐标(u,v)，即：( , ),u vr x ys x y 第61页/共71页 2 和 未知条件下的几何校正 假定基准图像像素的空间坐标(x,y)和被校正图像对应像素的空间坐标(u,v)之间的关系用二元多项式来表示。式中，N为多项式的次数，aij和bij为各项待定系数。( , )r x y( , )s x y1100( , )NNijijijur x ya x y1100( , )NNijijijvs x yb x y第62页/共71页对于线性失真：对于一般的（非线性）二次失真：001001( , )ur x yaa xa y001001( , )vs x ybb xb y00100111( , )ur x yaa xa ya xy第63页/共71页 利用“连接点”建立失真图像与校正图像之间其他像素空间位置的对应关系，而这些“连接点”在输入（失真）图像和输出（校正）图像中的位置是精确已知的。 第64页/共71页例6.9f=imread(fig620.jpg); %读取图像figure(1),imshow(f); %显示原图像k=0.7; %变换拉伸系数theta=pi/6; %变换旋转角度T=k*cos(theta) k*sin(theta) 0 -k*sin(theta) k*cos(theta) 0 0 0 1; %变换矩阵tform1=maketform(affine,T); %几何变换结构g1=imtransform(f,tform1,nearest); %以最近邻插值进行仿射变换figure(2),imshow(g1); %显示变换结果g2=imtransform(f,tform1,bilinear); %以双线性插值进行仿射变换figure(3),imshow(g2); %显示变换结果第65页/共71页第66页/共71页例6.10f=imread(fig622.jpg); %读取原图像figure(1),imshow(f); %显示原图像g=imread(fig622b.jpg); %读取几何失真图像figure(2),imshow(g); %显示几何失真图像 %利用cpselect(g, f)交互选择如下9对连接点base_points=61.7500 99.2500; 36.7500 149.2500; 86.7500 148.7500; 128.2500 109.2500; 168.7500 156.7500; 166.2500 117.7500; 211.7500 91.7500; 212.2500 133.2500; 211.7500 173.2500; %从f中选择的9个连接点的坐标矩阵 input_points=111.7500 98.7500; 110.2500 148.2500; 162.7500 148.7500; 184.7500 108.7500; 251.7500 159.7500; 224.2500 117.2500; 258.2500 91.7500; 279.7500 133.7500; 299.7500 173.7500; %从g中选择的9个连接点的坐标矩阵 第67页/共71页tform=cp2tform(input_points,base_points,projective); %由9对连接点坐标矩阵建立几何变换结构 gp=imtransform(g,tform,XData,1 256,YData,1 256); %由得到的几何变换结构对失真图像g进行校正 figure(3),imshow(gp); %显示校正后的图像第68页/共71页第69页/共71页第70页/共71页感谢您的观看！第71页/共71页