create `summarize_river.R` based on old river summary script to take a CSV as argument

[rburghol]

Desired actions for summarize_river.R

Script to use: https://github.com/HARPgroup/om/blob/master/R/summarize/waterSupplyModelNode.R

• Split into 2 separate scripts: summarize_values.R and summarize_figures.R

summarize_values.R Use: summarize_values.R [riverseg] [scenario] [hydr_path] [model_version] [json_dir(?)] Example: Rscript HARP-2022-Summer/AutomatedScripts/summarize_values.R JU4_7260_7380 hsp2_2022 /media/model/p6/out/river/hsp2_2022/hydr/JU4_7260_7380_hydrd_wy.csv cbp-6.0 /media/model/p6/out/river/hsp2_2022/json/

summarize_figures.R Use: summarize_figures.R [riverseg] [scenario] [hydr_path] [model_version] [image_dir] [json path] Example: Rscript HARP-2022-Summer/AutomatedScripts/summarize_figures.R JU4_7260_7380 hsp2_2022 /media/model/p6/out/river/hsp2_2022/hydr/JU4_7260_7380_hydrd_wy.csv cbp-6.0 /media/model/p6/out/river/hsp2_2022/images/ /media/model/p6/out/river/hsp2_2022/json/JU4_7260_7380_summ.json

To finalize this script:

