API Reference

grid_l2(config, grid_spec, time_spec; kwargs...)

Grid Level-2 satellite data with footprint-aware oversampling.

Reads input NetCDF files as specified by config, applies quality filters, subdivides each satellite footprint into sub-pixels for proper overlap handling, and accumulates gridded averages. Output is written to a NetCDF4 file.

Arguments

• config::DataSourceConfig: Data source configuration (from TOML/JSON)
• grid_spec::GridSpec: Output grid geometry
• time_spec::TimeSpec: Temporal binning parameters

Keyword Arguments

• n_oversample::Union{Nothing,Int}=nothing: Sub-pixel factor. nothing = auto-compute from footprint/grid ratio.
• compute_std::Bool=false: Also compute per-cell standard deviation
• outfile::String="gridded_output.nc": Output file path
• backend: KernelAbstractions backend (default: nothing for sequential Welford path). Use CPU() for threaded KA execution or CUDABackend() for GPU.

grid_center(config, grid_spec, time_spec; kwargs...)

Grid data using center coordinates only (no footprint bounds), suitable for MODIS-style data where each pixel maps to exactly one grid cell.

Keyword Arguments

• geo_table::Union{Nothing,String}=nothing: Path to a geolocation lookup table (NetCDF)
• veg_indices::Bool=false: Compute vegetation indices (EVI, NDVI, NIRv, NDWI)
• compute_std::Bool=false: Compute per-cell standard deviation
• outfile::String="gridded_output.nc": Output file path

load_config(config_path::String) -> DataSourceConfig

Load a configuration file (TOML or JSON) and return a DataSourceConfig.

The file must contain [basic] and [grid] sections. Filters are defined in an optional [filter] section using intuitive string expressions:

Example TOML

filePattern = "S5P_PAL__L2B_SIF____YYYYMMDD*.nc"
folder = "/path/to/data/"

[basic]
lat = "PRODUCT/latitude"
lon = "PRODUCT/longitude"
lat_bnd = "PRODUCT/SUPPORT_DATA/GEOLOCATIONS/latitude_bounds"
lon_bnd = "PRODUCT/SUPPORT_DATA/GEOLOCATIONS/longitude_bounds"

[grid]
sif_743 = "PRODUCT/SIF_743"
sif_735 = "PRODUCT/SIF_735"

[filter]
"PRODUCT/SUPPORT_DATA/GEOLOCATIONS/solar_zenith_angle" = "< 80"
"PRODUCT/qa_value" = "> 0.5"
"PRODUCT/methane" = "1600 < x < 2200"
"PRODUCT/mode" = "== 1"

The legacy JSON format with separate filter_gt/filter_lt/filter_eq sections is still supported for backward compatibility.

GridSpec{T<:AbstractFloat}

Specification of the output grid geometry.

Fields

• lat_min, lat_max: Latitude bounds (degrees)
• lon_min, lon_max: Longitude bounds (degrees)
• dlat, dlon: Grid cell size (degrees)
• lat, lon: Vectors of cell center coordinates

TimeSpec

Specification of the temporal gridding parameters.

Fields

• start_date, stop_date: Date range for processing
• time_step: Temporal bin size (Dates.Day or Dates.Month)
• oversample_temporal: Multiplier for the actual averaging window (>1 gives moving-average-like behavior)

DataSourceConfig

Configuration loaded from a TOML/JSON file that defines a satellite data source.

Fields

• basic: Maps internal keys ("lat", "lon", "lat_bnd", "lon_bnd") to variable paths in the NetCDF files
• grid_vars: Ordered mapping of output variable names to input variable paths (all will be gridded)
• filters: Quality filter rules parsed from [filter] section
• file_pattern: Glob pattern with YYYY/MM/DD/DOY placeholders for finding input files
• folder: Root folder for input data (may also contain YYYY/MM/DD placeholders)

FilterRule

A single filter condition for a NetCDF variable.

Fields

• variable: Path to the NetCDF variable (may include group separators /)
• op: Operation — :lt, :gt, :eq, :between
• lo: Lower bound (threshold for :gt/:lt/:eq, low end for :between)
• hi: Upper bound (only used for :between)

Config syntax (TOML [filter] section)

[filter]
"solar_zenith_angle" = "< 80"         # less than
"qa_value" = "> 0.5"                  # greater than
"mode" = "== 1"                       # equality
"methane" = "1600 < x < 2200"         # range (between)

sort_corners_ccw(lat_corners, lon_corners) -> (sorted_lat, sorted_lon)

Sort 4 footprint corners into counter-clockwise order by computing the angle from the centroid to each corner. This ensures edges 1→2 and 4→3 are always opposite edges, which is required by compute_subpixels!.

