tbjh_20080407_runprog.m

和基于欧标(DVB_T)的均衡器设计的实际数据验证的MATLAB程序

    delta_t(j)=Tao(j)-Td(j);                                                    %符号定时偏移小数部分：采样定时偏移
end 

% %--------------------------------------------------------------------------
% %% 精符号估计（利用连续导频） %%
% receiv_jtb=receiv_zpjz;
% Tao=zeros(1,symbols);
% for j=1:symbols
%     Tao(j)=0;
%     pilot_cons=CP+1;                                                            %连续导频位置向量
%     pilot_value=receiv_jtb(j,pilot_cons);                                       %连续导频位置值的向量
%     p=length(CP); 
%     for n=1:p-1
%         delta_Fi=angle(pilot_value(n+1)/pilot_value(n));
%         delta_k=pilot_cons(n+1)-pilot_cons(n);
%         Tao(j)=Tao(j)+delta_Fi./delta_k;
%     end
%     Tao(j)=N/(2*pi)*Tao(j)/(p-1);                                               %第j个符号定时同步误差
%     Td(j)=fix(Tao(j));                                                          %符号定时偏移整数部分：符号定时偏移
%     delta_t(j)=Tao(j)-Td(j);                                                    %符号定时偏移小数部分：采样钟偏移 
% end

%--------------------------------------------------------------------------
%% 符号定时偏差校正（反馈到FFT窗口） %%
delta_jtb=fix(mean(Tao));                                                        %符号定时偏差估计值
                              
num_jtb=num_jtb+1;
if num_jtb==10   
    break;   
end

end  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     %结束精符号同步循环

%--------------------------------------------------------------------------
%% 采样定时偏差校正（反馈到A/D采样前端） %%
delta_sample= mean(delta_t);                                                    %采样定时偏差估计值
delta_cyjz=round(Num_cy*delta_sample);

num_cytb=num_cytb+1;
if num_cytb==10
    break;
end

end  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     %结束采样定时同步循环
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%--------------------------------------------------------------------------
%% 粗符号同步作图 %%
figure(1);
title('粗符号同步');
T=1:3*N+2*G;
subplot(211);
plot(T,ctb_cor);
xlabel('采样点');
ylabel('相关值');
grid on;

subplot(212);
plot(T,e);
xlabel('采样点');
ylabel('小数频偏估计值');
grid on;

disp('粗符号同步偏移为：');
disp(num2str(start_place-1)); 

%--------------------------------------------------------------------------
%% 精符号同步作图 %%
figure(2);
title('精符号同步');
for i=1:length(T)                                                               %FFT窗移动的距离
    data1=rx_jtb(i:G+i-1);                                                      %取GI长个数据存data1
    data2=rx_jtb(N+i:N+G+i-1);                                                  %取相距IFFT_bin_length的GI长个数据存data2
    s_jtb(i)=data1*data2';                                                      %求互相关值    
    cor_jtb(i)=abs(s_jtb(i));           
end

subplot(211);
plot(T,cor_jtb);
xlabel('载波数');
ylabel('相关值');
grid on;

for i=1:length(T)
    deltaF_jtb(i)=1/(2*pi)*angle(conj(s_jtb(i)));
    r1=rx_jtb(i:G+i-1);
    r2=rx_jtb(i+N:G+i+N-1);
    e_jtb(i)=-1/(2*pi)*atan(imag(r1*r2')/(real(r1*r2')));
end

subplot(212);
plot(T,e_jtb);
xlabel('采样点');
ylabel('小数频偏估计值');
grid on;

disp('符号定时同步偏移为：');
disp(num2str(delta_symbol));  

disp('采样定时同步偏移为：');
disp(num2str(delta_sample));                          
                
%##########################################################################
%################################  均衡部分  ##############################
%##########################################################################
%% 均衡前判断分散导频所在符号位置 %%
receiv_use=receiv_jtb(:,(1:N_carrier)+172);                                     %从解调符号中舍去作OFDM所加零点位
s_cor=zeros(4,2*1705-1);                                                        %存放4次相关的相关值
p_max=zeros(1,4);                                                               %存放4次相关的最大值
for i=1:4
    symbols_first=receiv_use(1,:);                                              %选取要处理数据第一个符号
    pilot_select=flag1_cycle(i,:);                                              %在4个导频符号周期内依次选取导频所在符号
    cor=xcorr(abs(symbols_first),abs(pilot_select));                            %对上述两组数据做互相关
    s_cor(i,:)=abs(cor);
    [Max_cor p_max(i)]=max(s_cor(i,:));                                         %得到相关峰值位置               