• [x] Verify elfgen_huc() works correctly and shows up in VAHydro
• [x] Verify pid retrieval is correct @ end of script
• [x] Verify wd and ps plots are exported when wd and ps are present. For this, the wd and ps handling routine #541 must be run for the segment being used before using the summary script. Hydr csv with wd and demand but no ps:JL1_6770_6850_hydr.csv
• [x] Go through currently unused loop and make sure vars are correct
• [x] Implement changes for daily hydr csv into summary script once made in conversion script -- the daily csv will be the data used for the analysis (so that metrics are all correct since IHA expects daily data)
• [x] Separate back to metrics & figures -- desirable for modularity in the future Either use a df in the temp directory to store metrics needed by fig script, or keep script merged Use json object retrieval (in R) for this ^. If management of temp directory becomes difficult, json files can be stored in our export directory instead (/media/model/... Into hsp_hydr_conversion.R:
• [x] Export/save daily aggregated csv (NOT replacing original hydr file), name as daily, rename hourly hydr?

[juliabruneau]

@rburghol

• Is the goal for this script to combine all the Rscripts that we have done? (Extract csvs, convert values, analyze/push to VAHydro) OR
• Is it to export new values from ps/wd to VAHydro similarly to the hsp_hydr_analysis script?

I guess I'm just confused on what the goal for this script is even after looking at the old script linked to this issue.

[rburghol]

@juliabruneau The goal for the script is to create a version of this old script that behaves just like the new scripts you all have already made.

Essentially you'll be Making it behave exactly like this this: HARP-2022-Summer/AutomatedScripts/hsp_hydr_analysis.R

Now, one might argue that you could just move all of the things from the old summarize script to hsp_hydr_analysis.R -- and honestly I would not object to that, but we need to change the name of the script.

And I also believe that the new script that you all did may use some of the more sophisticated data source functions (adding properties with ds ) that I would like to preserve.

[glenncampagna]

Laying out Variables from watersupplyModelNode.R

• Arguments: Currently uses pid, elid, and runidwhich will be replaced by river_segment_name, scenario_name, input_file_path, image_directory_path, and model_version

• Variables exported to VAhydro: Note: some will exist, others will need placeholders (column of zeros) Note 2: names listed below are what the variables are named in the script, not their names in vahydro net_consumption_mgd = wd_cumulative_mgd - ps_cumulative_mgd wd_mgd wd_cumulative_mgd = set to 0 ps_mgd ps_cumulative_mgd = set to 0 Qout Qbaseline ps_nextdown_mgd = set to 0 consumptive_use_frac daily_consumptive_use_frac l90_Qout l90_year - year which the 90 day low flow occurs l30_Qout l30_year 7q10 ml8 mne9_10 unmet_demand_mgd - the difference between demand and withdrawal l90_Qout_trim l90_year_trim l30_Qout_trim l30_year_trim usable_pct_p0 usable_pct_p10 usable_pct_p50 remaining_days_p0 remaining_days_p10 remaining_days_p50

• Other variables: model_run_start - earliest/first date in index col of hydr file model_run_end - latest/last date in index col of hydr file wd_imp_child_mgd - the 'child' variables won't apply to hsp2 = set to 0

• Assumptions (we will implement in new script): imp_off = 1 (there will not be an impoundment) - set to 0 in modelNode script impoundment_Qin = 0 ^These will be added as columns to our hydr file before analysis (using a modified hsp_hydr_conversion.R script)

• Figures exported: fig.l90_imp_storage fig.l90_imp_storage.2yr fig.fdc.all fig.imp_storage.all elev90_imp_storage.all fig.l90_flows.2yr fig.flows.all fig.fdc fig.hydrograph_dry

[juliabruneau]

Update

We decided to base summarize_river.R heavily on how waterSupplyModelNode.R is set up, but we made some tweaks to make it clearer for us. The general steps are as follows:

1. Call in data with arguments
2. Set up data source and entity type in order to push to VAHydro
3. Perform analysis on Qout, wd, and ps 3.1 Push to VAHydro
4. Calculate metrics that require 'Zoo' (IHA) 4.1 Push to VAHydro
5. Calculate metrics with a trimmed climate change timescale 5.1 Push to VAHydro

• At this point we created another script called summarize_river_figures.R, that has all the code to generate the figures. We are still confused on this part, because all of the figures are generated within loops, and it will require us some more time to understand them/think if there is a better alternative.

All the changes so far are pushed to GitHUB, with information commented at each step.

@rburghol Could you explain why the timeseries is trimmed (step 5) to the "climate change scenario timescale"? Is there specific analysis that is done later on climate change effects on the river segments?

[glenncampagna]

1st version is ready for testing, I was stalled by a permission error. Someone else can try

Here is the example command: Rscript ~/HARParchive/HARP-2022-Summer/AutomatedScripts/summarize_river.R JL1_6770_6850 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_685_hydr.csv cbp-5.3.2 General Use: summarize_river.R [riverseg] [scenario] [hydr_path] [model_version] This script does not incorporate any figures or any 'trimmed' variables since we will have a separate figure analysis script and we are unsure about the trimmed values/metrics (for period 1990 - 2000). Are they still desired? This script can be used/tested on any hydr file that has been converted and has Qout, date cols, and recently implemented wd_mgd

[juliabruneau]

Successful Trial

I ran the river segment: PL3_5250_0001

Use: Rscript ~/HARParchive/HARP-2022-Summer/AutomatedScripts/summarize_river.R PL3_5250_0001 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/PL3_5250_0001_hydr.csv cbp-5.3.2

Results: http://deq1.bse.vt.edu/d.dh/om-model-info/6859600

• There are some values in there from running previous analysis scripts

NEXT STEPS

• Some values could be rounded/cleaned up
  • Instead of: -0.000000000000024961514336989 we should do - 2.49 10^-x
• We can go through once more the naming conventions for the pushed values to make sure they are clear/sufficient
• We should check the exported value calculations once more to check that all the units are correct

I tried including the trimmed low flow values, but when I test the if loop in R, I don't get any values (but no errors?). This is the error I got in the terminal:

Error in vahydro_post_metric_to_scenprop(model_scenario$pid, "om_class_Constant",  :
  object 'l90_Qout_trim' not found
Execution halted

This makes me believe that something within the if-loop is not working correctly.

Everything ran fine up until this point.

[juliabruneau]

Next steps

• We verified with Joey that the trimmed climate change scenario values were not needed for our analysis
• The summarize_river.R is now verified to run without errors

To do:

• [x] Make sure all values are calculated correctly
• [x] Why is l30 negative?
• [x] Modify summarize_river_figures.R
• [x] Push summarize_river_figures.R successfully to VAHydro

[juliabruneau]

Start on summarize_river_figures.R

The main loop of the figure generation script is set by imp_off:

• Since we set imp_off = 1 manually, the script that generates plots for us starts around row 380

After this point, these figures are generated:

• fig.l90_flows.2yr (zoom in on the Lowest 90 Day Flow Period)

• fig.flows.all (Qout, ps, and wd - results as zero since wd_cumulative and ps_cumulative are 0)
• fig.fdc (flow duration curve - with Qbaseline and Qout)
• fig.hydrograph_dry (plots the lowest 90 day flow period)

Reasons why some plots will not give meaningful results:

• wd_cumulative = 0
• ps_cumulative = 0
• Our Qbaseline = Qout, since Qbaseline = Qout + (wd - ps)

We set the filenames and url's, so that our next step is to test if these plots will be exported in the right place.

[glenncampagna]

• [x] Set cumulative variables equal to rates (cumulative does not refer to time in this case, hence the mgd)
• [x] Try to substitute window() function for filter() (we are supporting Rob and going anti-dplyr)
• [x] Re-run with newest hydr csv (OR1_7700_7980_hydr.csv), see how exports to VAhydro change

[juliabruneau]

Update (new script: summarize_river_values.R)

• OR1_7700_7980_hydr.csv pushed to VAHydro successfully
• Use: Rscript ~/HARParchive/HARP-2022-Summer/AutomatedScripts/summarize_river.R OR1_7700_7980 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/OR1_7700_7980_hydr.csv cbp-5.3.2
• Results: http://deq1.bse.vt.edu/d.dh/om-model-info/6853355

In order to start testing pushing the figures into VAhydro:

• A new function was created: summarize_river_values.R

  
  summarize_river_values <- function (hydr) {
  syear = as.integer(min(hydr$year), na.rm = TRUE)
  eyear = as.integer(max(hydr$year), na.rm = TRUE)
  model_run_start <- min(hydr$date)   
  model_run_end <- max(hydr$date)
  years <- seq(syear,eyear)

[...calculations]

values <- list(hydr$index, hydr$date, hydr$hour, hydr$day, hydr$month, hydr$year, hydr$Qout, hydr$Qbaseline, hydr$wd_mgd, hydr$wd_cumulative_mgd, hydr$ps_mgd, hydr$ps_cumulative_mgd, hydr$ps_nextdown_mgd, hydr$net_consumption_mgd, hydr$unmet_demand_mgd) names(values) <- c("index", "date", "hour", "day", "month", 'year', "Qout", "Qbaseline", "wd_mgd", "wd_cumulative_mgd", "ps_mgd", "ps_cumulative_mgd", "ps_nextdown_mgd", "net_consumption_mgd", "unmet_demand_mgd")

values2 <- list(l90_Qout, l90_year, l30_Qout, l30_year, imp_off) names(values2) <- c("l90_Qout", "l90_year", "l30_Qout", "l30_year","imp_off")

values3 <- list(values, values2)

return(values3) }

- This will output all the required calculations from the summarize_river script to be used to generate figures
- Output (as of now):

![image](https://user-images.githubusercontent.com/104520318/199523215-e81b26b4-29dc-44c0-9317-e41e2222ef4d.png)

 - Struggles so far was realizing that functions only return one "value", which is why the different length values are outputted as a list for now
  - Still have to figure out how to call the values from the list - should I convert them a data frame (and individual values) **OR** edit the script when calling for different values (hesitant on the latter since we don't know if this is a permanent solution)

NEXT STEPS:

- [x] Plot figures from the values that were generated with the function
- [x] Push figures to VAHydro

[glenncampagna]

Adding unmet_demand_mgd

hydr$unmet_demand_mgd = as.numeric(hydr$demand_mgd) - as.numeric(hydr$wd_mgd) We didn't previously have unmet demand as a column in our hydr table so I created this line to generate it.. is it correct?

[juliabruneau]

Update

Notable Changes:

• We were able to generate values and figures with zoo data by setting the data as ‘numeric’ afterwards
• The data was transformed into daily data in order to speed up the figure generation and to get the correct 90 and 30-day low flows
• The values and figures are now generated in the same script called: summarize_river_merge.R (latest version of changes in the script)
• Figure generation was successful

Next steps:

• Testing pushing the data into VAHydro
• We are encountering errors that haven’t shown up in our desktop R. This might be due to a different R version in the terminal.
• Continuing with this next week...

[glenncampagna]

Testing merged script (metrics and figures) on the deq server

• Testing successful in R studio but receive error in terminal: In zoo(hydr, order.by = hydr$index) : some methods for “zoo” objects do not work if the index entries in ‘order.by’ are not unique Error in rval[i, j, drop = drop., ...] : subscript out of bounds Calls: mean -> univar -> do.call -> $ -> $.zoo -> [ -> [.zoo Execution halted Note: changing the mode/class of the hydr zoo to numeric currently makes all dates/timestamps NA (index and date columns NA). This is overcome by replacing the 'index' col of the zoo with the previous index vector. This is done after the original hydr table is aggregated to daily data

• [x] Look for any duplicated index values None: > anyDuplicated(hydr$index) [1] 0 Also verified there are no duplicates in date col after aggregating, however I get the same error about non-unique values: > anyDuplicated(hydr$date) [1] 0 But get the same error when ordering the zoo by the date:

  
  In zoo(hydr, order.by = hydr$date) :
  some methods for “zoo” objects do not work if the index entries in ‘order.by’ are not unique
  Error in rval[i, j, drop = drop., ...] : subscript out of bounds
  Calls: mean -> univar -> do.call -> $ -> $.zoo -> [ -> [.zoo
  Execution halted

This does confirm that the error occurs when the original zoo for the whole hydr table is created:
`hydr <- zoo(hydr, order.by = hydr$date)`
Asking myself: why is the terminal R seeing duplicates while R studio does not?
Manual testing in terminal R results in successfully creating the zoo with the same hydr csv file ..

**Solution:** `hydr <- zoo(hydr, order.by = as.Date(hydr$Group.1))` , where `Group.1` is a date added by the aggregation
- Exports to VAhydro from testing seg OR1_7700_7980 located @ http://deq1.bse.vt.edu/d.dh/om-model-info/6853355 
  Note: 'empty' plots don't show up in VAhydro (that compared wd and ps which don't exist for this seg)

[juliabruneau]

Separating the script using JSON

The two scripts are now summarize_values.R and summarize_figures.R

• The txt file created with JSON is called: summarize_values_temp.R for now

  I tried using two different commands in summarize_values.R:

  1. toJSON()
     • didn't accept data as zoo; had to be converted to a data frame
     • then using rbind_pages() from the JSON package, data frames of different lengths can be listed and exported as 1 data frame (resulted in a frame with rows of NA's)
     • size of JSON file: 3.8 MB
  2. serializeJSON()
     • let's the user input any type of data
     • I created a list that had 2 elements (list of daily data, and list of the data with just 1 value)
     • created a txt file successfully with write(), with a size of only 1.8 MB

     • clearly the better option due to keeping the data as a zoo & the size

       In summarize_figures.R

• Using unserializeJSON() to access the data in the txt file

[glenncampagna]

• Julia and I did some work on the summary script(s) this morning and here are our updates:
  1. We are writing an updated hydr csv in the values script (1st part of summary scripts) with all columns needed by the figures script that will then be read by the figures script, but since this hydr file has been trimmed down to water year it will not replace the original hydr file (name = seg_hydr_summ.csv). This file probably won't need to be exported from our temp directory in the model
  2. We streamlined our json file to just contain the values/variables needed. For now these are imp_off and l90_year
  3. My attempt at testing gives an error I haven't seen before, and almost looks like it might be REST-related but I don't know:

     
     Rscript ~/HARParchive/HARP-2022-Summer/AutomatedScripts/summarize_values.R JL1_6770_6850 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_6850_hydr.csv cbp-6.0 /media/model/p532/out/river/hsp2_2022/hydr/
     [1] "REST AUTH INFO HAS BEEN SUPPLIED"
     [1] "RETRIEVING REST TOKEN"
     [1] "Login attempt **unsuccessful**"
     Warning message:
     no DISPLAY variable so Tk is not available 
     reading from http://deq1.bse.vt.edu/d.dh
     [1] "REST AUTH INFO HAS BEEN SUPPLIED"
     [1] "RETRIEVING REST TOKEN"
     [1] "Login attempt unsuccessful"
     Error in request(headers = c(..., .headers)) : 
     is.character(headers) is not TRUE
     Calls: <Anonymous> ... request_build -> as.request -> <Anonymous> -> request -> stopifnot
     In addition: Warning message:
     In fn_get_rest(entity_type, pk, config, self$site, private$token) :
     NAs introduced by coercion
     Execution halted

[juliabruneau]

Work Session after 12/9 Meeting

Here is the conclusion on the workflow we agreed on:

1. hsp_hydr_conversion.R
   • Generates all columns based on wd and ps (some will be moved from summarize_values.R)
   • Saves a normal hydr.csv with all additional columns
   • Trims to daily data
   • Trims to a water year with a buffer
   • Saves this into a csv: hydrd_wy.csv (?)
2. summarize_values.R
   • Performs summaries: l90, l30 …
   • Pushes values to VAhydro
   • Creates a json file with values needed for figure generation
   • Unique file name - [riverseg]_summ.json
3. summarize_figures.R
   • Calls in daily hydr file in the water year, and the json file
   • Generates and pushes figures to VAhydro

Resulting file from hsp_hydr_conversion: [riverseg]_hydrd_wy.csv

> head(hydrd)
           index               DEP       IVOL      O1 O2 O3        OVOL1 OVOL2 OVOL3     PRSUPY    RO        ROVOL     SAREA     TAU        USTAR     VOL       VOLEV      
1984-10-01 1984-10-01 11:30:00 1.2943616 5.8630603 0  0  53.200393 0     0     4.3632592 0.3510285 53.200393 4.3632592 175.80083 0.11918396 0.2479337 213.23919 0.033429171
1984-10-02 1984-10-02 11:30:00 1.3006436 4.0104952 0  0  53.633711 0     0     4.4425998 0         53.633711 4.4425998 175.90797 0.11972911 0.2485568 214.31405 0.030532922
1984-10-03 1984-10-03 11:30:00 1.2381018 3.5885918 0  0  47.862436 0     0     3.9644563 0         47.862436 3.9644563 174.84134 0.11446638 0.2430337 203.34863 0.032575537
1984-10-04 1984-10-04 11:30:00 1.1915302 3.3681695 0  0  43.587347 0     0     3.6083385 0         43.587347 3.6083385 174.04708 0.11052086 0.2388105 195.22596 0.039316895
1984-10-05 1984-10-05 11:30:00 1.1577942 3.2100393 0  0  41.062119 0     0     3.3966015 0         41.062119 3.3966015 173.47173 0.10764766 0.2356864 189.36554 0.044858515
1984-10-06 1984-10-06 11:30:00 1.1257915 3.0716329 0  0  39.394389 0     0     3.2585724 0         39.394389 3.2585724 172.92593 0.10491024 0.2326703 183.82426 0.041597132
           date       hour day month year divr_cfs diva_cfs ps_afd Qout      wd_mgd ps_mgd demand_mgd Qbaseline wd_imp_child_mgd wd_cumulative_mgd ps_cumulative_mgd ps_nextdown_mgd
1984-10-01 1984-10-01 11.5 1   10    1984 0        0        0      52.795436 0      0      0          52.795436 0                0                 0                 0              
1984-10-02 1984-10-02 11.5 2   10    1984 0        0        0      53.755458 0      0      0          53.755458 0                0                 0                 0              
1984-10-03 1984-10-03 11.5 3   10    1984 0        0        0      47.969921 0      0      0          47.969921 0                0                 0                 0              
1984-10-04 1984-10-04 11.5 4   10    1984 0        0        0      43.660896 0      0      0          43.660896 0                0                 0                 0              
1984-10-05 1984-10-05 11.5 5   10    1984 0        0        0      41.098878 0      0      0          41.098878 0                0                 0                 0              
1984-10-06 1984-10-06 11.5 6   10    1984 0        0        0      39.428726 0      0      0          39.428726 0                0                 0                 0              
> 

• This file has wd_imp_child,_mgd wd_cumulative_mgd, ps_cumulative_mgd, and ps_nextdown_mgd added to it from the summarize_values.R script
• Will double check that those are the only columns needed for the figure generation in summarize_figures.R

  • _l90year and _impoff are still exported as a json file from the summary

    Updated uses of the scripts

1. Use: hsp_hydr_conversion.R [filename] Example: Rscript hsp_hydr_conversion.R JL1_6770_6850_hydr.csv

   • Output: JL1_6770_6850_hydr.csv and JL1_6770_6850_hydrd_wy.csv

2. Use: summarize_values.R [riverseg] [scenario] [hydrd_wy.csv path] [model_version] [json_file] Example: Rscript summarize_values.R JL1_6770_6850 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_6850_hydrr_wy.csv cbp-6.0 /media/model/p532/out/river/hsp2_2022/json/

   • Output: JL1_6770_6850_summ.json and summary values pushed to VAHydro

3. Use: summarize_values.R [riverseg] [scenario] [hydrd_wy path] [model_version] [image_dir] [json_file] Example: Rscript summarize_figures.R JL1_6770_6850 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_6850_hydrd_wy.csv cbp-6.0 /media/model/p532/out/river/hsp2_2022/images/ /media/model/p532/out/river/hsp2_2022/json/JL1_6770_6850_summ.json

   • Output: figures pushed to VAHydro

     Next Steps

   • [x] Making sure the conversion script works in the terminal (couldn't test due to permission issue - Julia)
   • [x] Testing all 3 scripts back to back in the terminal; see if everything gets pushed to VAHydro
   • [ ] Testing all 3 scripts in the meta model framework

[juliabruneau]

Update - debugging

These are the working arguments for the 2 summary scripts:

summarize_values: Rscript ~/HARParchive/HARP-2022-Summer/AutomatedScripts/summarize_values.R JL1_6770_6850 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_6850_hydrd.csv cbp-5.3.2 /media/model/p532/out/river/hsp2_2022/hydr/

• Ran successfully

summarize_figures: Rscript ~/HARParchive/HARP-2022-Summer/AutomatedScripts/summarize_figures.R JL1_6770_6850 hsp2_2022 /media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_6850_hydrd.csv cbp-5.3.2 /media/model/p532/out/river/hsp2_2022/images/ /media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_6850_summ.json

• Read the json file successfully - has to be run in the directory: /media/model/p532/out/river/hsp2_2022/**hydr**

In order to read the json file, the file name had to be extracted. It is in a format so that it will select the last element of the list - this way it will work universally.

json_split <- strsplit(json_dir, split = '/')
last_element <- as.numeric(length(json_split[[1]]))
json_file <- json_split[[1]][[last_element]]  # selecting just the json file name
        Selects the last element of the list
values <- unserializeJSON(readLines(json_file))

New Error with figure script:

• This is the error code I'm receiving at the moment, but will fix later (stopping for now):

  Plotting critical flow periods
  Error in attr(x, "tsp") <- c(1, NROW(x), 1) :
  invalid time series parameters specified
  Calls: window -> window.default -> hasTsp
  Execution halted

[juliabruneau]

Debugging

Plotting critical flow periods Error in attr(x, "tsp") <- c(1, NROW(x), 1) : invalid time series parameters specified Calls: window -> window.default -> hasTsp Execution halted

The error occurs when trying to start the if-loop that will generate the figures.

message("Plotting critical flow periods")
# does this have an active impoundment sub-comp
if (imp_off == 0) {

imp_off is read from the new json file, which is a list:

imp_off <- as.numeric(values[[1]])
l90_year <- as.numeric(values[[2]])

• The error message seems to not recognize the value 1 of the imp_off?
• I tried changing the as.numeric to as.factor, but received the same error

When looking for the attributes of imp_off:

> attributes(imp_off)
NULL

• It seems like the json file is being extracted properly:

  
  > json_dir <- '/media/model/p532/out/river/hsp2_2022/hydr/JL1_6770_6850_summ.json'
  > json_split <- strsplit(json_dir, split = '/')
  > last_element <- as.numeric(length(json_split[[1]]))
  > json_file <- json_split[[1]][[last_element]]
  > json_file
  [1] "JL1_6770_6850_summ.json"
  > values <- unserializeJSON(readLines(json_file))
  > values
  $imp_off
  [1] 1

$l90_year [1] 2019

imp_off <- as.numeric(values[[1]]) l90_year <- as.numeric(values[[2]]) imp_off [1] 1

### Update

- I don't think the error is in recognizing the imp_off as "1", since my test in the terminal R worked:

if(imp_off == 0){

• print("zero")
• } else {
• print("one")
• } [1] "one"

### Found the location of the error:

hydrpd <- window(hydr, start = pdstart, end = pdend)

- This is causing the error
- Fixed by setting the data as 'zoo' again. For some reason, the data read from the directory is not automatically kept as a zoo.

[juliabruneau]

Error with ylim and cbind

When testing in local and terminal R, there are no errors. But when running the R script with arguments, there is an error with setting the ylim with cbind:

Error in plot.window(...) : need finite 'ylim' values
Calls: plot -> plot -> plot.default -> localWindow -> plot.window
In addition: Warning message:
In .Primitive("max")(numeric(0), na.rm = TRUE) :
  no non-missing arguments to max; returning -Inf
Execution halted

The line of script where the error occurs:

ymx <- as.numeric(max(cbind(hydrpd$Qbaseline, hydrpd$Qout)), na.rm = TRUE)
  plot(
    as.numeric(hydrpd$Qbaseline), ylim = c(0,ymx), xlim=c(xmn,xmx),  #Placeholders for xlim, come back to this and create xlim based on hydrpd
    ylab="Flow/WD/PS (cfs)",
    xlab=paste("Lowest 90 Day Flow Period",pdstart,"to",pdend)
  )

• Have tried debugging by making sure everything with the hydr dataset is like it was before - a zoo and set as a 'window'

[juliabruneau]

Progress Update

• All 3 scripts work when tested in the terminal
• River segment: JL1_6770_6850 that has withdrawals data was used in testing, and verified that all plots work

• the elfgen_huc function has not worked in the terminal yet (I found a document of the elfgen package, but this particular function is not mentioned in it)

  • Next steps include testing these scripts in the meta-model, and using the graphs in the summary markdown document.

[glenncampagna]

Test

Rscript HARP-2022-Summer/AutomatedScripts/summarize_values.R JU3_6650_7300 hsp2_2022 /media/model/p6/out/river/hsp2_2022/hydr/JU3_6650_7300_hydrd_wy.csv cbp-6.0 /media/model/p6/out/river/hsp2_2022/hydr/

Rscript HARP-2022-Summer/AutomatedScripts/summarize_figures.R JU3_6650_7300 hsp2_2022 /media/model/p6/out/river/hsp2_2022/hydr/JU3_6650_7300_hydrd_wy.csv cbp-6.0 /media/model/p6/out/river/hsp2_2022/images/ /media/model/p6/out/river/hsp2_2022/hydr/JU3_6650_7300_summ.json

Output:

Error in file(con, "r") : cannot open the connection
Calls: unserializeJSON -> unpack -> parseJSON -> readLines -> file
In addition: Warning message:
In file(con, "r") :
  cannot open file 'JU3_6650_7300_summ.json': No such file or directory
Execution halted

• Error means that no json file exists OR we aren't allowed to read it with unserializeJSON()

  Debugging

  While testing in local R - no errors occur, and json file is created/extracted fine.

Testing in terminal R (summarize_values.R):

river_seg <- 'JL1_6770_6850'
scenario_name <- 'hsp2_2022'
hydr_file_path <- '/media/model/p6/out/river/hsp2_2022/hydr/JL1_6770_6850_hydrd_wy.csv'
model_version <- 'cbp-6.0'
json_dir <- '/media/model/p6/out/river/hsp2_2022/hydr/'

Turns out to be a permission issue = we aren't allowed to write the json file into the folder:

> write(values_json, file= paste0(json_dir, river_seg, "_summ.json"), sep = ",")
Error in file(file, ifelse(append, "a", "w")) :
  cannot open the connection
In addition: Warning message:
In file(file, ifelse(append, "a", "w")) :
  cannot open file '/media/model/p6/out/river/hsp2_2022/hydr/JL1_6770_6850_summ.json': Permission denied

If @glenncampagna wants to test this to see if it is a permission error for him as well, it should verify the source of this error.

[jdkleiner]

Spent time looking into the elfgen_huc() item.

• elfgen_huc() (called at the bottom of summarize_figures.R)
  • First ensure you're loading the file rseg_elfgen.R at the top of your script, as this conatins the function definition for elfgen_huc()
    • source(paste(om_location,'R/summarize','rseg_elfgen.R',sep='/'))
  • To be able to generate the elfgen analysis, there must be a consumptive_use_frac metric on the model scenario container
    • This property should be getting generated in summarize_values.R but it currently resolves to NA
    • Note that consumptive_use_frac <- 1.0 - (Qout / Qbaseline)
  • This issue here lies in the calculation of Qbaseline
  • Refer back to the calculation in waterSupplyModelNode.R https://github.com/HARPgroup/om/blob/master/R/summarize/waterSupplyModelNode.R#L134
    • We should calculate it this same way if possible.

      dat$Qbaseline <- dat$Qout +
      (dat$wd_cumulative_mgd - dat$ps_cumulative_mgd ) * 1.547

    • So where is the hangup?
      • The hydr file is all 0s for wd_cumulative_mgd?
      • The hydr file doesnt have a column for ps_cumulative_mgd?
• I think the following should still work even if your scripts are not successfully setting a consumptive_use_frac value
  • Though richness_change_abs and richness_change_pct will be 0

    
    source(paste(om_location,'R/summarize','rseg_elfgen.R',sep='/'))

###############################################

RSEG ELFGEN

###############################################

GET RSEG HYDROID FROM RSEG MODEL PID

rseg <- RomProperty$new(ds, list(pid=model$pid), TRUE) rseg_hydroid<-rseg$featureid

huc_level <- 'huc8' dataset <- 'VAHydro-EDAS'

elfgen_huc(scenario_name, rseg_hydroid, huc_level, dataset, model_scenario, ds, image_dir, save_url, site)

- result will look something like: 
   - http://deq1.bse.vt.edu:81/d.dh/om-model-info/6972498
![image](https://user-images.githubusercontent.com/29379385/218560257-e3aad9d8-df0a-446e-ad6a-62dc064972aa.png)

[glenncampagna]

Thanks for getting back to us @jdkleiner. Errors related to plotting for us were fixed when the NAs for Qbaseline were fixed and now I can run the script till the end. I was able to get the elfgen container to show up in vahydro but did get an error at the very end of the script and the elfgen container was empty The error:

Error in (function (..., row.names = NULL, check.rows = FALSE, check.names = TRUE,  : 
  arguments imply differing number of rows: 1, 0
Calls: elfgen_huc ... as.data.frame -> as.data.frame.list -> do.call -> <Anonymous>

[jdkleiner]

@glenncampagna Nice progress with those fixes! What segment id are you testing with that results in the error? I'll try to replicate on my end

[glenncampagna]

@jdkleiner I used JU4_7260_7380

[jdkleiner]

@glenncampagna Thanks!

• It runs without a hitch locally for me (check it out: http://deq1.bse.vt.edu:81/d.dh/om-model-info/6973111)
• Testing running from the server now, but running into that same error you'd been having (even though we know that json file is indeed there). Did you have a way of getting around this error?

  
  jkleiner@deq2:/opt/model/HARParchive$  Rscript HARP-2022-Summer/AutomatedScripts/summarize_figures.R JU4_7260_7380 hsp2_2022 /media/model/p6/out/river/hsp2_2022/hydr/JU4_7260_7380_hydrd_wy.csv cbp-6.0 /media/model/p6/out/river/hsp2_2022/images/ /media/model/p6/out/river/hsp2_2022/hydr/JU4_7260_7380_summ.json

Error in file(con, "r") : cannot open the connection Calls: unserializeJSON -> unpack -> parseJSON -> readLines -> file In addition: Warning message: In file(con, "r") : cannot open file 'JU4_7260_7380_summ.json': No such file or directory Execution halted

[glenncampagna]

@jdkleiner that is a path-related error that comes from a bug that I haven't quite had the chance to fix today, and it can be overcome by running the summarize_figures.R script from the directory where the json file is located (should be [scenario]/json/)

• [ ] fix path-related error in summarize_figures.R

[jdkleiner]

@glenncampagna Thanks, that's just the workaround I needed to get it running.

• After testing your scripts, this ended up being a simple case of needing to pull the latest copy of the HARParchive master to deq2. This was needed to implement the following fix:
  • old: elfgen_huc(scenario_name, rseg_hydroid, huc_level, hydraset, scenprop, ds, image_dir, save_url, site)
  • new: elfgen_huc(scenario_name, rseg_hydroid, huc_level, Dataset, model_scenario, ds, image_dir, save_url, site)
• I ran it for your test id JU4_7260_7380 and looks successful to me, but would be good if you could verify: http://deq1.bse.vt.edu:81/d.dh/om-model-info/6973111

[jdkleiner]

Success using jdk account:

jdk@deq2:/opt/model/HARParchive$ Rscript HARP-2022-Summer/AutomatedScripts/summarize_figures.R JU4_7260_7380 hsp2_2022 /media/model/p6/out/river/hsp2_2022/hydr/JU4_7260_7380_hydrd_wy.csv cbp-6.0 /media/model/p6/out/river/hsp2_2022/images/ /media/model/p6/out/river/hsp2_2022/json/JU4_7260_7380_summ.json

[glenncampagna]

Moving/renaming final version(s)

summarize_values.R moved into hsp_hydr_analysis.R summarize_figures.R renamed as hsp_hydr_plots.R