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.
48 #include "gromacs/fileio/gmxfio.h"
49 #include "gromacs/math/units.h"
50 #include "gromacs/math/vec.h"
51 #include "gromacs/pbcutil/pbc.h"
52 #include "gromacs/topology/atomprop.h"
53 #include "gromacs/topology/ifunc.h"
54 #include "gromacs/topology/residuetypes.h"
55 #include "gromacs/topology/symtab.h"
56 #include "gromacs/topology/topology.h"
57 #include "gromacs/utility/coolstuff.h"
58 #include "gromacs/utility/cstringutil.h"
59 #include "gromacs/utility/fatalerror.h"
60 #include "gromacs/utility/futil.h"
61 #include "gromacs/utility/smalloc.h"
67 typedef struct gmx_conect_t
{
70 gmx_conection_t
*conect
;
73 static const char *pdbtp
[epdbNR
] = {
74 "ATOM ", "HETATM", "ANISOU", "CRYST1",
75 "COMPND", "MODEL", "ENDMDL", "TER", "HEADER", "TITLE", "REMARK",
80 #define REMARK_SIM_BOX "REMARK THIS IS A SIMULATION BOX"
82 static void xlate_atomname_pdb2gmx(char *name
)
87 length
= std::strlen(name
);
88 if (length
> 3 && std::isdigit(name
[0]))
91 for (i
= 1; i
< length
; i
++)
95 name
[length
-1] = temp
;
99 static void xlate_atomname_gmx2pdb(char *name
)
104 length
= std::strlen(name
);
105 if (length
> 3 && std::isdigit(name
[length
-1]))
107 temp
= name
[length
-1];
108 for (i
= length
-1; i
> 0; --i
)
117 void gmx_write_pdb_box(FILE *out
, int ePBC
, matrix box
)
119 real alpha
, beta
, gamma
;
123 ePBC
= guess_ePBC(box
);
126 if (ePBC
== epbcNONE
)
131 if (norm2(box
[YY
])*norm2(box
[ZZ
]) != 0)
133 alpha
= RAD2DEG
*std::acos(cos_angle_no_table(box
[YY
], box
[ZZ
]));
139 if (norm2(box
[XX
])*norm2(box
[ZZ
]) != 0)
141 beta
= RAD2DEG
*std::acos(cos_angle_no_table(box
[XX
], box
[ZZ
]));
147 if (norm2(box
[XX
])*norm2(box
[YY
]) != 0)
149 gamma
= RAD2DEG
*std::acos(cos_angle_no_table(box
[XX
], box
[YY
]));
155 fprintf(out
, "REMARK THIS IS A SIMULATION BOX\n");
156 if (ePBC
!= epbcSCREW
)
158 fprintf(out
, "CRYST1%9.3f%9.3f%9.3f%7.2f%7.2f%7.2f %-11s%4d\n",
159 10*norm(box
[XX
]), 10*norm(box
[YY
]), 10*norm(box
[ZZ
]),
160 alpha
, beta
, gamma
, "P 1", 1);
164 /* Double the a-vector length and write the correct space group */
165 fprintf(out
, "CRYST1%9.3f%9.3f%9.3f%7.2f%7.2f%7.2f %-11s%4d\n",
166 20*norm(box
[XX
]), 10*norm(box
[YY
]), 10*norm(box
[ZZ
]),
167 alpha
, beta
, gamma
, "P 21 1 1", 1);
172 static void read_cryst1(char *line
, int *ePBC
, matrix box
)
175 char sa
[12], sb
[12], sc
[12], sg
[SG_SIZE
+1], ident
;
176 double fa
, fb
, fc
, alpha
, beta
, gamma
, cosa
, cosb
, cosg
, sing
;
177 int syma
, symb
, symc
;
180 sscanf(line
, "%*s%s%s%s%lf%lf%lf", sa
, sb
, sc
, &alpha
, &beta
, &gamma
);
183 if (strlen(line
) >= 55)
185 strncpy(sg
, line
+55, SG_SIZE
);
191 sscanf(sg
, "%c %d %d %d", &ident
, &syma
, &symb
, &symc
);
192 if (ident
== 'P' && syma
== 1 && symb
<= 1 && symc
<= 1)
194 fc
= strtod(sc
, NULL
)*0.1;
195 ePBC_file
= (fc
> 0 ? epbcXYZ
: epbcXY
);
197 if (ident
== 'P' && syma
== 21 && symb
== 1 && symc
== 1)
199 ePBC_file
= epbcSCREW
;
209 fa
= strtod(sa
, NULL
)*0.1;
210 fb
= strtod(sb
, NULL
)*0.1;
211 fc
= strtod(sc
, NULL
)*0.1;
212 if (ePBC_file
== epbcSCREW
)
218 if ((alpha
!= 90.0) || (beta
!= 90.0) || (gamma
!= 90.0))
222 cosa
= std::cos(alpha
*DEG2RAD
);
230 cosb
= std::cos(beta
*DEG2RAD
);
238 cosg
= std::cos(gamma
*DEG2RAD
);
239 sing
= std::sin(gamma
*DEG2RAD
);
246 box
[YY
][XX
] = fb
*cosg
;
247 box
[YY
][YY
] = fb
*sing
;
248 box
[ZZ
][XX
] = fc
*cosb
;
249 box
[ZZ
][YY
] = fc
*(cosa
- cosb
*cosg
)/sing
;
250 box
[ZZ
][ZZ
] = std::sqrt(fc
*fc
251 - box
[ZZ
][XX
]*box
[ZZ
][XX
] - box
[ZZ
][YY
]*box
[ZZ
][YY
]);
261 void write_pdbfile_indexed(FILE *out
, const char *title
,
262 t_atoms
*atoms
, rvec x
[],
263 int ePBC
, matrix box
, char chainid
,
264 int model_nr
, int nindex
, const int index
[],
265 gmx_conect conect
, gmx_bool bTerSepChains
)
267 gmx_conect_t
*gc
= (gmx_conect_t
*)conect
;
268 char resnm
[6], nm
[6];
271 enum PDB_record type
;
272 unsigned char resic
, ch
;
276 int chainnum
, lastchainnum
;
277 gmx_residuetype_t
*rt
;
278 const char *p_restype
;
279 const char *p_lastrestype
;
281 gmx_residuetype_init(&rt
);
283 fprintf(out
, "TITLE %s\n", (title
&& title
[0]) ? title
: gmx::bromacs().c_str());
284 if (box
&& ( norm2(box
[XX
]) || norm2(box
[YY
]) || norm2(box
[ZZ
]) ) )
286 gmx_write_pdb_box(out
, ePBC
, box
);
290 /* Check whether any occupancies are set, in that case leave it as is,
291 * otherwise set them all to one
294 for (ii
= 0; (ii
< nindex
) && bOccup
; ii
++)
297 bOccup
= bOccup
&& (atoms
->pdbinfo
[i
].occup
== 0.0);
305 fprintf(out
, "MODEL %8d\n", model_nr
> 0 ? model_nr
: 1);
310 for (ii
= 0; ii
< nindex
; ii
++)
313 resind
= atoms
->atom
[i
].resind
;
314 chainnum
= atoms
->resinfo
[resind
].chainnum
;
315 p_lastrestype
= p_restype
;
316 gmx_residuetype_get_type(rt
, *atoms
->resinfo
[resind
].name
, &p_restype
);
318 /* Add a TER record if we changed chain, and if either the previous or this chain is protein/DNA/RNA. */
319 if (bTerSepChains
&& ii
> 0 && chainnum
!= lastchainnum
)
321 /* Only add TER if the previous chain contained protein/DNA/RNA. */
322 if (gmx_residuetype_is_protein(rt
, p_lastrestype
) || gmx_residuetype_is_dna(rt
, p_lastrestype
) || gmx_residuetype_is_rna(rt
, p_lastrestype
))
324 fprintf(out
, "TER\n");
326 lastchainnum
= chainnum
;
329 strncpy(resnm
, *atoms
->resinfo
[resind
].name
, sizeof(resnm
)-1);
330 resnm
[sizeof(resnm
)-1] = 0;
331 strncpy(nm
, *atoms
->atomname
[i
], sizeof(nm
)-1);
332 nm
[sizeof(nm
)-1] = 0;
334 /* rename HG12 to 2HG1, etc. */
335 xlate_atomname_gmx2pdb(nm
);
336 resnr
= atoms
->resinfo
[resind
].nr
;
337 resic
= atoms
->resinfo
[resind
].ic
;
344 ch
= atoms
->resinfo
[resind
].chainid
;
353 resnr
= resnr
% 10000;
357 type
= static_cast<enum PDB_record
>(atoms
->pdbinfo
[i
].type
);
358 altloc
= atoms
->pdbinfo
[i
].altloc
;
359 if (!isalnum(altloc
))
363 occup
= bOccup
? 1.0 : atoms
->pdbinfo
[i
].occup
;
364 bfac
= atoms
->pdbinfo
[i
].bfac
;
374 gmx_fprintf_pdb_atomline(out
,
383 10*x
[i
][XX
], 10*x
[i
][YY
], 10*x
[i
][ZZ
],
386 atoms
->atom
[i
].elem
);
388 if (atoms
->pdbinfo
&& atoms
->pdbinfo
[i
].bAnisotropic
)
390 fprintf(out
, "ANISOU%5d %-4.4s%4.4s%c%4d%c %7d%7d%7d%7d%7d%7d\n",
391 (i
+1)%100000, nm
, resnm
, ch
, resnr
,
392 (resic
== '\0') ? ' ' : resic
,
393 atoms
->pdbinfo
[i
].uij
[0], atoms
->pdbinfo
[i
].uij
[1],
394 atoms
->pdbinfo
[i
].uij
[2], atoms
->pdbinfo
[i
].uij
[3],
395 atoms
->pdbinfo
[i
].uij
[4], atoms
->pdbinfo
[i
].uij
[5]);
399 fprintf(out
, "TER\n");
400 fprintf(out
, "ENDMDL\n");
404 /* Write conect records */
405 for (i
= 0; (i
< gc
->nconect
); i
++)
407 fprintf(out
, "CONECT%5d%5d\n", gc
->conect
[i
].ai
+1, gc
->conect
[i
].aj
+1);
411 gmx_residuetype_destroy(rt
);
414 void write_pdbfile(FILE *out
, const char *title
, t_atoms
*atoms
, rvec x
[],
415 int ePBC
, matrix box
, char chainid
, int model_nr
, gmx_conect conect
, gmx_bool bTerSepChains
)
419 snew(index
, atoms
->nr
);
420 for (i
= 0; i
< atoms
->nr
; i
++)
424 write_pdbfile_indexed(out
, title
, atoms
, x
, ePBC
, box
, chainid
, model_nr
,
425 atoms
->nr
, index
, conect
, bTerSepChains
);
429 int line2type(char *line
)
434 for (k
= 0; (k
< 6); k
++)
440 for (k
= 0; (k
< epdbNR
); k
++)
442 if (std::strncmp(type
, pdbtp
[k
], strlen(pdbtp
[k
])) == 0)
451 static void read_anisou(char line
[], int natom
, t_atoms
*atoms
)
455 char anr
[12], anm
[12];
459 for (k
= 0; (k
< 5); k
++, j
++)
465 for (k
= 0; (k
< 4); k
++, j
++)
472 /* Strip off spaces */
475 /* Search backwards for number and name only */
476 atomnr
= std::strtol(anr
, NULL
, 10);
477 for (i
= natom
-1; (i
>= 0); i
--)
479 if ((std::strcmp(anm
, *(atoms
->atomname
[i
])) == 0) &&
480 (atomnr
== atoms
->pdbinfo
[i
].atomnr
))
487 fprintf(stderr
, "Skipping ANISOU record (atom %s %d not found)\n",
492 if (sscanf(line
+29, "%d%d%d%d%d%d",
493 &atoms
->pdbinfo
[i
].uij
[U11
], &atoms
->pdbinfo
[i
].uij
[U22
],
494 &atoms
->pdbinfo
[i
].uij
[U33
], &atoms
->pdbinfo
[i
].uij
[U12
],
495 &atoms
->pdbinfo
[i
].uij
[U13
], &atoms
->pdbinfo
[i
].uij
[U23
])
498 atoms
->pdbinfo
[i
].bAnisotropic
= TRUE
;
502 fprintf(stderr
, "Invalid ANISOU record for atom %d\n", i
);
503 atoms
->pdbinfo
[i
].bAnisotropic
= FALSE
;
508 void get_pdb_atomnumber(t_atoms
*atoms
, gmx_atomprop_t aps
)
510 int i
, atomnumber
, len
;
512 char anm
[6], anm_copy
[6], *ptr
;
518 gmx_incons("Trying to deduce atomnumbers when no pdb information is present");
520 for (i
= 0; (i
< atoms
->nr
); i
++)
522 std::strcpy(anm
, atoms
->pdbinfo
[i
].atomnm
);
523 std::strcpy(anm_copy
, atoms
->pdbinfo
[i
].atomnm
);
524 bool atomNumberSet
= false;
526 if ((anm
[0] != ' ') && ((len
<= 2) || ((len
> 2) && !std::isdigit(anm
[2]))))
529 if (gmx_atomprop_query(aps
, epropElement
, "???", anm_copy
, &eval
))
531 atomnumber
= std::round(eval
);
532 atomNumberSet
= true;
537 if (gmx_atomprop_query(aps
, epropElement
, "???", anm_copy
, &eval
))
539 atomnumber
= std::round(eval
);
540 atomNumberSet
= true;
547 while ((k
< std::strlen(anm
)) && (std::isspace(anm
[k
]) || std::isdigit(anm
[k
])))
551 anm_copy
[0] = anm
[k
];
553 if (gmx_atomprop_query(aps
, epropElement
, "???", anm_copy
, &eval
))
555 atomnumber
= std::round(eval
);
556 atomNumberSet
= true;
561 atoms
->atom
[i
].atomnumber
= atomnumber
;
562 ptr
= gmx_atomprop_element(aps
, atomnumber
);
565 fprintf(debug
, "Atomnumber for atom '%s' is %d\n",
573 std::strncpy(atoms
->atom
[i
].elem
, ptr
== NULL
? "" : ptr
, 4);
577 static int read_atom(t_symtab
*symtab
,
578 char line
[], int type
, int natom
,
579 t_atoms
*atoms
, rvec x
[], int chainnum
, gmx_bool bChange
)
584 char anr
[12], anm
[12], anm_copy
[12], altloc
, resnm
[12], rnr
[12], elem
[3];
585 char xc
[12], yc
[12], zc
[12], occup
[12], bfac
[12];
588 int resnr
, atomnumber
;
590 if (natom
>= atoms
->nr
)
592 gmx_fatal(FARGS
, "\nFound more atoms (%d) in pdb file than expected (%d)",
598 for (k
= 0; (k
< 5); k
++, j
++)
605 for (k
= 0; (k
< 4); k
++, j
++)
610 std::strcpy(anm_copy
, anm
);
616 for (k
= 0; (k
< 4); k
++, j
++)
626 for (k
= 0; (k
< 4); k
++, j
++)
632 resnr
= std::strtol(rnr
, NULL
, 10);
636 /* X,Y,Z Coordinate */
637 for (k
= 0; (k
< 8); k
++, j
++)
642 for (k
= 0; (k
< 8); k
++, j
++)
647 for (k
= 0; (k
< 8); k
++, j
++)
654 for (k
= 0; (k
< 6); k
++, j
++)
661 for (k
= 0; (k
< 7); k
++, j
++)
671 for (k
= 0; (k
< 2); k
++, j
++)
680 atomn
= &(atoms
->atom
[natom
]);
682 atoms
->resinfo
[atoms
->atom
[natom
-1].resind
].nr
!= resnr
||
683 atoms
->resinfo
[atoms
->atom
[natom
-1].resind
].ic
!= resic
||
684 (strcmp(*atoms
->resinfo
[atoms
->atom
[natom
-1].resind
].name
, resnm
) != 0))
692 atomn
->resind
= atoms
->atom
[natom
-1].resind
+ 1;
694 atoms
->nres
= atomn
->resind
+ 1;
695 t_atoms_set_resinfo(atoms
, natom
, symtab
, resnm
, resnr
, resic
, chainnum
, chainid
);
699 atomn
->resind
= atoms
->atom
[natom
-1].resind
;
703 xlate_atomname_pdb2gmx(anm
);
705 atoms
->atomname
[natom
] = put_symtab(symtab
, anm
);
708 atomn
->atomnumber
= atomnumber
;
709 strncpy(atomn
->elem
, elem
, 4);
711 x
[natom
][XX
] = strtod(xc
, NULL
)*0.1;
712 x
[natom
][YY
] = strtod(yc
, NULL
)*0.1;
713 x
[natom
][ZZ
] = strtod(zc
, NULL
)*0.1;
716 atoms
->pdbinfo
[natom
].type
= type
;
717 atoms
->pdbinfo
[natom
].atomnr
= strtol(anr
, NULL
, 10);
718 atoms
->pdbinfo
[natom
].altloc
= altloc
;
719 strcpy(atoms
->pdbinfo
[natom
].atomnm
, anm_copy
);
720 atoms
->pdbinfo
[natom
].bfac
= strtod(bfac
, NULL
);
721 atoms
->pdbinfo
[natom
].occup
= strtod(occup
, NULL
);
728 gmx_bool
is_hydrogen(const char *nm
)
732 std::strcpy(buf
, nm
);
739 else if ((std::isdigit(buf
[0])) && (buf
[1] == 'H'))
746 gmx_bool
is_dummymass(const char *nm
)
750 std::strcpy(buf
, nm
);
753 if ((buf
[0] == 'M') && std::isdigit(buf
[strlen(buf
)-1]))
761 static void gmx_conect_addline(gmx_conect_t
*con
, char *line
)
764 char format
[32], form2
[32];
766 sprintf(form2
, "%%*s");
767 sprintf(format
, "%s%%d", form2
);
768 if (sscanf(line
, format
, &ai
) == 1)
772 std::strcat(form2
, "%*s");
773 sprintf(format
, "%s%%d", form2
);
774 n
= sscanf(line
, format
, &aj
);
777 srenew(con
->conect
, ++con
->nconect
);
778 con
->conect
[con
->nconect
-1].ai
= ai
-1;
779 con
->conect
[con
->nconect
-1].aj
= aj
-1;
786 void gmx_conect_dump(FILE *fp
, gmx_conect conect
)
788 gmx_conect_t
*gc
= (gmx_conect_t
*)conect
;
791 for (i
= 0; (i
< gc
->nconect
); i
++)
793 fprintf(fp
, "%6s%5d%5d\n", "CONECT",
794 gc
->conect
[i
].ai
+1, gc
->conect
[i
].aj
+1);
798 gmx_conect
gmx_conect_init()
807 void gmx_conect_done(gmx_conect conect
)
809 gmx_conect_t
*gc
= conect
;
814 gmx_bool
gmx_conect_exist(gmx_conect conect
, int ai
, int aj
)
816 gmx_conect_t
*gc
= conect
;
822 for (i
= 0; (i
< gc
->nconect
); i
++)
824 if (((gc
->conect
[i
].ai
== ai
) &&
825 (gc
->conect
[i
].aj
== aj
)) ||
826 ((gc
->conect
[i
].aj
== ai
) &&
827 (gc
->conect
[i
].ai
== aj
)))
835 void gmx_conect_add(gmx_conect conect
, int ai
, int aj
)
837 gmx_conect_t
*gc
= (gmx_conect_t
*)conect
;
842 if (!gmx_conect_exist(conect
, ai
, aj
))
844 srenew(gc
->conect
, ++gc
->nconect
);
845 gc
->conect
[gc
->nconect
-1].ai
= ai
;
846 gc
->conect
[gc
->nconect
-1].aj
= aj
;
850 int read_pdbfile(FILE *in
, char *title
, int *model_nr
,
851 t_atoms
*atoms
, t_symtab
*symtab
, rvec x
[], int *ePBC
,
852 matrix box
, gmx_bool bChange
, gmx_conect conect
)
854 gmx_conect_t
*gc
= conect
;
856 gmx_bool bConnWarn
= FALSE
;
861 gmx_bool bStop
= FALSE
;
865 /* Only assume pbc when there is a CRYST1 entry */
877 while (!bStop
&& (fgets2(line
, STRLEN
, in
) != NULL
))
879 line_type
= line2type(line
);
885 natom
= read_atom(symtab
, line
, line_type
, natom
, atoms
, x
, chainnum
, bChange
);
891 read_anisou(line
, natom
, atoms
);
896 read_cryst1(line
, ePBC
, box
);
901 if (std::strlen(line
) > 6)
904 /* skip HEADER or TITLE and spaces */
913 /* truncate after title */
914 d
= std::strstr(c
, " ");
919 if (std::strlen(c
) > 0)
921 std::strcpy(title
, c
);
927 if ((!std::strstr(line
, ": ")) || (std::strstr(line
+6, "MOLECULE:")))
929 if (!(c
= std::strstr(line
+6, "MOLECULE:")) )
933 /* skip 'MOLECULE:' and spaces */
942 /* truncate after title */
946 while ( (d
[-1] == ';') && d
> c
)
956 std::strcat(title
, "; ");
957 std::strcat(title
, c
);
961 std::strcpy(title
, c
);
975 sscanf(line
, "%*s%d", model_nr
);
985 gmx_conect_addline(gc
, line
);
989 fprintf(stderr
, "WARNING: all CONECT records are ignored\n");
1002 void get_pdb_coordnum(FILE *in
, int *natoms
)
1007 while (fgets2(line
, STRLEN
, in
))
1009 if (std::strncmp(line
, "ENDMDL", 6) == 0)
1013 if ((std::strncmp(line
, "ATOM ", 6) == 0) || (std::strncmp(line
, "HETATM", 6) == 0))
1020 void gmx_pdb_read_conf(const char *infile
,
1021 t_topology
*top
, rvec x
[], int *ePBC
, matrix box
)
1023 FILE *in
= gmx_fio_fopen(infile
, "r");
1025 read_pdbfile(in
, title
, NULL
, &top
->atoms
, &top
->symtab
, x
, ePBC
, box
, TRUE
, NULL
);
1026 top
->name
= put_symtab(&top
->symtab
, title
);
1030 gmx_conect
gmx_conect_generate(t_topology
*top
)
1035 /* Fill the conect records */
1036 gc
= gmx_conect_init();
1038 for (f
= 0; (f
< F_NRE
); f
++)
1042 for (i
= 0; (i
< top
->idef
.il
[f
].nr
); i
+= interaction_function
[f
].nratoms
+1)
1044 gmx_conect_add(gc
, top
->idef
.il
[f
].iatoms
[i
+1],
1045 top
->idef
.il
[f
].iatoms
[i
+2]);
1053 gmx_fprintf_pdb_atomline(FILE * fp
,
1054 enum PDB_record record
,
1055 int atom_seq_number
,
1056 const char * atom_name
,
1057 char alternate_location
,
1058 const char * res_name
,
1061 char res_insertion_code
,
1067 const char * element
)
1069 char tmp_atomname
[6], tmp_resname
[6];
1070 gmx_bool start_name_in_col13
;
1073 if (record
!= epdbATOM
&& record
!= epdbHETATM
)
1075 gmx_fatal(FARGS
, "Can only print PDB atom lines as ATOM or HETATM records");
1078 /* Format atom name */
1079 if (atom_name
!= NULL
)
1081 /* If the atom name is an element name with two chars, it should start already in column 13.
1082 * Otherwise it should start in column 14, unless the name length is 4 chars.
1084 if ( (element
!= NULL
) && (std::strlen(element
) >= 2) && (gmx_strncasecmp(atom_name
, element
, 2) == 0) )
1086 start_name_in_col13
= TRUE
;
1090 start_name_in_col13
= (std::strlen(atom_name
) >= 4);
1092 sprintf(tmp_atomname
, start_name_in_col13
? "" : " ");
1093 std::strncat(tmp_atomname
, atom_name
, 4);
1094 tmp_atomname
[5] = '\0';
1098 tmp_atomname
[0] = '\0';
1101 /* Format residue name */
1102 std::strncpy(tmp_resname
, (res_name
!= NULL
) ? res_name
: "", 4);
1103 /* Make sure the string is terminated if strlen was > 4 */
1104 tmp_resname
[4] = '\0';
1105 /* String is properly terminated, so now we can use strcat. By adding a
1106 * space we can write it right-justified, and if the original name was
1107 * three characters or less there will be a space added on the right side.
1109 std::strcat(tmp_resname
, " ");
1111 /* Truncate integers so they fit */
1112 atom_seq_number
= atom_seq_number
% 100000;
1113 res_seq_number
= res_seq_number
% 10000;
1116 "%-6s%5d %-4.4s%c%4.4s%c%4d%c %8.3f%8.3f%8.3f%6.2f%6.2f %2s\n",
1128 (element
!= NULL
) ? element
: "");