Reading a COG from a URL Is Slow or Downloads Everything
On this page
Problem statement
A remote read that should transfer a megabyte transfers the whole file, or takes seconds when it should take milliseconds. There are five causes, and they are distinguishable.
Measured against a byte-counting server on a 4096 × 4096 Sentinel-2 band (Example 4):
whole scene at 1/16 resolution
file with overviews 3 requests 0.131 MB
file without overviews 5 requests 24.135 MB <- the whole file
four scattered windows from the COG
default GDAL settings 18 requests 4.996 MB
GDAL_DISABLE_READDIR_ON_OPEN 9 requests 4.996 MB
The first is a property of the file. The second is a property of your configuration. Diagnosing which you have is the whole job.
Quick answer
import rasterio
url = served("striped.tif")
with rasterio.open(url) as ds:
print(f"blocks {ds.block_shapes[0]}")
print(f"overviews {ds.overviews(1) or 'NONE'}")
print(f"compress {ds.profile.get('compress')}")blocks (1, 4096)
overviews NONE
compress deflate
A block height of 1 means the file is striped, not tiled. An empty overview list means every downsampled read costs full resolution. Either makes the file behave like an ordinary GeoTIFF, and no client setting fixes it.
No overviews, striping, probing, NumPy downsampling, scattered reads
1. Cause one: the file has no overviews
The largest single effect. Without overviews, GDAL produces a downsampled image by reading every full-resolution pixel and averaging.
whole scene at 1/16
with overviews 0.131 MB
without 24.135 MB
There is no client-side fix. Either rewrite the file with overviews, or accept full-resolution reads.
Where you cannot rewrite the original, an external .ovr sidecar alongside it works — but only if the server hosts it and CPL_VSIL_CURL_ALLOWED_EXTENSIONS permits GDAL to look.
2. Cause two: the file is striped
one 512 x 512 window
tiled 1.473 MB
striped 3.097 MB
A striped file stores whole rows, so a narrow window reads full-width strips. The waste grows with the raster's width: here a 512-pixel window on a 4096-wide file read 8 times the pixels it needed, and on a full 10,980-pixel Sentinel-2 band it would be 21 times.
Again, no client fix.
3. Cause three: sidecar probing
default settings 18 requests
GDAL_DISABLE_READDIR_ON_OPEN 9 requests
Nine requests were GDAL looking around the file: a directory listing and HEAD requests for .aux.xml, .aux, .xml and .msk sidecars, each in two spellings, that do not exist. Each 404 is a full round trip, and none of them transfers a byte of the image.
rasterio.Env(GDAL_DISABLE_READDIR_ON_OPEN="EMPTY_DIR",
CPL_VSIL_CURL_ALLOWED_EXTENSIONS=".tif")This is the highest-value client-side setting, and against a high-latency store it often dominates everything else.
4. Cause four: reading full resolution and downsampling in NumPy
data = ds.read(1)[::16, ::16] # reads everything
data = ds.read(1, out_shape=(h // 16, w // 16),
resampling=Resampling.average) # reads an overviewThe first line is the single most common way to accidentally download a whole file. It looks like a downsampled read and is not.
5. Cause five: an access pattern that is not tile-shaped
1 px 393.2 kB 393,216 bytes/px
256 px 393.2 kB 6.0 bytes/px
512 px 1,489.1 kB 5.7 bytes/px
One pixel costs a whole tile. Sampling a raster at a thousand scattered points can cost a thousand tiles; reading the block containing them costs a handful.
A diagnostic that names which cause it is
The data used on this page
A 4096 × 4096 crop of the 10 m near-infrared band of Sentinel-2 scene S2A_30UVD_20250922_0_L2A from Earth Search, written three ways with the same deflate compression: striped.tif (one row per strip, no overviews), tiled.tif (512-pixel tiles, no overviews) and cog.tif (512-pixel tiles and four overview levels, written by GDAL's COG driver so the layout is right too). A small HTTP server in a separate process answers byte-range requests for them and logs every request, including HEAD requests and 404s; cost(fn, path, **gdal) runs fn on a fresh URL under the given GDAL settings and returns the number of requests and the bytes sent. The same harness measures the COG explainer.
import atexit
import itertools
import os
import subprocess
import sys
import pystac_client
import rasterio
from rasterio.windows import Window
item = next(pystac_client.Client.open("https://earth-search.aws.element84.com/v1")
.search(ids=["S2A_30UVD_20250922_0_L2A"], collections=["sentinel-2-l2a"]).items())
if not os.path.exists("nir.tif"): # a 4096 x 4096 crop of the 10 m near-infrared band
with rasterio.open(item.assets["nir"].href) as src:
win = Window(3000, 3000, 4096, 4096)
profile = {"driver": "GTiff", "dtype": "uint16", "count": 1, "width": 4096, "height": 4096,
"crs": src.crs, "transform": src.window_transform(win)}
with rasterio.open("nir.tif", "w", **profile) as dst:
dst.write(src.read(1, window=win), 1)
SERVER = r"""
import http.server, os, sys
class Counting(http.server.SimpleHTTPRequestHandler):
def do_GET(self):
path = self.translate_path(self.path.split("?")[0])
if not os.path.isfile(path): # sidecar probes and directory listings
self.send_error(404)
with open("requests.log", "a") as log:
log.write("0\n")
return
size = os.path.getsize(path)
first, last = 0, size - 1
if "Range" in self.headers:
a, b = self.headers["Range"].split("=")[1].split("-")
first, last = int(a), min(int(b) if b else size - 1, size - 1)
self.send_response(206)
self.send_header("Content-Range", f"bytes {first}-{last}/{size}")
else:
self.send_response(200)
self.send_header("Content-Length", str(last - first + 1))
self.send_header("Accept-Ranges", "bytes")
self.end_headers()
with open(path, "rb") as f:
f.seek(first)
self.wfile.write(f.read(last - first + 1))
with open("requests.log", "a") as log: # one line per GET: bytes sent
log.write(f"{last - first + 1}\n")
def do_HEAD(self):
path = self.translate_path(self.path.split("?")[0])
if os.path.isfile(path):
self.send_response(200)
self.send_header("Content-Length", str(os.path.getsize(path)))
self.end_headers()
else:
self.send_error(404)
with open("requests.log", "a") as log:
log.write("0\n")
def log_message(self, *args):
pass
server = http.server.ThreadingHTTPServer(("127.0.0.1", 0), Counting)
print(server.server_port, flush=True)
server.serve_forever()
"""
# Its own process: GDAL holds the GIL while it waits for a response, so a server
# thread inside this interpreter would never get to answer.
proc = subprocess.Popen([sys.executable, "-c", SERVER], stdout=subprocess.PIPE, stderr=subprocess.DEVNULL, text=True)
atexit.register(proc.terminate)
BASE = f"http://127.0.0.1:{proc.stdout.readline().strip()}"
_fresh = itertools.count()
def served(path):
"""A URL no earlier read has cached: a fresh hard link to the file, so GDAL sees a new name."""
stem, ext = os.path.splitext(path)
link = f"{stem}.read{next(_fresh)}{ext}"
if os.path.exists(link):
os.remove(link)
os.link(path, link)
return f"{BASE}/{link}"
def cost(fn, path, **gdal):
"""Requests and bytes the server sent while fn(url) ran, under the given GDAL settings."""
open("requests.log", "w").close()
with rasterio.Env(**gdal):
fn(served(path))
sizes = [int(line) for line in open("requests.log")]
return len(sizes), sum(sizes)
from rasterio.enums import Resampling
from rasterio.shutil import copy as rio_copy
for name, tiled in (("striped.tif", False), ("tiled.tif", True)):
with rasterio.open("nir.tif") as src: # the same pixels, three layouts, same compression
profile = src.profile | {"compress": "deflate", "predictor": 2}
if tiled:
profile |= {"tiled": True, "blockxsize": 512, "blockysize": 512}
with rasterio.open(name, "w", **profile) as dst:
dst.write(src.read(1), 1)
rio_copy("nir.tif", "cog.tif", driver="COG", BLOCKSIZE=512, COMPRESS="DEFLATE", PREDICTOR=2,
OVERVIEW_RESAMPLING="AVERAGE", OVERVIEW_COUNT=4) # GDAL's COG driver: the layout, too
print({name: os.path.getsize(name) for name in ("striped.tif", "tiled.tif", "cog.tif")}){'striped.tif': 24865194, 'tiled.tif': 24135158, 'cog.tif': 33122273}
Example 1 — a diagnostic that names the cause
import time
import rasterio
from rasterio.enums import Resampling
def diagnose_remote(url, settings=None):
"""Identify which of the five causes applies."""
env = settings or dict(GDAL_DISABLE_READDIR_ON_OPEN="EMPTY_DIR",
CPL_VSIL_CURL_ALLOWED_EXTENSIONS=".tif")
problems = []
start = time.time()
with rasterio.Env(**env):
with rasterio.open(url) as ds:
open_ms = (time.time() - start) * 1000
block_h, block_w = ds.block_shapes[0]
overviews = ds.overviews(1)
print(f" {ds.width} x {ds.height} {ds.dtypes[0]}, "
f"{ds.count} band(s)")
print(f" blocks {(block_h, block_w)}, overviews "
f"{overviews or 'NONE'}, compress "
f"{ds.profile.get('compress')}")
print(f" open took {open_ms:.0f} ms")
if block_h == 1 or block_h == ds.height:
problems.append(
f"striped ({block_h} x {block_w}): a narrow window reads "
f"full-width rows — {ds.width / 512:.0f}x more pixels than "
"a 512-wide window needs")
if not overviews:
problems.append(
"no overviews: any downsampled read costs the full "
"resolution. Rewrite the file or add an .ovr sidecar")
elif max(overviews) * max(block_h, block_w) < max(ds.width, ds.height):
problems.append(
f"coarsest overview is {max(ds.width, ds.height) / max(overviews):.0f} px "
"— a full-extent view still needs several requests")
if open_ms > 500:
problems.append(
"slow open: check GDAL_DISABLE_READDIR_ON_OPEN and the "
"latency to the host")
for p in problems:
print(f" ! {p}")
return problems
for path in ("striped.tif", "cog.tif"):
print(path)
diagnose_remote(served(path))
print("the archive's own band")
remote = diagnose_remote(item.assets["nir"].href)striped.tif
4096 x 4096 uint16, 1 band(s)
blocks (1, 4096), overviews NONE, compress deflate
open took 1 ms
! striped (1 x 4096): a narrow window reads full-width rows — 8x more pixels than a 512-wide window needs
! no overviews: any downsampled read costs the full resolution. Rewrite the file or add an .ovr sidecar
cog.tif
4096 x 4096 uint16, 1 band(s)
blocks (512, 512), overviews [2, 4, 8, 16], compress deflate
open took 1 ms
the archive's own band
10980 x 10980 uint16, 1 band(s)
blocks (1024, 1024), overviews [2, 4, 8, 16], compress deflate
open took 308 ms
Example 2 — proving that the settings are the problem
import time
import rasterio
from rasterio.windows import Window
CONFIGS = {
"defaults": {},
"no readdir": {"GDAL_DISABLE_READDIR_ON_OPEN": "EMPTY_DIR"},
"no readdir + extensions": {
"GDAL_DISABLE_READDIR_ON_OPEN": "EMPTY_DIR",
"CPL_VSIL_CURL_ALLOWED_EXTENSIONS": ".tif",
},
"everything": {
"GDAL_DISABLE_READDIR_ON_OPEN": "EMPTY_DIR",
"CPL_VSIL_CURL_ALLOWED_EXTENSIONS": ".tif",
"VSI_CACHE": "TRUE",
"GDAL_HTTP_MERGE_CONSECUTIVE_RANGES": "YES",
},
}
def compare_settings(path, windows):
"""Count what the same reads cost under several GDAL configurations."""
def read(url):
with rasterio.open(url) as ds:
for window in windows:
ds.read(1, window=window)
for name, env in CONFIGS.items():
start = time.time()
requests, sent = cost(read, path, **env)
print(f" {name:26} {requests:3d} requests {sent / 1e6:6.3f} MB "
f"{(time.time() - start) * 1000:5.0f} ms")
scattered = [Window(0, 0, 256, 256), Window(3000, 500, 256, 256),
Window(1000, 3500, 256, 256), Window(3500, 3500, 256, 256)]
compare_settings("cog.tif", scattered) defaults 18 requests 4.996 MB 27 ms
no readdir 9 requests 4.996 MB 22 ms
no readdir + extensions 9 requests 4.996 MB 23 ms
everything 9 requests 5.012 MB 22 ms
The milliseconds here are to a server on the same machine, where a round trip costs almost nothing; against object storage each extra request is tens of milliseconds, which is why the request count is the number to watch. Measure each configuration on a URL nothing has read before — GDAL caches across Env blocks, so a second configuration measured on the same URL benefits from the first's cache. cost does that by handing out a fresh name every time.
Example 3 — restructuring an access pattern
import numpy as np
import rasterio
from rasterio.windows import Window
def sample_many(url, xs, ys, band=1, max_block_px=50_000_000):
"""One block read instead of a read per point."""
with rasterio.open(url) as ds:
rows, cols = rasterio.transform.rowcol(ds.transform, xs, ys)
rows, cols = np.asarray(rows), np.asarray(cols)
row0, row1 = max(rows.min(), 0), min(rows.max() + 1, ds.height)
col0, col1 = max(cols.min(), 0), min(cols.max() + 1, ds.width)
span = (row1 - row0) * (col1 - col0)
print(f" {len(xs)} points span {col1 - col0} x {row1 - row0} px "
f"({span / 1e6:.1f} Mpx)")
if span > max_block_px:
print(" ! too spread out for one block — cluster the points and "
"read one block per cluster")
return None
block = ds.read(band, window=Window(col0, row0, col1 - col0, row1 - row0))
return block[rows - row0, cols - col0]
rng = np.random.default_rng(0)
with rasterio.open("cog.tif") as ds:
left, top = ds.transform * (1024, 1024)
xs = left + rng.uniform(0, 20_480, 1000) # 1,000 points in a 2048-pixel square
ys = top - rng.uniform(0, 20_480, 1000)
values = []
requests, sent = cost(lambda url: values.append(sample_many(url, xs, ys)), "cog.tif",
GDAL_DISABLE_READDIR_ON_OPEN="EMPTY_DIR")
print(f" one block: {requests} requests, {sent / 1e6:.2f} MB")
def one_by_one(url):
with rasterio.open(url) as ds:
values.append(np.array([v[0] for v in ds.sample(zip(xs, ys))]))
requests, sent = cost(one_by_one, "cog.tif", GDAL_DISABLE_READDIR_ON_OPEN="EMPTY_DIR")
print(f" ds.sample: {requests} requests, {sent / 1e6:.2f} MB; same values: {np.array_equal(*values)}") 1000 points span 2047 x 2042 px (4.2 Mpx)
one block: 6 requests, 5.92 MB
ds.sample: 18 requests, 6.44 MB; same values: True
ds.sample() issues one read per point, and GDAL's block cache decides what that costs. Here the thousand points fell in a 2048-pixel square, 16 tiles, and the cache absorbed nearly all of it: 18 requests and 6.44 MB against the block read's 6 and 5.92 MB. Spread the same points over a whole scene and the cache stops helping; a block read, or one per cluster, is predictable either way.
Example 4 — the five causes, measured
def window(url, size=512, at=1792):
with rasterio.open(url) as ds:
ds.read(1, window=Window(at, at, size, size))
def thumbnail(url):
with rasterio.open(url) as ds:
ds.read(1, out_shape=(256, 256), resampling=Resampling.average)
def numpy_thumbnail(url):
with rasterio.open(url) as ds:
ds.read(1)[::16, ::16]
def four_windows(url):
with rasterio.open(url) as ds:
for w in scattered:
ds.read(1, window=w)
FAST = {"GDAL_DISABLE_READDIR_ON_OPEN": "EMPTY_DIR"}
rows = [
("1. thumbnail, file with overviews", thumbnail, "cog.tif", FAST),
("1. thumbnail, file without overviews", thumbnail, "tiled.tif", FAST),
("2. 512 px window, tiled", window, "tiled.tif", FAST),
("2. 512 px window, striped", window, "striped.tif", FAST),
("3. four windows, default settings", four_windows, "cog.tif", {}),
("3. four windows, no readdir", four_windows, "cog.tif", FAST),
("4. thumbnail by out_shape", thumbnail, "cog.tif", FAST),
("4. thumbnail by [::16, ::16]", numpy_thumbnail, "cog.tif", FAST),
]
for label, fn, path, env in rows:
requests, sent = cost(fn, path, **env)
print(f"{label:38} {requests:3d} requests {sent / 1e6:8.3f} MB")
for size in (1, 256, 512):
requests, sent = cost(lambda url: window(url, size), "cog.tif", **FAST)
print(f"5. {size:4d} px window {sent / 1e3:9,.1f} kB {sent / size ** 2:10,.1f} bytes/px")1. thumbnail, file with overviews 3 requests 0.131 MB
1. thumbnail, file without overviews 5 requests 24.135 MB
2. 512 px window, tiled 4 requests 1.473 MB
2. 512 px window, striped 5 requests 3.097 MB
3. four windows, default settings 18 requests 4.996 MB
3. four windows, no readdir 9 requests 4.996 MB
4. thumbnail by out_shape 3 requests 0.131 MB
4. thumbnail by [::16, ::16] 5 requests 24.160 MB
5. 1 px window 393.2 kB 393,216.0 bytes/px
5. 256 px window 393.2 kB 6.0 bytes/px
5. 512 px window 1,472.8 kB 5.6 bytes/px
Why some causes cannot be fixed from the client
Why a missing overview is unfixable from the client
Downsampling requires reading the values being averaged. Without overviews those values exist only at full resolution, so producing a 256 × 256 thumbnail from a 4096 × 4096 image genuinely requires reading 16.7 million pixels.
No client setting changes that arithmetic. The saving from overviews is that the averaging was done once, at write time, and stored — 0.131 MB instead of 24.135 MB.
The only remote workaround is an external .ovr file, which contains the same pre-computed pyramid alongside the original. It needs to be hosted next to the file and permitted by CPL_VSIL_CURL_ALLOWED_EXTENSIONS.
Why striping wastes so much
A strip spans the full width of the image. Reading a 512-pixel-wide window from a 4096-wide striped file transfers 4096 pixels per row and uses 512 — eight times the data.
On a Sentinel-2 band at 10,980 pixels wide the ratio is 21. Measured on the 4096-wide test raster, the striped file transferred 3.097 MB against the tiled file's 1.473 MB.
Striping is the default in many writers, which is why "it is a GeoTIFF on S3" is not the same as "it is a COG".
Why the sidecar probes matter so much
GDAL looks for auxiliary files when opening: statistics, external overviews, masks, world files. Locally these are stat calls costing microseconds.
Over HTTP each is a request. Four missing sidecars is four round trips — 400 ms against a store 100 ms away, before a single byte of data.
The measured effect was ten requests down to six. In wall-clock terms on a high-latency link, that is often the largest single improvement available, and it is one line of configuration.
Why scattered access is the worst pattern
Tiles are compressed independently, so the minimum transfer for any pixel is its whole tile — 393 kB on the file measured here.
Reading a thousand scattered pixels therefore costs up to a thousand tiles — here that bound would be 393 MB — to obtain a thousand values. Reading a block containing them costs a few tiles.
The general principle: make the access pattern match the storage layout. Where points are genuinely spread across a large raster, cluster them and read one block per cluster.
Overview lists, striped blocks and the GDAL cache
- Check
ds.overviews(1)first. An empty list is the most common cause. block_shapes[0][0] == 1means striped, and no client setting helps.GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIRremoved 9 of 18 requests here.ds.read(1)[::16, ::16]reads everything. Useout_shape.- A 1-pixel read costs a whole tile.
- Open once and reuse the handle; a fresh open costs two requests and clears the cache.
- GDAL caches across
Envblocks, so benchmark configurations in separate processes. - Requester-pays buckets need
AWS_REQUEST_PAYER=requester, or every read fails.
Internal links
- How to read a COG from a URL without downloading the whole file — the settings in full
- Cloud-optimised GeoTIFF explained — why tiles and overviews matter
- How to write a cloud-optimised GeoTIFF in Python — fixing the file
- GDAL cannot open an S3 or HTTPS path — when the open fails outright
- How to read spatial data from S3 and other object storage — credentials and requester-pays
- How to choose chunk and tile sizes that actually help — matching access to layout
- How to profile a slow Python GIS script and find the real bottleneck — general profiling
- Cloud-native geospatial explained — the wider architecture
FAQ
Why does reading a small window download the whole file?
Either the file is striped, so a window reads full-width rows, or you asked for a downsampled read from a file with no overviews.
How do I check whether a remote file is really a COG?
Open it and print ds.block_shapes[0] and ds.overviews(1). A block height of 1 means striped; an empty overview list means no pyramid.
What is the single most useful GDAL setting?
GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR. It removed four of ten requests in a benchmark, all of them probes for files that do not exist.
Why is my thumbnail read so slow?
Almost certainly ds.read(1)[::16, ::16], which reads everything before subsampling. Use out_shape with a resampling method.
Can I fix a file that has no overviews without rewriting it?
Only by adding an external .ovr sidecar next to it, hosted on the same server and permitted by CPL_VSIL_CURL_ALLOWED_EXTENSIONS.
Why is sampling points so slow?
Each point can cost a whole tile — 393 kB here. Read one block covering the points and index into it instead.
Does increasing the cache help?
Only if reads are clustered or repetitive. In a four-scattered-window benchmark, enabling VSI_CACHE changed nothing.