clc;
clear?all;
K=2;?????????%認(rèn)知用戶數(shù)
N=64;????????%系統(tǒng)總的子載波數(shù)
cmax=6;??????%每個(gè)子載波平均所能攜帶的最大比特?cái)?shù)目
cmin=120;????%滿足自身QoS需求所需要的最小信道容量
B=1e6;???????%總的信道帶寬
Pe=1e-7;?????%誤碼率
P_av=1;?
SNR_av=0;?
Noise=P_av./10.^（SNR_av./10）;?%噪聲功率
N0=Noise./（B/N）;??????????????%噪聲功率譜密度
%?N0=1e-6;
%?PL=B/N*N0*（10.^（SNR_av./10））;
PL=0.0078;%保證認(rèn)知用戶能夠正常通信的功率需求
Pt=1;?????????????????????????%所有認(rèn)知用戶的總功率上限
PT=0.5;???????????????????????%各認(rèn)知用戶的總功率上限
l=1:cmax;?????????????????????%子載波每次分配增加的比特?cái)?shù)

%?子載波分配
%?h?=?rand（KN）;??%信道系數(shù)
%?h?=?random（‘rayleigh‘1KN）;?%信道增益
%?H（k）=?h（k:）.^2;??????????????%信道增益的平方
Nk_initial=zeros（K1）;??????????%各用戶初始分配的子載波數(shù)
h=zeros（KN）;???????????????????%信道增益
for?k=1:K
????Nk_initial（k）=floor（cmin/cmax）;%每個(gè)用戶最少占用的子載波數(shù)目
????if?sum（Nk_initial）>N???????????%如果所有用戶需要的最少的子載波數(shù)目之和仍大于所有的子載波數(shù)目，則報(bào)錯(cuò)
????????error（‘There?are?two?many?users!\n‘）;
????end
????for?n=1:N
????????h（kn）=1/sqrt（2）*（randn（1）?+?1j*randn（1））;??%信道增益
????end
end
H=abs（h）;
%?sort（H2‘descend‘）;%對(duì)各用戶的信道增益進(jìn)行降序排列
Htemp?=?H.^2;?????????%信道系數(shù)的臨時(shí)記錄值，若該子載波被占用，則信道系數(shù)設(shè)置為0

P=zeros（K1）;????%各用戶的功率需求
omiga=cell（K1）;?%前20個(gè)子載波組成的集合，初始化為空集
f=zeros（KN）;????%用于存儲(chǔ)接收到c個(gè)比特時(shí)所需的功率
p=zeros（KN）;????%用于存儲(chǔ)認(rèn)知用戶k在第n個(gè)子載波上消耗的功率???
Htemp1=Htemp;
for?k=1:K
????m=1;
????while?m<=Nk_initial（k）
????????Hmax=max（max（Htemp1（k:）））;%找出最大的信道系數(shù)
????????[kh?index]=find（Htemp1（k:）==Hmax）;%找出該信道系數(shù)的坐標(biāo)
????????%此時(shí)并計(jì)算認(rèn)知用戶k的第n個(gè)子載波的功率
????????f（km）=N0/3*（qfuncinv（Pe/4）*qfuncinv（Pe/4））*（2^（cmin./Nk_initial（k））-1）;%MQAM
????????p（km）=f（km）/Htemp1（kindex）;
????????omiga{k}=union（omiga{k}index）;
????????Htemp1（kindex）=0;
????????m=m+1;
????end
????P（k）=sum（p（k:）2）;
end
%?P=sum（p2）;
Ptotal=sum（P）;??????????%總功率需求
Ntotal=sum（Nk_initial）;?%分配的子載波數(shù)總和

Pp=zeros（K1）;
PP=zeros（K1）;
ff=zeros（KN）;???
pp=zeros（KN）;
CP=zeros（K1）;%記錄各用戶每多分配一個(gè)子載波時(shí)所減小的功率
NK=Nk_initial;
Htemp2=Htemp;
while?Ptotal>Pt*（Ntotal/N）
%?????Htemp2=Htemp;
????for?k=1:K
????????m=1;
????????while?m<=NK（k）
????????????Hmax2=max（max（Htemp2（k:）））;%找出最大的信道系數(shù)
????????????[kh1?index1]=find（Htemp2（k:）==Hmax2）;%找出該信道系數(shù)的坐標(biāo)
%?????????????[Hmax?nm]=max（Htemp2（k:））;%找出最大的信道系數(shù)??
????????????%此時(shí)并計(jì)算認(rèn)知用戶k的第n個(gè)子載波的功率
????????????f（km）=N0/3*（qfuncinv（Pe/4）*qfuncinv（Pe/4））*（2^（cmin./NK（k））-1）;%MQAM
%?????????????p（km）=f（km）/Hmax;
????????????p（km）=f（km）/Htemp2（kindex1）;
????????????ff（km）=N0/3*（qfuncinv（Pe/4）*qfuncinv（Pe/4））*（2^（cmin./（NK（k）+1））-1）;%MQAM
%?????????????pp（km）=ff（km）/Hmax;
????????????pp（km）=ff（km）/Htemp2（kindex1）;
????????????omiga{k}=union（omiga{k}index1）;
????????????Htemp2（kindex1）=0;
%?????????????omiga{k}=union（omiga{k}nm）;
%?????????????Htemp2（knm）=0;
????????????m=m+1;
????????end
????????Pp（k）=sum（p（k:）2）;
????????PP（k）=sum（pp（k:）2）;
????end
????????Htemp4=Htemp2;
????????for?kk=1:K
????????????Hmax3=max（max（Htemp4（kk:）））;%找出最大的信道系數(shù)
????????????[kh2?index2]=find（Htemp4（kk:）==Hmax3）;%找出該信道系數(shù)的坐標(biāo)
????????????%此時(shí)并計(jì)算認(rèn)知用戶k的在原來基礎(chǔ)上增加一個(gè)子載波的功率
????????