src/ui/events/gesture_detection/velocity_tracker.cc

   1 // Copyright 2014 The Chromium Authors. All rights reserved.
   2 // Use of this source code is governed by a BSD-style license that can be
   3 // found in the LICENSE file.
   4
   5 #include "ui/events/gesture_detection/velocity_tracker.h"
   6
   7 #include <cmath>
   8
   9 #include "base/logging.h"
  10 #include "ui/events/gesture_detection/motion_event.h"
  11
  12 using base::TimeDelta;
  13 using base::TimeTicks;
  14
  15 namespace ui {
  16
  17 // Implements a particular velocity tracker algorithm.
  18 class VelocityTrackerStrategy {
  19  public:
  20   virtual ~VelocityTrackerStrategy() {}
  21
  22   virtual void Clear() = 0;
  23   virtual void ClearPointers(BitSet32 id_bits) = 0;
  24   virtual void AddMovement(const base::TimeTicks& event_time,
  25                            BitSet32 id_bits,
  26                            const Position* positions) = 0;
  27   virtual bool GetEstimator(uint32_t id, Estimator* out_estimator) const = 0;
  28
  29  protected:
  30   VelocityTrackerStrategy() {}
  31 };
  32
  33 namespace {
  34
  35 COMPILE_ASSERT(MotionEvent::MAX_POINTER_ID < 32, max_pointer_id_too_large);
  36
  37 // Threshold between ACTION_MOVE events for determining that a pointer has
  38 // stopped moving. Some input devices do not send ACTION_MOVE events in the case
  39 // where a pointer has stopped.  We need to detect this case so that we can
  40 // accurately predict the velocity after the pointer starts moving again.
  41 const int kAssumePointerMoveStoppedTimeMs = 40;
  42
  43 // Threshold between ACTION_MOVE and ACTION_{UP|POINTER_UP} events for
  44 // determining that a pointer has stopped moving. This is a larger threshold
  45 // than |kAssumePointerMoveStoppedTimeMs|, as some devices may delay synthesis
  46 // of ACTION_{UP|POINTER_UP} to reduce risk of noisy release.
  47 const int kAssumePointerUpStoppedTimeMs = 80;
  48
  49 struct Position {
  50   float x, y;
  51 };
  52
  53 struct Estimator {
  54   static const uint8_t kMaxDegree = 4;
  55
  56   // Estimator time base.
  57   TimeTicks time;
  58
  59   // Polynomial coefficients describing motion in X and Y.
  60   float xcoeff[kMaxDegree + 1], ycoeff[kMaxDegree + 1];
  61
  62   // Polynomial degree (number of coefficients), or zero if no information is
  63   // available.
  64   uint32_t degree;
  65
  66   // Confidence (coefficient of determination), between 0 (no fit)
  67   // and 1 (perfect fit).
  68   float confidence;
  69
  70   inline void Clear() {
  71     time = TimeTicks();
  72     degree = 0;
  73     confidence = 0;
  74     for (size_t i = 0; i <= kMaxDegree; i++) {
  75       xcoeff[i] = 0;
  76       ycoeff[i] = 0;
  77     }
  78   }
  79 };
  80
  81 float VectorDot(const float* a, const float* b, uint32_t m) {
  82   float r = 0;
  83   while (m--) {
  84     r += *(a++) * *(b++);
  85   }
  86   return r;
  87 }
  88
  89 float VectorNorm(const float* a, uint32_t m) {
  90   float r = 0;
  91   while (m--) {
  92     float t = *(a++);
  93     r += t * t;
  94   }
  95   return sqrtf(r);
  96 }
  97
  98 // Velocity tracker algorithm based on least-squares linear regression.
  99 class LeastSquaresVelocityTrackerStrategy : public VelocityTrackerStrategy {
 100  public:
 101   enum Weighting {
 102     // No weights applied.  All data points are equally reliable.
 103     WEIGHTING_NONE,
 104
 105     // Weight by time delta.  Data points clustered together are weighted less.
 106     WEIGHTING_DELTA,
 107
 108     // Weight such that points within a certain horizon are weighed more than
 109     // those outside of that horizon.
 110     WEIGHTING_CENTRAL,
 111
 112     // Weight such that points older than a certain amount are weighed less.
 113     WEIGHTING_RECENT,
 114   };
 115
 116   enum Restriction {
 117     // There's no restriction on the output of the velocity tracker.
 118     RESTRICTION_NONE,
 119
 120     // If the velocity determined by the tracker is in a sufficiently different
 121     // direction from the primary motion of the finger for the events being
 122     // considered for velocity calculation, return a velocity of 0.
 123     RESTRICTION_ALIGNED_DIRECTIONS
 124   };
 125
 126   // Number of samples to keep.
 127   static const uint8_t kHistorySize = 20;
 128
 129   // Degree must be no greater than Estimator::kMaxDegree.
 130   LeastSquaresVelocityTrackerStrategy(
 131       uint32_t degree,
 132       Weighting weighting,
 133       Restriction restriction);
 134   ~LeastSquaresVelocityTrackerStrategy() override;
 135
 136   void Clear() override;
 137   void ClearPointers(BitSet32 id_bits) override;
 138   void AddMovement(const TimeTicks& event_time,
 139                    BitSet32 id_bits,
 140                    const Position* positions) override;
 141   bool GetEstimator(uint32_t id, Estimator* out_estimator) const override;
 142
 143  private:
 144   // Sample horizon.
 145   // We don't use too much history by default since we want to react to quick
 146   // changes in direction.
 147   static const uint8_t kHorizonMS = 100;
 148
 149   struct Movement {
 150     TimeTicks event_time;
 151     BitSet32 id_bits;
 152     Position positions[VelocityTracker::MAX_POINTERS];
 153
 154     inline const Position& GetPosition(uint32_t id) const {
 155       return positions[id_bits.get_index_of_bit(id)];
 156     }
 157   };
 158
 159   float ChooseWeight(uint32_t index) const;
 160
 161   const uint32_t degree_;
 162   const Weighting weighting_;
 163   const Restriction restriction_;
 164   uint32_t index_;
 165   Movement movements_[kHistorySize];
 166 };
 167
 168 // Velocity tracker algorithm that uses an IIR filter.
 169 class IntegratingVelocityTrackerStrategy : public VelocityTrackerStrategy {
 170  public:
 171   // Degree must be 1 or 2.
 172   explicit IntegratingVelocityTrackerStrategy(uint32_t degree);
 173   ~IntegratingVelocityTrackerStrategy() override;
 174
 175   void Clear() override;
 176   void ClearPointers(BitSet32 id_bits) override;
 177   void AddMovement(const TimeTicks& event_time,
 178                    BitSet32 id_bits,
 179                    const Position* positions) override;
 180   bool GetEstimator(uint32_t id, Estimator* out_estimator) const override;
 181
 182  private:
 183   // Current state estimate for a particular pointer.
 184   struct State {
 185     TimeTicks update_time;
 186     uint32_t degree;
 187
 188     float xpos, xvel, xaccel;
 189     float ypos, yvel, yaccel;
 190   };
 191
 192   const uint32_t degree_;
 193   BitSet32 pointer_id_bits_;
 194   State mPointerState[MotionEvent::MAX_POINTER_ID + 1];
 195
 196   void InitState(State& state,
 197                  const TimeTicks& event_time,
 198                  float xpos,
 199                  float ypos) const;
 200   void UpdateState(State& state,
 201                    const TimeTicks& event_time,
 202                    float xpos,
 203                    float ypos) const;
 204   void PopulateEstimator(const State& state, Estimator* out_estimator) const;
 205 };
 206
 207 VelocityTrackerStrategy* CreateStrategy(VelocityTracker::Strategy strategy) {
 208   LeastSquaresVelocityTrackerStrategy::Weighting none =
 209       LeastSquaresVelocityTrackerStrategy::WEIGHTING_NONE;
 210   LeastSquaresVelocityTrackerStrategy::Restriction no_restriction =
 211       LeastSquaresVelocityTrackerStrategy::RESTRICTION_NONE;
 212   switch (strategy) {
 213     case VelocityTracker::LSQ1:
 214       return new LeastSquaresVelocityTrackerStrategy(1, none, no_restriction);
 215     case VelocityTracker::LSQ2:
 216       return new LeastSquaresVelocityTrackerStrategy(2, none, no_restriction);
 217     case VelocityTracker::LSQ2_RESTRICTED:
 218       return new LeastSquaresVelocityTrackerStrategy(
 219           2, LeastSquaresVelocityTrackerStrategy::WEIGHTING_NONE,
 220           LeastSquaresVelocityTrackerStrategy::RESTRICTION_ALIGNED_DIRECTIONS);
 221     case VelocityTracker::LSQ3:
 222       return new LeastSquaresVelocityTrackerStrategy(3, none, no_restriction);
 223     case VelocityTracker::WLSQ2_DELTA:
 224       return new LeastSquaresVelocityTrackerStrategy(
 225           2, LeastSquaresVelocityTrackerStrategy::WEIGHTING_DELTA,
 226           no_restriction);
 227     case VelocityTracker::WLSQ2_CENTRAL:
 228       return new LeastSquaresVelocityTrackerStrategy(
 229           2, LeastSquaresVelocityTrackerStrategy::WEIGHTING_CENTRAL,
 230           no_restriction);
 231     case VelocityTracker::WLSQ2_RECENT:
 232       return new LeastSquaresVelocityTrackerStrategy(
 233           2, LeastSquaresVelocityTrackerStrategy::WEIGHTING_RECENT,
 234           no_restriction);
 235     case VelocityTracker::INT1:
 236       return new IntegratingVelocityTrackerStrategy(1);
 237     case VelocityTracker::INT2:
 238       return new IntegratingVelocityTrackerStrategy(2);
 239   }
 240   NOTREACHED() << "Unrecognized velocity tracker strategy: " << strategy;
 241   // Quadratic regression is a safe default.
 242   return CreateStrategy(VelocityTracker::STRATEGY_DEFAULT);
 243 }
 244
 245 }  // namespace
 246
 247 // --- VelocityTracker ---
 248
 249 VelocityTracker::VelocityTracker(Strategy strategy)
 250     : current_pointer_id_bits_(0),
 251       active_pointer_id_(-1),
 252       strategy_(CreateStrategy(strategy)) {}
 253
 254 VelocityTracker::~VelocityTracker() {}
 255
 256 void VelocityTracker::Clear() {
 257   current_pointer_id_bits_.clear();
 258   active_pointer_id_ = -1;
 259   strategy_->Clear();
 260 }
 261
 262 void VelocityTracker::ClearPointers(BitSet32 id_bits) {
 263   BitSet32 remaining_id_bits(current_pointer_id_bits_.value & ~id_bits.value);
 264   current_pointer_id_bits_ = remaining_id_bits;
 265
 266   if (active_pointer_id_ >= 0 && id_bits.has_bit(active_pointer_id_)) {
 267     active_pointer_id_ = !remaining_id_bits.is_empty()
 268                              ? remaining_id_bits.first_marked_bit()
 269                              : -1;
 270   }
 271
 272   strategy_->ClearPointers(id_bits);
 273 }
 274
 275 void VelocityTracker::AddMovement(const TimeTicks& event_time,
 276                                   BitSet32 id_bits,
 277                                   const Position* positions) {
 278   while (id_bits.count() > MAX_POINTERS)
 279     id_bits.clear_last_marked_bit();
 280
 281   if ((current_pointer_id_bits_.value & id_bits.value) &&
 282       (event_time - last_event_time_) >=
 283           base::TimeDelta::FromMilliseconds(kAssumePointerMoveStoppedTimeMs)) {
 284     // We have not received any movements for too long. Assume that all pointers
 285     // have stopped.
 286     strategy_->Clear();
 287   }
 288   last_event_time_ = event_time;
 289
 290   current_pointer_id_bits_ = id_bits;
 291   if (active_pointer_id_ < 0 || !id_bits.has_bit(active_pointer_id_))
 292     active_pointer_id_ = id_bits.is_empty() ? -1 : id_bits.first_marked_bit();
 293
 294   strategy_->AddMovement(event_time, id_bits, positions);
 295 }
 296
 297 void VelocityTracker::AddMovement(const MotionEvent& event) {
 298   int32_t actionMasked = event.GetAction();
 299
 300   switch (actionMasked) {
 301     case MotionEvent::ACTION_DOWN:
 302       // case MotionEvent::HOVER_ENTER:
 303       // Clear all pointers on down before adding the new movement.
 304       Clear();
 305       break;
 306     case MotionEvent::ACTION_POINTER_DOWN: {
 307       // Start a new movement trace for a pointer that just went down.
 308       // We do this on down instead of on up because the client may want to
 309       // query the final velocity for a pointer that just went up.
 310       BitSet32 downIdBits;
 311       downIdBits.mark_bit(event.GetPointerId(event.GetActionIndex()));
 312       ClearPointers(downIdBits);
 313       break;
 314     }
 315     case MotionEvent::ACTION_MOVE:
 316       // case MotionEvent::ACTION_HOVER_MOVE:
 317       break;
 318     case MotionEvent::ACTION_UP:
 319     case MotionEvent::ACTION_POINTER_UP:
 320       // Note that ACTION_UP and ACTION_POINTER_UP always report the last known
 321       // position of the pointers that went up.  ACTION_POINTER_UP does include
 322       // the new position of pointers that remained down but we will also
 323       // receive an ACTION_MOVE with this information if any of them actually
 324       // moved.  Since we don't know how many pointers will be going up at once
 325       // it makes sense to just wait for the following ACTION_MOVE before adding
 326       // the movement. However, if the up event itself is delayed because of
 327       // (difficult albeit possible) prolonged stationary screen contact, assume
 328       // that motion has stopped.
 329       if ((event.GetEventTime() - last_event_time_) >=
 330           base::TimeDelta::FromMilliseconds(kAssumePointerUpStoppedTimeMs))
 331         strategy_->Clear();
 332       return;
 333     default:
 334       // Ignore all other actions because they do not convey any new information
 335       // about pointer movement.  We also want to preserve the last known
 336       // velocity of the pointers.
 337       return;
 338   }
 339
 340   size_t pointer_count = event.GetPointerCount();
 341   if (pointer_count > MAX_POINTERS) {
 342     pointer_count = MAX_POINTERS;
 343   }
 344
 345   BitSet32 id_bits;
 346   for (size_t i = 0; i < pointer_count; i++) {
 347     id_bits.mark_bit(event.GetPointerId(i));
 348   }
 349
 350   uint32_t pointer_index[MAX_POINTERS];
 351   for (size_t i = 0; i < pointer_count; i++) {
 352     pointer_index[i] = id_bits.get_index_of_bit(event.GetPointerId(i));
 353   }
 354
 355   Position positions[MAX_POINTERS];
 356   size_t historySize = event.GetHistorySize();
 357   for (size_t h = 0; h < historySize; h++) {
 358     for (size_t i = 0; i < pointer_count; i++) {
 359       uint32_t index = pointer_index[i];
 360       positions[index].x = event.GetHistoricalX(i, h);
 361       positions[index].y = event.GetHistoricalY(i, h);
 362     }
 363     AddMovement(event.GetHistoricalEventTime(h), id_bits, positions);
 364   }
 365
 366   for (size_t i = 0; i < pointer_count; i++) {
 367     uint32_t index = pointer_index[i];
 368     positions[index].x = event.GetX(i);
 369     positions[index].y = event.GetY(i);
 370   }
 371   AddMovement(event.GetEventTime(), id_bits, positions);
 372 }
 373
 374 bool VelocityTracker::GetVelocity(uint32_t id,
 375                                   float* out_vx,
 376                                   float* out_vy) const {
 377   Estimator estimator;
 378   if (GetEstimator(id, &estimator) && estimator.degree >= 1) {
 379     *out_vx = estimator.xcoeff[1];
 380     *out_vy = estimator.ycoeff[1];
 381     return true;
 382   }
 383   *out_vx = 0;
 384   *out_vy = 0;
 385   return false;
 386 }
 387
 388 void LeastSquaresVelocityTrackerStrategy::AddMovement(
 389     const TimeTicks& event_time,
 390     BitSet32 id_bits,
 391     const Position* positions) {
 392   if (++index_ == kHistorySize) {
 393     index_ = 0;
 394   }
 395
 396   Movement& movement = movements_[index_];
 397   movement.event_time = event_time;
 398   movement.id_bits = id_bits;
 399   uint32_t count = id_bits.count();
 400   for (uint32_t i = 0; i < count; i++) {
 401     movement.positions[i] = positions[i];
 402   }
 403 }
 404
 405 bool VelocityTracker::GetEstimator(uint32_t id,
 406                                    Estimator* out_estimator) const {
 407   return strategy_->GetEstimator(id, out_estimator);
 408 }
 409
 410 // --- LeastSquaresVelocityTrackerStrategy ---
 411
 412 LeastSquaresVelocityTrackerStrategy::LeastSquaresVelocityTrackerStrategy(
 413     uint32_t degree,
 414     Weighting weighting,
 415     Restriction restriction)
 416     : degree_(degree),
 417       weighting_(weighting),
 418       restriction_(restriction) {
 419   DCHECK_LT(degree_, static_cast<uint32_t>(Estimator::kMaxDegree));
 420   Clear();
 421 }
 422
 423 LeastSquaresVelocityTrackerStrategy::~LeastSquaresVelocityTrackerStrategy() {}
 424
 425 void LeastSquaresVelocityTrackerStrategy::Clear() {
 426   index_ = 0;
 427   movements_[0].id_bits.clear();
 428 }
 429
 430 /**
 431  * Solves a linear least squares problem to obtain a N degree polynomial that
 432  * fits the specified input data as nearly as possible.
 433  *
 434  * Returns true if a solution is found, false otherwise.
 435  *
 436  * The input consists of two vectors of data points X and Y with indices 0..m-1
 437  * along with a weight vector W of the same size.
 438  *
 439  * The output is a vector B with indices 0..n that describes a polynomial
 440  * that fits the data, such the sum of W[i] * W[i] * abs(Y[i] - (B[0] + B[1]
 441  * X[i] * + B[2] X[i]^2 ... B[n] X[i]^n)) for all i between 0 and m-1 is
 442  * minimized.
 443  *
 444  * Accordingly, the weight vector W should be initialized by the caller with the
 445  * reciprocal square root of the variance of the error in each input data point.
 446  * In other words, an ideal choice for W would be W[i] = 1 / var(Y[i]) = 1 /
 447  * stddev(Y[i]).
 448  * The weights express the relative importance of each data point.  If the
 449  * weights are* all 1, then the data points are considered to be of equal
 450  * importance when fitting the polynomial.  It is a good idea to choose weights
 451  * that diminish the importance of data points that may have higher than usual
 452  * error margins.
 453  *
 454  * Errors among data points are assumed to be independent.  W is represented
 455  * here as a vector although in the literature it is typically taken to be a
 456  * diagonal matrix.
 457  *
 458  * That is to say, the function that generated the input data can be
 459  * approximated by y(x) ~= B[0] + B[1] x + B[2] x^2 + ... + B[n] x^n.
 460  *
 461  * The coefficient of determination (R^2) is also returned to describe the
 462  * goodness of fit of the model for the given data.  It is a value between 0
 463  * and 1, where 1 indicates perfect correspondence.
 464  *
 465  * This function first expands the X vector to a m by n matrix A such that
 466  * A[i][0] = 1, A[i][1] = X[i], A[i][2] = X[i]^2, ..., A[i][n] = X[i]^n, then
 467  * multiplies it by w[i]./
 468  *
 469  * Then it calculates the QR decomposition of A yielding an m by m orthonormal
 470  * matrix Q and an m by n upper triangular matrix R.  Because R is upper
 471  * triangular (lower part is all zeroes), we can simplify the decomposition into
 472  * an m by n matrix Q1 and a n by n matrix R1 such that A = Q1 R1.
 473  *
 474  * Finally we solve the system of linear equations given by
 475  * R1 B = (Qtranspose W Y) to find B.
 476  *
 477  * For efficiency, we lay out A and Q column-wise in memory because we
 478  * frequently operate on the column vectors.  Conversely, we lay out R row-wise.
 479  *
 480  * http://en.wikipedia.org/wiki/Numerical_methods_for_linear_least_squares
 481  * http://en.wikipedia.org/wiki/Gram-Schmidt
 482  */
 483 static bool SolveLeastSquares(const float* x,
 484                               const float* y,
 485                               const float* w,
 486                               uint32_t m,
 487                               uint32_t n,
 488                               float* out_b,
 489                               float* out_det) {
 490   // MSVC does not support variable-length arrays (used by the original Android
 491   // implementation of this function).
 492 #if defined(COMPILER_MSVC)
 493   const uint32_t M_ARRAY_LENGTH =
 494       LeastSquaresVelocityTrackerStrategy::kHistorySize;
 495   const uint32_t N_ARRAY_LENGTH = Estimator::kMaxDegree;
 496   DCHECK_LE(m, M_ARRAY_LENGTH);
 497   DCHECK_LE(n, N_ARRAY_LENGTH);
 498 #else
 499   const uint32_t M_ARRAY_LENGTH = m;
 500   const uint32_t N_ARRAY_LENGTH = n;
 501 #endif
 502
 503   // Expand the X vector to a matrix A, pre-multiplied by the weights.
 504   float a[N_ARRAY_LENGTH][M_ARRAY_LENGTH];  // column-major order
 505   for (uint32_t h = 0; h < m; h++) {
 506     a[0][h] = w[h];
 507     for (uint32_t i = 1; i < n; i++) {
 508       a[i][h] = a[i - 1][h] * x[h];
 509     }
 510   }
 511
 512   // Apply the Gram-Schmidt process to A to obtain its QR decomposition.
 513
 514   // Orthonormal basis, column-major order.
 515   float q[N_ARRAY_LENGTH][M_ARRAY_LENGTH];
 516   // Upper triangular matrix, row-major order.
 517   float r[N_ARRAY_LENGTH][N_ARRAY_LENGTH];
 518   for (uint32_t j = 0; j < n; j++) {
 519     for (uint32_t h = 0; h < m; h++) {
 520       q[j][h] = a[j][h];
 521     }
 522     for (uint32_t i = 0; i < j; i++) {
 523       float dot = VectorDot(&q[j][0], &q[i][0], m);
 524       for (uint32_t h = 0; h < m; h++) {
 525         q[j][h] -= dot * q[i][h];
 526       }
 527     }
 528
 529     float norm = VectorNorm(&q[j][0], m);
 530     if (norm < 0.000001f) {
 531       // vectors are linearly dependent or zero so no solution
 532       return false;
 533     }
 534
 535     float invNorm = 1.0f / norm;
 536     for (uint32_t h = 0; h < m; h++) {
 537       q[j][h] *= invNorm;
 538     }
 539     for (uint32_t i = 0; i < n; i++) {
 540       r[j][i] = i < j ? 0 : VectorDot(&q[j][0], &a[i][0], m);
 541     }
 542   }
 543
 544   // Solve R B = Qt W Y to find B.  This is easy because R is upper triangular.
 545   // We just work from bottom-right to top-left calculating B's coefficients.
 546   float wy[M_ARRAY_LENGTH];
 547   for (uint32_t h = 0; h < m; h++) {
 548     wy[h] = y[h] * w[h];
 549   }
 550   for (uint32_t i = n; i-- != 0;) {
 551     out_b[i] = VectorDot(&q[i][0], wy, m);
 552     for (uint32_t j = n - 1; j > i; j--) {
 553       out_b[i] -= r[i][j] * out_b[j];
 554     }
 555     out_b[i] /= r[i][i];
 556   }
 557
 558   // Calculate the coefficient of determination as 1 - (SSerr / SStot) where
 559   // SSerr is the residual sum of squares (variance of the error),
 560   // and SStot is the total sum of squares (variance of the data) where each
 561   // has been weighted.
 562   float ymean = 0;
 563   for (uint32_t h = 0; h < m; h++) {
 564     ymean += y[h];
 565   }
 566   ymean /= m;
 567
 568   float sserr = 0;
 569   float sstot = 0;
 570   for (uint32_t h = 0; h < m; h++) {
 571     float err = y[h] - out_b[0];
 572     float term = 1;
 573     for (uint32_t i = 1; i < n; i++) {
 574       term *= x[h];
 575       err -= term * out_b[i];
 576     }
 577     sserr += w[h] * w[h] * err * err;
 578     float var = y[h] - ymean;
 579     sstot += w[h] * w[h] * var * var;
 580   }
 581   *out_det = sstot > 0.000001f ? 1.0f - (sserr / sstot) : 1;
 582   return true;
 583 }
 584
 585 void LeastSquaresVelocityTrackerStrategy::ClearPointers(BitSet32 id_bits) {
 586   BitSet32 remaining_id_bits(movements_[index_].id_bits.value & ~id_bits.value);
 587   movements_[index_].id_bits = remaining_id_bits;
 588 }
 589
 590 bool LeastSquaresVelocityTrackerStrategy::GetEstimator(
 591     uint32_t id,
 592     Estimator* out_estimator) const {
 593   out_estimator->Clear();
 594
 595   // Iterate over movement samples in reverse time order and collect samples.
 596   float x[kHistorySize];
 597   float y[kHistorySize];
 598   float w[kHistorySize];
 599   float time[kHistorySize];
 600   uint32_t m = 0;
 601   uint32_t index = index_;
 602   const base::TimeDelta horizon = base::TimeDelta::FromMilliseconds(kHorizonMS);
 603   const Movement& newest_movement = movements_[index_];
 604   const Movement* first_movement = nullptr;
 605
 606   do {
 607     const Movement& movement = movements_[index];
 608     if (!movement.id_bits.has_bit(id))
 609       break;
 610
 611     first_movement = &movement;
 612     TimeDelta age = newest_movement.event_time - movement.event_time;
 613     if (age > horizon)
 614       break;
 615
 616     const Position& position = movement.GetPosition(id);
 617     x[m] = position.x;
 618     y[m] = position.y;
 619     w[m] = ChooseWeight(index);
 620     time[m] = -static_cast<float>(age.InSecondsF());
 621     index = (index == 0 ? kHistorySize : index) - 1;
 622   } while (++m < kHistorySize);
 623
 624   if (m == 0)
 625     return false;  // no data
 626
 627   // Calculate a least squares polynomial fit.
 628   uint32_t degree = degree_;
 629   if (degree > m - 1)
 630     degree = m - 1;
 631
 632   if (degree >= 1) {
 633     float xdet, ydet;
 634     uint32_t n = degree + 1;
 635     if (SolveLeastSquares(time, x, w, m, n, out_estimator->xcoeff, &xdet) &&
 636         SolveLeastSquares(time, y, w, m, n, out_estimator->ycoeff, &ydet)) {
 637       if (restriction_ == RESTRICTION_ALIGNED_DIRECTIONS) {
 638         DCHECK(first_movement);
 639         float dx = newest_movement.GetPosition(id).x -
 640                    first_movement->GetPosition(id).x;
 641         float dy = newest_movement.GetPosition(id).y -
 642                    first_movement->GetPosition(id).y;
 643
 644         // If the velocity is in a sufficiently different direction from the
 645         // primary movement, ignore it.
 646         if (out_estimator->xcoeff[1] * dx + out_estimator->ycoeff[1] * dy < 0)
 647           return false;
 648       }
 649
 650       out_estimator->time = newest_movement.event_time;
 651       out_estimator->degree = degree;
 652       out_estimator->confidence = xdet * ydet;
 653       return true;
 654     }
 655   }
 656
 657   // No velocity data available for this pointer, but we do have its current
 658   // position.
 659   out_estimator->xcoeff[0] = x[0];
 660   out_estimator->ycoeff[0] = y[0];
 661   out_estimator->time = newest_movement.event_time;
 662   out_estimator->degree = 0;
 663   out_estimator->confidence = 1;
 664   return true;
 665 }
 666
 667 float LeastSquaresVelocityTrackerStrategy::ChooseWeight(uint32_t index) const {
 668   switch (weighting_) {
 669     case WEIGHTING_DELTA: {
 670       // Weight points based on how much time elapsed between them and the next
 671       // point so that points that "cover" a shorter time span are weighed less.
 672       //   delta  0ms: 0.5
 673       //   delta 10ms: 1.0
 674       if (index == index_) {
 675         return 1.0f;
 676       }
 677       uint32_t next_index = (index + 1) % kHistorySize;
 678       float delta_millis =
 679           static_cast<float>((movements_[next_index].event_time -
 680                               movements_[index].event_time).InMillisecondsF());
 681       if (delta_millis < 0)
 682         return 0.5f;
 683       if (delta_millis < 10)
 684         return 0.5f + delta_millis * 0.05f;
 685
 686       return 1.0f;
 687     }
 688
 689     case WEIGHTING_CENTRAL: {
 690       // Weight points based on their age, weighing very recent and very old
 691       // points less.
 692       //   age  0ms: 0.5
 693       //   age 10ms: 1.0
 694       //   age 50ms: 1.0
 695       //   age 60ms: 0.5
 696       float age_millis =
 697           static_cast<float>((movements_[index_].event_time -
 698                               movements_[index].event_time).InMillisecondsF());
 699       if (age_millis < 0)
 700         return 0.5f;
 701       if (age_millis < 10)
 702         return 0.5f + age_millis * 0.05f;
 703       if (age_millis < 50)
 704         return 1.0f;
 705       if (age_millis < 60)
 706         return 0.5f + (60 - age_millis) * 0.05f;
 707
 708       return 0.5f;
 709     }
 710
 711     case WEIGHTING_RECENT: {
 712       // Weight points based on their age, weighing older points less.
 713       //   age   0ms: 1.0
 714       //   age  50ms: 1.0
 715       //   age 100ms: 0.5
 716       float age_millis =
 717           static_cast<float>((movements_[index_].event_time -
 718                               movements_[index].event_time).InMillisecondsF());
 719       if (age_millis < 50) {
 720         return 1.0f;
 721       }
 722       if (age_millis < 100) {
 723         return 0.5f + (100 - age_millis) * 0.01f;
 724       }
 725       return 0.5f;
 726     }
 727
 728     case WEIGHTING_NONE:
 729     default:
 730       return 1.0f;
 731   }
 732 }
 733
 734 // --- IntegratingVelocityTrackerStrategy ---
 735
 736 IntegratingVelocityTrackerStrategy::IntegratingVelocityTrackerStrategy(
 737     uint32_t degree)
 738     : degree_(degree) {
 739   DCHECK_LT(degree_, static_cast<uint32_t>(Estimator::kMaxDegree));
 740 }
 741
 742 IntegratingVelocityTrackerStrategy::~IntegratingVelocityTrackerStrategy() {}
 743
 744 void IntegratingVelocityTrackerStrategy::Clear() { pointer_id_bits_.clear(); }
 745
 746 void IntegratingVelocityTrackerStrategy::ClearPointers(BitSet32 id_bits) {
 747   pointer_id_bits_.value &= ~id_bits.value;
 748 }
 749
 750 void IntegratingVelocityTrackerStrategy::AddMovement(
 751     const TimeTicks& event_time,
 752     BitSet32 id_bits,
 753     const Position* positions) {
 754   uint32_t index = 0;
 755   for (BitSet32 iter_id_bits(id_bits); !iter_id_bits.is_empty();) {
 756     uint32_t id = iter_id_bits.clear_first_marked_bit();
 757     State& state = mPointerState[id];
 758     const Position& position = positions[index++];
 759     if (pointer_id_bits_.has_bit(id))
 760       UpdateState(state, event_time, position.x, position.y);
 761     else
 762       InitState(state, event_time, position.x, position.y);
 763   }
 764
 765   pointer_id_bits_ = id_bits;
 766 }
 767
 768 bool IntegratingVelocityTrackerStrategy::GetEstimator(
 769     uint32_t id,
 770     Estimator* out_estimator) const {
 771   out_estimator->Clear();
 772
 773   if (pointer_id_bits_.has_bit(id)) {
 774     const State& state = mPointerState[id];
 775     PopulateEstimator(state, out_estimator);
 776     return true;
 777   }
 778
 779   return false;
 780 }
 781
 782 void IntegratingVelocityTrackerStrategy::InitState(State& state,
 783                                                    const TimeTicks& event_time,
 784                                                    float xpos,
 785                                                    float ypos) const {
 786   state.update_time = event_time;
 787   state.degree = 0;
 788   state.xpos = xpos;
 789   state.xvel = 0;
 790   state.xaccel = 0;
 791   state.ypos = ypos;
 792   state.yvel = 0;
 793   state.yaccel = 0;
 794 }
 795
 796 void IntegratingVelocityTrackerStrategy::UpdateState(
 797     State& state,
 798     const TimeTicks& event_time,
 799     float xpos,
 800     float ypos) const {
 801   const base::TimeDelta MIN_TIME_DELTA = TimeDelta::FromMicroseconds(2);
 802   const float FILTER_TIME_CONSTANT = 0.010f;  // 10 milliseconds
 803
 804   if (event_time <= state.update_time + MIN_TIME_DELTA)
 805     return;
 806
 807   float dt = static_cast<float>((event_time - state.update_time).InSecondsF());
 808   state.update_time = event_time;
 809
 810   float xvel = (xpos - state.xpos) / dt;
 811   float yvel = (ypos - state.ypos) / dt;
 812   if (state.degree == 0) {
 813     state.xvel = xvel;
 814     state.yvel = yvel;
 815     state.degree = 1;
 816   } else {
 817     float alpha = dt / (FILTER_TIME_CONSTANT + dt);
 818     if (degree_ == 1) {
 819       state.xvel += (xvel - state.xvel) * alpha;
 820       state.yvel += (yvel - state.yvel) * alpha;
 821     } else {
 822       float xaccel = (xvel - state.xvel) / dt;
 823       float yaccel = (yvel - state.yvel) / dt;
 824       if (state.degree == 1) {
 825         state.xaccel = xaccel;
 826         state.yaccel = yaccel;
 827         state.degree = 2;
 828       } else {
 829         state.xaccel += (xaccel - state.xaccel) * alpha;
 830         state.yaccel += (yaccel - state.yaccel) * alpha;
 831       }
 832       state.xvel += (state.xaccel * dt) * alpha;
 833       state.yvel += (state.yaccel * dt) * alpha;
 834     }
 835   }
 836   state.xpos = xpos;
 837   state.ypos = ypos;
 838 }
 839
 840 void IntegratingVelocityTrackerStrategy::PopulateEstimator(
 841     const State& state,
 842     Estimator* out_estimator) const {
 843   out_estimator->time = state.update_time;
 844   out_estimator->confidence = 1.0f;
 845   out_estimator->degree = state.degree;
 846   out_estimator->xcoeff[0] = state.xpos;
 847   out_estimator->xcoeff[1] = state.xvel;
 848   out_estimator->xcoeff[2] = state.xaccel / 2;
 849   out_estimator->ycoeff[0] = state.ypos;
 850   out_estimator->ycoeff[1] = state.yvel;
 851   out_estimator->ycoeff[2] = state.yaccel / 2;
 852 }
 853
 854 }  // namespace ui