例4.1 使用均匀随机数生成具有规定发布的随机数
matlab实现:
x=rand(1,1000); % 产生随机数x~U (0,1)
y=sqrt(-8.*log(1-x)); % 产生服从2
M=80; % 设置划分y 值的统计区间数目
[N,Y]=hist(y,M); % 统计落入每个区间的次数赋值给N ,区间坐标赋值给Y
N=N/1000/((max(y)-min(y))/M); % 由统计的区间内次数求PDF
bar(Y,N); % 绘制PDF 直方图
title('随机数y 的PDF 估计');
ylabel('p(x)');
xlabel('y');
例4.2 使用数imnoise2所生成的直方图:
matlab实现:
r = imnoise2('gaussian',100000,1,0,1);
hist(r,50)
%r = imnoise2('gaussian',100000,1,0,1);
%hist(r,50)
%以下显示所有类型的直方图,没有赋值的取默认值
R=imnoise2('gaussian',100000,1,0,1);
subplot(3,2,1),hist(R,50),title('高斯')%用hist显示直方图
R=imnoise2('salt & pepper',100000,10);
subplot(3,2,2),hist(R,50),title('焦盐')
R=imnoise2('lognormal',100000,10);
subplot(3,2,3),hist(R,50),title('对数正态')
R=imnoise2('rayleigh',100000,1);
subplot(3,2,4),hist(R,50),title('瑞利')
R=imnoise2('exponential',100000,1);
subplot(3,2,5),hist(R,50),title('指数')
R=imnoise2('erlang',100000,1);
subplot(3,2,6),hist(R,50),title('厄兰')
例4.3
C=[0 64;0 128;32 32;64 0;128 0;-32 32];
[r,R,S]=imnoise3(512,512,C);
subplot(3,2,1),imshow(S,[])%显示频谱
subplot(3,2,2),imshow(r,[])%显示空间正弦噪声模式
C=[0 32;0 64;16 16;32 0;64 0;-16 16];%改变冲击位置,观察频谱和空间正弦噪声模式变化
[r,R,S]=imnoise3(512,512,C);
subplot(3,2,3),imshow(S,[])
subplot(3,2,4),imshow(r,[])
C=[6 32;-2 2];%改变冲击位置,观察空间正弦噪声模式变化
[r,R,S]=imnoise3(512,512,C);
subplot(3,2,5),imshow(r,[])
A=[1 5];%使用非默认振幅向量A,观察空间正弦噪声模式变化
[r,R,S]=imnoise3(512,512,C,A);
subplot(3,2,6),imshow(r,[])
例4.5 使用函数spfilt
matlab实现:
% 胡椒噪声
f = imread('E:\图像处理\Digital-Image-Processing-main\matlab教材示例图片\dipum_images_ch05\5.tif');
[M,N] = size(f);
R = imnoise2('salt & pepper',M,N,0.1,0);
c = find(R==0);
gp = f;
gp(c) = 0;
% 盐粒噪声
R = imnoise2('salt & pepper',M,N,0,0.1);
c = find(R==1);
gs = f;
gs(c) = 255;
% 过渡胡椒噪声的较好方法是使用Q为正值的反调和滤波器
fp = spfilt(gp,'chmean',3,3,1.5);
% 过渡盐粒噪声的较好方法是使用Q为负值的反调和滤波器
fs = spfilt(gs,'chmean',3,3,-1.5);
% 使用最大和最小滤波器可以得到类似的结果
fpmax = spfilt(gp,'max',3,3);
fsmin = spfilt(gs,'min',3,3);
subplot(3,2,1);imshow(gp);title('仅被胡椒噪声污染');
subplot(3,2,2);imshow(gs);title('仅被盐粒噪声污染');
subplot(3,2,3);imshow(fp);title('Q=1.5反调和滤波器处理胡椒噪声');
subplot(3,2,4);imshow(fs);title('Q=-1.5反调和滤波器处理盐粒噪声');
subplot(3,2,5);imshow(fpmax);title('最大滤波器处理胡椒噪声');
subplot(3,2,6);imshow(fsmin);title('最小滤波器处理盐粒噪声');
例4.6
matlab实现:
f = imread('E:\图像处理\Digital-Image-Processing-main\matlab教材示例图片\dipum_images_ch05\5.tif');
g = imnoise(f,'salt & pepper',.25);
f1 = medfilt2(g,[7 7],'symmetric');
f2 = adpmedian(g,7);
figure,subplot(131),imshow(f);
subplot(132),imshow(f1);
subplot(133),imshow(f2);
例4.7 模糊的带噪图像的建模
matlab是实现:
%checkerboard产生测试板图像,第一个参数为每个正方形一边的像素数,第二个为参数行数,第三个为列数
f=checkerboard(8);%产生一个一边为8个正方形的测试板
subplot(2,2,1),imshow(f),title('原图像')
PSF=fspecial('motion',7,45);%用fspecial处理图像的运动模糊,PSF为空间滤波器
gb=imfilter(f,PSF,'circular');%减少边界效应
subplot(2,2,2),imshow(gb);title('运动模糊图像')
noise=imnoise(zeros(size(f)),'gaussian',0,0.001);%高斯噪声
subplot(2,2,3),imshow(noise);title('运动模糊图像')
g=gb+noise;%添加运动模糊的高斯噪声
subplot(2,2,4),imshow(pixeldup(g,8),[]);title('运动模糊+高斯噪声图像')%大图像运算过程过慢,选用小图像来节省时间
例4.8 使用函数deconvwnr复原模糊的带噪图像
matlab实现:
%生成模糊图像
f=checkerboard(8);
PSF=fspecial('motion',7,45);
gb=imfilter(f,PSF,'circular');
subplot(2,2,1),imshow(gb),title('模糊图像')
%生成逆滤波图像
noise=imnoise(zeros(size(f)),'gaussian',0,0.001);
g=gb+noise;
fr1=deconvwnr(g,PSF);%使用deconvwnr函数生成逆滤波图像
subplot(2,2,2),imshow(fr1),title('逆滤波图像')
Sn=abs(fft2(noise)).^2;%噪声功率谱
nA=sum(Sn(:))/prod(size(noise));%噪声平均功率
Sf=abs(fft2(f)).^2;%图像功率谱
fA=sum(Sf(:))/prod(size(f));%图像平均功率
R=nA/fA;%计算噪声和信号比例
fr2=deconvwnr(g,PSF,R);%使用deconvwnr函数生成常数比率的维纳滤波图像
subplot(2,2,3),imshow(fr2),title('常数比率的维纳滤波图像')
NCORR=fftshift(real(ifft2(Sn)));%噪声自相关函数
ICORR=fftshift(real(ifft2(Sf)));%原图像自相关函数
fr3=deconvwnr(g,PSF,NCORR,ICORR);%使用deconvwnr函数生成自相关函数的维纳滤波图像
subplot(2,2,4),imshow(fr3),title('自相关函数的维纳滤波图像')
例4.9 使用函数radon
matlab实现:
g1 = zeros(600,600);
g1(100:500,250:350) =1;
g2 = phantom('Modified Shepp-Logan',600);
figure,subplot(2,2,1),imshow(g1);
subplot(2,2,2),imshow(g2);
theta = 0:0.5:179.5;
[R1,xp1] = radon(g1,theta);
[R2,xp2] = radon(g2,theta);
R1 = flipud(R1');
R2 = flipud(R2');
subplot(2,2,3),imshow(R1,[],'XData',xp1([1 end]),'YData',[179.5 0]);
axis xy
axis on
xlabel('\rho'),ylabel('\theta');
subplot(2,2,4),imshow(R2,[],'XData',xp2([1 end]),'YData',[179.5 0]);
axis xy
axis on
xlabel('\rho'),ylabel('\theta');