Signed-off-by: Miroslav Lichvar <mlich...@redhat.com>
---
 linreg.c | 63 +++++++++++++++++++++++++++++++++++++++++++--------------------
 1 file changed, 43 insertions(+), 20 deletions(-)

diff --git a/linreg.c b/linreg.c
index 3f7fe9a..e2915c5 100644
--- a/linreg.c
+++ b/linreg.c
@@ -42,6 +42,7 @@
 struct point {
        uint64_t x;
        uint64_t y;
+       double w;
 };
 
 struct result {
@@ -70,6 +71,8 @@ struct linreg_servo {
        uint64_t last_update;
        /* Regression results for all sizes */
        struct result results[MAX_SIZE - MIN_SIZE + 1];
+       /* Selected size */
+       unsigned int size;
        /* Current frequency offset of the clock */
        double clock_freq;
        /* Expected interval between updates */
@@ -120,12 +123,13 @@ static void update_reference(struct linreg_servo *s, 
uint64_t local_ts)
        s->last_update = local_ts;
 }
 
-static void add_sample(struct linreg_servo *s, int64_t offset)
+static void add_sample(struct linreg_servo *s, int64_t offset, double weight)
 {
        s->last_point = (s->last_point + 1) % MAX_POINTS;
 
        s->points[s->last_point].x = s->reference.x;
        s->points[s->last_point].y = s->reference.y - offset;
+       s->points[s->last_point].w = weight;
 
        if (s->num_points < MAX_POINTS)
                s->num_points++;
@@ -133,11 +137,12 @@ static void add_sample(struct linreg_servo *s, int64_t 
offset)
 
 static void regress(struct linreg_servo *s)
 {
-       double x, y, y0, e, x_sum, y_sum, xy_sum, x2_sum;
+       double x, y, y0, e, x_sum, y_sum, xy_sum, x2_sum, w, w_sum, w2_sum;
        unsigned int i, l, n, size;
        struct result *res;
 
        x_sum = 0.0, y_sum = 0.0, xy_sum = 0.0, x2_sum = 0.0;
+       w_sum = 0.0; w2_sum = 0.0;
        i = 0;
 
        y0 = (int64_t)(s->points[s->last_point].y - s->reference.y);
@@ -169,27 +174,31 @@ static void regress(struct linreg_servo *s)
 
                        x = (int64_t)(s->points[l].x - s->reference.x);
                        y = (int64_t)(s->points[l].y - s->reference.y);
-
-                       x_sum += x;
-                       y_sum += y;
-                       xy_sum += x * y;
-                       x2_sum += x * x;
+                       w = s->points[l].w;
+
+                       x_sum += x * w;
+                       y_sum += y * w;
+                       xy_sum += x * y * w;
+                       x2_sum += x * x * w;
+                       w_sum += w;
+                       w2_sum += w * w;
                }
 
                /* Get new intercept and slope */
-               res->slope = (xy_sum - x_sum * y_sum / n) /
-                               (x2_sum - x_sum * x_sum / n);
-               res->intercept = (y_sum - res->slope * x_sum) / n;
+               res->slope = (xy_sum - x_sum * y_sum / w_sum) /
+                               (x2_sum - x_sum * x_sum / w_sum);
+               res->intercept = (y_sum - res->slope * x_sum) / w_sum;
        }
 }
 
-/* Return largest size with smallest prediction error */
-static int get_best_size(struct linreg_servo *s)
+static void update_size(struct linreg_servo *s)
 {
        struct result *res;
        double best_err;
        int size, best_size;
 
+       /* Find largest size with smallest prediction error */
+
        best_size = 0;
        best_err = 0.0;
 
@@ -203,7 +212,19 @@ static int get_best_size(struct linreg_servo *s)
                }
        }
 
-       return best_size;
+       s->size = best_size;
+}
+
+static int linreg_weight_samples(struct servo *servo)
+{
+       struct linreg_servo *s = container_of(servo, struct linreg_servo, 
servo);
+
+       /*
+        * When the minimum size is selected, assume it is still too large
+        * (frequency is changing too quickly) and weighted samples would make
+        * it even more difficult to track.
+        */
+       return s->size > MIN_SIZE;
 }
 
 static double linreg_sample(struct servo *servo,
@@ -214,7 +235,7 @@ static double linreg_sample(struct servo *servo,
 {
        struct linreg_servo *s = container_of(servo, struct linreg_servo, 
servo);
        struct result *res;
-       int size, corr_interval;
+       int corr_interval;
 
        /*
         * The current time and the time when will be the frequency of the
@@ -225,21 +246,21 @@ static double linreg_sample(struct servo *servo,
         */
 
        update_reference(s, local_ts);
-       add_sample(s, offset);
+       add_sample(s, offset, weight);
        regress(s);
 
-       size = get_best_size(s);
+       update_size(s);
 
-       if (size < MIN_SIZE) {
+       if (s->size < MIN_SIZE) {
                /* Not enough points, wait for more */
                *state = SERVO_UNLOCKED;
                return -s->clock_freq;
        }
 
-       res = &s->results[size - MIN_SIZE];
+       res = &s->results[s->size - MIN_SIZE];
 
        pr_debug("linreg: points %d slope %.9f intercept %.0f err %.0f",
-                       1 << size, res->slope, res->intercept, res->err);
+                1 << s->size, res->slope, res->intercept, res->err);
 
        if ((servo->first_update &&
             servo->first_step_threshold &&
@@ -265,7 +286,7 @@ static double linreg_sample(struct servo *servo,
         * correction slowing down the clock will result in an overshoot. With
         * the system clock's maximum adjustment of 10% that's acceptable.
         */
-       corr_interval = size <= 4 ? 1 : size / 2;
+       corr_interval = s->size <= 4 ? 1 : s->size / 2;
        s->clock_freq += res->intercept / s->update_interval / corr_interval;
 
        /* Clamp the frequency to the allowed maximum */
@@ -293,6 +314,7 @@ static void linreg_reset(struct servo *servo)
 
        s->num_points = 0;
        s->last_update = 0;
+       s->size = 0;
        s->frequency_ratio = 1.0;
 
        for (i = MIN_SIZE; i <= MAX_SIZE; i++) {
@@ -331,6 +353,7 @@ struct servo *linreg_servo_create(int fadj)
                return NULL;
 
        s->servo.destroy = linreg_destroy;
+       s->servo.weight_samples = linreg_weight_samples;
        s->servo.sample = linreg_sample;
        s->servo.sync_interval = linreg_sync_interval;
        s->servo.reset = linreg_reset;
-- 
2.1.0


------------------------------------------------------------------------------
Dive into the World of Parallel Programming. The Go Parallel Website,
sponsored by Intel and developed in partnership with Slashdot Media, is your
hub for all things parallel software development, from weekly thought
leadership blogs to news, videos, case studies, tutorials and more. Take a
look and join the conversation now. http://goparallel.sourceforge.net/
_______________________________________________
Linuxptp-devel mailing list
Linuxptp-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/linuxptp-devel

Reply via email to