%SOM神經網絡，20051210
%實驗類內距離和類內離差陣的關系
%======每5類DOA相差9，檢查對DOA的適應情況=DOA變化從3倍sita變化=== 識別率為   100 ========================

format long;
 P=[35 41 77 83 99;  
 9310 9416 9822 9728 9834;
 4.1 4.2 4.3 4.4 4.5]'; 
%  4.05 4.1 4.15 4.2 4.25]'; %  
   P=[65 71 77 83 89;  
 9810 9816 9822 9828 9834;
 4.1 4.2 3.3 4.4 4.5]'; 
D1=[];
Nu=length(P);%數據位數
Pro=[];
PP1=zeros(3,1);PPPP=[];
PP1=[];
TT=zeros(3,1);
TT=[]
for i=1:Nu
sigma_D=3;
sigma_F=3;
sigma_PW=0.1;
DOA=P(i,1);
Fa=P(i,2);
PW=P(i,3);
sita_MX1=sqrt(sigma_D.^2+sigma_F.^2+sigma_PW.^2);%門限
M1=4;%每類的個數
 T=[0 0 0 0 1;0 0 0 1 0;0 0 1 0 0;0 1 0 0 0;1 0 0 0 0]  ;
for j=1:M1 
    D=normrnd(DOA,sigma_D);%隨機產生高斯分布的數據
    F=normrnd(Fa,sigma_F);
    F1=normrnd(PW,sigma_PW);
    PP=[D F F1];
%     PP=[D F P(i,3) P(i,4)];
PP1=[PP1;PP];% 每種類型中有M1組數據
TT=[TT;T(i,:)];
end
for j=1:2*M1 
    D=normrnd(DOA,sigma_D);%隨機產生高斯分布的數據
    F=normrnd(Fa,sigma_F);
    F1=normrnd(PW,sigma_PW);
    PP=[D F F1];
%     PP=[D F P(i,3) P(i,4)];
PPPP=[PPPP;PP];% 每種類型中有M1組數據
end
end 

PP1=PPPP([1:M1 2*M1+1:3*M1 4*M1+1:5*M1 6*M1+1:7*M1 8*M1+1:9*M1],:);%訓練樣本
PP2=PPPP([1*M1+1:2*M1 3*M1+1:4*M1 5*M1+1:6*M1 7*M1+1:8*M1 9*M1+1:10*M1],:);%測試樣本
p{1,1}=PP1(1:4,:);p{2}=PP1(5:8,:);p{3}=PP1(9:12,:);p{4}=PP1(13:16,:);p{5}=PP1(17:20,:);
m1=mean(p1);m2=mean(p2);m3=mean(p3);m4=mean(p4);m5=mean(p5);m=mean(PP1);
j1=0;
for i=1:4
    j1=j1+(p1(i,:)-m1)*(p1(i,:)-m1)';
end
j2=0;
for i=1:4
    j2=j2+(p2(i,:)-m2)*(p2(i,:)-m2)';
end
j3=0;
for i=1:4
    j3=j3+(p3(i,:)-m3)*(p3(i,:)-m3)';
end
j4=0;
for i=1:4
    j4=j4+(p4(i,:)-m4)*(p4(i,:)-m4)';
end
j5=0;
for i=1:4
    j5=j5+(p5(i,:)-m5)*(p5(i,:)-m5)';
end
Jw=j1+j2+j3+j4+j5;
JB=(m1-m)*(m1-m)'+(m2-m)*(m2-m)'+(m3-m)*(m3-m)'+(m4-m)*(m4-m)'+(m5-m)*(m5-m)';
J1=JB/Jw
%%%%%%
j1=zeros(3,3);
for i=1:4
    j1=j1+(p1(i,:)-m1)'*(p1(i,:)-m1);
end
j2=zeros(3,3);
for i=1:4
    j2=j2+(p2(i,:)-m2)'*(p2(i,:)-m2);
end
j3=zeros(3,3);
for i=1:4
    j3=j3+(p3(i,:)-m3)'*(p3(i,:)-m3);
end
j4=zeros(3,3);
for i=1:4
    j4=j4+(p4(i,:)-m4)'*(p4(i,:)-m4);
end
j5=zeros(3,3);
for i=1:4
    j5=j5+(p5(i,:)-m5)'*(p5(i,:)-m5);
end
Sw=1/20*(j1+j2+j3+j4+j5);
SB=1/5*((m1-m)'*(m1-m)+(m2-m)'*(m2-m)+(m3-m)'*(m3-m)+(m4-m)'*(m4-m)+(m5-m)'*(m5-m));
l=[1,0,0;0,1,0;0,0,1];
J2=sum(sum(l.*(SB*inv(Sw))))