diff --git a/Shape_detection/examples/Shape_detection/region_growing_spheres_on_point_set_3.cpp b/Shape_detection/examples/Shape_detection/region_growing_spheres_on_point_set_3.cpp index fed431aa972..4851d3e199b 100644 --- a/Shape_detection/examples/Shape_detection/region_growing_spheres_on_point_set_3.cpp +++ b/Shape_detection/examples/Shape_detection/region_growing_spheres_on_point_set_3.cpp @@ -42,22 +42,24 @@ int main (int argc, char** argv) // Default parameters for data/spheres.ply const std::size_t k = 12; - const double tolerance = 0.01; + const double max_distance = 0.01; const double max_angle = 10.; const std::size_t min_region_size = 50; - // No constraint on radius - const double min_radius = 0.; - const double max_radius = std::numeric_limits::infinity(); - Neighbor_query neighbor_query( points, CGAL::parameters:: k_neighbors(k). point_map(points.point_map())); - Region_type region_type(points, tolerance, max_angle, min_region_size, - min_radius, max_radius, - points.point_map(), points.normal_map()); - Region_growing region_growing(points, neighbor_query, region_type); + Region_type region_type( + points, CGAL::parameters:: + maximum_distance(max_distance). + maximum_angle(max_angle). + minimum_region_size(min_region_size). + point_map(points.point_map()). + normal_map(points.normal_map())); + + Region_growing region_growing( + points, neighbor_query, region_type); // Add maps to get colored output Point_set::Property_map diff --git a/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_circle_fit_region.h b/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_circle_fit_region.h index 7612ababa26..1054db332d4 100644 --- a/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_circle_fit_region.h +++ b/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_circle_fit_region.h @@ -83,9 +83,6 @@ namespace Point_set { using Sqrt = typename Get_sqrt::Sqrt; public: - - /// @} - /// \name Initialization /// @{ @@ -202,7 +199,8 @@ namespace Point_set { CGAL_precondition(m_min_radius >= FT(0)); const FT m_max_radius = parameters::choose_parameter( - parameters::get_parameter(np, internal_np::maximum_radius), FT(std::numeric_limits::max())); + parameters::get_parameter(np, internal_np::maximum_radius), + FT(std::numeric_limits::max())); CGAL_precondition(m_max_radius >= m_min_radius); } @@ -288,7 +286,8 @@ namespace Point_set { \param query_index index of the query point - \param indices indices of the inliers of the region + \param indices + indices of the inliers of the region The first parameter is not used in this implementation. diff --git a/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_region.h b/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_region.h index 43cfde3c5fb..42c0b549288 100644 --- a/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_region.h +++ b/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_region.h @@ -16,299 +16,407 @@ #include -#include -#include - -#include -#include -#include - +// Internal includes. #include namespace CGAL { namespace Shape_detection { namespace Point_set { -/*! - \ingroup PkgShapeDetectionRGOnPoints - - \brief Region type based on the quality of the least squares sphere - fit applied to 3D points. - - This class fits a sphere to chunks of points in a 3D point set and - controls the quality of this fit. If all quality conditions are - satisfied, the chunk is accepted as a valid region, otherwise - rejected. - - \tparam GeomTraits - must be a model of `Kernel`. - - \tparam InputRange - must be a model of `ConstRange` whose iterator type is `RandomAccessIterator`. - - \tparam PointMap - must be an `LvaluePropertyMap` whose key type is the value type of the input - range and value type is `Kernel::Point_3`. - - \tparam NormalMap - must be an `LvaluePropertyMap` whose key type is the value type of the input - range and value type is `Kernel::Vector_3`. - - \cgalModels `RegionType` -*/ -template -class Least_squares_sphere_fit_region -{ - -public: - - /// \name Types - /// @{ - - /// \cond SKIP_IN_MANUAL - using Traits = GeomTraits; - using Input_range = InputRange; - using Point_map = PointMap; - using Normal_map = NormalMap; - /// \endcond - - /// Number type. - typedef typename GeomTraits::FT FT; - - /// \cond SKIP_IN_MANUAL - using Point_3 = typename Traits::Point_3; - using Vector_3 = typename Traits::Vector_3; - using Plane_3 = typename Traits::Plane_3; - - using Squared_length_3 = typename Traits::Compute_squared_length_3; - using Squared_distance_3 = typename Traits::Compute_squared_distance_3; - using Scalar_product_3 = typename Traits::Compute_scalar_product_3; - - using Get_sqrt = internal::Get_sqrt; - using Sqrt = typename Get_sqrt::Sqrt; - -private: - - const Input_range& m_input_range; - - const FT m_distance_threshold; - const FT m_normal_threshold; - const std::size_t m_min_region_size; - const FT m_min_radius; - const FT m_max_radius; - - const Point_map m_point_map; - const Normal_map m_normal_map; - - const Squared_length_3 m_squared_length_3; - const Squared_distance_3 m_squared_distance_3; - const Scalar_product_3 m_scalar_product_3; - const Sqrt m_sqrt; - - Point_3 m_center; - FT m_radius; - -public: - - /// \endcond - - /// @} - - /// \name Initialization - /// @{ - /*! - \brief initializes all internal data structures. + \ingroup PkgShapeDetectionRGOnPoints - \param input_range an instance of `InputRange` with 3D points and - corresponding 3D normal vectors + \brief Region type based on the quality of the least squares sphere + fit applied to 3D points. - \param distance_threshold the maximum distance from a point to a - sphere. %Default is 1. + This class fits a sphere to chunks of points in a 3D point set and + controls the quality of this fit. If all quality conditions are + satisfied, the chunk is accepted as a valid region, otherwise + rejected. - \param angle_threshold the maximum accepted angle in degrees - between the normal of a point and the radius of a sphere. %Default - is 25 degrees. + \tparam GeomTraits + a model of `Kernel` - \param min_region_size the minimum number of 3D points a region - must have. %Default is 3. + \tparam InputRange + a model of `ConstRange` whose iterator type is `RandomAccessIterator` - \param minimum_radius the radius below which an estimated sphere - is considered as invalid and discarded. %Default is 0 (no limit). + \tparam PointMap + a model of `ReadablePropertyMap` whose key type is the value type of the input + range and value type is `Kernel::Point_3` - \param maximum_radius the radius above which an estimated sphere - is considered as invalid and discarded. %Default is infinity (no - limit). + \tparam NormalMap + a model of `ReadablePropertyMap` whose key type is the value type of the input + range and value type is `Kernel::Vector_3` - \param point_map an instance of `PointMap` that maps an item from - `input_range` to `Kernel::Point_3` - - \param normal_map an instance of `NormalMap` that maps an item - from `input_range` to `Kernel::Vector_3` - - \param traits an instance of `GeomTraits`. - - \pre `input_range.size() > 0` - \pre `distance_threshold >= 0` - \pre `angle_threshold >= 0 && angle_threshold <= 90` - \pre `min_region_size > 0` + \cgalModels `RegionType` */ - Least_squares_sphere_fit_region (const InputRange& input_range, - const FT distance_threshold = FT(1), - const FT angle_threshold = FT(25), - const std::size_t min_region_size = 3, - const FT minimum_radius = FT(0), - const FT maximum_radius = std::numeric_limits::infinity(), - const PointMap point_map = PointMap(), - const NormalMap normal_map = NormalMap(), - const GeomTraits traits = GeomTraits()) - : m_input_range(input_range) - , m_distance_threshold(distance_threshold) - , m_normal_threshold(static_cast( - std::cos( - CGAL::to_double( - (angle_threshold * static_cast(CGAL_PI)) / FT(180))))) - , m_min_region_size(min_region_size) - , m_min_radius (minimum_radius) - , m_max_radius (maximum_radius) - , m_point_map(point_map) - , m_normal_map(normal_map) - , m_squared_length_3(traits.compute_squared_length_3_object()) - , m_squared_distance_3(traits.compute_squared_distance_3_object()) - , m_scalar_product_3(traits.compute_scalar_product_3_object()) - , m_sqrt(Get_sqrt::sqrt_object(traits)) - { - CGAL_precondition(input_range.size() > 0); - CGAL_precondition(distance_threshold >= FT(0)); - CGAL_precondition(angle_threshold >= FT(0) && angle_threshold <= FT(90)); - CGAL_precondition(min_region_size > 0); - CGAL_precondition(minimum_radius >= 0); - CGAL_precondition(maximum_radius > minimum_radius); - } + template< + typename GeomTraits, + typename InputRange, + typename PointMap, + typename NormalMap> + class Least_squares_sphere_fit_region { - /// @} + public: + /// \name Types + /// @{ - /// \name Access - /// @{ + /// \cond SKIP_IN_MANUAL + using Traits = GeomTraits; + using Input_range = InputRange; + using Point_map = PointMap; + using Normal_map = NormalMap; + /// \endcond - /*! - \brief implements `RegionType::is_part_of_region()`. + /// Number type. + typedef typename GeomTraits::FT FT; - This function controls if a point with the index `query_index` is - within the `distance_threshold` from the corresponding sphere and - if the angle between its normal and the sphere radius is within - the `angle_threshold`. If both conditions are satisfied, it - returns `true`, otherwise `false`. + /// @} - \param query_index index of the query point + private: + using Point_3 = typename Traits::Point_3; + using Vector_3 = typename Traits::Vector_3; + using Plane_3 = typename Traits::Plane_3; - \param indices indices of the inliers of the region + using Squared_length_3 = typename Traits::Compute_squared_length_3; + using Squared_distance_3 = typename Traits::Compute_squared_distance_3; + using Scalar_product_3 = typename Traits::Compute_scalar_product_3; - The first parameter is not used in this implementation. + using Get_sqrt = internal::Get_sqrt; + using Sqrt = typename Get_sqrt::Sqrt; - \return Boolean `true` or `false` + public: + /// \name Initialization + /// @{ - \pre `query_index >= 0 && query_index < input_range.size()` - */ - bool is_part_of_region (const std::size_t, - const std::size_t query_index, - const std::vector& indices) const - { - CGAL_precondition(query_index < m_input_range.size()); + /*! + \brief initializes all internal data structures. + + \tparam NamedParameters + a sequence of \ref bgl_namedparameters "Named Parameters" + + \param input_range + an instance of `InputRange` with 3D points and + corresponding 3D normal vectors + + \param np + a sequence of \ref bgl_namedparameters "Named Parameters" + among the ones listed below + + \cgalNamedParamsBegin + \cgalParamNBegin{maximum_distance} + \cgalParamDescription{the maximum distance from a point to a sphere} + \cgalParamType{`GeomTraits::FT`} + \cgalParamDefault{1} + \cgalParamNEnd + \cgalParamNBegin{maximum_angle} + \cgalParamDescription{the maximum angle in degrees between + the normal of a point and the radius of a sphere} + \cgalParamType{`GeomTraits::FT`} + \cgalParamDefault{25 degrees} + \cgalParamNEnd + \cgalParamNBegin{cosine_value} + \cgalParamDescription{the cos value computed as `cos(maximum_angle * PI / 180)`, + this parameter can be used instead of the `maximum_angle`} + \cgalParamType{`GeomTraits::FT`} + \cgalParamDefault{`cos(25 * PI / 180)`} + \cgalParamNEnd + \cgalParamNBegin{minimum_region_size} + \cgalParamDescription{the minimum number of 3D points a region must have} + \cgalParamType{`std::size_t`} + \cgalParamDefault{3} + \cgalParamNEnd + \cgalParamNBegin{minimum_radius} + \cgalParamDescription{the radius below which an estimated sphere + is considered as invalid and discarded} + \cgalParamType{`GeomTraits::FT`} + \cgalParamDefault{0, no limit} + \cgalParamNEnd + \cgalParamNBegin{maximum_radius} + \cgalParamDescription{the radius above which an estimated sphere + is considered as invalid and discarded.} + \cgalParamType{`GeomTraits::FT`} + \cgalParamDefault{+infinity, no limit} + \cgalParamNEnd + \cgalParamNBegin{point_map} + \cgalParamDescription{an instance of `PointMap` that maps an item from `input_range` + to `Kernel::Point_3`} + \cgalParamDefault{`PointMap()`} + \cgalParamNEnd + \cgalParamNBegin{normal_map} + \cgalParamDescription{ an instance of `NormalMap` that maps an item from `input_range` + to `Kernel::Vector_3`} + \cgalParamDefault{`NormalMap()`} + \cgalParamNEnd + \cgalParamNBegin{geom_traits} + \cgalParamDescription{an instance of `GeomTraits`} + \cgalParamDefault{`GeomTraits()`} + \cgalParamNEnd + \cgalNamedParamsEnd + + \pre `input_range.size() > 0` + \pre `maximum_distance >= 0` + \pre `maximum_angle >= 0 && maximum_angle <= 90` + \pre `cosine_value >= 0 && cosine_value <= 1` + \pre `minimum_region_size > 0` + \pre `minimum_radius >= 0` + \pre `maximum_radius >= minimum_radius` + */ + template + Least_squares_sphere_fit_region( + const InputRange& input_range, + const NamedParameters& np) : + m_input_range(input_range), + m_point_map(parameters::choose_parameter(parameters::get_parameter( + np, internal_np::point_map), PointMap())), + m_normal_map(parameters::choose_parameter(parameters::get_parameter( + np, internal_np::normal_map), NormalMap())), + m_traits(parameters::choose_parameter(parameters::get_parameter( + np, internal_np::geom_traits), GeomTraits())), + m_squared_length_3(m_traits.compute_squared_length_3_object()), + m_squared_distance_3(m_traits.compute_squared_distance_3_object()), + m_scalar_product_3(m_traits.compute_scalar_product_3_object()), + m_sqrt(Get_sqrt::sqrt_object(m_traits)) { + + CGAL_precondition(input_range.size() > 0); + const FT max_distance = parameters::choose_parameter( + parameters::get_parameter(np, internal_np::maximum_distance), FT(1)); + CGAL_precondition(max_distance >= FT(0)); + m_distance_threshold = max_distance; + + const FT max_angle = parameters::choose_parameter( + parameters::get_parameter(np, internal_np::maximum_angle), FT(25)); + CGAL_precondition(max_angle >= FT(0) && max_angle <= FT(90)); + + m_min_region_size = parameters::choose_parameter( + parameters::get_parameter(np, internal_np::minimum_region_size), 3); + CGAL_precondition(m_min_region_size > 0); + + const FT default_cos_value = static_cast(std::cos(CGAL::to_double( + (max_angle * static_cast(CGAL_PI)) / FT(180)))); + const FT cos_value = parameters::choose_parameter( + parameters::get_parameter(np, internal_np::cosine_value), default_cos_value); + CGAL_precondition(cos_value >= FT(0) && cos_value <= FT(1)); + m_cos_value_threshold = cos_value; + + const FT m_min_radius = parameters::choose_parameter( + parameters::get_parameter(np, internal_np::minimum_radius), FT(0)); + CGAL_precondition(m_min_radius >= FT(0)); + + const FT m_max_radius = parameters::choose_parameter( + parameters::get_parameter(np, internal_np::maximum_radius), + FT(std::numeric_limits::max())); + CGAL_precondition(m_max_radius >= m_min_radius); + } + + /*! + \brief initializes all internal data structures. + + \deprecated This constructor is deprecated since the version 5.5 of \cgal. + + \param input_range + an instance of `InputRange` with 3D points and + corresponding 3D normal vectors + + \param distance_threshold + the maximum distance from a point to a sphere. %Default is 1. + + \param angle_threshold + the maximum accepted angle in degrees between the normal of a point + and the radius of a sphere. %Default is 25 degrees. + + \param min_region_size + the minimum number of 3D points a region must have. %Default is 3. + + \param minimum_radius + the radius below which an estimated sphere is considered as invalid + and discarded. %Default is 0 (no limit). + + \param maximum_radius + the radius above which an estimated sphere is considered as invalid + and discarded. %Default is infinity (no limit). + + \param point_map + an instance of `PointMap` that maps an item from `input_range` to `Kernel::Point_3` + + \param normal_map + an instance of `NormalMap` that maps an item from `input_range` to `Kernel::Vector_3` + + \param traits + an instance of `GeomTraits` + + \pre `input_range.size() > 0` + \pre `distance_threshold >= 0` + \pre `angle_threshold >= 0 && angle_threshold <= 90` + \pre `min_region_size > 0` + \pre `minimum_radius >= 0` + \pre `maximum_radius >= minimum_radius` + */ + CGAL_DEPRECATED_MSG("This constructor is deprecated since the version 5.5 of CGAL!") + Least_squares_sphere_fit_region( + const InputRange& input_range, + const FT distance_threshold = FT(1), + const FT angle_threshold = FT(25), + const std::size_t min_region_size = 3, + const FT minimum_radius = FT(0), + const FT maximum_radius = FT(std::numeric_limits::max()), + const PointMap point_map = PointMap(), + const NormalMap normal_map = NormalMap(), + const GeomTraits traits = GeomTraits()) : + Least_squares_sphere_fit_region( + input_range, CGAL::parameters:: + maximum_distance(distance_threshold). + maximum_angle(angle_threshold). + minimum_region_size(min_region_size). + minimum_radius(minimum_radius). + maximum_radius(maximum_radius). + point_map(point_map). + normal_map(normal_map). + geom_traits(traits)) + { } + + /// @} + + /// \name Access + /// @{ + + /*! + \brief implements `RegionType::is_part_of_region()`. + + This function controls if a point with the index `query_index` is within + the `maximum_distance` from the corresponding sphere and if the angle between + its normal and the sphere radius is within the `maximum_angle`. If both conditions + are satisfied, it returns `true`, otherwise `false`. + + \param query_index + index of the query point + + \param indices + indices of the inliers of the region + + The first parameter is not used in this implementation. + + \return Boolean `true` or `false` + + \pre `query_index < input_range.size()` + */ + bool is_part_of_region( + const std::size_t, + const std::size_t query_index, + const std::vector& indices) const { + + CGAL_precondition(query_index < m_input_range.size()); + + // First, we need to integrate at least 6 points so that the + // computed sphere means something + if (indices.size() < 6) + return true; + + // TODO: Why do we get so many nan in this class? + if (std::isnan(m_radius)) + return false; + + // If radius is out of bound, nothing fits, early ending + if (m_radius < m_min_radius || m_radius > m_max_radius) + return false; + + const auto& key = *(m_input_range.begin() + query_index); + const Point_3& query_point = get(m_point_map, key); + Vector_3 normal = get(m_normal_map, key); + + const FT sq_dist = m_squared_distance_3(query_point, m_center); + if (std::isnan(sq_dist)) return false; + FT distance_to_center = m_sqrt (sq_dist); + FT distance_to_sphere = CGAL::abs (distance_to_center - m_radius); + + if (distance_to_sphere > m_distance_threshold) + return false; + + const FT sq_norm = normal * normal; + if (std::isnan(sq_norm)) return false; + normal = normal / m_sqrt (sq_norm); + + Vector_3 ray (m_center, query_point); + const FT sq_ray = ray * ray; + if (std::isnan(sq_ray)) return false; + ray = ray / m_sqrt (sq_ray); + + if (CGAL::abs (normal * ray) < m_cos_value_threshold) + return false; - // First, we need to integrate at least 6 points so that the - // computed sphere means something - if (indices.size() < 6) return true; + } - // TODO: Why do we get so many nan in this class? - if (std::isnan(m_radius)) - return false; + /*! + \brief implements `RegionType::is_valid_region()`. - // If radius is out of bound, nothing fits, early ending - if (m_radius < m_min_radius || m_radius > m_max_radius) - return false; + This function controls if the estimated radius is between `minimum_radius` + and `maximum_radius` and if the `region` contains at least `min_region_size` points. - const auto& key = *(m_input_range.begin() + query_index); - const Point_3& query_point = get(m_point_map, key); - Vector_3 normal = get(m_normal_map, key); + \param region + indices of points included in the region - const FT sq_dist = m_squared_distance_3(query_point, m_center); - if (std::isnan(sq_dist)) return false; - FT distance_to_center = m_sqrt (sq_dist); - FT distance_to_sphere = CGAL::abs (distance_to_center - m_radius); + \return Boolean `true` or `false` + */ + inline bool is_valid_region(const std::vector& region) const { + return ( + (m_min_radius <= m_radius && m_radius <= m_max_radius) && + (region.size() >= m_min_region_size) + ); + } - if (distance_to_sphere > m_distance_threshold) - return false; + /*! + \brief implements `RegionType::update()`. - const FT sq_norm = normal * normal; - if (std::isnan(sq_norm)) return false; - normal = normal / m_sqrt (sq_norm); + This function fits the least squares sphere to all points from the `region`. - Vector_3 ray (m_center, query_point); - const FT sq_ray = ray * ray; - if (std::isnan(sq_ray)) return false; - ray = ray / m_sqrt (sq_ray); + \param region + indices of points included in the region - if (CGAL::abs (normal * ray) < m_normal_threshold) - return false; + \return Boolean `true` if the sphere fitting succeeded and `false` otherwise - return true; - } + \pre `region.size() > 0` + */ + bool update(const std::vector& region) { - /*! - \brief implements `RegionType::is_valid_region()`. + CGAL_precondition(region.size() > 0); + auto unary_function + = [&](const std::size_t& idx) -> const Point_3& + { + return get (m_point_map, *(m_input_range.begin() + idx)); + }; - This function controls if the estimated radius is between - `minimum_radius` and `maximum_radius` and if the `region` contains - at least `min_region_size` points. + internal::sphere_fit + (make_range (boost::make_transform_iterator + (region.begin(), unary_function), + boost::make_transform_iterator + (region.end(), unary_function)), + m_sqrt, m_squared_distance_3, m_center, m_radius); + return true; + } - \param region - indices of points included in the region + /// @} - \return Boolean `true` or `false` - */ - inline bool is_valid_region(const std::vector& region) const - { - return ((m_min_radius <= m_radius && m_radius <= m_max_radius) - && (region.size() >= m_min_region_size)); - } + private: + const Input_range& m_input_range; + const Point_map m_point_map; + const Normal_map m_normal_map; + const Traits m_traits; - /*! - \brief implements `RegionType::update()`. + FT m_distance_threshold; + FT m_cos_value_threshold; + std::size_t m_min_region_size; + FT m_min_radius; + FT m_max_radius; - This function fits the least squares sphere to all points from the `region`. + const Squared_length_3 m_squared_length_3; + const Squared_distance_3 m_squared_distance_3; + const Scalar_product_3 m_scalar_product_3; + const Sqrt m_sqrt; - \param region - indices of points included in the region - - \pre `region.size() > 0` - */ - bool update(const std::vector& region) - { - CGAL_precondition(region.size() > 0); - - auto unary_function - = [&](const std::size_t& idx) -> const Point_3& - { - return get (m_point_map, *(m_input_range.begin() + idx)); - }; - - internal::sphere_fit - (make_range (boost::make_transform_iterator - (region.begin(), unary_function), - boost::make_transform_iterator - (region.end(), unary_function)), - m_sqrt, m_squared_distance_3, m_center, m_radius); - return true; - } - - /// @} - -}; + Point_3 m_center; + FT m_radius; + }; } // namespace Point_set } // namespace Shape_detection diff --git a/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_sorting.h b/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_sorting.h index a28f131be38..f9b7aa00b68 100644 --- a/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_sorting.h +++ b/Shape_detection/include/CGAL/Shape_detection/Region_growing/Region_growing_on_point_set/Least_squares_sphere_fit_sorting.h @@ -8,7 +8,7 @@ // SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial // // -// Author(s) : Simon Giraudot, Dmitry Anisimov +// Author(s) : Simon Giraudot // #ifndef CGAL_SHAPE_DETECTION_REGION_GROWING_POINT_SET_LEAST_SQUARES_SPHERE_FIT_SORTING_H @@ -16,202 +16,235 @@ #include -// STL includes. -#include -#include - -// CGAL includes. -#include -#include -#include - // Internal includes. -#include #include namespace CGAL { namespace Shape_detection { namespace Point_set { -/*! - \ingroup PkgShapeDetectionRGOnPoints - - \brief Sorting of 3D points with respect to the local sphere fit quality. - - Indices of 3D input points are sorted with respect to the quality of the - least squares sphere fit applied to the neighboring points of each point. - - \tparam GeomTraits - must be a model of `Kernel`. - - \tparam InputRange - must be a model of `ConstRange` whose iterator type is `RandomAccessIterator`. - - \tparam NeighborQuery - must be a model of `NeighborQuery`. - - \tparam PointMap - must be an `LvaluePropertyMap` whose key type is the value type of the input - range and value type is `Kernel::Point_3`. -*/ -template -class Least_squares_sphere_fit_sorting -{ - -public: - - /// \name Types - /// @{ - - /// \cond SKIP_IN_MANUAL - using Traits = GeomTraits; - using Input_range = InputRange; - using Neighbor_query = NeighborQuery; - using Point_map = PointMap; - using Seed_map = internal::Seed_property_map; - /// \endcond - -#ifdef DOXYGEN_RUNNING /*! - an `LvaluePropertyMap` whose key and value type is `std::size_t`. - This map provides an access to the ordered indices of input points. + \ingroup PkgShapeDetectionRGOnPoints + + \brief Sorting of 3D points with respect to the local sphere fit quality. + + Indices of 3D input points are sorted with respect to the quality of the + least squares sphere fit applied to the neighboring points of each point. + + \tparam GeomTraits + a model of `Kernel` + + \tparam InputRange + a model of `ConstRange` whose iterator type is `RandomAccessIterator` + + \tparam NeighborQuery + a model of `NeighborQuery` + + \tparam PointMap + a model of `ReadablePropertyMap` whose key type is the value type of the input + range and value type is `Kernel::Point_3` */ - typedef unspecified_type Seed_map; -#endif + template< + typename GeomTraits, + typename InputRange, + typename NeighborQuery, + typename PointMap> + class Least_squares_sphere_fit_sorting { - /// @} + public: - /// \name Initialization - /// @{ + /// \name Types + /// @{ - /*! - \brief initializes all internal data structures. + /// \cond SKIP_IN_MANUAL + using Traits = GeomTraits; + using Input_range = InputRange; + using Neighbor_query = NeighborQuery; + using Point_map = PointMap; + using Seed_map = internal::Seed_property_map; + /// \endcond - \param input_range - an instance of `InputRange` with 3D points + #ifdef DOXYGEN_RUNNING + /*! + a model of `ReadablePropertyMap` whose key and value type is `std::size_t`. + This map provides an access to the ordered indices of input points. + */ + typedef unspecified_type Seed_map; + #endif - \param neighbor_query - an instance of `NeighborQuery` that is used internally to - access point's neighbors + /// @} - \param point_map - an instance of `PointMap` that maps an item from `input_range` - to `Kernel::Point_3` + private: + // using FT = typename Traits::FT; + // using Compare_scores = internal::Compare_scores; - \pre `input_range.size() > 0` - */ - Least_squares_sphere_fit_sorting (const InputRange& input_range, - NeighborQuery& neighbor_query, - const PointMap point_map = PointMap()) - : m_input_range(input_range) - , m_neighbor_query(neighbor_query) - , m_point_map(point_map) - , m_to_local_converter() - { - CGAL_precondition(input_range.size() > 0); + using Local_traits = Exact_predicates_inexact_constructions_kernel; + using Local_FT = typename Local_traits::FT; + using Local_point_3 = typename Local_traits::Point_3; + using To_local_converter = Cartesian_converter; + using Compare_scores = internal::Compare_scores; - m_order.resize(m_input_range.size()); - for (std::size_t i = 0; i < m_input_range.size(); ++i) - m_order[i] = i; - m_scores.resize(m_input_range.size()); - } + public: + /// \name Initialization + /// @{ - /// @} + /*! + \brief initializes all internal data structures. - /// \name Sorting - /// @{ + \tparam NamedParameters + a sequence of \ref bgl_namedparameters "Named Parameters" - /*! - \brief sorts indices of input points. - */ - void sort() - { - compute_scores(); - CGAL_postcondition(m_scores.size() > 0); + \param input_range + an instance of `InputRange` with 3D points - Compare_scores cmp(m_scores); - std::sort(m_order.begin(), m_order.end(), cmp); - } + \param neighbor_query + an instance of `NeighborQuery` that is used internally to + access point's neighbors - /// @} + \param np + a sequence of \ref bgl_namedparameters "Named Parameters" + among the ones listed below - /// \name Access - /// @{ + \cgalNamedParamsBegin + \cgalParamNBegin{point_map} + \cgalParamDescription{an instance of `PointMap` that maps an item from `input_range` + to `Kernel::Point_3`} + \cgalParamDefault{`PointMap()`} + \cgalParamNEnd + \cgalParamNBegin{geom_traits} + \cgalParamDescription{an instance of `GeomTraits`} + \cgalParamDefault{`GeomTraits()`} + \cgalParamNEnd + \cgalNamedParamsEnd - /*! - \brief returns an instance of `Seed_map` to access the ordered indices - of input points. - */ - Seed_map seed_map() - { - return Seed_map(m_order); - } + \pre `input_range.size() > 0` + */ + template + Least_squares_sphere_fit_sorting( + const InputRange& input_range, + NeighborQuery& neighbor_query, + const NamedParameters& np) : + m_input_range(input_range), + m_neighbor_query(neighbor_query), + m_point_map(parameters::choose_parameter(parameters::get_parameter( + np, internal_np::point_map), PointMap())), + m_traits(parameters::choose_parameter(parameters::get_parameter( + np, internal_np::geom_traits), GeomTraits())), + m_to_local_converter() { - /// @} - -private: - - // Types. - using Local_traits = Exact_predicates_inexact_constructions_kernel; - using Local_FT = typename Local_traits::FT; - using Local_point_3 = typename Local_traits::Point_3; - using To_local_converter = Cartesian_converter; - using Compare_scores = internal::Compare_scores; - - // Functions. - void compute_scores() - { - std::vector neighbors; - std::vector points; - - typename internal::Get_sqrt::Sqrt sqrt; - typename Local_traits::Compute_squared_distance_3 squared_distance_3; - - for (std::size_t i = 0; i < m_input_range.size(); ++i) - { - - neighbors.clear(); - m_neighbor_query(i, neighbors); - neighbors.push_back(i); - - points.clear(); - for (std::size_t j = 0; j < neighbors.size(); ++j) - { - CGAL_precondition(neighbors[j] < m_input_range.size()); - - const auto& key = *(m_input_range.begin() + neighbors[j]); - points.push_back(m_to_local_converter(get(m_point_map, key))); - } - CGAL_postcondition(points.size() == neighbors.size()); - - Local_point_3 fitted_center; - Local_FT fitted_radius; - - if (internal::sphere_fit (points, sqrt, squared_distance_3, fitted_center, fitted_radius)) - { - // Score is min squared distance to sphere - m_scores[i] = Local_FT(0); - for (const Local_point_3& p : points) - m_scores[i] += abs (sqrt(squared_distance_3(p, fitted_center)) - fitted_radius); - } - else - m_scores[i] = Local_FT(std::numeric_limits::infinity()); + CGAL_precondition(input_range.size() > 0); + m_order.resize(m_input_range.size()); + std::iota(m_order.begin(), m_order.end(), 0); + m_scores.resize(m_input_range.size()); } - } - // Fields. - const Input_range& m_input_range; - Neighbor_query& m_neighbor_query; - const Point_map m_point_map; + /*! + \brief initializes all internal data structures. - std::vector m_order; - std::vector m_scores; + \deprecated This constructor is deprecated since the version 5.5 of \cgal. - const To_local_converter m_to_local_converter; -}; + \param input_range + an instance of `InputRange` with 3D points + + \param neighbor_query + an instance of `NeighborQuery` that is used internally to + access point's neighbors + + \param point_map + an instance of `PointMap` that maps an item from `input_range` + to `Kernel::Point_3` + + \pre `input_range.size() > 0` + */ + CGAL_DEPRECATED_MSG("This constructor is deprecated since the version 5.5 of CGAL!") + Least_squares_sphere_fit_sorting( + const InputRange& input_range, + NeighborQuery& neighbor_query, + const PointMap point_map = PointMap()) : + Least_squares_sphere_fit_sorting( + input_range, neighbor_query, CGAL::parameters:: + point_map(point_map)) + { } + + /// @} + + /// \name Sorting + /// @{ + + /*! + \brief sorts indices of input points. + */ + void sort() { + + compute_scores(); + CGAL_postcondition(m_scores.size() > 0); + Compare_scores cmp(m_scores); + std::sort(m_order.begin(), m_order.end(), cmp); + } + + /// @} + + /// \name Access + /// @{ + + /*! + \brief returns an instance of `Seed_map` to access the ordered indices + of input points. + */ + Seed_map seed_map() { + return Seed_map(m_order); + } + + /// @} + + private: + const Input_range& m_input_range; + Neighbor_query& m_neighbor_query; + const Point_map m_point_map; + const Traits m_traits; + std::vector m_order; + std::vector m_scores; + const To_local_converter m_to_local_converter; + + void compute_scores() { + std::vector neighbors; + std::vector points; + + typename internal::Get_sqrt::Sqrt sqrt; + typename Local_traits::Compute_squared_distance_3 squared_distance_3; + + for (std::size_t i = 0; i < m_input_range.size(); ++i) + { + + neighbors.clear(); + m_neighbor_query(i, neighbors); + neighbors.push_back(i); + + points.clear(); + for (std::size_t j = 0; j < neighbors.size(); ++j) + { + CGAL_precondition(neighbors[j] < m_input_range.size()); + + const auto& key = *(m_input_range.begin() + neighbors[j]); + points.push_back(m_to_local_converter(get(m_point_map, key))); + } + CGAL_postcondition(points.size() == neighbors.size()); + + Local_point_3 fitted_center; + Local_FT fitted_radius; + + if (internal::sphere_fit (points, sqrt, squared_distance_3, fitted_center, fitted_radius)) + { + // Score is min squared distance to sphere + m_scores[i] = Local_FT(0); + for (const Local_point_3& p : points) + m_scores[i] += abs (sqrt(squared_distance_3(p, fitted_center)) - fitted_radius); + } + else + m_scores[i] = Local_FT(std::numeric_limits::infinity()); + } + } + }; } // namespace Point_set } // namespace Shape_detection