<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Retrieve NetCDF and model gridded climate time-series for a watershed
Case study
Step1:  Establish a secure connection with HydroShare by instantiating the hydroshare class that is defined within hs_utils. In addition to connecting with HydroShare, this command also sets and prints environment variables for several parameters that will be useful for saving work back to HydroShare.
Step4:  If you are curious about where the data is being downloaded, click on the Jupyter Notebook dashboard icon to return to the File System view.  The homedir directory location printed above is where you can find the data and contents you will download to a HydroShare JupyterHub server.  At the end of this work session, you can migrate this data to the HydroShare iRods server as a Generic Resource. 
2. Get list of gridded climate points for the watershed
For visualization purposes, we will also remap the study site shapefile, which is stored in HydroShare at the following url
Step5:  Summarize the file availability from each watershed mapping file
Step6:  3.  Compare Hydrometeorology
This section performs computations and generates plots of the Livneh 2013 and Salathe 2014 mean temperature and mean total monthly precipitation in order to compare them with each other. The generated plots are automatically downloaded and saved as .png files within the "homedir" directory.
Let's compare the Livneh 2013 and Salathe 2014 using the period of overlapping history.
Step7:  NetCDF retrieval and clipping to a spatial extent
The function get_x_dailywrf_salathe2014 retrieves and clips NetCDF files archived within the UW Rocinante NNRP repository. This archive contains daily data from January 1970 through December 1979 (10 years). Each netcdf file is comprised of meteorologic and VIC hydrologic outputs for a calendar month. The expected number of files would be 360 files (12 months for 30 years). 
In the code chunk below, 40 parallel workers will be initialized to distribute file retrieval and spatial clipping tasks. For each worker, they will wget the requested file, clip the netcdf file to gridded cell centroids within the the provided bounding box, then return the location of the spatially clipped output files.
Provide the home and subdirectory where the cropped NetCDF files will be stored. Also provide the spatial bounds (in WGS84) to crop the NetCDF files upon download. Finally, provide the number of workers to carry out the download tasks, and the start and end date of the files of interest.
Step8:  Convert collection of NetCDF files into a collection of ASCII files
Provide the home and subdirectory where the ASCII files will be stored, the source_directory of netCDF files, and the mapping file to which the resulting ASCII files will be cataloged. Also, provide the Pandas Datetime code for the frequency of the time steps. Finally, provide the catalog label that will be used for the mapping file catalog and the metadata file label.
Step9:  Create a dictionary of climate variables for the long-term mean (ltm).
INPUT
Step10:  Compute the total monthly and yearly precipitation, as well as the mean values across time and across stations
INPUT
Step11:  Visualize the "average monthly total precipitations"
INPUT
Step12:  generate a raster
Step13:  Higher resolution children of gridded cells
get data from Lower resolution parent grid cells to the children
Step14:  Visualize the "average monthly total precipitation"
5. Save the results back into HydroShare
<a name="creation"></a>
Using the hs_utils library, the results of the Geoprocessing steps above can be saved back into HydroShare.  First, define all of the required metadata for resource creation, i.e. title, abstract, keywords, content files.  In addition, we must define the type of resource that will be created, in this case genericresource.  
Note
<ASSISTANT_TASK:>
Python Code:
#!conda install -c conda-forge ogh libgdal gdal pygraphviz ncurses matplotlib=2.2.3 --yes
# silencing warning
import warnings
warnings.filterwarnings("ignore")
# data processing
import os
import pandas as pd, numpy as np, dask
# data migration library
import ogh
import ogh_xarray_landlab as oxl
from utilities import hydroshare
from ecohydrology_model_functions import run_ecohydrology_model, plot_results
from landlab import imshow_grid, CLOSED_BOUNDARY
# modeling input params
InputFile = os.path.join(os.getcwd(),'ecohyd_inputs.yaml')
# plotting and shape libraries
import matplotlib.pyplot as plt
%matplotlib inline
InputFile
# initialize ogh_meta
meta_file = dict(ogh.ogh_meta())
sorted(meta_file.keys())
Explanation: Retrieve NetCDF and model gridded climate time-series for a watershed
Case study:  the Sauk-Suiattle Watershed
<img src="http://www.sauk-suiattle.com/images/Elliott.jpg" 
style="float:right;width:150px;padding:20px">
Use this Jupyter Notebook to:
1. HydroShare setup and preparation
2. Re-establish the paths to the mapping file
3. Compute daily, monthly, and annual temperature and precipitation statistics
4. Visualize precipitation results relative to the forcing data
5. Visualize the time-series trends
6. Save results back into HydroShare
<br/><br/><br/>
<img src="https://www.washington.edu/brand/files/2014/09/W-Logo_Purple_Hex.png"
style="float:right;width:150px;padding:20px">
<br/><br/>
This data is compiled to digitally observe the watersheds, powered by HydroShare. <br/>Provided by the Watershed Dynamics Group, Dept. of Civil and Environmental Engineering, University of Washington
1.  Prepare HydroShare Setup and Preparation
To run this notebook, we must import several libaries. These are listed in order of 1) Python standard libraries, 2) hs_utils library provides functions for interacting with HydroShare, including resource querying, dowloading and creation, and 3) the observatory_gridded_hydromet library that is downloaded with this notebook.
End of explanation
notebookdir = os.getcwd()
hs=hydroshare.hydroshare()
homedir = hs.getContentPath(os.environ["HS_RES_ID"])
os.chdir(homedir)
Explanation: Establish a secure connection with HydroShare by instantiating the hydroshare class that is defined within hs_utils. In addition to connecting with HydroShare, this command also sets and prints environment variables for several parameters that will be useful for saving work back to HydroShare.
End of explanation
1/16-degree Gridded cell centroids
# List of available data
hs.getResourceFromHydroShare('ef2d82bf960144b4bfb1bae6242bcc7f')
NAmer = hs.content['NAmer_dem_list.shp']
Sauk
# Watershed extent
hs.getResourceFromHydroShare('c532e0578e974201a0bc40a37ef2d284')
sauk = hs.content['wbdhub12_17110006_WGS84_Basin.shp']
# reproject the shapefile into WGS84
ogh.reprojShapefile(sourcepath=sauk)
Explanation: If you are curious about where the data is being downloaded, click on the Jupyter Notebook dashboard icon to return to the File System view.  The homedir directory location printed above is where you can find the data and contents you will download to a HydroShare JupyterHub server.  At the end of this work session, you can migrate this data to the HydroShare iRods server as a Generic Resource. 
2. Get list of gridded climate points for the watershed
For visualization purposes, we will also remap the study site shapefile, which is stored in HydroShare at the following url: https://www.hydroshare.org/resource/c532e0578e974201a0bc40a37ef2d284/. Since the shapefile was previously migrated, we can select 'N' for no overwriting.
In the usecase1 notebook, the treatgeoself function identified the gridded cell centroid coordinates that overlap with our study site. These coordinates were documented within the mapping file, which will be remapped here. In the usecase2 notebook, the downloaded files were cataloged within the mapping file, so we will use the mappingfileSummary function to characterize the files available for Sauk-Suiattle for each gridded data product.
End of explanation
%%time
# map the mappingfiles from usecase1
mappingfile1=ogh.treatgeoself(shapefile=sauk, NAmer=NAmer, buffer_distance=0.06,
                              mappingfile=os.path.join(homedir,'Sauk_mappingfile.csv'))
