1 /*
2 * Copyright (C) 2012 The Android Open Source Project
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16
17 #define LOG_TAG "VelocityTracker"
18 //#define LOG_NDEBUG 0
19
20 // Log debug messages about velocity tracking.
21 #define DEBUG_VELOCITY 0
22
23 // Log debug messages about the progress of the algorithm itself.
24 #define DEBUG_STRATEGY 0
25
26 #include <math.h>
27 #include <limits.h>
28
29 #include <cutils/properties.h>
30 #include <input/VelocityTracker.h>
31 #include <utils/BitSet.h>
32 #include <utils/String8.h>
33 #include <utils/Timers.h>
34
35 namespace android {
36
37 // Nanoseconds per milliseconds.
38 static const nsecs_t NANOS_PER_MS = 1000000;
39
40 // Threshold for determining that a pointer has stopped moving.
41 // Some input devices do not send ACTION_MOVE events in the case where a pointer has
42 // stopped. We need to detect this case so that we can accurately predict the
43 // velocity after the pointer starts moving again.
44 static const nsecs_t ASSUME_POINTER_STOPPED_TIME = 40 * NANOS_PER_MS;
45
46
vectorDot(const float * a,const float * b,uint32_t m)47 static float vectorDot(const float* a, const float* b, uint32_t m) {
48 float r = 0;
49 while (m--) {
50 r += *(a++) * *(b++);
51 }
52 return r;
53 }
54
vectorNorm(const float * a,uint32_t m)55 static float vectorNorm(const float* a, uint32_t m) {
56 float r = 0;
57 while (m--) {
58 float t = *(a++);
59 r += t * t;
60 }
61 return sqrtf(r);
62 }
63
64 #if DEBUG_STRATEGY || DEBUG_VELOCITY
vectorToString(const float * a,uint32_t m)65 static String8 vectorToString(const float* a, uint32_t m) {
66 String8 str;
67 str.append("[");
68 while (m--) {
69 str.appendFormat(" %f", *(a++));
70 if (m) {
71 str.append(",");
72 }
73 }
74 str.append(" ]");
75 return str;
76 }
77
matrixToString(const float * a,uint32_t m,uint32_t n,bool rowMajor)78 static String8 matrixToString(const float* a, uint32_t m, uint32_t n, bool rowMajor) {
79 String8 str;
80 str.append("[");
81 for (size_t i = 0; i < m; i++) {
82 if (i) {
83 str.append(",");
84 }
85 str.append(" [");
86 for (size_t j = 0; j < n; j++) {
87 if (j) {
88 str.append(",");
89 }
90 str.appendFormat(" %f", a[rowMajor ? i * n + j : j * m + i]);
91 }
92 str.append(" ]");
93 }
94 str.append(" ]");
95 return str;
96 }
97 #endif
98
99
100 // --- VelocityTracker ---
101
102 // The default velocity tracker strategy.
103 // Although other strategies are available for testing and comparison purposes,
104 // this is the strategy that applications will actually use. Be very careful
105 // when adjusting the default strategy because it can dramatically affect
106 // (often in a bad way) the user experience.
107 const char* VelocityTracker::DEFAULT_STRATEGY = "lsq2";
108
VelocityTracker(const char * strategy)109 VelocityTracker::VelocityTracker(const char* strategy) :
110 mLastEventTime(0), mCurrentPointerIdBits(0), mActivePointerId(-1) {
111 char value[PROPERTY_VALUE_MAX];
112
113 // Allow the default strategy to be overridden using a system property for debugging.
114 if (!strategy) {
115 int length = property_get("debug.velocitytracker.strategy", value, NULL);
116 if (length > 0) {
117 strategy = value;
118 } else {
119 strategy = DEFAULT_STRATEGY;
120 }
121 }
122
123 // Configure the strategy.
124 if (!configureStrategy(strategy)) {
125 ALOGD("Unrecognized velocity tracker strategy name '%s'.", strategy);
126 if (!configureStrategy(DEFAULT_STRATEGY)) {
127 LOG_ALWAYS_FATAL("Could not create the default velocity tracker strategy '%s'!",
128 strategy);
129 }
130 }
131 }
132
~VelocityTracker()133 VelocityTracker::~VelocityTracker() {
134 delete mStrategy;
135 }
136
configureStrategy(const char * strategy)137 bool VelocityTracker::configureStrategy(const char* strategy) {
138 mStrategy = createStrategy(strategy);
139 return mStrategy != NULL;
140 }
141
createStrategy(const char * strategy)142 VelocityTrackerStrategy* VelocityTracker::createStrategy(const char* strategy) {
143 if (!strcmp("lsq1", strategy)) {
144 // 1st order least squares. Quality: POOR.
145 // Frequently underfits the touch data especially when the finger accelerates
146 // or changes direction. Often underestimates velocity. The direction
147 // is overly influenced by historical touch points.
148 return new LeastSquaresVelocityTrackerStrategy(1);
149 }
150 if (!strcmp("lsq2", strategy)) {
151 // 2nd order least squares. Quality: VERY GOOD.
152 // Pretty much ideal, but can be confused by certain kinds of touch data,
153 // particularly if the panel has a tendency to generate delayed,
154 // duplicate or jittery touch coordinates when the finger is released.
155 return new LeastSquaresVelocityTrackerStrategy(2);
156 }
157 if (!strcmp("lsq3", strategy)) {
158 // 3rd order least squares. Quality: UNUSABLE.
159 // Frequently overfits the touch data yielding wildly divergent estimates
160 // of the velocity when the finger is released.
161 return new LeastSquaresVelocityTrackerStrategy(3);
162 }
163 if (!strcmp("wlsq2-delta", strategy)) {
164 // 2nd order weighted least squares, delta weighting. Quality: EXPERIMENTAL
165 return new LeastSquaresVelocityTrackerStrategy(2,
166 LeastSquaresVelocityTrackerStrategy::WEIGHTING_DELTA);
167 }
168 if (!strcmp("wlsq2-central", strategy)) {
169 // 2nd order weighted least squares, central weighting. Quality: EXPERIMENTAL
170 return new LeastSquaresVelocityTrackerStrategy(2,
171 LeastSquaresVelocityTrackerStrategy::WEIGHTING_CENTRAL);
172 }
173 if (!strcmp("wlsq2-recent", strategy)) {
174 // 2nd order weighted least squares, recent weighting. Quality: EXPERIMENTAL
175 return new LeastSquaresVelocityTrackerStrategy(2,
176 LeastSquaresVelocityTrackerStrategy::WEIGHTING_RECENT);
177 }
178 if (!strcmp("int1", strategy)) {
179 // 1st order integrating filter. Quality: GOOD.
180 // Not as good as 'lsq2' because it cannot estimate acceleration but it is
181 // more tolerant of errors. Like 'lsq1', this strategy tends to underestimate
182 // the velocity of a fling but this strategy tends to respond to changes in
183 // direction more quickly and accurately.
184 return new IntegratingVelocityTrackerStrategy(1);
185 }
186 if (!strcmp("int2", strategy)) {
187 // 2nd order integrating filter. Quality: EXPERIMENTAL.
188 // For comparison purposes only. Unlike 'int1' this strategy can compensate
189 // for acceleration but it typically overestimates the effect.
190 return new IntegratingVelocityTrackerStrategy(2);
191 }
192 if (!strcmp("legacy", strategy)) {
193 // Legacy velocity tracker algorithm. Quality: POOR.
194 // For comparison purposes only. This algorithm is strongly influenced by
195 // old data points, consistently underestimates velocity and takes a very long
196 // time to adjust to changes in direction.
197 return new LegacyVelocityTrackerStrategy();
198 }
199 return NULL;
200 }
201
clear()202 void VelocityTracker::clear() {
203 mCurrentPointerIdBits.clear();
204 mActivePointerId = -1;
205
206 mStrategy->clear();
207 }
208
clearPointers(BitSet32 idBits)209 void VelocityTracker::clearPointers(BitSet32 idBits) {
210 BitSet32 remainingIdBits(mCurrentPointerIdBits.value & ~idBits.value);
211 mCurrentPointerIdBits = remainingIdBits;
212
213 if (mActivePointerId >= 0 && idBits.hasBit(mActivePointerId)) {
214 mActivePointerId = !remainingIdBits.isEmpty() ? remainingIdBits.firstMarkedBit() : -1;
215 }
216
217 mStrategy->clearPointers(idBits);
218 }
219
addMovement(nsecs_t eventTime,BitSet32 idBits,const Position * positions)220 void VelocityTracker::addMovement(nsecs_t eventTime, BitSet32 idBits, const Position* positions) {
221 while (idBits.count() > MAX_POINTERS) {
222 idBits.clearLastMarkedBit();
223 }
224
225 if ((mCurrentPointerIdBits.value & idBits.value)
226 && eventTime >= mLastEventTime + ASSUME_POINTER_STOPPED_TIME) {
227 #if DEBUG_VELOCITY
228 ALOGD("VelocityTracker: stopped for %0.3f ms, clearing state.",
229 (eventTime - mLastEventTime) * 0.000001f);
230 #endif
231 // We have not received any movements for too long. Assume that all pointers
232 // have stopped.
233 mStrategy->clear();
234 }
235 mLastEventTime = eventTime;
236
237 mCurrentPointerIdBits = idBits;
238 if (mActivePointerId < 0 || !idBits.hasBit(mActivePointerId)) {
239 mActivePointerId = idBits.isEmpty() ? -1 : idBits.firstMarkedBit();
240 }
241
242 mStrategy->addMovement(eventTime, idBits, positions);
243
244 #if DEBUG_VELOCITY
245 ALOGD("VelocityTracker: addMovement eventTime=%lld, idBits=0x%08x, activePointerId=%d",
246 eventTime, idBits.value, mActivePointerId);
247 for (BitSet32 iterBits(idBits); !iterBits.isEmpty(); ) {
248 uint32_t id = iterBits.firstMarkedBit();
249 uint32_t index = idBits.getIndexOfBit(id);
250 iterBits.clearBit(id);
251 Estimator estimator;
252 getEstimator(id, &estimator);
253 ALOGD(" %d: position (%0.3f, %0.3f), "
254 "estimator (degree=%d, xCoeff=%s, yCoeff=%s, confidence=%f)",
255 id, positions[index].x, positions[index].y,
256 int(estimator.degree),
257 vectorToString(estimator.xCoeff, estimator.degree + 1).string(),
258 vectorToString(estimator.yCoeff, estimator.degree + 1).string(),
259 estimator.confidence);
260 }
261 #endif
262 }
263
addMovement(const MotionEvent * event)264 void VelocityTracker::addMovement(const MotionEvent* event) {
265 int32_t actionMasked = event->getActionMasked();
266
267 switch (actionMasked) {
268 case AMOTION_EVENT_ACTION_DOWN:
269 case AMOTION_EVENT_ACTION_HOVER_ENTER:
270 // Clear all pointers on down before adding the new movement.
271 clear();
272 break;
273 case AMOTION_EVENT_ACTION_POINTER_DOWN: {
274 // Start a new movement trace for a pointer that just went down.
275 // We do this on down instead of on up because the client may want to query the
276 // final velocity for a pointer that just went up.
277 BitSet32 downIdBits;
278 downIdBits.markBit(event->getPointerId(event->getActionIndex()));
279 clearPointers(downIdBits);
280 break;
281 }
282 case AMOTION_EVENT_ACTION_MOVE:
283 case AMOTION_EVENT_ACTION_HOVER_MOVE:
284 break;
285 default:
286 // Ignore all other actions because they do not convey any new information about
287 // pointer movement. We also want to preserve the last known velocity of the pointers.
288 // Note that ACTION_UP and ACTION_POINTER_UP always report the last known position
289 // of the pointers that went up. ACTION_POINTER_UP does include the new position of
290 // pointers that remained down but we will also receive an ACTION_MOVE with this
291 // information if any of them actually moved. Since we don't know how many pointers
292 // will be going up at once it makes sense to just wait for the following ACTION_MOVE
293 // before adding the movement.
294 return;
295 }
296
297 size_t pointerCount = event->getPointerCount();
298 if (pointerCount > MAX_POINTERS) {
299 pointerCount = MAX_POINTERS;
300 }
301
302 BitSet32 idBits;
303 for (size_t i = 0; i < pointerCount; i++) {
304 idBits.markBit(event->getPointerId(i));
305 }
306
307 uint32_t pointerIndex[MAX_POINTERS];
308 for (size_t i = 0; i < pointerCount; i++) {
309 pointerIndex[i] = idBits.getIndexOfBit(event->getPointerId(i));
310 }
311
312 nsecs_t eventTime;
313 Position positions[pointerCount];
314
315 size_t historySize = event->getHistorySize();
316 for (size_t h = 0; h < historySize; h++) {
317 eventTime = event->getHistoricalEventTime(h);
318 for (size_t i = 0; i < pointerCount; i++) {
319 uint32_t index = pointerIndex[i];
320 positions[index].x = event->getHistoricalX(i, h);
321 positions[index].y = event->getHistoricalY(i, h);
322 }
323 addMovement(eventTime, idBits, positions);
324 }
325
326 eventTime = event->getEventTime();
327 for (size_t i = 0; i < pointerCount; i++) {
328 uint32_t index = pointerIndex[i];
329 positions[index].x = event->getX(i);
330 positions[index].y = event->getY(i);
331 }
332 addMovement(eventTime, idBits, positions);
333 }
334
getVelocity(uint32_t id,float * outVx,float * outVy) const335 bool VelocityTracker::getVelocity(uint32_t id, float* outVx, float* outVy) const {
336 Estimator estimator;
337 if (getEstimator(id, &estimator) && estimator.degree >= 1) {
338 *outVx = estimator.xCoeff[1];
339 *outVy = estimator.yCoeff[1];
340 return true;
341 }
342 *outVx = 0;
343 *outVy = 0;
344 return false;
345 }
346
getEstimator(uint32_t id,Estimator * outEstimator) const347 bool VelocityTracker::getEstimator(uint32_t id, Estimator* outEstimator) const {
348 return mStrategy->getEstimator(id, outEstimator);
349 }
350
351
352 // --- LeastSquaresVelocityTrackerStrategy ---
353
354 const nsecs_t LeastSquaresVelocityTrackerStrategy::HORIZON;
355 const uint32_t LeastSquaresVelocityTrackerStrategy::HISTORY_SIZE;
356
LeastSquaresVelocityTrackerStrategy(uint32_t degree,Weighting weighting)357 LeastSquaresVelocityTrackerStrategy::LeastSquaresVelocityTrackerStrategy(
358 uint32_t degree, Weighting weighting) :
359 mDegree(degree), mWeighting(weighting) {
360 clear();
361 }
362
~LeastSquaresVelocityTrackerStrategy()363 LeastSquaresVelocityTrackerStrategy::~LeastSquaresVelocityTrackerStrategy() {
364 }
365
clear()366 void LeastSquaresVelocityTrackerStrategy::clear() {
367 mIndex = 0;
368 mMovements[0].idBits.clear();
369 }
370
clearPointers(BitSet32 idBits)371 void LeastSquaresVelocityTrackerStrategy::clearPointers(BitSet32 idBits) {
372 BitSet32 remainingIdBits(mMovements[mIndex].idBits.value & ~idBits.value);
373 mMovements[mIndex].idBits = remainingIdBits;
374 }
375
addMovement(nsecs_t eventTime,BitSet32 idBits,const VelocityTracker::Position * positions)376 void LeastSquaresVelocityTrackerStrategy::addMovement(nsecs_t eventTime, BitSet32 idBits,
377 const VelocityTracker::Position* positions) {
378 if (++mIndex == HISTORY_SIZE) {
379 mIndex = 0;
380 }
381
382 Movement& movement = mMovements[mIndex];
383 movement.eventTime = eventTime;
384 movement.idBits = idBits;
385 uint32_t count = idBits.count();
386 for (uint32_t i = 0; i < count; i++) {
387 movement.positions[i] = positions[i];
388 }
389 }
390
391 /**
392 * Solves a linear least squares problem to obtain a N degree polynomial that fits
393 * the specified input data as nearly as possible.
394 *
395 * Returns true if a solution is found, false otherwise.
396 *
397 * The input consists of two vectors of data points X and Y with indices 0..m-1
398 * along with a weight vector W of the same size.
399 *
400 * The output is a vector B with indices 0..n that describes a polynomial
401 * that fits the data, such the sum of W[i] * W[i] * abs(Y[i] - (B[0] + B[1] X[i]
402 * + B[2] X[i]^2 ... B[n] X[i]^n)) for all i between 0 and m-1 is minimized.
403 *
404 * Accordingly, the weight vector W should be initialized by the caller with the
405 * reciprocal square root of the variance of the error in each input data point.
406 * In other words, an ideal choice for W would be W[i] = 1 / var(Y[i]) = 1 / stddev(Y[i]).
407 * The weights express the relative importance of each data point. If the weights are
408 * all 1, then the data points are considered to be of equal importance when fitting
409 * the polynomial. It is a good idea to choose weights that diminish the importance
410 * of data points that may have higher than usual error margins.
411 *
412 * Errors among data points are assumed to be independent. W is represented here
413 * as a vector although in the literature it is typically taken to be a diagonal matrix.
414 *
415 * That is to say, the function that generated the input data can be approximated
416 * by y(x) ~= B[0] + B[1] x + B[2] x^2 + ... + B[n] x^n.
417 *
418 * The coefficient of determination (R^2) is also returned to describe the goodness
419 * of fit of the model for the given data. It is a value between 0 and 1, where 1
420 * indicates perfect correspondence.
421 *
422 * This function first expands the X vector to a m by n matrix A such that
423 * A[i][0] = 1, A[i][1] = X[i], A[i][2] = X[i]^2, ..., A[i][n] = X[i]^n, then
424 * multiplies it by w[i]./
425 *
426 * Then it calculates the QR decomposition of A yielding an m by m orthonormal matrix Q
427 * and an m by n upper triangular matrix R. Because R is upper triangular (lower
428 * part is all zeroes), we can simplify the decomposition into an m by n matrix
429 * Q1 and a n by n matrix R1 such that A = Q1 R1.
430 *
431 * Finally we solve the system of linear equations given by R1 B = (Qtranspose W Y)
432 * to find B.
433 *
434 * For efficiency, we lay out A and Q column-wise in memory because we frequently
435 * operate on the column vectors. Conversely, we lay out R row-wise.
436 *
437 * http://en.wikipedia.org/wiki/Numerical_methods_for_linear_least_squares
438 * http://en.wikipedia.org/wiki/Gram-Schmidt
439 */
solveLeastSquares(const float * x,const float * y,const float * w,uint32_t m,uint32_t n,float * outB,float * outDet)440 static bool solveLeastSquares(const float* x, const float* y,
441 const float* w, uint32_t m, uint32_t n, float* outB, float* outDet) {
442 #if DEBUG_STRATEGY
443 ALOGD("solveLeastSquares: m=%d, n=%d, x=%s, y=%s, w=%s", int(m), int(n),
444 vectorToString(x, m).string(), vectorToString(y, m).string(),
445 vectorToString(w, m).string());
446 #endif
447
448 // Expand the X vector to a matrix A, pre-multiplied by the weights.
449 float a[n][m]; // column-major order
450 for (uint32_t h = 0; h < m; h++) {
451 a[0][h] = w[h];
452 for (uint32_t i = 1; i < n; i++) {
453 a[i][h] = a[i - 1][h] * x[h];
454 }
455 }
456 #if DEBUG_STRATEGY
457 ALOGD(" - a=%s", matrixToString(&a[0][0], m, n, false /*rowMajor*/).string());
458 #endif
459
460 // Apply the Gram-Schmidt process to A to obtain its QR decomposition.
461 float q[n][m]; // orthonormal basis, column-major order
462 float r[n][n]; // upper triangular matrix, row-major order
463 for (uint32_t j = 0; j < n; j++) {
464 for (uint32_t h = 0; h < m; h++) {
465 q[j][h] = a[j][h];
466 }
467 for (uint32_t i = 0; i < j; i++) {
468 float dot = vectorDot(&q[j][0], &q[i][0], m);
469 for (uint32_t h = 0; h < m; h++) {
470 q[j][h] -= dot * q[i][h];
471 }
472 }
473
474 float norm = vectorNorm(&q[j][0], m);
475 if (norm < 0.000001f) {
476 // vectors are linearly dependent or zero so no solution
477 #if DEBUG_STRATEGY
478 ALOGD(" - no solution, norm=%f", norm);
479 #endif
480 return false;
481 }
482
483 float invNorm = 1.0f / norm;
484 for (uint32_t h = 0; h < m; h++) {
485 q[j][h] *= invNorm;
486 }
487 for (uint32_t i = 0; i < n; i++) {
488 r[j][i] = i < j ? 0 : vectorDot(&q[j][0], &a[i][0], m);
489 }
490 }
491 #if DEBUG_STRATEGY
492 ALOGD(" - q=%s", matrixToString(&q[0][0], m, n, false /*rowMajor*/).string());
493 ALOGD(" - r=%s", matrixToString(&r[0][0], n, n, true /*rowMajor*/).string());
494
495 // calculate QR, if we factored A correctly then QR should equal A
496 float qr[n][m];
497 for (uint32_t h = 0; h < m; h++) {
498 for (uint32_t i = 0; i < n; i++) {
499 qr[i][h] = 0;
500 for (uint32_t j = 0; j < n; j++) {
501 qr[i][h] += q[j][h] * r[j][i];
502 }
503 }
504 }
505 ALOGD(" - qr=%s", matrixToString(&qr[0][0], m, n, false /*rowMajor*/).string());
506 #endif
507
508 // Solve R B = Qt W Y to find B. This is easy because R is upper triangular.
509 // We just work from bottom-right to top-left calculating B's coefficients.
510 float wy[m];
511 for (uint32_t h = 0; h < m; h++) {
512 wy[h] = y[h] * w[h];
513 }
514 for (uint32_t i = n; i-- != 0; ) {
515 outB[i] = vectorDot(&q[i][0], wy, m);
516 for (uint32_t j = n - 1; j > i; j--) {
517 outB[i] -= r[i][j] * outB[j];
518 }
519 outB[i] /= r[i][i];
520 }
521 #if DEBUG_STRATEGY
522 ALOGD(" - b=%s", vectorToString(outB, n).string());
523 #endif
524
525 // Calculate the coefficient of determination as 1 - (SSerr / SStot) where
526 // SSerr is the residual sum of squares (variance of the error),
527 // and SStot is the total sum of squares (variance of the data) where each
528 // has been weighted.
529 float ymean = 0;
530 for (uint32_t h = 0; h < m; h++) {
531 ymean += y[h];
532 }
533 ymean /= m;
534
535 float sserr = 0;
536 float sstot = 0;
537 for (uint32_t h = 0; h < m; h++) {
538 float err = y[h] - outB[0];
539 float term = 1;
540 for (uint32_t i = 1; i < n; i++) {
541 term *= x[h];
542 err -= term * outB[i];
543 }
544 sserr += w[h] * w[h] * err * err;
545 float var = y[h] - ymean;
546 sstot += w[h] * w[h] * var * var;
547 }
548 *outDet = sstot > 0.000001f ? 1.0f - (sserr / sstot) : 1;
549 #if DEBUG_STRATEGY
550 ALOGD(" - sserr=%f", sserr);
551 ALOGD(" - sstot=%f", sstot);
552 ALOGD(" - det=%f", *outDet);
553 #endif
554 return true;
555 }
556
getEstimator(uint32_t id,VelocityTracker::Estimator * outEstimator) const557 bool LeastSquaresVelocityTrackerStrategy::getEstimator(uint32_t id,
558 VelocityTracker::Estimator* outEstimator) const {
559 outEstimator->clear();
560
561 // Iterate over movement samples in reverse time order and collect samples.
562 float x[HISTORY_SIZE];
563 float y[HISTORY_SIZE];
564 float w[HISTORY_SIZE];
565 float time[HISTORY_SIZE];
566 uint32_t m = 0;
567 uint32_t index = mIndex;
568 const Movement& newestMovement = mMovements[mIndex];
569 do {
570 const Movement& movement = mMovements[index];
571 if (!movement.idBits.hasBit(id)) {
572 break;
573 }
574
575 nsecs_t age = newestMovement.eventTime - movement.eventTime;
576 if (age > HORIZON) {
577 break;
578 }
579
580 const VelocityTracker::Position& position = movement.getPosition(id);
581 x[m] = position.x;
582 y[m] = position.y;
583 w[m] = chooseWeight(index);
584 time[m] = -age * 0.000000001f;
585 index = (index == 0 ? HISTORY_SIZE : index) - 1;
586 } while (++m < HISTORY_SIZE);
587
588 if (m == 0) {
589 return false; // no data
590 }
591
592 // Calculate a least squares polynomial fit.
593 uint32_t degree = mDegree;
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, outEstimator->xCoeff, &xdet)
601 && solveLeastSquares(time, y, w, m, n, outEstimator->yCoeff, &ydet)) {
602 outEstimator->time = newestMovement.eventTime;
603 outEstimator->degree = degree;
604 outEstimator->confidence = xdet * ydet;
605 #if DEBUG_STRATEGY
606 ALOGD("estimate: degree=%d, xCoeff=%s, yCoeff=%s, confidence=%f",
607 int(outEstimator->degree),
608 vectorToString(outEstimator->xCoeff, n).string(),
609 vectorToString(outEstimator->yCoeff, n).string(),
610 outEstimator->confidence);
611 #endif
612 return true;
613 }
614 }
615
616 // No velocity data available for this pointer, but we do have its current position.
617 outEstimator->xCoeff[0] = x[0];
618 outEstimator->yCoeff[0] = y[0];
619 outEstimator->time = newestMovement.eventTime;
620 outEstimator->degree = 0;
621 outEstimator->confidence = 1;
622 return true;
623 }
624
chooseWeight(uint32_t index) const625 float LeastSquaresVelocityTrackerStrategy::chooseWeight(uint32_t index) const {
626 switch (mWeighting) {
627 case WEIGHTING_DELTA: {
628 // Weight points based on how much time elapsed between them and the next
629 // point so that points that "cover" a shorter time span are weighed less.
630 // delta 0ms: 0.5
631 // delta 10ms: 1.0
632 if (index == mIndex) {
633 return 1.0f;
634 }
635 uint32_t nextIndex = (index + 1) % HISTORY_SIZE;
636 float deltaMillis = (mMovements[nextIndex].eventTime- mMovements[index].eventTime)
637 * 0.000001f;
638 if (deltaMillis < 0) {
639 return 0.5f;
640 }
641 if (deltaMillis < 10) {
642 return 0.5f + deltaMillis * 0.05;
643 }
644 return 1.0f;
645 }
646
647 case WEIGHTING_CENTRAL: {
648 // Weight points based on their age, weighing very recent and very old points less.
649 // age 0ms: 0.5
650 // age 10ms: 1.0
651 // age 50ms: 1.0
652 // age 60ms: 0.5
653 float ageMillis = (mMovements[mIndex].eventTime - mMovements[index].eventTime)
654 * 0.000001f;
655 if (ageMillis < 0) {
656 return 0.5f;
657 }
658 if (ageMillis < 10) {
659 return 0.5f + ageMillis * 0.05;
660 }
661 if (ageMillis < 50) {
662 return 1.0f;
663 }
664 if (ageMillis < 60) {
665 return 0.5f + (60 - ageMillis) * 0.05;
666 }
667 return 0.5f;
668 }
669
670 case WEIGHTING_RECENT: {
671 // Weight points based on their age, weighing older points less.
672 // age 0ms: 1.0
673 // age 50ms: 1.0
674 // age 100ms: 0.5
675 float ageMillis = (mMovements[mIndex].eventTime - mMovements[index].eventTime)
676 * 0.000001f;
677 if (ageMillis < 50) {
678 return 1.0f;
679 }
680 if (ageMillis < 100) {
681 return 0.5f + (100 - ageMillis) * 0.01f;
682 }
683 return 0.5f;
684 }
685
686 case WEIGHTING_NONE:
687 default:
688 return 1.0f;
689 }
690 }
691
692
693 // --- IntegratingVelocityTrackerStrategy ---
694
IntegratingVelocityTrackerStrategy(uint32_t degree)695 IntegratingVelocityTrackerStrategy::IntegratingVelocityTrackerStrategy(uint32_t degree) :
696 mDegree(degree) {
697 }
698
~IntegratingVelocityTrackerStrategy()699 IntegratingVelocityTrackerStrategy::~IntegratingVelocityTrackerStrategy() {
700 }
701
clear()702 void IntegratingVelocityTrackerStrategy::clear() {
703 mPointerIdBits.clear();
704 }
705
clearPointers(BitSet32 idBits)706 void IntegratingVelocityTrackerStrategy::clearPointers(BitSet32 idBits) {
707 mPointerIdBits.value &= ~idBits.value;
708 }
709
addMovement(nsecs_t eventTime,BitSet32 idBits,const VelocityTracker::Position * positions)710 void IntegratingVelocityTrackerStrategy::addMovement(nsecs_t eventTime, BitSet32 idBits,
711 const VelocityTracker::Position* positions) {
712 uint32_t index = 0;
713 for (BitSet32 iterIdBits(idBits); !iterIdBits.isEmpty();) {
714 uint32_t id = iterIdBits.clearFirstMarkedBit();
715 State& state = mPointerState[id];
716 const VelocityTracker::Position& position = positions[index++];
717 if (mPointerIdBits.hasBit(id)) {
718 updateState(state, eventTime, position.x, position.y);
719 } else {
720 initState(state, eventTime, position.x, position.y);
721 }
722 }
723
724 mPointerIdBits = idBits;
725 }
726
getEstimator(uint32_t id,VelocityTracker::Estimator * outEstimator) const727 bool IntegratingVelocityTrackerStrategy::getEstimator(uint32_t id,
728 VelocityTracker::Estimator* outEstimator) const {
729 outEstimator->clear();
730
731 if (mPointerIdBits.hasBit(id)) {
732 const State& state = mPointerState[id];
733 populateEstimator(state, outEstimator);
734 return true;
735 }
736
737 return false;
738 }
739
initState(State & state,nsecs_t eventTime,float xpos,float ypos) const740 void IntegratingVelocityTrackerStrategy::initState(State& state,
741 nsecs_t eventTime, float xpos, float ypos) const {
742 state.updateTime = eventTime;
743 state.degree = 0;
744
745 state.xpos = xpos;
746 state.xvel = 0;
747 state.xaccel = 0;
748 state.ypos = ypos;
749 state.yvel = 0;
750 state.yaccel = 0;
751 }
752
updateState(State & state,nsecs_t eventTime,float xpos,float ypos) const753 void IntegratingVelocityTrackerStrategy::updateState(State& state,
754 nsecs_t eventTime, float xpos, float ypos) const {
755 const nsecs_t MIN_TIME_DELTA = 2 * NANOS_PER_MS;
756 const float FILTER_TIME_CONSTANT = 0.010f; // 10 milliseconds
757
758 if (eventTime <= state.updateTime + MIN_TIME_DELTA) {
759 return;
760 }
761
762 float dt = (eventTime - state.updateTime) * 0.000000001f;
763 state.updateTime = eventTime;
764
765 float xvel = (xpos - state.xpos) / dt;
766 float yvel = (ypos - state.ypos) / dt;
767 if (state.degree == 0) {
768 state.xvel = xvel;
769 state.yvel = yvel;
770 state.degree = 1;
771 } else {
772 float alpha = dt / (FILTER_TIME_CONSTANT + dt);
773 if (mDegree == 1) {
774 state.xvel += (xvel - state.xvel) * alpha;
775 state.yvel += (yvel - state.yvel) * alpha;
776 } else {
777 float xaccel = (xvel - state.xvel) / dt;
778 float yaccel = (yvel - state.yvel) / dt;
779 if (state.degree == 1) {
780 state.xaccel = xaccel;
781 state.yaccel = yaccel;
782 state.degree = 2;
783 } else {
784 state.xaccel += (xaccel - state.xaccel) * alpha;
785 state.yaccel += (yaccel - state.yaccel) * alpha;
786 }
787 state.xvel += (state.xaccel * dt) * alpha;
788 state.yvel += (state.yaccel * dt) * alpha;
789 }
790 }
791 state.xpos = xpos;
792 state.ypos = ypos;
793 }
794
populateEstimator(const State & state,VelocityTracker::Estimator * outEstimator) const795 void IntegratingVelocityTrackerStrategy::populateEstimator(const State& state,
796 VelocityTracker::Estimator* outEstimator) const {
797 outEstimator->time = state.updateTime;
798 outEstimator->confidence = 1.0f;
799 outEstimator->degree = state.degree;
800 outEstimator->xCoeff[0] = state.xpos;
801 outEstimator->xCoeff[1] = state.xvel;
802 outEstimator->xCoeff[2] = state.xaccel / 2;
803 outEstimator->yCoeff[0] = state.ypos;
804 outEstimator->yCoeff[1] = state.yvel;
805 outEstimator->yCoeff[2] = state.yaccel / 2;
806 }
807
808
809 // --- LegacyVelocityTrackerStrategy ---
810
811 const nsecs_t LegacyVelocityTrackerStrategy::HORIZON;
812 const uint32_t LegacyVelocityTrackerStrategy::HISTORY_SIZE;
813 const nsecs_t LegacyVelocityTrackerStrategy::MIN_DURATION;
814
LegacyVelocityTrackerStrategy()815 LegacyVelocityTrackerStrategy::LegacyVelocityTrackerStrategy() {
816 clear();
817 }
818
~LegacyVelocityTrackerStrategy()819 LegacyVelocityTrackerStrategy::~LegacyVelocityTrackerStrategy() {
820 }
821
clear()822 void LegacyVelocityTrackerStrategy::clear() {
823 mIndex = 0;
824 mMovements[0].idBits.clear();
825 }
826
clearPointers(BitSet32 idBits)827 void LegacyVelocityTrackerStrategy::clearPointers(BitSet32 idBits) {
828 BitSet32 remainingIdBits(mMovements[mIndex].idBits.value & ~idBits.value);
829 mMovements[mIndex].idBits = remainingIdBits;
830 }
831
addMovement(nsecs_t eventTime,BitSet32 idBits,const VelocityTracker::Position * positions)832 void LegacyVelocityTrackerStrategy::addMovement(nsecs_t eventTime, BitSet32 idBits,
833 const VelocityTracker::Position* positions) {
834 if (++mIndex == HISTORY_SIZE) {
835 mIndex = 0;
836 }
837
838 Movement& movement = mMovements[mIndex];
839 movement.eventTime = eventTime;
840 movement.idBits = idBits;
841 uint32_t count = idBits.count();
842 for (uint32_t i = 0; i < count; i++) {
843 movement.positions[i] = positions[i];
844 }
845 }
846
getEstimator(uint32_t id,VelocityTracker::Estimator * outEstimator) const847 bool LegacyVelocityTrackerStrategy::getEstimator(uint32_t id,
848 VelocityTracker::Estimator* outEstimator) const {
849 outEstimator->clear();
850
851 const Movement& newestMovement = mMovements[mIndex];
852 if (!newestMovement.idBits.hasBit(id)) {
853 return false; // no data
854 }
855
856 // Find the oldest sample that contains the pointer and that is not older than HORIZON.
857 nsecs_t minTime = newestMovement.eventTime - HORIZON;
858 uint32_t oldestIndex = mIndex;
859 uint32_t numTouches = 1;
860 do {
861 uint32_t nextOldestIndex = (oldestIndex == 0 ? HISTORY_SIZE : oldestIndex) - 1;
862 const Movement& nextOldestMovement = mMovements[nextOldestIndex];
863 if (!nextOldestMovement.idBits.hasBit(id)
864 || nextOldestMovement.eventTime < minTime) {
865 break;
866 }
867 oldestIndex = nextOldestIndex;
868 } while (++numTouches < HISTORY_SIZE);
869
870 // Calculate an exponentially weighted moving average of the velocity estimate
871 // at different points in time measured relative to the oldest sample.
872 // This is essentially an IIR filter. Newer samples are weighted more heavily
873 // than older samples. Samples at equal time points are weighted more or less
874 // equally.
875 //
876 // One tricky problem is that the sample data may be poorly conditioned.
877 // Sometimes samples arrive very close together in time which can cause us to
878 // overestimate the velocity at that time point. Most samples might be measured
879 // 16ms apart but some consecutive samples could be only 0.5sm apart because
880 // the hardware or driver reports them irregularly or in bursts.
881 float accumVx = 0;
882 float accumVy = 0;
883 uint32_t index = oldestIndex;
884 uint32_t samplesUsed = 0;
885 const Movement& oldestMovement = mMovements[oldestIndex];
886 const VelocityTracker::Position& oldestPosition = oldestMovement.getPosition(id);
887 nsecs_t lastDuration = 0;
888
889 while (numTouches-- > 1) {
890 if (++index == HISTORY_SIZE) {
891 index = 0;
892 }
893 const Movement& movement = mMovements[index];
894 nsecs_t duration = movement.eventTime - oldestMovement.eventTime;
895
896 // If the duration between samples is small, we may significantly overestimate
897 // the velocity. Consequently, we impose a minimum duration constraint on the
898 // samples that we include in the calculation.
899 if (duration >= MIN_DURATION) {
900 const VelocityTracker::Position& position = movement.getPosition(id);
901 float scale = 1000000000.0f / duration; // one over time delta in seconds
902 float vx = (position.x - oldestPosition.x) * scale;
903 float vy = (position.y - oldestPosition.y) * scale;
904 accumVx = (accumVx * lastDuration + vx * duration) / (duration + lastDuration);
905 accumVy = (accumVy * lastDuration + vy * duration) / (duration + lastDuration);
906 lastDuration = duration;
907 samplesUsed += 1;
908 }
909 }
910
911 // Report velocity.
912 const VelocityTracker::Position& newestPosition = newestMovement.getPosition(id);
913 outEstimator->time = newestMovement.eventTime;
914 outEstimator->confidence = 1;
915 outEstimator->xCoeff[0] = newestPosition.x;
916 outEstimator->yCoeff[0] = newestPosition.y;
917 if (samplesUsed) {
918 outEstimator->xCoeff[1] = accumVx;
919 outEstimator->yCoeff[1] = accumVy;
920 outEstimator->degree = 1;
921 } else {
922 outEstimator->degree = 0;
923 }
924 return true;
925 }
926
927 } // namespace android
928