Camera tracking refactoring:

- Deduplicate patetrn sampling used in esm and lmicklt trackers and
  move SamplePattern to image/sample.h - Move computation of Pearson
  product-moment correlation into own function in new file image/correlation.h
  so all trackers can use it to check final correlation.
- Remove SAD tracker. It's almost the same as brute tracker, with only two differences:
  1. It does brute search of affine transformation which in some cases helps to track
     rotating features 2. It didn't use common tracker api which probably gave some
     speed advantage, but lead to a real headache to use it together with other
     trackers leading to duplicated code in blender side.
- Switch blenedr to use brute tracker instead of sad tracker which tracking made
  source code much more simple to follow.

11 files changed, 116 insertions, 390 deletions

diff --git a/extern/libmv/CMakeLists.txt b/extern/libmv/CMakeLists.txt
index 5239c7ff7ca..005bcdf0a3b 100644
--- a/extern/libmv/CMakeLists.txt
+++ b/extern/libmv/CMakeLists.txt
@@ -80,7 +80,6 @@ set(SRC
 	libmv/tracking/lmicklt_region_tracker.cc
 	libmv/tracking/pyramid_region_tracker.cc
 	libmv/tracking/retrack_region_tracker.cc
-	libmv/tracking/sad.cc
 	libmv/tracking/trklt_region_tracker.cc
 
 	third_party/fast/fast_10.c
@@ -104,6 +103,7 @@ set(SRC
 	libmv/image/array_nd.h
 	libmv/image/convolve.h
 	libmv/image/image.h
+	libmv/image/correlation.h
 	libmv/image/sample.h
 	libmv/image/tuple.h
 	libmv/logging/logging.h
@@ -137,7 +137,6 @@ set(SRC
 	libmv/tracking/pyramid_region_tracker.h
 	libmv/tracking/region_tracker.h
 	libmv/tracking/retrack_region_tracker.h
-	libmv/tracking/sad.h
 	libmv/tracking/trklt_region_tracker.h
 
 	third_party/fast/fast.h
diff --git a/extern/libmv/libmv-capi.cpp b/extern/libmv/libmv-capi.cpp
index 465bf519dbd..7fd4bfdd9d1 100644
--- a/extern/libmv/libmv-capi.cpp
+++ b/extern/libmv/libmv-capi.cpp
@@ -44,8 +44,6 @@
 #include "libmv/tracking/lmicklt_region_tracker.h"
 #include "libmv/tracking/pyramid_region_tracker.h"
 
-#include "libmv/tracking/sad.h"
-
 #include "libmv/simple_pipeline/callbacks.h"
 #include "libmv/simple_pipeline/tracks.h"
 #include "libmv/simple_pipeline/initialize_reconstruction.h"
@@ -137,12 +135,25 @@ libmv_RegionTracker *libmv_hybridRegionTrackerNew(int max_iterations, int half_w
 	libmv::BruteRegionTracker *brute_region_tracker = new libmv::BruteRegionTracker;
 	brute_region_tracker->half_window_size = half_window_size;
 
+	/* do not use correlation check for brute checker itself,
+	 * this check will happen in esm tracker */
+	brute_region_tracker->minimum_correlation = 0.0;
+
 	libmv::HybridRegionTracker *hybrid_region_tracker =
 		new libmv::HybridRegionTracker(brute_region_tracker, esm_region_tracker);
 
 	return (libmv_RegionTracker *)hybrid_region_tracker;
 }
 
+libmv_RegionTracker *libmv_bruteRegionTrackerNew(int half_window_size, double minimum_correlation)
+{
+	libmv::BruteRegionTracker *brute_region_tracker = new libmv::BruteRegionTracker;
+	brute_region_tracker->half_window_size = half_window_size;
+	brute_region_tracker->minimum_correlation = minimum_correlation;
+
+	return (libmv_RegionTracker *)brute_region_tracker;
+}
+
 static void floatBufToImage(const float *buf, int width, int height, libmv::FloatImage *image)
 {
 	int x, y, a = 0;
@@ -316,33 +327,6 @@ void libmv_regionTrackerDestroy(libmv_RegionTracker *libmv_tracker)
 
 /* ************ Tracks ************ */
 
-void libmv_SADSamplePattern(unsigned char *image, int stride,
-			float warp[3][2], unsigned char *pattern, int pattern_size)
-{
-	libmv::mat32 mat32;
-
-	memcpy(mat32.data, warp, sizeof(float)*3*2);
-
-	libmv::SamplePattern(image, stride, mat32, pattern, pattern_size);
-}
-
-float libmv_SADTrackerTrack(unsigned char *pattern, unsigned char *warped, int pattern_size, unsigned char *image, int stride,
-			int width, int height, float warp[3][2])
-{
-	float result;
-	libmv::mat32 mat32;
-
-	memcpy(mat32.data, warp, sizeof(float)*3*2);
-
-	result = libmv::Track(pattern, warped, pattern_size, image, stride, width, height, &mat32, 16, 16);
-
-	memcpy(warp, mat32.data, sizeof(float)*3*2);
-
-	return result;
-}
-
-/* ************ Tracks ************ */
-
 libmv_Tracks *libmv_tracksNew(void)
 {
 	libmv::Tracks *libmv_tracks = new libmv::Tracks();
diff --git a/extern/libmv/libmv-capi.h b/extern/libmv/libmv-capi.h
index e10d4ef842a..e5728188afb 100644
--- a/extern/libmv/libmv-capi.h
+++ b/extern/libmv/libmv-capi.h
@@ -45,16 +45,11 @@ void libmv_setLoggingVerbosity(int verbosity);
 /* RegionTracker */
 struct libmv_RegionTracker *libmv_pyramidRegionTrackerNew(int max_iterations, int pyramid_level, int half_window_size, double minimum_correlation);
 struct libmv_RegionTracker *libmv_hybridRegionTrackerNew(int max_iterations, int half_window_size, double minimum_correlation);
+struct libmv_RegionTracker *libmv_bruteRegionTrackerNew(int half_window_size, double minimum_correlation);
 int libmv_regionTrackerTrack(struct libmv_RegionTracker *libmv_tracker, const float *ima1, const float *ima2,
 			int width, int height, double  x1, double  y1, double *x2, double *y2);
 void libmv_regionTrackerDestroy(struct libmv_RegionTracker *libmv_tracker);
 
-/* SAD Tracker */
-void libmv_SADSamplePattern(unsigned char *image, int stride,
-			float warp[3][2], unsigned char *pattern, int pattern_size);
-float libmv_SADTrackerTrack(unsigned char *pattern, unsigned char *warped, int pattern_size, unsigned char *image,
-			int stride, int width, int height, float warp[3][2]);
-
 /* Tracks */
 struct libmv_Tracks *libmv_tracksNew(void);
 void libmv_tracksInsert(struct libmv_Tracks *libmv_tracks, int image, int track, double x, double y);
diff --git a/extern/libmv/libmv/image/correlation.h b/extern/libmv/libmv/image/correlation.h
new file mode 100644
index 00000000000..9d6aceecceb
--- /dev/null
+++ b/extern/libmv/libmv/image/correlation.h
@@ -0,0 +1,51 @@
+// Copyright (c) 2012 libmv authors.
+//
+// Permission is hereby granted, free of charge, to any person obtaining a copy
+// of this software and associated documentation files (the "Software"), to
+// deal in the Software without restriction, including without limitation the
+// rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
+// sell copies of the Software, and to permit persons to whom the Software is
+// furnished to do so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included in
+// all copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
+// IN THE SOFTWARE.
+
+#ifndef LIBMV_IMAGE_CORRELATION_H
+#define LIBMV_IMAGE_CORRELATION_H
+
+#include "libmv/image/image.h"
+
+namespace libmv {
+
+inline double PearsonProductMomentCorrelation(Array3Df image_and_gradient1_sampled,
+                                              Array3Df image_and_gradient2_sampled,
+                                              int width) {
+  double sX=0,sY=0,sXX=0,sYY=0,sXY=0;
+  for (int r = 0; r < width; ++r) {
+    for (int c = 0; c < width; ++c) {
+      double x = image_and_gradient1_sampled(r, c, 0);
+      double y = image_and_gradient2_sampled(r, c, 0);
+      sX += x;
+      sY += y;
+      sXX += x*x;
+      sYY += y*y;
+      sXY += x*y;
+    }
+  }
+  double N = width*width;
+  sX /= N, sY /= N, sXX /= N, sYY /= N, sXY /= N;
+  double correlation = (sXY-sX*sY)/sqrt(double((sXX-sX*sX)*(sYY-sY*sY)));
+  return correlation;
+}
+
+}  // namespace libmv
+
+#endif  // LIBMV_IMAGE_IMAGE_CORRELATION_H
diff --git a/extern/libmv/libmv/image/sample.h b/extern/libmv/libmv/image/sample.h
index cd361231b58..ea86edf117d 100644
--- a/extern/libmv/libmv/image/sample.h
+++ b/extern/libmv/libmv/image/sample.h
@@ -98,6 +98,24 @@ inline void DownsampleChannelsBy2(const Array3Df &in, Array3Df *out) {
 
 }
 
+// Sample a region centered at x,y in image with size extending by half_width
+// from x,y. Channels specifies the number of channels to sample from.
+inline void SamplePattern(const FloatImage &image,
+                   double x, double y,
+                   int half_width,
+                   int channels,
+                   FloatImage *sampled) {
+  sampled->Resize(2 * half_width + 1, 2 * half_width + 1, channels);
+  for (int r = -half_width; r <= half_width; ++r) {
+    for (int c = -half_width; c <= half_width; ++c) {
+      for (int i = 0; i < channels; ++i) {
+        (*sampled)(r + half_width, c + half_width, i) =
+            SampleLinear(image, y + r, x + c, i);
+      }
+    }
+  }
+}
+
 }  // namespace libmv
 
 #endif  // LIBMV_IMAGE_SAMPLE_H_
diff --git a/extern/libmv/libmv/tracking/brute_region_tracker.cc b/extern/libmv/libmv/tracking/brute_region_tracker.cc
index a416710ad2b..af7f673472e 100644
--- a/extern/libmv/libmv/tracking/brute_region_tracker.cc
+++ b/extern/libmv/libmv/tracking/brute_region_tracker.cc
@@ -44,6 +44,7 @@
 
 #include "libmv/image/image.h"
 #include "libmv/image/convolve.h"
+#include "libmv/image/correlation.h"
 #include "libmv/image/sample.h"
 #include "libmv/logging/logging.h"
 
@@ -332,6 +333,29 @@ bool BruteRegionTracker::Track(const FloatImage &image1,
   if (best_sad != INT_MAX) {
     *x2 = best_j + half_window_size;
     *y2 = best_i + half_window_size;
+
+    if (minimum_correlation > 0) {
+      Array3Df image_and_gradient1_sampled, image_and_gradient2_sampled;
+
+      SamplePattern(image_and_gradient1, x1, y1, half_window_size, 3,
+                    &image_and_gradient1_sampled);
+      SamplePattern(image_and_gradient2, *x2, *y2, half_window_size, 3,
+                    &image_and_gradient2_sampled);
+
+      // Compute the Pearson product-moment correlation coefficient to check
+      // for sanity.
+      double correlation = PearsonProductMomentCorrelation(image_and_gradient1_sampled,
+                                                           image_and_gradient2_sampled,
+                                                           pattern_width);
+      LG << "Final correlation: " << correlation;
+
+      if (correlation < minimum_correlation) {
+        LG << "Correlation " << correlation << " greater than "
+           << minimum_correlation << "; bailing.";
+        return false;
+      }
+    }
+
     return true;
   }
   return false;
diff --git a/extern/libmv/libmv/tracking/brute_region_tracker.h b/extern/libmv/libmv/tracking/brute_region_tracker.h
index 7a9b06d8a31..8ff0c1b7965 100644
--- a/extern/libmv/libmv/tracking/brute_region_tracker.h
+++ b/extern/libmv/libmv/tracking/brute_region_tracker.h
@@ -27,7 +27,9 @@
 namespace libmv {
 
 struct BruteRegionTracker : public RegionTracker {
-  BruteRegionTracker() : half_window_size(4) {}
+  BruteRegionTracker()
+      : half_window_size(4),
+      minimum_correlation(0.78) {}
   
   virtual ~BruteRegionTracker() {}
 
@@ -39,6 +41,7 @@ struct BruteRegionTracker : public RegionTracker {
 
   // No point in creating getters or setters.
   int half_window_size;
+  double minimum_correlation;
 };
 
 }  // namespace libmv
diff --git a/extern/libmv/libmv/tracking/esm_region_tracker.cc b/extern/libmv/libmv/tracking/esm_region_tracker.cc
index df9c89a46d4..e0b85f19943 100644
--- a/extern/libmv/libmv/tracking/esm_region_tracker.cc
+++ b/extern/libmv/libmv/tracking/esm_region_tracker.cc
@@ -26,6 +26,7 @@
 #include "libmv/logging/logging.h"
 #include "libmv/image/image.h"
 #include "libmv/image/convolve.h"
+#include "libmv/image/correlation.h"
 #include "libmv/image/sample.h"
 #include "libmv/numeric/numeric.h"
 
@@ -55,24 +56,6 @@ static bool RegionIsInBounds(const FloatImage &image1,
   return true;
 }
 
-// Sample a region centered at x,y in image with size extending by half_width
-// from x,y. Channels specifies the number of channels to sample from.
-static void SamplePattern(const FloatImage &image,
-                   double x, double y,
-                   int half_width,
-                   int channels,
-                   FloatImage *sampled) {
-  sampled->Resize(2 * half_width + 1, 2 * half_width + 1, channels);
-  for (int r = -half_width; r <= half_width; ++r) {
-    for (int c = -half_width; c <= half_width; ++c) {
-      for (int i = 0; i < channels; ++i) {
-        (*sampled)(r + half_width, c + half_width, i) =
-            SampleLinear(image, y + r, x + c, i);
-      }
-    }
-  }
-}
-
 // Estimate "reasonable" error by computing autocorrelation for a small shift.
 // TODO(keir): Add a facility for 
 static double EstimateReasonableError(const FloatImage &image,
@@ -276,22 +259,9 @@ bool EsmRegionTracker::Track(const FloatImage &image1,
     if (d.squaredNorm() < min_update_squared_distance) {
       // Compute the Pearson product-moment correlation coefficient to check
       // for sanity.
-      // TODO(keir): Put this somewhere smarter.
-      double sX=0,sY=0,sXX=0,sYY=0,sXY=0;
-      for (int r = 0; r < width; ++r) {
-        for (int c = 0; c < width; ++c) {
-          double x = image_and_gradient1_sampled(r, c, 0);
-          double y = image_and_gradient2_sampled[new_image](r, c, 0);
-          sX += x;
-          sY += y;
-          sXX += x*x;
-          sYY += y*y;
-          sXY += x*y;
-        }
-      }
-      double N = width*width;
-      sX /= N, sY /= N, sXX /= N, sYY /= N, sXY /= N;
-      double correlation = (sXY-sX*sY)/sqrt(double((sXX-sX*sX)*(sYY-sY*sY)));
+      double correlation = PearsonProductMomentCorrelation(image_and_gradient1_sampled,
+                                                           image_and_gradient2_sampled[new_image],
+                                                           width);
       LG << "Final correlation: " << correlation;
 
       if (correlation < minimum_correlation) {
diff --git a/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc b/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc
index c06a1d3302c..581e984b569 100644
--- a/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc
+++ b/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc
@@ -52,24 +52,6 @@ static bool RegionIsInBounds(const FloatImage &image1,
   return true;
 }
 
-// Sample a region centered at x,y in image with size extending by half_width
-// from x,y. Channels specifies the number of channels to sample from.
-static void SamplePattern(const FloatImage &image,
-                   double x, double y,
-                   int half_width,
-                   int channels,
-                   FloatImage *sampled) {
-  sampled->Resize(2 * half_width + 1, 2 * half_width + 1, channels);
-  for (int r = -half_width; r <= half_width; ++r) {
-    for (int c = -half_width; c <= half_width; ++c) {
-      for (int i = 0; i < channels; ++i) {
-        (*sampled)(r + half_width, c + half_width, i) =
-            SampleLinear(image, y + r, x + c, i);
-      }
-    }
-  }
-}
-
 // Estimate "reasonable" error by computing autocorrelation for a small shift.
 static double EstimateReasonableError(const FloatImage &image,
                                double x, double y,
diff --git a/extern/libmv/libmv/tracking/sad.cc b/extern/libmv/libmv/tracking/sad.cc
deleted file mode 100644
index 0876ef2fe73..00000000000
--- a/extern/libmv/libmv/tracking/sad.cc
+++ /dev/null
@@ -1,191 +0,0 @@
-/****************************************************************************
-**
-** Copyright (c) 2011 libmv authors.
-**
-** Permission is hereby granted, free of charge, to any person obtaining a copy
-** of this software and associated documentation files (the "Software"), to
-** deal in the Software without restriction, including without limitation the
-** rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
-** sell copies of the Software, and to permit persons to whom the Software is
-** furnished to do so, subject to the following conditions:
-**
-** The above copyright notice and this permission notice shall be included in
-** all copies or substantial portions of the Software.
-**
-** THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
-** IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
-** FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
-** AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
-** LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
-** FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
-** IN THE SOFTWARE.
-**
-****************************************************************************/
-
-#include "libmv/tracking/sad.h"
-#include <stdlib.h>
-#include <math.h>
-#include <stdio.h>
-
-namespace libmv {
-
-void LaplaceFilter(ubyte* src, ubyte* dst, int width, int height, int strength) {
-  for(int y=1; y<height-1; y++) for(int x=1; x<width-1; x++) {
-    const ubyte* s = &src[y*width+x];
-    int l = 128 +
-        s[-width-1] + s[-width] + s[-width+1] +
-        s[1]        - 8*s[0]    + s[1]        +
-        s[ width-1] + s[ width] + s[ width+1] ;
-    int d = ((256-strength)*s[0] + strength*l) / 256;
-    if(d < 0) d=0;
-    if(d > 255) d=255;
-    dst[y*width+x] = d;
-  }
-}
-
-struct vec2 {
-  float x,y;
-  inline vec2(float x, float y):x(x),y(y){}
-};
-inline vec2 operator*(mat32 m, vec2 v) {
-  return vec2(v.x*m(0,0)+v.y*m(0,1)+m(0,2),v.x*m(1,0)+v.y*m(1,1)+m(1,2));
-}
-
-//! fixed point bilinear sample with precision k
-template <int k> inline int sample(const ubyte* image,int stride, int x, int y, int u, int v) {
-  const ubyte* s = &image[y*stride+x];
-  return ((s[     0] * (k-u) + s[       1] * u) * (k-v)
-        + (s[stride] * (k-u) + s[stride+1] * u) * (  v) ) / (k*k);
-}
-
-#ifdef __SSE__
-#include <xmmintrin.h>
-int lround(float x) { return _mm_cvtss_si32(_mm_set_ss(x)); }
-#elif defined(_MSC_VER)
-int lround(float x) { return x+0.5; }
-#endif
-
-//TODO(MatthiasF): SSE optimization
-void SamplePattern(ubyte* image, int stride, mat32 warp, ubyte* pattern, int size) {
-  const int k = 256;
-  for (int i = 0; i < size; i++) for (int j = 0; j < size; j++) {
-    vec2 p = warp*vec2(j-size/2,i-size/2);
-    int fx = lround(p.x*k), fy = lround(p.y*k);
-    int ix = fx/k, iy = fy/k;
-    int u = fx%k, v = fy%k;
-    pattern[i*size+j] = sample<k>(image,stride,ix,iy,u,v);
-  }
-}
-
-#ifdef __SSE2__
-#include <emmintrin.h>
-  static uint SAD(/*const*/ ubyte* pattern, /*const*/ ubyte* image, int stride, int size) {
-  uint sad = 0;
-  __m128i a = _mm_setzero_si128();
-
-  for(int i = 0; i < size; i++) {
-    int j = 0;
-
-    for(j = 0; j < size/16; j++) {
-      if((i*size/16+j) % 32 == 0) {
-        sad += _mm_extract_epi16(a,0) + _mm_extract_epi16(a,4);
-        a = _mm_setzero_si128();
-      }
-
-      a = _mm_adds_epu16(a, _mm_sad_epu8( _mm_loadu_si128((__m128i*)(pattern+i*size+j*16)),
-                                          _mm_loadu_si128((__m128i*)(image+i*stride+j*16))));
-    }
-
-    for(j = j*16; j < size; j++) {
-      sad += abs((int)pattern[i*size+j] - image[i*stride+j]);
-    }
-  }
-
-  sad += _mm_extract_epi16(a,0) + _mm_extract_epi16(a,4);
-
-  return sad;
-}
-#else
-static uint SAD(const ubyte* pattern, const ubyte* image, int stride, int size) {
-  uint sad=0;
-  for(int i = 0; i < size; i++) {
-    for(int j = 0; j < size; j++) {
-      sad += abs((int)pattern[i*size+j] - image[i*stride+j]);
-    }
-  }
-  return sad;
-}
-#endif
-
-float sq(float x) { return x*x; }
-float Track(ubyte* reference, ubyte* warped, int size, ubyte* image, int stride, int w, int h, mat32* warp, float areaPenalty, float conditionPenalty) {
-  mat32 m=*warp;
-  uint min=-1;
-
-  // exhaustive search integer pixel translation
-  int ix = m(0,2), iy = m(1,2);
-  for(int y = size/2; y < h-size/2; y++) {
-    for(int x = size/2; x < w-size/2; x++) {
-      m(0,2) = x, m(1,2) = y;
-      uint sad = SAD(warped,&image[(y-size/2)*stride+(x-size/2)],stride,size);
-      // TODO: using chroma could help disambiguate some cases
-      if(sad < min) {
-        min = sad;
-        ix = x, iy = y;
-      }
-    }
-  }
-  m(0,2) = ix, m(1,2) = iy;
-  min=-1; //reset score since direct warped search match too well (but the wrong pattern).
-
-  // 6D coordinate descent to find affine transform
-  ubyte* match = new ubyte[size*size];
-  float step = 0.5;
-  for(int p = 0; p < 8; p++) { //foreach precision level
-    for(int i = 0; i < 2; i++) { // iterate twice per precision level
-      //TODO: other sweep pattern might converge better
-      for(int d=0; d < 6; d++) { // iterate dimension sequentially (cyclic descent)
-        for(float e = -step; e <= step; e+=step) { //solve subproblem (evaluate only along one coordinate)
-          mat32 t = m;
-          t.data[d] += e;
-          //TODO: better performance would also allow a more exhaustive search
-          SamplePattern(image,stride,t,match,size);
-          uint sad = SAD(reference,match,size,size);
-          // regularization: keep constant area and good condition
-          float area = t(0,0)*t(1,1)-t(0,1)*t(1,0);
-          float x = sq(t(0,0))+sq(t(0,1)), y = sq(t(1,0))+sq(t(1,1));
-          float condition = x>y ? x/y : y/x;
-          sad += size*size*( areaPenalty*sq(area-1) + conditionPenalty*sq(condition-1) );
-          if(sad < min) {
-            min = sad;
-            m = t;
-          }
-        }
-      }
-    }
-    step /= 2;
-  }
-  *warp = m;
-
-  // Compute Pearson product-moment correlation coefficient
-  uint sX=0,sY=0,sXX=0,sYY=0,sXY=0;
-  SamplePattern(image,stride,m,match,size);
-  SAD(reference,match,size,size);
-  for(int i = 0; i < size; i++) {
-    for(int j = 0; j < size; j++) {
-      int x = reference[i*size+j];
-      int y = match[i*size+j];
-      sX += x;
-      sY += y;
-      sXX += x*x;
-      sYY += y*y;
-      sXY += x*y;
-    }
-  }
-  delete[] match;
-  const int N = size*size;
-  sX /= N, sY /= N, sXX /= N, sYY /= N, sXY /= N;
-  return (sXY-sX*sY)/sqrt(double((sXX-sX*sX)*(sYY-sY*sY)));
-}
-
-}  // namespace libmv
diff --git a/extern/libmv/libmv/tracking/sad.h b/extern/libmv/libmv/tracking/sad.h
deleted file mode 100644
index 9fe323b74c4..00000000000
--- a/extern/libmv/libmv/tracking/sad.h
+++ /dev/null
@@ -1,109 +0,0 @@
-/****************************************************************************
-**
-** Copyright (c) 2011 libmv authors.
-**
-** Permission is hereby granted, free of charge, to any person obtaining a copy
-** of this software and associated documentation files (the "Software"), to
-** deal in the Software without restriction, including without limitation the
-** rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
-** sell copies of the Software, and to permit persons to whom the Software is
-** furnished to do so, subject to the following conditions:
-**
-** The above copyright notice and this permission notice shall be included in
-** all copies or substantial portions of the Software.
-**
-** THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
-** IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
-** FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
-** AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
-** LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
-** FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
-** IN THE SOFTWARE.
-**
-****************************************************************************/
-
-#ifndef LIBMV_TRACKING_SAD_H_
-#define LIBMV_TRACKING_SAD_H_
-
-#ifdef __cplusplus
-namespace libmv {
-#endif
-
-typedef unsigned char ubyte;
-typedef unsigned int uint;
-
-/*!
-    Convolve \a src into \a dst with the discrete laplacian operator.
-
-    \a src and \a dst should be \a width x \a height images.
-    \a strength is an interpolation coefficient (0-256) between original image and the laplacian.
-
-    \note Make sure the search region is filtered with the same strength as the pattern.
-*/
-void LaplaceFilter(ubyte* src, ubyte* dst, int width, int height, int strength);
-
-/// Affine transformation matrix in column major order.
-struct mat32 {
-  float data[3*2];
-#ifdef __cplusplus
-  inline mat32(int d=1) { for(int i=0;i<3*2;i++) data[i]=0; if(d!=0) for(int i=0;i<2;i++) m(i,i)=d; }
-  inline float m(int i, int j) const { return data[j*2+i]; }
-  inline float& m(int i, int j) { return data[j*2+i]; }
-  inline float operator()(int i, int j) const { return m(i,j); }
-  inline float& operator()(int i, int j) { return m(i,j); }
-  inline operator bool() const { for (int i=0; i<3*2; i++) if(data[i]!=0) return true; return false; }
-#endif
-};
-
-/*!
-    Sample \a pattern from \a image.
-
-    \a warp is the transformation to apply to \a image when sampling the \a pattern.
-*/
-void SamplePattern(ubyte* image, int stride, mat32 warp, ubyte* pattern, int size);
-
-/*!
-    Track \a pattern in \a image.
-
-    This template matcher computes the
-    \link http://en.wikipedia.org/wiki/Sum_of_absolute_differences Sum of Absolute Differences (SAD) \endlink
-    for each integer pixel position in the search region and then iteratively
-    refine subpixel position using a square search.
-    A similar method is used for motion estimation in video encoders.
-
-    \a reference is the pattern to track.
-    \a warped is a warped version of reference for fast unsampled integer search.
-       Best is to directly extract an already warped pattern from previous frame.
-    The \a size of the patterns should be aligned to 16.
-    \a image is a reference to the region to search.
-    \a stride is size of \a image lines.
-
-    On input, \a warp is the predicted affine transformation (e.g from previous frame)
-    On return, \a warp is the affine transformation which best match the reference \a pattern
-
-    \a areaPenalty and conditionPenalty control the regularization and need to be tweaked depending on the motion.
-       Setting them to 0 will allow any transformation (including unrealistic distortions and scaling).
-       Good values are between 0-32. 16 can be used as a realistic default.
-       areaPenalty control scaling (decrease to allow pull/zoom, increase to allow only 2D rotation).
-       a large conditionPenalty avoid a large ratio between the largest and smallest axices.
-       It need to be decreased for non-2D rotation (when pattern appears to scale along an axis).
-
-    \return Pearson product-moment correlation coefficient between reference and matched pattern.
-            This measure of the linear dependence between the patterns
-            ranges from −1 (negative correlation) to 1 (positive correlation).
-            A value of 0 implies that there is no linear correlation between the variables.
-
-    \note To track affine features:
-     - Sample reference pattern using estimated (e.g previous frame) warp.
-     -
-    \note \a stride allow you to reference your search region instead of copying.
-    \note For a 16x speedup, compile this tracker with SSE2 support.
-*/
-float Track(ubyte* reference, ubyte* warped, int size, ubyte* image, int stride, int width, int height, mat32* warp,
-            float areaPenalty, float conditionPenalty);
-
-#ifdef __cplusplus
-}  // namespace libmv
-#endif
-
-#endif  // LIBMV_TRACKING_SAD_H_