clear;
[H2,G2,n1,k1]=hammgen(3);
U=[1 0 0 1; 1 0 1 0; 0 0 1 1; 0 1 0 1];
C=mod(U*G2,2);
[C_m,C_n]=size(C);
for i=1:C_m
        C(i,C_n+1)=0;
end   
for i=1:C_m
    for j=1:C_n
        C(i,C_n+1)=mod((C(i,j)+C(i,C_n+1)),2);
    end
end
C=C';
C=mod(C*G2,2);
[C_m,C_n]=size(C);
for i=1:C_m
        C(i,C_n+1)=0;
end   
for i=1:C_m
    for j=1:C_n
        C(i,C_n+1)=mod((C(i,j)+C(i,C_n+1)),2);
    end
end
C=C';
%--------------------------------------------------------------------
[R_inital,noise_std]=modulate_demodulate(C,C_m);%%%%BPSK調制與解調
for i=1:length(R_inital)
    if R_inital(i)>0
        Y_inital(i)=1;
    else 
        Y_inital(i)=0;
    end
end
[C_m,C_n]=size(C);
n=1;
for i=1:C_m
    for j=1:C_n
        R(i,j)=R_inital(n);  %接收信號R
        Y(i,j)=Y_inital(n);  %R的硬判決值Y
        n=n+1;
    end
end
%---------------------------------------------------------------------
m_max=8;
p=3;
q=2^p;
for i=1:C_m-1
        rr_s(i,:)=R(i,1:(C_n-1));
        rr_h(i,:)=Y(i,1:(C_n-1));
end
%------------------------------------------------------------
%Rout=SISO(r_s,r_h,m_max,p,H2,noise_std);
 R_soft=rr_s;
 R_hard=rr_h;
 [M,N]=size(rr_s);
 %********************************************************
for m=1:m_max
%----------initialize---------------------------------
%  行譯碼
 for ii=1:M
     r_s=R_soft(ii,:);
     r_h=R_hard(ii,:);
     RR(ii,:)=r_s,
     B(m)=m/m_max,
    switch(m)
       case 1,
            RRi(ii,:)=RR(ii,:);
            W(ii,:)=zeros(1,N);,
            A(m)=0,
            D(ii,:)=chase(p,r_s,r_h,H2),
            W(ii,:)=D(ii,:).*B(m);
            RRout(ii,:)=RRi(ii,:)+W(ii,:);
       otherwise,
          A(m)=(m*(noise_std^2))/var(W(ii,:)),
           for j=1:N
           RRi(ii,j)=r_s(j)+W(ii,j).*A(m),
           end
           D_f(ii,:)=D(ii,:);
          D(ii,:)=chase(p,r_s,r_h,H2),
          xishu=0,
          for j=1:N
              xishu=xishu+((RRi(ii,j)-D_f(ii,j))^2-(RRi(ii,j)-D(ii,j))^2),
          end
           xishu=xishu/4;
         for i=1:N
             if D(ii,i)==D_f(ii,i)
                 W(ii,i)=B(m)*D(ii,i);
                 RRout(ii,i)=RRi(ii,i)+W(ii,i);
             else
                 W(ii,i)=xishu*D(ii,j)-RR(ii,i);
                 RRout(ii,i)=RRi(ii,i)+W(ii,i);
             end
         end
     end
 end
 %------------------------------------------------------------------------
%  解交織
R_soft=R_soft';
R_hard=R_hard';
W=W';
RR=RR';
RRi=RRi';
RRout=RRout';
D=D';
% D_f=D_f';
%-------------------------------------------------------------------------
%---------------------------列譯碼-----------------------
 for ii=1:N
     r_s=R_soft(ii,:);
     r_h=R_hard(ii,:);
     RR(ii,:)=r_s;
     B(m)=m/m_max;
    switch(m)
       case 1,
            RRi(ii,:)=RR(ii,:);
            W(ii,:)=zeros(1,M);,
            A(m)=0,
            D(ii,:)=chase(p,r_s,r_h,H2),
            W(ii,:)=D(ii,:).*B(m);
            RRout(ii,:)=RRi(ii,:)+W(ii,:);
       otherwise,
          A(m)=(m*(noise_std^2))/var(W(ii,:)),
           for j=1:M
           RRi(ii,j)=r_s(j)+W(ii,j).*A(m),
           end
           D_f(ii,:)=D(ii,:);
          D(ii,:)=chase(p,r_s,r_h,H2),
          xishu=0,
          for j=1:M
              xishu=xishu+((RRi(ii,j)-D_f(ii,j))^2-(RRi(ii,j)-D(ii,j))^2),
          end
           xishu=xishu/4;
         for i=1:M
             if D(ii,i)==D_f(ii,i)
                 W(ii,i)=B(m)*D(ii,i);
                 RRout(ii,i)=RRi(ii,i)+W(ii,i);
             else
                 W(ii,i)=xishu*D(ii,j)-RR(ii,i);
                 RRout(ii,i)=RRi(ii,i)+W(ii,i);
             end
         end
     end
 end
%-------------------------------------------------------------
 
%  解交織
R_soft=R_soft';
R_hard=R_hard';
W=W'
RR=RR';
RRi=RRi';
RRout=RRout';
D=D';
% D_f=D_f';
end
%------------------------------------------------------------ 
 for i=1:M
     for j=1:N
if(RRout(i,j)>0) Rout(i,j)=1;
               else Rout(i,j)=0;
               end
     end
 end