tming_test1.m

通信中位同步的gardner算法的matlab仿真

clear
clc
close all

c1 = 1/128;
c2 = 1/2048;
c3 = 0.0001;
% 误差值为0.0894
c_1 = inline('-1/6*u^3+1/6*u');
c_2 = inline('1/2*u^3+1/2*u^2-u');
c_3 = inline('-1/2*u^3-u^2+1/2*u+1');
c_4 = inline('1/6*u^3+1/2*u^2+1/3*u');

dt = 1/88;
pi2 = 2*pi;

fid = fopen('rand_data.m','r');
source =  fread(fid)';
fclose(fid);
source1 = [source zeros(1,100)];
Qpsk_sig(1,:) = 2*source1(1:2:end)-1;
Qpsk_sig(2,:) = 2*source1(2:2:end)-1;

Expand_sig = zeros(2,size(Qpsk_sig,2)*8);
Expand_sig(1,1:8:end) = Qpsk_sig(1,:);
Expand_sig(2,1:8:end) = Qpsk_sig(2,:);

coeff_5 = firrcos(28,5.5,0.25,88,'rolloff','sqrt');

Filter_sig(1,:) = conv(Expand_sig(1,:),coeff_5);
Filter_sig(2,:) = conv(Expand_sig(2,:),coeff_5);

t = [0:dt:(size(Filter_sig,2)-1)*dt];
Md_sig_c =  Filter_sig(1,:).*cos(pi2*22*t) + Filter_sig(2,:).*sin(pi2*22*t);
Md_sig_s = -Filter_sig(1,:).*sin(pi2*22*t) + Filter_sig(2,:).*cos(pi2*22*t);   % 与希尔伯特变换等效,

clear Filter_sig;
clear Expand_sig;
clear Qpsk_sig;

sig_stor_c = resample(Md_sig_c,8802,8800);
sig_stor_s = resample(Md_sig_s,8802,8800);
dt1 = 1/88.02;
t = [0:dt1:(length(sig_stor_s)-1)*dt1];

sig_ave_pow = (norm(sig_stor_c)^2+norm(sig_stor_s)^2)/(length(sig_stor_s));
dc_rig = zeros(1,29);
ds_rig = zeros(1,29);
decode = zeros(2,2);
count = -32;
count1 = 0;
decode_out = [];
store = [];

%-----------------Timing_parament--------------------------------------------
it_flag = 0;
count2 = 0;
count3= 0;
u = 0;
timing_err = 0;
match_count = 0;
tnco_out = 1;
interp_regist = zeros(2,4);  % 滤波器的存储器，将4个点存储。如果得到插值命令就用这4个值插值
timing_regist_i = zeros(2,2);
timing_regist_q = zeros(2,2);
%     Ts = 1/22; 
xt1 = zeros(1,2);
xt2 = zeros(1,2);
xt3 = zeros(1,2);
contr_word = [0.25,0];

