/*
Copyright (C) 2003 by Sean David Fleming

sean@ivec.org

This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License
as published by the Free Software Foundation; either version 2
of the License, or (at your option) any later version.

This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

You should have received a copy of the GNU General Public License
along with this program; if not, write to the Free Software
Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.

The GNU GPL can also be found at http://www.gnu.org
*/

#include <stdio.h>
#include <string.h>
#include <time.h>

#include "gdis.h"
#include "coords.h"
#include "model.h"
#include "file.h"
#include "parse.h"
#include "matrix.h"
#include "interface.h"

#define DEBUG_MORE 0
#define MAX_KEYS 15

/* main structures */
extern struct sysenv_pak sysenv;
extern struct elem_pak elements[];


/****************************/
/* write a coordinate block */
/****************************/
gint write_arc_frame(FILE *fp, struct model_pak *model)
{
gint m, flag;
gdouble x[3];
GSList *clist, *mlist;
struct core_pak *core;
struct mol_pak *mol;
time_t t1;

/* get & print the date (streamlined slightly by sxm) */
t1 = time(NULL);
fprintf(fp,"!DATE %s",ctime(&t1));

if (model->periodic)
  {
  if (model->periodic == 2)
    {
/* saving a surface as a 3D model - make depth large enough to fit everything */
    fprintf(fp,"PBC %9.4f %9.4f %9.4f %9.4f %9.4f %9.4f (P1)\n",
                model->pbc[0], model->pbc[1], 2.0*model->rmax,
                180.0*model->pbc[3]/PI,
                180.0*model->pbc[4]/PI,
                180.0*model->pbc[5]/PI);
    }
  else
    {
    fprintf(fp,"PBC %9.4f %9.4f %9.4f %9.4f %9.4f %9.4f (P1)\n",
                model->pbc[0], model->pbc[1], model->pbc[2],
                180.0*model->pbc[3]/PI,
                180.0*model->pbc[4]/PI,
                180.0*model->pbc[5]/PI);
    }
  }

m=0;
for (mlist=model->moles ; mlist ; mlist=g_slist_next(mlist))
  {
  mol = (struct mol_pak *) mlist->data;

/* flag - atom written to file? (ie if UNDELETED) */
  flag=0;
/* m - number of atoms in molecule */
  for (clist=mol->cores ; clist ; clist=g_slist_next(clist))
    {
    core = (struct core_pak *) clist->data;
    if (core->status & DELETED)
      continue;

/* everything is cartesian after latmat mult */
    ARR3SET(x, core->x);
    vecmat(model->latmat, x);
/* unique molecule numbering for BIOSYM files (>0) */
    if (core->atom_type)
      {
      fprintf(fp, "%-4s  %14.9f %14.9f %14.9f CORE %-6d %-7s %-2s ",
                  core->atom_label,x[0],x[1],x[2],m+1, core->atom_type,
                  elements[core->atom_code].symbol);
      }
    else
      {
      fprintf(fp, "%-4s  %14.9f %14.9f %14.9f CORE %-6d ?       %-2s ",
                  core->atom_label,x[0],x[1],x[2],m+1, 
                  elements[core->atom_code].symbol);
      }

  if (core->charge < 0.0)
    fprintf(fp, "%5.3f\n", core->charge);
  else
    fprintf(fp, " %4.3f\n", core->charge);

/* indicate we've written at least one atom */
    flag++;
    }
/* if at least one atom in the molecule was written - inc mol numbering */
  if (flag)
    {
    fprintf(fp,"end\n");
    m++;
    }
  }
fprintf(fp,"end\n");

return(0);
}

/***********************/
/* write biosym header */
/***********************/
gint write_arc_header(FILE *fp, struct model_pak *model)
{
/* print header */
fprintf(fp,"!BIOSYM archive 3\n");
if (model->periodic)
  fprintf(fp,"PBC=ON\n");
else
  fprintf(fp,"PBC=OFF\n");
gdis_blurb(fp);

return(0);
}

/******************/
/* Biosym writing */
/******************/
gint write_arc(gchar *filename, struct model_pak *data)
{
FILE *fp;

/* checks */
g_return_val_if_fail(data != NULL, 1);
g_return_val_if_fail(filename != NULL, 2);

/* open the file */
fp = fopen(filename, "wt");
if (!fp)
  return(3);

/* TODO - write multiple frames if it has them */
write_arc_header(fp, data);
write_arc_frame(fp, data);

fclose(fp);
return(0);
}

/****************************************************************/
/* determine if a BIOSYM fragment label was generated by MARVIN */
/****************************************************************/
gint is_marvin_label(const gchar *label)
{
if (label[0] != 'R' && label[0] != 'r')
  return(FALSE);

if (!g_ascii_isdigit(label[1]))
  return(FALSE);

if (label[2] != 'A' && label[2] != 'a' && label[2] != 'B' && label[2] != 'b')
  return(FALSE);

if (label[3] != 'C' && label[3] != 'c' && label[3] != 'S' && label[3] != 's')
  return(FALSE);

return(TRUE);
}

gint marvin_region(const gchar *label)
{
return(label[1] - '1');
}

gint marvin_core(const gchar *label)
{
if (label[3] == 'C' || label[3] == 'c')
  return(TRUE);
return(FALSE);
}