Different satellite products store corners in different orders (e.g., TROPOMI uses CCW order while OCO-2 uses a different convention). This function normalizes them.

sort_corners_ccw!(lat_mat, lon_mat)

Sort corners in-place for all soundings in N×4 matrices.

sort_corners_ccw_ka!(backend, lat_mat, lon_mat)

KA kernel that sorts footprint corners into CCW order for all soundings in parallel. Uses a 5-comparator sorting network (optimal for 4 elements, no allocations).

• lat_mat, lon_mat: N×4 matrices of corner coordinates

compute_subpixels!(points, vert_lat, vert_lon, n, lats_0, lons_0, lats_1, lons_1)

Fill the n×n×2 array points with sub-pixel positions inside the quadrilateral defined by 4 corner coordinates (vert_lat[1:4], vert_lon[1:4]).

The algorithm:

1. Subdivide edge 1→2 into n baseline points
2. Subdivide edge 4→3 into n baseline points
3. For each pair of corresponding baseline points, subdivide the connecting line into n interior points

This produces n² sub-pixels that fill the footprint quadrilateral. Temporary buffers lats_0, lons_0, lats_1, lons_1 (each length n) must be provided.

floor_indices!(ix, iy, points)

Convert the floating-point sub-pixel positions in points to integer grid cell indices via floor. Results stored in ix and iy (each length n²).

compute_n_oversample(footprint_extent, cell_size; n_min=2, n_max=20)

Estimate an appropriate oversampling factor n based on the ratio of footprint size to grid cell size. Aims for ~3 sub-pixels per grid cell width.

Returns an integer in [n_min, n_max].

compute_footprint_indices_ka!(backend, ix_out, iy_out, skip_flag,
                              lat_corners, lon_corners, n)

KA kernel that computes n×n sub-pixel grid indices for each footprint. Also sets skip_flag[fp] = true if the footprint's longitude extent >= n (too wide to oversample meaningfully).

For footprints where all 4 corners fall in the same cell, the kernel stores the single cell index in position 1 and marks all other positions with ix=-1 (sentinel). The accumulation kernel handles these fast-path footprints by checking skip_flag.

• ix_out: (n_fp, n²) Int32 — lat grid cell indices
• iy_out: (n_fp, n²) Int32 — lon grid cell indices
• skip_flag: (n_fp,) Int32 — 0=oversample, 1=fast-path (single cell), 2=skip (too wide)
• lat_corners, lon_corners: (n_fp, 4) — fractional grid indices of corners
• n: oversampling factor

accumulate_footprint!(grid_data, grid_std, grid_weights, compute_std,
                      lat_idx, lon_idx, values, n_pixels, n_vars, n_oversample,
                      grid_spec, points_buf, ix_buf, iy_buf,
                      lats0_buf, lons0_buf, lats1_buf, lons1_buf)

Accumulate a batch of satellite observations into the output grid using footprint-aware oversampling and Welford's online averaging algorithm.

For each pixel:

• Fast path: If all 4 corner coordinates fall into the same grid cell, the full observation value is assigned to that cell (weight=1).
• Oversample path: If the footprint spans multiple cells, it is subdivided into n_oversample × n_oversample sub-pixels. Each sub-pixel contributes weight 1/n² to whichever grid cell it falls in.

The Welford algorithm maintains a running mean and sum-of-squared-deviations (M2), enabling numerically stable computation of mean and variance in a single pass.

Arguments

• grid_data: Output array (n_lon, n_lat, n_vars) — running mean
• grid_std: Output array (n_lon, n_lat, n_vars) — running M2 (sum of squared deviations)
• grid_weights: Output array (n_lon, n_lat) — accumulated weights
• compute_std: Whether to track variance
• lat_idx, lon_idx: Pre-computed fractional grid indices for corners, size (n_pixels, 4)
• values: Input data, size (n_pixels, n_vars)
• n_pixels: Number of observations in this batch
• n_vars: Number of variables being gridded
• n_oversample: Sub-pixel subdivision factor (n)
• Remaining arguments are pre-allocated work buffers

accumulate_batch!(backend, grid_sum, grid_weights,
                  lat_corners, lon_corners, values, n;
                  compute_std=false, grid_std=nothing, grid_mean=nothing)

Accumulate observations into the grid using the KA kernel pipeline. This function works identically on CPU and GPU backends.

The pipeline:

1. Sort corners into CCW order
2. Compute n×n subpixel grid indices per footprint
3. Scatter-accumulate weighted sums into grid cells (using atomics)

After processing all files, call finalize_mean! to convert sums to means.

For compute_std=true, a two-pass approach is used:

1. First call accumulate_batch! without std to accumulate sums
2. Call finalize_mean! to get the mean
3. Call accumulate_std_pass! with the computed mean to accumulate variance

Arguments

• backend: KernelAbstractions backend (CPU() or CUDABackend())
• grid_sum: (nlon, nlat, n_vars) — weighted sum accumulator
• grid_weights: (nlon, nlat) — weight accumulator
• lat_corners, lon_corners: (n_fp, 4) — fractional grid indices of corners
• values: (nfp, nvars) — observation values
• n: oversampling factor

accumulate_center!(grid_data, grid_weights, lat, lon, values, grid_spec)

Accumulate observations using center coordinates only (no footprint bounds). Used for MODIS-style gridding where each pixel maps to exactly one grid cell. Simple sum accumulation — divide by grid_weights afterward to get the mean.

scatter_accumulate!(grid_sum, grid_weights, ix, iy, skip_flag, values, n, n_vars)

Sequential scatter-accumulate of weighted values into grid cells. Uses precomputed grid indices from compute_footprint_indices_ka!.

For each footprint:

• skip_flag == 1 (fast path): add weight=1, value to single cell
• skip_flag == 0 (oversample): add weight=1/n² to each subpixel's cell
• skip_flag == 2 (skip): do nothing

This function is sequential (no atomics needed) and works with any array type. For GPU backends, use the @kernel version with @atomic instead.

• grid_sum: (nlon, nlat, n_vars) — weighted sum accumulator
• grid_weights: (nlon, nlat) — weight accumulator
• ix, iy: (n_fp, n²) — precomputed grid indices
• skip_flag: (n_fp,) — 0=oversample, 1=fast-path, 2=skip
• values: (nfp, nvars) — observation values

scatter_accumulate_ka!(backend, grid_sum, grid_weights,
                       ix, iy, skip_flag, values, n, n_vars)

KA kernel version of scatter-accumulate using @atomic for thread safety. Required for GPU backends where multiple threads write to the same grid cells. For CPU backends, use scatter_accumulate! (sequential) instead.

finalize_mean!(grid_data, grid_weights)

Convert weighted-sum accumulator to mean: mean = sum / weight. Cells with near-zero weight are left as zero.

find_files(pattern::String, folder::String, date::DateTime) -> Vector{String}

Find input files matching the given pattern for a specific date.

Resolves YYYY, MM, DD, and DOY placeholders in both pattern and folder, then globs for matching files. Returns only non-empty files.

read_nc_variable(dataset, path::String; bounds=false)

Read a variable from a NetCDF dataset, handling nested groups via / separators.

If bounds=true, reshape the result so that footprint vertices are in the last dimension, returning an N×4 matrix (N soundings, 4 corners).

read_nc_attribute(dataset, path::String, attr::String)

Read an attribute from a possibly group-nested variable.

create_output_dataset(outfile, grid_spec, time_spec, grid_vars, compute_std) -> Dataset, Dict

Create the output NetCDF4 file with dimensions, coordinate variables, and data variable definitions. Returns the dataset and a Dict mapping variable names to NCDatasets variables.

apply_filters(dataset, config::DataSourceConfig, lat_bounds, lon_bounds,
              grid_spec::GridSpec) -> Vector{Int}

Apply spatial bounding box and quality filters, returning indices of pixels that pass all criteria.

The spatial filter requires all 4 corners to be within the grid bounds and the footprint to not cross the dateline (lon span < 50°). Quality filters from config.filters are applied as additional criteria.

compute_evi(red, nir, blue; G=2.5f0, C1=6.0f0, C2=7.5f0, L=10000.0f0)

Compute the Enhanced Vegetation Index (EVI).

EVI = G × (NIR - RED) / (NIR + C1×RED - C2×BLUE + L)

compute_ndvi(red, nir)

Compute the Normalized Difference Vegetation Index (NDVI).

NDVI = (NIR - RED) / (NIR + RED)

compute_nirv(red, nir)

Compute the Near-Infrared Reflectance of Vegetation (NIRv).

NIRv = NDVI × NIR

compute_ndwi(nir, swir)

Compute the Normalized Difference Water Index (NDWI).

NDWI = (NIR - SWIR) / (NIR + SWIR)

parse_l2_args(args=ARGS) -> Dict

Parse command-line arguments for L2 gridding. Matches the interface of the original gridL2_Dates.jl script.

parse_center_args(args=ARGS) -> Dict

Parse command-line arguments for center-coordinate gridding (MODIS-style).

args_to_grid_spec(args::Dict) -> GridSpec{Float32}

Convert parsed CLI arguments to a GridSpec.

args_to_time_spec(args::Dict) -> TimeSpec

Convert parsed CLI arguments to a TimeSpec.