Merge pull request #2690 from sgiraudot/Classification-ETHZ_random_forest-GF

Classification: ETHZ Random Forest

Algebraic_kernel_d/include/CGAL/Algebraic_kernel_d/Curve_pair_analysis_2.h

@@ -32,7 +32,7 @@
 #include <boost/optional.hpp>
 
 #include <CGAL/Handle_with_policy.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <CGAL/Polynomial.h>
 #include <CGAL/Polynomial_traits_d.h>

Algebraic_kernel_d/include/CGAL/Algebraic_kernel_d/shear.h

@@ -32,7 +32,7 @@
 #include <vector>
 #include <functional>
 #include <iterator>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <CGAL/boost/iterator/counting_iterator.hpp>
 
 namespace CGAL {

BGL/examples/BGL_LCC/transform_iterator_lcc.cpp

@@ -4,7 +4,7 @@
 #include <CGAL/boost/graph/iterator.h>
 #include <fstream>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <algorithm>
 
 typedef CGAL::Simple_cartesian<double>              Kernel;

BGL/examples/BGL_polyhedron_3/transform_iterator.cpp

@@ -2,7 +2,7 @@
 #include <CGAL/Polyhedron_3.h>
 
 #include <CGAL/boost/graph/iterator.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <fstream>
 #include <algorithm>

BGL/include/CGAL/boost/graph/Face_filtered_graph.h

@@ -30,7 +30,7 @@
 #include <boost/foreach.hpp>
 #include <boost/unordered_set.hpp>
 #include <boost/graph/graph_traits.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/foreach.hpp>
 #include <boost/iterator/filter_iterator.hpp>
 #include <boost/dynamic_bitset.hpp>

BGL/include/CGAL/boost/graph/Graph_with_descriptor_with_graph.h

@@ -25,7 +25,7 @@
 #include <CGAL/boost/graph/properties.h>
 #include <boost/graph/graph_traits.hpp>
 #include <CGAL/boost/graph/iterator.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <CGAL/boost/graph/Graph_with_descriptor_with_graph_fwd.h>
 

BGL/include/CGAL/boost/graph/graph_traits_OpenMesh.h

@@ -24,7 +24,7 @@
 #include <boost/graph/graph_traits.hpp>
 #include <boost/graph/properties.hpp>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <CGAL/boost/graph/internal/OM_iterator_from_circulator.h>
 #include <CGAL/boost/graph/iterator.h>

BGL/include/CGAL/boost/graph/graph_traits_Seam_mesh.h

@@ -24,7 +24,7 @@
 #include <boost/graph/graph_traits.hpp>
 #include <boost/graph/properties.hpp>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <CGAL/boost/graph/properties_Seam_mesh.h>
 #include <CGAL/boost/graph/iterator.h>

Bounding_volumes/include/CGAL/Min_annulus_d.h

@@ -40,7 +40,7 @@
 #include <CGAL/QP_solver/QP_full_filtered_pricing.h>
 #include <CGAL/QP_solver/QP_full_exact_pricing.h>
 #include <CGAL/boost/iterator/counting_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/functional.hpp>
 
 // here is how it works. We have d+2 variables: 

Classification/doc/Classification/Classification.txt

@@ -8,6 +8,10 @@ namespace CGAL {
 
 This component implements the algorithm described in \cgalCite{cgal:lm-clscm-12} (section 2), generalized to handle different types of data, multiple features and multiple labels. It classifies a data set into a user-defined set of labels, such as _ground_, _vegetation_ and _buildings_. A flexible API is provided so that users can classify any type of data which they can index and for which they can compute relevant features, compute their own local features on the input data set and define their own labels.
 
+\note This component depends on the Boost libraries
+[Serialization](http://www.boost.org/libs/serialization) and
+[IO Streams](http://www.boost.org/libs/iostreams) (compiled with the GZIP dependency).
+
 \section Classification_Organization Package Organization
 
 %Classification of data sets is achieved as follows (see Figure \cgalFigureRef{Classification_organization_fig}):
@@ -80,7 +84,9 @@ For more details about how these different features can help to identify one lab
 Finally, if the input data set has additional properties, these can also be used as features. For example, \cgal provides the following features:
 
 - [Echo_scatter](@ref CGAL::Classification::Feature::Echo_scatter) uses the number of returns (echo) provided by most LIDAR scanners if available;
-- [Hsv](@ref CGAL::Classification::Feature::Hsv) uses input color information if available.
+- [Hsv](@ref CGAL::Classification::Feature::Hsv) uses input color information if available;
+- [Simple_feature](@ref CGAL::Classification::Feature::Simple_feature) uses any property map applicable to the input range and whose value type is castable to `float` (useful if an additional property of the input set should be used as is, for example and `intensity` measurement).
+
 
 In the following code snippet, a subset of these features are instantiated. Note that all the predefined features can also be automatically generated in multiple scales (see \ref Classification_feature_generator).
 
@@ -117,7 +123,7 @@ The following snippet shows how to use the feature generator:
 
 %Classification relies on a classifier: this classifier is an object that, from the set of values taken by the features at an input item, computes the energy that measures the likelihood of an input item to belong to one label or another. A model of the concept `CGAL::Classification::Classifier` must take the index of an input item and store the energies associated to each label in a vector. If a classifier returns the value 0 for a pair of label and input item, it means that this item belongs to this label with certainty; large values mean that this item is not likely to belong to this label.
 
-\cgal provides two models for this concept, [Sum_of_weighted_features_classifier](@ref CGAL::Classification::Sum_of_weighted_features_classifier) and [Random_forest_classifier](@ref CGAL::Classification::Random_forest_classifier).
+\cgal provides three models for this concept, [Sum_of_weighted_features_classifier](@ref CGAL::Classification::Sum_of_weighted_features_classifier), [ETHZ_random_forest_classifier](@ref CGAL::Classification::ETHZ_random_forest_classifier) and [OpenCV_random_forest_classifier](@ref CGAL::Classification::OpenCV_random_forest_classifier).
 
 To perform classification based on these classifiers, please refer to \ref Classification_classification_functions.
 
@@ -194,22 +200,45 @@ Figure \cgalFigureRef{Classification_sowf_result_fig} shows an example of output
 Example of classification on a point set with medium noise and outliers (left: input, right: output). _Ground_ is orange, _roofs_ are pink, _vegetation_ is green. Outliers are classified with an additional label _outlier_ in black.
 \cgalFigureEnd
 
+\subsection Classification_ETHZ_random_forest ETHZ Random Forest
 
-\subsection Classification_random_forest Random Forest
+\cgal provides [ETHZ_random_forest_classifier](@ref CGAL::Classification::ETHZ_random_forest_classifier),
-
+a classifier based on the Random Forest Template Library developed by
-This second classifier is [Random_forest_classifier](@ref CGAL::Classification::Random_forest_classifier).
+Stefan Walk at ETH Zurich \cgalCite{cgal:w-erftl-14} (the library is
-It uses the \ref thirdpartyOpenCV library, more specifically the
+distributed under the MIT license and is included with the \cgal release,
-[Random Trees](http://docs.opencv.org/2.4/modules/ml/doc/random_trees.html)
+the user does not have to install anything more). This classifier uses
-package. This classifier uses a ground truth training set to construct several
+a ground truth training set to construct several decision trees that
-decision trees that are then used to assign a label to each input
+are then used to assign a label to each input item.
-item.
 
 This classifier cannot be set up by hand and requires a ground truth
-training set. The training algorithm is faster but usually requires a
+training set. The training algorithm is significantly faster but
-higher number of inliers than the previous classifier.
+usually requires a higher number of inliers than the previous
+classifier. The training algorithm uses more memory at runtime and the
+configuration files are larger than those produced by
+[Sum_of_weighted_features_classifier](@ref CGAL::Classification::Sum_of_weighted_features_classifier),
+but the output quality is usually significantly better, especially in
+the cases where many labels are used (more than five).
 
-An [example](\ref Classification_example_random_forest) shows how to
+An [example](\ref Classification_example_ethz_random_forest) shows how to
-use this classifier.  For more details about the method, please refer
+use this classifier. For more details about the algorithm, please refer
+to README provided in the [ETH Zurich's code archive](https://www.ethz.ch/content/dam/ethz/special-interest/baug/igp/photogrammetry-remote-sensing-dam/documents/sourcecode-and-datasets/Random%20Forest/rforest.zip).
+
+\subsection Classification_OpenCV_random_forest OpenCV Random Forest
+
+This last classifier is [OpenCV_random_forest_classifier](@ref CGAL::Classification::OpenCV_random_forest_classifier).
+It uses the \ref thirdpartyOpenCV library, more specifically the
+[Random Trees](http://docs.opencv.org/2.4/modules/ml/doc/random_trees.html)
+package.
+
+Note that this classifier usually produces results with a lower
+quality than [ETHZ_random_forest_classifier](@ref CGAL::Classification::ETHZ_random_forest_classifier).
+
+It is provided for the sake of completeness and for testing purposes,
+but if you are not sure what to use, we advise using the ETHZ Random
+Forest instead.
+
+An [example](\ref Classification_example_opencv_random_forest) shows how to
+use this classifier. For more details about the algorithm, please refer
 to [the official documentation](http://docs.opencv.org/2.4/modules/ml/doc/random_trees.html)
 of OpenCV.
 
@@ -373,11 +402,17 @@ The following example:
 
 \cgalExample{Classification/example_generation_and_training.cpp}
 
-\subsection Classification_example_random_forest Random Forest
+\subsection Classification_example_ethz_random_forest ETHZ Random Forest
 
-The following example shows how to use the classifier [Random_forest_classifier](@ref CGAL::Classification::Random_forest_classifier) using an input training set.
+The following example shows how to use the classifier [ETHZ_random_forest_classifier](@ref CGAL::Classification::ETHZ_random_forest_classifier) using an input training set.
 
-\cgalExample{Classification/example_random_forest.cpp}
+\cgalExample{Classification/example_ethz_random_forest.cpp}
+
+\subsection Classification_example_opencv_random_forest OpenCV Random Forest
+
+The following example shows how to use the classifier [OpenCV_random_forest_classifier](@ref CGAL::Classification::OpenCV_random_forest_classifier) using an input training set.
+
+\cgalExample{Classification/example_opencv_random_forest.cpp}
 
 \section Classification_history History
 

Classification/doc/Classification/Concepts/Classifier.h

@@ -13,7 +13,8 @@ Concept describing a classifier used by classification functions (see
 `CGAL::Classification::classify_with_graphcut()`).
 
 \cgalHasModel `CGAL::Classification::Sum_of_weighted_features_classifier`
-\cgalHasModel `CGAL::Classification::Random_forest_classifier`
+\cgalHasModel `CGAL::Classification::ETHZ_random_forest_classifier`
+\cgalHasModel `CGAL::Classification::OpenCV_random_forest_classifier`
 
 */
 class Classifier

Classification/doc/Classification/PackageDescription.txt

@@ -47,6 +47,8 @@ Features are defined as scalar fields that associates each input item with a spe
 
 \cgalPkgShortInfoBegin
 \cgalPkgSince{4.9}
+\cgalPkgDependsOn{\ref PkgSolverSummary, \ref PkgSpatialSearchingDSummary, [Boost Serialization](http://www.boost.org/libs/serialization) and
+[Boost IO Streams](http://www.boost.org/libs/iostreams)}
 \cgalPkgBib{cgal:lm-clscm-12}
 \cgalPkgLicense{\ref licensesGPL "GPL"}
 \cgalPkgDemo{Operations on Polyhedra,polyhedron_3.zip}
@@ -70,7 +72,8 @@ Features are defined as scalar fields that associates each input item with a spe
 ## Classifiers ##
 
 - `CGAL::Classification::Sum_of_weighted_features_classifier`
-- `CGAL::Classification::Random_forest_classifier`
+- `CGAL::Classification::ETHZ_random_forest_classifier`
+- `CGAL::Classification::OpenCV_random_forest_classifier`
 
 ## Data Structures ##
 
@@ -103,6 +106,7 @@ Features are defined as scalar fields that associates each input item with a spe
 - `CGAL::Classification::Feature::Linearity`
 - `CGAL::Classification::Feature::Omnivariance`
 - `CGAL::Classification::Feature::Planarity`
+- `CGAL::Classification::Feature::Simple_feature<InputRange, PropertyMap>`
 - `CGAL::Classification::Feature::Sphericity`
 - `CGAL::Classification::Feature::Sum_eigenvalues`
 - `CGAL::Classification::Feature::Surface_variation`

Classification/doc/Classification/examples.txt

@@ -2,5 +2,6 @@
 \example Classification/example_classification.cpp
 \example Classification/example_feature.cpp
 \example Classification/example_generation_and_training.cpp
-\example Classification/example_random_forest.cpp
+\example Classification/example_ethz_random_forest.cpp
+\example Classification/example_opencv_random_forest.cpp
 */

Classification/examples/Classification/CMakeLists.txt

@@ -10,10 +10,8 @@ cmake_minimum_required(VERSION 2.8.11)
 find_package( CGAL QUIET COMPONENTS  )
 
 if ( NOT CGAL_FOUND )
-
   message(STATUS "This project requires the CGAL library, and will not be compiled.")
   return()  
-
 endif()
 
 # include helper file
@@ -24,14 +22,20 @@ include( ${CGAL_USE_FILE} )
 find_package( Boost REQUIRED )
 
 if ( NOT Boost_FOUND )
-
   message(STATUS "This project requires the Boost library, and will not be compiled.")
-
   return()  
-
 endif()
 
-find_package( TBB )
+find_package( Boost OPTIONAL_COMPONENTS serialization iostreams )
+
+if( WIN32 )
+# to avoid a warning with old cmake 
+  set(_Boost_BZIP2_HEADERS             "boost/iostreams/filter/bzip2.hpp")
+  set(_Boost_ZLIB_HEADERS              "boost/iostreams/filter/zlib.hpp")
+  find_package( Boost OPTIONAL_COMPONENTS zlib)
+endif()
+
+find_package(TBB QUIET)
 
 find_package(OpenCV QUIET)
 
@@ -47,19 +51,67 @@ include_directories( BEFORE ../../include )
 
 include( CGAL_CreateSingleSourceCGALProgram )
 
+# Classification examples
+set(targets
+  example_classification
+  example_ethz_random_forest
+  example_feature
+  example_generation_and_training)
+
+# Classification requires some C++11 features
 set(needed_cxx_features cxx_rvalue_references cxx_variadic_templates)
-create_single_source_cgal_program( "example_classification.cpp" CXX_FEATURES ${needed_cxx_features} )
-if(TBB_FOUND)
-  CGAL_target_use_TBB(example_classification)
-endif()
-create_single_source_cgal_program( "example_generation_and_training.cpp" CXX_FEATURES ${needed_cxx_features} )
-create_single_source_cgal_program( "example_feature.cpp" CXX_FEATURES ${needed_cxx_features} )
 
-if( OpenCV_FOUND )
+# Libraries and flags
-  include_directories( ${OpenCV_INCLUDE_DIRS} )
+set(classification_linked_libraries)
-  create_single_source_cgal_program( "example_random_forest.cpp" CXX_FEATURES ${needed_cxx_features} )
+set(classification_compile_definitions)
-  target_link_libraries( example_random_forest PRIVATE ${OpenCV_LIBS} )
+
+
+if (Boost_SERIALIZATION_FOUND AND Boost_IOSTREAMS_FOUND)
+  set(classification_linked_libraries ${classification_linked_libraries}
+    ${Boost_SERIALIZATION_LIBRARY}
+    ${Boost_IOSTREAMS_LIBRARY})
 else()
-  message(STATUS "OpenCV not found, random forest example won't be compiled.")
+  message(STATUS "NOTICE: Boost Serialization or IO Streams not found, no example will be compiled.")
+  return()
 endif()
 
+if( WIN32 )
+  if (Boost_ZLIB_FOUND)
+    set(classification_linked_libraries ${classification_linked_libraries}
+      ${Boost_ZLIB_LIBRARY})
+  else()
+    message(STATUS "NOTICE: Boost ZLIB not found, no example will be compiled.")
+    return()
+  endif()
+endif()
+
+find_package(OpenCV QUIET)
+if (OpenCV_FOUND)
+  message(STATUS "Found OpenCV ${OpenCV_VERSION}")
+  include_directories(${OpenCV_INCLUDE_DIRS})
+  set(classification_linked_libraries ${classification_linked_libraries}
+    ${OpenCV_LIBS})
+  set(classification_compile_definitions ${classification_compile_definitions}
+    "-DCGAL_LINKED_WITH_OPENCV")
+
+  set(targets ${targets} example_opencv_random_forest)
+else()
+  message(STATUS "NOTICE: OpenCV was not found. OpenCV random forest predicate for classification won't be available.")
+endif()
+
+# Creating targets with correct libraries and flags
+foreach(target ${targets})
+  create_single_source_cgal_program( "${target}.cpp" CXX_FEATURES ${needed_cxx_features} )
+  if(TARGET ${target})
+    target_link_libraries(${target} PUBLIC ${classification_linked_libraries})
+    target_compile_definitions(${target} PUBLIC ${classification_compile_definitions})
+    if(TBB_FOUND)
+      CGAL_target_use_TBB(${target})
+    endif()
+  endif()
+endforeach()
+
+
+
+
+

Classification/examples/Classification/example_ethz_random_forest.cpp

@@ -10,7 +10,6 @@
 
 #include <CGAL/Simple_cartesian.h>
 #include <CGAL/Classification.h>
-#include <CGAL/Classification/Random_forest_classifier.h>
 #include <CGAL/Point_set_3.h>
 #include <CGAL/Point_set_3/IO.h>
 
@@ -32,14 +31,16 @@ typedef Classification::Feature_handle
 typedef Classification::Label_set                                               Label_set;
 typedef Classification::Feature_set                                             Feature_set;
 
-typedef Classification::Random_forest_classifier Classifier;
-
 typedef Classification::Point_set_feature_generator<Kernel, Point_set, Pmap>    Feature_generator;
 
 
 int main (int argc, char** argv)
 {
-  std::string filename (argc > 1 ? argv[1] : "data/b9_training.ply");
+  std::string filename = "data/b9_training.ply";
+  
+  if (argc > 1)
+    filename = argv[1];
+
   std::ifstream in (filename.c_str(), std::ios::binary);
   Point_set pts;
 
@@ -76,7 +77,10 @@ int main (int argc, char** argv)
   Label_handle vegetation = labels.add ("vegetation");
   Label_handle roof = labels.add ("roof");
 
-  Classifier classifier (labels, features);
+  std::vector<int> label_indices(pts.size(), -1);
+  
+  std::cerr << "Using ETHZ Random Forest Classifier" << std::endl;
+  Classification::ETHZ_random_forest_classifier classifier (labels, features);
   
   std::cerr << "Training" << std::endl;
   t.reset();
@@ -87,12 +91,12 @@ int main (int argc, char** argv)
 
   t.reset();
   t.start();
-  std::vector<int> label_indices(pts.size(), -1);
   Classification::classify_with_graphcut<CGAL::Sequential_tag>
     (pts, pts.point_map(), labels, classifier,
      generator.neighborhood().k_neighbor_query(12),
      0.2f, 1, label_indices);
   t.stop();
+  
   std::cerr << "Classification with graphcut done in " << t.time() << " second(s)" << std::endl;
 
   std::cerr << "Precision, recall, F1 scores and IoU:" << std::endl;

Classification/examples/Classification/example_opencv_random_forest.cpp

@@ -0,0 +1,151 @@
+#if defined (_MSC_VER) && !defined (_WIN64)
+#pragma warning(disable:4244) // boost::number_distance::distance()
+                              // converts 64 to 32 bits integers
+#endif
+
+#include <cstdlib>
+#include <fstream>
+#include <iostream>
+#include <string>
+
+#include <CGAL/Simple_cartesian.h>
+#include <CGAL/Classification.h>
+#include <CGAL/Point_set_3.h>
+#include <CGAL/Point_set_3/IO.h>
+
+#include <CGAL/Real_timer.h>
+
+typedef CGAL::Simple_cartesian<double> Kernel;
+typedef Kernel::Point_3 Point;
+typedef CGAL::Point_set_3<Point> Point_set;
+typedef Kernel::Iso_cuboid_3 Iso_cuboid_3;
+
+typedef Point_set::Point_map Pmap;
+typedef Point_set::Property_map<int> Imap;
+typedef Point_set::Property_map<unsigned char> UCmap;
+
+namespace Classification = CGAL::Classification;
+
+typedef Classification::Label_handle                                            Label_handle;
+typedef Classification::Feature_handle                                          Feature_handle;
+typedef Classification::Label_set                                               Label_set;
+typedef Classification::Feature_set                                             Feature_set;
+
+typedef Classification::Point_set_feature_generator<Kernel, Point_set, Pmap>    Feature_generator;
+
+
+int main (int argc, char** argv)
+{
+  std::string filename = "data/b9_training.ply";
+  
+  if (argc > 1)
+    filename = argv[1];
+
+  std::ifstream in (filename.c_str(), std::ios::binary);
+  Point_set pts;
+
+  std::cerr << "Reading input" << std::endl;
+  in >> pts;
+  
+  Imap label_map;
+  bool lm_found = false;
+  boost::tie (label_map, lm_found) = pts.property_map<int> ("label");
+  if (!lm_found)
+  {
+    std::cerr << "Error: \"label\" property not found in input file." << std::endl;
+    return EXIT_FAILURE;
+  }
+
+  std::vector<int> ground_truth;
+  ground_truth.reserve (pts.size());
+  std::copy (pts.range(label_map).begin(), pts.range(label_map).end(),
+             std::back_inserter (ground_truth));
+
+  Feature_set features;
+  
+  std::cerr << "Generating features" << std::endl;
+  CGAL::Real_timer t;
+  t.start();
+  Feature_generator generator (features, pts, pts.point_map(),
+                               5);  // using 5 scales
+  t.stop();
+  std::cerr << "Done in " << t.time() << " second(s)" << std::endl;
+
+  // Add types
+  Label_set labels;
+  Label_handle ground = labels.add ("ground");
+  Label_handle vegetation = labels.add ("vegetation");
+  Label_handle roof = labels.add ("roof");
+
+  std::vector<int> label_indices(pts.size(), -1);
+  
+  std::cerr << "Using OpenCV Random Forest Classifier" << std::endl;
+  Classification::OpenCV_random_forest_classifier classifier (labels, features);
+  
+  std::cerr << "Training" << std::endl;
+  t.reset();
+  t.start();
+  classifier.train (ground_truth);
+  t.stop();
+  std::cerr << "Done in " << t.time() << " second(s)" << std::endl;
+
+  t.reset();
+  t.start();
+  Classification::classify_with_graphcut<CGAL::Sequential_tag>
+    (pts, pts.point_map(), labels, classifier,
+     generator.neighborhood().k_neighbor_query(12),
+     0.2f, 1, label_indices);
+  t.stop();
+  
+  std::cerr << "Classification with graphcut done in " << t.time() << " second(s)" << std::endl;
+
+  std::cerr << "Precision, recall, F1 scores and IoU:" << std::endl;
+  Classification::Evaluation evaluation (labels, ground_truth, label_indices);
+  
+  for (std::size_t i = 0; i < labels.size(); ++ i)
+  {
+    std::cerr << " * " << labels[i]->name() << ": "
+              << evaluation.precision(labels[i]) << " ; "
+              << evaluation.recall(labels[i]) << " ; "
+              << evaluation.f1_score(labels[i]) << " ; "
+              << evaluation.intersection_over_union(labels[i]) << std::endl;
+  }
+
+  std::cerr << "Accuracy = " << evaluation.accuracy() << std::endl
+            << "Mean F1 score = " << evaluation.mean_f1_score() << std::endl
+            << "Mean IoU = " << evaluation.mean_intersection_over_union() << std::endl;
+
+  // Color point set according to class
+  UCmap red = pts.add_property_map<unsigned char>("red", 0).first;
+  UCmap green = pts.add_property_map<unsigned char>("green", 0).first;
+  UCmap blue = pts.add_property_map<unsigned char>("blue", 0).first;
+
+  for (std::size_t i = 0; i < label_indices.size(); ++ i)
+  {
+    label_map[i] = label_indices[i]; // update label map with computed classification
+    
+    Label_handle label = labels[label_indices[i]];
+
+    if (label == ground)
+    {
+      red[i] = 245; green[i] = 180; blue[i] =   0;
+    }
+    else if (label == vegetation)
+    {
+      red[i] =   0; green[i] = 255; blue[i] =  27;
+    }
+    else if (label == roof)
+    {
+      red[i] = 255; green[i] =   0; blue[i] = 170;
+    }
+  }
+
+  // Write result
+  std::ofstream f ("classification.ply");
+  f.precision(18);
+  f << pts;
+
+  std::cerr << "All done" << std::endl;
+  
+  return EXIT_SUCCESS;
+}

Classification/include/CGAL/Classification.h

@@ -25,9 +25,10 @@
 
 #include <CGAL/Classification/classify.h>
 #include <CGAL/Classification/Sum_of_weighted_features_classifier.h>
+#include <CGAL/Classification/ETHZ_random_forest_classifier.h>
 
 #ifdef CGAL_LINKED_WITH_OPENCV
-#include <CGAL/Classification/Random_forest_classifier.h>
+#include <CGAL/Classification/OpenCV_random_forest_classifier.h>
 #endif
 
 #include <CGAL/Classification/Color.h>
@@ -46,6 +47,7 @@
 #include <CGAL/Classification/Feature/Eigen.h>
 #include <CGAL/Classification/Feature/Elevation.h>
 #include <CGAL/Classification/Feature/Hsv.h>
+#include <CGAL/Classification/Feature/Simple_feature.h>
 #include <CGAL/Classification/Feature/Vertical_dispersion.h>
 #include <CGAL/Classification/Feature/Verticality.h>
 

Classification/include/CGAL/Classification/ETHZ_random_forest_classifier.h

@@ -0,0 +1,250 @@
+// Copyright (c) 2017 GeometryFactory Sarl (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(s)     : Simon Giraudot
+
+#ifndef CGAL_CLASSIFICATION_ETHZ_RANDOM_FOREST_CLASSIFIER_H
+#define CGAL_CLASSIFICATION_ETHZ_RANDOM_FOREST_CLASSIFIER_H
+
+#include <CGAL/license/Classification.h>
+
+#include <CGAL/Classification/Feature_set.h>
+#include <CGAL/Classification/Label_set.h>
+
+#ifdef CGAL_CLASSIFICATION_VERBOSE
+#define VERBOSE_TREE_PROGRESS 1
+#endif
+
+// Disable warnings from auxiliary library
+#ifdef BOOST_MSVC
+#  pragma warning(push)
+#  pragma warning(disable:4141)
+#  pragma warning(disable:4244)
+#  pragma warning(disable:4267)
+#  pragma warning(disable:4275)
+#  pragma warning(disable:4251)
+#endif
+
+#include <CGAL/Classification/internal/auxiliary/random-forest/node-gini.hpp>
+#include <CGAL/Classification/internal/auxiliary/random-forest/forest.hpp>
+
+#include <boost/archive/text_iarchive.hpp>
+#include <boost/archive/text_oarchive.hpp>
+#include <boost/iostreams/filtering_stream.hpp>
+#include <boost/iostreams/filter/gzip.hpp>
+
+#ifdef BOOST_MSVC
+#  pragma warning(pop)
+#endif
+
+
+namespace CGAL {
+
+namespace Classification {
+
+/*!
+  \ingroup PkgClassificationClassifiers
+
+  \brief %Classifier based on the ETH Zurich version of random forest algorithm \cgalCite{cgal:w-erftl-14}.
+
+  \note This classifier is distributed under the MIT license.
+
+  \cgalModels `CGAL::Classification::Classifier`
+*/
+class ETHZ_random_forest_classifier
+{
+  typedef CGAL::internal::liblearning::RandomForest::RandomForest
+  < CGAL::internal::liblearning::RandomForest::NodeGini
+    < CGAL::internal::liblearning::RandomForest::AxisAlignedSplitter> > Forest;
+  
+  const Label_set& m_labels;
+  const Feature_set& m_features;
+  Forest* m_rfc;
+
+public:
+  
+  /// \name Constructor
+  /// @{
+  
+  /*!
+    \brief Instantiate the classifier using the sets of `labels` and `features`.
+
+  */
+  ETHZ_random_forest_classifier (const Label_set& labels,
+                                 const Feature_set& features)
+    : m_labels (labels), m_features (features), m_rfc (NULL)
+  { }
+  
+  /// \cond SKIP_IN_MANUAL
+  ~ETHZ_random_forest_classifier ()
+  {
+    if (m_rfc != NULL)
+      delete m_rfc;
+  }
+  /// \endcond
+  
+  /// @}
+
+  /// \name Training
+  
+  /// @{
+  /*!
+    \brief Runs the training algorithm.
+
+    From the set of provided ground truth, this algorithm estimates
+    sets up the random trees that produce the most accurate result
+    with respect to this ground truth.
+
+    \pre At least one ground truth item should be assigned to each
+    label.
+
+    \param ground_truth vector of label indices. It should contain for
+    each input item, in the same order as the input set, the index of
+    the corresponding label in the `Label_set` provided in the
+    constructor. Input items that do not have a ground truth
+    information should be given the value `-1`.
+
+    \param reset_trees should be set to `false` if the users wants to
+    _add_ new trees to the existing forest, and kept to `true` if the
+    training should be recomputing from scratch (discarding the
+    current forest).
+
+    \param num_trees number of trees generated by the training
+    algorithm. Higher values may improve result at the cost of higher
+    computation times (in general, using a few dozens of trees is
+    enough).
+
+    \param max_depth maximum depth of the trees. Higher values will
+    improve how the forest fits the training set. A overly low value
+    will underfit the test data and conversely an overly high value
+    will likely overfit.
+  */
+  template <typename LabelIndexRange>
+  void train (const LabelIndexRange& ground_truth,
+              bool reset_trees = true,
+              std::size_t num_trees = 25,
+              std::size_t max_depth = 20)
+  {
+    CGAL::internal::liblearning::RandomForest::ForestParams params;
+    params.n_trees   = num_trees;
+    params.max_depth = max_depth;
+
+    std::vector<int> gt;
+    std::vector<float> ft;
+    
+    for (std::size_t i = 0; i < ground_truth.size(); ++ i)
+      if (ground_truth[i] != -1)
+      {
+        for (std::size_t f = 0; f < m_features.size(); ++ f)
+          ft.push_back(m_features[f]->value(i));
+        gt.push_back(ground_truth[i]);
+      }
+
+    CGAL::internal::liblearning::DataView2D<int> label_vector (&(gt[0]), gt.size(), 1);    
+    CGAL::internal::liblearning::DataView2D<float> feature_vector(&(ft[0]), gt.size(), ft.size() / gt.size());
+
+    if (m_rfc != NULL && reset_trees)
+    {
+      delete m_rfc;
+      m_rfc = NULL;
+    }
+    
+    if (m_rfc == NULL)
+      m_rfc = new Forest (params);
+
+    CGAL::internal::liblearning::RandomForest::AxisAlignedRandomSplitGenerator generator;
+    
+    m_rfc->train(feature_vector, label_vector, CGAL::internal::liblearning::DataView2D<int>(), generator, 0, false, reset_trees);
+  }
+
+  /// \cond SKIP_IN_MANUAL
+  void operator() (std::size_t item_index, std::vector<float>& out) const
+  {
+    out.resize (m_labels.size(), 0.);
+    
+    std::vector<float> ft;
+    ft.reserve (m_features.size());
+    for (std::size_t f = 0; f < m_features.size(); ++ f)
+      ft.push_back (m_features[f]->value(item_index));
+
+    std::vector<float> prob (m_labels.size());
+
+    m_rfc->evaluate (ft.data(), prob.data());
+    
+    for (std::size_t i = 0; i < out.size(); ++ i)
+    {
+      out[i] = - std::log (prob[i]);
+      if (out[i] < 0.f)
+        out[i] = -out[i];
+    }
+  }
+  /// \endcond
+  
+  /// @}
+
+  /// \name Input/Output
+  /// @{
+  
+  /*!
+    \brief Saves the current configuration in the stream `output`.
+
+    This allows to easily save and recover a specific classification
+    configuration.
+
+    The output file is written in an GZIP container that is readable
+    by the `load_configuration()` method.
+  */
+  void save_configuration (std::ostream& output)
+  {
+    boost::iostreams::filtering_ostream outs;
+    outs.push(boost::iostreams::gzip_compressor());
+    outs.push(output);
+    boost::archive::text_oarchive oas(outs);
+    oas << BOOST_SERIALIZATION_NVP(*m_rfc);
+  }
+
+  /*!
+    \brief Loads a configuration from the stream `input`.
+
+    The input file should be a GZIP container written by the
+    `save_configuration()` method. The feature set of the classifier
+    should contain the exact same features in the exact same order as
+    the ones present when the file was generated using
+    `save_configuration()`.
+  */
+  void load_configuration (std::istream& input)
+  {
+    CGAL::internal::liblearning::RandomForest::ForestParams params;
+    if (m_rfc != NULL)
+      delete m_rfc;
+    m_rfc = new Forest (params);
+    
+    boost::iostreams::filtering_istream ins;
+    ins.push(boost::iostreams::gzip_decompressor());
+    ins.push(input);
+    boost::archive::text_iarchive ias(ins);
+    ias >> BOOST_SERIALIZATION_NVP(*m_rfc);
+  }
+
+};
+
+}
+
+}
+
+#endif // CGAL_CLASSIFICATION_ETHZ_RANDOM_FOREST_CLASSIFIER_H

Classification/include/CGAL/Classification/Evaluation.h

@@ -107,19 +107,25 @@ public:
     m_mean_iou = 0.;
     m_mean_f1 = 0.;
 
+    std::size_t correct_labels = 0;
+    
     for (std::size_t j = 0; j < labels.size(); ++ j)
     {
       m_precision[j] = true_positives[j] / float(true_positives[j] + false_positives[j]);
       m_recall[j] = true_positives[j] / float(true_positives[j] + false_negatives[j]);
       m_iou[j] = true_positives[j] / float(true_positives[j] + false_positives[j] + false_negatives[j]);
 
+      if (std::isnan(m_iou[j]))
+        continue;
+
+      ++ correct_labels;
       m_mean_iou += m_iou[j];
       m_mean_f1 += 2.f * (m_precision[j] * m_recall[j])
         / (m_precision[j] + m_recall[j]);
     }
 
-    m_mean_iou /= labels.size();
+    m_mean_iou /= correct_labels;
-    m_mean_f1 /= labels.size();
+    m_mean_f1 /= correct_labels;
     m_accuracy = sum_true_positives / float(total);
   }
 

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++)
-        Scatter(i,j)=0;
+        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,43 +99,43 @@ 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;
-        typename Grid::iterator end = grid.indices_end(i,j);
-        for (typename Grid::iterator it = grid.indices_begin(i,j); it != end; ++ it)
         {
-          mean += float(get(point_map, *(input.begin()+(*it))).z());
+          float mean = 0.;
-          ++ nb;
+          std::size_t nb = 0;
+          typename Grid::iterator end = grid.indices_end(i,j);
+          for (typename Grid::iterator it = grid.indices_begin(i,j); it != end; ++ it)
+          {
+            mean += float(get(point_map, *(input.begin()+(*it))).z());
+            ++ nb;
+          }
+          if (nb == 0)
+            continue;
+          mean /= nb;
+          dem(i,j) = mean;
         }
-        if (nb == 0)
-          continue;
-        mean /= nb;
-        dem(i,j) = mean;
-      }
 
     std::size_t square = (std::size_t)(0.5 * radius_dtm / grid.resolution()) + 1;
     
     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);
+          std::size_t squareXmin = (i < square ? 0 : i - square);
+          std::size_t squareXmax = (std::min)(grid.width() - 1, i + square);
 
-        std::vector<float> z;
+          std::vector<float> z;
-        z.reserve(squareXmax - squareXmin +1 );
+          z.reserve(squareXmax - squareXmin +1 );
-        for(std::size_t k = squareXmin; k <= squareXmax; k++)
+          for(std::size_t k = squareXmin; k <= squareXmax; k++)
-          if (dem(k,j) != 0.)
+            if (dem(k,j) != 0.)
-            z.push_back (dem(k,j));
+              z.push_back (dem(k,j));
-        if (z.empty())
+          if (z.empty())
-          continue;
+            continue;
-        std::nth_element (z.begin(), z.begin() + (z.size() / 10), z.end());
+          std::nth_element (z.begin(), z.begin() + (z.size() / 10), z.end());
-        dtm_x(i,j) = z[z.size() / 10];
+          dtm_x(i,j) = z[z.size() / 10];
-      }
+        }
-    }
     dem.free();
 
 #ifdef CGAL_CLASSIFICATION_PRECOMPUTE_FEATURES
@@ -145,22 +145,21 @@ 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);
+          std::size_t squareYmin = (j < square ? 0 : j - square);
-        std::vector<float> z;
+          std::size_t squareYmax = (std::min)(grid.height() - 1, j + square);
-        z.reserve(squareYmax - squareYmin +1 );
+          std::vector<float> z;
-        for(std::size_t l = squareYmin; l <= squareYmax; l++)
+          z.reserve(squareYmax - squareYmin +1 );
-          if (dtm_x(i,l) != 0.)
+          for(std::size_t l = squareYmin; l <= squareYmax; l++)
-            z.push_back (dtm_x(i,l));
+            if (dtm_x(i,l) != 0.)
-        if (z.empty())
+              z.push_back (dtm_x(i,l));
-          continue;
+          if (z.empty())
-        std::nth_element (z.begin(), z.begin() + (z.size() / 10), z.end());
+            continue;
-        dtm(i,j) = z[z.size() / 10];
+          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/Simple_feature.h

@@ -0,0 +1,80 @@
+// Copyright (c) 2017 GeometryFactory Sarl (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(s)     : Simon Giraudot
+
+#ifndef CGAL_CLASSIFICATION_SIMPLE_FEATURE_H
+#define CGAL_CLASSIFICATION_SIMPLE_FEATURE_H
+
+#include <CGAL/license/Classification.h>
+
+#include <vector>
+
+namespace CGAL {
+
+namespace Classification {
+
+namespace Feature {
+
+  /*!
+    \ingroup PkgClassificationFeatures
+
+    %Feature based on a user-defined scalar field.
+
+    \tparam InputRange model of `ConstRange`. Its iterator type
+    is `RandomAccessIterator`.
+    \tparam PropertyMap model of `ReadablePropertyMap` whose key
+    type is the value type of the iterator of `PointRange` and value type
+    is statically castable to `float`.
+  */
+template <typename InputRange, typename PropertyMap>
+class Simple_feature : public Feature_base
+{
+  const InputRange& m_input;
+  PropertyMap m_pmap;
+  
+public:
+  /*!
+    \brief Constructs the feature using an input range and a property map.
+
+    \param input point range.
+    \param property_map property map to access scalar field.
+    \param name name of the property (no default name is given).
+  */
+  Simple_feature (const InputRange& input,
+                  PropertyMap property_map,
+                  const std::string& name)
+    : m_input (input), m_pmap (property_map)
+  {
+    this->set_name (name);
+  }
+  /// \cond SKIP_IN_MANUAL
+  virtual float value (std::size_t pt_index)
+  {
+    return static_cast<float> (get (m_pmap, *(m_input.begin()+pt_index)));
+  }
+  /// \endcond
+};
+
+} // namespace Feature
+
+} // namespace Classification
+
+} // namespace CGAL
+
+#endif // CGAL_CLASSIFICATION_SIMPLE_FEATURE_H

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

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

Classification/include/CGAL/Classification/Feature_set.h

@@ -48,6 +48,14 @@ class Feature_set
   typedef std::vector<Feature_handle> Base;
   Base m_features;
 
+  struct Compare_name
+  {
+    bool operator() (const Feature_handle& a, const Feature_handle& b) const
+    {
+      return a->name() < b->name();
+    }
+  };
+  
 #ifdef CGAL_LINKED_WITH_TBB
   tbb::mutex m_mutex;
   void mutex_lock() { m_mutex.lock(); }
@@ -136,6 +144,12 @@ public:
   {
     m_features[i] = Feature_handle();
   }
+
+  void sort_features_by_name()
+  {
+    std::sort (m_features.begin(), m_features.end(),
+               Compare_name());               
+  }
   /// \endcond
   
 };

Classification/include/CGAL/Classification/Image.h

@@ -23,6 +23,10 @@
 
 #include <CGAL/license/Classification.h>
 
+#include <boost/shared_ptr.hpp>
+
+#define CGAL_CLASSIFICATION_IMAGE_SIZE_LIMIT 100000000
+
 namespace CGAL {
 namespace Classification {
 
@@ -31,10 +35,19 @@ 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;
-  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&)
+  {
+  }
   
 public:
 
@@ -47,50 +60,30 @@ public:
       m_height (height)
   {
     if (m_width * m_height > 0)
-      m_raw = new Type[width * height]();
+    {
-    else
+      if (m_width * m_height < CGAL_CLASSIFICATION_IMAGE_SIZE_LIMIT)
-      m_raw = NULL;
+        m_raw = boost::shared_ptr<Vector> (new Vector(m_width * m_height));
+      else
+        m_sparse = boost::shared_ptr<Map> (new Map());
+    }
   }
   
   ~Image ()
   {
-    free();
   }
 
   void free()
   {
-    if (m_raw != NULL)
+    m_raw = boost::shared_ptr<Vector>();
-      delete[] m_raw;
+    m_sparse = boost::shared_ptr<Map>();
-    m_raw = NULL;
   }
 
-  Image (const Image& other)
-    : m_width (other.width()),
-      m_height (other.height())
-
-  {
-    if (m_width * m_height > 0)
-    {
-      m_raw = new Type[m_width * m_height];
-      std::copy (other.m_raw, other.m_raw + (m_width * m_height), this->m_raw);
-    }
-    else
-      m_raw = NULL;
-  }
   Image& operator= (const Image& other)
   {
-    if (m_raw != NULL)
+    m_raw = other.m_raw;
-      delete[] m_raw;
+    m_sparse = other.m_sparse;
-
-    m_raw = NULL;
     m_width = other.width();
     m_height = other.height();
-    if (m_width * m_height > 0)
-    {
-      m_raw = new Type[m_width * m_height];
-      std::copy (other.m_raw, other.m_raw + (m_width * m_height), this->m_raw);
-    }
-    
     return *this;
   }
   
@@ -99,13 +92,25 @@ 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
-    return m_raw[x * m_height + y];
+    {
+      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
-    return m_raw[x * m_height + y];
+    {
+      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];
   }
   
 

Classification/include/CGAL/Classification/OpenCV_random_forest_classifier.h

@@ -18,8 +18,8 @@
 //
 // Author(s)     : Simon Giraudot
 
-#ifndef CGAL_CLASSIFICATION_RANDOM_FOREST_CLASSIFIER_H
+#ifndef CGAL_CLASSIFICATION_OPENCV_RANDOM_FOREST_CLASSIFIER_H
-#define CGAL_CLASSIFICATION_RANDOM_FOREST_CLASSIFIER_H
+#define CGAL_CLASSIFICATION_OPENCV_RANDOM_FOREST_CLASSIFIER_H
 
 #include <CGAL/license/Classification.h>
 
@@ -36,13 +36,13 @@ namespace Classification {
 /*!
   \ingroup PkgClassificationClassifiers
 
-  \brief %Classifier based on a random forest algorithm.
+  \brief %Classifier based on the OpenCV version of random forest algorithm.
 
   \note This class requires the \ref thirdpartyOpenCV library.
 
   \cgalModels `CGAL::Classification::Classifier`
 */
-class Random_forest_classifier
+class OpenCV_random_forest_classifier
 {
   const Label_set& m_labels;
   const Feature_set& m_features;
@@ -60,6 +60,9 @@ class Random_forest_classifier
 
 public:
   
+  /// \name Constructor
+  /// @{
+  
 /*!
   \brief Instantiate the classifier using the sets of `labels` and `features`.
 
@@ -73,13 +76,13 @@ public:
   \param max_number_of_trees_in_the_forest The maximum number of trees in the forest (surprise, surprise). Typically the more trees you have the better the accuracy. However, the improvement in accuracy generally diminishes and asymptotes pass a certain number of trees. Also to keep in mind, the number of tree increases the prediction time linearly.
   \param forest_accuracy Sufficient accuracy (OOB error).
 */
-  Random_forest_classifier (const Label_set& labels,
+  OpenCV_random_forest_classifier (const Label_set& labels,
-                            const Feature_set& features,
+                                   const Feature_set& features,
-                            int max_depth = 20,
+                                   int max_depth = 20,
-                            int min_sample_count = 5,
+                                   int min_sample_count = 5,
-                            int max_categories = 15,
+                                   int max_categories = 15,
-                            int max_number_of_trees_in_the_forest = 100,
+                                   int max_number_of_trees_in_the_forest = 100,
-                            float forest_accuracy = 0.01f)
+                                   float forest_accuracy = 0.01f)
     : m_labels (labels), m_features (features),
       m_max_depth (max_depth), m_min_sample_count (min_sample_count),
       m_max_categories (max_categories),
@@ -91,7 +94,7 @@ public:
   {  }
 
   /// \cond SKIP_IN_MANUAL
-  ~Random_forest_classifier ()
+  ~OpenCV_random_forest_classifier ()
   {
 #if (CV_MAJOR_VERSION < 3)
     if (rtree != NULL)
@@ -100,6 +103,24 @@ public:
   }
   /// \endcond
   
+  /// @}
+
+  /// \name Parameters
+  /// @{
+  
+  void set_max_depth (int max_depth) { m_max_depth = max_depth; }
+  void set_min_sample_count (int min_sample_count) { m_min_sample_count = min_sample_count; }
+  void set_max_categories (int max_categories) { m_max_categories = max_categories; }
+  void set_max_number_of_trees_in_the_forest (int max_number_of_trees_in_the_forest)
+  { m_max_number_of_trees_in_the_forest = max_number_of_trees_in_the_forest; }
+  void set_forest_accuracy (float forest_accuracy) { m_forest_accuracy = forest_accuracy; }
+  
+  
+  /// @}
+
+  /// \name Training
+  /// @{
+  
   /*!
     \brief Runs the training algorithm.
 
@@ -191,12 +212,7 @@ public:
 
   }
 
-  void set_max_depth (int max_depth) { m_max_depth = max_depth; }
+  /// @}
-  void set_min_sample_count (int min_sample_count) { m_min_sample_count = min_sample_count; }
-  void set_max_categories (int max_categories) { m_max_categories = max_categories; }
-  void set_max_number_of_trees_in_the_forest (int max_number_of_trees_in_the_forest)
-  { m_max_number_of_trees_in_the_forest = max_number_of_trees_in_the_forest; }
-  void set_forest_accuracy (float forest_accuracy) { m_forest_accuracy = forest_accuracy; }
 
   /// \cond SKIP_IN_MANUAL
   void operator() (std::size_t item_index, std::vector<float>& out) const
@@ -232,12 +248,35 @@ public:
 
 #endif
   }
+  /// \endcond
 
+  /// \name Input/Output
+  /// @{
+
+
+  /*!
+    \brief Saves the current configuration in the file named `filename`.
+
+    This allows to easily save and recover a specific classification
+    configuration.
+
+    The output file is written in an XML format that is readable by
+    the `load_configuration()` method.
+  */
   void save_configuration (const char* filename)
   {
     rtree->save(filename);
   }
 
+  /*!
+    \brief Loads a configuration from the file named `filename`.
+
+    The input file should be in the XML format written by the
+    `save_configuration()` method. The feature set of the classifier
+    should contain the exact same features in the exact same order as
+    the ones present when the file was generated using
+    `save_configuration()`.
+  */
   void load_configuration (const char* filename)
   {
 #if (CV_MAJOR_VERSION < 3)
@@ -249,7 +288,7 @@ public:
     rtree = cv::ml::StatModel::load<cv::ml::RTrees> (filename);
 #endif
   }
-  /// \endcond
+
 
 };
 
@@ -257,4 +296,4 @@ public:
 
 }
 
-#endif // CGAL_CLASSIFICATION_RANDOM_FOREST_CLASSIFIER_H
+#endif // CGAL_CLASSIFICATION_OPENCV_RANDOM_FOREST_CLASSIFIER_H

Classification/include/CGAL/Classification/Point_set_feature_generator.h

@@ -343,6 +343,8 @@ public:
                             get_parameter<Vmap>(normal_map),
                             get_parameter<Cmap>(color_map),
                             get_parameter<Emap>(echo_map));
+
+    m_features->sort_features_by_name();
   }
 
   

Classification/include/CGAL/Classification/Point_set_neighborhood.h

@@ -25,7 +25,7 @@
 
 #include <vector>
 
-#include <boost/iterator/counting_iterator.hpp>
+#include <CGAL/boost/iterator/counting_iterator.hpp>
 
 #include <CGAL/Search_traits_3.h>
 #include <CGAL/Fuzzy_sphere.h>

Classification/include/CGAL/Classification/Sum_of_weighted_features_classifier.h

@@ -42,6 +42,7 @@
 #include <tbb/mutex.h>
 #endif // CGAL_LINKED_WITH_TBB
 
+#define CLASSIFICATION_TRAINING_QUICK_ESTIMATION
 
 namespace CGAL {
 
@@ -311,6 +312,11 @@ public:
         ++ nb_tot;
       }
 
+#ifdef CLASSIFICATION_TRAINING_QUICK_ESTIMATION
+    for (std::size_t i = 0; i < m_labels.size(); ++ i)
+      std::random_shuffle (training_sets[i].begin(), training_sets[i].end());
+#endif
+    
     CGAL_CLASSIFICATION_CERR << "Training using " << nb_tot << " inliers" << std::endl;
     
     for (std::size_t i = 0; i < m_labels.size(); ++ i)
@@ -841,12 +847,19 @@ private:
                                   
     for (std::size_t j = 0; j < m_labels.size(); ++ j)
     {
-      for (std::size_t k = 0; k < training_sets[j].size(); ++ k)
+#ifdef CLASSIFICATION_TRAINING_QUICK_ESTIMATION
+      std::size_t training_set_size = (std::min) (std::size_t(0.1 * training_sets[j].size()),
+                                                  std::size_t(10000));
+#else
+      std::size_t training_set_size = training_sets[j].size();
+#endif
+      
+      for (std::size_t k = 0; k < training_set_size; ++ k)
       {
         float val = normalized(feature, training_sets[j][k]);
         mean[j] += val;
       }
-      mean[j] /= training_sets[j].size();
+      mean[j] /= training_set_size;
     }
 
     std::vector<float> sd (m_labels.size(), 0.);
@@ -855,12 +868,18 @@ private:
     {
       Label_handle clabel = m_labels[j];
             
-      for (std::size_t k = 0; k < training_sets[j].size(); ++ k)
+#ifdef CLASSIFICATION_TRAINING_QUICK_ESTIMATION
+      std::size_t training_set_size = (std::min) (std::size_t(0.1 * training_sets[j].size()),
+                                                  std::size_t(10000));
+#else
+      std::size_t training_set_size = training_sets[j].size();
+#endif
+      for (std::size_t k = 0; k < training_set_size; ++ k)
       {
         float val = normalized(feature, training_sets[j][k]);
         sd[j] += (val - mean[j]) * (val - mean[j]);
       }
-      sd[j] = std::sqrt (sd[j] / training_sets[j].size());
+      sd[j] = std::sqrt (sd[j] / training_set_size);
       if (mean[j] - sd[j] > (2./3.))
         set_effect (j, feature, FAVORING);
       else if (mean[j] + sd[j] < (1./3.))

Classification/include/CGAL/Classification/internal/auxiliary/dataview.h

@@ -0,0 +1,178 @@
+// Copyright (c) 2014 Stefan Walk
+//
+// Permission is hereby granted, free of charge, to any person obtaining a copy of
+// this software and associated documentation files (the "Software"), to deal in
+// the Software without restriction, including without limitation the rights to
+// use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies
+// of the Software, and to permit persons to whom the Software is furnished to do
+// so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included in all
+// copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+// SOFTWARE.
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: MIT
+//
+// Author(s)     : Stefan Walk
+
+// Modifications from original library:
+//  * changed inclusion protection tag
+//  * moved to namespace CGAL::internal::
+
+
+#ifndef CGAL_INTERNAL_LIBLEARNING_DATAVIEW_H
+#define CGAL_INTERNAL_LIBLEARNING_DATAVIEW_H
+#include <vector>
+
+namespace CGAL { namespace internal {
+    
+namespace liblearning {
+//! \brief A view for one-dimensional data
+//
+// The data can be e.g. a row or column of a matrix or the elements of an
+// array or a std::vector - anything with a fixed step size between pointers
+template <typename ElementType>
+struct DataView1D {
+    //! \brief Element access
+    //
+    // \param idx the index of the element
+    ElementType& operator()(size_t idx) const
+    {
+        return *(data + step * idx);
+    }
+    //! \brief Construct using pointer, size and (optional) step size
+    DataView1D(ElementType* ptr, size_t size, ptrdiff_t step_size = 1) :
+        data(ptr), step(step_size), num_elements(size)
+    {
+    }
+    //! \brief Construct view from std::vector
+    //
+    // \param vec the vector, by-ref since we don't take ownership - no
+    //            temporaries allowed. maybe needs to change
+    DataView1D(std::vector<ElementType>& vec) :
+        data(&vec[0]), step(1), num_elements(vec.size())
+    {
+    }
+    
+    //! \brief Construct empty view
+    DataView1D() : data(0), step(1), num_elements(0)
+    {
+    }
+    ElementType* data; //!< Pointer to first element
+    ptrdiff_t step;    //!< Step size between elements - for a vector this is 1, no need to multiply by sizeof(ElementType)
+    size_t num_elements; //!< Number of elements in the view
+};
+
+template <typename ElementType>
+struct DataView2D {
+    //! \brief Construct empty view
+    DataView2D() : data(0), row_step(1), col_step(1), rows(0), cols(0)
+    {}
+    //! \brief Construct view from memory using given step sizes
+    DataView2D(ElementType* ptr, size_t rows_, size_t cols_, ptrdiff_t row_step_, ptrdiff_t col_step_) :
+        data(ptr), row_step(row_step_), col_step(col_step_), rows(rows_), cols(cols_)
+    {
+    }
+    //! \brief Construct view from a continuous block of memory in row-major order
+    DataView2D(ElementType* ptr, size_t rows_, size_t cols_) :
+        data(ptr), row_step(cols_), col_step(1), rows(rows_), cols(cols_)
+    {
+    }
+    //! \brief Element access
+    ElementType& operator()(size_t row_idx, size_t col_idx) const
+    {
+        return *(data + row_step * row_idx + col_step * col_idx);
+    }
+    //! \brief Return a 1D view of a row
+    DataView1D<ElementType> row(size_t row_idx)
+    {
+        return DataView1D<ElementType>(data + row_step * row_idx, cols, col_step);
+    }
+    //! \brief Return a 1D view of a column
+    DataView1D<ElementType> col(size_t col_idx)
+    {
+        return DataView1D<ElementType>(data + col_step * col_idx, rows, row_step);
+    }
+    //! \brief Return a new view, using a subset of rows
+    DataView2D row_range(size_t start_row, size_t end_row) const
+    {
+        DataView2D ret(*this);
+        ret.data = data + row_step * start_row;
+        ret.rows = end_row - start_row;
+        return ret;
+    }
+    //! \brief Return a new view, using a subset of columns
+    DataView2D col_range(size_t start_col, size_t end_col) const
+    {
+        DataView2D ret(*this);
+        ret.data = data + col_step * start_col;
+        ret.cols = end_col - start_col;
+        return ret;
+    }
+    //! \brief Transpose the matrix (actually done by swapping the steps, no
+    //                               copying)
+    DataView2D transpose() const
+    {
+        DataView2D ret(*this);
+        std::swap(ret.row_step, ret.col_step);
+        std::swap(ret.rows,     ret.cols);
+        return ret;
+    }
+    //! \brief Create a 2D view from a 1D one, represent as a column
+    static DataView2D column(DataView1D<ElementType> vec)
+    {
+        return DataView2D(vec.data, vec.num_elements, 1);
+    }
+    //! \brief Create a 2D view from a 1D one, represent as a row
+    static DataView2D row(DataView1D<ElementType> vec)
+    {
+        return DataView2D(vec.data, 1, vec.num_elements);
+    }
+    //! \brief Return the number of elements in this view
+    size_t num_elements() { return rows * cols; }
+    //! \brief Return true if the view is empty (no elements)
+    bool   empty() { return num_elements() == 0; }
+    //! \brief Return true if rows in this view are continuous in memory
+    bool row_continuous() {
+        return col_step == 1;
+    }
+    //! \brief Return true if columns in this view are continuous in memory
+    bool col_continuous() {
+        return row_step == 1;
+    }
+    //! \brief Get pointer to row (only valid if row-continuous)
+    ElementType* row_pointer(size_t row_idx) {
+        return data + row_idx * row_step;
+    }
+    //! \brief Get pointer to column (only valid if column-continuous)
+    ElementType* col_pointer(size_t col_idx) {
+        return data + col_idx * col_step;
+    }
+
+    ElementType* data; //!< Pointer to first element
+    ptrdiff_t row_step; //!< Pointer difference between an element and its right neighbor - 1 if row-continuous
+    ptrdiff_t col_step; //!< Pointer difference between an element and its bottom neighbor - 1 if column-continuous
+    size_t rows; //!< Number of rows in the view
+    size_t cols; //!< Number of columns in the view
+};
+
+//! \brief Determine if two 2D views have the same dimensions
+template <typename A, typename B>
+bool equal_dims(DataView2D<A> view_a, DataView2D<B> view_b)
+{
+    return view_a.rows == view_b.rows && view_a.cols == view_b.cols;
+}
+}
+
+}} // namespace CGAL::internal::
+
+#endif

Classification/include/CGAL/Classification/internal/auxiliary/random-forest/common-libraries.hpp

@@ -0,0 +1,340 @@
+// Copyright (c) 2014 Stefan Walk
+//
+// Permission is hereby granted, free of charge, to any person obtaining a copy of
+// this software and associated documentation files (the "Software"), to deal in
+// the Software without restriction, including without limitation the rights to
+// use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies
+// of the Software, and to permit persons to whom the Software is furnished to do
+// so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included in all
+// copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+// SOFTWARE.
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: MIT
+//
+// Author(s)     : Stefan Walk
+
+// Modifications from original library:
+//  * changed inclusion protection tag
+//  * moved to namespace CGAL::internal::
+
+#ifndef CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_COMMON_LIBRARIES_H
+#define CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_COMMON_LIBRARIES_H
+#include <algorithm>
+#include <numeric>
+#include <limits>
+#include <list>
+#include <boost/version.hpp>
+#include <boost/bind.hpp>
+#include <boost/foreach.hpp>
+#include <boost/ptr_container/ptr_vector.hpp>
+#include <boost/random/mersenne_twister.hpp>
+#if BOOST_VERSION >= 104700
+#  include <boost/random/uniform_int_distribution.hpp>
+#else 
+#  include <boost/random/uniform_int.hpp>
+#endif
+#include <boost/random/uniform_01.hpp>
+#include <boost/random/normal_distribution.hpp>
+#include <boost/serialization/vector.hpp>
+#include <boost/scoped_ptr.hpp>
+#include <boost/shared_ptr.hpp>
+#include <boost/make_shared.hpp>
+#include <boost/unordered_set.hpp>
+#include <iostream>
+#include <cstdio>
+
+#include "../dataview.h"
+
+namespace CGAL { namespace internal {
+
+namespace liblearning {
+namespace RandomForest {
+
+typedef std::vector< std::pair<float, int> > FeatureClassDataFloat;
+inline void init_feature_class_data(FeatureClassDataFloat& data, int /*n_classes*/, int n_samples)
+{
+    data.resize(n_samples);
+}
+typedef boost::unordered_set<int> FeatureSet;
+
+#if BOOST_VERSION >= 104700
+typedef boost::random::uniform_int_distribution<> UniformIntDist;
+typedef boost::random::normal_distribution<> NormalDist;
+typedef boost::random::mt19937 RandomGen;
+typedef boost::random::uniform_01<> UnitDist;
+#else 
+typedef boost::uniform_int<> UniformIntDist;
+typedef boost::normal_distribution<> NormalDist;
+typedef boost::uniform_01<> UnitDist;
+typedef boost::mt19937 RandomGen;
+#endif
+
+struct ForestParams { 
+    size_t n_classes;
+    size_t n_features;
+    size_t n_samples;
+    size_t n_in_bag_samples;
+    size_t max_depth;
+    size_t n_trees;
+    size_t min_samples_per_node;
+    float  sample_reduction;
+    ForestParams() :
+        n_classes(0),
+        n_features(0),
+        n_samples(0),
+        n_in_bag_samples(0),
+        max_depth(42),
+        n_trees(100),
+        min_samples_per_node(5),
+        sample_reduction(0)
+    {}
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /*version*/)
+    {
+        ar & BOOST_SERIALIZATION_NVP(n_classes);
+        ar & BOOST_SERIALIZATION_NVP(n_features);
+        ar & BOOST_SERIALIZATION_NVP(n_samples);
+        ar & BOOST_SERIALIZATION_NVP(n_in_bag_samples);
+        ar & BOOST_SERIALIZATION_NVP(max_depth);
+        ar & BOOST_SERIALIZATION_NVP(n_trees);
+        ar & BOOST_SERIALIZATION_NVP(min_samples_per_node);
+        ar & BOOST_SERIALIZATION_NVP(sample_reduction);
+    }
+};
+
+struct QuadraticSplitter {
+    typedef float FeatureType;
+    typedef FeatureClassDataFloat FeatureClassData;
+    int n_features;
+    std::vector<FeatureType> w;
+    FeatureType threshold;
+    QuadraticSplitter() : n_features(0) {}
+    QuadraticSplitter(int n_features, std::vector<FeatureType> const& w) :
+        n_features(n_features), w(w)
+    {}
+    void set_threshold(FeatureType new_threshold) {
+        threshold = new_threshold;
+    }
+    FeatureType map_sample(FeatureType const* v) const {
+        double result = 0.0;
+        int i_feature = 0;
+        for (; i_feature < n_features; ++i_feature) {
+            result += w[i_feature] * v[i_feature];
+        }
+        for (int i_1 = 0; i_1 < n_features; ++i_1) {
+            for (int i_2 = 0; i_2 < n_features; ++i_2) {
+                result += w[i_feature++] * v[i_1] * v[i_2];
+            }
+        }
+        return result;
+    }
+    bool classify_sample(FeatureType const* v) const {
+        return map_sample(v) > threshold;
+    }
+    void map_points(DataView2D<FeatureType> samples,
+                    DataView2D<int>   labels,
+                    int const*      sample_idxes,
+                    int             n_samples,
+                    FeatureClassData&        data_points) const
+    {
+        for (int i_sample = 0; i_sample < n_samples; ++i_sample) {
+            int sample_idx    = sample_idxes[i_sample];
+            int sample_class  = labels(sample_idx, 0);
+            FeatureType sample_fval = map_sample(samples.row_pointer(sample_idx));
+            data_points[i_sample] = std::make_pair(sample_fval, sample_class);
+        }
+    }
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /*version*/)
+    {
+        ar & BOOST_SERIALIZATION_NVP(n_features);
+        ar & BOOST_SERIALIZATION_NVP(w);
+        ar & BOOST_SERIALIZATION_NVP(threshold);
+    }
+};
+
+struct LinearSplitter {
+    typedef float FeatureType;
+    typedef FeatureClassDataFloat FeatureClassData;
+    std::vector<FeatureType> w;
+    FeatureType threshold;
+    LinearSplitter() {}
+    LinearSplitter(std::vector<FeatureType> const& w) :
+        w(w)
+    {}
+    void set_threshold(FeatureType new_threshold) {
+        threshold = new_threshold;
+    }
+    bool classify_sample(FeatureType const* v) const {
+        return std::inner_product(w.begin(), w.end(), v, 0.0f) > threshold;
+    }
+    void map_points(DataView2D<FeatureType> samples,
+                    DataView2D<int>   labels,
+                    int const*      sample_idxes,
+                    int             n_samples,
+                    FeatureClassData&        data_points) const
+    {
+        for (int i_sample = 0; i_sample < n_samples; ++i_sample) {
+            int sample_idx    = sample_idxes[i_sample];
+            int sample_class  = labels(sample_idx, 0);
+            FeatureType sample_fval = std::inner_product(w.begin(), w.end(), 
+                                                   samples.row_pointer(sample_idx), 0.0f);
+            data_points[i_sample] = std::make_pair(sample_fval, sample_class);
+        }
+    }
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /*version*/)
+    {
+        ar & BOOST_SERIALIZATION_NVP(w);
+        ar & BOOST_SERIALIZATION_NVP(threshold);
+    }
+};
+
+struct AxisAlignedSplitter {
+    typedef float FeatureType;
+    typedef FeatureClassDataFloat FeatureClassData;
+    int feature;
+    FeatureType threshold;
+    AxisAlignedSplitter() : feature(-1) {}
+    AxisAlignedSplitter(int feature) : 
+        feature(feature) 
+    {}
+    void set_threshold(FeatureType new_threshold) {
+        threshold = new_threshold;
+    }
+    bool classify_sample(FeatureType const* v) const {
+        return v[feature] > threshold;
+    }
+    void map_points(DataView2D<FeatureType> samples,
+                    DataView2D<int>   labels,
+                    int const*      sample_idxes,
+                    int             n_samples,
+                    FeatureClassData&        data_points) const
+    {
+        for (int i_sample = 0; i_sample < n_samples; ++i_sample) {
+            // determine index of this sample ...
+            int sample_idx    = sample_idxes[i_sample];
+            // determine class ...
+            int sample_class  = labels(sample_idx, 0);
+            // determine value of the selected feature for this sample
+            FeatureType sample_fval = samples(sample_idx, feature);
+            data_points[i_sample] = std::make_pair(sample_fval, sample_class);
+        }
+    }
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /*version*/)
+    {
+        ar & BOOST_SERIALIZATION_NVP(feature);
+        ar & BOOST_SERIALIZATION_NVP(threshold);
+    }
+};
+
+struct AxisAlignedRandomSplitGenerator {
+    typedef float FeatureType;
+    FeatureSet features;
+    FeatureSet::const_iterator it;
+
+    void init(DataView2D<FeatureType> samples,
+              DataView2D<int>   /*labels*/,
+              int*            /*sample_idxes*/,
+              int             /*n_samples*/,
+              size_t          /*n_classes*/,
+              RandomGen&      gen)
+    {
+        features.clear();
+        int n_features = samples.cols;
+        size_t n_used_features = std::sqrt(n_features);
+        UniformIntDist dist(0, n_features - 1);
+        // insert into set until required number of unique features is found
+        while (features.size() < n_used_features) {
+            features.insert(dist(gen));
+        }
+        it = features.begin();
+    }
+    AxisAlignedSplitter gen_proposal(RandomGen& /*gen*/)
+    {
+        if (it == features.end()) {
+            it = features.begin();
+        }
+        return AxisAlignedSplitter(*it++);
+    }
+    size_t num_proposals() const {
+        return features.size();
+    }
+};
+
+struct LinearSplitGenerator {
+    typedef float FeatureType;
+    size_t n_features;
+    size_t n_proposals;
+    LinearSplitGenerator(size_t n_proposals = 5) : 
+        n_proposals(n_proposals)
+    {}
+    void init(DataView2D<FeatureType> samples,
+              DataView2D<int>   /*labels*/,
+              int*            /*sample_idxes*/,
+              int             /*n_samples*/,
+              size_t          /*n_classes*/,
+              RandomGen&      /*gen*/)
+    {
+        n_features = samples.cols;
+    }
+    size_t num_proposals() const {
+        return n_proposals;
+    }
+    LinearSplitter gen_proposal(RandomGen& gen) {
+        NormalDist dist;
+        std::vector<FeatureType> weights(n_features);
+        for (size_t i_feature = 0; i_feature < n_features; ++i_feature) {
+            weights[i_feature] = dist(gen);
+        }
+        return LinearSplitter(weights);
+    }
+};
+
+struct QuadraticSplitGenerator {
+    typedef float FeatureType;
+    size_t n_features;
+    size_t n_proposals;
+    QuadraticSplitGenerator(size_t n_proposals = 5) : 
+        n_proposals(n_proposals)
+    {}
+    void init(DataView2D<FeatureType> samples,
+              DataView2D<int>   /*labels*/,
+              int*            /*sample_idxes*/,
+              int             /*n_samples*/,
+              size_t          /*n_classes*/,
+              RandomGen&      /*gen*/)
+    {
+        n_features = samples.cols;
+    }
+    size_t num_proposals() const {
+        return n_proposals;
+    }
+    QuadraticSplitter gen_proposal(RandomGen& gen) {
+        NormalDist dist;
+        std::vector<FeatureType> weights(n_features + n_features*n_features);
+        for (size_t i_feature = 0; i_feature < weights.size(); ++i_feature) {
+            weights[i_feature] = dist(gen);
+        }
+        return QuadraticSplitter(n_features, weights);
+    }
+};
+
+}
+}
+
+}} // namespace CGAL::internal::
+    
+#endif

