Merge branch 'master' into pygromacs

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

diff --git a/src/gromacs/gmxana/gmx_helixorient.c b/src/gromacs/gmxana/gmx_helixorient.cpp
similarity index 86%
rename from src/gromacs/gmxana/gmx_helixorient.c
rename to src/gromacs/gmxana/gmx_helixorient.cpp
index 0cb969d0533b3221a92d6ed8b4e0752c1fb808a8..ddd92064abd578fa46ba5e9c08181dac19da1045 100644 (file)
--- a/src/gromacs/gmxana/gmx_helixorient.c
+++ b/src/gromacs/gmxana/gmx_helixorient.cpp
@@ -36,7 +36,7 @@
  */
 #include "gmxpre.h"
 
-#include <math.h>
+#include <cmath>
 
 #include "gromacs/commandline/pargs.h"
 #include "gromacs/fileio/trxio.h"
@@ -77,20 +77,19 @@ int gmx_helixorient(int argc, char *argv[])
 
     t_topology      *top = NULL;
     real             t;
-    rvec            *x = NULL, dx;
+    rvec            *x = NULL;
     matrix           box;
     t_trxstatus     *status;
     int              natoms;
     real             theta1, theta2, theta3;
 
-    int              d, i, j, teller = 0;
+    int              i, j, teller = 0;
     int              iCA, iSC;
     atom_id         *ind_CA;
     atom_id         *ind_SC;
     char            *gn_CA;
     char            *gn_SC;
-    rvec             averageaxis;
-    rvec             v1, v2, p1, p2, vtmp, vproj;
+    rvec             v1, v2;
     rvec            *x_CA, *x_SC;
     rvec            *r12;
     rvec            *r23;
@@ -103,8 +102,6 @@ int gmx_helixorient(int argc, char *argv[])
     rvec            *residuevector;
     rvec            *sidechainvector;
 
-    rvec             axes_t0[3];
-    rvec             axes[3];
     rvec            *residuehelixaxis_t0;
     rvec            *residuevector_t0;
     rvec            *axis3_t0;
@@ -122,7 +119,7 @@ int gmx_helixorient(int argc, char *argv[])
     real            *rise, *residuerise;
     real            *residuebending;
 
-    real             tmp, rotangle;
+    real             tmp;
     real             weight[3];
     t_pbc            pbc;
     matrix           A;
@@ -284,9 +281,9 @@ int gmx_helixorient(int argc, char *argv[])
             svmul(1.0/norm(helixaxis[i]), helixaxis[i], helixaxis[i]);
 
             tmp       = cos_angle(diff13[i], diff24[i]);
-            twist[i]  = 180.0/M_PI * acos( tmp );
-            radius[i] = sqrt( norm(diff13[i])*norm(diff24[i]) ) / (2.0* (1.0-tmp) );
-            rise[i]   = fabs(iprod(r23[i], helixaxis[i]));
+            twist[i]  = 180.0/M_PI * std::acos( tmp );
+            radius[i] = std::sqrt( norm(diff13[i])*norm(diff24[i]) ) / (2.0* (1.0-tmp) );
+            rise[i]   = std::abs(iprod(r23[i], helixaxis[i]));
 
             svmul(radius[i]/norm(diff13[i]), diff13[i], v1);
             svmul(radius[i]/norm(diff24[i]), diff24[i], v2);
@@ -306,7 +303,7 @@ int gmx_helixorient(int argc, char *argv[])
             residueradius[i]  = 0.5*(radius[i-2]+radius[i-1]);
             residuetwist[i]   = 0.5*(twist[i-2]+twist[i-1]);
             residuerise[i]    = 0.5*(rise[i-2]+rise[i-1]);
-            residuebending[i] = 180.0/M_PI*acos( cos_angle(helixaxis[i-2], helixaxis[i-1]) );
+            residuebending[i] = 180.0/M_PI*std::acos( cos_angle(helixaxis[i-2], helixaxis[i-1]) );
         }
         residueradius[iCA-2]  = radius[iCA-4];
         residuetwist[iCA-2]   = twist[iCA-4];
