Distribution Statement
# # # This source code is subject to the license referenced at
# # # NRLMMD-GEOIPS.

Build a Reader: a hands-on walkthrough#

This tutorial builds two complete, class-based readers from scratch and runs them through Order-Based Processing. It complements the conceptual reader reference with a step-by-step example.

We build:

  1. A single-file reader — for data delivered in one file (as geostationary sensors often are).

  2. A multi-file reader — for data spread across many files (as polar-orbiter passes often are).

To keep the focus on the reader itself, both examples use a small fictional sensor whose data are stored in HDF5. Apply the same patterns to your own data.

What a reader does#

A reader does two jobs:

  • It reads a data file (or files) into the GeoIPS internal representation — a dictionary of xarray.Dataset objects.

  • It reports the metadata GeoIPS needs to decide whether and how to process the data (data source, start/end times, geolocation).

At minimum, a reader must return datasets containing:

  • latitude and longitude (2-D, matching the data),

  • your data variable(s), named however you like,

  • and the required attributes source_name, start_datetime, and end_datetime.

Important

Variable names and capitalization matter. The variable names your reader produces must match the names your product requests.

Where readers fit in#

In Order-Based Processing, a reader is a reader step in a workflow. GeoIPS may call your reader twice: first with metadata_only=True to read just the metadata (so it can decide whether the run is worth doing and set up the area definition), then again to read the full data. Always honor metadata_only — see challenge #1.

Note

A reader’s call method keeps the argument names fnames and chans. When a reader runs as an OBP workflow step, GeoIPS supplies these as filenames and variables and renames them for you — so author your reader with fnames / chans (as all GeoIPS readers do) and reference filenames / variables in the workflow YAML.

Know your data first#

Before writing a reader, understand your file’s contents — its variables, dimensions, and metadata attributes. Documentation, ncdump/HDFView (NetCDF4/HDF5), or wgrib (GRIB) all help. A quick throwaway script also works for NetCDF4/HDF5:

import xarray

xobj = xarray.open_dataset("myfile.h5")

print("FILE METADATA ATTRIBUTES")
for attr in xobj.attrs:
    print(attr, ":", xobj.attrs[attr])

print("\nVARIABLES + METADATA ATTRIBUTES")
for var in xobj:
    print(var, ":", xobj[var].attrs or "No attributes")

Reader #1 — all data in one file#

Create the reader in your package’s class-based readers directory:

cd $MY_PKG_DIR/$MY_PKG_NAME/plugins/classes/readers

Create a file named fl89_value_hdf5.py. The naming convention is <source>_<datatype>_<format> (recommended, not required).

Plugin header and class#

Start with a module docstring, imports, the logger, and the class attributes that tell GeoIPS what this plugin is:

"""Reader for the single-file FL-89 sensor (HDF5)."""

# Python standard libraries
import datetime
import logging

# Third-party libraries
import numpy
import xarray

from geoips.interfaces.class_based.readers import BaseReaderPlugin

LOG = logging.getLogger(__name__)


class GreatValueHDF5ReaderPlugin(BaseReaderPlugin):
    """Read FL-89 data from a single HDF5 file."""

    interface = "readers"
    family = "standard"
    name = "fl89_value_hdf5"
    source_names = ["fl89"]

interface/family/name identify the plugin; source_names declares which data sources it can read (products reference the same source name).

The call method#

call is what GeoIPS invokes to use your reader. Readers receive a standard set of arguments:

def call(self, fnames, metadata_only=False, chans=None, area_def=None,
         self_register=False):
    """Read FL-89 data from a list of filenames.

    Parameters
    ----------
    fnames : list of str
        Full paths to the files to read.
    metadata_only : bool, default=False
        Return after reading metadata only, if True.
    chans : list of str, default=None
        Channels/variables requested by the product; read all if None.
    area_def : pyresample.AreaDefinition, default=None
        Region to constrain the read (if the reader implements it).
    self_register : str or bool, default=False
        Return unprojected/self-registered data (if implemented).

    Returns
    -------
    dict of xarray.Dataset
        Keyed by descriptive dataset ids, plus a ``"METADATA"`` entry.

    See Also
    --------
    :ref:`xarray_standards`
    """

Note

area_def, self_register, chans, and any other argument only do something if your reader acts on them. They are part of the standard signature, but honoring them is up to you.

Read the metadata#

This sensor delivers a single file, so read it and set the required attributes:

fname = fnames[0]  # this sensor delivers a single file
LOG.interactive("Reading %s", fname)

# xarray reads HDF5/NetCDF4 without importing the backend directly
xobj = xarray.open_dataset(fname)

sdt = datetime.datetime.strptime(
    xobj.attrs["time_data_start"], "%Y-%m-%dT%H:%M:%SZ"
)
edt = datetime.datetime.strptime(
    xobj.attrs["time_data_end"], "%Y-%m-%dT%H:%M:%SZ"
)

