1.研究噪声特性的必要性
本文的内容主要介绍了常见噪声的分类与其特性。将噪声建模,然后用模型去实现各式各样的噪声。
实际生活中的各种照片的老化,都可以归结为以下老化模型。
这个模型很简单,也可以直接用以下公式来表达。
在频域内,用以下公式区表示。
根据以上式子,可以看出,老旧照片的复原,主要分为两个任务,一个是去噪;另一个是去卷积,或者称为逆滤波,也就是将老化滤波器做反处理。
本文首先由噪声类型与其建模。随后的博文,会介绍几种基础的去噪方法和基础的逆滤波方法。
2.噪声的实现
2.1 评价用图像与其直方图
![](https://i-blog.csdnimg.cn/blog_migrate/2119c4453f6548102d3fdf265e2f379b.jpeg)
2.2 高斯噪声
- a = 0;
- b = 0.08;
- n_gaussian = a + b .* randn(M,N);
![](https://i-blog.csdnimg.cn/blog_migrate/e4a61b2b29039cee4cb4e033b6f24a59.jpeg)
2.3 瑞利噪声
瑞利噪声相比高斯噪声而言,其形状向右歪斜,这对于拟合某些歪斜直方图噪声很有用。
瑞利噪声的实现可以借由平均噪声来实现。如下所示。
这里的表示均值为0,方差为1的均匀分布的噪声。Matlab里,使用函数rand(M,N)就可以产生一个均值为0,方差为1的均匀噪声。
- a = -0.2;
- b = 0.03;
- n_rayleigh = a + (-b .* log(1 - rand(M,N))).^0.5;
![](https://i-blog.csdnimg.cn/blog_migrate/dc25914dbb9327810f52be99a792f745.jpeg)
2.4 伽马噪声
伽马噪声的分布,服从了伽马曲线的分布。伽马噪声的实现,需要使用b个服从指数分布的噪声叠加而来。指数分布的噪声,可以使用均匀分布来实现。
使用若干个(这里用b表示)均匀分布叠加,就可以得到伽马噪声。
当然,当b=1的时候,就可以得到指数噪声了。
- a = 25;
- b = 3;
- n_Erlang = zeros(M,N);
- for j=1:b
- n_Erlang = n_Erlang + (-1/a)*log(1 - rand(M,N));
- end
2.5 均匀噪声
如同前面所示,均匀噪声可以由函数rand(M,N)直接产生。
- a = 0;
- b = 0.3;
- n_Uniform = a + (b-a)*rand(M,N);
![](https://i-blog.csdnimg.cn/blog_migrate/347d7b075fe9a34d8f04e4209cea33c5.jpeg)
2.6 椒盐噪声
椒盐噪声也成为双脉冲噪声。在早期的印刷电影胶片上,由于胶片化学性质的不稳定和播放时候的损伤,会使得胶片表面的感光材料和胶片的基底欠落,在播放时候,产生一些或白或黑的损伤。事实上,这也可以归结为特殊的椒盐噪声。椒盐噪声的实现,需要一些逻辑判断。这里我们的思路是,产生均匀噪声,然后将超过阈值的点设置为黑点,或白点。当然,如果需要拟合电影胶片的损伤的话,可以选用别的类型噪声去拟合。
- a = 0.05;
- b = 0.05;
- x = rand(M,N);
- g_sp = zeros(M,N);
- g_sp = f;
- g_sp(find(x<=a)) = 0;
- g_sp(find(x > a & x<(a+b))) = 1;
3.总结
本文,实现的几类较为基本的噪声。并给出了其实现的方法,代码在下面。下一篇博文,会进行几个常用去噪滤波器的比较。- close all;
- clear all;
- clc;
- f = imread('./original_pattern.tif');
- f = mat2gray(f,[0 255]);
- [M,N] = size(f);
- figure();
- subplot(1,2,1);
- imshow(f,[0 1]);
- xlabel('a).Original image');
- subplot(1,2,2);
- x = linspace(-0.2,1.2,358);
- h = hist(f,x)/(M*N);
- Histogram = zeros(358,1);
- for y = 1:256
- Histogram = Histogram + h(:,y);
- end
- bar(-0.2:1/255:1.2,Histogram);
- axis([-0.2 1.2 0 0.014]),grid;
- xlabel('b).The Histogram of a');
- ylabel('Number of pixels');
- %% ---------------gaussian-------------------
- a = 0;
- b = 0.08;
- n_gaussian = a + b .* randn(M,N);
- g_gaussian = f + n_gaussian;
- figure();
- subplot(1,2,1);
- imshow(g_gaussian,[0 1]);
- xlabel('a).Ruselt of Gaussian noise');
- subplot(1,2,2);
- x = linspace(-0.2,1.2,358);
- h = hist(g_gaussian,x)/(M*N);
- Histogram = zeros(358,1);
- for y = 1:256
- Histogram = Histogram + h(:,y);
- end
- bar(-0.2:1/255:1.2,Histogram);
- axis([-0.2 1.2 0 0.014]),grid;
- xlabel('b).The Histogram of a');
- ylabel('Number of pixels');
- %% ---------------rayleigh-------------------
- a = -0.2;
- b = 0.03;
- n_rayleigh = a + (-b .* log(1 - rand(M,N))).^0.5;
- g_rayleigh = f + n_rayleigh;
- figure();
- subplot(1,2,1);
- imshow(g_rayleigh,[0 1]);
- xlabel('a).Ruselt of Rayleigh noise');
- subplot(1,2,2);
- x = linspace(-0.2,1.2,358);
- h = hist(g_rayleigh,x)/(M*N);
- Histogram = zeros(358,1);
- for y = 1:256
- Histogram = Histogram + h(:,y);
- end
- bar(-0.2:1/255:1.2,Histogram);
- axis([-0.2 1.2 0 0.014]),grid;
- xlabel('b).The Histogram of a');
- ylabel('Number of pixels');
- %% ---------------Erlang-------------------
- a = 25;
- b = 3;
- n_Erlang = zeros(M,N);
- for j=1:b
- n_Erlang = n_Erlang + (-1/a)*log(1 - rand(M,N));
- end
- g_Erlang = f + n_Erlang;
- figure();
- subplot(1,2,1);
- imshow(g_Erlang,[0 1]);
- xlabel('a).Ruselt of Erlang noise');
- subplot(1,2,2);
- x = linspace(-0.2,1.2,358);
- h = hist(g_Erlang,x)/(M*N);
- Histogram = zeros(358,1);
- for y = 1:256
- Histogram = Histogram + h(:,y);
- end
- bar(-0.2:1/255:1.2,Histogram);
- axis([-0.2 1.2 0 0.014]),grid;
- xlabel('b).The Histogram of a');
- ylabel('Number of pixels');
- %% ---------------Exponential-------------------
- a = 9;
- n_Ex = (-1/a)*log(1 - rand(M,N));
- g_Ex = f + n_Ex;
- figure();
- subplot(1,2,1);
- imshow(g_Ex,[0 1]);
- xlabel('a).Ruselt of Exponential noise');
- subplot(1,2,2);
- x = linspace(-0.2,1.2,358);
- h = hist(g_Ex,x)/(M*N);
- Histogram = zeros(358,1);
- for y = 1:256
- Histogram = Histogram + h(:,y);
- end
- bar(-0.2:1/255:1.2,Histogram);
- axis([-0.2 1.2 0 0.014]),grid;
- xlabel('b).The Histogram of a');
- ylabel('Number of pixels');
- %% ---------------Uniform-------------------
- a = 0;
- b = 0.3;
- n_Uniform = a + (b-a)*rand(M,N);
- g_Uniform = f + n_Uniform;
- figure();
- subplot(1,2,1);
- imshow(g_Uniform,[0 1]);
- xlabel('a).Ruselt of Uniform noise');
- subplot(1,2,2);
- x = linspace(-0.2,1.2,358);
- h = hist(g_Uniform,x)/(M*N);
- Histogram = zeros(358,1);
- for y = 1:256
- Histogram = Histogram + h(:,y);
- end
- bar(-0.2:1/255:1.2,Histogram);
- axis([-0.2 1.2 0 0.014]),grid;
- xlabel('b).The Histogram of a');
- ylabel('Number of pixels');
- %% ---------------Salt & pepper-------------------
- a = 0.05;
- b = 0.05;
- x = rand(M,N);
- g_sp = zeros(M,N);
- g_sp = f;
- g_sp(find(x<=a)) = 0;
- g_sp(find(x > a & x<(a+b))) = 1;
- figure();
- subplot(1,2,1);
- imshow(g_sp,[0 1]);
- xlabel('a).Ruselt of Salt & pepper noise');
- subplot(1,2,2);
- x = linspace(-0.2,1.2,358);
- h = hist(g_sp,x)/(M*N);
- Histogram = zeros(358,1);
- for y = 1:256
- Histogram = Histogram + h(:,y);
- end
- bar(-0.2:1/255:1.2,Histogram);
- axis([-0.2 1.2 0 0.3]),grid;
- xlabel('b).The Histogram of a');
- ylabel('Number of pixels');