逆滤波和维纳滤波

简介

在图像的获取、传输以及记录保存过程中，由于各种因素，如成像设备与目标物体的相对运动，大气的湍流效应，光学系统的相差，成像系统的非线性畸变，环境的随机噪声等原因都会使图像产生一定程度的退化，图像退化的典型表现是图像出现模糊、失真，出现附加噪声等。由于图像的退化，使得最终获取的图像不再是原始图像，图像效果明显变差。为此，要较好地显示原始图像，必须对退化后的图像进行处理，恢复出真实的原始图像，这一过程就称为图像复原。

大气湍流退化

对经过大气湍流退化的图片实现全逆滤波，半径受限逆滤波以及维纳滤波，并对比。

三种滤波思想

直接全逆滤波、半径受限逆滤波和维纳滤波

直接全逆滤波

如果退化图像为 g(x, y)，原始图像为 f (x, y)，在不考虑噪声的情况下，由傅立叶变换的卷积定理可知有下式成立 G(u, v) = H(u, v) ∗ F (u, v)。由此可见，如果已知退化图像的傅里叶变换和系统的冲击响应函数 G(u, v)，则可以求出原图像的傅里叶变换 H(u, v)，之后使用傅里叶逆变换即可得到原图像。这就是全逆滤波的主要思想。

半径受限逆滤波

全逆滤波在有噪声的情况下：
F ( u , v ) = G ( u , v ) / H ( u , v ) F (u, v) = G(u, v) /H(u, v) F(u,v)=G(u,v)/H(u,v)
F ’ ( u , v ) = G ( u , v ) / H ( u , v ) − N ( u , v ) / H ( u , v ) F ’(u, v) =G(u, v)/H(u, v) −N (u, v) /H(u, v) F’(u,v)=G(u,v)/H(u,v)−N(u,v)/H(u,v)
N (u, v) 是噪声的傅里叶变换。如果此时在 N (u, v) 不为 0 的地方，H(u, v) 过零或者比较小，则会导致噪声被放大。因此对 F (u, v) 首先使用低通滤波器过滤掉高频的噪声，之后再除以H(u, v)，这样虽然有可能会导致图像的高频细节被忽略，但是相对于全逆滤波更加稳定。

维纳滤波

维纳滤波就是最小二乘滤波，它是使原始图像 f(x,y) 与其恢复图像 f’(x,y) 之间的均方误差最小的复原方法。具有较好的抑制噪音的能力。由于在日常使用中，无法确定图片的信噪比，因此在公式中使用参数 k 来代替信噪比。其公式为：
F ’ ( u , v ) = 1 H ( u , v ) ∗ ∣ H ( u , v ) ∣ 2 ∣ H ( u , v ) ∣ 2 + k F’(u,v) = \frac{1}{H(u,v)} * \frac{|H(u,v)|^{2}}{|H(u,v)|^{2} + k} F’(u,v)=H(u,v)1∗∣H(u,v)∣2+k∣H(u,v)∣2

代码1

主函数

clc;
clear;
close all;%% 课本图 5.28
% 读取图片
im = imread('demo-1.jpg');   % 原始图像 480x480 uint8%% 图像退化（大气湍流模型）
% Output（H：退化模型， im_f：退化后图片）
[H, im_f,im_F] = atmosph(im);        %% 全逆滤波，半径受限逆滤波
D0 = 60;
% Input（im_f：退化图片，H：退化模型，D0：半径）
% Output（im_inverse：全逆滤波结果，im_inverse_b：半径受限逆滤波）[im_inverse, im_inverse_b] = my_inverse(im_f, H, D0);
% [im_inverse, im_inverse_b] = my_inverse_F(im_f, H, D0,0.001); %% 维纳滤波
K = 0.0001;
% Input（im_f：退化图片，H：退化模型，K：维纳滤波常数）
im_wiener = my_wiener(im_f, H, K);%% 显示结果
figure(1);
subplot(231); imshow(im); title('原图'); axis on
subplot(232); imshow(im_f); title('大气湍流(k=0.0025)'); axis on
subplot(233); imshow(im_inverse); title('全逆滤波'); axis on
subplot(234); imshow(im_inverse_b); title('半径受限的逆滤波'); axis on
subplot(235); imshow(im_wiener); title('维纳滤波'); axis on

全逆滤波和半径受限的全逆滤波

