% common emitter amplifier, inductor applied to remove dc, influence of the
% inductor current ac ripple
% all voltages normalized by Vcc
% all currents normalized by Icq
% inductor normalized to l=(w L Icq/Vcc)
% optimal load resistor, R=Vcc/Icq

clear all
close all

%initialization
npoint=720; % reduce if the simulation is too slow
i=1:npoint;
wt=2*pi*(i-0.5)/npoint;
deg=wt*180/pi;
ps=sin(wt);
pc=cos(wt);

np=40; % number of points over l, reduce if the simulation is too slow
ndec=4; % number of decades to be covered
l0=10^(-ndec/2);
lstep=10^(ndec/np);

for k=1:np
    l0=l0*lstep;
    l(k)=l0; % actual normalized inductance
    p1=l0/sqrt(1+l0^2);
    Vm(k)=p1; % amplitude of the output voltage
    p2=1/sqrt(1+l0^2);
    mce=1-p1*ps;
    jl=1-p2*pc;
    jout=p1*ps;
    jc=jl+jout;
    mout=1-mce;
    pcc=1*jc; % 1* is for educational purposes, a symbol for Vcc*
    pd=mce.*jc;
    pout=mout.*jout;
    pl=jl.*mout;
    Pcc(k)=mean(pcc);
    Pd(k)=mean(pd);
    Pout(k)=mean(pout);
    Pl(k)=mean(pl);
    efficiency(k)=Pout(k)/Pcc(k)*100;
    figure(1)
    subplot(4,2,1)
    plot(deg,mce)
    xlabel('wt [deg]')
    ylabel('v_{CE}/V_{CC}')
    axis([0 360 0 2])
    set(gca,'xtick',[0 30 60 90 120 150 180 210 240 270 300 330 360])
    subplot(4,2,3)
    plot(deg,jc)
    xlabel('wt [deg]')
    ylabel('i_C/I_{CQ}')
    axis([0 360 0 2])
    set(gca,'xtick',[0 30 60 90 120 150 180 210 240 270 300 330 360])
    subplot(4,2,5)
    plot(deg,mout)
    xlabel('wt [deg]')
    ylabel('v_{OUT}/V_{CC}')
    axis([0 360 -1 1])
    set(gca,'xtick',[0 30 60 90 120 150 180 210 240 270 300 330 360])
    subplot(4,2,2)
    plot(deg,pcc)
    xlabel('wt [deg]')
    ylabel('p_{CC}/(V_{CC} I_{CQ})')
    axis([0 360 0 2])
    set(gca,'xtick',[0 30 60 90 120 150 180 210 240 270 300 330 360])
    subplot(4,2,4)
    plot(deg,pd)
    xlabel('wt [deg]')
    ylabel('p_D/(V_{CC} I_{CQ})')
    axis([0 360 0 2])
    set(gca,'xtick',[0 30 60 90 120 150 180 210 240 270 300 330 360])
    subplot(4,2,6)
    plot(deg,pout)
    xlabel('wt [deg]')
    ylabel('p_{OUT}/(V_{CC} I_{CQ})')
    axis([0 360 0 2])
    set(gca,'xtick',[0 30 60 90 120 150 180 210 240 270 300 330 360])
    subplot(4,2,7)
    plot(mce,jc,'.')
    xlabel('v_{CE}/V_{CC}')
    ylabel('i_C/I_{CQ}')
    axis([0 2 0 2])
    subplot(4,2,8)
    plot(deg,pl)
    xlabel('wt [deg]')
    ylabel('p_L/(V_{CC} I_{CQ})')
    axis([0 360 -1 1])
    set(gca,'xtick',[0 30 60 90 120 150 180 210 240 270 300 330 360])
    pause(1)
end

figure(2)
hold on
plot(l,Pcc,'g')
plot(l,Pd,'b')
plot(l,Pout,'r')
plot(l,Pl,'y')
xlabel('w*L*Icq/Vcc')
ylabel('P_{CC} [green], P_D [blue], P_{OUT} [red], P_L [yellow], all /(V_{CC} I_{CQ})')
axis([min(l) max(l) -0.1 1.1])

figure(3)
plot(l,efficiency)
xlabel('w L I_{CQ}/V_{CC}')
ylabel('efficiency [%]')
axis([min(l) max(l) -25 75])

figure(4)
loglog(l,Vm)
xlabel('w L I_{CQ}/V_{CC}')
ylabel('V_m/V_{CC}')