from_stac(groupby="solar_day") — group tiled optical imagery by acquisition date¶

When you build a time series from a STAC catalog of tiled optical satellite imagery (Sentinel-2, Landsat, HLS, MODIS, ...), a single pass of the satellite over your area of interest is almost never one file. The provider chops each overpass into many granules / tiles, so a STAC search returns N items for what is really one acquisition on one date.

DatasetCollection.from_stac(..., groupby="solar_day") collapses those tiles back into one timestep per acquisition date — an analysis-ready, one-step-per-date cube — instead of N tile-timesteps. The default groupby=None keeps one timestep per item.

How it works¶

1. Solar day of each item = its UTC timestamp shifted by centroid_longitude / 15 hours, reduced to a calendar date. 15° of longitude ≈ 1 hour of local solar time, so the shift converts UTC to local solar time. This keeps one overpass on a single date instead of splitting it across the UTC-midnight boundary (a tile imaged at 23:30 UTC far to the east belongs to the next solar day).
2. Group all items that share a solar day.
3. Mosaic each group with merge_rasters(method="first") — the first non-no-data pixel wins where neighbouring tiles overlap.
4. The resulting time_length is the number of distinct solar days, in chronological order.

Single-asset only (you mosaic one asset, e.g. one band or the visual product).

A tiny offline catalog¶

from_stac is fully duck-typed: a "STAC item" is any object (or plain dict) exposing properties["datetime"], a bbox, and assets[<key>]["href"]. So we can demo the whole thing offline with three small GeoTIFFs and three hand-built item dicts — no network, no pystac.

Scenario: two adjacent tiles of the same overpass on 2024-01-15, plus one tile from the next day's overpass on 2024-01-16.

In [1]:

%matplotlib inline

import os
import tempfile
from datetime import datetime, timedelta, timezone

import numpy as np

from pyramids.dataset import Dataset, DatasetCollection

work = tempfile.mkdtemp(prefix="solar_day_demo_")
work

%matplotlib inline import os import tempfile from datetime import datetime, timedelta, timezone import numpy as np from pyramids.dataset import Dataset, DatasetCollection work = tempfile.mkdtemp(prefix="solar_day_demo_") work

2026-07-11 14:39:12 | INFO | pyramids.base.config | Logging is configured.

Out[1]:

'/tmp/solar_day_demo_z7h3u8vt'

A small helper writes a 1-band GeoTIFF tile so we can build a few overlapping scenes.

In [2]:

def make_tile(name, top_left, value):
    """Write a small 1-band GeoTIFF tile and return its path."""
    path = os.path.join(work, name)
    arr = np.full((20, 20), value, dtype="float32")
    Dataset.create_from_array(
        arr, top_left_corner=top_left, cell_size=0.05, epsg=4326,
        no_data_value=-9999.0, driver_type="GTiff", path=path,
    )
    return path

def make_tile(name, top_left, value): """Write a small 1-band GeoTIFF tile and return its path.""" path = os.path.join(work, name) arr = np.full((20, 20), value, dtype="float32") Dataset.create_from_array( arr, top_left_corner=top_left, cell_size=0.05, epsg=4326, no_data_value=-9999.0, driver_type="GTiff", path=path, ) return path

Create two tiles from the same overpass (2024-01-15) plus one from the next day.

In [3]:

# Two adjacent tiles of the SAME overpass (2024-01-15) + one NEXT-day tile.
tile_a = make_tile("a.tif", top_left=(10.0, 41.0), value=1.0)  # day 1, west
tile_b = make_tile("b.tif", top_left=(11.0, 41.0), value=2.0)  # day 1, east
tile_c = make_tile("c.tif", top_left=(10.0, 41.0), value=3.0)  # day 2
[tile_a, tile_b, tile_c]

# Two adjacent tiles of the SAME overpass (2024-01-15) + one NEXT-day tile. tile_a = make_tile("a.tif", top_left=(10.0, 41.0), value=1.0) # day 1, west tile_b = make_tile("b.tif", top_left=(11.0, 41.0), value=2.0) # day 1, east tile_c = make_tile("c.tif", top_left=(10.0, 41.0), value=3.0) # day 2 [tile_a, tile_b, tile_c]

Out[3]:

