% Script for calculating the derivatives of [0, 1] periodic sawtooth 
% and swell functions, given in the m-files saw.m and swell.m

  N = 16;

  x = (0:(N-1))/N;
  M = [0:(N/2-1) (-N/2):(-1)];
  k = 2*pi*M;

  Nhi = 8*N; % Define hi-res grid for plotting
  xhi = (0:(Nhi-1))/Nhi;

%  Plot exact and DFT derivative for sawtooth

  [y,yhat] = saw(x,M);
  Y = fft(y);
  dydx = real(ifft(1i*k.*Y));
  dydx_ex = 2*sign(xhi - 0.5);
  subplot(2,2,1)
  plot(xhi,dydx_ex,'-',x,dydx,'x')
  xlabel('x')
  ylabel('dy/dx')
  title(['Sawtooth, N = ' num2str(N)])
  legend('Exact','DFT',0)

%  Plot exact and DFT derivative for swell.

  [y,yhat] = swell(x,M);
  Y = fft(y);
  dydx = real(ifft(1i*k.*Y));
  dydx_ex = -1.2*pi*sin(2*pi*xhi)./(1 - 0.6*cos(2*pi*xhi)).^2;
  subplot(2,2,2)
  plot(xhi,dydx_ex,'-',x,dydx,'x')
  xlabel('x')
  ylabel('dy/dx')
  title(['1/(1 - 0.6 cos(2\pix)), N = ' num2str(N)])
  legend('Exact','DFT',0)