
Data was accessed from the USGS The National Map website. To clarify data availability, pre-fire data from 2018 included a point cloud and a digital surface model for the extent of the Eaton fire. Post-fire, the point cloud itself was not available.
The basic methodology was to access the LAZ files, mosaic them, generate a DSM, create a hillshade, and compare pre- and post-fire data. Generating a DSM from the LAZ files for the pre-fire data was not necessary, as a DSM is available from USGS The National Map, but we include this process for reference.


We also uploaded these images into the King Lab Juxtapose software to allow for a more dynamic comparison
The R code is provided below for your reference.
# Call needed packages & Import files -------------------------------------
# packages
library(lidR)
library(raster)
library(terra)
# 2025 DSM file path
file_path_2025 <- "/Users/nathanroe/Downloads/post_fire_dsm_11SLT003945037845.tif"
# Create terra SpatRaster for 2025 DSM
post <- terra::rast(file_path_2025)
# 2018 pre-fire LAZ file paths
file_paths_laz <- c("/Users/nathanroe/Downloads/USGS_LPC_CA_LosAngeles_2016_L4_6519_1894a_LAS_2018.laz",
"/Users/nathanroe/Downloads/USGS_LPC_CA_LosAngeles_2016_L4_6519_1889b_LAS_2018.laz",
"/Users/nathanroe/Downloads/USGS_LPC_CA_LosAngeles_2016_L4_6514_1889c_LAS_2018.laz",
"/Users/nathanroe/Downloads/USGS_LPC_CA_LosAngeles_2016_L4_6514_1894d_LAS_2018.laz")
# Read LAZ files into list
las_files <- lapply(file_paths_laz, readLAS)
# Generate DSM from pre-fire LAZ ------------------------------------------
# Merge the LAS files using rbind
las_combined <- las_files[[1]]
for(i in 2:length(las_files)) {
las_combined <- rbind(las_combined, las_files[[i]])
}
# Generate DSM using non-ground points
pre <- rasterize_canopy(las_combined, res = 1, algorithm = dsmtin())
# Create hillshades -------------------------------------------------------
# Pre-fire hillshade generation
pre_slope <- terrain(pre, "slope", unit="radians")
pre_aspect <- terrain(pre, "aspect", unit="radians")
pre_shade <- shade(pre_slope, pre_aspect, 40, 270)
# Post-fire hillshade generation
post_slope <- terrain(post, "slope", unit="radians")
post_aspect <- terrain(post, "aspect", unit="radians")
post_shade <- shade(post_slope, post_aspect, 40, 270)
# Ensure the same coordinate reference system is used
pre_shade <- project(pre_shade, terra::crs(post_shade))
# Crop pre-fire hillshade to post-fire hillshade
pre_shade <- crop(pre_shade, ext(post_shade))
# Plot! -------------------------------------------------------------------
plot(pre_shade, main="Eaton Pre-Fire", col=gray.colors(100))
plot(post_shade,main="Eaton Post-Fire", col=gray.colors(100))
Thanks for reading to the end!
















