How do I add LiDAR data to JOSM?

Help and support
, , ,
1

I’m trying to add this LiDAR imagery that I found on a municipal website to JOSM, but I’m not exactly sure how to go about it.

2

Typical LIDAR Data is a 3D point cloud not just Imagery. I dont think you can add that to Josm.

You would need to postprocess the Lidar data to something useful, e.g. Building heights, road widths, distances depending on the observers position when the Lidar data was taken (From above or while driving on the roads).

2 Likes
3

Also, are you sure its even allowed to use this data for OSM? I did a quick read through Master Agreements - Products and Services Terms of Use and Changes | Esri Legal and i am not sure it is allowed.

After a quick search i found a nice article that implies it is Contributing GIS Data to OpenStreetMap and also explains steps that should be taken. Maybe worth a read? :slight_smile:

1 Like
4

It is licensed under “Open Government Licence - City of Abbotsford”:
https://opendata-abbotsford.hub.arcgis.com

https://osmfoundation.org/wiki/Licence/Licence_Compatibility#Open_Government_Licence_(OGL)_based_licences

3 Likes
5

That appears to be an index of lidar data, rather than the actual data. You will need to get the actual data and render it as an image, e.g. a hill shade. Perhaps the municipality has already done this. A hill shade can be very useful for precisely locating features such as streams.

1 Like
6

There are 547 las-files, approx. 2 TB in total:
https://opendata-abbotsford.hub.arcgis.com/datasets/1b0b439b14934da2b3625f0ae3c1f87f_0/explore

1 Like
7

You should be able to use qgis to convert las to DEM, and then create a hill shade from the DEM. Then use geoserver to create a wms which can be consumed by JOSM

1 Like
8

Could that all be done on-the-fly, or would I have to download the entire 2TB dataset first, then convert it in QGIS, then serve it to JOSM?

9

Not “on-the-fly” in the sense that when you pan to a new location in JOSM all of those steps are done automatically, but you should be able to manually go through the process a file at a time.

BTW, what types of features do you hope to map using this LiDAR data?

1 Like
10

Dang, that sounds rather tedious.

Trails mostly.

11

Yes, I think that assessment is correct.

Many trails will show up in elevation data derived from LiDAR, so this is a good idea. I often use it myself. Have you taken a look at the USGS 3D Elevation data? A hill shade of it is available in JOSM (main menu → imagery → imagery preferences → search “USGS” → select oppropriate entry → Activate). This might be derived from the same LiDAR data to which you are referring as the USGS often cooperates with local authorities to collect data like this. Even if that isn’t the case here, the USGS data might be good enough for your purpose.

1 Like
12

Unfortunately, this imagery was too coarse to be able to discern any detail fine enough for trail mapping the area that I’m interested in.

13

So you’re probably interested in the ground points. I have created the hillshade from level 8-19 (266078 tiles).
It is available for download here:

The tiles can be used in JOSM as TMS background image. After unpacking, in JOSM enter this URL at point 4 “Edit generated TMS URL (optional)” and adjust the path accordingly:
tms[8,19]:file:///path/tiles/{z}/{x}/{y}.png

2 Likes
14

This is amazing! Thank you so much!

How did you generate this? Did you use QGIS? I would love to know the process so that I can do it myself in the future!

15

Here is the script (slightly simplified, some steps may be unnecessary, but removing them could have unwanted side effects) that I created and run on a server with Ubuntu 22.04.
Be aware that some files no longer needed for processing will be removed immediately to save storage space, especially the large las files.

#!/bin/bash

# Please install the necessary packages first:
# apt install pdal gdal-bin parallel optipng aria2

# Base storage directory
BASE_DIR="/mnt/storage"

LAS_DIR="$BASE_DIR/las"
DEM_DIR="$BASE_DIR/dem"
TILES_DIR="$BASE_DIR/tiles"

mkdir -p "$LAS_DIR"
mkdir -p "$DEM_DIR"
mkdir -p "$TILES_DIR"

wget -qO- "https://hub.arcgis.com/api/v3/datasets/1b0b439b14934da2b3625f0ae3c1f87f_0/downloads/data?format=csv&spatialRefId=26910" | awk -F',' '/https/ {print $3}' > $BASE_DIR/las_list.txt

# Lines to process from las_list.txt (inclusive)
START_LINE=1
END_LINE=$(wc -l < $BASE_DIR/las_list.txt)

