Fallback mode using map if vector is too large (useful for sparse scenes)

Classification/include/CGAL/Classification/Feature/Echo_scatter.h

@@ -81,14 +81,15 @@ public:
     Image_float Scatter(grid.width(), grid.height());
     for (std::size_t j = 0; j < grid.height(); j++)
       for (std::size_t i = 0; i < grid.width(); i++)
+        if (grid.has_points(i,j))
           Scatter(i,j)=0;
 
     std::size_t square = (std::size_t)(0.5 * radius_neighbors / grid.resolution()) + 1;
 
-    for (std::size_t j = 0; j < grid.height(); j++){	
+    for (std::size_t j = 0; j < grid.height(); j++)
-      for (std::size_t i = 0; i < grid.width(); i++){
+      for (std::size_t i = 0; i < grid.width(); i++)
-						
+        if(grid.has_points(i,j))
-        if(grid.has_points(i,j)){
+        {
 
           std::size_t squareXmin = (i < square ? 0 : i - square);
           std::size_t squareXmax = (std::min) (grid.width()-1, i + square);
@@ -121,12 +122,8 @@ public:
           }
 					
           Scatter(i,j)=(float)NB_echo_sup/NB_echo_total;
-				
         }
   
-      }
-		
-    }
     for(std::size_t i = 0; i < input.size(); i++){
       std::size_t I= grid.x(i);
       std::size_t J= grid.y(i);

Classification/include/CGAL/Classification/Feature/Elevation.h

@@ -99,6 +99,7 @@ public:
 
     for (std::size_t j = 0; j < grid.height(); ++ j)
       for (std::size_t i = 0; i < grid.width(); ++ i)
+        if (grid.has_points(i,j))
         {
           float mean = 0.;
           std::size_t nb = 0;
@@ -119,8 +120,8 @@ public:
     Image_float dtm_x(grid.width(),grid.height());
     
     for (std::size_t j = 0; j < grid.height(); ++ j)
-    {
       for (std::size_t i = 0; i < grid.width(); ++ i)
+        if (grid.has_points(i,j))
         {
           std::size_t squareXmin = (i < square ? 0 : i - square);
           std::size_t squareXmax = (std::min)(grid.width() - 1, i + square);
@@ -135,7 +136,6 @@ public:
           std::nth_element (z.begin(), z.begin() + (z.size() / 10), z.end());
           dtm_x(i,j) = z[z.size() / 10];
         }
-    }
     dem.free();
 
 #ifdef CGAL_CLASSIFICATION_PRECOMPUTE_FEATURES
@@ -145,8 +145,8 @@ public:
 #endif
     
     for (std::size_t i = 0; i < grid.width(); ++ i)
-    {
       for (std::size_t j = 0; j < grid.height(); ++ j)
+        if (grid.has_points(i,j))
         {
           std::size_t squareYmin = (j < square ? 0 : j - square);
           std::size_t squareYmax = (std::min)(grid.height() - 1, j + square);
@@ -160,7 +160,6 @@ public:
           std::nth_element (z.begin(), z.begin() + (z.size() / 10), z.end());
           dtm(i,j) = z[z.size() / 10];
         }
-    }
     dtm_x.free();
 
 #ifdef CGAL_CLASSIFICATION_PRECOMPUTE_FEATURES

Classification/include/CGAL/Classification/Feature/Vertical_dispersion.h

@@ -99,6 +99,7 @@ public:
     
     for (std::size_t j = 0; j < grid.height(); j++)	
       for (std::size_t i = 0; i < grid.width(); i++)
+        if (grid.has_points(i,j))
           Dispersion(i,j)=0;
     
     std::size_t square = (std::size_t)(0.5 * radius_neighbors / grid.resolution()) + 1;

Classification/include/CGAL/Classification/Image.h

@@ -25,6 +25,8 @@
 
 #include <boost/shared_ptr.hpp>
 
+#define CGAL_CLASSIFICATION_IMAGE_SIZE_LIMIT 10000000
+
 namespace CGAL {
 namespace Classification {
 
@@ -33,9 +35,14 @@ namespace Classification {
 template <typename Type>
 class Image
 {
+  typedef std::vector<Type> Vector;
+  typedef std::map<std::size_t, Type> Map;
+  
   std::size_t m_width;
   std::size_t m_height;
-  boost::shared_ptr<std::vector<Type> > m_raw;
+  boost::shared_ptr<Vector> m_raw;
+  boost::shared_ptr<Map> m_sparse;
+  Type m_default;
 
   // Forbid using copy constructor
   Image (const Image&)
@@ -53,7 +60,12 @@ public:
       m_height (height)
   {
     if (m_width * m_height > 0)
-      m_raw = boost::shared_ptr<std::vector<Type> > (new std::vector<Type>(m_width * m_height));
+    {
+      if (m_width * m_height < CGAL_CLASSIFICATION_IMAGE_SIZE_LIMIT)
+        m_raw = boost::shared_ptr<Vector> (new Vector(m_width * m_height));
+      else
+        m_sparse = boost::shared_ptr<Map> (new Map());
+    }
   }
   
   ~Image ()
@@ -62,12 +74,14 @@ public:
 
   void free()
   {
-    m_raw = boost::shared_ptr<std::vector<Type> >();
+    m_raw = boost::shared_ptr<Vector>();
+    m_sparse = boost::shared_ptr<Map>();
   }
 
   Image& operator= (const Image& other)
   {
     m_raw = other.m_raw;
+    m_sparse = other.m_sparse;
     m_width = other.width();
     m_height = other.height();
     return *this;
@@ -78,12 +92,24 @@ public:
 
   Type& operator() (const std::size_t& x, const std::size_t& y)
   {
-    //    return m_raw[y * m_width + x];
+    if (m_raw == boost::shared_ptr<Vector>()) // sparse case
+    {
+      typename Map::iterator inserted = m_sparse->insert (std::make_pair (x * m_height + y, Type())).first;
+      return inserted->second;
+    }
+
     return (*m_raw)[x * m_height + y];
   }
   const Type& operator() (const std::size_t& x, const std::size_t& y) const
   {
-    //    return m_raw[y * m_width + x];
+    if (m_raw == boost::shared_ptr<Vector>()) // sparse case
+    {
+      typename Map::iterator found = m_sparse->find (x * m_height + y);
+      if (found != m_sparse->end())
+        return found->second;
+      return m_default;
+    }
+
     return (*m_raw)[x * m_height + y];
   }