Merge branch 'master' into pygromacs
[alexxy/gromacs.git] / src / gromacs / gmxana / gmx_mindist.cpp
similarity index 95%
rename from src/gromacs/gmxana/gmx_mindist.c
rename to src/gromacs/gmxana/gmx_mindist.cpp
index 214ba9ef4099c7dfc095bb48d21773c5dadc37f6..44f5d36af840cee470734d9551051289ca65c27f 100644 (file)
  */
 #include "gmxpre.h"
 
-#include <math.h>
-#include <stdlib.h>
-#include <string.h>
+#include <cmath>
+#include <cstdlib>
+#include <cstring>
+
+#include <algorithm>
 
 #include "gromacs/commandline/pargs.h"
 #include "gromacs/fileio/confio.h"
@@ -55,6 +57,7 @@
 #include "gromacs/topology/index.h"
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/futil.h"
+#include "gromacs/utility/gmxassert.h"
 #include "gromacs/utility/smalloc.h"
 
 
@@ -67,10 +70,10 @@ static void periodic_dist(int ePBC,
     real sqr_box, r2min, r2max, r2;
     rvec shift[NSHIFT_MAX], d0, d;
 
-    sqr_box = min(norm2(box[XX]), norm2(box[YY]));
+    sqr_box = std::min(norm2(box[XX]), norm2(box[YY]));
     if (ePBC == epbcXYZ)
     {
-        sqr_box = min(sqr_box, norm2(box[ZZ]));
+        sqr_box = std::min(sqr_box, norm2(box[ZZ]));
         nsz     = 1;
     }
     else if (ePBC == epbcXY)
@@ -131,8 +134,8 @@ static void periodic_dist(int ePBC,
         }
     }
 
-    *rmin = sqrt(r2min);
-    *rmax = sqrt(r2max);
+    *rmin = std::sqrt(r2min);
+    *rmax = std::sqrt(r2max);
 }
 
 static void periodic_mindist_plot(const char *trxfn, const char *outfn,
@@ -147,7 +150,7 @@ static void periodic_mindist_plot(const char *trxfn, const char *outfn,
     rvec        *x;
     matrix       box;
     int          natoms, ind_min[2] = {0, 0}, ind_mini = 0, ind_minj = 0;
-    real         r, rmin, rmax, rmint, tmint;
+    real         rmin, rmax, rmint, tmint;
     gmx_bool     bFirst;
     gmx_rmpbc_t  gpbc = NULL;
 
@@ -187,7 +190,7 @@ static void periodic_mindist_plot(const char *trxfn, const char *outfn,
             ind_mini = ind_min[0];
             ind_minj = ind_min[1];
         }
-        if (bSplit && !bFirst && fabs(t/output_env_get_time_factor(oenv)) < 1e-5)
+        if (bSplit && !bFirst && std::abs(t/output_env_get_time_factor(oenv)) < 1e-5)
         {
             fprintf(out, "%s\n", output_env_get_print_xvgr_codes(oenv) ? "&" : "");
         }
@@ -316,8 +319,8 @@ static void calc_dist(real rcut, gmx_bool bPBC, int ePBC, matrix box, rvec x[],
             *nmax += nmax_j;
         }
     }
-    *rmin = sqrt(rmin2);
-    *rmax = sqrt(rmax2);
+    *rmin = std::sqrt(rmin2);
+    *rmax = std::sqrt(rmax2);
 }
 
 void dist_plot(const char *fn, const char *afile, const char *dfile,
@@ -335,16 +338,17 @@ void dist_plot(const char *fn, const char *afile, const char *dfile,
     real             t, dmin, dmax, **mindres = NULL, **maxdres = NULL;
     int              nmin, nmax;
     t_trxstatus     *status;
-    int              i = -1, j, k, natoms;
-    int              min1, min2, max1, max2, min1r, min2r, max1r, max2r;
+    int              i = -1, j, k;
+    int              min2, max2, min1r, min2r, max1r, max2r;
+    int              min1 = 0;
+    int              max1 = 0;
     atom_id          oindex[2];
     rvec            *x0;
     matrix           box;
-    t_trxframe       frout;
     gmx_bool         bFirst;
     FILE            *respertime = NULL;
 
-    if ((natoms = read_first_x(oenv, &status, fn, &t, &x0, box)) == 0)
+    if (read_first_x(oenv, &status, fn, &t, &x0, box) == 0)
     {
         gmx_fatal(FARGS, "Could not read coordinates from statusfile\n");
     }
@@ -419,7 +423,6 @@ void dist_plot(const char *fn, const char *afile, const char *dfile,
 
     }
 
-    j = 0;
     if (nres)
     {
         snew(mindres, ng-1);
@@ -438,7 +441,7 @@ void dist_plot(const char *fn, const char *afile, const char *dfile,
     bFirst = TRUE;
     do
     {
-        if (bSplit && !bFirst && fabs(t/output_env_get_time_factor(oenv)) < 1e-5)
+        if (bSplit && !bFirst && std::abs(t/output_env_get_time_factor(oenv)) < 1e-5)
         {
             fprintf(dist, "%s\n", output_env_get_print_xvgr_codes(oenv) ? "&" : "");
             if (num)
@@ -503,8 +506,8 @@ void dist_plot(const char *fn, const char *afile, const char *dfile,
                         calc_dist(rcut, bPBC, ePBC, box, x0, residue[j+1]-residue[j], gnx[i],
                                   &(index[0][residue[j]]), index[i], bGroup,
                                   &dmin, &dmax, &nmin, &nmax, &min1r, &min2r, &max1r, &max2r);
-                        mindres[i-1][j] = min(mindres[i-1][j], dmin);
-                        maxdres[i-1][j] = max(maxdres[i-1][j], dmax);
+                        mindres[i-1][j] = std::min(mindres[i-1][j], dmin);
+                        maxdres[i-1][j] = std::max(maxdres[i-1][j], dmax);
                     }
                 }
             }
@@ -658,9 +661,6 @@ int gmx_mindist(int argc, char *argv[])
         "of the three box vectors.[PAR]",
         "Also [gmx-distance] and [gmx-pairdist] calculate distances."
     };
-    const char     *bugs[] = {
-        "The [TT]-pi[tt] option is very slow."
-    };
 
     static gmx_bool bMat             = FALSE, bPI = FALSE, bSplit = FALSE, bMax = FALSE, bPBC = TRUE;
     static gmx_bool bGroup           = FALSE;
@@ -693,13 +693,11 @@ int gmx_mindist(int argc, char *argv[])
     t_topology     *top  = NULL;
     int             ePBC = -1;
     char            title[256];
-    real            t;
     rvec           *x;
     matrix          box;
     gmx_bool        bTop = FALSE;
 
-    FILE           *atm;
-    int             i, j, nres = 0;
+    int             i, nres = 0;
     const char     *trxfnm, *tpsfnm, *ndxfnm, *distfnm, *numfnm, *atmfnm, *oxfnm, *resfnm;
     char          **grpname;
     int            *gnx;
@@ -790,8 +788,9 @@ int gmx_mindist(int argc, char *argv[])
 
     if (resfnm)
     {
-        nres = find_residues(top ? &(top->atoms) : NULL,
-                             gnx[0], index[0], &residues);
+        GMX_RELEASE_ASSERT(top != NULL, "top pointer cannot be NULL when finding residues");
+        nres = find_residues(&(top->atoms), gnx[0], index[0], &residues);
+
         if (debug)
         {
             dump_res(debug, nres, residues, index[0]);