function [im_inverse, im_inverse_b] = my_inverse(img, H, D0)[M,N] = size(img);%图像进行二位傅里叶变换，之后移动到中心点
im_double = mat2gray(img,[0 255]);
im_F = fftshift(fft2(im_double));      % 空域 > 频域% 10阶巴特沃斯低通滤波器
H2 = zeros(M,N);
for x = (-M/2):1:(M/2)-1for y = (-N/2):1:(N/2)-1D2 = (x^2 + y^2)^(1/2);H2(x+(M/2)+1,y+(N/2)+1) = 1/(1+(D2/D0)^20);%十阶   end
end
%滤波后图像
im_H2 = im_F .* H2; % 全逆滤波
% 抑制低幅值的频率的全逆滤波
% im_inverse = zeros(M,N);
% num_F = abs(im_F);
% for x = 1:M
%      for y = 1:N
%          if(num_F(x,y)>0.78)
%              im_inverse(x,y) = im_F(x,y) / H(x,y);
%          else
%              im_inverse(x,y) = im_F(x,y);
%          end
%      end
% end%直接全逆滤波
im_inverse = im_F ./ H;  im_inverse_double = real(ifft2(ifftshift(im_inverse)));    % 频域 > 空域
im_inverse = im2uint8(mat2gray(im_inverse_double));% 半径受限逆滤波
% 抑制低幅值的频率的全逆滤波
% im_inverse_b = zeros(M,N);
% num_H2 = abs(im_H2);
% for x = 1:M
%      for y = 1:N
%          if(num_H2(x,y)>0.78)
%              im_inverse_b(x,y) = im_H2(x,y) / H(x,y);
%          else
%              im_inverse_b(x,y) = im_H2(x,y);
%          end
%      end
% end
% 直接全逆滤波
im_inverse_b = im_H2 ./ H;
im_inverse_b_double = real(ifft2(ifftshift(im_inverse_b)));    % 频域 > 空域
im_inverse_b = im2uint8(mat2gray(im_inverse_b_double));
end

实验结果1

在本次实验中，首先对原图进行大气湍流退化（退化系数 k = 0.0025），之后分别对退化的图像进行全逆滤波，半径受限逆滤波（取半径 D0 = 60），维纳滤波（取信噪比 k = 0.0001）。

实验结果如下：

如上图所示，实验中全逆滤波的结果是一团噪声，但是在图像退化过程中未引入噪声。在实验中，经过反复调试，终于发现问题出在图像退化时的傅里叶变换的时候。在图像退化的反傅里叶变换时，存在一次对图像的取实部，这个操作没有问题，但是由于图像是离散的，当图像再次傅里叶变换，会在各个频段出现看似随机的噪声。这直接导致了由于 H(u,v) 过零导致的噪声放大。
为了说明全逆滤波对无噪声图像有着极好的恢复能力，重写给定退化函数，增加返还量im_F，不做反变换，用来使用全逆滤波。新的实验结果如下图：

可以看出此时全逆滤波的效果很明显。

实验分析1

从上述实验结果可以看出，全逆滤波对无噪声的图像，具有极好的恢复能力，但是及其不稳定，可能由于各种各样的噪声，导致图像不仅没有变得清晰，反而变的更加模糊。而半径受限的全逆滤波，尽管理想的效果不如全逆滤波，会丢失很多细节，但是对于噪声的忍耐能力更强，相对也更加稳定。从上面比较来看，维纳滤波是效果最好的，并且其对噪声的抑制能力最强，但是需要知道图像的信噪比，才能最好的滤波。

运动模糊和高斯噪声污染

运动模糊是图像中的对象由于平移运动，导致的多幅图像的线性堆积。高斯噪声是服从高斯分布的亮度随机变化，是图像中常见的噪声污染。

运动模糊和高斯噪声设计思路

运动模糊是由于物体运动导致的，高斯噪声是一种由于设备、环境导致的随机性误差。

代码2

主函数

