/*
 * paramxp.c
 *
 * Summary:  Reader for XPLOR force field parameter files.
 *
 * Author:  David Hardy
 */

#include <stdio.h>
#include <string.h>
#include <ctype.h>
#include <math.h>
#include "param.h"
#include "paramxp.h"
#include "pdata.h"
#include "ihash.h"
#include "error.h"


#undef  WRDLEN
#undef  MSGLEN

#define  WRDLEN   127
#define  MSGLEN   255
#define  BUFSIZE  512

#undef  BOND
#undef  ANGLE
#undef  DIHED
#undef  IMPR
#undef  NONB
#undef  NBFIX
#undef  UNUSED

#define  BOND    0
#define  ANGLE   1
#define  DIHED   2
#define  IMPR    3
#define  NONB    4
#define  NBFIX   5
#define  UNUSED  6

typedef struct Keywdval_Tag {
  const char *keywd;
  int val;
} Keywdval;

static const Keywdval KEYWD[] = {
  { "REMARK",     UNUSED },
  { "REMARKS",    UNUSED },
  { "SET",        UNUSED },
  { "END",        UNUSED },
  { "RESET",      UNUSED },
  { "BOND",       BOND   },
  { "ANGL",       ANGLE  },
  { "ANGLE",      ANGLE  },
  { "DIHE",       DIHED  },
  { "DIHEDRAL",   DIHED  },
  { "IMPR",       IMPR   },
  { "IMPROPER",   IMPR   },
  { "NONB",       NONB   },
  { "NONBONDED",  NONB   },
  { "NBF",        NBFIX  },
  { "NBFIX",      NBFIX  },
  { "NBON",       UNUSED },
  { "NBONDS",     UNUSED },
  { "HBON",       UNUSED },
  { "HBONDED",    UNUSED },
  { "HBONDS",     UNUSED },
  { "CUTN",       UNUSED },
  { "CUTNB",      UNUSED },
  { "GROU",       UNUSED },
  { "GROUP",      UNUSED },
  { "ATOM",       UNUSED },
  { "TOLE",       UNUSED },
  { "TOLERANCE",  UNUSED },
  { "WMIN",       UNUSED },
  { "CTOFNB",     UNUSED },
  { "CTONNB",     UNUSED },
  { "TRUN",       UNUSED },
  { "TRUNCATION", UNUSED },
  { "SWI",        UNUSED },
  { "SWITCH",     UNUSED },
  { "SHIF",       UNUSED },
  { "SHIFT",      UNUSED },
  { "CDIE",       UNUSED },
  { "RDIE",       UNUSED },
  { "E14F",       UNUSED },
  { "E14FAC",     UNUSED },
  { "EPS",        UNUSED },
  { "VSWI",       UNUSED },
  { "VSWITCH",    UNUSED },
  { "REPE1",      UNUSED },
  { "REXP",       UNUSED },
  { "REXPONENT",  UNUSED },
  { "IREX",       UNUSED },
  { "IREXPONENT", UNUSED },
  { "RCON",       UNUSED },
  { "RCONST",     UNUSED },
  { "NBXM",       UNUSED },
  { "NBXMOD",     UNUSED },
  { "INHI",       UNUSED },
  { "INHIBIT",    UNUSED },
  { "DCUT",       UNUSED },
  { "DOFF",       UNUSED },
  { "DON",        UNUSED },
  { "ACUT",       UNUSED },
  { "AOFF",       UNUSED },
  { "AON",        UNUSED },
  { "ACCE",       UNUSED },
  { "ACCEPTOR",   UNUSED },
  { "AEXP",       UNUSED },
  { "HAEX",       UNUSED },
  { "AAEX",       UNUSED }
};

static const int KEYWD_LEN = 
    (int) (sizeof(KEYWD) / sizeof(KEYWD[0]));


int paramxp_init(Paramxp *p, Pdata *db, FILE *f)
{
  int i;
  Ihash *h = &(p->table);

  p->status = 0;
  MDIO_ASSERT(f != NULL);
  p->f = f;
  MDIO_ASSERT(db != NULL);
  p->pdata = db;
  if (ihash_init(h, 2 * KEYWD_LEN) == IHASH_MEMALLOC) {
    mdio_error("ihash_init() failed, file %s, line %d\n", __FILE__, __LINE__);
    p->status |= MDIO_ERR_MEMALLOC;
    return -1;
  }
  for (i = 0;  i < KEYWD_LEN;  i++) {
    if (ihash_insert(h, KEYWD[i].keywd, KEYWD[i].val) == IHASH_MEMALLOC) {
      mdio_error("ihash_insert() failed, file %s, line %d\n",__FILE__,__LINE__);
      p->status |= MDIO_ERR_MEMALLOC;
      return -1;
    }
  }
  return 0;
}


