Merge branch 'master' into pygromacs

[alexxy/gromacs.git] / src / gromacs / gmxana / gmx_hydorder.cpp

diff --git a/src/gromacs/gmxana/gmx_hydorder.c b/src/gromacs/gmxana/gmx_hydorder.cpp
similarity index 95%
rename from src/gromacs/gmxana/gmx_hydorder.c
rename to src/gromacs/gmxana/gmx_hydorder.cpp
index 98e33e5c9a3fdd3a63a4bb3afe5fc07f5beb89ed..d1c27c207a878a7e645e045dd0fe9cf2bd60dcfa 100644 (file)
--- a/src/gromacs/gmxana/gmx_hydorder.c
+++ b/src/gromacs/gmxana/gmx_hydorder.cpp
@@ -35,8 +35,8 @@
 
 #include "gmxpre.h"
 
-#include <math.h>
-#include <string.h>
+#include <cmath>
+#include <cstring>
 
 #include "gromacs/commandline/pargs.h"
 #include "gromacs/fileio/confio.h"
@@ -44,6 +44,7 @@
 #include "gromacs/fileio/trxio.h"
 #include "gromacs/fileio/xvgr.h"
 #include "gromacs/gmxana/binsearch.h"
+#include "gromacs/gmxana/gmx_ana.h"
 #include "gromacs/gmxana/gstat.h"
 #include "gromacs/gmxana/powerspect.h"
 #include "gromacs/legacyheaders/macros.h"
@@ -55,6 +56,7 @@
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/fatalerror.h"
 #include "gromacs/utility/futil.h"
+#include "gromacs/utility/gmxassert.h"
 #include "gromacs/utility/smalloc.h"
 
 /* Print name of first atom in all groups in index file */
