Fischer比特功率分配算法及MATLAB实现

原理介绍

Fischer 算法的优化准则是在确定的系统总速率和发射功率的约束条件下，使系统的误比特率达到最小。Fischer 算法给出了比特分配的闭式解，它首先把各个子信道上的噪声功率值${\log }_2{N}_i$ 存储下来，接下来就只需进行一种加法和除数为整数的除法，因此它的复杂度较 Chow算法有了进一步的降低。Fischer 算法是目前效率比较高的算法之一。Fischer的表达式无论是对于实际应用还是理论分析，都具有指导作用

比特分配

（1）首先必须已知各个子信道上的噪声方差$N_i\forall i=1,\cdots ,N$ 也可以认为是信道增益平方分倒数。设置目标比特速率，也就是要分配的比特总数$R_T$ 。设$N^{\prime }$ 为已使用的子信道数，设 $N^{\prime }$的初始值$N^{\prime }\text{=}N$ 。记激活的子信道集合为$I$ ，设 $I$的初始值 I = { 1 , 2 , ⋯   , N } I\text{=}\left\{ 1,2,\cdots ,N \right\} I={1,2,⋯,N}
（2）计算各个子信道的$LD{N}_i={\log }_2{N}_{i}i=1,\cdots ,N$并把 这些值存储起来
（3）计算$I$ 中各个子信道可分配的比特数目
$$b(i)=\left(R_T+\sum _{l\in I}LD{N}_i\right)/R^{\prime }-LD{N}_i$$
若 $b(i)\le 0$且$i\in I$ ，那么$N^{\prime }\text{=}{N}^{\prime }-1$ ，把第 $i$个子信道从 $I$中剔除。然后转到步骤（2），继续下去，直到$b(i)>0$ 且 $i\in I$
（4） $\hat{b}(i)=round(b(i))，diff(i)=b(i)-\hat{b}(i)$
（5）计算 $R\text{=}sum\left(\left\{\hat{b}\left(i\right)\right\}\right)$
（6）若$R\text{=}{R}_T$ ，则 进行功率分配
（7）若$R>R_T$ ，找到最小的$diff\left(i\right)$ ，相应地 $\hat{b}(i)\text{=}\hat{b}(i)\text{-}1$， $R=R-1$ ，$diff\left(i\right)=diff\left(i\right)+1$ ，重复此步骤直至$R\text{=}{R}_T$ ，转到（6）
（8）若 $R<R_T$，找到最大的$diff\left(i\right)$ ，相应地$\hat{b}(i)\text{=}\hat{b}(i)\text{+}1$ ，$R=R\text{+}1$ ，$diff\left(i\right)=diff\left(i\right)\text{-}1$ ，重复此步骤直至$R\text{=}{R}_T$ ，转到（6）

功率分配

进行功率分配。每个激活的子信道上分配的发射功率按下式算。
$$P_i=\frac{P_T\cdot N_i\cdot 2^{b(i)}}{\sum _{l\in I}{N}_i\cdot 2^{b(i)}}$$

仿真参数设置

参数名称	参数值
平均信噪比	10
子载波数目	32
平均发射功率	1
比特总数上限	128

某次仿真结果如下：

从图中可以看出，经过Fischer算法分配之后，由于该算法使系统的误比特率达到最小，那么在信道状况比较差的信道就会分配较多的功率，从而满足误比特率最小，结果也与该算法的要求相一致。
代码如下：主函数

clear;
close all;
clc;
SNR_av=10;   %平均信噪比
Num_subc=32;   %载波数
P_av=1;    %平均功率
Pt=P_av*Num_subc;   %发射功率
Rt=128;       %每个OFDM符号所包含的比特数
noise=P_av./10.^(SNR_av./10);   %噪声功率
gain_subc=random('rayleigh',1,1,Num_subc);  %子载波的瑞利衰落
[bit_alloc,power_alloc]=Fischer(Num_subc,Rt,Pt,gain_subc); %Fisher资源分配
power_alloc=Pt.*(power_alloc./sum(power_alloc));    %功率分配
figure(1); 
subplot(2,1,1); 
plot(gain_subc,'-r'); 
legend('信道增益'); 
hold on; 
stem(bit_alloc); 
title('Fischer比特分配与功率分配'); 
ylabel('比特分配结果'); 
xlabel('子载波序号'); 
subplot(2,1,2); 
stem(power_alloc); 
ylabel('功率分配结果'); 
xlabel('子载波序号'); 

