Accessing 3-D point cloud data in IDL using the ENVI PointCloud API
13168 Rate this article:
No rating
Accessing 3-D point cloud data in IDL using the ENVI PointCloud API
Anonym Thursday, November 5, 2015
The following is an example of using the new PointCloud API included in ENVI 5.3. This allows easy access to 3-D data stored in LAS format. The data and metadata becomes accessible from IDL. The new API includes ability to read as well as write 3-D point cloud data in LAS format (there is also much more functionality included beyond that). There are 2 fundamentally different ways to get read access. One method is to simply read the points in the order they were stored in the LAS file. This would normally be the order the points were collected by the LiDAR sensor, although point cloud data could also be derived from other sources. The second way to access the point cloud data is by creating an ENVIPointCloud project, which means that the points will be spatially sorted such that they can be indexed quickly based on spatial tiles. The sorting process takes a little time initially, but once that is done, points can be retrieved quickly based on spatial querying methods.
This first example shows how to access point cloud data using the simple, non-sorted API calls.
ENVI> nv = envi(/current)
ENVI> pc=nv.QueryPointCloud('C:\Users\borsholm\Downloads\Lincoln.laz')
ENVI> print,pc
ENVIPOINTCLOUDQUERY <717772>
DATA_RANGE = 692204.97, 4519668.7, 340.79242, 695331.36, 4522642.4, 485.24673
NPOINTS = 9278073
SPATIALREF = !NULL
URI = 'C:\Users\borsholm\Downloads\Lincoln.laz'
ENVI> pc.metadata
ENVIPOINTCLOUDMETADATA <718218>
File Creation Day = 0
File Creation Year = 0
File Source ID = 0
Generating Software = 'Lidar Explorer by ProLogic, Inc.'
Global Encoding = 0
Max X = 695331.36
Max Y = 4522642.4
Max Z = 485.24673
Min X = 692204.97
Min Y = 4519668.7
Min Z = 340.79242
Number Of Point Records = 9278073
Number Of Variable Length = 3
Point Data Format = 0
Point Data Record Length = 20
Project ID GUID Data 1 = 0
Project ID GUID Data 2 = 0
Project ID GUID Data 3 = 0
Project ID GUID Data 4 = 0, 0, 0, 0, 0, 0, ...
System Identifier = ''
Version Major = 1
Version Minor = 0
X Scale Factor = 1.4558386e-006
X Offset = 692204.97
Y Scale Factor = 1.3847266e-006
Y Offset = 4519668.7
Z Scale Factor = 6.7266782e-008
Z Offset = 340.79242
ENVI> pts=pc.GetPointsInRange(0,pc.npoints,intensity=i)
ENVI> help,pts,i
PTS DOUBLE = Array[3, 9278073]
I UINT = Array[9278073]
ENVI> pc.Close
This shows how get access to metadata associated with the file, and to return the points as a 3xN array, in the original order that they were stored in the file.
This second example shows how to spatially sort the points in a new project, and access points in a given rectangle.
ENVI> pc = nv.OpenPointCloud('C:\Users\borsholm\Downloads\Lincoln.laz')
ENVI> pts = pc.GetPointsInRect(694000, 4520000, 695000, 4521000)
ENVI> help,pts
PTS DOUBLE = Array[3, 1000002]
ENVI> pc.Close
Automated Cloud and Cloud Shadow Detection Change Detection using Airbus WorldDEM™ and Elevation from Optical Imagery