Add a way to pass input normal to faces

Shape_detection/include/CGAL/Polygon_mesh_processing/region_growing.h

@@ -192,13 +192,13 @@ region_growing_of_planes_on_faces(const PolygonMesh& mesh,
                              Static_boolean_property_map<edge_descriptor, false>());
 
   using Neighbor_query = internal::One_ring_neighbor_query_with_constraints<PolygonMesh, ECM>;
-  using Region_type = RG_PM::Least_squares_plane_fit_region<Traits, PolygonMesh, VPM>;
-  using Sorting = RG_PM::Least_squares_plane_fit_sorting<Traits, PolygonMesh, Neighbor_query, VPM>;
+  using Region_type = RG_PM::Plane_face_region<Traits, PolygonMesh, VPM>;
+  using Sorting = RG_PM::Face_area_sorting<Traits, PolygonMesh, VPM>;
   using Region_growing = CGAL::Shape_detection::Region_growing<Neighbor_query, Region_type, RegionMap>;
 
   Neighbor_query neighbor_query(mesh, ecm);
   Region_type region_type(mesh, np);
-  Sorting sorting(mesh, neighbor_query, np);
+  Sorting sorting(mesh, np);
   sorting.sort();
 
   std::vector<typename Region_growing::Primitive_and_region> regions;

Shape_detection/include/CGAL/Shape_detection/Region_growing/Polygon_mesh.h

@@ -29,5 +29,7 @@
 
 #include <CGAL/Shape_detection/Region_growing/Polygon_mesh/Least_squares_plane_fit_region.h>
 #include <CGAL/Shape_detection/Region_growing/Polygon_mesh/Least_squares_plane_fit_sorting.h>
+#include <CGAL/Shape_detection/Region_growing/Polygon_mesh/Plane_face_region.h>
+#include <CGAL/Shape_detection/Region_growing/Polygon_mesh/Face_area_sorting.h>
 
 #endif // CGAL_SHAPE_DETECTION_REGION_GROWING_POLYGON_MESH_H

Shape_detection/include/CGAL/Shape_detection/Region_growing/Polygon_mesh/Face_area_sorting.h

