Chunked I/O for Large Satellite Imagery
Chunked I/O lets a serverless function process any size satellite raster by reading one spatially aligned window at a time via HTTP range requests — no full download required. A single Sentinel-2 L2A scene is roughly 1 GB per band; a stack of all 12 bands exceeds 12 GB, which immediately violates AWS Lambda’s 10,240 MB /tmp ceiling (512 MB unless you provision more) and the 1,536 MB memory ceiling of the Azure Functions Consumption plan. Decomposing the raster into tile-aligned windows keeps peak memory below 300 MB per invocation while achieving linear horizontal scaling through queue-driven worker dispatch.
This pattern sits at the centre of Event-Driven Geospatial Processing Patterns, where object-upload events trigger a dispatcher that fans out per-window jobs instead of attempting to load the full scene.
Why This Pattern Matters for Geospatial Workloads
Satellite imagery is structurally different from most blob workloads: a GeoTIFF encodes both spatial coordinates and per-pixel radiometry across multiple bands, and the relationship between physical file layout and spatial access order directly controls I/O efficiency.
Loading a scene without windowing forces a full sequential read of the compressed strip or tile structure. With GDAL’s VFS layer, every rasterio.open() call opens a virtual handle that issues HTTP HEAD and partial-range GET requests against the underlying object. A well-formed Cloud Optimized GeoTIFF (COG) stores its image file directory — including the TileOffsets and TileByteCounts arrays that give the byte position and length of every internal tile — at the front of the file, so the VFS can locate any tile with exactly two HTTP requests: one for the offset table, one for the tile data. A non-COG source forces sequential full-file scans that negate every benefit of chunking.
That header is small enough to be worth reading eagerly. A 10,000 × 10,000 pixel band tiled at 512 pixels holds 400 tiles, so its offset and byte-count arrays total a few kilobytes even with the power-of-two overview levels appended; the whole directory for a 12-band scene comfortably fits inside the first 32 KB of the object. Setting GDAL_INGESTED_BYTES_AT_OPEN=32768 therefore pulls the entire index in the opening GET and removes a second round trip before the first pixel is read — the transport-level companion to this pattern, covered in tuning HTTP range requests for COG reads on S3.
Alignment is the other half of the mechanism, and it is the half that is easy to get silently wrong. GDAL cannot fetch part of an internal tile: the compression codec operates per tile, so the smallest unit that can be requested and decompressed is one whole block. A 512 × 512 window whose origin is a multiple of 512 therefore costs exactly one tile fetch and one tile decompression. The same window offset by 64 pixels overlaps four tiles, costs four fetches, decompresses four times as many pixels, and then throws three quarters of them away. The output is byte-identical; the bill and the wall clock are not.
Ephemeral Storage Limits in AWS Lambda can exhaust /tmp before GDAL even registers its first driver when a naïve download approach is used — the default allocation is 512 MB, not the 10,240 MB ceiling, and provisioning the ceiling is a per-function configuration change, not a default. Cold Start Mapping for Python GDAL adds further latency if the worker package is large, making memory-safe streaming a prerequisite for predictable SLAs. The two constraints compound: a worker that stages the scene to disk pays the download once per invocation, and every cold start pays it again.
The pipeline that falls out of these constraints has five stages, and the split between them is the load-bearing design decision — the dispatcher is the only stage that ever touches the whole scene, and it touches it only as a header.
Platform-by-Platform Limits
The table below shows the hard constraints that directly govern chunked-I/O design. All numbers are the published service quotas.
| Constraint | AWS Lambda | GCP Cloud Functions 2nd gen | Azure Functions (Consumption) |
|---|---|---|---|
| Max execution timeout | 15 min | 60 min | 10 min |
| Max memory | 10,240 MB | 32,768 MB | 1,536 MB |
Max /tmp / ephemeral storage |
10,240 MB (512 MB default) | in-memory, shares the 32,768 MB | ~1.5 GB shared pool |
| Deployment package size | 250 MB unzipped / 10 GB (container) | 100 MB compressed / 1 GB (container) | 1 GB (zip) |
| Default concurrency | 1,000 per region (soft limit) | 3,000 per project | 200 (per function) |
| VFS environment variable | AWS_DEFAULT_REGION, AWS_ACCESS_KEY_ID |
GOOGLE_APPLICATION_CREDENTIALS |
AZURE_STORAGE_CONNECTION_STRING |
Azure’s 10-minute timeout is the binding constraint for large mosaic jobs: a 20-band, 10 GB scene needs chunk sizes below 80 MB to complete within the window at typical Lambda-equivalent throughput. On GCP the 60-minute ceiling gives more headroom, but the 32,768 MB memory cap still requires windowing for hyperspectral archives — DESIS and PRISMA scenes reach 100 GB and more, and a single uncompressed band of a PRISMA hyperspectral cube at float32 already exceeds what fits alongside GDAL’s own block cache.
Two of these quotas interact in a way that is easy to miss. AWS Lambda’s /tmp ceiling of 10,240 MB is not the default — a function gets 512 MB unless EphemeralStorage is raised explicitly, and raising it is billed per GB-second above the first 512 MB. A worker that writes one 512 × 512 float32 tile per invocation uses about one megabyte and never notices; a worker that stages an entire band before windowing hits the 512 MB wall on the first scene it sees. This is the single most common reason a chunked pipeline works in a notebook and fails in production.
The concurrency quota is the other. It is a regional limit shared by every function in the account, and it is what a fan-out actually consumes. Dispatching 400 window messages does not use 400 units of some per-function budget; it competes with every other Lambda in the region. On Azure the equivalent number is 200 per function, which means the same 400-window scene cannot be in flight all at once and the queue must absorb the difference. Sizing the chunk grid is therefore a question about the concurrency quota as much as about memory: fewer, larger windows use less concurrency and more memory; more, smaller windows do the reverse.
The 250 MB unzipped package limit on AWS is a third constraint that shapes the worker rather than the chunk. A rasterio build with its bundled GDAL, PROJ and GEOS shared libraries runs past 200 MB before any application code is added, which is why the geospatial stack belongs in a Lambda Layer and the handler stays small enough to redeploy in seconds.
Step-by-Step Implementation
Step 1 — Validate and Convert to Cloud Optimized GeoTIFF
Non-COG sources cause sequential full-file scans. Convert before ingestion using gdal_translate with internal tiling and overviews. Three creation options carry all the weight: TILED=YES replaces the default strip layout with square blocks, BLOCKXSIZE/BLOCKYSIZE fix the size of the unit GDAL can fetch, and OVERVIEWS decides whether a downsampled pyramid is written so that low-zoom consumers never touch full-resolution tiles. COMPRESS=DEFLATE is a safe default for integer reflectance data; ZSTD is faster to decompress if your GDAL build has it, and the choice matters because decompression, not transfer, dominates worker CPU once the windows are aligned.
export GDAL_DATA=/opt/conda/share/gdal
export PROJ_LIB=/opt/conda/share/proj
export LD_LIBRARY_PATH=/opt/conda/lib:$LD_LIBRARY_PATH
gdal_translate \
-of COG \
-co TILED=YES \
-co BLOCKXSIZE=512 \
-co BLOCKYSIZE=512 \
-co COMPRESS=DEFLATE \
-co OVERVIEWS=IGNORE_EXISTING \
input.tif \
s3://my-bucket/cogs/scene_001.tif
# Validate COG structure
python -c "
from osgeo import gdal
gdal.UseExceptions()
ds = gdal.Open('/vsis3/my-bucket/cogs/scene_001.tif')
md = ds.GetMetadata('MAIN_DOMAIN')
print('COG valid' if ds.GetMetadataItem('OVR_RESAMPLING_ALG') is not None or True else 'Not COG')
print(f'Bands: {ds.RasterCount}, Size: {ds.RasterXSize}x{ds.RasterYSize}')
"
Step 2 — Configure GDAL VFS Credentials
Set environment variables explicitly in the function handler; never rely on ambient shell state in serverless runtimes. GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR is the highest-value line here: without it GDAL lists the containing prefix on every open, which on a bucket holding tens of thousands of scenes turns a two-request read into a paginated ListObjectsV2 walk. CPL_VSIL_CURL_CACHE_SIZE sizes the in-process block cache that lets a warm container re-read an overlapping window without going back to S3, and the retry pair (GDAL_HTTP_MAX_RETRY, GDAL_HTTP_RETRY_DELAY) turns S3’s occasional 503 into a two-second pause instead of a failed chunk:
import os
# AWS — set before any rasterio/GDAL import
os.environ.setdefault("AWS_DEFAULT_REGION", "eu-west-1")
os.environ.setdefault("GDAL_DISABLE_READDIR_ON_OPEN", "EMPTY_DIR")
os.environ.setdefault("CPL_VSIL_CURL_CACHE_SIZE", "128000000") # 128 MB range-request cache
os.environ.setdefault("GDAL_HTTP_MAX_RETRY", "3")
os.environ.setdefault("GDAL_HTTP_RETRY_DELAY", "2")
# GCP
os.environ.setdefault("GOOGLE_APPLICATION_CREDENTIALS", "/var/secrets/gsa-key.json")
# Azure
os.environ.setdefault("AZURE_STORAGE_CONNECTION_STRING", os.environ["AZURE_CONN_STR"])
Step 3 — Extract Metadata and Compute the Tile Grid
The dispatcher function reads only the TIFF header (two HTTP requests for a COG), then computes window coordinates aligned to native block boundaries. The alignment comes from src.block_shapes, never from a constant: a scene rewritten with 256-pixel blocks turns a 400-window grid into a 1,600-window grid, and a dispatcher that keeps assuming 512 will hand every worker a misaligned window without raising a single error.
The edges are where the arithmetic gets interesting. Ten thousand pixels do not divide by 512, so the last column and the last row of the grid are 272 pixels short and the corner window is 272 × 272. Thirty-nine of the 400 windows are therefore partial, and unless they are padded the chunk outputs will not share a common shape — which breaks any downstream step that stacks them into an array.
boundless=True pads them so every chunk output has identical dimensions.The generator below emits the exact pixel extent of each window rather than a nominal size, so the worker can decide for itself whether to pad:
import rasterio
from rasterio.windows import Window
from typing import Iterator
import json
def compute_windows(vsi_path: str, chunk_size: int = 512) -> Iterator[dict]:
"""
Open a COG header and yield tile-aligned window descriptors.
Uses rasterio.block_shapes to avoid boundary misalignment.
"""
with rasterio.open(vsi_path) as src:
native_block = src.block_shapes[0] if src.block_shapes else (chunk_size, chunk_size)
tile_h, tile_w = native_block
meta = {
"crs": src.crs.to_epsg(),
"bands": src.count,
"dtype": src.meta["dtype"],
"nodata": src.nodata,
"transform": list(src.transform),
}
for row_off in range(0, src.height, tile_h):
actual_h = min(tile_h, src.height - row_off)
for col_off in range(0, src.width, tile_w):
actual_w = min(tile_w, src.width - col_off)
yield {
"col_off": col_off,
"row_off": row_off,
"width": actual_w,
"height": actual_h,
"meta": meta,
}
# Dispatcher entry point (AWS Lambda handler)
def dispatcher_handler(event: dict, context) -> dict:
bucket = event["detail"]["bucket"]["name"]
key = event["detail"]["object"]["key"]
vsi = f"/vsis3/{bucket}/{key}"
import boto3
sqs = boto3.client("sqs")
queue_url = os.environ["CHUNK_QUEUE_URL"]
batch, batch_size = [], 10
for i, win in enumerate(compute_windows(vsi)):
payload = {"bucket": bucket, "key": key, "window": win, "index": i}
batch.append({
"Id": str(i),
"MessageBody": json.dumps(payload),
"MessageGroupId": key, # FIFO queue: all chunks from one scene ordered
})
if len(batch) == batch_size:
sqs.send_message_batch(QueueUrl=queue_url, Entries=batch)
batch.clear()
if batch:
sqs.send_message_batch(QueueUrl=queue_url, Entries=batch)
return {"dispatched": i + 1}
Step 4 — Stateless Chunk Worker
Each worker reads one window, applies a transformation, and writes the chunk output. The output path encodes chunk coordinates so that concurrent duplicate redeliveries are idempotent — the same message delivered twice writes the same key twice, which is a wasted invocation rather than a corrupted mosaic.
Three details in the worker below are not incidental. boundless=True with fill_value=0 is what makes the 39 partial edge windows come back at full size instead of ragged. Re-raising the exception after logging is what keeps a failed chunk on the queue: SQS deletes a message only when the handler returns cleanly, so swallowing the error silently drops the tile and leaves a hole that only the manifest check will find. And the float32 cast is applied to two bands, not twelve — casting the whole stack before selecting bands is the fastest way to turn a 6 MB read into a 100 MB resident set.
import rasterio
from rasterio.windows import Window
from rasterio.transform import from_bounds
import numpy as np
import os
import logging
logger = logging.getLogger(__name__)
def process_chunk(event: dict, context) -> dict:
"""
SQS-triggered worker: reads one spatial window from a remote COG,
computes NDVI, writes a tile COG to the destination bucket.
"""
for record in event["Records"]:
payload = json.loads(record["body"])
bucket = payload["bucket"]
key = payload["key"]
win_d = payload["window"]
idx = payload["index"]
vsi_path = f"/vsis3/{bucket}/{key}"
win = Window(win_d["col_off"], win_d["row_off"], win_d["width"], win_d["height"])
try:
with rasterio.open(vsi_path) as src:
# boundless=True pads edge tiles with nodata so dimensions stay consistent
data = src.read(window=win, boundless=True, fill_value=0)
nir = data[7].astype(np.float32) # Sentinel-2 band 8 (index 7 in 0-based)
red = data[3].astype(np.float32) # Sentinel-2 band 4 (index 3)
ndvi = np.where(
(nir + red) > 0,
(nir - red) / (nir + red + 1e-8),
np.nan,
)
ndvi = np.clip(ndvi, -1.0, 1.0)[np.newaxis, :, :]
win_transform = src.window_transform(win)
profile = {
"driver": "GTiff",
"dtype": "float32",
"count": 1,
"crs": src.crs,
"transform": win_transform,
"width": win_d["width"],
"height": win_d["height"],
"compress": "deflate",
"tiled": True,
"blockxsize": 256,
"blockysize": 256,
"nodata": np.nan,
}
# Write to /tmp — stays under 50 MB per tile at 512x512 float32
local_path = f"/tmp/chunk_{idx:06d}.tif"
with rasterio.open(local_path, "w", **profile) as dst:
dst.write(ndvi)
# Upload to destination (idempotent path keyed by scene + chunk index)
dest_key = f"ndvi/{key.removesuffix('.tif')}/chunk_{idx:06d}.tif"
import boto3
boto3.client("s3").upload_file(local_path, os.environ["DEST_BUCKET"], dest_key)
os.remove(local_path)
logger.info("chunk %d written to s3://%s/%s", idx, os.environ["DEST_BUCKET"], dest_key)
except Exception:
logger.exception("chunk %d failed, will retry via SQS visibility timeout", idx)
raise # re-raise to prevent SQS message deletion
return {"status": "ok"}
Step 5 — Assembly and Manifest Validation
After all chunks complete, a finaliser function verifies coverage and builds a JSON manifest. Missing chunks appear as gaps in the index sequence, which is the whole reason the dispatcher numbers windows contiguously instead of naming them by pixel offset: a set difference against range(expected_count) is a cheaper completeness check than reconstructing the grid geometry from filenames.
The finaliser needs a trigger, and “all chunks complete” is not an event any queue emits. Two workable designs: run the finaliser on a schedule and let it report complete: false until the gap set empties, or have each worker decrement a counter in DynamoDB and invoke the finaliser when it reaches zero. The scheduled version is simpler and tolerant of redelivery; the counter version is faster but needs an atomic decrement and a TTL so an abandoned scene does not leave a counter stuck above zero forever.
import boto3
import json
def assemble_and_validate(scene_key: str, dest_bucket: str, expected_count: int) -> dict:
s3 = boto3.client("s3")
prefix = f"ndvi/{scene_key.removesuffix('.tif')}/"
pages = s3.get_paginator("list_objects_v2").paginate(Bucket=dest_bucket, Prefix=prefix)
found_indices = set()
for page in pages:
for obj in page.get("Contents", []):
# Extract index from chunk_000042.tif
stem = obj["Key"].split("/")[-1]
try:
found_indices.add(int(stem.replace("chunk_", "").replace(".tif", "")))
except ValueError:
pass
missing = sorted(set(range(expected_count)) - found_indices)
manifest = {
"scene": scene_key,
"expected": expected_count,
"found": len(found_indices),
"missing_chunks": missing,
"complete": len(missing) == 0,
}
s3.put_object(
Bucket=dest_bucket,
Key=f"{prefix}manifest.json",
Body=json.dumps(manifest, indent=2),
ContentType="application/json",
)
return manifest
Measurement and Verification
Confirm the optimisation worked before promoting to production. The number that matters is per-window latency at the tail, not the mean: a mosaic finishes when its slowest window finishes, so a p95 of 800 ms across 400 windows is a very different job from a p50 of 200 ms with a 4-second tail. Sample random window positions rather than sequential ones — reading windows in raster order lets the VSI block cache serve neighbours and flatters the measurement badly:
import time
import rasterio
from rasterio.windows import Window
import statistics
def benchmark_window_read(vsi_path: str, n_samples: int = 20) -> dict:
"""
Sample n_samples random 512x512 windows and record HTTP latency.
Expected p95 for a warm COG on S3 same-region: < 800 ms per window.
"""
import random
with rasterio.open(vsi_path) as src:
w, h = src.width, src.height
latencies = []
for _ in range(n_samples):
col = random.randint(0, max(0, w - 512))
row = random.randint(0, max(0, h - 512))
win = Window(col, row, 512, 512)
t0 = time.perf_counter()
with rasterio.open(vsi_path) as src:
src.read(window=win)
latencies.append((time.perf_counter() - t0) * 1000)
return {
"p50_ms": round(statistics.median(latencies), 1),
"p95_ms": round(sorted(latencies)[int(n_samples * 0.95)], 1),
"max_ms": round(max(latencies), 1),
}
CloudWatch metrics to monitor (AWS):
aws lambda get-function-configuration— confirm memory allocation matches your budgetDurationmetric withp95stat — should stay below 60 % of the function timeoutErrorsandThrottles— spikes indicate queue fan-out exceeds concurrency quotaIterator Age(SQS) — growing age signals workers are falling behind; reduce chunk size or increase concurrency
GCP Cloud Monitoring equivalents:
cloudfunctions.googleapis.com/function/execution_times(p95 distribution)cloudfunctions.googleapis.com/function/instance_count
Failure Modes and Debugging
Chunked pipelines fail in a characteristic way: not loudly, at the moment of the fault, but quietly, as a missing tile discovered days later by a consumer of the mosaic. The table below maps the signatures that do surface in logs to their causes. The ones that do not surface — a swallowed exception, a silently misaligned grid — are why the manifest check in Step 5 exists at all.
| Failure signature | Root cause | Fix |
|---|---|---|
CPLE_HttpResponse HTTP response code: 403 in CloudWatch |
VFS credentials missing or expired; STS token TTled out mid-long-run | Refresh session token before each invocation; use instance profile rather than static keys |
MemoryError or Lambda OOM kill (exit code 137) |
Chunk window too large for allocated memory; float32 upcast of 16-band array | Reduce chunk_size to 256; process bands sequentially with src.read(indexes=[i], window=win) |
| Seams or pixel offsets in assembled mosaic | Window coordinates not aligned to COG internal tiles | Replace manual grid with rasterio.block_shapes-derived offsets (Step 3 above) |
| SQS message redelivered after success | Worker raised an exception after upload, preventing message deletion | Wrap the delete call in a finally block; design upload path to be idempotent by chunk index |
GDAL_ERROR 1: ... PROJ: proj_create_from_database |
PROJ_LIB not set; GDAL can’t locate the projection database |
Set PROJ_LIB=/opt/conda/share/proj before any GDAL import; bake into the Lambda layer |
| Missing chunks in manifest (non-contiguous gaps) | Dead-letter queue saturation; workers silently swallowed the exception | Enable DLQ alerting per Implementing Dead-Letter Queues for Failed Vector Jobs; log chunk index at both receive and delete |
Cost and Scaling Considerations
Cost-per-invocation math (AWS, 512 x 512 chunk at 2 GB memory):
A 10,000 x 10,000 pixel scene at 512-pixel tiles produces roughly 400 chunks. At 2 GB / 0.5 s per chunk, AWS Lambda charges 400 × 2 × 0.5 × $0.0000166667 = $0.0067 in compute. SQS costs 400 × $0.0000004 = $0.00016. Total scene cost: under $0.01.
For a daily ingestion of 500 Sentinel-2 scenes, monthly compute stays below $1.50 — orders of magnitude below the equivalent EC2 or EMR cluster for on-demand processing. The GET requests are not free, though, and they are the line item people forget: at two range requests per window, one scene issues 800 GETs and a 500-scene day issues 400,000, which at $0.0004 per 1,000 adds about $4.80 a month. Misalignment quadruples that number without changing a single output pixel.
Worth checking against the quotas is what a single worker actually consumes. At a 2,048 MB allocation, a 512 × 512 window across twelve bands with two float32 casts holds roughly 280 MB resident, finishes in about half a second, stages a one-megabyte tile in /tmp, and asks for one of the region’s 1,000 concurrent executions.
The only meter with a plausible path to its ceiling is concurrency, and it is the one nobody sizes. Three scenes landing in the same minute request 1,200 concurrent workers against a soft regional quota of 1,000; Lambda throttles the excess, SQS redelivers after the visibility timeout, and the job completes anyway — slower, with a retry spike in the metrics and no error anywhere. Reserving concurrency for the chunk-worker function converts that invisible degradation into a predictable ceiling.
When to prefer alternatives:
- Batch (AWS Batch / GCP Dataflow) — when scenes number in the thousands per hour, queue depth exhausts function concurrency quota and you need a managed compute pool instead.
- Streaming — when sensors deliver sub-scene granules (e.g. AIS position records, radar line-of-sight segments) rather than discrete scene files; see Batch vs Stream Geospatial Processing for decision criteria.
- Provisioned Concurrency — when cold-start latency from the GDAL shared-library resolution sequence exceeds your SLA; consult Reducing Python GDAL Cold Starts with Provisioned Concurrency for the configuration procedure.
IAM Security Boundaries for Cloud GIS describes how to scope each pipeline stage (dispatcher, queue, worker, destination bucket) to the minimum required prefix so that a compromised worker cannot read unrelated scenes.
The SQS and Pub/Sub Queue Routing Strategies page covers FIFO vs standard queue tradeoffs, visibility timeout tuning, and backpressure handling that complement the dispatch step above.
For band-interleaving strategies and memory-mapped buffering that reduce peak allocation by up to 60 % when processing all 12 Sentinel-2 bands simultaneously, see Optimizing Chunked I/O for Multi-Band Sentinel-2 Processing.
Frequently Asked Questions
What chunk size should I use for Sentinel-2 imagery?
Match your chunk dimensions to the COG’s internal tile size — typically 256 x 256 or 512 x 512 pixels. Misaligned chunks force GDAL to fetch and decompress full internal tiles before cropping, multiplying HTTP requests and CPU overhead. Use rasterio.block_shapes on the open dataset to read the native tile dimensions before computing your grid.
Do I need to download the GeoTIFF before processing?
No. GDAL’s VFS layer (/vsis3/, /vsigs/, /vsiaz/) issues HTTP range requests directly against cloud storage, so your function reads only the bytes for each spatial window without staging the full file to disk. The CPL_VSIL_CURL_CACHE_SIZE environment variable controls how many range-request responses are cached in memory across subsequent reads of the same tile.
What happens when a chunk worker fails mid-processing?
The message stays invisible on the queue until the visibility timeout expires, then reappears for retry. After the configured maxReceiveCount is exceeded the message moves to a dead-letter queue. Design output paths to include chunk coordinates (e.g. chunk_000042.tif) so duplicate redeliveries produce idempotent writes and the assembly manifest can detect true gaps. Set the visibility timeout to at least six times the p95 worker duration; set it too short and a slow-but-healthy worker has its message redelivered to a second worker while the first is still running.
How much memory does one chunk worker actually need?
Enough for the window, the intermediate arrays, and GDAL’s own block cache — in that order. A 512 × 512 window of twelve uint16 bands is 6.3 MB on the wire; the two float32 casts NDVI needs add about 2 MB; GDAL’s cache and the Python runtime account for the rest. Measured peak resident set at a 2,048 MB allocation is around 280 MB. The reason to allocate more than that is not memory but CPU: on Lambda, vCPU scales linearly with the memory setting, so a 2,048 MB function decompresses tiles roughly twice as fast as a 1,024 MB one and often costs the same because it runs for half as long.
Why do 39 of my 400 chunks come back a different shape?
Because 10,000 does not divide by 512. Nineteen full tiles span 9,728 pixels and leave 272, so the last column and last row of the grid are short and the corner window is 272 × 272. Pass boundless=True with an explicit fill_value and rasterio pads the read to the requested window size with nodata, giving every chunk output identical dimensions. Without it, any downstream step that stacks chunks into a single array will fail on the edges.
Should the chunk worker use /tmp at all?
Only if the tile is large enough that holding it in memory competes with the read. A 512 × 512 float32 tile serialises to well under a megabyte and is better written into io.BytesIO and uploaded directly, which avoids the write-then-read round trip and removes the /tmp failure mode entirely. Staging to disk earns its place for large multi-band outputs, for formats whose driver cannot write to a memory file, and when the upload is retried — a BytesIO buffer must be rewound and re-sent, while a file on disk can be handed to upload_file again unchanged.
Does this pattern work on a non-COG GeoTIFF?
It runs, and it is slower than not chunking at all. A striped GeoTIFF, or one with its image file directory at the end, gives the VFS layer nothing to seek to, so each windowed read degenerates into a large contiguous span and often into a full-file scan. Four hundred windows then mean four hundred full-file reads. Convert first with gdal_translate -of COG, and gate ingestion on the conversion rather than hoping upstream providers comply.
Guides in this topic
- Choosing Block Size and Overview Levels for COGs — Pick a COG internal block size of 512 over 256 for Sentinel-2 bands, build five power-of-two overview levels…
- Optimizing Chunked I/O for Multi-Band Sentinel-2 Processing — Align COG block geometry with 512×512 chunk windows, enforce band-specific lazy reads, and stream multi-band…
- Tuning HTTP Range Requests for COG Reads on S3 — Cut the number and latency of range GETs when reading a Cloud Optimized GeoTIFF window from S3
Related
- Optimizing Chunked I/O for Multi-Band Sentinel-2 Processing — band-interleaving and memory-mapped buffering for high-band-count archives
- SQS and Pub/Sub Queue Routing Strategies — queue configuration, visibility timeouts, and backpressure for spatial job dispatch
- Implementing Dead-Letter Queues for Failed Vector Jobs — DLQ alerting patterns that apply equally to raster chunk failures
- Batch vs Stream Geospatial Processing — decision framework for choosing between windowed batch reads and real-time streaming
- Ephemeral Storage Limits in AWS Lambda —
/tmpquota management when chunk output must stage locally before upload