Merge branch 'master' into rotation

Conflicts:
	src/tools/Makefile.am

diff --git a/include/gmx_statistics.h b/include/gmx_statistics.h
index 5d28891b98faa401154d8f801dd8a0eb4480e17c..afc124894a4bac075ab8f17ec85e3f7ac16bdea9 100644 (file)
--- a/include/gmx_statistics.h
+++ b/include/gmx_statistics.h
@@ -60,6 +60,8 @@ int gmx_stats_done(gmx_stats_t stats);
    sigma. Level needs to be larger than one obviously. */
 int gmx_stats_remove_outliers(gmx_stats_t stats,double level);
 
+
+
 int gmx_stats_add_point(gmx_stats_t stats,double x,double y,
                               double dx,double dy);
 
@@ -72,9 +74,11 @@ int gmx_stats_add_points(gmx_stats_t stats,int n,real *x,real *y,
    more points, and returns estatsNOPOINTS when the last point has
    been returned. Should be used in a while loop. Variables for either
    pointer may be NULL, in which case the routine can be used as an
-   expensive point counter. */
+   expensive point counter. 
+   If level > 0 then the outliers outside level*sigma are reported
+   only. */
 int gmx_stats_get_point(gmx_stats_t stats,real *x,real *y,
-                              real *dx,real *dy);
+                       real *dx,real *dy,real level);
 
 /* Fit the data to y = ax + b, possibly weighted, if uncertainties
    have been input. Returns slope in *a and intercept in b, *return

diff --git a/src/gmxlib/statistics/gmx_statistics.c b/src/gmxlib/statistics/gmx_statistics.c
index eb723bed61e13c1f96c3338abbd16fc9f5e2b7b1..4fd7ae4d5320227f6f09c2e267aa458d20a7003b 100644 (file)
--- a/src/gmxlib/statistics/gmx_statistics.c
+++ b/src/gmxlib/statistics/gmx_statistics.c
@@ -129,20 +129,34 @@ int gmx_stats_add_point(gmx_stats_t gstats,double x,double y,
 }
 
 int gmx_stats_get_point(gmx_stats_t gstats,real *x,real *y,
-                        real *dx,real *dy)
+                        real *dx,real *dy,real level)
 {
     gmx_stats *stats = (gmx_stats *) gstats;
-  
-    if (stats->np_c < stats->np) 
+    int  ok,outlier;
+    real rmsd,r;
+    
+    if ((ok = gmx_stats_get_rmsd(gstats,&rmsd)) != estatsOK)
+    {
+        return ok;
+    }
+    outlier = 0;
+    while ((outlier == 0) && (stats->np_c < stats->np))
     {
-        if (NULL != x)  *x  = stats->x[stats->np_c];
-        if (NULL != y)  *y  = stats->y[stats->np_c];
-        if (NULL != dx) *dx = stats->dx[stats->np_c];
-        if (NULL != dy) *dy = stats->dy[stats->np_c];
+        r = fabs(stats->x[stats->np_c] - stats->y[stats->np_c]);
+        outlier = (r > rmsd*level);
+        if (outlier)
+        {
+            if (NULL != x)  *x  = stats->x[stats->np_c];
+            if (NULL != y)  *y  = stats->y[stats->np_c];
+            if (NULL != dx) *dx = stats->dx[stats->np_c];
+            if (NULL != dy) *dy = stats->dy[stats->np_c];
+        }
         stats->np_c++;
-    
-        return estatsOK;
+        
+        if (outlier)    
+            return estatsOK;
     }
+    
     stats->np_c = 0;
   
     return estatsNO_POINTS;

diff --git a/src/tools/Makefile.am b/src/tools/Makefile.am
index 7ce8dbc44a866d36da77bf5689d0ee877096dfd1..9241a65c96624533618cb5c476eb7505f9887d71 100644 (file)
--- a/src/tools/Makefile.am
+++ b/src/tools/Makefile.am
@@ -64,7 +64,7 @@ bin_PROGRAMS = \
        g_enemat        g_energy        g_lie           g_filter        \
        g_gyrate        g_h2order       g_hbond         g_helix         \
        g_mindist       g_msd           g_morph         g_nmeig         \
-       g_nmens         g_order         g_membed        \
+       g_nmens         g_order         g_kinetics      g_membed        \
        g_polystat      g_potential     g_rama          \
        g_rdf           g_rms           g_rmsdist       g_rmsf          \
        g_rotacf        g_rotmat        g_saltbr        g_sas           \

diff --git a/src/tools/gmx_msd.c b/src/tools/gmx_msd.c
index e46f6f40776781ff84199db56b4581e81cebe912..3665f390557c168fe22383a744ae6cd0831d0c04 100644 (file)
--- a/src/tools/gmx_msd.c
+++ b/src/tools/gmx_msd.c
@@ -508,7 +508,7 @@ void printmol(t_corr *curr,const char *fn,
     for(j=0; (j<curr->nrestart); j++) {
       real xx,yy,dx,dy;
       
-      while(gmx_stats_get_point(curr->lsq[j][i],&xx,&yy,&dx,&dy) == estatsOK)
+      while(gmx_stats_get_point(curr->lsq[j][i],&xx,&yy,&dx,&dy,0) == estatsOK)
           gmx_stats_add_point(lsq1,xx,yy,dx,dy);
     }
     gmx_stats_get_ab(lsq1,elsqWEIGHT_NONE,&a,&b,NULL,NULL,NULL,NULL);