decompositionviterbi.m

关于分解的viterbi译码的matlab算法

%%making use of the result of the reference<A Decomposition approach to CPM>
%rect
%L=1

clear;
clc;

Rb=2400;
Fc=0.5*Rb;
N_sample=4;
Fs=N_sample*Rb;
% N=1000;
M=2;%the modulation level
h=1/2;
p=2;%p is changing with the h
T=1/Rb;
L=1;
N_observe=4;%make a decision when observe N_observe bit
error_ratio=[];
Number=100;

for snr=0:16
snr
NUM=0;%NUM累加每次循环时所得的误码数


for number_index=1:Number
snr
number_index
    
    
N_data=1000;
N_back=N_data;
Ik=round(rand(1,N_data));%the Data for transition
Ik_back=Ik;%back up the Data
Ik=[Ik,zeros(1,N_observe-1)];
Ik=2*Ik-1;
N=length(Ik);

Ik_temp=(Ik+(M-1))/2;%chang the Data form the space [-(M-1),......,M-1] to the space[0, 1,.....,M-1]


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%    modulation
state=0;
for i=1:N
    for k=0:(N_sample-1)
        
%         faik((i-1)*N_sample+k)=rem((2*pi*h*state+4*pi*h*Ik(i)*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi);
        faik((i-1)*N_sample+k+1)=rem((2*pi*h*state+4*pi*h*Ik_temp(i)*k/(2*N_sample)),2*pi);
    end
    
    state=rem(state+Ik_temp(i),p);
end
modu_out=exp(j*faik);
t=0:T/N_sample:(N*T-T/N_sample);
% Modu=cos(2*pi*Fc*t+faik);
% r_signal=modu_out;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% adding noise
% snr=10;
sigma=sqrt(0.5*N_sample*10^(-snr/10));
r_real=real(modu_out)+sigma*randn(1,N_sample*N);
r_imag=imag(modu_out)+sigma*randn(1,N_sample*N);
r_signal=r_real+j*r_imag;





%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% detecting using vitebi algorithm

% state=0;
Data=[];%Restore the decoded Data

state_temp=0;
r_temp=r_signal(1:N_sample);
U0=0;
U1=1;%U1 standing for the transmitting data when the ending state is 1
tao=0:T/N_sample:(T-T/N_sample);
state_temp0=0;%standing for the oringinal state is 0
state_temp1=1;
s0=exp(j*rem((2*pi*h*state_temp+4*pi*h*U0*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));
s1=exp(j*rem((2*pi*h*state_temp+4*pi*h*U1*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));%s0 standing for the ending state is 0
% CM0=sum(real(r_temp.*s0));
% CM1=sum(real(r_temp.*s1));

CM0=sum(real(r_temp).*real(s0)+imag(r_temp).*imag(s0));
CM1=sum(real(r_temp).*real(s1)+imag(r_temp).*imag(s1));

Data0=[0];
Data1=[1];% Restore the temp trellis,Data1 restore the data constructing the trellis ending at the state 1,and Data0 restore the data constructing 
          %the trellis ending at the state 0

for i=1:N_back
    
    for k=1:(N_observe-1)
        r_temp=r_signal((i+k-1)*N_sample+1:(i+k)*N_sample);
        
        U0=0;
        U1=1;
        
        s0_0=exp(j*rem((2*pi*h*state_temp0+4*pi*h*U0*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));%s0_1 standing for the signal format if the 
        s0_1=exp(j*rem((2*pi*h*state_temp1+4*pi*h*U1*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));% transmitting signal is 1 ,and the ending state is 0
        s1_0=exp(j*rem((2*pi*h*state_temp1+4*pi*h*U0*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));%and the same  with others
        s1_1=exp(j*rem((2*pi*h*state_temp0+4*pi*h*U1*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));
        
%         CM0_0=sum(real(r_temp.*s0_0));%CM1_0 standing for the metric of state 1 which comes from the state 0 ,and the same with the other
%         CM0_1=sum(real(r_temp.*s0_1));
%         CM1_0=sum(real(r_temp.*s1_1));
%         CM1_1=sum(real(r_temp.*s1_0));
        
        CM0_0=sum(real(r_temp).*real(s0_0)+imag(r_temp).*imag(s0_0));
        CM0_1=sum(real(r_temp).*real(s0_1)+imag(r_temp).*imag(s0_1));
        CM1_0=sum(real(r_temp).*real(s1_1)+imag(r_temp).*imag(s1_1));
        CM1_1=sum(real(r_temp).*real(s1_0)+imag(r_temp).*imag(s1_0));
        
        
        CM0_temp0=CM0+CM0_0;%CM1_temp0 standing for the metrics of the state 1 whose trellis comes from the state 0
        CM0_temp1=CM1+CM0_1;
        CM1_temp0=CM0+CM1_0;
        CM1_temp1=CM1+CM1_1;
        

        if CM0_temp0>CM0_temp1
            CM0=CM0_temp0;
            Data0_temp=[Data0 0];%CM0 standing for the metrics of the state 0,CM1 standing for the state 1
        else
            CM0=CM0_temp1;
            Data0_temp=[Data1 1];
        end
        
        if CM1_temp0>CM1_temp1
            CM1=CM1_temp0;
            Data1_temp=[Data0 1];
        else
            CM1=CM1_temp1;
            Data1_temp=[Data1 0];
        end
        Data0=Data0_temp;
        Data1=Data1_temp;
        
                
    end
            
    if CM0>CM1
        Data=[Data Data0(1)];
    else
        Data=[Data Data1(1)];
    end
    
    state_temp=rem(state_temp+Data(i),2);
    
    
    if state_temp==0%setting the original parameter based on the first state
        U0=0;
        U1=1;
        Data0=[0];
        Data1=[1];
    else
        U0=1;
        U1=0;
        Data0=[1];
        Data1=[0];
    end
        
        
        
        
    
    r_temp=r_signal((i)*N_sample+1:(i+1)*N_sample);
    s0=exp(j*rem((2*pi*h*state_temp+4*pi*h*U0*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));
    s1=exp(j*rem((2*pi*h*state_temp+4*pi*h*U1*(tao/(2*T))+pi*h*tao/T-2*pi*h*(tao/(2*T))),2*pi));
%     CM0=sum(real(r_temp.*s0));
%     CM1=sum(real(r_temp.*s1));
    
    CM0=sum(real(r_temp).*real(s0)+imag(r_temp).*imag(s0));
    CM1=sum(real(r_temp).*real(s1)+imag(r_temp).*imag(s1));
    
end

% Data=rem(Data+1,2);
[NUM_temp,Ratio]=biterr(Ik_back,Data);
NUM=NUM+NUM_temp;
NUM_temp
end%number
error_ratio_temp=NUM/(Number*N_data);
error_ratio=[error_ratio,error_ratio_temp];

end%snr

snr=0:16;
semilogy(snr,error_ratio)