Fix PQR file output

PQR files from editconf were always written as fixed format PDB files
with just the field information added. As pointed out in the linked
redmine, this can violate the PQR file standard if the field lengths are
too long, even though the file would still be a valid PDB.

This adds a slightly different form of the writeout that has a flexible,
PQR conform format.

Fixes #2511

Change-Id: I626380b642a0214970753da289e9c969ce411ea7

diff --git a/docs/release-notes/2018/2018.2.rst b/docs/release-notes/2018/2018.2.rst
index c168c76bfd1674b2f976a610fb8607c99f7d1df3..139212485dfafe3ed33b5c82849e859e7dbdd074 100644 (file)
--- a/docs/release-notes/2018/2018.2.rst
+++ b/docs/release-notes/2018/2018.2.rst
@@ -40,6 +40,15 @@ Fixed enemat when the .edr file had no matching energy groups
 
 :issue:`2508`
 
+Fix PQR file output
+""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+PQR files from gmx editconf violated the standard for the format because
+they were always written in fixed format. This commit fixes the issue by
+introducing a different output method for PQR files that follows the
+standard.
+
+:issue:`2511`
+
 Fixes to improve portability
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 

diff --git a/src/gromacs/fileio/confio.cpp b/src/gromacs/fileio/confio.cpp
index b9015383b4d3410a56b6184d0cdbf57e8b56d8e6..cfb8d0fac3f0ac5aa2b9887e5f4129b4e99fe320 100644 (file)
--- a/src/gromacs/fileio/confio.cpp
+++ b/src/gromacs/fileio/confio.cpp
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -101,7 +101,7 @@ void write_sto_conf_indexed(const char *outfile, const char *title,
         case efENT:
         case efPQR:
             out = gmx_fio_fopen(outfile, "w");
-            write_pdbfile_indexed(out, title, atoms, x, ePBC, box, ' ', -1, nindex, index, nullptr, TRUE);
+            write_pdbfile_indexed(out, title, atoms, x, ePBC, box, ' ', -1, nindex, index, nullptr, TRUE, TRUE);
             gmx_fio_fclose(out);
             break;
         case efESP:

diff --git a/src/gromacs/fileio/pdbio.cpp b/src/gromacs/fileio/pdbio.cpp
index d6c7579199759b7517429e5fca7df086e2a1d568..da9893d6b4f7a2a6d38451c02d6f02caa718c014 100644 (file)
--- a/src/gromacs/fileio/pdbio.cpp
+++ b/src/gromacs/fileio/pdbio.cpp
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -259,11 +259,56 @@ static void read_cryst1(char *line, int *ePBC, matrix box)
     }
 }
 
