案例7夜间图像去噪

案例7夜间图像去噪崔丽案例7夜间图像去噪共61页，您现在浏览的是第1页!（1）.噪声的扫描变换现在图像系统的输入光电变换都是先把二维图像信号扫描变换成一维电信号再进行处理加工。最后再将一维电信号变成二维图像信号。噪声也存在着同样的变换方式。（2）.噪声与图像的相关性使用光导摄象管的摄像机，可以认为，信号幅度和噪声幅度无关。而使用超正析摄像机的信号和噪声相关，黑暗部分噪声大，明亮部分噪声小，在数字图像处理技术中量化噪声是肯定存在的，它和图像相位有关，如图像内容接近平坦时，量化噪声呈现伪轮廓，但在此时图像信号中的随机噪声就会因为颤噪效应反而使量化噪声变得不那么明显。（3）.噪声的迭加性在串联图像传输系统中，各部分窜入噪声若是同类噪声可以进行功率相加，依次信噪比要下降。若不是同类噪声应区别对待，而且要考虑视觉检出特性的影响。但是因为视觉检出特性中的许多问题还没有研究情趣，所以也只能进行一些主观的评价试验。如空间频率特性不同的噪声迭加要考虑到视觉空间频谱的带通特性。而时间特性不同的噪声迭加就要考虑视觉滞留和其闪烁的特性等等。亮度和色度噪声的迭加一定要清楚视觉的彩色特性。而以上的这些都因为视觉特性的未获解决而无法进行分析。图像系统噪声的特点案例7夜间图像去噪共61页，您现在浏览的是第2页!噪声图像显示案例7夜间图像去噪共61页，您现在浏览的是第3页!从图像中找到一块平滑区域，根据统计数字来获得噪声方差，计算其均值和方差。Matlab命令：选取小块——roipoly(f)灰度直方图——imhist(f)自编程序_冈萨雷斯《数字图像处理》：几何中心矩的计算——statement噪声分析——imnoise2噪声方差估计案例7夜间图像去噪共61页，您现在浏览的是第4页!Statement函数解读：function[v,unv]=statmoments(p,n)%STATMOMENTS计算图像直方图的中心矩，输入：图像直方图p，n为矩的阶数返回值：v为中心矩，unv为没有平移中心的矩案例7夜间图像去噪共61页，您现在浏览的是第5页!estimate_noise.m文件估计出噪声强度，通过矩来计算用imnoise2，把噪声重现查看属于哪一类噪声。估计噪声实现案例7夜间图像去噪共61页，您现在浏览的是第6页!Matlab命令deconvregJ=deconvreg(I,PSF)J=deconvreg(I,PSF,NOISEPOWER)J=deconvreg(I,PSF,NOISEPOWER,LRANGE)程序解读：I是结果图像，PSF是滤波器（比原图像要小），NOISEPOWER是噪声强度（实数）去噪函数见ConstrainedLeastsquaresFiltering（CLSF.m）模型是g=hf+η还有待研究。案例7夜间图像去噪共61页，您现在浏览的是第7页!2维小波变换Matlab命令[CA,CH,CV,CD]=DWT2(X,‘wname’)或者[CA,CH,CV,CD]=DWT2(X,Lo_D,Hi_D)；X=IDWT2(CA,CH,CV,CD,‘wname’)或者X=IDWT2(CA,CH,CV,CD,Lo_R,Hi_R)；单独重构可用[]代替不用的部分。案例7夜间图像去噪共61页，您现在浏览的是第8页!小波分析去噪原理对信号进行小波分解针对小波系数进行处理方法：模极大值去噪阈值法去噪相关性去噪重构信号衡量标准案例7夜间图像去噪共61页，您现在浏览的是第9页!案例7夜间图像去噪共61页，您现在浏览的是第10页!案例7夜间图像去噪共61页，您现在浏览的是第11页!阈值的估计σ