@@ -0,0 +1,193 @@
+// Copyright (c) 2018 INRIA Sophia-Antipolis (France).
+// All rights reserved.
+//
+// This file is part of CGAL (www.cgal.org).
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
+//
+//
+// Author(s)     : Florent Lafarge, Simon Giraudot, Thien Hoang, Dmitry Anisimov
+//
+
+#ifndef CGAL_SHAPE_DETECTION_REGION_GROWING_POLYGON_MESH_FACE_AREA_SORTING_H
+#define CGAL_SHAPE_DETECTION_REGION_GROWING_POLYGON_MESH_FACE_AREA_SORTING_H
+
+#include <CGAL/license/Shape_detection.h>
+#include <CGAL/Shape_detection/Region_growing/internal/utils.h>
+
+namespace CGAL {
+namespace Shape_detection {
+namespace Polygon_mesh {
+
+  /*!
+    \ingroup PkgShapeDetectionRGOnMesh
+
+    \brief Sorting of polygon mesh faces with respect to the local plane fit quality.
+
+    `Items` of faces in a polygon mesh are sorted with respect to the quality of the
+    least squares plane fit applied to the vertices of incident faces of each face.
+
+    \tparam GeomTraits
+    a model of `Kernel`
+
+    \tparam PolygonMesh
+    a model of `FaceListGraph`
+
+    \tparam VertexToPointMap
+    a model of `ReadablePropertyMap` whose key type is the vertex type of a polygon mesh and
+    value type is `Kernel::Point_3`
+  */
+  template<
+  typename GeomTraits,
+  typename PolygonMesh,
+  typename VertexToPointMap = typename boost::property_map<PolygonMesh, CGAL::vertex_point_t>::const_type>
+  class Face_area_sorting
+  {
+    using vertex_descriptor = typename boost::graph_traits<PolygonMesh>::vertex_descriptor;
+    using face_descriptor = typename boost::graph_traits<PolygonMesh>::face_descriptor;
+    using halfedge_descriptor = typename boost::graph_traits<PolygonMesh>::halfedge_descriptor;
+    using FT = typename GeomTraits::FT;
+
+  public:
+    /// \name Types
+    /// @{
+
+    /// Item type.
+    using Item = typename boost::graph_traits<PolygonMesh>::face_descriptor;
+
+    /// Seed range.
+    using Seed_range = std::vector<Item>;
+
+    /// @}
+
+    /// \name Initialization
+    /// @{
+
+    /*!
+      \brief initializes all internal data structures.
+
+      \tparam NamedParameters
+      a sequence of \ref bgl_namedparameters "Named Parameters"
+
+      \param pmesh
+      an instance of `PolygonMesh` that represents a polygon mesh
+
+      \param np
+      a sequence of \ref bgl_namedparameters "Named Parameters"
+      among the ones listed below
+
+      \cgalNamedParamsBegin
+        \cgalParamNBegin{vertex_point_map}
+          \cgalParamDescription{an instance of `VertexToPointMap` that maps a polygon mesh
+          vertex to `Kernel::Point_3`}
+          \cgalParamDefault{`boost::get(CGAL::vertex_point, pmesh)`}
+        \cgalParamNEnd
+        \cgalParamNBegin{geom_traits}
+          \cgalParamDescription{an instance of `GeomTraits`}
+          \cgalParamDefault{`GeomTraits()`}
+        \cgalParamNEnd
+      \cgalNamedParamsEnd
+
+      \pre `faces(pmesh).size() > 0`
+    */
+    template<typename CGAL_NP_TEMPLATE_PARAMETERS>
+    Face_area_sorting(
+      const PolygonMesh& pmesh,
+      const CGAL_NP_CLASS& np = parameters::default_values())
+        : m_face_graph(pmesh)
+        , m_vpm(parameters::choose_parameter(parameters::get_parameter(np, internal_np::vertex_point),
+                                                             get_const_property_map(CGAL::vertex_point, pmesh)))
+        , m_traits(parameters::choose_parameter<GeomTraits>(parameters::get_parameter(np, internal_np::geom_traits)))
+    {
+      CGAL_precondition(faces(pmesh).size() > 0);
+
+      m_ordered.resize(faces(pmesh).size());
+
+      std::size_t index = 0;
+      for (Item item : faces(pmesh))
+        m_ordered[index++] = item;
+      m_scores.resize(m_ordered.size(), 0.);
+    }
+
+    /// @}
+
+    /// \name Sorting
+    /// @{
+
+    /*!
+      \brief sorts `Items` of input faces.
+    */
+    void sort() {
+      compute_scores();
+      CGAL_precondition(m_scores.size() > 0);
+
+      auto cmp = [this](const std::size_t i, const std::size_t j)
+      {
+        CGAL_precondition(i < m_scores.size());
+        CGAL_precondition(j < m_scores.size());
+        return m_scores[i] > m_scores[j];
+      };
+      std::vector<std::size_t> order(m_ordered.size());
+      std::iota(order.begin(), order.end(), 0);
+      std::sort(order.begin(), order.end(), cmp);
+
+      std::vector<Item> tmp(m_ordered.size());
+      for (std::size_t i = 0; i < m_ordered.size(); i++)
+        tmp[i] = m_ordered[order[i]];
+
+      m_ordered.swap(tmp);
+    }
+    /// @}
+
+    /// \name Access
+    /// @{
+
+    /*!
+      \brief returns an instance of `Seed_range` to access the ordered `Items`
+      of input faces.
+    */
+    const Seed_range &ordered() {
+      return m_ordered;
+    }
+    /// @}
+
+  private:
+    const PolygonMesh& m_face_graph;
+    VertexToPointMap m_vpm;
+    const GeomTraits m_traits;
+    Seed_range m_ordered;
+    std::vector<FT> m_scores;
+
+    void compute_scores()
+    {
+      auto squared_area = m_traits.compute_squared_area_3_object();
+      std::size_t idx = 0;
+      for (Item item : m_ordered)
+      {
+        halfedge_descriptor hd = halfedge(item, m_face_graph);
+        std::vector<typename GeomTraits::Point_3> pts;
+
+        for (vertex_descriptor v : vertices_around_face(hd,m_face_graph))
+          pts.push_back( get(m_vpm, v) );
+
+        if (pts.size()==3)
+          m_scores[idx++] = approximate_sqrt(squared_area(pts[0], pts[1], pts[2]));
+        else
+        {
+          std::vector<typename GeomTraits::Triangle_3> triangles;
+          internal::triangulate_face<GeomTraits>(pts, triangles);
+          for (const typename GeomTraits::Triangle_3& tr : triangles)
+            m_scores[idx++] += approximate_sqrt(squared_area(tr[0], tr[1], tr[2]));
+        }
+      }
+    }
+  };
+
+} // namespace Polygon_mesh
+} // namespace Shape_detection
+} // namespace CGAL
+
+#endif // CGAL_SHAPE_DETECTION_REGION_GROWING_POLYGON_MESH_FACE_AREA_SORTING_H
+

