readme.md

** GENERAL **

All the test data are placed in the folder examples/data. The name of the folder indicates the type of data or its complexity. The complexity 1 is low while the complexity 6 is high.

Examples are the entry point to the package:

• kinetic_2d_example: 2D kinetic algorithm on a random set of segments
• kinetic_precomputed_shapes_example: 3D kinetic algorithm on a set of user-defined polygons read from a file
• kinetic_random_shapes_example: 3D kinetic algorithm on a set of random polygons, you can provide the number of input polygons to be generated and the number of vertices each polygon should have
• kinetic_reconstruction_example: given a ply with classified point cloud, it first fits polygons using Region Growing, then runs 3D kinetic, and finally filters out exterior volumes creating a reconstructed model

Tests:

• kinetic_2d_stress_test: tests the 2D kinetic algorithm
• kinetic_3d_test_all: tests the 3D kinetic algorithm on all data from the examples folder

The 3D algorithm is running on macOS with zero warnings. On Windows and Linux, it compiles and works, but can generate warnings related to conversions between std:size_t and int e.g. or similar warnings. The 2D algorithm should work on all platforms as well.

** CURRENT ISSUES **

-- EPICK versus EPECK --

By running: ./kinetic_precomputed_shapes_example data/stress-test-4/test-9-rnd-polygons-12-4.off

first with EPICK and then with EPECK shows a huge difference in runtime. This data set contains 12 random polygons, each having 4 vertices. The time for EPICK is milliseconds while for EPECK is about 3 minutes. It can also be noticed that most of the time is spent at the last iterations while the first iterations are very fast. It is even slower for Simple_cartesian<Gmpq>.

It is probably happening because at the latest iterations the computation involves all operations carried out from the first iteration. Evaluating these cascading operations is a slow process. E.g. computing an intersection between two lines at the first iterations takes seconds while at the higher iterations it takes minutes.

Another issue shows up when running the code with the following data set: ./kinetic_precomputed_shapes_example data/real-data-test/test-40-polygons.ply 6 that is with k = 6, where k is the number of allowed intersections between support planes of the input polygons. Here, at about 5300 iteration, assertions start failing but not because some parts of the algorithm are not implemented (as indicated in the TODO) but because the scaling is no longer uniform due to the accumulating errors. It leads to wrong configurations and hence to the wrong events, which in practice with the correct precision cannot happen.

-- Possible solutions --

• Use EPICK and compute always with respect to the input. E.g. intersections between lines should be done not between lines at the current and previous iterations but between lines at the current and first iterations. The reason for this optimization is that if we use simply EPICK, when the number of input polygons grow, we bump into an issue of the accumulating error that gets bigger and bigger with each iteration and at some point can break the results. It does not happen with EPECK but we lose speed there.

• Use EPECK only when it is absolutely necessary like when computing intersections and use EPICK for all other computations. This way we avoid accumulating errors and should preserve speed.

A few ideas shortly after the meeting with Sebastien:

• speed up kinetic with EPECK by optimizing the parts, which are slow, vtune it
• use different speed for each polygon edge during the uniform scaling, a similar trick to what is used in the straight skeleton now
• point_2() wrapper inside Support_plane.h: it can compute point not based on the point constructed at the previous event (which is cascading) but based on the interesting polygon edge and intersection graph edge
• avoid certain events in the queue or keep track only of the last 10 events instead of all of them

** INTERNALS **

-- File descriptions --

• KSR sub-folder, used in both KSR2 and KSR3:

  • conversions.h: the kinetic_traits class that performs intersections, all intersections in the code are called from this class, here the hybrid mode could be potentially implemented.

  • debug.h: a file with a bunch of functions which enable to export / dump to a file different intermediate steps, events, data structure, etc.

  • enum.h: enumerations used by the reconstruction code from the Reconstruction.h.

  • parameters.h: internal parameters used in the code.

  • property_map.h: property maps used by the reconstruction code from the Reconstruction.h.

  • utils.h: different internal utilities such as normals, distances, angles, etc. the most important one is tolerance(): this is a global tolerance used everywhere in the code

• KSR_2 sub-folder, 2D kinetic algorithm, works and tested, fully implemented by Simon.

