Gradient and Integration

xr-scipy wraps scipy.gradient and some of scipy.integrate functions. Let’s create a simple example DataArray:

In [1]: arr = xr.DataArray(np.sin(np.linspace(0, 6.28, 30)) ** 2,
   ...:                    dims=('x'), coords={'x': np.linspace(0, 5, 30)})
   ...: 

In [2]: arr
Out[2]: 
<xarray.DataArray (x: 30)>
array([0.00000000e+00, 4.61661813e-02, 1.76139460e-01, 3.65918356e-01,
       5.80457402e-01, 7.80138804e-01, 9.28088451e-01, 9.96985262e-01,
       9.74106426e-01, 8.63676857e-01, 6.86089001e-01, 4.74137072e-01,
       2.66961112e-01, 1.02819215e-01, 1.20225984e-02, 1.13381944e-02,
       1.00892388e-01, 2.64147680e-01, 4.70956575e-01, 6.83128766e-01,
       8.61483534e-01, 9.73085045e-01, 9.97324436e-01, 9.29725546e-01,
       7.82771507e-01, 5.83599545e-01, 3.68989697e-01, 1.78572829e-01,
       4.75122222e-02, 1.01461475e-05])
Coordinates:
  * x        (x) float64 0.0 0.1724 0.3448 0.5172 ... 4.483 4.655 4.828 5.0

Our gradient() takes an xarray object (possibly high dimensional) and a coordinate name which direction we compute the gradient of the array,

In [3]: grad = xrscipy.gradient(arr, 'x')

In [4]: grad
Out[4]: 
<xarray.DataArray (x: 30)>
array([ 0.26776385,  0.51080443,  0.92728131,  1.17252203,  1.2012393 ,
        1.00813004,  0.62885473,  0.13345213, -0.38659437, -0.83525053,
       -1.12966538, -1.21547088, -1.07682178, -0.73932169, -0.26529496,
        0.25772239,  0.73314751,  1.07318614,  1.21504515,  1.13252818,
        0.84087321,  0.39393861, -0.12574255, -0.62220349, -1.0037654 ,
       -1.19996725, -1.17457748, -0.93228468, -0.51783178, -0.27551204])
Coordinates:
  * x        (x) float64 0.0 0.1724 0.3448 0.5172 ... 4.483 4.655 4.828 5.0

The return type is also a DataArray with coordinates.

In [5]: arr.plot(label='arr')
Out[5]: [<matplotlib.lines.Line2D at 0x7fdc21f31450>]

In [6]: grad.plot(label='gradient')
Out[6]: [<matplotlib.lines.Line2D at 0x7fdc21dc1e50>]

In [7]: plt.legend()
Out[7]: <matplotlib.legend.Legend at 0x7fdc21efbd90>

In [8]: plt.show()
_images/grad.png

The other options (edge_order) are the same to numpy.gradient. See gradient().

Similar to gradient(), xr-scipy wraps some functions in scipy.integrate module. Our integration function also takes an xarray object and coordinate name along which the array to be integrated. The return type is also a DataArray,

# trapz computes definite integration
In [9]: xrscipy.integrate.trapz(arr, coord='x')
Out[9]: 
<xarray.DataArray ()>
array(2.50124814)

# cumurative integration returns a same shaped array
In [10]: integ = xrscipy.integrate.cumtrapz(arr, 'x')

In [11]: integ
Out[11]: 
<xarray.DataArray (x: 30)>
array([0.        , 0.00397984, 0.02314412, 0.06987324, 0.15145736,
       0.26875014, 0.41601111, 0.58196574, 0.75188744, 0.91031703,
       1.04391753, 1.14393702, 1.2078248 , 1.23970241, 1.24960257,
       1.25161643, 1.26129148, 1.29276045, 1.35613151, 1.45562162,
       1.58877786, 1.74693032, 1.91679321, 2.08291821, 2.23054726,
       2.34833787, 2.43045763, 2.4776613 , 2.49715139, 2.50124814])
Coordinates:
  * x        (x) float64 0.0 0.1724 0.3448 0.5172 ... 4.483 4.655 4.828 5.0

In [12]: arr.plot(label='arr')
Out[12]: [<matplotlib.lines.Line2D at 0x7fdc21d2d8d0>]

In [13]: integ.plot(label='integration')
Out[13]: [<matplotlib.lines.Line2D at 0x7fdc22088810>]

In [14]: plt.legend()
Out[14]: <matplotlib.legend.Legend at 0x7fdc21d5de50>

In [15]: plt.show()
_images/cumtrapz.png

See trapz() for other options.

Note

There are slight difference from the original implementations. Our gradient() does not accept multiple coordinates. Our cumtrapz() always assume initial=0.