@@ -184,7 +186,7 @@ static void find_tetra_order_grid(t_topology top, int ePBC,
         rmean = 0;
         for (j = 0; (j < 4); j++)
         {
-            r_nn[j][i] = sqrt(r_nn[j][i]);
+            r_nn[j][i] = std::sqrt(r_nn[j][i]);
             rmean     += r_nn[j][i];
         }
         rmean /= 4;
@@ -275,17 +277,15 @@ static void calc_tetra_order_interface(const char *fnNDX, const char *fnTPS, con
                                        output_env_t oenv)
 {
     FILE         *fpsg   = NULL, *fpsk = NULL;
-    char         *sgslfn = "sg_ang_mesh";  /* Hardcoded filenames for debugging*/
-    char         *skslfn = "sk_dist_mesh";
     t_topology    top;
     int           ePBC;
-    char          title[STRLEN], subtitle[STRLEN];
+    char          title[STRLEN];
     t_trxstatus  *status;
     int           natoms;
     real          t;
     rvec         *xtop, *x;
     matrix        box;
-    real          sg, sk, sgintf, pos;
+    real          sg, sk, sgintf;
     atom_id     **index   = NULL;
     char        **grpname = NULL;
     int           i, j, k, n, *isize, ng, nslicez, framenr;
@@ -300,9 +300,9 @@ static void calc_tetra_order_interface(const char *fnNDX, const char *fnTPS, con
 
     read_tps_conf(fnTPS, title, &top, &ePBC, &xtop, NULL, box, FALSE);
 
-    *nslicex = (int)(box[XX][XX]/binw + onehalf); /*Calculate slicenr from binwidth*/
-    *nslicey = (int)(box[YY][YY]/binw + onehalf);
-    nslicez  = (int)(box[ZZ][ZZ]/binw +  onehalf);
+    *nslicex = static_cast<int>(box[XX][XX]/binw + onehalf); /*Calculate slicenr from binwidth*/
+    *nslicey = static_cast<int>(box[YY][YY]/binw + onehalf);
+    nslicez  = static_cast<int>(box[ZZ][ZZ]/binw +  onehalf);
 
 
 
@@ -367,6 +367,7 @@ static void calc_tetra_order_interface(const char *fnNDX, const char *fnTPS, con
 
         find_tetra_order_grid(top, ePBC, natoms, box, x, isize[0], index[0],
                               &sg, &sk, *nslicex, *nslicey, nslicez, sg_grid, sk_grid);
+        GMX_RELEASE_ASSERT(sk_fravg != NULL, "Trying to dereference NULL sk_fravg pointer");
         for (i = 0; i < *nslicex; i++)
         {
             for (j = 0; j < *nslicey; j++)
@@ -383,6 +384,7 @@ static void calc_tetra_order_interface(const char *fnNDX, const char *fnTPS, con
 
         if (framenr%tblock == 0)
         {
+            GMX_RELEASE_ASSERT(sk_4d != NULL, "Trying to dereference NULL sk_4d pointer");
             sk_4d[*nframes] = sk_fravg;
             sg_4d[*nframes] = sg_fravg;
             (*nframes)++;
@@ -399,8 +401,8 @@ static void calc_tetra_order_interface(const char *fnNDX, const char *fnTPS, con
     /*Debugging for printing out the entire order parameter meshes.*/
     if (debug)
     {
-        fpsg = xvgropen(sgslfn, "S\\sg\\N Angle Order Parameter / Meshpoint", "(nm)", "S\\sg\\N", oenv);
-        fpsk = xvgropen(skslfn, "S\\sk\\N Distance Order Parameter / Meshpoint", "(nm)", "S\\sk\\N", oenv);
+        fpsg = xvgropen("sg_ang_mesh", "S\\sg\\N Angle Order Parameter / Meshpoint", "(nm)", "S\\sg\\N", oenv);
+        fpsk = xvgropen("sk_dist_mesh", "S\\sk\\N Distance Order Parameter / Meshpoint", "(nm)", "S\\sk\\N", oenv);
         for (n = 0; n < (*nframes); n++)
         {
             fprintf(fpsg, "%i\n", n);
@@ -466,13 +468,8 @@ static void calc_tetra_order_interface(const char *fnNDX, const char *fnTPS, con
     }
 
 
-    /*sfree(perm);*/
     sfree(sk_4d);
     sfree(sg_4d);
-    /*sfree(sg_grid);*/
-    /*sfree(sk_grid);*/
-
-
 }
 
 static void writesurftoxpms(real ***surf, int tblocks, int xbins, int ybins, real bw, char **outfiles, int maplevels )
@@ -608,7 +605,7 @@ int gmx_hydorder(int argc, char *argv[])
     static gmx_bool    bRawOut   = FALSE;
     int                frames, xslices, yslices; /* Dimensions of interface arrays*/
     real            ***intfpos;                  /* Interface arrays (intfnr,t,xy) -potentially large */
-    static char       *normal_axis[] = { NULL, "z", "x", "y", NULL };
+    static const char *normal_axis[] = { NULL, "z", "x", "y", NULL };
 
     t_pargs            pa[] = {
         { "-d",   FALSE, etENUM, {normal_axis},
@@ -659,15 +656,16 @@ int gmx_hydorder(int argc, char *argv[])
     trxfnm = ftp2fn(efTRX, NFILE, fnm);
 
     /* Calculate axis */
-    if (strcmp(normal_axis[0], "x") == 0)
+    GMX_RELEASE_ASSERT(normal_axis[0] != NULL, "Option setting inconsistency; normal_axis[0] is NULL");
+    if (std::strcmp(normal_axis[0], "x") == 0)
     {
         axis = XX;
     }
-    else if (strcmp(normal_axis[0], "y") == 0)
+    else if (std::strcmp(normal_axis[0], "y") == 0)
     {
         axis = YY;
     }
-    else if (strcmp(normal_axis[0], "z") == 0)
+    else if (std::strcmp(normal_axis[0], "z") == 0)
     {
         axis = ZZ;
     }