/* nag_partial_corr (g02byc) Example Program.
 *
 * Copyright 2000 Numerical Algorithms Group.
 *
 * Mark 6, 2000.
 */

#include <stdio.h>
#include <nag.h>
#include <nag_stdlib.h>
#include <nagg02.h>

int main(void)
{

  double *r=0, *std=0, *v, *x=0, *xbar=0, sw;
  Integer *sz=0, j, k, m, n, nx, ny;
  Integer exit_status=0;
  NagError fail;

#define X(I,J) x[((I)-1)*m + ((J)-1)]
#define R(I,J) r[((I)-1)*m + ((J)-1)]

  INIT_FAIL(fail);
  Vprintf("g02byc Example Program Results\n");

  /* Skip heading in data file */
  Vscanf("%*[^\n]");

  Vscanf("%ld %ld", &n, &m);
  if (!(r=NAG_ALLOC(m*m, double))
      || !(std=NAG_ALLOC(m, double))
      || !(v=NAG_ALLOC(m*m, double))
      || !(x=NAG_ALLOC(n*m, double))
      || !(xbar=NAG_ALLOC(m, double))
      || !(sz=NAG_ALLOC(m, Integer)))
    {
      Vprintf("Allocation failure\n");
      exit_status = -1;
      goto END;
    }
  for (j = 1; j <= n; ++j)
    for (k = 1; k <= m; ++k)
      Vscanf("%lf", &X(j,k));

  /* Calculate correlation matrix */
  g02bxc(n, m, x, m, 0, 0, &sw, xbar, std, r, m, v, m,
         &fail);
  if (fail.code == NE_NOERROR)
    {
      /* Print the correlation matrix */
      Vprintf("\nCorrelation Matrix\n\n");
      for (j=1; j<=m; j++)
        {
          for(k=1; k<=m; k++)
            if (j>k)
              Vprintf("%11s", "");
            else
              Vprintf("%7.4f%4s",  R(j,k),"");
          Vprintf("\n");
        }

      Vscanf("%ld %ld", &ny, &nx);
      for (j = 1; j <= m; ++j)
        Vscanf("%ld", &sz[j - 1]);

      /* Calculate partial correlation matrix */
      g02byc(m, ny, nx, sz, v, m, r, m, &fail);
      if (fail.code != NE_NOERROR)
        {
          Vprintf("Error from g02byc.\n%s\n", fail.message);
          exit_status = 1;
          goto END;
        }
      /* Print partial correlation matrix */
      Vprintf("\n");
      Vprintf("\nPartial Correlation Matrix\n\n");
      for (j=1; j<=ny; j++)
        {
          for(k=1; k<=ny; k++)
            {
              if (j>k)
                Vprintf("%11s", "");
              else if (j==k)
                Vprintf("%7.4f%4s",  1.0, "");
              else
                Vprintf("%7.4f%4s",  R(j,k), "");
            }
          Vprintf("\n");
        }

    }
  else
    {
      Vprintf("Error from g02bxc.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
 END:
  if (r) NAG_FREE(r);
  if (std) NAG_FREE(std);
  if (v) NAG_FREE(v);
  if (x) NAG_FREE(x);
  if (xbar) NAG_FREE(xbar);
  if (sz) NAG_FREE(sz);
  return exit_status;
}