['/tmp/solar_day_demo_z7h3u8vt/a.tif',
 '/tmp/solar_day_demo_z7h3u8vt/b.tif',
 '/tmp/solar_day_demo_z7h3u8vt/c.tif']

Wrap the tiles as duck-typed STAC items¶

from_stac never calls a STAC API — it only reads three fields off each item, so a plain dict is a valid "item". Below we hand-build one per tile: two share the 2024-01-15 overpass, the third is the next day. The bbox is what drives grouping — its centroid longitude is what shifts the UTC timestamp to a solar day.

Field	Used for
properties["datetime"]	The UTC acquisition time
bbox ([minx, miny, maxx, maxy])	Centroid longitude → solar-day shift
assets[<key>]["href"]	Path/URL of the raster to mosaic

In [4]:

def item(item_id, when_iso, bbox, href):
    """A minimal duck-typed STAC item dict."""
    return {
        "id": item_id,
        "properties": {"datetime": when_iso},
        "bbox": bbox,  # [minx, miny, maxx, maxy] -> centroid lon
        "assets": {"data": {"href": href}},
    }


items = [
    item("a", "2024-01-15T10:00:00Z", [10.0, 40.0, 11.0, 41.0], tile_a),
    item("b", "2024-01-15T10:01:00Z", [11.0, 40.0, 12.0, 41.0], tile_b),
    item("c", "2024-01-16T10:00:00Z", [10.0, 40.0, 11.0, 41.0], tile_c),
]
len(items)

