X-Git-Url: http://git.kpe.io/?p=ctsim.git;a=blobdiff_plain;f=libctsupport%2Fmathfuncs.cpp;h=2eaa55153c6f73c6325dd66e4335dad3c0fcd912;hp=bc3b38308bfde55d9b6b2f7002d1801b15661b98;hb=f7d2b7144f32a7bd157b7689022e62944b82fcc1;hpb=2d39e823ba389fc68e5317c422b55be006094252

diff --git a/libctsupport/mathfuncs.cpp b/libctsupport/mathfuncs.cpp
index bc3b383..2eaa551 100644
--- a/libctsupport/mathfuncs.cpp
+++ b/libctsupport/mathfuncs.cpp
@@ -2,7 +2,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
-**  $Id: mathfuncs.cpp,v 1.1 2000/06/22 10:17:28 kevin Exp $
+**  $Id: mathfuncs.cpp,v 1.4 2000/12/21 03:40:58 kevin Exp $
 **
 **  This program is free software; you can redistribute it and/or modify
 **  it under the terms of the GNU General Public License (version 2) as
@@ -63,10 +63,10 @@ integrateSimpson (const double xmin, const double xmax, const double *y, const i
 
 
 /* NAME
- *    norm_angle       Normalize angle to 0 to 2 * PI range
+ *    normalizeAngle       Normalize angle to 0 to 2 * PI range
  *
  * SYNOPSIS
- *    t = norm_angle (theta)
+ *    t = normalizeAngle (theta)
  *    double t	       Normalized angle
  *    double theta     Input angle
  */
@@ -81,3 +81,57 @@ normalizeAngle (double theta)
   
   return (theta);
 }
+
+
+void 
+vectorNumericStatistics (std::vector<double> vec, const int nPoints, double& min, double& max, double& mean, double& mode, double& median, double& stddev)
+{
+	if (nPoints <= 0)
+		return;
+
+    mean = 0;
+    min = vec[0];
+    max = vec[0];
+	int i;
+    for (i = 0; i < nPoints; i++) {
+		double v = vec[i];
+		if (v > max)
+			max = v;
+		if (v < min)
+			min = v;
+		mean += v;
+    }
+    mean /= nPoints;
+	
+    static const int nbin = 1024;
+    int hist[ nbin ] = {0};
+    double spread = max - min;
+    mode = 0;
+    stddev = 0;
+    for (i = 0; i < nPoints; i++) {
+		double v = vec[i];
+		int b = static_cast<int>((((v - min) / spread) * (nbin - 1)) + 0.5);
+		hist[b]++;
+		double diff = (v - mean);
+		stddev += diff * diff;
+    }
+    stddev = sqrt (stddev / nPoints);
+	
+    int max_binindex = 0;
+    int max_bin = -1;
+    for (int ibin = 0; ibin < nbin; ibin++) {
+		if (hist[ibin] > max_bin) {
+			max_bin = hist[ibin];
+			max_binindex = ibin;
+		}
+    }
+	
+    mode = (max_binindex * spread / (nbin - 1)) + min;
+	
+	std::sort(vec.begin(), vec.end());
+	
+	if (nPoints % 2)  // Odd
+		median = vec[((nPoints - 1) / 2)];
+	else        // Even
+		median = (vec[ (nPoints / 2) - 1 ] + vec[ nPoints / 2 ]) / 2;
+}