%% Headding
% Written by: William Brooks
% Date Created: 2016-4-6
% Date Last Modified: 2016-4-7
%
% MATLAB Script to Plot Transfer Curves Given H(jw)
% User Defined: Frequency Range and Plotting vs. Omega in [Hz] or [rad/sec]
% Plots Response in Both Linear and Logarithmic (Bode) Formats

%% Housekeeping
close all; clear; clc;

%% Plotting Parameters (USER INPUT REQUIRED)
startFreq = 10.^(-2); % Starting Frequency [rad/s]
endFreq = 10.^(8); % Ending Frequency [rad/s]
numPt = 5; % Log_{10}{Number of Points}
freqScale = 1; % Control the Reporting of Frequency
    % (1) [Hz]
    % (else) [rad/sec]
circ = 1; % Circuit Selection

%% Plotting Vectors (NO USER INPUT REQUIRED)
w = ones(2,10.^numPt + 1); % Pre-Allocate Omega Array of Positive Frequencies
    % Row 1 is a Linearly Spaced Vector
    % Row 2 is a Logarithmically Spaced Vector
w(1,:) = linspace(startFreq, endFreq, 10.^numPt + 1);
w(2,:) = logspace(log10(startFreq), log10(endFreq), 10.^numPt + 1);

%% Function Definition (USER INPUT REQUIRED)
% (Voltage) Transfer Function Modelling the Circuit
% MATLAB expects 1i to be used as the complex constant (j)

switch (circ)

case 1 % Circuit 1
strFunct = 'H(j\omega)=[1]/[1 - w^2lc + jwrc]';
% LaTeX Formatted Title
l1 = 0.47*10.^(-3);
c1 = 0.1.*10.^(-9);
r1 = 1.0*10.^(3);
transFunct = (1)./(1 - w.^2.*l1.*c1 + 1i.*w.*r1.*c1);

case 2 % Circuit 2
strFunct = 'H(j\omega)=[-j\omega C_1R_2]/[(1+j\omega C_1R_1)(1 + j\omega C_2R_2)]';
% LaTeX Formatted Title
r1 = 1.0*10.^(3); r2 = 100.0*10.^(3);
c1 = 0.1.*10^(-6); c2 = c1;
transFunct = (-1.*1i.*w.*c1.*r2)./((1 + 1i.*w.*c1.*r1).*(1 + 1i.*w.*c2.*r2));

end

%% Calculation Routines (NO USER INPUT REQUIRED)
rad2deg = 180./pi; % Conversion Factor from Radians to Degrees
transMag = zeros(size(w)); % Pre-Allocate Magnitude Array
transMag(1,:) = abs(transFunct(1,:)); % Linear Magnitude
transMag(2,:) = 20.*log10(abs(transFunct(2,:))); % Magnitude in dB [V/V]
transAngle = unwrap(angle(transFunct)).*rad2deg; % Phase Angle of H(jw)
if(freqScale == 1) % Custom Setting of X-Axis Units
    rad2hz = 1./(2.*pi);
    w = w.*rad2hz; % Convert w to units of [Hz]
else
    % w = w; % Plot Omega in Rad/Sec
end

%% Plotting Routines (NO USER INPUT REQUIRED)
figure('name','Linear Plot of Transfer Function');

% Linear Plot
subplot(2,1,1) % Magnitude Response
    ptLinMag = plot(w(1,:),transMag(1,:),'linewidth',1.5);
    grid on % Enable Grid
    title(['Linear Plot of ', strFunct]); % Title Sub-Plots With Function
    if(freqScale == 1) % Ensure Appropriate Labelling of X-Axis
        xlabel('\omega [Hz]'); % X-Axis Label
    else
        xlabel('\omega [rad/sec]'); % X-Axis Label
    end
    ylabel('Magnitude [V/V]'); % Y-Axis Label
    legend('|H(j\omega)|','Location','Best'); % Legend

subplot(2,1,2) % Phase Response
    ptLinPhase = plot(w(1,:),transAngle(1,:),'linewidth',1.5);
    grid on % Enable Grid
    if(freqScale == 1) % Ensure Appropriate Labelling of X-Axis
        xlabel('\omega [Hz]'); % X-Axis Label
    else
        xlabel('\omega [rad/sec]'); % X-Axis Label
    end
    ylabel('Phase Angle [degrees]'); % Y-Axis Label
    legend('Angle(H(j\omega))','Location','Best'); % Legend

% Bode Plot
figure('name','Logarithmic Plot of Transfer Function (Bode Plot)');

subplot(2,1,1) % Magnitude Response
    ptLinMag = semilogx(w(2,:),transMag(2,:),'linewidth',1.5);
    grid on; % Enable Grid
    title(['Bode Plot of ', strFunct]); % Title Sub-Plots With Function
    if(freqScale == 1) % Ensure Appropriate Labelling of X-Axis
        xlabel('\omega [Hz]'); % X-Axis Label
    else
        xlabel('\omega [rad/sec]'); % X-Axis Label
    end
    ylabel('Magnitude [dB]'); % Y-Axis Label
    legend('|H(j\omega)|','Location','Best'); % Legend

subplot(2,1,2) % Phase Response
    ptLinPhase = semilogx(w(2,:),transAngle(2,:),'linewidth',1.5);
    grid on; % Enable Grid
    if(freqScale == 1) % Ensure Appropriate Labelling of X-Axis
        xlabel('\omega [Hz]'); % X-Axis Label
    else
        xlabel('\omega [rad/sec]'); % X-Axis Label
    end
    ylabel('Phase Angle [degrees]'); % Y-Axis Label
    legend('Angle(H(j\omega))','Location','Best'); % Legend
