<jupyter_start><jupyter_text>3.a xarray![xarray logo](images/xarray_logo.png)https://xarray.pydata.org/en/stable/index.html**xarray** is a python package which allows us to handle multi-dimensional datasets in a simple way. It provides a huge set of functions for advanced analytics and visualization. It is part of the SciPy and Pangeo ecosystem.**xarray** data structure deals with scientific data by using labels, attributes, dimensions and coordinates, and extend the capabilities of **NumPy** and **pandas**. Data structures- DataArray- DataSet- Dimensions- CoordinatesDataArray:     N-dimensional array with dimensions. The objects add dimension names, coordinates, and attibutes to the underlying data structure (numpy and dask arrays).Dataset:     Dict-like collection of DataArray objects with aligned dimensions. Similar use of variables, dimensions, coordinates, and attributes like for DataArray. You can see an xarray Dataset as a netCDF file like object. Dimensions:     Named dimension axes, if missing the dimension names are dim_0, dim_1, ...Coordinates:     An array which labels a dimension. Two types are defined a) dimension coordinates - 1-dimensional coordinate array assigned to the DataArray with a name and dimension name. b) Non-dimensional coordinate - a coordinate array assigned to DataArray with the name assigned to the coordinates and not to the dimensions.<jupyter_code>%matplotlib inline

import numpy as np
import pandas as pd
import xarray as xr<jupyter_output><empty_output><jupyter_text>Working with DataArraysFirst, we create a random data array a with 20 values with numpy's ```random.rand()``` function.<jupyter_code>a = np.random.rand(20)

print(a)<jupyter_output><empty_output><jupyter_text>Make an xarray DataArray from the numpy array a with ```xarrays.DataArray()```.<jupyter_code>da_a = xr.DataArray(a)

print(da_a)<jupyter_output><empty_output><jupyter_text>As you can see a dimension ```dim_0``` is assed to the array.For n-dimensional arrays, a corresponding number of dimensions are used.E.g. 3D data array:<jupyter_code>data = np.random.rand(4,90,180)

print(data)
print(xr.DataArray(data))<jupyter_output><empty_output><jupyter_text>The dimensions have no names and we want to change it in the next step with the ```coords``` and ```dims``` parameters.<jupyter_code>time = pd.date_range("2020-01-01", periods=4)
lat = np.linspace( -90.0, 90.0,  90) 
lon = np.linspace(-180., 180.0, 180)

da = xr.DataArray(data, coords=[time,lat,lon], dims=['time','lat','lon'])

print(da)<jupyter_output><empty_output><jupyter_text>Available DataArray attributes are- values- data- coords- dims- sizes- name- attrs<jupyter_code>print(da.values)
print(da.data)
print(da.coords)
print(da.dims)
print(da.sizes)
print(da.name)

da.name = 'data'

print(da.name)
print(da.attrs)

da.attrs['units'] = 'data units'

print(da.attrs['units'])<jupyter_output><empty_output><jupyter_text>Like numpy xarray provides some array methods e.g. where, min, and max.<jupyter_code>A = xr.DataArray(np.arange(1, 26).reshape(5, 5), dims=('x', 'y'))
A      #-- if used without print notebooks will give additional informations when available<jupyter_output><empty_output><jupyter_text>```xr.DataArray.where()``` example<jupyter_code>print(A.where(A.x > 2))
print(A.where(A.x + A.y > 2))
print(A.where(A > 10))<jupyter_output><empty_output><jupyter_text>Use numpy's 'where' method for replacing values.<jupyter_code>print(np.where(A > 10, A, -9999.9))<jupyter_output><empty_output><jupyter_text>Working with DatasetsA Dataset can contain multiple variables with different dimensions and coordinates.Define two random data arrays, temp and prec of size (12,90,180).<jupyter_code>temp = np.random.uniform(low=265, high=310, size=(12,90,180)) 
prec = np.random.uniform(low=0.0001, high=0.001, size=(12,90,180))<jupyter_output><empty_output><jupyter_text>Now, we want to generate and add coordinate variables to the dataset.To create a time coordinate we use pandas ```date_range()```function. 12 time steps, 15th Jan to 15th Dec 2020.<jupyter_code>time = pd.date_range(start='2020-01-1', periods=12, freq='SM')

print(time)<jupyter_output><empty_output><jupyter_text>Create the coordinate variable arrays for longitude and latitude with numpy's ```linspace()``` function.<jupyter_code>lat = np.linspace(-90.0, 90.0, 90)
lon = np.linspace(-180.0, 180.0, 180)

print(lat)
print(lon)<jupyter_output><empty_output><jupyter_text>All we need is defined and we can create the dataset. The coordinate variables and the variable temp will be assigned to the dataset.<jupyter_code>ds = xr.Dataset(data_vars={'temperature':(['time','lat','lon'], temp),}, 
                coords={'time':('time', time), 
                        'lat':(['lat'], lat), 
                        'lon':(['lon'], lon)})

