% Plot the truncated trigonometric Fourier series of the square wave
% on page 196 of the Lathi textbook.
%
% ELEC 320, Rich Kozick
%

N = input('Number of harmonics: ');

clear a
t = -3*pi:6/1000:3*pi;  % Time values to evaluate function: -3 to 3, steps of 6/1000
T = 2*pi;		  % Period
w0 = 2*pi/T;	  % Fundamental frequency

a0 = 0.5;
xN = a0*ones(1,length(t));   % dc component

% Compute truncated exponential Fourier series

for n=1:N
  if (rem(n,2) == 0)	% Even harmonic
    a(n) = 0;
  else			% Odd harmonic
    a(n) = 2 / (pi*n) * sin(n*pi/2);
  end
  xN = xN + a(n)*cos(n*w0*t);
end

% Define the square wave
g = zeros(size(t));
z1 = find(abs(t) < pi/2);
z2 = find((abs(t) > 3*pi/2) & (abs(t) < 5*pi/2));
g(z1) = ones(size(z1));
g(z2) = ones(size(z2));

figure(1)
subplot(211)
stem((0:N)*w0, [a0, a])
xlabel('FREQUENCY (rad/sec)')
ylabel('a_n')
title('Fourier coefficients, a_n')
subplot(212)
plot(t,xN, 'b-', t, g, 'r-')
title(['Square wave with harmonics <= ',num2str(N)])
xlabel('Time (sec)')
ylabel(['x',num2str(N),'(t)'])

