python - Generating surfaces of different heights from LIDAR points

I am seeking a way to generate surfaces (polygons in a 3D space) of different heights from LIDAR point clouds, whose spacing is 4 points per square meter. The typical application is like roofs, which may include holes, due to different heights.

However, most tools, which I can find, only generate the boundaries or footprints.

1. 
   
   

   ArcGIS: the examples are like: http://desktop.arcgis.com/en/arcmap/10.3/manage-data/las-dataset/lidar-solutions-data-area-delineation-from-lidar-points.htm and Creating boundary polygon shapefile from set of LAS files using ArcGIS Desktop?

   
   


2. 
   

   LAStools lasboundary.exe at http://www.cs.unc.edu/~isenburg/lastools/

   
   


3. 
   

   Python mayavi package: the example is like: https://stackoverflow.com/questions/33376657/from-point-cloud-to-surface-using-python

   
   

My initial thought is to do conversion from grid to polygons. Las file can be easy to convert to TIN, and then to grid, but grids are missing height information. Thus, the way is difficult to delineate detailed roof structures.

I also have a 30-centimeter and a 1-meter resolution orthophotos, if these do help, but no stereo photos.

The concrete expectation is like the following picture. The multi-color points are LIDAR point cloud, and the pink lines, which are 'digitized manually', indicate our expectation. However, it is easy to get the outside building boundary (the selected blue-green line), once buildings are classified, but no roofs with different heights (smaller inner polygons) are delineated.

Answer

First of all, if you want to work with heights, normalize the LiDAR point cloud.

You have a dense point cloud (4 pts/m²) and also high resolution aerial photos, hence, as you said getting the outside roof perimeter is not an issue. Therefore, assuming you'll manage to have a shapefile with the outer roof boundaries, use it to horizontally clip the point cloud and so, isolating building-only points. See here, for example.

Then, from each new .las file generated in previous step retrieve the vertical profile to identify how many surfaces there are in the scene. Tweak the class interval in the histogram, in order to identify all major surfaces. Get the height breaks to be used in the next step.

Clip the point clouds vertically in order to isolate points within same vertical surface (for example, with Fusion use ClipData with switches zmin and zmax). Then, generate individual Digital Surface Models (DSMs; examples 1 and 2) and convert the outputs to polygons using tools of type Raster to Polygon. The polygons will have the surface height in the TOC.

Try to use a programming environment to automate all the processing steps (example).