Monday, February 26, 2018

Use PDAL to register a point cloud to control points

Sometimes, some 3D LiDAR point cloud data are not at the correct location or they are not aligned with other point clouds. When this happens, at least 3 control points from the point cloud and the reference need to be identified and a transformation matrix (usually a 4x4 homogeneous matrix) be calculated.

Once the matrix has been calculated, PDAL's transformation filter can then be used to perform the registration or transformation of the point cloud to the reference.

While PDAL can apply the transformation, it does not have a tool for calculating the matrix. This has to be done externally. This post uses a simple WebApp calculator at this site https://dominoc925-pages.appspot.com/webapp/calc_transf3d/default.html to calculate a rigid 4x4 transformation matrix given 2 corresponding sets of 3 or more control points.

Identifying 2 sets of corresponding control points

  1. Using a point cloud viewer such as FugroViewer, identify 3 or more source control points from the unregistered point cloud., e.g. p1, p2, p3 and p4.




  2. Using a GPS receiver or from a database, identify the corresponding points from the destination reference system.



Calculating the 4x4 transformation matrix
  1. Open up url https://dominoc925-pages.appspot.com/webapp/calc_transf3d/default.html in an Internet browser.
  2. In the Source 3D points field, type or paste in the source control point coordinates from the previous section.


  3. In the Destination 3D points field, type or paste in the corresponding control points from the destination reference system from the previous section.
  4. Click Calculate.

    The 4x4 transformation matrix and root mean square error is calculated.


    Note: if the RMSE is large, then there may be some error in the input coordinates.
  5. Using a text editor, create a PDAL pipeline JSON file which specifies the input, filter, and output. Copy and paste the 4x4 transformation matrix into the matrix field in the JSON file as shown below. Replace the commas with blanks. 
{
  "pipeline":[
    "input.laz",
    {
      "type":"filters.transformation",
      "matrix":"
0.21095277662308629 0.9774781497556848 -0.00594918357528593 23238.318780198064 
-0.9774332258657229 0.21100438075203337 0.01007175640974158 30100.055196603218 
0.011100225616467532 0.0036902647131527543 0.9999315811282344 32.53070599039722 
0 0 0 1"
    },
    {
      "type":"writers.las",
   "compression": "laszip",
      "filename":"output.laz"
    }
  ]
}

Running the transformation

  1. In Windows, open up the OSGeo4W command prompt.

    The Administrator: OSGeo4W Shell appears.


  2. In the Command Prompt, type in the following and press RETURN:

    C:> pdal pipeline transform_pipeline.json

    where transform_pipeline.json is the JSON file created in the previous section.

    The point cloud is transformed.

No comments: