% ----------------------------------------------------------------------
% *********************INITIAL SETUP********************* 

clear all;
format long;

% ----------------------------------------------------------------------
% *********************PLL PARAMETERS********************* 
f = linspace (10000, 1e9, 20000);
w = 2*pi*f;
% Input Clock Frequency
FIN = 5e6;
fs = 1/10e-12;

% Charge Pump parametere
ICP = 100e-6;
err_1 = 0;
err_2 = 0;
SR = 10e9; %Slew Rate is 1V in 100ps

% Loop Filter parameters
Rz = 15e3;
Cz = 42e-12;
Cp = 5.3e-12;


% KVCO parameters. 
Kv = 1e9; %In Hz/V for simulink
KVCO = 1e9*2*pi; %In rad/V
KvcoOffset = 0;
Fo = 2.3e9;


% Divider Parameters
N = 480;

% ----------------------------------------------------------------------
% *********************START SIMULATION********************* 

options=simset('MaxStep',1/fs,...
    'RelTol',3e-3/fs,'AbsTol',1e-4/fs, ...
    'Solver','ode45',...
    'ZeroCross','on' );

sim('PLL_ab',[0 5e-6], options)

% ----------------------------------------------------------------------
% *********************POST PROCESSING********************* 

        Kpd = ICP/(2*pi);
        Kfwd = Kpd*KVCO;
        Kloop = Kfwd/N ;
        a1 = Rz*Cz;
        b21 = Rz*Cz*Cp;
        b22 = (Cp+Cz);
       
        Aopenloop_num = Kloop*[a1 1];
        Aopenloop_den = [b21 b22 0 0];
        Aopenloop = freqs(Aopenloop_num,Aopenloop_den,w);
        Gain_Open_Loop = db(abs(Aopenloop));    
        Phase_Open_Loop = (180/pi)*angle(Aopenloop);
       
        Zero = abs(roots(Aopenloop_num)/(2*pi));
        Pole = abs(roots(Aopenloop_den)/(2*pi));
        BW = Kpd*KVCO*Rz/N*(Cz/(Cz+Cp))/(2*pi);
        DF = Rz/2*(Kpd*KVCO*Cz/N)^(1/2);
               
        [y x1] = min(abs(Gain_Open_Loop));
        Unity_Gain_Frequency = f(x1);
        Phase_Margin = Phase_Open_Loop(x1) + 180;
        Natural_freq = Unity_Gain_Frequency/(2*DF);      
         figure(1)
         subplot(2,1,1)
         semilogx(f,Gain_Open_Loop,'Color','b','LineWidth',2);
         title('The Gain Plot of the Open Loop Transfer Function');
         xlabel('frequency (Hz)');
         ylabel('Gain (dB)');
         grid on;
         
         subplot(2,1,2)
         semilogx(f,Phase_Open_Loop,'Color','b','LineWidth',2);
         title('The Phase Plot of the Open Loop Transfer Function');
         xlabel('frequency (Hz)');
         ylabel('Phase (degrees)');
         grid on;
         
 % Closed Loop Analysis of the VCO (Low Pass Transfer Function)
      
        Alowpass_num = Kfwd*N*[a1 1];
        Alowpass_den = [N*b21 N*b22 a1*Kfwd Kfwd];
        Alowpass = freqs(Alowpass_num,Alowpass_den,w);
        Gain_Low_Pass = db(abs(Alowpass));        
           
% Closed Loop Analysis of the VCO (High Pass Transfer Function from the VCO
% input to the output)

        Ahighpass_num = (N)*[b21 b22 0 0];
        Ahighpass_den = [N*b21 N*b22 a1*Kfwd Kfwd];
        Ahighpass = freqs(Ahighpass_num,Ahighpass_den,w);
        Gain_High_Pass = db(abs(Ahighpass));   
        
% Closed Loop Analysis of the VCO (Band Pass Transfer Function from the VCO
% output to the output)

        Abandpass_num =  N*KVCO*[b21 b22 0];
        Abandpass_den = [N*b21 N*b22 a1*Kfwd Kfwd];
        Abandpass = freqs(Abandpass_num,Abandpass_den,w);
        Gain_Band_Pass = db(abs(Abandpass));        
        
         figure
         subplot(3,1,1)
         semilogx(f,Gain_Low_Pass,'Color','g','LineWidth',2 );
         Title('Low Pass Transfer Function (for Reference Noise)');
         xlabel('frequency (Hz)');
         ylabel('Gain (dB)');
         grid on;
%         
         subplot(3,1,2)
         semilogx(f,Gain_High_Pass,'Color','r','LineWidth',2);
         Title('High Pass Transfer Function (for CLK Buffer Noise)');
         xlabel('frequency (Hz)');
         ylabel('Gain (dB)');
         grid on;
         
         subplot(3,1,3) 
         semilogx(f,Gain_Band_Pass,'Color','b','LineWidth',2);
         Title('Band Pass Transfer Function (for VCO input Noise)');
         xlabel('frequency (Hz)');
         ylabel('Gain (dB)');
         grid on;
         
%figure;

