diff --git a/extern/libmv/CMakeLists.txt b/extern/libmv/CMakeLists.txt index 80a212b64f8..ef7e49b4d81 100644 --- a/extern/libmv/CMakeLists.txt +++ b/extern/libmv/CMakeLists.txt @@ -103,6 +103,7 @@ if(WITH_LIBMV) libmv/image/array_nd.h libmv/image/convolve.h libmv/image/correlation.h + libmv/image/image_converter.h libmv/image/image.h libmv/image/sample.h libmv/image/tuple.h diff --git a/extern/libmv/ChangeLog b/extern/libmv/ChangeLog index c34ea066172..225138b0166 100644 --- a/extern/libmv/ChangeLog +++ b/extern/libmv/ChangeLog @@ -1,3 +1,30 @@ +commit b63b8d6989f460fda7d963a2c8b21e8ffa6f8066 +Author: Sergey Sharybin +Date: Fri Jan 24 19:32:34 2014 +0600 + + Rework detector API and implement Harris detector + + Switch the detector API to a single function which accepts + a float image and detector options. This makes usage of + feature detection more unified across different algorithms. + + Options structure is pretty much straightforward and contains + detector to be used and all the detector-specific settings. + + Also implemented Harris feature detection algorithm which + is not as fast as FAST one but is expected to detect more + robust feature points. + + Reviewers: keir + + Differential Revision: https://developer.blender.org/D258 + +commit 6458915f64fceba108c5279b7320ca8c76e8a742 +Author: Sergey Sharybin +Date: Fri Jan 24 19:14:18 2014 +0600 + + Add arcanist configuration file + commit 0a69fadadc5aacbd339f839ac5bd12c3571278b1 Author: Sergey Sharybin Date: Thu Jan 9 15:50:11 2014 +0600 @@ -659,29 +686,3 @@ Date: Sat Apr 6 00:38:40 2013 +0600 Also fixed variable shadowing which lead to wrong markers positions (were either zero or undefined). - -commit a597ca9291ff7c6e209e0500945568afa84d8b4e -Author: Sergey Sharybin -Date: Sat Apr 6 00:12:12 2013 +0600 - - Camera intrinsics unit tests fixes - - - Existing test ApplyIsInvertibleSimple was not - doing right thing - distortion model used there - was ininvertible. - - Tweaked parameters in a way model is invertible now. - - - Added some own tests which tests; - - * Principal point always maps from pixel space to - zero normalized position. - - * Some basic tests to check whether individual - apply/invert works correct. - -commit a9a3016c7beb02e62655221d46783b206010d04c -Author: Sergey Sharybin -Date: Thu Apr 4 02:59:58 2013 +0600 - - Suppress strict compiler warnings in glags/glog libraries diff --git a/extern/libmv/files.txt b/extern/libmv/files.txt index e5297c7cced..a639ea3afa5 100644 --- a/extern/libmv/files.txt +++ b/extern/libmv/files.txt @@ -7,6 +7,7 @@ libmv/image/array_nd.h libmv/image/convolve.cc libmv/image/convolve.h libmv/image/correlation.h +libmv/image/image_converter.h libmv/image/image.h libmv/image/sample.h libmv/image/tuple.h diff --git a/extern/libmv/libmv-capi.cc b/extern/libmv/libmv-capi.cc index 474501afbcc..9a536f718aa 100644 --- a/extern/libmv/libmv-capi.cc +++ b/extern/libmv/libmv-capi.cc @@ -58,6 +58,14 @@ #include "libmv/simple_pipeline/reconstruction_scale.h" #include "libmv/simple_pipeline/keyframe_selection.h" +#include "libmv/multiview/homography.h" + +#include "libmv/image/convolve.h" + +#ifdef _MSC_VER +# define snprintf _snprintf +#endif + struct libmv_Reconstruction { libmv::EuclideanReconstruction reconstruction; @@ -69,7 +77,7 @@ struct libmv_Reconstruction { }; struct libmv_Features { - int count, margin; + int count; libmv::Feature *features; }; @@ -107,6 +115,21 @@ void libmv_setLoggingVerbosity(int verbosity) /* ************ Utility ************ */ +static void byteBufToImage(const unsigned char *buf, int width, int height, int channels, libmv::FloatImage *image) +{ + int x, y, k, a = 0; + + image->Resize(height, width, channels); + + for (y = 0; y < height; y++) { + for (x = 0; x < width; x++) { + for (k = 0; k < channels; k++) { + (*image)(y, x, k) = (float)buf[a++] / 255.0f; + } + } + } +} + static void floatBufToImage(const float *buf, int width, int height, int channels, libmv::FloatImage *image) { int x, y, k, a = 0; @@ -833,65 +856,90 @@ struct libmv_CameraIntrinsics *libmv_reconstructionExtractIntrinsics(struct libm /* ************ Feature detector ************ */ -struct libmv_Features *libmv_detectFeaturesFAST(const unsigned char *data, - int width, int height, int stride, - int margin, int min_trackness, int min_distance) +static libmv_Features *libmv_featuresFromVector( + const libmv::vector &features) { - libmv::Feature *features = NULL; - std::vector v; struct libmv_Features *libmv_features = LIBMV_STRUCT_NEW(libmv_Features, 1); - int i = 0, count; - - if (margin) { - data += margin * stride+margin; - width -= 2 * margin; - height -= 2 * margin; - } - - v = libmv::DetectFAST(data, width, height, stride, min_trackness, min_distance); - - count = v.size(); - + int count = features.size(); if (count) { - features = LIBMV_STRUCT_NEW(libmv::Feature, count); + libmv_features->features = LIBMV_STRUCT_NEW(libmv::Feature, count); - for(std::vector::iterator it = v.begin(); it != v.end(); it++) { - features[i++] = *it; + for (int i = 0; i < count; i++) { + libmv_features->features[i] = features.at(i); } } - - libmv_features->features = features; - libmv_features->count = count; - libmv_features->margin = margin; - - return (struct libmv_Features *)libmv_features; -} - -struct libmv_Features *libmv_detectFeaturesMORAVEC(const unsigned char *data, - int width, int height, int stride, - int margin, int count, int min_distance) -{ - libmv::Feature *features = NULL; - struct libmv_Features *libmv_features = LIBMV_STRUCT_NEW(libmv_Features, 1); - - if (count) { - if (margin) { - data += margin * stride+margin; - width -= 2 * margin; - height -= 2 * margin; - } - - features = LIBMV_STRUCT_NEW(libmv::Feature, count); - libmv::DetectMORAVEC(data, stride, width, height, features, &count, min_distance, NULL); + else { + libmv_features->features = NULL; } libmv_features->count = count; - libmv_features->margin = margin; - libmv_features->features = features; return libmv_features; } +static void libmv_convertDetectorOptions(libmv_DetectOptions *options, + libmv::DetectOptions *detector_options) +{ + switch (options->detector) { +#define LIBMV_CONVERT(the_detector) \ + case LIBMV_DETECTOR_ ## the_detector: \ + detector_options->type = libmv::DetectOptions::the_detector; \ + break; + LIBMV_CONVERT(FAST) + LIBMV_CONVERT(MORAVEC) + LIBMV_CONVERT(HARRIS) +#undef LIBMV_CONVERT + } + detector_options->margin = options->margin; + detector_options->min_distance = options->min_distance; + detector_options->fast_min_trackness = options->fast_min_trackness; + detector_options->moravec_max_count = options->moravec_max_count; + detector_options->moravec_pattern = options->moravec_pattern; + detector_options->harris_threshold = options->harris_threshold; +} + +struct libmv_Features *libmv_detectFeaturesByte(const unsigned char *image_buffer, + int width, int height, int channels, + libmv_DetectOptions *options) +{ + // Prepare the image. + libmv::FloatImage image; + byteBufToImage(image_buffer, width, height, channels, &image); + + // Configure detector. + libmv::DetectOptions detector_options; + libmv_convertDetectorOptions(options, &detector_options); + + // Run the detector. + libmv::vector detected_features; + libmv::Detect(image, detector_options, &detected_features); + + // Convert result to C-API. + libmv_Features *result = libmv_featuresFromVector(detected_features); + return result; +} + +struct libmv_Features *libmv_detectFeaturesFloat(const float *image_buffer, + int width, int height, int channels, + libmv_DetectOptions *options) +{ + // Prepare the image. + libmv::FloatImage image; + floatBufToImage(image_buffer, width, height, channels, &image); + + // Configure detector. + libmv::DetectOptions detector_options; + libmv_convertDetectorOptions(options, &detector_options); + + // Run the detector. + libmv::vector detected_features; + libmv::Detect(image, detector_options, &detected_features); + + // Convert result to C-API. + libmv_Features *result = libmv_featuresFromVector(detected_features); + return result; +} + void libmv_featuresDestroy(struct libmv_Features *libmv_features) { if (libmv_features->features) { @@ -910,8 +958,8 @@ void libmv_getFeature(const struct libmv_Features *libmv_features, int number, d { libmv::Feature feature = libmv_features->features[number]; - *x = feature.x + libmv_features->margin; - *y = feature.y + libmv_features->margin; + *x = feature.x; + *y = feature.y; *score = feature.score; *size = feature.size; } diff --git a/extern/libmv/libmv-capi.h b/extern/libmv/libmv-capi.h index 13cc3ae8499..d183bc4cd41 100644 --- a/extern/libmv/libmv-capi.h +++ b/extern/libmv/libmv-capi.h @@ -118,10 +118,29 @@ double libmv_reprojectionError(const struct libmv_Reconstruction *libmv_reconstr struct libmv_CameraIntrinsics *libmv_reconstructionExtractIntrinsics(struct libmv_Reconstruction *libmv_Reconstruction); /* Feature detector */ -struct libmv_Features *libmv_detectFeaturesFAST(const unsigned char *data, int width, int height, int stride, - int margin, int min_trackness, int min_distance); -struct libmv_Features *libmv_detectFeaturesMORAVEC(const unsigned char *data, int width, int height, int stride, - int margin, int count, int min_distance); +enum { + LIBMV_DETECTOR_FAST, + LIBMV_DETECTOR_MORAVEC, + LIBMV_DETECTOR_HARRIS, +}; + +typedef struct libmv_DetectOptions { + int detector; + int margin; + int min_distance; + int fast_min_trackness; + int moravec_max_count; + unsigned char *moravec_pattern; + double harris_threshold; +} libmv_DetectOptions; + +struct libmv_Features *libmv_detectFeaturesByte(const unsigned char *image_buffer, + int width, int height, int channels, + libmv_DetectOptions *options); +struct libmv_Features *libmv_detectFeaturesFloat(const float *image_buffer, + int width, int height, int channels, + libmv_DetectOptions *options); + void libmv_featuresDestroy(struct libmv_Features *libmv_features); int libmv_countFeatures(const struct libmv_Features *libmv_features); void libmv_getFeature(const struct libmv_Features *libmv_features, int number, double *x, double *y, double *score, diff --git a/extern/libmv/libmv-capi_stub.cc b/extern/libmv/libmv-capi_stub.cc index bda9605b422..bd5d16c9077 100644 --- a/extern/libmv/libmv-capi_stub.cc +++ b/extern/libmv/libmv-capi_stub.cc @@ -146,15 +146,16 @@ void libmv_reconstructionDestroy(struct libmv_Reconstruction * /*libmv_reconstru /* ************ feature detector ************ */ -struct libmv_Features *libmv_detectFeaturesFAST(const unsigned char * /*data*/, int /*width*/, int /*height*/, - int /*stride*/, int /*margin*/, int /*min_trackness*/, - int /*min_distance*/) +struct libmv_Features *libmv_detectFeaturesByte(const unsigned char */*image_buffer*/, + int /*width*/, int /*height*/, int /*channels*/, + libmv_DetectOptions */*options*/) { return NULL; } -struct libmv_Features *libmv_detectFeaturesMORAVEC(const unsigned char * /*data*/, int /*width*/, int /*height*/, - int /*stride*/, int /*margin*/, int /*count*/, int /*min_distance*/) +struct libmv_Features *libmv_detectFeaturesFloat(const float */*image_buffer*/, + int /*width*/, int /*height*/, int /*channels*/, + libmv_DetectOptions */*options*/) { return NULL; } diff --git a/extern/libmv/libmv/image/array_nd.h b/extern/libmv/libmv/image/array_nd.h index f40fadcbde7..c5099f24d7b 100644 --- a/extern/libmv/libmv/image/array_nd.h +++ b/extern/libmv/libmv/image/array_nd.h @@ -471,6 +471,23 @@ void MultiplyElements(const ArrayND &a, } } +template +void MultiplyElements(const Array3D &a, + const Array3D &b, + Array3D *c) { + // Specialization for N==3 + c->ResizeLike(a); + assert(a.Shape(0) == b.Shape(0)); + assert(a.Shape(1) == b.Shape(1)); + assert(a.Shape(2) == b.Shape(2)); + for (int i = 0; i < a.Shape(0); ++i) { + for (int j = 0; j < a.Shape(1); ++j) { + for (int k = 0; k < a.Shape(2); ++k) { + (*c)(i, j, k) = TC(a(i, j, k) * b(i, j, k)); + } + } + } +} } // namespace libmv diff --git a/extern/libmv/libmv/image/convolve.cc b/extern/libmv/libmv/image/convolve.cc index 4682c619f04..464043581d2 100644 --- a/extern/libmv/libmv/image/convolve.cc +++ b/extern/libmv/libmv/image/convolve.cc @@ -182,6 +182,23 @@ void ConvolveGaussian(const Array3Df &in, ConvolveHorizontal(tmp, kernel, out_pointer); } +void ImageDerivatives(const Array3Df &in, + double sigma, + Array3Df *gradient_x, + Array3Df *gradient_y) { + Vec kernel, derivative; + ComputeGaussianKernel(sigma, &kernel, &derivative); + Array3Df tmp; + + // Compute first derivative in x. + ConvolveVertical(in, kernel, &tmp); + ConvolveHorizontal(tmp, derivative, gradient_x); + + // Compute first derivative in y. + ConvolveHorizontal(in, kernel, &tmp); + ConvolveVertical(tmp, derivative, gradient_y); +} + void BlurredImageAndDerivatives(const Array3Df &in, double sigma, Array3Df *blurred_image, diff --git a/extern/libmv/libmv/image/image_converter.h b/extern/libmv/libmv/image/image_converter.h new file mode 100644 index 00000000000..b3a3fa2bf8c --- /dev/null +++ b/extern/libmv/libmv/image/image_converter.h @@ -0,0 +1,81 @@ +// Copyright (c) 2009 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_IMAGE_CONVERTER_H +#define LIBMV_IMAGE_IMAGE_CONVERTER_H + +#include "libmv/image/array_nd.h" + +namespace libmv { + +// The factor comes from http://www.easyrgb.com/ +// RGB to XYZ : Y is the luminance channel +// var_R * 0.2126 + var_G * 0.7152 + var_B * 0.0722 +template +inline T RGB2GRAY(const T r, const T g, const T b) { + return static_cast(r * 0.2126 + g * 0.7152 + b * 0.0722); +} + +/* +// Specialization for the uchar type. (that do not want to work...) +template<> +inline unsigned char RGB2GRAY(const unsigned char r, + const unsigned char g, + const unsigned char b) { + return (unsigned char)(r * 0.2126 + g * 0.7152 + b * 0.0722 +0.5); +}*/ + +template +void Rgb2Gray(const ImageIn &imaIn, ImageOut *imaOut) { + // It is all fine to cnvert RGBA image here as well, + // all the additional channels will be nicely ignored. + assert(imaIn.Depth() >= 3); + + imaOut->Resize(imaIn.Height(), imaIn.Width(), 1); + // Convert each RGB pixel into Gray value (luminance) + + for (int j = 0; j < imaIn.Height(); ++j) { + for (int i = 0; i < imaIn.Width(); ++i) { + (*imaOut)(j, i) = RGB2GRAY(imaIn(j, i, 0) , imaIn(j, i, 1), imaIn(j, i, 2)); + } + } +} + +// Convert given float image to an unsigned char array of pixels. +template +unsigned char *FloatImageToUCharArray(const Image &image) { + unsigned char *buffer = + new unsigned char[image.Width() * image.Height() * image.Depth()]; + + for (int y = 0; y < image.Height(); y++) { + for (int x = 0; x < image.Width(); x++) { + for (int d = 0; d < image.Depth(); d++) { + int index = (y * image.Width() + x) * image.Depth() + d; + buffer[index] = 255.0 * image(y, x, d); + } + } + } + + return buffer; +} + +} // namespace libmv + +#endif // LIBMV_IMAGE_IMAGE_CONVERTER_H diff --git a/extern/libmv/libmv/simple_pipeline/detect.cc b/extern/libmv/libmv/simple_pipeline/detect.cc index 27639e958a0..1ddb23dcd28 100644 --- a/extern/libmv/libmv/simple_pipeline/detect.cc +++ b/extern/libmv/libmv/simple_pipeline/detect.cc @@ -22,96 +22,133 @@ ** ****************************************************************************/ -#include "libmv/simple_pipeline/detect.h" -#include #include #include +#include +#include + +#include "libmv/base/scoped_ptr.h" +#include "libmv/image/array_nd.h" +#include "libmv/image/image_converter.h" +#include "libmv/image/convolve.h" +#include "libmv/logging/logging.h" +#include "libmv/simple_pipeline/detect.h" #ifdef __SSE2__ -#include +# include #endif namespace libmv { -typedef unsigned int uint; +namespace { -static int featurecmp(const void *a_v, const void *b_v) { - Feature *a = (Feature*)a_v; - Feature *b = (Feature*)b_v; +class FeatureComparison { + public: + bool operator() (const Feature &left, const Feature &right) const { + return right.score > left.score; + } +}; - return b->score - a->score; +// Filter the features so there are no features closer than +// minimal distance to each other. +// This is a naive implementation with O(n^2) asymptotic. +void FilterFeaturesByDistance(const vector &all_features, + int min_distance, + vector *detected_features) { + const int min_distance_squared = min_distance * min_distance; + + // Use priority queue to sort the features by their score. + // + // Do this on copy of the input features to prevent possible + // distortion in callee function behavior. + std::priority_queue, + FeatureComparison> priority_features; + + for (int i = 0; i < all_features.size(); i++) { + priority_features.push(all_features.at(i)); + } + + while (!priority_features.empty()) { + bool ok = true; + Feature a = priority_features.top(); + + for (int i = 0; i < detected_features->size(); i++) { + Feature &b = detected_features->at(i); + if (Square(a.x - b.x) + Square(a.y - b.y) < min_distance_squared) { + ok = false; + break; + } + } + + if (ok) { + detected_features->push_back(a); + } + + priority_features.pop(); + } } -std::vector DetectFAST(const unsigned char* data, - int width, int height, - int stride, - int min_trackness, - int min_distance) { - std::vector features; +void DetectFAST(const FloatImage &grayscale_image, + const DetectOptions &options, + vector *detected_features) { + const int min_distance = options.min_distance; + const int min_trackness = options.fast_min_trackness; + const int margin = options.margin; + const int width = grayscale_image.Width() - 2 * margin; + const int height = grayscale_image.Width() - 2 * margin; + const int stride = grayscale_image.Width(); + + // TODO(sergey): Use scoped_array to guard de-allocation of the array. + // Same goes to `scores` and `all` arrays here. + unsigned char *byte_image = FloatImageToUCharArray(grayscale_image); + const int byte_image_offset = margin * stride + margin; + // TODO(MatthiasF): Support targetting a feature count (binary search trackness) int num_features; - xy* all = fast9_detect(data, width, height, - stride, min_trackness, &num_features); + xy *all = fast9_detect(byte_image + byte_image_offset, + width, + height, + stride, + min_trackness, + &num_features); if (num_features == 0) { free(all); - return features; + delete [] byte_image; + return; } - int* scores = fast9_score(data, stride, all, num_features, min_trackness); + int *scores = fast9_score(byte_image + byte_image_offset, + stride, + all, + num_features, + min_trackness); // TODO(MatthiasF): merge with close feature suppression - xy* nonmax = nonmax_suppression(all, scores, num_features, &num_features); + xy *nonmax = nonmax_suppression(all, scores, num_features, &num_features); free(all); + delete [] byte_image; // Remove too close features // TODO(MatthiasF): A resolution independent parameter would be better than // distance e.g. a coefficient going from 0 (no minimal distance) to 1 // (optimal circle packing) // FIXME(MatthiasF): this method will not necessarily give all maximum markers if (num_features) { - Feature *all_features = new Feature[num_features]; - + vector all_features; for (int i = 0; i < num_features; ++i) { - Feature a = { (float)nonmax[i].x, (float)nonmax[i].y, (float)scores[i], 0 }; - all_features[i] = a; + Feature new_feature = {(float)nonmax[i].x + margin, + (float)nonmax[i].y + margin, + (float)scores[i], + 0}; + all_features.push_back(new_feature); } - - qsort((void *)all_features, num_features, sizeof(Feature), featurecmp); - - features.reserve(num_features); - - int prev_score = all_features[0].score; - const int min_distance_squared = min_distance * min_distance; - for (int i = 0; i < num_features; ++i) { - bool ok = true; - Feature a = all_features[i]; - if (a.score>prev_score) - abort(); - prev_score = a.score; - - // compare each feature against filtered set - for (int j = 0; j < features.size(); j++) { - Feature& b = features[j]; - if ((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y) < min_distance_squared) { - // already a nearby feature - ok = false; - break; - } - } - - if (ok) { - // add the new feature - features.push_back(a); - } - } - - delete [] all_features; + FilterFeaturesByDistance(all_features, min_distance, detected_features); } free(scores); free(nonmax); - return features; } #ifdef __SSE2__ -static uint SAD(const ubyte* imageA, const ubyte* imageB, - int strideA, int strideB) { +static unsigned int SAD(const ubyte* imageA, const ubyte* imageB, + int strideA, int strideB) { __m128i a = _mm_setzero_si128(); for (int i = 0; i < 16; i++) { a = _mm_adds_epu16(a, @@ -121,9 +158,9 @@ static uint SAD(const ubyte* imageA, const ubyte* imageB, return _mm_extract_epi16(a, 0) + _mm_extract_epi16(a, 4); } #else -static uint SAD(const ubyte* imageA, const ubyte* imageB, - int strideA, int strideB) { - uint sad = 0; +static unsigned int SAD(const ubyte* imageA, const ubyte* imageB, + int strideA, int strideB) { + unsigned int sad = 0; for (int i = 0; i < 16; i++) { for (int j = 0; j < 16; j++) { sad += abs((int)imageA[i*strideA+j] - imageB[i*strideB+j]); @@ -133,11 +170,20 @@ static uint SAD(const ubyte* imageA, const ubyte* imageB, } #endif -void DetectMORAVEC(const ubyte* image, - int stride, int width, int height, - Feature* detected, int* count, - int distance, - ubyte* pattern) { +void DetectMORAVEC(const FloatImage &grayscale_image, + const DetectOptions &options, + vector *detected_features) { + const int distance = options.min_distance; + const int margin = options.margin; + const unsigned char *pattern = options.moravec_pattern; + const int count = options.moravec_max_count; + const int width = grayscale_image.Width() - 2 * margin; + const int height = grayscale_image.Width() - 2 * margin; + const int stride = grayscale_image.Width(); + + // TODO(sergey): Use scoped_array to guard de-allocation of the array. + unsigned char *byte_image = FloatImageToUCharArray(grayscale_image); + unsigned short histogram[256]; memset(histogram, 0, sizeof(histogram)); ubyte* scores = new ubyte[width*height]; @@ -145,7 +191,7 @@ void DetectMORAVEC(const ubyte* image, const int r = 1; // radius for self similarity comparison for (int y = distance; y < height-distance; y++) { for (int x = distance; x < width-distance; x++) { - const ubyte* s = &image[y*stride+x]; + const ubyte* s = &byte_image[y*stride+x]; int score = // low self-similarity with overlapping patterns // OPTI: load pattern once SAD(s, s-r*stride-r, stride, stride)+SAD(s, s-r*stride, stride, stride)+SAD(s, s-r*stride+r, stride, stride)+ @@ -182,19 +228,105 @@ void DetectMORAVEC(const ubyte* image, int min = 255, total = 0; for (; min > 0; min--) { int h = histogram[min]; - if (total+h > *count) break; + if (total+h > count) break; total += h; } - int i = 0; for (int y = 16; y < height-16; y++) { for (int x = 16; x < width-16; x++) { int s = scores[y*width+x]; Feature f = { (float)x+8.0f, (float)y+8.0f, (float)s, 16 }; - if (s > min) detected[i++] = f; + if (s > min) { + detected_features->push_back(f); + } } } - *count = i; delete[] scores; + delete[] byte_image; +} + +void DetectHarris(const FloatImage &grayscale_image, + const DetectOptions &options, + vector *detected_features) { + const double alpha = 0.06; + const double sigma = 0.9; + + const int min_distance = options.min_distance; + const int margin = options.margin; + const double threshold = options.harris_threshold; + + FloatImage gradient_x, gradient_y; + ImageDerivatives(grayscale_image, sigma, &gradient_x, &gradient_y); + + FloatImage gradient_xx, gradient_yy, gradient_xy; + MultiplyElements(gradient_x, gradient_x, &gradient_xx); + MultiplyElements(gradient_y, gradient_y, &gradient_yy); + MultiplyElements(gradient_x, gradient_y, &gradient_xy); + + FloatImage gradient_xx_blurred, + gradient_yy_blurred, + gradient_xy_blurred; + ConvolveGaussian(gradient_xx, sigma, &gradient_xx_blurred); + ConvolveGaussian(gradient_yy, sigma, &gradient_yy_blurred); + ConvolveGaussian(gradient_xy, sigma, &gradient_xy_blurred); + + vector all_features; + for (int y = margin; y < gradient_xx_blurred.Height() - margin; ++y) { + for (int x = margin; x < gradient_xx_blurred.Width() - margin; ++x) { + // Construct matrix + // + // A = [ Ix^2 Ix*Iy ] + // [ Ix*Iy Iy^2 ] + Mat2 A; + A(0, 0) = gradient_xx_blurred(y, x); + A(1, 1) = gradient_yy_blurred(y, x); + A(0, 1) = A(1, 0) = gradient_xy_blurred(y, x); + + double detA = A.determinant(); + double traceA = A.trace(); + double harris_function = detA - alpha * traceA * traceA; + if (harris_function > threshold) { + Feature new_feature = {(float)x, (float)y, (float)harris_function, 0.0f}; + all_features.push_back(new_feature); + } + } + } + + FilterFeaturesByDistance(all_features, min_distance, detected_features); +} + +} // namespace + +DetectOptions::DetectOptions() + : type(DetectOptions::HARRIS), + margin(0), + min_distance(120), + fast_min_trackness(128), + moravec_max_count(0), + moravec_pattern(NULL), + harris_threshold(0.0) {} + +void Detect(const FloatImage &image, + const DetectOptions &options, + vector *detected_features) { + // Currently all the detectors requires image to be grayscale. + // Do it here to avoid code duplication. + FloatImage grayscale_image; + if (image.Depth() != 1) { + Rgb2Gray(image, &grayscale_image); + } else { + // TODO(sergey): Find a way to avoid such image duplication/ + grayscale_image = image; + } + + if (options.type == DetectOptions::FAST) { + DetectFAST(grayscale_image, options, detected_features); + } else if (options.type == DetectOptions::MORAVEC) { + DetectMORAVEC(grayscale_image, options, detected_features); + } else if (options.type == DetectOptions::HARRIS) { + DetectHarris(grayscale_image, options, detected_features); + } else { + LOG(FATAL) << "Unknown detector has been passed to featur detection"; + } } } // namespace libmv diff --git a/extern/libmv/libmv/simple_pipeline/detect.h b/extern/libmv/libmv/simple_pipeline/detect.h index f5bcd07afd7..2c782294d1b 100644 --- a/extern/libmv/libmv/simple_pipeline/detect.h +++ b/extern/libmv/libmv/simple_pipeline/detect.h @@ -27,6 +27,9 @@ #include +#include "libmv/base/vector.h" +#include "libmv/image/image.h" + namespace libmv { typedef unsigned char ubyte; @@ -50,47 +53,51 @@ struct Feature { float size; }; -/*! - Detect features in an image. +struct DetectOptions { + DetectOptions(); - You need to input a single channel 8-bit image using pointer to image \a data, - \a width, \a height and \a stride (i.e bytes per line). + // TODO(sergey): Write descriptions to each of algorithms. + enum DetectorType { + FAST, + MORAVEC, + HARRIS, + }; + DetectorType type; - You can tweak the count of detected features using \a min_trackness, which is - the minimum score to add a feature, and \a min_distance which is the minimal - distance accepted between two featuress. + // Margin in pixels from the image boundary. + // No features will be detected within the margin. + int margin; - \note You can binary search over \a min_trackness to get a given feature count. + // Minimal distance between detected features. + int min_distance; - \note a way to get an uniform distribution of a given feature count is: - \a min_distance = \a width * \a height / desired_feature_count ^ 2 + // Minimum score to add a feature. Only used by FAST detector. + // + // TODO(sergey): Write more detailed documentation about which + // units this value is measured in and so on. + int fast_min_trackness; - \return All detected feartures matching given parameters -*/ -std::vector DetectFAST(const unsigned char* data, int width, int height, - int stride, int min_trackness = 128, - int min_distance = 120); + // Maximum count to detect. Only used by MORAVEC detector. + int moravec_max_count; -/*! - Detect features in an image. + // Find only features similar to this pattern. Only used by MORAVEC detector. + // + // This is an image patch denoted in byte array with dimensions of 16px by 16px + // used to filter features by similarity to this patch. + unsigned char *moravec_pattern; - \a image is a single channel 8-bit image of size \a width x \a height + // Threshold value of the Harris function to add new featrue + // to the result. + double harris_threshold; +}; - \a detected is an array with space to hold \a *count features. - \a *count is the maximum count to detect on input and the actual - detected count on output. - - \a distance is the minimal distance between detected features. - - if \a pattern is null all good features will be found. - if \a pattern is not null only features similar to \a pattern will be found. - - \note \a You can crop the image (to avoid detecting markers near the borders) without copying: - image += marginY*stride+marginX, width -= 2*marginX, height -= 2*marginY; -*/ -void DetectMORAVEC(const ubyte* image, int stride, int width, int height, - Feature* detected, int* count, int distance /*=32*/, - ubyte* pattern /*=0*/); +// Detect features on a given image using given detector options. +// +// Image could have 1-4 channels, it'll be converted to a grayscale +// by the detector function if needed. +void Detect(const FloatImage &image, + const DetectOptions &options, + vector *detected_features); } // namespace libmv diff --git a/source/blender/blenkernel/BKE_tracking.h b/source/blender/blenkernel/BKE_tracking.h index 9936416847b..62149bbaf88 100644 --- a/source/blender/blenkernel/BKE_tracking.h +++ b/source/blender/blenkernel/BKE_tracking.h @@ -242,6 +242,10 @@ void BKE_tracking_detect_fast(struct MovieTracking *tracking, struct ListBase *t int framenr, int margin, int min_trackness, int min_distance, struct bGPDlayer *layer, bool place_outside_layer); +void BKE_tracking_detect_harris(struct MovieTracking *tracking, struct ListBase *tracksbase, struct ImBuf *ibuf, + int framenr, int margin, float threshold, int min_distance, struct bGPDlayer *layer, + bool place_outside_layer); + /* **** 2D stabilization **** */ void BKE_tracking_stabilization_data_get(struct MovieTracking *tracking, int framenr, int width, int height, float translation[2], float *scale, float *angle); diff --git a/source/blender/blenkernel/intern/tracking_detect.c b/source/blender/blenkernel/intern/tracking_detect.c index b262844d40f..65fd997ae6c 100644 --- a/source/blender/blenkernel/intern/tracking_detect.c +++ b/source/blender/blenkernel/intern/tracking_detect.c @@ -113,8 +113,11 @@ static void detect_retrieve_libmv_features(MovieTracking *tracking, ListBase *tr libmv_getFeature(features, a, &x, &y, &score, &size); - xu = x / width; - yu = y / height; + /* In Libmv integer coordinate points to pixel center, in blender + * it's not. Need to add 0.5px offset to center. + */ + xu = (x + 0.5f) / width; + yu = (y + 0.5f) / height; if (layer) ok = check_point_in_layer(layer, xu, yu) != place_outside_layer; @@ -128,36 +131,26 @@ static void detect_retrieve_libmv_features(MovieTracking *tracking, ListBase *tr } } -/* Get a gray-scale unsigned char buffer from given image buffer - * wich will be used for feature detection. - */ -static unsigned char *detect_get_frame_ucharbuf(ImBuf *ibuf) +static void run_configured_detector(MovieTracking *tracking, ListBase *tracksbase, + ImBuf *ibuf, int framenr, bGPDlayer *layer, bool place_outside_layer, + libmv_DetectOptions *options) { - int x, y; - unsigned char *pixels, *cp; + struct libmv_Features *features = NULL; - cp = pixels = MEM_callocN(ibuf->x * ibuf->y * sizeof(unsigned char), "tracking ucharBuf"); - for (y = 0; y < ibuf->y; y++) { - for (x = 0; x < ibuf->x; x++) { - int pixel = ibuf->x * y + x; - - if (ibuf->rect_float) { - const float *rrgbf = ibuf->rect_float + pixel * 4; - const float gray_f = 0.2126f * rrgbf[0] + 0.7152f * rrgbf[1] + 0.0722f * rrgbf[2]; - - *cp = FTOCHAR(gray_f); - } - else { - const unsigned char *rrgb = (unsigned char *)ibuf->rect + pixel * 4; - - *cp = 0.2126f * rrgb[0] + 0.7152f * rrgb[1] + 0.0722f * rrgb[2]; - } - - cp++; - } + if (ibuf->rect_float) { + features = libmv_detectFeaturesFloat(ibuf->rect_float, ibuf->x, ibuf->y, 4, options); + } + else if (ibuf->rect) { + features = libmv_detectFeaturesByte((unsigned char *) ibuf->rect, ibuf->x, ibuf->y, 4, options); } - return pixels; + if (features != NULL) { + detect_retrieve_libmv_features(tracking, tracksbase, features, + framenr, ibuf->x, ibuf->y, layer, + place_outside_layer); + + libmv_featuresDestroy(features); + } } /* Detect features using FAST detector */ @@ -165,17 +158,29 @@ void BKE_tracking_detect_fast(MovieTracking *tracking, ListBase *tracksbase, ImB int framenr, int margin, int min_trackness, int min_distance, bGPDlayer *layer, bool place_outside_layer) { - struct libmv_Features *features; - unsigned char *pixels = detect_get_frame_ucharbuf(ibuf); + libmv_DetectOptions options = {0}; - features = libmv_detectFeaturesFAST(pixels, ibuf->x, ibuf->y, ibuf->x, - margin, min_trackness, min_distance); + options.detector = LIBMV_DETECTOR_FAST; + options.margin = margin; + options.min_distance = min_distance; + options.fast_min_trackness = min_trackness; - MEM_freeN(pixels); - - detect_retrieve_libmv_features(tracking, tracksbase, features, - framenr, ibuf->x, ibuf->y, layer, - place_outside_layer); - - libmv_featuresDestroy(features); + run_configured_detector(tracking, tracksbase, ibuf, framenr, layer, + place_outside_layer, &options); +} + +/* Detect features using Harris detector */ +void BKE_tracking_detect_harris(MovieTracking *tracking, ListBase *tracksbase, ImBuf *ibuf, + int framenr, int margin, float threshold, int min_distance, bGPDlayer *layer, + bool place_outside_layer) +{ + libmv_DetectOptions options = {0}; + + options.detector = LIBMV_DETECTOR_HARRIS; + options.margin = margin; + options.min_distance = min_distance; + options.harris_threshold = threshold; + + run_configured_detector(tracking, tracksbase, ibuf, framenr, layer, + place_outside_layer, &options); } diff --git a/source/blender/editors/space_clip/tracking_ops.c b/source/blender/editors/space_clip/tracking_ops.c index 1269f9536b6..59cfe88e8c9 100644 --- a/source/blender/editors/space_clip/tracking_ops.c +++ b/source/blender/editors/space_clip/tracking_ops.c @@ -2913,8 +2913,8 @@ static int detect_features_exec(bContext *C, wmOperator *op) MovieTrackingTrack *track = tracksbase->first; int placement = RNA_enum_get(op->ptr, "placement"); int margin = RNA_int_get(op->ptr, "margin"); - int min_trackability = RNA_int_get(op->ptr, "min_trackability"); int min_distance = RNA_int_get(op->ptr, "min_distance"); + float threshold = RNA_float_get(op->ptr, "threshold"); int place_outside_layer = 0; int framenr = ED_space_clip_get_clip_frame_number(sc); bGPDlayer *layer = NULL; @@ -2938,8 +2938,8 @@ static int detect_features_exec(bContext *C, wmOperator *op) track = track->next; } - BKE_tracking_detect_fast(tracking, tracksbase, ibuf, framenr, margin, - min_trackability, min_distance, layer, place_outside_layer); + BKE_tracking_detect_harris(tracking, tracksbase, ibuf, framenr, margin, + threshold / 100000.0f, min_distance, layer, place_outside_layer); IMB_freeImBuf(ibuf); @@ -2972,9 +2972,9 @@ void CLIP_OT_detect_features(wmOperatorType *ot) /* properties */ RNA_def_enum(ot->srna, "placement", placement_items, 0, "Placement", "Placement for detected features"); - RNA_def_int(ot->srna, "margin", 16, 0, INT_MAX, "Margin", "Only corners further than margin pixels from the image edges are considered", 0, 300); - RNA_def_int(ot->srna, "min_trackability", 16, 0, INT_MAX, "Trackability", "Minimum trackability score to add a corner", 0, 300); - RNA_def_int(ot->srna, "min_distance", 120, 0, INT_MAX, "Distance", "Minimal distance accepted between two corners", 0, 300); + RNA_def_int(ot->srna, "margin", 16, 0, INT_MAX, "Margin", "Only features further than margin pixels from the image edges are considered", 0, 300); + RNA_def_float(ot->srna, "threshold", 1.0f, 0.0001f, FLT_MAX, "Threshold", "Threshold level to consider feature good enough for tracking", 0.0001f, FLT_MAX); + RNA_def_int(ot->srna, "min_distance", 120, 0, INT_MAX, "Distance", "Minimal distance accepted between two features", 0, 300); } /********************** frame jump operator *********************/