为噪声方差，N为信号长度，T值估计应考虑平稳性与信噪比（自适应的选取）案例7夜间图像去噪共61页，您现在浏览的是第12页!带噪的一维信号nblocr1.mat；nblocr2.mat（带噪的方波信号）nbumpr1.mat；nbumpr2.mat；nbumpr3.mat；ndoppr1.mat；nelec.matnoischir.mat;noisbloc.mat；noisbump.mat；noisdopp.mat；noismima.mat；noissin.mat（高斯噪声）案例7夜间图像去噪共61页，您现在浏览的是第13页!如何人为加噪声1.读入信号或者图像X2.加入随机噪声，

随机噪声＝e*rand（M，N）;3.带有噪声的图像为

X＋随机噪声案例7夜间图像去噪共61页，您现在浏览的是第14页!进行软硬阈值函数运算命令Y=WTHRESH(X,SORH,T)输入参数：X为输入矩阵，SORH为's'(软阈值函数)，'h'(硬阈值函数)，T是阈值T返回值：用软(硬)阈值处理后的矩阵。[XD,CXD,LXD]=WDEN(X,TPTR,SORH,SCAL,N,'wname')1-D信号去噪函数输入参数：X为原始带噪信号，TPTR为阈值选择准则的字符串：见THSELECT命令

