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