p_u = 0;
for i = 1:length(sig_stor_c)  
    count = count+1;
    Md_sig_c = sig_stor_c ;
    Md_sig_s = sig_stor_s;
    Md_sig = Md_sig_c + Md_sig_s;
    
    dc = Md_sig_c(i)*cos(pi2*22*t(i)) - Md_sig_s(i)*sin(pi2*22*t(i));
    ds = Md_sig_c(i)*sin(pi2*22*t(i)) + Md_sig_s(i)*cos(pi2*22*t(i));
    dc_rig(2:end) = dc_rig(1:end-1);
    ds_rig(2:end) = ds_rig(1:end-1);
    dc_rig(1) = dc;
    ds_rig(1) = ds;
    
    d_convert_c = dc_rig*coeff_5';
    d_convert_s = ds_rig*coeff_5';
    
    %-------------- timing---------------------------------------------
    if count>0
        %     store = [store d_convert_c];
    end
    count2 = count2+1;                                % 调整计数器，反映从滤波器出来的数据个数，（包括插值的数值）
    interp_regist(:,2:end) = interp_regist(:,1:end-1);  % 滤波器的存储器，将4个点存储。如果得到插值命令就用这4个值插值
    interp_regist(1,1) = d_convert_c;
    interp_regist(2,1) = d_convert_s;
    
    if it_flag == 1      
        if contr_word(1)==0
            u=0;
        else        
            u = (p_u)/contr_word(1);    
        end 
        interp_val_i = c_1(u)*interp_regist(1,4)+c_2(u)*interp_regist(1,3)... % 得到nco的指令，要求插值
            +c_3(u)*interp_regist(1,2)+c_4(u)*interp_regist(1,1);
        interp_val_q = c_1(u)*interp_regist(2,4)+c_2(u)*interp_regist(2,3)... % 得到nco的指令，要求插值
            +c_3(u)*interp_regist(2,2)+c_4(u)*interp_regist(2,1);
        %         store = [store interp_val_i];
    end
    if mod(count2,4)==1% & count>0                     % 数据滤波，将数据率下降到原来的1/4
        count3 = mod(count3+1,2);                     % count3,寻找准点采样值.           
        if count3 == 0
            timing_regist_i(1,2) = timing_regist_i(1,1);
            timing_regist_i(1,1) = d_convert_c;     
            timing_regist_q(1,2) = timing_regist_q(1,1);
            timing_regist_q(1,1) = d_convert_s;
            
            %             timing_err = (timing_regist_i(1,1))/(abs(timing_regist_i(1,2))+abs(timing_regist_i(1,1)));
        elseif count3 ==1  
            timing_regist_i(2,2) = timing_regist_i(2,1);
            timing_regist_i(2,1) = d_convert_c;
            timing_regist_q(2,2) = timing_regist_q(2,1);
            timing_regist_q(2,1) = d_convert_s;
            
            timing_err = (timing_regist_i(2,2)*(timing_regist_i(1,2)-timing_regist_i(1,1))+...
                timing_regist_q(2,2)*(timing_regist_q(1,2)-timing_regist_q(1,1)));  
            store = [store timing_err];
        end     
    end
    
    xt1(2) = xt1(1);
    xt1(1) = c3*timing_err;  % 相位误差
    xt2(2) = xt2(1);
    xt3(2) = xt3(1);
    contr_word(2)=contr_word(1);
    
    xt2(1) = ((c1+c2)*xt1(1)+c1*xt1(2));  % 环路滤波器的输出
    contr_word(1) = contr_word(2)-xt2(1); % 得到控制字。W(m)
    xt3(1) = xt2(1)+xt3(2)-0.25;                     % nco的输出
    p_u = tnco_out;
    tnco_out = tnco_out-(contr_word(1));   
    if tnco_out<0
        it_flag =1;
        tnco_out = mod(tnco_out,1);
        continue;
    end
    
    %-----------------------------------------------------------------------
    
    if it_flag ==1;
        it_flag=0;
        count1 = count1+1;
        decode(1,2) = decode(1,1);
        decode(2,2) = decode(2,1);
        decode(1,1) = interp_val_i;
        decode(2,1) = interp_val_q;
        
        if mod(count1,2)==0 & count1>6
            tempt1 = decode(1,1) + decode(1,2);
            tempt2 = decode(2,1) + decode(2,2);
            
            if tempt1>0 & tempt2>0
                decode_out = [decode_out 1 1];
            elseif tempt1>0 & tempt2<0
                decode_out = [decode_out 1 0];
            elseif tempt1<0 & tempt2>0
                decode_out = [decode_out 0 1];
            else 
                decode_out = [decode_out 0 0];
            end
            %             store = [store 0.5*sign(tempt1)];
        end
        
    end
    
end

err_num = 0;
err_bit = zeros(1,length(source));
for i = 1:length(source)
    if decode_out(i)~=source(i)
        err_num = err_num+1;
        err_bit(i) = 1;
    end
end
err_rate = err_num/length(source)
plot(store)