/* nag_stat_summary_onevar_combine (g01auc) Example Program.
 *
 * Copyright 2022 Numerical Algorithms Group.
 *
 * Mark 28.3, 2022.
 */

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

#define MRCOMM(I, J) mrcomm[(J)*20 + (I)]

int main(void) {
  /* Integer scalar and array declarations */
  Integer b, i, iwt, j, nb, pn;
  Integer exit_status = 0;

  /* NAG structures and types */
  NagError fail;

  /* Double scalar and array declarations */
  double xkurt, xmax, xmean, xmin, xsd, xskew;
  double rcomm[20];
  double *mrcomm = 0, *wt = 0, *x = 0;

  /* Initialize the error structure */
  INIT_FAIL(fail);

  printf("nag_stat_summary_onevar_combine (g01auc) Example Program Results"
         "\n\n");

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

  /* Read in the number of block of data we have */
  scanf("%" NAG_IFMT "%*[^\n] ", &b);

  if (!(mrcomm = NAG_ALLOC(20 * b, double))) {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }

  /* Loop over each block of data */
  for (j = 0; j < b; j++) {
    /* Read in the number of observations in this block and a flag indicating
       whether weights have been supplied (iwt = 1) or not (iwt = 0) */
    scanf("%" NAG_IFMT "%" NAG_IFMT "%*[^\n] ", &nb, &iwt);

    /* Reallocate X to the required size */
    NAG_FREE(x);
    if (!(x = NAG_ALLOC(nb, double))) {
      printf("Allocation failure\n");
      exit_status = -2;
      goto END;
    }

    /* Read in the data for this block */
    if (iwt) {
      /* Weights supplied, so reallocate WT to the required size */
      NAG_FREE(wt);
      if (!(wt = NAG_ALLOC(nb, double))) {
        printf("Allocation failure\n");
        exit_status = -3;
        goto END;
      }
      for (i = 0; i < nb; i++)
        scanf("%lf%lf", &x[i], &wt[i]);
    } else {
      /* No weights */
      NAG_FREE(wt);
      wt = 0;

      for (i = 0; i < nb; i++)
        scanf("%lf", &x[i]);
    }
    scanf("%*[^\n] ");

    /* Call nag_stat_summary_onevar (g01atc) to summarise
       this block of data */
    pn = 0;
    nag_stat_summary_onevar(nb, x, wt, &pn, &xmean, &xsd, &xskew, &xkurt, &xmin,
                            &xmax, rcomm, &fail);
    if (fail.code != NE_NOERROR && fail.code != NE_CASES_ONE &&
        fail.code != NE_ZERO_VARIANCE && fail.code != NE_CASES_ZERO) {
      printf("Error from nag_stat_summary_onevar (g01atc).\n%s\n",
             fail.message);
      exit_status = 1;
      goto END;
    }

    /* Save RCOMM for future reference */
    for (i = 0; i < 20; i++)
      MRCOMM(i, j) = rcomm[i];

    /* Display the results for this block */
    printf(" Summary for block %" NAG_IFMT "\n", j + 1);
    if (fail.code == NE_CASES_ZERO)
      printf(" No valid observations supplied. All weights are zero.\n");
    else {
      printf(" %" NAG_IFMT " valid observations\n", pn);
      printf("   Mean           %13.2f\n", xmean);
      if (fail.code == NE_CASES_ONE) {
        printf("   Unable to calculate the standard deviation,");
        printf("skewness or kurtosis\n");
      } else {
        printf("   Std devn       %13.2f\n", xsd);
        if (fail.code == NE_ZERO_VARIANCE)
          printf("   Unable to calculate the skewness and kurtosis\n");
        else {
          printf("   Skewness       %13.2f\n", xskew);
          printf("   Kurtosis       %13.2f\n", xkurt);
        }
      }
      printf("   Minimum        %13.2f\n", xmin);
      printf("   Maximum        %13.2f\n\n", xmax);
    }
  }

  /* Call nag_stat_summary_onevar_combine (g01auc) to combine the
     summaries across all the blocks */
  nag_stat_summary_onevar_combine(b, mrcomm, &pn, &xmean, &xsd, &xskew, &xkurt,
                                  &xmin, &xmax, rcomm, &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_stat_summary_onevar_combine (g01auc).\n%s\n",
           fail.message);
    exit_status = 2;
    goto END;
  }

  /* Display the combined results */
  printf(" Summary for the combined data\n");
  if (fail.code == NE_CASES_ZERO)
    printf(" No valid observations supplied. All weights are zero.\n");
  else {
    printf(" %" NAG_IFMT " valid observations\n", pn);
    printf("   Mean           %13.2f\n", xmean);
    if (fail.code == NE_CASES_ONE) {
      printf("   Unable to calculate the standard deviation,");
      printf(" skewness or kurtosis\n");
    } else {
      printf("   Std devn       %13.2f\n", xsd);
      if (fail.code == NE_ZERO_VARIANCE)
        printf("   Unable to calculate the skewness and kurtosis\n");
      else {
        printf("   Skewness       %13.2f\n", xskew);
        printf("   Kurtosis       %13.2f\n", xkurt);
      }
    }
    printf("   Minimum        %13.2f\n", xmin);
    printf("   Maximum        %13.2f\n", xmax);
  }

END:
  NAG_FREE(x);
  NAG_FREE(wt);
  NAG_FREE(mrcomm);

  return (exit_status);
}