From 898ba93a12117de413a9a5cd7b3fdf8c0d39ee7b Mon Sep 17 00:00:00 2001 From: Sergey Sharybin Date: Fri, 5 Apr 2013 09:23:20 +0000 Subject: Update libmv to current upstream version - Solves some strict compilation warning - Style/code cleanup --- extern/libmv/CMakeLists.txt | 2 +- extern/libmv/ChangeLog | 431 ++++++++------------- extern/libmv/files.txt | 2 +- extern/libmv/libmv/base/vector.h | 2 +- extern/libmv/libmv/base/vector_utils.h | 2 +- extern/libmv/libmv/image/array_nd.cc | 10 +- extern/libmv/libmv/image/array_nd.h | 44 ++- extern/libmv/libmv/image/convolve.cc | 60 +-- extern/libmv/libmv/image/convolve.h | 6 +- extern/libmv/libmv/image/correlation.h | 7 +- extern/libmv/libmv/image/image.h | 17 +- extern/libmv/libmv/image/sample.h | 8 +- extern/libmv/libmv/image/tuple.h | 8 +- extern/libmv/libmv/multiview/conditioning.cc | 6 +- extern/libmv/libmv/multiview/conditioning.h | 5 +- .../libmv/libmv/multiview/euclidean_resection.cc | 126 +++--- extern/libmv/libmv/multiview/euclidean_resection.h | 14 +- extern/libmv/libmv/multiview/fundamental.cc | 4 +- extern/libmv/libmv/multiview/homography.cc | 189 ++++----- extern/libmv/libmv/multiview/homography.h | 12 +- .../libmv/multiview/homography_parameterization.h | 12 +- extern/libmv/libmv/multiview/nviewtriangulation.h | 8 +- extern/libmv/libmv/multiview/panography.cc | 125 ++++++ extern/libmv/libmv/multiview/panography.h | 94 +---- extern/libmv/libmv/multiview/projection.cc | 59 +-- extern/libmv/libmv/multiview/projection.h | 24 +- extern/libmv/libmv/multiview/resection.h | 16 +- extern/libmv/libmv/multiview/triangulation.cc | 17 +- extern/libmv/libmv/multiview/triangulation.h | 10 +- extern/libmv/libmv/numeric/dogleg.h | 8 +- extern/libmv/libmv/numeric/function_derivative.h | 4 +- extern/libmv/libmv/numeric/levenberg_marquardt.h | 14 +- extern/libmv/libmv/numeric/numeric.cc | 6 +- extern/libmv/libmv/numeric/numeric.h | 69 ++-- extern/libmv/libmv/numeric/poly.cc | 6 +- extern/libmv/libmv/numeric/poly.h | 28 +- extern/libmv/libmv/simple_pipeline/bundle.cc | 15 +- extern/libmv/libmv/simple_pipeline/bundle.h | 13 +- extern/libmv/libmv/simple_pipeline/callbacks.cc | 29 -- extern/libmv/libmv/simple_pipeline/callbacks.h | 3 +- .../libmv/simple_pipeline/camera_intrinsics.cc | 151 +++++--- .../libmv/simple_pipeline/camera_intrinsics.h | 5 +- extern/libmv/libmv/simple_pipeline/detect.cc | 121 +++--- extern/libmv/libmv/simple_pipeline/detect.h | 6 +- .../simple_pipeline/initialize_reconstruction.cc | 14 +- .../simple_pipeline/initialize_reconstruction.h | 6 +- extern/libmv/libmv/simple_pipeline/intersect.cc | 7 +- extern/libmv/libmv/simple_pipeline/intersect.h | 6 +- extern/libmv/libmv/simple_pipeline/modal_solver.cc | 12 +- extern/libmv/libmv/simple_pipeline/modal_solver.h | 2 +- extern/libmv/libmv/simple_pipeline/pipeline.cc | 35 +- extern/libmv/libmv/simple_pipeline/pipeline.h | 9 +- .../libmv/libmv/simple_pipeline/reconstruction.h | 16 +- extern/libmv/libmv/simple_pipeline/resect.cc | 23 +- extern/libmv/libmv/simple_pipeline/resect.h | 6 +- extern/libmv/libmv/simple_pipeline/tracks.cc | 8 +- .../libmv/libmv/tracking/brute_region_tracker.cc | 19 +- extern/libmv/libmv/tracking/brute_region_tracker.h | 2 +- extern/libmv/libmv/tracking/esm_region_tracker.cc | 18 +- extern/libmv/libmv/tracking/esm_region_tracker.h | 2 +- .../libmv/libmv/tracking/hybrid_region_tracker.h | 2 +- extern/libmv/libmv/tracking/klt_region_tracker.cc | 6 +- extern/libmv/libmv/tracking/klt_region_tracker.h | 2 +- .../libmv/libmv/tracking/lmicklt_region_tracker.cc | 7 +- .../libmv/libmv/tracking/lmicklt_region_tracker.h | 2 +- .../libmv/libmv/tracking/pyramid_region_tracker.cc | 3 +- .../libmv/libmv/tracking/retrack_region_tracker.cc | 4 +- extern/libmv/libmv/tracking/track_region.cc | 70 ++-- .../libmv/libmv/tracking/trklt_region_tracker.cc | 3 +- extern/libmv/third_party/fast/fast_10.c | 3 +- extern/libmv/third_party/fast/fast_11.c | 3 +- extern/libmv/third_party/fast/fast_12.c | 3 +- extern/libmv/third_party/fast/fast_9.c | 3 +- .../libmv/third_party/gflags/gflags_completions.cc | 1 + extern/libmv/third_party/glog/src/logging.cc | 4 + extern/libmv/third_party/glog/src/utilities.cc | 2 +- extern/libmv/third_party/glog/src/vlog_is_on.cc | 6 + 77 files changed, 1067 insertions(+), 1040 deletions(-) create mode 100644 extern/libmv/libmv/multiview/panography.cc delete mode 100644 extern/libmv/libmv/simple_pipeline/callbacks.cc (limited to 'extern') diff --git a/extern/libmv/CMakeLists.txt b/extern/libmv/CMakeLists.txt index a4910d94a37..025481f9954 100644 --- a/extern/libmv/CMakeLists.txt +++ b/extern/libmv/CMakeLists.txt @@ -45,12 +45,12 @@ set(SRC libmv/multiview/euclidean_resection.cc libmv/multiview/fundamental.cc libmv/multiview/homography.cc + libmv/multiview/panography.cc libmv/multiview/projection.cc libmv/multiview/triangulation.cc libmv/numeric/numeric.cc libmv/numeric/poly.cc libmv/simple_pipeline/bundle.cc - libmv/simple_pipeline/callbacks.cc libmv/simple_pipeline/camera_intrinsics.cc libmv/simple_pipeline/detect.cc libmv/simple_pipeline/initialize_reconstruction.cc diff --git a/extern/libmv/ChangeLog b/extern/libmv/ChangeLog index a4f25f83a1c..a3d39d46f65 100644 --- a/extern/libmv/ChangeLog +++ b/extern/libmv/ChangeLog @@ -1,3 +1,169 @@ +commit e3b2bccba6145290738a6677c14f7369ec7a38cd +Author: Sergey Sharybin +Date: Thu Apr 4 02:59:58 2013 +0600 + + Suppress strict compiler warnings in glags/glog libraries + +commit 5fca459adcf0a3419fa9cd8d983dc2c952d02647 +Author: Sergey Sharybin +Date: Thu Apr 4 01:20:18 2013 +0600 + + Lint cleanup, mostly white space and line width. + + Also moved own includes to the top of files. + + Should be no functional changes :) + +commit 9a9dd458a622928b91dbd3c79900577923283838 +Author: Sergey Sharybin +Date: Fri Mar 29 00:20:29 2013 +0600 + + Fix for TransformTracks in uncalibrated pipeline + + Transformation matrix was completely ignored by + TransformTracks() and final marker coordinate + exactly matched it's source coordinates. + + Seems to be just a typo in vector usage: need to + use "b" (which is transformed one) instead of "a" + when converting projective coordinate to 2D space. + +commit d35766cc9901609e32f4d80faba715695bea3c40 +Author: Sergey Sharybin +Date: Fri Mar 29 00:19:11 2013 +0600 + + Revert part of e2eb58c4230f94ef0c72fb4005e5434088d52e80 + + That commit included one change which shall have been + go as separate commit with more detailed description. + +commit e8d71b4e96fd78eb60773b6557d66da672e65753 +Author: Sergey Sharybin +Date: Wed Mar 27 20:37:05 2013 +0600 + + Silenced more warnings + + - Added includes of own header to fast implementation files. + + - Camera intrinsics wouldn't complain about unknown pragma when + building without OpenMP support. + + TODO: Make it a CMake option to build libmv with OpenMP support. + Currently multi-threaded intrinsics only available when + using custom CMake rules for bundled libmv version + (as it's done in Blender). + +commit ad442812654f270dc088394410fda1b81b8dc450 +Author: Sergey Sharybin +Date: Wed Mar 27 20:18:51 2013 +0600 + + Multithreaded camera intrinsics + + Implemented multithreaded buffer (un)distortion + for camera intrinsics using OpenMP. + + By default, (un)distortion is single-threaded, + but it is possible to as CameraIntrinsics to + use more threads by calling SetThreads method. + +commit c88b4881096174a16a9f9e6fc2c9dcad3e255b25 +Author: Sergey Sharybin +Date: Wed Mar 27 18:45:09 2013 +0600 + + Movie functions implementation from panography header + into own CC implementation file. + + Before this all panography functions were declared as + static, which is not so much useful from re-useability + point of view. + +commit 2d2faf9104bc035722cff6775e1b8e7c93143aba +Author: Sergey Sharybin +Date: Wed Mar 27 18:37:36 2013 +0600 + + Build shared Ceres library only if BUILD_SHARED_LIBS is enabled + +commit daa3ddd3260ccaf2bf9c72eadb89213d91e549ec +Author: Sergey Sharybin +Date: Wed Mar 27 18:21:52 2013 +0600 + + Update Ceres to upstream version 1.5.0 + +commit cf5dc678878345ea3f221ce50cb2b9e539c2ab38 +Author: Sergey Sharybin +Date: Wed Mar 27 15:06:24 2013 +0600 + + Code cleanup: removed more deprecated FFmpeg API usage + + This time in qt-tracker application. + +commit e2eb58c4230f94ef0c72fb4005e5434088d52e80 +Author: Sergey Sharybin +Date: Tue Mar 26 17:19:51 2013 +0600 + + Code cleanup: silent unused variables warnings + +commit af89bb24667e39b7e655173ea807fdcfbeef4422 +Author: Sergey Sharybin +Date: Tue Mar 26 16:54:14 2013 +0600 + + Code cleanup: no need to declare empty body for ProgressUpdateCallback:invoke + + Make force this method to be overridden by derivative classes. + Also removed currently unneeded callbacks.cc. + +commit 0441d4ee06fad0219256a5704f931eec916a3868 +Author: Sergey Sharybin +Date: Tue Mar 26 16:37:27 2013 +0600 + + Code cleanup: silent type narrowing in qt-tracker + +commit cd4b61c976448d0fdedefb3ed4b21d70e078f94b +Author: Sergey Sharybin +Date: Tue Mar 26 16:26:39 2013 +0600 + + Changes to unit testing + + - Move ceres test binaries to ${LIBMV_TESTS_OUTPUT_DIR}/ceres, + so they don't mess with libmv's application binaries and + tests. + + - Removed ceres_ prefix from ceres unit tests, only use this + prefix for targets (targets need to be unique name). + + - Added unit tests data for ceres, otherwise system_test fails. + + - Restored "test" makefile target. + +commit cf704ada08acc8b26167e7bfb3e1e88fd278de23 +Author: Sergey Sharybin +Date: Tue Mar 26 15:01:15 2013 +0600 + + Code cleanup: use rw-rw-r-- mode for source files + +commit 64b31e3e43acb52aaf6f591b9d1c2449bf6ef3bd +Author: Sergey Sharybin +Date: Tue Mar 26 14:57:46 2013 +0600 + + Code cleanup: don't use deprecated FFmpeg API functions + +commit 2a3676499548ad5dba5a5c5eadf3bb71e640b612 +Author: Sergey Sharybin +Date: Tue Mar 5 17:40:52 2013 +0600 + + Switch from DENSE_NORMAL_CHOLESKY to DENSE_QR + + DENSE_QR is better behaved numerically and after recent + changes from Sameer there's no big difference in speed. + +commit bcb920df02133da5b7e55fbc74edb9222004eecc +Author: Sergey Sharybin +Date: Tue Mar 5 17:15:43 2013 +0600 + + Update Ceres to 1.5RC3 + + It brings optimization of DENSE_QR and DENSE_SCHUR solvers. + commit 473996468a4e67e7c860169181a4ff31ce9b8c80 Author: Sergey Sharybin Date: Fri Mar 1 17:44:54 2013 +0600 @@ -533,268 +699,3 @@ Date: Fri Jun 8 17:42:17 2012 +0000 aborts immediately. The workaround for now is to disable the correlation checking, and examine your tracks carefully. A fix will get added shortly. - -commit 81d028f13738ebe2304287dfce90e91bc782e2cf -Author: Keir Mierle -Date: Fri May 18 20:04:43 2012 +0000 - - Remove an unnecessary template<> line in libmv. Convert debug logs to LG. - -commit 238aaba241ef99995d254aadc974db719da04b96 -Author: Keir Mierle -Date: Fri May 18 12:05:10 2012 +0000 - - Support normalization in the tracking prepass - - The last tracker commit added normalized tracking. This makes - tracking patches undergoing uniform illumination change easier. - However, the prepass which computes a quick translation-only - estimate of the warp did not take this into account. This commit - fixes that. - - This works reasonably well but in some examples the brute - initialization fails. I suspect this is due to the warped template - estimate in the current frame being too different from the - original, so there are multiple peaks in the normalized-SAD - correlation function. - - The solution is to use the previous frame for the brute - initialization and the keyframe for refinement, but that requires - architecture changes. - -commit 981ca4f6a679cd9ac3d086eae3cd946ce72ca8a5 -Author: Keir Mierle -Date: Fri May 18 02:12:47 2012 +0000 - - Add light-normalized tracking to the planar tracker - - This commit adds the ability to normalize patterns by their - average value while tracking, to make them invariant to global - illumination changes. - - To see this in action, check out the "Lobby" scene from Hollywood - VFX. If you track the markers that are shadowed by the actress, - previously they would not track. With the scale adaption on, the - tracker would shrink the area to compensate for the changed - illumination, losing the track. With "Normalize" turned on, the - patch is correctly tracked and scale is maintained. - - A remaining problem is that only the Ceres cost function is - updated to handle the normalization. The brute translation search - does not take this into account. Perhaps "Prepass" (see below) - should get disabled if normalization is enabled until I fix the - prepass to normalize as well. - - There are a few other changes: - - - Bail out of the sampling loop early if the mask is zero; this - saves expensive samples of the image derivatives. - - - Fix a bug where the mask was ignored when sampling in the cost - functor. - -commit e9384b15fb2a6a5b81346d5758fa136f0911e945 -Author: Keir Mierle -Date: Thu May 17 23:53:32 2012 +0000 - - Implement support for affine tracking in the planar tracker; cleanups. - -commit 021d41eed8b4ce6a4e37786ccd357ed5dc83a13f -Author: Keir Mierle -Date: Thu May 17 21:26:06 2012 +0000 - - For the planar tracker, initialize the warp from the four correspondences - after brute force translation search. - -commit 003d1bf6145cfd30938b35f6e10d43708dbf916c -Author: Sergey Sharybin -Date: Thu Dec 6 16:56:01 2012 +0600 - - Correction to region tracker options initialization. - - Based on patch from Keir to Blender: - https://svn.blender.org/svnroot/bf-blender/branches/soc-2011-tomato@46743 - -commit 6af47b218cfdf5219f0ebb3cb95459817cf9abf2 -Author: Keir Mierle -Date: Thu May 17 02:31:52 2012 +0000 - - Add new planar tracker features and use the new planar API - - This commit removes the use of the legacy RegionTracker API from - Blender, and replaces it with the new TrackRegion API. This also - adds several features to the planar tracker in libmv: - - - Do a brute-force initialization of tracking similar to "Hybrid" - mode in the stable release, but using all floats. This is slower - but more accurate. It is still necessary to evaluate if the - performance loss is worth it. In particular, this change is - necessary to support high bit depth imagery. - - - Add support for masks over the search window. This is a step - towards supporting user-defined tracker masks. The tracker masks - will make it easy for users to make a mask for e.g. a ball. - - - Add Pearson product moment correlation coefficient checking (aka - "Correlation" in the UI. This causes tracking failure if the - tracked patch is not linearly related to the template. - - - Add support for warping a few points in addition to the supplied - points. This is useful because the tracking code deliberately - does not expose the underlying warp representation. Instead, - warps are specified in an aparametric way via the correspondences. - - - Remove the "num_samples_xy" concept and replace it with - automatic determination of the number of samples. This makes the - API easier for users. - - - Fix various bugs in the parameterizations. - - There remains a bug with subpixel precision tracking when in - "keyframe" mode; this will get fixed shortly. - -commit 16a46db104468cec80bd31ca9d5f8bffbe3e003e -Author: Keir Mierle -Date: Mon May 14 12:15:38 2012 +0000 - - "Efficient Second-order Minimization" for the planar tracker - - This implements the "Efficient Second-order Minimization" - scheme, as supported by the existing translation tracker. - This increases the amount of per-iteration work, but - decreases the number of iterations required to converge and - also increases the size of the basin of attraction for the - optimization. - -commit 23243b1b1f3e1ab3ef862b47bca06ee876ac2cf4 -Author: Keir Mierle -Date: Sun May 13 23:08:56 2012 +0000 - - Add a planar tracking implementation to libmv - - This adds a new planar tracking implementation to libmv. The - tracker is based on Ceres[1], the new nonlinear minimizer that - myself and Sameer released from Google as open source. Since - the motion model is more involved, the interface is - different than the RegionTracker interface used previously - in Blender. - - The ESM tracker, also known as the KLT tracker in the UI, is - temporarily changed to use the new Ceres-based planar - tracker in translation-only mode. Currently it is a bit - slower than ESM and also doesn't have all the bells and - whistles implemented. Those will come soon. Longer term, - both trackers will remain since Ceres is unlikely to be as - fast as ESM for pure translation solving, due to its - generality. - - The next step is to implement a new tracking UI. The current - UI assumes a translational motion model; the new one must - support arbitrary perspective transforms of the pattern - regions. - - [1] http://code.google.com/p/ceres-solver - -commit 52be92b53eb4decb1a316690b162196f227cc441 -Author: Sergey Sharybin -Date: Thu Dec 6 16:06:08 2012 +0600 - - Initial Ceres integration - - Currently only put sources to src/third_party/ceres and made sure they're - not giving compilation issues. - - Used Ceres upstream version 1.3.0. - - Needed to make some modifications to it's CMakeLists.txt also to glog and - fglags. They're described in README.libmv of this libraries. - - Basically: - - - Added -fPIC to glog/gflags, so shared ceres library could be linked - statically against this libraries. - - - Tweaked Ceres's build rules to use needed libraries from libmv's - third_party folder. - -commit b13f9d13122e091cb85855c2094386ccdef6e5a4 -Author: Sergey Sharybin -Date: Wed Dec 5 19:05:34 2012 +0600 - - Update Eigen to version 3.1.2 - - Mainly because of lots of warnings generating by gcc-4.7 which are - resolved in newer eigen version. - -commit 1f0dd94e8e37d3fe2df89282ec16a6a685fdde0b -Author: Sergey Sharybin -Date: Fri May 25 16:36:44 2012 +0600 - - - Added avutil to qt-tracker linking when building with FFmpeg support. - On some platforms it seems to be required - - Synchronized QT Creator project for qt-tracker with changes in sources, - so no it might be compiled from QT Creator. - -commit b813dbe3f46bbbc7e73ac791d4665622e4fc7ba5 -Author: Sergey Sharybin -Date: Wed May 9 19:01:10 2012 +0600 - - Modal solver: Detect rigid transformation between initial frame and current - instead of detecting it between two neighbour frames. - - This prevents accumulation of error and seems to be working better in footages i've tested. - -commit 9254621c76daaf239ec1f535e197ca792eea97b6 -Author: Sergey Sharybin -Date: Wed May 9 18:57:00 2012 +0600 - - Backport changes made by Keir in Blender: - - - Enhance logging in libmv's trackers. - - Cleanups in brute_region_tracker.cc. - -commit d9c56b9d3c63f886d83129ca0ebed1e76d9c93d7 -Author: Sergey Sharybin -Date: Fri Apr 27 16:20:41 2012 +0600 - - Fixes for MinGW64 support by Caleb Joseph with slight modifications by Antony Riakiotakis - - - Functions snprintf and sincos shouldn't be redefined for MinGW64 - - Type pid_t shouldn't be re-defined for MinGW64 - -commit e1902b6938676011607ac99986b8b140bdbf090e -Author: Sergey Sharybin -Date: Fri Apr 27 16:04:19 2012 +0600 - - Fixes for Qt calibration tool - - - Passing directory with images via command line argument now isn't - required -- it there's no such directory specified standard open - dialog might be used for this (before application used to abort - due to accessing to non-existing list element). - - Conversion of source images to grayscale now happens correct. - It was needed to build grayscale palette for 8bit indexed buffer. - -commit 05f1a0a78ad8ff6646d1e8da97e6f7575b891536 -Author: Sergey Sharybin -Date: Sat Apr 14 17:21:29 2012 +0600 - - Make QtTracker compilable again porting it to recent API change and code cleanup: - - - It was using SAD tracker with own API, now it's using standard RegionTracker API - which should make it easier to switch between different trackers. - - Restored LaplaceFilter from old SAD module which convolves images with the - discrete laplacian operator. - -commit a44312a7beb2963b8e3bf8015c516d2eff40cc3d -Author: Sergey Sharybin -Date: Thu Apr 12 13:56:02 2012 +0600 - - Added solver for modal camera motion, currently supports only tripod solving - - This solver is intended to deal with such camera motions as tripod and panning, - where it's impossible to reconstruct exact position of markers in 3d view. - - It projects markers onto sphere and uses rigid registration of rotation to - find rotation angles which makes bundles from previous and current frame be - as closest as it's possible. diff --git a/extern/libmv/files.txt b/extern/libmv/files.txt index 16afdb36371..d77d0ea8147 100644 --- a/extern/libmv/files.txt +++ b/extern/libmv/files.txt @@ -21,6 +21,7 @@ libmv/multiview/homography.cc libmv/multiview/homography.h libmv/multiview/homography_parameterization.h libmv/multiview/nviewtriangulation.h +libmv/multiview/panography.cc libmv/multiview/panography.h libmv/multiview/projection.cc libmv/multiview/projection.h @@ -36,7 +37,6 @@ libmv/numeric/poly.cc libmv/numeric/poly.h libmv/simple_pipeline/bundle.cc libmv/simple_pipeline/bundle.h -libmv/simple_pipeline/callbacks.cc libmv/simple_pipeline/callbacks.h libmv/simple_pipeline/camera_intrinsics.cc libmv/simple_pipeline/camera_intrinsics.h diff --git a/extern/libmv/libmv/base/vector.h b/extern/libmv/libmv/base/vector.h index c2a3298087a..067db4ba7f7 100644 --- a/extern/libmv/libmv/base/vector.h +++ b/extern/libmv/libmv/base/vector.h @@ -39,7 +39,7 @@ namespace libmv { // - doesn't support iterators. // - impede compatibility with code using STL. // - the STL already provide support for custom allocators -// it could be replaced with a simple +// it could be replaced with a simple // template class vector : std::vector {} declaration // provided it doesn't break code relying on libmv::vector specific behavior template clear(); } -#endif // LIBMV_BASE_VECTOR_UTILS_H_ +#endif // LIBMV_BASE_VECTOR_UTILS_H_ diff --git a/extern/libmv/libmv/image/array_nd.cc b/extern/libmv/libmv/image/array_nd.cc index 3a77e3e4881..469a19aabf1 100644 --- a/extern/libmv/libmv/image/array_nd.cc +++ b/extern/libmv/libmv/image/array_nd.cc @@ -35,8 +35,8 @@ void FloatArrayToScaledByteArray(const Array3Df &float_array, for (int i = 0; i < float_array.Height(); ++i) { for (int j = 0; j < float_array.Width(); ++j) { for (int k = 0; k < float_array.Depth(); ++k) { - minval = std::min(minval, float_array(i,j,k)); - maxval = std::max(maxval, float_array(i,j,k)); + minval = std::min(minval, float_array(i, j, k)); + maxval = std::max(maxval, float_array(i, j, k)); } } } @@ -47,8 +47,8 @@ void FloatArrayToScaledByteArray(const Array3Df &float_array, for (int i = 0; i < float_array.Height(); ++i) { for (int j = 0; j < float_array.Width(); ++j) { for (int k = 0; k < float_array.Depth(); ++k) { - float unscaled = (float_array(i,j,k) - minval) / (maxval - minval); - (*byte_array)(i,j,k) = (unsigned char)(255 * unscaled); + float unscaled = (float_array(i, j, k) - minval) / (maxval - minval); + (*byte_array)(i, j, k) = (unsigned char)(255 * unscaled); } } } @@ -60,7 +60,7 @@ void ByteArrayToScaledFloatArray(const Array3Du &byte_array, for (int i = 0; i < byte_array.Height(); ++i) { for (int j = 0; j < byte_array.Width(); ++j) { for (int k = 0; k < byte_array.Depth(); ++k) { - (*float_array)(i,j,k) = float(byte_array(i,j,k)) / 255.0f; + (*float_array)(i, j, k) = float(byte_array(i, j, k)) / 255.0f; } } } diff --git a/extern/libmv/libmv/image/array_nd.h b/extern/libmv/libmv/image/array_nd.h index 6d7570cda9b..f40fadcbde7 100644 --- a/extern/libmv/libmv/image/array_nd.h +++ b/extern/libmv/libmv/image/array_nd.h @@ -57,9 +57,13 @@ class ArrayND : public BaseArray { ArrayND(int s0) : data_(NULL), own_data(true) { Resize(s0); } ArrayND(int s0, int s1) : data_(NULL), own_data(true) { Resize(s0, s1); } - ArrayND(int s0, int s1, int s2) : data_(NULL), own_data(true) { Resize(s0, s1, s2); } + ArrayND(int s0, int s1, int s2) : data_(NULL), own_data(true) { + Resize(s0, s1, s2); + } - ArrayND(T* data, int s0, int s1, int s2) : data_(data), own_data(false) { Resize(s0, s1, s2); } + ArrayND(T* data, int s0, int s1, int s2) : data_(data), own_data(false) { + Resize(s0, s1, s2); + } /// Destructor deletes pixel data. ~ArrayND() { @@ -93,7 +97,7 @@ class ArrayND : public BaseArray { for (int i = N - 1; i > 0; --i) { strides_(i - 1) = strides_(i) * shape_(i); } - if(own_data) { + if (own_data) { delete [] data_; data_ = NULL; if (Size() > 0) { @@ -103,7 +107,7 @@ class ArrayND : public BaseArray { } template - void ResizeLike(const ArrayND &other) { + void ResizeLike(const ArrayND &other) { Resize(other.Shape()); } @@ -136,12 +140,12 @@ class ArrayND : public BaseArray { /// Resize a 3D array to shape (s0,s1,s2). void Resize(int s0, int s1, int s2) { assert(N == 3); - int shape[] = {s0,s1,s2}; + int shape[] = {s0, s1, s2}; Resize(shape); } template - void CopyFrom(const ArrayND &other) { + void CopyFrom(const ArrayND &other) { ResizeLike(other); T *data = Data(); const D *other_data = other.Data(); @@ -238,7 +242,7 @@ class ArrayND : public BaseArray { T &operator()(int i0, int i1) { assert(0 <= i0 && i0 < Shape(0)); assert(0 <= i1 && i1 < Shape(1)); - return *( Data() + Offset(i0,i1) ); + return *(Data() + Offset(i0, i1)); } /// 3D specialization. @@ -246,29 +250,29 @@ class ArrayND : public BaseArray { assert(0 <= i0 && i0 < Shape(0)); assert(0 <= i1 && i1 < Shape(1)); assert(0 <= i2 && i2 < Shape(2)); - return *( Data() + Offset(i0,i1,i2) ); + return *(Data() + Offset(i0, i1, i2)); } /// Return a constant reference to the element at position index. const T &operator()(const Index &index) const { - return *( Data() + Offset(index) ); + return *(Data() + Offset(index)); } /// 1D specialization. const T &operator()(int i0) const { - return *( Data() + Offset(i0) ); + return *(Data() + Offset(i0)); } /// 2D specialization. const T &operator()(int i0, int i1) const { assert(0 <= i0 && i0 < Shape(0)); assert(0 <= i1 && i1 < Shape(1)); - return *( Data() + Offset(i0,i1) ); + return *(Data() + Offset(i0, i1)); } /// 3D specialization. const T &operator()(int i0, int i1, int i2) const { - return *( Data() + Offset(i0,i1,i2) ); + return *(Data() + Offset(i0, i1, i2)); } /// True if index is inside array. @@ -343,14 +347,14 @@ class Array3D : public ArrayND { Array3D() : Base() { } - Array3D(int height, int width, int depth=1) + Array3D(int height, int width, int depth = 1) : Base(height, width, depth) { } - Array3D(T* data, int height, int width, int depth=1) + Array3D(T* data, int height, int width, int depth = 1) : Base(data, height, width, depth) { } - void Resize(int height, int width, int depth=1) { + void Resize(int height, int width, int depth = 1) { Base::Resize(height, width, depth); } @@ -376,12 +380,12 @@ class Array3D : public ArrayND { T &operator()(int i0, int i1, int i2 = 0) { assert(0 <= i0 && i0 < Height()); assert(0 <= i1 && i1 < Width()); - return Base::operator()(i0,i1,i2); + return Base::operator()(i0, i1, i2); } const T &operator()(int i0, int i1, int i2 = 0) const { assert(0 <= i0 && i0 < Height()); assert(0 <= i1 && i1 < Width()); - return Base::operator()(i0,i1,i2); + return Base::operator()(i0, i1, i2); } }; @@ -413,13 +417,13 @@ void ByteArrayToScaledFloatArray(const Array3Du &byte_array, Array3Df *float_array); template -void MultiplyElements( const AArrayType &a, +void MultiplyElements(const AArrayType &a, const BArrayType &b, - CArrayType *c ) { + CArrayType *c) { // This function does an element-wise multiply between // the two Arrays A and B, and stores the result in C. // A and B must have the same dimensions. - assert( a.Shape() == b.Shape() ); + assert(a.Shape() == b.Shape()); c->ResizeLike(a); // To perform the multiplcation, a "current" index into the N-dimensions of diff --git a/extern/libmv/libmv/image/convolve.cc b/extern/libmv/libmv/image/convolve.cc index 63ff7d6d8ff..4682c619f04 100644 --- a/extern/libmv/libmv/image/convolve.cc +++ b/extern/libmv/libmv/image/convolve.cc @@ -18,10 +18,11 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/image/convolve.h" + #include #include "libmv/image/image.h" -#include "libmv/image/convolve.h" namespace libmv { @@ -118,7 +119,7 @@ void Convolve(const Array3Df &in, // fast path. int half_width = kernel.size() / 2; switch (half_width) { -#define static_convolution( size ) case size: \ +#define static_convolution(size) case size: \ FastConvolve(kernel, width, height, src, src_stride, \ src_line_stride, dst, dst_stride); break; static_convolution(1) @@ -136,13 +137,15 @@ void Convolve(const Array3Df &in, double sum = 0; // Slow path: this loop cannot be unrolled. for (int k = -dynamic_size; k <= dynamic_size; ++k) { - if(vertical) { + if (vertical) { if (y + k >= 0 && y + k < height) { - sum += src[k * src_line_stride] * kernel(2 * dynamic_size - (k + dynamic_size)); + sum += src[k * src_line_stride] * + kernel(2 * dynamic_size - (k + dynamic_size)); } } else { if (x + k >= 0 && x + k < width) { - sum += src[k * src_stride] * kernel(2 * dynamic_size - (k + dynamic_size)); + sum += src[k * src_stride] * + kernel(2 * dynamic_size - (k + dynamic_size)); } } } @@ -234,25 +237,25 @@ void BoxFilterHorizontal(const Array3Df &in, int half_width = (window_size - 1) / 2; for (int k = 0; k < in.Depth(); ++k) { - for (int i=0; i 255) d=255; - dst[y*width+x] = d; - } +void LaplaceFilter(unsigned char* src, + unsigned char* dst, + int width, + int height, + int strength) { + for (int y = 1; y < height-1; y++) + for (int x = 1; x < width-1; x++) { + const unsigned char* s = &src[y*width+x]; + int l = 128 + + s[-width-1] + s[-width] + s[-width+1] + + s[1] - 8*s[0] + s[1] + + s[ width-1] + s[ width] + s[ width+1]; + int d = ((256-strength)*s[0] + strength*l) / 256; + if (d < 0) d=0; + if (d > 255) d=255; + dst[y*width+x] = d; + } } } // namespace libmv diff --git a/extern/libmv/libmv/image/convolve.h b/extern/libmv/libmv/image/convolve.h index a005dc31f10..d3b6da9794b 100644 --- a/extern/libmv/libmv/image/convolve.h +++ b/extern/libmv/libmv/image/convolve.h @@ -95,7 +95,11 @@ void BoxFilter(const FloatImage &in, \note Make sure the search region is filtered with the same strength as the pattern. */ -void LaplaceFilter(unsigned char* src, unsigned char* dst, int width, int height, int strength); +void LaplaceFilter(unsigned char* src, + unsigned char* dst, + int width, + int height, + int strength); } // namespace libmv diff --git a/extern/libmv/libmv/image/correlation.h b/extern/libmv/libmv/image/correlation.h index ba64951a167..ac1f8edab49 100644 --- a/extern/libmv/libmv/image/correlation.h +++ b/extern/libmv/libmv/image/correlation.h @@ -26,9 +26,10 @@ namespace libmv { -inline double PearsonProductMomentCorrelation(Array3Df image_and_gradient1_sampled, - Array3Df image_and_gradient2_sampled, - int width) { +inline double PearsonProductMomentCorrelation( + Array3Df image_and_gradient1_sampled, + Array3Df image_and_gradient2_sampled, + int width) { double sX = 0, sY = 0, sXX = 0, sYY = 0, sXY = 0; for (int r = 0; r < width; ++r) { diff --git a/extern/libmv/libmv/image/image.h b/extern/libmv/libmv/image/image.h index d158b0e0397..e0f200a4c93 100644 --- a/extern/libmv/libmv/image/image.h +++ b/extern/libmv/libmv/image/image.h @@ -60,8 +60,7 @@ class Image { // Size in bytes that the image takes in memory. int MemorySizeInBytes() { int size; - switch (array_type_) - { + switch (array_type_) { case BYTE: size = reinterpret_cast(array_)->MemorySizeInBytes(); break; @@ -83,8 +82,7 @@ class Image { } ~Image() { - switch (array_type_) - { + switch (array_type_) { case BYTE: delete reinterpret_cast(array_); @@ -109,23 +107,22 @@ class Image { Image& operator= (const Image& f) { if (this != &f) { array_type_ = f.array_type_; - switch (array_type_) - { + switch (array_type_) { case BYTE: delete reinterpret_cast(array_); - array_ = new Array3Du( *(Array3Du *)f.array_); + array_ = new Array3Du(*(Array3Du *)f.array_); break; case FLOAT: delete reinterpret_cast(array_); - array_ = new Array3Df( *(Array3Df *)f.array_); + array_ = new Array3Df(*(Array3Df *)f.array_); break; case INT: delete reinterpret_cast(array_); - array_ = new Array3Di( *(Array3Di *)f.array_); + array_ = new Array3Di(*(Array3Di *)f.array_); break; case SHORT: delete reinterpret_cast(array_); - array_ = new Array3Ds( *(Array3Ds *)f.array_); + array_ = new Array3Ds(*(Array3Ds *)f.array_); break; default: assert(0); diff --git a/extern/libmv/libmv/image/sample.h b/extern/libmv/libmv/image/sample.h index a8850effeab..a95d08a9815 100644 --- a/extern/libmv/libmv/image/sample.h +++ b/extern/libmv/libmv/image/sample.h @@ -71,8 +71,8 @@ inline T SampleLinear(const Array3D &image, float y, float x, int v = 0) { const T im21 = image(y2, x1, v); const T im22 = image(y2, x2, v); - return T( dy * ( dx * im11 + (1.0 - dx) * im12 ) + - (1 - dy) * ( dx * im21 + (1.0 - dx) * im22 )); + return T( dy * (dx * im11 + (1.0 - dx) * im12) + + (1 - dy) * (dx * im21 + (1.0 - dx) * im22)); } /// Linear interpolation, of all channels. The sample is assumed to have the @@ -95,8 +95,8 @@ inline void SampleLinear(const Array3D &image, float y, float x, T *sample) { const T im21 = image(y2, x1, i); const T im22 = image(y2, x2, i); - sample[i] = T( dy * ( dx * im11 + (1.0 - dx) * im12 ) + - (1 - dy) * ( dx * im21 + (1.0 - dx) * im22 )); + sample[i] = T( dy * (dx * im11 + (1.0 - dx) * im12) + + (1 - dy) * (dx * im21 + (1.0 - dx) * im22)); } } diff --git a/extern/libmv/libmv/image/tuple.h b/extern/libmv/libmv/image/tuple.h index 79acc9579d0..c8dc36f2e18 100644 --- a/extern/libmv/libmv/image/tuple.h +++ b/extern/libmv/libmv/image/tuple.h @@ -37,10 +37,10 @@ class Tuple { Tuple(D *values) { Reset(values); } template - Tuple(const Tuple &b) { Reset(b); } + Tuple(const Tuple &b) { Reset(b); } template - Tuple& operator=(const Tuple& b) { + Tuple& operator=(const Tuple& b) { Reset(b); return *this; } @@ -50,14 +50,14 @@ class Tuple { template void Reset(D *values) { - for(int i=0;iresize(2,n); + transformed_points->resize(2, n); Mat3X p(3, n); EuclideanToHomogeneous(points, &p); p = T * p; @@ -96,4 +96,4 @@ void UnnormalizerI::Unnormalize(const Mat3 &T1, const Mat3 &T2, Mat3 *H) { *H = T2.inverse() * (*H) * T1; } -} // namespace libmv +} // namespace libmv diff --git a/extern/libmv/libmv/multiview/conditioning.h b/extern/libmv/libmv/multiview/conditioning.h index 181d7485374..8f3e3a76070 100644 --- a/extern/libmv/libmv/multiview/conditioning.h +++ b/extern/libmv/libmv/multiview/conditioning.h @@ -54,7 +54,6 @@ struct UnnormalizerT { static void Unnormalize(const Mat3 &T1, const Mat3 &T2, Mat3 *H); }; -} //namespace libmv +} // namespace libmv - -#endif // LIBMV_MULTIVIEW_CONDITIONNING_H_ +#endif // LIBMV_MULTIVIEW_CONDITIONNING_H_ diff --git a/extern/libmv/libmv/multiview/euclidean_resection.cc b/extern/libmv/libmv/multiview/euclidean_resection.cc index 062c8b9067c..10c7330770f 100644 --- a/extern/libmv/libmv/multiview/euclidean_resection.cc +++ b/extern/libmv/libmv/multiview/euclidean_resection.cc @@ -18,6 +18,8 @@ // 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 #include @@ -26,7 +28,6 @@ #include "libmv/base/vector.h" #include "libmv/logging/logging.h" -#include "libmv/multiview/euclidean_resection.h" #include "libmv/multiview/projection.h" namespace libmv { @@ -34,7 +35,7 @@ namespace euclidean_resection { typedef unsigned int uint; -bool EuclideanResection(const Mat2X &x_camera, +bool EuclideanResection(const Mat2X &x_camera, const Mat3X &X_world, Mat3 *R, Vec3 *t, ResectionMethod method, @@ -52,7 +53,7 @@ bool EuclideanResection(const Mat2X &x_camera, return false; } -bool EuclideanResection(const Mat &x_image, +bool EuclideanResection(const Mat &x_image, const Mat3X &X_world, const Mat3 &K, Mat3 *R, Vec3 *t, @@ -80,11 +81,11 @@ void AbsoluteOrientation(const Mat3X &X, // Normalize the two point sets. Mat3X Xn(3, num_points), Xpn(3, num_points); - for( int i = 0; i < num_points; ++i ){ + for (int i = 0; i < num_points; ++i) { Xn.col(i) = X.col(i) - C; Xpn.col(i) = Xp.col(i) - Cp; } - + // Construct the N matrix (pg. 635). double Sxx = Xn.row(0).dot(Xpn.row(0)); double Syy = Xn.row(1).dot(Xpn.row(1)); @@ -101,7 +102,7 @@ void AbsoluteOrientation(const Mat3X &X, Syz - Szy, Sxx - Syy - Szz, Sxy + Syx, Szx + Sxz, Szx - Sxz, Sxy + Syx, -Sxx + Syy - Szz, Syz + Szy, Sxy - Syx, Szx + Sxz, Syz + Szy, -Sxx - Syy + Szz; - + // Find the unit quaternion q that maximizes qNq. It is the eigenvector // corresponding to the lagest eigenvalue. Vec4 q = N.jacobiSvd(Eigen::ComputeFullU).matrixU().col(0); @@ -180,12 +181,12 @@ static Vec MatrixToConstraint(const Mat &A, C.setZero(); int idx = 0; for (int i = 0; i < num_lambda; ++i) { - for( int j = i; j < num_lambda; ++j) { + for (int j = i; j < num_lambda; ++j) { C(idx) = A(i, j); - if (i != j){ + if (i != j) { C(idx) += A(j, i); } - ++ idx; + ++idx; } } return C; @@ -194,14 +195,14 @@ static Vec MatrixToConstraint(const Mat &A, // Normalizes the columns of vectors. static void NormalizeColumnVectors(Mat3X *vectors) { int num_columns = vectors->cols(); - for (int i = 0; i < num_columns; ++i){ + for (int i = 0; i < num_columns; ++i) { vectors->col(i).normalize(); } } -void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, - const Mat3X &X_world, - Mat3 *R, +void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, + const Mat3X &X_world, + Mat3 *R, Vec3 *t) { CHECK(x_camera.cols() == X_world.cols()); CHECK(x_camera.cols() > 3); @@ -214,7 +215,7 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, x_camera_unit.block(0, 0, 2, num_points) = x_camera; x_camera_unit.row(2).setOnes(); NormalizeColumnVectors(&x_camera_unit); - + int num_m_rows = num_points * (num_points - 1) / 2; int num_tt_variables = num_points * (num_points + 1) / 2; int num_m_columns = num_tt_variables + 1; @@ -224,7 +225,7 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, // Create the constraint equations for the t_ij variables (7) and arrange // them into the M matrix (8). Also store the initial (i, j) indices. - int row=0; + int row = 0; for (int i = 0; i < num_points; ++i) { for (int j = i+1; j < num_points; ++j) { M(row, row) = -2 * x_camera_unit.col(i).dot(x_camera_unit.col(j)); @@ -243,7 +244,7 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, } int num_lambda = num_points + 1; // Dimension of the null space of M. - Mat V = M.jacobiSvd(Eigen::ComputeFullV).matrixV().block(0, + Mat V = M.jacobiSvd(Eigen::ComputeFullV).matrixV().block(0, num_m_rows, num_m_columns, num_lambda); @@ -263,12 +264,11 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, int counter_k_row = 0; for (int idx1 = num_m_rows; idx1 < num_tt_variables; ++idx1) { for (int idx2 = 0; idx2 < num_m_rows; ++idx2) { - unsigned int i = ij_index(idx1, 0); unsigned int j = ij_index(idx2, 0); unsigned int k = ij_index(idx2, 1); - if( i != j && i != k ){ + if (i != j && i != k) { int idx3 = IJToPointIndex(i, j, num_points); int idx4 = IJToPointIndex(i, k, num_points); @@ -318,13 +318,12 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, L_sq = L_sq * Sign(L_sq(IJToIndex(max_L_sq_index, max_L_sq_index, num_lambda))); - - + Vec L(num_lambda); L(max_L_sq_index) = sqrt(L_sq(IJToIndex(max_L_sq_index, max_L_sq_index, num_lambda))); - + for (int i = 0; i < num_lambda; ++i) { if (i != max_L_sq_index) { L(i) = L_sq(IJToIndex(max_L_sq_index, i, num_lambda)) / L(max_L_sq_index); @@ -334,7 +333,7 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, // Correct the scale using the fact that the last constraint is equal to 1. L = L / (V.row(num_m_columns - 1).dot(L)); Vec X = V * L; - + // Recover the distances from the camera center to the 3D points Q. Vec d(num_points); d.setZero(); @@ -344,7 +343,7 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, // Create the 3D points in the camera system. Mat X_cam(3, num_points); - for (int c_point = 0; c_point < num_points; ++c_point ) { + for (int c_point = 0; c_point < num_points; ++c_point) { X_cam.col(c_point) = d(c_point) * x_camera_unit.col(c_point); } // Recover the camera translation and rotation. @@ -363,17 +362,17 @@ static void SelectControlPoints(const Mat3X &X_world, X_control_points->col(0) = mean; // Computes PCA - X_centered->resize (3, num_points); + X_centered->resize(3, num_points); for (size_t c = 0; c < num_points; c++) { - X_centered->col(c) = X_world.col (c) - mean; + X_centered->col(c) = X_world.col(c) - mean; } Mat3 X_centered_sq = (*X_centered) * X_centered->transpose(); Eigen::JacobiSVD X_centered_sq_svd(X_centered_sq, Eigen::ComputeFullU); Vec3 w = X_centered_sq_svd.singularValues(); Mat3 u = X_centered_sq_svd.matrixU(); for (size_t c = 0; c < 3; c++) { - double k = sqrt (w (c) / num_points); - X_control_points->col (c + 1) = mean + k * u.col (c); + double k = sqrt(w(c) / num_points); + X_control_points->col(c + 1) = mean + k * u.col(c); } } @@ -382,7 +381,7 @@ static void ComputeBarycentricCoordinates(const Mat3X &X_world_centered, const Mat34 &X_control_points, Mat4X *alphas) { size_t num_points = X_world_centered.cols(); - Mat3 C2 ; + Mat3 C2; for (size_t c = 1; c < 4; c++) { C2.col(c-1) = X_control_points.col(c) - X_control_points.col(0); } @@ -400,7 +399,7 @@ static void ComputeBarycentricCoordinates(const Mat3X &X_world_centered, // Estimates the coordinates of all real points in the camera coordinate frame static void ComputePointsCoordinatesInCameraFrame( - const Mat4X &alphas, + const Mat4X &alphas, const Vec4 &betas, const Eigen::Matrix &U, Mat3X *X_camera) { @@ -423,7 +422,7 @@ static void ComputePointsCoordinatesInCameraFrame( // Check the sign of the z coordinate of the points (should be positive) uint num_z_neg = 0; for (size_t i = 0; i < X_camera->cols(); ++i) { - if ((*X_camera)(2,i) < 0) { + if ((*X_camera)(2, i) < 0) { num_z_neg++; } } @@ -432,7 +431,7 @@ static void ComputePointsCoordinatesInCameraFrame( if (num_z_neg > 0.5 * X_camera->cols()) { C2b = -C2b; *X_camera = -(*X_camera); - } + } } bool EuclideanResectionEPnP(const Mat2X &x_camera, @@ -442,16 +441,16 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, CHECK(x_camera.cols() == X_world.cols()); CHECK(x_camera.cols() > 3); size_t num_points = X_world.cols(); - + // Select the control points. Mat34 X_control_points; Mat X_centered; SelectControlPoints(X_world, &X_centered, &X_control_points); - + // Compute the barycentric coordinates. Mat4X alphas(4, num_points); ComputeBarycentricCoordinates(X_centered, X_control_points, &alphas); - + // Estimates the M matrix with the barycentric coordinates Mat M(2 * num_points, 12); Eigen::Matrix sub_M; @@ -462,17 +461,17 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, double a3 = alphas(3, c); double ui = x_camera(0, c); double vi = x_camera(1, c); - M.block(2*c, 0, 2, 12) << a0, 0, + M.block(2*c, 0, 2, 12) << a0, 0, a0*(-ui), a1, 0, - a1*(-ui), a2, 0, + a1*(-ui), a2, 0, a2*(-ui), a3, 0, - a3*(-ui), 0, + a3*(-ui), 0, a0, a0*(-vi), 0, a1, a1*(-vi), 0, a2, a2*(-vi), 0, a3, a3*(-vi); } - + // TODO(julien): Avoid the transpose by rewriting the u2.block() calls. Eigen::JacobiSVD MtMsvd(M.transpose()*M, Eigen::ComputeFullU); Eigen::Matrix u2 = MtMsvd.matrixU().transpose(); @@ -495,18 +494,18 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, dv2.row(3) = u2.block(10, 3, 1, 3) - u2.block(10, 6, 1, 3); dv2.row(4) = u2.block(10, 3, 1, 3) - u2.block(10, 9, 1, 3); dv2.row(5) = u2.block(10, 6, 1, 3) - u2.block(10, 9, 1, 3); - dv3.row(0) = u2.block( 9, 0, 1, 3) - u2.block( 9, 3, 1, 3); - dv3.row(1) = u2.block( 9, 0, 1, 3) - u2.block( 9, 6, 1, 3); - dv3.row(2) = u2.block( 9, 0, 1, 3) - u2.block( 9, 9, 1, 3); - dv3.row(3) = u2.block( 9, 3, 1, 3) - u2.block( 9, 6, 1, 3); - dv3.row(4) = u2.block( 9, 3, 1, 3) - u2.block( 9, 9, 1, 3); - dv3.row(5) = u2.block( 9, 6, 1, 3) - u2.block( 9, 9, 1, 3); - dv4.row(0) = u2.block( 8, 0, 1, 3) - u2.block( 8, 3, 1, 3); - dv4.row(1) = u2.block( 8, 0, 1, 3) - u2.block( 8, 6, 1, 3); - dv4.row(2) = u2.block( 8, 0, 1, 3) - u2.block( 8, 9, 1, 3); - dv4.row(3) = u2.block( 8, 3, 1, 3) - u2.block( 8, 6, 1, 3); - dv4.row(4) = u2.block( 8, 3, 1, 3) - u2.block( 8, 9, 1, 3); - dv4.row(5) = u2.block( 8, 6, 1, 3) - u2.block( 8, 9, 1, 3); + dv3.row(0) = u2.block(9, 0, 1, 3) - u2.block(9, 3, 1, 3); + dv3.row(1) = u2.block(9, 0, 1, 3) - u2.block(9, 6, 1, 3); + dv3.row(2) = u2.block(9, 0, 1, 3) - u2.block(9, 9, 1, 3); + dv3.row(3) = u2.block(9, 3, 1, 3) - u2.block(9, 6, 1, 3); + dv3.row(4) = u2.block(9, 3, 1, 3) - u2.block(9, 9, 1, 3); + dv3.row(5) = u2.block(9, 6, 1, 3) - u2.block(9, 9, 1, 3); + dv4.row(0) = u2.block(8, 0, 1, 3) - u2.block(8, 3, 1, 3); + dv4.row(1) = u2.block(8, 0, 1, 3) - u2.block(8, 6, 1, 3); + dv4.row(2) = u2.block(8, 0, 1, 3) - u2.block(8, 9, 1, 3); + dv4.row(3) = u2.block(8, 3, 1, 3) - u2.block(8, 6, 1, 3); + dv4.row(4) = u2.block(8, 3, 1, 3) - u2.block(8, 9, 1, 3); + dv4.row(5) = u2.block(8, 6, 1, 3) - u2.block(8, 9, 1, 3); Eigen::Matrix L; for (size_t r = 0; r < 6; r++) { @@ -520,7 +519,7 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, 2.0 * dv2.row(r).dot(dv4.row(r)), 2.0 * dv3.row(r).dot(dv4.row(r)), dv4.row(r).dot(dv4.row(r)); - } + } Vec6 rho; rho << (X_control_points.col(0) - X_control_points.col(1)).squaredNorm(), (X_control_points.col(0) - X_control_points.col(2)).squaredNorm(), @@ -528,7 +527,7 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, (X_control_points.col(1) - X_control_points.col(2)).squaredNorm(), (X_control_points.col(1) - X_control_points.col(3)).squaredNorm(), (X_control_points.col(2) - X_control_points.col(3)).squaredNorm(); - + // There are three possible solutions based on the three approximations of L // (betas). Below, each one is solved for then the best one is chosen. Mat3X X_camera; @@ -548,7 +547,8 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, // check, is the right approach. So far this seems the case. // // TODO(sergey): Made it an option for now, in some cases it makes sense to - // still fallback to reprojection solution (see bug [#32765] from Blender bug tracker) + // still fallback to reprojection solution + // For details see bug [#32765] from Blender bug tracker // double kSuccessThreshold = std::numeric_limits::max(); double kSuccessThreshold = success_threshold; @@ -559,15 +559,15 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, Vec4 betas = Vec4::Zero(); Eigen::Matrix l_6x4; for (size_t r = 0; r < 6; r++) { - l_6x4.row(r) << L(r, 0), L(r, 1), L(r, 3), L(r, 6); + l_6x4.row(r) << L(r, 0), L(r, 1), L(r, 3), L(r, 6); } - Eigen::JacobiSVD svd_of_l4(l_6x4, + Eigen::JacobiSVD svd_of_l4(l_6x4, Eigen::ComputeFullU | Eigen::ComputeFullV); Vec4 b4 = svd_of_l4.solve(rho); if ((l_6x4 * b4).isApprox(rho, kSuccessThreshold)) { if (b4(0) < 0) { b4 = -b4; - } + } b4(0) = std::sqrt(b4(0)); betas << b4(0), b4(1) / b4(0), b4(2) / b4(0), b4(3) / b4(0); ComputePointsCoordinatesInCameraFrame(alphas, betas, u2, &X_camera); @@ -578,14 +578,14 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, ts[0].setZero(); rmse(0) = std::numeric_limits::max(); } - + // Find the second possible solution for R, t corresponding to: // Betas = [b00 b01 b11 b02 b12 b22 b03 b13 b23 b33] // Betas_approx_2 = [b00 b01 b11] betas.setZero(); Eigen::Matrix l_6x3; l_6x3 = L.block(0, 0, 6, 3); - Eigen::JacobiSVD svdOfL3(l_6x3, + Eigen::JacobiSVD svdOfL3(l_6x3, Eigen::ComputeFullU | Eigen::ComputeFullV); Vec3 b3 = svdOfL3.solve(rho); VLOG(2) << " rho = " << rho; @@ -611,14 +611,14 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, ts[1].setZero(); rmse(1) = std::numeric_limits::max(); } - + // Find the third possible solution for R, t corresponding to: // Betas = [b00 b01 b11 b02 b12 b22 b03 b13 b23 b33] // Betas_approx_3 = [b00 b01 b11 b02 b12] betas.setZero(); Eigen::Matrix l_6x5; l_6x5 = L.block(0, 0, 6, 5); - Eigen::JacobiSVD svdOfL5(l_6x5, + Eigen::JacobiSVD svdOfL5(l_6x5, Eigen::ComputeFullU | Eigen::ComputeFullV); Vec5 b5 = svdOfL5.solve(rho); if ((l_6x5 * b5).isApprox(rho, kSuccessThreshold)) { @@ -650,7 +650,7 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, ts[2].setZero(); rmse(2) = std::numeric_limits::max(); } - + // Finally, with all three solutions, select the (R, t) with the best RMSE. VLOG(2) << "RMSE for solution 0: " << rmse(0); VLOG(2) << "RMSE for solution 1: " << rmse(1); @@ -675,5 +675,5 @@ bool EuclideanResectionEPnP(const Mat2X &x_camera, return true; } -} // namespace resection -} // namespace libmv +} // namespace resection +} // namespace libmv diff --git a/extern/libmv/libmv/multiview/euclidean_resection.h b/extern/libmv/libmv/multiview/euclidean_resection.h index b0428ec61fd..1a329702c2a 100644 --- a/extern/libmv/libmv/multiview/euclidean_resection.h +++ b/extern/libmv/libmv/multiview/euclidean_resection.h @@ -26,7 +26,7 @@ namespace libmv { namespace euclidean_resection { - + enum ResectionMethod { RESECTION_ANSAR_DANIILIDIS, @@ -48,7 +48,7 @@ enum ResectionMethod { * \param success_threshold Threshold of an error which is still considered a success * (currently used by EPnP algorithm only) */ -bool EuclideanResection(const Mat2X &x_camera, +bool EuclideanResection(const Mat2X &x_camera, const Mat3X &X_world, Mat3 *R, Vec3 *t, ResectionMethod method = RESECTION_EPNP, @@ -67,7 +67,7 @@ bool EuclideanResection(const Mat2X &x_camera, * \param t Solution for the camera translation vector * \param method Resection method */ -bool EuclideanResection(const Mat &x_image, +bool EuclideanResection(const Mat &x_image, const Mat3X &X_world, const Mat3 &K, Mat3 *R, Vec3 *t, @@ -101,7 +101,7 @@ void AbsoluteOrientation(const Mat3X &X, * This is the algorithm described in: "Linear Pose Estimation from Points or * Lines", by Ansar, A. and Daniilidis, PAMI 2003. vol. 25, no. 5. */ -void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, +void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, const Mat3X &X_world, Mat3 *R, Vec3 *t); /** @@ -121,12 +121,12 @@ void EuclideanResectionAnsarDaniilidis(const Mat2X &x_camera, * \note: the non-linear optimization is not implemented here. */ bool EuclideanResectionEPnP(const Mat2X &x_camera, - const Mat3X &X_world, + const Mat3X &X_world, Mat3 *R, Vec3 *t, double success_threshold = 1e-3); -} // namespace euclidean_resection -} // namespace libmv +} // namespace euclidean_resection +} // namespace libmv #endif /* LIBMV_MULTIVIEW_EUCLIDEAN_RESECTION_H_ */ diff --git a/extern/libmv/libmv/multiview/fundamental.cc b/extern/libmv/libmv/multiview/fundamental.cc index 80f155e804d..96576724150 100644 --- a/extern/libmv/libmv/multiview/fundamental.cc +++ b/extern/libmv/libmv/multiview/fundamental.cc @@ -18,13 +18,14 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/multiview/fundamental.h" + #include "libmv/logging/logging.h" #include "libmv/numeric/numeric.h" #include "libmv/numeric/poly.h" #include "libmv/multiview/conditioning.h" #include "libmv/multiview/projection.h" #include "libmv/multiview/triangulation.h" -#include "libmv/multiview/fundamental.h" namespace libmv { @@ -304,7 +305,6 @@ void MotionFromEssential(const Mat3 &E, std::vector *ts) { Eigen::JacobiSVD USV(E, Eigen::ComputeFullU | Eigen::ComputeFullV); Mat3 U = USV.matrixU(); - Vec3 d = USV.singularValues(); Mat3 Vt = USV.matrixV().transpose(); // Last column of U is undetermined since d = (a a 0). diff --git a/extern/libmv/libmv/multiview/homography.cc b/extern/libmv/libmv/multiview/homography.cc index 538c62598c0..317eb3f85c6 100644 --- a/extern/libmv/libmv/multiview/homography.cc +++ b/extern/libmv/libmv/multiview/homography.cc @@ -18,8 +18,9 @@ // 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/homography.h" + +#include "libmv/logging/logging.h" #include "libmv/multiview/homography_parameterization.h" namespace libmv { @@ -55,30 +56,30 @@ static bool Homography2DFromCorrespondencesLinearEuc( Mat b = Mat::Zero(n * 3, 1); for (int i = 0; i < n; ++i) { int j = 3 * i; - L(j, 0) = x1(0, i); // a - L(j, 1) = x1(1, i); // b - L(j, 2) = 1.0; // c - L(j, 6) = -x2(0, i) * x1(0, i); // g - L(j, 7) = -x2(0, i) * x1(1, i); // h - b(j, 0) = x2(0, i); // i + L(j, 0) = x1(0, i); // a + L(j, 1) = x1(1, i); // b + L(j, 2) = 1.0; // c + L(j, 6) = -x2(0, i) * x1(0, i); // g + L(j, 7) = -x2(0, i) * x1(1, i); // h + b(j, 0) = x2(0, i); // i ++j; - L(j, 3) = x1(0, i); // d - L(j, 4) = x1(1, i); // e - L(j, 5) = 1.0; // f - L(j, 6) = -x2(1, i) * x1(0, i); // g - L(j, 7) = -x2(1, i) * x1(1, i); // h - b(j, 0) = x2(1, i); // i - + L(j, 3) = x1(0, i); // d + L(j, 4) = x1(1, i); // e + L(j, 5) = 1.0; // f + L(j, 6) = -x2(1, i) * x1(0, i); // g + L(j, 7) = -x2(1, i) * x1(1, i); // h + b(j, 0) = x2(1, i); // i + // This ensures better stability // TODO(julien) make a lite version without this 3rd set ++j; - L(j, 0) = x2(1, i) * x1(0, i); // a - L(j, 1) = x2(1, i) * x1(1, i); // b - L(j, 2) = x2(1, i); // c - L(j, 3) = -x2(0, i) * x1(0, i); // d - L(j, 4) = -x2(0, i) * x1(1, i); // e - L(j, 5) = -x2(0, i) ; // f + L(j, 0) = x2(1, i) * x1(0, i); // a + L(j, 1) = x2(1, i) * x1(1, i); // b + L(j, 2) = x2(1, i); // c + L(j, 3) = -x2(0, i) * x1(0, i); // d + L(j, 4) = -x2(0, i) * x1(1, i); // e + L(j, 5) = -x2(0, i); // f } // Solve Lx=B Vec h = L.fullPivLu().solve(b); @@ -103,7 +104,7 @@ bool Homography2DFromCorrespondencesLinear(const Mat &x1, Mat3 *H, double expected_precision) { if (x1.rows() == 2) { - return Homography2DFromCorrespondencesLinearEuc(x1, x2, H, + return Homography2DFromCorrespondencesLinearEuc(x1, x2, H, expected_precision); } assert(3 == x1.rows()); @@ -119,29 +120,29 @@ bool Homography2DFromCorrespondencesLinear(const Mat &x1, Mat b = Mat::Zero(n * 3, 1); for (int i = 0; i < n; ++i) { int j = 3 * i; - L(j, 0) = x2(w, i) * x1(x, i);//a - L(j, 1) = x2(w, i) * x1(y, i);//b - L(j, 2) = x2(w, i) * x1(w, i);//c - L(j, 6) = -x2(x, i) * x1(x, i);//g - L(j, 7) = -x2(x, i) * x1(y, i);//h + L(j, 0) = x2(w, i) * x1(x, i); // a + L(j, 1) = x2(w, i) * x1(y, i); // b + L(j, 2) = x2(w, i) * x1(w, i); // c + L(j, 6) = -x2(x, i) * x1(x, i); // g + L(j, 7) = -x2(x, i) * x1(y, i); // h b(j, 0) = x2(x, i) * x1(w, i); ++j; - L(j, 3) = x2(w, i) * x1(x, i);//d - L(j, 4) = x2(w, i) * x1(y, i);//e - L(j, 5) = x2(w, i) * x1(w, i);//f - L(j, 6) = -x2(y, i) * x1(x, i);//g - L(j, 7) = -x2(y, i) * x1(y, i);//h + L(j, 3) = x2(w, i) * x1(x, i); // d + L(j, 4) = x2(w, i) * x1(y, i); // e + L(j, 5) = x2(w, i) * x1(w, i); // f + L(j, 6) = -x2(y, i) * x1(x, i); // g + L(j, 7) = -x2(y, i) * x1(y, i); // h b(j, 0) = x2(y, i) * x1(w, i); // This ensures better stability ++j; - L(j, 0) = x2(y, i) * x1(x, i);//a - L(j, 1) = x2(y, i) * x1(y, i);//b - L(j, 2) = x2(y, i) * x1(w, i);//c - L(j, 3) = -x2(x, i) * x1(x, i);//d - L(j, 4) = -x2(x, i) * x1(y, i);//e - L(j, 5) = -x2(x, i) * x1(w, i);//f + L(j, 0) = x2(y, i) * x1(x, i); // a + L(j, 1) = x2(y, i) * x1(y, i); // b + L(j, 2) = x2(y, i) * x1(w, i); // c + L(j, 3) = -x2(x, i) * x1(x, i); // d + L(j, 4) = -x2(x, i) * x1(y, i); // e + L(j, 5) = -x2(x, i) * x1(w, i); // f } // Solve Lx=B Vec h = L.fullPivLu().solve(b); @@ -165,11 +166,11 @@ bool Homography2DFromCorrespondencesLinear(const Mat &x1, * |0 0 0 1| -> |-x2w| |-1 0 0 0| -> | x2x| |0 0 0 0| -> | 0 | * |0 0 0 0| | 0 | | 0 0 0 0| | 0 | |0 0 0 1| |-x2w| * |0-1 0 0| | x2y| | 0 0 0 0| | 0 | |0 0-1 0| | x2z| - * |a b c d| - * A = |e f g h| + * |a b c d| + * A = |e f g h| * |i j k l| - * |m n o 1| - * + * |m n o 1| + * * x2^t * H1 * A *x1 = (-x2w*a +x2x*m )*x1x + (-x2w*b +x2x*n )*x1y + (-x2w*c +x2x*o )*x1z + (-x2w*d +x2x*1 )*x1w = 0 * x2^t * H2 * A *x1 = (-x2z*e +x2y*i )*x1x + (-x2z*f +x2y*j )*x1y + (-x2z*g +x2y*k )*x1z + (-x2z*h +x2y*l )*x1w = 0 * x2^t * H3 * A *x1 = (-x2z*a +x2x*i )*x1x + (-x2z*b +x2x*j )*x1y + (-x2z*c +x2x*k )*x1z + (-x2z*d +x2x*l )*x1w = 0 @@ -196,64 +197,64 @@ bool Homography3DFromCorrespondencesLinear(const Mat &x1, Mat b = Mat::Zero(n * 6, 1); for (int i = 0; i < n; ++i) { int j = 6 * i; - L(j, 0) = -x2(w, i) * x1(x, i);//a - L(j, 1) = -x2(w, i) * x1(y, i);//b - L(j, 2) = -x2(w, i) * x1(z, i);//c - L(j, 3) = -x2(w, i) * x1(w, i);//d - L(j,12) = x2(x, i) * x1(x, i);//m - L(j,13) = x2(x, i) * x1(y, i);//n - L(j,14) = x2(x, i) * x1(z, i);//o - b(j, 0) = -x2(x, i) * x1(w, i); - + L(j, 0) = -x2(w, i) * x1(x, i); // a + L(j, 1) = -x2(w, i) * x1(y, i); // b + L(j, 2) = -x2(w, i) * x1(z, i); // c + L(j, 3) = -x2(w, i) * x1(w, i); // d + L(j, 12) = x2(x, i) * x1(x, i); // m + L(j, 13) = x2(x, i) * x1(y, i); // n + L(j, 14) = x2(x, i) * x1(z, i); // o + b(j, 0) = -x2(x, i) * x1(w, i); + ++j; - L(j, 4) = -x2(z, i) * x1(x, i);//e - L(j, 5) = -x2(z, i) * x1(y, i);//f - L(j, 6) = -x2(z, i) * x1(z, i);//g - L(j, 7) = -x2(z, i) * x1(w, i);//h - L(j, 8) = x2(y, i) * x1(x, i);//i - L(j, 9) = x2(y, i) * x1(y, i);//j - L(j,10) = x2(y, i) * x1(z, i);//k - L(j,11) = x2(y, i) * x1(w, i);//l - + L(j, 4) = -x2(z, i) * x1(x, i); // e + L(j, 5) = -x2(z, i) * x1(y, i); // f + L(j, 6) = -x2(z, i) * x1(z, i); // g + L(j, 7) = -x2(z, i) * x1(w, i); // h + L(j, 8) = x2(y, i) * x1(x, i); // i + L(j, 9) = x2(y, i) * x1(y, i); // j + L(j, 10) = x2(y, i) * x1(z, i); // k + L(j, 11) = x2(y, i) * x1(w, i); // l + ++j; - L(j, 0) = -x2(z, i) * x1(x, i);//a - L(j, 1) = -x2(z, i) * x1(y, i);//b - L(j, 2) = -x2(z, i) * x1(z, i);//c - L(j, 3) = -x2(z, i) * x1(w, i);//d - L(j, 8) = x2(x, i) * x1(x, i);//i - L(j, 9) = x2(x, i) * x1(y, i);//j - L(j,10) = x2(x, i) * x1(z, i);//k - L(j,11) = x2(x, i) * x1(w, i);//l - + L(j, 0) = -x2(z, i) * x1(x, i); // a + L(j, 1) = -x2(z, i) * x1(y, i); // b + L(j, 2) = -x2(z, i) * x1(z, i); // c + L(j, 3) = -x2(z, i) * x1(w, i); // d + L(j, 8) = x2(x, i) * x1(x, i); // i + L(j, 9) = x2(x, i) * x1(y, i); // j + L(j, 10) = x2(x, i) * x1(z, i); // k + L(j, 11) = x2(x, i) * x1(w, i); // l + ++j; - L(j, 4) = -x2(w, i) * x1(x, i);//e - L(j, 5) = -x2(w, i) * x1(y, i);//f - L(j, 6) = -x2(w, i) * x1(z, i);//g - L(j, 7) = -x2(w, i) * x1(w, i);//h - L(j,12) = x2(y, i) * x1(x, i);//m - L(j,13) = x2(y, i) * x1(y, i);//n - L(j,14) = x2(y, i) * x1(z, i);//o - b(j, 0) = -x2(y, i) * x1(w, i); - + L(j, 4) = -x2(w, i) * x1(x, i); // e + L(j, 5) = -x2(w, i) * x1(y, i); // f + L(j, 6) = -x2(w, i) * x1(z, i); // g + L(j, 7) = -x2(w, i) * x1(w, i); // h + L(j, 12) = x2(y, i) * x1(x, i); // m + L(j, 13) = x2(y, i) * x1(y, i); // n + L(j, 14) = x2(y, i) * x1(z, i); // o + b(j, 0) = -x2(y, i) * x1(w, i); + ++j; - L(j, 0) = -x2(y, i) * x1(x, i);//a - L(j, 1) = -x2(y, i) * x1(y, i);//b - L(j, 2) = -x2(y, i) * x1(z, i);//c - L(j, 3) = -x2(y, i) * x1(w, i);//d - L(j, 4) = x2(x, i) * x1(x, i);//e - L(j, 5) = x2(x, i) * x1(y, i);//f - L(j, 6) = x2(x, i) * x1(z, i);//g - L(j, 7) = x2(x, i) * x1(w, i);//h - + L(j, 0) = -x2(y, i) * x1(x, i); // a + L(j, 1) = -x2(y, i) * x1(y, i); // b + L(j, 2) = -x2(y, i) * x1(z, i); // c + L(j, 3) = -x2(y, i) * x1(w, i); // d + L(j, 4) = x2(x, i) * x1(x, i); // e + L(j, 5) = x2(x, i) * x1(y, i); // f + L(j, 6) = x2(x, i) * x1(z, i); // g + L(j, 7) = x2(x, i) * x1(w, i); // h + ++j; - L(j, 8) = -x2(w, i) * x1(x, i);//i - L(j, 9) = -x2(w, i) * x1(y, i);//j - L(j,10) = -x2(w, i) * x1(z, i);//k - L(j,11) = -x2(w, i) * x1(w, i);//l - L(j,12) = x2(z, i) * x1(x, i);//m - L(j,13) = x2(z, i) * x1(y, i);//n - L(j,14) = x2(z, i) * x1(z, i);//o - b(j, 0) = -x2(z, i) * x1(w, i); + L(j, 8) = -x2(w, i) * x1(x, i); // i + L(j, 9) = -x2(w, i) * x1(y, i); // j + L(j, 10) = -x2(w, i) * x1(z, i); // k + L(j, 11) = -x2(w, i) * x1(w, i); // l + L(j, 12) = x2(z, i) * x1(x, i); // m + L(j, 13) = x2(z, i) * x1(y, i); // n + L(j, 14) = x2(z, i) * x1(z, i); // o + b(j, 0) = -x2(z, i) * x1(w, i); } // Solve Lx=B Vec h = L.fullPivLu().solve(b); diff --git a/extern/libmv/libmv/multiview/homography.h b/extern/libmv/libmv/multiview/homography.h index a295c4366b6..8d2dff930eb 100644 --- a/extern/libmv/libmv/multiview/homography.h +++ b/extern/libmv/libmv/multiview/homography.h @@ -50,12 +50,12 @@ namespace libmv { bool Homography2DFromCorrespondencesLinear(const Mat &x1, const Mat &x2, Mat3 *H, - double expected_precision = + double expected_precision = EigenDouble::dummy_precision()); /** * 3D Homography transformation estimation. - * + * * This function can be used in order to estimate the homography transformation * from a list of 3D correspondences. * @@ -68,16 +68,16 @@ bool Homography2DFromCorrespondencesLinear(const Mat &x1, * |m n o 1| * \param[in] expected_precision The expected precision in order for instance * to accept almost homography matrices. - * + * * \return true if the transformation estimation has succeeded - * + * * \note Need at least 5 non coplanar points * \note Points coordinates must be in homogeneous coordinates */ bool Homography3DFromCorrespondencesLinear(const Mat &x1, const Mat &x2, Mat4 *H, - double expected_precision = + double expected_precision = EigenDouble::dummy_precision()); /** @@ -87,6 +87,6 @@ bool Homography3DFromCorrespondencesLinear(const Mat &x1, */ double SymmetricGeometricDistance(Mat3 &H, Vec2 &x1, Vec2 &x2); -} // namespace libmv +} // namespace libmv #endif // LIBMV_MULTIVIEW_HOMOGRAPHY_H_ diff --git a/extern/libmv/libmv/multiview/homography_parameterization.h b/extern/libmv/libmv/multiview/homography_parameterization.h index b31642eea15..ca8fbd8066e 100644 --- a/extern/libmv/libmv/multiview/homography_parameterization.h +++ b/extern/libmv/libmv/multiview/homography_parameterization.h @@ -1,15 +1,15 @@ // 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 @@ -40,7 +40,7 @@ class Homography2DNormalizedParameterization { typedef Eigen::Matrix Parameterized; // H /// Convert from the 8 parameters to a H matrix. - static void To(const Parameters &p, Parameterized *h) { + static void To(const Parameters &p, Parameterized *h) { *h << p(0), p(1), p(2), p(3), p(4), p(5), p(6), p(7), 1.0; @@ -70,7 +70,7 @@ class Homography3DNormalizedParameterization { typedef Eigen::Matrix Parameterized; // H /// Convert from the 15 parameters to a H matrix. - static void To(const Parameters &p, Parameterized *h) { + static void To(const Parameters &p, Parameterized *h) { *h << p(0), p(1), p(2), p(3), p(4), p(5), p(6), p(7), p(8), p(9), p(10), p(11), @@ -86,6 +86,6 @@ class Homography3DNormalizedParameterization { } }; -} // namespace libmv +} // namespace libmv #endif // LIBMV_MULTIVIEW_HOMOGRAPHY_PARAMETERIZATION_H_ diff --git a/extern/libmv/libmv/multiview/nviewtriangulation.h b/extern/libmv/libmv/multiview/nviewtriangulation.h index b4f521f185d..f4614ab1a5c 100644 --- a/extern/libmv/libmv/multiview/nviewtriangulation.h +++ b/extern/libmv/libmv/multiview/nviewtriangulation.h @@ -1,15 +1,15 @@ // 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 @@ -17,7 +17,7 @@ // 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. -// +// // Compute a 3D position of a point from several images of it. In particular, // compute the projective point X in R^4 such that x = PX. // diff --git a/extern/libmv/libmv/multiview/panography.cc b/extern/libmv/libmv/multiview/panography.cc new file mode 100644 index 00000000000..b62802948c4 --- /dev/null +++ b/extern/libmv/libmv/multiview/panography.cc @@ -0,0 +1,125 @@ +// 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/panography.h" + +namespace libmv { + +static bool Build_Minimal2Point_PolynomialFactor( + const Mat & x1, const Mat & x2, + double * P) { // P must be a double[4] + assert(2 == x1.rows()); + assert(2 == x1.cols()); + assert(x1.rows() == x2.rows()); + assert(x1.cols() == x2.cols()); + + // Setup the variable of the input problem: + Vec xx1 = (x1.col(0)).transpose(); + Vec yx1 = (x1.col(1)).transpose(); + + double a12 = xx1.dot(yx1); + Vec xx2 = (x2.col(0)).transpose(); + Vec yx2 = (x2.col(1)).transpose(); + double b12 = xx2.dot(yx2); + + double a1 = xx1.squaredNorm(); + double a2 = yx1.squaredNorm(); + + double b1 = xx2.squaredNorm(); + double b2 = yx2.squaredNorm(); + + // Build the 3rd degre polynomial in F^2. + // + // f^6 * p + f^4 * q + f^2* r + s = 0; + // + // Coefficients in ascending powers of alpha, i.e. P[N]*x^N. + // Run panography_coeffs.py to get the below coefficients. + P[0] = b1*b2*a12*a12-a1*a2*b12*b12; + P[1] = -2*a1*a2*b12+2*a12*b1*b2+b1*a12*a12+b2*a12*a12-a1*b12*b12-a2*b12*b12; + P[2] = b1*b2-a1*a2-2*a1*b12-2*a2*b12+2*a12*b1+2*a12*b2+a12*a12-b12*b12; + P[3] = b1+b2-2*b12-a1-a2+2*a12; + + // If P[3] equal to 0 we get ill conditionned data + return (P[3] != 0.0); +} + +// This implements a minimal solution (2 points) for panoramic stitching: +// +// http://www.cs.ubc.ca/~mbrown/minimal/minimal.html +// +// [1] M. Brown and R. Hartley and D. Nister. Minimal Solutions for Panoramic +// Stitching. CVPR07. +void F_FromCorrespondance_2points(const Mat &x1, const Mat &x2, + vector *fs) { + // Build Polynomial factor to get squared focal value. + double P[4]; + Build_Minimal2Point_PolynomialFactor(x1, x2, &P[0]); + + // Solve it by using F = f^2 and a Cubic polynomial solver + // + // F^3 * p + F^2 * q + F^1 * r + s = 0 + // + double roots[3]; + int num_roots = SolveCubicPolynomial(P, roots); + for (int i = 0; i < num_roots; ++i) { + if (roots[i] > 0.0) { + fs->push_back(sqrt(roots[i])); + } + } +} + +// Compute the 3x3 rotation matrix that fits two 3D point clouds in the least +// square sense. The method is from: +// +// K. Arun,T. Huand and D. Blostein. Least-squares fitting of 2 3-D point +// sets. IEEE Transactions on Pattern Analysis and Machine Intelligence, +// 9:698-700, 1987. +void GetR_FixedCameraCenter(const Mat &x1, const Mat &x2, + const double focal, + Mat3 *R) { + assert(3 == x1.rows()); + assert(2 <= x1.cols()); + assert(x1.rows() == x2.rows()); + assert(x1.cols() == x2.cols()); + + // Build simplified K matrix + Mat3 K(Mat3::Identity() * 1.0/focal); + K(2, 2)= 1.0; + + // Build the correlation matrix; equation (22) in [1]. + Mat3 C = Mat3::Zero(); + for (int i = 0; i < x1.cols(); ++i) { + Mat r1i = (K * x1.col(i)).normalized(); + Mat r2i = (K * x2.col(i)).normalized(); + C += r2i * r1i.transpose(); + } + + // Solve for rotation. Equations (24) and (25) in [1]. + Eigen::JacobiSVD svd(C, Eigen::ComputeThinU | Eigen::ComputeThinV); + Mat3 scale = Mat3::Identity(); + scale(2, 2) = ((svd.matrixU() * svd.matrixV().transpose()).determinant() > 0.0) + ? 1.0 + : -1.0; + + (*R) = svd.matrixU() * scale * svd.matrixV().transpose(); +} + +} // namespace libmv diff --git a/extern/libmv/libmv/multiview/panography.h b/extern/libmv/libmv/multiview/panography.h index 6f01beb0ffd..6e87bd71304 100644 --- a/extern/libmv/libmv/multiview/panography.h +++ b/extern/libmv/libmv/multiview/panography.h @@ -28,45 +28,6 @@ namespace libmv { -static bool Build_Minimal2Point_PolynomialFactor( - const Mat & x1, const Mat & x2, - double * P) //P must be a double[4] -{ - assert(2 == x1.rows()); - assert(2 == x1.cols()); - assert(x1.rows() == x2.rows()); - assert(x1.cols() == x2.cols()); - - // Setup the variable of the input problem: - Vec xx1 = (x1.col(0)).transpose(); - Vec yx1 = (x1.col(1)).transpose(); - - double a12 = xx1.dot(yx1); - Vec xx2 = (x2.col(0)).transpose(); - Vec yx2 = (x2.col(1)).transpose(); - double b12 = xx2.dot(yx2); - - double a1 = xx1.squaredNorm(); - double a2 = yx1.squaredNorm(); - - double b1 = xx2.squaredNorm(); - double b2 = yx2.squaredNorm(); - - // Build the 3rd degre polynomial in F^2. - // - // f^6 * p + f^4 * q + f^2* r + s = 0; - // - // Coefficients in ascending powers of alpha, i.e. P[N]*x^N. - // Run panography_coeffs.py to get the below coefficients. - P[0] = b1*b2*a12*a12-a1*a2*b12*b12; - P[1] = -2*a1*a2*b12+2*a12*b1*b2+b1*a12*a12+b2*a12*a12-a1*b12*b12-a2*b12*b12; - P[2] = b1*b2-a1*a2-2*a1*b12-2*a2*b12+2*a12*b1+2*a12*b2+a12*a12-b12*b12; - P[3] = b1+b2-2*b12-a1-a2+2*a12; - - // If P[3] equal to 0 we get ill conditionned data - return (P[3] != 0.0); -} - // This implements a minimal solution (2 points) for panoramic stitching: // // http://www.cs.ubc.ca/~mbrown/minimal/minimal.html @@ -92,25 +53,8 @@ static bool Build_Minimal2Point_PolynomialFactor( // K = [0 f 0] // [0 0 1] // -static void F_FromCorrespondance_2points(const Mat &x1, const Mat &x2, - vector *fs) { - - // Build Polynomial factor to get squared focal value. - double P[4]; - Build_Minimal2Point_PolynomialFactor(x1, x2, &P[0]); - - // Solve it by using F = f^2 and a Cubic polynomial solver - // - // F^3 * p + F^2 * q + F^1 * r + s = 0 - // - double roots[3]; - int num_roots = SolveCubicPolynomial(P, roots); - for (int i = 0; i < num_roots; ++i) { - if (roots[i] > 0.0) { - fs->push_back(sqrt(roots[i])); - } - } -} +void F_FromCorrespondance_2points(const Mat &x1, const Mat &x2, + vector *fs); // Compute the 3x3 rotation matrix that fits two 3D point clouds in the least // square sense. The method is from: @@ -146,36 +90,10 @@ static void F_FromCorrespondance_2points(const Mat &x1, const Mat &x2, // // R = arg min || X2 - R * x1 || // -static void GetR_FixedCameraCenter(const Mat &x1, const Mat &x2, - const double focal, - Mat3 *R) { - assert(3 == x1.rows()); - assert(2 <= x1.cols()); - assert(x1.rows() == x2.rows()); - assert(x1.cols() == x2.cols()); - - // Build simplified K matrix - Mat3 K( Mat3::Identity() * 1.0/focal ); - K(2,2)= 1.0; - - // Build the correlation matrix; equation (22) in [1]. - Mat3 C = Mat3::Zero(); - for(int i = 0; i < x1.cols(); ++i) { - Mat r1i = (K * x1.col(i)).normalized(); - Mat r2i = (K * x2.col(i)).normalized(); - C += r2i * r1i.transpose(); - } - - // Solve for rotation. Equations (24) and (25) in [1]. - Eigen::JacobiSVD svd(C, Eigen::ComputeThinU | Eigen::ComputeThinV); - Mat3 scale = Mat3::Identity(); - scale(2,2) = ((svd.matrixU() * svd.matrixV().transpose()).determinant() > 0.0) - ? 1.0 - : -1.0; - - (*R) = svd.matrixU() * scale * svd.matrixV().transpose(); -} +void GetR_FixedCameraCenter(const Mat &x1, const Mat &x2, + const double focal, + Mat3 *R); } // namespace libmv -#endif // LIBMV_MULTIVIEW_PANOGRAPHY_H +#endif // LIBMV_MULTIVIEW_PANOGRAPHY_H diff --git a/extern/libmv/libmv/multiview/projection.cc b/extern/libmv/libmv/multiview/projection.cc index a7d1a058e9c..f8bece3de68 100644 --- a/extern/libmv/libmv/multiview/projection.cc +++ b/extern/libmv/libmv/multiview/projection.cc @@ -1,15 +1,15 @@ // 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 @@ -38,10 +38,11 @@ void KRt_From_P(const Mat34 &P, Mat3 *Kp, Mat3 *Rp, Vec3 *tp) { // Set K(2,1) to zero. if (K(2, 1) != 0) { - double c = -K(2,2); - double s = K(2,1); + double c = -K(2, 2); + double s = K(2, 1); double l = sqrt(c * c + s * s); - c /= l; s /= l; + c /= l; + s /= l; Mat3 Qx; Qx << 1, 0, 0, 0, c, -s, @@ -54,7 +55,8 @@ void KRt_From_P(const Mat34 &P, Mat3 *Kp, Mat3 *Rp, Vec3 *tp) { double c = K(2, 2); double s = K(2, 0); double l = sqrt(c * c + s * s); - c /= l; s /= l; + c /= l; + s /= l; Mat3 Qy; Qy << c, 0, s, 0, 1, 0, @@ -67,11 +69,12 @@ void KRt_From_P(const Mat34 &P, Mat3 *Kp, Mat3 *Rp, Vec3 *tp) { double c = -K(1, 1); double s = K(1, 0); double l = sqrt(c * c + s * s); - c /= l; s /= l; + c /= l; + s /= l; Mat3 Qz; - Qz << c,-s, 0, - s, c, 0, - 0, 0, 1; + Qz << c, -s, 0, + s, c, 0, + 0, 0, 1; K = K * Qz; Q = Qz.transpose() * Q; } @@ -83,19 +86,19 @@ void KRt_From_P(const Mat34 &P, Mat3 *Kp, Mat3 *Rp, Vec3 *tp) { // - K(0,0) > 0 // - K(2,2) = 1 // - det(R) = 1 - if (K(2,2) < 0) { + if (K(2, 2) < 0) { K = -K; R = -R; } - if (K(1,1) < 0) { + if (K(1, 1) < 0) { Mat3 S; - S << 1, 0, 0, - 0,-1, 0, - 0, 0, 1; + S << 1, 0, 0, + 0, -1, 0, + 0, 0, 1; K = K * S; R = S * R; } - if (K(0,0) < 0) { + if (K(0, 0) < 0) { Mat3 S; S << -1, 0, 0, 0, 1, 0, @@ -106,13 +109,13 @@ void KRt_From_P(const Mat34 &P, Mat3 *Kp, Mat3 *Rp, Vec3 *tp) { // Compute translation. Vec p(3); - p << P(0,3), P(1,3), P(2,3); + p << P(0, 3), P(1, 3), P(2, 3); // TODO(pau) This sould be done by a SolveLinearSystem(A, b, &x) call. // TODO(keir) use the eigen LU solver syntax... Vec3 t = K.inverse() * p; // scale K so that K(2,2) = 1 - K = K / K(2,2); + K = K / K(2, 2); *Kp = K; *Rp = R; @@ -140,10 +143,10 @@ void ProjectionChangeAspectRatio(const Mat34 &P, 0, aspect_ratio_new / aspect_ratio, 0, 0, 0, 1; Mat34 P_temp; - - ProjectionShiftPrincipalPoint(P, principal_point, Vec2(0,0), &P_temp); + + ProjectionShiftPrincipalPoint(P, principal_point, Vec2(0, 0), &P_temp); P_temp = T * P_temp; - ProjectionShiftPrincipalPoint(P_temp, Vec2(0,0), principal_point, P_new); + ProjectionShiftPrincipalPoint(P_temp, Vec2(0, 0), principal_point, P_new); } void HomogeneousToEuclidean(const Mat &H, Mat *X) { @@ -198,15 +201,15 @@ void EuclideanToHomogeneous(const Vec3 &X, Vec4 *H) { // TODO(julien) Call conditioning.h/ApplyTransformationToPoints ? void EuclideanToNormalizedCamera(const Mat2X &x, const Mat3 &K, Mat2X *n) { - Mat3X x_image_h; - EuclideanToHomogeneous(x, &x_image_h); - Mat3X x_camera_h = K.inverse() * x_image_h; - HomogeneousToEuclidean(x_camera_h, n); + Mat3X x_image_h; + EuclideanToHomogeneous(x, &x_image_h); + Mat3X x_camera_h = K.inverse() * x_image_h; + HomogeneousToEuclidean(x_camera_h, n); } void HomogeneousToNormalizedCamera(const Mat3X &x, const Mat3 &K, Mat2X *n) { - Mat3X x_camera_h = K.inverse() * x; - HomogeneousToEuclidean(x_camera_h, n); + Mat3X x_camera_h = K.inverse() * x; + HomogeneousToEuclidean(x_camera_h, n); } double Depth(const Mat3 &R, const Vec3 &t, const Vec3 &X) { diff --git a/extern/libmv/libmv/multiview/projection.h b/extern/libmv/libmv/multiview/projection.h index bc353b64c08..3220bc2dbbc 100644 --- a/extern/libmv/libmv/multiview/projection.h +++ b/extern/libmv/libmv/multiview/projection.h @@ -94,7 +94,7 @@ inline Vec3 EuclideanToHomogeneous(const Vec2 &x) { inline Vec4 EuclideanToHomogeneous(const Vec3 &x) { return Vec4(x(0), x(1), x(2), 1); } -// Conversion from image coordinates to normalized camera coordinates +// Conversion from image coordinates to normalized camera coordinates void EuclideanToNormalizedCamera(const Mat2X &x, const Mat3 &K, Mat2X *n); void HomogeneousToNormalizedCamera(const Mat3X &x, const Mat3 &K, Mat2X *n); @@ -180,16 +180,16 @@ double Depth(const Mat3 &R, const Vec3 &t, const Vec4 &X); * Returns true if the homogenious 3D point X is in front of * the camera P. */ -inline bool isInFrontOfCamera(const Mat34 &P, const Vec4 &X){ +inline bool isInFrontOfCamera(const Mat34 &P, const Vec4 &X) { double condition_1 = P.row(2).dot(X) * X[3]; double condition_2 = X[2] * X[3]; - if( condition_1 > 0 && condition_2 > 0 ) + if (condition_1 > 0 && condition_2 > 0) return true; else return false; } -inline bool isInFrontOfCamera(const Mat34 &P, const Vec3 &X){ +inline bool isInFrontOfCamera(const Mat34 &P, const Vec3 &X) { Vec4 X_homo; X_homo.segment<3>(0) = X; X_homo(3) = 1; @@ -200,14 +200,14 @@ inline bool isInFrontOfCamera(const Mat34 &P, const Vec3 &X){ * Transforms a 2D point from pixel image coordinates to a 2D point in * normalized image coordinates. */ -inline Vec2 ImageToNormImageCoordinates(Mat3 &Kinverse, Vec2 &x){ +inline Vec2 ImageToNormImageCoordinates(Mat3 &Kinverse, Vec2 &x) { Vec3 x_h = Kinverse*EuclideanToHomogeneous(x); return HomogeneousToEuclidean( x_h ); } /// Estimates the root mean square error (2D) -inline double RootMeanSquareError(const Mat2X &x_image, - const Mat4X &X_world, +inline double RootMeanSquareError(const Mat2X &x_image, + const Mat4X &X_world, const Mat34 &P) { size_t num_points = x_image.cols(); Mat2X dx = Project(P, X_world) - x_image; @@ -215,10 +215,10 @@ inline double RootMeanSquareError(const Mat2X &x_image, } /// Estimates the root mean square error (2D) -inline double RootMeanSquareError(const Mat2X &x_image, - const Mat3X &X_world, - const Mat3 &K, - const Mat3 &R, +inline double RootMeanSquareError(const Mat2X &x_image, + const Mat3X &X_world, + const Mat3 &K, + const Mat3 &R, const Vec3 &t) { Mat34 P; P_From_KRt(K, R, t, &P); @@ -226,6 +226,6 @@ inline double RootMeanSquareError(const Mat2X &x_image, Mat2X dx = Project(P, X_world) - x_image; return dx.norm() / num_points; } -} // namespace libmv +} // namespace libmv #endif // LIBMV_MULTIVIEW_PROJECTION_H_ diff --git a/extern/libmv/libmv/multiview/resection.h b/extern/libmv/libmv/multiview/resection.h index e4623899147..c142d6deeb2 100644 --- a/extern/libmv/libmv/multiview/resection.h +++ b/extern/libmv/libmv/multiview/resection.h @@ -1,15 +1,15 @@ // 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 @@ -17,7 +17,7 @@ // 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. -// +// // Compute the projection matrix from a set of 3D points X and their // projections x = PX in 2D. This is useful if a point cloud is reconstructed. // @@ -46,10 +46,10 @@ void Resection(const Matrix &x, T xi = x(0, i); T yi = x(1, i); // See equation (7.2) on page 179 of H&Z. - design.template block<1,4>(2*i, 4) = -X.col(i).transpose(); - design.template block<1,4>(2*i, 8) = yi*X.col(i).transpose(); - design.template block<1,4>(2*i + 1, 0) = X.col(i).transpose(); - design.template block<1,4>(2*i + 1, 8) = -xi*X.col(i).transpose(); + design.template block<1, 4>(2*i, 4) = -X.col(i).transpose(); + design.template block<1, 4>(2*i, 8) = yi*X.col(i).transpose(); + design.template block<1, 4>(2*i + 1, 0) = X.col(i).transpose(); + design.template block<1, 4>(2*i + 1, 8) = -xi*X.col(i).transpose(); } Matrix p; Nullspace(&design, &p); diff --git a/extern/libmv/libmv/multiview/triangulation.cc b/extern/libmv/libmv/multiview/triangulation.cc index b9d38cef936..4d146c8f21b 100644 --- a/extern/libmv/libmv/multiview/triangulation.cc +++ b/extern/libmv/libmv/multiview/triangulation.cc @@ -1,15 +1,15 @@ // 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 @@ -18,9 +18,10 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/multiview/triangulation.h" + #include "libmv/numeric/numeric.h" #include "libmv/multiview/projection.h" -#include "libmv/multiview/triangulation.h" namespace libmv { @@ -30,10 +31,10 @@ void TriangulateDLT(const Mat34 &P1, const Vec2 &x1, Vec4 *X_homogeneous) { Mat4 design; for (int i = 0; i < 4; ++i) { - design(0,i) = x1(0) * P1(2,i) - P1(0,i); - design(1,i) = x1(1) * P1(2,i) - P1(1,i); - design(2,i) = x2(0) * P2(2,i) - P2(0,i); - design(3,i) = x2(1) * P2(2,i) - P2(1,i); + design(0, i) = x1(0) * P1(2, i) - P1(0, i); + design(1, i) = x1(1) * P1(2, i) - P1(1, i); + design(2, i) = x2(0) * P2(2, i) - P2(0, i); + design(3, i) = x2(1) * P2(2, i) - P2(1, i); } Nullspace(&design, X_homogeneous); } diff --git a/extern/libmv/libmv/multiview/triangulation.h b/extern/libmv/libmv/multiview/triangulation.h index c35774d9e1b..be878890242 100644 --- a/extern/libmv/libmv/multiview/triangulation.h +++ b/extern/libmv/libmv/multiview/triangulation.h @@ -1,15 +1,15 @@ // 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 @@ -32,7 +32,7 @@ void TriangulateDLT(const Mat34 &P1, const Vec2 &x1, void TriangulateDLT(const Mat34 &P1, const Vec2 &x1, const Mat34 &P2, const Vec2 &x2, Vec3 *X_euclidean); - -} // namespace libmv + +} // namespace libmv #endif // LIBMV_MULTIVIEW_TRIANGULATION_H_ diff --git a/extern/libmv/libmv/numeric/dogleg.h b/extern/libmv/libmv/numeric/dogleg.h index f05882c1191..bf6f996ddaf 100644 --- a/extern/libmv/libmv/numeric/dogleg.h +++ b/extern/libmv/libmv/numeric/dogleg.h @@ -41,7 +41,7 @@ namespace libmv { template, typename Solver = Eigen::PartialPivLU< - Matrix > > class Dogleg { @@ -52,8 +52,8 @@ class Dogleg { typedef Matrix JMatrixType; - typedef Matrix AMatrixType; enum Status { @@ -75,7 +75,7 @@ class Dogleg { : f_(f), df_(f) {} struct SolverParameters { - SolverParameters() + SolverParameters() : gradient_threshold(1e-16), relative_step_threshold(1e-16), error_threshold(1e-16), diff --git a/extern/libmv/libmv/numeric/function_derivative.h b/extern/libmv/libmv/numeric/function_derivative.h index d7bc437b2e0..9820885f04e 100644 --- a/extern/libmv/libmv/numeric/function_derivative.h +++ b/extern/libmv/libmv/numeric/function_derivative.h @@ -36,7 +36,7 @@ enum NumericJacobianMode { FORWARD, }; -template +template class NumericJacobian { public: typedef typename Function::XMatrixType Parameters; @@ -58,7 +58,7 @@ class NumericJacobian { XScalar mean_eps = eps.sum() / eps.rows(); if (mean_eps == XScalar(0)) { // TODO(keir): Do something better here. - mean_eps = 1e-8; // ~sqrt(machine precision). + mean_eps = 1e-8; // ~sqrt(machine precision). } // TODO(keir): Elimininate this needless function evaluation for the // central difference case. diff --git a/extern/libmv/libmv/numeric/levenberg_marquardt.h b/extern/libmv/libmv/numeric/levenberg_marquardt.h index a7877e0270b..2af9a62cf7b 100644 --- a/extern/libmv/libmv/numeric/levenberg_marquardt.h +++ b/extern/libmv/libmv/numeric/levenberg_marquardt.h @@ -40,7 +40,7 @@ namespace libmv { template, typename Solver = Eigen::PartialPivLU< - Matrix > > class LevenbergMarquardt { @@ -51,8 +51,8 @@ class LevenbergMarquardt { typedef Matrix JMatrixType; - typedef Matrix AMatrixType; // TODO(keir): Some of these knobs can be derived from each other and @@ -69,7 +69,7 @@ class LevenbergMarquardt { : f_(f), df_(f) {} struct SolverParameters { - SolverParameters() + SolverParameters() : gradient_threshold(1e-16), relative_step_threshold(1e-16), error_threshold(1e-16), @@ -140,7 +140,7 @@ class LevenbergMarquardt { if (solved && dx.norm() <= params.relative_step_threshold * x.norm()) { results.status = RELATIVE_STEP_SIZE_TOO_SMALL; break; - } + } if (solved) { x_new = x + dx; // Rho is the ratio of the actual reduction in error to the reduction @@ -156,8 +156,8 @@ class LevenbergMarquardt { u = u*std::max(1/3., 1 - (tmp*tmp*tmp)); v = 2; continue; - } - } + } + } // Reject the update because either the normal equations failed to solve // or the local linear model was not good (rho < 0). Instead, increase u // to move closer to gradient descent. diff --git a/extern/libmv/libmv/numeric/numeric.cc b/extern/libmv/libmv/numeric/numeric.cc index 0ca3a64e4f4..9007663c8e2 100644 --- a/extern/libmv/libmv/numeric/numeric.cc +++ b/extern/libmv/libmv/numeric/numeric.cc @@ -121,15 +121,15 @@ void HorizontalStack(const Mat &left, const Mat &right, Mat *stacked) { void MatrixColumn(const Mat &A, int i, Vec2 *v) { assert(A.rows() == 2); - *v << A(0,i), A(1,i); + *v << A(0, i), A(1, i); } void MatrixColumn(const Mat &A, int i, Vec3 *v) { assert(A.rows() == 3); - *v << A(0,i), A(1,i), A(2,i); + *v << A(0, i), A(1, i), A(2, i); } void MatrixColumn(const Mat &A, int i, Vec4 *v) { assert(A.rows() == 4); - *v << A(0,i), A(1,i), A(2,i), A(3,i); + *v << A(0, i), A(1, i), A(2, i), A(3, i); } } // namespace libmv diff --git a/extern/libmv/libmv/numeric/numeric.h b/extern/libmv/libmv/numeric/numeric.h index d2cadf4b579..a95723d59cf 100644 --- a/extern/libmv/libmv/numeric/numeric.h +++ b/extern/libmv/libmv/numeric/numeric.h @@ -34,11 +34,11 @@ #include #if (defined(_WIN32) || defined(__APPLE__) || defined(__FreeBSD__)) && !defined(__MINGW64__) - void static sincos (double x, double *sinx, double *cosx) { + static void sincos(double x, double *sinx, double *cosx) { *sinx = sin(x); *cosx = cos(x); } -#endif //_WIN32 || __APPLE__ +#endif // _WIN32 || __APPLE__ #if (defined(WIN32) || defined(WIN64)) && !defined(__MINGW32__) inline long lround(double d) { @@ -48,7 +48,7 @@ return (d>0) ? int(d+0.5) : int(d-0.5); } typedef unsigned int uint; -#endif //_WIN32 +#endif // _WIN32 namespace libmv { @@ -86,16 +86,16 @@ typedef Eigen::Matrix RMat4; typedef Eigen::Matrix Mat2X; typedef Eigen::Matrix Mat3X; typedef Eigen::Matrix Mat4X; -typedef Eigen::Matrix MatX2; -typedef Eigen::Matrix MatX3; -typedef Eigen::Matrix MatX4; -typedef Eigen::Matrix MatX5; -typedef Eigen::Matrix MatX6; -typedef Eigen::Matrix MatX7; -typedef Eigen::Matrix MatX8; -typedef Eigen::Matrix MatX9; -typedef Eigen::Matrix MatX15; -typedef Eigen::Matrix MatX16; +typedef Eigen::Matrix MatX2; +typedef Eigen::Matrix MatX3; +typedef Eigen::Matrix MatX4; +typedef Eigen::Matrix MatX5; +typedef Eigen::Matrix MatX6; +typedef Eigen::Matrix MatX7; +typedef Eigen::Matrix MatX8; +typedef Eigen::Matrix MatX9; +typedef Eigen::Matrix MatX15; +typedef Eigen::Matrix MatX16; typedef Eigen::Vector2d Vec2; typedef Eigen::Vector3d Vec3; @@ -286,18 +286,21 @@ void MeanAndVarianceAlongRows(const Mat &A, #if _WIN32 // TODO(bomboze): un-#if this for both platforms once tested under Windows - /* This solution was extensively discussed here http://forum.kde.org/viewtopic.php?f=74&t=61940 */ - #define SUM_OR_DYNAMIC(x,y) (x==Eigen::Dynamic||y==Eigen::Dynamic)?Eigen::Dynamic:(x+y) + /* This solution was extensively discussed here + http://forum.kde.org/viewtopic.php?f=74&t=61940 */ + #define SUM_OR_DYNAMIC(x, y) (x == Eigen::Dynamic || y == Eigen::Dynamic) ? Eigen::Dynamic : (x+y) template struct hstack_return { typedef typename Derived1::Scalar Scalar; enum { RowsAtCompileTime = Derived1::RowsAtCompileTime, - ColsAtCompileTime = SUM_OR_DYNAMIC(Derived1::ColsAtCompileTime, Derived2::ColsAtCompileTime), + ColsAtCompileTime = SUM_OR_DYNAMIC(Derived1::ColsAtCompileTime, + Derived2::ColsAtCompileTime), Options = Derived1::Flags&Eigen::RowMajorBit ? Eigen::RowMajor : 0, MaxRowsAtCompileTime = Derived1::MaxRowsAtCompileTime, - MaxColsAtCompileTime = SUM_OR_DYNAMIC(Derived1::MaxColsAtCompileTime, Derived2::MaxColsAtCompileTime) + MaxColsAtCompileTime = SUM_OR_DYNAMIC(Derived1::MaxColsAtCompileTime, + Derived2::MaxColsAtCompileTime) }; typedef Eigen::Matrix - typename hstack_return::type - HStack (const Eigen::MatrixBase& lhs, const Eigen::MatrixBase& rhs) { - typename hstack_return::type res; + typename hstack_return::type + HStack(const Eigen::MatrixBase& lhs, + const Eigen::MatrixBase& rhs) { + typename hstack_return::type res; res.resize(lhs.rows(), lhs.cols()+rhs.cols()); res << lhs, rhs; return res; @@ -321,10 +325,12 @@ void MeanAndVarianceAlongRows(const Mat &A, struct vstack_return { typedef typename Derived1::Scalar Scalar; enum { - RowsAtCompileTime = SUM_OR_DYNAMIC(Derived1::RowsAtCompileTime, Derived2::RowsAtCompileTime), + RowsAtCompileTime = SUM_OR_DYNAMIC(Derived1::RowsAtCompileTime, + Derived2::RowsAtCompileTime), ColsAtCompileTime = Derived1::ColsAtCompileTime, Options = Derived1::Flags&Eigen::RowMajorBit ? Eigen::RowMajor : 0, - MaxRowsAtCompileTime = SUM_OR_DYNAMIC(Derived1::MaxRowsAtCompileTime, Derived2::MaxRowsAtCompileTime), + MaxRowsAtCompileTime = SUM_OR_DYNAMIC(Derived1::MaxRowsAtCompileTime, + Derived2::MaxRowsAtCompileTime), MaxColsAtCompileTime = Derived1::MaxColsAtCompileTime }; typedef Eigen::Matrix - typename vstack_return::type - VStack (const Eigen::MatrixBase& lhs, const Eigen::MatrixBase& rhs) { - typename vstack_return::type res; + typename vstack_return::type + VStack(const Eigen::MatrixBase& lhs, + const Eigen::MatrixBase& rhs) { + typename vstack_return::type res; res.resize(lhs.rows()+rhs.rows(), lhs.cols()); res << lhs, rhs; return res; }; -#else //_WIN32 +#else // _WIN32 // Since it is not possible to typedef privately here, use a macro. // Always take dynamic columns if either side is dynamic. @@ -400,7 +407,7 @@ void MeanAndVarianceAlongRows(const Mat &A, } #undef COLS #undef ROWS -#endif //_WIN32 +#endif // _WIN32 @@ -454,7 +461,7 @@ inline bool isnan(double i) { /// and negative values template FloatType ceil0(const FloatType& value) { - FloatType result = std::ceil( std::fabs( value ) ); + FloatType result = std::ceil(std::fabs(value)); return (value < 0.0) ? -result : result; } @@ -470,8 +477,8 @@ inline Mat3 SkewMat(const Vec3 &x) { /// the first two (independent) lines inline Mat23 SkewMatMinimal(const Vec2 &x) { Mat23 skew; - skew << 0,-1, x(1), - 1, 0, -x(0); + skew << 0, -1, x(1), + 1, 0, -x(0); return skew; } @@ -485,6 +492,6 @@ inline Mat3 RotationFromEulerVector(Vec3 euler_vector) { Mat3 w_hat = CrossProductMatrix(w); return Mat3::Identity() + w_hat*sin(theta) + w_hat*w_hat*(1 - cos(theta)); } -} // namespace libmv +} // namespace libmv #endif // LIBMV_NUMERIC_NUMERIC_H diff --git a/extern/libmv/libmv/numeric/poly.cc b/extern/libmv/libmv/numeric/poly.cc index d96e3c104c7..376403616c3 100644 --- a/extern/libmv/libmv/numeric/poly.cc +++ b/extern/libmv/libmv/numeric/poly.cc @@ -1,15 +1,15 @@ // 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 diff --git a/extern/libmv/libmv/numeric/poly.h b/extern/libmv/libmv/numeric/poly.h index cb1d65b32c4..76ba062d475 100644 --- a/extern/libmv/libmv/numeric/poly.h +++ b/extern/libmv/libmv/numeric/poly.h @@ -1,15 +1,15 @@ // 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 @@ -36,7 +36,7 @@ namespace libmv { // // The GSL cubic solver was used as a reference for this routine. template -int SolveCubicPolynomial(Real a, Real b, Real c, +int SolveCubicPolynomial(Real a, Real b, Real c, Real *x0, Real *x1, Real *x2) { Real q = a * a - 3 * b; Real r = 2 * a * a * a - 9 * a * b + 27 * c; @@ -52,7 +52,7 @@ int SolveCubicPolynomial(Real a, Real b, Real c, if (R == 0 && Q == 0) { // Tripple root in one place. - *x0 = *x1 = *x2 = -a / 3 ; + *x0 = *x1 = *x2 = -a / 3; return 3; } else if (CR2 == CQ3) { @@ -62,7 +62,7 @@ int SolveCubicPolynomial(Real a, Real b, Real c, // Due to finite precision some double roots may be missed, and considered // to be a pair of complex roots z = x +/- epsilon i close to the real // axis. - Real sqrtQ = sqrt (Q); + Real sqrtQ = sqrt(Q); if (R > 0) { *x0 = -2 * sqrtQ - a / 3; *x1 = sqrtQ - a / 3; @@ -76,13 +76,13 @@ int SolveCubicPolynomial(Real a, Real b, Real c, } else if (CR2 < CQ3) { // This case is equivalent to R2 < Q3. - Real sqrtQ = sqrt (Q); + Real sqrtQ = sqrt(Q); Real sqrtQ3 = sqrtQ * sqrtQ * sqrtQ; - Real theta = acos (R / sqrtQ3); + Real theta = acos(R / sqrtQ3); Real norm = -2 * sqrtQ; - *x0 = norm * cos (theta / 3) - a / 3; - *x1 = norm * cos ((theta + 2.0 * M_PI) / 3) - a / 3; - *x2 = norm * cos ((theta - 2.0 * M_PI) / 3) - a / 3; + *x0 = norm * cos(theta / 3) - a / 3; + *x1 = norm * cos((theta + 2.0 * M_PI) / 3) - a / 3; + *x2 = norm * cos((theta - 2.0 * M_PI) / 3) - a / 3; // Put the roots in ascending order. if (*x0 > *x1) { @@ -95,10 +95,10 @@ int SolveCubicPolynomial(Real a, Real b, Real c, } } return 3; - } + } Real sgnR = (R >= 0 ? 1 : -1); - Real A = -sgnR * pow (fabs (R) + sqrt (R2 - Q3), 1.0/3.0); - Real B = Q / A ; + Real A = -sgnR * pow(fabs(R) + sqrt(R2 - Q3), 1.0/3.0); + Real B = Q / A; *x0 = A + B - a / 3; return 1; } diff --git a/extern/libmv/libmv/simple_pipeline/bundle.cc b/extern/libmv/libmv/simple_pipeline/bundle.cc index 3e36a78c9df..d84f7eb8215 100644 --- a/extern/libmv/libmv/simple_pipeline/bundle.cc +++ b/extern/libmv/libmv/simple_pipeline/bundle.cc @@ -1,15 +1,15 @@ // Copyright (c) 2011, 2012, 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 @@ -18,6 +18,8 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/simple_pipeline/bundle.h" + #include #include "ceres/ceres.h" @@ -29,7 +31,6 @@ #include "libmv/multiview/projection.h" #include "libmv/numeric/numeric.h" #include "libmv/simple_pipeline/camera_intrinsics.h" -#include "libmv/simple_pipeline/bundle.h" #include "libmv/simple_pipeline/reconstruction.h" #include "libmv/simple_pipeline/tracks.h" @@ -125,8 +126,8 @@ struct OpenCVReprojectionError { // single parameter block. struct RotationMatrixPlus { template - bool operator()(const T* R_array, // Rotation 3x3 col-major. - const T* delta, // Angle-axis delta + bool operator()(const T* R_array, // Rotation 3x3 col-major. + const T* delta, // Angle-axis delta T* R_plus_delta_array) const { T angle_axis[3]; @@ -293,7 +294,7 @@ void EuclideanBundleCommonIntrinsics(const Tracks &tracks, } LG << "Number of residuals: " << num_residuals; - if(!num_residuals) { + if (!num_residuals) { LG << "Skipping running minimizer with zero residuals"; return; } diff --git a/extern/libmv/libmv/simple_pipeline/bundle.h b/extern/libmv/libmv/simple_pipeline/bundle.h index 573bcf4cc21..6235b46f1e7 100644 --- a/extern/libmv/libmv/simple_pipeline/bundle.h +++ b/extern/libmv/libmv/simple_pipeline/bundle.h @@ -1,15 +1,15 @@ // 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 @@ -89,14 +89,15 @@ enum BundleIntrinsics { BUNDLE_TANGENTIAL = 48, }; enum BundleConstraints { - BUNDLE_NO_CONSTRAINTS = 0, - BUNDLE_NO_TRANSLATION = 1, + BUNDLE_NO_CONSTRAINTS = 0, + BUNDLE_NO_TRANSLATION = 1, }; void EuclideanBundleCommonIntrinsics(const Tracks &tracks, int bundle_intrinsics, EuclideanReconstruction *reconstruction, CameraIntrinsics *intrinsics, - int bundle_constraints = BUNDLE_NO_CONSTRAINTS); + int bundle_constraints = + BUNDLE_NO_CONSTRAINTS); /*! Refine camera poses and 3D coordinates using bundle adjustment. diff --git a/extern/libmv/libmv/simple_pipeline/callbacks.cc b/extern/libmv/libmv/simple_pipeline/callbacks.cc deleted file mode 100644 index 7e4bca6efc7..00000000000 --- a/extern/libmv/libmv/simple_pipeline/callbacks.cc +++ /dev/null @@ -1,29 +0,0 @@ -// Copyright (c) 2011 libmv authors. -// -// Permission is hereby granted, free of charge, to any person obtaining a copy -// of this software and associated documentation files (the "Software"), to -// deal in the Software without restriction, including without limitation the -// rights to use, copy, modify, merge, publish, distribute, sublicense, and/or -// sell copies of the Software, and to permit persons to whom the Software is -// furnished to do so, subject to the following conditions: -// -// The above copyright notice and this permission notice shall be included in -// all copies or substantial portions of the Software. -// -// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING -// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS -// IN THE SOFTWARE. - -#include "libmv/simple_pipeline/callbacks.h" - -namespace libmv { - -void ProgressUpdateCallback::invoke(double progress, const char* message) -{ -} - -} // namespace libmv diff --git a/extern/libmv/libmv/simple_pipeline/callbacks.h b/extern/libmv/libmv/simple_pipeline/callbacks.h index 675f73c1cab..58f7b0d3cc9 100644 --- a/extern/libmv/libmv/simple_pipeline/callbacks.h +++ b/extern/libmv/libmv/simple_pipeline/callbacks.h @@ -25,7 +25,8 @@ namespace libmv { class ProgressUpdateCallback { public: - virtual void invoke(double progress, const char *message); + virtual ~ProgressUpdateCallback() {} + virtual void invoke(double progress, const char *message) = 0; }; } // namespace libmv diff --git a/extern/libmv/libmv/simple_pipeline/camera_intrinsics.cc b/extern/libmv/libmv/simple_pipeline/camera_intrinsics.cc index 543fddedd0c..aee2ecb5882 100644 --- a/extern/libmv/libmv/simple_pipeline/camera_intrinsics.cc +++ b/extern/libmv/libmv/simple_pipeline/camera_intrinsics.cc @@ -25,7 +25,7 @@ namespace libmv { struct Offset { short ix, iy; - unsigned char fx,fy; + unsigned char fx, fy; }; struct Grid { @@ -34,8 +34,7 @@ struct Grid { double overscan; }; -static struct Grid *copyGrid(struct Grid *from) -{ +static struct Grid *copyGrid(struct Grid *from) { struct Grid *to = NULL; if (from) { @@ -74,8 +73,7 @@ CameraIntrinsics::CameraIntrinsics(const CameraIntrinsics &from) k3_(from.k3_), p1_(from.p1_), p2_(from.p2_), - threads_(from.threads_) -{ + threads_(from.threads_) { distort_ = copyGrid(from.distort_); undistort_ = copyGrid(from.undistort_); } @@ -122,9 +120,8 @@ void CameraIntrinsics::SetTangentialDistortion(double p1, double p2) { FreeLookupGrid(); } -void CameraIntrinsics::SetThreads(int threads) -{ - threads_ = threads; +void CameraIntrinsics::SetThreads(int threads) { + threads_ = threads; } void CameraIntrinsics::ApplyIntrinsics(double normalized_x, @@ -189,68 +186,76 @@ void CameraIntrinsics::InvertIntrinsics(double image_x, // TODO(MatthiasF): downsample lookup template -void CameraIntrinsics::ComputeLookupGrid(Grid* grid, int width, int height, double overscan) { +void CameraIntrinsics::ComputeLookupGrid(Grid* grid, int width, int height, + double overscan) { double w = (double)width / (1 + overscan); double h = (double)height / (1 + overscan); double aspx = (double)w / image_width_; double aspy = (double)h / image_height_; - +#if defined(_OPENMP) #pragma omp parallel for schedule(dynamic) num_threads(threads_) if (threads_ > 1 && height > 100) +#endif for (int y = 0; y < height; y++) { for (int x = 0; x < width; x++) { - double src_x = (x - 0.5 * overscan * w) / aspx, src_y = (y - 0.5 * overscan * h) / aspy; + double src_x = (x - 0.5 * overscan * w) / aspx, + src_y = (y - 0.5 * overscan * h) / aspy; double warp_x, warp_y; - WarpFunction(this,src_x,src_y,&warp_x,&warp_y); + WarpFunction(this, src_x, src_y, &warp_x, &warp_y); warp_x = warp_x*aspx + 0.5 * overscan * w; warp_y = warp_y*aspy + 0.5 * overscan * h; int ix = int(warp_x), iy = int(warp_y); int fx = round((warp_x-ix)*256), fy = round((warp_y-iy)*256); - if(fx == 256) { fx=0; ix++; } - if(fy == 256) { fy=0; iy++; } + if (fx == 256) { fx = 0; ix++; } // NOLINT + if (fy == 256) { fy = 0; iy++; } // NOLINT // Use nearest border pixel - if( ix < 0 ) { ix = 0, fx = 0; } - if( iy < 0 ) { iy = 0, fy = 0; } - if( ix >= width-2 ) ix = width-2; - if( iy >= height-2 ) iy = height-2; + if (ix < 0) { ix = 0, fx = 0; } // NOLINT + if (iy < 0) { iy = 0, fy = 0; } // NOLINT + if (ix >= width - 2) ix = width-2; + if (iy >= height - 2) iy = height-2; - Offset offset = { (short)(ix-x), (short)(iy-y), (unsigned char)fx, (unsigned char)fy }; + Offset offset = { (short)(ix-x), (short)(iy-y), + (unsigned char)fx, (unsigned char)fy }; grid->offset[y*width+x] = offset; } } } // TODO(MatthiasF): cubic B-Spline image sampling, bilinear lookup -template +template static void Warp(const Grid* grid, const T* src, T* dst, int width, int height, int threads) { + (void) threads; // Ignored if OpenMP is disabled +#if defined(_OPENMP) #pragma omp parallel for schedule(dynamic) num_threads(threads) if (threads > 1 && height > 100) +#endif for (int y = 0; y < height; y++) { for (int x = 0; x < width; x++) { Offset offset = grid->offset[y*width+x]; const T* s = &src[((y+offset.iy)*width+(x+offset.ix))*N]; for (int i = 0; i < N; i++) { - dst[(y*width+x)*N+i] = ((s[ i] * (256-offset.fx) + s[ N+i] * offset.fx) * (256-offset.fy) - +(s[width*N+i] * (256-offset.fx) + s[width*N+N+i] * offset.fx) * offset.fy) / (256*256); + dst[(y*width+x)*N+i] = ((s[ i] * (256-offset.fx) + s[ N+i] * offset.fx) * (256-offset.fy) // NOLINT + +(s[width*N+i] * (256-offset.fx) + s[width*N+N+i] * offset.fx) * offset.fy) / (256*256); // NOLINT } } } } void CameraIntrinsics::FreeLookupGrid() { - if(distort_) { + if (distort_) { delete distort_->offset; delete distort_; distort_ = NULL; } - if(undistort_) { + if (undistort_) { delete undistort_->offset; delete undistort_; undistort_ = NULL; } } -// FIXME: C++ templates limitations makes thing complicated, but maybe there is a simpler method. +// FIXME: C++ templates limitations makes thing complicated, +// but maybe there is a simpler method. struct ApplyIntrinsicsFunction { ApplyIntrinsicsFunction(CameraIntrinsics* intrinsics, double x, double y, double *warp_x, double *warp_y) { @@ -263,16 +268,21 @@ struct ApplyIntrinsicsFunction { struct InvertIntrinsicsFunction { InvertIntrinsicsFunction(CameraIntrinsics* intrinsics, double x, double y, double *warp_x, double *warp_y) { - intrinsics->InvertIntrinsics(x,y,warp_x,warp_y); - *warp_x = *warp_x*intrinsics->focal_length_x()+intrinsics->principal_point_x(); - *warp_y = *warp_y*intrinsics->focal_length_y()+intrinsics->principal_point_y(); + intrinsics->InvertIntrinsics(x, y, warp_x, warp_y); + + *warp_x = *warp_x * intrinsics->focal_length_x() + + intrinsics->principal_point_x(); + + *warp_y = *warp_y * intrinsics->focal_length_y() + + intrinsics->principal_point_y(); } }; -void CameraIntrinsics::CheckDistortLookupGrid(int width, int height, double overscan) -{ - if(distort_) { - if(distort_->width != width || distort_->height != height || distort_->overscan != overscan) { +void CameraIntrinsics::CheckDistortLookupGrid(int width, int height, + double overscan) { + if (distort_) { + if (distort_->width != width || distort_->height != height || + distort_->overscan != overscan) { delete [] distort_->offset; distort_->offset = NULL; } @@ -281,9 +291,10 @@ void CameraIntrinsics::CheckDistortLookupGrid(int width, int height, double over distort_->offset = NULL; } - if(!distort_->offset) { - distort_->offset = new Offset[width*height]; - ComputeLookupGrid(distort_,width,height,overscan); + if (!distort_->offset) { + distort_->offset = new Offset[width * height]; + ComputeLookupGrid(distort_, width, + height, overscan); } distort_->width = width; @@ -291,10 +302,11 @@ void CameraIntrinsics::CheckDistortLookupGrid(int width, int height, double over distort_->overscan = overscan; } -void CameraIntrinsics::CheckUndistortLookupGrid(int width, int height, double overscan) -{ - if(undistort_) { - if(undistort_->width != width || undistort_->height != height || undistort_->overscan != overscan) { +void CameraIntrinsics::CheckUndistortLookupGrid(int width, int height, + double overscan) { + if (undistort_) { + if (undistort_->width != width || undistort_->height != height || + undistort_->overscan != overscan) { delete [] undistort_->offset; undistort_->offset = NULL; } @@ -303,9 +315,10 @@ void CameraIntrinsics::CheckUndistortLookupGrid(int width, int height, double ov undistort_->offset = NULL; } - if(!undistort_->offset) { - undistort_->offset = new Offset[width*height]; - ComputeLookupGrid(undistort_,width,height,overscan); + if (!undistort_->offset) { + undistort_->offset = new Offset[width * height]; + ComputeLookupGrid(undistort_, width, + height, overscan); } undistort_->width = width; @@ -313,39 +326,53 @@ void CameraIntrinsics::CheckUndistortLookupGrid(int width, int height, double ov undistort_->overscan = overscan; } -void CameraIntrinsics::Distort(const float* src, float* dst, int width, int height, double overscan, int channels) { +void CameraIntrinsics::Distort(const float* src, float* dst, + int width, int height, + double overscan, + int channels) { CheckDistortLookupGrid(width, height, overscan); - if(channels==1) Warp(distort_,src,dst,width,height,threads_); - else if(channels==2) Warp(distort_,src,dst,width,height,threads_); - else if(channels==3) Warp(distort_,src,dst,width,height,threads_); - else if(channels==4) Warp(distort_,src,dst,width,height,threads_); + if (channels==1) Warp(distort_, src, dst, width, height, threads_); // NOLINT + else if (channels==2) Warp(distort_, src, dst, width, height, threads_); // NOLINT + else if (channels==3) Warp(distort_, src, dst, width, height, threads_); // NOLINT + else if (channels==4) Warp(distort_, src, dst, width, height, threads_); // NOLINT //else assert("channels must be between 1 and 4"); } -void CameraIntrinsics::Distort(const unsigned char* src, unsigned char* dst, int width, int height, double overscan, int channels) { +void CameraIntrinsics::Distort(const unsigned char* src, + unsigned char* dst, + int width, int height, + double overscan, + int channels) { CheckDistortLookupGrid(width, height, overscan); - if(channels==1) Warp(distort_,src,dst,width,height,threads_); - else if(channels==2) Warp(distort_,src,dst,width,height,threads_); - else if(channels==3) Warp(distort_,src,dst,width,height,threads_); - else if(channels==4) Warp(distort_,src,dst,width,height,threads_); + if (channels == 1) Warp(distort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 2) Warp(distort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 3) Warp(distort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 4) Warp(distort_, src, dst, width, height, threads_); // NOLINT //else assert("channels must be between 1 and 4"); } -void CameraIntrinsics::Undistort(const float* src, float* dst, int width, int height, double overscan, int channels) { +void CameraIntrinsics::Undistort(const float* src, float* dst, + int width, int height, + double overscan, + int channels) { CheckUndistortLookupGrid(width, height, overscan); - if(channels==1) Warp(undistort_,src,dst,width,height,threads_); - else if(channels==2) Warp(undistort_,src,dst,width,height,threads_); - else if(channels==3) Warp(undistort_,src,dst,width,height,threads_); - else if(channels==4) Warp(undistort_,src,dst,width,height,threads_); + if (channels == 1) Warp(undistort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 2) Warp(undistort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 3) Warp(undistort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 4) Warp(undistort_, src, dst, width, height, threads_); // NOLINT //else assert("channels must be between 1 and 4"); } -void CameraIntrinsics::Undistort(const unsigned char* src, unsigned char* dst, int width, int height, double overscan, int channels) { +void CameraIntrinsics::Undistort(const unsigned char* src, + unsigned char* dst, + int width, int height, + double overscan, + int channels) { CheckUndistortLookupGrid(width, height, overscan); - if(channels==1) Warp(undistort_,src,dst,width,height,threads_); - else if(channels==2) Warp(undistort_,src,dst,width,height,threads_); - else if(channels==3) Warp(undistort_,src,dst,width,height,threads_); - else if(channels==4) Warp(undistort_,src,dst,width,height,threads_); + if (channels == 1) Warp(undistort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 2) Warp(undistort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 3) Warp(undistort_, src, dst, width, height, threads_); // NOLINT + else if (channels == 4) Warp(undistort_, src, dst, width, height, threads_); // NOLINT //else assert("channels must be between 1 and 4"); } diff --git a/extern/libmv/libmv/simple_pipeline/camera_intrinsics.h b/extern/libmv/libmv/simple_pipeline/camera_intrinsics.h index 632922c38ff..29d26c7f9ab 100644 --- a/extern/libmv/libmv/simple_pipeline/camera_intrinsics.h +++ b/extern/libmv/libmv/simple_pipeline/camera_intrinsics.h @@ -135,7 +135,10 @@ class CameraIntrinsics { int width, int height, double overscan, int channels); private: - template void ComputeLookupGrid(struct Grid* grid, int width, int height, double overscan); + template void ComputeLookupGrid(struct Grid* grid, + int width, + int height, + double overscan); void CheckUndistortLookupGrid(int width, int height, double overscan); void CheckDistortLookupGrid(int width, int height, double overscan); void FreeLookupGrid(); diff --git a/extern/libmv/libmv/simple_pipeline/detect.cc b/extern/libmv/libmv/simple_pipeline/detect.cc index 9e3edf32d71..8451e47d374 100644 --- a/extern/libmv/libmv/simple_pipeline/detect.cc +++ b/extern/libmv/libmv/simple_pipeline/detect.cc @@ -35,37 +35,40 @@ namespace libmv { typedef unsigned int uint; -static int featurecmp(const void *a_v, const void *b_v) -{ +static int featurecmp(const void *a_v, const void *b_v) { Feature *a = (Feature*)a_v; Feature *b = (Feature*)b_v; return b->score - a->score; } -std::vector DetectFAST(const unsigned char* data, int width, int height, int stride, - int min_trackness, int min_distance) { +std::vector DetectFAST(const unsigned char* data, + int width, int height, + int stride, + int min_trackness, + int min_distance) { std::vector features; // TODO(MatthiasF): Support targetting a feature count (binary search trackness) int num_features; xy* all = fast9_detect(data, width, height, stride, min_trackness, &num_features); - if(num_features == 0) { + if (num_features == 0) { free(all); return features; } int* scores = fast9_score(data, stride, all, num_features, min_trackness); - // TODO: merge with close feature suppression + // TODO(MatthiasF): merge with close feature suppression xy* nonmax = nonmax_suppression(all, scores, num_features, &num_features); free(all); // Remove too close features - // TODO(MatthiasF): A resolution independent parameter would be better than distance - // e.g. a coefficient going from 0 (no minimal distance) to 1 (optimal circle packing) + // TODO(MatthiasF): A resolution independent parameter would be better than + // distance e.g. a coefficient going from 0 (no minimal distance) to 1 + // (optimal circle packing) // FIXME(MatthiasF): this method will not necessarily give all maximum markers - if(num_features) { + if (num_features) { Feature *all_features = new Feature[num_features]; - for(int i = 0; i < num_features; ++i) { + for (int i = 0; i < num_features; ++i) { Feature a = { (float)nonmax[i].x, (float)nonmax[i].y, (float)scores[i], 0 }; all_features[i] = a; } @@ -75,15 +78,15 @@ std::vector DetectFAST(const unsigned char* data, int width, int height features.reserve(num_features); int prev_score = all_features[0].score; - for(int i = 0; i < num_features; ++i) { + for (int i = 0; i < num_features; ++i) { bool ok = true; Feature a = all_features[i]; - if(a.score>prev_score) + if (a.score>prev_score) abort(); prev_score = a.score; // compare each feature against filtered set - for(int j = 0; j < features.size(); j++) { + for (int j = 0; j < features.size(); j++) { Feature& b = features[j]; if ( (a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y) < min_distance*min_distance ) { // already a nearby feature @@ -92,7 +95,7 @@ std::vector DetectFAST(const unsigned char* data, int width, int height } } - if(ok) { + if (ok) { // add the new feature features.push_back(a); } @@ -106,19 +109,22 @@ std::vector DetectFAST(const unsigned char* data, int width, int height } #ifdef __SSE2__ -static uint SAD(const ubyte* imageA, const ubyte* imageB, int strideA, int strideB) { +static uint SAD(const ubyte* imageA, const ubyte* imageB, + int strideA, int strideB) { __m128i a = _mm_setzero_si128(); - for(int i = 0; i < 16; i++) { - a = _mm_adds_epu16(a, _mm_sad_epu8( _mm_loadu_si128((__m128i*)(imageA+i*strideA)), - _mm_loadu_si128((__m128i*)(imageB+i*strideB)))); + for (int i = 0; i < 16; i++) { + a = _mm_adds_epu16(a, + _mm_sad_epu8(_mm_loadu_si128((__m128i*)(imageA+i*strideA)), + _mm_loadu_si128((__m128i*)(imageB+i*strideB)))); } - return _mm_extract_epi16(a,0) + _mm_extract_epi16(a,4); + return _mm_extract_epi16(a, 0) + _mm_extract_epi16(a, 4); } #else -static uint SAD(const ubyte* imageA, const ubyte* imageB, int strideA, int strideB) { - uint sad=0; - for(int i = 0; i < 16; i++) { - for(int j = 0; j < 16; j++) { +static uint SAD(const ubyte* imageA, const ubyte* imageB, + int strideA, int strideB) { + uint sad = 0; + for (int i = 0; i < 16; i++) { + for (int j = 0; j < 16; j++) { sad += abs((int)imageA[i*strideA+j] - imageB[i*strideB+j]); } } @@ -126,59 +132,68 @@ static uint SAD(const ubyte* imageA, const ubyte* imageB, int strideA, int strid } #endif -void DetectMORAVEC(ubyte* image, int stride, int width, int height, Feature* detected, int* count, int distance, ubyte* pattern) { +void DetectMORAVEC(ubyte* image, + int stride, int width, int height, + Feature* detected, int* count, + int distance, + ubyte* pattern) { unsigned short histogram[256]; - memset(histogram,0,sizeof(histogram)); + memset(histogram, 0, sizeof(histogram)); ubyte* scores = new ubyte[width*height]; - memset(scores,0,width*height); - const int r = 1; //radius for self similarity comparison - for(int y=distance; y255) score=255; // clip + score /= 256; // normalize + if (pattern) // find only features similar to pattern + score -= SAD(s, pattern, stride, 16); + if (score <= 16) continue; // filter very self-similar features + score -= 16; // translate to score/histogram values + if (score>255) score=255; // clip ubyte* c = &scores[y*width+x]; - for(int i=-distance; i<0; i++) { - for(int j=-distance; j= score) goto nonmax; - c[i*width+j]=0, histogram[s]--; + if (s == 0) continue; + if (s >= score) goto nonmax; + c[i*width+j] = 0; + histogram[s]--; } } - for(int i=0, j=-distance; j<0; j++) { + for (int i = 0, j = -distance; j < 0; j++) { int s = c[i*width+j]; - if(s == 0) continue; - if(s >= score) goto nonmax; - c[i*width+j]=0, histogram[s]--; + if (s == 0) continue; + if (s >= score) goto nonmax; + c[i*width+j] = 0; + histogram[s]--; } c[0] = score, histogram[score]++; - nonmax:; + nonmax: + { } // Do nothing. } } - int min=255, total=0; - for(; min>0; min--) { + int min = 255, total = 0; + for (; min > 0; min--) { int h = histogram[min]; - if(total+h > *count) break; + if (total+h > *count) break; total += h; } - int i=0; - for(int y=16; ymin) detected[i++] = f; + if (s > min) detected[i++] = f; } } *count = i; delete[] scores; } -} +} // namespace libmv diff --git a/extern/libmv/libmv/simple_pipeline/detect.h b/extern/libmv/libmv/simple_pipeline/detect.h index bbe7aed784c..6d0f55a082f 100644 --- a/extern/libmv/libmv/simple_pipeline/detect.h +++ b/extern/libmv/libmv/simple_pipeline/detect.h @@ -88,8 +88,10 @@ std::vector DetectFAST(const unsigned char* data, int width, int height \note \a You can crop the image (to avoid detecting markers near the borders) without copying: image += marginY*stride+marginX, width -= 2*marginX, height -= 2*marginY; */ -void DetectMORAVEC(ubyte* image, int stride, int width, int height, Feature* detected, int* count, int distance /*=32*/, ubyte* pattern /*=0*/); +void DetectMORAVEC(ubyte* image, int stride, int width, int height, + Feature* detected, int* count, int distance /*=32*/, + ubyte* pattern /*=0*/); -} +} // namespace libmv #endif diff --git a/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.cc b/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.cc index 84d143b77d6..246828a644f 100644 --- a/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.cc +++ b/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.cc @@ -1,15 +1,15 @@ // 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 @@ -18,13 +18,14 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/simple_pipeline/initialize_reconstruction.h" + #include "libmv/base/vector.h" #include "libmv/logging/logging.h" #include "libmv/multiview/fundamental.h" #include "libmv/multiview/projection.h" #include "libmv/numeric/levenberg_marquardt.h" #include "libmv/numeric/numeric.h" -#include "libmv/simple_pipeline/initialize_reconstruction.h" #include "libmv/simple_pipeline/reconstruction.h" #include "libmv/simple_pipeline/tracks.h" @@ -109,7 +110,8 @@ namespace { Mat3 DecodeF(const Vec9 &encoded_F) { // Decode F and force it to be rank 2. Map full_rank_F(encoded_F.data(), 3, 3); - Eigen::JacobiSVD svd(full_rank_F, Eigen::ComputeFullU | Eigen::ComputeFullV); + Eigen::JacobiSVD svd(full_rank_F, + Eigen::ComputeFullU | Eigen::ComputeFullV); Vec3 diagonal = svd.singularValues(); diagonal(2) = 0; Mat3 F = svd.matrixU() * diagonal.asDiagonal() * svd.matrixV().transpose(); @@ -130,7 +132,7 @@ struct FundamentalSampsonCostFunction { Vec operator()(const Vec9 &encoded_F) const { // Decode F and force it to be rank 2. Mat3 F = DecodeF(encoded_F); - + Vec residuals(markers.size() / 2); residuals.setZero(); for (int i = 0; i < markers.size() / 2; ++i) { diff --git a/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.h b/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.h index f512c9a3439..744436246b0 100644 --- a/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.h +++ b/extern/libmv/libmv/simple_pipeline/initialize_reconstruction.h @@ -1,15 +1,15 @@ // 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 diff --git a/extern/libmv/libmv/simple_pipeline/intersect.cc b/extern/libmv/libmv/simple_pipeline/intersect.cc index d0f139a991a..d2dde4b36af 100644 --- a/extern/libmv/libmv/simple_pipeline/intersect.cc +++ b/extern/libmv/libmv/simple_pipeline/intersect.cc @@ -1,15 +1,15 @@ // 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 @@ -25,7 +25,6 @@ #include "libmv/multiview/projection.h" #include "libmv/multiview/triangulation.h" #include "libmv/multiview/nviewtriangulation.h" -#include "libmv/multiview/projection.h" #include "libmv/numeric/numeric.h" #include "libmv/numeric/levenberg_marquardt.h" #include "libmv/simple_pipeline/reconstruction.h" diff --git a/extern/libmv/libmv/simple_pipeline/intersect.h b/extern/libmv/libmv/simple_pipeline/intersect.h index edbf4a0335b..3a0ffa7418b 100644 --- a/extern/libmv/libmv/simple_pipeline/intersect.h +++ b/extern/libmv/libmv/simple_pipeline/intersect.h @@ -1,15 +1,15 @@ // 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 diff --git a/extern/libmv/libmv/simple_pipeline/modal_solver.cc b/extern/libmv/libmv/simple_pipeline/modal_solver.cc index d79c71508cc..90dfde15660 100644 --- a/extern/libmv/libmv/simple_pipeline/modal_solver.cc +++ b/extern/libmv/libmv/simple_pipeline/modal_solver.cc @@ -34,7 +34,7 @@ namespace libmv { namespace { -void ProjectMarkerOnSphere(Marker &marker, Vec3 &X) { +void ProjectMarkerOnSphere(const Marker &marker, Vec3 &X) { X(0) = marker.x; X(1) = marker.y; X(2) = 1.0; @@ -43,8 +43,7 @@ void ProjectMarkerOnSphere(Marker &marker, Vec3 &X) { } void ModalSolverLogProress(ProgressUpdateCallback *update_callback, - double progress) -{ + double progress) { if (update_callback) { char message[256]; @@ -56,13 +55,14 @@ void ModalSolverLogProress(ProgressUpdateCallback *update_callback, } struct ModalReprojectionError { - ModalReprojectionError(double observed_x, double observed_y, Vec3 &bundle) + ModalReprojectionError(double observed_x, + double observed_y, + const Vec3 &bundle) : observed_x(observed_x), observed_y(observed_y), bundle(bundle) { } template bool operator()(const T* quaternion, // Rotation quaternion T* residuals) const { - T R[9]; ceres::QuaternionToRotation(quaternion, R); @@ -96,7 +96,7 @@ struct ModalReprojectionError { }; } // namespace -void ModalSolver(Tracks &tracks, +void ModalSolver(const Tracks &tracks, EuclideanReconstruction *reconstruction, ProgressUpdateCallback *update_callback) { int max_image = tracks.MaxImage(); diff --git a/extern/libmv/libmv/simple_pipeline/modal_solver.h b/extern/libmv/libmv/simple_pipeline/modal_solver.h index 560b37c2987..9801fd21d81 100644 --- a/extern/libmv/libmv/simple_pipeline/modal_solver.h +++ b/extern/libmv/libmv/simple_pipeline/modal_solver.h @@ -39,7 +39,7 @@ namespace libmv { Reconstructed cameras and projected bundles would be added to reconstruction object. */ -void ModalSolver(Tracks &tracks, +void ModalSolver(const Tracks &tracks, EuclideanReconstruction *reconstruction, ProgressUpdateCallback *update_callback = NULL); diff --git a/extern/libmv/libmv/simple_pipeline/pipeline.cc b/extern/libmv/libmv/simple_pipeline/pipeline.cc index 463738e42bb..f6013d71867 100644 --- a/extern/libmv/libmv/simple_pipeline/pipeline.cc +++ b/extern/libmv/libmv/simple_pipeline/pipeline.cc @@ -18,10 +18,11 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/simple_pipeline/pipeline.h" + #include #include "libmv/logging/logging.h" -#include "libmv/simple_pipeline/pipeline.h" #include "libmv/simple_pipeline/bundle.h" #include "libmv/simple_pipeline/intersect.h" #include "libmv/simple_pipeline/resect.h" @@ -92,6 +93,9 @@ struct ProjectivePipelineRoutines { static bool Resect(const ReconstructionOptions &options, const vector &markers, ProjectiveReconstruction *reconstruction, bool final_pass) { + (void) options; // Ignored. + (void) final_pass; // Ignored. + return ProjectiveResect(markers, reconstruction); } @@ -120,17 +124,19 @@ struct ProjectivePipelineRoutines { } // namespace -static void CompleteReconstructionLogProress(ProgressUpdateCallback *update_callback, +static void CompleteReconstructionLogProress( + ProgressUpdateCallback *update_callback, double progress, - const char *step = NULL) -{ - if(update_callback) { + const char *step = NULL) { + if (update_callback) { char message[256]; - if(step) - snprintf(message, sizeof(message), "Completing solution %d%% | %s", (int)(progress*100), step); + if (step) + snprintf(message, sizeof(message), "Completing solution %d%% | %s", + (int)(progress*100), step); else - snprintf(message, sizeof(message), "Completing solution %d%%", (int)(progress*100)); + snprintf(message, sizeof(message), "Completing solution %d%%", + (int)(progress*100)); update_callback->invoke(progress, message); } @@ -207,7 +213,8 @@ void InternalCompleteReconstruction( if (reconstructed_markers.size() >= 5) { CompleteReconstructionLogProress(update_callback, (double)tot_resects/(max_image)); - if (PipelineRoutines::Resect(options, reconstructed_markers, reconstruction, false)) { + if (PipelineRoutines::Resect(options, reconstructed_markers, + reconstruction, false)) { num_resects++; tot_resects++; LG << "Ran Resect() for image " << image; @@ -243,7 +250,8 @@ void InternalCompleteReconstruction( if (reconstructed_markers.size() >= 5) { CompleteReconstructionLogProress(update_callback, (double)tot_resects/(max_image)); - if (PipelineRoutines::Resect(options, reconstructed_markers, reconstruction, true)) { + if (PipelineRoutines::Resect(options, reconstructed_markers, + reconstruction, true)) { num_resects++; LG << "Ran final Resect() for image " << image; } else { @@ -260,9 +268,10 @@ void InternalCompleteReconstruction( } template -double InternalReprojectionError(const Tracks &image_tracks, - const typename PipelineRoutines::Reconstruction &reconstruction, - const CameraIntrinsics &intrinsics) { +double InternalReprojectionError( + const Tracks &image_tracks, + const typename PipelineRoutines::Reconstruction &reconstruction, + const CameraIntrinsics &intrinsics) { int num_skipped = 0; int num_reprojected = 0; double total_error = 0.0; diff --git a/extern/libmv/libmv/simple_pipeline/pipeline.h b/extern/libmv/libmv/simple_pipeline/pipeline.h index 11c11297d78..d8489012b95 100644 --- a/extern/libmv/libmv/simple_pipeline/pipeline.h +++ b/extern/libmv/libmv/simple_pipeline/pipeline.h @@ -49,10 +49,11 @@ namespace libmv { \sa EuclideanResect, EuclideanIntersect, EuclideanBundle */ -void EuclideanCompleteReconstruction(const ReconstructionOptions &options, - const Tracks &tracks, - EuclideanReconstruction *reconstruction, - ProgressUpdateCallback *update_callback = NULL); +void EuclideanCompleteReconstruction( + const ReconstructionOptions &options, + const Tracks &tracks, + EuclideanReconstruction *reconstruction, + ProgressUpdateCallback *update_callback = NULL); /*! Estimate camera matrices and homogeneous 3D coordinates for all frames and diff --git a/extern/libmv/libmv/simple_pipeline/reconstruction.h b/extern/libmv/libmv/simple_pipeline/reconstruction.h index 71789e3a245..1610029b9d2 100644 --- a/extern/libmv/libmv/simple_pipeline/reconstruction.h +++ b/extern/libmv/libmv/simple_pipeline/reconstruction.h @@ -27,14 +27,14 @@ namespace libmv { struct ReconstructionOptions { - // threshold value of reconstruction error which is still considered successful - // if reconstruction error bigger than this value, fallback reconstruction - // algorithm would be used (if enabled) - double success_threshold; - - // use fallback reconstruction algorithm in cases main reconstruction algorithm - // failed to reconstruct - bool use_fallback_reconstruction; + // threshold value of reconstruction error which is still considered successful + // if reconstruction error bigger than this value, fallback reconstruction + // algorithm would be used (if enabled) + double success_threshold; + + // use fallback reconstruction algorithm in cases main reconstruction algorithm + // failed to reconstruct + bool use_fallback_reconstruction; }; /*! diff --git a/extern/libmv/libmv/simple_pipeline/resect.cc b/extern/libmv/libmv/simple_pipeline/resect.cc index 3929271e66f..941c95cf237 100644 --- a/extern/libmv/libmv/simple_pipeline/resect.cc +++ b/extern/libmv/libmv/simple_pipeline/resect.cc @@ -1,15 +1,15 @@ // 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 @@ -18,6 +18,8 @@ // 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 #include "libmv/base/vector.h" @@ -27,7 +29,6 @@ #include "libmv/multiview/projection.h" #include "libmv/numeric/numeric.h" #include "libmv/numeric/levenberg_marquardt.h" -#include "libmv/simple_pipeline/resect.h" #include "libmv/simple_pipeline/reconstruction.h" #include "libmv/simple_pipeline/tracks.h" @@ -108,13 +109,13 @@ bool EuclideanResect(const ReconstructionOptions &options, double success_threshold = std::numeric_limits::max(); - if(options.use_fallback_reconstruction) + if (options.use_fallback_reconstruction) success_threshold = options.success_threshold; - if (0 || !euclidean_resection::EuclideanResection(points_2d, points_3d, &R, &t, - euclidean_resection::RESECTION_EPNP, - success_threshold)) - { + if (0 || !euclidean_resection::EuclideanResection( + points_2d, points_3d, &R, &t, + euclidean_resection::RESECTION_EPNP, + success_threshold)) { // printf("Resection for image %d failed\n", markers[0].image); LG << "Resection for image " << markers[0].image << " failed;" << " trying fallback projective resection."; @@ -172,7 +173,7 @@ bool EuclideanResect(const ReconstructionOptions &options, Solver solver(resect_cost); Solver::SolverParameters params; - Solver::Results results = solver.minimize(params, &dRt); + /* Solver::Results results = */ solver.minimize(params, &dRt); LG << "LM found incremental rotation: " << dRt.head<3>().transpose(); // TODO(keir): Check results to ensure clean termination. @@ -264,7 +265,7 @@ bool ProjectiveResect(const vector &markers, Solver solver(resect_cost); Solver::SolverParameters params; - Solver::Results results = solver.minimize(params, &vector_P); + /* Solver::Results results = */ solver.minimize(params, &vector_P); // TODO(keir): Check results to ensure clean termination. // Unpack the projection matrix. diff --git a/extern/libmv/libmv/simple_pipeline/resect.h b/extern/libmv/libmv/simple_pipeline/resect.h index 1691e7ee245..47a6c6b60ea 100644 --- a/extern/libmv/libmv/simple_pipeline/resect.h +++ b/extern/libmv/libmv/simple_pipeline/resect.h @@ -1,15 +1,15 @@ // 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 diff --git a/extern/libmv/libmv/simple_pipeline/tracks.cc b/extern/libmv/libmv/simple_pipeline/tracks.cc index 620f6fc880a..f9e50d20af9 100644 --- a/extern/libmv/libmv/simple_pipeline/tracks.cc +++ b/extern/libmv/libmv/simple_pipeline/tracks.cc @@ -18,12 +18,13 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/simple_pipeline/tracks.h" + #include #include #include #include "libmv/numeric/numeric.h" -#include "libmv/simple_pipeline/tracks.h" namespace libmv { @@ -82,7 +83,8 @@ vector Tracks::MarkersInBothImages(int image1, int image2) const { return markers; } -vector Tracks::MarkersForTracksInBothImages(int image1, int image2) const { +vector Tracks::MarkersForTracksInBothImages(int image1, + int image2) const { std::vector image1_tracks; std::vector image2_tracks; @@ -106,7 +108,7 @@ vector Tracks::MarkersForTracksInBothImages(int image1, int image2) cons vector markers; for (int i = 0; i < markers_.size(); ++i) { if ((markers_[i].image == image1 || markers_[i].image == image2) && - std::binary_search(intersection.begin(),intersection.end(), + std::binary_search(intersection.begin(), intersection.end(), markers_[i].track)) { markers.push_back(markers_[i]); } diff --git a/extern/libmv/libmv/tracking/brute_region_tracker.cc b/extern/libmv/libmv/tracking/brute_region_tracker.cc index ef36dffb56d..01279a04347 100644 --- a/extern/libmv/libmv/tracking/brute_region_tracker.cc +++ b/extern/libmv/libmv/tracking/brute_region_tracker.cc @@ -63,13 +63,13 @@ void *aligned_malloc(int size, int alignment) { #elif __FreeBSD__ void *result; - if(posix_memalign(&result, alignment, size)) { + if (posix_memalign(&result, alignment, size)) { // non-zero means allocation error // either no allocation or bad alignment value return NULL; } return result; -#else // This is for Linux. +#else // This is for Linux. return memalign(alignment, size); #endif } @@ -107,7 +107,7 @@ bool RegionIsInBounds(const FloatImage &image1, #ifdef __SSE2__ -// Compute the sub of absolute differences between the arrays "a" and "b". +// Compute the sub of absolute differences between the arrays "a" and "b". // The array "a" is assumed to be 16-byte aligned, while "b" is not. The // result is returned as the first and third elements of __m128i if // interpreted as a 4-element 32-bit integer array. The SAD is the sum of the @@ -122,7 +122,7 @@ inline static __m128i SumOfAbsoluteDifferencesContiguousSSE( unsigned int size, __m128i sad) { // Do the bulk of the work as 16-way integer operations. - for(unsigned int j = 0; j < size / 16; j++) { + for (unsigned int j = 0; j < size / 16; j++) { sad = _mm_add_epi32(sad, _mm_sad_epu8( _mm_load_si128 ((__m128i*)(a + 16 * j)), _mm_loadu_si128((__m128i*)(b + 16 * j)))); } @@ -307,7 +307,8 @@ bool BruteRegionTracker::Track(const FloatImage &image1, // Convert the search area directly to bytes without sampling. unsigned char *search_area; int search_area_stride; - FloatArrayToByteArrayWithPadding(image_and_gradient2, &search_area, &search_area_stride); + FloatArrayToByteArrayWithPadding(image_and_gradient2, &search_area, + &search_area_stride); // Try all possible locations inside the search area. Yes, everywhere. int best_i = -1, best_j = -1, best_sad = INT_MAX; @@ -361,9 +362,11 @@ bool BruteRegionTracker::Track(const FloatImage &image1, // Compute the Pearson product-moment correlation coefficient to check // for sanity. - double correlation = PearsonProductMomentCorrelation(image_and_gradient1_sampled, - image_and_gradient2_sampled, - pattern_width); + double correlation = PearsonProductMomentCorrelation( + image_and_gradient1_sampled, + image_and_gradient2_sampled, + pattern_width); + LG << "Final correlation: " << correlation; if (correlation < minimum_correlation) { diff --git a/extern/libmv/libmv/tracking/brute_region_tracker.h b/extern/libmv/libmv/tracking/brute_region_tracker.h index 8ff0c1b7965..a699c42ee92 100644 --- a/extern/libmv/libmv/tracking/brute_region_tracker.h +++ b/extern/libmv/libmv/tracking/brute_region_tracker.h @@ -30,7 +30,7 @@ struct BruteRegionTracker : public RegionTracker { BruteRegionTracker() : half_window_size(4), minimum_correlation(0.78) {} - + virtual ~BruteRegionTracker() {} // Tracker interface. diff --git a/extern/libmv/libmv/tracking/esm_region_tracker.cc b/extern/libmv/libmv/tracking/esm_region_tracker.cc index a8dc46d439b..1e44012e208 100644 --- a/extern/libmv/libmv/tracking/esm_region_tracker.cc +++ b/extern/libmv/libmv/tracking/esm_region_tracker.cc @@ -72,7 +72,7 @@ bool EsmRegionTracker::Track(const FloatImage &image1, << ", hw=" << half_window_size << "."; return false; } - + // XXX // TODO(keir): Delete the block between the XXX's once the planar tracker is // integrated into blender. @@ -101,7 +101,7 @@ bool EsmRegionTracker::Track(const FloatImage &image1, options.max_iterations = 20; options.sigma = sigma; options.use_esm = true; - + TrackRegionResult result; TrackRegion(image1, image2, xx1, yy1, options, xx2, yy2, &result); @@ -129,7 +129,7 @@ bool EsmRegionTracker::Track(const FloatImage &image1, &image_and_gradient1_sampled); // Step 0: Initialize delta = 0.01. - // + // // Ignored for my "normal" LM loop. // Step 1: Warp I with W(x, p) to compute I(W(x; p). @@ -179,8 +179,7 @@ bool EsmRegionTracker::Track(const FloatImage &image1, } } - double tau = 1e-4, eps1, eps2, eps3; - eps1 = eps2 = eps3 = 1e-15; + double tau = 1e-4; double mu = tau * std::max(H_image1(0, 0), H_image1(1, 1)); double nu = M_E; @@ -201,7 +200,7 @@ bool EsmRegionTracker::Track(const FloatImage &image1, image_and_gradient1_sampled(r, c, 2)); Vec2 g2(image_and_gradient2_sampled[current_image](r, c, 1), image_and_gradient2_sampled[current_image](r, c, 2)); - Vec2 g = g1 + g2; // Should be / 2.0, but do that outside the loop. + Vec2 g = g1 + g2; // Should be / 2.0, but do that outside the loop. H += g * g.transpose(); } } @@ -282,9 +281,10 @@ bool EsmRegionTracker::Track(const FloatImage &image1, if (d.squaredNorm() < min_update_squared_distance) { // Compute the Pearson product-moment correlation coefficient to check // for sanity. - double correlation = PearsonProductMomentCorrelation(image_and_gradient1_sampled, - image_and_gradient2_sampled[new_image], - width); + double correlation = PearsonProductMomentCorrelation( + image_and_gradient1_sampled, + image_and_gradient2_sampled[new_image], + width); LG << "Final correlation: " << correlation; // Note: Do the comparison here to handle nan's correctly (since all diff --git a/extern/libmv/libmv/tracking/esm_region_tracker.h b/extern/libmv/libmv/tracking/esm_region_tracker.h index 602b32f1675..97b8e27dba2 100644 --- a/extern/libmv/libmv/tracking/esm_region_tracker.h +++ b/extern/libmv/libmv/tracking/esm_region_tracker.h @@ -40,7 +40,7 @@ struct EsmRegionTracker : public RegionTracker { min_update_squared_distance(1e-4), sigma(0.9), minimum_correlation(0.78) {} - + virtual ~EsmRegionTracker() {} // Tracker interface. diff --git a/extern/libmv/libmv/tracking/hybrid_region_tracker.h b/extern/libmv/libmv/tracking/hybrid_region_tracker.h index a18341f0ce4..967d2afd1e5 100644 --- a/extern/libmv/libmv/tracking/hybrid_region_tracker.h +++ b/extern/libmv/libmv/tracking/hybrid_region_tracker.h @@ -34,7 +34,7 @@ class HybridRegionTracker : public RegionTracker { RegionTracker *fine_tracker) : coarse_tracker_(coarse_tracker), fine_tracker_(fine_tracker) {} - + virtual ~HybridRegionTracker() {} // Tracker interface. diff --git a/extern/libmv/libmv/tracking/klt_region_tracker.cc b/extern/libmv/libmv/tracking/klt_region_tracker.cc index c8e605de572..dbbf9f0b996 100644 --- a/extern/libmv/libmv/tracking/klt_region_tracker.cc +++ b/extern/libmv/libmv/tracking/klt_region_tracker.cc @@ -18,8 +18,9 @@ // 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/tracking/klt_region_tracker.h" + +#include "libmv/logging/logging.h" #include "libmv/image/image.h" #include "libmv/image/convolve.h" #include "libmv/image/sample.h" @@ -144,7 +145,8 @@ bool KltRegionTracker::Track(const FloatImage &image1, LG << "Determinant " << determinant << " is too small; failing tracking."; return false; } - LG << "x=" << *x2 << ", y=" << *y2 << ", dx=" << dx << ", dy=" << dy << ", det=" << determinant; + LG << "x=" << *x2 << ", y=" << *y2 << ", dx=" << dx << ", dy=" << dy + << ", det=" << determinant; // If the update is small, then we probably found the target. if (dx * dx + dy * dy < min_update_squared_distance) { diff --git a/extern/libmv/libmv/tracking/klt_region_tracker.h b/extern/libmv/libmv/tracking/klt_region_tracker.h index 2b2d8a9a49d..43977757084 100644 --- a/extern/libmv/libmv/tracking/klt_region_tracker.h +++ b/extern/libmv/libmv/tracking/klt_region_tracker.h @@ -33,7 +33,7 @@ struct KltRegionTracker : public RegionTracker { min_determinant(1e-6), min_update_squared_distance(1e-6), sigma(0.9) {} - + virtual ~KltRegionTracker() {} // Tracker interface. diff --git a/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc b/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc index 581e984b569..9e7446c7d08 100644 --- a/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc +++ b/extern/libmv/libmv/tracking/lmicklt_region_tracker.cc @@ -79,7 +79,7 @@ bool LmickltRegionTracker::Track(const FloatImage &image1, << ", hw=" << half_window_size << "."; return false; } - + int width = 2 * half_window_size + 1; // TODO(keir): Avoid recomputing gradients for e.g. the pyramid tracker. @@ -100,7 +100,7 @@ bool LmickltRegionTracker::Track(const FloatImage &image1, &image_and_gradient1_sampled); // Step 0: Initialize delta = 0.01. - // + // // Ignored for my "normal" LM loop. double reasonable_error = @@ -152,8 +152,7 @@ bool LmickltRegionTracker::Track(const FloatImage &image1, } } - double tau = 1e-3, eps1, eps2, eps3; - eps1 = eps2 = eps3 = 1e-15; + double tau = 1e-3; double mu = tau * std::max(H(0, 0), H(1, 1)); double nu = 2.0; diff --git a/extern/libmv/libmv/tracking/lmicklt_region_tracker.h b/extern/libmv/libmv/tracking/lmicklt_region_tracker.h index 744f2dfe277..95f4ba5f9b5 100644 --- a/extern/libmv/libmv/tracking/lmicklt_region_tracker.h +++ b/extern/libmv/libmv/tracking/lmicklt_region_tracker.h @@ -40,7 +40,7 @@ struct LmickltRegionTracker : public RegionTracker { min_determinant(1e-6), min_update_squared_distance(1e-6), sigma(0.9) {} - + virtual ~LmickltRegionTracker() {} // Tracker interface. diff --git a/extern/libmv/libmv/tracking/pyramid_region_tracker.cc b/extern/libmv/libmv/tracking/pyramid_region_tracker.cc index c177f9c5a83..4db501050f3 100644 --- a/extern/libmv/libmv/tracking/pyramid_region_tracker.cc +++ b/extern/libmv/libmv/tracking/pyramid_region_tracker.cc @@ -18,13 +18,14 @@ // 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 #include "libmv/image/convolve.h" #include "libmv/image/image.h" #include "libmv/image/sample.h" #include "libmv/logging/logging.h" -#include "libmv/tracking/pyramid_region_tracker.h" namespace libmv { diff --git a/extern/libmv/libmv/tracking/retrack_region_tracker.cc b/extern/libmv/libmv/tracking/retrack_region_tracker.cc index b3230b1b173..4d230086d28 100644 --- a/extern/libmv/libmv/tracking/retrack_region_tracker.cc +++ b/extern/libmv/libmv/tracking/retrack_region_tracker.cc @@ -18,11 +18,11 @@ // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS // IN THE SOFTWARE. +#include "libmv/tracking/retrack_region_tracker.h" + #include #include -#include "libmv/tracking/retrack_region_tracker.h" - namespace libmv { bool RetrackRegionTracker::Track(const FloatImage &image1, diff --git a/extern/libmv/libmv/tracking/track_region.cc b/extern/libmv/libmv/tracking/track_region.cc index 036b7c94966..6e9fb95654d 100644 --- a/extern/libmv/libmv/tracking/track_region.cc +++ b/extern/libmv/libmv/tracking/track_region.cc @@ -60,6 +60,8 @@ struct JetOps { } static void ScaleDerivative(double scale_by, T *value) { // For double, there is no derivative to scale. + (void) scale_by; // Ignored. + (void) value; // Ignored. } }; @@ -85,7 +87,9 @@ struct Chain { const FunctionType dfdx[kNumArgs], const ArgumentType x[kNumArgs]) { // In the default case of scalars, there's nothing to do since there are no - // derivatives to propagate. + // derivatives to propagate. + (void) dfdx; // Ignored. + (void) x; // Ignored. return f; } }; @@ -276,8 +280,8 @@ class PixelDifferenceCostFunctor { int num_samples_y, const Warp &warp) : options_(options), - image_and_gradient1_(image_and_gradient1), - image_and_gradient2_(image_and_gradient2), + image_and_gradient1_(image_and_gradient1), + image_and_gradient2_(image_and_gradient2), canonical_to_image1_(canonical_to_image1), num_samples_x_(num_samples_x), num_samples_y_(num_samples_y), @@ -348,7 +352,7 @@ class PixelDifferenceCostFunctor { // // Note that partial masks are not short circuited. To see why short // circuiting produces bitwise-exact same results, consider that the - // residual for each pixel is + // residual for each pixel is // // residual = mask * (src - dst) , // @@ -431,11 +435,10 @@ class PixelDifferenceCostFunctor { return true; } - // For normalized matching, the average and + // For normalized matching, the average and template void ComputeNormalizingCoefficient(const T *warp_parameters, T *dst_mean) const { - *dst_mean = T(0.0); double num_samples = 0.0; for (int r = 0; r < num_samples_y_; ++r) { @@ -443,7 +446,7 @@ class PixelDifferenceCostFunctor { // Use the pre-computed image1 position. Vec2 image1_position(pattern_positions_(r, c, 0), pattern_positions_(r, c, 1)); - + // Sample the mask early; if it's zero, this pixel has no effect. This // allows early bailout from the expensive sampling that happens below. double mask_value = 1.0; @@ -484,28 +487,28 @@ class PixelDifferenceCostFunctor { LG << "Normalization for dst:" << *dst_mean; } - // TODO(keir): Consider also computing the cost here. - double PearsonProductMomentCorrelationCoefficient( - const double *warp_parameters) const { - for (int i = 0; i < Warp::NUM_PARAMETERS; ++i) { - VLOG(2) << "Correlation warp_parameters[" << i << "]: " - << warp_parameters[i]; - } + // TODO(keir): Consider also computing the cost here. + double PearsonProductMomentCorrelationCoefficient( + const double *warp_parameters) const { + for (int i = 0; i < Warp::NUM_PARAMETERS; ++i) { + VLOG(2) << "Correlation warp_parameters[" << i << "]: " + << warp_parameters[i]; + } - // The single-pass PMCC computation is somewhat numerically unstable, but - // it's sufficient for the tracker. - double sX = 0, sY = 0, sXX = 0, sYY = 0, sXY = 0; + // The single-pass PMCC computation is somewhat numerically unstable, but + // it's sufficient for the tracker. + double sX = 0, sY = 0, sXX = 0, sYY = 0, sXY = 0; - // Due to masking, it's important to account for fractional samples. - // For example, samples with a 50% mask are counted as a half sample. - double num_samples = 0; + // Due to masking, it's important to account for fractional samples. + // For example, samples with a 50% mask are counted as a half sample. + double num_samples = 0; - for (int r = 0; r < num_samples_y_; ++r) { - for (int c = 0; c < num_samples_x_; ++c) { + for (int r = 0; r < num_samples_y_; ++r) { + for (int c = 0; c < num_samples_x_; ++c) { // Use the pre-computed image1 position. Vec2 image1_position(pattern_positions_(r, c, 0), pattern_positions_(r, c, 1)); - + double mask_value = 1.0; if (options_.image1_mask != NULL) { mask_value = pattern_mask_(r, c); @@ -795,7 +798,7 @@ struct TranslationRotationWarp { // Obtain the rotation via orthorgonal procrustes. Mat2 correlation_matrix; for (int i = 0; i < 4; ++i) { - correlation_matrix += q1.CornerRelativeToCentroid(i) * + correlation_matrix += q1.CornerRelativeToCentroid(i) * q2.CornerRelativeToCentroid(i).transpose(); } Mat2 R = OrthogonalProcrustes(correlation_matrix); @@ -863,7 +866,7 @@ struct TranslationRotationScaleWarp { // Obtain the rotation via orthorgonal procrustes. Mat2 correlation_matrix; for (int i = 0; i < 4; ++i) { - correlation_matrix += q1.CornerRelativeToCentroid(i) * + correlation_matrix += q1.CornerRelativeToCentroid(i) * q2.CornerRelativeToCentroid(i).transpose(); } Mat2 R = OrthogonalProcrustes(correlation_matrix); @@ -938,7 +941,7 @@ struct AffineWarp { Vec2 v1 = q1.CornerRelativeToCentroid(i); Vec2 v2 = q2.CornerRelativeToCentroid(i); - Q1.row(2 * i + 0) << v1[0], v1[1], 0, 0 ; + Q1.row(2 * i + 0) << v1[0], v1[1], 0, 0; Q1.row(2 * i + 1) << 0, 0, v1[0], v1[1]; Q2(2 * i + 0) = v2[0]; @@ -1039,6 +1042,9 @@ struct HomographyWarp { void PickSampling(const double *x1, const double *y1, const double *x2, const double *y2, int *num_samples_x, int *num_samples_y) { + (void) x2; // Ignored. + (void) y2; // Ignored. + Vec2 a0(x1[0], y1[0]); Vec2 a1(x1[1], y1[1]); Vec2 a2(x1[2], y1[2]); @@ -1074,6 +1080,10 @@ void PickSampling(const double *x1, const double *y1, bool SearchAreaTooBigForDescent(const FloatImage &image2, const double *x2, const double *y2) { // TODO(keir): Check the bounds and enable only when it makes sense. + (void) image2; // Ignored. + (void) x2; // Ignored. + (void) y2; // Ignored. + return true; } @@ -1150,12 +1160,12 @@ void CreateBrutePattern(const double *x1, const double *y1, inverse_warp.Forward(inverse_warp.parameters, dst_x, dst_y, &src_x, &src_y); - + if (PointInQuad(x1, y1, src_x, src_y)) { (*pattern)(i, j) = SampleLinear(image1, src_y, src_x); (*mask)(i, j) = 1.0; if (image1_mask) { - (*mask)(i, j) = SampleLinear(*image1_mask, src_y, src_x);; + (*mask)(i, j) = SampleLinear(*image1_mask, src_y, src_x); } } else { (*pattern)(i, j) = 0.0; @@ -1218,7 +1228,7 @@ bool BruteTranslationOnlyInitialize(const FloatImage &image1, // // TODO(keir): There are a number of possible optimizations here. One choice // is to make a grid and only try one out of every N possible samples. - // + // // Another, slightly more clever idea, is to compute some sort of spatial // frequency distribution of the pattern patch. If the spatial resolution is // high (e.g. a grating pattern or fine lines) then checking every possible @@ -1375,7 +1385,7 @@ void TemplatedTrackRegion(const FloatImage &image1, num_samples_x, num_samples_y, warp); - problem.AddResidualBlock( + problem.AddResidualBlock( new ceres::AutoDiffCostFunction< PixelDifferenceCostFunctor, ceres::DYNAMIC, diff --git a/extern/libmv/libmv/tracking/trklt_region_tracker.cc b/extern/libmv/libmv/tracking/trklt_region_tracker.cc index f19315b6b11..05ef3d1d272 100644 --- a/extern/libmv/libmv/tracking/trklt_region_tracker.cc +++ b/extern/libmv/libmv/tracking/trklt_region_tracker.cc @@ -158,7 +158,8 @@ bool TrkltRegionTracker::Track(const FloatImage &image1, LG << "Determinant " << determinant << " is too small; failing tracking."; return false; } - LG << "x=" << *x2 << ", y=" << *y2 << ", dx=" << d[0] << ", dy=" << d[1] << ", det=" << determinant; + LG << "x=" << *x2 << ", y=" << *y2 << ", dx=" << d[0] << ", dy=" << d[1] + << ", det=" << determinant; // If the update is small, then we probably found the target. diff --git a/extern/libmv/third_party/fast/fast_10.c b/extern/libmv/third_party/fast/fast_10.c index 3af63869478..7929a653a12 100644 --- a/extern/libmv/third_party/fast/fast_10.c +++ b/extern/libmv/third_party/fast/fast_10.c @@ -1,8 +1,7 @@ /*This is mechanically generated code*/ #include -typedef struct { int x, y; } xy; -typedef unsigned char byte; +#include "fast.h" int fast10_corner_score(const byte* p, const int pixel[], int bstart) { diff --git a/extern/libmv/third_party/fast/fast_11.c b/extern/libmv/third_party/fast/fast_11.c index b4af4309521..589ddb169f1 100644 --- a/extern/libmv/third_party/fast/fast_11.c +++ b/extern/libmv/third_party/fast/fast_11.c @@ -1,8 +1,7 @@ /*This is mechanically generated code*/ #include -typedef struct { int x, y; } xy; -typedef unsigned char byte; +#include "fast.h" int fast11_corner_score(const byte* p, const int pixel[], int bstart) { diff --git a/extern/libmv/third_party/fast/fast_12.c b/extern/libmv/third_party/fast/fast_12.c index f73f68dd043..27426f8c0d8 100644 --- a/extern/libmv/third_party/fast/fast_12.c +++ b/extern/libmv/third_party/fast/fast_12.c @@ -1,8 +1,7 @@ /*This is mechanically generated code*/ #include -typedef struct { int x, y; } xy; -typedef unsigned char byte; +#include "fast.h" int fast12_corner_score(const byte* p, const int pixel[], int bstart) { diff --git a/extern/libmv/third_party/fast/fast_9.c b/extern/libmv/third_party/fast/fast_9.c index 6d33daeffbb..c22c1c0e6a0 100644 --- a/extern/libmv/third_party/fast/fast_9.c +++ b/extern/libmv/third_party/fast/fast_9.c @@ -1,8 +1,7 @@ /*This is mechanically generated code*/ #include -typedef struct { int x, y; } xy; -typedef unsigned char byte; +#include "fast.h" int fast9_corner_score(const byte* p, const int pixel[], int bstart) { diff --git a/extern/libmv/third_party/gflags/gflags_completions.cc b/extern/libmv/third_party/gflags/gflags_completions.cc index 75e9236ef63..f46a2e0814e 100644 --- a/extern/libmv/third_party/gflags/gflags_completions.cc +++ b/extern/libmv/third_party/gflags/gflags_completions.cc @@ -57,6 +57,7 @@ #include #include +#include "gflags/gflags_completions.h" #include "gflags/gflags.h" #include "util.h" diff --git a/extern/libmv/third_party/glog/src/logging.cc b/extern/libmv/third_party/glog/src/logging.cc index bf4d85be5b3..f36e14d2ce4 100644 --- a/extern/libmv/third_party/glog/src/logging.cc +++ b/extern/libmv/third_party/glog/src/logging.cc @@ -1464,6 +1464,10 @@ void LogToStderr() { namespace base { namespace internal { +/* Put prototypes here to suppress strict compiler warnings */ +bool GetExitOnDFatal(); +void SetExitOnDFatal(bool value); + bool GetExitOnDFatal() { MutexLock l(&log_mutex); return exit_on_dfatal; diff --git a/extern/libmv/third_party/glog/src/utilities.cc b/extern/libmv/third_party/glog/src/utilities.cc index c9db2b7f7cd..159b094c170 100644 --- a/extern/libmv/third_party/glog/src/utilities.cc +++ b/extern/libmv/third_party/glog/src/utilities.cc @@ -84,7 +84,7 @@ static void DebugWriteToStderr(const char* data, void *) { } } -void DebugWriteToString(const char* data, void *arg) { +static void DebugWriteToString(const char* data, void *arg) { reinterpret_cast(arg)->append(data); } diff --git a/extern/libmv/third_party/glog/src/vlog_is_on.cc b/extern/libmv/third_party/glog/src/vlog_is_on.cc index 5eefc96324a..4396aa8d066 100644 --- a/extern/libmv/third_party/glog/src/vlog_is_on.cc +++ b/extern/libmv/third_party/glog/src/vlog_is_on.cc @@ -62,6 +62,12 @@ _START_GOOGLE_NAMESPACE_ namespace glog_internal_namespace_ { +// Put protytype here to suppress strict compiler flags +GOOGLE_GLOG_DLL_DECL bool SafeFNMatch_(const char* pattern, + size_t patt_len, + const char* str, + size_t str_len); + // Implementation of fnmatch that does not need 0-termination // of arguments and does not allocate any memory, // but we only support "*" and "?" wildcards, not the "[...]" patterns. -- cgit v1.2.3