Classification/include/CGAL/Classification/internal/auxiliary/random-forest/forest.hpp

@@ -0,0 +1,187 @@
+// Copyright (c) 2014 Stefan Walk
+//
+// Permission is hereby granted, free of charge, to any person obtaining a copy of
+// this software and associated documentation files (the "Software"), to deal in
+// the Software without restriction, including without limitation the rights to
+// use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies
+// of the Software, and to permit persons to whom the Software is furnished to do
+// so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included in all
+// copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+// SOFTWARE.
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: MIT
+//
+// Author(s)     : Stefan Walk
+
+// Modifications from original library:
+//  * changed inclusion protection tag
+//  * moved to namespace CGAL::internal::
+//  * add parameter "reset_trees" to train() to be able to construct
+//    forest with several iterations
+
+#ifndef CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_FOREST_H
+#define CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_FOREST_H
+#include "common-libraries.hpp"
+#include "tree.hpp"
+#include <boost/ptr_container/serialize_ptr_vector.hpp>
+#if VERBOSE_TREE_PROGRESS
+#include <cstdio>
+#endif
+
+namespace CGAL { namespace internal {
+
+namespace liblearning {
+namespace RandomForest {
+
+template <typename NodeT>
+class RandomForest {
+public:
+    typedef typename NodeT::ParamType ParamType;
+    typedef typename NodeT::FeatureType FeatureType;
+    typedef Tree<NodeT> TreeType;
+    ParamType params;
+
+    std::vector<uint8_t> was_oob_data;
+    DataView2D<uint8_t> was_oob;
+
+    boost::ptr_vector< Tree<NodeT> > trees;
+
+    RandomForest() {}
+    RandomForest(ParamType const& params) : params(params) {}
+
+    template<typename SplitGenerator>
+    void train(DataView2D<FeatureType> samples, 
+               DataView2D<int> labels, 
+               DataView2D<int> train_sample_idxes, 
+               SplitGenerator const& split_generator,
+               size_t seed_start = 1,
+               bool register_oob = true,
+               bool reset_trees = true
+               ) 
+    {
+        if (reset_trees)
+          trees.clear();
+
+        params.n_classes  = *std::max_element(&labels(0,0), &labels(0,0)+labels.num_elements()) + 1;
+        params.n_features = samples.cols;
+        params.n_samples  = samples.rows;
+
+        std::vector<int> sample_idxes;
+
+        if (train_sample_idxes.empty()) {
+            // no indexes were passed, generate vector with all indexes
+            sample_idxes.resize(params.n_samples);
+            for (size_t i_sample = 0; i_sample < params.n_samples; ++i_sample) {
+                sample_idxes[i_sample] = i_sample;
+            }
+        } else {
+            // copy indexes
+            sample_idxes.assign(&train_sample_idxes(0,0), &train_sample_idxes(0,0)+train_sample_idxes.num_elements());
+        }
+
+        size_t n_idxes = sample_idxes.size();
+        params.n_in_bag_samples = n_idxes * (1 - params.sample_reduction);
+
+        // Random distribution over indexes
+        UniformIntDist dist(0, n_idxes - 1);
+
+        // Store for each sample and each tree if sample was used for tree
+        if (register_oob) {
+            was_oob_data.assign(n_idxes*params.n_trees, 1);
+            was_oob = DataView2D<uint8_t>(&was_oob_data[0], n_idxes, params.n_trees);
+        }
+
+        std::size_t nb_trees = trees.size();
+        for (size_t i_tree = nb_trees; i_tree < nb_trees + params.n_trees; ++i_tree) {
+#if VERBOSE_TREE_PROGRESS
+            std::printf("Training tree %zu/%zu, max depth %zu\n", i_tree+1, nb_trees + params.n_trees, params.max_depth);
+#endif
+            // new tree
+            trees.push_back(new TreeType(&params));
+            // initialize random generator with sequential seeds (one for each
+            // tree)
+            RandomGen gen(seed_start + i_tree);
+            // Bagging: draw random sample indexes used for this tree
+            std::vector<int> in_bag_samples(params.n_in_bag_samples);
+            for (size_t i_sample = 0; i_sample < in_bag_samples.size(); ++i_sample) {
+                int random_idx = dist(gen);
+                in_bag_samples[i_sample] = sample_idxes[random_idx];
+                if (register_oob && was_oob(random_idx, i_tree)) {
+                    was_oob(random_idx, i_tree) = 0;
+                }
+            }
+#ifdef TREE_GRAPHVIZ_STREAM
+            TREE_GRAPHVIZ_STREAM << "digraph Tree {" << std::endl;
+#endif
+            // Train the tree
+            trees.back().train(samples, labels, &in_bag_samples[0], in_bag_samples.size(), split_generator, gen);
+#ifdef TREE_GRAPHVIZ_STREAM
+            TREE_GRAPHVIZ_STREAM << "}" << std::endl << std::endl;
+#endif
+        }
+    }
+    int evaluate(FeatureType const* sample, float* results) {
+        // initialize output probabilities to 0
+        std::fill_n(results, params.n_classes, 0.0f);
+        // accumulate votes of the trees
+        for (size_t i_tree = 0; i_tree < trees.size(); ++i_tree) {
+            float const* tree_result = trees[i_tree].evaluate(sample);
+            for (size_t i_cls = 0; i_cls < params.n_classes; ++i_cls) {
+                results[i_cls] += tree_result[i_cls];
+            }
+        }
+        float best_val   = 0.0;
+        int   best_class = 0;
+        float scale      = 1.0 / trees.size();
+        for (size_t i_cls = 0; i_cls < params.n_classes; ++i_cls) {
+            // divide by number of trees to normalize probability
+            results[i_cls] *= scale;
+            // determine best class
+            if (results[i_cls] > best_val) {
+                best_val = results[i_cls];
+                best_class = i_cls;
+            }
+        }
+        return best_class;
+    }
+#if 0
+    float similarity_endnode(float const* sample_1, float const* sample_2) {
+        double sum = 0.0;
+        for (size_t i_tree = 0; i_tree < trees.size(); ++i_tree) {
+            sum += trees[i_tree].similarity_endnode(sample_1, sample_2);
+        }
+        return sum/trees.size();
+    }
+    float similarity_path(float const* sample_1, float const* sample_2) {
+        double sum = 0.0;
+        for (size_t i_tree = 0; i_tree < trees.size(); ++i_tree) {
+            sum += trees[i_tree].similarity_path(sample_1, sample_2);
+        }
+        return sum/trees.size();
+    }
+#endif
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /* version */)
+    {
+        ar & BOOST_SERIALIZATION_NVP(params);
+        ar & BOOST_SERIALIZATION_NVP(trees);
+    }
+};
+
+}
+}
+
+}} // namespace CGAL::internal::
+
+#endif

