//---------------------------------------------------------------------------
//  fftplot.r

//  Syntax:	fftplot ( Y )
//		fftplot ( Y , T )

//  Fftplot compute and plot the magnitude in dB and phase in degrees
//  of the FFT data in Y
//
//  T = sampling period, if this parameter is not given it is assumed
//      to be equal to one.

//  The magnitude and phase are returned in list elements `mag' and
//  `ph' respectively. The values for the frequency scale are returned
//  in list element `freq'. The magnitude and phase are truncated
//  above the Nyquist frequency. If plotf exists, then a plot of the
//  magnitude of Y (in dB) versus frequency is produced.
//
//---------------------------------------------------------------------------

fftplot = function ( Y, T, plotf )
{

  lab = ["DIGITAL FREQUENCY / pi" , "FREQUENCY HERTZ"];
  ly = log (Y[:]);
  mag = 20 * real(ly) / log(10);
  ph = imag (ly);
  n = max (size (Y))/2;
  f = ((0:n-1)/n)';
  if (!exist (T)) 
  { 
    labid = 1;
  else
    f = 0.5*f/T;
    labid = 2;
  }

  if (exist (plotf))
  {
    plot ("set title 'MAGNITUDE'");
    plot ("set xlabel '"+lab[labid]+"'");
    plot ("set ylabel 'dB'");

    plot ( [f, mag[1:n] ] );
  }

  return << freq = f; mag = mag[1:n]; ph = ph[1:n] >>
};