@@ -368,29 +365,6 @@ int gmx_helixorient(int argc, char *argv[])
                 fprintf(fpradius, "%15.12g  ", residueradius[i]);
                 fprintf(fptwist, "%15.12g  ", residuetwist[i]);
                 fprintf(fpbending, "%15.12g  ", residuebending[i]);
-
-                /* angle with local vector? */
-                /*
-                   printf("res[%2d]:  axis: %g %g %g    origin: %g %g %g   vector: %g %g %g   angle: %g\n",i,
-                       residuehelixaxis[i][0],
-                       residuehelixaxis[i][1],
-                       residuehelixaxis[i][2],
-                       residueorigin[i][0],
-                       residueorigin[i][1],
-                       residueorigin[i][2],
-                       residuevector[i][0],
-                       residuevector[i][1],
-                       residuevector[i][2],
-                       180.0/M_PI*acos( cos_angle(residuevector[i],residuehelixaxis[i]) ));
-                 */
-                /*      fprintf(fp,"%15.12g %15.12g %15.12g   %15.12g %15.12g %15.12g\n",
-                              residuehelixaxis[i][0],
-                              residuehelixaxis[i][1],
-                              residuehelixaxis[i][2],
-                              residuevector[i][0],
-                              residuevector[i][1],
-                              residuevector[i][2]);
-                 */
             }
         }
         fprintf(fprise, "\n");
@@ -443,17 +417,6 @@ int gmx_helixorient(int argc, char *argv[])
                     copy_rvec(residuevector[i], newaxes[1]);
                     copy_rvec(axis3[i], newaxes[2]);
 
-                    /*
-                       printf("frame %d, i=%d:\n  old: %g %g %g , %g %g %g , %g %g %g\n  new:  %g %g %g , %g %g %g , %g %g %g\n",
-                       teller,i,
-                       refaxes[0][0],refaxes[0][1],refaxes[0][2],
-                               refaxes[1][0],refaxes[1][1],refaxes[1][2],
-                               refaxes[2][0],refaxes[2][1],refaxes[2][2],
-                               newaxes[0][0],newaxes[0][1],newaxes[0][2],
-                               newaxes[1][0],newaxes[1][1],newaxes[1][2],
-                               newaxes[2][0],newaxes[2][1],newaxes[2][2]);
-                     */
-
                     /* rotate reference frame onto unit axes */
                     calc_fit_R(3, 3, weight, unitaxes, refaxes, A);
                     for (j = 0; j < 3; j++)
@@ -470,16 +433,11 @@ int gmx_helixorient(int argc, char *argv[])
                      * A contains rotation column vectors.
                      */
 
-                    /*
-                       printf("frame %d, i=%d, A: %g %g %g  , %g %g %g , %g %g %g\n",
-                       teller,i,A[0][0],A[0][1],A[0][2],A[1][0],A[1][1],A[1][2],A[2][0],A[2][1],A[2][2]);
-                     */
-
-                    theta1 = 180.0/M_PI*atan2(A[0][2], A[0][0]);
-                    theta2 = 180.0/M_PI*asin(-A[0][1]);
-                    theta3 = 180.0/M_PI*atan2(A[2][1], A[1][1]);
+                    theta1 = 180.0/M_PI*std::atan2(A[0][2], A[0][0]);
+                    theta2 = 180.0/M_PI*std::asin(-A[0][1]);
+                    theta3 = 180.0/M_PI*std::atan2(A[2][1], A[1][1]);
 
-                    tilt     = sqrt(theta1*theta1+theta2*theta2);
+                    tilt     = std::sqrt(theta1*theta1+theta2*theta2);
                     rotation = theta3;
                     fprintf(fptheta1, "%15.12g  ", theta1);
                     fprintf(fptheta2, "%15.12g  ", theta2);