Explanation: Summarize the file availability from each watershed mapping file
End of explanation
help(ogh.getDailyMET_livneh2013)
help(oxl.get_x_dailymet_Livneh2013_raw)
Explanation: 3.  Compare Hydrometeorology
This section performs computations and generates plots of the Livneh 2013 and Salathe 2014 mean temperature and mean total monthly precipitation in order to compare them with each other. The generated plots are automatically downloaded and saved as .png files within the "homedir" directory.
Let's compare the Livneh 2013 and Salathe 2014 using the period of overlapping history.
End of explanation
maptable, nstations = ogh.mappingfileToDF(mappingfile1)
spatialbounds = {'minx':maptable.LONG_.min(), 'maxx':maptable.LONG_.max(),
                 'miny':maptable.LAT.min(), 'maxy':maptable.LAT.max()}
outputfiles = oxl.get_x_dailymet_Livneh2013_raw(homedir=homedir,
                                                subdir='livneh2013/Daily_MET_1970_1970/raw_netcdf',
                                                spatialbounds=spatialbounds,
                                                nworkers=6,
                                                start_date='1970-01-01', end_date='1970-12-31')
Explanation: NetCDF retrieval and clipping to a spatial extent
The function get_x_dailywrf_salathe2014 retrieves and clips NetCDF files archived within the UW Rocinante NNRP repository. This archive contains daily data from January 1970 through December 1979 (10 years). Each netcdf file is comprised of meteorologic and VIC hydrologic outputs for a calendar month. The expected number of files would be 360 files (12 months for 30 years). 
In the code chunk below, 40 parallel workers will be initialized to distribute file retrieval and spatial clipping tasks. For each worker, they will wget the requested file, clip the netcdf file to gridded cell centroids within the the provided bounding box, then return the location of the spatially clipped output files.
Provide the home and subdirectory where the cropped NetCDF files will be stored. Also provide the spatial bounds (in WGS84) to crop the NetCDF files upon download. Finally, provide the number of workers to carry out the download tasks, and the start and end date of the files of interest.
End of explanation
%%time
# convert the netCDF files into daily ascii time-series files for each gridded location
outfilelist = oxl.netcdf_to_ascii(homedir=homedir, 
                                  subdir='livneh2013/Daily_MET_1970_1970/raw_ascii', 
                                  source_directory=os.path.join(homedir, 'livneh2013/Daily_MET_1970_1970/raw_netcdf'),
                                  mappingfile=mappingfile1,
                                  temporal_resolution='D',
                                  meta_file=meta_file,
                                  catalog_label='sp_dailymet_livneh_1970_1970')
