=== versions === python 3.12.0 | dask 2026.8.0 | distributed 2026.8.0 | xarray 2026.9.0 | rioxarray 0.23.0 | rasterio 1.5.2 | numpy 2.5.4 machine: 16 cores / 32 threads, 68 GB RAM === Step 1: eight low-cloud Sentinel-2 scenes of tile 30UXC (2025) === 2025-03-06 S2B_MSIL2A_20250306T110829_R137_T30UXC_20250306T133711 tile cloud 3.7% no data 0.0002% 2025-04-02 S2A_MSIL2A_20250402T110651_R137_T30UXC_20250402T201716 tile cloud 1.0% no data 0.0001% 2025-04-05 S2B_MSIL2A_20250405T110619_R137_T30UXC_20250405T132414 tile cloud 0.3% no data 0.0000% 2025-04-30 S2C_MSIL2A_20250430T110641_R137_T30UXC_20250430T151423 tile cloud 6.5% no data 0.0000% 2025-05-02 S2A_MSIL2A_20250502T110701_R137_T30UXC_20250502T181413 tile cloud 2.1% no data 0.0000% 2025-06-19 S2C_MSIL2A_20250619T110641_R137_T30UXC_20250619T151116 tile cloud 0.2% no data 0.0000% 2025-07-11 S2A_MSIL2A_20250711T110651_R137_T30UXC_20250711T153113 tile cloud 5.2% no data 0.0000% 2025-11-21 S2B_MSIL2A_20251121T111259_R137_T30UXC_20251121T120801 tile cloud 4.9% no data 0.0005% 16 files on disk, 4.58 GB compressed (0.00 GB downloaded this run, 0 s) each file: 10980 x 10980 px, uint16, deflate, internal tiles (512, 512), overviews [2, 4, 8, 16] uncompressed: 8 dates x 2 bands x 10980 x 10980 uint16 = 3.86 GB (7.72 GB as float32) === Step 2: open lazily with rioxarray and choose a chunk size === opened 16 files in 0.10 s without reading any pixels Size: 2GB dask.array Coordinates: * time (time) datetime64[ns] 64B 2025-03-06 2025-04-02 ... 2025-11-21 * y (y) float64 88kB 5.8e+06 5.8e+06 5.8e+06 ... 5.69e+06 5.69e+06 * x (x) float64 88kB 6e+05 6e+05 6e+05 ... 7.098e+05 7.098e+05 spatial_ref int64 8B 0 Attributes: AREA_OR_POINT: Area _FillValue: 0 scale_factor: 1.0 add_offset: 0.0 chunk size options for a 10980 x 10980 band: 256 px: 1849 chunks per band, 0.3 MB per float32 chunk, 2 MB per 8-date block for the median 512 px: 484 chunks per band, 1.0 MB per float32 chunk, 8 MB per 8-date block for the median 1024 px: 121 chunks per band, 4.2 MB per float32 chunk, 34 MB per 8-date block for the median 2048 px: 36 chunks per band, 16.8 MB per float32 chunk, 134 MB per 8-date block for the median 4096 px: 9 chunks per band, 67.1 MB per float32 chunk, 537 MB per 8-date block for the median === Step 3: build the median NDVI as a lazy task graph === Size: 482MB dask.array Coordinates: * y (y) float64 88kB 5.8e+06 5.8e+06 5.8e+06 ... 5.69e+06 5.69e+06 * x (x) float64 88kB 6e+05 6e+05 6e+05 ... 7.098e+05 7.098e+05 spatial_ref int64 8B 0 2048 px chunks: 6136 tasks in 81 layers; output 482 MB in 36 chunks 256 px chunks: 314346 tasks, 1849 output chunks 512 px chunks: 82296 tasks, 484 output chunks 1024 px chunks: 20586 tasks, 121 output chunks 2048 px chunks: 6136 tasks, 36 output chunks 4096 px chunks: 1546 tasks, 9 output chunks two dates in 2 x 2 chunks: 180 tasks as written, 52 after dask fuses chains of element-wise steps graph as written: 180 tasks, 228 dependencies, 14 tasks on the longest path task types: concatenate 16, elementwise 104, median 4, open file 4, read chunk 48, rechunk 4 === Step 4: compute with dask diagnostics (8 threads, 2048 px chunks) === [ ] | 0% Completed | 294.90 us [####################### ] | 59% Completed | 10.01 s [########################################] | 100% Completed | 20.01 s computed in 20.1 s: 1513 tasks; process memory 0.35 GB at start, 6.64 GB at peak; mean CPU 694% (800% = 8 busy threads) median NDVI: mean 0.584; 100.00% of pixels valid; 6,874 km² above 0.6 and 192 km² below 0.1 === Step 5: peak memory, eager NumPy against dask (fresh process per run) === eager NumPy, 1 date(s): 18.9 s, peak 8.95 GB (mean NDVI 0.4903) eager NumPy, 2 date(s): 35.3 s, peak 9.92 GB (mean NDVI 0.5125) eager NumPy, 4 date(s): 62.9 s, peak 12.37 GB (mean NDVI 0.5475) dask, 1 date(s): 3.9 s, peak 2.12 GB (mean NDVI 0.4903) dask, 2 date(s): 6.1 s, peak 2.50 GB (mean NDVI 0.5125) dask, 4 date(s): 10.4 s, peak 4.12 GB (mean NDVI 0.5475) dask, 8 date(s): 19.0 s, peak 6.36 GB (mean NDVI 0.5836) eager peak grows by 1.15 GB per date: about 17 GB for 8 dates and 88 GB for the 70 or so scenes a year of this tile === Step 6: time chunk sizes and thread counts (8 dates, fresh process per run) === 256 px chunks, 2 threads: 125.0 s, peak 4.94 GB 256 px chunks, 4 threads: 88.8 s, peak 4.92 GB 256 px chunks, 8 threads: 60.3 s, peak 5.04 GB 256 px chunks, 16 threads: 55.4 s, peak 5.23 GB 512 px chunks, 1 threads: 136.5 s, peak 4.58 GB 512 px chunks, 2 threads: 71.6 s, peak 4.57 GB 512 px chunks, 4 threads: 37.2 s, peak 4.69 GB 512 px chunks, 8 threads: 20.4 s, peak 4.80 GB 512 px chunks, 16 threads: 18.5 s, peak 4.88 GB 1024 px chunks, 2 threads: 60.0 s, peak 4.45 GB 1024 px chunks, 4 threads: 34.0 s, peak 4.47 GB 1024 px chunks, 8 threads: 18.9 s, peak 4.70 GB 1024 px chunks, 16 threads: 12.8 s, peak 5.39 GB 2048 px chunks, 1 threads: 122.3 s, peak 5.01 GB 2048 px chunks, 2 threads: 67.1 s, peak 5.07 GB 2048 px chunks, 4 threads: 35.9 s, peak 5.33 GB 2048 px chunks, 8 threads: 21.4 s, peak 6.54 GB 2048 px chunks, 16 threads: 14.2 s, peak 9.00 GB 4096 px chunks, 2 threads: 67.9 s, peak 8.22 GB 4096 px chunks, 4 threads: 40.0 s, peak 10.44 GB 4096 px chunks, 8 threads: 23.1 s, peak 13.94 GB 4096 px chunks, 16 threads: 18.0 s, peak 20.92 GB 2048 px chunks, 8 threads, lock=True: 108.2 s, peak 5.03 GB LocalCluster, 4 processes x 2 threads, 2048 px chunks: 21.2 s (memory is in the worker processes) fastest: 1024 px chunks with 16 threads, 12.8 s; 2048 px: 1 thread 122 s -> 8 threads 21 s (5.7x) === Step 7: write the result to disk chunk by chunk === wrote ndvi_median_2025_30UXC.tif: 370 MB on disk in 21.1 s; memory rose by at most 2.49 GB for a 482 MB result check: 10980 x 10980 float32, EPSG:32630, deflate, 484 blocks of (512, 512) finished in 25.4 min