/*************************************************/
/* read a car/arc block into the model structure */
/*************************************************/
/* NB: assumes fp is at a !DATE line */
void read_arc_block(FILE *fp, struct model_pak *data)
{
gint region, core_flag, end_count, num_tokens;
gchar **buff;
GSList *clist, *slist;
struct core_pak *core;
struct shel_pak *shel;

g_assert(fp != NULL);
g_assert(data != NULL);

/* FIXME - this breaks measurement updates */
/*
free_core_list(data);
*/

/* init */
clist = data->cores;
slist = data->shels;
end_count=0;

/* loop while there's data */
for (;;)
  {
  buff = get_tokenized_line(fp, &num_tokens);
  if (!buff)
    break;

/* increment/reset end counter according to input line */
  if (g_ascii_strncasecmp("end", *buff, 3) == 0)
    end_count++;
  else
    end_count = 0;

/* skip single end, terminate on double end */
  if (end_count == 1)
    {
    g_strfreev(buff);
    continue;
    }
  if (end_count == 2)
    break;

/* cell dimension search */
  if (g_ascii_strncasecmp("pbc", *buff, 3) == 0)
    {
    if (num_tokens > 6)
      {
      data->pbc[0] = str_to_float(*(buff+1));
      data->pbc[1] = str_to_float(*(buff+2));
      data->pbc[2] = str_to_float(*(buff+3));
      data->pbc[3] = PI*str_to_float(*(buff+4))/180.0;
      data->pbc[4] = PI*str_to_float(*(buff+5))/180.0;
      data->pbc[5] = PI*str_to_float(*(buff+6))/180.0;
      }
/* cope with 2D */
    if (data->pbc[2] == 0)
      data->periodic = 2;
    continue;
    }

/* process a single coord line */
  if (num_tokens < 8) 
    {
    gui_text_show(ERROR, "unexpected end of file reading arc file\n");
    g_strfreev(buff);
    return;
    }

/* process core/shell candidates */
  if (num_tokens > 8)
  if (elem_symbol_test(*(buff+7)))
    {
    core_flag=TRUE;
    region = 0;

/* MARVIN core/shell/region parse */ 
    if (is_marvin_label(*(buff+4)))
      {
      if (!marvin_core(*(buff+4)))
        core_flag = FALSE;
      region = marvin_region(*(buff+4));
      data->region_empty[region] = FALSE;
      }

/* overwrite existing core (if any) */
    if (core_flag)
      {
      if (clist)
        {
        core = (struct core_pak *) clist->data;
        clist = g_slist_next(clist);
        }
      else
        {
/* otherwise, add a new core */
/* NB: must append (not prepend) since we may overwrite & the order must be the same */
/* TODO - we could prepend cores (it's faster) and then reverse the list at the end */
/* but it's more complicated as we should ONLY reverse newly added cores */
        core = new_core(*(buff+7), data);
        data->cores = g_slist_append(data->cores, core);
        }

/* set the proper label */
g_free(core->atom_label);
core->atom_label = g_strdup(*buff);

/* set values */
      core->x[0] = str_to_float(*(buff+1));
      core->x[1] = str_to_float(*(buff+2));
      core->x[2] = str_to_float(*(buff+3));
      core->charge = str_to_float(*(buff+8));

g_free(core->atom_type);
core->atom_type = g_strdup(*(buff+6));

      core->region = region;
      core->lookup_charge = FALSE;
      }
    else
      {
/* overwrite existing shell (if any) */
      if (slist)
        {
        shel = (struct shel_pak *) slist->data;
        slist = g_slist_next(slist);
        }
      else
        {
/* otherwise, add a new shell */
        shel = new_shell(*(buff+7), data);
        data->shels = g_slist_append(data->shels, shel);
        }

/* set the proper label */
g_free(shel->shell_label);
shel->shell_label = g_strdup(*buff);

/* set values */
      shel->x[0] = str_to_float(*(buff+1));
      shel->x[1] = str_to_float(*(buff+2));
      shel->x[2] = str_to_float(*(buff+3));
      shel->charge = str_to_float(*(buff+8));
      shel->region = region;
      shel->lookup_charge = FALSE;
      }
    }
  g_strfreev(buff);
  }

/* convert input cartesian coords to fractional */
if( data->periodic)
  {
  matrix_lattice_init(data);
  coords_make_fractional(data);
  }
data->fractional = TRUE;

g_strfreev(buff);
}

/*********************/
/* biosym frame read */
/*********************/
gint read_arc_frame(FILE *fp, struct model_pak *data)
{
/* frame overwrite */
read_arc_block(fp, data);
return(0);
}

/******************/
/* Biosym reading */
/******************/
#define DEBUG_READ_ARC 0
gint read_arc(gchar *filename, struct model_pak *data)
{
gint frame=0, num_tokens;
gchar **buff;
FILE *fp;

/* checks */
g_return_val_if_fail(data != NULL, 1);
g_return_val_if_fail(filename != NULL, 2);

fp = fopen(filename, "rt");
if (!fp)
  return(3);

/* loop while there's data */
for (;;)
  {
  buff = get_tokenized_line(fp, &num_tokens);
  if (!buff)
    {
    g_strfreev(buff);
    break;
    }

/* pbc on/off search */
  if (g_ascii_strncasecmp("pbc=on", *buff, 6) == 0)
    data->periodic = 3;
  if (g_ascii_strncasecmp("pbc=off", *buff, 7) == 0)
    data->periodic = 0;

/* coords search */
  if (g_ascii_strncasecmp("!date", *buff, 5) == 0)
    {
/* NEW */
    add_frame_offset(fp, data);

/* only load the required frame */
    if (frame == data->cur_frame)
      read_arc_block(fp, data);

/* increment counter */
    frame++;
    }
  g_strfreev(buff);
  }

/* got everything */
data->num_frames = frame;

/* get rid of frame list if only one frame */
if (data->num_frames == 1)
  {
  free_list(data->frame_list);
  data->frame_list = NULL;
  }

/* model setup */
strcpy(data->filename, filename);
g_free(data->basename);
data->basename = strdup_basename(filename);

model_prep(data);

return(0);
}