
%% R-C Low Pass Filter
R = 30;
C = 50e-12;

%% Circuit Frequency Response
j = 1i;
w = logspace(5, 12, 1e4);
H = -j./w/C./(R-j./w/C);

%% Plot (linear) Frequency Response
figure(1);
subplot(2,1,1)
plot(w/2/pi/1e3,abs(H), 'LineWidth',3)
xlabel('frequency [kHz]');
ylabel('|H|');

subplot(2,1,2)
plot(w/2/pi/1e3,angle(H)*180/pi, 'LineWidth',3)
xlabel('frequency [kHz]');
ylabel('\angle{H} [\circ]')

%% Plot (log-scale) Frequency Response
figure(2);
subplot(2,1,1)
Bode = loglog(w/2/pi,abs(H), 'LineWidth',3);
xlabel('frequency [Hz] (Log Scale)');
ylabel('|H| (Log Scale)');

subplot(2,1,2)
semilogx(w/2/pi,angle(H)*180/pi, 'LineWidth',3);
xlabel('frequency [Hz] (Log Scale)');
ylabel('\angle{H} [\circ]')

%% Digital Signal Fourier Series
Vcomm = 1.2;        % Supply Voltage
f = 4e9;            % clock frequency
wo = 2*pi*f;

Nh = 1000;          % number of harmonics to consider

B0 = Vcomm/2;       % Average (DC) value

n = 1:2:Nh;         % only odd harmonics of sin signals are non-zero
Bn = Vcomm*2./n/pi; % Fourier series coefficients (as solved in class)

t = linspace(0, 5/f, 1e4);  % time vector (5 periods of signal)

%% recreate the signal from the first Nh harmonics
% fhat is our approximation to the squarewave signal
fhat = B0*ones(size(t));    % start with the average value
for i = 1:length(n)         % ...then add each harmonic one-by-one
    fhat = fhat + Bn(i)*sin(2*pi*f*n(i)*t);
end

%% Plot signal approximation trough first Nh harmonics
figure(3);
plot(t*1e9,B0+B0*sign(sin(wo*t)), ':k', 'Linewidth', 2);
hold on;
plot(t*1e9,fhat, 'Linewidth', 2);
xlabel('time [ns]');
ylabel('V_{in} [V]');
legend('f(t)', ['\Sigma_{k=1}^{' num2str(Nh) '}b_ksin(k\omega_0 t)']);

%% Find output signal by multiplication of fhat and frequency response
% yhat is the output signal
yhat = B0*ones(size(t));        %start with dc solution
for i = 1:length(n)             %for each frequency
    w = 2*pi*f*n(i);
    H(i) = -j./w/C./(R-j./w/C);    %find H(jw)
    
    %multiply magnitudes and add phases
    yhat = yhat + Bn(i)*sin(w*t + angle(H(i)))*abs(H(i));
end

%% Plot the predicted output from first Nh terms of Fourier Series
figure(4);
plot(t*1e9,yhat, 'Linewidth', 2);
xlabel('time [ns]');
ylabel('V_{out} [V]');

%% Plot Fourier Series of input and output along with frequency response on the same axes
figure(2);
subplot(2,1,1);
hold on;
for i = 1:15
   Fhat(i) = plot([f*n(i) f*n(i)], [1e-9, Bn(i)], 'r', 'Linewidth', 2);
   Yhat(i) = plot([f*n(i) f*n(i)], [1e-9, abs(Bn(i)*H(i))], 'k', 'Linewidth', 2);
end

subplot(2,1,1)
xlim([1e6, 1e11])
ylim([1e-9, 1e1])
legend([Bode, Fhat(1), Yhat(1)],'|H(jw)|', 'F[V_{in}]', 'F[V_{in}]*|H(j\omega)|', 'Location', 'southwest');

subplot(2,1,2)
xlim([1e6, 1e11])







