simu_ppgj_daopin.m

DTTB（中国数字地面电视标准）数字电视均衡仿真分析

%######## 基于国家数字电视地面广播（DTTB）标准,单载波调制系统的频偏估计 #######%
%##############  仿真条件：PN595帧头，16QAM调制，多径衰落信道模型  ############%
%#############################  利用导频估计频偏  #############################%

clear all;
close all;

len_head=595;                                                                   %帧头长度（帧头模式2）
len_carriers=3744;                                                              %有效信息符号长度
len_sysmes=36;                                                                  %系统信息符号长度
len_body=len_carriers+len_sysmes;                                               %帧体长度
len_frame=len_head+len_body;                                                    %信号帧长度
symbols=40;
bits=2;                                                                         %帧体每符号上的比特数
fd=2000;                                                                        %采样间隔2000Hz
Rs=7.56*1e6;                                                                    %符号率
ts=(1/7.56)*1e-6;                                                               %基带符号周期（1/7.56us）
fs=2*Rs;                                                                        %采样率

%----------------------------------TRANSMISSION---------------------------%
%--------------------------------------------------------------------------
%% 发送数据调制 %%
N=symbols*len_carriers*bits;                                                    %发送有效数据包含比特数
rand('state',0);                                                                %确定随机发生器种子
sendbits=round(rand(1,N));                                                      %产生数据位

%% 符号映射 %%
data_qam=fun_GB_QAM(sendbits,bits);                                             %调用16QAM调制函数
trans_data=reshape(data_qam,len_carriers,symbols).';                            %串并变换，每一行代表一个信息符号

%--------------------------------------------------------------------------
%% 系统信息调制 %%
sysmes=zeros(1,len_sysmes);                                                     %36位系统信息
sysmes(1:4)=[0 0 0 0];                                                          %4位系统信息：子载波个数C=1模式
bits_sysmes='00010001101000010100011110111100';                                 %32位系统信息：16QAM以及LDPC码率1，交织模式1
symmes(5:36)=bits_sysmes-48;                                                    
sysmes_bits=zeros(1,2*len_sysmes);                                              %进行I,Q路数据相同的4QAM调制
sysmes_bits(1:2:end)=sysmes;
sysmes_bits(2:2:end)=sysmes;
sysmes_qam=fun_GB_QAM(sysmes_bits,2);                                           %进行I,Q路数据相同的4QAM调制
trans_sysmes=repmat(sysmes_qam,symbols,1);                                      %表示symbols个符号的系统信息的矩阵

%--------------------------------------------------------------------------
%% 帧体数据处理 %%
trans_body=[trans_sysmes trans_data];
frame_day=[zeros(symbols,len_head) trans_body];                         
len=size(frame_day);                                                    
frame_day_serial=reshape(frame_day.',1,len(1)*len(2));                          %并串变换
frame_day_serial(1:2:end)=frame_day_serial(1:2:end)-1;                          %偶数符号，实部-1（第一个符号看作第0个符号，为偶数）
frame_day_serial(2:2:end)=frame_day_serial(2:2:end)+1;                          %奇数符号，实部+1（第二个符号看作第1个符号，为奇数）
frame_day_parallel=reshape(frame_day_serial.',len(2),len(1)).';                 %串并变换                         
body_process=frame_day_parallel(:,len_head+1:end);                              %处理完成后的帧体数据

%--------------------------------------------------------------------------
%% 已知帧头信号 %%
head_PN595=zeros(symbols,len_head);                                             %PN595序列生成
g=bin2dec('10010000001');                                                       %伪随机二进制序列生成多项式:1+x3+x10 
state=bin2dec('0000000001');                                                    %伪随机二进制序列生成寄存器初始状态
N=2^10-1;                                                                       %生成二进制序列长度
m=mgen(g,state,N);                                                              %完成帧头0到+1，1到-1的映射
PN595=m(1:len_head);                                                            %帧头数据取m序列的前595个数值
for i=1:len_head                                                                %调用mgen函数，生成1023位二进制随机序列
    if PN595(i)==0
        PN595(i)=-4.5*(1+1j);
    else
        PN595(i)=4.5*(1+1j);
    end
end
head_PN595=repmat(PN595,symbols,1);                                             %组成symbols个帧头
body_power=mean(mean(abs(body_process).^2,2))/(4.5^2*2);                        %帧体数据平均功率
head_power=body_power;                                                          %帧头数据的平均功率和帧体相同  
head_PN595=head_PN595.*sqrt(head_power);                                        %帧头发送数据

%--------------------------------------------------------------------------
%% 组帧 %%
trans_signal=[head_PN595 body_process];                                         %将每个符号的帧头和帧体组成一个完整信号帧    

%--------------------------------------------------------------------------
%% 并串变换 %%
trans_frame=reshape(trans_signal.',1,symbols*len_frame);                        %并串变换，将变换为串行数据

%--------------------------------------------------------------------------
%% 加入频偏 %%
data_srrc=trans_frame;
len_tx=length(data_srrc);
delta_f=6000;
t = (1:len_tx)*ts;
Tx_pp = data_srrc.*exp(1j*2*pi*t*delta_f);

%--------------------------------------------------------------------------
%% 加高斯白噪声 %%
SNR=15;                                                    
frame_passchannel = awgn(Tx_pp,SNR,'measured');

%--------------------------------------------------------------------------
%% 利用导频进行频偏估计 %%
Tx_sample=upsample(frame_passchannel,2);
rx=Tx_sample(1:10*7560);                                                        
len_pp=length(rx);
move_f=Rs/2;
t=(1:len_pp)*(1/fs);

data_mix=rx.*exp(-1j*2*pi*move_f*t);                                            %频谱搬移

figure(1);
subplot(211);
ff=(-len_pp/2:1:len_pp/2-1)*fs/len_pp;
plot(ff,abs(fftshift(fft(rx))));
title('频谱搬移前');
xlabel('频率/(Hz)');
ylabel('幅度');
subplot(212);
plot(ff,abs(fftshift(fft(data_mix))));
title('频谱搬移后');
xlabel('频率/(Hz)');
ylabel('幅度');

figure(2);
ff=(-len_pp/2:1:len_pp/2-1)*fs/len_pp;
plot(ff,abs(fftshift(fft(data_mix))));
axis([-50*fd 50*fd 0 4e4]);
grid on;
title('频偏估计');
xlabel('频率/(Hz)');
ylabel('幅度');