clc;
clear;
close all;%% 课本图 5.29
% 读取图片
im = imread('demo-2.jpg');   % 原始图像 688x688 uint8%% 图像退化（运动模糊+高斯噪声）
% [~, im1_f] = motionblur(im, 0.01);        % 噪声方差=0.01
% [~, im2_f] = motionblur(im, 0.001);
% [H, im3_f] = motionblur(im, 0.0000001);        % H：退化模型[~, im1_f,im1_F] = motionblur(im, 0.01);        % 噪声方差=0.01[~, im2_f,im2_F] = motionblur(im, 0.001);        [H, im3_f,im3_F] = motionblur(im, 0.0000001);        % H：退化模型%% 全逆滤波，半径受限逆滤波
D0 = 33;[~, im1_inverse] = my_inverse(im1_f, H, D0);[~, im2_inverse] = my_inverse(im2_f, H, D0);[~, im3_inverse] = my_inverse(im3_f, H, D0);
%  [im1_inverse, ~] = my_inverse(im1_f, H, D0);
%  [im2_inverse, ~] = my_inverse(im2_f, H, D0);
%  [im3_inverse, ~] = my_inverse(im3_f, H, D0);
%  [~, im1_inverse] = my_inverse_F(im1_f, H, 33);
%  [~, im2_inverse] = my_inverse_F(im2_f, H, 33);
%  [~, im3_inverse] = my_inverse_F(im3_f, H, 33);
%  [im1_inverse, ~] = my_inverse_F(im1_f, H, D0 ,0.95);
%  [im2_inverse, ~] = my_inverse_F(im2_f, H, D0 ,0.31);
%  [im3_inverse, ~] = my_inverse_F(im3_f, H, D0 ,0.005);%% 维纳滤波
%K = 0.0001;%0.1 0.008 0.0001
im1_wiener = my_wiener(im1_f, H, 0.1);
im2_wiener = my_wiener(im2_f, H, 0.008);
im3_wiener = my_wiener(im3_f, H, 0.00001);%% 显示结果
figure(1);
subplot(331); imshow(im1_f); title('运动模糊+加性噪声(sigma)'); axis on
subplot(332); imshow(im1_inverse); title('逆滤波结果'); axis on
subplot(333); imshow(im1_wiener); title('维纳滤波结果'); axis on
subplot(334); imshow(im2_f); title('运动模糊+加性噪声(sigma*0.1)'); axis on
subplot(335); imshow(im2_inverse); title('逆滤波结果'); axis on
subplot(336); imshow(im2_wiener); title('维纳滤波结果'); axis on
subplot(337); imshow(im3_f); title('运动模糊+加性噪声(sigma*0.00001)'); axis on
subplot(338); imshow(im3_inverse); title('逆滤波结果'); axis on
subplot(339); imshow(im3_wiener); title('维纳滤波结果'); axis on

运动模糊

function [H,im_blured, im_blured_F] = motionblur(img, sigma)[M,N] = size(img);
%图像进行二位傅里叶变换，之后移动到中心点
im_double = mat2gray(img,[0 255]);
im_F = fftshift(fft2(im_double));      % 空域 > 频域H = zeros(M,N);
gauss = 0 + sqrt(sigma) * randn(M,N);%二维高斯分布矩阵 0是均值
% 运动模糊
a = 0.1;
b = 0.1;
T = 1;for u = -M/2:1:M/2-1for v = -N/2:1:N/2-1L = ((a*u)+(b*v))*pi;if (L == 0)H(u + M/2 + 1,v + N/2 + 1) = 1;elseH(u + M/2 + 1,v + N/2 + 1) = T * sin(L) * exp( -L * 1i )/L;endend
endim_G = im_F .* H;
%im_double = real(ifft2(ifftshift(im_G)));    % 频域 > 空域% noise
im_gauss = fftshift(fft2(gauss));      % 空域 > 频域im_blured_F = im_G + im_gauss;
%im_blured_F = im_G;
im_blured = real(ifft2(ifftshift(im_blured_F)));    % 频域 > 空域
end

维纳滤波

function im_wiener = my_wiener(img, H, K)
[M,N] = size(img);
H2 = zeros(M,N);
Mo = zeros(M,N);
Me = zeros(M,N);
% 空域 > 频域
im_double = mat2gray(img,[0 255]);
im_G = fftshift(fft2(im_double));%得到系数H2
for u = 1 : Mfor v = 1 : NMo(u,v) = (abs(H(u,v)))^2;Me(u,v) = (H(u,v)*(Mo(u,v) + K));H2(u,v) = Mo(u,v)/Me(u,v);end
endim_F = im_G.* H2;
im_wiener_double = real(ifft2(ifftshift(im_F)));    % 频域 > 空域
im_wiener = im2uint8(mat2gray(im_wiener_double));
end