Fischer.m

function [bit_alloc,power_alloc]=Fischer(Num_subc,Rt,Pt,gain_subc) 
%------------------bit allocation ----------- 
N_H=1./gain_subc.^2; %信道增益平方的倒数作为已知子载波的噪声方差
LD=log2(N_H);            %计算并存储，不用下次计算
index_use0=1:Num_subc;     %激活的子载波数
num_use=Num_subc;          %使用的子载波数目
bit_temp0=(Rt+sum(LD))./Num_subc-LD; %计算各个子载波可分配的比特数目
bit_temp=zeros(1,Num_subc);  %暂存比特数目
flag=1;  
while flag 
    id_remove=find(bit_temp0(index_use0)<=0); %找到没有用的载波数
    index_remove=index_use0(id_remove); %将其移出激活的子载波集合
    if(length(index_remove)==0)   %所有都激活，则直接跳出
        break; 
    end     
    index_use0=setdiff(index_use0,index_remove); %统计已经激活的子载波，即将上述需要移出的子载波进行移出
    num_use=length(index_use0);   %可用的子载波数目
    bit_temp0(index_use0)=(Rt+sum(LD(index_use0)))./num_use-LD(index_use0); %计算各个子载波可分配的比特数目
    flag=1; 
end 
index_use=index_use0;   %可用的子载波
bit_temp(index_use)=bit_temp0(index_use);  
%-----------------bit round--------------------------- 
bit_round=zeros(1,Num_subc);   
bit_round(index_use)=round(bit_temp(index_use)); %四舍五入取整
bit_diff(index_use)=bit_temp(index_use)-bit_round(index_use); %两者差值
bit_total=sum(bit_round(index_use));  %计算包含的比特数
%------------------bit alteration-------------------- 
while(bit_total>Rt)                            %如果 bit_total>Rt，找到最小的bit_diff(i)
    id_alter=find(bit_round(index_use)>0);     %相应地bit_round减1，bit_total=bit_total-1
    index_alter=index_use(id_alter);           %bit_diff(i)加1，直至bit_total=Rt
    min_diff=min(bit_diff(index_alter));  
    id_min=find(bit_diff(index_alter)==min_diff,1);
    index_min=index_alter(id_min);
    bit_round(index_min)=bit_round(index_min)-1; 
    bit_total=bit_total-1; 
    bit_diff(index_min)=bit_diff(index_min)+1; 
end 
while(bit_total<Rt)                                 %如果 bit_total<Rt，找到最大的bit_diff(i)
    id_alter=find(bit_round(index_use)>0);          %相应地bit_round加1，bit_total=bit_total+1
    index_alter=index_use(id_alter);                %bit_diff(i)减1，直至bit_total=Rt
    max_diff=max(bit_diff(index_alter)); 
    id_max=find(bit_diff(index_alter)==max_diff,1);
    index_max=index_alter(id_max);
    bit_round(index_max)=bit_round(index_max)+1; 
    bit_total=bit_total+1; 
    bit_diff(index_max)=bit_diff(index_max)-1; 
end 
%----------------bit allocation------------------ 
bit_alloc=bit_round; 
%------------------power allocation----------------------- 
power_alloc=zeros(1,Num_subc);    
%---------------power allocation----------------------- 
index_use2=find(bit_alloc>0); 
power_alloc(index_use2)=Pt.*(N_H(index_use2)).*2.^bit_round(index_use2)./... 
    (sum(N_H(index_use2)).*2.^bit_round(index_use2));   
end