From 3f7d4c5afffb5975d9aa92e27511a01e56058934 Mon Sep 17 00:00:00 2001 From: steven_f Date: Wed, 26 Jun 2024 11:49:32 -0600 Subject: [PATCH 1/3] Add use of JSON configuration file instead of changing scripts for each run. --- README.md | 32 +++++---- scripts/01_PrepareDataForFusion.py | 72 +++++-------------- .../01_PrepareDataForFusion_MultiProjects.py | 54 ++++++-------- scripts/02_CreateAPSettingsPRP.R | 50 +++++-------- scripts/03_CreateGriddedMetrics.py | 26 +++---- scripts/config.json | 16 +++++ 6 files changed, 105 insertions(+), 145 deletions(-) create mode 100644 scripts/config.json diff --git a/README.md b/README.md index f3d24da..f3d6e45 100644 --- a/README.md +++ b/README.md @@ -13,27 +13,29 @@ Verify that the file structure of your data is compatible with these scripts. ## lidar processing scripts ### `scripts/AP` -The FUSION Processing Scripts. Edit these at your own risk. +The FUSION Processing Scripts. Edit these at your own risk. + +### `scripts/config.json` +Edit this JSON to define the paths and processing parameters to use across the three processing scripts. The specific parameters and how they're used in the script details below. ### `scripts/01_PrepareDataForFusion.py` This script calls PDAL and FUSION. -User needs to edit the following: +User needs to edit the following in the config.json: - `project` - the name of the lidar project -- `dirLidarOriginal` - file path to the input lidar data -- `dirFUSION` - file path to FUSION executables +- `dirInputBase` - base folder path to the input lidar data. All LAZ files should be in `[dirInputBase]\\project\\Points\\LAZ` +- `dirFUSION` - folder path to FUSION executables - `nCoresMax` - maximum number of processing cores available -- `dirBase` - main output directory +- `dirBase` - main output directory of preparations +- `srsIn` - string with input horizontal and optionally vertical EPSG code(s) of data, e.g. "6428" or "6428+6360". Use *null* if the spatial reference is correctly defined in the LAZ header. ### `scripts/02_CreateAPSettingsPRP.R` This script calls FUSION. -User needs to edit the following: -- `SRS.Lidar` - Spatial refernce system of the lidar files +User needs to edit the following in the config.json: - `dirFUSION` - file path to FUSION executables - `project` - the name of the lidar project - `CELLSIZE` - Spatial resolution of the FUSION products -- `NCORESMAX` - maximum number of processing cores available -- `DIR_BASE` - main output directory -- `DIR_LIDAR` - directory of the lidar files projected in `01_PrepareDataForFusion.py` +- `nCoresMax` - maximum number of processing cores available +- `dirBase` - main output directory - `DIRSCRIPTS` - directory of the FUSION AreaProcessor scripts (these scripts are found in `scripts/AP`) ### `scripts/03_CreateGriddedMetrics.py` @@ -41,15 +43,17 @@ This script calls FUSION. User needs to edit the following: - `project` - the name of the lidar project - `dirBase` - main output directory -- `dirFinalProducts` - directory where the final products should be saved +- `dirFinalProducts` - base directory where the final products should be saved ## Usage Setup workflow. -Run `scripts/01_PrepareDataForFusion.py` and `scripts/02_CreateAPSettingsPRP.R`. -The output of `scripts/02_CreateAPSettingsPRP.R` is a PRP file that is used to set up the FUSION processing run. +Edit the config.json with the correct paths and parameters. Call all scripts from the command line with the path to the config.json as the first argument (e.g., `python 01_PrepareDataForFusion.py C:\\path\\to\\config.json` or `Rscript.exe 02_CreateAPSettingsPRP.R C:\\path\\to\\config.json`). +Run `scripts/01_PrepareDataForFusion.py`, which will output a copy of the files in EPSG:5070 along with a QAQC of the data and file list to be used by the next script. +Run `scripts/02_CreateAPSettingsPRP.R` which outputs is a PRP file that is used to set up the FUSION processing run. Open the FUSION program `AreaProcessor.exe` and load the PRP file. Create the processing layout. Create the processing scripts. Run `scripts/03_CreateGriddedMetrics.py`. This script runs the batch file created in `[DIR_BASE]/[studyArea]/Processing/AP/APFusion.bat`, cleans the FUSION grids, and copies various products to a user-specified directory. -Note: There is an alternative script `scripts/01_PrepareDataForFusion_MultiProjects.py` that is designed to loop over multiple lidar projects. The advantage of this script is FUSION QAQC is run in parallel with one project per job. Users are welcome to alter the other scripts such that they loop over multiple projects. +Note: There is an alternative script `scripts/01_PrepareDataForFusion_MultiProjects.py` that is designed to loop over multiple lidar projects. The advantage of this script is FUSION QAQC is run in parallel with one project per job. To use this script edit the config.json to provide a list of projects (e.g. ["CO_TellerCo_2016", "CO_Aspen_2020"], and a corresponding list of srsIn can also be provided (e.g. ["6428+6360", null]), but if null then the SRS will be read from the LAZ header for all projects. +Users are welcome to alter the other scripts such that they loop over multiple projects. diff --git a/scripts/01_PrepareDataForFusion.py b/scripts/01_PrepareDataForFusion.py index 78fa2d3..7762631 100644 --- a/scripts/01_PrepareDataForFusion.py +++ b/scripts/01_PrepareDataForFusion.py @@ -22,11 +22,13 @@ # ----------------------------------------------------------------------------- import os import shutil -import pdal +import sys import json import time -from joblib import Parallel, delayed import subprocess +import pdal +from types import SimpleNamespace +from joblib import Parallel, delayed # ----------------------------------------------------------------------------- @@ -36,60 +38,26 @@ # ----------------------------------------------------------------------------- start = time.time() -# Assign a project -project = "CO_ARRA_ParkCo_2010" +# Load config data into namespace +config_path = sys.argv[1] +with open(config_path, 'r') as f: + config_dict = json.load(f) +config = SimpleNamespace(**config_dict) -print(project) +print(config.project) # Directory of lidar data needing to be processed (e.g., external HHD) -dirLidarOriginal = os.path.join(r"L:\Lidar", project, "Points", "LAZ") - -# FUSION directory -dirFUSION = r"C:\Fusion" - -# Maximum number of processing cores -nCoresMax = 26 +dirLidarOriginal = os.path.join(config.dirInputBase, config.project, "Points", "LAZ") # main output directory -dirBase = r"D:\LidarProcessing" -if not os.path.exists(dirBase): - os.mkdir(dirBase) +if not os.path.exists(config.dirBase): + os.mkdir(config.dirBase) # directory for the specific lidar project; HOME_FOLDER in FUSION scripts -dirHomeFolder = os.path.join(dirBase, project) +dirHomeFolder = os.path.join(config.dirBase, config.project) if not os.path.exists(dirHomeFolder): os.mkdir(dirHomeFolder) -# ----------------------------------------------------------------------------- -# ----------------------------------------------------------------------------- -# Assign Spatial Reference Systems -# ----------------------------------------------------------------------------- -# ----------------------------------------------------------------------------- - -# PDAL can read SRS directly from the LAS file. Some of the older LAS files -# do not have embeded SRS info. If you want to use the SRS in the lidar file, -# then do not include the key:pair (studyArea:epgsCode) in the dictionary. -# If the SRS is missing from the lidar file, or you want to be explicity, -# include the SRS in the dictionary - -# Spatial Reference Systems of different project areas -dictSRS = { - "CO_CheesmanLake_2004": 26913, - "CO_DenverDNC_2008": 26913, - "CO_Denver_2008": 26913, - "CO_ElPasoCoCentral_2018": "6428+8228", - "CO_FremontCo_2016": "6428+8228", - "CO_HuerfanoCo_2018": "6432+8228", - "CO_LarimerCo_GlenHaven_2013": "2231+8228", - "CO_LovelandE_2016": "6430+8228", - "CO_LovelandW_2016": "6430+8228", - "CO_MesaCo_QL2_2015": 6341, - "CO_PitkinCo_2016": "6428+8228", - "CO_RouttCo_2016": "6430+8228", - "WY_Casper_2010": "", - "WY_NRCS_LIDAR_2006": 26913, -} - # ----------------------------------------------------------------------------- # ----------------------------------------------------------------------------- # Define Functions @@ -185,7 +153,6 @@ def calcNCores(x, nCoresMax): lidarFilesOriginal = os.listdir(dirLidarOriginal) lidarFilesOriginal.sort() - # ---------------------------------------------------------------------------- # Process Point Data # ---------------------------------------------------------------------------- @@ -211,11 +178,6 @@ def calcNCores(x, nCoresMax): # project to EPSG 5070 print("\tProjecting Lidar Files") -srsNeedsDefining = project in dictSRS -if srsNeedsDefining: - srsIn = dictSRS[project] -else: - srsIn = None dirLAZ5070 = os.path.join(dirPoints, "LAZ5070") if not os.path.exists(dirLAZ5070): @@ -224,9 +186,9 @@ def calcNCores(x, nCoresMax): lidarFilesCopy = os.listdir(dirLidarCopy) lidarFilesCopy.sort() -nCores = calcNCores(lidarFilesCopy, nCoresMax) +nCores = calcNCores(lidarFilesCopy, config.nCoresMax) Parallel(n_jobs=nCores)( - delayed(parallelProjectFunc)(lidarFile, dirLidarCopy, dirLAZ5070, srsIn) + delayed(parallelProjectFunc)(lidarFile, dirLidarCopy, dirLAZ5070, config.srsIn) for lidarFile in lidarFilesCopy ) del nCores @@ -267,7 +229,7 @@ def calcNCores(x, nCoresMax): del lidarFiles5070 # Run Catalog -exeCatalog = os.path.join(dirFUSION, "Catalog.exe") +exeCatalog = os.path.join(config.dirFUSION, "Catalog.exe") fpQAQCOut = os.path.join(dirQAQC, "QAQC.csv") cmdCatalog = ( diff --git a/scripts/01_PrepareDataForFusion_MultiProjects.py b/scripts/01_PrepareDataForFusion_MultiProjects.py index a8e3ca3..c53128e 100644 --- a/scripts/01_PrepareDataForFusion_MultiProjects.py +++ b/scripts/01_PrepareDataForFusion_MultiProjects.py @@ -22,9 +22,10 @@ import pdal import json import time -from joblib import Parallel, delayed +import sys import subprocess - +from types import SimpleNamespace +from joblib import Parallel, delayed # ----------------------------------------------------------------------------- # ----------------------------------------------------------------------------- @@ -33,8 +34,11 @@ # ----------------------------------------------------------------------------- start = time.time() -# Assign a project -projects = ["CO_ARRA_ParkCo_2010", "CO_ARRA_GrandCo_2010"] +# Load config data into namespace +config_path = sys.argv[1] +with open(config_path, 'r') as f: + config_dict = json.load(f) +config = SimpleNamespace(**config_dict) # ----------------------------------------------------------------------------- @@ -50,22 +54,7 @@ # include the SRS in the dictionary # Spatial Reference Systems of different project areas -dictSRS = { - "CO_CheesmanLake_2004": 26913, - "CO_DenverDNC_2008": 26913, - "CO_Denver_2008": 26913, - "CO_ElPasoCoCentral_2018": "6428+8228", - "CO_FremontCo_2016": "6428+8228", - "CO_HuerfanoCo_2018": "6432+8228", - "CO_LarimerCo_GlenHaven_2013": "2231+8228", - "CO_LovelandE_2016": "6430+8228", - "CO_LovelandW_2016": "6430+8228", - "CO_MesaCo_QL2_2015": 6341, - "CO_PitkinCo_2016": "6428+8228", - "CO_RouttCo_2016": "6430+8228", - "WY_Casper_2010": "", - "WY_NRCS_LIDAR_2006": 26913, -} +dictSRS = dict(zip(config.project, config.srsIn)) # ----------------------------------------------------------------------------- # ----------------------------------------------------------------------------- @@ -194,28 +183,28 @@ def calcNCores(x, nCoresMax): # ----------------------------------------------------------------------------- # ----------------------------------------------------------------------------- -for project in projects: +for project in config.project: # ------------------------------------------------------------------------- # Setup # ------------------------------------------------------------------------- print(project) # Directory of lidar data needing to be processed (e.g., external HHD) - dirLidarOriginal = os.path.join(r"L:\Lidar", project, "Points", "LAZ") + dirLidarOriginal = os.path.join(config.dirInputBase, project, "Points", "LAZ") # FUSION directory - dirFUSION = r"C:\Fusion" +# dirFUSION = r"C:\Fusion" # Maximum number of processing cores - nCoresMax = 26 +# nCoresMax = 26 # main output directory - dirBase = r"D:\LidarProcessing" - if not os.path.exists(dirBase): - os.mkdir(dirBase) +# dirBase = r"D:\LidarProcessing" + if not os.path.exists(config.dirBase): + os.mkdir(config.dirBase) # directory for the specific lidar project; HOME_FOLDER in FUSION scripts - dirHomeFolder = os.path.join(dirBase, project) + dirHomeFolder = os.path.join(config.dirBase, project) if not os.path.exists(dirHomeFolder): os.mkdir(dirHomeFolder) @@ -247,8 +236,7 @@ def calcNCores(x, nCoresMax): # project to EPSG 5070 print("\tProjecting Lidar Files") - srsNeedsDefining = project in dictSRS - if srsNeedsDefining: + if config.srsIn: srsIn = dictSRS[project] else: srsIn = None @@ -260,7 +248,7 @@ def calcNCores(x, nCoresMax): lidarFilesCopy = os.listdir(dirLidarCopy) lidarFilesCopy.sort() - nCores = calcNCores(lidarFilesCopy, nCoresMax) + nCores = calcNCores(lidarFilesCopy, config.nCoresMax) Parallel(n_jobs=nCores)( delayed(parallelProjectFunc)(lidarFile, dirLidarCopy, dirLAZ5070, srsIn) for lidarFile in lidarFilesCopy @@ -310,9 +298,9 @@ def calcNCores(x, nCoresMax): # ---------------------------------------------------------------------------- print("\nRunning FUSION Catalog\n") -nCores = calcNCores(projects, nCoresMax) +nCores = calcNCores(config.project, config.nCoresMax) Parallel(n_jobs=nCores)( - delayed(parallelRunQAQC)(project, dirBase, dirFUSION) for project in projects + delayed(parallelRunQAQC)(project, config.dirBase, config.dirFUSION) for project in config.project ) del nCores diff --git a/scripts/02_CreateAPSettingsPRP.R b/scripts/02_CreateAPSettingsPRP.R index f7d18ce..49e5184 100644 --- a/scripts/02_CreateAPSettingsPRP.R +++ b/scripts/02_CreateAPSettingsPRP.R @@ -18,6 +18,7 @@ rm(list=ls()) # ---------------------------------------------------------------------------- library(rgdal) +library(rjson) # ---------------------------------------------------------------------------- # ---------------------------------------------------------------------------- @@ -25,29 +26,15 @@ library(rgdal) # ---------------------------------------------------------------------------- # ---------------------------------------------------------------------------- -#Spatial reference system of lidar data (e.g., EPSG 5070) -SRS.Lidar = CRS("+init=epsg:5070") - -#File path to Fusion Executuables -dirFUSION <- "C:\\FUSION\\" - -#Lidar Project Name -project <- "CO_ARRA_ParkCo_2010" - -#Raster resolution -CELLSIZE <- 30 +# Load variables from JSON +args = commandArgs(trailingOnly=TRUE) +config_path = args[1] +config <- fromJSON(file=config_path) -#Maximum number of processing cores -NCORESMAX <- 26 - -#Working directory (same as dirWD in script 01_PrepareDataForFusion.py) -DIR_BASE <- paste0("D:\\LidarProcessing\\") - -#Direcotry of lidar data -DIR_LAZ5070 <- paste0(DIR_BASE, project, "\\", "Points\\", "LAZ5070\\") +SRS.Lidar = CRS("+init=epsg:5070") +CELLSIZE <- config$CELLSIZE +DIR_LAZ5070 <- paste0(config$dirBase, "\\", config$project, "\\", "Points\\", "LAZ5070\\") -#Directory of AP scripts -DIRSCRIPTS <- "C:\\Users\\pafekety\\Desktop\\CMS2LidarProcessing\\scripts\\AP" # Do not include the trailing \\ # ---------------------------------------------------------------------------- @@ -113,7 +100,7 @@ RunDtmDescribe <- function(DtmDirectory, OutputDirectory, OutputName, ALL=TRUE, #Run Fusion DTMDescribe shell( - paste0(dirFUSION, "DTMDescribe.exe ", DtmFileTxt, " ", OutputDirectory, OutputName, "_Summary.csv") + paste0(config$dirFUSION, "\\DTMDescribe.exe ", DtmFileTxt, " ", OutputDirectory, OutputName, "_Summary.csv") ) } @@ -387,6 +374,7 @@ writeSection6 <- function(DIR_LIDAR, ext = ""){ getLidarExtent <- function(DIR_LIDAR, ext = ""){ #file path for lidar files + print(DIR_LIDAR) ext = toupper(ext) if(!is.element(ext, c("LAS", "LAZ"))) stop ("ext must be LAS or LAZ") if(ext == "LAS"){ @@ -767,11 +755,11 @@ writePRP <- function(PROJECT, LATITUDE, HOMEFOLDER, PROCESSINGHOME, DIRSCRIPTS, #Ground data if(length(dir(DTMSPEC)) > 0){ - writeSection7(DTMSPEC = DTMSPEC, PROCESSINGHOME = PROCESSINGHOME, dirFUSION = dirFUSION) + writeSection7(DTMSPEC = DTMSPEC, PROCESSINGHOME = PROCESSINGHOME, dirFUSION = config$dirFUSION) } #Density data - writeSection8(PRODUCTHOME = PRODUCTHOME, PROCESSINGHOME = PROCESSINGHOME, dirFUSION = dirFUSION) + writeSection8(PRODUCTHOME = PRODUCTHOME, PROCESSINGHOME = PROCESSINGHOME, dirFUSION = config$dirFUSION) #Mask Layers writeSection9() @@ -803,7 +791,7 @@ createPRP <- function(project, cellSize, nCores, DIR_BASE, DIRSCRIPTS, DIR_LIDAR print(paste0("Creating PRP file for ", project));flush.console() #These folders are following the LTK naming scheme from setup.bat - HOMEFOLDER <- paste0(DIR_BASE, project, '\\') + HOMEFOLDER <- paste0(DIR_BASE, "\\", project, '\\') PROCESSINGHOME <- paste0(HOMEFOLDER, 'Processing\\AP\\') DTMSPEC <- paste0(HOMEFOLDER, 'Deliverables\\DTM\\') PRODUCTHOME <- paste0(HOMEFOLDER, "Products\\") @@ -889,12 +877,12 @@ createPRP <- function(project, cellSize, nCores, DIR_BASE, DIRSCRIPTS, DIR_LIDAR # ---------------------------------------------------------------------------- createPRP( - project = project, - cellSize = CELLSIZE, - nCores = NCORESMAX, - DIR_BASE = DIR_BASE, - DIRSCRIPTS = DIRSCRIPTS, + project = config$project, + cellSize = config$CELLSIZE, + nCores = config$nCoresMax, + DIR_BASE = config$dirBase, + DIRSCRIPTS = config$DIRSCRIPTS, DIR_LIDAR = DIR_LAZ5070, - dirFUSION = dirFUSION + dirFUSION = config$dirFUSION ) diff --git a/scripts/03_CreateGriddedMetrics.py b/scripts/03_CreateGriddedMetrics.py index 2e535c3..e3ecf04 100644 --- a/scripts/03_CreateGriddedMetrics.py +++ b/scripts/03_CreateGriddedMetrics.py @@ -21,6 +21,9 @@ import shutil import subprocess import time +import json +import sys +from types import SimpleNamespace import rasterio as rio # from rasterio.plot import show @@ -33,17 +36,16 @@ # ----------------------------------------------------------------------------- start = time.time() -# Assign a project -project = "CO_ARRA_ParkCo_2010" -print(project) +# Load config data into namespace +config_path = sys.argv[1] +with open(config_path, 'r') as f: + config_dict = json.load(f) +config = SimpleNamespace(**config_dict) -# main output directory -dirBase = r"D:\LidarProcessing" +print(config.project) -# directory where the final products should be saved -dirFinalProducts = r"G:\FusionRuns" -if not os.path.exists(dirFinalProducts): - os.mkdir(dirFinalProducts) +if not os.path.exists(config.dirFinalProducts): + os.mkdir(config.dirFinalProducts) # ----------------------------------------------------------------------------- @@ -87,7 +89,7 @@ def cleanGrids(inFiles, outDir, fpElev): # ----------------------------------------------------------------------------- # ----------------------------------------------------------------------------- # directory for the specific lidar project; HOME_FOLDER in FUSION scripts -dirHomeFolder = os.path.join(dirBase, project) +dirHomeFolder = os.path.join(config.dirBase, config.project) # Directory of FUSION products dirFusionProducts = os.path.join(dirHomeFolder, "Products") @@ -98,7 +100,7 @@ def cleanGrids(inFiles, outDir, fpElev): # Put in a check to see if the script has run if not os.path.exists(os.path.join(dirFusionProducts, "complete.txt")): - print("Running FUSION for " + project) + print("Running FUSION for " + config.project) subprocess.run(cmdFusionMainBat, shell=True) while not os.path.exists(os.path.join(dirFusionProducts, "complete.txt")): @@ -111,7 +113,7 @@ def cleanGrids(inFiles, outDir, fpElev): # ----------------------------------------------------------------------------- # Project-specific directory for saved outputs -dirOutProject = os.path.join(dirFinalProducts, project) +dirOutProject = os.path.join(config.dirFinalProducts, config.project) if not os.path.exists(dirOutProject): os.mkdir(dirOutProject) diff --git a/scripts/config.json b/scripts/config.json new file mode 100644 index 0000000..090609a --- /dev/null +++ b/scripts/config.json @@ -0,0 +1,16 @@ +{ + "project": "CO_Aspen_2020", + "dirInputBase": "F:\\Lidar", + "dirFUSION": "C:\\Fusion", + "dirBase": "C:\\scratch\\LidarProcessing", + "nCoresMax": 6, + "srsIn": null, + + + "SRS.Lidar": "+init=epsg:5070", + "CELLSIZE": 30, + "DIRSCRIPTS": "J:\\users\\stevenf\\code\\fork\\CMS2LidarProcessing\\scripts\\AP", + + "dirFinalProducts": "C:\\scratch\\FusionRuns" +} + \ No newline at end of file From 6fe3c914822dc1b81a202d3586bf16c3961cd48c Mon Sep 17 00:00:00 2001 From: Steven Filippelli <37000538+FevetS@users.noreply.github.com> Date: Wed, 26 Jun 2024 11:55:10 -0600 Subject: [PATCH 2/3] Update config.json --- scripts/config.json | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/scripts/config.json b/scripts/config.json index 090609a..c50f7cd 100644 --- a/scripts/config.json +++ b/scripts/config.json @@ -6,11 +6,9 @@ "nCoresMax": 6, "srsIn": null, - - "SRS.Lidar": "+init=epsg:5070", "CELLSIZE": 30, "DIRSCRIPTS": "J:\\users\\stevenf\\code\\fork\\CMS2LidarProcessing\\scripts\\AP", "dirFinalProducts": "C:\\scratch\\FusionRuns" } - \ No newline at end of file + From e02192cc286c0990b2e11df839df890f6cd4ab4b Mon Sep 17 00:00:00 2001 From: Steven Filippelli <37000538+FevetS@users.noreply.github.com> Date: Wed, 26 Jun 2024 11:58:29 -0600 Subject: [PATCH 3/3] Update 01_PrepareDataForFusion_MultiProjects.py --- scripts/01_PrepareDataForFusion_MultiProjects.py | 10 +--------- 1 file changed, 1 insertion(+), 9 deletions(-) diff --git a/scripts/01_PrepareDataForFusion_MultiProjects.py b/scripts/01_PrepareDataForFusion_MultiProjects.py index c53128e..6357d1a 100644 --- a/scripts/01_PrepareDataForFusion_MultiProjects.py +++ b/scripts/01_PrepareDataForFusion_MultiProjects.py @@ -191,15 +191,7 @@ def calcNCores(x, nCoresMax): # Directory of lidar data needing to be processed (e.g., external HHD) dirLidarOriginal = os.path.join(config.dirInputBase, project, "Points", "LAZ") - - # FUSION directory -# dirFUSION = r"C:\Fusion" - - # Maximum number of processing cores -# nCoresMax = 26 - - # main output directory -# dirBase = r"D:\LidarProcessing" + if not os.path.exists(config.dirBase): os.mkdir(config.dirBase)