USGS-R / drb-gw-hw-model-prep

Code repo to prepare groundwater and headwater-related datasets for modeling river temperature in the Delaware River Basin
Other
0 stars 3 forks source link

Fetch coarse sediment #36

Closed msleckman closed 2 years ago

msleckman commented 2 years ago

This PR adds a target p2_coarse_sediment_area_reaches_along_nhm description th relative area of coarse statified sediment surface material by reach buffers.

head(p2_coarse_sediment_area_reaches_along_nhm)

# A tibble: 6 × 5
  PRMS_segid cs_area_km2 total_reach_buffer_area_km2 cs_area_proportion                     geometry
  <chr>           [km^2]                      [km^2]              <dbl>                <POLYGON [m]>
1 1_1              0                            6.94              0     ((1739484 2336890, 1739490 …
2 10_1             0                            1.24              0     ((1712258 2309595, 1712296 …
3 11_1             0                            1.10              0     ((1730997 2318387, 1731006 …
4 111_1            1.48                         5.48              0.269 ((1730444 2347353, 1730488 …
5 112_1            0.199                       14.4               0.014 ((1713217 2332271, 1713220 …
6 113_1            0                            1.07              0     ((1713687 2331813, 1713687 

here is a snapshot of the buffered reaches, with the legend discribing coarse sediment proportion value within give nhm reach buffer:

image image

I added 1 function for this target - coarse_sediment_area_calc() in 2_process/src/coarse_stratified_sediment_processing.R

Additional comments I've added in line.

Please review and let me know if this is an adequate output for this input gw variable. Initially I expected to calc an area weighted mean but there are values attributed to the source geospatial dataset from Soller et al., 2009, simple an deptiction of relative area within our AOI.

msleckman commented 2 years ago

Ah, yeah, that makes sense! Either way works, so if it would take much time to add the area to p2_buffered_nhd_reaches_along_nhm, it's OK to leave here

@lekoenig I went ahead and moved the total_reach_buffer_area_km2 col creation to the p2_buffered_nhd_reaches_along_nhm target in https://github.com/USGS-R/drb-gw-hw-model-prep/pull/36/commits/12e0e8db686383fe4f4b5fe92944aaa7be304611

msleckman commented 2 years ago

I really like how you've used units here! I often do things like st_area(x)/(1000*1000) to get km2 and this approach to explicitly define the units is great - I'll have to remember for next time 🤩 (and FWIW the use of units:: notation can often help make the code more interpretable/reproducible IMO, even if the library is loaded within the pipeline.)

There is also an argument to be made that it's a bit annoying to have these units displayed. and instead keeping track on the user end and with colnames that are explicit (e.g. I write _km2) at end of colname. It's also confusing because when you view the df in Rstudio: View(p2_coarse_sediment_area_reaches_along_nhm %>% head), it shows teh tibble in a df and the units in the cell. image

However, in the console the tibble shows those units where they ought to be seen below the colname:

> p2_coarse_sediment_area_reaches_along_nhm %>% head
Simple feature collection with 6 features and 4 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 1712174 ymin: 2309159 xmax: 1743545 ymax: 2348925
Projected CRS: NAD83 / Conus Albers
# A tibble: 6 × 5
  PRMS_segid total_reach_buffer_area_km2 cs_area_km2 cs_area_proportion                                       geometry
  <chr>                           [km^2]      [km^2]              <dbl>                                  <POLYGON [m]>
1 1_1                               6.94       0                  0     ((1739484 2336890, 1739490 2336919, 1739498 2…
2 10_1                              1.24       0                  0     ((1712258 2309595, 1712296 2309652, 1712303 2…

and if I want to pull just hte value from the tibble, unit is shown in the output, and it's considered a obj of class 'units' however I can run any calculation on it as a simple numeric value.

> p2_coarse_sediment_area_reaches_along_nhm[[4,3]]
1.477084 [km^2]
> class(p2_coarse_sediment_area_reaches_along_nhm[[4,3]])
[1] "units"
> p2_coarse_sediment_area_reaches_along_nhm[[4,3]]*150
221.5626 [km^2]

Some things to keep in mind for downstream processing.

msleckman commented 2 years ago

@lekoenig I've addressed all your suggestions. Realize I did not request approval for merge. But let me know if these changes work for you and if I can merge!