dhardestylewis / dem2basin Goto Github PK

Preprocessing scripts to produce watershed-delineated input files from standard source data

License: GNU General Public License v3.0

Python 4.93% Shell 0.04% Jupyter Notebook 95.02%

python digital-elevation-model raster vector national-hydrography-dataset watershed lidar

dem2basin's Introduction

DEM to basin preprocessing

This preprocessing library takes source 1-meter digital elevation model (DEM) data and splits, crops, buffers, and reprojects it to individual hydrologic basins (identified by their unique identifier, the "HUC12" ID).

This preprocessing script also produces ancillary data products corresponding to each new HUC12 DEM raster to describe their sub-basins (ie "catchments"), their streams (ie "flowlines"), and the roughness of each streambed.

Taken together, these are the major inputs needed to run GeoFlood, which creates short-term flood projections.

Alpha-release Python library

The underlying Python library is available on PyPi and can be installed by:

pip install dem2basin

Main Python script

The recommended way to run dem2basin.py:

python3 dem2basin.py \
    --shapefile study_area_polygon.shp \
    --huc12 WBD-HUC12s.shp \
    --nhd NHD_catchments_and_flowlines.gdb/ \
    --raster TNRIS-LIDAR-Datasets/ \
    --availability TNRIS-LIDAR-Dataset_availability.shp \
    --directory HUC12-DEM_outputs/ \
    --restart dem2basin-study_area.pickle

Required source data inputs

There are 5 required inputs.

--shapefile Study area polygon vector GIS file: a vector GIS file of a single polygon which defines the study area
--huc12 USGS Watershed Boundary Dataset (WBD): a HUC12 vector GIS file from USGS's WBD or a subset of the WBD
--nhd NHD Medium Resolution (MR): a GeoDataBase of NHD MR catchments and flowlines
--raster TNRIS Lidar: the parent directory to a collection of TNRIS Lidar datasets
--availability TNRIS Lidar availability: a vector GIS file of TNRIS Lidar availability, provided by TNRIS here

Optional parameters

--directory Outputs directory: a directory to store outputs, which will each be sorted by HUC12
--restart Restart file: a Python Pickle file from which you can restart the preprocessing if it's interrupted
--overwrite Overwrite flag: optional flag to overwrite all files found in output directory
--overwrite_rasters Overwrite rasters flag: optional flag to overwrite just the raster outputs
--overwrite_flowlines Overwrite flowlines flag: optional flag to overwrite just the flowline outputs
--overwrite_catchments Overwrite catchments flag: optional flag to overwrite just the catchment outputs
--overwrite_roughnesses Overwrite roughness table flag: optional flag to overwrite the roughness table
--log Log file: a file to store runtime log

Description of outputs

There are 4 outputs per HUC12.

Cropped & buffered DEM:
- buffered 500m
- cropped to each HUC12 intersecting the study area
- at least 1m resolution
- mosaicked with preference for lowest resolution tiles
- reprojected to the study area's projections
corresponding NHD MR flowlines:
- subset of NHD MR flowlines
- each flowline's median point along the line lies within the HUC12
- reprojected to the study area's projections
corresponding NHD MR catchments:
- subset of NHD MR catchments
- correspond with the NHD MR flowlines above
- reprojected to the study area's projections
Manning's n roughness table:
- organized by flowline using their ComIDs
- vary by stream order

Here is an example of these outputs, originally visualized by Prof David Maidment.

Already preprocessed DEMs

Already preprocessed DEMs are now available for the vast majority of Texas's HUC12s if you are a TACC user. You can request a TACC account here.

Notes about preprocessed DEMs

The DEMs are not provided for any HUC12s that have any gap in 1m resolution data.
All of the DEMS are reprojected to WGS 84 / UTM 14N, even if the HUC12 is outside of UTM 14.

Where to find them

The DEMs are located on Stampede2 at /scratch/projects/tnris/dhl-flood-modelling/TX-HUC12-DEM_outputs.

If you run into trouble

Please submit a ticket if you have trouble accessing this data. You may also contact me directly at @dhardestylewis or [email protected]

Available preprocessed HUC12s

These HUC12 DEMs are available right now on Stampede2.

Confirmed successfully preprocessed HUC12s

These HUC12 DEMs have been successfully preprocessed in the past, and will soon be available once again on Stampede2. If you need any of these right now, please contact me.

Preprocessing workflow

If you would like an understanding of the preprocessing workflow, I provide a simplified but representative example in this Jupyter notebook. This Jupyter notebook was presented at the inaugural TACC Institute on Planet Texas 2050 Cyberecosystem Tools in August, 2020. Please contact me if you would like a recording.

dem2basin's People

Contributors

Stargazers

Watchers

Forkers

edgargarza7 fernando-aristizabal khushbooagarwalhub

dem2basin's Issues

Draft workflow diagram of this data pipeline

Manually re-write 2 remaining missing dirnames in coverage file

For remaining 2 mis-written dirnames in Lidar availability file, consider manually rewriting just these 2

Draft Jupyter notebooks to demonstrate basic use-cases of library

Basic use-cases include:

GeoFlood workflow preparation
Non-HUC-oriented workflow preparation
HEC-RAS preparation excluding road network / bridge segment aspects
10m
3m
1m
Backfill 1m with 3m or 10m
QA/QC existing DEMs for common issues including nodata gaps, missing source data, and wrong/bad projections
Determine accuracy of resulting products

Consider using crosswalk table to associate NHD MR and WBD

Consider using crosswalk table to associate NHD MR and WBD:
https://www.sciencebase.gov/catalog/item/5c86a747e4b09388244b3da1

