2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
44 #include "gromacs/fileio/espio.h"
45 #include "gromacs/fileio/filenm.h"
46 #include "gromacs/fileio/g96io.h"
47 #include "gromacs/fileio/gmxfio.h"
48 #include "gromacs/fileio/groio.h"
49 #include "gromacs/fileio/pdbio.h"
50 #include "gromacs/fileio/tpxio.h"
51 #include "gromacs/fileio/trx.h"
52 #include "gromacs/fileio/trxio.h"
53 #include "gromacs/math/vec.h"
54 #include "gromacs/topology/atomprop.h"
55 #include "gromacs/topology/atoms.h"
56 #include "gromacs/topology/block.h"
57 #include "gromacs/topology/mtop_util.h"
58 #include "gromacs/topology/symtab.h"
59 #include "gromacs/topology/topology.h"
60 #include "gromacs/utility/cstringutil.h"
61 #include "gromacs/utility/fatalerror.h"
62 #include "gromacs/utility/smalloc.h"
64 void write_sto_conf_indexed(const char *outfile, const char *title,
66 rvec x[], rvec *v, int ePBC, matrix box,
67 atom_id nindex, atom_id index[])
73 ftp = fn2ftp(outfile);
77 out = gmx_fio_fopen(outfile, "w");
78 write_hconf_indexed_p(out, title, atoms, nindex, index, 3, x, v, box);
82 clear_trxframe(&fr, TRUE);
85 fr.natoms = atoms->nr;
96 copy_mat(box, fr.box);
97 out = gmx_fio_fopen(outfile, "w");
98 write_g96_conf(out, &fr, nindex, index);
105 out = gmx_fio_fopen(outfile, "w");
106 write_pdbfile_indexed(out, title, atoms, x, ePBC, box, ' ', -1, nindex, index, NULL, TRUE);
110 out = gmx_fio_fopen(outfile, "w");
111 write_espresso_conf_indexed(out, title, atoms, nindex, index, x, v, box);
115 gmx_fatal(FARGS, "Sorry, can not write a topology to %s", outfile);
118 gmx_incons("Not supported in write_sto_conf_indexed");
122 void write_sto_conf(const char *outfile, const char *title, t_atoms *atoms,
123 rvec x[], rvec *v, int ePBC, matrix box)
129 ftp = fn2ftp(outfile);
133 write_conf_p(outfile, title, atoms, 3, x, v, box);
136 clear_trxframe(&fr, TRUE);
139 fr.natoms = atoms->nr;
150 copy_mat(box, fr.box);
151 out = gmx_fio_fopen(outfile, "w");
152 write_g96_conf(out, &fr, -1, NULL);
158 out = gmx_fio_fopen(outfile, "w");
159 write_pdbfile(out, title, atoms, x, ePBC, box, ' ', -1, NULL, TRUE);
163 out = gmx_fio_fopen(outfile, "w");
164 write_espresso_conf_indexed(out, title, atoms, atoms->nr, NULL, x, v, box);
168 gmx_fatal(FARGS, "Sorry, can not write a topology to %s", outfile);
171 gmx_incons("Not supported in write_sto_conf");
175 void write_sto_conf_mtop(const char *outfile, const char *title,
177 rvec x[], rvec *v, int ePBC, matrix box)
183 ftp = fn2ftp(outfile);
187 out = gmx_fio_fopen(outfile, "w");
188 write_hconf_mtop(out, title, mtop, 3, x, v, box);
192 /* This is a brute force approach which requires a lot of memory.
193 * We should implement mtop versions of all writing routines.
195 atoms = gmx_mtop_global_atoms(mtop);
197 write_sto_conf(outfile, title, &atoms, x, v, ePBC, box);
204 void get_stx_coordnum(const char *infile, int *natoms)
207 int ftp, tpxver, tpxgen;
209 char g96_line[STRLEN+1];
211 ftp = fn2ftp(infile);
212 range_check(ftp, 0, efNR);
216 get_coordnum(infile, natoms);
219 in = gmx_fio_fopen(infile, "r");
226 *natoms = read_g96_conf(in, infile, &fr, g96_line);
232 in = gmx_fio_fopen(infile, "r");
233 get_pdb_coordnum(in, natoms);
237 *natoms = get_espresso_coordnum(infile);
243 read_tpxheader(infile, &tpx, TRUE, &tpxver, &tpxgen);
244 *natoms = tpx.natoms;
248 gmx_fatal(FARGS, "File type %s not supported in get_stx_coordnum",
253 static void tpx_make_chain_identifiers(t_atoms *atoms, t_block *mols)
259 /* We always assign a new chain number, but save the chain id characters
260 * for larger molecules.
262 #define CHAIN_MIN_ATOMS 15
266 for (m = 0; m < mols->nr; m++)
269 a1 = mols->index[m+1];
270 if ((a1-a0 >= CHAIN_MIN_ATOMS) && (chainid <= 'Z'))
279 for (a = a0; a < a1; a++)
281 atoms->resinfo[atoms->atom[a].resind].chainnum = chainnum;
282 atoms->resinfo[atoms->atom[a].resind].chainid = c;
287 /* Blank out the chain id if there was only one chain */
290 for (r = 0; r < atoms->nres; r++)
292 atoms->resinfo[r].chainid = ' ';
297 void read_stx_conf(const char *infile, char *title, t_atoms *atoms,
298 rvec x[], rvec *v, int *ePBC, matrix box)
305 char g96_line[STRLEN+1];
309 fprintf(stderr, "Warning: Number of atoms in %s is 0\n", infile);
311 else if (atoms->atom == NULL)
313 gmx_mem("Uninitialized array atom");
321 ftp = fn2ftp(infile);
325 read_whole_conf(infile, title, atoms, x, v, box);
329 fr.natoms = atoms->nr;
334 in = gmx_fio_fopen(infile, "r");
335 read_g96_conf(in, infile, &fr, g96_line);
337 copy_mat(fr.box, box);
338 std::strncpy(title, fr.title, STRLEN);
343 read_pdb_conf(infile, title, atoms, x, ePBC, box, TRUE, NULL);
346 read_espresso_conf(infile, title, atoms, x, v, box);
350 i = read_tpx(infile, NULL, box, &natoms, x, v, NULL, mtop);
356 strcpy(title, *(mtop->name));
358 /* Free possibly allocated memory */
361 *atoms = gmx_mtop_global_atoms(mtop);
362 top = gmx_mtop_t_to_t_topology(mtop);
363 tpx_make_chain_identifiers(atoms, &top.mols);
366 /* The strings in the symtab are still in use in the returned t_atoms
367 * structure, so we should not free them. But there is no place to put the
368 * symbols; the only choice is to leak the memory...
369 * So we clear the symbol table before freeing the topology structure. */
370 free_symtab(&top.symtab);
375 gmx_incons("Not supported in read_stx_conf");
379 static void done_gmx_groups_t(gmx_groups_t *g)
383 for (i = 0; (i < egcNR); i++)
385 if (NULL != g->grps[i].nm_ind)
387 sfree(g->grps[i].nm_ind);
388 g->grps[i].nm_ind = NULL;
390 if (NULL != g->grpnr[i])
396 /* The contents of this array is in symtab, don't free it here */
400 gmx_bool read_tps_conf(const char *infile, char *title, t_topology *top, int *ePBC,
401 rvec **x, rvec **v, matrix box, gmx_bool bMass)
404 int natoms, i, version, generation;
405 gmx_bool bTop, bXNULL = FALSE;
409 bTop = fn2bTPX(infile);
413 read_tpxheader(infile, &header, TRUE, &version, &generation);
416 snew(*x, header.natoms);
420 snew(*v, header.natoms);
423 *ePBC = read_tpx(infile, NULL, box, &natoms,
424 (x == NULL) ? NULL : *x, (v == NULL) ? NULL : *v, NULL, mtop);
425 *top = gmx_mtop_t_to_t_topology(mtop);
426 /* In this case we need to throw away the group data too */
427 done_gmx_groups_t(&mtop->groups);
429 std::strcpy(title, *top->name);
430 tpx_make_chain_identifiers(&top->atoms, &top->mols);
434 get_stx_coordnum(infile, &natoms);
435 init_t_atoms(&top->atoms, natoms, (fn2ftp(infile) == efPDB));
446 read_stx_conf(infile, title, &top->atoms, *x, (v == NULL) ? NULL : *v, ePBC, box);
454 aps = gmx_atomprop_init();
455 for (i = 0; (i < natoms); i++)
457 if (!gmx_atomprop_query(aps, epropMass,
458 *top->atoms.resinfo[top->atoms.atom[i].resind].name,
459 *top->atoms.atomname[i],
460 &(top->atoms.atom[i].m)))
464 fprintf(debug, "Can not find mass for atom %s %d %s, setting to 1\n",
465 *top->atoms.resinfo[top->atoms.atom[i].resind].name,
466 top->atoms.resinfo[top->atoms.atom[i].resind].nr,
467 *top->atoms.atomname[i]);
471 gmx_atomprop_destroy(aps);
473 top->idef.ntypes = -1;