Classification/include/CGAL/Classification/internal/auxiliary/random-forest/node-gini.hpp

@@ -0,0 +1,116 @@
+// Copyright (c) 2014 Stefan Walk
+//
+// Permission is hereby granted, free of charge, to any person obtaining a copy of
+// this software and associated documentation files (the "Software"), to deal in
+// the Software without restriction, including without limitation the rights to
+// use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies
+// of the Software, and to permit persons to whom the Software is furnished to do
+// so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included in all
+// copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+// SOFTWARE.
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: MIT
+//
+// Author(s)     : Stefan Walk
+
+// Modifications from original library:
+//  * changed inclusion protection tag
+//  * moved to namespace CGAL::internal::
+
+#ifndef CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_NODE_GINI_H
+#define CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_NODE_GINI_H
+#include "node.hpp"
+#include "common-libraries.hpp"
+
+namespace CGAL { namespace internal {
+
+namespace liblearning {
+namespace RandomForest {
+
+/*
+template <typename T>
+class X : Y<X> {}
+-> http://en.wikipedia.org/wiki/Curiously_recurring_template_pattern#Static_polymorphism
+*/
+
+template <typename Splitter>
+class NodeGini : public Node< NodeGini<Splitter>, ForestParams, Splitter > {
+public:
+    typedef typename Node< NodeGini<Splitter>, ForestParams, Splitter>::ParamType ParamType;
+    typedef typename Splitter::FeatureType FeatureType;
+    typedef typename Splitter::FeatureClassData FeatureClassData;
+    using Node< NodeGini<Splitter>, ForestParams, Splitter>::params;
+    NodeGini() {}
+    NodeGini(int depth, ParamType const* params) :
+        Node< NodeGini<Splitter>, ForestParams, Splitter>(depth, params)
+    { 
+    }
+
+    uint64_t gini_square_term(std::vector<uint64_t> const& frequencies) const 
+    {
+        return std::inner_product( frequencies.begin(), frequencies.end(), frequencies.begin(), uint64_t(0));
+    }
+    std::pair<FeatureType, float> determine_best_threshold(FeatureClassData& data_points,
+                                                     std::vector<uint64_t>& classes_l,
+                                                     std::vector<uint64_t>& classes_r,
+                                                     RandomGen&        gen)
+    {
+        double best_loss = std::numeric_limits<double>::infinity();
+        float best_thresh = 0;
+
+        UnitDist fraction_dist;
+        classes_l.assign(params->n_classes, 0);
+        classes_r.assign(params->n_classes, 0);
+        double n_l = 0;
+        double n_r = 0;
+        for (size_t i_sample = 0; i_sample < data_points.size(); ++i_sample) {
+            classes_r[data_points[i_sample].second]++;
+            n_r += 1;
+        }
+        // sort data so thresholding is easy based on position in array
+        std::sort(data_points.begin(), data_points.end());
+        // loop over data, update class distributions left&right
+        for (size_t i_point = 1; i_point < data_points.size(); ++i_point) {
+            int cls = data_points[i_point-1].second;
+            classes_l[cls]++; // sample with class cls gets moved to left ...
+            classes_r[cls]--; // from right
+            n_l += 1;
+            n_r -= 1;
+            // don't split here if values are the same
+            if (data_points[i_point-1].first == data_points[i_point].first)
+                continue;
+            // weighted average
+            double gini = n_l - gini_square_term(classes_l) * 1.0 / n_l + n_r - gini_square_term(classes_r) * 1.0 / n_r;
+            if (gini < best_loss) {
+                best_loss = gini;
+                double fraction = fraction_dist(gen);
+                best_thresh = fraction * data_points[i_point-1].first + (1-fraction) * data_points[i_point].first;
+            }
+        }
+        return std::make_pair(best_thresh, float(best_loss));
+    }
+
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /* version */)
+    {
+        ar & boost::serialization::make_nvp("base",  boost::serialization::base_object< Node< NodeGini<Splitter>, ForestParams, Splitter > >(*this));
+    }
+};
+
+}
+}
+
+}} // namespace CGAL::internal::
+
+#endif