void paramxp_destroy(Paramxp *p)
{
  ihash_destroy(&(p->table));
  p->pdata = NULL;
  p->f = NULL;
  p->status = 0;
}


int paramxp_status(Paramxp *p)
{
  return p->status;
}


void paramxp_clear_status(Paramxp *p)
{
  if (ferror(p->f) != 0) {
    clearerr(p->f);
  }
  p->status = 0;
}


static int get_cooked_line(char *buf, int buflen,
                           FILE *f, Ihash *h, int *status)
{
  int c, i, match = -1;
  char s[WRDLEN+1];

  /* loop until we match a recognized param word */
  c = fgetc(f);
  while ( c != EOF  &&  match == -1 ) {

    /* skip any comments and leading space */
    while ( isspace(c)  ||  c == '!'  ||  c == '{' ) {
      if (c == '!') {
        while ( (c = fgetc(f)) != EOF  &&  c != '\n' )  ;
      }
      else if (c == '{') {
        while ( (c = fgetc(f)) != EOF  &&  c != '}' )  ;
        if (c == EOF) {
          mdio_warn("unmatched comment found in parameter file, " \
              "unexpected EOF; file %s, line %d\n", __FILE__, __LINE__);
          *status |= MDIO_WARN_COMMENT;
        }
        else  c = fgetc(f);  /* skip the '}' character */
      }
      else  c = fgetc(f);
    }
    if (c == EOF) {
      buf[0] = '\0';
      return EOF;
    }

    /* get the first word */
    i = 0;
    s[i++] = (char) toupper(c);
    while ( (c = fgetc(f)) != EOF && isalnum(c) && i < WRDLEN ) {
      s[i++] = (char) toupper(c);
    }
    s[i] = '\0';

    /* match word from hash table */
    match = ihash_lookup(h, s);

    /* if match, store rest of line up to buflen characters */
    if (match == IHASH_FAIL) {
      mdio_warn("unrecognized word found in parameter file: %s ; " \
          "file %s, line %d\n", s, __FILE__, __LINE__);
      match = -1;
      *status |= MDIO_WARN_UNKNOWN;
      while ( c != EOF  &&  c != '\n' )  c = fgetc(f);
    }
    else if (match >= UNUSED) {
      match = -1;
      while ( c != EOF  &&  c != '\n' )  c = fgetc(f);
    }
    else {

      /* skip leading spaces */
      while ( isspace(c) )  c = fgetc(f);
      i = 0;

      /* fill the buffer up with anything significant in the line */
      while ( c != EOF  &&  c != '\n'  &&  i < buflen ) {

        /* place in buffer only if not a comment */
        if ( c != '{'  &&  c != '!' ) {
          buf[i++] = (char) c;
          c = fgetc(f);
        }

        /* skip any embedded comment (even if it spans a line) */
        else if (c == '{') {
          while ( (c = fgetc(f)) != EOF  &&  c != '}' )  ;
          if (c == EOF) {
            mdio_warn("unmatched comment found in parameter file, " \
                "unexpected EOF; file %s, line %d\n", __FILE__, __LINE__);
            *status |= MDIO_WARN_COMMENT;
          }
          else  c = fgetc(f);  /* skip the '}' character */
        }

        /* clean up rest of line for any EOL comment */
        else {
          while ( (c = fgetc(f)) != EOF  &&  c != '\n' )  ;
        }
      }
      if (i < buflen) {
        buf[i] = '\0';
      }
      else {
        buf[--i] = '\0';
        mdio_warn("parameter file line too long for input buffer; " \
            "file %s, line %d\n", __FILE__, __LINE__);
        *status |= MDIO_WARN_TRUNC;
      }
    }

  } /* end outer while loop */

  return (match != -1 ? match : EOF);
}


