HARPgroup / cbp_wsm

1 stars 0 forks source link

Land Use Subdivide for subwatershed creation #96

Open rburghol opened 1 year ago

rburghol commented 1 year ago

Will need a function that takes 1 to n source land/river seg entries in a CSV and subtracts the area ratio proportion from the receiving watershed, and creates new lines for each land/riverseg combo using the new riverseg, and each land seg flowing into the receiving watershed.

Code 1: Small alteration based on code in pull request: https://github.com/HARPgroup/HARParchive/pull/719 . This preserves the land segment names, so there can be no confusion on that front.

new_landuses <- receiving_landuses
new_landuses[-1:-2] <- new_landuses[-1:-2] * propor
new_landuses$riverseg <- subshed

Code 2: Update recieving segment LU.

receiving_landuses[-1:-2] <- receiving_landuses[-1:-2] * (1.0 - propor)
rburghol commented 1 year ago

@megpritch @juliabruneau Super exciting outcome on this script. I did a couple tweaks in #723, but nothing big.

Better new: Tried this same script on another file that had multiple LRsegs in it, config/catalog/geo/vahydro/land_water_area.csv, and it friggin worked out of the box! The only hitch was that land_water_area.csv has funny column headers with parentheses "(" and the 3-asterisks at the end.
Those get stripped and replaced with ".", so the header goes from: riverseg,landseg,area(acres)*** to riverseg,landseg,area.acres....

Did you all figure out a clean way to restore the files header line that would work without knowing what the header line was supposed to be?

If we can figure out the header glitch, this script should work for any land/river data set that needs to be proportioned, so a lot of work gets done for us...

rburghol commented 1 year ago

Testing water balance:

basin=JL1_6560_6440
segs=`cbp get_riversegs $basin`
scenario="p6ss"
model="hspf_cbp6"

for i in $segs; do
  deps=`cbp mm_river_deps $scenario $i`
  echo "Segment $i Found deps= $deps"
  job_name=`cbp mm_job_name $scenario $i`
  echo "sbatch --job-name=$job_name $deps /opt/model/meta_model/run_model $model $scenario $i  auto river"
  sbatch --job-name=$job_name $deps /opt/model/meta_model/run_model $model $scenario $i auto river
done
rburghol commented 1 year ago
# use R to query and sum land use 
library("data.table")
library("sqldf")

lu <- read.csv("/opt/model/p6/vadeq/input/scenario/river/land_use/land_use_2013VAHYDRO2018615.csv")
# compare to:
# lu <- read.csv("/opt/model/p6/vadeq/input/scenario/river/land_use/land_use_2013BASE20180615.csv")
lwa <- read.csv("/opt/model/p6/vadeq/config/catalog/geo/vahydro/land_water_area.csv")

# isolate watersheds
lu_mr <- sqldf("select * from lu where riverseg = 'JL1_6560_6440'")
lu_bc <- sqldf("select * from lu where riverseg = 'JL1_6562_6560'")
lwa_mr <- sqldf("select * from lwa where riverseg = 'JL1_6560_6440'")
lwa_bc <- sqldf("select * from lwa where riverseg = 'JL1_6562_6560'")

# verify LU just the acre columns:
sum(lu_mr[1,3:length(lu_mr[1,])])
sum(lu_bc[1,3:length(lu_bc[1,])])

# verify lwa
sum(lwa_mr[1,3:length(lwa_mr[1,])])
sum(lwa_bc[1,3:length(lwa_bc[1,])])
rburghol commented 1 year ago

Check the WDMs for eos and stream

rburghol commented 1 year ago

Hourly files:

rburghol commented 1 year ago