Classification/include/CGAL/Classification/internal/auxiliary/random-forest/node.hpp

@@ -0,0 +1,262 @@
+// Copyright (c) 2014 Stefan Walk
+//
+// Permission is hereby granted, free of charge, to any person obtaining a copy of
+// this software and associated documentation files (the "Software"), to deal in
+// the Software without restriction, including without limitation the rights to
+// use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies
+// of the Software, and to permit persons to whom the Software is furnished to do
+// so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included in all
+// copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+// SOFTWARE.
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: MIT
+//
+// Author(s)     : Stefan Walk
+
+// Modifications from original library:
+//  * changed inclusion protection tag
+//  * moved to namespace CGAL::internal::
+//  * fix computation of node_dist[label] so that results are always <= 1.0
+
+#ifndef CGAL_INTERNAL_LIBLEARNING_RANDOMFORESTS_NODE_H
+#define CGAL_INTERNAL_LIBLEARNING_RANDOMFORESTS_NODE_H
+#include "../dataview.h"
+#include "common-libraries.hpp"
+#include <boost/serialization/scoped_ptr.hpp>
+#include <boost/serialization/vector.hpp>
+#if VERBOSE_NODE_LEARNING
+#include <cstdio>
+#endif
+
+namespace CGAL { namespace internal {
+
+namespace liblearning {
+namespace RandomForest {
+
+template <typename Derived, typename ParamT, typename Splitter>
+class Node {
+public:
+    typedef typename Splitter::FeatureType FeatureType;
+    bool is_leaf;
+    size_t n_samples;
+    size_t depth;
+    typedef ParamT ParamType;
+    ParamType const* params;
+    Splitter splitter;
+    
+    boost::scoped_ptr<Derived> left;
+    boost::scoped_ptr<Derived> right;
+    std::vector<float> node_dist;
+
+    Node() : is_leaf(true), n_samples(0), depth(-1), params(0) {}
+    Node(size_t depth, ParamType const* params) :
+        is_leaf(true), n_samples(0), depth(depth), params(params)
+    {}
+
+    bool pure(DataView2D<int> labels, int* sample_idxes) const {
+        if (n_samples < 2) 
+            return true; // an empty node is by definition pure
+        int first_sample_idx = sample_idxes[0];
+        int seen_class = labels(first_sample_idx, 0);
+        // check if all classes are equal to the first class
+        for (size_t i_sample = 1; i_sample < n_samples; ++i_sample) {
+            int sample_idx = sample_idxes[i_sample];
+            if (labels(sample_idx, 0) != seen_class)
+                return false;
+        }
+        return true;
+    }
+
+    float const* votes() const {
+        return (float const*)&node_dist[0];
+    }
+
+    int partition_samples(DataView2D<FeatureType> samples, int* sample_idxes) {
+        // sort samples in bag so that left-samples precede right-samples
+        // works like std::partition 
+        int low  = 0;
+        int high = n_samples;
+
+        while (true) {
+            while (true) {
+                if (low == high) {
+                    return low;
+                } else if (!splitter.classify_sample(samples.row_pointer(sample_idxes[low]))) {
+                    ++low;
+                } else {
+                    break;
+                }
+            }
+            --high;
+            while (true) {
+                if (low == high) {
+                    return low;
+                } else if (splitter.classify_sample(samples.row_pointer(sample_idxes[high]))) {
+                    --high;
+                } else {
+                    break;
+                }
+            }
+            std::swap(sample_idxes[low], sample_idxes[high]);
+            ++low;
+        }
+    }
+
+    Derived const* split (FeatureType const* sample) const {
+        if (splitter.classify_sample(sample)) {
+            return right.get();
+        } else {
+            return left.get();
+        }
+    }
+
+    typedef std::list<Derived const*> NodeList;
+
+    NodeList get_all_childs() {
+        NodeList ret;
+        ret.push_back(this);
+        if (!is_leaf) {
+            NodeList left_childs  = left->get_all_childs();
+            ret.splice(ret.end(), left_childs);
+            NodeList right_childs = right->get_all_childs();
+            ret.splice(ret.end(), right_childs);
+        }
+        return ret;
+    }
+
+    template<typename SplitGenerator>
+    void determine_best_split(DataView2D<FeatureType> samples,
+                              DataView2D<int>         labels,
+                              int*                  sample_idxes,
+                              SplitGenerator        split_generator,
+                              RandomGen&            gen
+                              )
+    {
+        typename Splitter::FeatureClassData data_points;
+        init_feature_class_data(data_points, params->n_classes, n_samples);
+        float best_loss = std::numeric_limits<float>::infinity();
+
+        std::vector<uint64_t> classes_l;
+        std::vector<uint64_t> classes_r;
+
+        // pass information about data to split generator
+        split_generator.init(samples, 
+                             labels, 
+                             sample_idxes,
+                             n_samples,
+                             params->n_classes,
+                             gen);
+
+        size_t n_proposals = split_generator.num_proposals();
+
+        std::pair<FeatureType, float> results; // (threshold, loss)
+        
+        for (size_t i_proposal = 0; i_proposal < n_proposals; ++i_proposal) {
+            // generate proposal
+            Splitter split = split_generator.gen_proposal(gen);
+            // map samples to numbers using proposal
+            split.map_points(samples, labels, sample_idxes, n_samples, data_points);
+            // check best loss using this proposal
+            results = static_cast<Derived*>(this)->determine_best_threshold(data_points, classes_l, classes_r, gen);
+            if (results.second < best_loss) {
+                // Proposal resulted into new optimum
+                best_loss = results.second;
+                split.set_threshold(results.first);
+                splitter = split;
+            }
+        }
+    }
+
+    template<typename SplitGenerator>
+    void train(DataView2D<FeatureType> samples, 
+               DataView2D<int> labels, 
+               int* sample_idxes, 
+               size_t n_samples_, 
+               SplitGenerator const& split_generator,
+               RandomGen& gen
+               ) 
+    {
+        n_samples = n_samples_;
+        node_dist.resize(params->n_classes, 0.0f);
+        for (size_t i_sample = 0; i_sample < n_samples; ++i_sample) {
+            int label = labels(sample_idxes[i_sample], 0);
+            node_dist[label] += 1.0f;
+        }
+
+        if (n_samples != 0)
+          for (std::size_t i = 0; i < node_dist.size(); ++ i)
+            node_dist[i] /= n_samples;
+
+        bool do_split = // Only split if ...
+            (n_samples >= params->min_samples_per_node) && // enough samples are available
+            !pure(labels, sample_idxes) && // this node is not already pure
+            (depth < params->max_depth); // we did not reach max depth
+        if (!do_split) {
+            splitter.threshold = 0.0;
+            return;
+        }
+
+        is_leaf = false; 
+
+#if VERBOSE_NODE_LEARNING
+        std::printf("Determining the best split at depth %zu/%zu\n", depth, params->max_depth);
+#endif
+        determine_best_split(samples, labels, sample_idxes, split_generator, gen);
+
+        left.reset(new Derived(depth + 1, params));
+        right.reset(new Derived(depth + 1, params));
+
+        // sort samples in bag so that left-samples precede right-samples
+        int low = partition_samples(samples, sample_idxes);
+        int n_samples_left  = low;
+        int n_samples_right = n_samples - low;
+        int offset_left     = 0;
+        int offset_right    = low;
+#ifdef TREE_GRAPHVIZ_STREAM
+        if (depth <= TREE_GRAPHVIZ_MAX_DEPTH) {
+            TREE_GRAPHVIZ_STREAM << "p" << std::hex << (unsigned long)this       
+                << " -> " 
+                << "p" << std::hex << (unsigned long)left.get() 
+                << std::dec << " [label=\"" << n_samples_left << "\"];" <<  std::endl;
+            TREE_GRAPHVIZ_STREAM << "p" << std::hex << (unsigned long)this       
+                << " -> " 
+                << "p" << std::hex << (unsigned long)right.get() 
+                << std::dec << " [label=\"" << n_samples_right << "\"];" << std::endl;
+        }
+#endif
+        // train left and right side of split
+        left->train (samples, labels, sample_idxes + offset_left,  n_samples_left,  split_generator, gen);
+        right->train(samples, labels, sample_idxes + offset_right, n_samples_right, split_generator, gen);
+    } 
+
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /*version*/)
+    {
+        ar & BOOST_SERIALIZATION_NVP(is_leaf);
+        ar & BOOST_SERIALIZATION_NVP(n_samples);
+        ar & BOOST_SERIALIZATION_NVP(depth);
+        ar & BOOST_SERIALIZATION_NVP(params);
+        ar & BOOST_SERIALIZATION_NVP(splitter);
+        ar & BOOST_SERIALIZATION_NVP(node_dist);
+        ar & BOOST_SERIALIZATION_NVP(left);
+        ar & BOOST_SERIALIZATION_NVP(right);
+    }
+};
+
+}
+}
+
+}} // namespace CGAL::internal::
+
+#endif

