/* nag_tsa_spectrum_univar (g13cbc) Example Program.

 *

 * Copyright 1996 Numerical Algorithms Group.

 *

 * Mark 4, 1996.

 * Mark 8 revised, 2004.

 *

 */



#include <nag.h>

#include <nag_stdlib.h>

#include <stdio.h>

#include <nagg13.h>



#define KCMAX 400

#define NXMAX KCMAX/2



int main(void)

{



  Integer exit_status=0, i, kc, l, mw, ng, nx;

  NagError fail;

  double *g, pw, px, *stats=0, *x=0, *xh=0;



  INIT_FAIL(fail);



  Vprintf("nag_tsa_spectrum_univar (g13cbc) Example Program Results\n");

  /* Skip heading in data file */

  Vscanf("%*[^\n] ");



  Vscanf("%ld ",&nx);

  if (nx>=1 && nx <= NXMAX)

    {

      if ( !( stats = NAG_ALLOC(4, double)) ||

           !( x = NAG_ALLOC(KCMAX, double)) ||

           !( xh = NAG_ALLOC(NXMAX, double))

           )

        {

          Vprintf("Allocation failure\n");

          exit_status = -1;

          goto END;

        }

    }

  else

    {

      Vprintf("Invalid nx.\n");

      exit_status = 1;

      return exit_status;

    }

  for (i = 1; i <= nx; ++i)

    Vscanf("%lf ", &xh[i - 1]);

  px = 0.2;

  mw = nx;

  pw = 0.5;

  kc = KCMAX;

  l = 100;

  while ((scanf("%ld ", &mw)) != EOF)

    {

      if (mw > 0 && mw <= nx)

        {

          for (i = 1; i <= nx; ++i)

            x[i - 1] = xh[i - 1];



          /* nag_tsa_spectrum_univar (g13cbc).

           * Univariate time series, smoothed sample spectrum using

           * spectral smoothing by the trapezium frequency (Daniell)

           * window

           */

          nag_tsa_spectrum_univar(nx, Nag_Mean, px, mw, pw, l, kc, Nag_Logged,

                                  x, &g, &ng, stats, &fail);

          if (fail.code != NE_NOERROR)

            {

              Vprintf("Error from nag_tsa_spectrum_univar (g13cbc).\n%s\n",

                      fail.message);

              exit_status = 1;

              goto END;

            }



          if (mw == nx)

            Vprintf("\n No smoothing\n\n");

          else

            Vprintf("\n Frequency width of smoothing window = 1/%ld"

                    "\n\n", mw);

          Vprintf(" Degrees of freedom =%4.1f      Bandwidth =%7.4f\n\n",

                  stats[0],stats[3]);

          Vprintf(" 95 percent confidence limits -     Lower =%7.4f  "

                  "Upper =%7.4f\n\n", stats[1], stats[2]);

          Vprintf("       Spectrum       Spectrum       Spectrum"

                  "       Spectrum\n");

          Vprintf("       estimate       estimate       estimate"

                  "       estimate\n\n");

          for (i = 1; i <= ng; ++i)

            Vprintf("%5ld%10.4f%s",i,g[i - 1], (i%4==0 ? "\n": ""));

          Vprintf("\n");

          if (g)

            NAG_FREE(g);

        }

    }

 END:

  if (stats) NAG_FREE(stats);

  if (x) NAG_FREE(x);

  if (xh) NAG_FREE(xh);

  return exit_status;

}