+static int
+gmx_fprintf_pqr_atomline(FILE *            fp,
+                         enum PDB_record   record,
+                         int               atom_seq_number,
+                         const char *      atom_name,
+                         const char *      res_name,
+                         char              chain_id,
+                         int               res_seq_number,
+                         real              x,
+                         real              y,
+                         real              z,
+                         real              occupancy,
+                         real              b_factor)
+{
+    GMX_RELEASE_ASSERT(record == epdbATOM || record == epdbHETATM,
+                       "Can only print PQR atom lines as ATOM or HETATM records");
+
+    /* Check atom name */
+    GMX_RELEASE_ASSERT(atom_name != nullptr,
+                       "Need atom information to print pqr");
+
+    /* Check residue name */
+    GMX_RELEASE_ASSERT(res_name != nullptr,
+                       "Need residue information to print pqr");
+
+    /* Truncate integers so they fit */
+    atom_seq_number = atom_seq_number % 100000;
+    res_seq_number  = res_seq_number % 10000;
+
+    int n = fprintf(fp,
+                    "%s %d %s %s %c %d %8.3f %8.3f %8.3f %6.2f %6.2f\n",
+                    pdbtp[record],
+                    atom_seq_number,
+                    atom_name,
+                    res_name,
+                    chain_id,
+                    res_seq_number,
+                    x, y, z,
+                    occupancy,
+                    b_factor);
+
+    return n;
+}
+
 void write_pdbfile_indexed(FILE *out, const char *title,
                            const t_atoms *atoms, const rvec x[],
                            int ePBC, const matrix box, char chainid,
                            int model_nr, int nindex, const int index[],
-                           gmx_conect conect, gmx_bool bTerSepChains)
+                           gmx_conect conect, gmx_bool bTerSepChains,
+                           bool usePqrFormat)
 {
     gmx_conect_t     *gc = (gmx_conect_t *)conect;
     char              resnm[6], nm[6];
@@ -370,29 +415,44 @@ void write_pdbfile_indexed(FILE *out, const char *title,
         }
         occup = bOccup ? 1.0 : pdbinfo.occup;
         bfac  = pdbinfo.bfac;
-
-        gmx_fprintf_pdb_atomline(out,
-                                 type,
-                                 i+1,
-                                 nm,
-                                 altloc,
-                                 resnm,
-                                 ch,
-                                 resnr,
-                                 resic,
-                                 10*x[i][XX], 10*x[i][YY], 10*x[i][ZZ],
-                                 occup,
-                                 bfac,
-                                 atoms->atom[i].elem);
-
-        if (atoms->pdbinfo && atoms->pdbinfo[i].bAnisotropic)
+        if (!usePqrFormat)
+        {
+            gmx_fprintf_pdb_atomline(out,
+                                     type,
+                                     i+1,
+                                     nm,
+                                     altloc,
+                                     resnm,
+                                     ch,
+                                     resnr,
+                                     resic,
+                                     10*x[i][XX], 10*x[i][YY], 10*x[i][ZZ],
+                                     occup,
+                                     bfac,
+                                     atoms->atom[i].elem);
+
+            if (atoms->pdbinfo && atoms->pdbinfo[i].bAnisotropic)
+            {
+                fprintf(out, "ANISOU%5d  %-4.4s%4.4s%c%4d%c %7d%7d%7d%7d%7d%7d\n",
+                        (i+1)%100000, nm, resnm, ch, resnr,
+                        (resic == '\0') ? ' ' : resic,
+                        atoms->pdbinfo[i].uij[0], atoms->pdbinfo[i].uij[1],
+                        atoms->pdbinfo[i].uij[2], atoms->pdbinfo[i].uij[3],
+                        atoms->pdbinfo[i].uij[4], atoms->pdbinfo[i].uij[5]);
+            }
+        }
+        else
         {
-            fprintf(out, "ANISOU%5d  %-4.4s%4.4s%c%4d%c %7d%7d%7d%7d%7d%7d\n",
-                    (i+1)%100000, nm, resnm, ch, resnr,
-                    (resic == '\0') ? ' ' : resic,
-                    atoms->pdbinfo[i].uij[0], atoms->pdbinfo[i].uij[1],
-                    atoms->pdbinfo[i].uij[2], atoms->pdbinfo[i].uij[3],
-                    atoms->pdbinfo[i].uij[4], atoms->pdbinfo[i].uij[5]);
+            gmx_fprintf_pqr_atomline(out,
+                                     type,
+                                     i+1,
+                                     nm,
+                                     resnm,
+                                     ch,
+                                     resnr,
+                                     10*x[i][XX], 10*x[i][YY], 10*x[i][ZZ],
+                                     occup,
+                                     bfac);
         }
     }
 
@@ -422,7 +482,7 @@ void write_pdbfile(FILE *out, const char *title, const t_atoms *atoms, const rve
         index[i] = i;
     }
     write_pdbfile_indexed(out, title, atoms, x, ePBC, box, chainid, model_nr,
-                          atoms->nr, index, conect, bTerSepChains);
+                          atoms->nr, index, conect, bTerSepChains, false);
     sfree(index);
 }
 