SORH:选择's'或者'h'，'h'是硬阈值，'s'为软阈值。见WTHRESH命令。SCAL定义多种方法的阈值重新调节:'one'：没有重新调节。'sln'基于层系数进行噪声估计来重新调节'mln'层层独立进行噪声估计来重新调节N为小波分解层数，'wname'为小波函数返回值：XD为去噪后的信号；CXD,LXD去噪后信号的小波分解的系数。另一种表达：[XD,CXD,LXD]=WDEN(C,L,TPTR,SORH,SCAL,N,'wname')与种不同之处在于用分解系数第N层的小波分解系数[C,L]代替直接的带噪信号X案例7夜间图像去噪共61页，您现在浏览的是第15页!案例7夜间图像去噪共61页，您现在浏览的是第16页!案例7夜间图像去噪共61页，您现在浏览的是第17页!案例7夜间图像去噪共61页，您现在浏览的是第18页!function[v,unv]=statmoments(p,n)%STATMOMENTS计算图像灰度直方图的统计中心矩。%参数说明：输入p为图像直方图，n为矩的阶数；%返回值V(1)为均值,V(2)为方差,V(3)为三阶中心矩，以此类推V(n)为第n阶中心矩。%返回值UNV为像素值取值再[0,255]非归一化的时候的矩的值。Lp=length(p);G=Lp-1;%直方图归一化p=p/sum(p);%灰度值随机变量及其归一化取值范围到[0,1].z=0:G;z=z./G;%均值m=z*p;%平移中心z=z-m;%计算中心矩v=zeros(1,n);v(1)=m;forj=2:nv(j)=(z.^j)*p;endifnargout>1%计算非中心矩unv=zeros(1,n);unv(1)=m.*G;forj=2:nunv(j)=((z*G).^j)*p;endend********************************************************************************************估计噪声参数案例7夜间图像去噪共61页，您现在浏览的是第19页!案例7夜间图像去噪共61页，您现在浏览的是第20页!Fourier变换去噪>>loadnoisbloc;>>X=noisbloc;>>Y=fft(X);>>plot(abs(Y));>>T=abs(Y);>>tt=find(T>1000);>>Y(tt)=0;>>IY=ifft(Y);>>subplot(2,1,1),plot(X),title('原始信号');subplot(2,1,2),plot(IY),title('Fourier变换去噪后信号');案例7夜间图像去噪共61页，您现在浏览的是第21页!子函数：用模极大值重构信号，采用交错投影法，见李世雄《小波变换域应用》functionint=Py(int,len);%Py投影functionint=Pg(int,lev,sr);%Pg投影functionw2=Pyg(w1,wp,sr);Pg和Py投影去噪函数functions=mden(f,lev,n,wf)模极大值去噪法案例7夜间图像去噪共61页，您现在浏览的是第22页!%该函数对一个区间进行Pgama投影%输入参数：int为极值点处差、lev为层号、sr为采样率%返回值为修正的区间functionint=Pg(int,lev,sr);T=length(int);ifT==2%只有两个点int=int;%不动elset=linspace(0,(T-1)/sr,T);%在0和(T-1)/sr中均匀插点，使最后点的个数为T。para=(([1,1;exp(2^(-lev)*t(T)),exp(-2^(-lev)*t(T))])\[int(1),int(T)]')';%求解对应Euler方程的解的参数，参看参考文献李世雄书中的（5.24）式ap=para(1);bt=para(2);int=ap.*exp(2^(-lev).*t)+bt.*exp(-2^(-lev).*t);%修正的区间end*****************************************************************案例7夜间图像去噪共61页，您现在浏览的是第23页!%模极大值法去噪：%参数C，ep，sr可以调节。%输入参数为：f为带去噪信号，lev为小波分解层数，n去噪重构算法的迭代次数，wf为小波函数名称functions=mden(f,lev,n,wf)pp=size(f);pp=pp(2);%所处理数据的长度

sr=360;%抽样率[Lo_D,Hi_D,Lo_R,Hi_R]=wfilters(wf);%二进小波变换[swa,swd]=swt(f,lev,Lo_D,Hi_D);

%%%%求二进小波变换的模极大值及其位置案例7夜间图像去噪共61页，您现在浏览的是第24页!wp=dw.*swd;%处理模极大值：最高层的模极大值点开始Dwp(:,lev)=wp(lev,:)';M=max(Dwp(:,lev));%模极大值的阈值用去噪函数ddencmp来获得。[Thr,sorh,keepapp]=ddencmp('den','wv',f);C=0.8;%阈值参数Thr=C*M/lev;%阈值计算，见参考文献”夏敏华，3mm波段脉冲雷达系统研究和小波去噪分析”。%将大于阈值Thr的最后一层的模极大值留下Dwp(:,lev)=Dwp(:,lev).*(abs(Dwp(:,lev))>Thr);%处理最后一层%模极大值的处理方式：%在尺度j上极大值点位置，构造一个搜索区域，%在尺度j-1中，极大值点落在该区域的点保留，其他的置0；ep=3;%该参数确定邻域forj=1:lev-2Dwp(:,lev-j)=wp(lev-j,:)';Dp(:,lev-j+1)=(Dwp(:,lev-j+1)~=0);DD=Dp(:,lev-j+1);forPd=ep:(length(Dp(:,lev-j+1))-ep);ifDD(Pd)==1fori=1:ep-1;DD(Pd-i)=1;endendendDp(:,lev-j+1)=DD;Dwp(:,lev-j)=Dwp(:,lev-j).*Dp(:,lev-j+1);Dp(:,lev-j)=(Dwp(:,lev-j)~=0);end案例7夜间图像去噪共61页，您现在浏览的是第25页!>>loadsumsin;%读入带噪信号>>f=[sumsinzeros(1,24)];%补零使得信号变成2的幂次，这里是信号长度变为1024。>>lev=4;%模极大值提取时候的分解层数>>n=20;%模极大值重构信号时候的迭代次数>>wf='db3';%采用的小波函数>>s=mden(f,lev,n,wf);%模极大值去噪算法主函数：案例7夜间图像去噪共61页，您现在浏览的是第26页!直接用ddencmp阈值，采用不同的传递方式，可以得到下面的图形案例7夜间图像去噪共61页，您现在浏览的是第27页!案例7夜间图像去噪共61页，您现在浏览的是第28页!mask_max=ones(1,length(T));T=[T((1:(n-1)),:);mask_max];%最后一层系数处理矩阵全置1。s1=s1.*T;f_n=iswt(swa,s1,wf);%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%画图：figure;%信号滤波前后比较。subplot(2,1,1);plot(f,'r');gridon;title('原始信号');subplot(2,1,2);plot(f_n);gridon;title('相关性去噪后信号')****************************************************************

案例7夜间图像去噪共61页，您现在浏览的是第29页!注意：相关性去噪也可以采用先估计噪声方差，进而阈值处理一下，再进行相关性去噪另外，相关性去噪的改进算法，还可以继续研究。案例7夜间图像去噪共61页，您现在浏览的是第30页!案例7夜间图像去噪共61页，您现在浏览的是第31页!傻瓜相机效果案例7夜间图像去噪共61页，您现在浏览的是第32页!几何矩定义:中心矩中心坐标：中心矩：

向量的几何矩案例7夜间图像去噪共61页，您现在浏览的是第33页!Imnoise2函数解读functionR=imnoise2(type,M,N,a,b)输入：type表示什么样的噪声，M，N为产生噪声的矩阵大小，a，b为相应的噪声参数；输出：R为大小为M*N矩阵Type型噪声，该噪声带有参数为a，b。案例7夜间图像去噪共61页，您现在浏览的是第34页!通过imnoise,imnoise2函数实现。具体程序见fnoise.m模糊化处理程序EstimateH，如何找到这个模糊核？产生带不同噪声的图像案例7夜间图像去噪共61页，您现在浏览的是第35页!在共享的文件夹里面，还有若干函数，可以供大家学习各种滤波。如spatialfilter.m(空间滤波)课上可完成学习inversefilter.m(反滤波和Wiener滤波)Inversewiener.m(Wiener滤波)Notchfilter.m(Notch滤波)OptimalNotch.m(最优Notch滤波)推荐有兴趣同学看书——冈萨雷斯《数字图像处理（Matlab版）》案例7夜间图像去噪共61页，您现在浏览的是第36页!Wavemenu——wavelet2D案例7夜间图像去噪共61页，您现在浏览的是第37页!阈值函数案例7夜间图像去噪共61页，您现在浏览的是第38页!缺点：硬阈值函数在w=T时是不连续的，用它重构信号时会产生振荡；软阈值函数虽然连续性好，但|w|>T时，它与w存在恒定的偏差，直接影响重构信号的性质。案例7夜间图像去噪共61页，您现在浏览的是第39页!案例7夜间图像去噪共61页，您现在浏览的是第40页!去噪阈值的获取ddencmp，thselect，wbmpen，wdcbm2实现小波去噪的函数有wden，wdencmp，wpdencmp，wthresh，wpthcoef和wthcoef2Wavemenu工具箱一维二维小波变换对应下的denoise。案例7夜间图像去噪共61页，您现在浏览的是第41页!带噪的二维图像noiswom.mat；noissi2d.mat;nbarb1.mat案例7夜间图像去噪共61页，您现在浏览的是第42页!阈值去噪Matlab命令：不同阈值选择方法下的阈值THR=THSELECT(X,TPTR)用TPTR字符串定义的阈值选取方法选择与输入信号X有关的阈值。输入参数：X为带噪信号TPTR为阈值选择准则的字符串：'rigrsure'，用Stein's无偏估计的原理进行的自适应阈值；'heursure'是启发式方法；'sqtwolog'通用阈值方法，阈值为；'minimaxi'最大最小阈值法。阈值选取方法的噪声模型为y=f(t)+e其中e是满足N(0,1)的白噪声处理非白噪声或者非校正的噪声可以用SCAL参数来重新调整，见wden命令。返回值：THR为阈值。案例7夜间图像去噪共61页，您现在浏览的是第43页!Wavemenu工具箱：wavelet1D

de-noise按键案例7夜间图像去噪共61页，您现在浏览的是第44页!Wavelet2Dde—noise按键案例7夜间图像去噪共61页，您现在浏览的是第45页!%在原始图像中加入Gaussian白噪声，椒盐噪声和乘性噪声。>>f=imread('original_pattern.tif');>>m=64/255;>>var=400/255^2;>>g_gauss=imnoise(f,'gaussian',m,var);>>d=0.05;%d表示噪声强度>>g_salt=imnoise(f,'salt&pepper',d);>>v=0.06;>>g_speckle=imnoise(f,'speckle',v);

%显示加噪图像>>figure(1)subplot(2,2,1),imshow(f),title('原始图像');subplot(2,2,2),imshow(g_gauss),title('Guassian噪声图像');subplot(2,2,3),imshow(g_salt),title('椒盐噪声图像');subplot(2,2,4),imshow(g_speckle),title('speckle噪声图像');

>>figure(2)title('直方图对比');subplot(2,2,1),imhist(f);%原始图像的灰度直方图subplot(2,2,2),imhist(g_gauss);%Guassian噪声图像的直方图subplot(2,2,3),imhist(g_salt);%椒盐噪声图像的直方图subplot(2,2,4),imhist(g_speckle);%speckle噪声图像的直方图产生模拟噪声案例7夜间图像去噪共61页，您现在浏览的是第46页!案例7夜间图像去噪共61页，您现在浏览的是第47页!>>g_Gaussian=imread('gaussian-noise.tif');%取出图像的一块>>[B,c,r]=roipoly(g_Gaussian);%画出原始图像f在B这块的直方图p>>p=imhist(g_Gaussian(B));%求和>>npix=sum(B(:));%计算2阶中心矩>>[v,unv]=statmoments(p,2)>>v

v=

0.78680.0033

>>unv

unv=

200.6315214.4707

%生成Gaussian噪声，均值为unv(1)，方差为unv(2)>>X=unv(1)+sqrt(unv(2))*randn(npix,1);>>M=max(p)

M=

86>>subplot(2,2,1),imshow(g_Gaussian),title('噪声图像');subplot(2,2,2),imshow(B),title('交互模式选取的区域');subplot(2,2,3),bar(p,1);subplot(2,2,4),hist(noise,130),axis([0,300,0,150])案例7夜间图像去噪共61页，您现在浏览的是第48页!************************************************************************************%计算两个信号f和g的峰值信噪比tfunctiont=PSNR(f,g)[m,n]=size(f);[m1,n1]=size(g);if(m1~=m|n1~=n)error('Signalsarenotthesamesize.',msg)endt=10*log(255^2/(1/(m*n)*sum((f-g).^2)));*************************************************************************************************************************************************************************%计算两个信号f和g的均方根误差tfunctiont=mse(f,g)[m,n]=size(f);[m1,n1]=size(g);if(m1~=m|n1~=n)error('Signalsarenotthesamesize.',msg)endt=1/(m*n)*sum((f-g).^2);**********************************************************峰值信噪比和均方误差的函数案例7夜间图像去噪共61页，您现在浏览的是第49页!案例7夜间图像去噪共61页，您现在浏览的是第50页!%单区间Py投影：对区间进行裁减即Py投影，返回裁剪后的区间信号%输入参数：int为单区间的点、len为区间中点的个数%返回值int为Py投影后的区间functionint=Py(int,len);ifsign(int(1))==sign(int(len))%若区间左右端点同号int=int.*(sign(int)==sign(int(1)));%区间上，只保留本身符号与左端点符号相同的那些inte=interp1([1,len],[int(1),int(len)],(1:len),'linear');%只用区间端点处的值，在区间上进行线性插值，得到区间的差值后的点，与原来区间点个数相同int=sign(int(1))*(abs(inte)-(abs(inte)-abs(int)).*((abs(inte)-abs(int))>0));%只保留区间上的点的模值比线性插值对应点模值小的那些点，其余的为0，符号与左端点一致。elsesgn=sign(int(len)-int(1));%若区间端点不同号，则取相减的符号。两极值点异号，中间有单调性intem=max([int(1),int(len)]);将端点中最大最小点分别赋值inten=min([int(1),int(len)]);%从区间端点开始，循环找寻该区间的极大值点和极小值点fori=1:len-2ifsign(int(i+1)-int(i))~=sgn%若差商符号不是sgnint(i+1)=int(i);%令二者相等，保证两个极值点之间的单调性endifint(i+1)>intem%比端点值大int(i+1)=intem;%则是端点值endifint(i+1)<intenint(i+1)=inten;%比端点值小，则是端点值endendend**********************************************************************************案例7夜间图像去噪共61页，您现在浏览的是第51页!******************************************************************************************%用于进行Pgama和Py投影，%输入参数：w1为加入细节信息改进的分解信号，wp为小波多层分解后的信号的极值点，sr为采样率functionw2=Pyg(w1,wp,sr);err=wp-w1.*(wp~=0);%极值点位置处原来与现在的差w2=zeros(size(wp));[r,c]=size(wp);%对每一层小波分别进行处理form=1:r%处理每一层fr=find(wp(m,:));%模极值点位置num_int=length(fr)-1;%用模极值点将区间分段，区间个数为极值点个数-1%先找到以模极大划分的区间,然后对每一区间进行Py投影forj=1:num_intint=w1(m,fr(j):fr(j+1));%极值点划分的区间中的点len=length(int);%区间中点的个数iflen>2w1(m,fr(j):fr(j+1))=Py(int,len);%区间内多于两个点，即中间有点，则进行Py投影endend%再逐一区间进行Pgama投影forj=1:num_intint=err(m,fr(j):fr(j+1));%在两极值点之间的原来与现在的差err(m,fr(j):fr(j+1))=Pg(int,m,sr);endw2(m,:)=w1(m,:)+err(m,:);%投影最终结果，对应参考文献【3】中的（5.28）式。end***************************************************************************************案例7夜间图像去噪共61页，您现在浏览的是第52页!

%%初始结果%dw为局部极大位置，wp为局部极大序列。dw=zeros(size(swd));zdw=dw;fdw=dw;%要找模极大值，把小波系数中大于0，和小于0的分别考虑。%小波系数大于0的赋值给zwzw=swd.*(swd>0);%留下左边元素比右边元素小的值的位置记为1，存入zdw中zdw=((zw(:,1:pp-1)-zw(:,2:pp))<0);%再计算1,0间隔的点，即找到模极大值的点的位置zdw(:,2:pp-1)=((zdw(:,1:pp-2)-zdw(:,2:pp-1))>0);%小波系数小于0的赋值给fwfw=swd.*(swd<0);%留下左边元素比右边元素大的值的位置记为1，存入fdw中fdw=((fw(:,1:pp-1)-fw(:,2:pp))>0);%再计算1,0间隔的点，即找到模极大值的点的位置fdw(:,2:pp-1)=((fdw(:,1:pp-2)-fdw(:,2:pp-1))>0);%将zdw和fdw中的模极大值点位置合并dw=zdw|fdw;%列最后一列保留dw(:,1)=1;dw(:,pp)=1;

%wp存放模极大值点的值案例7夜间图像去噪共61页，您现在浏览的是第53页!%层单独处理，在第二层极大值点位置上，保留层相应极大值点：Dwp(:,1)=wp(1,:);Dwp(:,1)=Dwp(:,1).*Dp(:,2);wp=Dwp;wp=wp';

%重构信号：s=swa(lev,:);%s:为待重建的信号w0=zeros(1,pp);[a,d]=swt(w0,lev,Lo_D,Hi_D);%从0开始w2=d;%重建高频forj=1:nw2=Pyg(d,wp,sr);%Py投影和Pgama投影w0=iswt(s,w2,Lo_R,Hi_R);%Pv投影[a,d]=swt(w0,lev,Lo_D,Hi_D);ends=iswt(swa(lev,:),w2,Lo_R,Hi_R);%计算重建信号figure,subplot(2,1,1)plot(f,'r');gridon;title('原始信号');subplot(2,1,2);plot(s);gridon;title('模极值去噪后信号')**********************************************************************************案例7夜间图像去噪共61页，您现在浏览的是第54页!案例7夜间图像去噪共61页，您现在浏览的是第55页!阈值法去噪法loadsumsin;f=[sumsinzeros(1,24)];lev=4;n=5;wf='db3';%软硬阈值去噪s_s=wden(f,'heursure','s','one',lev,wf);s_h=wden(f,'heursure','h','one',lev,wf);%计算PSNR值P1=PSNR(f,s_s)P2=PSNR(f,s_h)m1=mse(f,s_s)m2=ms