Raster

Background

This module extends xarray and rioxarray to provide additional functionality for raster processing and post-processing. These are used throughout the datasources to process various raster data, but are also useful directly for the end user, and is thus available for direct use. The extension is based on the guidance for how to extend xarray.

Usage

Accessor

To make the extension accessible, you only need to import ochanticipy directly. Then, extension methods and properties are simply accessed using da.oap.method() or da.oap.property on a data array or dataset. As a simple example, we can create a data array that has incorrect coordinates (ascending latitude and descending longitude). Then we can use the extension to invert them to the correct ordering.

import xarray
import numpy

import ochanticipy

da = xarray.DataArray(
  numpy.arange(16).reshape(4,4),
  coords={"lat":numpy.array([87, 88, 89, 90]),
          "lon":numpy.array([70, 69, 68, 67])}
)
da_inv = da.oap.invert_coordinates()

# check they have inverted
da_inv.get_index("lon")
#> Int64Index([67, 68, 69, 70], dtype='int64', name='lon')
da_inv.get_index("lat")
#> Int64Index([90, 89, 88, 87], dtype='int64', name='lat')

Raster statistics

One of the more useful functionalities of the module is providing an easy ability to calcuate raster statistics across a specified set of geometries defined in a geodataframe. Paired with the administrative boundaries and raster data sources available through this library, we can easily calculate raster statistics.

A full snippet of example code is available below.

import ochanticipy
import datetime

# load the administrative boundaries
country_config = ochanticipy.create_country_config(iso3="eth")
codab = ochanticipy.CodAB(country_config)
codab.download()
codab_eth = codab.load(admin_level=2)

# get geobounding box for CHIRPS downloading
geo_bounding_box = ochanticipy.GeoBoundingBox.from_shape(codab_eth)

# load CHIRPs data for processing
start_date = datetime.date(year=2001, month=2, day=1)
end_date = datetime.date(year=2006, month=3, day=31)

chirps_monthly = ChirpsMonthly(
    country_config=country_config,
    geo_bounding_box=geo_bounding_box,
    start_date=start_date,
    end_date=end_date
    )

chirps_monthly.download()
chirps_monthly.process()
chirps_monthly_data = chirps_monthly.load()

# compute raster statistics

chirps_monthly_data.oap.compute_raster_stats(
  gdf=codab_eth,
  feature_col="ADM2_PCODE"
)

Properties

The user should be careful when accessing attributes of a data array when using the raster module. This module builds on rioxarray <https://corteva.github.io/rioxarray>_ extensions, and thus methods and attributes accessible via da.rio.method() or da.rio.property are also accessible using da.oap.method() or da.oap.property. However, original rioxarray properties should be accessed using da.rio.property.

Let’s create a simple data array where we want to specify the spatial dimensions explicitly because the coordinate names are not automatically detected.

import xarray
import numpy
import ochanticipy

da = xarray.DataArray(
    numpy.arange(16).reshape(4,4),
    coords={"a":numpy.array([90, 89, 88, 87]),
            "b":numpy.array([70, 69, 68, 67])}
)

We can set the spatial dimensions using da.rio.set_spatial_dims() or call it directly from da.oap.

da_new = da.oap.set_spatial_dims(
  x_dim="a",
  y_dim="b"
)

However, even though we can set the dimensions using either accessor, we have to be careful accessing the properties.

da_new.rio.x_dim
#> 'a'

da_new.oap.x_dim
#> MissingSpatialDimensionError: x dimension not found.
#> 'rio.set_spatial_dims()' or using 'rename()' to change
#> the dimension name to 'x' can address this.

Even though the method was called using oap, the property is not accessible through it. Users need to be careful about accessing rioxarray properties using the oap accessor.

For best practice, rioxarray methods and properties should all be accessed using rio. These properties are rio.x_dim, rio.y_dim, rio.shape, rio.width, rio.height, and rio.crs. This module’s methods and properties should be accessed using the oap accessor. These properties are oap.t_dim and oap.longitude_range.