diff --git a/src/gromacs/fileio/pdbio.h b/src/gromacs/fileio/pdbio.h
index 27d5106d53b1faec01042e8392240b4e4a1cd886..c88a5c2834cb5ce45f8c94f9857258fff277d588 100644 (file)
--- a/src/gromacs/fileio/pdbio.h
+++ b/src/gromacs/fileio/pdbio.h
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -95,7 +95,8 @@ void gmx_write_pdb_box(FILE *out, int ePBC, const matrix box);
 void write_pdbfile_indexed(FILE *out, const char *title, const t_atoms *atoms,
                            const rvec x[], int ePBC, const matrix box, char chain,
                            int model_nr, int nindex, const int index[],
-                           gmx_conect conect, gmx_bool bTerSepChains);
+                           gmx_conect conect, gmx_bool bTerSepChains,
+                           bool usePqrFormat);
 /* REALLY low level */
 
 void write_pdbfile(FILE *out, const char *title, const t_atoms *atoms,

diff --git a/src/gromacs/fileio/trxio.cpp b/src/gromacs/fileio/trxio.cpp
index a9932699c4a49fbc5e2bb5df409408b34a3a02bb..2ed5c41cc1da40227882008a0b326e367f22bb14 100644 (file)
--- a/src/gromacs/fileio/trxio.cpp
+++ b/src/gromacs/fileio/trxio.cpp
@@ -427,7 +427,7 @@ int write_trxframe_indexed(t_trxstatus *status, const t_trxframe *fr, int nind,
             else
             {
                 write_pdbfile_indexed(gmx_fio_getfp(status->fio), title, fr->atoms,
-                                      fr->x, -1, fr->box, ' ', fr->step, nind, ind, gc, TRUE);
+                                      fr->x, -1, fr->box, ' ', fr->step, nind, ind, gc, TRUE, FALSE);
             }
             break;
         case efG96:

diff --git a/src/gromacs/gmxana/gmx_do_dssp.cpp b/src/gromacs/gmxana/gmx_do_dssp.cpp
index a8aab1bc21c1734db7b7b0a47c0d119ffbc53d00..039e96388168dd3c75c787a40aa2398f02c1ead3 100644 (file)
--- a/src/gromacs/gmxana/gmx_do_dssp.cpp
+++ b/src/gromacs/gmxana/gmx_do_dssp.cpp
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2012,2013,2014,2015,2017, by the GROMACS development team, led by
+ * Copyright (c) 2012,2013,2014,2015,2017,2018, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -660,7 +660,7 @@ int gmx_do_dssp(int argc, char *argv[])
         }
         gmx_rmpbc(gpbc, natoms, box, x);
         tapein = gmx_ffopen(pdbfile, "w");
-        write_pdbfile_indexed(tapein, nullptr, atoms, x, ePBC, box, ' ', -1, gnx, index, nullptr, TRUE);
+        write_pdbfile_indexed(tapein, nullptr, atoms, x, ePBC, box, ' ', -1, gnx, index, nullptr, TRUE, FALSE);
         gmx_ffclose(tapein);
 
         if (0 != system(dssp))

diff --git a/src/gromacs/gmxana/gmx_editconf.cpp b/src/gromacs/gmxana/gmx_editconf.cpp
index 7abb509a260f6f31239802534e3306fa322b7119..f26bd749cb7346bb0e53de011b1c314effed5efe 100644 (file)
--- a/src/gromacs/gmxana/gmx_editconf.cpp
+++ b/src/gromacs/gmxana/gmx_editconf.cpp
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -1265,7 +1265,7 @@ int gmx_editconf(int argc, char *argv[])
         if (outftp == efPDB)
         {
             out = gmx_ffopen(outfile, "w");
-            write_pdbfile_indexed(out, *top_tmp->name, &atoms, x, ePBC, box, ' ', 1, isize, index, conect, TRUE);
+            write_pdbfile_indexed(out, *top_tmp->name, &atoms, x, ePBC, box, ' ', 1, isize, index, conect, TRUE, FALSE);
             gmx_ffclose(out);
         }
         else
@@ -1311,7 +1311,18 @@ int gmx_editconf(int argc, char *argv[])
                     atoms.resinfo[atoms.atom[i].resind].chainid = label[0];
                 }
             }
-            write_pdbfile(out, *top_tmp->name, &atoms, x, ePBC, box, ' ', -1, conect, TRUE);
+            /* Need to bypass the regular write_pdbfile because I don't want to change
+             * all instances to include the boolean flag for writing out PQR files.
+             */
+            int *index;
+            snew(index, atoms.nr);
+            for (int i = 0; i < atoms.nr; i++)
+            {
+                index[i] = i;
+            }
+            write_pdbfile_indexed(out, *top_tmp->name, &atoms, x, ePBC, box, ' ', -1, atoms.nr, index, conect,
+                                  TRUE, outftp == efPQR ? true : false);
+            sfree(index);
             if (bLegend)
             {
                 pdb_legend(out, atoms.nr, atoms.nres, &atoms, x);