t1 = ogh.mappingfileSummary(listofmappingfiles = [mappingfile1], 
                            listofwatershednames = ['Sauk-Suiattle river'],
                            meta_file=meta_file)
t1
# Save the metadata
ogh.saveDictOfDf(dictionaryObject=meta_file, outfilepath='test.json')
Explanation: Convert collection of NetCDF files into a collection of ASCII files
Provide the home and subdirectory where the ASCII files will be stored, the source_directory of netCDF files, and the mapping file to which the resulting ASCII files will be cataloged. Also, provide the Pandas Datetime code for the frequency of the time steps. Finally, provide the catalog label that will be used for the mapping file catalog and the metadata file label.
End of explanation
meta_file['sp_dailymet_livneh_1970_1970']['variable_list']
%%time
ltm = ogh.gridclim_dict(mappingfile=mappingfile1,
                        metadata=meta_file,
                        dataset='sp_dailymet_livneh_1970_1970',
                        variable_list=['Prec','Tmax','Tmin'])
sorted(ltm.keys())
Explanation: Create a dictionary of climate variables for the long-term mean (ltm).
INPUT: gridded meteorology ASCII files located from the Sauk-Suiattle Mapping file. The inputs to gridclim_dict() include the folder location and name of the hydrometeorology data, the file start and end, the analysis start and end, and the elevation band to be included in the analsyis (max and min elevation). <br/>OUTPUT: dictionary of dataframes where rows are temporal summaries and columns are spatial summaries
End of explanation
# extract metadata
dr = meta_file['sp_dailymet_livneh_1970_1970']
# compute sums and mean monthly an yearly sums
ltm = ogh.aggregate_space_time_sum(df_dict=ltm,
                                   suffix='Prec_sp_dailymet_livneh_1970_1970',
                                   start_date=dr['start_date'],
                                   end_date=dr['end_date'])
