from IPython.display import Image
Image(url='http://xray.readthedocs.org/en/latest/_images/dataset-diagram.png', embed=True, width=950, height=300)

"""
Explanation: Examinging NetCDF files using xray
The simplest way that I have found for opening and exploring NetCDF files in python, depends on the python package called xray. Here is a little graphical representation of the way to think about this data. For clarification on how multidimensional data are represented in xray, or to figure out how to download the package, visit: http://xray.readthedocs.org/en/latest/
End of explanation
"""


import os
import posixpath      # similar to os, but less dependant on operating system
import numpy as np
import pandas as pd
import xray

NETCDF_DIR = os.getcwd().replace('\\','/')+'/raw_netCDF_output/'
datafile = 'soil'

nc_file = os.listdir(NETCDF_DIR+datafile)[-1]
nc_path = posixpath.join(NETCDF_DIR, datafile, nc_file)
ds = xray.open_dataset(nc_path)
ds

"""
Explanation: Loading a NetCDF file into a dataset
The first step when you receive a NetCDF file is to open it up and see what it contains.
End of explanation
"""


ds.sel(site='Open_W').coords

"""
Explanation: Inspecting and selecting from dataset
To inspect the coordinates at a specific site (for example, 'Open_W') we just write:
End of explanation
"""


print ds.VW_05cm_Avg.sel(site='Open_W', time='2015-06-02T06:10:00')

"""
Explanation: Now if we are only interested in soil moisture at the upper depth at a specific time (don't forget that the time is in UTC unless the timezone is explicit), we can pull out just that one data point:
End of explanation
"""


print ds.VW_05cm_Avg.sel(site='Open_W', time='2015-06-02T06:10:00').values

"""
Explanation: And if we are only interested in the actual value rather than all the attributes:
End of explanation
"""


ds_dict = {}
nc_files = os.listdir(NETCDF_DIR+datafile)
for nc_file in nc_files:
    nc_path = posixpath.join(NETCDF_DIR, datafile, nc_file)
    ds = xray.open_dataset(nc_path)
    date = nc_file.split('Tower_')[1].split('.')[0]
    ds_dict.update({date: ds})

ds_dict.keys()

xray.concat(ds_dict.values(), dim='time')

"""
Explanation: Concatenate files
End of explanation
"""


df = ds.to_dataframe()

for i in range(len(df.index.levels)):
    print 'df.index.levels[{i}]\n{index}\n'.format(i=i, index=df.index.levels[i])

"""
Explanation: Convert to pandas dataframe
Some analyses are no doubt easier to carry out in pandas, but luckily xray makes it very easy to move back and forth between the two packages. To converat an xray Dataset object to a pandas MultiIndex DataFrame object, just run:
End of explanation
"""