print(ds)<jupyter_output><empty_output><jupyter_text>Instead of using the print function, the info method of xarray Datasets can be used. The result looks very similar to the output of ncdump.<jupyter_code>ds.info()<jupyter_output><empty_output><jupyter_text>Read data from fileThe function ```open_dataset()``` of xarray is used to read the content of the file.<jupyter_code>import xarray as xr
import numpy as np

fname = './data/tsurf.nc'

ds = xr.open_dataset(fname)

ds.info()<jupyter_output><empty_output><jupyter_text>Printing the dataset content gives you an overview of the dimension and variable names, their sizes, and the global file attributes. Show variable names and coordinatesIt is always good to have a closer look at your data, and this can be done very easily.Ok, show me the variables stored in that file (ups - just one :D) and the coordinate variables, too.<jupyter_code>coords    = ds.coords
variables = ds.variables

print('--> coords:    \n\n', coords)
print('--> variables: \n\n', variables)<jupyter_output><empty_output><jupyter_text>Ah, that's better. Here we can see the time displayed in a readable way, because xarray use the datetime64 module under the hood. Also the variable and coordinate attributes are shown.  Select variable and coordinate variablesAt the moment, we only have created a dataset respectively a file object containing the coordinate variables and variable data. Now, we want to select the variable **tsurf** and the coordinate variables **lat** and **lon**.<jupyter_code>tsurf = ds.tsurf
lat   = tsurf.lat
lon   = tsurf.lon

print('Variable tsurf:            \n', tsurf.data)
print('\nCoordinate variable lat: \n', lat.data)
print('\nCoordinate variable lon: \n', lon.data)<jupyter_output><empty_output><jupyter_text>The variable types have the type ```xr.DataArray()```.<jupyter_code>print(type(tsurf))<jupyter_output><empty_output><jupyter_text>Dimensions, shape and sizeTo get more informations about the dimension, shape and size of a variable we can use the approbriate attributes.<jupyter_code>dimensions = ds.dims
shape = tsurf.shape
size  = tsurf.size
rank  = len(shape)

print('dimensions: ', dimensions)
print('shape:      ', shape)
print('size:       ', size)
print('rank:       ', rank)<jupyter_output><empty_output><jupyter_text>Variable attributesVariable attributes are very important to work in a correct manor with the data.<jupyter_code>attributes = list(tsurf.attrs)

print('attributes: ', attributes)<jupyter_output><empty_output><jupyter_text>Let's see how we can get the content of an attribute.<jupyter_code>long_name = tsurf.long_name
units = tsurf.units

print('long_name: ', long_name)
print('units:     ', units)<jupyter_output><empty_output><jupyter_text>TimeXarray is able to convert the time values to readable times using the internally datetime64 module.<jupyter_code>time = ds.time.data

print('timestep 0: ', time[0])<jupyter_output><empty_output><jupyter_text>Read a GRIB fileTo read a GRIB file xarray needs an additional module ```cfgrib```, which is used as an so called _engine_.<jupyter_code>import cfgrib

ds2 = xr.open_dataset('./data/MET9_IR108_cosmode_0909210000.grb2', engine='cfgrib')

variables2 = ds2.variables

print('--> variables2: \n\n', variables2)<jupyter_output><empty_output><jupyter_text>ReshapingThere are different ways to swap the dimensions of an array from (x,y) to (y,x).<jupyter_code>B = xr.DataArray(np.arange(1, 31).reshape(6, 5), dims=('x', 'y'))

print(B)
print(B.transpose())
print(B.T)<jupyter_output><empty_output><jupyter_text>Computations (xarray methods)Xarray includes the scientific libraries of Python stack, Numpy and pandas. This means we can use Numpy's arithmetic functions for computations. Arithmetic computations with arrays<jupyter_code>C = xr.DataArray(np.random.uniform(low=0, high=100, size=(10,20)) , dims=('x','y'))
D = xr.DataArray(np.random.uniform(low=0, high=100, size=(10,20)) , dims=('x','y'))

print(C)
print(D)<jupyter_output><empty_output><jupyter_text>Addition of constant value<jupyter_code>C_add_value = C + 5.

print('Original value: %f  new value: %f' % (C[3,0], C_add_value[3,0]))<jupyter_output><empty_output><jupyter_text>Addition of two arrays same size<jupyter_code>CD_add = C + D

print('C[3,0]: %f  D[3,0]: %f  Added: %f' % (C[3,0], D[3,0], CD_add[3,0]))<jupyter_output><empty_output><jupyter_text>Basic methods<jupyter_code>print('Minimum value = ', C.min().values)
print('Maximum value = ', C.max().values)
print('Sum           = ', C.sum().values)<jupyter_output><empty_output><jupyter_text>Advanced methods<jupyter_code>print(C.mean(dim='x'))<jupyter_output><empty_output><jupyter_text>Working with missing valuesSet the value -9999 in array to missing value. Numpy's np.nan method is used to define a missing value.<jupyter_code>tarray = xr.DataArray(data=[0, 1, -9999, 3, 4, 5, 6, 7, -9999, 9, 10], dims='x')

