Welcome to mirror list, hosted at ThFree Co, Russian Federation.

git.blender.org/blender.git - Unnamed repository; edit this file 'description' to name the repository.
summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorSergey Sharybin <sergey.vfx@gmail.com>2014-06-18 16:49:17 +0400
committerCampbell Barton <ideasman42@gmail.com>2014-06-18 20:09:16 +0400
commit306cbb82ecf0d7c1ba4fb0a1240175b1976bd25b (patch)
treee9eac65bf57126ac233f22fe88cc00df137ec4e5 /extern/libmv
parent47ec0394ca3d03e07c07a67e8f8d1625aedd39dd (diff)
GTest unit testing framework
Currently covers only small set of functionality.
Diffstat (limited to 'extern/libmv')
-rw-r--r--extern/libmv/CMakeLists.txt112
-rw-r--r--extern/libmv/SConscript2
-rwxr-xr-xextern/libmv/bundle.sh87
-rw-r--r--extern/libmv/files.txt36
-rw-r--r--extern/libmv/libmv/base/scoped_ptr_test.cc79
-rw-r--r--extern/libmv/libmv/base/vector_test.cc223
-rw-r--r--extern/libmv/libmv/image/array_nd_test.cc324
-rw-r--r--extern/libmv/libmv/image/convolve_test.cc110
-rw-r--r--extern/libmv/libmv/image/image_drawing.h285
-rw-r--r--extern/libmv/libmv/image/image_test.cc45
-rw-r--r--extern/libmv/libmv/image/sample_test.cc89
-rw-r--r--extern/libmv/libmv/image/tuple_test.cc83
-rw-r--r--extern/libmv/libmv/multiview/euclidean_resection_test.cc237
-rw-r--r--extern/libmv/libmv/multiview/fundamental_test.cc162
-rw-r--r--extern/libmv/libmv/multiview/homography_error.h248
-rw-r--r--extern/libmv/libmv/multiview/homography_test.cc261
-rw-r--r--extern/libmv/libmv/multiview/nviewtriangulation_test.cc94
-rw-r--r--extern/libmv/libmv/multiview/panography_kernel.cc51
-rw-r--r--extern/libmv/libmv/multiview/panography_kernel.h54
-rw-r--r--extern/libmv/libmv/multiview/panography_test.cc144
-rw-r--r--extern/libmv/libmv/multiview/projection_test.cc115
-rw-r--r--extern/libmv/libmv/multiview/resection_test.cc61
-rw-r--r--extern/libmv/libmv/multiview/test_data_sets.cc196
-rw-r--r--extern/libmv/libmv/multiview/test_data_sets.h105
-rw-r--r--extern/libmv/libmv/multiview/triangulation_test.cc47
-rw-r--r--extern/libmv/libmv/multiview/two_view_kernel.h137
-rw-r--r--extern/libmv/libmv/numeric/dogleg_test.cc95
-rw-r--r--extern/libmv/libmv/numeric/function_derivative_test.cc57
-rw-r--r--extern/libmv/libmv/numeric/levenberg_marquardt_test.cc56
-rw-r--r--extern/libmv/libmv/numeric/numeric_test.cc439
-rw-r--r--extern/libmv/libmv/numeric/poly_test.cc98
-rw-r--r--extern/libmv/libmv/simple_pipeline/camera_intrinsics_test.cc239
-rw-r--r--extern/libmv/libmv/simple_pipeline/detect_test.cc230
-rw-r--r--extern/libmv/libmv/simple_pipeline/intersect_test.cc81
-rw-r--r--extern/libmv/libmv/simple_pipeline/keyframe_selection_test.cc307
-rw-r--r--extern/libmv/libmv/simple_pipeline/modal_solver_test.cc79
-rw-r--r--extern/libmv/libmv/simple_pipeline/resect_test.cc234
-rw-r--r--extern/libmv/libmv/tracking/brute_region_tracker_test.cc51
-rw-r--r--extern/libmv/libmv/tracking/klt_region_tracker_test.cc51
-rw-r--r--extern/libmv/libmv/tracking/pyramid_region_tracker_test.cc80
40 files changed, 5434 insertions, 50 deletions
diff --git a/extern/libmv/CMakeLists.txt b/extern/libmv/CMakeLists.txt
index c3c5143cc16..62b387e8968 100644
--- a/extern/libmv/CMakeLists.txt
+++ b/extern/libmv/CMakeLists.txt
@@ -30,6 +30,9 @@ set(INC
.
)
+set(INC_SYS
+)
+
set(SRC
libmv-capi.h
libmv-capi_intern.h
@@ -67,6 +70,7 @@ if(WITH_LIBMV)
libmv/multiview/fundamental.cc
libmv/multiview/homography.cc
libmv/multiview/panography.cc
+ libmv/multiview/panography_kernel.cc
libmv/multiview/projection.cc
libmv/multiview/triangulation.cc
libmv/numeric/numeric.cc
@@ -92,10 +96,6 @@ if(WITH_LIBMV)
libmv/tracking/track_region.cc
libmv/tracking/trklt_region_tracker.cc
- third_party/gflags/gflags.cc
- third_party/gflags/gflags_completions.cc
- third_party/gflags/gflags_reporting.cc
-
libmv-util.h
libmv/base/aligned_malloc.h
libmv/base/id_generator.h
@@ -106,6 +106,7 @@ if(WITH_LIBMV)
libmv/image/convolve.h
libmv/image/correlation.h
libmv/image/image_converter.h
+ libmv/image/image_drawing.h
libmv/image/image.h
libmv/image/sample.h
libmv/image/tuple.h
@@ -113,13 +114,16 @@ if(WITH_LIBMV)
libmv/multiview/conditioning.h
libmv/multiview/euclidean_resection.h
libmv/multiview/fundamental.h
+ libmv/multiview/homography_error.h
libmv/multiview/homography.h
libmv/multiview/homography_parameterization.h
libmv/multiview/nviewtriangulation.h
libmv/multiview/panography.h
+ libmv/multiview/panography_kernel.h
libmv/multiview/projection.h
libmv/multiview/resection.h
libmv/multiview/triangulation.h
+ libmv/multiview/two_view_kernel.h
libmv/numeric/dogleg.h
libmv/numeric/function_derivative.h
libmv/numeric/levenberg_marquardt.h
@@ -149,18 +153,86 @@ if(WITH_LIBMV)
libmv/tracking/track_region.h
libmv/tracking/trklt_region_tracker.h
+ third_party/msinttypes/inttypes.h
+ third_party/msinttypes/stdint.h
+ )
+
+ if(WIN32)
+ list(APPEND INC
+ third_party/glog/src/windows
+ )
+
+ if(NOT MINGW)
+ list(APPEND INC
+ third_party/msinttypes
+ )
+ endif()
+ endif()
+
+ if(WITH_TESTS)
+ blender_add_lib(libmv_test_dataset "./libmv/multiview/test_data_sets.cc" "" "")
+ endif()
+
+ BLENDER_SRC_GTEST("libmv_scoped_ptr" "./libmv/base/scoped_ptr_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_vector" "./libmv/base/vector_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_array_nd" "./libmv/image/array_nd_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_convolve" "./libmv/image/convolve_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_image" "./libmv/image/image_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_sample" "./libmv/image/sample_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_tuple" "./libmv/image/tuple_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_euclidean_resection" "./libmv/multiview/euclidean_resection_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_fundamental" "./libmv/multiview/fundamental_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_homography" "./libmv/multiview/homography_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_nviewtriangulation" "./libmv/multiview/nviewtriangulation_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_panography" "./libmv/multiview/panography_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_projection" "./libmv/multiview/projection_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_resection" "./libmv/multiview/resection_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_triangulation" "./libmv/multiview/triangulation_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_dogleg" "./libmv/numeric/dogleg_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_function_derivative" "./libmv/numeric/function_derivative_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_levenberg_marquardt" "./libmv/numeric/levenberg_marquardt_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_numeric" "./libmv/numeric/numeric_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_poly" "./libmv/numeric/poly_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_camera_intrinsics" "./libmv/simple_pipeline/camera_intrinsics_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_detect" "./libmv/simple_pipeline/detect_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_intersect" "./libmv/simple_pipeline/intersect_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_keyframe_selection" "./libmv/simple_pipeline/keyframe_selection_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_modal_solver" "./libmv/simple_pipeline/modal_solver_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_resect" "./libmv/simple_pipeline/resect_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_brute_region_tracker" "./libmv/tracking/brute_region_tracker_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_klt_region_tracker" "./libmv/tracking/klt_region_tracker_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+ BLENDER_SRC_GTEST("libmv_pyramid_region_tracker" "./libmv/tracking/pyramid_region_tracker_test.cc" "libmv_test_dataset;extern_libmv;extern_ceres")
+else()
+ list(APPEND SRC
+ libmv-capi_stub.cc
+ )
+endif()
+
+blender_add_lib(extern_libmv "${SRC}" "${INC}" "${INC_SYS}")
+
+if(WITH_LIBMV)
+ add_subdirectory(third_party)
+endif()
+
+# make GLog a separate target, so it can be used for gtest as well.
+if(WITH_LIBMV OR WITH_TESTS)
+ # We compile GLog together with GFlag so we don't worry about
+ # adding extra lib to linker.
+ set(GLOG_SRC
+ third_party/gflags/gflags.cc
+ third_party/gflags/gflags_completions.cc
+ third_party/gflags/gflags_reporting.cc
+
third_party/gflags/config.h
third_party/gflags/gflags/gflags_completions.h
third_party/gflags/gflags/gflags_declare.h
third_party/gflags/gflags/gflags.h
third_party/gflags/mutex.h
third_party/gflags/util.h
- third_party/msinttypes/inttypes.h
- third_party/msinttypes/stdint.h
)
if(WIN32)
- list(APPEND SRC
+ list(APPEND GLOG_SRC
third_party/glog/src/logging.cc
third_party/glog/src/raw_logging.cc
third_party/glog/src/utilities.cc
@@ -185,18 +257,8 @@ if(WITH_LIBMV)
third_party/glog/src/windows/port.h
third_party/glog/src/windows/config.h
)
-
- list(APPEND INC
- third_party/glog/src/windows
- )
-
- if(NOT MINGW)
- list(APPEND INC
- third_party/msinttypes
- )
- endif()
else()
- list(APPEND SRC
+ list(APPEND GLOG_SRC
third_party/glog/src/demangle.cc
third_party/glog/src/logging.cc
third_party/glog/src/raw_logging.cc
@@ -228,14 +290,14 @@ if(WITH_LIBMV)
third_party/glog/src/utilities.h
)
endif()
-else()
- list(APPEND SRC
- libmv-capi_stub.cc
+
+ set(GLOG_INC
+ third_party/gflags
+ third_party/glog/src
)
-endif()
-blender_add_lib(extern_libmv "${SRC}" "${INC}" "${INC_SYS}")
+ set(GLOG_INC_SYS
+ )
-if(WITH_LIBMV)
- add_subdirectory(third_party)
+ blender_add_lib(extern_glog "${GLOG_SRC}" "${GLOG_INC}" "${GLOG_INC_SYS}")
endif()
diff --git a/extern/libmv/SConscript b/extern/libmv/SConscript
index dc129501847..a267c96520e 100644
--- a/extern/libmv/SConscript
+++ b/extern/libmv/SConscript
@@ -44,6 +44,8 @@ if env['WITH_BF_LIBMV']:
else:
src = env.Glob("libmv-capi_stub.cc")
+src = [src for src in src if src.find('_test.cc') == -1]
+
env.BlenderLib ( libname = 'extern_libmv', sources=src, includes=Split(incs), defines=defs, libtype=['extern', 'player'], priority=[20,137] )
if env['WITH_BF_LIBMV']:
diff --git a/extern/libmv/bundle.sh b/extern/libmv/bundle.sh
index e8e698a89f2..36df51770c5 100755
--- a/extern/libmv/bundle.sh
+++ b/extern/libmv/bundle.sh
@@ -30,15 +30,20 @@ rm -rf $tmp
chmod 664 ./third_party/glog/src/windows/*.cc ./third_party/glog/src/windows/*.h ./third_party/glog/src/windows/glog/*.h
-sources=`find ./libmv -type f -iname '*.cc' -or -iname '*.cpp' -or -iname '*.c' | sed -r 's/^\.\//\t\t/' | sort -d`
-headers=`find ./libmv -type f -iname '*.h' | sed -r 's/^\.\//\t\t/' | sort -d`
+sources=`find ./libmv -type f -iname '*.cc' -or -iname '*.cpp' -or -iname '*.c' | grep -v _test.cc | grep -v test_data_sets | sed -r 's/^\.\//\t\t/' | sort -d`
+headers=`find ./libmv -type f -iname '*.h' | grep -v test_data_sets | sed -r 's/^\.\//\t\t/' | sort -d`
-third_sources=`find ./third_party -type f -iname '*.cc' -or -iname '*.cpp' -or -iname '*.c' | grep -v glog | grep -v ceres | sed -r 's/^\.\//\t\t/' | sort -d`
-third_headers=`find ./third_party -type f -iname '*.h' | grep -v glog | grep -v ceres | sed -r 's/^\.\//\t\t/' | sort -d`
+third_sources=`find ./third_party -type f -iname '*.cc' -or -iname '*.cpp' -or -iname '*.c' | grep -v glog | grep -v gflags | grep -v ceres | sed -r 's/^\.\//\t\t/' | sort -d`
+third_headers=`find ./third_party -type f -iname '*.h' | grep -v glog | grep -v gflags | grep -v ceres | sed -r 's/^\.\//\t\t/' | sort -d`
third_glog_sources=`find ./third_party -type f -iname '*.cc' -or -iname '*.cpp' -or -iname '*.c' | grep glog | grep -v windows | sed -r 's/^\.\//\t\t\t/' | sort -d`
third_glog_headers=`find ./third_party -type f -iname '*.h' | grep glog | grep -v windows | sed -r 's/^\.\//\t\t\t/' | sort -d`
+third_gflags_sources=`find ./third_party -type f -iname '*.cc' -or -iname '*.cpp' -or -iname '*.c' | grep gflags | grep -v windows | sed -r 's/^\.\//\t\t/' | sort -d`
+third_gflags_headers=`find ./third_party -type f -iname '*.h' | grep gflags | grep -v windows | sed -r 's/^\.\//\t\t/' | sort -d`
+
+tests=`find ./libmv -type f -iname '*_test.cc' | sort -d | awk ' { name=gensub(".*/([A-Za-z_]+)_test.cc", "\\\\1", $1); printf("\tBLENDER_SRC_GTEST(\"libmv_%s\" \"%s\" \"libmv_test_dataset;extern_libmv;extern_ceres\")\n", name, $1) } '`
+
src_dir=`find ./libmv -type f -iname '*.cc' -exec dirname {} \; -or -iname '*.cpp' -exec dirname {} \; -or -iname '*.c' -exec dirname {} \; | sed -r 's/^\.\//\t\t/' | sort -d | uniq`
src_third_dir=`find ./third_party -type f -iname '*.cc' -exec dirname {} \; -or -iname '*.cpp' -exec dirname {} \; -or -iname '*.c' -exec dirname {} \; | grep -v ceres | sed -r 's/^\.\//\t\t/' | sort -d | uniq`
src=""
@@ -118,6 +123,9 @@ set(INC
.
)
+set(INC_SYS
+)
+
set(SRC
libmv-capi.h
libmv-capi_intern.h
@@ -148,9 +156,7 @@ if(WITH_LIBMV)
libmv-capi.cc
libmv-util.cc
${sources}
-
${third_sources}
-
libmv-util.h
${headers}
@@ -158,7 +164,46 @@ ${third_headers}
)
if(WIN32)
- list(APPEND SRC
+ list(APPEND INC
+ third_party/glog/src/windows
+ )
+
+ if(NOT MINGW)
+ list(APPEND INC
+ third_party/msinttypes
+ )
+ endif()
+ endif()
+
+ if(WITH_TESTS)
+ blender_add_lib(libmv_test_dataset "./libmv/multiview/test_data_sets.cc" "${INC}" "${INC_SYS}")
+ endif()
+
+${tests}
+else()
+ list(APPEND SRC
+ libmv-capi_stub.cc
+ )
+endif()
+
+blender_add_lib(extern_libmv "\${SRC}" "\${INC}" "\${INC_SYS}")
+
+if(WITH_LIBMV)
+ add_subdirectory(third_party)
+endif()
+
+# make GLog a separate target, so it can be used for gtest as well.
+if(WITH_LIBMV OR WITH_TESTS)
+ # We compile GLog together with GFlag so we don't worry about
+ # adding extra lib to linker.
+ set(GLOG_SRC
+${third_gflags_sources}
+
+${third_gflags_headers}
+ )
+
+ if(WIN32)
+ list(APPEND GLOG_SRC
third_party/glog/src/logging.cc
third_party/glog/src/raw_logging.cc
third_party/glog/src/utilities.cc
@@ -183,33 +228,23 @@ ${third_headers}
third_party/glog/src/windows/port.h
third_party/glog/src/windows/config.h
)
-
- list(APPEND INC
- third_party/glog/src/windows
- )
-
- if(NOT MINGW)
- list(APPEND INC
- third_party/msinttypes
- )
- endif()
else()
- list(APPEND SRC
+ list(APPEND GLOG_SRC
${third_glog_sources}
${third_glog_headers}
)
endif()
-else()
- list(APPEND SRC
- libmv-capi_stub.cc
+
+ set(GLOG_INC
+ third_party/gflags
+ third_party/glog/src
)
-endif()
-blender_add_lib(extern_libmv "\${SRC}" "\${INC}" "\${INC_SYS}")
+ set(GLOG_INC_SYS
+ )
-if(WITH_LIBMV)
- add_subdirectory(third_party)
+ blender_add_lib(extern_glog "\${GLOG_SRC}" "\${GLOG_INC}" "\${GLOG_INC_SYS}")
endif()
EOF
@@ -254,6 +289,8 @@ ${win_src}
else:
src = env.Glob("libmv-capi_stub.cc")
+src = [src for src in src if src.find('_test.cc') == -1]
+
env.BlenderLib ( libname = 'extern_libmv', sources=src, includes=Split(incs), defines=defs, libtype=['extern', 'player'], priority=[20,137] )
if env['WITH_BF_LIBMV']:
diff --git a/extern/libmv/files.txt b/extern/libmv/files.txt
index 60a99307835..f1eb9580ff5 100644
--- a/extern/libmv/files.txt
+++ b/extern/libmv/files.txt
@@ -2,60 +2,92 @@ libmv/base/aligned_malloc.cc
libmv/base/aligned_malloc.h
libmv/base/id_generator.h
libmv/base/scoped_ptr.h
+libmv/base/scoped_ptr_test.cc
libmv/base/vector.h
+libmv/base/vector_test.cc
libmv/base/vector_utils.h
libmv/image/array_nd.cc
libmv/image/array_nd.h
+libmv/image/array_nd_test.cc
libmv/image/convolve.cc
libmv/image/convolve.h
+libmv/image/convolve_test.cc
libmv/image/correlation.h
libmv/image/image_converter.h
+libmv/image/image_drawing.h
libmv/image/image.h
+libmv/image/image_test.cc
libmv/image/sample.h
+libmv/image/sample_test.cc
libmv/image/tuple.h
+libmv/image/tuple_test.cc
libmv/logging/logging.h
libmv/multiview/conditioning.cc
libmv/multiview/conditioning.h
libmv/multiview/euclidean_resection.cc
libmv/multiview/euclidean_resection.h
+libmv/multiview/euclidean_resection_test.cc
libmv/multiview/fundamental.cc
libmv/multiview/fundamental.h
+libmv/multiview/fundamental_test.cc
libmv/multiview/homography.cc
+libmv/multiview/homography_error.h
libmv/multiview/homography.h
libmv/multiview/homography_parameterization.h
+libmv/multiview/homography_test.cc
libmv/multiview/nviewtriangulation.h
+libmv/multiview/nviewtriangulation_test.cc
libmv/multiview/panography.cc
libmv/multiview/panography.h
+libmv/multiview/panography_kernel.cc
+libmv/multiview/panography_kernel.h
+libmv/multiview/panography_test.cc
libmv/multiview/projection.cc
libmv/multiview/projection.h
+libmv/multiview/projection_test.cc
libmv/multiview/resection.h
+libmv/multiview/resection_test.cc
+libmv/multiview/test_data_sets.cc
+libmv/multiview/test_data_sets.h
libmv/multiview/triangulation.cc
libmv/multiview/triangulation.h
+libmv/multiview/triangulation_test.cc
+libmv/multiview/two_view_kernel.h
libmv/numeric/dogleg.h
+libmv/numeric/dogleg_test.cc
libmv/numeric/function_derivative.h
+libmv/numeric/function_derivative_test.cc
libmv/numeric/levenberg_marquardt.h
+libmv/numeric/levenberg_marquardt_test.cc
libmv/numeric/numeric.cc
libmv/numeric/numeric.h
+libmv/numeric/numeric_test.cc
libmv/numeric/poly.cc
libmv/numeric/poly.h
+libmv/numeric/poly_test.cc
libmv/simple_pipeline/bundle.cc
libmv/simple_pipeline/bundle.h
libmv/simple_pipeline/callbacks.h
libmv/simple_pipeline/camera_intrinsics.cc
libmv/simple_pipeline/camera_intrinsics.h
libmv/simple_pipeline/camera_intrinsics_impl.h
+libmv/simple_pipeline/camera_intrinsics_test.cc
libmv/simple_pipeline/detect.cc
libmv/simple_pipeline/detect.h
+libmv/simple_pipeline/detect_test.cc
libmv/simple_pipeline/distortion_models.cc
libmv/simple_pipeline/distortion_models.h
libmv/simple_pipeline/initialize_reconstruction.cc
libmv/simple_pipeline/initialize_reconstruction.h
libmv/simple_pipeline/intersect.cc
libmv/simple_pipeline/intersect.h
+libmv/simple_pipeline/intersect_test.cc
libmv/simple_pipeline/keyframe_selection.cc
libmv/simple_pipeline/keyframe_selection.h
+libmv/simple_pipeline/keyframe_selection_test.cc
libmv/simple_pipeline/modal_solver.cc
libmv/simple_pipeline/modal_solver.h
+libmv/simple_pipeline/modal_solver_test.cc
libmv/simple_pipeline/pipeline.cc
libmv/simple_pipeline/pipeline.h
libmv/simple_pipeline/reconstruction.cc
@@ -64,16 +96,20 @@ libmv/simple_pipeline/reconstruction_scale.cc
libmv/simple_pipeline/reconstruction_scale.h
libmv/simple_pipeline/resect.cc
libmv/simple_pipeline/resect.h
+libmv/simple_pipeline/resect_test.cc
libmv/simple_pipeline/tracks.cc
libmv/simple_pipeline/tracks.h
libmv/tracking/brute_region_tracker.cc
libmv/tracking/brute_region_tracker.h
+libmv/tracking/brute_region_tracker_test.cc
libmv/tracking/hybrid_region_tracker.cc
libmv/tracking/hybrid_region_tracker.h
libmv/tracking/klt_region_tracker.cc
libmv/tracking/klt_region_tracker.h
+libmv/tracking/klt_region_tracker_test.cc
libmv/tracking/pyramid_region_tracker.cc
libmv/tracking/pyramid_region_tracker.h
+libmv/tracking/pyramid_region_tracker_test.cc
libmv/tracking/region_tracker.h
libmv/tracking/retrack_region_tracker.cc
libmv/tracking/retrack_region_tracker.h
diff --git a/extern/libmv/libmv/base/scoped_ptr_test.cc b/extern/libmv/libmv/base/scoped_ptr_test.cc
new file mode 100644
index 00000000000..ce1d56b500a
--- /dev/null
+++ b/extern/libmv/libmv/base/scoped_ptr_test.cc
@@ -0,0 +1,79 @@
+// 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.
+
+#include "libmv/base/scoped_ptr.h"
+#include "testing/testing.h"
+
+namespace libmv {
+namespace {
+
+struct FreeMe {
+ FreeMe(int *freed) : freed(freed) {}
+ ~FreeMe() { (*freed)++; }
+ int *freed;
+};
+
+TEST(ScopedPtr, NullDoesNothing) {
+ scoped_ptr<FreeMe> x(NULL);
+ // Does nothing.
+}
+
+TEST(ScopedPtr, FreesWhenOutOfScope) {
+ int frees = 0;
+ {
+ scoped_ptr<FreeMe> scoped(new FreeMe(&frees));
+ EXPECT_EQ(0, frees);
+ }
+ EXPECT_EQ(1, frees);
+}
+
+TEST(ScopedPtr, Operators) {
+ int tag = 123;
+ scoped_ptr<FreeMe> scoped(new FreeMe(&tag));
+ EXPECT_EQ(123, *(scoped->freed));
+ EXPECT_EQ(123, *((*scoped).freed));
+}
+
+TEST(ScopedPtr, Reset) {
+ int frees = 0;
+ scoped_ptr<FreeMe> scoped(new FreeMe(&frees));
+ EXPECT_EQ(0, frees);
+ scoped.reset(new FreeMe(&frees));
+ EXPECT_EQ(1, frees);
+}
+
+TEST(ScopedPtr, ReleaseAndGet) {
+ int frees = 0;
+ FreeMe *allocated = new FreeMe(&frees);
+ FreeMe *released = NULL;
+ {
+ scoped_ptr<FreeMe> scoped(allocated);
+ EXPECT_EQ(0, frees);
+ EXPECT_EQ(allocated, scoped.get());
+ released = scoped.release();
+ EXPECT_EQ(0, frees);
+ EXPECT_EQ(released, allocated);
+ }
+ EXPECT_EQ(0, frees);
+ delete released;
+}
+
+} // namespace
+} // namespace libmv
diff --git a/extern/libmv/libmv/base/vector_test.cc b/extern/libmv/libmv/base/vector_test.cc
new file mode 100644
index 00000000000..f17718c3926
--- /dev/null
+++ b/extern/libmv/libmv/base/vector_test.cc
@@ -0,0 +1,223 @@
+// 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.
+
+#include "libmv/base/vector.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+#include <algorithm>
+
+namespace {
+using namespace libmv;
+
+// This uses a Vec2d which is a fixed-size vectorizable Eigen type. It is
+// necessary to test vectorizable types to ensure that the alignment asserts
+// trigger if the alignment is not correct.
+TEST(VectorAlignmentTest, PushBack) {
+ Vec2 x1, x2;
+ x1 << 1, 2;
+ x2 << 3, 4;
+
+ vector<Vec2> vs;
+ vs.push_back(x1);
+ EXPECT_EQ(1, vs.size());
+ EXPECT_EQ(1, vs.capacity());
+
+ vs.push_back(x2);
+ EXPECT_EQ(2, vs.size());
+ EXPECT_EQ(2, vs.capacity());
+
+ // The following is necessary because of some bug in gtest; the expected
+ // parameter can't be a fixed size vectorizable type with alignment
+ // requirements.
+ Vec x1r = x1;
+ Vec x2r = x2;
+ EXPECT_EQ(x1r, vs[0]);
+ EXPECT_EQ(x2r, vs[1]);
+
+ vs.push_back(x2);
+ vs.push_back(x2);
+ vs.push_back(x2);
+ EXPECT_EQ(5, vs.size());
+ EXPECT_EQ(8, vs.capacity());
+}
+
+// Count the number of destruct calls to test that the destructor gets called.
+int foo_construct_calls = 0;
+int foo_destruct_calls = 0;
+
+struct Foo {
+ public:
+ Foo() : value(5) { foo_construct_calls++; }
+ ~Foo() { foo_destruct_calls++; }
+ int value;
+};
+
+struct VectorTest : public testing::Test {
+ VectorTest() {
+ foo_construct_calls = 0;
+ foo_destruct_calls = 0;
+ }
+};
+
+TEST_F(VectorTest, EmptyVectorDoesNotConstruct) {
+ {
+ vector<Foo> v;
+ EXPECT_EQ(0, v.size());
+ EXPECT_EQ(0, v.capacity());
+ }
+ EXPECT_EQ(0, foo_construct_calls);
+ EXPECT_EQ(0, foo_destruct_calls);
+}
+
+TEST_F(VectorTest, DestructorGetsCalled) {
+ {
+ vector<Foo> v;
+ v.resize(5);
+ }
+ EXPECT_EQ(5, foo_construct_calls);
+ EXPECT_EQ(5, foo_destruct_calls);
+}
+
+TEST_F(VectorTest, ReserveDoesNotCallConstructorsOrDestructors) {
+ vector<Foo> v;
+ EXPECT_EQ(0, v.size());
+ EXPECT_EQ(0, v.capacity());
+ EXPECT_EQ(0, foo_construct_calls);
+ EXPECT_EQ(0, foo_destruct_calls);
+
+ v.reserve(5);
+ EXPECT_EQ(0, v.size());
+ EXPECT_EQ(5, v.capacity());
+ EXPECT_EQ(0, foo_construct_calls);
+ EXPECT_EQ(0, foo_destruct_calls);
+}
+
+TEST_F(VectorTest, ResizeConstructsAndDestructsAsExpected) {
+ vector<Foo> v;
+
+ // Create one object.
+ v.resize(1);
+ EXPECT_EQ(1, v.size());
+ EXPECT_EQ(1, v.capacity());
+ EXPECT_EQ(1, foo_construct_calls);
+ EXPECT_EQ(0, foo_destruct_calls);
+ EXPECT_EQ(5, v[0].value);
+
+ // Create two more.
+ v.resize(3);
+ EXPECT_EQ(3, v.size());
+ EXPECT_EQ(3, v.capacity());
+ EXPECT_EQ(3, foo_construct_calls);
+ EXPECT_EQ(0, foo_destruct_calls);
+
+ // Delete the last one.
+ v.resize(2);
+ EXPECT_EQ(2, v.size());
+ EXPECT_EQ(3, v.capacity());
+ EXPECT_EQ(3, foo_construct_calls);
+ EXPECT_EQ(1, foo_destruct_calls);
+
+ // Delete the remaining two.
+ v.resize(0);
+ EXPECT_EQ(0, v.size());
+ EXPECT_EQ(3, v.capacity());
+ EXPECT_EQ(3, foo_construct_calls);
+ EXPECT_EQ(3, foo_destruct_calls);
+}
+
+TEST_F(VectorTest, PushPopBack) {
+ vector<Foo> v;
+
+ Foo foo;
+ foo.value = 10;
+ v.push_back(foo);
+ EXPECT_EQ(1, v.size());
+ EXPECT_EQ(10, v.back().value);
+
+ v.pop_back();
+ EXPECT_EQ(0, v.size());
+ EXPECT_EQ(1, foo_construct_calls);
+ EXPECT_EQ(1, foo_destruct_calls);
+}
+
+TEST_F(VectorTest, CopyConstructor) {
+ vector<int> a;
+ a.push_back(1);
+ a.push_back(5);
+ a.push_back(3);
+
+ vector<int> b(a);
+ EXPECT_EQ(a.size(), b.size());
+ //EXPECT_EQ(a.capacity(), b.capacity());
+ for (int i = 0; i < a.size(); ++i) {
+ EXPECT_EQ(a[i], b[i]);
+ }
+}
+
+TEST_F(VectorTest, OperatorEquals) {
+ vector<int> a, b;
+ a.push_back(1);
+ a.push_back(5);
+ a.push_back(3);
+
+ b = a;
+
+ EXPECT_EQ(a.size(), b.size());
+ //EXPECT_EQ(a.capacity(), b.capacity());
+ for (int i = 0; i < a.size(); ++i) {
+ EXPECT_EQ(a[i], b[i]);
+ }
+}
+
+TEST_F(VectorTest, STLFind) {
+ vector<int> a;
+ a.push_back(1);
+ a.push_back(5);
+ a.push_back(3);
+
+ // Find return an int *
+ EXPECT_EQ(std::find(&a[0], &a[2], 1) == &a[0], true);
+ EXPECT_EQ(std::find(&a[0], &a[2], 5) == &a[1], true);
+ EXPECT_EQ(std::find(&a[0], &a[2], 3) == &a[2], true);
+
+ // Find return a const int *
+ EXPECT_EQ(std::find(a.begin(), a.end(), 1) == &a[0], true);
+ EXPECT_EQ(std::find(a.begin(), a.end(), 5) == &a[1], true);
+ EXPECT_EQ(std::find(a.begin(), a.end(), 3) == &a[2], true);
+
+ // Search value that are not in the vector
+ EXPECT_EQ(std::find(a.begin(), a.end(), 0) == a.end(), true);
+ EXPECT_EQ(std::find(a.begin(), a.end(), 52) == a.end(), true);
+}
+
+TEST(Vector, swap) {
+ vector<int> a, b;
+ a.push_back(1);
+ a.push_back(2);
+ b.push_back(3);
+ a.swap(b);
+ EXPECT_EQ(1, a.size());
+ EXPECT_EQ(3, a[0]);
+ EXPECT_EQ(2, b.size());
+ EXPECT_EQ(1, b[0]);
+ EXPECT_EQ(2, b[1]);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/image/array_nd_test.cc b/extern/libmv/libmv/image/array_nd_test.cc
new file mode 100644
index 00000000000..313f21b60e9
--- /dev/null
+++ b/extern/libmv/libmv/image/array_nd_test.cc
@@ -0,0 +1,324 @@
+// Copyright (c) 2007, 2008 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/image/array_nd.h"
+#include "testing/testing.h"
+
+using libmv::ArrayND;
+using libmv::Array3D;
+using libmv::Array3Df;
+
+namespace {
+
+TEST(ArrayND, EmptyConstructor) {
+ ArrayND<int, 2> a;
+
+ EXPECT_EQ(0, a.Shape(0));
+ EXPECT_EQ(0, a.Shape(1));
+}
+
+TEST(ArrayND, IndexConstructor) {
+ int s[] = {1, 2, 3};
+ ArrayND<int, 3>::Index shape(s);
+ ArrayND<int, 3> a(shape);
+
+ EXPECT_EQ(1, a.Shape(0));
+ EXPECT_EQ(2, a.Shape(1));
+ EXPECT_EQ(3, a.Shape(2));
+}
+
+TEST(ArrayND, PointerConstructor) {
+ int s[] = {1, 2, 3};
+ ArrayND<int, 3> a(s);
+
+ EXPECT_EQ(1, a.Shape(0));
+ EXPECT_EQ(2, a.Shape(1));
+ EXPECT_EQ(3, a.Shape(2));
+}
+
+TEST(ArrayND, CopyConstructor) {
+ int s[] = {1, 2, 3};
+ ArrayND<int, 3> a(s);
+ a(0, 1, 1) = 3;
+ a(0, 1, 2) = 3;
+ ArrayND<int, 3> b(a);
+ EXPECT_EQ(1, b.Shape(0));
+ EXPECT_EQ(2, b.Shape(1));
+ EXPECT_EQ(3, b.Shape(2));
+ EXPECT_EQ(3, b(0, 1, 1));
+ b(0, 1, 2) = 2;
+ EXPECT_EQ(3, a(0, 1, 2));
+}
+
+TEST(ArrayND, Assignation) {
+ int s[] = {1, 2, 3};
+ ArrayND<int, 3> a(s);
+ a(0, 1, 1) = 3;
+ a(0, 1, 2) = 3;
+ ArrayND<int, 3> b;
+ b = a;
+ EXPECT_EQ(1, b.Shape(0));
+ EXPECT_EQ(2, b.Shape(1));
+ EXPECT_EQ(3, b.Shape(2));
+ EXPECT_EQ(3, b(0, 1, 1));
+ b(0, 1, 2) = 2;
+ EXPECT_EQ(3, a(0, 1, 2));
+}
+
+TEST(ArrayND, Fill) {
+ int s[] = {2, 2};
+ ArrayND<int, 2> a(s);
+ a.Fill(42);
+ EXPECT_EQ(42, a(0, 0));
+ EXPECT_EQ(42, a(0, 1));
+ EXPECT_EQ(42, a(1, 0));
+ EXPECT_EQ(42, a(1, 1));
+}
+
+TEST(ArrayND, Size) {
+ int s[] = {1, 2, 3};
+ ArrayND<int, 3>::Index shape(s);
+ ArrayND<int, 3> a(shape);
+
+ int l[] = {0, 1, 2};
+ ArrayND<int, 3>::Index last(l);
+
+ EXPECT_EQ(a.Size(), a.Offset(last)+1);
+ EXPECT_TRUE(a.Contains(last));
+ EXPECT_FALSE(a.Contains(shape));
+}
+
+TEST(ArrayND, MemorySizeInBytes) {
+ int s[] = {2, 3};
+ ArrayND<int, 2>::Index shape(s);
+ ArrayND<int, 2> a(shape);
+
+ int size = 24 + sizeof(a);
+ EXPECT_EQ(size, a.MemorySizeInBytes());
+}
+
+TEST(ArrayND, Parenthesis) {
+ typedef ArrayND<int, 2>::Index Index;
+
+ int s[] = {3, 3};
+ ArrayND<int, 2> a(s);
+
+ *(a.Data()+0) = 0;
+ *(a.Data()+5) = 5;
+
+ int i1[] = {0, 0};
+ EXPECT_EQ(0, a(Index(i1)));
+ int i2[] = {1, 2};
+ EXPECT_EQ(5, a(Index(i2)));
+}
+
+TEST(ArrayND, 1DConstructor) {
+ ArrayND<int, 1> a(3);
+
+ EXPECT_EQ(3, a.Shape(0));
+}
+
+TEST(ArrayND, 2DConstructor) {
+ ArrayND<int, 2> a(1, 2);
+
+ EXPECT_EQ(1, a.Shape(0));
+ EXPECT_EQ(2, a.Shape(1));
+}
+
+TEST(ArrayND, 3DConstructor) {
+ ArrayND<int, 3> a(1, 2, 3);
+
+ EXPECT_EQ(1, a.Shape(0));
+ EXPECT_EQ(2, a.Shape(1));
+ EXPECT_EQ(3, a.Shape(2));
+}
+
+TEST(ArrayND, 1DAccessor) {
+ ArrayND<int, 1> a(3);
+ a(0) = 1;
+ a(1) = 2;
+
+ EXPECT_EQ(1, a(0));
+ EXPECT_EQ(2, a(1));
+ EXPECT_EQ(1, *(a.Data()));
+ EXPECT_EQ(2, *(a.Data() + a.Stride(0)));
+}
+
+TEST(ArrayND, 2DAccessor) {
+ ArrayND<int, 2> a(3, 3);
+ a(0, 0) = 1;
+ a(1, 1) = 2;
+
+ EXPECT_EQ(1, a(0, 0));
+ EXPECT_EQ(2, a(1, 1));
+ EXPECT_EQ(1, *(a.Data()));
+ EXPECT_EQ(2, *(a.Data() + a.Stride(0) + a.Stride(1)));
+}
+
+TEST(ArrayND, 3DAccessor) {
+ ArrayND<int, 3> a(3, 3, 3);
+ a(0, 0, 0) = 1;
+ a(1, 1, 1) = 2;
+
+ EXPECT_EQ(1, a(0, 0, 0));
+ EXPECT_EQ(2, a(1, 1, 1));
+ EXPECT_EQ(1, *(a.Data()));
+ EXPECT_EQ(2, *(a.Data() + a.Stride(0) + a.Stride(1) + a.Stride(2)));
+}
+
+TEST(ArrayND, CopyFrom) {
+ ArrayND<int, 3> a(2, 2, 1);
+ a(0, 0, 0) = 1;
+ a(0, 1, 0) = 2;
+ a(1, 0, 0) = 3;
+ a(1, 1, 0) = 4;
+ ArrayND<float, 3> b;
+ b.CopyFrom(a);
+ EXPECT_FLOAT_EQ(1.f, b(0, 0, 0));
+ EXPECT_FLOAT_EQ(2.f, b(0, 1, 0));
+ EXPECT_FLOAT_EQ(3.f, b(1, 0, 0));
+ EXPECT_FLOAT_EQ(4.f, b(1, 1, 0));
+}
+
+TEST(ArrayND, MultiplyElements) {
+ ArrayND<int, 3> a(2, 2, 1);
+ a(0, 0, 0) = 1;
+ a(0, 1, 0) = 2;
+ a(1, 0, 0) = 3;
+ a(1, 1, 0) = 4;
+ ArrayND<int, 3> b(2, 2, 1);
+ b(0, 0, 0) = 6;
+ b(0, 1, 0) = 5;
+ b(1, 0, 0) = 4;
+ b(1, 1, 0) = 3;
+ ArrayND<int, 3> c;
+ MultiplyElements(a, b, &c);
+ EXPECT_FLOAT_EQ(6, c(0, 0, 0));
+ EXPECT_FLOAT_EQ(10, c(0, 1, 0));
+ EXPECT_FLOAT_EQ(12, c(1, 0, 0));
+ EXPECT_FLOAT_EQ(12, c(1, 1, 0));
+}
+
+TEST(ArrayND, IsEqualOperator) {
+ ArrayND<int, 3> a(2, 2, 1);
+ a(0, 0, 0) = 1;
+ a(0, 1, 0) = 2;
+ a(1, 0, 0) = 3;
+ a(1, 1, 0) = 4;
+ ArrayND<int, 3> b(2, 2, 1);
+ b(0, 0, 0) = 1;
+ b(0, 1, 0) = 2;
+ b(1, 0, 0) = 3;
+ b(1, 1, 0) = 4;
+ EXPECT_TRUE(a == b);
+ EXPECT_FALSE(a != b);
+ b(1, 1, 0) = 5;
+ EXPECT_TRUE(a != b);
+ EXPECT_FALSE(a == b);
+}
+
+TEST(Array3D, Sizes) {
+ Array3D<int> array;
+ EXPECT_EQ(array.Height(), 0);
+ EXPECT_EQ(array.Width(), 0);
+ EXPECT_EQ(array.Depth(), 0);
+ EXPECT_EQ(array.Shape(0), 0);
+}
+
+TEST(Array3D, CopyConstructor) {
+ Array3D<int> array(10, 10);
+ array(0, 0) = 1;
+ array(0, 1) = 1;
+ Array3D<int> copy(array);
+ EXPECT_EQ(copy.Height(), 10);
+ EXPECT_EQ(copy.Width(), 10);
+ EXPECT_EQ(copy.Depth(), 1);
+ EXPECT_EQ(copy(0, 0), 1);
+ copy(0, 1) = 2;
+ EXPECT_EQ(array(0, 1), 1);
+}
+
+TEST(Array3D, Assignation) {
+ Array3D<int> array(10, 10);
+ array(0, 0) = 1;
+ array(0, 1) = 1;
+ Array3D<int> copy;
+ copy = array;
+ EXPECT_EQ(copy.Height(), 10);
+ EXPECT_EQ(copy.Width(), 10);
+ EXPECT_EQ(copy.Depth(), 1);
+ EXPECT_EQ(copy(0, 0), 1);
+ copy(0, 1) = 2;
+ EXPECT_EQ(array(0, 1), 1);
+}
+
+TEST(Array3D, Parenthesis) {
+ Array3D<int> array(1, 2, 3);
+ array(0, 1, 0) = 3;
+ EXPECT_EQ(array(0, 1), 3);
+}
+
+TEST(Array3Df, SplitChannels) {
+ Array3Df array(1, 2, 3);
+ array(0, 0, 0) = 1;
+ array(0, 1, 0) = 1;
+ array(0, 0, 1) = 2;
+ array(0, 1, 1) = 2;
+ array(0, 0, 2) = 3;
+ array(0, 1, 2) = 3;
+ Array3Df c0, c1, c2;
+ SplitChannels(array, &c0, &c1, &c2);
+ for (int row = 0; row < 1; ++row) {
+ for (int column = 0; column < 2; ++column) {
+ EXPECT_EQ(array(row, column, 0), c0(row, column));
+ EXPECT_EQ(array(row, column, 1), c1(row, column));
+ EXPECT_EQ(array(row, column, 2), c2(row, column));
+ }
+ }
+}
+
+TEST(ArrayND, MultiplyElementsGeneric) {
+ ArrayND<double, 5> A;
+ ArrayND<int, 5> B;
+ ArrayND<double, 5> C;
+ int shape[] = {1, 3, 5, 7, 1};
+ A.Resize(shape);
+ B.Resize(shape);
+
+ A.Fill(1.1);
+ B.Fill(2);
+ MultiplyElements(A, B, &C);
+
+ ArrayND<double, 5>::Index cIndex;
+ for (int d0 = 0; d0 < shape[0]; ++d0)
+ for (int d1 = 0; d1 < shape[1]; ++d1)
+ for (int d2 = 0; d2 < shape[2]; ++d2)
+ for (int d3 = 0; d3 < shape[3]; ++d3)
+ for (int d4 = 0; d4 < shape[4]; ++d4) {
+ cIndex(0) = d0;
+ cIndex(1) = d1;
+ cIndex(2) = d2;
+ cIndex(3) = d3;
+ cIndex(4) = d4;
+ EXPECT_EQ(2.2, C(cIndex));
+ }
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/image/convolve_test.cc b/extern/libmv/libmv/image/convolve_test.cc
new file mode 100644
index 00000000000..0cdef8e1e72
--- /dev/null
+++ b/extern/libmv/libmv/image/convolve_test.cc
@@ -0,0 +1,110 @@
+// Copyright (c) 2007, 2008 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 <iostream>
+
+#include "libmv/image/convolve.h"
+#include "libmv/image/image.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+
+using namespace libmv;
+
+namespace {
+
+TEST(Convolve, ComputeGaussianKernel) {
+ Vec kernel, derivative;
+ ComputeGaussianKernel(1, &kernel, &derivative);
+ EXPECT_EQ(7, kernel.size());
+ // TODO(keir): Put in a more thorough test here!
+}
+
+TEST(Convolve, ConvolveGaussian) {
+ FloatImage im(10, 10);
+ im.Fill(1);
+ FloatImage blured;
+ ConvolveGaussian(im, 3, &blured);
+ EXPECT_NEAR(im(5, 5), 1, 1e-7);
+}
+
+TEST(Convolve, BoxFilterHorizontal) {
+ FloatImage im(10, 10), convolved, filtered;
+ im.Fill(1);
+ BoxFilterHorizontal(im, 3, &filtered);
+ Vec kernel(3);
+ kernel.setConstant(1.);
+ ConvolveHorizontal(im, kernel, &convolved);
+ EXPECT_EQ(filtered(5, 5), 3);
+ EXPECT_TRUE(filtered == convolved);
+}
+
+TEST(Convolve, BoxFilter) {
+ FloatImage image(5, 5), filtered;
+ // A single 1.0 inside a 5x5 image should expand to a 3x3 square.
+ image.Fill(0);
+ image(2, 2) = 1.0;
+ BoxFilter(image, 3, &filtered);
+ for (int j = 0; j < 5; j++) {
+ for (int i = 0; i < 5; i++) {
+ if (i == 0 || i == 4 || j == 0 || j == 4) {
+ EXPECT_EQ(0.0, filtered(j, i));
+ } else {
+ EXPECT_EQ(1.0, filtered(j, i));
+ }
+ }
+ }
+}
+
+TEST(Convolve, BlurredImageAndDerivativesChannelsFlat) {
+ FloatImage im(10, 10), blurred_and_derivatives;
+ im.Fill(1);
+ BlurredImageAndDerivativesChannels(im, 1.0, &blurred_and_derivatives);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 0), 1.0, 1e-7);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 1), 0.0, 1e-7);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 2), 0.0, 1e-7);
+}
+
+TEST(Convolve, BlurredImageAndDerivativesChannelsHorizontalSlope) {
+ FloatImage image(10, 10), blurred_and_derivatives;
+ for (int j = 0; j < 10; ++j) {
+ for (int i = 0; i < 10; ++i) {
+ image(j, i) = 2*i;
+ }
+ }
+ BlurredImageAndDerivativesChannels(image, 0.9, &blurred_and_derivatives);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 0), 10.0, 1e-7);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 1), 2.0, 1e-7);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 2), 0.0, 1e-7);
+}
+
+TEST(Convolve, BlurredImageAndDerivativesChannelsVerticalSlope) {
+ FloatImage image(10, 10), blurred_and_derivatives;
+ for (int j = 0; j < 10; ++j) {
+ for (int i = 0; i < 10; ++i) {
+ image(j, i) = 2*j;
+ }
+ }
+ BlurredImageAndDerivativesChannels(image, 0.9, &blurred_and_derivatives);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 0), 10.0, 1e-7);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 1), 0.0, 1e-7);
+ EXPECT_NEAR(blurred_and_derivatives(5, 5, 2), 2.0, 1e-7);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/image/image_drawing.h b/extern/libmv/libmv/image/image_drawing.h
new file mode 100644
index 00000000000..f50e48b75a3
--- /dev/null
+++ b/extern/libmv/libmv/image/image_drawing.h
@@ -0,0 +1,285 @@
+// 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.
+
+// Generic Image Processing Algorithm (GIPA)
+// Use an ImageModel class that must implement the following :
+//
+// ::Contains(int y, int x) <= Tell if a point is inside or not the image
+// ::operator(int y,int x) <= Modification accessor over the pixel (y,x)
+// ::Width()
+// ::Height()
+
+#ifndef LIBMV_IMAGE_IMAGE_DRAWING_H
+#define LIBMV_IMAGE_IMAGE_DRAWING_H
+
+namespace libmv {
+
+/// Put the pixel in the image to the given color only if the point (xc,yc)
+/// is inside the image.
+template <class Image, class Color>
+inline void safePutPixel(int yc, int xc, const Color & col, Image *pim) {
+ if (!pim)
+ return;
+ if (pim->Contains(yc, xc)) {
+ (*pim)(yc, xc) = col;
+ }
+}
+/// Put the pixel in the image to the given color only if the point (xc,yc)
+/// is inside the image. This function support multi-channel color
+/// \note The color pointer col must have size as the image depth
+template <class Image, class Color>
+inline void safePutPixel(int yc, int xc, const Color *col, Image *pim) {
+ if (!pim)
+ return;
+ if (pim->Contains(yc, xc)) {
+ for (int i = 0; i < pim->Depth(); ++i)
+ (*pim)(yc, xc, i) = *(col + i);
+ }
+}
+
+// Bresenham approach to draw ellipse.
+// http://raphaello.univ-fcomte.fr/ig/algorithme/ExemplesGLUt/BresenhamEllipse.htm
+// Add the rotation of the ellipse.
+// As the algo. use symmetry we must use 4 rotations.
+template <class Image, class Color>
+void DrawEllipse(int xc, int yc, int radiusA, int radiusB,
+ const Color &col, Image *pim, double angle = 0.0) {
+ int a = radiusA;
+ int b = radiusB;
+
+ // Counter Clockwise rotation matrix.
+ double matXY[4] = { cos(angle), sin(angle),
+ -sin(angle), cos(angle)};
+ int x, y;
+ double d1, d2;
+ x = 0;
+ y = b;
+ d1 = b*b - a*a*b + a*a/4;
+
+ float rotX = (matXY[0] * x + matXY[1] * y);
+ float rotY = (matXY[2] * x + matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (matXY[0] * x - matXY[1] * y);
+ rotY = (matXY[2] * x - matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (-matXY[0] * x - matXY[1] * y);
+ rotY = (-matXY[2] * x - matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (-matXY[0] * x + matXY[1] * y);
+ rotY = (-matXY[2] * x + matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+
+ while (a*a*(y-.5) > b*b*(x+1)) {
+ if (d1 < 0) {
+ d1 += b*b*(2*x+3);
+ ++x;
+ } else {
+ d1 += b*b*(2*x+3) + a*a*(-2*y+2);
+ ++x;
+ --y;
+ }
+ rotX = (matXY[0] * x + matXY[1] * y);
+ rotY = (matXY[2] * x + matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (matXY[0] * x - matXY[1] * y);
+ rotY = (matXY[2] * x - matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (-matXY[0] * x - matXY[1] * y);
+ rotY = (-matXY[2] * x - matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (-matXY[0] * x + matXY[1] * y);
+ rotY = (-matXY[2] * x + matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ }
+ d2 = b*b*(x+.5)*(x+.5) + a*a*(y-1)*(y-1) - a*a*b*b;
+ while (y > 0) {
+ if (d2 < 0) {
+ d2 += b*b*(2*x+2) + a*a*(-2*y+3);
+ --y;
+ ++x;
+ } else {
+ d2 += a*a*(-2*y+3);
+ --y;
+ }
+ rotX = (matXY[0] * x + matXY[1] * y);
+ rotY = (matXY[2] * x + matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (matXY[0] * x - matXY[1] * y);
+ rotY = (matXY[2] * x - matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (-matXY[0] * x - matXY[1] * y);
+ rotY = (-matXY[2] * x - matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ rotX = (-matXY[0] * x + matXY[1] * y);
+ rotY = (-matXY[2] * x + matXY[3] * y);
+ safePutPixel(yc + rotY, xc + rotX, col, pim);
+ }
+}
+
+// Bresenham approach do not allow to draw concentric circle without holes.
+// So it's better the use the Andres method.
+// http://fr.wikipedia.org/wiki/Algorithme_de_tracé_de_cercle_d'Andres.
+template <class Image, class Color>
+void DrawCircle(int x, int y, int radius, const Color &col, Image *pim) {
+ Image &im = *pim;
+ if ( im.Contains(y + radius, x + radius)
+ || im.Contains(y + radius, x - radius)
+ || im.Contains(y - radius, x + radius)
+ || im.Contains(y - radius, x - radius)) {
+ int x1 = 0;
+ int y1 = radius;
+ int d = radius - 1;
+ while (y1 >= x1) {
+ // Draw the point for each octant.
+ safePutPixel( y1 + y, x1 + x, col, pim);
+ safePutPixel( x1 + y, y1 + x, col, pim);
+ safePutPixel( y1 + y, -x1 + x, col, pim);
+ safePutPixel( x1 + y, -y1 + x, col, pim);
+ safePutPixel(-y1 + y, x1 + x, col, pim);
+ safePutPixel(-x1 + y, y1 + x, col, pim);
+ safePutPixel(-y1 + y, -x1 + x, col, pim);
+ safePutPixel(-x1 + y, -y1 + x, col, pim);
+ if (d >= 2 * x1) {
+ d = d - 2 * x1 - 1;
+ x1 += 1;
+ } else {
+ if (d <= 2 * (radius - y1)) {
+ d = d + 2 * y1 - 1;
+ y1 -= 1;
+ } else {
+ d = d + 2 * (y1 - x1 - 1);
+ y1 -= 1;
+ x1 += 1;
+ }
+ }
+ }
+ }
+}
+
+// Bresenham algorithm
+template <class Image, class Color>
+void DrawLine(int xa, int ya, int xb, int yb, const Color &col, Image *pim) {
+ Image &im = *pim;
+
+ // If one point is outside the image
+ // Replace the outside point by the intersection of the line and
+ // the limit (either x=width or y=height).
+ if (!im.Contains(ya, xa) || !im.Contains(yb, xb)) {
+ int width = pim->Width();
+ int height = pim->Height();
+ const bool xdir = xa < xb, ydir = ya < yb;
+ float nx0 = xa, nx1 = xb, ny0 = ya, ny1 = yb,
+ &xleft = xdir?nx0:nx1, &yleft = xdir?ny0:ny1,
+ &xright = xdir?nx1:nx0, &yright = xdir?ny1:ny0,
+ &xup = ydir?nx0:nx1, &yup = ydir?ny0:ny1,
+ &xdown = ydir?nx1:nx0, &ydown = ydir?ny1:ny0;
+
+ if (xright < 0 || xleft >= width) return;
+ if (xleft < 0) {
+ yleft -= xleft*(yright - yleft)/(xright - xleft);
+ xleft = 0;
+ }
+ if (xright >= width) {
+ yright -= (xright - width)*(yright - yleft)/(xright - xleft);
+ xright = width - 1;
+ }
+ if (ydown < 0 || yup >= height) return;
+ if (yup < 0) {
+ xup -= yup*(xdown - xup)/(ydown - yup);
+ yup = 0;
+ }
+ if (ydown >= height) {
+ xdown -= (ydown - height)*(xdown - xup)/(ydown - yup);
+ ydown = height - 1;
+ }
+
+ xa = (int) xleft;
+ xb = (int) xright;
+ ya = (int) yleft;
+ yb = (int) yright;
+ }
+
+ int xbas, xhaut, ybas, yhaut;
+ // Check the condition ybas < yhaut.
+ if (ya <= yb) {
+ xbas = xa;
+ ybas = ya;
+ xhaut = xb;
+ yhaut = yb;
+ } else {
+ xbas = xb;
+ ybas = yb;
+ xhaut = xa;
+ yhaut = ya;
+ }
+ // Initialize slope.
+ int x, y, dx, dy, incrmX, incrmY, dp, N, S;
+ dx = xhaut - xbas;
+ dy = yhaut - ybas;
+ if (dx > 0) { // If xhaut > xbas we will increment X.
+ incrmX = 1;
+ } else {
+ incrmX = -1; // else we will decrement X.
+ dx *= -1;
+ }
+ if (dy > 0) { // Positive slope will increment X.
+ incrmY = 1;
+ } else { // Negative slope.
+ incrmY = -1;
+ }
+ if (dx >= dy) {
+ dp = 2 * dy - dx;
+ S = 2 * dy;
+ N = 2 * (dy - dx);
+ y = ybas;
+ x = xbas;
+ while (x != xhaut) {
+ safePutPixel(y, x, col, pim);
+ x += incrmX;
+ if (dp <= 0) { // Go in direction of the South Pixel.
+ dp += S;
+ } else { // Go to the North.
+ dp += N;
+ y+=incrmY;
+ }
+ }
+ } else {
+ dp = 2 * dx - dy;
+ S = 2 * dx;
+ N = 2 * (dx - dy);
+ x = xbas;
+ y = ybas;
+ while (y < yhaut) {
+ safePutPixel(y, x, col, pim);
+ y += incrmY;
+ if (dp <= 0) { // Go in direction of the South Pixel.
+ dp += S;
+ } else { // Go to the North.
+ dp += N;
+ x += incrmX;
+ }
+ }
+ }
+ safePutPixel(y, x, col, pim);
+}
+
+} // namespace libmv
+
+#endif // LIBMV_IMAGE_IMAGE_DRAWING_H
diff --git a/extern/libmv/libmv/image/image_test.cc b/extern/libmv/libmv/image/image_test.cc
new file mode 100644
index 00000000000..241f49f2244
--- /dev/null
+++ b/extern/libmv/libmv/image/image_test.cc
@@ -0,0 +1,45 @@
+// Copyright (c) 2007, 2008 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 <iostream>
+
+#include "libmv/image/image.h"
+#include "testing/testing.h"
+
+using libmv::Image;
+using libmv::Array3Df;
+
+namespace {
+
+TEST(Image, SimpleImageAccessors) {
+ Array3Df *array = new Array3Df(2, 3);
+ Image image(array);
+ EXPECT_EQ(array, image.AsArray3Df());
+ EXPECT_TRUE(NULL == image.AsArray3Du());
+}
+
+TEST(Image, MemorySizeInBytes) {
+ Array3Df *array = new Array3Df(2, 3);
+ Image image(array);
+ int size = sizeof(image) + array->MemorySizeInBytes();
+ EXPECT_EQ(size, image.MemorySizeInBytes());
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/image/sample_test.cc b/extern/libmv/libmv/image/sample_test.cc
new file mode 100644
index 00000000000..c8a0ce470c2
--- /dev/null
+++ b/extern/libmv/libmv/image/sample_test.cc
@@ -0,0 +1,89 @@
+// Copyright (c) 2007, 2008 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/image/sample.h"
+#include "testing/testing.h"
+
+using namespace libmv;
+
+namespace {
+
+TEST(Image, Nearest) {
+ Array3Du image(2, 2);
+ image(0, 0) = 0;
+ image(0, 1) = 1;
+ image(1, 0) = 2;
+ image(1, 1) = 3;
+ EXPECT_EQ(0, SampleNearest(image, -0.4f, -0.4f));
+ EXPECT_EQ(0, SampleNearest(image, 0.4f, 0.4f));
+ EXPECT_EQ(3, SampleNearest(image, 0.6f, 0.6f));
+ EXPECT_EQ(3, SampleNearest(image, 1.4f, 1.4f));
+}
+
+TEST(Image, Linear) {
+ Array3Df image(2, 2);
+ image(0, 0) = 0;
+ image(0, 1) = 1;
+ image(1, 0) = 2;
+ image(1, 1) = 3;
+ EXPECT_EQ(1.5, SampleLinear(image, 0.5, 0.5));
+}
+
+TEST(Image, DownsampleBy2) {
+ Array3Df image(2, 2);
+ image(0, 0) = 0;
+ image(0, 1) = 1;
+ image(1, 0) = 2;
+ image(1, 1) = 3;
+ Array3Df resampled_image;
+ DownsampleChannelsBy2(image, &resampled_image);
+ ASSERT_EQ(1, resampled_image.Height());
+ ASSERT_EQ(1, resampled_image.Width());
+ ASSERT_EQ(1, resampled_image.Depth());
+ EXPECT_FLOAT_EQ(6./4., resampled_image(0, 0));
+}
+
+TEST(Image, DownsampleBy2MultiChannel) {
+ Array3Df image(2, 2, 3);
+ image(0, 0, 0) = 0;
+ image(0, 1, 0) = 1;
+ image(1, 0, 0) = 2;
+ image(1, 1, 0) = 3;
+
+ image(0, 0, 1) = 5;
+ image(0, 1, 1) = 6;
+ image(1, 0, 1) = 7;
+ image(1, 1, 1) = 8;
+
+ image(0, 0, 2) = 9;
+ image(0, 1, 2) = 10;
+ image(1, 0, 2) = 11;
+ image(1, 1, 2) = 12;
+
+ Array3Df resampled_image;
+ DownsampleChannelsBy2(image, &resampled_image);
+ ASSERT_EQ(1, resampled_image.Height());
+ ASSERT_EQ(1, resampled_image.Width());
+ ASSERT_EQ(3, resampled_image.Depth());
+ EXPECT_FLOAT_EQ((0+1+2+3)/4., resampled_image(0, 0, 0));
+ EXPECT_FLOAT_EQ((5+6+7+8)/4., resampled_image(0, 0, 1));
+ EXPECT_FLOAT_EQ((9+10+11+12)/4., resampled_image(0, 0, 2));
+}
+} // namespace
diff --git a/extern/libmv/libmv/image/tuple_test.cc b/extern/libmv/libmv/image/tuple_test.cc
new file mode 100644
index 00000000000..df44e5638b5
--- /dev/null
+++ b/extern/libmv/libmv/image/tuple_test.cc
@@ -0,0 +1,83 @@
+// Copyright (c) 2007, 2008 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/image/tuple.h"
+#include "testing/testing.h"
+
+#include <algorithm>
+
+using libmv::Tuple;
+
+namespace {
+
+TEST(Tuple, InitConstantValue) {
+ Tuple<int, 3> t(5);
+ EXPECT_EQ(t(0), 5);
+ EXPECT_EQ(t(0), 5);
+ EXPECT_EQ(t(0), 5);
+}
+
+TEST(Tuple, InitFromPointer) {
+ float vals[3] = {1.0f, 2.0f, 3.0f};
+
+ Tuple<float, 3> t(vals);
+ for (int i = 0; i < 3; i++)
+ EXPECT_EQ(t(i), vals[i]);
+
+ Tuple<int, 3> b(t);
+ EXPECT_EQ(b(0), int(vals[0]));
+ EXPECT_EQ(b(1), int(vals[1]));
+ EXPECT_EQ(b(2), int(vals[2]));
+}
+
+TEST(Tuple, Swap) {
+ unsigned char vala[3] = {1, 2, 3};
+ unsigned char valb[3] = {4, 5, 6};
+
+ Tuple<unsigned char, 3> a(vala);
+ Tuple<unsigned char, 3> b(valb);
+
+ std::swap(a, b);
+
+ EXPECT_EQ(a(0), int(valb[0]));
+ EXPECT_EQ(a(1), int(valb[1]));
+ EXPECT_EQ(a(2), int(valb[2]));
+ EXPECT_EQ(b(0), int(vala[0]));
+ EXPECT_EQ(b(1), int(vala[1]));
+ EXPECT_EQ(b(2), int(vala[2]));
+}
+
+TEST(Tuple, IsEqualOperator) {
+ Tuple<int, 3> a;
+ a(0) = 1;
+ a(1) = 2;
+ a(2) = 3;
+ Tuple<int, 3> b;
+ b(0) = 1;
+ b(1) = 2;
+ b(2) = 3;
+ EXPECT_TRUE(a == b);
+ EXPECT_FALSE(a != b);
+ b(1) = 5;
+ EXPECT_TRUE(a != b);
+ EXPECT_FALSE(a == b);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/multiview/euclidean_resection_test.cc b/extern/libmv/libmv/multiview/euclidean_resection_test.cc
new file mode 100644
index 00000000000..5ec9dbda3cf
--- /dev/null
+++ b/extern/libmv/libmv/multiview/euclidean_resection_test.cc
@@ -0,0 +1,237 @@
+// 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.
+
+#include "libmv/multiview/euclidean_resection.h"
+#include "libmv/numeric/numeric.h"
+#include "libmv/logging/logging.h"
+#include "libmv/multiview/projection.h"
+#include "testing/testing.h"
+
+using namespace libmv::euclidean_resection;
+using namespace libmv;
+
+// Generates all necessary inputs and expected outputs for EuclideanResection.
+void CreateCameraSystem(const Mat3& KK,
+ const Mat3X& x_image,
+ const Vec& X_distances,
+ const Mat3& R_input,
+ const Vec3& T_input,
+ Mat2X *x_camera,
+ Mat3X *X_world,
+ Mat3 *R_expected,
+ Vec3 *T_expected) {
+ int num_points = x_image.cols();
+
+ Mat3X x_unit_cam(3, num_points);
+ x_unit_cam = KK.inverse() * x_image;
+
+ // Create normalized camera coordinates to be used as an input to the PnP
+ // function, instead of using NormalizeColumnVectors(&x_unit_cam).
+ *x_camera = x_unit_cam.block(0, 0, 2, num_points);
+ for (int i = 0; i < num_points; ++i) {
+ x_unit_cam.col(i).normalize();
+ }
+
+ // Create the 3D points in the camera system.
+ Mat X_camera(3, num_points);
+ for (int i = 0; i < num_points; ++i) {
+ X_camera.col(i) = X_distances(i) * x_unit_cam.col(i);
+ }
+
+ // Apply the transformation to the camera 3D points
+ Mat translation_matrix(3, num_points);
+ translation_matrix.row(0).setConstant(T_input(0));
+ translation_matrix.row(1).setConstant(T_input(1));
+ translation_matrix.row(2).setConstant(T_input(2));
+
+ *X_world = R_input * X_camera + translation_matrix;
+
+ // Create the expected result for comparison.
+ *R_expected = R_input.transpose();
+ *T_expected = *R_expected * (-T_input);
+};
+
+TEST(AbsoluteOrientation, QuaternionSolution) {
+ int num_points = 4;
+ Mat X;
+ Mat Xp;
+ X = 100 * Mat::Random(3, num_points);
+
+ // Create a random translation and rotation.
+ Mat3 R_input;
+ R_input = Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitY())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ());
+
+ Vec3 t_input;
+ t_input.setRandom();
+ t_input = 100 * t_input;
+
+ Mat translation_matrix(3, num_points);
+ translation_matrix.row(0).setConstant(t_input(0));
+ translation_matrix.row(1).setConstant(t_input(1));
+ translation_matrix.row(2).setConstant(t_input(2));
+
+ // Create the transformed 3D points Xp as Xp = R * X + t.
+ Xp = R_input * X + translation_matrix;
+
+ // Output variables.
+ Mat3 R;
+ Vec3 t;
+
+ AbsoluteOrientation(X, Xp, &R, &t);
+
+ EXPECT_MATRIX_NEAR(t, t_input, 1e-6);
+ EXPECT_MATRIX_NEAR(R, R_input, 1e-8);
+}
+
+TEST(EuclideanResection, Points4KnownImagePointsRandomTranslationRotation) {
+ // In this test only the translation and rotation are random. The image
+ // points are selected from a real case and are well conditioned.
+ Vec2i image_dimensions;
+ image_dimensions << 1600, 1200;
+
+ Mat3 KK;
+ KK << 2796, 0, 804,
+ 0 , 2796, 641,
+ 0, 0, 1;
+
+ // The real image points.
+ int num_points = 4;
+ Mat3X x_image(3, num_points);
+ x_image << 1164.06, 734.948, 749.599, 430.727,
+ 681.386, 844.59, 496.315, 580.775,
+ 1, 1, 1, 1;
+
+
+ // A vector of the 4 distances to the 3D points.
+ Vec X_distances = 100 * Vec::Random(num_points).array().abs();
+
+ // Create the random camera motion R and t that resection should recover.
+ Mat3 R_input;
+ R_input = Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitY())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ());
+
+ Vec3 T_input;
+ T_input.setRandom();
+ T_input = 100 * T_input;
+
+ // Create the camera system, also getting the expected result of the
+ // transformation.
+ Mat3 R_expected;
+ Vec3 T_expected;
+ Mat3X X_world;
+ Mat2X x_camera;
+ CreateCameraSystem(KK, x_image, X_distances, R_input, T_input,
+ &x_camera, &X_world, &R_expected, &T_expected);
+
+ // Finally, run the code under test.
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_ANSAR_DANIILIDIS);
+
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+
+ // For now, the EPnP doesn't have a non-linear optimization step and so is
+ // not precise enough with only 4 points.
+ //
+ // TODO(jmichot): Reenable this test when there is nonlinear refinement.
+#if 0
+ R_output.setIdentity();
+ T_output.setZero();
+
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_EPNP);
+
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);*/
+#endif
+}
+
+// TODO(jmichot): Reduce the code duplication here with the code above.
+TEST(EuclideanResection, Points6AllRandomInput) {
+ Mat3 KK;
+ KK << 2796, 0, 804,
+ 0 , 2796, 641,
+ 0, 0, 1;
+
+ // Create random image points for a 1600x1200 image.
+ int w = 1600;
+ int h = 1200;
+ int num_points = 6;
+ Mat3X x_image(3, num_points);
+ x_image.row(0) = w * Vec::Random(num_points).array().abs();
+ x_image.row(1) = h * Vec::Random(num_points).array().abs();
+ x_image.row(2).setOnes();
+
+ // Normalized camera coordinates to be used as an input to the PnP function.
+ Mat2X x_camera;
+ Vec X_distances = 100 * Vec::Random(num_points).array().abs();
+
+ // Create the random camera motion R and t that resection should recover.
+ Mat3 R_input;
+ R_input = Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitY())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ());
+
+ Vec3 T_input;
+ T_input.setRandom();
+ T_input = 100 * T_input;
+
+ // Create the camera system.
+ Mat3 R_expected;
+ Vec3 T_expected;
+ Mat3X X_world;
+ CreateCameraSystem(KK, x_image, X_distances, R_input, T_input,
+ &x_camera, &X_world, &R_expected, &T_expected);
+
+ // Test each of the resection methods.
+ {
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_ANSAR_DANIILIDIS);
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+ }
+ {
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_EPNP);
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+ }
+ {
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_image, X_world, KK,
+ &R_output, &T_output);
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+ }
+}
diff --git a/extern/libmv/libmv/multiview/fundamental_test.cc b/extern/libmv/libmv/multiview/fundamental_test.cc
new file mode 100644
index 00000000000..da0eb449b8f
--- /dev/null
+++ b/extern/libmv/libmv/multiview/fundamental_test.cc
@@ -0,0 +1,162 @@
+// Copyright (c) 2007, 2008 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 <iostream>
+
+#include "libmv/logging/logging.h"
+#include "libmv/multiview/conditioning.h"
+#include "libmv/multiview/fundamental.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/multiview/test_data_sets.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+
+namespace {
+
+using namespace libmv;
+
+TEST(Fundamental, FundamentalFromProjections) {
+ Mat34 P1_gt, P2_gt;
+ P1_gt << 1, 0, 0, 0,
+ 0, 1, 0, 0,
+ 0, 0, 1, 0;
+ P2_gt << 1, 1, 1, 3,
+ 0, 2, 0, 3,
+ 0, 1, 1, 0;
+ Mat3 F_gt;
+ FundamentalFromProjections(P1_gt, P2_gt, &F_gt);
+
+ Mat34 P1, P2;
+ ProjectionsFromFundamental(F_gt, &P1, &P2);
+
+ Mat3 F;
+ FundamentalFromProjections(P1, P2, &F);
+
+ EXPECT_MATRIX_PROP(F_gt, F, 1e-6);
+}
+
+TEST(Fundamental, PreconditionerFromPoints) {
+ int n = 4;
+ Mat points(2, n);
+ points << 0, 0, 1, 1,
+ 0, 2, 1, 3;
+
+ Mat3 T;
+ PreconditionerFromPoints(points, &T);
+
+ Mat normalized_points;
+ ApplyTransformationToPoints(points, T, &normalized_points);
+
+ Vec mean, variance;
+ MeanAndVarianceAlongRows(normalized_points, &mean, &variance);
+
+ EXPECT_NEAR(0, mean(0), 1e-8);
+ EXPECT_NEAR(0, mean(1), 1e-8);
+ EXPECT_NEAR(2, variance(0), 1e-8);
+ EXPECT_NEAR(2, variance(1), 1e-8);
+}
+
+TEST(Fundamental, EssentialFromFundamental) {
+ TwoViewDataSet d = TwoRealisticCameras();
+
+ Mat3 E_from_Rt;
+ EssentialFromRt(d.R1, d.t1, d.R2, d.t2, &E_from_Rt);
+
+ Mat3 E_from_F;
+ EssentialFromFundamental(d.F, d.K1, d.K2, &E_from_F);
+
+ EXPECT_MATRIX_PROP(E_from_Rt, E_from_F, 1e-6);
+}
+
+TEST(Fundamental, MotionFromEssential) {
+ TwoViewDataSet d = TwoRealisticCameras();
+
+ Mat3 E;
+ EssentialFromRt(d.R1, d.t1, d.R2, d.t2, &E);
+
+ Mat3 R;
+ Vec3 t;
+ RelativeCameraMotion(d.R1, d.t1, d.R2, d.t2, &R, &t);
+ NormalizeL2(&t);
+
+ std::vector<Mat3> Rs;
+ std::vector<Vec3> ts;
+ MotionFromEssential(E, &Rs, &ts);
+ bool one_solution_is_correct = false;
+ for (size_t i = 0; i < Rs.size(); ++i) {
+ if (FrobeniusDistance(Rs[i], R) < 1e-8 && DistanceL2(ts[i], t) < 1e-8) {
+ one_solution_is_correct = true;
+ break;
+ }
+ }
+ EXPECT_TRUE(one_solution_is_correct);
+}
+
+TEST(Fundamental, MotionFromEssentialChooseSolution) {
+ TwoViewDataSet d = TwoRealisticCameras();
+
+ Mat3 E;
+ EssentialFromRt(d.R1, d.t1, d.R2, d.t2, &E);
+
+ Mat3 R;
+ Vec3 t;
+ RelativeCameraMotion(d.R1, d.t1, d.R2, d.t2, &R, &t);
+ NormalizeL2(&t);
+
+ std::vector<Mat3> Rs;
+ std::vector<Vec3> ts;
+ MotionFromEssential(E, &Rs, &ts);
+
+ Vec2 x1, x2;
+ MatrixColumn(d.x1, 0, &x1);
+ MatrixColumn(d.x2, 0, &x2);
+ int solution = MotionFromEssentialChooseSolution(Rs, ts, d.K1, x1, d.K2, x2);
+
+ EXPECT_LE(0, solution);
+ EXPECT_LE(solution, 3);
+ EXPECT_LE(FrobeniusDistance(Rs[solution], R), 1e-8);
+ EXPECT_LE(DistanceL2(ts[solution], t), 1e-8);
+}
+
+TEST(Fundamental, MotionFromEssentialAndCorrespondence) {
+ TwoViewDataSet d = TwoRealisticCameras();
+
+ Mat3 E;
+ EssentialFromRt(d.R1, d.t1, d.R2, d.t2, &E);
+
+ Mat3 R;
+ Vec3 t;
+ RelativeCameraMotion(d.R1, d.t1, d.R2, d.t2, &R, &t);
+ NormalizeL2(&t);
+
+ Vec2 x1, x2;
+ MatrixColumn(d.x1, 0, &x1);
+ MatrixColumn(d.x2, 0, &x2);
+
+ Mat3 R_estimated;
+ Vec3 t_estimated;
+ MotionFromEssentialAndCorrespondence(E, d.K1, x1, d.K2, x2,
+ &R_estimated, &t_estimated);
+
+ EXPECT_LE(FrobeniusDistance(R_estimated, R), 1e-8);
+ EXPECT_LE(DistanceL2(t_estimated, t), 1e-8);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/multiview/homography_error.h b/extern/libmv/libmv/multiview/homography_error.h
new file mode 100644
index 00000000000..f8b9d45e73c
--- /dev/null
+++ b/extern/libmv/libmv/multiview/homography_error.h
@@ -0,0 +1,248 @@
+// 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_MULTIVIEW_HOMOGRAPHY_ERRORS_H_
+#define LIBMV_MULTIVIEW_HOMOGRAPHY_ERRORS_H_
+
+#include "libmv/multiview/projection.h"
+
+namespace libmv {
+namespace homography {
+namespace homography2D {
+
+ /**
+ * Structure for estimating the asymmetric error between a vector x2 and the
+ * transformed x1 such that
+ * Error = ||x2 - Psi(H * x1)||^2
+ * where Psi is the function that transforms homogeneous to euclidean coords.
+ * \note It should be distributed as Chi-squared with k = 2.
+ */
+struct AsymmetricError {
+ /**
+ * Computes the asymmetric residuals between a set of 2D points x2 and the
+ * transformed 2D point set x1 such that
+ * Residuals_i = x2_i - Psi(H * x1_i)
+ * where Psi is the function that transforms homogeneous to euclidean coords.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \param[in] x2 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \param[out] dx A 2xN matrix of column vectors of residuals errors
+ */
+ static void Residuals(const Mat &H, const Mat &x1,
+ const Mat &x2, Mat2X *dx) {
+ dx->resize(2, x1.cols());
+ Mat3X x2h_est;
+ if (x1.rows() == 2)
+ x2h_est = H * EuclideanToHomogeneous(static_cast<Mat2X>(x1));
+ else
+ x2h_est = H * x1;
+ dx->row(0) = x2h_est.row(0).array() / x2h_est.row(2).array();
+ dx->row(1) = x2h_est.row(1).array() / x2h_est.row(2).array();
+ if (x2.rows() == 2)
+ *dx = x2 - *dx;
+ else
+ *dx = HomogeneousToEuclidean(static_cast<Mat3X>(x2)) - *dx;
+ }
+ /**
+ * Computes the asymmetric residuals between a 2D point x2 and the transformed
+ * 2D point x1 such that
+ * Residuals = x2 - Psi(H * x1)
+ * where Psi is the function that transforms homogeneous to euclidean coords.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \param[in] x2 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \param[out] dx A vector of size 2 of the residual error
+ */
+ static void Residuals(const Mat &H, const Vec &x1,
+ const Vec &x2, Vec2 *dx) {
+ Vec3 x2h_est;
+ if (x1.rows() == 2)
+ x2h_est = H * EuclideanToHomogeneous(static_cast<Vec2>(x1));
+ else
+ x2h_est = H * x1;
+ if (x2.rows() == 2)
+ *dx = x2 - x2h_est.head<2>() / x2h_est[2];
+ else
+ *dx = HomogeneousToEuclidean(static_cast<Vec3>(x2)) -
+ x2h_est.head<2>() / x2h_est[2];
+ }
+ /**
+ * Computes the squared norm of the residuals between a set of 2D points x2
+ * and the transformed 2D point set x1 such that
+ * Error = || x2 - Psi(H * x1) ||^2
+ * where Psi is the function that transforms homogeneous to euclidean coords.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \param[in] x2 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \return The squared norm of the asymmetric residuals errors
+ */
+ static double Error(const Mat &H, const Mat &x1, const Mat &x2) {
+ Mat2X dx;
+ Residuals(H, x1, x2, &dx);
+ return dx.squaredNorm();
+ }
+ /**
+ * Computes the squared norm of the residuals between a 2D point x2 and the
+ * transformed 2D point x1 such that rms = || x2 - Psi(H * x1) ||^2
+ * where Psi is the function that transforms homogeneous to euclidean coords.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \param[in] x2 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \return The squared norm of the asymmetric residual error
+ */
+ static double Error(const Mat &H, const Vec &x1, const Vec &x2) {
+ Vec2 dx;
+ Residuals(H, x1, x2, &dx);
+ return dx.squaredNorm();
+ }
+};
+
+ /**
+ * Structure for estimating the symmetric error
+ * between a vector x2 and the transformed x1 such that
+ * Error = ||x2 - Psi(H * x1)||^2 + ||x1 - Psi(H^-1 * x2)||^2
+ * where Psi is the function that transforms homogeneous to euclidean coords.
+ * \note It should be distributed as Chi-squared with k = 4.
+ */
+struct SymmetricError {
+ /**
+ * Computes the squared norm of the residuals between x2 and the
+ * transformed x1 such that
+ * Error = ||x2 - Psi(H * x1)||^2 + ||x1 - Psi(H^-1 * x2)||^2
+ * where Psi is the function that transforms homogeneous to euclidean coords.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \param[in] x2 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \return The squared norm of the symmetric residuals errors
+ */
+ static double Error(const Mat &H, const Vec &x1, const Vec &x2) {
+ // TODO(keir): This is awesomely inefficient because it does a 3x3
+ // inversion for each evaluation.
+ Mat3 Hinv = H.inverse();
+ return AsymmetricError::Error(H, x1, x2) +
+ AsymmetricError::Error(Hinv, x2, x1);
+ }
+ // TODO(julien) Add residuals function \see AsymmetricError
+};
+ /**
+ * Structure for estimating the algebraic error (cross product)
+ * between a vector x2 and the transformed x1 such that
+ * Error = ||[x2] * H * x1||^^2
+ * where [x2] is the skew matrix of x2.
+ */
+struct AlgebraicError {
+ // TODO(julien) Make an AlgebraicError2Rows and AlgebraicError3Rows
+
+ /**
+ * Computes the algebraic residuals (cross product) between a set of 2D
+ * points x2 and the transformed 2D point set x1 such that
+ * [x2] * H * x1 where [x2] is the skew matrix of x2.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \param[in] x2 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \param[out] dx A 3xN matrix of column vectors of residuals errors
+ */
+ static void Residuals(const Mat &H, const Mat &x1,
+ const Mat &x2, Mat3X *dx) {
+ dx->resize(3, x1.cols());
+ Vec3 col;
+ for (int i = 0; i < x1.cols(); ++i) {
+ Residuals(H, x1.col(i), x2.col(i), &col);
+ dx->col(i) = col;
+ }
+ }
+ /**
+ * Computes the algebraic residuals (cross product) between a 2D point x2
+ * and the transformed 2D point x1 such that
+ * [x2] * H * x1 where [x2] is the skew matrix of x2.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \param[in] x2 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \param[out] dx A vector of size 3 of the residual error
+ */
+ static void Residuals(const Mat &H, const Vec &x1,
+ const Vec &x2, Vec3 *dx) {
+ Vec3 x2h_est;
+ if (x1.rows() == 2)
+ x2h_est = H * EuclideanToHomogeneous(static_cast<Vec2>(x1));
+ else
+ x2h_est = H * x1;
+ if (x2.rows() == 2)
+ *dx = SkewMat(EuclideanToHomogeneous(static_cast<Vec2>(x2))) * x2h_est;
+ else
+ *dx = SkewMat(x2) * x2h_est;
+ // TODO(julien) This is inefficient since it creates an
+ // identical 3x3 skew matrix for each evaluation.
+ }
+ /**
+ * Computes the squared norm of the algebraic residuals between a set of 2D
+ * points x2 and the transformed 2D point set x1 such that
+ * [x2] * H * x1 where [x2] is the skew matrix of x2.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \param[in] x2 A set of 2D points (2xN or 3xN matrix of column vectors).
+ * \return The squared norm of the asymmetric residuals errors
+ */
+ static double Error(const Mat &H, const Mat &x1, const Mat &x2) {
+ Mat3X dx;
+ Residuals(H, x1, x2, &dx);
+ return dx.squaredNorm();
+ }
+ /**
+ * Computes the squared norm of the algebraic residuals between a 2D point x2
+ * and the transformed 2D point x1 such that
+ * [x2] * H * x1 where [x2] is the skew matrix of x2.
+ *
+ * \param[in] H The 3x3 homography matrix.
+ * The estimated homography should approximatelly hold the condition y = H x.
+ * \param[in] x1 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \param[in] x2 A 2D point (vector of size 2 or 3 (euclidean/homogeneous))
+ * \return The squared norm of the asymmetric residual error
+ */
+ static double Error(const Mat &H, const Vec &x1, const Vec &x2) {
+ Vec3 dx;
+ Residuals(H, x1, x2, &dx);
+ return dx.squaredNorm();
+ }
+};
+// TODO(keir): Add error based on ideal points.
+
+} // namespace homography2D
+// TODO(julien) add homography3D errors
+} // namespace homography
+} // namespace libmv
+
+#endif // LIBMV_MULTIVIEW_HOMOGRAPHY_ERRORS_H_
diff --git a/extern/libmv/libmv/multiview/homography_test.cc b/extern/libmv/libmv/multiview/homography_test.cc
new file mode 100644
index 00000000000..8d7266e3d11
--- /dev/null
+++ b/extern/libmv/libmv/multiview/homography_test.cc
@@ -0,0 +1,261 @@
+// 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 "testing/testing.h"
+#include "libmv/logging/logging.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/multiview/homography.h"
+
+namespace {
+using namespace libmv;
+
+namespace {
+
+// Check whether homography transform M actually transforms
+// given vectors x1 to x2. Used to check validness of a reconstructed
+// homography matrix.
+// TODO(sergey): Consider using this in all tests since possible homography
+// matrix is not fixed to a single value and different-looking matrices
+// might actually crrespond to the same exact transform.
+void CheckHomography2DTransform(const Mat3 &H,
+ const Mat &x1,
+ const Mat &x2) {
+ for (int i = 0; i < x2.cols(); ++i) {
+ Vec3 x2_expected = x2.col(i);
+ Vec3 x2_observed = H * x1.col(i);
+ x2_observed /= x2_observed(2);
+ EXPECT_MATRIX_NEAR(x2_expected, x2_observed, 1e-8);
+ }
+}
+
+} // namespace
+
+TEST(Homography2DTest, Rotation45AndTranslationXY) {
+ Mat x1(3, 4);
+ x1 << 0, 1, 0, 5,
+ 0, 0, 2, 3,
+ 1, 1, 1, 1;
+
+ double angle = 45.0;
+ Mat3 m;
+ m << cos(angle), -sin(angle), -2,
+ sin(angle), cos(angle), 5,
+ 0, 0, 1;
+
+ Mat x2 = x1;
+ // Transform point from ground truth matrix
+ for (int i = 0; i < x2.cols(); ++i)
+ x2.col(i) = m * x1.col(i);
+
+ Mat3 homography_mat;
+ EXPECT_TRUE(Homography2DFromCorrespondencesLinear(x1, x2, &homography_mat));
+ VLOG(1) << "Mat Homography2D ";
+ VLOG(1) << homography_mat;
+ VLOG(1) << "Mat GT ";
+ VLOG(1) << m;
+ EXPECT_MATRIX_NEAR(homography_mat, m, 1e-8);
+}
+
+TEST(Homography2DTest, AffineGeneral4) {
+ // TODO(julien) find why it doesn't work with 4 points!!!
+ Mat x1(3, 4);
+ x1 << 0, 1, 0, 2,
+ 0, 0, 1, 2,
+ 1, 1, 1, 1;
+ Mat3 m;
+ m << 3, -1, 4,
+ 6, -2, -3,
+ 0, 0, 1;
+
+ Mat x2 = x1;
+ for (int i = 0; i < x2.cols(); ++i) {
+ x2.col(i) = m * x1.col(i);
+ }
+
+ Mat3 homography_mat;
+ EXPECT_TRUE(Homography2DFromCorrespondencesLinear(x1, x2, &homography_mat));
+ VLOG(1) << "Mat Homography2D";
+ VLOG(1) << homography_mat;
+ CheckHomography2DTransform(homography_mat, x1, x2);
+
+ // Test with euclidean coordinates
+ Mat eX1, eX2;
+ HomogeneousToEuclidean(x1, &eX1);
+ HomogeneousToEuclidean(x2, &eX2);
+ homography_mat.setIdentity();
+ EXPECT_TRUE(Homography2DFromCorrespondencesLinear(eX1, eX2, &homography_mat));
+
+ VLOG(1) << "Mat Homography2D ";
+ VLOG(1) << homography_mat;
+ CheckHomography2DTransform(homography_mat, x1, x2);
+}
+
+TEST(Homography2DTest, AffineGeneral5) {
+ Mat x1(3, 5);
+ x1 << 0, 1, 0, 2, 5,
+ 0, 0, 1, 2, 2,
+ 1, 1, 1, 1, 1;
+ Mat3 m;
+ m << 3, -1, 4,
+ 6, -2, -3,
+ 0, 0, 1;
+
+ Mat x2 = x1;
+ for (int i = 0; i < x2.cols(); ++i)
+ x2.col(i) = m * x1.col(i);
+
+ Mat3 homography_mat;
+ EXPECT_TRUE(Homography2DFromCorrespondencesLinear(x1, x2, &homography_mat));
+
+ VLOG(1) << "Mat Homography2D ";
+ VLOG(1) << homography_mat;
+ EXPECT_MATRIX_NEAR(homography_mat, m, 1e-8);
+
+ // Test with euclidean coordinates
+ Mat eX1, eX2;
+ HomogeneousToEuclidean(x1, &eX1);
+ HomogeneousToEuclidean(x2, &eX2);
+ homography_mat.setIdentity();
+ EXPECT_TRUE(Homography2DFromCorrespondencesLinear(eX1, eX2, &homography_mat));
+
+ VLOG(1) << "Mat Homography2D ";
+ VLOG(1) << homography_mat;
+ EXPECT_MATRIX_NEAR(homography_mat, m, 1e-8);
+}
+
+TEST(Homography2DTest, HomographyGeneral) {
+ Mat x1(3, 4);
+ x1 << 0, 1, 0, 5,
+ 0, 0, 2, 3,
+ 1, 1, 1, 1;
+ Mat3 m;
+ m << 3, -1, 4,
+ 6, -2, -3,
+ 1, -3, 1;
+
+ Mat x2 = x1;
+ for (int i = 0; i < x2.cols(); ++i)
+ x2.col(i) = m * x1.col(i);
+
+ Mat3 homography_mat;
+ EXPECT_TRUE(Homography2DFromCorrespondencesLinear(x1, x2, &homography_mat));
+
+ VLOG(1) << "Mat Homography2D ";
+ VLOG(1) << homography_mat;
+ EXPECT_MATRIX_NEAR(homography_mat, m, 1e-8);
+}
+
+TEST(Homography3DTest, RotationAndTranslationXYZ) {
+ Mat x1(4, 5);
+ x1 << 0, 0, 1, 5, 2,
+ 0, 1, 2, 3, 5,
+ 0, 2, 0, 1, 5,
+ 1, 1, 1, 1, 1;
+ Mat4 M;
+ M.setIdentity();
+ /*
+ M = AngleAxisd(45.0, Vector3f::UnitZ())
+ * AngleAxisd(25.0, Vector3f::UnitX())
+ * AngleAxisd(5.0, Vector3f::UnitZ());*/
+
+ // Rotation on x + translation
+ double angle = 45.0;
+ Mat4 rot;
+ rot << 1, 0, 0, 1,
+ 0, cos(angle), -sin(angle), 3,
+ 0, sin(angle), cos(angle), -2,
+ 0, 0, 0, 1;
+ M *= rot;
+ // Rotation on y
+ angle = 25.0;
+ rot << cos(angle), 0, sin(angle), 0,
+ 0, 1, 0, 0,
+ -sin(angle), 0, cos(angle), 0,
+ 0, 0, 0, 1;
+ M *= rot;
+ // Rotation on z
+ angle = 5.0;
+ rot << cos(angle), -sin(angle), 0, 0,
+ sin(angle), cos(angle), 0, 0,
+ 0, 0, 1, 0,
+ 0, 0, 0, 1;
+ M *= rot;
+ Mat x2 = x1;
+ for (int i = 0; i < x2.cols(); ++i) {
+ x2.col(i) = M * x1.col(i);
+ }
+
+ Mat4 homography_mat;
+ EXPECT_TRUE(Homography3DFromCorrespondencesLinear(x1, x2, &homography_mat));
+
+ VLOG(1) << "Mat Homography3D " << homography_mat;
+ VLOG(1) << "Mat GT " << M;
+ EXPECT_MATRIX_NEAR(homography_mat, M, 1e-8);
+}
+
+TEST(Homography3DTest, AffineGeneral) {
+ Mat x1(4, 5);
+ x1 << 0, 0, 1, 5, 2,
+ 0, 1, 2, 3, 5,
+ 0, 2, 0, 1, 5,
+ 1, 1, 1, 1, 1;
+ Mat4 m;
+ m << 3, -1, 4, 1,
+ 6, -2, -3, -6,
+ 1, 0, 1, 2,
+ 0, 0, 0, 1;
+
+ Mat x2 = x1;
+ for (int i = 0; i < x2.cols(); ++i) {
+ x2.col(i) = m * x1.col(i);
+ }
+
+ Mat4 homography_mat;
+ EXPECT_TRUE(Homography3DFromCorrespondencesLinear(x1, x2, &homography_mat));
+ VLOG(1) << "Mat Homography3D ";
+ VLOG(1) << homography_mat;
+ EXPECT_MATRIX_NEAR(homography_mat, m, 1e-8);
+}
+
+TEST(Homography3DTest, HomographyGeneral) {
+ Mat x1(4, 5);
+ x1 << 0, 0, 1, 5, 2,
+ 0, 1, 2, 3, 5,
+ 0, 2, 0, 1, 5,
+ 1, 1, 1, 1, 1;
+ Mat4 m;
+ m << 3, -1, 4, 1,
+ 6, -2, -3, -6,
+ 1, 0, 1, 2,
+ -3, 1, 0, 1;
+
+ Mat x2 = x1;
+ for (int i = 0; i < x2.cols(); ++i) {
+ x2.col(i) = m * x1.col(i);
+ }
+
+ Mat4 homography_mat;
+ EXPECT_TRUE(Homography3DFromCorrespondencesLinear(x1, x2, &homography_mat));
+ VLOG(1) << "Mat Homography3D";
+ VLOG(1) << homography_mat;
+ EXPECT_MATRIX_NEAR(homography_mat, m, 1e-8);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/multiview/nviewtriangulation_test.cc b/extern/libmv/libmv/multiview/nviewtriangulation_test.cc
new file mode 100644
index 00000000000..5a4d8499753
--- /dev/null
+++ b/extern/libmv/libmv/multiview/nviewtriangulation_test.cc
@@ -0,0 +1,94 @@
+// 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.
+
+#include <iostream>
+
+#include "libmv/base/vector.h"
+#include "libmv/logging/logging.h"
+#include "libmv/multiview/nviewtriangulation.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/multiview/test_data_sets.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+
+namespace {
+
+using namespace libmv;
+
+TEST(NViewTriangulate, FiveViews) {
+ int nviews = 5;
+ int npoints = 6;
+ NViewDataSet d = NRealisticCamerasFull(nviews, npoints);
+
+ // Collect P matrices together.
+ vector<Mat34> Ps(nviews);
+ for (int j = 0; j < nviews; ++j) {
+ Ps[j] = d.P(j);
+ }
+
+ for (int i = 0; i < npoints; ++i) {
+ // Collect the image of point i in each frame.
+ Mat2X xs(2, nviews);
+ for (int j = 0; j < nviews; ++j) {
+ xs.col(j) = d.x[j].col(i);
+ }
+ Vec4 X;
+ NViewTriangulate(xs, Ps, &X);
+
+ // Check reprojection error. Should be nearly zero.
+ for (int j = 0; j < nviews; ++j) {
+ Vec3 x_reprojected = Ps[j]*X;
+ x_reprojected /= x_reprojected(2);
+ double error = (x_reprojected.head(2) - xs.col(j)).norm();
+ EXPECT_NEAR(error, 0.0, 1e-9);
+ }
+ }
+}
+
+TEST(NViewTriangulateAlgebraic, FiveViews) {
+ int nviews = 5;
+ int npoints = 6;
+ NViewDataSet d = NRealisticCamerasFull(nviews, npoints);
+
+ // Collect P matrices together.
+ vector<Mat34> Ps(nviews);
+ for (int j = 0; j < nviews; ++j) {
+ Ps[j] = d.P(j);
+ }
+
+ for (int i = 0; i < npoints; ++i) {
+ // Collect the image of point i in each frame.
+ Mat2X xs(2, nviews);
+ for (int j = 0; j < nviews; ++j) {
+ xs.col(j) = d.x[j].col(i);
+ }
+ Vec4 X;
+ NViewTriangulate(xs, Ps, &X);
+
+ // Check reprojection error. Should be nearly zero.
+ for (int j = 0; j < nviews; ++j) {
+ Vec3 x_reprojected = Ps[j]*X;
+ x_reprojected /= x_reprojected(2);
+ double error = (x_reprojected.head<2>() - xs.col(j)).norm();
+ EXPECT_NEAR(error, 0.0, 1e-9);
+ }
+ }
+}
+} // namespace
diff --git a/extern/libmv/libmv/multiview/panography_kernel.cc b/extern/libmv/libmv/multiview/panography_kernel.cc
new file mode 100644
index 00000000000..8fdc9e79aed
--- /dev/null
+++ b/extern/libmv/libmv/multiview/panography_kernel.cc
@@ -0,0 +1,51 @@
+// Copyright (c) 2008, 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.
+
+#include "libmv/multiview/panography_kernel.h"
+#include "libmv/multiview/panography.h"
+
+namespace libmv {
+namespace panography {
+namespace kernel {
+
+void TwoPointSolver::Solve(const Mat &x1, const Mat &x2, vector<Mat3> *Hs) {
+ // Solve for the focal lengths.
+ vector<double> fs;
+ F_FromCorrespondance_2points(x1, x2, &fs);
+
+ // Then solve for the rotations and homographies.
+ Mat x1h, x2h;
+ EuclideanToHomogeneous(x1, &x1h);
+ EuclideanToHomogeneous(x2, &x2h);
+ for (int i = 0; i < fs.size(); ++i) {
+ Mat3 K1 = Mat3::Identity() * fs[i];
+ K1(2, 2) = 1.0;
+
+ Mat3 R;
+ GetR_FixedCameraCenter(x1h, x2h, fs[i], &R);
+ R /= R(2, 2);
+
+ (*Hs).push_back(K1 * R * K1.inverse());
+ }
+}
+
+} // namespace kernel
+} // namespace panography
+} // namespace libmv
diff --git a/extern/libmv/libmv/multiview/panography_kernel.h b/extern/libmv/libmv/multiview/panography_kernel.h
new file mode 100644
index 00000000000..a6adbd54b20
--- /dev/null
+++ b/extern/libmv/libmv/multiview/panography_kernel.h
@@ -0,0 +1,54 @@
+// 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_MULTIVIEW_PANOGRAPHY_KERNEL_H
+#define LIBMV_MULTIVIEW_PANOGRAPHY_KERNEL_H
+
+#include "libmv/base/vector.h"
+#include "libmv/multiview/conditioning.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/multiview/two_view_kernel.h"
+#include "libmv/multiview/homography_error.h"
+#include "libmv/numeric/numeric.h"
+
+namespace libmv {
+namespace panography {
+namespace kernel {
+
+struct TwoPointSolver {
+ enum { MINIMUM_SAMPLES = 2 };
+ static void Solve(const Mat &x1, const Mat &x2, vector<Mat3> *Hs);
+};
+
+typedef two_view::kernel::Kernel<
+ TwoPointSolver, homography::homography2D::AsymmetricError, Mat3>
+ UnnormalizedKernel;
+
+typedef two_view::kernel::Kernel<
+ two_view::kernel::NormalizedSolver<TwoPointSolver, UnnormalizerI>,
+ homography::homography2D::AsymmetricError,
+ Mat3>
+ Kernel;
+
+} // namespace kernel
+} // namespace panography
+} // namespace libmv
+
+#endif // LIBMV_MULTIVIEW_PANOGRAPHY_KERNEL_H
diff --git a/extern/libmv/libmv/multiview/panography_test.cc b/extern/libmv/libmv/multiview/panography_test.cc
new file mode 100644
index 00000000000..f6faf0f6022
--- /dev/null
+++ b/extern/libmv/libmv/multiview/panography_test.cc
@@ -0,0 +1,144 @@
+// 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.
+
+#include "libmv/logging/logging.h"
+#include "libmv/multiview/panography.h"
+#include "libmv/multiview/panography_kernel.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+
+namespace libmv {
+namespace {
+
+TEST(Panography, PrintSomeSharedFocalEstimationValues) {
+ Mat x1(2, 2), x2(2, 2);
+ x1<< 158, 78,
+ 124, 113;
+ x2<< 300, 214,
+ 125, 114;
+
+ // Normalize data (set principal point 0,0 and image border to 1.0).
+ x1.block<1, 2>(0, 0) /= 320;
+ x1.block<1, 2>(1, 0) /= 240;
+ x2.block<1, 2>(0, 0) /= 320;
+ x2.block<1, 2>(1, 0) /= 240;
+ x1+=Mat2::Constant(0.5);
+ x2+=Mat2::Constant(0.5);
+
+ vector<double> fs;
+ F_FromCorrespondance_2points(x1, x2, &fs);
+
+ // Assert we found a valid solution.
+ EXPECT_EQ(1, fs.size());
+ EXPECT_NEAR(1.01667, fs[1], 1e-3);
+}
+
+TEST(Panography, GetR_FixedCameraCenterWithIdentity) {
+ Mat x1(3, 3);
+ x1 << 0.5, 0.6, 0.7,
+ 0.5, 0.5, 0.4,
+ 10.0, 10.0, 10.0;
+
+ Mat3 R;
+ GetR_FixedCameraCenter(x1, x1, 1.0, &R);
+ R /= R(2, 2);
+ EXPECT_MATRIX_NEAR(Mat3::Identity(), R, 1e-8);
+ LOG(INFO) << "R \n" << R;
+}
+
+TEST(Panography, Homography_GetR_Test_PitchY30) {
+ int n = 3;
+
+ Mat x1(3, n);
+ x1 << 0.5, 0.6, 0.7,
+ 0.5, 0.5, 0.4,
+ 10, 10, 10;
+
+ Mat x2 = x1;
+ const double alpha = 30.0 * M_PI / 180.0;
+ Mat3 rotY;
+ rotY << cos(alpha), 0, -sin(alpha),
+ 0, 1, 0,
+ sin(alpha), 0, cos(alpha);
+
+ for (int i = 0; i < n; ++i) {
+ x2.block<3, 1>(0, i) = rotY * x1.col(i);
+ }
+
+ Mat3 R;
+ GetR_FixedCameraCenter(x1, x2, 1.0, &R);
+
+ // Assert that residuals are small enough
+ for (int i = 0; i < n; ++i) {
+ Vec residuals = (R * x1.col(i)) - x2.col(i);
+ EXPECT_NEAR(0, residuals.norm(), 1e-6);
+ }
+
+ // Check that the rotation angle along Y is the expected one.
+ // Use the euler approximation to recover the angle.
+ double pitch_y = asin(R(2, 0)) * 180.0 / M_PI;
+ EXPECT_NEAR(30, pitch_y, 1e-4);
+}
+
+TEST(MinimalPanoramic, Real_Case_Kernel) {
+ const int n = 2;
+ Mat x1(2, n); // From image 0.jpg
+ x1<< 158, 78,
+ 124, 113;
+
+ Mat x2(2, n); // From image 3.jpg
+ x2<< 300, 214,
+ 125, 114;
+
+ Mat3 Ground_TruthHomography;
+ Ground_TruthHomography<< 1, 0.02, 129.83,
+ -0.02, 1.012, 0.07823,
+ 0, 0, 1;
+
+ vector<Mat3> Hs;
+
+ libmv::panography::kernel::TwoPointSolver::Solve(x1, x2, &Hs);
+
+ LOG(INFO) << "Got " << Hs.size() << " solutions.";
+ for (int j = 0; j < Hs.size(); ++j) {
+ Mat3 H = Hs[j];
+
+ EXPECT_MATRIX_NEAR(H, Ground_TruthHomography, 1e-1);
+
+ Mat x1h, x2h;
+ EuclideanToHomogeneous(x1, &x1h);
+ EuclideanToHomogeneous(x2, &x2h);
+
+ // Assert that residuals are small enough
+ for (int i = 0; i < n; ++i) {
+ Vec x1p = H * x1h.col(i);
+ Vec residuals = x1p/x1p(2) - x2h.col(i);
+ EXPECT_MATRIX_NEAR_ZERO(residuals, 1e-5);
+ }
+ }
+}
+
+} // namespace
+} // namespace libmv
+
+// TODO(pmoulon): Add a real test case based on images.
+// TODO(pmoulon): Add a check for the actual f value for the real images.
+// TODO(pmoulon): Add a test that has some inliers and outliers.
diff --git a/extern/libmv/libmv/multiview/projection_test.cc b/extern/libmv/libmv/multiview/projection_test.cc
new file mode 100644
index 00000000000..c060bfb0681
--- /dev/null
+++ b/extern/libmv/libmv/multiview/projection_test.cc
@@ -0,0 +1,115 @@
+// Copyright (c) 2007, 2008 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 <iostream>
+
+#include "libmv/multiview/projection.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+
+namespace {
+using namespace libmv;
+
+TEST(Projection, P_From_KRt) {
+ Mat3 K, Kp;
+ K << 10, 1, 30,
+ 0, 20, 40,
+ 0, 0, 1;
+
+ Mat3 R, Rp;
+ R << 1, 0, 0,
+ 0, 1, 0,
+ 0, 0, 1;
+
+ Vec3 t, tp;
+ t << 1, 2, 3;
+
+ Mat34 P;
+ P_From_KRt(K, R, t, &P);
+ KRt_From_P(P, &Kp, &Rp, &tp);
+
+ EXPECT_MATRIX_NEAR(K, Kp, 1e-8);
+ EXPECT_MATRIX_NEAR(R, Rp, 1e-8);
+ EXPECT_MATRIX_NEAR(t, tp, 1e-8);
+
+ // TODO(keir): Change the code to ensure det(R) == 1, which is not currently
+ // the case. Also add a test for that here.
+}
+
+Vec4 GetRandomPoint() {
+ Vec4 X;
+ X.setRandom();
+ X(3) = 1;
+ return X;
+}
+
+TEST(Projection, isInFrontOfCamera) {
+ Mat34 P;
+ P << 1, 0, 0, 0,
+ 0, 1, 0, 0,
+ 0, 0, 1, 0;
+
+ Vec4 X_front = GetRandomPoint();
+ Vec4 X_back = GetRandomPoint();
+ X_front(2) = 10; // Any point in the positive Z direction
+ // where Z > 1 is infront of the camera.
+ X_back(2) = -10; // Any point int he negative Z dirstaion
+ // is behind the camera.
+
+ bool res_front = isInFrontOfCamera(P, X_front);
+ bool res_back = isInFrontOfCamera(P, X_back);
+
+ EXPECT_EQ(true, res_front);
+ EXPECT_EQ(false, res_back);
+}
+
+TEST(AutoCalibration, ProjectionShiftPrincipalPoint) {
+ Mat34 P1, P2;
+ P1 << 1, 0, 0, 0,
+ 0, 1, 0, 0,
+ 0, 0, 1, 0;
+ P2 << 1, 0, 3, 0,
+ 0, 1, 4, 0,
+ 0, 0, 1, 0;
+ Mat34 P1_computed, P2_computed;
+ ProjectionShiftPrincipalPoint(P1, Vec2(0, 0), Vec2(3, 4), &P2_computed);
+ ProjectionShiftPrincipalPoint(P2, Vec2(3, 4), Vec2(0, 0), &P1_computed);
+
+ EXPECT_MATRIX_EQ(P1, P1_computed);
+ EXPECT_MATRIX_EQ(P2, P2_computed);
+}
+
+TEST(AutoCalibration, ProjectionChangeAspectRatio) {
+ Mat34 P1, P2;
+ P1 << 1, 0, 3, 0,
+ 0, 1, 4, 0,
+ 0, 0, 1, 0;
+ P2 << 1, 0, 3, 0,
+ 0, 2, 4, 0,
+ 0, 0, 1, 0;
+ Mat34 P1_computed, P2_computed;
+ ProjectionChangeAspectRatio(P1, Vec2(3, 4), 1, 2, &P2_computed);
+ ProjectionChangeAspectRatio(P2, Vec2(3, 4), 2, 1, &P1_computed);
+
+ EXPECT_MATRIX_EQ(P1, P1_computed);
+ EXPECT_MATRIX_EQ(P2, P2_computed);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/multiview/resection_test.cc b/extern/libmv/libmv/multiview/resection_test.cc
new file mode 100644
index 00000000000..368e2281cfa
--- /dev/null
+++ b/extern/libmv/libmv/multiview/resection_test.cc
@@ -0,0 +1,61 @@
+// 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.
+
+#include <iostream>
+
+#include "libmv/multiview/projection.h"
+#include "libmv/multiview/resection.h"
+#include "libmv/multiview/test_data_sets.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+#include "libmv/logging/logging.h"
+
+namespace {
+
+using namespace libmv;
+using namespace libmv::resection;
+
+TEST(Resection, ThreeViews) {
+ int nviews = 5;
+ int npoints = 6;
+ NViewDataSet d = NRealisticCamerasFull(nviews, npoints);
+ for (int i = 0; i < nviews; ++i) {
+ Mat4X X(4, npoints);
+ X.block(0, 0, 3, npoints) = d.X;
+ X.row(3).setOnes();
+ const Mat2X &x = d.x[i];
+ Mat34 P;
+ Resection(x, X, &P);
+ Mat34 P_expected = d.P(i);
+
+ // Because the P matrices are homogeneous, it is necessary to be tricky
+ // about the scale factor to make them match.
+ P_expected *= 1/P_expected.array().abs().sum();
+ P *= 1/P.array().abs().sum();
+ if (!((P(0, 0) > 0 && P_expected(0, 0) > 0) ||
+ (P(0, 0) < 0 && P_expected(0, 0) < 0))) {
+ P *= -1;
+ }
+
+ EXPECT_MATRIX_NEAR(P_expected, P, 1e-9);
+ }
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/multiview/test_data_sets.cc b/extern/libmv/libmv/multiview/test_data_sets.cc
new file mode 100644
index 00000000000..110bde6f762
--- /dev/null
+++ b/extern/libmv/libmv/multiview/test_data_sets.cc
@@ -0,0 +1,196 @@
+// Copyright (c) 2007, 2008 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/multiview/test_data_sets.h"
+
+#include <cmath>
+
+#include "libmv/numeric/numeric.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/multiview/fundamental.h"
+
+namespace libmv {
+
+TwoViewDataSet TwoRealisticCameras(bool same_K) {
+ TwoViewDataSet d;
+
+ d.K1 << 320, 0, 160,
+ 0, 320, 120,
+ 0, 0, 1;
+ if (same_K) {
+ d.K2 = d.K1;
+ } else {
+ d.K2 << 360, 0, 170,
+ 0, 360, 110,
+ 0, 0, 1;
+ }
+ d.R1 = RotationAroundZ(-0.1);
+ d.R2 = RotationAroundX(-0.1);
+ d.t1 << 1, 1, 10;
+ d.t2 << -2, -1, 10;
+ P_From_KRt(d.K1, d.R1, d.t1, &d.P1);
+ P_From_KRt(d.K2, d.R2, d.t2, &d.P2);
+
+ FundamentalFromProjections(d.P1, d.P2, &d.F);
+
+ d.X.resize(3, 30);
+ d.X.setRandom();
+
+ Project(d.P1, d.X, &d.x1);
+ Project(d.P2, d.X, &d.x2);
+
+ return d;
+}
+
+nViewDatasetConfigator::nViewDatasetConfigator(int fx , int fy,
+ int cx, int cy,
+ double distance,
+ double jitter_amount) {
+ _fx = fx;
+ _fy = fy;
+ _cx = cx;
+ _cy = cy;
+ _dist = distance;
+ _jitter_amount = jitter_amount;
+}
+
+NViewDataSet NRealisticCamerasFull(int nviews, int npoints,
+ const nViewDatasetConfigator config) {
+ NViewDataSet d;
+ d.n = nviews;
+ d.K.resize(nviews);
+ d.R.resize(nviews);
+ d.t.resize(nviews);
+ d.C.resize(nviews);
+ d.x.resize(nviews);
+ d.x_ids.resize(nviews);
+
+ d.X.resize(3, npoints);
+ d.X.setRandom();
+ d.X *= 0.6;
+
+ Vecu all_point_ids(npoints);
+ for (size_t j = 0; j < npoints; ++j)
+ all_point_ids[j] = j;
+
+ for (size_t i = 0; i < nviews; ++i) {
+ Vec3 camera_center, t, jitter, lookdir;
+
+ double theta = i * 2 * M_PI / nviews;
+ camera_center << sin(theta), 0.0, cos(theta);
+ camera_center *= config._dist;
+ d.C[i] = camera_center;
+
+ jitter.setRandom();
+ jitter *= config._jitter_amount / camera_center.norm();
+ lookdir = -camera_center + jitter;
+
+ d.K[i] << config._fx, 0, config._cx,
+ 0, config._fy, config._cy,
+ 0, 0, 1;
+ d.R[i] = LookAt(lookdir);
+ d.t[i] = -d.R[i] * camera_center;
+ d.x[i] = Project(d.P(i), d.X);
+ d.x_ids[i] = all_point_ids;
+ }
+ return d;
+}
+
+
+NViewDataSet NRealisticCamerasSparse(int nviews, int npoints,
+ float view_ratio, unsigned min_projections,
+ const nViewDatasetConfigator config) {
+ assert(view_ratio <= 1.0);
+ assert(view_ratio > 0.0);
+ assert(min_projections <= npoints);
+ NViewDataSet d;
+ d.n = nviews;
+ d.K.resize(nviews);
+ d.R.resize(nviews);
+ d.t.resize(nviews);
+ d.C.resize(nviews);
+ d.x.resize(nviews);
+ d.x_ids.resize(nviews);
+
+ d.X.resize(3, npoints);
+ d.X.setRandom();
+ d.X *= 0.6;
+
+ Mat visibility(nviews, npoints);
+ visibility.setZero();
+ Mat randoms(nviews, npoints);
+ randoms.setRandom();
+ randoms = (randoms.array() + 1)/2.0;
+ unsigned num_visibles = 0;
+ for (size_t i = 0; i < nviews; ++i) {
+ num_visibles = 0;
+ for (size_t j = 0; j < npoints; j++) {
+ if (randoms(i, j) <= view_ratio) {
+ visibility(i, j) = true;
+ num_visibles++;
+ }
+ }
+ if (num_visibles < min_projections) {
+ unsigned num_projections_to_add = min_projections - num_visibles;
+ for (size_t j = 0; j < npoints && num_projections_to_add > 0; ++j) {
+ if (!visibility(i, j)) {
+ num_projections_to_add--;
+ }
+ }
+ num_visibles += num_projections_to_add;
+ }
+ d.x_ids[i].resize(num_visibles);
+ d.x[i].resize(2, num_visibles);
+ }
+
+ size_t j_visible = 0;
+ Vec3 X;
+ for (size_t i = 0; i < nviews; ++i) {
+ Vec3 camera_center, t, jitter, lookdir;
+
+ double theta = i * 2 * M_PI / nviews;
+ camera_center << sin(theta), 0.0, cos(theta);
+ camera_center *= config._dist;
+ d.C[i] = camera_center;
+
+ jitter.setRandom();
+ jitter *= config._jitter_amount / camera_center.norm();
+ lookdir = -camera_center + jitter;
+
+ d.K[i] << config._fx, 0, config._cx,
+ 0, config._fy, config._cy,
+ 0, 0, 1;
+ d.R[i] = LookAt(lookdir);
+ d.t[i] = -d.R[i] * camera_center;
+ j_visible = 0;
+ for (size_t j = 0; j < npoints; j++) {
+ if (visibility(i, j)) {
+ X = d.X.col(j);
+ d.x[i].col(j_visible) = Project(d.P(i), X);
+ d.x_ids[i][j_visible] = j;
+ j_visible++;
+ }
+ }
+ }
+ return d;
+}
+
+
+} // namespace libmv
diff --git a/extern/libmv/libmv/multiview/test_data_sets.h b/extern/libmv/libmv/multiview/test_data_sets.h
new file mode 100644
index 00000000000..cf01663ca02
--- /dev/null
+++ b/extern/libmv/libmv/multiview/test_data_sets.h
@@ -0,0 +1,105 @@
+// Copyright (c) 2007, 2008 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_MULTIVIEW_TEST_DATA_SETS_H_
+#define LIBMV_MULTIVIEW_TEST_DATA_SETS_H_
+
+#include "libmv/base/vector.h"
+#include "libmv/multiview/fundamental.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/numeric/numeric.h"
+
+namespace libmv {
+
+struct TwoViewDataSet {
+ Mat3 K1, K2; // Internal parameters.
+ Mat3 R1, R2; // Rotation.
+ Vec3 t1, t2; // Translation.
+ Mat34 P1, P2; // Projection matrix, P = K(R|t)
+ Mat3 F; // Fundamental matrix.
+ Mat3X X; // 3D points.
+ Mat2X x1, x2; // Projected points.
+};
+
+// Two cameras at (-1,-1,-10) and (2,1,-10) looking approximately towards z+.
+TwoViewDataSet TwoRealisticCameras(bool same_K = false);
+
+// An N-view metric dataset . An important difference between this
+// and the other reconstruction data types is that all points are seen by all
+// cameras.
+struct NViewDataSet {
+ vector<Mat3> K; // Internal parameters (fx, fy, etc).
+ vector<Mat3> R; // Rotation.
+ vector<Vec3> t; // Translation.
+ vector<Vec3> C; // Camera centers.
+ Mat3X X; // 3D points.
+ vector<Mat2X> x; // Projected points; may have noise added.
+ vector<Vecu> x_ids; // Indexes of points corresponding to the projections
+
+ int n; // Actual number of cameras.
+
+ Mat34 P(int i) {
+ assert(i < n);
+ return K[i] * HStack(R[i], t[i]);
+ }
+ Mat3 F(int i, int j) {
+ Mat3 F_;
+ FundamentalFromProjections(P(i), P(j), &F_);
+ return F_;
+ }
+ void Reproject() {
+ for (int i = 0; i < n; ++i) {
+ x[i] = Project(P(i), X);
+ }
+ }
+ // TODO(keir): Add gaussian jitter functions.
+};
+
+struct nViewDatasetConfigator {
+ /// Internal camera parameters
+ int _fx;
+ int _fy;
+ int _cx;
+ int _cy;
+
+ /// Camera random position parameters
+ double _dist;
+ double _jitter_amount;
+
+ nViewDatasetConfigator(int fx = 1000, int fy = 1000,
+ int cx = 500, int cy = 500,
+ double distance = 1.5,
+ double jitter_amount = 0.01);
+};
+
+NViewDataSet NRealisticCamerasFull(int nviews, int npoints,
+ const nViewDatasetConfigator
+ config = nViewDatasetConfigator());
+
+// Generates sparse projections (not all points are projected)
+NViewDataSet NRealisticCamerasSparse(int nviews, int npoints,
+ float view_ratio = 0.6,
+ unsigned min_projections = 3,
+ const nViewDatasetConfigator
+ config = nViewDatasetConfigator());
+
+} // namespace libmv
+
+#endif // LIBMV_MULTIVIEW_TEST_DATA_SETS_H_
diff --git a/extern/libmv/libmv/multiview/triangulation_test.cc b/extern/libmv/libmv/multiview/triangulation_test.cc
new file mode 100644
index 00000000000..66d1ee25a62
--- /dev/null
+++ b/extern/libmv/libmv/multiview/triangulation_test.cc
@@ -0,0 +1,47 @@
+// Copyright (c) 2007, 2008 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 <iostream>
+
+#include "libmv/multiview/triangulation.h"
+#include "libmv/multiview/fundamental.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/multiview/test_data_sets.h"
+#include "libmv/numeric/numeric.h"
+#include "testing/testing.h"
+
+namespace {
+using namespace libmv;
+
+TEST(Triangulation, TriangulateDLT) {
+ TwoViewDataSet d = TwoRealisticCameras();
+
+ for (int i = 0; i < d.X.cols(); ++i) {
+ Vec2 x1, x2;
+ MatrixColumn(d.x1, i, &x1);
+ MatrixColumn(d.x2, i, &x2);
+ Vec3 X_estimated, X_gt;
+ MatrixColumn(d.X, i, &X_gt);
+ TriangulateDLT(d.P1, x1, d.P2, x2, &X_estimated);
+ EXPECT_NEAR(0, DistanceLInfinity(X_estimated, X_gt), 1e-8);
+ }
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/multiview/two_view_kernel.h b/extern/libmv/libmv/multiview/two_view_kernel.h
new file mode 100644
index 00000000000..7af0ed5ddab
--- /dev/null
+++ b/extern/libmv/libmv/multiview/two_view_kernel.h
@@ -0,0 +1,137 @@
+// 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_MULTIVIEW_TWO_VIEW_KERNEL_H_
+#define LIBMV_MULTIVIEW_TWO_VIEW_KERNEL_H_
+
+#include "libmv/base/vector.h"
+#include "libmv/logging/logging.h"
+#include "libmv/multiview/conditioning.h"
+#include "libmv/numeric/numeric.h"
+
+namespace libmv {
+namespace two_view {
+namespace kernel {
+
+template<typename Solver, typename Unnormalizer>
+struct NormalizedSolver {
+ enum { MINIMUM_SAMPLES = Solver::MINIMUM_SAMPLES };
+ static void Solve(const Mat &x1, const Mat &x2, vector<Mat3> *models) {
+ assert(2 == x1.rows());
+ assert(MINIMUM_SAMPLES <= x1.cols());
+ assert(x1.rows() == x2.rows());
+ assert(x1.cols() == x2.cols());
+
+ // Normalize the data.
+ Mat3 T1, T2;
+ Mat x1_normalized, x2_normalized;
+ NormalizePoints(x1, &x1_normalized, &T1);
+ NormalizePoints(x2, &x2_normalized, &T2);
+
+ Solver::Solve(x1_normalized, x2_normalized, models);
+
+ for (int i = 0; i < models->size(); ++i) {
+ Unnormalizer::Unnormalize(T1, T2, &(*models)[i]);
+ }
+ }
+};
+
+template<typename Solver, typename Unnormalizer>
+struct IsotropicNormalizedSolver {
+ enum { MINIMUM_SAMPLES = Solver::MINIMUM_SAMPLES };
+ static void Solve(const Mat &x1, const Mat &x2, vector<Mat3> *models) {
+ assert(2 == x1.rows());
+ assert(MINIMUM_SAMPLES <= x1.cols());
+ assert(x1.rows() == x2.rows());
+ assert(x1.cols() == x2.cols());
+
+ // Normalize the data.
+ Mat3 T1, T2;
+ Mat x1_normalized, x2_normalized;
+ NormalizeIsotropicPoints(x1, &x1_normalized, &T1);
+ NormalizeIsotropicPoints(x2, &x2_normalized, &T2);
+
+ Solver::Solve(x1_normalized, x2_normalized, models);
+
+ for (int i = 0; i < models->size(); ++i) {
+ Unnormalizer::Unnormalize(T1, T2, &(*models)[i]);
+ }
+ }
+};
+// This is one example (targeted at solvers that operate on correspondences
+// between two views) that shows the "kernel" part of a robust fitting
+// problem:
+//
+// 1. The model; Mat3 in the case of the F or H matrix.
+// 2. The minimum number of samples needed to fit; 7 or 8 (or 4).
+// 3. A way to convert samples to a model.
+// 4. A way to convert a sample and a model to an error.
+//
+// Of particular note is that the kernel does not expose what the samples are.
+// All the robust fitting algorithm sees is that there is some number of
+// samples; it is able to fit subsets of them (via the kernel) and check their
+// error, but can never access the samples themselves.
+//
+// The Kernel objects must follow the following concept so that the robust
+// fitting alogrithm can fit this type of relation:
+//
+// 1. Kernel::Model
+// 2. Kernel::MINIMUM_SAMPLES
+// 3. Kernel::Fit(vector<int>, vector<Kernel::Model> *)
+// 4. Kernel::Error(int, Model) -> error
+//
+// The fit routine must not clear existing entries in the vector of models; it
+// should append new solutions to the end.
+template<typename SolverArg,
+ typename ErrorArg,
+ typename ModelArg = Mat3>
+class Kernel {
+ public:
+ Kernel(const Mat &x1, const Mat &x2) : x1_(x1), x2_(x2) {}
+ typedef SolverArg Solver;
+ typedef ModelArg Model;
+ enum { MINIMUM_SAMPLES = Solver::MINIMUM_SAMPLES };
+ void Fit(const vector<int> &samples, vector<Model> *models) const {
+ Mat x1 = ExtractColumns(x1_, samples);
+ Mat x2 = ExtractColumns(x2_, samples);
+ Solver::Solve(x1, x2, models);
+ }
+ double Error(int sample, const Model &model) const {
+ return ErrorArg::Error(model,
+ static_cast<Vec>(x1_.col(sample)),
+ static_cast<Vec>(x2_.col(sample)));
+ }
+ int NumSamples() const {
+ return x1_.cols();
+ }
+ static void Solve(const Mat &x1, const Mat &x2, vector<Model> *models) {
+ // By offering this, Kernel types can be passed to templates.
+ Solver::Solve(x1, x2, models);
+ }
+ protected:
+ const Mat &x1_;
+ const Mat &x2_;
+};
+
+} // namespace kernel
+} // namespace two_view
+} // namespace libmv
+
+#endif // LIBMV_MULTIVIEW_TWO_VIEW_KERNEL_H_
diff --git a/extern/libmv/libmv/numeric/dogleg_test.cc b/extern/libmv/libmv/numeric/dogleg_test.cc
new file mode 100644
index 00000000000..90c46c31672
--- /dev/null
+++ b/extern/libmv/libmv/numeric/dogleg_test.cc
@@ -0,0 +1,95 @@
+// 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.
+
+#include "testing/testing.h"
+#include "libmv/numeric/dogleg.h"
+
+using namespace libmv;
+
+namespace {
+
+class F {
+ public:
+ typedef Vec4 FMatrixType;
+ typedef Vec3 XMatrixType;
+ Vec4 operator()(const Vec3 &x) const {
+ double x1 = x.x() - 2;
+ double y1 = x.y() - 5;
+ double z1 = x.z();
+ Vec4 fx; fx << x1*x1 + z1*z1,
+ y1*y1 + z1*z1,
+ z1*z1,
+ x1*x1;
+ return fx;
+ }
+};
+
+TEST(Dogleg, SimpleCase) {
+ Vec3 x; x << 0.76026643, -30.01799744, 0.55192142;
+ F f;
+ Dogleg<F>::SolverParameters params;
+ Dogleg<F> lm(f);
+ /* TODO(sergey): Better error handling. */
+ /* Dogleg<F>::Results results = */ lm.minimize(params, &x);
+ Vec3 expected_min_x; expected_min_x << 2, 5, 0;
+
+ EXPECT_MATRIX_NEAR(expected_min_x, x, 1e-5);
+}
+
+// Example 3.2 from [1]; page 11 of the pdf, 20 of the document. This is a
+// tricky problem because of the singluar Jacobian near the origin.
+class F32 {
+ public:
+ typedef Vec2 FMatrixType;
+ typedef Vec2 XMatrixType;
+ Vec2 operator()(const Vec2 &x) const {
+ double x1 = x(0);
+ double x2 = 10*x(0)/(x(0) + 0.1) + 2*x(1)*x(1);
+ Vec2 fx; fx << x1, x2;
+ return fx;
+ }
+};
+
+class JF32 {
+ public:
+ JF32(const F32 &f) { (void) f; }
+ Mat2 operator()(const Vec2 &x) {
+ Mat2 J; J << 1, 0,
+ 1./pow(x(0) + 0.1, 2), 4*x(1)*x(1);
+ return J;
+ }
+};
+
+// TODO(keir): Re-enable this when the dogleg code properly handles singular
+// normal equations.
+/*
+TEST(Dogleg, Example32) {
+ Vec2 x; x << 3, 1;
+ F32 f;
+ CheckJacobian<F32, JF32>(f, x);
+ Dogleg<F32, JF32> dogleg(f);
+ Dogleg<F32, JF32>::Results results = dogleg.minimize(&x);
+ Vec2 expected_min_x; expected_min_x << 0, 0;
+
+ EXPECT_MATRIX_NEAR(expected_min_x, x, 1e-5);
+}
+*/
+
+} // namespace
diff --git a/extern/libmv/libmv/numeric/function_derivative_test.cc b/extern/libmv/libmv/numeric/function_derivative_test.cc
new file mode 100644
index 00000000000..8d976d3e9a0
--- /dev/null
+++ b/extern/libmv/libmv/numeric/function_derivative_test.cc
@@ -0,0 +1,57 @@
+// Copyright (c) 2007, 2008, 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.
+
+#include "testing/testing.h"
+#include "libmv/numeric/numeric.h"
+#include "libmv/numeric/function_derivative.h"
+
+using namespace libmv;
+
+namespace {
+
+class F {
+ public:
+ typedef Vec2 FMatrixType;
+ typedef Vec3 XMatrixType;
+ Vec2 operator()(const Vec3 &x) const {
+ Vec2 fx;
+ fx << 0.19*x(0) + 0.19*x(1)*x(1) + x(2),
+ 3*sin(x(0)) + 2*cos(x(1));
+ return fx;
+ }
+ Mat23 J(const Vec3 &x) const {
+ Mat23 jacobian;
+ jacobian << 0.19, 2*0.19*x(1), 1.0,
+ 3*cos(x(0)), -2*sin(x(1)), 0;
+ return jacobian;
+ }
+};
+
+TEST(FunctionDerivative, SimpleCase) {
+ Vec3 x; x << 0.76026643, 0.01799744, 0.55192142;
+ F f;
+ NumericJacobian<F, CENTRAL> J(f);
+ EXPECT_MATRIX_NEAR(f.J(x), J(x), 1e-8);
+ NumericJacobian<F, FORWARD> J_forward(f);
+ // Forward difference is very inaccurate.
+ EXPECT_MATRIX_NEAR(f.J(x), J_forward(x), 1e-5);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/numeric/levenberg_marquardt_test.cc b/extern/libmv/libmv/numeric/levenberg_marquardt_test.cc
new file mode 100644
index 00000000000..fc3f9ebbb29
--- /dev/null
+++ b/extern/libmv/libmv/numeric/levenberg_marquardt_test.cc
@@ -0,0 +1,56 @@
+// Copyright (c) 2007, 2008, 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.
+
+#include "testing/testing.h"
+#include "libmv/numeric/levenberg_marquardt.h"
+
+using namespace libmv;
+
+namespace {
+
+class F {
+ public:
+ typedef Vec4 FMatrixType;
+ typedef Vec3 XMatrixType;
+ Vec4 operator()(const Vec3 &x) const {
+ double x1 = x.x() - 2;
+ double y1 = x.y() - 5;
+ double z1 = x.z();
+ Vec4 fx; fx << x1*x1 + z1*z1,
+ y1*y1 + z1*z1,
+ z1*z1,
+ x1*x1;
+ return fx;
+ }
+};
+
+TEST(LevenbergMarquardt, SimpleCase) {
+ Vec3 x(0.76026643, -30.01799744, 0.55192142);
+ F f;
+ LevenbergMarquardt<F>::SolverParameters params;
+ LevenbergMarquardt<F> lm(f);
+ /* TODO(sergey): Better error handling. */
+ /* LevenbergMarquardt<F>::Results results = */ lm.minimize(params, &x);
+ Vec3 expected_min_x(2, 5, 0);
+
+ EXPECT_MATRIX_NEAR(expected_min_x, x, 1e-5);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/numeric/numeric_test.cc b/extern/libmv/libmv/numeric/numeric_test.cc
new file mode 100644
index 00000000000..0cdfaf33ab2
--- /dev/null
+++ b/extern/libmv/libmv/numeric/numeric_test.cc
@@ -0,0 +1,439 @@
+// Copyright (c) 2007, 2008 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/numeric/numeric.h"
+#include "testing/testing.h"
+
+using namespace libmv;
+
+namespace {
+
+TEST(Numeric, DynamicSizedNullspace) {
+ Mat A(3, 4);
+ A << 0.76026643, 0.01799744, 0.55192142, 0.8699745,
+ 0.42016166, 0.97863392, 0.33711682, 0.14479271,
+ 0.51016811, 0.66528302, 0.54395496, 0.57794893;
+ Vec x;
+ double s = Nullspace(&A, &x);
+ EXPECT_NEAR(0.0, s, 1e-15);
+ EXPECT_NEAR(0.0, (A * x).norm(), 1e-15);
+ EXPECT_NEAR(1.0, x.norm(), 1e-15);
+}
+
+TEST(Numeric, FixedSizeMatrixNullspace) {
+ Mat34 A;
+ A << 0.76026643, 0.01799744, 0.55192142, 0.8699745,
+ 0.42016166, 0.97863392, 0.33711682, 0.14479271,
+ 0.51016811, 0.66528302, 0.54395496, 0.57794893;
+ Vec x;
+ double s = Nullspace(&A, &x);
+ EXPECT_NEAR(0.0, s, 1e-15);
+ EXPECT_NEAR(0.0, (A * x).norm(), 1e-15);
+ EXPECT_NEAR(1.0, x.norm(), 1e-15);
+}
+
+TEST(Numeric, NullspaceMatchesLapackSVD) {
+ Mat43 A;
+ A << 0.76026643, 0.01799744, 0.55192142,
+ 0.8699745, 0.42016166, 0.97863392,
+ 0.33711682, 0.14479271, 0.51016811,
+ 0.66528302, 0.54395496, 0.57794893;
+ Vec x;
+ double s = Nullspace(&A, &x);
+ EXPECT_NEAR(1.0, x.norm(), 1e-15);
+ EXPECT_NEAR(0.206694992663, s, 1e-9);
+ EXPECT_NEAR(0.206694992663, (A * x).norm(), 1e-9);
+
+ EXPECT_NEAR(-0.64999717, x(0), 1e-8);
+ EXPECT_NEAR(-0.18452646, x(1), 1e-8);
+ EXPECT_NEAR(0.7371931, x(2), 1e-8);
+}
+
+TEST(Numeric, Nullspace2) {
+ Mat43 A;
+ A << 0.76026643, 0.01799744, 0.55192142,
+ 0.8699745, 0.42016166, 0.97863392,
+ 0.33711682, 0.14479271, 0.51016811,
+ 0.66528302, 0.54395496, 0.57794893;
+ Vec3 x1, x2;
+ double s = Nullspace2(&A, &x1, &x2);
+ EXPECT_NEAR(1.0, x1.norm(), 1e-15);
+ EXPECT_NEAR(0.206694992663, s, 1e-9);
+ EXPECT_NEAR(0.206694992663, (A * x1).norm(), 1e-9);
+
+ EXPECT_NEAR(-0.64999717, x1(0), 1e-8);
+ EXPECT_NEAR(-0.18452646, x1(1), 1e-8);
+ EXPECT_NEAR( 0.7371931, x1(2), 1e-8);
+
+ if (x2(0) < 0) {
+ x2 *= -1;
+ }
+ EXPECT_NEAR( 0.34679618, x2(0), 1e-8);
+ EXPECT_NEAR(-0.93519689, x2(1), 1e-8);
+ EXPECT_NEAR( 0.07168809, x2(2), 1e-8);
+}
+
+TEST(Numeric, TinyMatrixSquareTranspose) {
+ Mat2 A;
+ A << 1.0, 2.0, 3.0, 4.0;
+ libmv::TransposeInPlace(&A);
+ EXPECT_EQ(1.0, A(0, 0));
+ EXPECT_EQ(3.0, A(0, 1));
+ EXPECT_EQ(2.0, A(1, 0));
+ EXPECT_EQ(4.0, A(1, 1));
+}
+
+TEST(Numeric, NormalizeL1) {
+ Vec2 x;
+ x << 1, 2;
+ double l1 = NormalizeL1(&x);
+ EXPECT_DOUBLE_EQ(3., l1);
+ EXPECT_DOUBLE_EQ(1./3., x(0));
+ EXPECT_DOUBLE_EQ(2./3., x(1));
+}
+
+TEST(Numeric, NormalizeL2) {
+ Vec2 x;
+ x << 1, 2;
+ double l2 = NormalizeL2(&x);
+ EXPECT_DOUBLE_EQ(sqrt(5.0), l2);
+ EXPECT_DOUBLE_EQ(1./sqrt(5.), x(0));
+ EXPECT_DOUBLE_EQ(2./sqrt(5.), x(1));
+}
+
+TEST(Numeric, Diag) {
+ Vec x(2);
+ x << 1, 2;
+ Mat D = Diag(x);
+ EXPECT_EQ(1, D(0, 0));
+ EXPECT_EQ(0, D(0, 1));
+ EXPECT_EQ(0, D(1, 0));
+ EXPECT_EQ(2, D(1, 1));
+}
+
+TEST(Numeric, Determinant) {
+ Mat A(2, 2);
+ A << 1, 2,
+ -1, 3;
+ double detA = A.determinant();
+ EXPECT_NEAR(5, detA, 1e-8);
+
+ Mat B(4, 4);
+ B << 0, 1, 2, 3,
+ 4, 5, 6, 7,
+ 8, 9, 10, 11,
+ 12, 13, 14, 15;
+ double detB = B.determinant();
+ EXPECT_NEAR(0, detB, 1e-8);
+
+ Mat3 C;
+ C << 0, 1, 2,
+ 3, 4, 5,
+ 6, 7, 1;
+ double detC = C.determinant();
+ EXPECT_NEAR(21, detC, 1e-8);
+}
+
+TEST(Numeric, Inverse) {
+ Mat A(2, 2), A1;
+ A << 1, 2,
+ -1, 3;
+ Mat I = A * A.inverse();
+
+ EXPECT_NEAR(1, I(0, 0), 1e-8);
+ EXPECT_NEAR(0, I(0, 1), 1e-8);
+ EXPECT_NEAR(0, I(1, 0), 1e-8);
+ EXPECT_NEAR(1, I(1, 1), 1e-8);
+
+ Mat B(4, 4), B1;
+ B << 0, 1, 2, 3,
+ 4, 5, 6, 7,
+ 8, 9, 2, 11,
+ 12, 13, 14, 4;
+ Mat I2 = B * B.inverse();
+ EXPECT_NEAR(1, I2(0, 0), 1e-8);
+ EXPECT_NEAR(0, I2(0, 1), 1e-8);
+ EXPECT_NEAR(0, I2(0, 2), 1e-8);
+ EXPECT_NEAR(0, I2(1, 0), 1e-8);
+ EXPECT_NEAR(1, I2(1, 1), 1e-8);
+ EXPECT_NEAR(0, I2(1, 2), 1e-8);
+ EXPECT_NEAR(0, I2(2, 0), 1e-8);
+ EXPECT_NEAR(0, I2(2, 1), 1e-8);
+ EXPECT_NEAR(1, I2(2, 2), 1e-8);
+}
+
+TEST(Numeric, MeanAndVarianceAlongRows) {
+ int n = 4;
+ Mat points(2, n);
+ points << 0, 0, 1, 1,
+ 0, 2, 1, 3;
+
+ Vec mean, variance;
+ MeanAndVarianceAlongRows(points, &mean, &variance);
+
+ EXPECT_NEAR(0.5, mean(0), 1e-8);
+ EXPECT_NEAR(1.5, mean(1), 1e-8);
+ EXPECT_NEAR(0.25, variance(0), 1e-8);
+ EXPECT_NEAR(1.25, variance(1), 1e-8);
+}
+
+TEST(Numeric, HorizontalStack) {
+ Mat x(2, 1), y(2, 1), z;
+ x << 1, 2;
+ y << 3, 4;
+
+ HorizontalStack(x, y, &z);
+
+ EXPECT_EQ(2, z.cols());
+ EXPECT_EQ(2, z.rows());
+ EXPECT_EQ(1, z(0, 0));
+ EXPECT_EQ(2, z(1, 0));
+ EXPECT_EQ(3, z(0, 1));
+ EXPECT_EQ(4, z(1, 1));
+}
+
+TEST(Numeric, HStack) {
+ Mat x(2, 1), y(2, 1), z(2, 2);
+ x << 1, 2;
+ y << 3, 4;
+ z << 1, 3,
+ 2, 4;
+ Vec2 xC = x, yC = y;
+
+ Mat2 xy = HStack(x, y);
+ EXPECT_MATRIX_EQ(z, xy);
+
+ EXPECT_MATRIX_EQ(z, HStack(x, y));
+ EXPECT_MATRIX_EQ(z, HStack(x, yC));
+ EXPECT_MATRIX_EQ(z, HStack(xC, y));
+ EXPECT_MATRIX_EQ(z, HStack(xC, yC));
+}
+
+// TODO(keir): Need some way of verifying that the compile time types of the
+// resulting stacked matrices properly propagate the fixed dimensions.
+TEST(Numeric, VStack) {
+ Mat x(2, 2), y(2, 2), z(4, 2);
+ x << 1, 2,
+ 3, 4;
+ y << 10, 20,
+ 30, 40;
+ z << 1, 2,
+ 3, 4,
+ 10, 20,
+ 30, 40;
+ Mat2 xC = x, yC = y;
+
+ Mat xy = VStack(x, y);
+ EXPECT_MATRIX_EQ(z, xy);
+
+ EXPECT_MATRIX_EQ(z, VStack(x, y));
+ EXPECT_MATRIX_EQ(z, VStack(x, yC));
+ EXPECT_MATRIX_EQ(z, VStack(xC, y));
+ EXPECT_MATRIX_EQ(z, VStack(xC, yC));
+}
+
+TEST(Numeric, VerticalStack) {
+ Mat x(1, 2), y(1, 2), z;
+ x << 1, 2;
+ y << 3, 4;
+ VerticalStack(x, y, &z);
+
+ EXPECT_EQ(2, z.cols());
+ EXPECT_EQ(2, z.rows());
+ EXPECT_EQ(1, z(0, 0));
+ EXPECT_EQ(2, z(0, 1));
+ EXPECT_EQ(3, z(1, 0));
+ EXPECT_EQ(4, z(1, 1));
+}
+
+TEST(Numeric, CrossProduct) {
+ Vec3 x, y, z;
+ x << 1, 0, 0;
+ y << 0, 1, 0;
+ z << 0, 0, 1;
+ Vec3 xy = CrossProduct(x, y);
+ Vec3 yz = CrossProduct(y, z);
+ Vec3 zx = CrossProduct(z, x);
+ EXPECT_NEAR(0, DistanceLInfinity(xy, z), 1e-8);
+ EXPECT_NEAR(0, DistanceLInfinity(yz, x), 1e-8);
+ EXPECT_NEAR(0, DistanceLInfinity(zx, y), 1e-8);
+}
+
+TEST(Numeric, CrossProductMatrix) {
+ Vec3 x, y;
+ x << 1, 2, 3;
+ y << 2, 3, 4;
+ Vec3 xy = CrossProduct(x, y);
+ Vec3 yx = CrossProduct(y, x);
+ Mat3 X = CrossProductMatrix(x);
+ Vec3 Xy, Xty;
+ Xy = X * y;
+ Xty = X.transpose() * y;
+ EXPECT_NEAR(0, DistanceLInfinity(xy, Xy), 1e-8);
+ EXPECT_NEAR(0, DistanceLInfinity(yx, Xty), 1e-8);
+}
+
+TEST(Numeric, MatrixColumn) {
+ Mat A2(2, 3);
+ Vec2 v2;
+ A2 << 1, 2, 3,
+ 4, 5, 6;
+ MatrixColumn(A2, 1, &v2);
+ EXPECT_EQ(2, v2(0));
+ EXPECT_EQ(5, v2(1));
+
+ Mat A3(3, 3);
+ Vec3 v3;
+ A3 << 1, 2, 3,
+ 4, 5, 6,
+ 7, 8, 9;
+ MatrixColumn(A3, 1, &v3);
+ EXPECT_EQ(2, v3(0));
+ EXPECT_EQ(5, v3(1));
+ EXPECT_EQ(8, v3(2));
+
+ Mat A4(4, 3);
+ Vec4 v4;
+ A4 << 1, 2, 3,
+ 4, 5, 6,
+ 7, 8, 9,
+ 10, 11, 12;
+ MatrixColumn(A4, 1, &v4);
+ EXPECT_EQ( 2, v4(0));
+ EXPECT_EQ( 5, v4(1));
+ EXPECT_EQ( 8, v4(2));
+ EXPECT_EQ(11, v4(3));
+}
+
+// This used to give a compile error with FLENS.
+TEST(Numeric, TinyMatrixView) {
+ Mat34 P;
+ Mat K = P.block(0, 0, 3, 3);
+}
+
+// This gives a compile error.
+TEST(Numeric, Mat3MatProduct) {
+ Mat3 A;
+ Mat3 B;
+ Mat C = A * B;
+}
+
+// This gives a compile error.
+TEST(Numeric, Vec3Negative) {
+ Vec3 y; y << 1, 2, 3;
+ Vec3 x = -y;
+ EXPECT_EQ(-1, x(0));
+ EXPECT_EQ(-2, x(1));
+ EXPECT_EQ(-3, x(2));
+}
+
+// This gives a compile error.
+TEST(Numeric, Vec3VecInteroperability) {
+ Vec y(3);
+ y << 1, 2, 3;
+ Vec3 x = y + y;
+ EXPECT_EQ(2, x(0));
+ EXPECT_EQ(4, x(1));
+ EXPECT_EQ(6, x(2));
+}
+
+// This segfaults inside lapack.
+TEST(Numeric, DeterminantLU7) {
+ Mat A(5, 5);
+ A << 1, 0, 0, 0, 0,
+ 0, 1, 0, 0, 0,
+ 0, 0, 1, 0, 0,
+ 0, 0, 0, 1, 0,
+ 0, 0, 0, 0, 1;
+ EXPECT_NEAR(1, A.determinant(), 1e-8);
+}
+
+// This segfaults inside lapack.
+TEST(Numeric, DeterminantLU) {
+ Mat A(2, 2);
+ A << 1, 2,
+ -1, 3;
+ EXPECT_NEAR(5, A.determinant(), 1e-8);
+}
+
+// This does unexpected things.
+// Keir: Not with eigen2!
+TEST(Numeric, InplaceProduct) {
+ Mat2 K, S;
+ K << 1, 0,
+ 0, 1;
+ S << 1, 0,
+ 0, 1;
+ K = K * S;
+ EXPECT_MATRIX_NEAR(Mat2::Identity(), K, 1e-8);
+}
+
+TEST(Numeric, ExtractColumns) {
+ Mat2X A(2, 5);
+ A << 1, 2, 3, 4, 5,
+ 6, 7, 8, 9, 10;
+ Vec2i columns; columns << 0, 2;
+ Mat2X extracted = ExtractColumns(A, columns);
+ EXPECT_NEAR(1, extracted(0, 0), 1e-15);
+ EXPECT_NEAR(3, extracted(0, 1), 1e-15);
+ EXPECT_NEAR(6, extracted(1, 0), 1e-15);
+ EXPECT_NEAR(8, extracted(1, 1), 1e-15);
+}
+
+TEST(Numeric, RotationRodrigues) {
+ Vec3 x, y, z;
+ x << 1, 0, 0;
+ y << 0, 1, 0;
+ z << 0, 0, 1;
+
+ Mat3 rodrigues_x = RotationRodrigues(x);
+ Mat3 rodrigues_y = RotationRodrigues(y);
+ Mat3 rodrigues_z = RotationRodrigues(z);
+
+ Mat3 Rx = RotationAroundX(1);
+ Mat3 Ry = RotationAroundY(1);
+ Mat3 Rz = RotationAroundZ(1);
+
+ EXPECT_MATRIX_NEAR(Rx, rodrigues_x, 1e-15);
+ EXPECT_MATRIX_NEAR(Ry, rodrigues_y, 1e-15);
+ EXPECT_MATRIX_NEAR(Rz, rodrigues_z, 1e-15);
+}
+
+TEST(Numeric, LookAt) {
+ // Simple orthogonality check.
+ Vec3 e; e << 1, 2, 3;
+ Mat3 R = LookAt(e), I = Mat3::Identity();
+ Mat3 RRT = R*R.transpose();
+ Mat3 RTR = R.transpose()*R;
+
+ EXPECT_MATRIX_NEAR(I, RRT, 1e-15);
+ EXPECT_MATRIX_NEAR(I, RTR, 1e-15);
+}
+
+TEST(Numeric, Reshape) {
+ Vec4 x; x << 1, 2, 3, 4;
+ Mat2 M, M_expected;
+ reshape(x, 2, 2, &M);
+ M_expected << 1, 2,
+ 3, 4;
+ EXPECT_MATRIX_NEAR(M_expected, M, 1e-15);
+}
+
+} // namespace
diff --git a/extern/libmv/libmv/numeric/poly_test.cc b/extern/libmv/libmv/numeric/poly_test.cc
new file mode 100644
index 00000000000..ea50383190f
--- /dev/null
+++ b/extern/libmv/libmv/numeric/poly_test.cc
@@ -0,0 +1,98 @@
+// Copyright (c) 2007, 2008 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/numeric/numeric.h"
+#include "libmv/numeric/poly.h"
+#include "testing/testing.h"
+
+using namespace libmv;
+
+namespace {
+
+// Find the polynomial coefficients of x in the equation
+//
+// (x - a)(x - b)(x - c) == 0
+//
+// by expanding to
+//
+// x^3 - (c+b+a) * x^2 + (a*b+(b+a)*c) * x - a*b*c = 0.
+// = p = q = r
+void CoeffsForCubicZeros(double a, double b, double c,
+ double *p, double *q, double *r) {
+ *p = -(c + b + a);
+ *q = (a * b + (b + a) * c);
+ *r = -a * b * c;
+}
+// Find the polynomial coefficients of x in the equation
+//
+// (x - a)(x - b)(x - c)(x - d) == 0
+//
+// by expanding to
+//
+// x^4 - (d+c+b+a) * x^3 + (d*(c+b+a) + a*b+(b+a)*c) * x^2
+// - (d*(a*b+(b+a)*c)+a*b*c) * x + a*b*c*d = 0.
+void CoeffsForQuarticZeros(double a, double b, double c, double d,
+ double *p, double *q, double *r, double *s) {
+ *p = -(d + c + b + a);
+ *q = (d * (c + b + a) + a * b + (b + a) * c);
+ *r = -(d * (a * b + (b + a) * c) + a * b * c);
+ *s = a * b * c *d;
+}
+
+TEST(Poly, SolveCubicPolynomial) {
+ double a, b, c, aa, bb, cc;
+ double p, q, r;
+
+ a = 1; b = 2; c = 3;
+ CoeffsForCubicZeros(a, b, c, &p, &q, &r);
+ ASSERT_EQ(3, SolveCubicPolynomial(p, q, r, &aa, &bb, &cc));
+ EXPECT_NEAR(a, aa, 1e-10);
+ EXPECT_NEAR(b, bb, 1e-10);
+ EXPECT_NEAR(c, cc, 1e-10);
+
+ a = 0; b = 1; c = 3;
+ CoeffsForCubicZeros(a, b, c, &p, &q, &r);
+ ASSERT_EQ(3, SolveCubicPolynomial(p, q, r, &aa, &bb, &cc));
+ EXPECT_NEAR(a, aa, 1e-10);
+ EXPECT_NEAR(b, bb, 1e-10);
+ EXPECT_NEAR(c, cc, 1e-10);
+
+ a = -10; b = 0; c = 1;
+ CoeffsForCubicZeros(a, b, c, &p, &q, &r);
+ ASSERT_EQ(3, SolveCubicPolynomial(p, q, r, &aa, &bb, &cc));
+ EXPECT_NEAR(a, aa, 1e-10);
+ EXPECT_NEAR(b, bb, 1e-10);
+ EXPECT_NEAR(c, cc, 1e-10);
+
+ a = -8; b = 1; c = 3;
+ CoeffsForCubicZeros(a, b, c, &p, &q, &r);
+ ASSERT_EQ(3, SolveCubicPolynomial(p, q, r, &aa, &bb, &cc));
+ EXPECT_NEAR(a, aa, 1e-10);
+ EXPECT_NEAR(b, bb, 1e-10);
+ EXPECT_NEAR(c, cc, 1e-10);
+
+ a = 28; b = 28; c = 105;
+ CoeffsForCubicZeros(a, b, c, &p, &q, &r);
+ ASSERT_EQ(3, SolveCubicPolynomial(p, q, r, &aa, &bb, &cc));
+ EXPECT_NEAR(a, aa, 1e-10);
+ EXPECT_NEAR(b, bb, 1e-10);
+ EXPECT_NEAR(c, cc, 1e-10);
+}
+} // namespace
diff --git a/extern/libmv/libmv/simple_pipeline/camera_intrinsics_test.cc b/extern/libmv/libmv/simple_pipeline/camera_intrinsics_test.cc
new file mode 100644
index 00000000000..96d35a29ef8
--- /dev/null
+++ b/extern/libmv/libmv/simple_pipeline/camera_intrinsics_test.cc
@@ -0,0 +1,239 @@
+// 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/simple_pipeline/camera_intrinsics.h"
+
+#include <iostream>
+
+#include "testing/testing.h"
+#include "libmv/image/image.h"
+#include "libmv/image/image_drawing.h"
+#include "libmv/logging/logging.h"
+
+namespace libmv {
+
+TEST(PolynomialCameraIntrinsics2, ApplyOnFocalCenter) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1300.0, 1300.0);
+ intrinsics.SetPrincipalPoint(640.0, 540.0);
+ intrinsics.SetRadialDistortion(-0.2, -0.1, -0.05);
+
+ double distorted_x, distorted_y;
+ intrinsics.ApplyIntrinsics(0.0, 0.0, &distorted_x, &distorted_y);
+
+ EXPECT_NEAR(640.0, distorted_x, 1e-8);
+ EXPECT_NEAR(540.0, distorted_y, 1e-8);
+}
+
+TEST(PolynomialCameraIntrinsics, InvertOnFocalCenter) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1300.0, 1300.0);
+ intrinsics.SetPrincipalPoint(640.0, 540.0);
+ intrinsics.SetRadialDistortion(-0.2, -0.1, -0.05);
+
+ double normalized_x, normalized_y;
+ intrinsics.InvertIntrinsics(640.0, 540.0, &normalized_x, &normalized_y);
+
+ EXPECT_NEAR(0.0, normalized_x, 1e-8);
+ EXPECT_NEAR(0.0, normalized_y, 1e-8);
+}
+
+TEST(PolynomialCameraIntrinsics, ApplyIntrinsics) {
+ const int N = 5;
+
+ double expected[N][N][2] = {
+ { {75.312500, -24.687500}, {338.982239, -62.035522},
+ {640.000000, -72.929688}, {941.017761, -62.035522},
+ {1204.687500, -24.687500}},
+
+ { {37.964478, 238.982239}, {323.664551, 223.664551},
+ {640.000000, 219.193420}, {956.335449, 223.664551},
+ {1242.035522, 238.982239}},
+
+ { {27.070312, 540.000000}, {319.193420, 540.000000},
+ {640.000000, 540.000000}, {960.806580, 540.000000},
+ {1252.929688, 540.000000}},
+
+ { {37.964478, 841.017761}, {323.664551, 856.335449},
+ {640.000000, 860.806580}, {956.335449, 856.335449},
+ {1242.035522, 841.017761}},
+
+ { {75.312500, 1104.687500}, {338.982239, 1142.035522},
+ {640.000000, 1152.929688}, {941.017761, 1142.035522},
+ {1204.687500, 1104.687500}}
+ };
+
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1300.0, 1300.0);
+ intrinsics.SetPrincipalPoint(640.0, 540.0);
+ intrinsics.SetRadialDistortion(-0.2, -0.1, -0.05);
+
+ double step = 1.0 / (N - 1);
+
+ for (int i = 0; i < N; i++) {
+ for (int j = 0; j < N; j++) {
+ double normalized_x = j * step - 0.5,
+ normalized_y = i * step - 0.5;
+
+ double distorted_x, distorted_y;
+ intrinsics.ApplyIntrinsics(normalized_x, normalized_y,
+ &distorted_x, &distorted_y);
+
+ EXPECT_NEAR(expected[i][j][0], distorted_x, 1e-6);
+ EXPECT_NEAR(expected[i][j][1], distorted_y, 1e-6);
+ }
+ }
+}
+
+TEST(PolynomialCameraIntrinsics, InvertIntrinsics) {
+ const int N = 5;
+
+ double expected[N][N][2] = {
+ { {-0.524482, -0.437069}, {-0.226237, -0.403994},
+ { 0.031876, -0.398446}, { 0.293917, -0.408218},
+ { 0.632438, -0.465028}},
+
+ { {-0.493496, -0.189173}, {-0.219052, -0.179936},
+ { 0.030975, -0.178107}, { 0.283742, -0.181280},
+ { 0.574557, -0.194335}},
+
+ { {-0.488013, 0.032534}, {-0.217537, 0.031077},
+ { 0.030781, 0.030781}, { 0.281635, 0.031293},
+ { 0.566344, 0.033314}},
+
+ { {-0.498696, 0.257660}, {-0.220424, 0.244041},
+ { 0.031150, 0.241409}, { 0.285660, 0.245985},
+ { 0.582670, 0.265629}},
+
+ { {-0.550617, 0.532263}, {-0.230399, 0.477255},
+ { 0.032380, 0.469510}, { 0.299986, 0.483311},
+ { 0.684740, 0.584043}}
+ };
+
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1300.0, 1300.0);
+ intrinsics.SetPrincipalPoint(600.0, 500.0);
+ intrinsics.SetRadialDistortion(-0.2, -0.1, -0.05);
+
+ double step_x = 1280.0 / (N - 1),
+ step_y = 1080.0 / (N - 1);
+
+ for (int i = 0; i < N; i++) {
+ for (int j = 0; j < N; j++) {
+ double distorted_x = j * step_x,
+ distorted_y = i * step_y;
+
+ double normalized_x, normalized_y;
+ intrinsics.InvertIntrinsics(distorted_x, distorted_y,
+ &normalized_x, &normalized_y);
+
+ EXPECT_NEAR(expected[i][j][0], normalized_x, 1e-6);
+ EXPECT_NEAR(expected[i][j][1], normalized_y, 1e-6);
+ }
+ }
+}
+
+TEST(PolynomialCameraIntrinsics, ApplyIsInvertibleSimple) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1300.0, 1300.0);
+ intrinsics.SetPrincipalPoint(640.0, 540.0);
+ intrinsics.SetRadialDistortion(-0.2, -0.1, -0.05);
+
+ // Scan over image coordinates, invert the intrinsics, then re-apply them to
+ // make sure the cycle gets back where it started.
+ for (double y = 0; y < 1000; y += 100) {
+ for (double x = 0; x < 1000; x += 100) {
+ double normalized_x, normalized_y;
+ intrinsics.InvertIntrinsics(x, y, &normalized_x, &normalized_y);
+
+ double xp, yp;
+ intrinsics.ApplyIntrinsics(normalized_x, normalized_y, &xp, &yp);
+
+ EXPECT_NEAR(x, xp, 1e-8) << "y: " << y;
+ EXPECT_NEAR(y, yp, 1e-8) << "x: " << x;
+ LG << "Error x: " << (x - xp);
+ LG << "Error y: " << (y - yp);
+ }
+ }
+}
+
+TEST(PolynomialCameraIntrinsics, IdentityDistortBuffer) {
+ const int w = 101, h = 101;
+ FloatImage image(h, w);
+ image.Fill(0);
+
+ DrawLine(0.0, h / 2.0, w - 1, h / 2.0, 1.0, &image);
+ DrawLine(0.0, h / 4.0, w - 1, h / 4.0, 1.0, &image);
+ DrawLine(0.0, h / 4.0 * 3.0, w - 1.0, h / 4.0 * 3.0, 1.0, &image);
+ DrawLine(w / 2.0, 0.0, w / 2.0, h - 1.0, 1.0, &image);
+ DrawLine(w / 4.0, 0.0, w / 4.0, h - 1.0, 1.0, &image);
+ DrawLine(w / 4.0 * 3.0, 0.0, w / 4.0 * 3.0, h - 1.0, 1.0, &image);
+
+ PolynomialCameraIntrinsics intrinsics;
+ FloatImage distorted_image(h, w);
+ intrinsics.SetImageSize(w, h);
+ intrinsics.SetFocalLength(10.0, 10.0);
+ intrinsics.SetPrincipalPoint((double) w / 2.0, (double) h / 2.0);
+ intrinsics.SetRadialDistortion(0.0, 0.0, 0.0);
+ intrinsics.DistortBuffer(image.Data(),
+ image.Width(), image.Height(),
+ 0.0,
+ image.Depth(),
+ distorted_image.Data());
+
+ for (int x = 0; x < image.Width(); ++x) {
+ for (int y = 0; y < image.Height(); ++y) {
+ EXPECT_EQ(image(y, x), distorted_image(y, x));
+ }
+ }
+}
+
+TEST(PolynomialCameraIntrinsics, IdentityUndistortBuffer) {
+ const int w = 101, h = 101;
+ FloatImage image(h, w);
+ image.Fill(0);
+
+ DrawLine(0.0, h / 2.0, w - 1, h / 2.0, 1.0, &image);
+ DrawLine(0.0, h / 4.0, w - 1, h / 4.0, 1.0, &image);
+ DrawLine(0.0, h / 4.0 * 3.0, w - 1.0, h / 4.0 * 3.0, 1.0, &image);
+ DrawLine(w / 2.0, 0.0, w / 2.0, h - 1.0, 1.0, &image);
+ DrawLine(w / 4.0, 0.0, w / 4.0, h - 1.0, 1.0, &image);
+ DrawLine(w / 4.0 * 3.0, 0.0, w / 4.0 * 3.0, h - 1.0, 1.0, &image);
+
+ PolynomialCameraIntrinsics intrinsics;
+ FloatImage distorted_image(h, w);
+ intrinsics.SetImageSize(w, h);
+ intrinsics.SetFocalLength(10.0, 10.0);
+ intrinsics.SetPrincipalPoint((double) w / 2.0, (double) h / 2.0);
+ intrinsics.SetRadialDistortion(0.0, 0.0, 0.0);
+ intrinsics.UndistortBuffer(image.Data(),
+ image.Width(), image.Height(),
+ 0.0,
+ image.Depth(),
+ distorted_image.Data());
+
+ for (int x = 0; x < image.Width(); ++x) {
+ for (int y = 0; y < image.Height(); ++y) {
+ EXPECT_EQ(image(y, x), distorted_image(y, x));
+ }
+ }
+}
+
+} // namespace libmv
diff --git a/extern/libmv/libmv/simple_pipeline/detect_test.cc b/extern/libmv/libmv/simple_pipeline/detect_test.cc
new file mode 100644
index 00000000000..fe57e3d04a2
--- /dev/null
+++ b/extern/libmv/libmv/simple_pipeline/detect_test.cc
@@ -0,0 +1,230 @@
+// Copyright (c) 2014 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/simple_pipeline/detect.h"
+
+#include "testing/testing.h"
+#include "libmv/logging/logging.h"
+
+namespace libmv {
+
+namespace {
+
+void PreformSinglePointTest(const DetectOptions &options) {
+ // Prepare the image.
+ FloatImage image(15, 15);
+ image.fill(1.0);
+ image(7, 7) = 0.0;
+
+ // Run the detector.
+ vector<Feature> detected_features;
+ Detect(image, options, &detected_features);
+
+ // Check detected features matches our expectations.
+ EXPECT_EQ(1, detected_features.size());
+ if (detected_features.size() == 1) {
+ Feature &feature = detected_features[0];
+ EXPECT_EQ(7, feature.x);
+ EXPECT_EQ(7, feature.y);
+ }
+}
+
+void PreformCheckerBoardTest(const DetectOptions &options) {
+ // Prepare the image.
+ FloatImage image(30, 30);
+ for (int y = 0; y < image.Height(); ++y) {
+ for (int x = 0; x < image.Width(); ++x) {
+ image(y, x) = (x / 10 + y / 10) % 2 ? 1.0 : 0.0;
+ }
+ }
+
+ // Run the detector.
+ vector<Feature> detected_features;
+ Detect(image, options, &detected_features);
+
+ // Check detected features matches our expectations.
+
+ // We expect here only corners of a center square to be
+ // considered a feature points.
+ EXPECT_EQ(4, detected_features.size());
+
+ // We don't know which side of the corner detector will choose,
+ // so what we're checking here is that detected feature is from
+ // any side of the corner.
+ //
+ // This doesn't check whether there're multiple features which
+ // are placed on different sides of the same corner. The way we
+ // deal with this is requiring min_distance to be greater than 2px.
+ for (int i = 0; i < detected_features.size(); ++i) {
+ Feature &feature = detected_features[i];
+ int rounded_x = ((feature.x + 1) / 10) * 10,
+ rounded_y = ((feature.y + 1) / 10) * 10;
+ EXPECT_LE(1, std::abs(feature.x - rounded_x));
+ EXPECT_LE(1, std::abs(feature.y - rounded_y));
+ }
+}
+
+void CheckExpectedFeatures(const vector<Feature> &detected_features,
+ const vector<Feature> &expected_features) {
+ EXPECT_EQ(expected_features.size(), detected_features.size());
+
+ // That's unsafe to iterate over vectors when their lengths
+ // doesn't match. And it doesn't make any sense actually since
+ // the test will already be considered failed here.
+ if (expected_features.size() != detected_features.size()) {
+ return;
+ }
+
+ for (int i = 0; i < expected_features.size(); ++i) {
+ const Feature &extected_feature = expected_features[i];
+ bool found = false;
+ for (int j = 0; j < detected_features.size(); ++j) {
+ const Feature &detected_feature = detected_features[j];
+ if (extected_feature.x == detected_feature.x &&
+ extected_feature.y == detected_feature.y) {
+ found = true;
+ break;
+ }
+ }
+ EXPECT_TRUE(found);
+ }
+}
+
+void PreformSingleTriangleTest(const DetectOptions &options) {
+ // Prepare the image.
+ FloatImage image(15, 21);
+ image.fill(1.0);
+
+ int vertex_x = 10, vertex_y = 5;
+ for (int i = 0; i < 6; ++i) {
+ int current_x = vertex_x - i,
+ current_y = vertex_y + i;
+ for (int j = 0; j < i * 2 + 1; ++j, ++current_x) {
+ image(current_y, current_x) = 0.0;
+ }
+ }
+
+ // Run the detector.
+ vector<Feature> detected_features;
+ Detect(image, options, &detected_features);
+
+ // Check detected features matches our expectations.
+ vector<Feature> expected_features;
+ expected_features.push_back(Feature(6, 10));
+ expected_features.push_back(Feature(14, 10));
+ expected_features.push_back(Feature(10, 6));
+
+ CheckExpectedFeatures(detected_features, expected_features);
+}
+
+} // namespace
+
+#ifndef LIBMV_NO_FAST_DETECTOR
+TEST(Detect, FASTSinglePointTest) {
+ DetectOptions options;
+ options.type = DetectOptions::FAST;
+ options.min_distance = 0;
+ options.fast_min_trackness = 1;
+
+ PreformSinglePointTest(options);
+}
+#endif // LIBMV_NO_FAST_DETECTOR
+
+#if 0
+// TODO(sergey): FAST doesn't detect checker board corners, but should it?
+TEST(Detect, FASTCheckerBoardTest) {
+ DetectOptions options;
+ options.type = DetectOptions::FAST;
+ options.min_distance = 0;
+ options.fast_min_trackness = 1;
+
+ PreformCheckerBoardTest(options);
+}
+#endif
+
+#if 0
+// TODO(sergey): FAST doesn't detect triangle corners!
+TEST(Detect, FASTSingleTriangleTest) {
+ DetectOptions options;
+ options.type = DetectOptions::FAST;
+ options.margin = 3;
+ options.min_distance = 0;
+ options.fast_min_trackness = 2;
+
+ PreformSingleTriangleTest(options);
+}
+#endif
+
+#if 0
+// TODO(sergey): This doesn't actually detect single point,
+// but should it or it's expected that Moravec wouldn't consider
+// single point as feature?
+//
+// Uncomment this or remove as soon as we know answer for the
+// question.
+TEST(Detect, MoravecSinglePointTest) {
+ DetectOptions options;
+ options.type = DetectOptions::MORAVEC;
+ options.min_distance = 0;
+ options.moravec_max_count = 10;
+
+ PreformSinglePointTest(options);
+}
+
+// TODO(sergey): Moravec doesn't detect checker board corners, but should it?
+TEST(Detect, MoravecCheckerBoardTest) {
+ DetectOptions options;
+ options.type = DetectOptions::MORAVEC;
+ options.min_distance = 0;
+ options.moravec_max_count = 10;
+
+ PreformCheckerBoardTest(options);
+}
+#endif
+
+TEST(Detect, HarrisSinglePointTest) {
+ DetectOptions options;
+ options.type = DetectOptions::HARRIS;
+
+ // Set this to non-zero so image corners are not considered
+ // a feature points and avoid center point neighbors to be
+ // considered a features as well.
+ options.margin = 3;
+ options.min_distance = 3;
+
+ PreformSinglePointTest(options);
+}
+
+TEST(Detect, HarrisSingleTriangleTest) {
+ DetectOptions options;
+ options.type = DetectOptions::HARRIS;
+
+ options.margin = 3;
+ options.min_distance = 2;
+ options.harris_threshold = 1e-3;
+
+ PreformSingleTriangleTest(options);
+}
+
+// TODO(sergey): Add tests for margin option.
+
+// TODO(sergey): Add tests for min_distance option.
+
+} // namespace libmv
diff --git a/extern/libmv/libmv/simple_pipeline/intersect_test.cc b/extern/libmv/libmv/simple_pipeline/intersect_test.cc
new file mode 100644
index 00000000000..dd4fdc789af
--- /dev/null
+++ b/extern/libmv/libmv/simple_pipeline/intersect_test.cc
@@ -0,0 +1,81 @@
+// Copyright (c) 2007, 2008 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/simple_pipeline/intersect.h"
+
+#include <iostream>
+
+#include "testing/testing.h"
+#include "libmv/multiview/projection.h"
+#include "libmv/numeric/numeric.h"
+#include "libmv/logging/logging.h"
+
+namespace libmv {
+
+TEST(Intersect, EuclideanIntersect) {
+ Mat3 K1 = Mat3::Identity();
+ // K1 << 320, 0, 160,
+ // 0, 320, 120,
+ // 0, 0, 1;
+ Mat3 K2 = Mat3::Identity();
+ // K2 << 360, 0, 170,
+ // 0, 360, 110,
+ // 0, 0, 1;
+ Mat3 R1 = RotationAroundZ(-0.1);
+ Mat3 R2 = RotationAroundX(-0.1);
+ Vec3 t1; t1 << 1, 1, 10;
+ Vec3 t2; t2 << -2, -1, 10;
+ Mat34 P1, P2;
+ P_From_KRt(K1, R1, t1, &P1);
+ P_From_KRt(K2, R2, t2, &P2);
+
+ //Mat3 F; FundamentalFromProjections(P1, P2, &F);
+
+ Mat3X X;
+ X.resize(3, 30);
+ X.setRandom();
+
+ Mat2X X1, X2;
+ Project(P1, X, &X1);
+ Project(P2, X, &X2);
+
+ for (int i = 0; i < X.cols(); ++i) {
+ Vec2 x1, x2;
+ MatrixColumn(X1, i, &x1);
+ MatrixColumn(X2, i, &x2);
+ Vec3 expected;
+ MatrixColumn(X, i, &expected);
+
+ EuclideanReconstruction reconstruction;
+ reconstruction.InsertCamera(1, R1, t1);
+ reconstruction.InsertCamera(2, R2, t2);
+
+ vector<Marker> markers;
+ Marker a = { 1, 0, x1.x(), x1.y(), 1.0 };
+ markers.push_back(a);
+ Marker b = { 2, 0, x2.x(), x2.y(), 1.0 };
+ markers.push_back(b);
+
+ EuclideanIntersect(markers, &reconstruction);
+ Vec3 estimated = reconstruction.PointForTrack(0)->X;
+ EXPECT_NEAR(0, DistanceLInfinity(estimated, expected), 1e-8);
+ }
+}
+} // namespace
diff --git a/extern/libmv/libmv/simple_pipeline/keyframe_selection_test.cc b/extern/libmv/libmv/simple_pipeline/keyframe_selection_test.cc
new file mode 100644
index 00000000000..9d88362cc88
--- /dev/null
+++ b/extern/libmv/libmv/simple_pipeline/keyframe_selection_test.cc
@@ -0,0 +1,307 @@
+// 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/simple_pipeline/keyframe_selection.h"
+
+#include "testing/testing.h"
+#include "libmv/simple_pipeline/camera_intrinsics.h"
+#include "libmv/logging/logging.h"
+
+namespace libmv {
+
+// Synthetic test, N markers with the same translation
+// Should not be keyframe
+TEST(KeyframeSelection, SyntheticNeighborFrame) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(900.0,900.0);
+ intrinsics.SetPrincipalPoint(640.0, 540.0);
+ intrinsics.SetRadialDistortion(0.0, 0.0, 0.0);
+
+ Tracks tracks;
+ const int markers_per_size = 15;
+
+ // Fill in tracks for homography estimation
+ for (int x = 0; x < markers_per_size; x++) {
+ for (int y = 0; y < markers_per_size; y++) {
+ double current_x = 10 + x * 40, current_y = 10 + y * 40;
+ double next_x = current_x + 10, next_y = current_y + 10;
+
+ intrinsics.InvertIntrinsics(current_x, current_y, &current_x, &current_y);
+ intrinsics.InvertIntrinsics(next_x, next_y, &next_x, &next_y);
+
+ tracks.Insert(1, y * markers_per_size + x, current_x, current_y);
+ tracks.Insert(2, y * markers_per_size + x, next_x, next_y);
+ }
+ }
+
+ vector<int> keyframes;
+ SelectKeyframesBasedOnGRICAndVariance(tracks, intrinsics, keyframes);
+
+ // Synthetic second frame shouldn't be considered a keyframe
+ EXPECT_EQ(0, keyframes.size());
+}
+
+// Frames 1 and 2 of FabrikEingang footage
+// Only one wall is tracked, should not be keyframes
+TEST(KeyframeSelection, FabrikEingangNeighborFrames) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1605.797, 1605.797);
+ intrinsics.SetPrincipalPoint(960.000, 544.000);
+ intrinsics.SetRadialDistortion(0.0, 0.0, 0.0);
+
+ Marker markers[] = {
+ {1, 0, 737.599983, 646.397594, 1.0}, {2, 0, 737.906628, 648.113327, 1.0}, {1, 1, 863.045425, 646.081905, 1.0},
+ {2, 1, 863.339767, 647.650040, 1.0}, {1, 2, 736.959972, 574.080151, 1.0}, {2, 2, 737.217350, 575.604900, 1.0},
+ {1, 3, 864.097424, 573.374908, 1.0}, {2, 3, 864.383469, 574.900307, 1.0}, {1, 4, 789.429073, 631.677521, 1.0},
+ {2, 4, 789.893131, 633.124451, 1.0}, {1, 5, 791.051960, 573.442028, 1.0}, {2, 5, 791.336575, 575.088890, 1.0},
+ {1, 6, 738.973961, 485.130308, 1.0}, {2, 6, 739.435501, 486.734207, 1.0}, {1, 7, 862.403240, 514.866074, 1.0},
+ {2, 7, 862.660618, 516.413261, 1.0}, {1, 8, 802.240162, 485.759838, 1.0}, {2, 8, 802.602253, 487.432899, 1.0},
+ {1, 9, 754.340630, 500.624559, 1.0}, {2, 9, 754.559956, 502.079920, 1.0}, {1, 10, 849.398689, 484.480545, 1.0},
+ {2, 10, 849.599934, 486.079937, 1.0}, {1, 11, 788.803768, 515.924391, 1.0}, {2, 11, 789.119911, 517.439932, 1.0},
+ {1, 12, 838.733940, 558.212688, 1.0}, {2, 12, 839.039898, 559.679916, 1.0}, {1, 13, 760.014782, 575.194466, 1.0},
+ {2, 13, 760.319881, 576.639904, 1.0}, {1, 14, 765.321636, 616.015957, 1.0}, {2, 14, 765.759945, 617.599915, 1.0},
+ {1, 15, 800.963230, 660.032082, 1.0}, {2, 15, 801.279945, 661.759876, 1.0}, {1, 16, 846.321087, 602.313053, 1.0},
+ {2, 16, 846.719913, 603.839878, 1.0}, {1, 17, 864.288311, 616.790524, 1.0}, {2, 17, 864.639931, 618.239918, 1.0},
+ {1, 18, 800.006790, 602.573425, 1.0}, {2, 18, 800.319958, 604.159912, 1.0}, {1, 19, 739.026890, 617.944138, 1.0},
+ {2, 19, 739.199924, 619.519924, 1.0}, {1, 20, 801.987419, 544.134888, 1.0}, {2, 20, 802.239933, 545.599911, 1.0},
+ {1, 21, 753.619823, 542.961300, 1.0}, {2, 21, 753.919945, 544.639874, 1.0}, {1, 22, 787.921257, 499.910206, 1.0},
+ {2, 22, 788.159924, 501.439917, 1.0}, {1, 23, 839.095459, 529.287903, 1.0}, {2, 23, 839.359932, 530.879934, 1.0},
+ {1, 24, 811.760330, 630.732269, 1.0}, {2, 24, 812.159901, 632.319859, 1.0}
+ };
+ int num_markers = sizeof(markers) / sizeof(Marker);
+
+ Tracks tracks;
+ for (int i = 0; i < num_markers; i++) {
+ double x = markers[i].x, y = markers[i].y;
+ intrinsics.InvertIntrinsics(x, y, &x, &y);
+ tracks.Insert(markers[i].image, markers[i].track, x, y);
+ }
+
+ vector<int> keyframes;
+ SelectKeyframesBasedOnGRICAndVariance(tracks, intrinsics, keyframes);
+
+ EXPECT_EQ(0, keyframes.size());
+}
+
+// Frames 120 and 200 from FabrikEingang footage
+// Should be enough of parallax for keyframing
+TEST(KeyframeSelection, FabrikEingangFarFrames) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1605.797, 1605.797);
+ intrinsics.SetPrincipalPoint(960.000, 544.000);
+ intrinsics.SetRadialDistortion(0.0, 0.0, 0.0);
+
+ Marker markers[] = {
+ {1, 0, 369.459200, 619.315258, 1.0}, {2, 0, 279.677496, 722.086842, 1.0}, {1, 1, 376.831970, 370.278397, 1.0},
+ {2, 1, 221.695247, 460.065418, 1.0}, {1, 2, 1209.139023, 567.705605, 1.0}, {2, 2, 1080.760117, 659.230083, 1.0},
+ {1, 3, 1643.495750, 903.620453, 1.0}, {2, 3, 1618.405037, 1015.374908, 1.0}, {1, 4, 1494.849815, 425.302460, 1.0},
+ {2, 4, 1457.467575, 514.727587, 1.0}, {1, 5, 1794.637299, 328.728609, 1.0}, {2, 5, 1742.161446, 408.988636, 1.0},
+ {1, 6, 1672.822723, 102.240358, 1.0}, {2, 6, 1539.287224, 153.536892, 1.0}, {1, 7, 1550.843925, 53.424943, 1.0},
+ {2, 7, 1385.579109, 96.450085, 1.0}, {1, 8, 852.953281, 465.399578, 1.0}, {2, 8, 779.404564, 560.091843, 1.0},
+ {1, 9, 906.853752, 299.827040, 1.0}, {2, 9, 786.923218, 385.570770, 1.0}, {1, 10, 406.322966, 87.556041, 1.0},
+ {2, 10, 140.339413, 150.877481, 1.0}, {1, 11, 254.811573, 851.296478, 1.0}, {2, 11, 94.478302, 969.350189, 1.0},
+ {1, 12, 729.087868, 806.092758, 1.0}, {2, 12, 606.212139, 919.876560, 1.0}, {1, 13, 1525.719452, 920.398083, 1.0},
+ {2, 13, 1495.579720, 1031.971218, 1.0}
+ };
+ int num_markers = sizeof(markers) / sizeof(Marker);
+
+ Tracks tracks;
+ for (int i = 0; i < num_markers; i++) {
+ double x = markers[i].x, y = markers[i].y;
+ intrinsics.InvertIntrinsics(x, y, &x, &y);
+ tracks.Insert(markers[i].image, markers[i].track, x, y);
+ }
+
+ vector<int> keyframes;
+ SelectKeyframesBasedOnGRICAndVariance(tracks, intrinsics, keyframes);
+
+ EXPECT_EQ(2, keyframes.size());
+}
+
+// Manually selected keyframes from copter footage from Sebastian
+// Keyframes were 167 and 237
+TEST(KeyframeSelection, CopterManualKeyFrames) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1155.043, 1155.043);
+ intrinsics.SetPrincipalPoint(640.000, 360.000);
+ intrinsics.SetRadialDistortion(-0.08590, 0.0, 0.0);
+
+ Marker markers[] = {
+ {1, 0, 645.792694, 403.115931, 1.0}, {2, 0, 630.641174, 307.996409, 1.0}, {1, 1, 783.469086, 403.904328, 1.0},
+ {2, 1, 766.001129, 308.998225, 1.0}, {1, 2, 650.000000, 160.000001, 1.0}, {1, 3, 785.225906, 158.619039, 1.0},
+ {2, 3, 767.526474, 70.449695, 1.0}, {1, 4, 290.640526, 382.335634, 1.0}, {2, 4, 273.001728, 86.993319, 1.0},
+ {1, 5, 291.162739, 410.602684, 1.0}, {2, 5, 273.287849, 111.937487, 1.0}, {1, 6, 136.919317, 349.895797, 1.0},
+ {1, 7, 490.844345, 47.572222, 1.0}, {1, 8, 454.406433, 488.935761, 1.0}, {1, 9, 378.655815, 618.522248, 1.0},
+ {2, 9, 357.061806, 372.265077, 1.0}, {1, 10, 496.011391, 372.668824, 1.0}, {2, 10, 477.979164, 222.986112, 1.0},
+ {1, 11, 680.060272, 256.103625, 1.0}, {2, 11, 670.587540, 204.830453, 1.0}, {1, 12, 1070.817108, 218.775322, 1.0},
+ {2, 12, 1046.129913, 128.969783, 1.0}, {1, 14, 242.516403, 596.048512, 1.0}, {2, 14, 224.182606, 248.272154, 1.0},
+ {1, 15, 613.936272, 287.519073, 1.0}, {2, 15, 600.467644, 196.085722, 1.0}, {1, 31, 844.637451, 256.354315, 1.0},
+ {2, 31, 823.200150, 165.714952, 1.0},
+ };
+ int num_markers = sizeof(markers) / sizeof(Marker);
+
+ Tracks tracks;
+ for (int i = 0; i < num_markers; i++) {
+ double x = markers[i].x, y = markers[i].y;
+ intrinsics.InvertIntrinsics(x, y, &x, &y);
+ tracks.Insert(markers[i].image, markers[i].track, x, y);
+ }
+
+ vector<int> keyframes;
+ SelectKeyframesBasedOnGRICAndVariance(tracks, intrinsics, keyframes);
+
+ EXPECT_EQ(2, keyframes.size());
+}
+
+// Used old friend elevator scene MMI_2366 with automatic feature selection
+// and manual outlier elimination and manual keyframe selection
+// Selected keyframes were 29 and 41
+TEST(KeyframeSelection, ElevatorManualKeyframesFrames) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1380.000, 1380.000);
+ intrinsics.SetPrincipalPoint(960.000, 540.000);
+ intrinsics.SetRadialDistortion(-0.034, 0.0, 0.0);
+
+ Marker markers[] = {
+ {1, 2, 1139.861412, 1034.634984, 1.0}, {2, 2, 1143.512192, 1065.355718, 1.0}, {1, 3, 1760.821953, 644.658036, 1.0},
+ {2, 3, 1770.901108, 697.899928, 1.0}, {1, 4, 858.071823, 1068.520746, 1.0}, {1, 6, 1633.952408, 797.050145, 1.0},
+ {2, 6, 1642.508469, 849.157140, 1.0}, {1, 8, 1716.695824, 451.805491, 1.0}, {2, 8, 1726.513939, 502.095687, 1.0},
+ {1, 9, 269.577627, 724.986935, 1.0}, {2, 9, 269.424820, 764.154246, 1.0}, {1, 10, 1891.321907, 706.948843, 1.0},
+ {2, 10, 1903.338547, 766.068377, 1.0}, {1, 12, 1806.227074, 956.089604, 1.0}, {2, 12, 1816.619568, 1013.767376, 1.0},
+ {1, 14, 269.544153, 1002.333570, 1.0}, {2, 14, 269.367542, 1043.509254, 1.0}, {1, 15, 1402.772141, 281.392962, 1.0},
+ {2, 15, 1409.089165, 318.731629, 1.0}, {1, 16, 195.877233, 919.454341, 1.0}, {2, 16, 192.531109, 997.367899, 1.0},
+ {1, 17, 1789.584045, 120.036661, 1.0}, {2, 17, 1800.391846, 167.822964, 1.0}, {1, 18, 999.363213, 765.004807, 1.0},
+ {2, 18, 1002.345772, 790.560122, 1.0}, {1, 19, 647.342491, 1044.805727, 1.0}, {2, 19, 649.328041, 1058.682940, 1.0},
+ {1, 20, 1365.486832, 440.901829, 1.0}, {2, 20, 1371.413040, 477.888730, 1.0}, {1, 21, 1787.125282, 301.431606, 1.0},
+ {2, 21, 1798.527260, 355.224531, 1.0}, {1, 22, 1257.805824, 932.797258, 1.0}, {2, 22, 1263.017578, 969.376774, 1.0},
+ {1, 23, 961.969528, 843.148112, 1.0}, {2, 23, 964.869461, 868.587620, 1.0}, {1, 24, 158.076110, 1052.643592, 1.0},
+ {1, 25, 1072.884521, 1005.296981, 1.0}, {2, 25, 1076.091156, 1032.776856, 1.0}, {1, 26, 1107.656937, 526.577228, 1.0},
+ {2, 26, 1111.618423, 555.524454, 1.0}, {1, 27, 1416.410751, 529.857581, 1.0}, {2, 27, 1422.663574, 570.025957, 1.0},
+ {1, 28, 1498.673630, 1005.453086, 1.0}, {2, 28, 1505.381813, 1051.827149, 1.0}, {1, 29, 1428.647804, 652.473629, 1.0},
+ {2, 29, 1434.898224, 692.715390, 1.0}, {1, 30, 1332.318764, 503.673599, 1.0}, {2, 30, 1338.000069, 540.507967, 1.0},
+ {1, 32, 1358.642693, 709.837904, 1.0}, {2, 32, 1364.231529, 748.678265, 1.0}, {1, 33, 1850.911560, 460.475668, 1.0},
+ {2, 33, 1862.221413, 512.797347, 1.0}, {1, 34, 1226.117821, 607.053959, 1.0}, {2, 34, 1230.736084, 641.091449, 1.0},
+ {1, 35, 619.598236, 523.341744, 1.0}, {2, 35, 621.601124, 554.453287, 1.0}, {1, 36, 956.591492, 958.223183, 1.0},
+ {2, 36, 959.289265, 983.289263, 1.0}, {1, 37, 1249.922218, 419.095856, 1.0}, {2, 37, 1255.005913, 452.556177, 1.0},
+ {1, 39, 1300.528450, 386.251166, 1.0}, {2, 39, 1305.957413, 420.185595, 1.0}, {1, 40, 1128.689919, 972.558346, 1.0},
+ {2, 40, 1132.413712, 1003.984737, 1.0}, {1, 41, 503.304749, 1053.504388, 1.0}, {2, 41, 505.019703, 1069.175613, 1.0},
+ {1, 42, 1197.352982, 472.681564, 1.0}, {2, 42, 1201.910706, 503.459880, 1.0}, {1, 43, 1794.391022, 383.911400, 1.0},
+ {2, 43, 1805.324135, 436.116468, 1.0}, {1, 44, 789.641418, 1058.045647, 1.0}, {1, 45, 1376.575241, 928.714979, 1.0},
+ {2, 45, 1381.995850, 969.511957, 1.0}, {1, 46, 1598.023567, 93.975592, 1.0}, {2, 46, 1606.937141, 136.827035, 1.0},
+ {1, 47, 1455.550232, 762.128685, 1.0}, {2, 47, 1462.014313, 805.164878, 1.0}, {1, 48, 1357.123489, 354.460326, 1.0},
+ {2, 48, 1363.071899, 390.363121, 1.0}, {1, 49, 939.792652, 781.549895, 1.0}, {2, 49, 942.802620, 806.164012, 1.0},
+ {1, 50, 1380.251083, 805.948620, 1.0}, {2, 50, 1385.637932, 845.592098, 1.0}, {1, 51, 1021.769943, 1049.802361, 1.0},
+ {1, 52, 1065.634918, 608.099055, 1.0}, {2, 52, 1069.142189, 635.361736, 1.0}, {1, 53, 624.324188, 463.202863, 1.0},
+ {2, 53, 626.395454, 494.994088, 1.0}, {1, 54, 1451.459885, 881.557624, 1.0}, {2, 54, 1457.679634, 924.345531, 1.0},
+ {1, 55, 1201.885986, 1057.079022, 1.0}, {1, 56, 581.157532, 947.661438, 1.0}, {2, 56, 583.242359, 960.831449, 1.0},
+ {1, 58, 513.593102, 954.175858, 1.0}, {2, 58, 515.470047, 971.309574, 1.0}, {1, 59, 928.069038, 901.774421, 1.0},
+ {2, 59, 930.847950, 925.613744, 1.0}, {1, 60, 1065.860023, 740.395389, 1.0}, {2, 60, 1069.484253, 768.971086, 1.0},
+ {1, 61, 990.479393, 906.264632, 1.0}, {2, 61, 993.217506, 933.088803, 1.0}, {1, 62, 1776.196747, 776.278453, 1.0},
+ {2, 62, 1786.292496, 831.136880, 1.0}, {1, 63, 834.454365, 1012.449725, 1.0}, {2, 63, 836.868324, 1033.451807, 1.0},
+ {1, 64, 1355.190697, 869.184809, 1.0}, {2, 64, 1360.736618, 909.773347, 1.0}, {1, 65, 702.072487, 897.519686, 1.0},
+ {2, 65, 704.203377, 911.931131, 1.0}, {1, 66, 1214.022903, 856.199934, 1.0}, {2, 66, 1218.109016, 890.753052, 1.0},
+ {1, 67, 327.676048, 236.814036, 1.0}, {2, 67, 328.335285, 277.251878, 1.0}, {1, 68, 289.064083, 454.793912, 1.0},
+ {2, 68, 288.651924, 498.882444, 1.0}, {1, 69, 1626.240692, 278.374350, 1.0}, {2, 69, 1634.131508, 315.853672, 1.0},
+ {1, 70, 1245.375710, 734.862142, 1.0}, {2, 70, 1250.047417, 769.670885, 1.0}, {1, 71, 497.015305, 510.718904, 1.0},
+ {2, 71, 498.682308, 541.070201, 1.0}, {1, 72, 1280.542030, 153.939185, 1.0}, {2, 72, 1286.993637, 198.436196, 1.0},
+ {1, 73, 1534.748840, 138.601043, 1.0}, {2, 73, 1542.961349, 180.170819, 1.0}, {1, 74, 1477.412682, 200.608061, 1.0},
+ {2, 74, 1484.683914, 240.413260, 1.0}, {1, 76, 450.637321, 407.279642, 1.0}, {2, 76, 451.695642, 441.666291, 1.0},
+ {1, 78, 246.981239, 220.786298, 1.0}, {2, 78, 244.524879, 290.016564, 1.0}, {1, 79, 36.696489, 420.023407, 1.0},
+ {2, 79, 21.364746, 591.245492, 1.0},
+ };
+ int num_markers = sizeof(markers) / sizeof(Marker);
+
+ Tracks tracks;
+ for (int i = 0; i < num_markers; i++) {
+ double x = markers[i].x, y = markers[i].y;
+ intrinsics.InvertIntrinsics(x, y, &x, &y);
+ tracks.Insert(markers[i].image, markers[i].track, x, y);
+ }
+
+ vector<int> keyframes;
+ SelectKeyframesBasedOnGRICAndVariance(tracks, intrinsics, keyframes);
+
+ EXPECT_EQ(2, keyframes.size());
+}
+
+// Elevator scene MMI_2366 with manual tracks, frames 1, 2, 3, 5 and 27
+TEST(KeyframeSelection, ElevatorReconstructionVarianceTest) {
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(1380.000, 1380.000);
+ intrinsics.SetPrincipalPoint(960.000, 540.000);
+ intrinsics.SetRadialDistortion(-0.034, 0.0, 0.0);
+
+ Marker markers[] = {
+ {1, 0, 182.999997, 1047.000010, 1.0}, {2, 0, 181.475730, 1052.091079, 1.0}, {3, 0, 181.741562, 1057.893341, 1.0},
+ {4, 0, 183.190498, 1068.310440, 1.0}, {1, 1, 271.000013, 666.000009, 1.0}, {2, 1, 270.596180, 668.665760, 1.0},
+ {3, 1, 270.523510, 671.559069, 1.0}, {4, 1, 271.856518, 676.818151, 1.0}, {5, 1, 268.989000, 727.051570, 1.0},
+ {1, 2, 264.999990, 1018.000031, 1.0}, {2, 2, 264.020061, 1021.157591, 1.0}, {3, 2, 264.606056, 1024.823506, 1.0},
+ {4, 2, 266.200933, 1031.168690, 1.0}, {1, 3, 270.000000, 938.000014, 1.0}, {2, 3, 269.022617, 941.153390, 1.0},
+ {3, 3, 269.605579, 944.454954, 1.0}, {4, 3, 271.281366, 949.452167, 1.0}, {5, 3, 268.963480, 1004.417453, 1.0},
+ {1, 4, 200.999994, 799.000003, 1.0}, {2, 4, 199.841366, 803.891838, 1.0}, {3, 4, 200.262208, 809.323246, 1.0},
+ {4, 4, 201.456513, 819.271195, 1.0}, {5, 4, 195.026493, 924.363234, 1.0}, {1, 5, 1775.000038, 49.999998, 1.0},
+ {2, 5, 1775.255127, 53.718264, 1.0}, {3, 5, 1776.449890, 55.951670, 1.0}, {4, 5, 1778.815727, 61.923309, 1.0},
+ {5, 5, 1790.274124, 123.074923, 1.0}, {1, 6, 164.000001, 927.999988, 1.0}, {2, 6, 162.665462, 933.169527, 1.0},
+ {3, 6, 163.067923, 938.577182, 1.0}, {4, 6, 164.370360, 948.840945, 1.0}, {5, 6, 157.199407, 1057.762341, 1.0},
+ {1, 7, 618.000011, 477.999998, 1.0}, {2, 7, 617.583504, 480.124243, 1.0}, {3, 7, 618.356495, 482.441897, 1.0},
+ {4, 7, 619.792500, 486.428132, 1.0}, {5, 7, 619.546051, 525.222627, 1.0}, {1, 8, 499.999981, 1036.999984, 1.0},
+ {2, 8, 499.080162, 1038.720160, 1.0}, {3, 8, 499.949398, 1039.014344, 1.0}, {4, 8, 501.828003, 1041.286647, 1.0},
+ {5, 8, 502.777576, 1055.196369, 1.0}, {1, 9, 1587.000046, 31.999999, 1.0}, {2, 9, 1586.988373, 34.635853, 1.0},
+ {3, 9, 1588.155899, 37.444186, 1.0}, {4, 9, 1589.973106, 42.492081, 1.0}, {5, 9, 1598.683205, 96.526332, 1.0},
+ {1, 10, 622.999992, 416.999999, 1.0}, {2, 10, 622.449017, 419.233485, 1.0}, {3, 10, 623.283234, 421.500703, 1.0},
+ {4, 10, 624.620132, 425.537406, 1.0}, {5, 10, 624.290829, 465.078338, 1.0}, {1, 11, 577.999992, 931.999998, 1.0},
+ {2, 11, 577.042294, 932.872703, 1.0}, {3, 11, 577.832451, 934.045451, 1.0}, {4, 11, 579.729137, 935.735435, 1.0},
+ {5, 11, 580.691242, 948.396256, 1.0}, {1, 12, 510.999985, 931.999998, 1.0}, {2, 12, 510.111237, 933.152146, 1.0},
+ {3, 12, 510.797081, 934.454219, 1.0}, {4, 12, 512.647362, 936.595910, 1.0}, {5, 12, 513.247204, 955.144157, 1.0},
+ {1, 13, 330.459995, 177.059993, 1.0}, {2, 13, 329.876347, 179.615586, 1.0}, {3, 13, 330.681696, 182.757810, 1.0},
+ {4, 13, 331.345053, 187.903853, 1.0}, {5, 13, 327.824135, 239.611639, 1.0}, {1, 14, 291.813097, 388.516195, 1.0},
+ {2, 14, 290.984058, 391.382725, 1.0}, {3, 14, 291.526737, 394.778595, 1.0}, {4, 14, 292.763815, 400.310973, 1.0},
+ {5, 14, 288.714552, 457.548015, 1.0}, {1, 15, 496.491680, 466.534005, 1.0}, {2, 15, 495.909519, 468.518561, 1.0},
+ {3, 15, 496.588383, 470.853596, 1.0}, {4, 15, 497.976780, 474.731458, 1.0}, {5, 15, 496.998882, 512.568694, 1.0},
+ {1, 16, 1273.000031, 89.000000, 1.0}, {2, 16, 1272.951965, 92.003637, 1.0}, {3, 16, 1273.934784, 94.972191, 1.0},
+ {4, 16, 1275.493584, 100.139952, 1.0}, {5, 16, 1281.003571, 156.880163, 1.0}, {1, 17, 1524.713173, 78.852922, 1.0},
+ {2, 17, 1524.782066, 81.427142, 1.0}, {3, 17, 1525.759048, 84.057939, 1.0}, {4, 17, 1527.579689, 88.966550, 1.0},
+ {5, 17, 1535.262451, 141.186054, 1.0}, {1, 18, 1509.425011, 94.371824, 1.0}, {1, 19, 451.000013, 357.000003, 1.0},
+ {2, 19, 450.354881, 359.312410, 1.0}, {3, 19, 451.107473, 361.837088, 1.0}, {4, 19, 452.186537, 366.318061, 1.0},
+ {5, 19, 450.507660, 409.257599, 1.0}, {1, 20, 254.004936, 114.784185, 1.0}, {2, 20, 253.291512, 119.288486, 1.0},
+ {3, 20, 253.745584, 124.114957, 1.0}, {4, 20, 254.453287, 132.795120, 1.0}, {5, 20, 246.772242, 225.165337, 1.0},
+ {1, 21, 65.262880, 147.889409, 1.0}, {2, 21, 63.634465, 157.656807, 1.0}, {3, 21, 63.306799, 169.067053, 1.0},
+ {4, 21, 62.462311, 189.724241, 1.0}, {5, 21, 35.396615, 430.308380, 1.0},
+ };
+ int num_markers = sizeof(markers) / sizeof(Marker);
+
+ Tracks tracks;
+ for (int i = 0; i < num_markers; i++) {
+ double x = markers[i].x, y = markers[i].y;
+ intrinsics.InvertIntrinsics(x, y, &x, &y);
+ tracks.Insert(markers[i].image, markers[i].track, x, y);
+ }
+
+ vector<int> keyframes;
+ SelectKeyframesBasedOnGRICAndVariance(tracks, intrinsics, keyframes);
+
+ EXPECT_EQ(2, keyframes.size());
+ if (keyframes.size() == 2) {
+ EXPECT_EQ(1, keyframes[0]);
+ EXPECT_EQ(5, keyframes[1]);
+ }
+}
+
+} // namespace libmv
diff --git a/extern/libmv/libmv/simple_pipeline/modal_solver_test.cc b/extern/libmv/libmv/simple_pipeline/modal_solver_test.cc
new file mode 100644
index 00000000000..8b87acd95bb
--- /dev/null
+++ b/extern/libmv/libmv/simple_pipeline/modal_solver_test.cc
@@ -0,0 +1,79 @@
+// Copyright (c) 2013 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/simple_pipeline/modal_solver.h"
+
+#include "testing/testing.h"
+#include "libmv/logging/logging.h"
+#include "libmv/simple_pipeline/bundle.h"
+#include "libmv/simple_pipeline/camera_intrinsics.h"
+
+#include <stdio.h>
+
+namespace libmv {
+
+TEST(ModalSolver, SyntheticCubeSceneMotion) {
+ double kTolerance = 1e-8;
+
+ PolynomialCameraIntrinsics intrinsics;
+ intrinsics.SetFocalLength(658.286, 658.286);
+ intrinsics.SetPrincipalPoint(480.0, 270.0);
+ intrinsics.SetRadialDistortion(0.0, 0.0, 0.0);
+
+ Marker markers[] = {
+ {1, 0, 212.172775, 354.713538, 1.0}, {2, 0, 773.468399, 358.735306, 1.0},
+ {1, 1, 62.415197, 287.905354, 1.0}, {2, 1, 619.103336, 324.402537, 1.0},
+ {1, 2, 206.847939, 237.567925, 1.0}, {2, 2, 737.496986, 247.881383, 1.0},
+ {1, 3, 351.743889, 316.415906, 1.0}, {2, 3, 908.779621, 290.703617, 1.0},
+ {1, 4, 232.941413, 54.265443, 1.0}, {2, 4, 719.444847, 63.062531, 1.0},
+ {1, 5, 96.391611, 119.283537, 1.0}, {2, 5, 611.413136, 160.890715, 1.0},
+ {1, 6, 363.444958, 150.838144, 1.0}, {2, 6, 876.374531, 114.916206, 1.0},
+ };
+ int num_markers = sizeof(markers) / sizeof(Marker);
+
+ Tracks tracks;
+ for (int i = 0; i < num_markers; i++) {
+ double x = markers[i].x, y = markers[i].y;
+ intrinsics.InvertIntrinsics(x, y, &x, &y);
+ tracks.Insert(markers[i].image, markers[i].track, x, y);
+ }
+
+ EuclideanReconstruction reconstruction;
+ ModalSolver(tracks, &reconstruction);
+ EuclideanBundleCommonIntrinsics(tracks,
+ BUNDLE_NO_INTRINSICS,
+ BUNDLE_NO_TRANSLATION,
+ &reconstruction,
+ &intrinsics,
+ NULL);
+
+ Mat3 expected_rotation;
+ expected_rotation << 0.98215101299251, 0.17798357184544, 0.06083778292258,
+ -0.16875286001759, 0.97665299913606, -0.13293378620359,
+ -0.08307743323957, 0.12029450291547, 0.98925596922871;
+
+ Mat3 &first_camera_R = reconstruction.CameraForImage(1)->R;
+ Mat3 &second_camera_R = reconstruction.CameraForImage(2)->R;
+
+ EXPECT_TRUE(Mat3::Identity().isApprox(first_camera_R, kTolerance));
+ EXPECT_TRUE(expected_rotation.isApprox(second_camera_R, kTolerance));
+}
+
+} // namespace libmv
diff --git a/extern/libmv/libmv/simple_pipeline/resect_test.cc b/extern/libmv/libmv/simple_pipeline/resect_test.cc
new file mode 100644
index 00000000000..811edd282d8
--- /dev/null
+++ b/extern/libmv/libmv/simple_pipeline/resect_test.cc
@@ -0,0 +1,234 @@
+// 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.
+
+#include "libmv/simple_pipeline/resect.h"
+#include "libmv/logging/logging.h"
+#include "testing/testing.h"
+
+#if 0
+// Generates all necessary inputs and expected outputs for EuclideanResection.
+void CreateCameraSystem(const Mat3& KK,
+ const Mat3X& x_image,
+ const Vec& X_distances,
+ const Mat3& R_input,
+ const Vec3& T_input,
+ Mat2X *x_camera,
+ Mat3X *X_world,
+ Mat3 *R_expected,
+ Vec3 *T_expected) {
+ int num_points = x_image.cols();
+
+ Mat3X x_unit_cam(3, num_points);
+ x_unit_cam = KK.inverse() * x_image;
+
+ // Create normalized camera coordinates to be used as an input to the PnP
+ // function, instead of using NormalizeColumnVectors(&x_unit_cam).
+ *x_camera = x_unit_cam.block(0, 0, 2, num_points);
+ for (int i = 0; i < num_points; ++i) {
+ x_unit_cam.col(i).normalize();
+ }
+
+ // Create the 3D points in the camera system.
+ Mat X_camera(3, num_points);
+ for (int i = 0; i < num_points; ++i) {
+ X_camera.col(i) = X_distances(i) * x_unit_cam.col(i);
+ }
+
+ // Apply the transformation to the camera 3D points
+ Mat translation_matrix(3, num_points);
+ translation_matrix.row(0).setConstant(T_input(0));
+ translation_matrix.row(1).setConstant(T_input(1));
+ translation_matrix.row(2).setConstant(T_input(2));
+
+ *X_world = R_input * X_camera + translation_matrix;
+
+ // Create the expected result for comparison.
+ *R_expected = R_input.transpose();
+ *T_expected = *R_expected * (-T_input);
+};
+
+TEST(AbsoluteOrientation, QuaternionSolution) {
+ int num_points = 4;
+ Mat X;
+ Mat Xp;
+ X = 100 * Mat::Random(3, num_points);
+
+ // Create a random translation and rotation.
+ Mat3 R_input;
+ R_input = Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitY())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ());
+
+ Vec3 t_input;
+ t_input.setRandom();
+ t_input = 100 * t_input;
+
+ Mat translation_matrix(3, num_points);
+ translation_matrix.row(0).setConstant(t_input(0));
+ translation_matrix.row(1).setConstant(t_input(1));
+ translation_matrix.row(2).setConstant(t_input(2));
+
+ // Create the transformed 3D points Xp as Xp = R * X + t.
+ Xp = R_input * X + translation_matrix;
+
+ // Output variables.
+ Mat3 R;
+ Vec3 t;
+
+ AbsoluteOrientation(X, Xp, &R, &t);
+
+ EXPECT_MATRIX_NEAR(t, t_input, 1e-6);
+ EXPECT_MATRIX_NEAR(R, R_input, 1e-8);
+}
+
+TEST(EuclideanResection, Points4KnownImagePointsRandomTranslationRotation) {
+ // In this test only the translation and rotation are random. The image
+ // points are selected from a real case and are well conditioned.
+ Vec2i image_dimensions;
+ image_dimensions << 1600, 1200;
+
+ Mat3 KK;
+ KK << 2796, 0, 804,
+ 0 , 2796, 641,
+ 0, 0, 1;
+
+ // The real image points.
+ int num_points = 4;
+ Mat3X x_image(3, num_points);
+ x_image << 1164.06, 734.948, 749.599, 430.727,
+ 681.386, 844.59, 496.315, 580.775,
+ 1, 1, 1, 1;
+
+
+ // A vector of the 4 distances to the 3D points.
+ Vec X_distances = 100 * Vec::Random(num_points).array().abs();
+
+ // Create the random camera motion R and t that resection should recover.
+ Mat3 R_input;
+ R_input = Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitY())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ());
+
+ Vec3 T_input;
+ T_input.setRandom();
+ T_input = 100 * T_input;
+
+ // Create the camera system, also getting the expected result of the
+ // transformation.
+ Mat3 R_expected;
+ Vec3 T_expected;
+ Mat3X X_world;
+ Mat2X x_camera;
+ CreateCameraSystem(KK, x_image, X_distances, R_input, T_input,
+ &x_camera, &X_world, &R_expected, &T_expected);
+
+ // Finally, run the code under test.
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_ANSAR_DANIILIDIS);
+
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+
+ // For now, the EPnP doesn't have a non-linear optimization step and so is
+ // not precise enough with only 4 points.
+ //
+ // TODO(jmichot): Reenable this test when there is nonlinear refinement.
+#if 0
+ R_output.setIdentity();
+ T_output.setZero();
+
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_EPNP);
+
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);*/
+#endif
+}
+
+// TODO(jmichot): Reduce the code duplication here with the code above.
+TEST(EuclideanResection, Points6AllRandomInput) {
+ Mat3 KK;
+ KK << 2796, 0, 804,
+ 0 , 2796, 641,
+ 0, 0, 1;
+
+ // Create random image points for a 1600x1200 image.
+ int w = 1600;
+ int h = 1200;
+ int num_points = 6;
+ Mat3X x_image(3, num_points);
+ x_image.row(0) = w * Vec::Random(num_points).array().abs();
+ x_image.row(1) = h * Vec::Random(num_points).array().abs();
+ x_image.row(2).setOnes();
+
+ // Normalized camera coordinates to be used as an input to the PnP function.
+ Mat2X x_camera;
+ Vec X_distances = 100 * Vec::Random(num_points).array().abs();
+
+ // Create the random camera motion R and t that resection should recover.
+ Mat3 R_input;
+ R_input = Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitY())
+ * Eigen::AngleAxisd(rand(), Eigen::Vector3d::UnitZ());
+
+ Vec3 T_input;
+ T_input.setRandom();
+ T_input = 100 * T_input;
+
+ // Create the camera system.
+ Mat3 R_expected;
+ Vec3 T_expected;
+ Mat3X X_world;
+ CreateCameraSystem(KK, x_image, X_distances, R_input, T_input,
+ &x_camera, &X_world, &R_expected, &T_expected);
+
+ // Test each of the resection methods.
+ {
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_ANSAR_DANIILIDIS);
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+ }
+ {
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_camera, X_world,
+ &R_output, &T_output,
+ RESECTION_EPNP);
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+ }
+ {
+ Mat3 R_output;
+ Vec3 T_output;
+ EuclideanResection(x_image, X_world, KK,
+ &R_output, &T_output);
+ EXPECT_MATRIX_NEAR(T_output, T_expected, 1e-5);
+ EXPECT_MATRIX_NEAR(R_output, R_expected, 1e-7);
+ }
+}
+#endif
diff --git a/extern/libmv/libmv/tracking/brute_region_tracker_test.cc b/extern/libmv/libmv/tracking/brute_region_tracker_test.cc
new file mode 100644
index 00000000000..9014797c7cf
--- /dev/null
+++ b/extern/libmv/libmv/tracking/brute_region_tracker_test.cc
@@ -0,0 +1,51 @@
+// 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/brute_region_tracker.h"
+#include "libmv/image/image.h"
+#include "libmv/logging/logging.h"
+#include "testing/testing.h"
+
+namespace libmv {
+namespace {
+
+TEST(KltRegionTracker, Track) {
+ Array3Df image1(51, 51);
+ image1.Fill(0);
+
+ Array3Df image2(image1);
+
+ int x0 = 25, y0 = 25;
+ int dx = 3, dy = 2;
+ image1(y0, x0) = 1.0f;
+ image2(y0 + dy, x0 + dx) = 1.0;
+
+ double x1 = x0;
+ double y1 = y0;
+
+ BruteRegionTracker tracker;
+ EXPECT_TRUE(tracker.Track(image1, image2, x0, y0, &x1, &y1));
+
+ EXPECT_NEAR(x1, x0 + dx, 0.001);
+ EXPECT_NEAR(y1, y0 + dy, 0.001);
+}
+
+} // namespace
+} // namespace libmv
diff --git a/extern/libmv/libmv/tracking/klt_region_tracker_test.cc b/extern/libmv/libmv/tracking/klt_region_tracker_test.cc
new file mode 100644
index 00000000000..07d5d6500e3
--- /dev/null
+++ b/extern/libmv/libmv/tracking/klt_region_tracker_test.cc
@@ -0,0 +1,51 @@
+// 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/klt_region_tracker.h"
+#include "libmv/image/image.h"
+#include "testing/testing.h"
+
+namespace libmv {
+namespace {
+
+TEST(KltRegionTracker, Track) {
+ Array3Df image1(51, 51);
+ image1.Fill(0);
+
+ Array3Df image2(image1);
+
+ int x0 = 25, y0 = 25;
+ int dx = 3, dy = 2;
+ image1(y0, x0) = 1.0f;
+ image2(y0 + dy, x0 + dx) = 1.0f;
+
+ double x1 = x0;
+ double y1 = y0;
+
+ KltRegionTracker tracker;
+ tracker.half_window_size = 6;
+ EXPECT_TRUE(tracker.Track(image1, image2, x0, y0, &x1, &y1));
+
+ EXPECT_NEAR(x1, x0 + dx, 0.001);
+ EXPECT_NEAR(y1, y0 + dy, 0.001);
+}
+
+} // namespace
+} // namespace libmv
diff --git a/extern/libmv/libmv/tracking/pyramid_region_tracker_test.cc b/extern/libmv/libmv/tracking/pyramid_region_tracker_test.cc
new file mode 100644
index 00000000000..d90a1012237
--- /dev/null
+++ b/extern/libmv/libmv/tracking/pyramid_region_tracker_test.cc
@@ -0,0 +1,80 @@
+// 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/pyramid_region_tracker.h"
+#include "libmv/tracking/klt_region_tracker.h"
+#include "libmv/image/image.h"
+#include "testing/testing.h"
+
+namespace libmv {
+namespace {
+
+TEST(PyramidKltRegionTracker, Track) {
+ Array3Df image1(100, 100);
+ image1.Fill(0);
+
+ Array3Df image2(image1);
+
+ int x1 = 25, y1 = 25;
+ image1(y1 + 0, x1 + 0) = 1.0f;
+ image1(y1 + 0, x1 + 1) = 1.0f;
+ image1(y1 + 1, x1 + 0) = 1.0f;
+ image1(y1 + 1, x1 + 1) = 1.0f;
+
+ // Make the displacement too large for a single-level KLT.
+ int x2 = x1 + 6, y2 = y1 + 5;
+ image2(y2 + 0, x2 + 0) = 1.0f;
+ image2(y2 + 0, x2 + 1) = 1.0f;
+ image2(y2 + 1, x2 + 0) = 1.0f;
+ image2(y2 + 1, x2 + 1) = 1.0f;
+
+ // Use a small 5x5 tracking region.
+ int half_window_size = 3;
+
+ // Ensure that the track doesn't work with one level of KLT.
+ {
+ double x2_actual = x1;
+ double y2_actual = y1;
+
+ KltRegionTracker tracker;
+ tracker.half_window_size = half_window_size;
+ EXPECT_FALSE(tracker.Track(image1, image2, x1, y1,
+ &x2_actual, &y2_actual));
+ }
+
+ // Verify that it works with the pyramid tracker.
+ {
+ double x2_actual = x1;
+ double y2_actual = y1;
+
+ KltRegionTracker *klt_tracker = new KltRegionTracker;
+ klt_tracker->half_window_size = half_window_size;
+
+ PyramidRegionTracker tracker(klt_tracker, 3);
+ EXPECT_TRUE(tracker.Track(image1, image2, x1, y1,
+ &x2_actual, &y2_actual));
+
+ EXPECT_NEAR(x2_actual, x2, 0.001);
+ EXPECT_NEAR(y2_actual, y2, 0.001);
+ }
+}
+
+} // namespace
+} // namespace libmv