% ECE2311D10 
% Phasor analysis example
% DRB 5-Apr-2010
% =================================================
% USER PARAMETERS
% =================================================
N = 100;                        % number of Fourier coefficients
w0 = pi/3;                      % fundamental frequency of input
t = -10:0.001:10;               % time vector (could do this with linspace too)
a0 = 1/2;                       % DC term
L = 1;                          % inductor
R = 1;                          % resistor
% =================================================

j = sqrt(-1);                   % unit imaginary number
x = zeros(1,length(t));         % allocate space for input array
y = zeros(1,length(t));         % allocate space for output array

for n=0:N,
    
    omega = n*w0;
    H = 1/(1+j*omega*L/R);
    
    if n==0,
        x = x+a0;
        y = y+abs(H)*a0;
    else
        % input coefficients
        an = 6/(pi^2*n^2)*(cos(n*w0)-cos(2*n*w0));
        bn = 0;
        Cn = sqrt(an^2+bn^2);
        if an>=0,
            thetan = 0;
        else
            thetan = pi;
        end
        x = x+Cn*cos(n*w0*t+thetan);
        
        % output coefficients
        Cnprime = abs(H)*Cn;
        thetanprime = thetan+angle(H);
        y = y+Cnprime*cos(n*w0*t+thetanprime);        
    end
    
end

% plot
plot(t,x,t,y);
xlabel('time');
ylabel('input and output signals');
legend('x(t)','y(t)');
grid on