实验结果2

首先是按照实验模板做了第一次实验，在实验中给三种情况下的维纳滤波分别调整的不同的信噪比系数 k，分别为 0.1、0.008、0.00001。但是针对逆滤波，无论是半径受限，还是未使用半径受限，都是一团噪声，无法产生需要的结果。如下图所示，（图 1 是全逆滤波，图 2是半径受限滤波）

实验思考

如上面所示，全滤波存在问题，但是按照理论知识，在噪声很小的情况下，逆滤波应该表现的相对清晰，因此，采用和大气湍流退化相同的处理方法，即对退化后的图像不进行反变换，直接全逆滤波处理，最终在半径受限的时候，得到了理想的处理结果，如下图所示：

此时，笔者在反复分析图像处理过程时，发现对于这种程度的高斯噪声，其影响的频率点的幅值不会特别高。因此，笔者写了一个低幅值抑制的函数先对图像进行处理，然后对图像使用逆滤波，最后得到了细节保留优于维纳滤波的处理算法。下面是处理的函数和算法。本次实验的下限 k 分别是0.95，0.31，0.005。

function [im_inverse, im_inverse_b] = my_inverse_F(img, H, D0, k)[M,N] = size(img);
%图像进行二位傅里叶变换，之后移动到中心点im_double = mat2gray(img,[0 255]);im_F = fftshift(fft2(im_double));      % 空域 > 频域
% im_F = img;
% 10阶巴特沃斯低通滤波器
H2 = zeros(M,N);
for x = (-M/2):1:(M/2)-1for y = (-N/2):1:(N/2)-1D2 = (x^2 + y^2)^(1/2);H2(x+(M/2)+1,y+(N/2)+1) = 1/(1+(D2/D0)^20);%十阶   end
end
%滤波后图像
im_H2 = im_F .* H2; % 全逆滤波
% 抑制低幅值的频率的全逆滤波
im_inverse = zeros(M,N);
num_F = abs(im_F);
for x = 1:Mfor y = 1:Nif(num_F(x,y)>k)im_inverse(x,y) = im_F(x,y) / H(x,y);elseim_inverse(x,y) = im_F(x,y);end end
end%直接全逆滤波
%im_inverse = im_F ./ H;  im_inverse_double = real(ifft2(ifftshift(im_inverse)));    % 频域 > 空域
im_inverse = im2uint8(mat2gray(im_inverse_double));% 半径受限逆滤波
% 抑制低幅值的频率的全逆滤波
im_inverse_b = zeros(M,N);
num_H2 = abs(im_H2);
for x = 1:Mfor y = 1:Nif(num_H2(x,y)>k)im_inverse_b(x,y) = im_H2(x,y) / H(x,y);elseim_inverse_b(x,y) = im_H2(x,y);end end
end% 直接全逆滤波
% im_inverse_b = im_H2 ./ H;  im_inverse_b_double = real(ifft2(ifftshift(im_inverse_b)));    % 频域 > 空域
im_inverse_b = im2uint8(mat2gray(im_inverse_b_double));
end

i f ( n u m F ( x , y ) > k ) 如果幅值大于下限 k if (num_F (x, y) > k) 如果幅值大于下限 k if(numF(x,y)>k)如果幅值大于下限k i n v e r s e ( x , y ) = i m F ( x , y ) / H ( x , y ) ; 正常使用逆滤波 inverse(x, y) = im_F (x, y)/H(x, y);正常使用逆滤波 inverse(x,y)=imF(x,y)/H(x,y);正常使用逆滤波 e l s e 否则，说明其信息重要度低 else 否则，说明其信息重要度低 else否则，说明其信息重要度低 i n v e r s e ( x , y ) = i m F ( x , y ) ; 保持原值 inverse(x, y) = im_F (x, y); 保持原值 inverse(x,y)=imF(x,y);保持原值

实验分析2

这次实验发现，使用维纳滤波交互性的选择信噪比 K，可以得到较好的图像显示效果。如期望的，全逆滤波产生的图像质量非常差，即使使用了半径限制的逆滤波依旧会产生较大的噪声干扰。这是因为 H 的小值较多，会使噪声放大。

