/* nag_tsa_uni_arma_roots (g13dxc) Example Program.
 *
 * Copyright 2023 Numerical Algorithms Group.
 *
 * Mark 29.3, 2023.
 */

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

int main(void) {
  /* Scalars */
  Integer exit_status, i, ip, k, npar;
  NagError fail;

  /* Arrays */
  double *par = 0, *ri = 0, *rmod = 0, *rr = 0;

  INIT_FAIL(fail);

  exit_status = 0;

  printf("nag_tsa_uni_arma_roots (g13dxc) Example Program Results\n");

  /* Skip heading in data file */
  scanf("%*[^\n] ");
  scanf("%" NAG_IFMT "%" NAG_IFMT "%*[^\n] ", &k, &ip);

  if (k > 0 && ip > 0) {
    /* Allocate arrays */
    if (!(par = NAG_ALLOC(k * k * ip, double)) ||
        !(ri = NAG_ALLOC(k * ip, double)) ||
        !(rmod = NAG_ALLOC(k * ip, double)) ||
        !(rr = NAG_ALLOC(k * ip, double))) {
      printf("Allocation failure\n");
      exit_status = -1;
      goto END;
    }

    /* Read the AR (or MA) parameters */
    npar = ip * k * k;
    for (i = 1; i <= npar; ++i)
      scanf("%lf", &par[i - 1]);
    scanf("%*[^\n] ");

    /* nag_tsa_uni_arma_roots (g13dxc).
     * Calculates the zeros of a vector autoregressive (or
     * moving average) operator
     */
    nag_tsa_uni_arma_roots(k, ip, par, rr, ri, rmod, &fail);
    if (fail.code != NE_NOERROR) {
      printf("Error from nag_tsa_uni_arma_roots (g13dxc).\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
    printf("\n");
    printf("        Eigenvalues       Moduli\n");
    printf("        -----------       ------\n");
    for (i = 1; i <= k * ip; ++i) {
      if (ri[i - 1] >= 0.0)
        printf("%10.3f  + %6.3f i  %8.3f\n", rr[i - 1], ri[i - 1], rmod[i - 1]);
      else
        printf("%10.3f  - %6.3f i  %8.3f\n", rr[i - 1], -ri[i - 1],
               rmod[i - 1]);
    }
  } else
    printf(" Either k or ip is out of range\n");

END:
  NAG_FREE(par);
  NAG_FREE(ri);
  NAG_FREE(rmod);
  NAG_FREE(rr);

  return exit_status;
}