New undocumented AABB primitive: cells of a 3D triangulation

The `Scene_c3t3_item` from our 3D demos now uses an AABB tree of
cells, instead of an AABB tree of triangles. That divides the number
of primitives by 4, and save times:

Before:
```
C3t3 facets AABB tree built in 22.010006904602051 wall-clock seconds
Scene_c3t3_item_priv::computeIntersections in 0.50893402099609375 wall-clock seconds
```

After:
```
C3t3 cells AABB tree built in 13.072829008102417 wall-clock seconds
Scene_c3t3_item_priv::computeIntersections in 0.41458892822265625 wall-clock seconds
```

The gain seems mostly in the construction of the tree (40% saved), and
the gain on the queries seems lower (about 20%).

AABB_tree/include/CGAL/AABB_triangulation_3_cell_primitive.h

@@ -0,0 +1,103 @@
+// Copyright (c) 2017 GeometryFactory (France).
+// All rights reserved.
+//
+// This file is part of CGAL (www.cgal.org).
+// You can redistribute it and/or modify it under the terms of the GNU
+// General Public License as published by the Free Software Foundation,
+// either version 3 of the License, or (at your option) any later version.
+//
+// Licensees holding a valid commercial license may use this file in
+// accordance with the commercial license agreement provided with the software.
+//
+// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
+// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: GPL-3.0+
+//
+//
+// Author        : Jane Tournois
+//
+
+#ifndef CGAL_AABB_TRIANGULATION_3_CELL_PRIMITIVE_H_
+#define CGAL_AABB_TRIANGULATION_3_CELL_PRIMITIVE_H_
+
+#include <CGAL/license/AABB_tree.h>
+
+
+#include <CGAL/AABB_primitive.h>
+#include <CGAL/result_of.h>
+#include <iterator>
+
+namespace CGAL
+{
+  namespace internal
+  {
+    template <class GeomTraits, class Iterator>
+    struct Point_from_cell_iterator_proprety_map
+    {
+      //classical typedefs
+      typedef Iterator key_type;
+      typedef typename GeomTraits::Point_3 value_type;
+      typedef typename cpp11::result_of<
+        typename GeomTraits::Construct_vertex_3(typename GeomTraits::Tetrahedron_3, int)
+      >::type reference;
+      typedef boost::readable_property_map_tag category;
+
+      inline friend
+        typename Point_from_cell_iterator_proprety_map<GeomTraits, Iterator>::reference
+        get(Point_from_cell_iterator_proprety_map<GeomTraits, Iterator>, Iterator it)
+      {
+        return it->vertex(1)->point().point();
+      }
+    };
+
+    template <class GeomTraits, class Iterator>
+    struct Tet_from_cell_iterator_proprety_map
+    {
+      //classical typedefs
+      typedef Iterator                           key_type;
+      typedef typename GeomTraits::Tetrahedron_3 value_type;
+      typedef value_type                         reference;
+      typedef boost::readable_property_map_tag category;
+
+      inline friend
+        reference
+        get(Tet_from_cell_iterator_proprety_map<GeomTraits, Iterator>, key_type it)
+      {
+        return value_type(it->vertex(0)->point().point(),
+                          it->vertex(1)->point().point(),
+                          it->vertex(2)->point().point(),
+                          it->vertex(3)->point().point());
+      }
+    };
+
+  }//namespace internal
+
+
+  template < class GeomTraits,
+             class Tr,
+             class CacheDatum = Tag_false,
+             class Handle = typename Tr::Cell_handle>
+  class AABB_triangulation_3_cell_primitive
+#ifndef DOXYGEN_RUNNING
+    : public AABB_primitive<  Handle,
+          internal::Tet_from_cell_iterator_proprety_map<GeomTraits, Handle>,
+          internal::Point_from_cell_iterator_proprety_map<GeomTraits, Handle>,
+          Tag_false,
+          CacheDatum >
+#endif
+  {
+    typedef AABB_primitive< Handle,
+          internal::Tet_from_cell_iterator_proprety_map<GeomTraits, Handle>,
+          internal::Point_from_cell_iterator_proprety_map<GeomTraits, Handle>,
+          Tag_false,
+          CacheDatum > Base;
+  public:
+    AABB_triangulation_3_cell_primitive(Handle h) : Base(h){}  };
+
+}  // end namespace CGAL
+
+
+#endif // CGAL_AABB_TRIANGULATION_3_CELL_PRIMITIVE_H_

Polyhedron/demo/Polyhedron/Scene_c3t3_item.cpp

@@ -28,13 +28,14 @@
 
 #include <CGAL/AABB_tree.h>
 #include <CGAL/AABB_traits.h>
-#include <CGAL/AABB_triangulation_3_triangle_primitive.h>
+#include <CGAL/AABB_triangulation_3_cell_primitive.h>
 #include <CGAL/IO/facets_in_complex_3_to_triangle_mesh.h>
 
 #include "Scene_polygon_soup_item.h"
 
 
-typedef CGAL::AABB_triangulation_3_triangle_primitive<EPICK,C3t3> Primitive;
+typedef CGAL::AABB_triangulation_3_cell_primitive<EPICK,
+                                                  C3t3::Triangulation> Primitive;
 typedef CGAL::AABB_traits<EPICK, Primitive> Traits;
 typedef CGAL::AABB_tree<Traits> Tree;
 typedef Tree::Point_and_primitive_id Point_and_primitive_id;