static int get_raw_line(char *buf, int buflen, FILE *f, int *status)
{
  int c, i = 0;

  /* make sure we don't return a "blank" line once comments are removed */
  while (i == 0) {

    /* skip any leading white space */
    while ( (c = fgetc(f)) != EOF  &&  c != '\n'  &&  isspace(c) )  ;

    /* put line into buffer */
    while ( c != EOF  &&  c != '\n'  &&  i < buflen ) {

      /* place in buffer only if not a comment */
      if ( c != '{'  &&  c != '!' ) {
        buf[i++] = (char) c;
        c = fgetc(f);
      }

      /* skip any embedded comment (even if it spans a line) */
      else if (c == '{') {
        while ( (c = fgetc(f)) != EOF  &&  c != '}' )  ;
        if (c == EOF) {
          mdio_warn("unmatched comment found in parameter file, " \
              "unexpected EOF; file %s, line %d\n", __FILE__, __LINE__);
          *status |= MDIO_WARN_COMMENT;
        }
        else  c = fgetc(f);  /* skip the '}' character */
      }

      /* clean up rest of line for any EOL comment */
      else {
        while ( (c = fgetc(f)) != EOF  &&  c != '\n' )  ;
      }
    }
  }
  if (i < buflen) {
    buf[i] = '\0';
  }
  else {
    buf[--i] = '\0';
    mdio_warn("parameter file line too long for input buffer; " \
        "file %s, line %d\n", __FILE__, __LINE__);
    *status |= MDIO_WARN_TRUNC;
  }
  return ( i == 0  &&  c == EOF  ?  EOF  :  0 );
}


static int read_bond(Pdata *db, const char *buf, int *status)
{
  MD_Bond_Param bond;
  float k, r0;
  int cnt, retval;
  char atom1[5], atom2[5];

  cnt = sscanf(buf, "%4s %4s %f %f", atom1, atom2, &k, &r0);
  if (cnt < 4) {
    mdio_error("bond format error in parameter file: %s\n", buf);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }
  bond.k = k;
  bond.r0 = r0;
  retval = pdata_add_bond(db, &bond, atom1, atom2);
  if (retval < 0) {
    mdio_error("can't enter bond into database: %s ; " \
        "pdata_add_bond() failed, file %s, line %d\n",buf,__FILE__,__LINE__);
    *status |= pdata_status(db);
    pdata_clear_status(db);
    return -1;
  }
  return retval;
}


static int read_angle(Pdata *db, const char *buf, int *status)
{
  MD_Angle_Param angle;
  float k_theta, theta0;
  int cnt, retval;
  char a1[5], a2[5], a3[5];

  cnt = sscanf(buf, "%4s %4s %4s %f %f", a1, a2, a3, &k_theta, &theta0);
  if (cnt < 5) {
    mdio_error("angle format error in parameter file: %s\n", buf);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }
  angle.k_theta = k_theta;
  angle.theta0 = theta0 * M_PI / 180.0;   /* convert angle to radians */
  angle.k_ub = 0.0;
  angle.r_ub = 0.0;
  retval = pdata_add_angle(db, &angle, a1, a2, a3);
  if (retval < 0) {
    mdio_error("can\'t enter angle into database: %s ; " \
        "pdata_add_angle() failed, file %s, line %d\n",buf,__FILE__,__LINE__);
    *status |= pdata_status(db);
    pdata_clear_status(db);
    return -1;
  }
  return retval;
}


static int read_dihed(Pdata *db, const char *buf, FILE *f,
                      int *status)
{
  MD_Dihed_Param dihed;
  float k, phi;
  int n, mult, cnt, len, i, retval;
  char buffer[BUFSIZE];
  char a1[5], a2[5], a3[5], a4[5];

  cnt = sscanf(buf, "%4s %4s %4s %4s MULTIPLE=%d %f %d %f", 
      a1, a2, a3, a4, &mult, &k, &n, &phi);
  if (cnt == 4) {
    cnt = sscanf(buf, "%*s %*s %*s %*s %f %d %f", &k, &n, &phi);
    if (cnt < 3) {
      mdio_error("dihedral format error in parameter file: %s\n", buf);
      *status |= MDIO_ERR_FORMAT;
      return -1;
    }
    dihed.k[0] = k;
    dihed.phi[0] = phi * M_PI / 180.0;   /* convert angle to radians */
    dihed.n[0] = n;
    dihed.mult = 1;
    len = 1;
  }
  else if (cnt == 8) {
    dihed.k[0] = k;
    dihed.phi[0] = phi * M_PI / 180.0;   /* convert angle to radians */
    dihed.n[0] = n;
    if (mult > MD_MULT_SIZE) {
      mdio_error("dihedral multiplicity exceeded in parameter file, " \
          "only %d lines used, other lines ignored: %s\n", 
          MD_MULT_SIZE, buf);
      *status |= MDIO_ERR_FORMAT;
      len = MD_MULT_SIZE;
    }
    else {
      len = mult;
    }
    dihed.mult = len;

    for (i = 1;  i < len;  i++) {
      if (get_raw_line(buffer, BUFSIZE, f, status) == EOF) {
        mdio_error("dihedral format error in parameter file, " \
            "unexpected EOF\n");
        *status |= MDIO_ERR_FORMAT;
        return -1;
      }
      cnt = sscanf(buffer, "%f %d %f", &k, &n, &phi);
      if (cnt < 3) {
        mdio_error("dihedral format error in parameter file: %s\n", buffer);
        *status |= MDIO_ERR_FORMAT;
        return -1;
      }
      dihed.k[i] = k;
      dihed.phi[i] = phi * M_PI / 180.0;   /* convert angle to radians */
      dihed.n[i] = n;
    }
    for ( ;  i < mult;  i++) {
      get_raw_line(buffer, BUFSIZE, f, status);  /* skip extra lines */
    }
  }
  else {
    mdio_error("dihedral format error in parameter file: %s\n", buf);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }
  retval = pdata_add_dihed(db, &dihed, a1, a2, a3, a4);
  if (retval < 0) {
    mdio_error("can\'t enter dihedral into database: %s ; " \
        "pdata_add_dihed() failed, file %s, line %d\n",buf,__FILE__,__LINE__);
    *status |= pdata_status(db);
    pdata_clear_status(db);
    return -1;
  }
  return retval;
}


static int read_impr(Pdata *db, const char *buf, FILE *f, int *status)
{
  MD_Impr_Param impr;
  float k, phi;
  int n, mult, cnt, len, i, retval;
  char buffer[BUFSIZE];
  char a1[5], a2[5], a3[5], a4[5];

  cnt = sscanf(buf, "%4s %4s %4s %4s MULTIPLE=%d %f %d %f", 
      a1, a2, a3, a4, &mult, &k, &n, &phi);
  if (cnt == 4) {
    cnt = sscanf(buf, "%*s %*s %*s %*s %f %d %f", &k, &n, &phi);
    if (cnt < 3) {
      mdio_error("improper format error in parameter file: %s\n", buf);
      *status |= MDIO_ERR_FORMAT;
      return -1;
    }
    impr.k[0] = k;
    impr.phi[0] = phi * M_PI / 180.0;   /* convert angle to radians */
    impr.n[0] = n;
    impr.mult = 1;
    len = 1;
  }
  else if (cnt == 8) {
    impr.k[0] = k;
    impr.phi[0] = phi * M_PI / 180.0;   /* convert angle to radians */
    impr.n[0] = n;
    if (mult > MD_MULT_SIZE) {
      mdio_error("improper multiplicity exceeded in parameter file, " \
          "only %d lines used, others ignored skipped: %s\n", 
          MD_MULT_SIZE, buf);
      *status |= MDIO_ERR_FORMAT;
      len = MD_MULT_SIZE;
    }
    else {
      len = mult;
    }
    impr.mult = len;

    for (i = 1;  i < len;  i++) {
      if (get_raw_line(buffer, BUFSIZE, f, status) == EOF) {
        mdio_error("improper format error in parameter file, " \
            "unexpected EOF\n");
        *status |= MDIO_ERR_FORMAT;
        return -1;
      }
      cnt = sscanf(buffer, "%f %d %f", &k, &n, &phi);
      if (cnt < 3) {
        mdio_error("improper format error in parameter file: %s\n", buffer);
        *status |= MDIO_ERR_FORMAT;
        return -1;
      }
      impr.k[i] = k;
      impr.phi[i] = phi * M_PI / 180.0;   /* convert angle to radians */
      impr.n[i] = n;
    }
    for ( ;  i < mult;  i++) {
      get_raw_line(buffer, BUFSIZE, f, status);  /* skip extra lines */
    }
  }
  else {
    mdio_error("improper format error in parameter file: %s\n", buf);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }
  retval = pdata_add_impr(db, &impr, a1, a2, a3, a4);
  if (retval < 0) {
    mdio_error("can\'t enter improper into database: %s ; " \
        "pdata_add_impr() failed, file %s, line %d\n",buf,__FILE__,__LINE__);
    *status |= pdata_status(db);
    pdata_clear_status(db);
    return -1;
  }
  return retval;
}