Classification/include/CGAL/Classification/internal/auxiliary/random-forest/tree.hpp

@@ -0,0 +1,145 @@
+// Copyright (c) 2014 Stefan Walk
+//
+// Permission is hereby granted, free of charge, to any person obtaining a copy of
+// this software and associated documentation files (the "Software"), to deal in
+// the Software without restriction, including without limitation the rights to
+// use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies
+// of the Software, and to permit persons to whom the Software is furnished to do
+// so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included in all
+// copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+// SOFTWARE.
+//
+// $URL$
+// $Id$
+// SPDX-License-Identifier: MIT
+//
+// Author(s)     : Stefan Walk
+
+// Modifications from original library:
+//  * changed inclusion protection tag
+//  * moved to namespace CGAL::internal::
+
+#ifndef CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_TREE_H
+#define CGAL_INTERNAL_LIBLEARNING_RANDOMFOREST_TREE_H
+#include "../dataview.h"
+#include "common-libraries.hpp"
+#include <boost/serialization/scoped_ptr.hpp>
+
+namespace CGAL { namespace internal {
+
+namespace liblearning {
+namespace RandomForest {
+
+template <typename NodeT>
+class Tree {
+public:
+    typedef NodeT NodeType;
+    typedef typename NodeType::ParamType ParamType;
+    typedef typename NodeType::FeatureType FeatureType;
+    boost::scoped_ptr<NodeT> root_node;
+    ParamType const* params;
+
+    Tree() : params(0) {}
+    Tree(ParamType const* params) :
+        params(params)
+    { }
+
+    template<typename SplitGenerator>
+    void train(DataView2D<FeatureType> samples, 
+               DataView2D<int> labels, 
+               int* sample_idxes, 
+               size_t n_samples, 
+               SplitGenerator split_generator,
+               RandomGen const& gen
+               ) 
+    {
+        // copy generator
+        RandomGen my_gen = gen;
+        // initialize root node
+        root_node.reset(new NodeT(0, params));
+        // train root node (other notes get trained recursively)
+        root_node->train(samples, labels, sample_idxes, n_samples, split_generator, my_gen);
+    }
+    float const* evaluate(FeatureType const* sample) {
+        // start with root
+        NodeT const* node = root_node.get();
+        // split until leaf
+        while (node && !node->is_leaf) {
+            node = node->split(sample);
+        }
+        if (!node) {
+            // bug if this is reached
+            return 0;
+        } else {
+            return node->votes();
+        }
+    }
+    std::list<NodeT const*> all_nodes() {
+        return root_node->get_all_childs();
+    }
+
+#if 0
+    float similarity_endnode(float const* sample_1, float const* sample_2) {
+        NodeT const* node_1 = root_node.get();
+        NodeT const* node_2 = root_node.get();
+        while (!node_1->is_leaf) {
+            node_1 = node_1->split(sample_1);
+        }
+        while (!node_2->is_leaf) {
+            node_2 = node_2->split(sample_2);
+        }
+        if (node_1 == node_2) {
+            return 1.0;
+        } else {
+            return 0.0;
+        }
+    }
+    float similarity_path(float const* sample_1, float const* sample_2) {
+        NodeT const* node_1 = root_node.get();
+        NodeT const* node_2 = root_node.get();
+        int n_common = 0;
+        int n_1 = 1;
+        int n_2 = 1;
+
+        while (node_1 == node_2) {
+            if (node_1->is_leaf) {
+                return 1;
+            }
+            n_common++;
+            node_1 = node_1->split(sample_1);
+            node_2 = node_2->split(sample_2);
+        }
+        while (!node_1->is_leaf) {
+            node_1 = node_1->split(sample_1);
+            n_1++;
+        }
+        while (!node_2->is_leaf) {
+            node_2 = node_2->split(sample_2);
+            n_2++;
+        }
+        return n_common/sqrt((n_common + n_1)*(n_common + n_2));
+    }
+#endif
+    template <typename Archive>
+    void serialize(Archive& ar, unsigned /*version*/)
+    {
+        ar & BOOST_SERIALIZATION_NVP(params);
+        ar & BOOST_SERIALIZATION_NVP(root_node);
+    }
+};
+
+}
+}
+
+}} // namespace CGAL::internal::
+
+#endif

Classification/package_info/Classification/license.txt

@@ -1 +1,2 @@
 GPL (v3 or later)
+MIT/X11 (BSD like)

Classification/test/Classification/CMakeLists.txt

@@ -23,20 +23,26 @@ include( ${CGAL_USE_FILE} )
 # Boost and its components
 find_package( Boost REQUIRED )
 
