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,2016, 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.
38 #ifndef GMX_FILEIO_PDBIO_H
39 #define GMX_FILEIO_PDBIO_H
43 #include "gromacs/math/vectypes.h"
44 #include "gromacs/topology/atoms.h"
45 #include "gromacs/utility/basedefinitions.h"
46 #include "gromacs/utility/real.h"
57 typedef struct gmx_conect_t
*gmx_conect
;
59 /* Write a PDB line with an ATOM or HETATM record directly to a file pointer.
61 * Returns the number of characters printed.
64 gmx_fprintf_pdb_atomline(FILE * fp
,
65 enum PDB_record record
,
67 const char * atom_name
,
68 char alternate_location
,
69 const char * res_name
,
72 char res_insertion_code
,
78 const char * element
);
80 /* Enumerated value for indexing an uij entry (anisotropic temperature factors) */
82 U11
, U22
, U33
, U12
, U13
, U23
85 void pdb_use_ter(gmx_bool bSet
);
86 /* set read_pdbatoms to read upto 'TER' or 'ENDMDL' (default, bSet=FALSE).
87 This function is fundamentally broken as far as thread-safety is concerned.*/
89 void gmx_write_pdb_box(FILE *out
, int ePBC
, const matrix box
);
90 /* write the box in the CRYST1 record,
91 * with ePBC=-1 the pbc is guessed from the box
92 * This function is fundamentally broken as far as thread-safety is concerned.
95 void write_pdbfile_indexed(FILE *out
, const char *title
, const t_atoms
*atoms
,
96 const rvec x
[], int ePBC
, const matrix box
, char chain
,
97 int model_nr
, int nindex
, const int index
[],
98 gmx_conect conect
, gmx_bool bTerSepChains
);
99 /* REALLY low level */
101 void write_pdbfile(FILE *out
, const char *title
, const t_atoms
*atoms
,
102 const rvec x
[], int ePBC
, const matrix box
, char chain
,
103 int model_nr
, gmx_conect conect
, gmx_bool bTerSepChains
);
104 /* Low level pdb file writing routine.
106 * ONLY FOR SPECIAL PURPOSES,
108 * USE write_sto_conf WHEN YOU CAN.
110 * override chain-identifiers with chain when chain>0
111 * write ENDMDL when bEndmodel is TRUE.
113 * If the gmx_conect structure is not NULL its content is dumped as CONECT records
114 * which may be useful for visualization purposes.
117 void get_pdb_atomnumber(const t_atoms
*atoms
, struct gmx_atomprop
*aps
);
118 /* Routine to extract atomic numbers from the atom names */
120 int read_pdbfile(FILE *in
, char *title
, int *model_nr
,
121 struct t_atoms
*atoms
, struct t_symtab
*symtab
,
122 rvec x
[], int *ePBC
, matrix box
,
123 gmx_bool bChange
, gmx_conect conect
);
124 /* Function returns number of atoms found.
125 * ePBC and gmx_conect structure may be NULL.
128 void gmx_pdb_read_conf(const char *infile
,
129 t_symtab
*symtab
, char ***name
, t_atoms
*atoms
,
130 rvec x
[], int *ePBC
, matrix box
);
131 /* Read a pdb file and extract ATOM and HETATM fields.
132 * Read a box from the CRYST1 line, return 0 box when no CRYST1 is found.
136 void get_pdb_coordnum(FILE *in
, int *natoms
);
137 /* Read a pdb file and count the ATOM and HETATM fields. */
139 gmx_bool
is_hydrogen(const char *nm
);
140 /* Return whether atom nm is a hydrogen */
142 gmx_bool
is_dummymass(const char *nm
);
143 /* Return whether atom nm is a dummy mass */
145 /* Routines to handle CONECT records if they have been read in */
146 void gmx_conect_dump(FILE *fp
, gmx_conect conect
);
148 gmx_bool
gmx_conect_exist(gmx_conect conect
, int ai
, int aj
);
149 /* Return TRUE if there is a conection between the atoms */
151 void gmx_conect_add(gmx_conect conect
, int ai
, int aj
);
152 /* Add a connection between ai and aj (numbered from 0 to natom-1) */
154 gmx_conect
gmx_conect_generate(const t_topology
*top
);
155 /* Generate a conect structure from a topology */
157 gmx_conect
gmx_conect_init(void);
158 /* Initiate data structure */
160 void gmx_conect_done(gmx_conect gc
);
167 #endif /* GMX_FILEIO_PDBIO_H */