%triangular wave filtering
clear all
hold off
j=sqrt(-1);
f0=500;     %fundamental freq of input triangular wave
T0 = 1/f0;  %period 
ts
tep = 0.005*T0;
no_sample = 3*T0/tstep + 1; %no. of samples  within  3*T0
no_sample1 = T0/tstep + 1; %no. of samples  within  T0
%tt = -0.5*T0:tstep:0.5*T0;
tt = -1.5*T0:tstep:1.5*T0;

tt1 = -0.5*T0:tstep:0.5*T0; % time vector for the period -0.5T0 to 0.5T0
gp1 = tt1/T0; %input - triangular wave in the period -0.5T0 to 0.5T0
gp_in = [gp1 gp1(2:no_sample1-1) gp1]; %3 cycles of the triangular wave  
Hp1 = plot(tt,gp_in)
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
title('input - time domain')
pause

      
      
K=1/(2*pi);
N=100; %no. of harmonics
nvec = -N:N;
for n = nvec
    m = n+N+1;
    absn = abs(n);
    if(n > 0)
      c_in(m) = j*K*((-1)^n)/n;
    elseif (n == 0)
      c_in(m) = 0.0;
    else 
      c_in(m) = -j*K*((-1)^absn)/absn;
    end
end
f = nvec*f0; %frequency vector
Hp1=plot(f,abs(c_in))
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
title('magnitude spectrum of input')
pause
Hp1=plot(f,angle(c_in))
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
axis([-0.1e4 0.1e4 -pi pi])
title('phase spectrum of input')
pause
R=3.3e3;
C=0.1e-6;
fc=1/(2*pi*R*C)     %cutoff freq of filter
%fc = 5000;
pause

Hf = 1 ./(1+j*f/fc) ;%filter transfer function
c_out = c_in .* Hf; %Fourier coefficients of the filter output
Hp1=plot(f,abs(c_out))
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
title('magnitude spectrum of filter output')
pause
hold on
Hp1=plot(f,abs(c_in),'r')
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
pause
hold off
Hp1=plot(f,angle(c_out))
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
title('phase spectrum of output')
axis([-0.1e4 0.1e4 -pi pi])
pause
hold on
Hp1=plot(f,angle(c_in),'r')
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
pause
hold off

%c_out = c_in;

A = zeros(2*N+1,ceil(no_sample));
for n = nvec
    m=n+N+1;
    A(m,:) = c_out(m) .* exp(j*2*pi*n*f0*tt);
end
gp_out = sum(A);
Hp1 = plot(tt,real(gp_out))
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)
title('filter output-time domain')
pause
hold on
Hp1 = plot(tt,gp_in,'r')
set(Hp1,'LineWidth',2)
Ha = gca;
set(Ha,'Fontsize',16)