+if ( NOT Boost_FOUND )
+  message(STATUS "This project requires the Boost library, and will not be compiled.")
+  return()  
+endif()
+
+find_package( Boost OPTIONAL_COMPONENTS serialization iostreams)
+
+if( WIN32 )
+# to avoid a warning with old cmake 
+  set(_Boost_BZIP2_HEADERS             "boost/iostreams/filter/bzip2.hpp")
+  set(_Boost_ZLIB_HEADERS              "boost/iostreams/filter/zlib.hpp")
+  find_package( Boost OPTIONAL_COMPONENTS zlib)
+endif()
+
 # Use Eigen
 find_package(Eigen3 3.1.0) #(requires 3.1.0 or greater)
 if (EIGEN3_FOUND)
   include( ${EIGEN3_USE_FILE} )
 endif()
 
-if ( NOT Boost_FOUND )
-
-  message(STATUS "This project requires the Boost library, and will not be compiled.")
-
-  return()  
-
-endif()
-
 find_package( TBB )
 
 # include for local directory
@@ -53,9 +59,36 @@ include( CGAL_CreateSingleSourceCGALProgram )
 
 #add_definitions("-DCGAL_DO_NOT_USE_BOYKOV_KOLMOGOROV_MAXFLOW_SOFTWARE")
 
+# Classification requires some C++11 features
 set(needed_cxx_features cxx_rvalue_references cxx_variadic_templates)
+
+# Libraries and flags
+set(classification_linked_libraries)
+
+if (Boost_SERIALIZATION_FOUND AND Boost_IOSTREAMS_FOUND)
+  set(classification_linked_libraries ${classification_linked_libraries}
+    ${Boost_SERIALIZATION_LIBRARY}
+    ${Boost_IOSTREAMS_LIBRARY})
+else()
+  message(STATUS "NOTICE: Boost Serialization or IO Streams not found, no test will be compiled.")
+  return()
+endif()
+
+if( WIN32 )
+  if (Boost_ZLIB_FOUND)
+    set(classification_linked_libraries ${classification_linked_libraries}
+      ${Boost_ZLIB_LIBRARY})
+  else()
+    message(STATUS "NOTICE: Boost ZLIB not found, no example will be compiled.")
+    return()
+  endif()
+endif()
+
 create_single_source_cgal_program( "test_classification_point_set.cpp" CXX_FEATURES ${needed_cxx_features} )
 
-if(TBB_FOUND)
+if(TARGET test_classification_point_set)
-  CGAL_target_use_TBB( test_classification_point_set )
+  target_link_libraries(test_classification_point_set PUBLIC ${classification_linked_libraries})
+  if (TBB_FOUND)
+    CGAL_target_use_TBB( test_classification_point_set )
+  endif()
 endif()

Documentation/doc/biblio/cgal_manual.bib

@@ -2125,6 +2125,16 @@ location    = {Salt Lake City, Utah, USA}
   ,update       = "97.04 kettner"
 }
 
+@Misc{ cgal:w-erftl-14,
+  key           = {ETHZRandomForest},
+  title         = {{ETH Zurich Random Forest Template Library}},
+  howpublished  = {Stefan Walk (ETH Zurich, Department of Civil,
+                  Environmental and Geomatic Engineering, Institute of
+                  Geodesy and Photogrammetry)},
+  url           = {http://www.prs.igp.ethz.ch/research/Source_code_and_datasets.html},
+  year          = 2014
+}
+
 @book{cgal:ww-smgd-02
    ,author = "Joe Warren and Henrik Weimer"
    ,title = "Subdivision Methods for Geometric Design"

HalfedgeDS/include/CGAL/boost/graph/graph_traits_HalfedgeDS.h

@@ -28,7 +28,7 @@
 
 #include <boost/config.hpp>
 #include <boost/iterator/iterator_adaptor.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/type_traits/remove_const.hpp>
 
 #include <boost/graph/graph_traits.hpp>

Kinetic_data_structures/archive/include/CGAL/Kinetic/Active_objects_set.h

@@ -31,7 +31,7 @@
 #include <map>
 #include <vector>
 #include <sstream>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <utility>
 
 namespace CGAL { namespace Kinetic {;

Kinetic_data_structures/include/CGAL/Kinetic/Active_objects_vector.h

@@ -32,7 +32,7 @@
 #include <vector>
 #include <sstream>
 #include <CGAL/iterator.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #if defined(BOOST_MSVC)
 #  pragma warning(push)

Mesh_2/include/CGAL/Mesh_2/Clusters.h

@@ -30,7 +30,7 @@
 #include <CGAL/Unique_hash_map.h>
 
 #include <utility>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 

Mesh_2/include/CGAL/Mesh_2/Refine_edges.h

@@ -34,7 +34,7 @@
 #include <utility>
 #include <iterator>
 #include <boost/iterator/filter_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 

Mesh_2/include/CGAL/Mesh_2/Refine_faces.h

@@ -27,7 +27,7 @@
 
 #include <CGAL/Mesh_2/Face_badness.h>
 #include <CGAL/Double_map.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <string>
 #include <sstream>

Mesh_3/benchmark/Mesh_3/StdAfx.h

@@ -29,7 +29,7 @@
 #include <boost/iterator/iterator_categories.hpp>
 #include <boost/iterator/iterator_facade.hpp>
 #include <boost/iterator/iterator_traits.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/iterator/zip_iterator.hpp>
 #include <boost/mpl/always.hpp>
 #include <boost/mpl/and.hpp>

Mesh_3/include/CGAL/Mesh_3/C3T3_helpers.h

@@ -47,7 +47,7 @@
 #include <boost/range/begin.hpp>
 #include <boost/range/end.hpp>
 #include <boost/optional.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/function_output_iterator.hpp>
 #include <boost/type_traits/is_same.hpp>
 #include <boost/type_traits/is_convertible.hpp>

Mesh_3/include/CGAL/Mesh_3/Protect_edges_sizing_field.h

@@ -78,7 +78,7 @@ namespace internal {
 #include <vector>
 #include <stack>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/optional.hpp>
 
 #include <CGAL/Mesh_3/utilities.h>

Mesh_3/include/CGAL/Mesh_3/Slivers_exuder.h

@@ -39,7 +39,7 @@
 #include <algorithm>
 #include <iomanip> // std::setprecision
 #include <iostream> // std::cerr/cout
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/function_output_iterator.hpp>
 #include <boost/bind.hpp>
 #include <boost/format.hpp>

Mesh_3/include/CGAL/Mesh_3/Triangulation_sizing_field.h

@@ -35,7 +35,7 @@
 #include <CGAL/Regular_triangulation_3.h>
 #include <CGAL/Triangulation_vertex_base_with_info_3.h>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 

Mesh_3/include/CGAL/Mesh_complex_3_in_triangulation_3.h

@@ -39,7 +39,7 @@
 #include <map>
 #include <boost/bimap/bimap.hpp>
 #include <boost/bimap/multiset_of.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/iterator/iterator_adaptor.hpp>
 #include <boost/mpl/if.hpp>
 #include <CGAL/internal/Mesh_3/Boundary_of_subdomain_of_complex_3_in_triangulation_3_to_off.h>

Mesh_3/include/CGAL/Polyhedral_complex_mesh_domain_3.h

@@ -48,7 +48,7 @@
 #include <boost/dynamic_bitset.hpp>
 #include <CGAL/boost/graph/split_graph_into_polylines.h>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/foreach.hpp>
 
 #include <string>

Mesh_3/include/CGAL/Polyhedral_mesh_domain_with_features_3.h

@@ -52,7 +52,7 @@
 #include <CGAL/boost/graph/split_graph_into_polylines.h>
 #include <CGAL/Default.h>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/foreach.hpp>
 
 #include <string>

NewKernel_d/include/CGAL/NewKernel_d/Types/Hyperplane.h

@@ -23,7 +23,7 @@
 #include <CGAL/enum.h>
 #include <CGAL/number_utils.h>
 #include <CGAL/NewKernel_d/store_kernel.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <CGAL/boost/iterator/counting_iterator.hpp>
 namespace CGAL {
 template <class R_> class Hyperplane {

Number_types/include/CGAL/Counted_number.h

@@ -29,7 +29,7 @@
 
 #include <CGAL/number_type_basic.h>
 #include <CGAL/atomic.h>
-#include <boost/iterator/transform_iterator.hpp> // for Root_of_selector
+#include <CGAL/boost/iterator/transform_iterator.hpp> // for Root_of_selector
 #include <iostream>
 
 namespace CGAL {

Number_types/include/CGAL/Lazy_exact_nt.h

@@ -30,7 +30,7 @@
 #include <CGAL/number_type_basic.h>
 #include <CGAL/assertions.h>
 
-#include <boost/iterator/transform_iterator.hpp> // for Root_of functor
+#include <CGAL/boost/iterator/transform_iterator.hpp> // for Root_of functor
 #include <boost/operators.hpp>
 #include <boost/type_traits/is_same.hpp>
 #include <boost/type_traits/is_arithmetic.hpp>

Optimal_transportation_reconstruction_2/include/CGAL/Optimal_transportation_reconstruction_2.h

@@ -42,7 +42,7 @@
 #include <boost/multi_index/mem_fun.hpp>
 #include <boost/multi_index/ordered_index.hpp>
 #include <boost/multi_index/identity.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/type_traits/is_float.hpp>
 
 namespace CGAL {

Point_set_shape_detection_3/include/CGAL/Shape_detection_3/Region_growing.h

@@ -34,8 +34,8 @@
 #include <CGAL/squared_distance_3.h>
 #include <CGAL/function.h>
 
-#include <boost/iterator/counting_iterator.hpp>
+#include <CGAL/boost/iterator/counting_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 namespace Shape_detection_3 {

Polygon_mesh_processing/benchmark/Polygon_mesh_processing/polygon_mesh_slicer.cpp

@@ -19,7 +19,7 @@
 
 #include <functional>
 #include <boost/foreach.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <fstream>
 

Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/Hole_filling/Triangulate_hole_polyline.h

@@ -37,7 +37,7 @@
 #include <stack>
 #include <map>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/unordered_set.hpp>
 
 namespace CGAL {

Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/intersection.h

@@ -29,7 +29,7 @@
 #include <boost/type_traits/is_same.hpp>
 #include <boost/utility/enable_if.hpp>
 #include <CGAL/Polygon_mesh_processing/bbox.h>
-#include <boost/iterator/counting_iterator.hpp>
+#include <CGAL/boost/iterator/counting_iterator.hpp>
 #include <boost/mpl/if.hpp>
 #include <CGAL/Polygon_mesh_processing/bbox.h>
 #include <CGAL/Polygon_mesh_processing/internal/named_params_helper.h>

Polyhedron/demo/Polyhedron/Plugins/Classification/CMakeLists.txt

@@ -2,23 +2,49 @@ include( polyhedron_demo_macros )
 
 if(EIGEN3_FOUND)
 
-  qt5_wrap_ui( classificationUI_FILES Classification_widget.ui)
+  set(classification_linked_libraries)
-  polyhedron_demo_plugin(classification_plugin Classification_plugin Point_set_item_classification ${classificationUI_FILES})
+  set(classification_compile_definitions)
-  find_package(OpenCV QUIET)
+  
-  if(OpenCV_FOUND)
+  find_package(Boost OPTIONAL_COMPONENTS serialization iostreams)
-    message(STATUS "Found OpenCV ${OpenCV_VERSION}")
+
-    include_directories(${OpenCV_INCLUDE_DIRS})
+  if( WIN32 )
-    target_link_libraries(classification_plugin PUBLIC scene_points_with_normal_item scene_polylines_item scene_polygon_soup_item scene_color_ramp ${OpenCV_LIBS})
+    find_package( Boost OPTIONAL_COMPONENTS zlib)
-    target_compile_definitions(classification_plugin PUBLIC "-DCGAL_LINKED_WITH_OPENCV")
+    if( Boost_ZLIB_FOUND )
+      set(classification_linked_libraries ${classification_linked_libraries}
+        ${Boost_ZLIB_LIBRARY})
+    endif()
+  endif()
+
+  if (Boost_SERIALIZATION_FOUND AND Boost_IOSTREAMS_FOUND AND (NOT WIN32 OR Boost_ZLIB_FOUND))
+    qt5_wrap_ui( classificationUI_FILES Classification_widget.ui)
+    polyhedron_demo_plugin(classification_plugin Classification_plugin Point_set_item_classification ${classificationUI_FILES})
+    
+    set(classification_linked_libraries ${classification_linked_libraries}
+      ${Boost_SERIALIZATION_LIBRARY}
+      ${Boost_IOSTREAMS_LIBRARY})
+    
+    find_package(OpenCV QUIET)
+    if (OpenCV_FOUND)
+      message(STATUS "Found OpenCV ${OpenCV_VERSION}")
+      include_directories(${OpenCV_INCLUDE_DIRS})
+      set(classification_linked_libraries ${classification_linked_libraries}
+        ${OpenCV_LIBS})
+      set(classification_compile_definitions ${classification_compile_definitions}
+        "-DCGAL_LINKED_WITH_OPENCV")
+    else()
+      message(STATUS "NOTICE: OpenCV was not found. OpenCV random forest predicate for classification won't be available.")
+    endif()
+    
+    target_link_libraries(classification_plugin PUBLIC scene_points_with_normal_item
+      scene_polylines_item scene_polygon_soup_item scene_color_ramp ${classification_linked_libraries})
+    target_compile_definitions(classification_plugin PUBLIC ${classification_compile_definitions})
   else()
-    target_link_libraries(classification_plugin PUBLIC scene_points_with_normal_item scene_polylines_item scene_polygon_soup_item scene_color_ramp)
+    message(STATUS "NOTICE: Boost Serialization or IO Streams or ZLIB not found. Classification plugin won't be available.")
-    message(STATUS "NOTICE: OpenCV was not found. Random forest predicate for classification will not be available.")
-  endif()
-  if(TBB_FOUND)
-    CGAL_target_use_TBB(classification_plugin)
   endif()
 
+
 else(EIGEN3_FOUND)
   message(STATUS "NOTICE: Eigen 3.1 (or greater) was not found. Classification plugin will not be available.")
 endif()
 
+

Polyhedron/demo/Polyhedron/Plugins/Classification/Classification_plugin.cpp

@@ -15,6 +15,9 @@
 #include <CGAL/Three/Polyhedron_demo_plugin_helper.h>
 
 #include <CGAL/Random.h>
+#include <CGAL/Real_timer.h>
+
+#include <QMultipleInputDialog.h>
 
 #include "ui_Classification_widget.h"
 
@@ -24,10 +27,8 @@
 #include <QCheckBox>
 #include <QInputDialog>
 #include <QMessageBox>
-#include <QFormLayout>
 #include <QSpinBox>
 #include <QDoubleSpinBox>
-#include <QDialogButtonBox>
 #include <QSlider>
 
 #include <map>
@@ -140,8 +141,12 @@ public:
     ui_widget.setupUi(dock_widget);
     addDockWidget(dock_widget);
 
-#ifndef CGAL_LINKED_WITH_OPENCV
+    ui_widget.classifier->addItem (tr("Random Forest (ETHZ)"));
-    ui_widget.classifier->removeItem(1);
+    
+#ifdef CGAL_LINKED_WITH_OPENCV
+    ui_widget.classifier->addItem (tr("Random Forest (OpenCV %1.%2)")
+                                   .arg(CV_MAJOR_VERSION)
+                                   .arg(CV_MINOR_VERSION));
 #endif
 
     color_att = QColor (75, 75, 77);
@@ -448,10 +453,31 @@ public Q_SLOTS:
         return; 
       }
 
-    QString filename = QFileDialog::getSaveFileName(mw,
+    QString filename;
-                                                    tr("Save classification configuration"),
+
-                                                    QString("config.xml"),
+    if (ui_widget.classifier->currentIndex() == 0)
-                                                    "Config file (*.xml);;");
+      filename = QFileDialog::getSaveFileName(mw,
+                                              tr("Save classification configuration"),
+                                              tr("%1 (CGAL classif config).xml").arg(classif->item()->name()),
+                                              "CGAL classification configuration (*.xml);;");
+    else if (ui_widget.classifier->currentIndex() == 1)
+      filename = QFileDialog::getSaveFileName(mw,
+                                              tr("Save classification configuration"),
+                                              tr("%1 (ETHZ random forest config).gz").arg(classif->item()->name()),
+                                              "Compressed ETHZ random forest configuration (*.gz);;");
+#ifdef CGAL_LINKED_WITH_OPENCV
+    else if (ui_widget.classifier->currentIndex() == 2)
+      filename = QFileDialog::getSaveFileName(mw,
+                                              tr("Save classification configuration"),
+                                              tr("%1 (OpenCV %2.%3 random forest config).xml")
+                                              .arg(classif->item()->name())
+                                              .arg(CV_MAJOR_VERSION)
+                                              .arg(CV_MINOR_VERSION),
+                                              tr("OpenCV %2.%3 random forest configuration (*.xml);;")
+                                              .arg(CV_MAJOR_VERSION)
+                                              .arg(CV_MINOR_VERSION));
+#endif
+    
     if (filename == QString())
       return;
 
@@ -474,10 +500,29 @@ public Q_SLOTS:
         print_message("Error: there is no point set classification item!");
         return; 
       }
-    QString filename = QFileDialog::getOpenFileName(mw,
+    QString filename;
-                                                    tr("Open classification configuration"),
+
-                                                    ".",
+    if (ui_widget.classifier->currentIndex() == 0)
-                                                    "Config file (*.xml);;All Files (*)");
+      filename = QFileDialog::getOpenFileName(mw,
+                                              tr("Open CGAL classification configuration"),
+                                              ".",
+                                              "CGAL classification configuration (*.xml);;All Files (*)");
+    else if (ui_widget.classifier->currentIndex() == 1)
+      filename = QFileDialog::getOpenFileName(mw,
+                                              tr("Open ETHZ random forest configuration"),
+                                              ".",
+                                              "Compressed ETHZ random forest configuration (*.gz);;All Files (*)");
+#ifdef CGAL_LINKED_WITH_OPENCV
+    else if (ui_widget.classifier->currentIndex() == 2)
+      filename = QFileDialog::getOpenFileName(mw,
+                                              tr("Open OpenCV %2.%3 random forest configuration")
+                                              .arg(CV_MAJOR_VERSION)
+                                              .arg(CV_MINOR_VERSION),
+                                              ".",
+                                              tr("OpenCV %2.%3 random forest configuration (*.xml);;All Files (*)")
+                                              .arg(CV_MAJOR_VERSION)
+                                              .arg(CV_MINOR_VERSION));
+#endif
 
     if (filename == QString())
       return;
@@ -515,7 +560,11 @@ public Q_SLOTS:
         return; 
       }
     QApplication::setOverrideCursor(Qt::WaitCursor);
+    CGAL::Real_timer t;
+    t.start();
     run (classif, 0);
+    t.stop();
+    std::cerr << "Raw classification computed in " << t.time() << " second(s)" << std::endl;
     QApplication::restoreOverrideCursor();
     item_changed(classif->item());
   }
@@ -531,10 +580,11 @@ public Q_SLOTS:
       }
 
     QApplication::setOverrideCursor(Qt::WaitCursor);
-    QTime time;
+    CGAL::Real_timer t;
-    time.start();
+    t.start();
     run (classif, 1);
-    std::cerr << "Smoothed classification computed in " << time.elapsed() / 1000 << " second(s)" << std::endl;
+    t.stop();
+    std::cerr << "Smoothed classification computed in " << t.time() << " second(s)" << std::endl;
     QApplication::restoreOverrideCursor();
     item_changed(classif->item());
   }
@@ -549,36 +599,25 @@ public Q_SLOTS:
         return; 
       }
 
-    QDialog dialog(mw);
+    QMultipleInputDialog dialog ("Classify with Graph Cut", mw);
-    dialog.setWindowTitle ("Classify with Graph Cut");
+    QSpinBox* subdivisions = dialog.add<QSpinBox> ("Number of subdivisons: ");
-    QFormLayout form (&dialog);
+    subdivisions->setRange (1, 9999);
-    QString label_sub = QString("Number of subdivisions: ");
+    subdivisions->setValue (16);
-    QSpinBox subdivisions (&dialog);
-    subdivisions.setRange (1, 9999);
-    subdivisions.setValue (16);
-    form.addRow(label_sub, &subdivisions);
 
-    QString label_smooth = QString("Regularization weight: ");
+    QDoubleSpinBox* smoothing = dialog.add<QDoubleSpinBox> ("Regularization weight: ");
-    QDoubleSpinBox smoothing (&dialog);
+    smoothing->setRange (0.0, 100.0);
-    smoothing.setRange (0.0, 100.0);
+    smoothing->setValue (0.5);
-    smoothing.setValue (0.5);
+    smoothing->setSingleStep (0.1);
-    smoothing.setSingleStep (0.1);
-    form.addRow(label_smooth, &smoothing);
-
-    QDialogButtonBox oknotok (QDialogButtonBox::Ok | QDialogButtonBox::Cancel,
-                              Qt::Horizontal, &dialog);
-    form.addRow (&oknotok);
-    QObject::connect (&oknotok, SIGNAL(accepted()), &dialog, SLOT(accept()));
-    QObject::connect (&oknotok, SIGNAL(rejected()), &dialog, SLOT(reject()));
 
     if (dialog.exec() != QDialog::Accepted)
       return;
     
     QApplication::setOverrideCursor(Qt::WaitCursor);
-    QTime time;
+    CGAL::Real_timer t;
-    time.start();
+    t.start();
-    run (classif, 2, std::size_t(subdivisions.value()), smoothing.value());
+    run (classif, 2, std::size_t(subdivisions->value()), smoothing->value());
-    std::cerr << "Graphcut classification computed in " << time.elapsed() / 1000 << " second(s)" << std::endl;
+    t.stop();
+    std::cerr << "Graph Cut classification computed in " << t.time() << " second(s)" << std::endl;
     QApplication::restoreOverrideCursor();
     item_changed(classif->item());
   }
@@ -742,6 +781,20 @@ public Q_SLOTS:
     classif->validate_selection();
   }
 
+  void on_select_random_region_clicked()
+  {
+    Item_classification_base* classif
+      = get_classification();
+    if(!classif)
+      {
+        print_message("Error: there is no point set classification item!");
+        return; 
+      }
+
+    classif->select_random_region();
+    item_changed(classif->item());
+  }
+
   void on_train_clicked()
   {
     Item_classification_base* classif
@@ -762,6 +815,8 @@ public Q_SLOTS:
       }
 
     int nb_trials = 0;
+    int num_trees = 0;
+    int max_depth = 0;
 
     if (ui_widget.classifier->currentIndex() == 0)
     {
@@ -773,9 +828,29 @@ public Q_SLOTS:
       if (!ok)
         return;
     }
+    else
+    {
+      QMultipleInputDialog dialog ("Train Random Forest Classifier", mw);
+      QSpinBox* trees = dialog.add<QSpinBox> ("Number of trees: ");
+      trees->setRange (1, 9999);
+      trees->setValue (25);
+      QSpinBox* depth = dialog.add<QSpinBox> ("Maximum depth of tree: ");
+      depth->setRange (1, 9999);
+      depth->setValue (20);
+
+      if (dialog.exec() != QDialog::Accepted)
+        return;
+      num_trees = trees->value();
+      max_depth = depth->value();
+    }
     
     QApplication::setOverrideCursor(Qt::WaitCursor);
-    classif->train(ui_widget.classifier->currentIndex(), nb_trials);
+    CGAL::Real_timer t;
+    t.start();
+    classif->train(ui_widget.classifier->currentIndex(), nb_trials,
+                   num_trees, max_depth);
+    t.stop();
+    std::cerr << "Done in " << t.time() << " second(s)" << std::endl;
     QApplication::restoreOverrideCursor();
     update_plugin_from_item(classif);
   }