end
[Min pilot_start]=min(abs((p_max-1705)));                                       %取4个峰值位置离理想值1705最近的点

%--------------------------------------------------------------------------
%% 导频位置校正（位置对应） %% 
flag1_jz=[flag1_cycle(pilot_start:4,:);flag1_rep(1:symbols-(5-pilot_start),:)];
flag_jz=flag1_jz;                                                               %flag为信息位与所有导频位标志数组，1代表导频（分散和连续导频），0代表信息位
for frame_no=0:symbols-1
    flag_jz(frame_no+1,CP+1)=1;                                                 %所生成的flag_jz位symbols行1705列，与一个周期的符号相对应
end

%--------------------------------------------------------------------------
%% 二维内插 %%
pilot_receiv=receiv_use.*flag1_jz;                                              %从解调符号中提出分散导频位
pilot_h=pilot_receiv./pilot;                                                    %接收导频除以本地导频，得到导频位的信道频域信息H(k)
receiv_h=[pilot_h(:,1) zeros(symbols,Kmax-1) pilot_h(:,N_carrier)];
for j=4:3:Kmax                                                                  %时间方向内插
    interp_time=pilot_h(:,j).';                                                
    st=find(interp_time~=0);                                                    %找时间内插起始位置点
    if length(st)>1                                                             %当含两个以上导频点时内插
        interp_time=interp1(st(1):4:symbols,interp_time(st(1):4:symbols),1:symbols);
    end
    receiv_h(:,j)=interp_time.';
end

for i=4:symbols-3                                                               %频率方向内插（只取中间完整内插符号）
    interp_freq=receiv_h(i,:);
    interp_freq=interp1(1:3:N_carrier,interp_freq(1:3:N_carrier),1:N_carrier);
    receiv_h(i,:)=interp_freq;
end

%--------------------------------------------------------------------------
%% 低通滤波 %%
IFFT_H=ifft(receiv_h,N_carrier,2);
h_filter=[IFFT_H(:,1:GI),zeros(symbols,N_carrier-GI)];
H_ce=fft(h_filter,N_carrier,2);
receiv_h=H_ce(:,1:N_carrier);

%--------------------------------------------------------------------------
%% 频域均衡 %%
receiv_equa=receiv_use./receiv_h;                                               %对数据作均衡（二维内插）
use_mean=mean(mean(abs(receiv_use(4:symbols-3,:)),2));         
equa_mean=mean(mean(abs(receiv_equa(4:symbols-3,:)),2));
multiple_equa=use_mean/equa_mean;                                               %数据幅值放大倍数
receiv_equa=receiv_equa.*multiple_equa;                                        

%--------------------------------------------------------------------------
%% IFFT变换 %%    
IFFT_CE=zeros((symbols-6),FS);
IFFT_CE(:,1+172:N_carrier+172)=receiv_equa(4:symbols-3,:);                      %有效数据子载波前有172个空子载波
Y1=ifft(IFFT_CE,FS,2);
    
%--------------------------------------------------------------------------
%% 加循环前缀保护间隔 %%
for i=1:FS
    Y0(:,i+GI)=Y1(:,i);
end
for i=1:GI
    Y0(:,i)=Y1(:,i+FS-GI);    
end

%--------------------------------------------------------------------------
%% 并串变换 %%
Y=reshape(Y0.',1,(symbols-6)*Ns);

%--------------------------------------------------------------------------
figure(3);                                                                      %均衡前信号自相关图
subplot(211);
Nc=FS;
y=carriers_jtb(4:symbols-3,:);                                                  
y=reshape(y.',1,(symbols-6)*Ns);                                                %均衡前时域串行数据
cyy1=xcorr(y,Nc);                                                               %生成长度为2*Nc+1的自相关序列
tt=(-Nc*ts:ts:Nc*ts)*1e6;
plot(tt,20*log10(abs(cyy1)/abs(cyy1(Nc+1))));
axis([-250 250 -80 0]);
xlabel('时间（us)');
ylabel('幅度（dB）');
title('均衡前自相关');

subplot(212);                                                                   %均衡前信号自相关图
cyy2=xcorr(Y,Nc);                                                               %生成长度为2*Nc+1的自相关序列
tt=(-Nc*ts:ts:Nc*ts)*1e6;
plot(tt,20*log10(abs(cyy2)/abs(cyy2(Nc+1))));
axis([-250 250 -80 0]);
xlabel('时间（us)');
ylabel('幅度（dB）');
title('均衡后自相关');