echo "start: $START_LINE"
echo "end: $END_LINE"

# Extract only the lines from START_LINE to END_LINE from las_list.txt and process them
while IFS= read -r url; do
    BASENAME_URL=$(basename "$url")
    LAS_FILE="$LAS_DIR/$BASENAME_URL"

    echo "Downloading LAS file $LAS_FILE with aria2c (4 segments)..."
    aria2c -x 4 -s 4 -j 1 -d "$LAS_DIR" -o "$BASENAME_URL" "$url"

    ls -lh "$LAS_FILE"

    if ! [ -f "$LAS_FILE" ]; then
      echo "File $LAS_FILE not found after download, skipping."
      continue
    fi

    BASENAME=$(basename "$LAS_FILE" .las)
    DEM_FILE="$DEM_DIR/${BASENAME}_dem.tif"
    
    RESOLUTION=0.25
    
    echo "Convert LAS $LAS_FILE to DEM $DEM_FILE"

    pdal pipeline --stdin <<EOF
{
  "pipeline": [
    {
      "type": "readers.las",
      "filename": "$LAS_FILE"
    },
    {
      "type": "filters.range",
      "limits": "Classification[2:2]"
    },
    {
      "type": "writers.gdal",
      "filename": "$DEM_FILE",
      "gdaldriver": "GTiff",
      "gdalopts": "compress=DEFLATE",
      "output_type": "min",
      "resolution": $RESOLUTION,
      "data_type": "float"
    }
  ]
}
EOF

    DEM_FILE_FILLED_TMP="$DEM_DIR/${BASENAME}_dem_filled_tmp.tif"
    DEM_FILE_FILLED="$DEM_DIR/${BASENAME}_dem_filled.tif"

    echo "Fill $DEM_FILE"
    gdal_fillnodata.py -md 100 "$DEM_FILE" "$DEM_FILE_FILLED_TMP"
    
    echo "Optimize $DEM_FILE_FILLED"
    gdal_translate -of GTiff -co COMPRESS=DEFLATE -co "PREDICTOR=3" -co "ZLEVEL=9" "$DEM_FILE_FILLED_TMP" "$DEM_FILE_FILLED"
    ls -lh "$DEM_FILE_FILLED"
    
    echo "Remove files"
    rm "$DEM_FILE_FILLED_TMP"
    rm "$DEM_FILE"
    rm "$LAS_FILE"
done < <(sed -n "${START_LINE},${END_LINE}p" "$BASE_DIR/las_list.txt")

echo "Generate $BASE_DIR/hillshade.tif"
gdalbuildvrt "$DEM_DIR/merged_dem.vrt" "$DEM_DIR"/*.tif
gdaldem hillshade "$DEM_DIR/merged_dem.vrt" "$BASE_DIR/hillshade.tif" -compute_edges -z 1.0 -s 1.0 -az 315 -alt 45

echo "Optimize $BASE_DIR/hillshade.tif"
gdal_translate -of GTiff -co BIGTIFF=YES -co TILED=YES -co COMPRESS=DEFLATE -co ZLEVEL=9 -co PREDICTOR=2 "$BASE_DIR/hillshade.tif" "$BASE_DIR/hillshade_compressed.tif"

rm "$BASE_DIR/hillshade.tif"

echo "Generate tiles"
gdal2tiles.py --processes=$(nproc) -x -z 8-19 -r bilinear --xyz "$BASE_DIR/hillshade_compressed.tif" "$TILES_DIR"

echo "Optimize tiles"
find "$TILES_DIR" -name "*.png" > "$BASE_DIR/png_list.txt"
cat "$BASE_DIR/png_list.txt" | parallel -j$(nproc) "optipng -o3 '{}' > /dev/null 2>&1"

echo "Delete aux files"
find $TILES_DIR -type f -name "*.png.aux.xml" -delete

process_file() {
    local FILE="$1"
    local channels
    channels=$(identify -format "%[channels]" "$FILE")
    
    if [[ "$channels" == *graya* || "$channels" == *srgba* ]]; then
        convert "$FILE" -background white -alpha remove -alpha off "$FILE"
    fi
}

export -f process_file

echo "Remove transparency"

find $TILES_DIR -type f -name "*.png" | parallel -j $(nproc) process_file {}

echo "Done."
3 Likes
16

@whb Great work! Thanks for doing this, and for providing your script. I am going to look for lidar data in my local area and attempt to create a hill shade from it using your approach.