二维 最大类间方差.docVIP

clear; clc; tic [name path]=uigetfile(*.jpg); a=imread([path name]); %Read a image %a=imread(d:\image\clock.tif,tif); %a=imread(d:\image\zhuangjia.tiff,tiff); %a=a(:,:,1); % % a=M; % %load(source_N0.02.mat); % %a=X; %a=imread(syn1-g2.gif); %load(装甲车N_0.01.mat); %a=b; %a=noise_h; subplot(1,2,1); imshow(a); [m,n]=size(a); %b=imnoise(a,salt pepper,0.003); %b=imnoise(b,gaussian,0,0.0015); %b = IMNOISE(a,speckle,0.09); %b=a; a0=double(a); h=1; a1=zeros(m,n); % 计算平均领域灰度的一维灰度直方图 for i=1:m for j=1:n for k=-h:h for w=-h:h; p=i+k; q=j+w; if (p=0)|( pm) p=i; end if (q=0)|(qn) q=j; end a1(i,j)=a0(p,q)+a1(i,j); end end a2(i,j)=uint8(1/9*a1(i,j)); end end fxy=zeros(256,256); % 计算二维直方图 for i=1:m for j=1:n c=a0(i,j); d=double(a2(i,j)); fxy(c+1,d+1)=fxy(c+1,d+1)+1; end end % figure, % mesh(fxy); % title(二维灰度直方图); Pxy=fxy/m/n; P0=zeros(256,256); Ui=zeros(256,256); Uj=zeros(256,256); P0(1,1)=Pxy(1,1); for i=2:256 P0(1,i)=P0(1,i-1)+Pxy(1,i); end for i=2:256 P0(i,1)=P0(i-1,1)+Pxy(i,1); end for i=2:256 for j=2:256 P0(i,j)=P0(i,j-1)+P0(i-1,j)-P0(i-1,j-1)+Pxy(i,j); end end P1=ones(256,256)-P0; Ui(1,1)=0; for i=2:256 Ui(1,i)=Ui(1,i-1)+(1-1)*Pxy(1,i); end for i=2:256 Ui(i,1)=Ui(i-1,1)+(i-1)*Pxy(i,1); end for i=2:256 for j=2:256 Ui(i,j)=Ui(i,j-1)+Ui(i-1,j)-Ui(i-1,j-1)+(i-1)*Pxy(i,j); end end Uj(1,1)=0; for i=2:256 Uj(1,i)=Uj(1,i-1)+(i-1)*Pxy(1,i); end for i=2:256 Uj(i,1)=Uj(i-1,1)+(1-1)*Pxy(i,1); end for i=2:256 for j=2:256 Uj(i,j)=Uj(i,j-1)+Uj(i-1,j)-Uj(i-1,j-1)+(j-1)*Pxy(i,j); end end uti=0; utj=0; for i=1:256 for j=1:256 uti=uti+(i-1)*Pxy(i,j); utj=utj+(j-1)*Pxy(i,j); end end % for i=1:256 % for j=1:256 %