def item(item_id, when_iso, bbox, href): """A minimal duck-typed STAC item dict.""" return { "id": item_id, "properties": {"datetime": when_iso}, "bbox": bbox, # [minx, miny, maxx, maxy] -> centroid lon "assets": {"data": {"href": href}}, } items = [ item("a", "2024-01-15T10:00:00Z", [10.0, 40.0, 11.0, 41.0], tile_a), item("b", "2024-01-15T10:01:00Z", [11.0, 40.0, 12.0, 41.0], tile_b), item("c", "2024-01-16T10:00:00Z", [10.0, 40.0, 11.0, 41.0], tile_c), ] len(items)

Out[4]:

3

Without grouping — one timestep per item¶

Three items in, three timesteps out. The two same-day tiles stay as two separate steps.

In [5]:

cube = DatasetCollection.from_stac(items, "data")
cube.time_length  # -> 3

cube = DatasetCollection.from_stac(items, "data") cube.time_length # -> 3

Out[5]:

3

With groupby="solar_day" — one timestep per acquisition date¶

The two 2024-01-15 tiles fuse into a single mosaicked timestep; 2024-01-16 is its own step. Three items in, two timesteps out.

In [6]:

daily = DatasetCollection.from_stac(items, "data", groupby="solar_day")
daily.time_length  # -> 2 (2024-01-15, 2024-01-16)

daily = DatasetCollection.from_stac(items, "data", groupby="solar_day") daily.time_length # -> 2 (2024-01-15, 2024-01-16)

Out[6]:

2

Plot the first timestep of the grouped daily collection — the two adjacent 2024-01-15 tiles fused into one solar-day mosaic (west tile = 1, east tile = 2).

In [7]:

daily.iloc(0).plot(band=0, title="Solar-day 1 mosaic (tiles a + b)")

daily.iloc(0).plot(band=0, title="Solar-day 1 mosaic (tiles a + b)")

Out[7]:

<cleopatra.array_glyph.ArrayGlyph at 0x7f5619b1b380>

Why shift by longitude? The UTC-midnight boundary¶

A satellite's overpass happens at roughly the same local time everywhere, but the UTC timestamp of a tile depends on its longitude. Near the date line a tile can be stamped late in one UTC day while its neighbour on the same pass is stamped early the next UTC day — a naive date() would split one overpass across two timesteps. Shifting by centroid_longitude / 15 hours puts both on the same solar day. The cell below reproduces the rule pyramids uses internally:

In [8]:

def solar_day(when_utc, centroid_lon):
    """The rule pyramids uses: UTC shifted to local solar time, then the date."""
    shifted = when_utc + timedelta(hours=centroid_lon / 15.0)
    return shifted.date().isoformat()


# A tile imaged at 23:30 UTC, far to the east (lon +150): its UTC day is the 15th,
# but its LOCAL solar day is the 16th.
late = datetime(2024, 1, 15, 23, 30, tzinfo=timezone.utc)
print("UTC date  :", late.date().isoformat())  # 2024-01-15
print("solar day :", solar_day(late, centroid_lon=150.0))  # 2024-01-16

def solar_day(when_utc, centroid_lon): """The rule pyramids uses: UTC shifted to local solar time, then the date.""" shifted = when_utc + timedelta(hours=centroid_lon / 15.0) return shifted.date().isoformat() # A tile imaged at 23:30 UTC, far to the east (lon +150): its UTC day is the 15th, # but its LOCAL solar day is the 16th. late = datetime(2024, 1, 15, 23, 30, tzinfo=timezone.utc) print("UTC date :", late.date().isoformat()) # 2024-01-15 print("solar day :", solar_day(late, centroid_lon=150.0)) # 2024-01-16

UTC date  : 2024-01-15
solar day : 2024-01-16

Where to use it — and where not¶

Use groupby="solar_day" when:

• The catalog is tiled optical EO (Sentinel-2 / Landsat / HLS / MODIS) and your AOI spans several tiles/granules per overpass.
• You want one analysis-ready timestep per acquisition date (e.g. a daily NDVI/true-colour stack) rather than many tile-steps.
• Same-overpass tiles overlap and a first-valid mosaic is acceptable.

Do not use it when:

• The data is not per-overpass imagery — climate-model output, reanalysis, or already-mosaicked / single-granule products. There groupby=None is correct.
• You need a reducer other than first-valid (median compositing, cloud-score weighting, BRDF normalisation): do that in your EO pipeline — see the recipe below.
• Items lack a datetime or bbox (the solar-day computation needs both).

Under the hood — the equivalent manual recipe¶

groupby="solar_day" is a convenience over two generic pieces: group the item hrefs by solar day, then merge_rasters each group. If you need a different reducer (median compositing, cloud-score weighting) just adapt this — group with the same key and swap the merge step.

In [9]:

from collections import defaultdict

from pyramids.dataset.merge import merge_rasters

def group_by_solar_day(items, asset, out_dir):
    """One mosaicked GeoTIFF per solar day — the building blocks groupby uses."""
    groups = defaultdict(list)
    for it in items:
        when = datetime.fromisoformat(it['properties']['datetime'].replace('Z', '+00:00'))
        minx, _, maxx, _ = it['bbox']
        day = solar_day(when, centroid_lon=(minx + maxx) / 2.0)
        groups[day].append(it['assets'][asset]['href'])
    paths = []
    for day in sorted(groups):
        out = os.path.join(out_dir, f'{day}.tif')
        merge_rasters(groups[day], out, method='first')  # swap for your own reducer
        paths.append(out)
    return paths

from collections import defaultdict from pyramids.dataset.merge import merge_rasters def group_by_solar_day(items, asset, out_dir): """One mosaicked GeoTIFF per solar day — the building blocks groupby uses.""" groups = defaultdict(list) for it in items: when = datetime.fromisoformat(it['properties']['datetime'].replace('Z', '+00:00')) minx, _, maxx, _ = it['bbox'] day = solar_day(when, centroid_lon=(minx + maxx) / 2.0) groups[day].append(it['assets'][asset]['href']) paths = [] for day in sorted(groups): out = os.path.join(out_dir, f'{day}.tif') merge_rasters(groups[day], out, method='first') # swap for your own reducer paths.append(out) return paths

Run the recipe: group the items by solar day, mosaic each day, and read the result back.

In [10]:

out_dir = tempfile.mkdtemp(prefix="solar_day_manual_")
day_paths = group_by_solar_day(items, "data", out_dir)
manual = DatasetCollection.from_files(day_paths)
manual.time_length  # -> 2, same as groupby='solar_day'

out_dir = tempfile.mkdtemp(prefix="solar_day_manual_") day_paths = group_by_solar_day(items, "data", out_dir) manual = DatasetCollection.from_files(day_paths) manual.time_length # -> 2, same as groupby='solar_day'

Out[10]:

2

The manual recipe yields time_length == 2 — identical to groupby="solar_day". That confirms the built-in flag is just this group-then-merge_rasters pattern behind one keyword; swapping the method="first" merge for your own reducer (median compositing, cloud-score weighting) is the extension point when first-valid mosaicking is not enough.