Postprocessing¶
Merge RAPID Output¶
-
class
RAPIDpy.postprocess.
ConvertRAPIDOutputToCF
(rapid_output_file, start_datetime, time_step, qinit_file='', comid_lat_lon_z_file='', rapid_connect_file='', project_name='Default RAPID Project', output_id_dim_name='rivid', output_flow_var_name='Qout', print_debug=False)[source]¶ Class to convert RAPID output to be CF compliant. You can also use this to combine consecutive RAPID output files into one file.
-
rapid_output_file
¶ str or list – Path to a single RAPID Qout file or a list of RAPID Qout files.
-
start_datetime
¶ datetime – Datetime object with the time of the start of the simulation.
-
time_step
¶ int or list – Time step of simulation in seconds if single Qout file or a list of time steps corresponding to each Qout file in the rapid_output_file.
-
qinit_file
¶ Optional[str] – Path to the Qinit file for the simulation. If used, it will use the values in the file for the flow at simulation time zero.
-
comid_lat_lon_z_file
¶ Optional[str] – Path to comid_lat_lon_z file. If included, the spatial information will be added to the output NetCDF file.
-
rapid_connect_file
¶ Optional[str] – Path to RAPID connect file. This is required if qinit_file is added.
-
project_name
¶ Optional[str] – Name of your project in the output file. Default is “Default RAPID Project”.
-
output_id_dim_name
¶ Optional[str] – Name of the output river ID dimension name. Default is ‘rivid’.
-
output_flow_var_name
¶ Optional[str] – Name of streamflow variable in output file, typically ‘Qout’ or ‘m3_riv’. Default is ‘Qout’.
-
print_debug
¶ Optional[bool] – If True, the debug output will be printed to the console. Default is False.
Warning
This code replaces the first file with the combined output and deletes the second file. BACK UP YOUR FILES!!!!
Example:
import datetime from RAPIDpy.postprocess import ConvertRAPIDOutputToCF file1 = "/path/to/Qout_1980to1981.nc" file2 = "/path/to/Qout_1981to1982.nc" cv = ConvertRAPIDOutputToCF(rapid_output_file=[file1, file2], start_datetime=datetime.datetime(2005,1,1), time_step=[3*3600, 3*3600], project_name="NLDAS(VIC)-RAPID historical flows by US Army ERDC", ) cv.convert()
-
Generate qinit from past qout¶
RAPIDpy also creates a qinit file from a RAPID qout file. This example shows how.
-
RAPID.
generate_qinit_from_past_qout
(qinit_file, time_index=-1, out_datetime=None)[source] Generate qinit from a RAPID qout file
Parameters: - qinit_file (str) – Path to output qinit_file.
- time_index (Optional[int]) – Index of simulation to generate initial flow file. Default is the end.
- out_datetime (Optional[datetime]) – Datetime object containing time of initialization.
Example:
from RAPIDpy import RAPID rapid_manager = RAPID(Qout_file='/output_mississippi-nfie/Qout_k2v1_2005to2009.nc', rapid_connect_file='/input_mississippi_nfie/rapid_connect_ECMWF.csv' ) rapid_manager.generate_qinit_from_past_qout(qinit_file='/input_mississippi_nfie/Qinit_2008_flood.csv', time_index=10162)
Generate seasonal qinit from past qout¶
-
RAPID.
generate_seasonal_intitialization
(qinit_file, datetime_start_initialization=datetime.datetime(2017, 7, 6, 18, 35, 48, 418889))[source] This creates a seasonal qinit file from a RAPID qout file. This requires a simulation Qout file with a longer time period of record and to be CF compliant. It takes the average of the current date +- 3 days and goes back as far as possible.
Parameters: - qinit_file (str) – Path to output qinit_file.
- datetime_start_initialization (Optional[datetime]) – Datetime object with date of simulation to go back through the years and get a running average to generate streamflow initialization. Default is utcnow.
This example shows how to use it:
from RAPIDpy.rapid import RAPID rapid_manager = RAPID(Qout_file='/output_mississippi-nfie/Qout_2000to2015.nc', rapid_connect_file='/input_mississippi_nfie/rapid_connect_ECMWF.csv' ) rapid_manager.generate_seasonal_intitialization(qinit_file='/input_mississippi_nfie/Qinit_seasonal_avg_jan_1.csv')
Goodness of Fit¶
To check how well your simulation performed versus observations, these functions can help you.
-
RAPIDpy.postprocess.
find_goodness_of_fit_csv
(observed_simulated_file, out_file=None)[source]¶ Finds the goodness of fit comparing observed and simulated flows In the file, the first column is the observed flows and the second column is the simulated flows.
Example:
33.5, 77.2 34.7, 73.0
Parameters: - observed_simulated_file (str) – Path to the csv file with the observed and simulated flows.
- out_file (Optional[str]) – Path to output file. If not provided, it will print to console.
Example:
from RAPIDpy.postprocess import find_goodness_of_fit_csv find_goodness_of_fit_csv('/united_kingdom-thames/flows_kingston_gage_noah.csv')
-
RAPIDpy.postprocess.
find_goodness_of_fit
(rapid_qout_file, reach_id_file, observed_file, out_analysis_file, daily=False, steps_per_group=1)[source]¶ Finds the goodness of fit comparing observed streamflow in a rapid Qout file with simulated flows in a csv file.
Parameters: - rapid_qout_file (str) – Path to the RAPID Qout file.
- reach_id_file (str) – Path to file with river reach ID’s associate with the RAPID Qout file. It is in the format of the RAPID observed flows reach ID file.
- observed_file (str) – Path to input csv with with observed flows corresponding to the RAPID Qout. It is in the format of the RAPID observed flows file.
- out_analysis_file (str) – Path to the analysis output csv file.
- daily (Optional[bool]) – If True and the file is CF-Compliant, it will compare the observed_file with daily average flow from Qout. Default is False.
Example with CF-Compliant RAPID Qout file:
import os from RAPIDpy.postprocess import find_goodness_of_fit INPUT_DATA_PATH = '/path/to/data' reach_id_file = os.path.join(INPUT_DATA_PATH, 'obs_reach_id.csv') observed_file = os.path.join(INPUT_DATA_PATH, 'obs_flow.csv') cf_input_qout_file = os.path.join(COMPARE_DATA_PATH, 'Qout_nasa_lis_3hr_20020830_CF.nc') cf_out_analysis_file = os.path.join(OUTPUT_DATA_PATH, 'cf_goodness_of_fit_results-daily.csv') find_goodness_of_fit(cf_input_qout_file, reach_id_file, observed_file, cf_out_analysis_file, daily=True)