Draft source-to-visualization pipeline

Include rough measures of completion of each aspect, basic features, links to various repositories, etc

Significantly revise get_bounding_boxes_by_project to center around project-level GDF objects

As soon as possible, revise the following code to operate on project_coverage GeoDataFrames, with a lambda function to find the files.

Also, dynamically compile a list of separate GeoDataFrames for each project, each row within these GDFs would represent each different DEM tile within that project's GDF, use geopandas.GeoDataFrame.to_crs to reproject these all to the same CRS, before concatenating them to a single master coverage file, which is the final output

    dem_tilenames = []
    filetypes = ('.img', '.dem', '.tif', '.jp2')
    for filetype in filetypes:
#        print(filetype)
        for project in projects:
#            print(project)
            for root, dirs, filenames in os.walk(str(project)):
#                print(root,dirs,filenames)
                for filename in filenames:
#                    print(filename)
                    if filename.endswith(filetype):
#                        print(True)
                        dem_tilenames.append(os.path.join(root,filename))

#    print(dem_tilenames)

    dem_tile_bounds = []
    for dem_tilename in dem_tilenames:
        dem_tile_bounds.append((
            dem_tilename,
            rasterio.open(dem_tilename).bounds
        ))

#    print(dem_tile_bounds)

    dem_tiles = gpd.GeoDataFrame(
        dem_tile_bounds,
        columns = ['lidar_file','bounds']
    )

#    print(dem_tiles)

Download all relevant datasets

Troubleshoot minor GDAL issues that arise during mosaicking at scale

Include following information in TNRIS Lidar metadata database

File name of DEM as found on Stampede2
Actual file type of each DEM (.dem,.img,.tif)
CRS
Geographic CRS
EPSG
Bounding box
Accuracy calculations for reprojections from source CRS to NAD83(2011) / UTM zone 14N & NAD83(2011) Texas Centric Albers Equal Area
Polygonizing the rasters
File size

Convert de facto standard datasets to existing TNRIS Lidar data model

Build basic functions to convert these datasets to existing TNRIS Lidar availability data model

Determine accuracy of general reprojection & mosaicking workflow

Especially whether inaccuracies introduced through so many reprojections are significant enough to matter

Check for bugs associated with not using select_utm option in various functions

Check for bugs associated with not using select_utm option in various functions

Restore functionality for using dem2basin as a script

Gather boundary data for TNRIS Lidar point-clouds

potentially could use some of the same methodologies developed in DEM2basin for DEMs (just a derivation of those functions)

Support for point-cloud to DEMs

Consider adding support for point-cloud to DEMs, but:

record timings
compare against DEM reprojection and mosaicking technique
estimate large-scale node-hours usages

Verify geometric association algorithm at scale across Texas

later consider verifying nationally

Include support for backfilling

Current approach: Merging TNRIS Lidar data modelled GDFs and mosaicking on highest resolution

Check to make sure prepare_and_output_geoflood_gis_inputs_by_huc works by default against National WBD

Slightly modify _get_mosaic_and_output_raster_simple for more general use-case

make dependence on having a coverage file optional
include option to explicitly specify a mask object (whether a polygon or a file)

Update Lidar index

Run regular scrape of latest Lidar data using tools developed in tnris-lidardata-update

Get GeoFlood Fathom3m timings

Get GeoFlood Fathom3m timings back to Mark Wang & Paola Passalacqua

Run pycallgraph against prepare_and_output_geoflood with/out new coverage file

Run pycallgraph against prepare_and_output_geoflood with/with-out new coverage file to determine why the new coverage file is not yet being created

Continue function-by-function documentation

Prioritize by functions appearing in Jupyter notebooks

Revise documentation & Jupyter notebook to reflect new directory structure

Provide HUC-by-HUC selection of TNRIS Lidar tiles

Provide HUC-by-HUC selection of TNRIS Lidar tiles, in a similar manner to the following web-portal:

https://www.coast.noaa.gov/inventory/

Validate general mosaicking algo from each pair of 2 adjacent TNRIS Lidar projects

Check area calculation difference before and after reproj

Check area calculation difference before and after reprojection to an equal area CRS

Include timing function

Include timing wrapper function:

def timethis(f):
    start_time = time.time()
    ret = f
    end_time = time.time()
    elapsed_time = end_time - start_time
    return((ret,elapsed_time))

Time large-scale use of pathlib vs os.path

Time for UnparsedPath in the following callgraph of DEM2basin indicates that the pathlib library may be too heavy to use at scale:
https://github.com/Khushbooagarwalhub/DEM2basin/blob/dev-kag/references/pycallgraph.png

Upload latest version to PyPi successfully

Add function to fill in 0-value in mosaicked rasters for the coast

double-check the datum

Missing catchments

Catch catchments not assigned to any HUC12 by algorithm. An example of such missing catchments is present in this screenshot

Consider the following back-ends for performance

Consider the following back-ends for performance:

In time-expensive functions (see pycallgraph.png), consider replacing Shapely-dependent objects with GDAL/OGR objects
This would also mean using either:

NetCDF, or
PostgreSQL

for I/O

.tif
.dem
.img
.jp2

Check under both */dem/* and */tiles/*

Make Python libraries' developers aware of developed extensions

Such libraries include:

Geopandas
Pandas
GDAL/OGR
Fiona
Shutil
etc

Revise README.md for use as a library of functions

Break up code to speed up each individual import

Adopt the following organization:

General GeoPandas extensions
Pandas extensions
Function wrappers
Misc extensions of other Python libraries
Hydro vector processing functions
Hydro raster processing functions