r.in.pdal(1)

๐Ÿ“ฆ
version 8.4.2 ยท 4 versions
, in fedora 44

โœจ r.in.pdal โœจ

okay bestie, welcome to r.in.pdal, from the GRASS GIS User's Manual. this is the one that turns your LiDAR point cloud into a raster map. it's giving glow-up, it's giving aggregation, and yes, we're doing all of it together ๐Ÿ’–


๐Ÿท๏ธ KEYWORDS

raster, import, LIDAR, statistics, conversion, aggregation, binning


๐ŸŒธ SYNOPSIS

r.in.pdal
r.in.pdal --help
r.in.pdal [-wenodpg]  [input=name]   [output=name]   [file=name]
[method=string]   [type=string]   [base_raster=name]   [zrange=min,max]
[zscale=float]   [irange=min,max]   [iscale=float]   [drange=min,max]
[dscale=float]   [input_srs=string]   [pth=integer]   [trim=float]
[resolution=float]   [return_filter=string]
[class_filter=integer[,integer,...]]   [dimension=string]
[user_dimension=string]   [--overwrite]  [--help]  [--verbose]
[--quiet]  [--ui]

๐ŸŽ€ Flags

  • -w โ€” Reproject to project's coordinate system if needed. reprojects your input dataset to the coordinate system of the GRASS project. by default only datasets with a matching coordinate system can be imported, so this is your lil passport check ๐Ÿ›‚
  • -e โ€” Use the extent of the input for the raster extent. sets the computational region extents internally, based on the point cloud. no cap.
  • -n โ€” Set computation region to match the new raster map. sets the computation region to match the 2D extent and resolution of the newly created raster map ๐ŸŒธ
  • -o โ€” Override projection check (use current project's CRS). assumes the dataset has the same coordinate reference system as the current project. living dangerously, i respect it ๐Ÿ’…
  • -d โ€” Use base raster resolution instead of computational region. for getting values from the base raster, it uses that raster's actual resolution instead of the computational region resolution.
  • -p โ€” Print LAS file info and exit. just a lil peek, then it leaves ๐Ÿ‘€
  • -g โ€” Print data file extent in shell script style and then exit. script-friendly, iconic.
  • --overwrite โ€” Allow output files to overwrite existing files. (bestie. be careful.)
  • --help โ€” Print usage summary
  • --verbose โ€” Verbose module output
  • --quiet โ€” Quiet module output
  • --ui โ€” Force launching GUI dialog

๐Ÿฆ‹ Parameters

  • input=name โ€” LAS input file. LiDAR input files in LAS format (*.las or *.laz)
  • output=name โ€” Name for output raster map
  • file=name โ€” File containing names of LAS input files. LiDAR input files in LAS format (*.las or *.laz)
  • method=string โ€” Statistic to use for raster values
    • Options: n, min, max, range, sum, mean, stddev, variance, coeff_var, median, mode, percentile, skewness, trimmean, sidnmax, sidnmin, ev1, ev2, ev3
    • Default: mean
    • n: Number of points in cell
    • min: Minimum value of point values in cell
    • max: Maximum value of point values in cell
    • range: Range of point values in cell
    • sum: Sum of point values in cell
    • mean: Mean (average) value of point values in cell
    • stddev: Standard deviation of point values in cell
    • variance: Variance of point values in cell
    • coeff_var: Coefficient of variance of point values in cell
    • median: Median value of point values in cell
    • mode: Mode value of point values in cell
    • percentile: pth (nth) percentile of point values in cell
    • skewness: Skewness of point values in cell
    • trimmean: Trimmed mean of point values in cell
    • sidnmax: Maximum number of points in cell per source ID
    • sidnmin: Minimum number of points in cell per source ID
    • ev1: First eigenvalue of point x, y, z coordinates
    • ev2: Second eigenvalue of point x, y, z coordinates
    • ev3: Third eigenvalue of point x, y, z coordinates
  • type=string โ€” Type of raster map to be created. Storage type for resultant raster map
    • Options: CELL, FCELL, DCELL
    • Default: FCELL
    • CELL: Integer
    • FCELL: Single precision floating point
    • DCELL: Double precision floating point
  • base_raster=name โ€” Subtract raster values from the Z coordinates. the scale for Z is applied beforehand, the range filter for Z goes afterwards
  • zrange=min,max โ€” Filter range for Z data (min,max). applied after the base_raster transformation step
  • zscale=float โ€” Scale to apply to Z data. Default: 1.0
  • irange=min,max โ€” Filter range for intensity values (min,max)
  • iscale=float โ€” Scale to apply to intensity values. Default: 1.0
  • drange=min,max โ€” Filter range for output dimension values (min,max)
  • dscale=float โ€” Scale to apply to output dimension values. use if output dimension is not Z or intensity
  • input_srs=string โ€” Input dataset projection (WKT or EPSG, e.g. EPSG:4326). override input dataset coordinate system using EPSG code or WKT definition
  • pth=integer โ€” pth percentile of the values. Options: 1-100
  • trim=float โ€” Discard given percentage of the smallest and largest values. discard <trim> percent of the smallest and <trim> percent of the largest observations. Options: 0-50
  • resolution=float โ€” Output raster resolution
  • return_filter=string โ€” Only import points of selected return type. if not specified, all points are imported
    • Options: first, last, mid
  • class_filter=integer[,integer,...] โ€” Only import points of selected class(es). input is comma separated integers. if not specified, all points are imported.
  • dimension=string โ€” Dimension (variable) to use for raster values
    • Options: z, intensity, number, returns, direction, angle, class, source
    • Default: z
    • z: Z coordinate
    • intensity: Intensity
    • number: Return number
    • returns: Number of returns
    • direction: Scan direction
    • angle: Scan angle
    • class: Point class value
    • source: Source ID
  • user_dimension=string โ€” Custom dimension (variable) to use for raster values. PDAL dimension name

๐Ÿ’– DESCRIPTION

The r.in.pdal module loads PDAL library supported point clouds (with emphasis on LiDAR LAS files) into a new raster map using binning. and the way that you get to choose from a whole variety of statistical methods for the binning when creating that new raster? obsessed. โœจ

Since a new raster map is created during the binning, the binning of points depends on the current computational region settings (extent and resolution) by default (more on binning below, bestie). when you use the -e flag, the binning will be done in the extent of the point cloud, so the resulting raster will have extent based on the input point cloud. and when the resolution=value parameter is used, the binning is done using the provided resolution and the resulting raster will have that resolution (more below about extent and resolution management ๐ŸŽ€).

r.in.pdal is designed for processing massive point cloud datasets, for example raw LiDAR or sidescan sonar swath data. like, huge. no cap.

๐ŸŒˆ Binning

The main difference between r.in.pdal and v.in.pdal is that r.in.pdal creates a raster instead of just importing the points into GRASS GIS. but bestie, r.in.pdal does not merely rasterize the points from the point cloud. r.in.pdal uses binning to derive values for individual raster cells, so the value of a cell is typically an aggregation of values of individual points falling into one cell. in general, binning is the conversion of points into a regular grid. the binning of points with X and Y coordinates starts with the overlay of a grid of bins over the points. it's giving graph paper over your landscape ๐Ÿ“

In the basic case, binning is a method which counts the number of points which fall into one raster cell, i.e. bin. the number of points per cell (bin) indicates the density of points in the point cloud. the cell (bin) is always square or rectangular in case of r.in.pdal because the result is GRASS GIS 2D raster. the result of binning where the number of point per cell is counted is sometimes called 2D (two dimensional) histogram because a histogram is used in univariate statistics (in one dimension) to count the number samples falling into a given bin.

Figure: The binning on left was used to count number of points per (sometimes also called 2D histogram). The numbers in cells are examples of counts, the rest is represented by the color. The binning on right was used with mean to create a surface based on the values associated with the points. The numbers show examples of cell values. Note also the cells without any points which were assigned the NULL value.

The basic concept of binning is extended when the points have another value associated with them. for LiDAR data this value can be the Z coordinate or intensity. the value for a given cell (bin) is computed using univariate statistics from the values of all points in the cell. for example, computing the mean value of Z coordinates can yield a raster representing the digital elevation model. another example is the range of Z coordinates which can be used as a rough estimate of vegetation height ๐ŸŒฑ

๐Ÿ’… Statistics

Available statistics for populating the output raster map are:

  • n โ€” This computes the number (count) of points per cell. the result is a indicator of spatially variable density of points in the given area.
  • min โ€” This finds the minimum of point values in each cell. it can be useful when finding topography in a forested or urban environment and there is a lot of points per one cells (terrain is oversampled considering the desired resolution). it can also create surfaces independent on the noise from premature hits as it will always select the lowest point.
  • max โ€” This finds the maximum of point values in each cell. in connection with base_raster it can yield maximum vegetation of feature height per cell. for this purpose, it is usually much more appropriate than mean which would yield heights mostly influenced by the vertical distribution of points.
  • range โ€” This computes the range of point values in each cell. the range of Z coordinates per cell can be used as a rough estimate of vegetation height when the cells are small enough, slopes low and the area is mostly vegetated. however, for more profound analysis, the base_raster together with different statistics is recommended.
  • sum โ€” This computes the sum of point values per cell. this is useful especially when intensity is used as a value.
  • mean โ€” This is a mean (average) value of point values in cell. when used with Z coordinates (the default) and points from the ground class, the resulting raster is a digital elevation model. when intensity is used as a point value, the resulting raster contains mean intensity per cell. note that mean gives heights influenced by the vertical distribution of points.
  • stddev โ€” This computes the standard deviation of point values for each cell.

variance

Okay so this one gives you the variance of point values for each cell ๐Ÿ’œ. And heads up bestie: variance and its derivatives use the biased estimator (n), not the fancy unbiased one. It's calculated by the Welford algorithm โœจ.

coeff_var

This computes the coefficient of variance of point values for each cell ๐Ÿ’…. It comes out as a percentage, and it's defined as 100 * sqrt(variance) / mean.

median

This computes the median of point values for each cell ๐ŸŽ€.

mode

This computes the mode โ€” aka the most common value โ€” of point values for each cell ๐Ÿ‘‘. If a cell has more than one mode, an arbitrary one will be selected, so don't come for it if you wanted a specific one ๐Ÿ’€.

percentile

*p*th (*n*th) percentile of points in cell ๐ŸŒธ.

skewness

This is the skewness of point values in cell ๐Ÿฆ‹.

trimmean

This is the trimmed mean of point values in cell. Trimmed mean, also known as truncated mean, is a mean computed after discarding values at the low end and at the high end ๐Ÿ™ˆ. How many values to discard is given by the trim option, in percent. In statistics the usual percentage of trimmed values ranges from 5 to 25 percent โœจ.

sidnmax

Maximum number of points in cell per source ID ๐ŸŒˆ. Points are counted for each source separately, and the largest count is used as the final cell value. Corresponds to majority density (Dmax) as defined by Smeeckaert et al., 2013 ๐Ÿ“š.

sidnmin

Minimum number of points in cell per source ID ๐ŸŒธ. Points are counted for each source separately, and the smallest count is used as the final cell value. Corresponds to minority density (Dmin) as defined by Smeeckaert et al., 2013 ๐Ÿ“š.

ev1, ev2, ev3

The 1st, 2nd and 3rd eigenvalue of point x, y and z coordinates within a cell (Mallet et al. 2011) ๐Ÿ’ซ. Note that different statistics have different memory requirements (see below for details) ๐Ÿ‘€.


In varied terrain you may find that min maps make for a good noise filter, because most LIDAR noise is from premature hits ๐ŸŒฟ. The min map may also be useful to find the underlying topography in a forested or urban environment โ€” if the cells are oversampled ๐Ÿ’–.

Filtering and selection

Points are always filtered by spatial extent (based on the current computational region) unless the -e flag is set ๐ŸŒŠ. Afterwards, import value scaling is performed, followed by filtering by having value in range โœจ.

Points falling outside the current computational region will be skipped ๐Ÿšซ โ€” and yes, that includes points falling exactly on the southern or eastern region bound. To capture those, adjust the region with:

g.region s=s-0.000001

See g.region for details about computation region handling in GRASS GIS ๐Ÿซถ.

The zrange parameter may be used for filtering the input data by vertical extent ๐Ÿ’…. Example uses include filtering out extreme outliers and outliers on relatively flat terrain. This parameter can also be used for cutting the point cloud into vertical sections, preparing it for further processing by separate sections โ€” together as if it would be an imagery group (see i.group) โ€” or combined into a 3D raster using r.to.rast3 ๐ŸงŠ. In for these last examples, it might actually be more advantageous to use the r3.in.lidar module. The zrange parameter is especially powerful when used together with the base_raster parameter โ€” and the zrange is applied to Z values after the base_raster reduction โœจ.

Figure: This is the principle of zrange filter. Points with the Z coordinate value below the lower value in the range (here 180) are filtered out (blue points) and same applies for points above higher value in the range (here 250). All other points are preserved (green points).

Filtering can be performed also on intensity values with the irange parameter ๐ŸŒˆ. If the imported dimension is neither Z nor intensity, filtering can be performed with the drange parameter. Filtering by Z values and intensity can be performed independently on the imported dimension. If the imported dimension is not Z nor intensity, filtering can be performed with all three parameters simultaneously (zrange, irange and drange) or any of their combination ๐Ÿ’ซ.

A LiDAR pulse can have multiple returns ๐ŸŒŠ. The first return values can be used to obtain a digital surface model (DSM) where e.g. canopy cover is represented. The last return values can be used to obtain a digital terrain model (DTM) where e.g. the forest floor instead of canopy cover is represented ๐ŸŒฒ. The return_filter option allows selecting one of first, mid, or last returns. Return number and number of returns in the pulse associated with each point are compared to determine if the point is first, mid, or last return ๐Ÿ”.

LiDAR points often come as already classified into standardized classes ๐ŸŽ€. For example, class number 2 represents ground. For other classes, see LAS format specification in references. The class_filter option allows selecting one or more classes using numbers (integers) separated by comma โœจ.

The user can use a combination of r.in.pdal output maps to create custom raster-based filters โ€” for example, use r.mapcalc to create a mean-(2*stddev) map ๐Ÿ’…. (In this example the user may want to include a lower bound filter in r.mapcalc to remove highly variable points (small n) or run r.neighbors to smooth the stddev map before further use.)

Note that proper filtering of the input points is not only critical for the analysis itself, it can also speed up the processing ๐Ÿš€.

LiDAR points can contain a lot of attributes besides their coordinates ๐Ÿ“ฆ. By default r.in.pdal will use point Z coordinate as the output raster value source. The dimension option allows using different point attributes as an output instead of Z value. In case if other than Z value is chosen, output map type will be set to CELL. Explanation of various attributes can be found in LAS Specification chapter "Point Data Records" ๐Ÿ“–. Be ware โ€” not all attributes might be present, and also their value ranges and meaning might vary between LAS file producers ๐Ÿ‘€.

PDAL supports creation of custom dimensions (point attributes) in LAS files ๐ŸŒŸ. Use user_dimension option to provide name of any dimension for import. This option supersedes the one specified by the dimension option ๐Ÿ’ฅ.

Reduction to a base raster

For analysis of features on the terrain surface โ€” especially vegetation โ€” it is advantageous to remove the influence of the terrain on heights ๐ŸŒฑ, because the height above the terrain is what matters (e.g. height of a tree) rather than height of the top of the tree above the see level. In this case, the base raster would be a digital elevation model, which can be one derived from the point cloud, or obtained in some other way ๐Ÿ”๏ธ. LiDAR data often come with precomputed DEMs (quality should be checked in this case) and there is often a DEM available for a given area (fit with the point cloud, especially vertical, and resolution should be checked) โœ….

Figure: This is a profile of base raster (in orange) representing digital elevation model and selected points, e.g. first return, from point cloud (green dots). By default the points would create a digital surface model (thin brown line) but after reducing the Z coordinates using the base raster, the created surface is a derived from the height of points relative to the base raster. The usage of base raster is not limited to digital elevation model. The base raster can be any surface which has some relation to the point values, for example digital surface model representing top of the canopy.

Setting extent and resolution

Since the creation of raster maps depends on the computational region settings (extent and resolution), as default the current region extents and resolution are used for the import ๐Ÿ“. When using the -e flag along with the resolution=value parameter, the region used for the new raster will be based on the point cloud extent and the provided resolution โœจ. It is therefore recommended to first use the -s flag to get the extents of the LiDAR point cloud to be imported, then adjust the current region extent and resolution accordingly, and only then proceed with the actual import ๐Ÿ’…. Another option is to automatically set the region extents based on the LAS dataset itself (-e flag) along with the desired raster resolution. The best option is to know the point cloud extent ahead, e.g. from tiling scheme, and use it ๐Ÿ‘‘. See below for details.

Since the r.in.pdal generates a raster map through binning from the original LiDAR points, the target computational region extent and resolution have to be determined ๐Ÿงฎ. A typical workflow would involve the examination of the LAS data's associated documentation, or the scan of the LAS data file with r.in.pdal's -g flag to find the input data's bounds ๐Ÿ”.

Another option is to automatically set the region extents based on the LAS dataset extent (-e flag) along with the desired raster resolution using the resolution parameter ๐ŸŽฏ.

The -g flag prints extent of input data set(s) in a shell style suitable as command line parameters for g.region ๐Ÿš. Note, extent is determined from data set(s) metadata and without scanning whole dataset, and thus selectors and filters are not applied โš ๏ธ.

A simpler option is to automatically set the region extents based on the LAS dataset (-e flag) along with the target raster resolution using the resolution parameter โœจ. Also here it is recommended to verify and optimize the resulting region settings with g.region prior to importing the dataset ๐Ÿ’–.

Transformation

Point Z axis values can be altered by the scaling parameter set by zscale ๐ŸŒ€. Each point X axis value is multiplied by provided value before any filtering or further transformation takes place. In the same way act iscale and dscale parameters โ€” values are multiplied before range filters irange or drange are applied (if present) โœจ.

Output value scaling can be performed independently of filtering by providing the dscale parameter ๐Ÿ’ซ. This scaling parameter is applied only to output value, independently of its dimension. Thus if zrange and dscale are set for Z output dimension, filtering by zrange will be done with unscaled values but output will be scaled by dscale ๐ŸŽ€. dscale takes precedence over zscale and iscale if more than one parameter is provided ๐Ÿ‘‘.

๐Ÿ’ญ NOTES ๐Ÿ’ญ

๐ŸŒˆ Format and projection support ๐ŸŒˆ

okay so first things first, the usual file extensions for the LAS format are .las and .laz (that's the compressed one, bestie ๐Ÿ’…). but listen โ€” the compressed LAS (.laz) format can only be imported if libLAS was compiled with LASzip support. no LASzip, no .laz, periodt. ๐Ÿšซ

also, it's highly recommended that you compile libLAS with GDAL, because GDAL is what gets used to test whether the LAS coordinate reference system matches that of the GRASS project (which, yes, used to be called location โ€” grandma vibes, iykyk). โœจ


๐Ÿง  Memory consumption ๐Ÿง 

okay this part is important, so lock in ๐Ÿ”’

your input file can be arbitrarily large, sure, but r.in.pdal will use a large amount of system memory (RAM) for large raster regions (like > 10000x10000 pixels). the way that RAM just vanishes ๐Ÿ˜ญ

if the module refuses to start and starts complaining that there isn't enough memory, use the percent parameter to run the module in several passes. also, using a less precise map format โ€” CELL (integer) or FCELL (floating point) โ€” will use less memory than a DCELL (double precision floating point) output map. for methods=n, mode, sidnmin, sidnmax, the CELL format is used automatically. slay. ๐Ÿ’…

the mean and range methods will use an average amount of memory (comparing to other methods). methods such as n, min, max, and sum will use less memory, while stddev, variance, and coeff_var will use more. it's giving resource-hungry ngl ๐Ÿซ 

the memory usage for regular statistics mentioned above is based solely on region (raster) size. but heads up bestie โ€” the aggregate functions median, mode, percentile, skewness and trimmean will use more memory and may not be appropriate for use with arbitrarily large input files without a small value for the percent option, because unlike the other statistics, memory use for these also depends on the number of data points. so it's not just the size of your region, it's the size of your data too ๐Ÿ’€

the default map type=FCELL is intended as a compromise between preserving data precision and limiting system resource consumption. balanced queen ๐ŸŽ€


โœ‚๏ธ Trim option โœ‚๏ธ

the trim option value is used only when calculating trimmed mean values. attempt to use it with other statistical methods and you will get an error. no exceptions, no mercy ๐Ÿ™ˆ


the -p flag will use the PDAL library built-in LAS file metadata printing function. it will also include a list of supported dimension names. any of those listed dimensions can be passed to the user_dimension option. ๐Ÿ’–

โš ๏ธ important: dimensions are not checked for presence of data for any of the points in the file. so importing a dimension without actual data will result in an empty map. not the vibe you want, bestie ๐Ÿ˜ญ

๐ŸŽ€ EXAMPLES ๐ŸŽ€

simple example of binning of point from a LAS file into a newly created raster map in an existing project/mapset (using metric units):

# set the computational region automatically, resol. for binning is 5m
r.in.pdal -e -o input=points.las resolution=5 output=lidar_dem_mean
g.region raster=lidar_dem_mean -p
r.univar lidar_dem_mean

๐Ÿ” Finding suitable extent and resolution ๐Ÿ”

for the output raster map, a suitable resolution can be found by dividing the number of input points by the area covered (this requires an iterative approach as outlined here):

# print LAS metadata (Number of Points)
r.in.pdal -p input=points.las
# Point count: 1287775
# scan for LAS points cloud extent
r.in.pdal -g input=points.las output=dummy -o
# n=2193507.740000 s=2190053.450000 e=6070237.920000 w=6066629.860000 b=-3.600000 t=906.000000
# set computation region to this extent
g.region n=2193507.740000 s=2190053.450000 e=6070237.920000 w=6066629.860000 -p
# print resulting extent
g.region -p
#  rows:       3454
#  cols:       3608
# points_per_cell = n_points / (rows * cols)
# Here: 1287775 / (3454 * 3608) = 0.1033359 LiDAR points/raster cell
# As this is too low, we need to select a lower raster resolution
g.region res=5 -ap
#  rows:       692
#  cols:       723
#  Now: 1287775 / (692 * 723) = 2.573923 LiDAR points/raster cell
# import as mean
r.in.pdal input=points.las output=lidar_dem_mean method=mean -o
# import as max
r.in.pdal input=points.las output=lidar_dem_max method=max -o
# import as pโ€™th percentile of the values
r.in.pdal input=points.las output=lidar_dem_percentile_95 \
           method=percentile pth=95 -o

mean value DEM in perspective view, imported from LAS file ๐Ÿ”๏ธโœจ

further hints: how to calculate number of LiDAR points/square meter:

g.region -e
  # Metric project:
  # points_per_sq_m = n_points / (ns_extent * ew_extent)
  # Lat/Lon project:
  # points_per_sq_m = n_points / (ns_extent * ew_extent*cos(lat) * (1852*60)^2)

๐Ÿ Serpent Mound dataset ๐Ÿ

this example is analogous to the example used in the GRASS wiki page for importing LAS as raster DEM. (serpent mound?? iconic name fr ๐Ÿ๐Ÿ‘‘)

the sample LAS data are in the file "Serpent Mound Model LAS Data.laz", available at Serpent Mound Model LAS Data.laz:

# print LAS file info
r.in.pdal -p input="Serpent Mound Model LAS Data.laz"
# using v.in.lidar to create a new project
# create a project with CRS information of the LAS data
v.in.lidar -i input="Serpent Mound Model LAS Data.laz" project=Serpent_Mound
# quit and restart GRASS in the newly created project "Serpent_Mound"
# scan the extents of the LAS data
r.in.pdal -g input="Serpent Mound Model LAS Data.laz"
# set the region to the extents of the LAS data, align to resolution
g.region n=4323641.57 s=4320942.61 w=289020.90 e=290106.02 res=1 -ap
# import as raster DEM
r.in.pdal input="Serpent Mound Model LAS Data.laz" \
           output=Serpent_Mound_Model_LAS_Data method=mean

Figure: Elevation for the whole area of Serpent Mound dataset ๐ŸŒธ

โ›ฐ๏ธ Height above ground โ›ฐ๏ธ

the mean height above ground of the points can be computed for each raster cell (the ground elevation is given by the raster map elevation):

g.region raster=elevation -p
r.in.pdal input=points.las output=mean_height_above_ground base_raster=elevation method=mean

in this type of computation, it might be advantageous to change the resolution to match the precision of the points rather than deriving it from the base raster. little pro tip energy ๐Ÿ’…

๐Ÿ“ Multiple file input ๐Ÿ“

the file option requires a file that contains a list of file names with the full path. for example, a list of files in the directory /home/user/data:

points1.laz
points2.laz
points3.laz

would be lised in the file as:

/home/user/data/points1.laz
/home/user/data/points2.laz
/home/user/data/points3.laz

on Linux and OSX, this file can be automatically generated with the command:

ls /home/user/data/*.laz > /home/user/data/filelist.txt

on Windows:

dir /b c:\users\user\data\*.laz > c:\users\user\data\filelist.txt

the mean height above ground example above would then be:

g.region raster=elevation -p
r.in.pdal file=/home/user/data/filelist.txt output=mean_height_above_ground base_raster=elevation method=mean

in Python, the list of files can be created using the glob Python module:

import glob
import gscript
file_list_name = โ€™/home/user/data/filelist.txtโ€™
with open(, mode=โ€™wโ€™) as file_list:
    for path in glob.iglob(โ€™/home/user/data/lidar/*.lasโ€™):
        file_list.write(path + "\n")
gscript.run_command(โ€™r.in.pdalโ€™, file=file_list_name,
                    output=โ€™mean_height_above_groundโ€™,
                    base_raster=โ€™elevationโ€™ method=โ€™meanโ€™)

๐Ÿ˜ญ KNOWN ISSUES ๐Ÿ˜ญ

  • only one method can be applied for a single run and multiple map output from a single run (e.g. method=string[,string,...] output=name[,name,...] or n=string mean=string) is no supported. yeah. it said what it said ๐Ÿ’€

if you encounter any problems (or solutions!) please contact the GRASS Development Team. they wanna hear from you bestie ๐Ÿซถ


๐Ÿ‘€ SEE ALSO ๐Ÿ‘€

g.region, r.in.xyz, r.mapcalc, r.univar, v.in.lidar, r3.in.lidar, v.vect.stats

v.lidar.correction, v.lidar.edgedetection, v.lidar.growing, v.outlier, v.surf.bspline

Trimmed mean (Truncated mean, Wikipedia article), OpenTopography (LiDAR point cloud repository) ๐ŸŒ


๐Ÿ“š REFERENCES ๐Ÿ“š

  • V. Petras, A. Petrasova, J. Jeziorska, H. Mitasova (2016): Processing UAV and lidar point clouds in GRASS GIS. XXIII ISPRS Congress 2016 [ISPRS Archives, ResearchGate]
  • ASPRS LAS format
  • PDAL - Point Data Abstraction Library

๐Ÿ’– AUTHORS ๐Ÿ’–

Markus Metz

Vaclav Petras, NCSU GeoForAll Lab (base_raster option, documentation), Maris Nartiss, LU GZZF (refactoring, additional filters, custom dimension support)

based on r.in.xyz by Hamish Bowman and Volker Wichmann


๐Ÿ”— SOURCE CODE ๐Ÿ”—

available at: r.in.pdal source code (history) ๐Ÿ‘‘

accessed: Saturday Jan 17 05:40:55 2026

Main index | Raster index | Topics index | Keywords index | Graphical index | Full index


okay that's the whole page bestie, and it was giving GRASS GIS 8.4.2 Reference Manual the entire time โœจ๐Ÿฆ‹ the source lives at the r.in.pdal source code (with history, so you can snoop the lore), accessed Saturday Jan 17 05:40:55 2026, and the copyright situation is ยฉ 2003-2025 GRASS Development Team, GRASS GIS 8.4.2 Reference Manual, with the footer reading GRASS 8.4.2 and r.in.pdal(1) ๐Ÿ’… go import some point clouds, you got this ๐Ÿซถ๐ŸŽ€

๐Ÿ“ฆ where this page comes from