%semilogx(f,Gain_Low_Pass,'Color','g','LineWidth',2);
%hold on; grid on;
%semilogx(f,Gain_High_Pass,'Color','r','LineWidth',2);
%hold on; grid on;
%semilogx(f,Gain_Band_Pass,'Color','b','LineWidth',2);
%legend ('Gain Low Pass','Gain High Pass','Gain Band Pass',4);
%Title('Noise Transfer function, N=480' );
%         xlabel('frequency (Hz)');
%         ylabel('Noise Transfer function (dB)');
%         grid on;


% For step response
%wz = 2*pi*Zero ;
%wp = 2*pi*max(Pole) ;
%wn = 2*pi*Natural_freq ;
%wd = wn*(1-DF^2)^0.5;

%figure;
%s = tf('s') ;
%HO = wn^2 * (1 + s/wz) / (s^2 * (1 + s/wp)) ;
%bode(HO) ;
%HC = wn^2 * (1 + s/wz) / (s^2 + ((wn^2)/wz)*s + wn^2) ;
%ltiview(HC) ;

%Jitter due to input ref Noise in an Ideal Second-Order PLL


        figure;
        semilogx(f,Gain_Low_Pass,'Color','b','LineWidth',2 );
         Title('Low Pass Transfer Function (for Reference Noise)');
         xlabel('frequency (Hz)');
         ylabel('Gain (dB)');
         hold on;
         grid on;

         wn_l = 2.3e6;
         df_l = 0.5;
         wd_l = wn_l*(1-df_l^2)^0.5;
         Num_l = N*[2*df_l*wn_l wn_l^2];
         den_l = [1 2*df_l*wn_l wn_l^2];
         NTF_LP=freqs(Num_l,den_l,w);
         Gain_LP=db(abs(NTF_LP));
         semilogx(f,Gain_LP,'Color','r','LineWidth',2);
         legend ('3rd Order PLL','2rd Order PLL',3); 

fout = 2.4e9;
wout = 2*pi*fout;
cycles = 10000;
T = 1/fout;
d_T = cycles*T;
Nclk_in = 1e-6;
k_l = (4*pi^2*Nclk_in*N/(wout)^2)^0.5;

% Jitter due by 1/f^2 noise
RMS_jitter_ref2 = (k_l^2*d_T*(1+1/(2*df_l*wn_l*d_T)+exp(-df_l*wn_l*d_T)/d_T*(sin(wd_l*d_T+pi/6)/(2*(1-df_l^2)*wn_l)-cos(wd_l*d_T)/(2*(1-df_l^2)*df_l*wn_l)-2*sin(wd_l*d_T)/wd_l)))^0.5

% Jitter due by white noise and T => 00
RMS_jitter_ref=(N*2/(wout)^2*Nclk_in*(wn_l*(4*df_l^2+1)/(4*df_l)))^0.5



%Jitter due to VCO Input Noise in an Ideal Second-Order PLL

figure;
semilogx(f,Gain_Band_Pass,'Color','b','LineWidth',2);
         Title('Band Pass Transfer Function (for VCO input Noise)');
         xlabel('frequency (Hz)');
         ylabel('Gain (dB)');
         hold on;
         grid on;

         
         wn_s = 2.3e6;
         df_s = 0.5;
         wd_s = wn_s*(1-df_s^2)^0.5;
         Num = [KVCO 0];
         den = [1 2*df_s*wn_s wn_s^2];
         second_band=freqs(Num,den,w);
         Gain_second_band=db(abs(second_band));
         semilogx(f,Gain_second_band, 'Color','r','LineWidth',2);
         legend ('3rd Order PLL','2rd Order PLL',3); 


NVCO=1e-3;
k = (4*pi^2*NVCO/(wout)^2)^0.5;
Jitter_VCO=(k^2*(1/(2*df_s*wn_s)+(exp(-df_s*wn_s*d_T)/(2*(1-df_s^2)))*(sin(wd_s*d_T+pi/2)/wn_s-cos(wd_s*d_T)/(df_s*wn_s))))^(1/2)

%Jitter due to VCOoutput(CLK Buff) Noise in an Ideal Second-Order PLL

figure; 
semilogx(f,Gain_High_Pass,'Color','b','LineWidth',2);
         Title('High Pass Transfer Function (for CLK Buffer Noise)');
         xlabel('frequency (Hz)');
         ylabel('Gain (dB)');
         hold on;
         grid on;
 
         wn_b = 2.8e6;
         df_b = 0.44;
         wd_b = wn_b*(1-df_b^2)^0.5;
         Num_b = [1 0 0];
         den_b = [1 2*df_b*wn_b wn_b^2];
         second_high=freqs(Num_b,den_b,w);
         Gain_second_high=db(abs(second_high));
         semilogx(f,Gain_second_high, 'Color','r','LineWidth',2);
         legend ('3rd Order PLL','2rd Order PLL',4); 
wbuf = 0;
Nbuf = 1e-6;
Jitter_Buff=(Nbuf/wout^2*(wbuf+wn_b*(1-12*df_b^2)/(2*df_b)-exp(-wbuf*d_T)*(wbuf-4*df_b*wn_b)-exp(-df_b*wn_b*d_T)*(-wn_b*sin(wd_b*d_T+3*pi/6.9-pi)/(2*(1-df_b^2))+wn_b*cos(wd_b*d_T)/(2*(1-df_b^2)*df_b)-2*wn_b*sin(wd_b*d_T+2*pi/6.9)/((1-df_b^2)^0.5))))^0.5







