Split out a function.

[foosball] / foorank.cpp

diff --git a/foorank.cpp b/foorank.cpp
index b7d9011194e338df5aed675437f6885144a69e70..080798caec02676ecb649982b36595c4e28ef15f 100644 (file)
--- a/foorank.cpp
+++ b/foorank.cpp
@@ -237,24 +237,34 @@ void solve3x3(double *A, double *x, double *B)
 // of statistical moments.
 void estimate_musigma(vector<pair<double, double> > &curve, double &mu_result, double &sigma_result)
 {
-       double sum_area = 0.0;
-       double ex = 0.0;
-       double ex2 = 0.0;
-
-       for (unsigned i = 1; i < curve.size(); ++i) {
-               double x1 = curve[i].first;
-               double x0 = curve[i-1].first;
-               double y1 = curve[i].second;
-               double y0 = curve[i-1].second;
-               double xm = 0.5 * (x0 + x1);
-               double ym = 0.5 * (y0 + y1);
-               sum_area += (x1-x0) * ym;
-               ex += (x1-x0) * xm * ym;
-               ex2 += (x1-x0) * xm * xm * ym;
+       double h = (curve.back().first - curve.front().first) / (curve.size() - 1);
+
+       double area = curve.front().second;
+       double ex = curve.front().first * curve.front().second;
+       double ex2 = curve.front().first * curve.front().first * curve.front().second;
+
+       for (unsigned i = 1; i < curve.size() - 1; i += 2) {
+               double x = curve[i].first;
+               double y = curve[i].second;
+               area += 4.0 * y;
+               ex += 4.0 * x * y;
+               ex2 += 4.0 * x * x * y;
        }
+       for (unsigned i = 2; i < curve.size() - 1; i += 2) {
+               double x = curve[i].first;
+               double y = curve[i].second;
+               area += 2.0 * y;
+               ex += 2.0 * x * y;
+               ex2 += 2.0 * x * x * y;
+       }
+       
+       area += curve.back().second;
+       ex += curve.back().first * curve.back().second;
+       ex2 += curve.back().first * curve.back().first * curve.back().second;
 
-       ex /= sum_area;
-       ex2 /= sum_area;
+       area = (h/3.0) * area;
+       ex = (h/3.0) * ex / area;
+       ex2 = (h/3.0) * ex2 / area;
 
        mu_result = ex;
        sigma_result = sqrt(ex2 - ex * ex);
@@ -333,14 +343,8 @@ void least_squares(vector<pair<double, double> > &curve, double mu1, double sigm
        sigma_result = sigma;
 }
 
-int main(int argc, char **argv)
+void compute_new_rating(double mu1, double sigma1, double mu2, double sigma2, int score1, int score2, double &mu, double &sigma)
 {
-       double mu1 = atof(argv[1]);
-       double sigma1 = atof(argv[2]);
-       double mu2 = atof(argv[3]);
-       double sigma2 = atof(argv[4]);
-       int score1 = atoi(argv[5]);
-       int score2 = atoi(argv[6]);
        vector<pair<double, double> > curve;
 
        if (score1 == 10) {
@@ -357,9 +361,22 @@ int main(int argc, char **argv)
                }
        }
 
-       double mu_est, sigma_est, mu, sigma;
+       double mu_est, sigma_est;
        normalize(curve);
        estimate_musigma(curve, mu_est, sigma_est);
        least_squares(curve, mu_est, sigma_est, mu, sigma);
+}
+
+int main(int argc, char **argv)
+{
+       double mu1 = atof(argv[1]);
+       double sigma1 = atof(argv[2]);
+       double mu2 = atof(argv[3]);
+       double sigma2 = atof(argv[4]);
+       int score1 = atoi(argv[5]);
+       int score2 = atoi(argv[6]);
+       double mu, sigma;
+       compute_new_rating(mu1, sigma1, mu2, sigma2, score1, score2, mu, sigma);
        printf("%f %f\n", mu, sigma);
 }
+