?? simu_whole_tbjh.m
字號:
%## 基于國家數字電視地面廣播（DTTB）標準,單載波調制系統的信道估計的原始模型 ##%
%################## 單載波應用模式：幀頭模式2，幀體模式C=1 ###################%
%############ 仿真條件：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=20;
bits=4; %幀體每符號上的比特數
fd=2000; %采樣間隔2000Hz
Rs=7.56*1e6; %符號率（7.56MHz）
ts=(1/7.56)*1e-6; %基帶符號周期（1/7.56us）
%----------------------------------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個幀頭
%--------------------------------------------------------------------------
%% 組幀 %%
trans_signal=[head_PN595 body_process]; %將每個符號的幀頭和幀體組成一個完整信號幀
%--------------------------------------------------------------------------
%% 并串變換 %%
trans_frame=reshape(trans_signal.',1,symbols*len_frame); %并串變換，將變換為串行數據
%--------------------------------------------------------------------------
%% 加入頻偏 %%
len_tx=length(trans_frame);
delta_f=6000;
t = (1:len_tx)*ts;
Tx_pp = trans_frame.*exp(j*2*pi*t*delta_f);
%% 信道模型:多徑衰落信道 %%
X7=Tx_pp;
r=3; %多徑數
a=[0.2 0.3 0.4]; %多徑的幅度
d=[30 50 80]; %各徑的延遲
channel1=zeros(size(X7));
channel1(1+d(1):length(X7))=a(1)*X7(1:length(X7)-d(1));
channel2=zeros(size(X7));
channel2(1+d(2):length(X7))=a(2)*X7(1:length(X7)-d(2));
channel3=zeros(size(X7));
channel3(1+d(3):length(X7))=a(3)*X7(1:length(X7)-d(3));
Tx_data=X7+channel1+channel2+channel3;
%--------------------------------------------------------------------------
%% 加高斯白噪聲 %%
SNR=20; %信噪比
frame_passchannel = awgn(Tx_data,SNR,'measured');
data_simu=frame_passchannel(1401:end); %幀起始位置偏移
%-----------------------------------RECEPTION-----------------------------%
%-------------------------------------------------------------------------%
%--------------------------------------------------------------------------
%% 選取接收數據 %%
symbols=8; %取8個信號幀長度的數據
rx=data_simu(1:symbols*len_frame);
%--------------------------------------------------------------------------
%% 符號估計 %%
N=len_frame;
G=len_head;
GI=len_head+len_sysmes; %各信號幀相同數據塊長度
T=1:3*N+2*G;
for i=1:length(T) %FFT窗移動的距離
data1=rx(i:GI+i-1); %data1對應一個符號的前GI個數據
data2=rx(N+i:N+GI+i-1); %data2對應一個符號的最后GI個數據
s(i)=data2*data1'; %求滑動相關值
ctb_cor(i)=abs(s(i)); %相關幅值
end
%--------------------------------------------------------------------------
%% 符號同步圖 %%
figure(1);
plot(T,ctb_cor);
xlabel('符號數');
ylabel('相關值');
grid on;
%% 符號校正 %%
[Max_cor start_place]=max(ctb_cor(1:N)); %start_place對應一個符號的第一個數據
rx_ctb=data_simu(start_place:symbols*len_frame+start_place-1); %以start_place所在處為數據起點
% disp(['符號起始位置：',num2str(start_place)]); %顯示粗同步偏移
% %--------------------------------------------------------------------------
% %% 幀頭PN序列相關圖 %%
% figure(2);
% rr_head=data_simu(start_place:start_place-1+595);
% plot(abs(xcorr(rr_head,PN595)));
% xlabel('PN序列位數');
% ylabel('相關值');
% grid on;
%--------------------------------------------------------------------------
%% 串并變換 %%
receiv_frame=reshape(rx_ctb,len_frame,symbols); %并串變換，將變換為串行數據
receiv_frame=receiv_frame.';
%--------------------------------------------------------------------------
%% 去幀頭PN序列 %%
receiv_head=receiv_frame(:,1:len_head); %接收到的幀頭序列
receiv_body=receiv_frame(:,len_head+1:end); %接收到的幀體數據
%--------------------------------------------------------------------------
%% 頻偏估計 %%
head_mean=mean(abs(receiv_head),2);
M=head_mean(2);
re_head=receiv_head(2,:);
L=10;
corr_head = re_head.*conj(PN595)./M;
corr_FFT=fft(corr_head,L*len_body);
%--------------------------------------------------------------------------
%% 頻偏估計作圖 %%
figure(3);
ff=(-len_body*L/2:1:len_body*L/2-1)*fd/L;
pp_cor=fftshift(abs(corr_FFT));
plot(ff,pp_cor);
title('頻偏估計');
xlabel('頻偏/Hz');
ylabel('幅值');
axis([-2e5 2e5 0 4000]);
%--------------------------------------------------------------------------
%% 頻偏校正 %%
pp_cor=fftshift(abs(corr_FFT));
[Max_corr x] = max(pp_cor);
xpp=(x-1-len_body*L/2)*fd/L; %估計出的頻偏值
tt=(1:length(rx_ctb))*ts;
rx_ppjz=rx_ctb.*exp(-1j*2*pi*xpp*tt);
%--------------------------------------------------------------------------
%% FFT %%
receiv_frame=reshape(rx_ppjz,len_frame,symbols); %并串變換，變換為串行數據
receiv_frame=receiv_frame.';
receiv_head=receiv_frame(:,1:len_head); %接收到的幀頭序列
receiv_body=receiv_frame(:,len_head+1:end); %接收到的幀體數據
head_FFT=fft(receiv_head,len_head,2);
body_FFT=fft(receiv_body,len_body,2);
%--------------------------------------------------------------------------
%% 信道估計 %%
PN595_FFT=fft(head_PN595,len_head,2);
PN595_FFT=PN595_FFT(1:symbols,:);
H_head=head_FFT./PN595_FFT; %由幀頭的PN序列估計信道
h_head=ifft(H_head,len_head,2);
H_body=fft(h_head,len_body,2);
% %--------------------------------------------------------------------------
% %% ZF頻域均衡 %%
% W_zf=1./H_body;
% body_equa_zf=body_FFT.*W_zf; %迫零均衡（ZF）
% body_IFFT_zf=ifft(body_equa_zf,len_body,2); %IFFT
%
% body_power=mean(mean(abs(body_IFFT_zf).^2,2)); %幀體數據平均功率
% head_power=body_power/(4.5^2*2); %幀頭數據的平均功率和幀體相同
% head_jh=head_PN595.*sqrt(head_power);
% frame_IFFT_zf=[head_jh body_IFFT_zf]; %加幀頭
% data_out_zf=reshape(frame_IFFT_zf.',1,symbols*len_frame);
%--------------------------------------------------------------------------
%% MMSE頻域均衡 %%
H_mmse_fm=abs(H_body).^2+1/10^(SNR/10);
W_mmse=conj(H_body)./H_mmse_fm;
body_equa_mmse=body_FFT.*W_mmse; %最小均方誤差均衡（MMSE）
body_IFFT_mmse=ifft(body_equa_mmse,len_body,2); %IFFT
%--------------------------------------------------------------------------
%% 均衡前后自相關圖 %%
figure(4);
Nc=len_frame;
subplot(211);
cyy=xcorr(receiv_body(1,:),Nc); %生成長度為2*Nc+1的自相關序列
tt=(-Nc*ts:ts:Nc*ts)*1e6;
plot(tt,20*log10(abs(cyy)/abs(cyy(Nc+1)))); %接收數據自相關圖
axis([-20 20 -40 0]);
zoom on;
xlabel('時間（us)');
ylabel('幅度（dB）');
title('通過信道后數據自相關');
subplot(212);
CYY=xcorr(body_IFFT_mmse(1,:),Nc); %生成長度為2*Nc+1的自相關序列
tt=(-Nc*ts:ts:Nc*ts)*1e6;
plot(tt,20*log10(abs(CYY)/abs(CYY(Nc+1)))); %均衡后數據自相關圖
axis([-20 20 -40 0]);
zoom on;
xlabel('時間（us)');
ylabel('幅度（dB）');
title('MMSE均衡后數據自相關');
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -