Eagerly attempt to refine a track before doing a brute search

Before the refinement phase of tracking, a brute force SAD search
is run across the search area. This works well but is slow;
especially if the guess for the track's location is accurate.

This patch runs a refinement phase before running a brute force
search, hoping that the guessed position (in x2, y2) is close to
the best answer. If it is, then no brute search is done. If it is
not, then a normal brute force search followed by refinement is
done.

In some cases this may produce worse tracks than before; the
regressions will need investigation. The predictive motion model
(to be implemented) will reduce the probability of that happening.

3 files changed, 58 insertions, 4 deletions

diff --git a/extern/libmv/libmv-capi.cc b/extern/libmv/libmv-capi.cc
index 0585bd3e8ac..f02ca0e98e0 100644
--- a/extern/libmv/libmv-capi.cc
+++ b/extern/libmv/libmv-capi.cc
@@ -306,6 +306,9 @@ int libmv_trackRegion(const libmv_TrackRegionOptions *options,
 	track_region_options.num_extra_points = 1;
 	track_region_options.image1_mask = NULL;
 	track_region_options.use_brute_initialization = options->use_brute;
+	/* TODO(keir): This will make some cases better, but may be a regression until
+	 * the motion model is in. Since this is on trunk, enable it for now. */
+	track_region_options.attempt_refine_before_brute = true;
 	track_region_options.use_normalized_intensities = options->use_normalization;
 
 	if (options->image1_mask) {
diff --git a/extern/libmv/libmv/tracking/track_region.cc b/extern/libmv/libmv/tracking/track_region.cc
index 786f16bf219..9e0bf9a77da 100644
--- a/extern/libmv/libmv/tracking/track_region.cc
+++ b/extern/libmv/libmv/tracking/track_region.cc
@@ -1289,6 +1289,15 @@ bool BruteTranslationOnlyInitialize(const FloatImage &image1,
   return true;
 }
 
+void CopyQuad(double *src_x, double *src_y,
+              double *dst_x, double *dst_y,
+              int num_extra_points) {
+  for (int i = 0; i < 4 + num_extra_points; ++i) {
+    dst_x[i] = src_x[i];
+    dst_y[i] = src_y[i];
+  }
+}
+
 }  // namespace
 
 template<typename Warp>
@@ -1298,11 +1307,46 @@ void TemplatedTrackRegion(const FloatImage &image1,
                           const TrackRegionOptions &options,
                           double *x2, double *y2,
                           TrackRegionResult *result) {
-  for (int i = 0; i < 4; ++i) {
+  for (int i = 0; i < 4 + options.num_extra_points; ++i) {
     LG << "P" << i << ": (" << x1[i] << ", " << y1[i] << "); guess ("
        << x2[i] << ", " << y2[i] << "); (dx, dy): (" << (x2[i] - x1[i]) << ", "
        << (y2[i] - y1[i]) << ").";
   }
+
+  // Since (x2, y2) contains a prediction for where the tracked point should
+  // go, try a refinement immediately in the hope that the prediction is close
+  // enough.
+  if (options.attempt_refine_before_brute) {
+    TrackRegionOptions modified_options = options;
+    modified_options.use_brute_initialization = false;
+    modified_options.attempt_refine_before_brute = false;
+
+    double x2_first_try[5];
+    double y2_first_try[5];
+    CopyQuad(x2, y2, x2_first_try, y2_first_try, options.num_extra_points);
+
+    TemplatedTrackRegion<Warp>(image1, image2,
+                               x1, y1, modified_options,
+                               x2_first_try, y2_first_try, result);
+
+    // Of the things that can happen in the first pass, don't try the brute
+    // pass (and second attempt) if the error is one of the terminations below.
+    if (result->termination == TrackRegionResult::PARAMETER_TOLERANCE ||
+        result->termination == TrackRegionResult::FUNCTION_TOLERANCE ||
+        result->termination == TrackRegionResult::GRADIENT_TOLERANCE ||
+        result->termination == TrackRegionResult::SOURCE_OUT_OF_BOUNDS ||
+        result->termination == TrackRegionResult::DESTINATION_OUT_OF_BOUNDS ||
+        result->termination == TrackRegionResult::INSUFFICIENT_PATTERN_AREA) {
+      LG << "Terminated with first try at refinement; no brute needed.";
+      // TODO(keir): Also check correlation?
+      CopyQuad(x2_first_try, y2_first_try, x2, y2, options.num_extra_points);
+      LG << "Early termination correlation: " << result->correlation;
+      return;
+    } else {
+      LG << "Initial eager-refinement failed; retrying normally.";
+    }
+  }
+
   if (options.use_normalized_intensities) {
     LG << "Using normalized intensities.";
   }
@@ -1368,7 +1412,6 @@ void TemplatedTrackRegion(const FloatImage &image1,
   int num_samples_y;
   PickSampling(x1, y1, x2, y2, &num_samples_x, &num_samples_y);
 
-
   // Compute the warp from rectangular coordinates.
   Mat3 canonical_homography = ComputeCanonicalHomography(x1, y1,
                                                          num_samples_x,
diff --git a/extern/libmv/libmv/tracking/track_region.h b/extern/libmv/libmv/tracking/track_region.h
index 6c7218f11d1..58742cab36c 100644
--- a/extern/libmv/libmv/tracking/track_region.h
+++ b/extern/libmv/libmv/tracking/track_region.h
@@ -55,6 +55,13 @@ struct TrackRegionOptions {
   // that the nearby minima is correct, or the search area is too small.
   bool use_brute_initialization;
 
+  // If true and brute initialization is enabled, first try refining with the
+  // initial guess instead of starting with the brute initialization. If the
+  // initial refinement fails, then a normal brute search followed by
+  // refinement is attempted. If the initial refinement succeeds, then the
+  // result is returned as is (skipping a costly brute search).
+  bool attempt_refine_before_brute;
+
   // If true, normalize the image patches by their mean before doing the sum of
   // squared error calculation. This is reasonable since the effect of
   // increasing light intensity is multiplicative on the pixel intensities.
@@ -67,8 +74,9 @@ struct TrackRegionOptions {
   // take the image derivative.
   double sigma;
 
-  // Extra points that should get transformed by the warp. This is useful
-  // because the actual warp parameters are not exposed.
+  // Extra points that should get transformed by the warp. These points are
+  // appended to the x and y arrays. This is useful because the actual warp
+  // parameters are not exposed.
   int num_extra_points;
 
   // For motion models other than translation, the optimizer sometimes has