Shape_detection/include/CGAL/Shape_detection/Region_growing/Polygon_mesh/Least_squares_plane_fit_region.h

@@ -137,6 +137,13 @@ namespace Polygon_mesh {
           vertex to `Kernel::Point_3`}
           \cgalParamDefault{`boost::get(CGAL::vertex_point, pmesh)`}
         \cgalParamNEnd
+        \cgalParamNBegin{face_normal_map}
+          \cgalParamDescription{a property map associating normal vectors to the faces of `pmesh`.}
+          \cgalParamType{a class model of `ReadablePropertyMap` with `boost::graph_traits<PolygonMesh>::%face_descriptor`
+                         as key type and `GT::Vector_3` as value type, `GT` being the type of the parameter `geom_traits`.}
+          \cgalParamDefault{If this parameter is omitted, face normals will be estimated using crossproducts of vectors created
+                            from consecutive vertices of the face.}
+        \cgalParamNEnd
         \cgalParamNBegin{geom_traits}
           \cgalParamDescription{an instance of `GeomTraits`}
           \cgalParamDefault{`GeomTraits()`}
@@ -161,10 +168,13 @@ namespace Polygon_mesh {
     m_squared_distance_3(m_traits.compute_squared_distance_3_object()),
     m_scalar_product_3(m_traits.compute_scalar_product_3_object()),
     m_cross_product_3(m_traits.construct_cross_product_vector_3_object()),
-    m_face_normals( get(CGAL::dynamic_face_property_t<Vector_3>(), pmesh) ),
+    m_face_normals(get(CGAL::dynamic_face_property_t<Vector_3>(), pmesh)),
     m_face_triangulations( get(CGAL::dynamic_face_property_t<std::vector<Triangle_3>>(), pmesh) )
     {
 
+      static constexpr bool use_input_face_normal =
+        !parameters::is_default_parameter<CGAL_NP_CLASS, internal_np::face_normal_t>::value;
+
 #ifdef CGAL_SD_RG_USE_PMP
     auto get_face_normal = [this](Item face, const PolygonMesh& pmesh)
     {
@@ -193,8 +203,19 @@ namespace Polygon_mesh {
     };
 #endif
 
+      if constexpr (!use_input_face_normal)
+      {
+        for (const Item &i : faces(pmesh))
+          put(m_face_normals, i, get_face_normal(i, pmesh));
+      }
+      else
+      {
+        auto fnm = parameters::get_parameter(np, internal_np::face_normal);
+        for (const Item &i : faces(pmesh))
+          put(m_face_normals, i, get(fnm, i));
+      }
+
       for (const Item &i : faces(pmesh)) {
-        put(m_face_normals, i, get_face_normal(i, pmesh));
         std::vector<Point_3> pts;
         auto h = halfedge(i, pmesh);
         auto s = h;