@@ -399,25 +400,19 @@ struct Scene_c3t3_item_priv {
     CGAL::Real_timer timer;
     timer.start();
     tree.clear();
-    for (Tr::Finite_facets_iterator
+    for (Tr::Finite_cells_iterator
-           fit = c3t3.triangulation().finite_facets_begin(),
+           cit = c3t3.triangulation().finite_cells_begin(),
-           end = c3t3.triangulation().finite_facets_end();
+           end = c3t3.triangulation().finite_cells_end();
-         fit != end; ++fit)
+         cit != end; ++cit)
     {
-      Tr::Cell_handle ch = fit->first, nh =ch->neighbor(fit->second);
+      Tr::Cell_handle ch = cit;
 
-      if( (!c3t3.is_in_complex(ch)) &&  (!c3t3.is_in_complex(nh)) )
+      if(!c3t3.is_in_complex(ch)) continue;
-        continue;
 
-      if(c3t3.is_in_complex(ch)){
+      tree.insert(Primitive(cit));
-        tree.insert(Primitive(fit));
-      } else{
-        int ni = nh->index(ch);
-        tree.insert(Primitive(Tr::Facet(nh,ni)));
-      }
     }
     tree.build();
-    std::cerr << "C3t3 facets AABB tree built in " << timer.time()
+    std::cerr << "C3t3 cells AABB tree built in " << timer.time()
               << " wall-clock seconds\n";
 
     is_aabb_tree_built = true;
@@ -517,7 +512,6 @@ struct Scene_c3t3_item_priv {
   typedef std::set<int> Indices;
   Indices surface_patch_indices_;
   Indices subdomain_indices_;
-  std::set<Tr::Cell_handle> intersected_cells;
   QSlider* tet_Slider;
 
   //!Allows OpenGL 2.0 context to get access to glDrawArraysInstanced.
@@ -1481,52 +1475,26 @@ void Scene_c3t3_item_priv::initializeBuffers(CGAL::Three::Viewer_interface *view
 
 
 
-void Scene_c3t3_item_priv::computeIntersection(const Primitive& facet)
+void Scene_c3t3_item_priv::computeIntersection(const Primitive& cell)
 {
   Geom_traits::Construct_point_3 wp2p
     = c3t3.triangulation().geom_traits().construct_point_3_object();
 
   typedef unsigned char UC;
-  Tr::Cell_handle ch = facet.id().first;
+  Tr::Cell_handle ch = cell.id();
-  if(intersected_cells.find(ch) == intersected_cells.end())
+  QColor c = this->colors_subdomains[ch->subdomain_index()].light(50);
-  {
-    QColor c = this->colors_subdomains[ch->subdomain_index()].light(50);
 
-    const Tr::Bare_point& pa = wp2p(ch->vertex(0)->point());
+  const Tr::Bare_point& pa = wp2p(ch->vertex(0)->point());
-    const Tr::Bare_point& pb = wp2p(ch->vertex(1)->point());
+  const Tr::Bare_point& pb = wp2p(ch->vertex(1)->point());
-    const Tr::Bare_point& pc = wp2p(ch->vertex(2)->point());
+  const Tr::Bare_point& pc = wp2p(ch->vertex(2)->point());
-    const Tr::Bare_point& pd = wp2p(ch->vertex(3)->point());
+  const Tr::Bare_point& pd = wp2p(ch->vertex(3)->point());
 
-    CGAL::Color color(UC(c.red()), UC(c.green()), UC(c.blue()));
+  CGAL::Color color(UC(c.red()), UC(c.green()), UC(c.blue()));
 
-    intersection->addTriangle(pb, pa, pc, color);
+  intersection->addTriangle(pb, pa, pc, color);
-    intersection->addTriangle(pa, pb, pd, color);
+  intersection->addTriangle(pa, pb, pd, color);
-    intersection->addTriangle(pa, pd, pc, color);
+  intersection->addTriangle(pa, pd, pc, color);
-    intersection->addTriangle(pb, pc, pd, color);
+  intersection->addTriangle(pb, pc, pd, color);
-    intersected_cells.insert(ch);
-  }
-  {
-    Tr::Cell_handle nh = ch->neighbor(facet.id().second);
-    if(c3t3.is_in_complex(nh)){
-      if(intersected_cells.find(nh) == intersected_cells.end())
-      {
-        const Tr::Bare_point& pa = wp2p(nh->vertex(0)->point());
-        const Tr::Bare_point& pb = wp2p(nh->vertex(1)->point());
-        const Tr::Bare_point& pc = wp2p(nh->vertex(2)->point());
-        const Tr::Bare_point& pd = wp2p(nh->vertex(3)->point());
-
-        QColor c = this->colors_subdomains[nh->subdomain_index()].light(50);
-
-        CGAL::Color color(UC(c.red()), UC(c.green()), UC(c.blue()));
-
-        intersection->addTriangle(pb, pa, pc, color);
-        intersection->addTriangle(pa, pb, pd, color);
-        intersection->addTriangle(pa, pd, pc, color);
-        intersection->addTriangle(pb, pc, pd, color);
-        intersected_cells.insert(nh);
-      }
-    }
-  }
 }
 
 struct ComputeIntersection {
@@ -1555,7 +1523,6 @@ void Scene_c3t3_item_priv::computeIntersections()
   const Geom_traits::Plane_3& plane = item->plane(offset);
   tree.all_intersected_primitives(plane,
         boost::make_function_output_iterator(ComputeIntersection(*this)));
-  intersected_cells.clear();
 }
 
 void Scene_c3t3_item_priv::computeSpheres()