%%load the relevant circuit parameters
clear all
close all

set(0,'DefaultLineLineWidth',3)
set(0,'DefaultAxesFontSize',20)
set(0,'DefaultTextFontSize',20)

HRDESIGN

m = zeros(N-1,N-1,Nvals);
Current_Source = [1,zeros(1,N-2)]';
%'

mm = sum(y,1);

for ii = 1:N-1

      m(ii,ii,:) = mm(1,ii,:);

      for jj = 1:N-1
            if  ii ~= jj 
                 m(ii,jj,:) = -y(jj,ii,:);
            end
      end
end


Volts = zeros(N-1,Nvals);
for ii = 1:Nvals
       Volts(:,ii) = inv(m(:,:,ii))*Current_Source;
end

%Plot out the voltage at node 1. Since this is the node where 1 Ampere is injected, we are in effect calculating the input impedance...

%Convert to Thevenin equivalent...factor of 2. as we have 6 dB into a 50 Ohm load
factor = 2./50.;

figure(1)
plot(frequency,20*log10(abs(factor.*Volts(6,:))),'r')
axis([min(frequency) max(frequency)  -60 10]);
xlabel('Freq (Hz)')
  title(cat(2,'transfer function, C_c= ',num2str(Cc)))
%axis([40e6 60e6  -60 10]);
figure(2)
Zin = Volts(1,:);
Zin = Zin*Zo./(Zo-Zin);
plot(frequency, real(Zin), 'k')
hold on
plot(frequency, imag(Zin), 'r')
hold off

xlabel('Freq (Hz)')
  title(cat(2,'Impedance, C_c= ',num2str(Cc)))