# print the name of the analytical dataframes and values within ltm
sorted(ltm.keys())
# initialize list of outputs
files=[]
# create the destination path for the dictionary of dataframes
ltm_sauk=os.path.join(homedir, 'ltm_1970_1970_sauk.json')
ogh.saveDictOfDf(dictionaryObject=ltm, outfilepath=ltm_sauk)
files.append(ltm_sauk)
# append the mapping file for Sauk-Suiattle gridded cell centroids
files.append(mappingfile1)
Explanation: Compute the total monthly and yearly precipitation, as well as the mean values across time and across stations
INPUT: daily precipitation for each station from the long-term mean dictionary (ltm) <br/>OUTPUT: Append the computed dataframes and values into the ltm dictionary
End of explanation
# # two lowest elevation locations
lowE_ref = ogh.findCentroidCode(mappingfile=mappingfile1, colvar='ELEV', colvalue=164)
# one highest elevation location
highE_ref = ogh.findCentroidCode(mappingfile=mappingfile1, colvar='ELEV', colvalue=2216)
# combine references together
reference_lines = highE_ref + lowE_ref
reference_lines
ogh.renderValueInBoxplot(ltm['meanbymonthsum_Prec_sp_dailymet_livneh_1970_1970'],
                         outfilepath='totalMonthlyRainfall.png', 
                         plottitle='Total monthly rainfall',
                         time_steps='month',
                         wateryear=True,
                         reference_lines=reference_lines,
                         ref_legend=True,
                         value_name='Total daily precipitation (mm)',
                         cmap='seismic_r',
                         figsize=(6,6))
ogh.renderValuesInPoints(ltm['meanbymonthsum_Prec_sp_dailymet_livneh_1970_1970'], 
                         vardf_dateindex=12, 
                         shapefile=sauk, 
                         cmap='seismic_r',
                         outfilepath='test.png', 
                         plottitle='December total rainfall',
                         colorbar_label='Total monthly rainfall (mm)', 
                         figsize=(1.5,1.5))
minx2, miny2, maxx2, maxy2 = oxl.calculateUTMbounds(mappingfile=mappingfile1,
                                                    mappingfile_crs={'init':'epsg:4326'},
                                                    spatial_resolution=0.06250)
print(minx2, miny2, maxx2, maxy2)
Explanation: Visualize the "average monthly total precipitations"
INPUT: dataframe with each month as a row and each station as a column. <br/>OUTPUT: A png file that represents the distribution across stations (in Wateryear order)
End of explanation
help(oxl.rasterDimensions)
# generate a raster
raster, row_list, col_list = oxl.rasterDimensions(minx=minx2, miny=miny2, maxx=maxx2, maxy=maxy2, dx=1000, dy=1000)
raster.shape
Explanation: generate a raster
End of explanation
help(oxl.mappingfileToRaster)
%%time
# landlab raster node crossmap to gridded cell id
nodeXmap, raster, m = oxl.mappingfileToRaster(mappingfile=mappingfile1, spatial_resolution=0.06250, 
                                           minx=minx2, miny=miny2, maxx=maxx2, maxy=maxy2, dx=1000, dy=1000)
# print the raster dimensions
raster.shape
%%time
nodeXmap.plot(column='ELEV', figsize=(10,10), legend=True)
# generate vector array of December monthly precipitation
prec_vector = ogh.rasterVector(vardf=ltm['meanbymonthsum_Prec_sp_dailymet_livneh_1970_1970'],
                          vardf_dateindex=12,
                          crossmap=nodeXmap,
                          nodata=-9999)
