% The script file pwoprs03 calculates the Fourier series
% coefficients for damped free vibrations looks at the influence of the
% number of points. The values of coefficients change but the amplitudes
% are similar.

g=9.81; T1=0.1; f=1/T1; omega1=2*pi*f; amax=5; dt=T1/100; 
eta=1/4/T1;

t=0:dt:4*T1; t=t'; N=4*T1/dt;
AM=exp(-eta*t); AME=amax*AM;
ZH=amax*exp(-i*omega1*t); XH=real(ZH); 
ZD=AM.*ZH; RZD=real(ZD);
plot(t,RZD,'k',t,AME,'r',t,-AME,'r',t,XH,'b')
title('Displacements in damped vibrations envelops and periodic component')
pause
close

a0= mean(RZD); RZD=RZD-a0*ones(size(RZD));
FS=fft(RZD);
   Ak=2*real(FS(2:21))/N; Bk=-2*imag(FS(2:21))/N;
   Nu=1:20; Nu=Nu';
   plot(Nu,Ak,'r*',Nu,Bk,'k. ')
   title('Coefficients of Fourier series x=0')
   pause
   close
   
   AS=abs(Ak+i*Bk);
   plot(Nu,AS,'k*')
   title('Amplitudes of Fourier Coefficients of Fourier series x=0')
   pause
   close

for s=1:10;
   RZDS=RZD(1:N-s);
   a0= mean(RZDS); RZDS=RZDS-a0*ones(size(RZDS));
   FS=fft(RZDS);
   Ak=2*real(FS(2:21))/N; Bk=-2*imag(FS(2:21))/N;
   Nu=1:20; Nu=Nu';
   plot(Nu,Ak,'r*',Nu,Bk,'k. ')
   title('Coefficients of Fourier series s=1:10')
   pause
   close
   AS=abs(Ak+i*Bk);
   plot(Nu,AS,'k*')
   title('Amplitudes of Fourier Coefficients s=1:10')
   pause
   close
end

