Read, Process, and Write Lidar Point Cloud Data

This example shows how to read a point cloud into the workspace, select a desired set of points, and then write the selected points to a point cloud file format.

Step 1: Read and Display Point Cloud

Read data from a .las file into the workspace by using the lasFileReader function. Display the properties stored in the output lasFileReader object.

fileName = fullfile(toolboxdir("lidar"),"lidardata","las","aerialLidarData.laz");
lasReader = lasFileReader(fileName)

lasReader = 
  lasFileReader with properties:

                 FileName: '/mathworks/devel/bat/filer/batfs2566-0/Bdoc24b.2725827/build/runnable/matlab/toolbox/lidar/lidardata/las/aerialLidarData.laz'
                    Count: 1018047
               LasVersion: '1.0'
                  XLimits: [4.2975e+05 4.3015e+05]
                  YLimits: [3.6798e+06 3.6801e+06]
                  ZLimits: [72.7900 125.8200]
            GPSTimeLimits: [3.3355e+05 sec    3.3443e+05 sec]
               NumReturns: 4
               NumClasses: 10
         SystemIdentifier: 'LAStools (c) by rapidlasso GmbH'
       GeneratingSoftware: 'TerraScan + OT'
         FileCreationDate: 28-Apr-2020
             FileSourceID: 0
                ProjectID: '0-0-0-00000000'
          PointDataFormat: 1
       ClassificationInfo: [6x3 table]
          LaserReturnInfo: [4x2 table]
    VariableLengthRecords: [3x3 table]

Read the point cloud from the .las file.

ptCloud = readPointCloud(lasReader);

Display the point cloud.

fig = figure(Position=[0 0 800 400]);
hPanel = uipanel(fig);
hPlot = axes(hPanel);
pcshow(ptCloud.Location,Parent=hPlot)

Step 2: Select Desired Set of Points

You can select a desired set of points in the input point cloud by specifying the classification value for the object classes and indices of the points within a region-of-interest (ROI).

Select points by specifying classification value

• To select the points by specifying the classification value, read information about the object classes in the input point cloud by using the ClassificationInfo property of the lasFileReader object.

disp(lasReader.ClassificationInfo)

    Classification Value        Class Name         Number of Points by Class
    ____________________    ___________________    _________________________

             1              "Unclassified"                  114842          
             2              "Ground"                        646632          
             4              "Medium Vegetation"             210101          
             6              "Building"                       45699          
             8              "Reserved(8)"                      751          
             9              "Water"                             22          

• Specify the classification value for the object class to read from the input point cloud by using the readpointCloud function. To read the points corresponding to medium vegetation region, set the value for Classification name-value argument to 4.

ptCloudB = readPointCloud(lasReader,Classification=4);

Display the point cloud.

fig1 = figure(Position=[0 0 800 400]);
hPanel1 = uipanel(fig1);
hPlot1 = axes(hPanel1);
pcshow(ptCloudB.Location,Parent=hPlot1)

Select points by specifying indices

Define a cuboid ROI within the range of the x, y and z coordinates of the input point cloud.

roi = [lasReader.XLimits(1)+200, lasReader.XLimits(2), ...
    lasReader.YLimits(1), lasReader.YLimits(2), lasReader.ZLimits(1), lasReader.ZLimits(2)];

Find the indices of the points that lie within the cuboid ROI.

indices = findPointsInROI(ptCloudB,roi);

Select the points that lie within the cuboid ROI and store as a point cloud object.

ptCloudC = select(ptCloudB,indices);

Display the point cloud.

fig2 = figure(Position=[0 0 800 400]);
hPanel2 = uipanel(fig2);
hPlot2 = axes(hPanel2);
pcshow(ptCloudC.Location,Parent=hPlot2)

Step 3: Write Selected Points to .las File format

Specify the name for the .las file and create a lasFileWriter object.

newfileName = "aerialvegetation.las";
lasWriter = lasFileWriter(newfileName);

Write the selected points to the .las file by using the writePointCloud function. The function creates the new file in the current working directory.

writePointCloud(lasWriter,ptCloudC);

Step 4: Check Properties of the Newly Written File

newlasReader = lasFileReader(newfileName)

newlasReader = 
  lasFileReader with properties:

                 FileName: '/tmp/Bdoc24b_2725827_1720804/tpd57d1183/lidar-ex04737654/aerialvegetation.las'
                    Count: 116598
               LasVersion: '1.2'
                  XLimits: [4.2995e+05 4.3015e+05]
                  YLimits: [3.6798e+06 3.6801e+06]
                  ZLimits: [84.9500 123.1100]
            GPSTimeLimits: [0 sec    0 sec]
               NumReturns: 1
               NumClasses: 1
         SystemIdentifier: 'MATLAB'
       GeneratingSoftware: 'LasFileWriter v1.0'
         FileCreationDate: 04-Sep-2024
             FileSourceID: 0
                ProjectID: '0-0-0-00000000'
          PointDataFormat: 3
       ClassificationInfo: [1x3 table]
          LaserReturnInfo: [1x2 table]
    VariableLengthRecords: [1x3 table]