/* nag_wav_dim2_sngl_fwd (c09eac) Example Program.
 *
 * Copyright 2024 Numerical Algorithms Group.
 *
 * Mark 30.2, 2024.
 */

#include <nag.h>
#include <string.h>

int main(void) {
  /* Scalars */
  Integer exit_status = 0;
  Integer i, j, m, n, nf, nwcm, nwcn, nwct, nwl, pda, pdb, pdc;
  /* Arrays */
  char mode[24], wavnam[20], title[50];
  double *a = 0, *b = 0, *ca = 0, *cd = 0, *ch = 0, *cv = 0;
  Integer icomm[(180)];
  /* NAG types */
  Nag_Wavelet wavnamenum;
  Nag_WaveletMode modenum;
  Nag_MatrixType matrix = Nag_GeneralMatrix;
  Nag_OrderType order = Nag_ColMajor;
  Nag_DiagType diag = Nag_NonUnitDiag;
  NagError fail;

  INIT_FAIL(fail);

  printf("nag_wav_dim2_sngl_fwd (c09eac) Example Program Results\n\n");

  /* Skip heading in data file and read problem parameters */
  scanf("%*[^\n] %" NAG_IFMT "%" NAG_IFMT "%*[^\n]", &m, &n);
  pda = m;
  pdb = m;
  scanf("%19s%23s%*[^\n]\n", wavnam, mode);
  if (!(a = NAG_ALLOC((pda) * (n), double)) ||
      !(b = NAG_ALLOC((pdb) * (n), double))) {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }
  printf(" Parameters read from file :: \n");
  printf(" DWT :: Wavelet : %s\n", wavnam);
  printf("        End mode: %s\n", mode);
  fflush(stdout);

  /*
   * nag_enum_name_to_value (x04nac).
   * Converts NAG enum member name to value
   */
  wavnamenum = (Nag_Wavelet)nag_enum_name_to_value(wavnam);
  modenum = (Nag_WaveletMode)nag_enum_name_to_value(mode);

  /* Read data array */
#define A(I, J) a[(J - 1) * pda + I - 1]
  for (i = 1; i <= m; i++)
    for (j = 1; j <= n; j++)
      scanf("%lf", &A(i, j));
  scanf("%*[^\n] ");

  printf("\n");
  fflush(stdout);
  nag_file_print_matrix_real_gen_comp(
      order, matrix, diag, m, n, a, pda, "%8.4f",
      "Input Data     A :", Nag_NoLabels, 0, Nag_NoLabels, 0, 80, 0, 0, &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_file_print_matrix_real_gen_comp (x04cbc).\n%s\n",
           fail.message);
    exit_status = 1;
    goto END;
  }
  printf("\n");

  /* nag_wav_dim2_init (c09abc).
   * Two-dimensional wavelet filter initialization
   */
  nag_wav_dim2_init(wavnamenum, Nag_SingleLevel, modenum, m, n, &nwl, &nf,
                    &nwct, &nwcn, icomm, &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_wav_dim2_init (c09abc).\n%s\n", fail.message);
    exit_status = 2;
    goto END;
  }
  nwcm = nwct / (4 * nwcn);
  if (!(ca = NAG_ALLOC((nwcm) * (nwcn), double)) ||
      !(cd = NAG_ALLOC((nwcm) * (nwcn), double)) ||
      !(cv = NAG_ALLOC((nwcm) * (nwcn), double)) ||
      !(ch = NAG_ALLOC((nwcm) * (nwcn), double))) {
    printf("Allocation failure\n");
    exit_status = -1;
    goto END;
  }

  pdc = nwcm;
  /* nag_wav_dim2_sngl_fwd (c09eac).
   * Two-dimensional discrete wavelet transform
   */
  nag_wav_dim2_sngl_fwd(m, n, a, pda, ca, pdc, ch, pdc, cv, pdc, cd, pdc, icomm,
                        &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_wav_dim2_sngl_fwd (c09eac).\n%s\n", fail.message);
    exit_status = 1;
    goto END;
  }
  fflush(stdout);

  /* Print decomposition */
  strcpy(title, "Approximation coefficients   CA :");
  nag_file_print_matrix_real_gen_comp(order, matrix, diag, nwcm, nwcn, ca, pdc,
                                      "%8.4f", title, Nag_NoLabels, 0,
                                      Nag_NoLabels, 0, 80, 0, 0, &fail);
  printf("\n");
  fflush(stdout);
  strcpy(title, "Diagonal coefficients        CD :");
  nag_file_print_matrix_real_gen_comp(order, matrix, diag, nwcm, nwcn, cd, pdc,
                                      "%8.4f", title, Nag_NoLabels, 0,
                                      Nag_NoLabels, 0, 80, 0, 0, &fail);
  printf("\n");
  fflush(stdout);
  strcpy(title, "Horizontal coefficients      CH :");
  nag_file_print_matrix_real_gen_comp(order, matrix, diag, nwcm, nwcn, ch, pdc,
                                      "%8.4f", title, Nag_NoLabels, 0,
                                      Nag_NoLabels, 0, 80, 0, 0, &fail);
  printf("\n");
  fflush(stdout);
  strcpy(title, "Vertical coefficients        CV :");
  nag_file_print_matrix_real_gen_comp(order, matrix, diag, nwcm, nwcn, cv, pdc,
                                      "%8.4f", title, Nag_NoLabels, 0,
                                      Nag_NoLabels, 0, 80, 0, 0, &fail);
  printf("\n");

  /* nag_wav_dim2_sngl_inv (c09ebc).
   * Two-dimensional inverse discrete wavelet transform
   */
  nag_wav_dim2_sngl_inv(m, n, ca, pdc, ch, pdc, cv, pdc, cd, pdc, b, pdb, icomm,
                        &fail);
  if (fail.code != NE_NOERROR) {
    printf("Error from nag_wav_dim2_sngl_inv (c09ebc).\n%s\n", fail.message);
    exit_status = 3;
    goto END;
  }
  fflush(stdout);

  /* Print reconstruction */
  strcpy(title, "Reconstruction                B :");
  nag_file_print_matrix_real_gen_comp(order, matrix, diag, m, n, b, pdb,
                                      "%8.4f", title, Nag_NoLabels, 0,
                                      Nag_NoLabels, 0, 80, 0, 0, &fail);
END:
  NAG_FREE(a);
  NAG_FREE(b);
  NAG_FREE(ca);
  NAG_FREE(cd);
  NAG_FREE(ch);
  NAG_FREE(cv);
  return exit_status;
}