clear all
global log10pe;
global Pt;
global gamma;
global Ptoptfinal2 Rch Enode Rdata Rdec thermalNoise boltzmann T W R maxdeg

Enode = 10e-12;
Rdata = 1.5e9;
Rdec = Rdata;
boltzmann = 1.38e-23;
T = 300;
W = 3e9;
thermalNoise = boltzmann*T*W;
pathlossexp = 3;
r = 10;
lambda = (1/200);
maxdeg =4;

Rch = Rdata/(2*W);
R =Rch;
pathloss = (r/lambda)^pathlossexp;

log10pevec =    [-1:-0.1:-5 -5.25:-0.25:-10 -11:-1:-32 -34:-2:-64];
gamma = Enode*Rdec/ (Rch*thermalNoise*pathloss);


for i = 1:length(log10pevec)
    

    
    log10pe = log10pevec(i)
    
    [PtoptNormalized(i), PtotalminNormalized(i)] = findptopt();
    pe = 10^log10pe;
    ShannonLimit(i) = qfuncinv( hbinv(1-Rch*(1-pe))  )* pathloss * thermalNoise;
    Ptopt(i) = PtoptNormalized(i) * pathloss * thermalNoise;
%    Ptopt(i) = max(Ptopt(i),ShannonLimit(i));
    Ptotalmin(i) = PtotalminNormalized(i) * pathloss * thermalNoise;
%    Ptotalmin(i) = max(Ptotalmin(i),ShannonLimit(i));

    Optimal_transmit_power = Ptopt(i)
    Total_power = Ptotalmin(i)
    
end
figure;
semilogy(Ptotalmin,log10pevec);
hold on;
semilogy(Ptopt,log10pevec,'r');
semilogy(ShannonLimit,log10pevec,'k');
beep
