/*
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);
}
syntax highlighted by Code2HTML, v. 0.9.1