Polyhedron/demo/Polyhedron/Plugins/Classification/Classification_widget.ui

@@ -105,11 +105,6 @@
           <string>Sum of Weighted Features</string>
          </property>
         </item>
-        <item>
-         <property name="text">
-          <string>Random Forest</string>
-         </property>
-        </item>
        </widget>
       </item>
      </layout>

Polyhedron/demo/Polyhedron/Plugins/Classification/Item_classification_base.h

@@ -8,9 +8,10 @@
 #include <CGAL/Classification/Feature_set.h>
 #include <CGAL/Classification/Label_set.h>
 #include <CGAL/Classification/Sum_of_weighted_features_classifier.h>
+#include <CGAL/Classification/ETHZ_random_forest_classifier.h>
 
 #ifdef CGAL_LINKED_WITH_OPENCV
-#include <CGAL/Classification/Random_forest_classifier.h>
+#include <CGAL/Classification/OpenCV_random_forest_classifier.h>
 #endif
 
 class Item_classification_base
@@ -21,9 +22,10 @@ public:
   typedef CGAL::Classification::Label_set   Label_set;
   typedef CGAL::Classification::Feature_set Feature_set;
   typedef CGAL::Classification::Sum_of_weighted_features_classifier Sum_of_weighted_features;
+  typedef CGAL::Classification::ETHZ_random_forest_classifier ETHZ_random_forest;
 
 #ifdef CGAL_LINKED_WITH_OPENCV
-  typedef CGAL::Classification::Random_forest_classifier Random_forest;
+  typedef CGAL::Classification::OpenCV_random_forest_classifier Random_forest;
 #endif
   
 public:
@@ -40,8 +42,10 @@ public:
   virtual void reset_training_set (const char* name) = 0;
   virtual void reset_training_sets() = 0;
 
+  virtual void select_random_region() = 0;
   virtual void validate_selection () = 0;
-  virtual void train(int classifier, unsigned int nb_trials) = 0;
+  virtual void train(int classifier, unsigned int nb_trials,
+                     std::size_t num_trees, std::size_t max_depth) = 0;
   virtual bool run (int method, int classifier, std::size_t subdivisions, double smoothing) = 0;
   
   virtual void update_color () = 0;
@@ -66,9 +70,13 @@ public:
   {
     m_labels.add(name);
     m_label_colors.push_back (color);
+    
     delete m_sowf;
     m_sowf = new Sum_of_weighted_features (m_labels, m_features);
 
+    delete m_ethz;
+    m_ethz = new ETHZ_random_forest (m_labels, m_features);
+
 #ifdef CGAL_LINKED_WITH_OPENCV
     delete m_random_forest;
     m_random_forest = new Random_forest (m_labels, m_features);
@@ -86,6 +94,9 @@ public:
     delete m_sowf;
     m_sowf = new Sum_of_weighted_features (m_labels, m_features);
 
+    delete m_ethz;
+    m_ethz = new ETHZ_random_forest (m_labels, m_features);
+
 #ifdef CGAL_LINKED_WITH_OPENCV
     delete m_random_forest;
     m_random_forest = new Random_forest (m_labels, m_features);
@@ -118,6 +129,11 @@ public:
       std::ofstream f (filename);
       m_sowf->save_configuration (f);
     }
+    else if (classifier == 1)
+    {
+      std::ofstream f (filename, std::ios_base::out | std::ios_base::binary);
+      m_ethz->save_configuration (f);
+    }
     else
     {
 #ifdef CGAL_LINKED_WITH_OPENCV
@@ -138,6 +154,11 @@ public:
       std::ifstream f (filename);
       m_sowf->load_configuration (f, true);
     }
