Source code for pycif.plugins.datastreams.fluxes.GridFED.get_domain

import numpy as np
from logging import debug
import os
import xarray as xr
import itertools

from .....utils.classes.setup import Setup
from .....utils.classes.domains import Domain
from .....utils.check.errclass import CifError
from .....utils.hdf5 import _hdf5_lock




def get_domain(ref_dir, ref_file, input_interval, target_dir, tracer=None):
    """Read information from the reference file
    to define the data horizontal and, if relevant, vertical domain.

    Args:
        ref_dir (str): the path to the input files
        ref_file (str): format of the input files
        input_interval (list): simulation interval (start and end dates)
        target_dir (str): where to copy
        tracer: the tracer Plugin, corresponding to the paragraph
            :bash:`datavect/components/fluxes/parameters/my_species` in the
            configuration yaml; can be needed to fetch extra information
            given by the user

    Return:
        domain (Domain): a domain class object, with the definition of the center grid
            cells coordinates, as well as corners

    """

    ref_file = list(itertools.chain.from_iterable(tracer.input_files.values()))[0]

    if not os.path.isfile(ref_file):
        raise CifError(f"Could not initialize the domain as no reference file is available. Expecting the following file: {ref_file}")

    # Grid cell centers
    with _hdf5_lock:
        coords = xr.open_dataset(ref_file)[["lat", "lon"]]
        lat = coords["lat"].values
        lon = coords["lon"].values

    nlat = len(lat)
    nlon = len(lon)

    zlon, zlat = np.meshgrid(lon, lat)

    # Grid cell corners
    dlat = np.mean(np.unique(np.diff(lat)))
    dlon = np.mean(np.unique(np.diff(lon)))

    latc = np.append(lat - dlat / 2, lat[-1] + dlat / 2)
    lonc = np.append(lon - dlon / 2, lon[-1] + dlon / 2)
    zlonc, zlatc = np.meshgrid(lonc, latc)

    # Vertical definition (surface level)
    pressure_unit = "Pa"
    nlev = 1
    sigma_a_mid = np.array([0])
    sigma_b_mid = np.array([1])

    # Put it to a domain Plugin
    domain = Domain(nlon=nlon, nlat=nlat,
                    zlon=zlon, zlat=zlat,
                    zlonc=zlonc, zlatc=zlatc,
                    nlev=nlev, pressure_unit=pressure_unit,
                    sigma_b_mid=sigma_b_mid, sigma_a_mid=sigma_a_mid,
                    lon_cyclic=True)

    return domain