# close-off areas without data
raster.status_at_node[prec_vector==-9999] = CLOSED_BOUNDARY
fig =plt.figure(figsize=(10,10))
imshow_grid(raster, 
            prec_vector,
            var_name='Monthly precipitation',
            var_units=meta_file['sp_dailymet_livneh_1970_1970']['variable_info']['Prec'].attrs['units'],
            color_for_closed='black', 
            cmap='seismic_r')
tmax_vector = ogh.rasterVector(vardf=ltm['meanbymonth_Tmax_sp_dailymet_livneh_1970_1970'],
                          vardf_dateindex=12,
                          crossmap=nodeXmap,
                          nodata=-9999)
fig = plt.figure(figsize=(10,10))
imshow_grid(raster, 
            tmax_vector,
            var_name='Daily maximum temperature',
            var_units=meta_file['sp_dailymet_livneh_1970_1970']['variable_info']['Tmax'].attrs['units'],
            color_for_closed='black', symmetric_cbar=False, cmap='magma')
tmin_vector = ogh.rasterVector(vardf=ltm['meanbymonth_Tmin_sp_dailymet_livneh_1970_1970'],
                          vardf_dateindex=12,
                          crossmap=nodeXmap,
                          nodata=-9999)
fig = plt.figure(figsize=(10,10))
imshow_grid(raster, 
            tmin_vector,
            var_name='Daily minimum temperature',
            var_units=meta_file['sp_dailymet_livneh_1970_1970']['variable_info']['Tmin'].attrs['units'],
            color_for_closed='black', symmetric_cbar=False, cmap='magma')
# convert a raster vector back to geospatial presentation
t4, t5 = oxl.rasterVectorToWGS(prec_vector, nodeXmap=nodeXmap, UTM_transformer=m)
t4.plot(column='value', figsize=(10,10), legend=True)
# this is one decade
inputvectors = {'precip_met': np.tile(ltm['meandaily_Prec_sp_dailymet_livneh_1970_1970'], 15000),
                'Tmax_met': np.tile(ltm['meandaily_Tmax_sp_dailymet_livneh_1970_1970'], 15000),
                'Tmin_met': np.tile(ltm['meandaily_Tmin_sp_dailymet_livneh_1970_1970'], 15000)}
%%time
(VegType_low, yrs_low, debug_low) = run_ecohydrology_model(raster,
                                                           input_data=inputvectors,
                                                           input_file=InputFile,
                                                           synthetic_storms=False,
                                                           number_of_storms=100000,
                                                           pet_method='PriestleyTaylor')
plot_results(raster, VegType_low, yrs_low, yr_step=yrs_low-1)
plt.show()
plt.savefig('grid_low.png')
Explanation: Higher resolution children of gridded cells
get data from Lower resolution parent grid cells to the children
End of explanation
len(files)
# for each file downloaded onto the server folder, move to a new HydroShare Generic Resource
title = 'Computed spatial-temporal summaries of two gridded data product data sets for Sauk-Suiattle'
abstract = 'This resource contains the computed summaries for the Meteorology data from Livneh et al. 2013 and the WRF data from Salathe et al. 2014.'
keywords = ['Sauk-Suiattle', 'Livneh 2013', 'Salathe 2014','climate','hydromet','watershed', 'visualizations and summaries'] 
rtype = 'genericresource'
# create the new resource
resource_id = hs.createHydroShareResource(abstract, 
                                          title,
                                          keywords=keywords, 
                                          resource_type=rtype, 
                                          content_files=files, 
                                          public=False)
Explanation: Visualize the "average monthly total precipitation"
5. Save the results back into HydroShare
<a name="creation"></a>
Using the hs_utils library, the results of the Geoprocessing steps above can be saved back into HydroShare.  First, define all of the required metadata for resource creation, i.e. title, abstract, keywords, content files.  In addition, we must define the type of resource that will be created, in this case genericresource.  
Note: Make sure you save the notebook at this point, so that all notebook changes will be saved into the new HydroShare resource.
Total files and image to migrate
End of explanation
<END_TASK>