• KSR_3 sub-folder, 3D kinetic algorithm + 3D reconstruction algorithm:

  • Data_structure.h: a class with all conversions and all basic operations, which are performed on the polygons and during events, it also stores the main intersection graph and all support planes.

  • Event_queue.h: a wrapper around boost multi_index_container that represents a queue of all events, which can be sorted by time or vertex type. It is based on doubles so even when we use EPECK globally, we need to convert to double here because boost multi_index_container fails to work with EPECK due to the fact that values represented by EPECK can vary when converted to_double() or to_interval().

  • Event.h: represents different event types used in the propagation.

  • Intersection_graph.h: a boost adjacency_list graph that stores an initial intersection graph fully in 3D that is used as constraints for kinetic propagation.

  • Support_plane.h: a wrapper around a support plane for each input polygon that stores a 2D surface mesh that represents an input polygon. So, initially this mesh has only 1 face that is input polygon while during propagation it is updated and new vertices, faces are added. We preserve the uniform scaling of the input polygon and validity of this mesh after each event. If any of these two fails, the whole algorithm fails. At the end of the support plane header, there is an overload of the operator()== that is used to compare if two planes are equal. It heavily depends on the tolerance parameters and in case two planes are found to be equal, they are merged at the initialization step.

  • Initializer.h: a class that gets input polygons, creates an enlarged bounding box, removes all near-collinear vertices of each polygon, merges all near-coplanar support planes, then inserts each support plane and intersects it with the bounding box. At the final step, it calls Polygon_splitter class in order to intersect all polygons within the bounding box, one support plane at a time.

  • Polygon_splitter.h: given a support plane, its original polygon, and a set of intersection edges created by intersecting all support planes with the bounding box (see Initializer.h), it builds a CDT, inserts all these edges in the CDT, marks all interior and exterior faces and creates a proper 2D surface mesh. This class can be parameterized by any kernel independently of the input kernel.

  • Propagation.h: is called after Initializer.h has done its work, this class creates an initial queue of events for all vertices of all polygons and handles these events one by one. This is the most time-consuming part of the total algorithm. The queue is updated and rebuilt until no valid events show up. This class can be parameterized by any kernel independently of the input kernel.

  • Finalizer.h: is called after Propagation.h, this class receives a set of 2D surface meshes from the propagation and an initial intersection graph and does three things: first it checks the validity of the results, it then searches for dangling / hanging faces if any and fills the gaps; and finally creates 3D volumes bounded by all faces of the 2D surface meshes. This class can be parameterized by any kernel independently of the input kernel.

  • Reconstruction.h: this class first gets a point cloud with or without semantic labels. If the labels are missing, some parts of the code have to be still finished. However, if labels are present, it assumes they come from urban areas that is they are walls, roofs, ground, and trees. It removes all tree points, fits a polygon to the ground points, runs region growing separately for wall and roof points and detects convex planar polygons, which approximate the input point cloud. Next, it runs 3D kinetic on these input polygons + ground polygon and gets a set of volumes, which partition the 3D space. It then estimates a probability of each volume to be inside or outside the point cloud using Visibility.h and runs the graph cut algorithm on these estimations using Graphcut.h. At the end, it removes all volumes, which have been labeled as exterior volumes and returns a reconstructed model: 3D surface mesh that approximates the input point cloud (or buildings in case labels are urban related).

  • Visibility.h: estimates a probability of each partition volume to be inside or outside a point cloud by using normals of the input points and sampling each volume.

  • Graphcut.h: takes initially estimated probabilities for each volume to be interior or exterior with respect to the point cloud, computed by Visibility.h, creates a graph where each node is a volume and each edge connects volumes to their neighboring volumes, and runs the mincut - maxflow Boykov optimization algorithm to define which volumes are inside and outside the point cloud. All edges are weighted by the face area that adjacent to two incident volumes and all nodes are weighted by the volume itself.

• Kinetic_shape_reconstruction_2.h: is an entry point to the 2D kinetic propagation.

• Kinetic_shape_reconstruction_3.h: is an entry point to the 3D kinetic propagation and 3D kinetic reconstruction algorithms.

-- Epsilon usage inside the code --

Note that epsilon tolerance is used throughout the code. It is defined inside utils.h in the function tolerance(). It is used everywhere where we expect a result of certain precision but it may not be the case. We check if we outside this tolerance and apply different sub-algorithms in order to be sure we that will generate the right results. This is mostly used for EPICK. When using EPECK, the precision is higher and the tolerance should be satisfied until there is a bug.

-- Important parts of the code --

• point_2() wrapper from the Support_plane.h. Currently all original points are 3D, this wrapper takes a point, and constructs the corresponding 2D point for the chosen 2D support plane and at the chosen time step. That means all these points are constructed that is a weak point when using EPECK because it is slow.

• operator==() at the bottom of the Support_plane.h: controls if two near-coplanar planes should be merged. If we merge too many planes because our input parameters are very weak, we fail to create a valid partition. If we do not merge at all, the near-coplanar support planes may lead to intricated errors in the kinetic propagation due to numerical instabilities.

• parameters: all all explained in the parameters.h. The parameters used in the Reconstruction.h are defined in the file examples/include/Parameters.h.

• FUTURE POINTS AND DIRECTIONS section at the bottom of the Data_structure.h, this is where the future points are computed and this is a part of the code that leads to multiple precision issues, identifying a future point from the previous event is hard, so instead we simply translate the lines and intersect them at the next time step to get the point at which our current vertex will be later, but these intersections are imprecise. If we lose precision here, we fail to scale polygon uniformly so our point can end up behind the bounding box or in the direction opposite to the one we need, especially if the lines that we intersect are near parallel.