LiDAR Processing
Open Topography Work Flow
Log into the OpenTopography myOpenTopo option
Select LiDAR Point Cloud & Processing option in the Data tab.
Select one of the CZO LiDAR sites, e.g. Jemez River Basin, Boulder Creek, Southern Sierra.
1a. Select a Region from which to calculate surfaces and extract data.
1b. Choose a return classification for the surface model generation.
- For Bare Earth surface modelling use only the 'Ground' classes checked on.
- For Canopy surface modelling use both the 'Ground' and 'Unclassified' classes.
- Open topography only allows 200 million points to be selected at a single time. If using only ground points you can select a larger area than if using the full point cloud.
2. Select the type of point cloud data you want to download
- Pick either LAS (uncompressed) or LAZ (compressed files)
3a. DEM generation (Local Gridding)
- OpenTopography has two options for generating DEM surfaces, the first option is local gridding which takes a statistic of the point cloud at a defined length scale square (e.g. 1m^2, 10m^2, etc.).
- For bare earth DEM gneation using the Local Gridding Method select the 'Calculate Zmean grid' and the 'Ground' only return classifications from 1b.
- For Canopy surface model generation select the 'Calculate Zmax grid'.
- For Gridding Parameters, set the Grid Resolution to the desired size.
- For this analysis I generated DEMs at 1m, 3m, 5m, 10m, 15m, and 30m resolution.
- For the Canopy surface models I set the resolution at 0.333m
- The Grid Format can be whatever you prefer - Arc ASCII Grid for ArcGIS users, or GeoTiff for everyone else.
- Set the 'Null Filling' option to 3x3 moving window.
3b. DEM Generation (TIN)
- OpenTopography's second option is using a TIN (triangular irregular network).
- Click on the 'Gridding Method' to 'Calculate TIN'
- Set the gridding parameters to the same option as for 3a.
- Leave the Max triangel size as 50 units (default)
- Select the Grid Format you like
4. Derivative Products
- If you are interested you can generate derivative products which include hillshades and slope
5. Visualization
- Optional, generate Google Earth files
6. TauDEM - Hydraulic Terrain Analysis Products
- OpenTopography has the option to generate TauDEM products directly.
- This is a great feature that can save time later on if you are interested in hydraulically corrected DEMs, and calculating flow direction, catchment area, and Topographic Wetness Index (TWI).
After you have selected your features you can execute the run and watch it go.
When OpenTopo is complete you will get an Email with the LiDAR Job Report with options to download the Job Results. There will also be a meta-data file that contains the full job description:
Dataset Information: Dataset Name: Jemez River Basin Snow-off LiDAR Survey(JRB_10_Jul) Dataset Acknowledgement: LiDAR data acquisition and processing were completed by the National Center for Airborne Laser Mapping (NCALM), funded by the National Science Foundation (NSF) Award EAR-0922307, Valles Caldera National Preserve (VCNP), Bandelier National Monument/National Park Service (BNM/NPS) and United States Geological Survey (USGS), and coordinated by Qinghua Guo for the Jemez River Basin and Santa Catalina Mountains Critical Zone Observatory funded by the National Science Foundation Award EAR-0724958. Data Access Acknowledgement: This material is based on [data, processing] services provided by the OpenTopography Facility with support from the National Science Foundation under NSF Award Numbers 0930731 & 0930643 Horizontal Coordinates: UTM Zone 13N NAD83 (CORS96) [EPSG: 26913] Vertical Coordinates: NAVD88 (GEOID03) [EPSG: 5703] Job Description: User:tswetnam@email.arizona.edu Job ID: 1434668225132265012811 Title: Redondo 0.333m Description:Redondo canopy surface 0.333m Job Processing Result: Submission Time: 2015-06-18 15:57:05 Completion Time: 2015-06-18 16:03:27 Duration: 382 seconds Number of Points: 52,812,660 Status: Done Data Selection Coordinates: Xmin: 358646.05 Xmax: 360472.918 Ymin: 3969729.673 Ymax: 3972072.639 Classifications: all Point Cloud Results: Download point cloud data in ASCII format: No Download point cloud data in LAS format: Yes Download point cloud data in LAZ format: No DEM Generation (Local Gridding): Calculate Zmin grid: No Calculate Zmax grid: Yes Calculate Zmean grid: No Calculate Zidw grid: No Calculate point count grid: No Grid Format: GTiff Resolution: 0.3333 meter Radius: 0.3333 meter Null Fill: 3 DEM Generation (TIN): Calculate TIN: Yes Grid Format: GTiff Resolution: 0.3333 Max Triangle Size: 50 Derivative Products: Generate Hillshades & Slopes: none Visualization: Images & Google Earth KML: No
LiDAR Data Processing
Download USFS FUSION
On a Windows OSmachine, open a command prompt in the search window type 'CMD'.
Example:
Microsoft Windows [Version 10.0.10586] (c) 2015 Microsoft Corporation. All rights reserved. C:\Users\tswetnam>e: E:\>cd 4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns>dir /b /s *.las>filelist.txt E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns>
Navigate to the root folder where the *.las or *.laz files are stored and type the following:
dir /b /s *.las>filelist.txt
This will create a file list text file in the root folder which lists the directory location and file name of all *.las files in the directory:
E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns\34110A7101.las
E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns\34110A7102.las
E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns\34110A7103.las
...
E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns\34111G3425.las- Next, prepare two batch files for processing the *.las of the entire directory.
- Save these files to the same folder where FUSION is installed.
- The first batch file, which you can name 'call_create_bare_surface.bat' calls the GridSurfaceCreate process with a *.BAT executable and the file list – this script ‘calls’ each file independently from the filelist.txt in a sequential order
Type all command line code on a single line, without carriage return:
for /F "eol=; tokens=1* delims=,. " %%i in (E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns\filelist.txt) do call create_bare_surface %%
Create an output folder in your desired location, i.e.: E:\4FRI_Phase_I\bare_earth\
Next write a second *.bat file called 'create_bare_surface.bat', also saved to the FUSION directory (You can put these files in another directory if you choose - but make sure to include their drive location, e.g. E:\folder_name\
GridSurfaceCreate /class:2 /median:3 /smooth:3 E:\4FRI_Phase_I\bare_earth\%~n1.dtm 1 m m 1 12 2 2 %1.las
Note: see the FUSION Manual for initializing the resolution and projection settings. Here I am creating a bare earth surface model using only 'class:2' which are the classified ground returns, a 3x3 median filter with a 3x3 smooth. The projection and resolution settings are: '1 m m 1 12 2 2' for: a 1 meter output surface in horizontal meters and vertical meters units in UTM zone 12 NAD83 NAVD88 projection. The %1.las calls the *.las file in that iteration of the file list.
Now that you have a bare earth model you can create a Canopy Height Model.
Create another 'call_create_canopy_height_model.bat' script which will use the same *.las file list as before:
for /F "eol=; tokens=1* delims=,. " %%i in (E:\4FRI_Phase_I\4FRI_LiDAR_as_delivered_2014_01_06\Points\All_Returns\filelist.txt) do call create_canopy_height_model %%
Create an output folder in your desired location, i.e.: E:\4FRI_Phase_I\canopy_model\
Next write a second *.bat file called 'create_bare_surface.bat', also saved to the FUSION directory
CanopyModel /ground:E:\4FRI_Phase_I\bare_earth\%~n1.dtm /outlier:-1,55 E:\4FRI_Phase_I\canopy_model\%~n1.dtm 0.333 m m 1 12 2 2 %1.las
Both Canopy Height Models and Bare Earth Models are in a *.DTM file format. These data are converted to ASCII file format with the DTM2ASCII FUSION command:
Generated a filelist.txt file using Windows CMD prompt, after navigating to the file folder where the *.DTM data are stored (See 1a above), change file type to *.dtm from *.las.
Called Process: dtm2ascii %1.dtm
Batch Process Command: for /F "eol=; tokens=1* delims=,. " %%i in (I:\BC2010_snoff\DTM\filelist.txt) do call d2a_bat %%i
A GridMetrics routine was run at 10 m and 30 m cell size.
The 2010 Snow-off Boulder Creek dataset contains 1,369 file tiles.
There are 40 different statistics for the GridMetrics.
ArcGIS Sub-Routines
GridMetrics commands were executed at 10 meters and 30 meters squared area. A minimum height of 2 meters above surface level was used as the minimum for estimating the canopy metrics.
These tiles were mosaicked in ArcGIS ‘Mosaic to New Raster’ tool with the ‘max’ setting for overlapping cells to avoid NoData pixel effects.
Mean Height > 2 m
Max Height > 2 m
Standard Deviation Height > 2 m
Kurtosis of Height > 2 m
Skew of Height > 2 m
Variance of Height > 2 m
All of the 30 m outputs were shifted with the ‘Shift’ tool in ArcGIS by -15 m (y-axis) and -15 m (x-axis) to align properly over the Canopy Height Model.
All of the 10m outputs were shifted by -5 m (y-axis) and -5 m (x-axis) to align properly over the Canopy Height Model.
The shifted layers were saved into the Grid_Metrics directory folder as GRID format.
MATLAB Sub-Routines
The individual tree list was generated by importing the ASCII CHM outputs into MATLAB and execute the Variable-area Local Maxima algorithm:
Execute the asc_import.m file in MATLAB:
tile=asc_import(I:\BC2010_snoff\CHM\ ot_450500_4433000.asc');
Run the VLM algorithm: tile_v=vlmt1(tile);
Project the stems with the tile UTM coordinates: tile_exp=export_utm(tile_v,tile);
Write the output to an ASCII *.CSV file: dlmwrite('I:\BC2010_snoff\STEMMAPS\36112E3301.csv',tile_exp,'delimiter',',','precision',10);
Open an ArcGIS ArcMAP instance and Add Data:
Add the Canopy Height Model *.ASC and project in NAD83 UTM Zone 13
Add the *.CSV file and right click to Display XY data.
You can open the *.CSV file in TextPad and add headers before you import into ArcMap if you choose, the five fields that are saved as columns are: OID, UTME, UTMN, HT, CR
Right Click on the *.CSV Events in the Table of Contents and select Properties
Go to ‘Symbology’ Tab
Select Quantities
Proportional Symbols
Fields: Value: Field 5 (or CR)
Set Units to Meters
Data represents Radius (1/2 width)
File Directory Explanations
The Catchment_GIS Directory is the Root folder for all of the data
The directory has five folders:
3 catchments
A Grid Metrics Directory of the entire Boulder Creek LiDAR flight
A Shape_files Directory with the 3 catchment polygons
The 3 catchments are identified as:
Betasso
Gordon_Gulch
Como_Creek
Within each catchment folder there are three sub folders:
Bare_Earth – 1 m surface model (minimum height hit return) for each of the tiles that over-lap the shape file of the catchment
Canopy_Height – 0.333 m surface model (high hit return) for each of the tiles that over-lap the shape file of the catcment
Stem_Map – the exported MATLAB output for local maxima over the Canopy_Height Models
The Grid Metrics folder contains both the 10 m and 30 m rasters that were mosaicked in ArcGIS:
Mean Height > 2 m
Max Height > 2 m
Standard Deviation Height > 2 m
Kurtosis of Height > 2 m
Skew of Height > 2 m
Variance of Height > 2 m