From 0875fa17ffaf43b229b9fde02a72e5ac4bfbc423 Mon Sep 17 00:00:00 2001 From: Simon Giraudot Date: Mon, 31 Aug 2015 09:34:17 +0200 Subject: [PATCH] Wrapper to automatically use Eigen_diagonalize if Eigen is available or Internal_diagonalize otherwise --- .../Approximate_min_ellipsoid_d_impl.h | 21 ++----- .../include/CGAL/Monge_via_jet_fitting.h | 49 +++++++++------- .../CGAL/linear_least_squares_fitting_2.h | 21 +------ .../CGAL/linear_least_squares_fitting_3.h | 19 +----- .../linear_least_squares_fitting_points_2.h | 7 --- ..._linear_least_squares_fitting_points_2.cpp | 21 ++----- ..._linear_least_squares_fitting_points_3.cpp | 20 ++----- ...inear_least_squares_fitting_segments_3.cpp | 23 ++------ ...linear_least_squares_fitting_spheres_3.cpp | 14 ++--- .../include/CGAL/Default_diagonalize_traits.h | 58 +++++++++++++++++++ 10 files changed, 112 insertions(+), 141 deletions(-) create mode 100644 Solver_interface/include/CGAL/Default_diagonalize_traits.h diff --git a/Bounding_volumes/include/CGAL/Approximate_min_ellipsoid_d/Approximate_min_ellipsoid_d_impl.h b/Bounding_volumes/include/CGAL/Approximate_min_ellipsoid_d/Approximate_min_ellipsoid_d_impl.h index f740211bf95..6b90635b5e9 100644 --- a/Bounding_volumes/include/CGAL/Approximate_min_ellipsoid_d/Approximate_min_ellipsoid_d_impl.h +++ b/Bounding_volumes/include/CGAL/Approximate_min_ellipsoid_d/Approximate_min_ellipsoid_d_impl.h @@ -19,11 +19,7 @@ // Author(s) : Kaspar Fischer -#ifdef CGAL_EIGEN3_ENABLED -#include -#else -#include -#endif +#include #include @@ -161,13 +157,8 @@ namespace CGAL { CGAL::cpp11::array eigenvectors; // Note: not neces. normalized. CGAL::cpp11::array eigenvalues; // Note: sorted ascendent. -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix + CGAL::Default_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix (matrix, eigenvalues, eigenvectors); -#else - CGAL::Internal_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix - (matrix, eigenvalues, eigenvectors); -#endif // normalize eigenvectors: double l1=1.0/std::sqrt(eigenvectors[2]*eigenvectors[2]+ @@ -208,13 +199,9 @@ namespace CGAL { CGAL::cpp11::array eigenvectors; // Note: not necessarily normalized. CGAL::cpp11::array eigenvalues; // Note: sorted ascendent. -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix + + CGAL::Default_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix (matrix, eigenvalues, eigenvectors); -#else - CGAL::Internal_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix - (matrix, eigenvalues, eigenvectors); -#endif // normalize eigenvectors: double l1 = 1.0/std::sqrt(eigenvectors[0] * eigenvectors[0]+ // x^2 diff --git a/Jet_fitting_3/include/CGAL/Monge_via_jet_fitting.h b/Jet_fitting_3/include/CGAL/Monge_via_jet_fitting.h index da1502f5bc6..597acffea5c 100644 --- a/Jet_fitting_3/include/CGAL/Monge_via_jet_fitting.h +++ b/Jet_fitting_3/include/CGAL/Monge_via_jet_fitting.h @@ -22,8 +22,8 @@ #include #include #include -#include #include +#include #include #include #ifdef CGAL_EIGEN3_ENABLED @@ -362,27 +362,30 @@ compute_PCA(InputIterator begin, InputIterator end) // assemble covariance matrix as a // semi-definite matrix. // Matrix numbering: - // 0 - // 1 2 - // 3 4 5 - FT covariance[6] = {xx,xy,yy,xz,yz,zz}; - FT eigen_values[3]; - FT eigen_vectors[9]; + // 0 1 2 + // 3 4 + // 5 + CGAL::cpp11::array covariance = {{ xx,xy,xz,yy,yz,zz }}; + CGAL::cpp11::array eigen_values = {{ 0., 0., 0. }}; + CGAL::cpp11::array eigen_vectors = {{ 0., 0., 0. }}; // solve for eigenvalues and eigenvectors. - // eigen values are sorted in descending order, + // eigen values are sorted in ascending order, // eigen vectors are sorted in accordance. - CGAL::internal::eigen_symmetric(covariance,3,eigen_vectors,eigen_values); + CGAL::Default_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix + (covariance, eigen_values, eigen_vectors); + //store in m_pca_basis for (int i=0; i<3; i++) { - m_pca_basis[i].first = eigen_values[i]; + m_pca_basis[i].first = eigen_values[2-i]; } - Vector_3 v1(eigen_vectors[0],eigen_vectors[1],eigen_vectors[2]); + + Vector_3 v1(eigen_vectors[6],eigen_vectors[7],eigen_vectors[8]); m_pca_basis[0].second = v1; Vector_3 v2(eigen_vectors[3],eigen_vectors[4],eigen_vectors[5]); m_pca_basis[1].second = v2; - Vector_3 v3(eigen_vectors[6],eigen_vectors[7],eigen_vectors[8]); + Vector_3 v3(eigen_vectors[0],eigen_vectors[1],eigen_vectors[2]); m_pca_basis[2].second = v3; switch_to_direct_orientation(m_pca_basis[0].second, m_pca_basis[1].second, @@ -523,15 +526,17 @@ compute_Monge_basis(const FT* A, Monge_form& monge_form) //in the new orthonormal basis (Y,Z) of the tangent plane : weingarten = inv *(1/det) * weingarten * change_XuXv2YZ; - //switch to eigen_symmetric algo for diagonalization of weingarten - FT W[3] = {weingarten(0,0), weingarten(1,0), weingarten(1,1)}; - FT eval[2]; - FT evec[4]; - //eval in decreasing order - CGAL::internal::eigen_symmetric(W,2,evec,eval); + // diagonalization of weingarten + CGAL::cpp11::array W = {{ weingarten(0,0), weingarten(1,0), weingarten(1,1) }}; + CGAL::cpp11::array eval = {{ 0., 0. }}; + CGAL::cpp11::array evec = {{ 0., 0., 0., 0. }}; - Vector_3 d_max = evec[0]*Y + evec[1]*Z, - d_min = evec[2]*Y + evec[3]*Z; + //eval in increasing order + CGAL::Default_diagonalize_traits::diagonalize_selfadjoint_covariance_matrix + (W, eval, evec); + + Vector_3 d_max = evec[2]*Y + evec[3]*Z, + d_min = evec[0]*Y + evec[1]*Z; switch_to_direct_orientation(d_max, d_min, normal); Aff_transformation change_basis (d_max[0], d_max[1], d_max[2], @@ -547,8 +552,8 @@ compute_Monge_basis(const FT* A, Monge_form& monge_form) monge_form.maximal_principal_direction() = L2D_converter(this->change_world2fitting.inverse()(d_max)); monge_form.minimal_principal_direction() = L2D_converter(this->change_world2fitting.inverse()(d_min)); monge_form.normal_direction() = L2D_converter(this->change_world2fitting.inverse()(normal)); - monge_form.coefficients()[0] = L2D_NTconverter()(eval[0]); - monge_form.coefficients()[1] = L2D_NTconverter()(eval[1]); + monge_form.coefficients()[0] = L2D_NTconverter()(eval[1]); + monge_form.coefficients()[1] = L2D_NTconverter()(eval[0]); } //end else } diff --git a/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_2.h b/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_2.h index 30b28a9e1eb..5d2bd42c6c8 100644 --- a/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_2.h +++ b/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_2.h @@ -29,12 +29,7 @@ #include #include #include -#include - -#ifdef CGAL_EIGEN3_ENABLED -#include -#endif - +#include #include namespace CGAL { @@ -77,12 +72,7 @@ linear_least_squares_fitting_2(InputIterator first, typedef typename Kernel_traits::Kernel Kernel; return CGAL::linear_least_squares_fitting_2 (first,beyond,line,centroid,tag,Kernel(), -#ifdef CGAL_EIGEN3_ENABLED - Eigen_diagonalize_traits() -#else - Internal_diagonalize_traits() -#endif - ); + Default_diagonalize_traits()); } template < typename InputIterator, @@ -99,12 +89,7 @@ linear_least_squares_fitting_2(InputIterator first, typedef typename Kernel_traits::Kernel Kernel; typename Kernel::Point_2 centroid; // unused return CGAL::linear_least_squares_fitting_2(first,beyond,line,centroid,tag,Kernel(), -#ifdef CGAL_EIGEN3_ENABLED - Eigen_diagonalize_traits() -#else - Internal_diagonalize_traits() -#endif - ); + Default_diagonalize_traits()); } diff --git a/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_3.h b/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_3.h index 37f33fa180c..448bc024725 100644 --- a/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_3.h +++ b/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_3.h @@ -29,10 +29,7 @@ #include #include -#include -#ifdef CGAL_EIGEN3_ENABLED -#include -#endif +#include #include @@ -80,12 +77,7 @@ linear_least_squares_fitting_3(InputIterator first, typedef typename std::iterator_traits::value_type Value_type; typedef typename Kernel_traits::Kernel Kernel; return CGAL::linear_least_squares_fitting_3(first,beyond,object,centroid,tag,Kernel(), -#ifdef CGAL_EIGEN3_ENABLED - Eigen_diagonalize_traits() -#else - Internal_diagonalize_traits() -#endif - ); + Default_diagonalize_traits()); } @@ -105,12 +97,7 @@ linear_least_squares_fitting_3(InputIterator first, typedef typename Kernel_traits::Kernel Kernel; typename Kernel::Point_3 centroid; // not used by caller return CGAL::linear_least_squares_fitting_3(first,beyond,object,centroid,tag,Kernel(), -#ifdef CGAL_EIGEN3_ENABLED - Eigen_diagonalize_traits() -#else - Internal_diagonalize_traits() -#endif - ); + Default_diagonalize_traits()); } diff --git a/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_2.h b/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_2.h index 4578d3afe99..e65a1321153 100644 --- a/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_2.h +++ b/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_2.h @@ -22,13 +22,6 @@ #include #include - -#include - -#ifdef CGAL_EIGEN3_ENABLED -#include -#endif - #include #include diff --git a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_2.cpp b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_2.cpp index 5ad47866509..b8aa1207cc3 100644 --- a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_2.cpp +++ b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_2.cpp @@ -3,11 +3,7 @@ #include #include #include -#include - -#ifdef CGAL_EIGEN3_ENABLED -#include -#endif +#include #include #include @@ -37,12 +33,7 @@ void test_2D() quality = linear_least_squares_fitting_2(points.begin(),points.end(),line,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_2(points.begin(),points.end(),line,centroid,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_2(points.begin(),points.end(),line,centroid,CGAL::Dimension_tag<0>(),k, -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits() -#else - CGAL::Internal_diagonalize_traits() -#endif - ); + CGAL::Default_diagonalize_traits()); std::cout << "done (quality: " << quality << ")" << std::endl; @@ -84,12 +75,8 @@ void test_2D_point_set(const unsigned int nb_points) quality = linear_least_squares_fitting_2(points.begin(),points.end(),line,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_2(points.begin(),points.end(),line,centroid,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_2(points.begin(),points.end(),line,centroid,CGAL::Dimension_tag<0>(),k, -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits() -#else - CGAL::Internal_diagonalize_traits() -#endif - ); + CGAL::Default_diagonalize_traits()); + std::cout << "done (quality: " << quality << ")" << std::endl; if(!line.is_horizontal()) diff --git a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_3.cpp b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_3.cpp index 6a1399e6204..bd94792cf91 100644 --- a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_3.cpp +++ b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_points_3.cpp @@ -3,10 +3,7 @@ #include #include -#include -#ifdef CGAL_EIGEN3_ENABLED -#include -#endif +#include #include @@ -49,12 +46,7 @@ void fit_point_set(std::list& points, quality = linear_least_squares_fitting_3(points.begin(),points.end(),line,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_3(points.begin(),points.end(),line,centroid,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_3(points.begin(),points.end(),line,centroid,CGAL::Dimension_tag<0>(),kernel, -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits() -#else - CGAL::Internal_diagonalize_traits() -#endif - ); + CGAL::Default_diagonalize_traits()); std::cout << "done (quality: " << quality << ")" << std::endl; @@ -62,12 +54,8 @@ void fit_point_set(std::list& points, quality = linear_least_squares_fitting_3(points.begin(),points.end(),plane,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_3(points.begin(),points.end(),plane,centroid,CGAL::Dimension_tag<0>()); quality = linear_least_squares_fitting_3(points.begin(),points.end(),plane,centroid,CGAL::Dimension_tag<0>(),kernel, -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits() -#else - CGAL::Internal_diagonalize_traits() -#endif - ); + CGAL::Default_diagonalize_traits()); + std::cout << "done (quality: " << quality << ")" << std::endl; } diff --git a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_segments_3.cpp b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_segments_3.cpp index 54df04e3094..3804156aedd 100644 --- a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_segments_3.cpp +++ b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_segments_3.cpp @@ -5,10 +5,7 @@ #include #include #include -#include -#ifdef CGAL_EIGEN3_ENABLED -#include -#endif +#include #include #include @@ -39,26 +36,14 @@ FT fit_set(std::list& segments, quality = linear_least_squares_fitting_3(segments.begin(),segments.end(),line,CGAL::Dimension_tag<1>()); quality = linear_least_squares_fitting_3(segments.begin(),segments.end(),line,centroid,CGAL::Dimension_tag<1>()); quality = linear_least_squares_fitting_3(segments.begin(),segments.end(),line,centroid,CGAL::Dimension_tag<1>(),kernel, -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits() -#else - CGAL::Internal_diagonalize_traits() -#endif - ); - + CGAL::Default_diagonalize_traits()); quality = linear_least_squares_fitting_3(segments.begin(),segments.end(),plane,CGAL::Dimension_tag<1>()); quality = linear_least_squares_fitting_3(segments.begin(),segments.end(),plane,centroid,CGAL::Dimension_tag<1>()); quality = linear_least_squares_fitting_3(segments.begin(),segments.end(),plane,centroid,CGAL::Dimension_tag<1>(),kernel, -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits() -#else - CGAL::Internal_diagonalize_traits() -#endif - ); + CGAL::Default_diagonalize_traits()); - - return quality; + return quality; } diff --git a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_spheres_3.cpp b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_spheres_3.cpp index 6cbdacb4740..d064d2034c5 100644 --- a/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_spheres_3.cpp +++ b/Principal_component_analysis/test/Principal_component_analysis/test_linear_least_squares_fitting_spheres_3.cpp @@ -2,6 +2,7 @@ #include #include +#include #include @@ -34,28 +35,23 @@ int main(void) Kernel kernel; Point centroid; -#ifdef CGAL_EIGEN3_ENABLED - CGAL::Eigen_diagonalize_traits diagonalize_traits; -#else - CGAL::Internal_diagonalize_traits diagonalize_traits; -#endif linear_least_squares_fitting_3(spheres.begin(),spheres.end(),line,CGAL::Dimension_tag<3>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),line,CGAL::Dimension_tag<2>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),line,centroid,CGAL::Dimension_tag<3>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),line,centroid,CGAL::Dimension_tag<2>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),line,centroid,CGAL::Dimension_tag<3>(),kernel, - diagonalize_traits); + CGAL::Default_diagonalize_traits()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),line,centroid,CGAL::Dimension_tag<2>(),kernel, - diagonalize_traits); + CGAL::Default_diagonalize_traits()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),plane,CGAL::Dimension_tag<3>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),plane,CGAL::Dimension_tag<2>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),plane,centroid,CGAL::Dimension_tag<3>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),plane,centroid,CGAL::Dimension_tag<2>()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),plane,centroid,CGAL::Dimension_tag<3>(),kernel, - diagonalize_traits); + CGAL::Default_diagonalize_traits()); linear_least_squares_fitting_3(spheres.begin(),spheres.end(),plane,centroid,CGAL::Dimension_tag<2>(),kernel, - diagonalize_traits); + CGAL::Default_diagonalize_traits()); return 0; } diff --git a/Solver_interface/include/CGAL/Default_diagonalize_traits.h b/Solver_interface/include/CGAL/Default_diagonalize_traits.h new file mode 100644 index 00000000000..93ead75f01b --- /dev/null +++ b/Solver_interface/include/CGAL/Default_diagonalize_traits.h @@ -0,0 +1,58 @@ +#ifndef CGAL_DEFAULT_DIAGONALIZE_TRAITS_H +#define CGAL_DEFAULT_DIAGONALIZE_TRAITS_H + +#ifdef CGAL_EIGEN3_ENABLED +#include +#else +#include +#endif + + +namespace CGAL { + + + // Wrapper class designed to automatically use + // Eigen_diagonalize_traits if Eigen is available and otherwise use + // the fallback Internal_diagonalize_traits class. + +template +class Default_diagonalize_traits{ + +#ifdef CGAL_EIGEN3_ENABLED + typedef Eigen_diagonalize_traits Base; +#else + typedef Internal_diagonalize_traits Base; +#endif + +public: + static bool + diagonalize_selfadjoint_covariance_matrix( + const cpp11::array& cov, + cpp11::array& eigenvalues) + { + return Base::diagonalize_selfadjoint_covariance_matrix (cov, eigenvalues); + } + + static bool + diagonalize_selfadjoint_covariance_matrix( + const cpp11::array& cov, + cpp11::array& eigenvalues, + cpp11::array& eigenvectors) + { + return Base::diagonalize_selfadjoint_covariance_matrix (cov, eigenvalues, eigenvectors); + } + + // Extract the eigenvector associated to the largest eigenvalue + static bool + extract_largest_eigenvector_of_covariance_matrix ( + const cpp11::array& cov, + cpp11::array &normal) + { + return Base::extract_largest_eigenvector_of_covariance_matrix (cov, normal); + } +}; + + +} // namespace CGAL + +#endif // CGAL_DEFAULT_DIAGONALIZE_TRAITS_H