static int read_nonb(Pdata *db, const char *buf, int *status)
{
  MD_Atom_Param atom;
  float eps, sig, eps14, sig14;
  double sixth_root_of_two = pow(2.0, 1.0/6.0);
  int cnt, retval;
  char atomname[5];

  cnt = sscanf(buf, "%4s %f %f %f %f", atomname, &eps, &sig, &eps14, &sig14);
  if (cnt < 5) {
    mdio_error("nonb format error in parameter file: %s\n", buf);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }
  atom.emin = -eps;
  atom.rmin = sig * sixth_root_of_two;
  atom.emin14 = -eps14;
  atom.rmin14 = sig14 * sixth_root_of_two;
  strcpy(atom.type, atomname);
  retval = pdata_add_atom(db, &atom, atomname);
  if (retval < 0) {
    mdio_error("can\'t enter nonb (atom) into database: %s ; " \
        "pdata_add_atom() failed, file %s, line %d\n",buf,__FILE__,__LINE__);
    *status |= pdata_status(db);
    pdata_clear_status(db);
    return -1;
  }
  return retval;
}


static int read_nbfix(Pdata *db, const char *buf, int *status)
{
  MD_Nbfix_Param nbfix;
  float a, b, a14, b14;
  static double one_sixth = 1.0/6.0;
  int cnt, i, retval;
  char atom1[5], atom2[5];

  cnt = sscanf(buf, "%4s %4s %f %f %f %f", atom1, atom2, &a, &b, &a14, &b14);
  if (cnt < 6) {
    mdio_error("nbfix format error in parameter file: %s\n", buf);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }

  /* assume atom (nonb) params for these atom types have already been read */
  if ((i = pdata_find_atom(db, atom1)) < 0) {
    mdio_error("atom (nonb) _%s_ params not found for _%s_%s_ nbfix, " \
        "pdata_find_atom() failed, file %s, line %d\n", 
        atom1, atom1, atom2, __FILE__, __LINE__);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }
  nbfix.param[0] = i;
  if ((i = pdata_find_atom(db, atom2)) < 0) {
    mdio_error("atom (nonb) _%s_ params not found for _%s_%s_ nbfix, " \
        "pdata_find_atom() failed, file %s, line %d\n",
        atom2, atom1, atom2, __FILE__, __LINE__);
    *status |= MDIO_ERR_FORMAT;
    return -1;
  }
  nbfix.param[1] = i;

  /* set other fields */
  nbfix.emin = -0.25 * b * b / a;
  nbfix.rmin = pow(2.0 * a / b, one_sixth);
  nbfix.emin14 = -0.25 * b14 * b14 / a14;
  nbfix.rmin14 = pow(2.0 * a14 / b14, one_sixth);
  retval = pdata_add_nbfix(db, &nbfix, atom1, atom2);
  if (retval < 0) {
    mdio_error("can\'t enter nbfix into database: %s ; " \
        "pdata_add_nbfix() failed, file %s, line %d\n",buf,__FILE__,__LINE__);
    *status |= pdata_status(db);
    pdata_clear_status(db);
    return -1;
  }
  return retval;
}


int paramxp_read(Paramxp *p)
{
  int keywd, notdone = 1, retval = 0;
  Pdata *pdata = p->pdata;
  FILE *f = p->f;
  Ihash *h = &(p->table);
  int *status = &(p->status);
  char buf[BUFSIZE];

  while (notdone && (retval == 0)) {
    keywd = get_cooked_line(buf, BUFSIZE, f, h, status);
    switch (keywd) {
    case BOND:
      retval = read_bond(pdata, buf, status);
      break;
    case ANGLE:
      retval = read_angle(pdata, buf, status);
      break;
    case DIHED:
      retval = read_dihed(pdata, buf, f, status);
      break;
    case IMPR:
      retval = read_impr(pdata, buf, f, status);
      break;
    case NONB:
      retval = read_nonb(pdata, buf, status);
      break;
    case NBFIX:
      retval = read_nbfix(pdata, buf, status);
      break;
    default:  /*  == EOF  */
      notdone = 0;
    }
  }
  if (retval < 0  &&  ferror(f)) {
    mdio_error("error occured while reading parameter file\n");
    *status |= MDIO_ERR_READ;
  }
  return retval;
}