#include "confio.h"
-#include <cerrno>
-#include <cmath>
#include <cstdio>
-#include <cstdlib>
#include <cstring>
-#include <algorithm>
-
+#include "gromacs/fileio/espio.h"
#include "gromacs/fileio/filenm.h"
#include "gromacs/fileio/g96io.h"
#include "gromacs/fileio/gmxfio.h"
#include "gromacs/fileio/trx.h"
#include "gromacs/fileio/trxio.h"
#include "gromacs/math/vec.h"
-#include "gromacs/pbcutil/pbc.h"
#include "gromacs/topology/atomprop.h"
#include "gromacs/topology/atoms.h"
#include "gromacs/topology/block.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/smalloc.h"
-static int get_espresso_word(FILE *fp, char word[])
-{
- int ret, nc, i;
-
- ret = 0;
- nc = 0;
-
- do
- {
- i = fgetc(fp);
- if (i != EOF)
- {
- if (i == ' ' || i == '\n' || i == '\t')
- {
- if (nc > 0)
- {
- ret = 1;
- }
- }
- else if (i == '{')
- {
- if (nc == 0)
- {
- word[nc++] = '{';
- }
- ret = 2;
- }
- else if (i == '}')
- {
- if (nc == 0)
- {
- word[nc++] = '}';
- }
- ret = 3;
- }
- else
- {
- word[nc++] = (char)i;
- }
- }
- }
- while (i != EOF && ret == 0);
-
- word[nc] = '\0';
-
- /* printf("'%s'\n",word); */
-
- return ret;
-}
-
-static int check_open_parenthesis(FILE *fp, int r,
- const char *infile, const char *keyword)
-{
- int level_inc;
- char word[STRLEN];
-
- level_inc = 0;
- if (r == 2)
- {
- level_inc++;
- }
- else
- {
- r = get_espresso_word(fp, word);
- if (r == 2)
- {
- level_inc++;
- }
- else
- {
- gmx_fatal(FARGS, "Expected '{' after '%s' in file '%s'",
- keyword, infile);
- }
- }
-
- return level_inc;
-}
-
-static int check_close_parenthesis(FILE *fp, int r,
- const char *infile, const char *keyword)
-{
- int level_inc;
- char word[STRLEN];
-
- level_inc = 0;
- if (r == 3)
- {
- level_inc--;
- }
- else
- {
- r = get_espresso_word(fp, word);
- if (r == 3)
- {
- level_inc--;
- }
- else
- {
- gmx_fatal(FARGS, "Expected '}' after section '%s' in file '%s'",
- keyword, infile);
- }
- }
-
- return level_inc;
-}
-
-enum {
- espID, espPOS, espTYPE, espQ, espV, espF, espMOLECULE, espNR
-};
-const char *esp_prop[espNR] = {
- "id", "pos", "type", "q", "v", "f",
- "molecule"
-};
-
-static void read_espresso_conf(const char *infile, char *title,
- t_atoms *atoms, rvec x[], rvec *v, matrix box)
-{
- t_symtab *symtab = NULL;
- FILE *fp;
- char word[STRLEN], buf[STRLEN];
- int level, r, nprop, p, i, m, molnr;
- int prop[32];
- double d;
- gmx_bool bFoundParticles, bFoundProp, bFoundVariable, bMol;
-
- if (!symtab)
- {
- snew(symtab, 1);
- open_symtab(symtab);
- }
- // TODO: The code does not understand titles it writes...
- title[0] = '\0';
-
- clear_mat(box);
-
- fp = gmx_fio_fopen(infile, "r");
-
- bFoundParticles = FALSE;
- bFoundVariable = FALSE;
- bMol = FALSE;
- level = 0;
- while ((r = get_espresso_word(fp, word)))
- {
- if (level == 1 && std::strcmp(word, "particles") == 0 && !bFoundParticles)
- {
- bFoundParticles = TRUE;
- level += check_open_parenthesis(fp, r, infile, "particles");
- nprop = 0;
- while (level == 2 && (r = get_espresso_word(fp, word)))
- {
- bFoundProp = FALSE;
- for (p = 0; p < espNR; p++)
- {
- if (strcmp(word, esp_prop[p]) == 0)
- {
- bFoundProp = TRUE;
- prop[nprop++] = p;
- /* printf(" prop[%d] = %s\n",nprop-1,esp_prop[prop[nprop-1]]); */
- }
- }
- if (!bFoundProp && word[0] != '}')
- {
- gmx_fatal(FARGS, "Can not read Espresso files with particle property '%s'", word);
- }
- if (bFoundProp && p == espMOLECULE)
- {
- bMol = TRUE;
- }
- if (r == 3)
- {
- level--;
- }
- }
-
- i = 0;
- while (level > 0 && (r = get_espresso_word(fp, word)))
- {
- if (r == 2)
- {
- level++;
- }
- else if (r == 3)
- {
- level--;
- }
- if (level == 2)
- {
- for (p = 0; p < nprop; p++)
- {
- switch (prop[p])
- {
- case espID:
- r = get_espresso_word(fp, word);
- /* Not used */
- break;
- case espPOS:
- for (m = 0; m < 3; m++)
- {
- r = get_espresso_word(fp, word);
- sscanf(word, "%lf", &d);
- x[i][m] = d;
- }
- break;
- case espTYPE:
- r = get_espresso_word(fp, word);
- atoms->atom[i].type = std::strtol(word, NULL, 10);
- break;
- case espQ:
- r = get_espresso_word(fp, word);
- sscanf(word, "%lf", &d);
- atoms->atom[i].q = d;
- break;
- case espV:
- for (m = 0; m < 3; m++)
- {
- r = get_espresso_word(fp, word);
- sscanf(word, "%lf", &d);
- v[i][m] = d;
- }
- break;
- case espF:
- for (m = 0; m < 3; m++)
- {
- r = get_espresso_word(fp, word);
- /* not used */
- }
- break;
- case espMOLECULE:
- r = get_espresso_word(fp, word);
- molnr = std::strtol(word, NULL, 10);
- if (i == 0 ||
- atoms->resinfo[atoms->atom[i-1].resind].nr != molnr)
- {
- atoms->atom[i].resind =
- (i == 0 ? 0 : atoms->atom[i-1].resind+1);
- atoms->resinfo[atoms->atom[i].resind].nr = molnr;
- atoms->resinfo[atoms->atom[i].resind].ic = ' ';
- atoms->resinfo[atoms->atom[i].resind].chainid = ' ';
- atoms->resinfo[atoms->atom[i].resind].chainnum = molnr; /* Not sure if this is right? */
- }
- else
- {
- atoms->atom[i].resind = atoms->atom[i-1].resind;
- }
- break;
- }
- }
- /* Generate an atom name from the particle type */
- sprintf(buf, "T%d", atoms->atom[i].type);
- atoms->atomname[i] = put_symtab(symtab, buf);
- if (bMol)
- {
- if (i == 0 || atoms->atom[i].resind != atoms->atom[i-1].resind)
- {
- atoms->resinfo[atoms->atom[i].resind].name =
- put_symtab(symtab, "MOL");
- }
- }
- else
- {
- /* Residue number is the atom number */
- atoms->atom[i].resind = i;
- /* Generate an residue name from the particle type */
- if (atoms->atom[i].type < 26)
- {
- sprintf(buf, "T%c", 'A'+atoms->atom[i].type);
- }
- else
- {
- sprintf(buf, "T%c%c",
- 'A'+atoms->atom[i].type/26, 'A'+atoms->atom[i].type%26);
- }
- t_atoms_set_resinfo(atoms, i, symtab, buf, i, ' ', 0, ' ');
- }
-
- if (r == 3)
- {
- level--;
- }
- i++;
- }
- }
- atoms->nres = atoms->nr;
-
- if (i != atoms->nr)
- {
- gmx_fatal(FARGS, "Internal inconsistency in Espresso routines, read %d atoms, expected %d atoms", i, atoms->nr);
- }
- }
- else if (level == 1 && std::strcmp(word, "variable") == 0 && !bFoundVariable)
- {
- bFoundVariable = TRUE;
- level += check_open_parenthesis(fp, r, infile, "variable");
- while (level == 2 && (r = get_espresso_word(fp, word)))
- {
- if (level == 2 && std::strcmp(word, "box_l") == 0)
- {
- for (m = 0; m < 3; m++)
- {
- r = get_espresso_word(fp, word);
- sscanf(word, "%lf", &d);
- box[m][m] = d;
- }
- level += check_close_parenthesis(fp, r, infile, "box_l");
- }
- }
- }
- else if (r == 2)
- {
- level++;
- }
- else if (r == 3)
- {
- level--;
- }
- }
-
- if (!bFoundParticles)
- {
- fprintf(stderr, "Did not find a particles section in Espresso file '%s'\n",
- infile);
- }
-
- gmx_fio_fclose(fp);
-}
-
-static int get_espresso_coordnum(const char *infile)
-{
- FILE *fp;
- char word[STRLEN];
- int natoms, level, r;
- gmx_bool bFoundParticles;
-
- natoms = 0;
-
- fp = gmx_fio_fopen(infile, "r");
-
- bFoundParticles = FALSE;
- level = 0;
- while ((r = get_espresso_word(fp, word)) && !bFoundParticles)
- {
- if (level == 1 && strcmp(word, "particles") == 0 && !bFoundParticles)
- {
- bFoundParticles = TRUE;
- level += check_open_parenthesis(fp, r, infile, "particles");
- while (level > 0 && (r = get_espresso_word(fp, word)))
- {
- if (r == 2)
- {
- level++;
- if (level == 2)
- {
- natoms++;
- }
- }
- else if (r == 3)
- {
- level--;
- }
- }
- }
- else if (r == 2)
- {
- level++;
- }
- else if (r == 3)
- {
- level--;
- }
- }
- if (!bFoundParticles)
- {
- fprintf(stderr, "Did not find a particles section in Espresso file '%s'\n",
- infile);
- }
-
- gmx_fio_fclose(fp);
-
- return natoms;
-}
-
-static void write_espresso_conf_indexed(FILE *out, const char *title,
- t_atoms *atoms, int nx, atom_id *index,
- rvec *x, rvec *v, matrix box)
-{
- int i, j;
-
- fprintf(out, "# %s\n", title);
- if (TRICLINIC(box))
- {
- gmx_warning("The Espresso format does not support triclinic unit-cells");
- }
- fprintf(out, "{variable {box_l %f %f %f}}\n", box[0][0], box[1][1], box[2][2]);
-
- fprintf(out, "{particles {id pos type q%s}\n", v ? " v" : "");
- for (i = 0; i < nx; i++)
- {
- if (index)
- {
- j = index[i];
- }
- else
- {
- j = i;
- }
- fprintf(out, "\t{%d %f %f %f %d %g",
- j, x[j][XX], x[j][YY], x[j][ZZ],
- atoms->atom[j].type, atoms->atom[j].q);
- if (v)
- {
- fprintf(out, " %f %f %f", v[j][XX], v[j][YY], v[j][ZZ]);
- }
- fprintf(out, "}\n");
- }
- fprintf(out, "}\n");
-}
-
void write_sto_conf_indexed(const char *outfile, const char *title,
t_atoms *atoms,
rvec x[], rvec *v, int ePBC, matrix box,