Add DSSP v4 support to do_dssp

[alexxy/gromacs.git] / src / gromacs / fileio / pdbio.cpp

diff --git a/src/gromacs/fileio/pdbio.cpp b/src/gromacs/fileio/pdbio.cpp
index 3de2fddbcd1b28b7346a751bbf882b95fa727924..c3c3995039774070b923690aa6330de820617a77 100644 (file)
--- a/src/gromacs/fileio/pdbio.cpp
+++ b/src/gromacs/fileio/pdbio.cpp
@@ -54,7 +54,6 @@
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/topology/atomprop.h"
 #include "gromacs/topology/ifunc.h"
-#include "gromacs/topology/residuetypes.h"
 #include "gromacs/topology/symtab.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/coolstuff.h"
@@ -87,8 +86,6 @@ const char* enumValueToString(PdbRecordType enumValue)
     return pdbRecordTypeName[enumValue];
 }
 
-#define REMARK_SIM_BOX "REMARK    THIS IS A SIMULATION BOX"
-
 void gmx_write_pdb_box(FILE* out, PbcType pbcType, const matrix box)
 {
     real alpha, beta, gamma;
@@ -299,7 +296,8 @@ void write_pdbfile_indexed(FILE*          out,
                            int            nindex,
                            const int      index[],
                            gmx_conect     conect,
-                           bool           usePqrFormat)
+                           bool           usePqrFormat,
+                           bool           standardCompliantMode)
 {
     gmx_conect_t* gc = static_cast<gmx_conect_t*>(conect);
     PdbRecordType type;
@@ -307,8 +305,27 @@ void write_pdbfile_indexed(FILE*          out,
     real          occup, bfac;
     gmx_bool      bOccup;
 
+    if (standardCompliantMode)
+    {
+        fprintf(out, "HEADER    GROMACS SIMULATION BOX                  01-JAN-00   0000\n");
+        fprintf(out,
+                "TITLE     Gromacs simulation box\n"
+                "COMPND    MOL_ID:  1;                                                           \n"
+                "COMPND   2 MOLECULE:  GROMACS SIMULATION BOX;                                   \n"
+                "COMPND   3 CHAIN: A;  \n"
+                "SOURCE    MOL_ID: 1;\n"
+                "SOURCE   2 SYNTHETIC\n"
+                "KEYWDS    GROMACS\n"
+                "EXPDTA    PURE PRODUCT OF COMPUTER SIMULATION\n"
+                "AUTHOR    GROMACS\n"
+                "REVDAT   1   01-JAN-00 0000    0\n");
+    }
+    else
+    {
+        fprintf(out, "TITLE     %s\n", (title && title[0]) ? title : gmx::bromacs().c_str());
+    }
+
 
-    fprintf(out, "TITLE     %s\n", (title && title[0]) ? title : gmx::bromacs().c_str());
     if (box && ((norm2(box[XX]) != 0.0F) || (norm2(box[YY]) != 0.0F) || (norm2(box[ZZ]) != 0.0F)))
     {
         gmx_write_pdb_box(out, pbcType, box);
@@ -332,7 +349,11 @@ void write_pdbfile_indexed(FILE*          out,
 
     fprintf(out, "MODEL %8d\n", model_nr > 0 ? model_nr : 1);
 
-    ResidueType rt;
+    // Collect last printed values for TER record
+    int         lastAtomNumber = 0;
+    std::string lastResName;
+    int         lastResNr = 0;
+
     for (int ii = 0; ii < nindex; ii++)
     {
         int         i      = index[ii];
@@ -395,6 +416,10 @@ void write_pdbfile_indexed(FILE*          out,
                                      bfac,
                                      atoms->atom[i].elem);
 
+            lastAtomNumber = i + 1;
+            lastResName    = resnm;
+            lastResNr      = resnr;
+
             if (atoms->pdbinfo && atoms->pdbinfo[i].bAnisotropic)
             {
                 fprintf(out,
@@ -430,7 +455,15 @@ void write_pdbfile_indexed(FILE*          out,
         }
     }
 
-    fprintf(out, "TER\n");
+    if (standardCompliantMode)
+    {
+        fprintf(out, "TER   %5d      %4.4s%c%4d\n", lastAtomNumber, lastResName.c_str(), chainid, lastResNr);
+    }
+    else
+    {
+        fprintf(out, "TER\n");
+    }
+
     fprintf(out, "ENDMDL\n");
 
     if (nullptr != gc)
@@ -575,18 +608,20 @@ void get_pdb_atomnumber(const t_atoms* atoms, AtomProperties* aps)
                 atomNumberSet = true;
             }
         }
-        std::string buf;
+        static constexpr size_t sc_maxElementNameLength = 3;
+        static_assert(sizeof(atoms->atom[i].elem) >= sc_maxElementNameLength + 1);
+        std::string element;
         if (atomNumberSet)
         {
             atoms->atom[i].atomnumber = atomnumber;
-            buf                       = aps->elementFromAtomNumber(atomnumber);
+            element                   = aps->elementFromAtomNumber(atomnumber);
             if (debug)
             {
                 fprintf(debug, "Atomnumber for atom '%s' is %d\n", anm, atomnumber);
             }
         }
-        buf.resize(3);
-        std::strncpy(atoms->atom[i].elem, buf.c_str(), 4);
+        element.resize(sc_maxElementNameLength);
+        std::strcpy(atoms->atom[i].elem, element.c_str());
     }
 }