Fixed g_rmsdist NOE calculation
authorCarlo Camilloni <carlo.camilloni@gmail.com>
Wed, 9 Apr 2014 22:38:08 +0000 (23:38 +0100)
committerCarlo Camilloni <carlo.camilloni@gmail.com>
Wed, 9 Apr 2014 22:38:08 +0000 (23:38 +0100)
This fix solve a series of bugs in NOE calculations:
1) incorrect number of frames in the calculation of the average distance
2) wrong selection of equivalent atoms
3) -equiv was wrongly documented

Change-Id: Icc6bbed0e1ba65774ad470c3edbd7a6a96e63ee6

src/tools/gmx_rmsdist.c

index 92e00cbff1df67eca147571d6e104e0e6fe5a78e..f5a9b58539445f1434aec869697865f84b1d6900 100644 (file)
@@ -327,7 +327,7 @@ static int analyze_noe_equivalent(const char *eq_fn,
                     if (bEquiv)
                     {
                         /* set index for matching atom */
-                        noe_index[j] = groupnr;
+                        noe_index[i] = groupnr;
                         /* skip matching atom */
                         i = j;
                     }
@@ -345,7 +345,7 @@ static int analyze_noe_equivalent(const char *eq_fn,
                    This is supposed to cover all CH3 groups and the like */
                 anmi   = *atoms->atomname[index[i]];
                 anmil  = strlen(anmi);
-                bMatch = i < isize-3 && anmi[anmil-1] == '1';
+                bMatch = i <= isize-3 && anmi[anmil-1] == '1';
                 if (bMatch)
                 {
                     for (j = 1; j < 3; j++)
@@ -648,7 +648,7 @@ int gmx_rmsdist (int argc, char *argv[])
         "equivalent atoms can be supplied ([TT]-equiv[tt]), each line containing",
         "a set of equivalent atoms specified as residue number and name and",
         "atom name; e.g.:[PAR]",
-        "[TT]3 SER  HB1 3 SER  HB2[tt][PAR]",
+        "[TT]HB* 3 SER  HB1 3 SER  HB2[tt][PAR]",
         "Residue and atom names must exactly match those in the structure",
         "file, including case. Specifying non-sequential atoms is undefined."
 
@@ -788,6 +788,7 @@ int gmx_rmsdist (int argc, char *argv[])
 
     /*do a first step*/
     natom = read_first_x(oenv, &status, ftp2fn(efTRX, NFILE, fnm), &t, &x, box);
+    teller=0;
 
     do
     {
@@ -795,14 +796,13 @@ int gmx_rmsdist (int argc, char *argv[])
 
         rmsnow = rms_diff(isize, d, d_r);
         fprintf(fp, "%g  %g\n", t, rmsnow);
+        teller++;
     }
     while (read_next_x(oenv, status, &t, natom, x, box));
     fprintf(stderr, "\n");
 
     ffclose(fp);
 
-    teller = nframes_read(status);
-
     close_trj(status);
 
     calc_rms(isize, teller, dtot, dtot2, mean, &meanmax, rms, &rmsmax, rmsc, &rmscmax);