[1pt] PR: v4.5.0.0 bug fixes

.github/PULL_REQUEST_TEMPLATE.md

@@ -12,11 +12,11 @@
 
 _How does the changes affect the product?_
 - [ ] Code only?
-- [ ] Require new or adjusted data inputs. Does it have start, end and duration code (in UTC)?
-- [ ] If new or updated data updates, has FIM related code been updated and tested?
+- [ ] Require new or adjusted data inputs? Does it have start, end and duration code (in UTC)?
+- [ ] If new or updated data sets, has the FIM code been updated and tested with the new/adjusted data (subset is fine, but must be a subset of the new data)?
 - [ ] Require new pre-clip set?
 - [ ] Has new or updated python packages?
-- [ ] If applicable, is the deployment plan be calculated and work with Deployment person/team?
+- [ ] If applicable, has a deployment plan be created with the deployment person/team?
 
 
 ### Issuer Checklist (For developer use)

docs/CHANGELOG.md

@@ -1,6 +1,22 @@
 All notable changes to this project will be documented in this file.
 We follow the [Semantic Versioning 2.0.0](http://semver.org/) format.
 
+## v4.5.0.1 - 2024-05-07 - [PR#1150](https://github.com/NOAA-OWP/inundation-mapping/pull/1150)
+
+Fixes two bugs discovered in v4.5.0.0:
+1. `echo` missing in bash command
+2. raster resolution of `dem_meters.tif` has now been explicitly set in `gdalwarp`.
+
+### Changes
+
+- `src/`
+    - `bash_variables.env`: added `res` environment variable for default raster cell size
+    - `delineate_hydros_and_produce_HAND.sh`: added missing `echo`
+    - `heal_bridges_osm.py`: fixed raster resolution and number of rows/columns
+    - `run_unit_wb.sh`: added `-tr` to gdalwarp when generating `dem_meters.tif`; removed extraneous `Tcount`
+
+<br/><br/>
+
 ## v4.5.0.0 - 2024-05-06 - [PR#1122](https://github.com/NOAA-OWP/inundation-mapping/pull/1122)
 
 This PR includes 2 scripts to add Open Street Map bridge data into the HAND process: a script that pulls data from OSM and a script that heals those bridges in the HAND grids. Both scripts should be run as part of a pre-processing step for FIM runs. They only need to be run if we think OSM data has changed a lot or for any new FIM versions.

src/add_crosswalk.py

@@ -244,12 +244,32 @@ def add_crosswalk(
                         (output_flows.From_Node == to_node) & (output_flows['order_'] == max_order)
                     ]['HydroID'].item()
 
+                # output_flows has a higher order than the max_order
+                elif output_flows.loc[(output_flows.From_Node == to_node), 'order_'].max() > max_order:
+                    update_id = output_flows.loc[
+                        (output_flows.From_Node == to_node)
+                        & (
+                            output_flows['order_']
+                            == output_flows.loc[(output_flows.From_Node == to_node), 'order_'].max()
+                        )
+                    ]['HydroID'].values[0]
+
                 # Get the first one
                 # Same stream order, without drainage area info it is hard to know which is the main channel.
                 else:
-                    update_id = output_flows.loc[
-                        (output_flows.From_Node == to_node) & (output_flows['order_'] == max_order)
-                    ]['HydroID'].values[0]
+                    if max_order in output_flows.loc[output_flows.From_Node == to_node, 'order_'].values:
+                        update_id = output_flows.loc[
+                            (output_flows.From_Node == to_node) & (output_flows['order_'] == max_order)
+                        ]['HydroID'].values[0]
+
+                    else:
+                        update_id = output_flows.loc[
+                            (output_flows.From_Node == to_node)
+                            & (
+                                output_flows['order_']
+                                == output_flows.loc[output_flows.From_Node == to_node, 'order_'].max()
+                            )
+                        ]['HydroID'].values[0]
 
             # no upstream segments; single downstream segment
             elif len(output_flows.loc[output_flows.From_Node == to_node]['HydroID']) == 1:

src/bash_variables.env

@@ -1,8 +1,11 @@
 ## Define inputs
 
-# NOTE: $inputsDir is defined in Dockerfile
+export res=10 # Default raster cell size in meters
+
 export             DEFAULT_FIM_PROJECTION_CRS=EPSG:5070
 export                             ALASKA_CRS=EPSG:3338
+
+# NOTE: $inputsDir is defined in Dockerfile
 export                              input_DEM=${inputsDir}/3dep_dems/10m_5070/fim_seamless_3dep_dem_10m_5070.vrt
 export                       input_DEM_Alaska=${inputsDir}/3dep_dems/10m_South_Alaska/23_11_07/FIM_3dep_dem_South_Alask_10m.vrt
 export                       input_DEM_domain=${inputsDir}/3dep_dems/10m_5070/HUC6_dem_domain.gpkg

src/delineate_hydros_and_produce_HAND.sh

@@ -240,9 +240,9 @@ if  [ -f $tempHucDataDir/osm_bridges_subset.gpkg ]; then
                 -s $tempHucDataDir/osm_bridges_subset.gpkg \
                 -p $tempCurrentBranchDataDir/bridge_lines_raster_$current_branch_id.tif \
                 -t $tempCurrentBranchDataDir/rem_zeroed_masked_$current_branch_id.tif \
-                -r 10
+                -r $res
 else
-    -e $startDiv"No applicable bridge data for $hucNumber"
+    echo -e $startDiv"No applicable bridge data for $hucNumber"
 fi
 
 

src/heal_bridges_osm.py

@@ -52,8 +52,8 @@ def process_bridges_in_huc(
     ########### setup new raster to save bridge max hand values #############
     bbox = hand_grid.bounds
     xmin, ymin, xmax, ymax = bbox
-    w = (xmax - xmin) // resolution
-    h = (ymax - ymin) // resolution
+    w = hand_grid.shape[1]
+    h = hand_grid.shape[0]
 
     out_meta = {
         "driver": "GTiff",
@@ -74,7 +74,7 @@ def process_bridges_in_huc(
         shapes = ((geom, value) for geom, value in zip(osm_gdf.geometry, osm_gdf.max_hand))
         # burn in values to any pixel that's touched by polygon and add nodata fill value
         burned = features.rasterize(
-            shapes=shapes, fill=-999999, out=out_arr, transform=out.transform, all_touched=True
+            shapes=shapes, fill=-999999, out=out_arr, transform=hand_grid.transform, all_touched=True
         )
         out.write_band(1, burned)
 
@@ -138,7 +138,7 @@ def process_bridges_in_huc(
         help='OPTIONAL: Resolution of HAND grid. Default value is 10m',
         required=False,
         default=10,
-        type=int,
+        type=float,
     )
 
     args = vars(parser.parse_args())

src/run_unit_wb.sh

@@ -119,13 +119,12 @@ echo -e $startDiv"Clipping rasters to branches $hucNumber $branch_zero_id"
 [ ! -f $tempCurrentBranchDataDir/dem_meters.tif ] && \
 gdalwarp -cutline $tempHucDataDir/wbd_buffered.gpkg -crop_to_cutline -ot Float32 -r bilinear -of "GTiff" \
     -overwrite -co "BLOCKXSIZE=512" -co "BLOCKYSIZE=512" -co "TILED=YES" -co "COMPRESS=LZW" \
-    -co "BIGTIFF=YES" -t_srs $huc_CRS $input_DEM $tempHucDataDir/dem_meters.tif
+    -co "BIGTIFF=YES" -t_srs $huc_CRS -tr $res $res $input_DEM $tempHucDataDir/dem_meters.tif
 
-Tcount
 
 ## GET RASTER METADATA
 echo -e $startDiv"Get DEM Metadata $hucNumber $branch_zero_id"
-read fsize ncols nrows ndv xmin ymin xmax ymax cellsize_resx cellsize_resy\
+read fsize ncols nrows ndv xmin ymin xmax ymax cellsize_resx cellsize_resy \
     <<<$($srcDir/getRasterInfoNative.py $tempHucDataDir/dem_meters.tif)
 
 ## RASTERIZE NLD MULTILINES ##