在实验中，通过限制低幅度值的频谱值，可以得到一定程度上优于维纳滤波的限制性逆滤波算法。但是这些滤波算法都没办法完全克服噪声的污染，随着加性的高斯噪声越来越大，结果也越来越差。

实验小结

本次实验，研究了全逆滤波，半径受限的逆滤波，维纳滤波，以及最后通过图像的频谱图设计出的频域幅度限制逆滤波算法。全逆滤波的抗干扰能力很差，但是如果是完全没噪声的图像，可以给出最好的修正效果。半径受限的逆滤波会损失高频的细节，但是可以有更好的抑制干扰的能力。维纳滤波要明显优于逆滤波，其通过修正信噪比 k 可以给出图像极好的恢复结果。

最后，通过图像的频谱图设计出的频域幅度限制逆滤波算法，可以通过对频谱幅值小的噪声有很好的过滤效果，尽管也会损失细节，但是可以通过设置噪声的下限来进行改善。在一定程度上，获得了比维纳滤波更好的结果。更换图片，重试发现算法有一定的普适性，具有不错的鲁棒性。

数字图像处理4：逆滤波及其变形和维纳滤波相关推荐

[数字图像处理]图像复原--逆滤波
1.逆滤波的问题点图像的老化,可以视为以下这样的一个过程.一个是退化函数的影响(致使图片模糊,褪色等),一个可加性噪声的影响. 用算式表示为前几篇博文,主要是介绍可加性噪声的 ...
数字图像处理--傅里叶(逆)变换
数字图像处理–傅里叶(逆)变换主要内容 (1)对一副图像进行缩放,显示原始图像和缩放后的图像,分别对其进行傅里叶变换,显示变换后结果: (2)对一副图像进行旋转,显示原始图像和旋转后的图像,分别对其 ...
领域平均法matlab代码实验,数字图像处理邻域平均法滤波实验报告matlab实现.doc...
数字图像处理邻域平均法滤波实验报告matlab实现数字图像处理实验报告实验三邻域平均法滤波学号姓名实验三邻域平均法滤波一.实验内容选取噪声较明显的图像,分别采用3*3.5*5.7* ...
数字图像处理之平滑滤波
数字图像处理之平滑滤波 by方阳版权声明:本文为博主原创文章,转载请标出转载地址 http://www.cnblog ...
matlab 数字图像滤波,数字图像处理 (基于Matlab) 滤波
<数字图像处理> 实验报告一.实验目的(不少于200字) 一.第一个实验用的是各种空间域的方式来滤波,也就是直接把图像和空间滤波器的模板做卷积,当然图像处理很重要的一个部分还有频域的处 ...
数字图像处理-频率域滤波原理
from:https://blog.csdn.net/forrest02/article/details/55510711?locationNum=15&fps=1 写在前面的话作者是一名在 ...
数字图像处理之频域滤波
前段时间看了很多的概念和知识,发现因为是走马观花的过了一遍,所以看得稀里糊涂的,然后许多地方混淆了概念,特别是关于图像频率域的部分的理解(包括图像频率域滤波之类的),所以下面总结一下这段时间重新看&l ...
系统学习数字图像处理之频域滤波
最近在看模板匹配,虽然很简单,但还是想认真过下基础,因此把信号处理频域相关的内容,接着图像处理再过一遍. 理论上,对连续变量t的连续函数f(t)的傅里叶变换为F(u),利用f(t)取样后的函数重建f( ...
数字图像处理——中值滤波及其改进算法
一.算法介绍中值滤波器是非线性滤波器的一个例子,它在保留图像特征方面非常有效. 但是,滤波器的窗口大小直接影响中值滤波器的性能. 较小的窗口保留了特征,但会导致噪声抑制的减少. 在较大窗口的情况下, ...

数字图像处理4：逆滤波及其变形和维纳滤波

逆滤波和维纳滤波

简介

大气湍流退化

三种滤波思想

直接全逆滤波

半径受限逆滤波

维纳滤波

代码1

主函数

全逆滤波和半径受限的全逆滤波

实验结果1

实验分析1

运动模糊和高斯噪声污染

运动模糊和高斯噪声设计思路

代码2

主函数

运动模糊

维纳滤波

实验结果2

实验思考

实验分析2

实验小结

数字图像处理4：逆滤波及其变形和维纳滤波相关推荐

最新文章

热门文章