xobj.attrs["start_datetime"] = sdt
xobj.attrs["end_datetime"] = edt
xobj.attrs["source_name"] = "fl89"
xobj.attrs["platform_name"] = "fl89"
xobj.attrs["data_provider"] = "Frito-Lay"

Return early for the metadata-only pass#

Challenge #1 — readers are called twice

Because GeoIPS calls your reader once with metadata_only=True and again with metadata_only=False, a large dataset (ABI, VIIRS, …) will be read twice if you don’t return early — doubling the memory and time cost. Always return after the metadata pass:

if metadata_only is True:
    LOG.info("metadata_only requested. Returning")
    return {"METADATA": xobj}

Build 2-D geolocation#

Challenge #2 — 1-D geolocation

Some datasets store 1-D latitude/longitude to save space, but GeoIPS expects latitude/longitude to be 2-D, matching the data. Build a 2-D grid with numpy.meshgrid:

LOG.info("Obtaining lat/lon from xarray")
lat = xobj.variables["Latitude"][...]
lon = xobj.variables["Longitude"][...]

LOG.info("Creating 2D lat/lon grid")
lon_final, lat_final = numpy.meshgrid(lon, lat)

xobj["latitude"] = xarray.DataArray(
    numpy.ma.array(lat_final), dims=("lat", "lon")
)
xobj["longitude"] = xarray.DataArray(
    numpy.ma.array(lon_final), dims=("lat", "lon")
)

# Drop the original 1-D coordinates now that we have 2-D versions
xobj = xobj.drop_vars(["Latitude", "Longitude"])

Return the data#

Return a dictionary of xarray datasets — one (or more) named datasets plus a METADATA entry — and declare PLUGIN_CLASS so the registry can find your reader:

        return {"FL89": xobj, "METADATA": xobj[[]]}


PLUGIN_CLASS = GreatValueHDF5ReaderPlugin

Your first reader is complete. Check your indentation: everything from class … is indented three spaces, everything from def call… an additional three, and PLUGIN_CLASS has no indentation.

Reader #2 — data spread across many files#

Polar-orbiter passes are often delivered as many files. This reader reads one or many files and stitches them together. Create fl98_value_hdf5.py in the same directory, with the same header pattern:

"""Reader for the multi-file FL-98 sensor (HDF5)."""

import datetime
import logging

import numpy
import xarray

from geoips.interfaces.class_based.readers import BaseReaderPlugin

LOG = logging.getLogger(__name__)


class Fl98ReaderPlugin(BaseReaderPlugin):
    """Read FL-98 data from one or more HDF5 files."""

    interface = "readers"
    family = "standard"
    name = "fl98_value_hdf5"
    source_names = ["fl98"]

Metadata from the first and last files#

With multiple files, derive the start/end times from the first and last files. Read only what you need for the metadata pass:

def call(self, fnames, metadata_only=False, chans=None, area_def=None,
         self_register=False):
    """Read FL-98 data from a list of filenames."""
    xobj_metadata = xarray.Dataset()

    xobj_temp = xarray.open_dataset(fnames[0])
    sdt = datetime.datetime.strptime(
        xobj_temp.attrs["time_data_start"], "%Y-%m-%dT%H:%M:%SZ"
    )
    xobj_temp.close()

    xobj_temp = xarray.open_dataset(fnames[-1])
    edt = datetime.datetime.strptime(
        xobj_temp.attrs["time_data_end"], "%Y-%m-%dT%H:%M:%SZ"
    )
    xobj_temp.close()

    xobj_metadata.attrs["start_datetime"] = sdt
    xobj_metadata.attrs["end_datetime"] = edt
    xobj_metadata.attrs["source_name"] = "fl98"
    xobj_metadata.attrs["platform_name"] = "fl98"
    xobj_metadata.attrs["data_provider"] = "Frito-Lay"

    if metadata_only is True:
        return {"METADATA": xobj_metadata}

Note

Using fnames[0] and fnames[-1] assumes the file list is sorted. Verify the ordering for your own data rather than relying on it.

A dedicated read function#

There are many valid ways to structure a reader. Here we factor the per-file read into a helper method that reads geolocation and only the requested channels. Add it to the class, above call:

def read_fl98_data(self, fname, chans=None):
    """Read one FL-98 file; return a dict of xarray datasets."""
    import h5py

    datain = h5py.File(fname)
    return_dict = {}
    temp_xarray = xarray.Dataset()

    # Geolocation → 2-D (see Challenge #2)
    lat = datain["Latitude"][:]
    lon = datain["Longitude"][:]
    coordict = dict(
        num_lats=range(lat.shape[0]), num_lons=range(lon.shape[0])
    )
    dims = ("num_lats", "num_lons")
    lon_final, lat_final = numpy.meshgrid(lon, lat)
    temp_xarray["latitude"] = xarray.DataArray(
        lat_final, coords=coordict, dims=dims
    )
    temp_xarray["longitude"] = xarray.DataArray(
        lon_final, coords=coordict, dims=dims
    )

    # Read ONLY the requested channels
    for chan in chans:
        temp_xarray[chan] = xarray.DataArray(
            datain[chan][:, :], coords=coordict, dims=dims
        )

    datain.close()
    return_dict["DATA"] = temp_xarray
    return return_dict

