// g05pm Example Program Text
// C# version, NAG Copyright 2008
using System;
using NagLibrary;
using System.Globalization;
using System.IO;
namespace NagDotNetExamples
{
  public class G05PME
  {
    static string datafile = "ExampleData/g05pme.d";
    const double one=1.00e0;
    const double two=2.00e0;
    static void Main(String[] args)
    {
      if (args.Length == 1)
      {
        datafile = args[0];
      }
      StartExample();
    }
    public static void StartExample()
    {
      try
      {
        DataReader sr = new DataReader(datafile);
        double ad, alpha, dv, tmp, var, z;
        int en, genid, i, itype, ival, j, k=0, mode, n, nf, nsim, p, smode, subid;
        int[] seed = new int[1];
        int ifail;
        Console.WriteLine("g05pm Example Program Results");
        //      Initialise the seed to a repeatable sequence
        seed[0] = 1762543;
        //      genid and subid identify the base generator.
        genid = 1;
        subid = 1;
        //      Initialise the generator to a repeatable sequence
        G05.G05State g05State = new G05.G05State(genid, subid, seed, out ifail);
        if (ifail != 0)
        {
          Console.WriteLine("   ** Generator initialisation returned with ifail = {0,5}",ifail);
          goto L140;
        }
        //      Skip headings in data file
        sr.Reset();
        //      Read in the initial arguments and check array sizes
        sr.Reset();
        mode = int.Parse(sr.Next());
        itype = int.Parse(sr.Next());
        n = int.Parse(sr.Next());
        nf = int.Parse(sr.Next());
        nsim = int.Parse(sr.Next());
        alpha = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
        double[,] blim = new double[2, nf];
        double[,] bsim = new double[nsim, nf];
        double[] e = new double[1];
        double[] fse = new double[nf];
        double[] fv = new double[nf];
        double[,] glim = new double[2, nf];
        double[,] gsim = new double[nsim, nf];
        double[] param = new double[4];
        double[] q = new double[2];
        double[] res = new double[n];
        double[] tsim1 = new double[n];
        double[] tsim2 = new double[n];
        double[] y = new double[n];
        double[] yhat = new double[n];
        if (n < 1  || nsim < 1)
        {
          Console.WriteLine("   ** Problem size is too small.");
          goto L140;
        }
        sr.Reset();
        for (i = 1 ; i <= n ; i++)
        {
          y[i - 1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
        }
        //      Read in the itype dependent arguments (skipping headings)
        if (itype == 1)
        {
          sr.Reset();
          param[0] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          p = 0;
          ival = 1;
        }
        else if (itype == 2)
        {
          sr.Reset();
          param[0] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          param[1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          p = 0;
          ival = 2;
        }
        else if (itype == 3)
        {
          sr.Reset();
          param[0] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          param[1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          param[2] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          p = 0;
          ival = 2;
        }
        else
        {
          sr.Reset();
          param[0] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          param[1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          param[2] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          param[3] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          p = int.Parse(sr.Next());
          ival = p + 2;
        }
        double[] init = new double[ival];
        double[] r = new double[p+13];
        //      Read in the mode dependent arguments (skipping headings)
        if (mode == 0)
        {
          //         User supplied initial values
          sr.Reset();
          for (i = 1 ; i <= ival ; i++)
          {
            init[i - 1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          }
        }
        else if (mode == 1)
        {
          //         Continuing from a previously saved R
          sr.Reset();
          for (i = 1 ; i <= p + 13 ; i++)
          {
            r[i - 1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
          }
        }
        else if (mode == 2)
        {
          //         Initial values calculated from first K observations
          sr.Reset();
          k = int.Parse(sr.Next());
        }
        // 
        //      Fit a smoothing model (parameter r in G05.g05PMF and state in G13.g13amf
        //      are in the same format)
        G13.g13am(mode, itype, p, param, n, y, k, init, nf, fv, fse, yhat, res, out dv,out ad, r,
        out ifail);
        if (ifail != 0)
        {
          Console.WriteLine("   ** g13am returned with ifail = {0,5}",ifail);
          goto L140;
        }
        // 
        //      Simulate forecast values from the model, and don't update R
        smode = 2;
        var = dv * dv;
        en = (int)n;
        //      Simulate nsim forecasts
        for (i = 1 ; i <= nsim ; i++)
        {
          //         Simulations assuming gaussian errors
          G05.g05pm(smode, nf, itype, p, param, init, var, r, g05State, e, 0, tsim1, out ifail);
          if (ifail != 0)
          {
            Console.WriteLine("   ** G05.g05pm returned with ifail = {0,5}",ifail);
            goto L140;
          }
          //         Bootstrapping errors
          G05.g05pm(smode, nf, itype, p, param, init, 0.00e0, r, g05State, res, en, tsim2,
          out ifail);
          if (ifail != 0)
          {
            Console.WriteLine("   ** G05.g05pm returned with ifail = {0,5}",ifail);
            goto L140;
          }
          //         Copy and transpose the simulated values
          for (j = 1 ; j <= nf ; j++)
          {
            gsim[i - 1 , j - 1] = tsim1[j - 1];
            bsim[i - 1 , j - 1] = tsim2[j - 1];
          }
        }
        // 
        //      Calculate ci based on the quantiles for each simulated forecast
        q[0] = alpha / two;
        q[1] = one - q[0];
        double [] gsim_vec = new double[nsim];
        double [] glim_vec = new double[nsim];
        double [] bsim_vec = new double[nsim];
        double [] blim_vec = new double[nsim];
        for (i = 1 ; i <= nf ; i++)
        {
          for (j = 1; j<= nsim; j++)
          {
            gsim_vec[j-1] = gsim[j-1, i-1];
            bsim_vec[j-1] = bsim[j-1, i-1];
          }
          G01.g01am(nsim, gsim_vec, 2, q, glim_vec, out ifail);
          if (ifail != 0)
          {
            Console.WriteLine("   ** g01am returned with ifail = {0,5}",ifail);
            goto L140;
          }
          G01.g01am(nsim, bsim_vec, 2, q, blim_vec, out ifail);
          if (ifail != 0)
          {
            Console.WriteLine("   ** g01am returned with ifail = {0,5}",ifail);
            goto L140;
          }
          for (j = 1; j<= 2; j++)
          {
            glim[j-1, i-1] = glim_vec[j-1];
            blim[j-1, i-1] = blim_vec[j-1];
          }
        }
        // 
        //      Display the forecast values and associated prediction intervals
        Console.WriteLine("");
        Console.WriteLine(" {0}","Initial values used:");
        for (i = 1 ; i <= ival ; i++)
        {
          Console.WriteLine("{0,4}  {1,12:f3}  ",i,init[i - 1]);
        }
        Console.WriteLine("");
        Console.WriteLine("{0}{1,12:e4}","Mean Deviation     = ",dv);
        Console.WriteLine("{0}{1,12:e4}","Absolute Deviation = ",ad);
        Console.WriteLine("");
        Console.WriteLine(" {0}","         Observed      1-Step");
        Console.WriteLine(" {0}"," Period   Values      Forecast      Residual");
        Console.WriteLine("");
        for (i = 1 ; i <= n ; i++)
        {
          Console.WriteLine("{0,4}  {1,12:f3}  {2,12:f3}  {3,12:f3}",i,y[i - 1],yhat[i - 1],res[i - 1]);
        }
        Console.WriteLine("");
        Console.WriteLine(" {0}",("") + ("      Simulated CI         Simulated CI"));
        Console.WriteLine(" {0}",("Obs.  Forecast      Estimated CI   ") + ("     (Gaussian Errors)    (Bootstrap Errors)"));
        z = G01.g01fa("L", q[1], out ifail);
        for (i = 1 ; i <= nf ; i++)
        {
          tmp = z * fse[i - 1];
          Console.WriteLine("{0,3}  {1,10:f3}  {2,10:f3}  {3,10:f3}  {4,10:f3}  {5,10:f3}  {6,10:f3}  {7,10:f3}",n + i,fv[i - 1],fv[i - 1] - tmp,fv[i - 1] + tmp,glim[0 , i - 1],glim[1 , i - 1],blim[0 , i - 1],blim[1 , i - 1]);
        }
        Console.WriteLine("  {0,5:f1}{1}",100.00e0 * (one - alpha),"% CIs were produced");
        // 
        L140: ;
        // 
      }
      catch (Exception e)
      {
        Console.WriteLine(e.Message);
        Console.WriteLine( "Exception Raised");
      }
    }
  }
}