tarray = tarray.where(tarray != -9999, np.nan)

print(tarray)<jupyter_output><empty_output><jupyter_text>Check if missing values exist. It returns a mask array of True/False elements.<jupyter_code>print(tarray.isnull())<jupyter_output><empty_output><jupyter_text>Now, create a mask array where the values are not missing values.<jupyter_code>print(tarray.notnull())<jupyter_output><empty_output><jupyter_text>Count value that are not missing values.<jupyter_code>print(tarray.count())<jupyter_output><empty_output><jupyter_text>Return all array elements that are not missing values.<jupyter_code>print(tarray.dropna(dim='x'))<jupyter_output><empty_output><jupyter_text>Set missing value to a constant number.<jupyter_code>print(tarray.fillna(0))<jupyter_output><empty_output><jupyter_text>Interpolation methodsInterpolation on 1D array.Define 1D-array data1D:    y-axis values in range -1.0 to 1.0 with increment 0.1    x-axis values in range 0.0 to 1 with increment 0.05<jupyter_code>data1D = xr.DataArray(data=[0., 3., 2.1, 2., 1.7, 5., 5.2, 3.3, 2.5, 4.], 
                                  dims='x', 
                                  coords={"x": np.linspace(0, 10, 10)})
print(data1D)

data1D.plot()<jupyter_output><empty_output><jupyter_text>Interpolate 1D values from 21 elements to 101 elements:<jupyter_code>data1D_interp = data1D.interp(x=np.linspace(0, 10, 30))

print(data1D_interp)

data1D_interp.plot()<jupyter_output><empty_output><jupyter_text>Interpolation on nD array.Define nD-array data1D:<jupyter_code>data2D = xr.DataArray(np.random.uniform(low=0, high=2, size=(15,20)), 
                                        dims=('y','x'), 
                                        coords={'y':range(15), 'x':range(20)})
print(data2D)
data2D.plot()<jupyter_output><empty_output><jupyter_text>Increase the resolution of the grid.<jupyter_code>data2D_linear = data2D.interp(x=np.linspace(0, 20, 100), y=np.linspace(0, 15, 100))

data2D_linear.plot()<jupyter_output><empty_output><jupyter_text>Write a Dataset to netCDF file<jupyter_code>temp = np.random.uniform(low=265, high=310, size=(12,90,180)) 
prec = np.random.uniform(low=0.0001, high=0.001, size=(12,90,180))

time = pd.date_range(start='2020-01-1', periods=12, freq='SM')
lat = np.linspace(-90.0, 90.0, 90)
lon = np.linspace(-180.0, 180.0, 180)

ds = xr.Dataset(data_vars={'temperature':(['time','lat','lon'], temp),}, 
                coords={'time':('time', time), 
                        'lat':(['lat'], lat), 
                        'lon':(['lon'], lon)})

ds.to_netcdf("my_data.nc")<jupyter_output><empty_output><jupyter_text>$ ncdump -h my_data.nc```netcdf my_data {dimensions:	time = 12 ;	lat = 90 ;	lon = 180 ;variables:	double temperature(time, lat, lon) ;		temperature:_FillValue = NaN ;	int64 time(time) ;		time:units = "days since 2020-01-15 00:00:00" ;		time:calendar = "proleptic_gregorian" ;	double lat(lat) ;		lat:_FillValue = NaN ;	double lon(lon) ;		lon:_FillValue = NaN ;}```  PlottingSome additional examples how to use the plot method.<jupyter_code>import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(6,4))

da = xr.DataArray(np.sin(np.linspace(0, 2 * np.pi, 10)), dims="x", coords={"x": np.linspace(0, 1, 10)})

da.plot.line('o', label='original')
da.interp(x=np.linspace(0, 1, 100)).plot.line(label='linear (default)')
da.interp(x=np.linspace(0, 1, 100), method='cubic').plot.line(label='cubic')
plt.legend()
data2D.plot()
xr.open_dataset('./data/tsurf.nc').tsurf[0,:,:].plot.contourf(levels=20)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1343679031037093, 'fraction_numerical': 0.03088569265707797, 'mean_word_length': 5.546580773042616, 'pattern_counts': {'":': 2, '<': 165, '<?xml version=': 0, '>': 172, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '41', 'dir': 'jupyter-structured-clean-dedup', 'id': '515267', 'license': 'MIT', 'n_tokens_mistral': 4809, 'n_tokens_neox': 4511, 'n_words': 1167, 'path': 'Tutorial/03a_xarray_introduction.ipynb', 'repo_name': '1271756664/-xESMF', 'size': '13210'}