Challenge #3 — awkward variable names

Some sensors name variables in ways that are unwieldy as GeoIPS variables (e.g. AMSR2’s Brightness_Temperature_89_GHz_AH). Map source names to GeoIPS names as you read:

varnames = {"Brightness_Temperature_89_GHz_AH": "tb89hA", ...}
for src_name, geoips_name in varnames.items():
    sub_xarray[geoips_name] = full_xarray[src_name]

Tip

Reading only the channels the product requested (chans) keeps memory and I/O down — don’t read variables nobody asked for.

Loop, stack, and concatenate#

Back in call (after the metadata_only return), read every file, stack each variable along a new dimension, then concatenate:

        xarray_fl98 = {}  # dict to be returned

        for fname in fnames:
            temp_dict = self.read_fl98_data(fname, chans)
            for key in temp_dict:
                stacked = temp_dict[key].stack(zed=temp_dict[key].dims)
                if key in xarray_fl98:
                    xarray_fl98[key].append(stacked)
                else:
                    xarray_fl98[key] = [stacked]

        for key in xarray_fl98:
            xarray_fl98[key] = xarray.concat(xarray_fl98[key], dim="zed")

        for dsname, curr_xarray in xarray_fl98.items():
            curr_xarray.attrs["start_datetime"] = sdt
            curr_xarray.attrs["end_datetime"] = edt
            curr_xarray.attrs["source_name"] = "fl98"
            curr_xarray.attrs["platform_name"] = "fl98"
            curr_xarray.attrs["data_provider"] = "Frito-Lay"

        xarray_fl98["METADATA"] = list(xarray_fl98.values())[0][[]]
        return xarray_fl98


PLUGIN_CLASS = Fl98ReaderPlugin

Challenge #4 — overlapping / masked data

stack uses the most recent file’s mask, which can silently drop data where files overlap with masked space. Be careful with geographically overlapping files.

Testing your reader#

To produce imagery, GeoIPS needs a product (what to make from the data) and a workflow (the ordered steps that make it).

1. Create a product#

cd $MY_PKG_DIR/$MY_PKG_NAME/plugins/yaml/products

Create fl89.yaml. The variable names must match what your reader produced:

interface: products
family: list
name: fl89
docstring: The fl89 products configuration.
spec:
  products:
    - name: coupon
      source_names: ["fl89"]
      docstring: FL-89 coupon RGB.
      product_defaults: RGB-Default
      spec:
        variables: ["ChanRed", "ChanGrn", "ChanBlu"]
        scale_factor: 0.003

See the products tutorial for more on product plugins.

2. Create a workflow#

cd $MY_PKG_DIR/$MY_PKG_NAME/plugins/yaml/workflows

Create fl89_coupon.yaml — a workflow that reads the data with your new reader and produces the coupon product:

apiVersion: geoips/v1
interface: workflows
family: order_based
name: fl89_coupon
docstring: FL-89 coupon RGB imagery.
package: <your_package>
spec:
  globals:
    product_name: coupon
  steps:
    sector:
      kind: sector
      name: global_cylindrical
    reader:
      kind: reader
      name: fl89_value_hdf5
      depends_on: [sector]
    fl89_coupon:
      kind: product
      name: [fl89, coupon]
      depends_on: [reader, sector]
    output_formatter:
      kind: output_formatter
      name: imagery_annotated
      depends_on:
        - fl89_coupon.algorithm
        - fl89_coupon.colormapper
        - sector

3. Register and run#

Install the package (only needed when you add new Python classes/modules) and rebuild the registries (needed whenever you add a class/module or edit a YAML plugin):

cd $GEOIPS_PACKAGES_DIR
pip install -e $MY_PKG_NAME
geoips config create-registries

Then run the workflow on your data file:

geoips run order_based fl89_coupon /path/to/your/fl89_data.h5

On success the log ends with Return Value 0 and prints the path to the output image. The FL-98 reader is tested the same way — create an fl98 product and an fl98_coupon workflow, then run it against your set of FL-98 files.

Tip

When do you need to reinstall vs. rebuild registries?

  • pip install -e — after adding a new Python class or module.

  • geoips config create-registries — after adding a class/module or editing a YAML plugin.

  • Neither — when editing the body of an already-installed, already-registered plugin in place.

More challenges to keep in mind#

  • Multiple file formats for one sensor. Some sensors ship the same data in different formats (e.g. VIIRS from NASA as HDF5 vs. NOAA as NetCDF), often structured completely differently — you may need more than one reader.

  • L1–L3 products. You can write one “super reader” or split across smaller readers per processing level.

  • Close your files. Prefer context management (with h5py.File(...) as f:); if you open a file another way, make sure you close() it.

In conclusion#

Writing a reader is as much art as science — there are many valid approaches, and the methods here won’t fit every dataset. Choose what works for your data, and be mindful of your users’ computational resources as you go.

See also#