Merge branch 'master' into pygromacs
[alexxy/gromacs.git] / src / gromacs / gmxana / gmx_nmtraj.cpp
similarity index 90%
rename from src/gromacs/gmxana/gmx_nmtraj.c
rename to src/gromacs/gmxana/gmx_nmtraj.cpp
index 7e98f3de6d82b546f2774e30873d6caaf1f2ca62..f319773d8445ae416c17cfd35659152dbf7d9188 100644 (file)
  */
 #include "gmxpre.h"
 
-#include <math.h>
-#include <stdlib.h>
-#include <string.h>
+#include <cctype>
+#include <cmath>
+#include <cstdlib>
+#include <cstring>
 
 #include "gromacs/commandline/pargs.h"
 #include "gromacs/fileio/confio.h"
@@ -55,7 +56,7 @@
 #include "gromacs/utility/fatalerror.h"
 #include "gromacs/utility/futil.h"
 #include "gromacs/utility/smalloc.h"
-
+#include "gromacs/utility/stringutil.h"
 
 int gmx_nmtraj(int argc, char *argv[])
 {
@@ -100,31 +101,24 @@ int gmx_nmtraj(int argc, char *argv[])
     t_atoms          *atoms;
     rvec             *xtop, *xref, *xav, *xout;
     int               nvec, *eignr = NULL;
-    int              *eigvalnr;
     rvec            **eigvec = NULL;
     matrix            box;
     int               natoms;
-    int               i, j, k, kmode, d, s, v;
+    int               i, j, k, kmode, d;
     gmx_bool          bDMR, bDMA, bFit;
-    char        *     indexfile;
 
-    char        *     grpname;
     real        *     eigval;
-    int               neigval;
     int        *      dummy;
     real        *     invsqrtm;
     char              title[STRLEN];
     real              fraction;
     int              *out_eigidx;
-    real             *out_eigval;
     rvec        *     this_eigvec;
-    real              omega, Ekin, sum, m, vel;
-    gmx_bool          found;
+    real              omega, Ekin, m, vel;
     int               nmodes, nphases;
     int              *imodes;
     real             *amplitude;
     real             *phases;
-    real              dum;
     const char       *p;
     char             *pe;
     output_env_t      oenv;
@@ -152,33 +146,20 @@ int gmx_nmtraj(int argc, char *argv[])
     /* Find vectors and phases */
 
     /* first find number of args in string */
-    nmodes = 0;
-    p      = eignrvec;
-    while (*p != 0)
-    {
-        dum = strtod(p, &pe);
-        p   = pe;
-        nmodes++;
-    }
+    nmodes = gmx::countWords(eignrvec);
 
     snew(imodes, nmodes);
     p = eignrvec;
     for (i = 0; i < nmodes; i++)
     {
         /* C indices start on 0 */
-        imodes[i] = strtol(p, &pe, 10)-1;
+        imodes[i] = std::strtol(p, &pe, 10)-1;
         p         = pe;
     }
 
     /* Now read phases */
-    nphases = 0;
-    p       = phasevec;
-    while (*p != 0)
-    {
-        dum = strtod(p, &pe);
-        p   = pe;
-        nphases++;
-    }
+    nphases = gmx::countWords(phasevec);
+
     if (nphases > nmodes)
     {
         gmx_fatal(FARGS, "More phases than eigenvector indices specified.\n");
@@ -274,7 +255,7 @@ int gmx_nmtraj(int argc, char *argv[])
             /* Derive amplitude from temperature and eigenvalue if we can */
 
             /* Convert eigenvalue to angular frequency, in units s^(-1) */
-            omega = sqrt(eigval[kmode]*1.0E21/(AVOGADRO*AMU));
+            omega = std::sqrt(eigval[kmode]*1.0E21/(AVOGADRO*AMU));
             /* Harmonic motion will be x=x0 + A*sin(omega*t)*eigenvec.
              * The velocity is thus:
              *
@@ -307,7 +288,7 @@ int gmx_nmtraj(int argc, char *argv[])
             Ekin *= AMU*1E-18;
 
             /* Set the amplitude so the energy is kT/2 */
-            amplitude[i] = sqrt(0.5*BOLTZMANN*temp/Ekin);
+            amplitude[i] = std::sqrt(0.5*BOLTZMANN*temp/Ekin);
         }
         else
         {
@@ -323,7 +304,7 @@ int gmx_nmtraj(int argc, char *argv[])
 
     for (i = 0; i < nframes; i++)
     {
-        fraction = (real)i/(real)nframes;
+        fraction = static_cast<real>(i)/nframes;
         for (j = 0; j < natoms; j++)
         {
             copy_rvec(xav[j], xout[j]);
@@ -338,11 +319,11 @@ int gmx_nmtraj(int argc, char *argv[])
             {
                 for (d = 0; d < DIM; d++)
                 {
-                    xout[j][d] += amplitude[k]*sin(2*M_PI*(fraction+phases[k]/360.0))*this_eigvec[j][d];
+                    xout[j][d] += amplitude[k]*std::sin(2*M_PI*(fraction+phases[k]/360.0))*this_eigvec[j][d];
                 }
             }
         }
-        write_trx(out, natoms, dummy, atoms, i, (real)i/(real)nframes, box, xout, NULL, NULL);
+        write_trx(out, natoms, dummy, atoms, i, static_cast<real>(i)/nframes, box, xout, NULL, NULL);
     }
 
     fprintf(stderr, "\n");