基于MATLAB,QPSK,16QAM,64QAM,软判决和硬判决对比

小心豪猪

已于 2024-03-18 14:50:43 修改

阅读量2.8k

点赞数 32

文章标签：信号处理 matlab 算法概率论矩阵

于 2024-03-18 13:48:23 首次发布

本文链接：https://blog.csdn.net/m0_72270392/article/details/136806526

版权

实验原理：

实验步骤：

确定实验比特数：要求BER到1e-4量级，那么实验比特数要达到10^6，采用分组设定1000轮，每轮1000个比特,规一化发送符号功率，确定N0。

%基本参数

M=[2 4 6];

Nmb_Symbol=1000;

Block=1000;

Es_N0=1:2:30;

n0=1./(10.^(Es_N0./10));

编码代码逻辑：归一化发送功率，以16qam举例，假设最近的4个点相隔为2.则能量计算为：

E=（1*1+1*1）*4+（3*3+1*1）*8+（3*3+3*3）*4=160

平均能量=E/16=10,平均功率=sqrt（10）

以gray_64pam举例，6位bit表示1个符号。前三位表示实部（即X坐标），第一位表示正负，后两位表示幅度（假设为正，则幅度有1,3,5,7），那么可以看做是4加减1,3，其中加减用第2位比特表示，第3位比特表示幅值大小.4+sign（bit（2）)*(2+bit(4)）。

%对不同进制下符号进行格雷编码

Bit_Code_2=2*randi([0 1],1,2*Nmb_Symbol)-1;

Bit_Gray_Code_4=sqrt(1/10)*Bit_Code_4(:,1:4:4*Nmb_Symbol).*(2+Bit_Code_4(:,2:4:4*Nmb_Symbol))+1i*sqrt(0.1)*Bit_Code_4(:,3:4:4*Nmb_Symbol).*(2+Bit_Code_4(:,4:4:4*Nmb_Symbol)) ;

%通过awgn信道的两种方法

% White_Noise=sqrt(n0(i)/2)*randn(1,Nmb_Symbol)+1i*sqrt(n0(i)/2)*randn(1,Nmb_Symbol);

% Bit_Gray_Code_2_Add_Noise=Bit_Gray_Code_2+White_Noise;

Bit_Gray_Code_2_Add_Noise=awgn(Bit_Gray_Code_2,Es_N0(i));

硬判决逻辑判断星座点正负，取绝对值，将半区图像移在中心点，再次判断正负，重复此过程。

%硬判决

Dec_Code_2=zeros(1,2*Nmb_Symbol);

Dec_Code_6(:,1:6:6*Nmb_Symbol)=sign(real(Bit_Gray_Code_6_Add_Noise));

Dec_Code_6(:,2:6:6*Nmb_Symbol)=sign(abs(real(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42));

Dec_Code_6(:,3:6:6*Nmb_Symbol)=sign(abs(abs(real(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42))-2*sqrt(1/42));

Dec_Code_6(:,4:6:6*Nmb_Symbol)=sign(imag(Bit_Gray_Code_6_Add_Noise));

Dec_Code_6(:,5:6:6*Nmb_Symbol)=sign(abs(imag(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42));

Dec_Code_6(:,6:6:6*Nmb_Symbol)=sign(abs(abs(imag(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42))-2*sqrt(1/42));

软解调原理：简单来说，软判决是对应每一位比特来说的。每一位比特的对数似比（LLR）就是软判决值。通过LLR值与0比较，得到对应比特位是bit 0还是bit 1。每一位比特LLR的计算可以通过一系列计算简化为，星座图上欧式距离的差。

首先构造信号集合的划分，存储所有的分组取值：

for divide_rounds=1:1:4

Collection_Symbols_4=[ones(1,2^(2-divide_rounds)),-ones(1,2^(2-divide_rounds))];

Constellation_Bits_4(divide_rounds,:)=kron(ones(1,2^(divide_rounds-1)),Collection_Symbols_4);

clear Collection_Symbols_4

End

以M=4为例，目标是获取可以判别4次比特（b1b2b3b4）正负的集合，代码逻辑：kron（a，b）函数可以对矩阵a（mn），矩阵b（pq）的所有可能取值进行组合，生成一个大小为(mn)*(nq)的矩阵。

计算信号中每个符号与星座点的欧式距离：

for loop_llr_det=1:1:Nmb_Symbol

for loop_xzd=1:1:2^2

Euclidean_Distance_2(loop_llr_det,loop_xzd)=-(abs(Bit_Gray_Code_2_A dd_Noise(loop_llr_det)-Constellation_Symbols_2(loop_xzd)))^2/n0(i);

end

for loop_llr_bit=1:1:2

ind_temp_plus_2 = find(Constellation_Bits_2(loop_llr_bit,:)==1);

ind_temp_minus_2=find(Constellation_Bits_2(loop_llr_bit,:)==-1);

a_temp_plus_2 = Euclidean_Distance_2(loop_llr_det,ind_temp_plus_2);

a_temp_minus_2=Euclidean_Distance_2(loop_llr_det,ind_temp_minus_2);

对数似然比计算公式

log(sum(exp(a_temp_plus_2)))-log(sum(exp(a_temp_minus_2)))

llr_det_2(1,(loop_llr_det-1)*2+loop_llr_bit)=max(a_temp_plus_2)+log(sum(exp(a_temp_plus_2-max(a_temp_plus_2))))-max(a_temp_minus_2)-log(sum(exp(a_temp_minus_2-max(a_temp_minus_2))));

llr_det_2(1,(loop_llr_det-1)*2+loop_llr_bit)=min(a_temp_minus_2)-min(a_temp_plus_2);

最终可化简为：

%%llr_det_2(1,(loop_llr_det-1)*2+loop_llr_bit)=max(a_temp_plus_2)-max(a_temp_minus_2);

计算误码：
%记录每轮误码个数，两种

err_num_QPSK_Hard =err_num_QPSK_Hard +length(find(Dec_Hard_Code_2~=Bit_Code_2));

err_num_QPSK_Soft =err_num_QPSK_Soft +length(find(Dec_Soft_Code_2~=Bit_Code_2));

end

BER_QPSK_Hard =[BER_QPSK_Hard err_num_QPSK_Hard/(Nmb_Symbol*2*Block)];

BER_QPSK_Soft =[BER_QPSK_Soft err_num_QPSK_Soft/(Nmb_Symbol*2*Block)];

%没改进

% for indx = 1:length(Dec_Hard_Code_2)

% if Dec_Hard_Code_2(indx)~=Bit_Code_2(indx)

% err_num_round_QPSK_Hard=err_num_round_QPSK_Hard+1;

绘图：

subplot(311)

semilogy(Es_N0,BER_QPSK_Hard,'-.X',Es_N0,BER_QPSK_Soft,'-.o');

title('QPSK调制下，软判决和硬判决BER和Es_N0的关系')

xlabel('Es_N0/dB')

ylabel('BER')

grid on

现象：软判决和硬判决并无差异。

解释：软判决的后一级经常跟着软输入软输出(SISO，Soft Input Soft Output)的模块，以利用比特似然比中的软信息，常见的卷积码或是Turbo码的译码器都有相应的软输入译码算法。反之，若经过软判决后，没有后续的模块，而直接是源信息，则软判决相对硬判决不能带来任何增益。

完整代码：

clear all
close all

%基本参数
M=[2 4 6];
Nmb_Symbol=10000;
Block=100;
Es_N0=1:1:30;
n0=1./(10.^(Es_N0./10));

Euclidean_Distance_2=zeros(1000,2^2);
Euclidean_Distance_4=zeros(1000,2^4);
Euclidean_Distance_6=zeros(1000,2^6);

% Constellation_Symbols_2=zeros(1,2^2);
% Constellation_Symbols_4=zeros(1,2^4);
% Constellation_Symbols_6=zeros(1,2^6);

Constellation_Bits_2=zeros(2,2^2);
Constellation_Bits_4=zeros(4,2^4);
Constellation_Bits_6=zeros(6,2^6);

% Constellation_P_2=zeros(1,2^2);
% Constellation_P_4=zeros(1,2^4);
% Constellation_P_6=zeros(1,2^6);

llr_det_2=zeros(1,2*Nmb_Symbol);
llr_det_4=zeros(1,4*Nmb_Symbol);
llr_det_6=zeros(1,6*Nmb_Symbol);

BER_QPSK_Hard=[];
BER_16QAM_Hard=[];
BER_64QAM_Hard=[];
BER_QPSK_Soft=[];
BER_16QAM_Soft=[];
BER_64QAM_Soft=[];

%软判决，符号集合的划分
for divide_rounds=1:1:2
    Collection_Symbols_2=[ones(1,2^(2-divide_rounds)),-ones(1,2^(2-divide_rounds))];
    Constellation_Bits_2(divide_rounds,:)=kron(ones(1,2^(divide_rounds-1)),Collection_Symbols_2);
    clear Collection_Symbols_2
end
for divide_rounds=1:1:4
    Collection_Symbols_4=[ones(1,2^(4-divide_rounds)),-ones(1,2^(4-divide_rounds))];
    Constellation_Bits_4(divide_rounds,:)=kron(ones(1,2^(divide_rounds-1)),Collection_Symbols_4);
    clear Collection_Symbols_4
end
for divide_rounds=1:1:6
    Collection_Symbols_6=[ones(1,2^(6-divide_rounds)),-ones(1,2^(6-divide_rounds))];
    Constellation_Bits_6(divide_rounds,:)=kron(ones(1,2^(divide_rounds-1)),Collection_Symbols_6);
    clear Collection_Symbols_6
end

%记录星座点及其功率
Constellation_Symbols_2=sqrt(0.5)*Constellation_Bits_2(1,:)+1i*sqrt(0.5)*Constellation_Bits_2(2,:);
Constellation_Symbols_4=sqrt(1/10)*(Constellation_Bits_4(1,:).*(2+Constellation_Bits_4(2,:)))+1i*sqrt(1/10)*(Constellation_Bits_4(3,:).*(2+Constellation_Bits_4(4,:)));
Constellation_Symbols_6=sqrt(1/42)*(Constellation_Bits_6(1,:).*(4+Constellation_Bits_6(2,:).*(2+Constellation_Bits_6(3,:))))+1i*sqrt(1/42)*(Constellation_Bits_6(4,:).*(4+Constellation_Bits_6(5,:).*(2+Constellation_Bits_6(6,:))));

Constellation_P_2=abs(Constellation_Symbols_2).^2;
Constellation_P_4=abs(Constellation_Symbols_4).^2;
Constellation_P_6=abs(Constellation_Symbols_6).^2;
                                        
%在同一信噪比下，进行1000轮，每次生成1000符号
for i=1:1:length(n0)
    err_num_QPSK_Hard=0;
    err_num_16QAM_Hard=0;
    err_num_64QAM_Hard=0;

    err_num_QPSK_Soft=0;
    err_num_16QAM_Soft=0;
    err_num_64QAM_Soft=0;
    for round=1:1:Block

        %对不同进制下符号进行格雷编码
        Bit_Code_2=2*randi([0 1],1,2*Nmb_Symbol)-1;
        Bit_Code_4=2*randi([0 1],1,4*Nmb_Symbol)-1;
        Bit_Code_6=2*randi([0 1],1,6*Nmb_Symbol)-1;

        Bit_Gray_Code_2=sqrt(1/2)*Bit_Code_2(:,1:2:2*Nmb_Symbol)+1i*sqrt(0.5)*Bit_Code_2(:,2:2:2*Nmb_Symbol) ;
        Bit_Gray_Code_4=sqrt(1/10)*Bit_Code_4(:,1:4:4*Nmb_Symbol).*(2+Bit_Code_4(:,2:4:4*Nmb_Symbol))+1i*sqrt(0.1)*Bit_Code_4(:,3:4:4*Nmb_Symbol).*(2+Bit_Code_4(:,4:4:4*Nmb_Symbol)) ;
        Bit_Gray_Code_6=sqrt(1/42)*(Bit_Code_6(:,1:6:6*Nmb_Symbol).*(4+Bit_Code_6(:,2:6:6*Nmb_Symbol).*(2+Bit_Code_6(:,3:6:6*Nmb_Symbol))))+1i*sqrt(1/42)* (Bit_Code_6(:,4:6:6*Nmb_Symbol).*(4+Bit_Code_6(:,5:6:6*Nmb_Symbol).*(2+Bit_Code_6(:,6:6:6*Nmb_Symbol)))) ;
        
%       M进制下星座图展示
%       for i=1:400
%           subplot(311)
%           plot(real(Bit_Gray_Code_2(i)),imag(Bit_Gray_Code_2(i)),'-ko');
%           title('QPSK星座图')
%           hold on
%           subplot(312)
%           plot(real(Bit_Gray_Code_4(i)),imag(Bit_Gray_Code_4(i)),'-ko');
%           title('16QAM星座图')
%           axis([-1.5,1.5,-1.5,1.5]);
%           hold on
%           subplot(313)
%           plot(real(Bit_Gray_Code_6(i)),imag(Bit_Gray_Code_6(i)),'-ko');
%           axis([-1.5,1.5,-1.5,1.5]);
%           title('64QAM星座图')
%           hold on
%       end
        %通过awgn信道
%         White_Noise=sqrt(n0(i)/2)*randn(1,Nmb_Symbol)+1i*sqrt(n0(i)/2)*randn(1,Nmb_Symbol);
%         Bit_Gray_Code_2_Add_Noise=Bit_Gray_Code_2+White_Noise;
%         Bit_Gray_Code_4_Add_Noise=Bit_Gray_Code_4+White_Noise;
%         Bit_Gray_Code_6_Add_Noise=Bit_Gray_Code_6+White_Noise;

        Bit_Gray_Code_2_Add_Noise=awgn(Bit_Gray_Code_2,Es_N0(i));
        Bit_Gray_Code_4_Add_Noise=awgn(Bit_Gray_Code_4,Es_N0(i));
        Bit_Gray_Code_6_Add_Noise=awgn(Bit_Gray_Code_6,Es_N0(i));

        %硬判决
        Dec_Hard_Code_2=zeros(1,2*Nmb_Symbol);
        Dec_Hard_Code_4=zeros(1,4*Nmb_Symbol);
        Dec_Hard_Code_6=zeros(1,6*Nmb_Symbol);

        Dec_Hard_Code_2(:,1:2:2*Nmb_Symbol)=sign(real(Bit_Gray_Code_2_Add_Noise));
        Dec_Hard_Code_2(:,2:2:2*Nmb_Symbol)=sign(imag(Bit_Gray_Code_2_Add_Noise));

        Dec_Hard_Code_4(:,1:4:4*Nmb_Symbol)=sign(real(Bit_Gray_Code_4_Add_Noise));
        Dec_Hard_Code_4(:,2:4:4*Nmb_Symbol)=sign(abs(real(Bit_Gray_Code_4_Add_Noise))-2*sqrt(1/10));
        Dec_Hard_Code_4(:,3:4:4*Nmb_Symbol)=sign(imag(Bit_Gray_Code_4_Add_Noise));
        Dec_Hard_Code_4(:,4:4:4*Nmb_Symbol)=sign(abs(imag(Bit_Gray_Code_4_Add_Noise))-2*sqrt(1/10));

        Dec_Hard_Code_6(:,1:6:6*Nmb_Symbol)=sign(real(Bit_Gray_Code_6_Add_Noise));
        Dec_Hard_Code_6(:,2:6:6*Nmb_Symbol)=sign(abs(real(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42));
        Dec_Hard_Code_6(:,3:6:6*Nmb_Symbol)=sign(abs(abs(real(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42))-2*sqrt(1/42));
        Dec_Hard_Code_6(:,4:6:6*Nmb_Symbol)=sign(imag(Bit_Gray_Code_6_Add_Noise));
        Dec_Hard_Code_6(:,5:6:6*Nmb_Symbol)=sign(abs(imag(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42));
        Dec_Hard_Code_6(:,6:6:6*Nmb_Symbol)=sign(abs(abs(imag(Bit_Gray_Code_6_Add_Noise))-4*sqrt(1/42))-2*sqrt(1/42));
        
        %软判决,先定义接收信号减去星座点，2^M次方的星座点
        for loop_llr_det=1:1:Nmb_Symbol
            for loop_xzd=1:1:2^2
                Euclidean_Distance_2(loop_llr_det,loop_xzd)=-(abs(Bit_Gray_Code_2_Add_Noise(loop_llr_det)-Constellation_Symbols_2(loop_xzd)))^2/n0(i);            
            end
            for loop_llr_bit=1:1:2
                ind_temp_plus_2 = find(Constellation_Bits_2(loop_llr_bit,:)==1);             
                ind_temp_minus_2=find(Constellation_Bits_2(loop_llr_bit,:)==-1);            
                a_temp_plus_2 = Euclidean_Distance_2(loop_llr_det,ind_temp_plus_2);
                a_temp_minus_2=Euclidean_Distance_2(loop_llr_det,ind_temp_minus_2);
                llr_det_2(1,(loop_llr_det-1)*2+loop_llr_bit)=log(sum(exp(a_temp_plus_2)))-log(sum(exp(a_temp_minus_2)));
               % llr_det_2(1,(loop_llr_det-1)*2+loop_llr_bit)=min(a_temp_minus_2)-min(a_temp_plus_2);
            end   

            for loop_xzd=1:1:2^4
                Euclidean_Distance_4(loop_llr_det,loop_xzd)=-(abs(Bit_Gray_Code_4_Add_Noise(loop_llr_det)-Constellation_Symbols_4(loop_xzd)))^2/n0(i);            
            end
            for loop_llr_bit=1:1:4
                ind_temp_plus_4 = find(Constellation_Bits_4(loop_llr_bit,:)==1);             
                ind_temp_minus_4=find(Constellation_Bits_4(loop_llr_bit,:)==-1);            
                a_temp_plus_4 = Euclidean_Distance_4(loop_llr_det,ind_temp_plus_4);
                a_temp_minus_4=Euclidean_Distance_4(loop_llr_det,ind_temp_minus_4);
                llr_det_4(1,(loop_llr_det-1)*4+loop_llr_bit)=max(a_temp_plus_4)+log(sum(exp(a_temp_plus_4-max(a_temp_plus_4))))-max(a_temp_minus_4)-log(sum(exp(a_temp_minus_4-max(a_temp_minus_4))));
            end 

            for loop_xzd=1:1:2^6
                Euclidean_Distance_6(loop_llr_det,loop_xzd)=-(abs(Bit_Gray_Code_6_Add_Noise(loop_llr_det)-Constellation_Symbols_6(loop_xzd)))^2/n0(i);            
            end
            for loop_llr_bit=1:1:6
                ind_temp_plus_6 = find(Constellation_Bits_6(loop_llr_bit,:)==1);            
                ind_temp_minus_6=find(Constellation_Bits_6(loop_llr_bit,:)==-1);            
                a_temp_plus_6 = Euclidean_Distance_6(loop_llr_det,ind_temp_plus_6);
                a_temp_minus_6=Euclidean_Distance_6(loop_llr_det,ind_temp_minus_6);
                llr_det_6(1,(loop_llr_det-1)*6+loop_llr_bit)=max(a_temp_plus_6)+log(sum(exp(a_temp_plus_6-max(a_temp_plus_6))))-max(a_temp_minus_6)-log(sum(exp(a_temp_minus_6-max(a_temp_minus_6))));
            end             
        end

        Dec_Soft_Code_2=sign(llr_det_2);
        Dec_Soft_Code_4=sign(llr_det_4);
        Dec_Soft_Code_6=sign(llr_det_6); 
        %记录每轮误码个数
%         err_num_round_QPSK_Hard=0;
%         err_num_round_16QAM_Hard=0;
%         err_num_round_64QAM_Hard=0;        
%         for indx = 1:length(Dec_Hard_Code_2)
%             if Dec_Hard_Code_2(indx)~=Bit_Code_2(indx)
%                 err_num_round_QPSK_Hard=err_num_round_QPSK_Hard+1;
%             end
%         end
%         for indx = 1:length(Dec_Hard_Code_4)
%             if Dec_Hard_Code_4(indx)~=Bit_Code_4(indx)
%                 err_num_round_16QAM_Hard=err_num_round_16QAM_Hard+1;
%             end
%         end
%         for indx = 1:length(Dec_Hard_Code_6)
%             if Dec_Hard_Code_6(indx)~=Bit_Code_6(indx)
%                 err_num_round_64QAM_Hard=err_num_round_64QAM_Hard+1;
%             end
%         end
%         err_num_QPSK_Hard=err_num_QPSK_Hard+err_num_round_QPSK_Hard;
%         err_num_16QAM_Hard=err_num_16QAM_Hard+err_num_round_16QAM_Hard;
%         err_num_64QAM_Hard=err_num_64QAM_Hard+err_num_round_64QAM_Hard;
        
        %记录每轮误码个数
        err_num_QPSK_Hard =err_num_QPSK_Hard +length(find(Dec_Hard_Code_2~=Bit_Code_2));
        err_num_16QAM_Hard=err_num_16QAM_Hard+length(find(Dec_Hard_Code_4~=Bit_Code_4));
        err_num_64QAM_Hard=err_num_64QAM_Hard+length(find(Dec_Hard_Code_6~=Bit_Code_6));

        err_num_QPSK_Soft =err_num_QPSK_Soft +length(find(Dec_Soft_Code_2~=Bit_Code_2));
        err_num_16QAM_Soft=err_num_16QAM_Soft+length(find(Dec_Soft_Code_4~=Bit_Code_4));
        err_num_64QAM_Soft=err_num_64QAM_Soft+length(find(Dec_Soft_Code_6~=Bit_Code_6));
    end
    BER_QPSK_Hard =[BER_QPSK_Hard err_num_QPSK_Hard/(Nmb_Symbol*2*Block)];
    BER_16QAM_Hard=[BER_16QAM_Hard err_num_16QAM_Hard/(Nmb_Symbol*4*Block)];
    BER_64QAM_Hard=[BER_64QAM_Hard err_num_64QAM_Hard/(Nmb_Symbol*6*Block)];   

    BER_QPSK_Soft =[BER_QPSK_Soft err_num_QPSK_Soft/(Nmb_Symbol*2*Block)];
    BER_16QAM_Soft=[BER_16QAM_Soft err_num_16QAM_Soft/(Nmb_Symbol*4*Block)];
    BER_64QAM_Soft=[BER_64QAM_Soft err_num_64QAM_Soft/(Nmb_Symbol*6*Block)];      
end
%绘制图像

semilogy(Es_N0,BER_QPSK_Hard,'-.X',Es_N0,BER_16QAM_Hard,Es_N0,BER_64QAM_Hard,'-.X',Es_N0,BER_QPSK_Soft,'-.o',Es_N0,BER_16QAM_Soft,'-.o',Es_N0,BER_64QAM_Soft,'-.o');
title('软判决和硬判决BER和Es_N0的关系')
xlabel('Es_N0/dB')
ylabel('BER')
grid on
% subplot(312)
% semilogy(Es_N0,BER_16QAM_Hard,'-.X',Es_N0,BER_16QAM_Soft,'-.o');
% title('16QAM调制下，软判决和硬判决BER和Es_N0的关系')
% xlabel('Es_N0/dB')
% ylabel('BER')
% grid on
% subplot(313)
% semilogy(Es_N0,BER_64QAM_Hard,'-.X',Es_N0,BER_64QAM_Soft,'-.o');
% title('64QAM调制下，软判决和硬判决BER和Es_N0的关系')
% xlabel('Es_N0/dB')
% ylabel('BER')
% grid on

文章参考：移动通信基础（7）软判决-CSDN博客