%用函数imnoise2 生成具有表5.1 中的CDF 的随机数function R=imnoise2(type,M,N,a,b)if nargin ==1 a=0;b=1; M=1;N=1;elseif nargin ==3 a=0;b=1;endswitch lower(type) case 'uniform' R=a+(b-a)*rand(M,N); case 'gaussian' R=a+b*randn(M,N); case 'salt & pepper' if nargin <=3 a=0.05;b=0.05; end if (a+b)>1; error('The sum Pa+Pb must not exceed 1.') end R(1:M, 1:N) = 0.5; X=rand(M,N); c=find(X<=a); R(c)=0; u=a+b; c=find(X>a & X<=u); R(c)=1; case 'rayleigh' R=a+(-b*log(1-rand(M,N))).^0.5; case 'exponential' if nargin <=3; a=1; end if a<=0 error('Parameter a must be positive for exponential type.') end k=-1/a; R=k*log(1-rand(M,N)); case 'erlang' if nargin<=3 a=2;b=5; end if (b~=round(b)|b<=0) error('Parameter b must be a positive integer for Erlang') end k=-1/a; R=zeros(M,N); for j=1:b R=R+k*log(1-rand(M,N)); end otherwise error('unknown distribution type.')end
function image=changeclass(class,varargin)switch class case 'uint8' image=im2uint8(varargin{:}); case 'uint16' image=im2uint16(varargin{:}); case 'double' image=im2double(varargin{:}); otherwise error('Unsupported IPT data class.');end%%%%% spfilt 函数与表中列出的任何滤波器在空间域执行滤波。function f = spfilt(g,type,m,n,parameter)if nargin ==2 m=3;n=3;Q=1.5;d=2;elseif nargin == 5 Q=parameter;d=parameter;elseif nargin== 4 Q=1.5; d=2;else error ('wrong number of inputs');endswitch type case 'amean' w=fspecial('average',[m,n]); f=imfilter(g,w, 'replicate'); case 'gmean' f=gmean(g,m,n); case 'hmean' f=harmean(g,m,n); case 'chmean' %f=charmean(g,m,n,Q); f=charmean(g,m,n,Q); case 'median' f=medfilt2(g,[m n], 'symmetric'); case 'max' f=ordfilt2(g,m*n,ones(m,n),'symmetric'); case 'min' f=ordfilt2(g,1,ones(m,n), 'symmetric'); case 'midpoint' f1=ordfilt2(g,1,ones(m,n), 'symmetric'); f2=ordfilt2(g,m*n,ones(m,n), 'symmetric'); f=imlincomb(0.5,f1,0.5,f2); case 'atrimmed' if(d<0)|(d/2~=round(d/2)) error('d must be a nonnegative, even integer.') end f=alphatrim(g,m,n,d); otherwise error('Unknown filter type.')end
function f=gmean(g,m,n)inclass =class (g);g=im2double(g);warning off;f=exp(imfilter(log(g),ones(m,n),'replicate')).^(1/m/n);warning on;f=changeclass(inclass, f);
function f=harmean(g,m,n)inclass=class(g);g=im2double(g);f=m*n./imfilter(1./(g+eps),ones(m,n),'replicate');f=changeclass(inclass,f);
function f=charmean(g,m,n,q)inclass=class(g);g=im2double(g);f= imfilter(g.^(q+1),ones(m,n),'replicate');f=f./ (imfilter(g.^q,ones(m,n),'replicate')+eps);f=changeclass(inclass,f);
function f=alphatrim(g,m,n,d)inclass = class(g);g=im2double(g);f=imfilter(g,ones(m,n),'symmetric');for k=1:d/2 f=imsubtract(f,ordfilt2(g,k,ones(m,n),'symmetric'));endfor k=(m*n – (d/2)+1):m*n f=imsubtract(f,ordfilt2(g,k,ones(m,n),'symmetric'));endf=f/(m*n-d);f=changeclass(inclass,f);
%使用函数spfiltclear allclcf=imread('D:\pic\DIP3E_CH04\FigP0438(left).tif');[M,N]=size(f);R=imnoise2('salt & pepper',M,N,0.1,0);%被概率只有0.1 的胡椒噪声污染c=find(R==0);gp=f;gp(c)=0;figure, imshow(gp);
R=imnoise2('salt & pepper',M,N,0,0.1);c=find(R==1);gs=f;gs(c)=255;figure,imshow(gs)
fp=spfilt(gp,'chmean',3,3,1.5);%使用Q 为正值的反调和滤波器figure, imshow(gp);fs=spfilt(gs,'chmean',3,3,-1.5);figure, imshow(gs);
fpmax=spfilt(gp,'max',3,3); %使用最大最小滤波器figure, imshow(gp);fsmin=spfilt(gs,'min',3,3);figure, imshow(gs);
评论