ch-图象复原ppt课件.ppt

逆滤波比较简单，但没有清楚地说明如何处逆滤波比较简单，但没有清楚地说明如何处理噪声，而维纳滤波综合了退化函数和噪声理噪声，而维纳滤波综合了退化函数和噪声统计特性两个方面进行复原处理统计特性两个方面进行复原处理.5 维纳滤波器维纳滤波器逆滤波方法不能完全恢复原始信号逆滤波方法不能完全恢复原始信号f(x,y)，而只能求出，而只能求出f(x,y)的一个估计值的一个估计值希望找到一种方法，在有噪声条件下，从退化图像希望找到一种方法，在有噪声条件下，从退化图像g(x,y)复原出复原出f(x,y)的估计值，该估计值符合一定的准则。的估计值，该估计值符合一定的准则。维纳滤波维纳滤波(Wiener filtering)=最小均方差滤波最小均方差滤波l维纳滤波是最常用的波是最常用的图像恢复方法像恢复方法l基于基于维纳滤波的波的图像恢复方法是像恢复方法是1967年提出的年提出的C.W.Helstrom,“Image restoration by the method of lest sqaures,”Journal of the Optical Scoiety of America,vol.57,no.3,pp.297-303,1967.C.W.Helstrom,This weeks citation classic,19821967-1982年年SCI引用超引用超过125次次.N.Wiener,“The extrapolation,interpolation and smoothing of stationary time series”,New York:Wiely,1949.在均方误差值最小的准则下得到的在均方误差值最小的准则下得到的 称为对称为对f(x,y)的的最小二乘方估计。最小二乘方估计。按照该准则得到的滤波器叫维纳滤波器。按照该准则得到的滤波器叫维纳滤波器。目标：目标：使得复原后图像使得复原后图像 与原始图像与原始图像 的均方的均方误差最小：误差最小：因此维纳滤波器又称为最小均方差滤波器因此维纳滤波器又称为最小均方差滤波器线性滤波：寻找点扩散函数hw(x,y)，使得则有则有由由Andrews和和Hunt推导满足这一要求的传递函数为：推导满足这一要求的传递函数为：这里，这里，是成像系统传递函数是成像系统传递函数H(u,v)H(u,v)的复共轭；的复共轭；S Sn(u,v)是噪声功率谱；是噪声功率谱；S Sf(u,v)是输入图像的功率谱。是输入图像的功率谱。或者：计算退化算退化图像像g(x,y)的二的二维Fourier变换G(u,v)Wiener滤波的波的过程程:计算理想算理想图像的像的频谱估估计计算点算点扩展函数展函数h(x,y)的二的二维Fourier变换H(u,v)计算退化算退化 图像和噪声的功率像和噪声的功率谱Sf(u,v),Sn(u,v)计算算滤波器波器HW(u,v)求反求反Fourier变换这一方法有如下特点：一方法有如下特点：（1）当当H（u,v）0或或幅幅值很很小小时，分分母母不不为零零，不不会会造造成成严重重的运算的运算误差。差。（2）在信噪比高的）在信噪比高的图像中，即像中，即Sn(u,v)Sf(u,v)维纳滤波复原法特点：如果没有噪声，就成为逆滤波如果没有噪声，就成为逆滤波(3)当理想当理想图像功率像功率谱Sf(u,v)0)时 ,表明我表明我们不可能从全是噪声的不可能从全是噪声的图像中恢复出任何有意像中恢复出任何有意义的信号的信号 （4）往往未退化图像的功率谱）往往未退化图像的功率谱Sf(u,v)难以知道，用下式近似表示：难以知道，用下式近似表示：维纳滤波函数：测试维纳滤波效果：逆滤波和维纳滤波的比较l维纳滤波的结果非常接近原始图像,比逆滤波要好全逆滤波的结果半径受限的逆滤波结果维纳滤波的结果(交互选择K)逆滤波和维纳滤波的比较l(a)运动模糊及均值为0方差为650的加性高斯噪声污染的图像 l(b)逆滤波的结果l(c)维纳滤波的结果l(d)-(f)噪声幅度的方差比(a)小1个数量级l(g)-(i)噪声幅度的方差比(a)小5个数量级逆滤波逆滤波(inverse filtering)线性代数恢复线性代数恢复图图像像恢恢复复方方法法频域图像恢复频域图像恢复维纳滤波维纳滤波(Wiener filtering)代数逆滤波代数逆滤波 奇异值伪逆滤波。奇异值伪逆滤波。线性代数恢复在性代数恢复在1970s提出提出B.R.Hunt,“A matrix theory proof the discrete convolution theorem,”IEEE Trans.Audio Electroacoust,vol.19,no.4,pp.285-288,1971.B.R.Hunt,“The application of constrained least squares estimation to image restoration by digital computer,”IEEE Trans.Computers,vol.22,no.9,pp.805-812,1973.约束最小二乘滤波约束最小二乘滤波 投影法。投影法。6 代数逆代数逆滤波波无无约束代数逆束代数逆滤波波(unconstrained restoration)已知退化已知退化图像的向量形式像的向量形式g和退化矩和退化矩阵H,则无无约束逆束逆滤波波恢复的的恢复的的图像像为:结果结果设恢复的恢复的图像像为 ，如果不考，如果不考虑噪声，用噪声，用 表示表示恢复恢复误差。差。则代数逆代数逆滤波的目的是最小化目波的目的是最小化目标函数函数即即推导推导:使使导数数为零零无无约束代数逆束代数逆滤波波(unconstrained restoration)注意注意Wiener滤波的差波的差别约束最小二乘束最小二乘滤波波(constrained least squares restoration)相相对于无于无约束最小二乘束最小二乘滤波，波，约束最小二乘束最小二乘滤在最小化目在最小化目标函函数数时，利用了估，利用了估计出的噪声的信息作出的噪声的信息作为约束条件。而我束条件。而我们已已经知道噪声模型等信息可以知道噪声模型等信息可以较容易地从容易地从图像中估像中估计出来。通常把出来。通常把噪声的范数作噪声的范数作为约束，即束，即基本思想基本思想令令Q为f的的线性矩性矩阵算子，算子，约束最小二乘恢复束最小二乘恢复问题就是在就是在满足足约束条件束条件 下，使下，使 最小化的最小化的问题。例如选择例如选择Q为有限差分矩阵，使得二阶差分的能量最小为有限差分矩阵，使得二阶差分的能量最小约束最小二乘束最小二乘滤波波(constrained least squares restoration)f(x,y)在在(x,y)处的二阶微分处的二阶微分用拉格朗日法求用拉格朗日法求微分，微分，可以用来调节以满足约束条件。可以用来调节以满足约束条件。约束最小二乘束最小二乘滤波波(constrained least squares restoration)推导推导clearF=checkerboard(8);subplot(2,2,1),imshow(F,);title(a)原始图像,FontSize,12)PSF=fspecial(motion,7,45);MF=imfilter(F,PSF,circular);noise=imnoise(zeros(size(F),gaussian,0,0.001);MEN=double(MF)+noise;subplot(2,2,2),imshow(MEN,);title(b)运动模糊+高斯噪声后的图像,FontSize,10)NSR=sum(noise(:).2)/sum(MEN(:).2);subplot(2,2,3),imshow(deconvwnr(MEN,PSF),);%matlab自带维纳滤波函数，若噪声功率谱设为0，则为逆滤波title(c)逆滤波后的图像,FontSize,10)subplot(2,2,4),imshow(deconvwnr(MEN,PSF,NSR),);%维纳滤波title(d)维纳滤波后的图像,FontSize,10)Matlab约束滤波函数：deconvreg 奇异奇异值（SVD）伪逆逆滤波波(SVD=Singular Value Decomposition)把退化矩把退化矩阵H进行行SVD分解，并得到分解，并得到H的广的广义逆矩逆矩阵。从而可以。从而可以把逆把逆滤波表示成迭代形式，通波表示成迭代形式，通过控制迭代次数，可以在恢复效控制迭代次数，可以在恢复效果果较好好时停止迭代，提供了一种人机交互机会。停止迭代，提供了一种人机交互机会。基本思想基本思想其中其中U的列向量是的列向量是AAT的特征特征值，V的列向量是的列向量是ATA的特征的特征值；是是对角矩角矩阵，对角角线上的上的值 i称称为H的奇异的奇异值。退化矩退化矩阵H的的SVD分解分解为于是恢复于是恢复图像像为则退化矩退化矩阵H的广的广义逆矩逆矩阵为写成迭代形式写成迭代形式7 7 7 7 图象灰度校正和几何畸变校正图象灰度校正和几何畸变校正图象灰度校正和几何畸变校正图象灰度校正和几何畸变校正一、灰度校正一、灰度校正 假设理想图象为假设理想图象为f（x,y），由于灰度失真因），由于灰度失真因子子D（x,y）的影响，实际得到的图象为）的影响，实际得到的图象为g（x,y），），g（x,y）=D（x,y）f（x,y）灰度校正是要从畸变的图象灰度校正是要从畸变的图象g（x,y）中复原原）中复原原 始图象始图象f（x,y）。一种最直接的方法是用光密度计）。一种最直接的方法是用光密度计 测量出被拍摄景物中的某一部分区域测量出被拍摄景物中的某一部分区域S内真实的灰内真实的灰 （亮）度数值（亮）度数值C，而对应的图象灰度为，而对应的图象灰度为gc（x,y），），则则 代入图象复原方程则有代入图象复原方程则有 三、几何畸变校正三、几何畸变校正 在图象获取或显示过程中可能产生图象的几何失真，在图象获取或显示过程中可能产生图象的几何失真，如下图所示如下图所示 从上图可以看出，几何畸变是将无失真坐标系中函从上图可以看出，几何畸变是将无失真坐标系中函 数数f（x,y）变换到另外一个坐标上，例如，原先在（）变换到另外一个坐标上，例如，原先在（x,y）点）点上的象素（灰度）变化到（上的象素（灰度）变化到（u,v），在图象上），在图象上 反映有些位置反映有些位置被挤压，而另一些位置被扩张。我们希望被挤压，而另一些位置被扩张。我们希望 找到这两个坐标找到这两个坐标系之间的关系。系之间的关系。几何基准图像的坐标系统用几何基准图像的坐标系统用(x,y)来表示来表示需要校正的图像的坐标系统用需要校正的图像的坐标系统用(x,y)表示表示设两个图像坐标系统之间的关系用解析式表示设两个图像坐标系统之间的关系用解析式表示通常通常h1(x,y)和和h2(x,y)用多项式来表示：用多项式来表示：通常用线性畸变来近似较小的几何畸变通常用线性畸变来近似较小的几何畸变更精确一些可以用二次型来近似更精确一些可以用二次型来近似 若基准图像为若基准图像为f(x,y)，畸变图像为，畸变图像为g(x,y)，对于景物上的同一，对于景物上的同一个点，假定其灰度不变，则个点，假定其灰度不变，则v 若已知若已知h1(x,y)和和h2(x,y)通常用线性畸变来近似较小的几何畸变通常用线性畸变来近似较小的几何畸变 零级内插零级内插 双线性插值双线性插值则：则：v h1(x,y)和和h2(x,y)未知未知通常用已知的多对对应点来确定系数通常用已知的多对对应点来确定系数a,b 线性畸变线性畸变可由基准图找出三个点可由基准图找出三个点(r1,s1),(r2,s2),(r3,s3)与畸变图像上三与畸变图像上三个点个点(u1,v1),(u2,v2),(u3,v3)一一对应。一一对应。将对应点代入，有：将对应点代入，有：解联立方程组，得出解联立方程组，得出6个系数。个系数。此时此时 二次型畸变二次型畸变有有12个未知量，需要个未知量，需要6对已知对应点对已知对应点代入上式代入上式,同样有同样有解方程组，得到解方程组，得到ai,bi.12个系数。个系数。记做向量矩阵形式记做向量矩阵形式:通常实际应用中，会取多余的对应点对，这时通常实际应用中，会取多余的对应点对，这时A不是方阵不是方阵而是高矩阵，这时矩阵的逆用广义逆矩阵来求解。高矩阵的而是高矩阵，这时矩阵的逆用广义逆矩阵来求解。高矩阵的广义逆矩阵为广义逆矩阵为在广告制作和计算机动画中常常要使物体变形。在广告制作和计算机动画中常常要使物体变形。内插方法除了零级内插和双线性内插外，还有内插方法除了零级内插和双线性内插外，还有B样条插值样条插值和和sinc插值函数内插插值函数内插 几何畸变复原的一套方法也可以用于使图像畸变的工作中。几何畸变复原的一套方法也可以用于使图像畸变的工作中。clearf=checkerboard(24);subplot(2,2,1),imshow(f)title(a)原始图像,FontSize,10)s=0.7;theta=pi/6;T=s*cos(theta)s*sin(theta)0;-s*sin(theta)s*cos(theta)0;0 0 1;tform=maketform(affine,T);g1=imtransform(f,tform,nearest);subplot(2,2,2),imshow(g1)title(b)最邻近插值变换,FontSize,10)g2=imtransform(f,tform);subplot(2,2,3),imshow(g2)title(c)双线性插值变换,FontSize,10)g3=imtransform(f,tform,Fillvalue,0.5);subplot(2,2,4),imshow(g3)title(e)修改c的背景为灰色,FontSize,10)clearf=imread(Fig0515(a)(base-with-control-points).tif);subplot(2,2,1),imshow(f)title(a)原始图像,FontSize,10)inputpoints=83 81;450 56;43 293;249 392;436 442;outputpoints=68 66;375 47;42 286;275 434;523 532;tform=cp2tform(inputpoints,outputpoints,projective);g1=imtransform(f,tform);subplot(2,2,2),imshow(g1)title(b)几何失真图像,FontSize,10)tform=cp2tform(outputpoints,inputpoints,projective);g2=imtransform(g1,tform,Xdata,1 500,Ydata,5 502);subplot(2,2,3),imshow(g2)title(c)复原图像,FontSize,10)附：一维卷积的离散形式：附：一维卷积的离散形式：若若f(x)和和h(x)都是周期为都是周期为N的序列，则的序列，则g(x)也具有周期也具有周期N.若若f(x)和和h(x)不具备周期性，则可以用延拓的方法使其成为周不具备周期性，则可以用延拓的方法使其成为周期函数。令周期期函数。令周期MA+B-1,一个非周期的卷积变成周期卷积后，可以用快速一个非周期的卷积变成周期卷积后，可以用快速fft来实来实现了。现了。写成矩阵的形式，有：写成矩阵的形式，有：其中，其中，H是是MM阶矩阵。阶矩阵。所以又可以写成：所以又可以写成：H为一个循环矩阵。为一个循环矩阵。