+    else if (classifier == 1)
+    {
+      std::ifstream f (filename, std::ios_base::in | std::ios_base::binary);
+      m_ethz->load_configuration (f);
+    }
     else
     {
 #ifdef CGAL_LINKED_WITH_OPENCV
@@ -165,6 +186,7 @@ protected:
   Feature_set m_features;
   std::vector<QColor> m_label_colors;
   Sum_of_weighted_features* m_sowf;
+  ETHZ_random_forest* m_ethz;
 #ifdef CGAL_LINKED_WITH_OPENCV
   Random_forest* m_random_forest;
 #endif

Polyhedron/demo/Polyhedron/Plugins/Classification/Point_set_item_classification.cpp

@@ -33,9 +33,31 @@ Point_set_item_classification::Point_set_item_classification(Scene_points_with_n
   training_found = !training_found; // add_property_map returns false if 
   classif_found = !classif_found;   // property was already there
 
+  bool las_found = false;
+
+  if (!classif_found)
+  {
+    Point_set::Property_map<unsigned char> las_classif;
+    boost::tie (las_classif, las_found) = m_points->point_set()->property_map<unsigned char>("classification");
+    if (las_found)
+    {
+      for (Point_set::const_iterator it = m_points->point_set()->begin();
+           it != m_points->point_set()->first_selected(); ++ it)
+      {
+        unsigned char uc = las_classif[*it];
+        m_classif[*it] = int(uc);
+        if (!training_found)
+          m_training[*it] = int(uc);
+      }
+      m_points->point_set()->remove_property_map (las_classif);
+      classif_found = true;
+      training_found = true;
+    }
+  }
+    
   if (training_found || classif_found)
   {
-    int lmax = 0;
+    std::vector<int> used_indices;
     
     for (Point_set::const_iterator it = m_points->point_set()->begin();
          it != m_points->point_set()->first_selected(); ++ it)
@@ -43,52 +65,165 @@ Point_set_item_classification::Point_set_item_classification(Scene_points_with_n
       if (training_found)
       {
         int l = m_training[*it];
-        lmax = (std::max)(l, lmax);
+        if (l >= 0)
+        {
+          if (std::size_t(l) >= used_indices.size())
+            used_indices.resize(std::size_t(l + 1), -1);
+          used_indices[std::size_t(l)] = 0;
+        }
       }
       if (classif_found)
       {
         int l = m_classif[*it];
-        lmax = (std::max)(l, lmax);
+        if (l >= 0)
+        {
+          if (std::size_t(l) >= used_indices.size())
+            used_indices.resize(std::size_t(l + 1), -1);
+          used_indices[std::size_t(l)] = 0;
+        }
       }
     }
 
-    // Try to recover label names from PLY comments
+    // map indices to filtered indices
-    std::map<int, std::string> label_names;
+    int current_idx = 0;
-    const std::string& comments = m_points->comments();
+    for (std::size_t i = 0; i < used_indices.size(); ++ i)
-    std::istringstream stream (comments);
-    std::string line;
-    while (getline(stream, line))
     {
-      std::stringstream iss (line);
+      if (las_found && (i < 2))
-      std::string tag;
-      if (iss >> tag && tag == "label")
       {
-        int idx;
+        used_indices[i] = -1;
-        std::string name;
+        continue;
-        if (iss >> idx >> name)
+      }
-          label_names.insert (std::make_pair (idx, name));
+      if (used_indices[i] == -1)
+        continue;
+
+      used_indices[i] = current_idx;
+      ++ current_idx;
+    }
+
+    if (std::size_t(current_idx) != used_indices.size()) // Empty indices -> reorder indices in point set
+    {
+      for (Point_set::const_iterator it = m_points->point_set()->begin();
+           it != m_points->point_set()->first_selected(); ++ it)
+      {
+        if (training_found)
+        {
+          if (las_found && (m_training[*it] == 0 || m_training[*it] == 1)) // Unclassified class in LAS
+            m_training[*it] = -1;
+          else if (m_training[*it] != -1)
+            m_training[*it] = used_indices[std::size_t(m_training[*it])];
+        }
+        if (classif_found)
+        {
+          if (las_found && (m_classif[*it] == 0 || m_classif[*it] == 1)) // Unclassified class in LAS
+            m_classif[*it] = -1;
+          else if (m_classif[*it] != -1)
+            m_classif[*it] = used_indices[std::size_t(m_classif[*it])];
+        }
       }
     }
     
-    for (int i = 0; i <= lmax; ++ i)
+    std::map<int, std::string> label_names;
+    if (las_found) // Use LAS standard
     {
+      // label_names.insert (std::make_pair (0, std::string("never_clfied")));
+      // label_names.insert (std::make_pair (1, std::string("unclassified")));
+      label_names.insert (std::make_pair (2, std::string("ground")));
+      label_names.insert (std::make_pair (3, std::string("low_veget")));
+      label_names.insert (std::make_pair (4, std::string("med_veget")));
+      label_names.insert (std::make_pair (5, std::string("high_veget")));
+      label_names.insert (std::make_pair (6, std::string("building")));
+      label_names.insert (std::make_pair (7, std::string("noise")));
+      label_names.insert (std::make_pair (8, std::string("reserved")));
+      label_names.insert (std::make_pair (9, std::string("water")));
+      label_names.insert (std::make_pair (10, std::string("rail")));
+      label_names.insert (std::make_pair (11, std::string("road_surface")));
+      label_names.insert (std::make_pair (12, std::string("reserved_2")));
+      label_names.insert (std::make_pair (13, std::string("wire_guard")));
+      label_names.insert (std::make_pair (14, std::string("wire_conduct")));
+      label_names.insert (std::make_pair (15, std::string("trans_tower")));
+      label_names.insert (std::make_pair (16, std::string("wire_connect")));
+      label_names.insert (std::make_pair (17, std::string("bridge_deck")));
+      label_names.insert (std::make_pair (18, std::string("high_noise")));
+    }
+    else // Try to recover label names from PLY comments
+    {
+
+      const std::string& comments = m_points->comments();
+      std::istringstream stream (comments);
+      std::string line;
+      while (getline(stream, line))
+      {
+        std::stringstream iss (line);
+        std::string tag;
+        if (iss >> tag && tag == "label")
+        {
+          int idx;
+          std::string name;
+          if (iss >> idx >> name)
+            label_names.insert (std::make_pair (idx, name));
+        }
+      }
+    }
+    
+    for (std::size_t i = 0; i < used_indices.size(); ++ i)
+    {
+      if (used_indices[i] == -1)
+        continue;
+
+      Label_handle new_label;
       std::map<int, std::string>::iterator found
-        = label_names.find (i);
+        = label_names.find (int(i));
       if (found != label_names.end())
-        m_labels.add(found->second.c_str());
+        new_label = m_labels.add(found->second.c_str());
       else
       {
         std::ostringstream oss;
         oss << "label_" << i;
-        m_labels.add(oss.str().c_str());
+        new_label = m_labels.add(oss.str().c_str());
       }
 
-      CGAL::Classification::HSV_Color hsv;
+      QColor color (64 + rand() % 192,
-      hsv[0] = 360. * (i / double(lmax + 1));
+                    64 + rand() % 192,
-      hsv[1] = 76.;
+                    64 + rand() % 192);
-      hsv[2] = 85.;
+      
-      Color rgb = CGAL::Classification::hsv_to_rgb(hsv);
+      if (new_label->name() == "ground")
-      m_label_colors.push_back (QColor(rgb[0], rgb[1], rgb[2]));
+        color = QColor (186, 189, 182);
+      else if (new_label->name() == "low_veget")
+        color = QColor (78, 154, 6);
+      else if (new_label->name() == "med_veget"
+               || new_label->name() == "vegetation")
+        color = QColor (138, 226, 52);
+      else if (new_label->name() == "high_veget")
+        color = QColor (204, 255, 201);
+      else if (new_label->name() == "building"
+               || new_label->name() == "roof")
+        color = QColor (245, 121, 0);
+      else if (new_label->name() == "noise")
+        color = QColor (0, 0, 0);
+      else if (new_label->name() == "reserved")
+        color = QColor (233, 185, 110);
+      else if (new_label->name() == "water")
+        color = QColor (114, 159, 207);
+      else if (new_label->name() == "rail")
+        color = QColor (136, 46, 25);
+      else if (new_label->name() == "road_surface")
+        color = QColor (56, 56, 56);
+      else if (new_label->name() == "reserved_2")
+        color = QColor (193, 138, 51);
+      else if (new_label->name() == "wire_guard")
+        color = QColor (37, 61, 136);
+      else if (new_label->name() == "wire_conduct")
+        color = QColor (173, 127, 168);
+      else if (new_label->name() == "trans_tower")
+        color = QColor (136, 138, 133);
+      else if (new_label->name() == "wire_connect")
+        color = QColor (145, 64, 236);
+      else if (new_label->name() == "bridge_deck")
+        color = QColor (213, 93, 93);
+      else if (new_label->name() == "high_noise")
+        color = QColor (255, 0, 0);
+      
+      m_label_colors.push_back (color);
     }
   }
   else
@@ -98,15 +233,16 @@ Point_set_item_classification::Point_set_item_classification(Scene_points_with_n
     m_labels.add("roof");
     m_labels.add("facade");
 
-    m_label_colors.push_back (QColor(245, 180, 0));
+    m_label_colors.push_back (QColor(186, 189, 182));
-    m_label_colors.push_back (QColor(0, 255, 27));
+    m_label_colors.push_back (QColor(138, 226, 52));
-    m_label_colors.push_back (QColor(255, 0, 170));
+    m_label_colors.push_back (QColor(245, 121, 0));
-    m_label_colors.push_back (QColor(100, 0, 255));
+    m_label_colors.push_back (QColor(233, 185, 110));
   }
   
   update_comments_of_point_set_item();
 
   m_sowf = new Sum_of_weighted_features (m_labels, m_features);
+  m_ethz = new ETHZ_random_forest (m_labels, m_features);
 #ifdef CGAL_LINKED_WITH_OPENCV
   m_random_forest = new Random_forest (m_labels, m_features);
 #endif
@@ -117,6 +253,8 @@ Point_set_item_classification::~Point_set_item_classification()
 {
   if (m_sowf != NULL)
     delete m_sowf;
+  if (m_ethz != NULL)
+    delete m_ethz;
 #ifdef CGAL_LINKED_WITH_OPENCV
   if (m_random_forest != NULL)
     delete m_random_forest;
@@ -335,6 +473,10 @@ void Point_set_item_classification::compute_features (std::size_t nb_scales)
   Point_set::Property_map<boost::uint8_t> echo_map;
   bool echo;
   boost::tie (echo_map, echo) = m_points->point_set()->template property_map<boost::uint8_t>("echo");
+  if (!echo)
+    boost::tie (echo_map, echo) = m_points->point_set()->template property_map<boost::uint8_t>("number_of_returns");
+
+  add_remaining_point_set_properties_as_features();
 
   if (!normals && !colors && !echo)
     m_generator = new Generator (m_features, *(m_points->point_set()), m_points->point_set()->point_map(), nb_scales);
@@ -362,6 +504,8 @@ void Point_set_item_classification::compute_features (std::size_t nb_scales)
 
   delete m_sowf;
   m_sowf = new Sum_of_weighted_features (m_labels, m_features);
+  delete m_ethz;
+  m_ethz = new ETHZ_random_forest (m_labels, m_features);
 #ifdef CGAL_LINKED_WITH_OPENCV
   delete m_random_forest;
   m_random_forest = new Random_forest (m_labels, m_features);
@@ -369,9 +513,96 @@ void Point_set_item_classification::compute_features (std::size_t nb_scales)
   std::cerr << "Features = " << m_features.size() << std::endl;
 }
 
+void Point_set_item_classification::select_random_region()
+{
+  m_points->point_set()->reset_indices();
   
+  std::size_t scale = (rand() % m_generator->number_of_scales());
 
-void Point_set_item_classification::train(int classifier, unsigned int nb_trials)
+  bool use_grid = (rand() % 2);
+
+  std::vector<std::size_t> selected;
+  
+  if (use_grid)
+  {
+    std::size_t x = (rand() % m_generator->grid(scale).width());
+    std::size_t y = (rand() % m_generator->grid(scale).height());
+    std::copy (m_generator->grid(scale).indices_begin(x,y),
+               m_generator->grid(scale).indices_end(x,y),
+               std::back_inserter (selected));
+  }
+  else
+  {
+    m_generator->neighborhood(0).sphere_neighbor_query (m_generator->radius_neighbors(scale))
+      (*(m_points->point_set()->points().begin() + (rand() % m_points->point_set()->size())),
+       std::back_inserter (selected));
+  }
+
+  if (selected.empty())
+    return;
+
+  std::sort (selected.begin(), selected.end());
+  std::size_t current_idx = 0;
+
+  std::vector<std::size_t> unselected;
+
+  for (Point_set::const_iterator it = m_points->point_set()->begin();
+       it != m_points->point_set()->end(); ++ it)
+    if (std::size_t(*it) == selected[current_idx])
+      current_idx ++;
+    else
+      unselected.push_back (*it);
+  
+  for (std::size_t i = 0; i < unselected.size(); ++ i)
+    *(m_points->point_set()->begin() + i) = unselected[i];
+  for (std::size_t i = 0; i < selected.size(); ++ i)
+    *(m_points->point_set()->begin() + (unselected.size() + i)) = selected[i];
+
+  m_points->point_set()->set_first_selected
+    (m_points->point_set()->begin() + unselected.size());
+
+}
+
+void Point_set_item_classification::add_remaining_point_set_properties_as_features()
+{
+  const std::vector<std::string>& prop = m_points->point_set()->base().properties();
+  
+  for (std::size_t i = 0; i < prop.size(); ++ i)
+  {
+    if (prop[i] == "index" ||
+        prop[i] == "point" ||
+        prop[i] == "normal" ||
+        prop[i] == "echo" ||
+        prop[i] == "number_of_returns" ||
+        prop[i] == "training" ||
+        prop[i] == "label" ||
+        prop[i] == "classification" ||
+        prop[i] == "real_color" ||
+        prop[i] == "red" || prop[i] == "green" || prop[i] == "blue" ||
+        prop[i] == "r" || prop[i] == "g" || prop[i] == "b")
+      continue;
+
+    if (try_adding_simple_feature<boost::int8_t>(prop[i]))
+      continue;
+    if (try_adding_simple_feature<boost::uint8_t>(prop[i]))
+      continue;
+    if (try_adding_simple_feature<boost::int16_t>(prop[i]))
+      continue;
+    if (try_adding_simple_feature<boost::uint16_t>(prop[i]))
+      continue;
+    if (try_adding_simple_feature<boost::int32_t>(prop[i]))
+      continue;
+    if (try_adding_simple_feature<boost::uint32_t>(prop[i]))
+      continue;
+    if (try_adding_simple_feature<float>(prop[i]))
+      continue;
+    if (try_adding_simple_feature<double>(prop[i]))
+      continue;
+  }
+}
+
+void Point_set_item_classification::train(int classifier, unsigned int nb_trials,
+                                          std::size_t num_trees, std::size_t max_depth)
 {
   if (m_features.size() == 0)
     {
@@ -394,9 +625,21 @@ void Point_set_item_classification::train(int classifier, unsigned int nb_trials
                                                        m_labels, *m_sowf,
                                                        indices);
     }
+  else if (classifier == 1)
+  {
+    m_ethz->train(training, true, num_trees, max_depth);
+    CGAL::Classification::classify<Concurrency_tag> (*(m_points->point_set()),
+                                                     m_labels, *m_ethz,
+                                                     indices);
+  }
   else
     {
 #ifdef CGAL_LINKED_WITH_OPENCV
+      if (m_random_forest != NULL)
+        delete m_random_forest;
+      m_random_forest = new Random_forest (m_labels, m_features,
+                                           int(max_depth), 5, 15,
+                                           int(num_trees));
       m_random_forest->train (training);
       CGAL::Classification::classify<Concurrency_tag> (*(m_points->point_set()),
                                                        m_labels, *m_random_forest,
@@ -424,6 +667,8 @@ bool Point_set_item_classification::run (int method, int classifier,
 
   if (classifier == 0)
     run (method, *m_sowf, subdivisions, smoothing);
+  else if (classifier == 1)
+    run (method, *m_ethz, subdivisions, smoothing);
 #ifdef CGAL_LINKED_WITH_OPENCV
   else
     run (method, *m_random_forest, subdivisions, smoothing);

Polyhedron/demo/Polyhedron/Plugins/Classification/Point_set_item_classification.h

@@ -46,6 +46,23 @@ class Point_set_item_classification : public Item_classification_base
   void erase_item() { m_points = NULL; }
 
   void compute_features (std::size_t nb_scales);
+  void add_remaining_point_set_properties_as_features();
+  
+  void select_random_region();
+
+  template <typename Type>
+  bool try_adding_simple_feature (const std::string& name)
+  {
+    typedef typename Point_set::template Property_map<Type> Pmap;
+    bool okay = false;
+    Pmap pmap;
+    boost::tie (pmap, okay) = m_points->point_set()->template property_map<Type>(name.c_str());
+    if (okay)
+      m_features.template add<CGAL::Classification::Feature::Simple_feature <Point_set, Pmap> >
+        (*(m_points->point_set()), pmap, name.c_str());
+
+    return okay;
+  }
   
   void add_selection_to_training_set (const char* name)
   {
@@ -91,7 +108,8 @@ class Point_set_item_classification : public Item_classification_base
     if (m_index_color == 1 || m_index_color == 2)
       change_color (m_index_color);
   }
-  void train(int classifier, unsigned int nb_trials);
+  void train(int classifier, unsigned int nb_trials,
+             std::size_t num_trees, std::size_t max_depth);
   bool run (int method, int classifier, std::size_t subdivisions, double smoothing);
 
   void update_color () { change_color (m_index_color); }

Polyhedron/demo/Polyhedron/Plugins/Convex_hull/Convex_hull_plugin.cpp

@@ -16,7 +16,7 @@
 
 #include <CGAL/convex_hull_3.h>
 #include <CGAL/boost/graph/copy_face_graph.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 using namespace CGAL::Three;
 class Polyhedron_demo_convex_hull_plugin : 
   public QObject,

Polyhedron/demo/Polyhedron/Plugins/PMP/Hole_filling_polyline_plugin.cpp

@@ -26,7 +26,7 @@
 #include <algorithm>
 
 #include <boost/function_output_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 
 struct Face : public CGAL::cpp11::array<int,3>

Polyhedron/demo/Polyhedron/Plugins/PMP/Point_inside_polyhedron_plugin.cpp

@@ -24,7 +24,7 @@
 #include <vector>
 #include <algorithm>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/optional/optional.hpp>
 
 using namespace CGAL::Three;

Polyhedron/demo/Polyhedron/include/QMultipleInputDialog.h

@@ -0,0 +1,53 @@
+// Note: this structure is inspired from the QInputDialog that allows
+// the user to easily create a dialog to get one value (integer,
+// double or string). This structure allows the user to easily create
+// a form inside a QDialog to get as manu values as needed.
+
+#ifndef CGAL_QMULTIPLEINPUTDIALOG_H
+#define CGAL_QMULTIPLEINPUTDIALOG_H
+
+#include <QDialog>
+#include <QFormLayout>
+#include <QDialogButtonBox>
+
+class QMultipleInputDialog
+{
+  QDialog* dialog;
+  QFormLayout* form;
+public:
+  QMultipleInputDialog (const char* name, QWidget* parent)
+  {
+    dialog = new QDialog (parent);
+    dialog->setWindowTitle (name);
+    form = new QFormLayout(dialog);
+  }
+  ~QMultipleInputDialog ()
+  {
+    delete dialog;
+  }
+
+  template <typename QObjectType>
+  QObjectType* add (const char* name)
+  {
+    QObjectType* out = new QObjectType (dialog);
+    form->addRow (QString(name), out);
+    return out;
+  }
+
+  int exec()
+  {
+    QDialogButtonBox* oknotok = new QDialogButtonBox
+      (QDialogButtonBox::Ok | QDialogButtonBox::Cancel,
+       Qt::Horizontal, dialog);
+      
+    form->addRow (oknotok);
+    QObject::connect (oknotok, SIGNAL(accepted()), dialog, SLOT(accept()));
+    QObject::connect (oknotok, SIGNAL(rejected()), dialog, SLOT(reject()));
+
+    return dialog->exec();
+  }
+};
+
+
+
+#endif // CGAL_QMULTIPLEINPUTDIALOG_H

Polynomial/include/CGAL/Polynomial.h

@@ -40,7 +40,7 @@
 #include <iostream>
 #include <CGAL/basic.h>
 #include <CGAL/Handle_with_policy.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/mpl/if.hpp>
 #include <CGAL/Flattening_iterator.h>
 

Polynomial/include/CGAL/Polynomial_traits_d.h

@@ -43,7 +43,7 @@
 #include <CGAL/Polynomial/subresultants.h>
 #include <CGAL/Polynomial/sturm_habicht_sequence.h>
 
-#include <boost/iterator/transform_iterator.hpp> 
+#include <CGAL/boost/iterator/transform_iterator.hpp> 
 
 
 #define CGAL_POLYNOMIAL_TRAITS_D_BASE_TYPEDEFS                          \

Principal_component_analysis/demo/Principal_component_analysis/types.h

@@ -3,7 +3,7 @@
 
 #include <CGAL/Simple_cartesian.h>
 #include <CGAL/Polyhedron_3.h>
-#include <boost/iterator/transform_iterator.hpp> 
+#include <CGAL/boost/iterator/transform_iterator.hpp> 
 
 typedef CGAL::Simple_cartesian<double> Kernel;
 typedef Kernel::FT FT;

QP_solver/examples/QP_solver/solve_convex_hull_containment_lp2.h

@@ -1,6 +1,6 @@
 // example: function to check whether a point is in the convex 
 // hull of other points; this version uses a maker
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <CGAL/Kernel_traits.h>
 #include <CGAL/QP_models.h>
 #include <CGAL/QP_functions.h>

QP_solver/examples/QP_solver/solve_convex_hull_containment_lp3.h

@@ -1,6 +1,6 @@
 // example: function to check whether a point is in the convex 
 // hull of other points; this version uses a maker
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <CGAL/Kernel_traits.h>
 #include <CGAL/QP_options.h>
 #include <CGAL/QP_models.h>

QP_solver/include/CGAL/QP_models.h

@@ -38,7 +38,7 @@
 #include <sstream>
 
 #include <CGAL/boost/iterator/counting_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 // this file defines the following models:
 // - Quadratic_program_from_iterators

QP_solver/include/CGAL/QP_solution.h

@@ -39,7 +39,7 @@
 #include <boost/function.hpp>
 
 #include <CGAL/boost/iterator/counting_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 

QP_solver/include/CGAL/QP_solver/QP_solver.h

@@ -48,7 +48,7 @@
 #include <boost/bind.hpp>
 #include <boost/function.hpp>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <vector>
 #include <numeric>

QP_solver/test/QP_solver/master_mps_to_derivatives.cpp

@@ -25,7 +25,7 @@
 #include <string>
 #include <cstdlib>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/iterator/zip_iterator.hpp>
 #include <boost/shared_ptr.hpp>
 

Scripts/developer_scripts/detect_packages_licenses

@@ -15,8 +15,8 @@ for p in ^*build*(/); do
       elif [ -z "`echo $l | grep -Ev '^L?GPL \(v3 or later\) *$'`" ]; then 
           echo 'GPL (v3 or later)' > "$licFile"
       else
-          echo "MULTIPLE!" > "$licFile"; 
+#          echo "MULTIPLE!" > "$licFile";
-          echo $l >> "$licFile"
+          echo $l > "$licFile"
       fi
   fi
 done

Segment_Delaunay_graph_2/include/CGAL/Segment_Delaunay_graph_2.h

@@ -57,7 +57,7 @@
 #include <CGAL/assertions.h>
 
 #include <CGAL/boost/iterator/counting_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 /*
   Conventions:

Spatial_searching/include/CGAL/internal/K_neighbor_search.h

@@ -33,7 +33,7 @@
 #include <CGAL/Euclidean_distance.h>
 #include <CGAL/Splitters.h>
 #include <CGAL/internal/bounded_priority_queue.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 namespace internal{

Spatial_searching/test/Spatial_searching/Point_with_info.h

@@ -1,4 +1,4 @@
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 template <class Pt>
 struct My_point_with_info

Surface_mesh/include/CGAL/boost/graph/graph_traits_Surface_mesh.h

@@ -31,7 +31,7 @@
 #include <boost/graph/graph_traits.hpp>
 #include <boost/graph/properties.hpp>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <CGAL/boost/graph/properties_Surface_mesh.h>
 #include <CGAL/boost/graph/iterator.h>

Surface_mesh_segmentation/include/CGAL/internal/Surface_mesh_segmentation/Alpha_expansion_graph_cut.h

@@ -296,7 +296,7 @@ public:
         cut_time += timer.time();
         #endif
 
-        if(min_cut - flow < flow * tolerance) {
+        if(min_cut - flow <= flow * tolerance) {
           continue;
         }
         min_cut = flow;
@@ -534,7 +534,7 @@ public:
         cut_time += timer.time();
         #endif
 
-        if(min_cut - flow < flow * tolerance) {
+        if(min_cut - flow <= flow * tolerance) {
           continue;
         }
         min_cut = flow;
@@ -664,7 +664,7 @@ public:
         cut_time += timer.time();
         #endif
 
-        if(min_cut - flow < flow * tolerance) {
+        if(min_cut - flow <= flow * tolerance) {
           continue;
         }
 

Surface_mesh_skeletonization/benchmark/Surface_mesh_skeletonization/solver_benchmark.cpp

@@ -11,7 +11,7 @@
 #include <boost/property_map/property_map.hpp>
 #include <boost/graph/graph_traits.hpp>
 #include <boost/graph/adjacency_list.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <Eigen/SparseLU>
 #include <Eigen/Sparse>

Surface_mesh_skeletonization/include/CGAL/Mean_curvature_flow_skeletonization.h

@@ -39,7 +39,7 @@
 #include <boost/unordered_map.hpp>
 #include <boost/property_map/property_map.hpp>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/foreach.hpp>
 
 #include <CGAL/boost/graph/iterator.h>

Surface_mesher/include/CGAL/Surface_mesher/Intersection_data_structure_3.h

@@ -31,7 +31,7 @@
 #include <boost/format.hpp>
 #include <boost/bind.hpp>
 #include <boost/ref.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <CGAL/intersections.h>
 

Triangulation/include/CGAL/Delaunay_triangulation.h

@@ -29,7 +29,7 @@
 #include <CGAL/Dimension.h>
 #include <CGAL/Default.h>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #include <algorithm>
 

Triangulation/include/CGAL/Regular_triangulation_traits_adapter.h

@@ -25,7 +25,7 @@
 
 #include <CGAL/basic.h>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 

Triangulation/include/CGAL/Triangulation.h

@@ -36,7 +36,7 @@
 #include <CGAL/Random.h>
 
 #include <boost/iterator/filter_iterator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {
 

Triangulation_2/benchmark/Triangulation_2/CDT_with_intersection_2.cpp

@@ -6,7 +6,7 @@
 #include <CGAL/Random.h>
 #include <CGAL/Timer.h>
 
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <boost/random.hpp>
 
 

Triangulation_2/examples/Triangulation_2/info_insert_with_transform_iterator_2.cpp

@@ -1,7 +1,7 @@
 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
 #include <CGAL/Delaunay_triangulation_2.h>
 #include <CGAL/Triangulation_vertex_base_with_info_2.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <vector>
 
 typedef CGAL::Exact_predicates_inexact_constructions_kernel         K;

Triangulation_2/include/CGAL/Triangulation_2.h

@@ -50,7 +50,7 @@
 #include <boost/random/linear_congruential.hpp>
 #include <boost/random/uniform_smallint.hpp>
 #include <boost/random/variate_generator.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 #ifndef CGAL_NO_STRUCTURAL_FILTERING
 #include <CGAL/internal/Static_filters/tools.h>

Triangulation_3/examples/Triangulation_3/info_insert_with_transform_iterator.cpp

@@ -1,7 +1,7 @@
 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
 #include <CGAL/Delaunay_triangulation_3.h>
 #include <CGAL/Triangulation_vertex_base_with_info_3.h>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <vector>
 
 typedef CGAL::Exact_predicates_inexact_constructions_kernel         K;

Visibility_2/include/CGAL/Triangular_expansion_visibility_2.h

@@ -28,7 +28,7 @@
 
 #include <CGAL/Arrangement_2.h>
 #include <boost/shared_ptr.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 #include <CGAL/Constrained_Delaunay_triangulation_2.h>
 #include <CGAL/Arr_observer.h>
 #include <CGAL/assertions.h>

Voronoi_diagram_2/include/CGAL/Voronoi_diagram_2.h

@@ -51,7 +51,7 @@
 #include <CGAL/Identity_policy_2.h>
 
 #include <boost/variant.hpp>
-#include <boost/iterator/transform_iterator.hpp>
+#include <CGAL/boost/iterator/transform_iterator.hpp>
 
 namespace CGAL {