Merge & stack

Free functions for combining multiple rasters into one.

Hold "Ctrl" to enable pan & zoom

flowchart LR
    M["<b>pyramids.dataset.merge</b>"]
    M --> MR["<b>merge_rasters</b><br/>mosaic overlapping rasters<br/>into one (reproject on the fly)"]
    M --> SB["<b>stack_bands</b><br/>stack single-band rasters<br/>into one multi-band Dataset"]

• merge_rasters — mosaic several (overlapping or adjacent) rasters into a single raster covering their union.
• stack_bands — stack several single-band rasters into one multi-band raster.

See the Mosaic & merge notebook for runnable examples.

pyramids.dataset.merge.merge_rasters(src, dst, no_data_value=INHERIT_NO_DATA, init='nan', n='nan', method='last', dst_crs=None, resampling=DEFAULT_RESAMPLING, signer=None, *, bbox=None, bbox_crs=None)

Merge a group of rasters into one raster, resolving overlaps by method.

The overlap-resolution method selects how overlapping pixels are combined:

• "last" (default) / "first" — z-order compositing: the last (or first) source covering a pixel wins. Implemented cheaply with :func:gdal.BuildVRT + :func:gdal.Translate.
• "min" / "max" / "sum" / "count" / "mean" — per-pixel reduction across every source overlapping that pixel, ignoring no-data. Each source is aligned onto the union grid and folded into a NaN-aware accumulator one strip at a time, so peak memory is O(strip) whichever reduction is asked for. "count" reports how many sources contributed a valid value to each cell, and "mean" divides the total by that same count -- never by the number of sources merely overlapping the cell.

"median" is recognised but not implemented: it would have to hold every overlapping source's value for a cell at once, which the one-accumulator strip reduction cannot do, and asking for it raises :class:NotImplementedError rather than returning an approximation.

Parameters:

Name	Type	Description	Default
src	Sequence[str | Path]	Paths to all input rasters.	required
dst	str | Path	Path to the output raster. Its extension alone selects the output driver (.tif -> GTiff, .nc -> netCDF, …) — the same resolution every other write path in the package uses — so one dst yields the same format for every method. COMPRESS=LZW is a GTiff creation option and is applied only when the extension resolves to GTiff; other formats are written with their driver defaults. A write-by-copy-only format is refused for every method -- that is .png, .jpg / .jpeg, .jp2 / .j2k and .asc, and .vrt on top (a VRT writes a reference, not a raster). The z-order path could produce several of them, since it writes via gdal.Translate, and the reduction path could not; letting method decide what dst may be is the same defect as letting it decide the format, so both take the stricter answer. .asc is the one this costs: it was writable before, through the z-order path only. Write a GTiff and convert.	required
no_data_value	float | int | str | None	Stamped on the output bands as the nodata marker. Omitted means inherit from the sources: the first value they declare wins, and a disagreement warns. Passing a value explicitly overrides that -- but note that any real cell holding it becomes unreadable, which is why 0 is a poor choice for an elevation, bathymetry, anomaly or difference raster (#1086). Passing None explicitly asks for no marker at all, on either path, and the string "inherit" -- the word the default renders as -- asks for the default. A wrapper that forwards an "unset" of its own can import the sentinel itself as from pyramids.dataset.merge import INHERIT_NO_DATA. Whichever marker is settled on fills the pixels no source covers too, on both write paths and however it was settled, so the gaps really hold what the band declares. init keeps those pixels only when you name it yourself: on Int32 tiles no_data_value=-1 alone writes gaps of -1 under a declared -1, while no_data_value=-1, init=7 writes gaps of 7 under a declared -1, which masks nothing (both measured). The exception is a value the output band cannot store — it warns, GDAL drops it, and init keeps the gaps because there is no storable marker to fill them with. The value is read from each source's band 1, and one marker is stamped on every band of the output. It is read from every source, including one BuildVRT goes on to drop for a band type that disagrees with the first's -- so on such a mosaic a source that contributes no pixel can still decide the marker, or raise the disagreement warning. That is a symptom of the drop, which GDAL reports separately, rather than of the inheritance. A multi-band merge therefore keeps band 1's answer throughout: where band 1 declares nothing and band 2 declares a marker, nothing is inherited, and where the bands declare different markers, band 1's is stamped over band 2's. For a GeoTIFF that costs nothing, because the format has no per-band marker to lose -- TIFFTAG_GDAL_NODATA holds one value for the whole dataset, and GDAL warns that setting a second "will be used for all bands on re-opening". It is observable only for sources that do carry one per band, such as a VRT. When no source declares one, a marker is chosen wherever there is something for it to mark: a mosaic generally has pixels no source covers, and leaving those undeclared makes them read as real data -- a literal 0 on an integer mosaic, which then poisons :meth:~pyramids.dataset.Dataset.stats. The reduction methods write Float64 and take NaN, which no real cell can hold, whatever the sources' footprints — except method="count", whose uncovered cells hold a real 0 and which therefore declares no marker at all, since any value a tally can hold would mask real counts (see method). The z-order methods take the value those pixels would otherwise hold, init — so NaN too on a floating mosaic, the common case. Where the band cannot store NaN, the footprints decide. Sources that tile their whole area leave no pixel to mark, so nothing is declared and no pixel is read to settle it: two contiguous Int32 tiles declaring nothing come back declaring nothing (measured). Sources that leave a gap earn a sentinel the dtype can store and the mosaic's own cells do not use -- a numeric init first, then -9999 where the dtype holds it, then the dtype's own extremes, an unsigned band offering its maximum before its 0, since 0 is the likelier real observation. A gapped signed integer mosaic of sources declaring nothing therefore comes back declaring -9999 and an unsigned one its maximum -- 65535 for UInt16 -- with the gaps holding it either way. Only a gapped mosaic whose data uses every value its dtype could spare is written without a marker, and that warns. An inherited value the output band cannot store -- a NaN from integer sources, which GDAL both accepts and reports back -- would be dropped by :func:gdal.Translate and leave the mosaic unmarked, so it warns and gives way to a storable one.	INHERIT_NO_DATA
init	float | int | str	Reported value for pixels with no source coverage in the VRT (z-order methods only — the reduction methods never read it). Maps to :func:gdal.BuildVRTOptions VRTNodata. It is the fallback, not the last word: whatever marker gets settled on is filled into the gaps instead, so init is what survives only when that marker turns out to be init itself. Sources declaring -9999 and an init=5 yield gaps holding -9999 (measured), and so does an explicit no_data_value=-9999. The default "nan" is therefore what a floating mosaic inheriting nothing ends up declaring, and on an integer mosaic — which cannot store it — the storable sentinel described under no_data_value takes its place. init has the gaps to itself only when you name it, which is the case that needs it to agree with the marker. A numeric init is also the sentinel search's preferred candidate, so an integer mosaic inheriting nothing declares init itself when its own cells do not hold that value. One they do hold is passed over -- init=0 over tiles containing real zeros is 1086 through a different argument -- and that warns. It must name a number, or be a spelling of NaN; anything else is refused up front, since whether it would have reached GDAL at all depends on the output's data type and on whether the sources leave a gap.	'nan'
n	float | int | str	Source pixels matching this value are ignored — both when building the VRT mosaic (z-order, as srcNodata) and when warping onto the union grid (the reduction methods, likewise). It is an override, applied to every source alike, and it replaces each source's own declaration rather than adding to it: over tiles declaring -9999 and -32768, n=-9999 composites the second tile's holes as real -32768 measurements (measured). The default "nan" therefore means "no override" — the spelling is matched, so "NaN" and float("nan") mean it too — and each source's own declared value marks that source's own holes, which is what lets the gaps and the inherited marker agree.	'nan'
method	str	Overlap-resolution rule: one of "first", "last" (default), "min", "max", "sum", "count", "mean". "count" writes the number of sources that contributed a valid value to each cell -- a source's own no-data holes are not counted, so a cell covered by three sources where one declares a hole counts 2. Its uncovered cells hold a real 0 rather than the no-data fill every other reduction writes there: zero contributions is the honest tally, no covered cell can hold it, and the alternative would be to report -9999 contributions. A count mosaic therefore declares no no-data value at all, whatever was inherited or passed: every cell holds a tally in 0..len(src_paths), and a marker inside that range masks the cells really holding that number. The case that forces it is 0 — the standard no-data of Sentinel-2, Landsat Collection 2 and most integer EO products, and the very tally every uncovered cell holds, so inheriting it would mask exactly the gaps this method reports honestly. An explicitly passed no_data_value is dropped here too, with a warning. "mean" divides each cell's total by that same valid count, so a cell covered by two of three sources averages the two values it holds. Its uncovered cells take the no-data fill like "min", "max" and "sum", there being no mean of nothing. "median" raises :class:NotImplementedError: see the summary above.	'last'
dst_crs	int | str | None	Target CRS for the mosaic, as an EPSG code (32632) or any GDAL-parseable CRS string ("EPSG:32632", a WKT, a PROJ string). Each source whose CRS differs from the target is reprojected onto it (via :func:gdal.Warp) before compositing, so tiles in different CRSs — e.g. a Sentinel-2 AOI straddling two UTM zones — mosaic correctly. None (default) keeps the previous behaviour: sources are assumed to share a CRS and are composited as-is, except when they are found to disagree, in which case they are reprojected onto the first source's CRS. Reprojection always happens before the BuildVRT/Warp compositing step because that step has no reprojection capability and assumes a single shared grid.	None
resampling	str	Resampling method used when a source is reprojected to dst_crs (or to the common CRS on auto-detect). Case-insensitive; any key of :data:pyramids.base._utils.INTERPOLATION_METHODS: "nearest" (alias "nearest neighbor", the default), "bilinear", "cubic", "cubic_spline", "lanczos", "average", "mode", "max", "min", "med", "q1", "q3", "sum", and "rms". Prefer "bilinear"/"cubic" for continuous data (reflectance, DEM) to avoid the blockiness nearest introduces across reprojection. Ignored when no source is reprojected.	DEFAULT_RESAMPLING
signer	Any	Optional signer exposing sign_href(str) -> str and gdal_env() -> dict[str, str] (e.g. a :class:pyramids.stac.signers.Signer). When given, both hooks are applied — exactly as :func:pyramids.stac.load_asset does: signer.sign_href rewrites every source path first (e.g. grafting a SAS token onto a blob URL), then signer.gdal_env() is installed via :class:~pyramids.base.remote.CloudConfig for the duration of the merge. This means URL-signing signers (Planetary Computer SAS, whose credential rides the href and whose gdal_env() is empty) and env-based signers (Requester-Pays, bearer) both authenticate without wrapping the call in a with CloudConfig(...) block. None (default) leaves source hrefs untouched and installs no extra config.	None
bbox	Sequence[float] | None	Optional (west, south, east, north) window to restrict the merge to. None (default) merges the full extent of every source, which is what the function has always done. This is not a convenience for cropping afterwards: without it GDAL is given no reason to read less, so a mosaic of remote sources pulls the entire source extent through /vsicurl even when the caller wants a fraction of it. A full Sentinel-2 tile is 10980x10980 px (verified against a public Earth Search COG), so an area of interest covering a fraction of one tile still costs the whole tile without a window. How much that saves depends on the ratio between the source footprint and the window. The window is resolved once, onto the mosaic's own pixel grid, and used by both methods: z-order passes it to :func:gdal.Translate as projWin, so only the byte ranges the window touches are requested; the reduction methods clip the union grid to it. Snapping outward onto whole pixels keeps the result a strict sub-grid, and resolving it once means both methods return the same grid for the same arguments. A bbox must be ordered west < east and south < north. A window crossing the antimeridian (west > east) is rejected rather than silently reinterpreted as its own complement; split it and merge the two halves. :meth:pyramids.dataset.Dataset.crop handles the seam directly, but a mosaic is composited on one grid, which a seam-crossing window would not have. For a lon/lat mosaic the window is rewritten into the mosaic's own longitude convention first, so a -180..180 bbox reads correctly against a 0..360 grid (and the reverse). A window that ends up spanning that convention's seam is rejected on the same grounds as an antimeridian one. A window extending past the mosaic is clipped to it; one that misses it entirely raises.	None
bbox_crs	int | str | None	CRS that bbox is expressed in — an EPSG code (4326), an authority string ("EPSG:4326"), a WKT, or anything :meth:pyproj.CRS.from_user_input accepts. None (default) means bbox is already in the mosaic's own CRS — which, when dst_crs is given, is dst_crs, since sources are reprojected before the window is applied. Ignored when bbox is None. Named to match :meth:pyramids.dataset.engines.cog.COG.read_part.	None

Returns:

Type	Description
None	None

Note

The z-order methods ("first"/"last") preserve the data type via BuildVRT + Translate — specifically that of the first source in z-order, which for method="first" is the last one you passed, since the list is reversed so that it wins. That is the source BuildVRT takes its band type from; one whose band type disagrees is skipped entirely, with a GDAL warning naming it, so the output then covers the remaining sources' extent rather than the union. The reduction methods ("min"/"max"/"sum"/"count"/"mean") align every source onto the union grid with gdal.Warp (nearest resampling — exact for already-aligned tiles) and write a single-precision-safe Float64 output regardless of the source dtype, so they may differ in dtype from a z-order merge of the same integer inputs. "count" is Float64 too, rather than an integer band of its own: one reduction path writes one dtype, and a tally is exact in Float64 far past any plausible source count. The warp itself is asked for Float64 as well, so that the NaN marking a strip's uncovered pixels survives it: rounded into an integer destination it would come back as a real 0 (or the type's minimum for Int32) and be folded into the result as data.

Raises:

Type	Description
NotImplementedError	method="median" — recognised, deliberately not implemented, and refused rather than approximated.
TypeError	resampling is not a string, or bbox is not four numbers (a string, a scalar, or a sequence holding a non-numeric element).
ValueError	method/resampling is not a supported value, dst_crs cannot be parsed as a CRS, a source carries no CRS, no_data_value, init or n names no number, or bbox is malformed (wrong length, non-finite, inverted, zero-area), crosses the antimeridian or the mosaic's longitude seam once reprojected, selects no whole pixel, does not overlap the mosaic, or cannot be projected into its CRS.
RuntimeError	GDAL failed to open a source, reproject it, or build the source mosaic. When a source is at fault the message names it and its position in src, and chains GDAL's own error; any credential in a signed URL is redacted.
DriverNotExistError	dst has no extension, or one the driver catalog does not know.
FileFormatNotSupportedError	dst's extension maps to a write-by-copy-only format, for any method.

Warns:

Type

Description

UserWarning

The sources declare more than one distinct no-data value (the first wins, and all of them are named); or the value they declare is one the output's data type cannot store, so a storable one is chosen in its place; or a numeric init was passed over as the marker because the mosaic's own cells hold that value, so the gaps hold the chosen sentinel instead; or the mosaic's cells use every value its dtype could spare as a marker, so it is written without one. Those four arise only while the marker is being inherited. Passing no_data_value explicitly silences them and raises one of its own where the value does not fit the output band: GDAL then drops it and the mosaic carries no marker at all. An explicit no_data_value under method="count" warns too, and is likewise dropped — see method for why a count mosaic declares nothing.

Examples:

• Mosaic two tiles, keeping the larger value wherever they overlap:

  >>> from pyramids.dataset.merge import merge_rasters
  >>> merge_rasters(  # doctest: +SKIP
  ...     ["tile_a.tif", "tile_b.tif"],
  ...     "mosaic_max.tif",
  ...     no_data_value=-9999.0,
  ...     method="max",
  ... )
  

• Composite the average of every source covering a pixel, and write the matching per-pixel tally of how many of them contributed it:

  >>> from pyramids.dataset.merge import merge_rasters
  >>> merge_rasters(  # doctest: +SKIP
  ...     ["scene_a.tif", "scene_b.tif", "scene_c.tif"],
  ...     "mosaic_mean.tif",
  ...     method="mean",
  ... )
  >>> merge_rasters(  # doctest: +SKIP
  ...     ["scene_a.tif", "scene_b.tif", "scene_c.tif"],
  ...     "mosaic_count.tif",
  ...     method="count",
  ... )
  

• Default last-wins compositing. The overlap rule is unchanged, but the mosaic's no-data marker is now inherited from the sources rather than set to 0 (#1086):

  >>> merge_rasters(["tile_a.tif", "tile_b.tif"], "mosaic.tif")  # doctest: +SKIP
  

• Override the inherited marker on a z-order merge of integer tiles. The pixels no tile covers are filled with -1 as well, so the marker the band declares is the one they hold:

  >>> merge_rasters(  # doctest: +SKIP
  ...     ["tile_a.tif", "tile_b.tif"],
  ...     "mosaic.tif",
  ...     no_data_value=-1,
  ... )
  

• Keep a different value in the uncovered pixels by naming init. The band still declares -1, but the gaps hold 0 and stay readable:

  >>> merge_rasters(  # doctest: +SKIP
  ...     ["tile_a.tif", "tile_b.tif"],
  ...     "mosaic.tif",
  ...     no_data_value=-1,
  ...     init=0,
  ... )
  

• Ask for no marker at all, so every value in the mosaic stays readable — including a real 0 m in an elevation model:

  >>> merge_rasters(  # doctest: +SKIP
  ...     ["dem_a.tif", "dem_b.tif"],
  ...     "dem.tif",
  ...     no_data_value=None,
  ... )
  

• Mosaic tiles from two UTM zones into a single CRS:

  >>> merge_rasters(  # doctest: +SKIP
  ...     ["utm32_tile.tif", "utm33_tile.tif"],
  ...     "mosaic_utm32.tif",
  ...     dst_crs=32632,
  ... )
  

• Mosaic Requester-Pays S3 tiles by passing a signer (no with block):

  >>> from pyramids.stac import AWSRequesterPaysSigner  # doctest: +SKIP
  >>> merge_rasters(  # doctest: +SKIP
  ...     ["s3://bucket/a.tif", "s3://bucket/b.tif"],
  ...     "mosaic.tif",
  ...     signer=AWSRequesterPaysSigner(region="us-west-2"),
  ... )
  

Source code in src/pyramids/dataset/merge.py

def merge_rasters(
    src: Sequence[str | Path],
    dst: str | Path,
    no_data_value: Any = INHERIT_NO_DATA,
    init: float | int | str = "nan",
    n: float | int | str = "nan",
    method: str = "last",
    dst_crs: int | str | None = None,
    resampling: str = DEFAULT_RESAMPLING,
    signer: Any = None,
    *,
    bbox: Sequence[float] | None = None,
    bbox_crs: int | str | None = None,
) -> None:
    """Merge a group of rasters into one raster, resolving overlaps by ``method``.

    The overlap-resolution ``method`` selects how overlapping pixels are
    combined:

    * ``"last"`` (default) / ``"first"`` — z-order compositing: the last (or
      first) source covering a pixel wins. Implemented cheaply with
      :func:`gdal.BuildVRT` + :func:`gdal.Translate`.
    * ``"min"`` / ``"max"`` / ``"sum"`` / ``"count"`` / ``"mean"`` — per-pixel
      reduction across every source overlapping that pixel, ignoring no-data.
      Each source is aligned onto the union grid and folded into a NaN-aware
      accumulator one strip at a time, so peak memory is ``O(strip)`` whichever
      reduction is asked for. ``"count"`` reports how many sources contributed a
      *valid* value to each cell, and ``"mean"`` divides the total by that same
      count -- never by the number of sources merely overlapping the cell.

    ``"median"`` is recognised but not implemented: it would have to hold every
    overlapping source's value for a cell at once, which the one-accumulator strip
    reduction cannot do, and asking for it raises :class:`NotImplementedError`
    rather than returning an approximation.

    Args:
        src (Sequence[str | Path]):
            Paths to all input rasters.
        dst (str | Path):
            Path to the output raster. Its extension alone selects the output
            driver (`.tif` -> GTiff, `.nc` -> netCDF, …) — the same
            resolution every other write path in the package uses — so one
            `dst` yields the same format for every `method`. `COMPRESS=LZW`
            is a GTiff creation option and is applied only when the extension
            resolves to GTiff; other formats are written with their driver
            defaults. A write-by-copy-only format is refused for every
            `method` -- that is `.png`, `.jpg` / `.jpeg`, `.jp2` / `.j2k` and
            `.asc`, and `.vrt` on top (a VRT writes a reference, not a
            raster). The z-order path could produce several of them, since it
            writes via `gdal.Translate`, and the reduction path could not;
            letting `method` decide what `dst` may be is the same defect as
            letting it decide the format, so both take the stricter answer.
            `.asc` is the one this costs: it was writable before, through the
            z-order path only. Write a GTiff and convert.
        no_data_value (float | int | str | None):
            Stamped on the output bands as the nodata marker. Omitted means
            **inherit from the sources**: the first value they declare wins, and
            a disagreement warns. Passing a value explicitly overrides that --
            but note that any real cell holding it becomes unreadable, which is
            why 0 is a poor choice for an elevation, bathymetry, anomaly or
            difference raster (#1086). Passing ``None`` explicitly asks for no
            marker at all, on either path, and the string ``"inherit"`` -- the
            word the default renders as -- asks for the default. A wrapper that
            forwards an "unset" of its own can import the sentinel itself as
            ``from pyramids.dataset.merge import INHERIT_NO_DATA``.

            **Whichever marker is settled on fills the pixels no source covers
            too**, on both write paths and however it was settled, so the gaps
            really hold what the band declares. `init` keeps those pixels only
            when you name it yourself: on Int32 tiles ``no_data_value=-1`` alone
            writes gaps of ``-1`` under a declared ``-1``, while
            ``no_data_value=-1, init=7`` writes gaps of ``7`` under a declared
            ``-1``, which masks nothing (both measured). The exception is a
            value the output band cannot store — it warns, GDAL drops it, and
            `init` keeps the gaps because there is no storable marker to fill
            them with.

            The value is read from each source's **band 1**, and one marker is
            stamped on every band of the output. It is read from *every* source,
            including one ``BuildVRT`` goes on to drop for a band type that
            disagrees with the first's -- so on such a mosaic a source that
            contributes no pixel can still decide the marker, or raise the
            disagreement warning. That is a symptom of the drop, which GDAL
            reports separately, rather than of the inheritance. A multi-band merge therefore
            keeps band 1's answer throughout: where band 1 declares nothing and
            band 2 declares a marker, nothing is inherited, and where the bands
            declare different markers, band 1's is stamped over band 2's. For a
            GeoTIFF that costs nothing, because the format has no per-band marker
            to lose -- ``TIFFTAG_GDAL_NODATA`` holds one value for the whole
            dataset, and GDAL warns that setting a second "will be used for all
            bands on re-opening". It is observable only for sources that do carry
            one per band, such as a VRT.

            When **no source declares one**, a marker is chosen wherever there
            is something for it to mark: a mosaic generally has pixels no source
            covers, and leaving those undeclared makes them read as real data --
            a literal `0` on an integer mosaic, which then poisons
            :meth:`~pyramids.dataset.Dataset.stats`. The reduction methods
            write Float64 and take ``NaN``, which no real cell can hold,
            whatever the sources' footprints — except ``method="count"``, whose
            uncovered cells hold a real ``0`` and which therefore declares no
            marker at all, since any value a tally can hold would mask real
            counts (see `method`). The z-order methods take the value
            those pixels would otherwise hold, `init` — so ``NaN`` too on a
            floating mosaic, the common case. Where the band cannot store
            ``NaN``, the footprints decide. Sources that **tile their whole
            area** leave no pixel to mark, so nothing is declared and no pixel
            is read to settle it: two contiguous Int32 tiles declaring nothing
            come back declaring nothing (measured). Sources that **leave a gap**
            earn a sentinel the dtype can store and the mosaic's own cells do
            not use -- a numeric `init` first, then ``-9999`` where the dtype
            holds it, then the dtype's own extremes, an unsigned band offering
            its maximum before its `0`, since `0` is the likelier real
            observation. A gapped signed integer mosaic of sources declaring
            nothing therefore comes back declaring ``-9999`` and an unsigned one
            its maximum -- ``65535`` for ``UInt16`` -- with the gaps holding it
            either way. Only a gapped mosaic whose data uses every value its
            dtype could spare is written without a marker, and that warns.

            An inherited value the output band cannot store -- a ``NaN`` from
            integer sources, which GDAL both accepts and reports back -- would
            be dropped by :func:`gdal.Translate` and leave the mosaic unmarked,
            so it warns and gives way to a storable one.
        init (float | int | str):
            Reported value for pixels with no source coverage in the VRT (z-order
            methods only — the reduction methods never read it). Maps to
            :func:`gdal.BuildVRTOptions` ``VRTNodata``.

            It is the fallback, not the last word: whatever marker gets settled
            on is filled into the gaps instead, so `init` is what survives only
            when that marker turns out to be `init` itself. Sources declaring
            ``-9999`` and an ``init=5`` yield gaps holding ``-9999`` (measured),
            and so does an explicit ``no_data_value=-9999``. The default
            ``"nan"`` is therefore what a floating mosaic inheriting nothing
            ends up declaring, and on an integer mosaic — which cannot store it
            — the storable sentinel described under `no_data_value` takes its
            place. `init` has the gaps to itself only when you name it, which is
            the case that needs it to agree with the marker.

            A **numeric** `init` is also the sentinel search's preferred
            candidate, so an integer mosaic inheriting nothing declares `init`
            itself when its own cells do not hold that value. One they do hold
            is passed over -- ``init=0`` over tiles containing real zeros is
            #1086 through a different argument -- and that warns. It must name a
            number, or be a spelling of ``NaN``; anything else is refused up
            front, since whether it would have reached GDAL at all depends on
            the output's data type and on whether the sources leave a gap.
        n (float | int | str):
            Source pixels matching this value are ignored — both when building
            the VRT mosaic (z-order, as ``srcNodata``) and when warping onto the
            union grid (the reduction methods, likewise). It is an **override**,
            applied to every source alike, and it *replaces* each source's own
            declaration rather than adding to it: over tiles declaring ``-9999``
            and ``-32768``, ``n=-9999`` composites the second tile's holes as
            real ``-32768`` measurements (measured). The default ``"nan"``
            therefore means "no override" — the spelling is matched, so
            ``"NaN"`` and ``float("nan")`` mean it too — and each source's own
            declared value marks that source's own holes, which is what lets the
            gaps and the inherited marker agree.
        method (str):
            Overlap-resolution rule: one of ``"first"``, ``"last"`` (default),
            ``"min"``, ``"max"``, ``"sum"``, ``"count"``, ``"mean"``.

            ``"count"`` writes the number of sources that contributed a valid
            value to each cell -- a source's own no-data holes are not counted,
            so a cell covered by three sources where one declares a hole counts
            ``2``. Its **uncovered cells hold a real ``0``** rather than the
            no-data fill every other reduction writes there: zero contributions
            is the honest tally, no covered cell can hold it, and the alternative
            would be to report ``-9999`` contributions. A count mosaic therefore
            **declares no no-data value at all**, whatever was inherited or
            passed: every cell holds a tally in ``0..len(src_paths)``, and a
            marker inside that range masks the cells really holding that number.
            The case that forces it is ``0`` — the standard no-data of
            Sentinel-2, Landsat Collection 2 and most integer EO products, and
            the very tally every uncovered cell holds, so inheriting it would
            mask exactly the gaps this method reports honestly. An explicitly
            passed ``no_data_value`` is dropped here too, with a warning.

            ``"mean"`` divides each cell's total by that same valid count, so a
            cell covered by two of three sources averages the two values it
            holds. Its uncovered cells take the no-data fill like ``"min"``,
            ``"max"`` and ``"sum"``, there being no mean of nothing.

            ``"median"`` raises :class:`NotImplementedError`: see the summary
            above.
        dst_crs (int | str | None):
            Target CRS for the mosaic, as an EPSG code (``32632``) or any
            GDAL-parseable CRS string (``"EPSG:32632"``, a WKT, a PROJ string).
            Each source whose CRS differs from the target is reprojected onto it
            (via :func:`gdal.Warp`) **before** compositing, so tiles in different
            CRSs — e.g. a Sentinel-2 AOI straddling two UTM zones — mosaic
            correctly. ``None`` (default) keeps the previous behaviour: sources
            are assumed to share a CRS and are composited as-is, *except* when
            they are found to disagree, in which case they are reprojected onto
            the first source's CRS. Reprojection always happens before the
            ``BuildVRT``/``Warp`` compositing step because that step has no
            reprojection capability and assumes a single shared grid.
        resampling (str):
            Resampling method used when a source is reprojected to ``dst_crs``
            (or to the common CRS on auto-detect). Case-insensitive; any key
            of :data:`pyramids.base._utils.INTERPOLATION_METHODS`: ``"nearest"``
            (alias ``"nearest neighbor"``, the default), ``"bilinear"``,
            ``"cubic"``, ``"cubic_spline"``, ``"lanczos"``, ``"average"``,
            ``"mode"``, ``"max"``, ``"min"``, ``"med"``, ``"q1"``, ``"q3"``,
            ``"sum"``, and ``"rms"``.
            Prefer ``"bilinear"``/``"cubic"`` for continuous data (reflectance,
            DEM) to avoid the blockiness nearest introduces across reprojection.
            Ignored when no source is reprojected.
        signer (Any):
            Optional signer exposing ``sign_href(str) -> str`` and
            ``gdal_env() -> dict[str, str]`` (e.g. a
            :class:`pyramids.stac.signers.Signer`). When given, **both** hooks
            are applied — exactly as :func:`pyramids.stac.load_asset` does:
            ``signer.sign_href`` rewrites every source path first (e.g. grafting
            a SAS token onto a blob URL), then ``signer.gdal_env()`` is installed
            via :class:`~pyramids.base.remote.CloudConfig` for the duration of
            the merge. This means URL-signing signers (Planetary Computer SAS,
            whose credential rides the href and whose ``gdal_env()`` is empty)
            and env-based signers (Requester-Pays, bearer) both authenticate
            without wrapping the call in a ``with CloudConfig(...)`` block.
            ``None`` (default) leaves source hrefs untouched and installs no
            extra config.
        bbox (Sequence[float] | None):
            Optional ``(west, south, east, north)`` window to restrict the merge
            to. ``None`` (default) merges the full extent of every source, which
            is what the function has always done.

            This is not a convenience for cropping afterwards: without it GDAL is
            given no reason to read less, so a mosaic of remote sources pulls the
            **entire** source extent through ``/vsicurl`` even when the caller
            wants a fraction of it. A full Sentinel-2 tile is 10980x10980 px
            (verified against a public Earth Search COG), so an area of interest
            covering a fraction of one tile still costs the whole tile without a
            window. How much that saves depends on the ratio between the source
            footprint and the window.

            The window is resolved once, onto the mosaic's own pixel grid, and
            used by both methods: z-order passes it to :func:`gdal.Translate` as
            ``projWin``, so only the byte ranges the window touches are requested;
            the reduction methods clip the union grid to it. Snapping outward onto
            whole pixels keeps the result a strict sub-grid, and resolving it once
            means both methods return the same grid for the same arguments.

            A bbox must be ordered ``west < east`` and ``south < north``. A window
            crossing the antimeridian (``west > east``) is rejected rather than
            silently reinterpreted as its own complement; split it and merge the
            two halves. :meth:`pyramids.dataset.Dataset.crop` handles the seam
            directly, but a mosaic is composited on one grid, which a seam-crossing
            window would not have.

            For a lon/lat mosaic the window is rewritten into the mosaic's own
            longitude convention first, so a ``-180..180`` bbox reads correctly
            against a ``0..360`` grid (and the reverse). A window that ends up
            spanning that convention's seam is rejected on the same grounds as an
            antimeridian one. A window extending past the mosaic is clipped to it;
            one that misses it entirely raises.
        bbox_crs (int | str | None):
            CRS that ``bbox`` is expressed in — an EPSG code (``4326``), an
            authority string (``"EPSG:4326"``), a WKT, or anything
            :meth:`pyproj.CRS.from_user_input` accepts. ``None`` (default) means
            ``bbox`` is already in the mosaic's own CRS — which, when ``dst_crs``
            is given, is ``dst_crs``, since sources are reprojected before the
            window is applied. Ignored when ``bbox`` is ``None``. Named to match
            :meth:`pyramids.dataset.engines.cog.COG.read_part`.

    Returns:
        None

    Note:
        The z-order methods (``"first"``/``"last"``) preserve the data type via
        ``BuildVRT`` + ``Translate`` — specifically that of the first source **in
        z-order**, which for ``method="first"`` is the *last* one you passed,
        since the list is reversed so that it wins. That is the source
        ``BuildVRT`` takes its band type from; one whose band type disagrees is
        skipped entirely, with a GDAL warning naming it, so the output then
        covers the remaining sources' extent rather than the union. The
        reduction methods (``"min"``/``"max"``/``"sum"``/``"count"``/``"mean"``)
        align every source onto
        the union grid with ``gdal.Warp`` (nearest resampling — exact for
        already-aligned tiles) and write a single-precision-safe **Float64**
        output regardless of the source dtype, so they may differ in dtype from a
        z-order merge of the same integer inputs. ``"count"`` is Float64 too,
        rather than an integer band of its own: one reduction path writes one
        dtype, and a tally is exact in Float64 far past any plausible source
        count. The warp itself is asked for
        Float64 as well, so that the ``NaN`` marking a strip's uncovered pixels
        survives it: rounded into an integer destination it would come back as a
        real ``0`` (or the type's minimum for Int32) and be folded into the
        result as data.

    Raises:
        NotImplementedError: ``method="median"`` — recognised, deliberately not
            implemented, and refused rather than approximated.
        TypeError: ``resampling`` is not a string, or ``bbox`` is not four numbers
            (a string, a scalar, or a sequence holding a non-numeric element).
        ValueError: ``method``/``resampling`` is not a supported value,
            ``dst_crs`` cannot be parsed as a CRS, a source carries no CRS,
            ``no_data_value``, ``init`` or ``n`` names no number, or ``bbox`` is
            malformed
            (wrong length, non-finite, inverted, zero-area), crosses the
            antimeridian or the mosaic's longitude seam once reprojected, selects
            no whole pixel, does not overlap the mosaic, or cannot be projected
            into its CRS.
        RuntimeError: GDAL failed to open a source, reproject it, or build the
            source mosaic. When a source is at fault the message names it and
            its position in `src`, and chains GDAL's own error; any credential
            in a signed URL is redacted.
        DriverNotExistError: `dst` has no extension, or one the driver catalog
            does not know.
        FileFormatNotSupportedError: `dst`'s extension maps to a
            write-by-copy-only format, for any `method`.

    Warns:
        UserWarning: The sources declare more than one distinct no-data value
            (the first wins, and all of them are named); or the value they
            declare is one the output's data type cannot store, so a storable
            one is chosen in its place; or a numeric `init` was passed over as
            the marker because the mosaic's own cells hold that value, so the
            gaps hold the chosen sentinel instead; or the mosaic's cells use
            every value its dtype could spare as a marker, so it is written
            without one. Those four arise only while the marker is being
            inherited. Passing `no_data_value` explicitly silences them and
            raises one of its own where the value does not fit the output band:
            GDAL then drops it and the mosaic carries no marker at all. An
            explicit `no_data_value` under ``method="count"`` warns too, and is
            likewise dropped — see `method` for why a count mosaic declares
            nothing.

    Examples:
        - Mosaic two tiles, keeping the larger value wherever they overlap:
            ```python
            >>> from pyramids.dataset.merge import merge_rasters
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["tile_a.tif", "tile_b.tif"],
            ...     "mosaic_max.tif",
            ...     no_data_value=-9999.0,
            ...     method="max",
            ... )

            ```
        - Composite the average of every source covering a pixel, and write the
          matching per-pixel tally of how many of them contributed it:
            ```python
            >>> from pyramids.dataset.merge import merge_rasters
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["scene_a.tif", "scene_b.tif", "scene_c.tif"],
            ...     "mosaic_mean.tif",
            ...     method="mean",
            ... )
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["scene_a.tif", "scene_b.tif", "scene_c.tif"],
            ...     "mosaic_count.tif",
            ...     method="count",
            ... )

            ```
        - Default last-wins compositing. The overlap rule is unchanged, but the
          mosaic's no-data marker is now inherited from the sources rather than
          set to 0 (#1086):
            ```python
            >>> merge_rasters(["tile_a.tif", "tile_b.tif"], "mosaic.tif")  # doctest: +SKIP

            ```
        - Override the inherited marker on a z-order merge of integer tiles.
          The pixels no tile covers are filled with ``-1`` as well, so the
          marker the band declares is the one they hold:
            ```python
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["tile_a.tif", "tile_b.tif"],
            ...     "mosaic.tif",
            ...     no_data_value=-1,
            ... )

            ```
        - Keep a different value in the uncovered pixels by naming `init`. The
          band still declares ``-1``, but the gaps hold ``0`` and stay readable:
            ```python
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["tile_a.tif", "tile_b.tif"],
            ...     "mosaic.tif",
            ...     no_data_value=-1,
            ...     init=0,
            ... )

            ```
        - Ask for no marker at all, so every value in the mosaic stays readable
          — including a real 0 m in an elevation model:
            ```python
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["dem_a.tif", "dem_b.tif"],
            ...     "dem.tif",
            ...     no_data_value=None,
            ... )

            ```
        - Mosaic tiles from two UTM zones into a single CRS:
            ```python
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["utm32_tile.tif", "utm33_tile.tif"],
            ...     "mosaic_utm32.tif",
            ...     dst_crs=32632,
            ... )

            ```
        - Mosaic Requester-Pays S3 tiles by passing a signer (no ``with`` block):
            ```python
            >>> from pyramids.stac import AWSRequesterPaysSigner  # doctest: +SKIP
            >>> merge_rasters(  # doctest: +SKIP
            ...     ["s3://bucket/a.tif", "s3://bucket/b.tif"],
            ...     "mosaic.tif",
            ...     signer=AWSRequesterPaysSigner(region="us-west-2"),
            ... )

            ```
    """
    if method in _DEFERRED_METHODS:
        # A name this function knows and deliberately does not answer, separated from
        # a name it does not know at all: "median is not implemented, and here is what
        # it would cost" is a different message from "median is not a method".
        raise NotImplementedError(
            f"method={method!r} is not implemented by merge_rasters: "
            f"{_DEFERRED_METHODS[method]}"
        )
    if method not in _MERGE_METHODS:
        raise ValueError(
            f"method must be one of {list(_MERGE_METHODS)}, got {method!r}."
        )
    # Resolve `dst` here, before anything is opened. Both write paths resolve it
    # again where they need the driver name, but a destination the catalog
    # cannot answer for is a pure argument error -- and it used to be reported
    # only after every source had been opened and possibly reprojected, which
    # for /vsicurl/ inputs is network work spent to reach a typo. It also
    # reported at two different points depending on `method`.
    # Resolved here so a destination the catalog cannot answer for fails
    # before any source is opened, and so both write paths below agree. They
    # call it again where they need the name; it is a cached catalog lookup, so
    # the repeat is free and keeps each path readable on its own.
    resolve_output_driver(dst)

    # `init` and `n` both default to the string `"nan"`, which round-trips
    # through GDAL as float NaN -- a value no integer band can store. `n` never
    # reaches one as itself: the default means "no override" and is not passed on
    # at all (see `_source_nodata`). `init` gives way to a storable marker
    # whenever one is settled on, which is while the marker is being inherited
    # (see `_storable_marker`) and, now, when the caller named it and left `init`
    # alone (see `_explicit_fill`). It still reaches the band as NaN when the
    # caller names both. The spellings are kept for backwards-compat with the
    # previous gdal_merge.main-based signature.
    src_paths = [str(p) for p in src]
    if signer is not None:
        # Apply the signer's href rewrite to every source (e.g. graft a SAS
        # token onto a blob URL) so URL-signing signers authenticate. A no-op
        # for signers that authenticate via gdal_env() only — the base
        # sign_href returns the href unchanged. This mirrors load_asset, which
        # applies BOTH signer hooks (sign_href + gdal_env); applying only the
        # env half here would silently read URL-signed sources unauthenticated.
        src_paths = [signer.sign_href(p) for p in src_paths]

    # Every read of the *sources* runs under the signer's cloud config (a no-op
    # when signer is None) so authenticated remote sources open with the right
    # credentials for the whole merge. A remote `dst` is not scoped by either
    # context -- writing a mosaic back to a bucket still needs an explicit
    # `path=` or ambient credentials that already cover it. `ambient_store_config` does the same for
    # credentials the caller supplied through a pathless `with CloudConfig(...)`
    # around this call: it re-scopes them onto these sources' buckets, which is
    # what makes them visible to the VSICurl workers that fetch the byte-ranges
    # (#983). Both wrap the whole merge, not just the open, because those GETs
    # happen during the mosaic build long after each source is opened.
    # The signer goes first so it claims the credential family: an ambient block
    # then fills only what identity does not cover (an endpoint, a region) and
    # cannot contribute half an identity -- an ambient access key signed
    # alongside a signer's session token.
    with _cloud_config(signer, path=src_paths), ambient_store_config(src_paths):
        # Put every source on one CRS before compositing. The BuildVRT/Warp
        # mosaic below cannot reproject — it stitches pixel grids assuming a
        # shared CRS — so mismatched sources must be warped first or they would
        # mis-align silently. `_keepalive` holds the in-memory warped VRTs so
        # GDAL does not free them while the mosaic is built.
        requested = _requested_no_data(no_data_value)
        _validated_init(init)
        inheriting = requested is INHERIT_NO_DATA
        sources, _keepalive = _prepare_sources(src_paths, dst_crs, resampling)
        # Under its own name rather than rebound onto the parameter, so that
        # `no_data_value` keeps meaning "what the caller asked for" throughout --
        # including the sentinel, which `inheriting` above is the only reading of.
        resolved_no_data = requested
        if inheriting:
            # Read it off the handles _prepare_sources already opened rather
            # than reopening: each open is billable under Requester-Pays.
            resolved_no_data = inherit_no_data(
                [handle.GetRasterBand(1).GetNoDataValue() for handle in sources]
            )

        if method in _REDUCE_METHODS:
            # Pair each handle with its path so a failure on the reduce path
            # names the source, not a SWIG proxy address (#1107).
            labelled = [
                _Source(f"{index + 1}/{len(src_paths)} {path!r}", handle)
                for index, (path, handle) in enumerate(zip(src_paths, sources))
            ]
            if inheriting and resolved_no_data is None:
                # Nothing to inherit. The reduction writes Float64, so NaN is
                # storable and no real cell can hold it -- the same rule the
                # z-order path applies, answered by the dtype it writes.
                resolved_no_data = float("nan")
            if (
                method == "count"
                and not inheriting
                and resolved_no_data is not None
                and not np.isnan(float(resolved_no_data))
            ):
                # Dropped by `_merge_reduce` because a count mosaic declares no
                # marker; say so rather than ignoring an explicit argument.
                warnings.warn(
                    f"no_data_value={resolved_no_data!r} is dropped on "
                    "method='count': every cell of a count mosaic holds a real "
                    "tally, so a declared marker would mask the cells holding that "
                    "number -- starting with the uncovered cells' honest 0. The "
                    "mosaic is written with no no-data value.",
                    stacklevel=2,
                )
            _merge_reduce(
                labelled, str(dst), method, resolved_no_data, n, bbox, bbox_crs
            )
            return

        # z-order: "last" keeps natural order (last source wins); "first"
        # reverses so the original first source is placed last in the VRT and
        # therefore wins.
        ordered = list(reversed(sources)) if method == "first" else sources
        resolved_no_data, vrt_fill = _marker_and_fill(
            ordered, src_paths, init, resolved_no_data, inheriting, bbox, bbox_crs
        )
        vrt_opts = gdal.BuildVRTOptions(
            srcNodata=_source_nodata(n),
            VRTNodata=vrt_fill,
        )
        vrt_ds = run_gdal_op(
            partial(gdal.BuildVRT, "", ordered, options=vrt_opts),
            error=RuntimeError,
            action="building the source mosaic",
            subject=f"sources {src_paths!r}",
            hint=_READABLE_RASTERS_HINT,
        )
        proj_win = None
        if bbox is not None:
            # Resolve the window here, in the mosaic's own CRS, and hand Translate a
            # `projWin` that already lies on the mosaic's pixel boundaries -- rather
            # than passing the caller's bbox with `projWinSRS` and letting GDAL
            # reproject it. Two independent reprojections of the same window round to
            # opposite sides of a pixel edge, and the two merge paths then return
            # different rasters for identical arguments (measured 3x3 at x=111364.8
            # against 4x3 at x=0.0 for `dst_crs=3857`). One resolver for both paths
            # makes them agree by construction, and `_restrict_grid` also rejects a
            # disjoint window, which GDAL does not: it writes a 1x1 no-data raster at
            # the window's origin, which reads back as a successful merge of nothing.
            clipped, window_x, window_y = _restrict_grid(
                vrt_ds.GetGeoTransform(),
                vrt_ds.RasterXSize,
                vrt_ds.RasterYSize,
                vrt_ds.GetProjection(),
                bbox,
                bbox_crs,
            )
            # (ulx, uly, lrx, lry) — already snapped, so GDAL's own rounding is a
            # no-op and the output grid matches the reduction path's exactly.
            proj_win = [
                clipped[0],
                clipped[3],
                clipped[0] + window_x * clipped[1],
                clipped[3] + window_y * clipped[5],
            ]

        # `projWin` is what stops the read at the window: without it GDAL has no
        # reason to restrict what it pulls through /vsicurl and materialises the
        # whole mosaic extent.
        # Hand gdal.Translate the driver the catalog resolved rather than
        # letting it re-infer from the extension. The two tables disagree: the
        # catalog knows `.nc4` as a netCDF alias, GDAL's netCDF driver
        # advertises only `nc`, so the same `dst` wrote a netCDF through the
        # reduction path and died with "Could not identify an output driver"
        # here -- the format depending on `method` again, which is exactly what
        # resolving once was meant to stop.
        # Strict gate, not `for_copy`: this method has two write paths and the
        # other builds with `Create`. Relaxing only this one would make `.png`
        # legal for method="last" and illegal for method="min" -- the very
        # asymmetry the shared resolution above exists to remove.
        out_driver = resolve_output_driver(dst)
        # LZW is a GTiff creation option; other drivers reject it.
        # "none" *removes* the marker rather than omitting the option: the VRT
        # above carries `init` (NaN by default) as VRTNodata, so merely leaving
        # noData unset would stamp NaN -- invalid on an integer band, which is
        # exactly the "sentinel the band cannot store" defect. It is reached
        # only when the caller asked for no marker by passing `None`, or when
        # `_storable_marker` found the data using every value it could have
        # chosen (and said so); inheriting nothing otherwise yields a marker
        # rather than removing one.
        translate_opts = gdal.TranslateOptions(
            format=out_driver,
            creationOptions=["COMPRESS=LZW"] if out_driver == "GTiff" else [],
            projWin=proj_win,
            noData=("none" if resolved_no_data is None else str(resolved_no_data)),
        )
        out_ds = run_gdal_op(
            partial(gdal.Translate, str(dst), vrt_ds, options=translate_opts),
            error=RuntimeError,
            action="writing the mosaic",
            subject=f"destination {str(dst)!r}",
            outcome="produced no output",
        )
        out_ds.FlushCache()
        out_ds = None
        vrt_ds = None

pyramids.dataset.merge.stack_bands(files, *, band_names=None, align=False, resampling=None, no_data_value=INHERIT_NO_DATA, path=None, signer=None)

Stack N single-band rasters into one multi-band :class:Dataset.

Free-function alias for :meth:pyramids.dataset.Dataset.from_band_files — see that method for the full contract, edge cases, and examples.

Parameters:

Name	Type	Description	Default
files	list[str | Path]	Single-band raster paths/URLs to stack (order = band order).	required
band_names	list[str] | None	Explicit per-band names; None derives them from the file names.	None
align	bool	When True, resample mismatched inputs onto files[0]'s grid instead of raising :class:~pyramids.base._errors.AlignmentError.	False
resampling	str | dict[str, str] | None	Which algorithm align uses. None (default) keeps nearest neighbour; a method name applies to every band; a {band name: method} mapping sets it per band, so a categorical band can stay nearest while a continuous one interpolates. Requires align=True — passing it otherwise raises, since it could not be honoured.	None
no_data_value	Any	No-data value stamped on the output bands. Omitted means inherit from the source rasters: the first file that declares one wins, a disagreement warns, and if none declares one the output declares none either. Pass a value to override that, or None for "declare no sentinel at all". The inheritance is the same rule :func:merge_rasters follows, resolved by the same :func:~pyramids.base._domain.inherit_no_data -- but what happens when nothing is inherited is not, and deliberately so. A stack covers one grid, so it has no uncovered pixel and declares nothing; a mosaic generally has them, and settles on a marker rather than leaving them to read as data.	INHERIT_NO_DATA
path	str | Path | None	Output path, whose extension selects the driver (.tif -> GTiff, .nc -> netCDF, …); None keeps the result in memory. COMPRESS=LZW is applied only when the extension resolves to GTiff. A write-by-copy-only format such as PNG is refused — see :meth:pyramids.dataset.Dataset.from_band_files for why both of its write paths answer alike.	None
signer	Any	Optional signer exposing sign_href(str) -> str and gdal_env() -> dict[str, str] (e.g. a :class:pyramids.stac.signers.Signer). When given, both hooks are applied (as in :func:pyramids.stac.load_asset): every input href is rewritten through signer.sign_href first, then signer.gdal_env() is installed via :class:~pyramids.base.remote.CloudConfig for the duration of the stack, so authenticated cloud inputs (URL-signed or env-credentialed) read with the right credentials. None (default) leaves behaviour unchanged.	None

Returns:

Name	Type	Description
Dataset	Dataset	A multi-band dataset, one band per input file.

Raises:

Type	Description
DriverNotExistError	path has no extension, or one the driver catalog does not know.
FileFormatNotSupportedError	path's extension maps to a write-by-copy-only format, whichever write path the inputs take.

Warns:

Type	Description
UserWarning	The inputs declare more than one distinct no-data value, while no_data_value is being inherited. The first file's wins.

Examples:

• Stack a Landsat scene's per-band GeoTIFFs; the band names come from the file names:

  >>> from pyramids.dataset.merge import stack_bands
  >>> scene = stack_bands(  # doctest: +SKIP
  ...     ["scene.B2.tif", "scene.B3.tif", "scene.B4.tif"]
  ... )
  >>> scene.band_names  # doctest: +SKIP
  ['B2', 'B3', 'B4']
  

• Name the bands and write the stack straight to disk:

  >>> from pyramids.dataset.merge import stack_bands
  >>> stack_bands(  # doctest: +SKIP
  ...     ["b2.tif", "b3.tif", "b4.tif"],
  ...     band_names=["blue", "green", "red"],
  ...     path="stack.tif",
  ... )
  

• Resample mismatched inputs onto the first file's grid, and declare no sentinel on the result:

  >>> from pyramids.dataset.merge import stack_bands
  >>> stack_bands(  # doctest: +SKIP
  ...     ["b2.tif", "coarse_b11.tif"],
  ...     align=True,
  ...     no_data_value=None,
  ... )
  

See Also

• :meth:pyramids.dataset.Dataset.from_band_files: The method this delegates to, with the full contract and runnable examples.
• :func:merge_rasters: Mosaics rasters that tile the same area, rather than stacking rasters that share a grid.

Source code in src/pyramids/dataset/merge.py

def stack_bands(
    files: list[str | Path],
    *,
    band_names: list[str] | None = None,
    align: bool = False,
    resampling: str | dict[str, str] | None = None,
    no_data_value: Any = INHERIT_NO_DATA,
    path: str | Path | None = None,
    signer: Any = None,
) -> Dataset:
    """Stack N single-band rasters into one multi-band :class:`Dataset`.

    Free-function alias for :meth:`pyramids.dataset.Dataset.from_band_files`
    — see that method for the full contract, edge cases, and examples.

    Args:
        files: Single-band raster paths/URLs to stack (order = band order).
        band_names: Explicit per-band names; ``None`` derives them from the
            file names.
        align: When ``True``, resample mismatched inputs onto ``files[0]``'s
            grid instead of raising :class:`~pyramids.base._errors.AlignmentError`.
        resampling: Which algorithm ``align`` uses. ``None`` (default) keeps
            nearest neighbour; a method name applies to every band; a
            ``{band name: method}`` mapping sets it per band, so a categorical
            band can stay nearest while a continuous one interpolates. Requires
            ``align=True`` — passing it otherwise raises, since it could not be
            honoured.
        no_data_value: No-data value stamped on the output bands. Omitted means
            **inherit from the source rasters**: the first file that declares
            one wins, a disagreement warns, and if none declares one the output
            declares none either. Pass a value to override that, or ``None`` for
            "declare no sentinel at all". The *inheritance* is the same rule
            :func:`merge_rasters` follows, resolved by the same
            :func:`~pyramids.base._domain.inherit_no_data` -- but what happens
            when nothing is inherited is not, and deliberately so. A stack
            covers one grid, so it has no uncovered pixel and declares nothing;
            a mosaic generally has them, and settles on a marker rather than
            leaving them to read as data.
        path: Output path, whose extension selects the driver (``.tif`` ->
            GTiff, ``.nc`` -> netCDF, …); ``None`` keeps the result in memory.
            `COMPRESS=LZW` is applied only when the extension resolves to
            GTiff. A write-by-copy-only format such as PNG is refused — see
            :meth:`pyramids.dataset.Dataset.from_band_files` for why both of
            its write paths answer alike.
        signer: Optional signer exposing ``sign_href(str) -> str`` and
            ``gdal_env() -> dict[str, str]`` (e.g. a
            :class:`pyramids.stac.signers.Signer`). When given, **both** hooks
            are applied (as in :func:`pyramids.stac.load_asset`): every input
            href is rewritten through ``signer.sign_href`` first, then
            ``signer.gdal_env()`` is installed via
            :class:`~pyramids.base.remote.CloudConfig` for the duration of the
            stack, so authenticated cloud inputs (URL-signed or env-credentialed)
            read with the right credentials. ``None`` (default) leaves behaviour
            unchanged.

    Returns:
        Dataset: A multi-band dataset, one band per input file.

    Raises:
        DriverNotExistError: `path` has no extension, or one the driver
            catalog does not know.
        FileFormatNotSupportedError: `path`'s extension maps to a
            write-by-copy-only format, whichever write path the inputs take.

    Warns:
        UserWarning: The inputs declare more than one distinct no-data value,
            while `no_data_value` is being inherited. The first file's wins.

    Examples:
        - Stack a Landsat scene's per-band GeoTIFFs; the band names come from
          the file names:
            ```python
            >>> from pyramids.dataset.merge import stack_bands
            >>> scene = stack_bands(  # doctest: +SKIP
            ...     ["scene.B2.tif", "scene.B3.tif", "scene.B4.tif"]
            ... )
            >>> scene.band_names  # doctest: +SKIP
            ['B2', 'B3', 'B4']

            ```
        - Name the bands and write the stack straight to disk:
            ```python
            >>> from pyramids.dataset.merge import stack_bands
            >>> stack_bands(  # doctest: +SKIP
            ...     ["b2.tif", "b3.tif", "b4.tif"],
            ...     band_names=["blue", "green", "red"],
            ...     path="stack.tif",
            ... )

            ```
        - Resample mismatched inputs onto the first file's grid, and declare no
          sentinel on the result:
            ```python
            >>> from pyramids.dataset.merge import stack_bands
            >>> stack_bands(  # doctest: +SKIP
            ...     ["b2.tif", "coarse_b11.tif"],
            ...     align=True,
            ...     no_data_value=None,
            ... )

            ```

    See Also:
        - :meth:`pyramids.dataset.Dataset.from_band_files`: The method this
          delegates to, with the full contract and runnable examples.
        - :func:`merge_rasters`: Mosaics rasters that tile the same area,
          rather than stacking rasters that share a grid.
    """
    if signer is not None:
        files = [signer.sign_href(str(f)) for f in files]
    # Same pairing as `merge_rasters`: the signer's config, plus a re-scope of
    # any credentials the caller supplied through a pathless `with
    # CloudConfig(...)` around this call, so both reach the VSICurl workers
    # (#983). The signer goes first so it claims the credential family.
    band_paths = [str(f) for f in files]
    with _cloud_config(signer, path=band_paths), ambient_store_config(band_paths):
        result = Dataset.from_band_files(
            files,
            band_names=band_names,
            align=align,
            resampling=resampling,
            no_data_value=no_data_value,
            path=path,
        )
    return result