Dockerized Sentinel-2 NDVI Pipeline

Build a containerized pipeline that maps vegetation health from satellite data.

Introduction

30 Second Summary

Every vegetation map you see in a weather app or climate report starts with raw satellite data. Scientists and farmers use these maps to track crop health, detect drought, and monitor deforestation from thousands of kilometers away.

In this project, you will build a containerized Python pipeline that fetches real Sentinel-2 satellite imagery from space, computes a vegetation health index, and outputs a color-coded map. You will do all of this inside a portable container that runs anywhere.

What You'll Build

A Docker-compatible container that takes any location on Earth and produces an NDVI vegetation health map from real satellite data.

By the end of this project, you'll have:

  • A Dockerfile that builds a geospatial Python environment with all satellite processing dependencies pre-installed.
  • A Python pipeline that discovers satellite imagery using the STAC protocol, downloads spectral bands from the cloud, and computes vegetation health.
  • A reusable container you can point at any coordinates and date range to generate a classified GeoTIFF and a color-mapped PNG visualization.
  • Secret Mission: Add a multi-date comparison that shows where vegetation has grown or declined over time.

Are there any prerequisites?

This project assumes you are comfortable writing Python scripts and navigating the Linux command line. You also need Podman installed with the docker compatibility layer so that docker commands work on your system.

Before We Start

Before We Start...

Set Up the Container Environment

Processing satellite imagery requires specialized geospatial libraries like GDAL. Installing GDAL natively is notoriously painful because it depends on dozens of system-level C libraries.

Containers solve this problem. In this step, you'll build a Podman container from a pre-built GDAL image, install your Python dependencies inside it, and establish a development workflow where you edit code in Zed and execute it inside the container.

In this step, get ready to:
  • Create the project directory, Dockerfile, and requirements file.
  • Build the container image and verify it works.
  • Set up the development workflow with volume mounts.
Create the project directory and Dockerfile

You need a folder to hold your project files and a Dockerfile that describes the container environment. The Dockerfile starts from a base image that already has GDAL compiled, saving you from installing it manually.

  • Press Ctrl+Alt+T to open a terminal.
  • Move to your Desktop by running this command:
cd ~/Desktop
  • Create your project folder and move into it by running these commands:
mkdir ndvi-pipeline
cd ndvi-pipeline
  • Confirm the folder was created by running:
pwd

You should see /home/<your-username>/Desktop/ndvi-pipeline printed in the output.

  • Open the project folder in Zed by running:
zed .

Why a container instead of a virtual environment?

GDAL requires compiled C libraries that pip alone cannot install. The osgeo/gdal:ubuntu-small-latest image ships with GDAL pre-compiled on Ubuntu 24.04 with Python 3.12. You get a working geospatial stack in one line of your Dockerfile.

  • In Zed, create a new file called Dockerfile (no extension) in the ndvi-pipeline folder.
  • Paste in the following code:
FROM osgeo/gdal:ubuntu-small-latest

RUN apt-get update && apt-get -y install python3-pip --fix-missing && \
    rm -rf /var/lib/apt/lists/*

WORKDIR /app

COPY requirements.txt /app/
RUN pip install --no-cache-dir --break-system-packages -r requirements.txt

ENTRYPOINT ["python3", "/app/pipeline.py"]

What does this Dockerfile do?

  • FROM osgeo/gdal:ubuntu-small-latest starts from the OSGeo GDAL image, which already has GDAL and Python 3.12.
  • RUN apt-get ... install python3-pip installs pip so you can add Python packages.
  • --break-system-packages is required on Ubuntu 24.04 because pip normally refuses to install system-wide due to PEP 668. Inside a container, this is safe.
  • COPY requirements.txt /app/ brings your dependency list into the container before installing.
  • ENTRYPOINT sets the default command to run your pipeline script when the container starts.
  • Save Dockerfile.

Now create the requirements file that lists your Python dependencies.

  • In Zed, create a new file called requirements.txt in the ndvi-pipeline folder.
  • Paste in the following dependencies:
pystac-client==0.8.5
rioxarray==0.22.0
matplotlib==3.10.3
numpy==2.2.6

What are these packages?

  • pystac-client connects to STAC APIs to search for satellite imagery.
  • rioxarray reads raster files (GeoTIFFs) directly from URLs into xarray DataArrays with coordinate and CRS metadata.
  • matplotlib generates the color-mapped PNG visualization.
  • numpy handles the array arithmetic for NDVI calculation.
  • Save requirements.txt.
  • Confirm both files exist by running this in your terminal:
ls

You should see Dockerfile and requirements.txt listed in the output.

✔️ Awesome, I've got everything!

Great. Double check you've saved both files before moving on.

ⓧ I'd like to double check the full code

Here is the full Dockerfile:

FROM osgeo/gdal:ubuntu-small-latest

RUN apt-get update && apt-get -y install python3-pip --fix-missing && \
    rm -rf /var/lib/apt/lists/*

WORKDIR /app

COPY requirements.txt /app/
RUN pip install --no-cache-dir --break-system-packages -r requirements.txt

ENTRYPOINT ["python3", "/app/pipeline.py"]

And here is requirements.txt:

pystac-client==0.8.5
rioxarray==0.22.0
matplotlib==3.10.3
numpy==2.2.6
Build the container image and verify it works

With your Dockerfile and requirements file in place, you can build the container image. This pulls the base GDAL image and installs all your Python packages inside it.

  • In your terminal, make sure you are inside the ndvi-pipeline folder.
  • Build the container image by running:
podman build -t ndvi-pipeline .

What does this command do?

podman build reads your Dockerfile and executes each instruction to create a container image. The -t ndvi-pipeline flag tags the image with a name so you can reference it easily later. The . tells Podman to use the current directory as the build context.

The build takes a few minutes. You will see output as each layer is pulled and each RUN instruction executes.

  • Once the build finishes, verify the image works by running this import check inside the container:
podman run --rm --entrypoint python3 ndvi-pipeline -c "import rioxarray; print('ready')"

What does this verification command do?

--rm removes the container after it exits so you don't accumulate stopped containers. --entrypoint python3 overrides the Dockerfile's default entrypoint (which expects pipeline.py to exist) and runs Python directly instead.

The -c flag passes a one-line script to Python that imports rioxarray and prints ready if everything installed correctly.

You should see ready printed in your terminal. This confirms that GDAL, rasterio, and all Python dependencies installed correctly inside the container.

Build failed or import error?

  • If you see a network error during the build, check your internet connection. Podman needs to pull the base image and download packages from PyPI.
  • If you see a file-not-found error about requirements.txt, make sure the file is saved in the same directory as your Dockerfile (not inside a subfolder).
  • If the import check prints an error instead of ready, try rebuilding the image with podman build --no-cache -t ndvi-pipeline . to start fresh.

Help me debug my Podman build or import error.

Set up the development workflow with volume mounts

Right now, your Python code lives on the host machine and the runtime lives inside the container. A volume mount bridges these two worlds by making a host directory visible inside the container at runtime.

  • Test the volume mount workflow by running:
podman run --rm -v $(pwd):/app --entrypoint ls ndvi-pipeline /app

How does the volume mount work?

-v $(pwd):/app maps your current directory on the host to /app inside the container. Any file you create or edit in Zed appears instantly inside the container without rebuilding the image.

This means your development loop is: edit in Zed, run inside the container. No rebuild needed between code changes.

You should see Dockerfile and requirements.txt listed in the output. These are your host files, visible inside the container at /app.

Not seeing your files listed?

  • Make sure you are running the command from inside the ndvi-pipeline directory. $(pwd) expands to your current working directory.
  • If you see a permission denied error, check that your user has read access to the project files.

Help me troubleshoot my Podman volume mount.

Your container environment is ready. Next up, you'll write a Python script that connects to a satellite data catalog and finds real Sentinel-2 imagery.

Search for Satellite Imagery

Your container environment is built and ready to run Python code. Now it's time to point it at real satellite data.

In this step, you'll write a function that searches for Sentinel-2 satellite imagery using the STAC protocol. STAC (SpatioTemporal Asset Catalog) is the industry standard for discovering and accessing geospatial data in the cloud.

In this step, get ready to:
  • Write a STAC search function that finds the least-cloudy satellite scene over any location.
  • Extract and verify the asset URLs for Red and NIR spectral bands.
Write the STAC search function

Earth Search is a free, public STAC API hosted by Element84 that indexes satellite data collections including Sentinel-2. Your function will connect to it, filter by location and cloud cover, and return the best available scene.

  • Create a new file called pipeline.py in your project directory (the same folder as your Dockerfile).
  • In Zed, add the following imports and constants to the top of pipeline.py:
import argparse
import sys

import matplotlib.pyplot as plt
import numpy as np
import rioxarray
from pystac_client import Client

# Earth Search STAC API endpoint (free, no auth required)
STAC_API_URL = "https://earth-search.aws.element84.com/v1"
# Sentinel-2 Level-2A surface reflectance collection
COLLECTION = "sentinel-2-l2a"

What does this code do?

The imports bring in the libraries you installed via requirements.txt.

  • pystac_client handles the STAC search requests to Earth Search.
  • rioxarray will load raster data from remote URLs in the next step.
  • STAC_API_URL points to Element84's free Earth Search endpoint. No API key needed.
  • COLLECTION targets Sentinel-2 Level-2A data, which is atmospherically corrected and ready for analysis.
  • Save pipeline.py.

Now add the search function that does the heavy lifting.

  • Below the COLLECTION constant, add the search_scene function:
def search_scene(bbox, date_range, max_cloud_cover=20):
    """Search for the least-cloudy Sentinel-2 scene."""
    client = Client.open(STAC_API_URL)
    search = client.search(
        collections=[COLLECTION],
        bbox=bbox,
        datetime=date_range,
        query={"eo:cloud_cover": {"lt": max_cloud_cover}},
        sortby=[{"field": "properties.eo:cloud_cover", "direction": "asc"}],
        max_items=1,
    )
    items = list(search.items())
    if not items:
        print("No scenes found. Try a wider date range or higher cloud cover limit.")
        sys.exit(1)
    return items[0]

What does this function do?

This function queries Earth Search for Sentinel-2 imagery matching your criteria.

  • bbox is a bounding box (west, south, east, north) that defines the geographic area.
  • query={"eo:cloud_cover": {"lt": max_cloud_cover}} filters out scenes with too many clouds.
  • sortby orders results by cloud cover ascending, so the first result is the clearest scene.
  • max_items=1 means you only get back the single best match.

Why filter by cloud cover?

Satellites photograph Earth on a fixed schedule regardless of weather. A scene with 80% cloud cover is useless for vegetation analysis because you can't see the ground.

The default threshold of 20% ensures you get mostly clear imagery. If your date range has no clear scenes, the function exits with a helpful message.

  • Save pipeline.py.

Seeing a syntax error?

  • Check that the search_scene function starts at column 0 (no extra indentation).
  • Verify you have two blank lines between the COLLECTION constant and the def line.
  • Make sure the closing parenthesis for client.search() is on its own line, indented to match the opening.

Help me debug a syntax error in my search_scene function.

Test the search and extract asset URLs

Each STAC item has assets containing URLs to the actual data files.

Sentinel-2 items in Earth Search use the keys red (Band 4) and nir (Band 8) for the spectral bands you need for NDVI.

  • Test your function by running this command inside the container, overriding the entrypoint to launch an interactive Python script:
podman run --rm -v $(pwd):/app --entrypoint python3 ndvi-pipeline -c "
from pipeline import search_scene
item = search_scene([12.4, 41.8, 12.6, 42.0], '2025-06-01/2025-06-30')
print('Found:', item.id)
print('Date:', item.properties['datetime'][:10])
print('Cloud cover:', str(item.properties['eo:cloud_cover']) + '%')
print('Red:', item.assets['red'].href)
print('NIR:', item.assets['nir'].href)
"

What does this command do?

  • --entrypoint python3 overrides the Dockerfile's default entrypoint so you can run a one-off script instead of pipeline.py directly.
  • The bounding box [12.4, 41.8, 12.6, 42.0] covers an area near Rome, Italy. You can swap in any coordinates you like.
  • item.assets['red'].href and item.assets['nir'].href extract the HTTPS URLs pointing to Cloud-Optimized GeoTIFFs stored in AWS S3.

You should see output containing a scene ID, a capture date, a cloud cover percentage, and two long HTTPS URLs ending in .tif.

What are Cloud-Optimized GeoTIFFs?

A Cloud-Optimized GeoTIFF (COG) is a regular GeoTIFF file organized internally so that software can read just the portion it needs over HTTP.

Instead of downloading a 500 MB file to read a small region, your pipeline will request only the relevant tiles. This is what makes remote satellite processing practical.

No scenes found?

  • Widen your date range. A single month might not have a cloud-free scene for some locations. Try a 3-month window.
  • Check your bounding box coordinates. The order is west, south, east, north. West must be less than east, and south must be less than north.
  • Try increasing max_cloud_cover to 40 or 50 if your area is frequently cloudy.

Help me troubleshoot why my STAC search returns no scenes.

Seeing a ModuleNotFoundError?

  • Make sure you included the -v $(pwd):/app flag. Without the volume mount, the container can't see your local pipeline.py file.
  • Confirm you're running the command from your project directory (the folder containing pipeline.py).

Help me fix a ModuleNotFoundError in my container.

✔️ Awesome, I've got everything!

Great. Double check you've saved pipeline.py in Zed before moving on.

ⓧ I'd like to double check the full code

import argparse
import sys

import matplotlib.pyplot as plt
import numpy as np
import rioxarray
from pystac_client import Client

STAC_API_URL = "https://earth-search.aws.element84.com/v1"
COLLECTION = "sentinel-2-l2a"


def search_scene(bbox, date_range, max_cloud_cover=20):
    """Search for the least-cloudy Sentinel-2 scene."""
    client = Client.open(STAC_API_URL)
    search = client.search(
        collections=[COLLECTION],
        bbox=bbox,
        datetime=date_range,
        query={"eo:cloud_cover": {"lt": max_cloud_cover}},
        sortby=[{"field": "properties.eo:cloud_cover", "direction": "asc"}],
        max_items=1,
    )
    items = list(search.items())
    if not items:
        print("No scenes found. Try a wider date range or higher cloud cover limit.")
        sys.exit(1)
    return items[0]

Your pipeline can now discover satellite imagery anywhere on Earth with a single function call. Next up, you'll load those band files and compute vegetation health from their pixel values.

Compute NDVI

Your pipeline can now find the best satellite scene for any location and date range. The next challenge is turning that raw satellite data into something meaningful.

Vegetation reflects near-infrared light strongly and absorbs red light for photosynthesis. By comparing these two bands, you can calculate NDVI (Normalized Difference Vegetation Index), a value from -1 to +1 that reveals how much healthy vegetation covers an area.

In this step, get ready to:
  • Load the Red and NIR spectral bands directly from their cloud URLs.
  • Compute NDVI using band arithmetic.
  • Classify the result into vegetation categories and save it as a GeoTIFF.
Load the Red and NIR bands to compute NDVI

The NDVI formula is (NIR - Red) / (NIR + Red). To evaluate it, you need the Red (Band 4) and Near-Infrared (Band 8) spectral bands from your Sentinel-2 scene.

Both bands are stored as Cloud-Optimized GeoTIFFs on AWS, so rioxarray can stream them directly from their HTTPS URLs without downloading the entire file first.

  • In Zed, open pipeline.py and find the closing return items[0] line at the bottom of the search_scene() function.
  • Below search_scene(), add the following compute_ndvi() function:
def compute_ndvi(item):
    """Load Red and NIR bands and compute NDVI."""
    red_url = item.assets["red"].href
    nir_url = item.assets["nir"].href

    print(f"Loading Red band: {red_url}")
    red = rioxarray.open_rasterio(red_url, masked=True)

    print(f"Loading NIR band: {nir_url}")
    nir = rioxarray.open_rasterio(nir_url, masked=True)

    # Align NIR to match Red's CRS and extent
    nir_matched = nir.rio.reproject_match(red)

    # Convert to float for division
    red_float = red.astype("float32")
    nir_float = nir_matched.astype("float32")

    # Apply the NDVI formula
    ndvi = (nir_float - red_float) / (nir_float + red_float)
    return ndvi

What does this code do?

This function takes a STAC item (from your search) and pulls the two bands needed for NDVI directly from remote URLs.

  • item.assets["red"].href and item.assets["nir"].href extract the download URLs for Band 4 (Red) and Band 8 (NIR) from the STAC metadata.
  • rioxarray.open_rasterio() streams the raster data from the URL using GDAL's virtual filesystem. It returns an xarray DataArray with coordinates and CRS attached.
  • reproject_match() aligns the NIR raster to match the Red raster's coordinate system and pixel grid, handling any edge-case misalignment between bands.
  • The final three lines cast both bands to float32 (so division works correctly) and apply the NDVI formula.
  • Save pipeline.py.

You won't see output yet. The next substep adds the classification and save logic that produces the final file.

Classify the NDVI result and save it as a GeoTIFF

Raw NDVI values are continuous numbers between -1 and +1. To make the output actionable, you will bin these values into four categories: water/bare soil, sparse vegetation, moderate vegetation, and dense vegetation.

  • Below the compute_ndvi() function you just added, paste the classify_ndvi() and save_geotiff() functions:
def classify_ndvi(ndvi):
    """Classify NDVI into vegetation categories."""
    # Bins: <0 water/bare, 0-0.2 sparse, 0.2-0.4 moderate, >0.4 dense
    bins = [-1.0, 0.0, 0.2, 0.4, 1.0]
    classified = ndvi.copy()
    classified.values = np.digitize(ndvi.values, bins).astype("int32")
    return classified


def save_geotiff(classified, output_path):
    """Save classified NDVI as a GeoTIFF."""
    classified.rio.to_raster(output_path, dtype="int32")
    print(f"Saved classified GeoTIFF: {output_path}")

What does this code do?

These two functions handle the final processing stage of the pipeline.

  • numpy.digitize() takes the continuous NDVI values and assigns each pixel to a bin (1 through 4) based on where it falls in the bins thresholds.
  • The categories map to: 1 = water/bare soil (NDVI < 0), 2 = sparse vegetation (0 to 0.2), 3 = moderate vegetation (0.2 to 0.4), 4 = dense vegetation (> 0.4).
  • rio.to_raster() writes the classified array as a GeoTIFF file, preserving the CRS and spatial coordinates from the original bands.
  • Save pipeline.py.

Now test the full processing chain by running the functions inside the container.

  • Create an output directory on your host by running:
mkdir -p output
  • Run the pipeline functions inside the container by running this command:
podman run --rm -v $(pwd):/app --entrypoint python3 ndvi-pipeline -c "
from pipeline import search_scene, compute_ndvi, classify_ndvi, save_geotiff
import os
os.makedirs('/app/output', exist_ok=True)
item = search_scene([12.4, 41.8, 12.6, 42.0], '2025-06-01/2025-06-30')
ndvi = compute_ndvi(item)
classified = classify_ndvi(ndvi)
save_geotiff(classified, '/app/output/ndvi_classified.tif')
"

What does this command do?

This overrides the container's default entrypoint with python3 and passes an inline script via -c. The script imports your functions, searches for a scene over Rome (bounding box 12.4 41.8 12.6 42.0), computes NDVI, classifies it, and saves the result.

This will take 30 to 60 seconds. The container is streaming band data from AWS as Cloud-Optimized GeoTIFFs, so it only downloads the pixels it needs rather than the entire file.

You should see output like:

Loading Red band: https://...
Loading NIR band: https://...
Saved classified GeoTIFF: /app/output/ndvi_classified.tif
  • Verify the output is a valid raster by running gdalinfo inside the container:
podman run --rm -v $(pwd):/app --entrypoint gdalinfo ndvi-pipeline /app/output/ndvi_classified.tif

You should see metadata including a Coordinate System (the CRS) and a Corner Coordinates section showing the bounding box of your raster.

Seeing "No scenes found" or a timeout?

  • If no scenes are found, the date range might not have clear imagery for that area. Try widening it (e.g. 2025-04-01/2025-06-30) or increasing the cloud cover tolerance by editing the max_cloud_cover default in search_scene().
  • If the command hangs, your network may be slow downloading the bands. Wait up to 2 minutes before cancelling.
  • If you see an import error, double check that all four functions are defined in pipeline.py and the file is saved.

help me debug my NDVI pipeline. I'm getting an error when running the compute and classify functions inside my Podman container.

✔️ Awesome, I've got everything!

Great. Double check you have saved pipeline.py and that output/ndvi_classified.tif exists on your host machine.

ⓧ I'd like to double check the full code

Here is the full pipeline.py at this point:

import argparse
import sys

import matplotlib.pyplot as plt
import numpy as np
import rioxarray
from pystac_client import Client

STAC_API_URL = "https://earth-search.aws.element84.com/v1"
COLLECTION = "sentinel-2-l2a"


def search_scene(bbox, date_range, max_cloud_cover=20):
    """Search for the least-cloudy Sentinel-2 scene."""
    client = Client.open(STAC_API_URL)
    search = client.search(
        collections=[COLLECTION],
        bbox=bbox,
        datetime=date_range,
        query={"eo:cloud_cover": {"lt": max_cloud_cover}},
        sortby=[{"field": "properties.eo:cloud_cover", "direction": "asc"}],
        max_items=1,
    )
    items = list(search.items())
    if not items:
        print("No scenes found. Try a wider date range or higher cloud cover limit.")
        sys.exit(1)
    return items[0]


def compute_ndvi(item):
    """Load Red and NIR bands and compute NDVI."""
    red_url = item.assets["red"].href
    nir_url = item.assets["nir"].href

    print(f"Loading Red band: {red_url}")
    red = rioxarray.open_rasterio(red_url, masked=True)

    print(f"Loading NIR band: {nir_url}")
    nir = rioxarray.open_rasterio(nir_url, masked=True)

    nir_matched = nir.rio.reproject_match(red)

    red_float = red.astype("float32")
    nir_float = nir_matched.astype("float32")

    ndvi = (nir_float - red_float) / (nir_float + red_float)
    return ndvi


def classify_ndvi(ndvi):
    """Classify NDVI into vegetation categories."""
    bins = [-1.0, 0.0, 0.2, 0.4, 1.0]
    classified = ndvi.copy()
    classified.values = np.digitize(ndvi.values, bins).astype("int32")
    return classified


def save_geotiff(classified, output_path):
    """Save classified NDVI as a GeoTIFF."""
    classified.rio.to_raster(output_path, dtype="int32")
    print(f"Saved classified GeoTIFF: {output_path}")

Your pipeline now transforms raw satellite reflectance data into a meaningful vegetation map. Next up, you will turn this classified raster into a beautiful color-mapped PNG and wire up command-line arguments so the container works as a reusable tool.

Visualize the Result

Your pipeline can already discover satellite imagery, compute NDVI, and save a classified GeoTIFF. That GeoTIFF is machine-readable, but humans need something visual to interpret vegetation patterns at a glance.

In this step, you will add a color-mapped PNG visualization using matplotlib, then wire the entire pipeline into a single command-line tool. By the end, your container takes a bounding box and date range as arguments and produces a complete vegetation health map.

In this step, get ready to:
  • Create a color-mapped PNG visualization with a legend and title.
  • Wire up the full pipeline with argument parsing so it runs from a single command.
Create the PNG visualization function

A GeoTIFF stores precise numeric data, but it requires GIS software to view. A color-mapped PNG turns those numbers into an image anyone can open and understand immediately.

  • In pipeline.py, add the following function below save_geotiff():
def save_png(ndvi, scene_date, output_path):
    """Save a color-mapped PNG visualization."""
    fig, ax = plt.subplots(1, 1, figsize=(10, 8))
    ndvi_plot = ndvi.squeeze()
    im = ax.imshow(
        ndvi_plot.values,
        cmap="RdYlGn",
        vmin=-0.2,
        vmax=0.8,
    )
    ax.set_title(f"NDVI - {scene_date}")
    ax.set_axis_off()
    fig.colorbar(im, ax=ax, label="NDVI", shrink=0.8)
    plt.tight_layout()
    plt.savefig(output_path, dpi=150, bbox_inches="tight")
    plt.close()
    print(f"Saved PNG visualization: {output_path}")

What does this code do?

  • ndvi.squeeze() removes the band dimension so the data is a 2D array that matplotlib can plot as an image.
  • cmap="RdYlGn" applies a Red-Yellow-Green diverging colormap. Red means low NDVI (bare soil or water), green means high NDVI (dense vegetation).
  • vmin=-0.2 and vmax=0.8 set the color range so it focuses on the meaningful NDVI values instead of stretching to the full -1 to +1 range.
  • fig.colorbar() adds a legend so viewers can map colors back to NDVI values.
  • plt.close() frees memory after saving. Inside a container, matplotlib has no display to show, so you write directly to a file.
  • Save pipeline.py.

You won't see a change yet. The next chunk wires up the main function that calls this visualization.

Wire up the full pipeline with argparse

Right now your functions exist but nothing calls them. A main() function with argument parsing turns the script into a command-line tool that accepts a location and date range.

  • In pipeline.py, add the following code below the save_png() function:
def main():
    parser = argparse.ArgumentParser(description="Sentinel-2 NDVI Pipeline")
    parser.add_argument(
        "--bbox",
        type=float,
        nargs=4,
        required=True,
        metavar=("WEST", "SOUTH", "EAST", "NORTH"),
        help="Bounding box as west south east north",
    )
    parser.add_argument(
        "--date",
        type=str,
        required=True,
        help="Date range as START/END (e.g. 2025-06-01/2025-06-30)",
    )
    parser.add_argument(
        "--output-dir",
        type=str,
        default="/app/output",
        help="Directory for output files",
    )
    args = parser.parse_args()

    print("=== Sentinel-2 NDVI Pipeline ===")
    print(f"Bounding box: {args.bbox}")
    print(f"Date range: {args.date}")

    print("\n[1/4] Searching for satellite imagery...")
    item = search_scene(args.bbox, args.date)
    scene_date = item.properties["datetime"][:10]
    cloud_cover = item.properties["eo:cloud_cover"]
    print(f"Found: {item.id}")
    print(f"  Date: {scene_date}")
    print(f"  Cloud cover: {cloud_cover}%")

    print("\n[2/4] Computing NDVI...")
    ndvi = compute_ndvi(item)

    print("\n[3/4] Classifying vegetation...")
    classified = classify_ndvi(ndvi)
    geotiff_path = f"{args.output_dir}/ndvi_classified.tif"
    save_geotiff(classified, geotiff_path)

    print("\n[4/4] Generating visualization...")
    png_path = f"{args.output_dir}/ndvi_map.png"
    save_png(ndvi, scene_date, png_path)

    print("\n=== Pipeline complete ===")
    print(f"Outputs in: {args.output_dir}")


if __name__ == "__main__":
    main()

What does this code do?

  • --bbox accepts exactly 4 floats (west, south, east, north) defining the geographic area to search.
  • --date accepts a date range string in the format START/END that the STAC API uses to filter scenes.
  • --output-dir defaults to /app/output, which maps to the output/ folder on your host through the volume mount.
  • The pipeline runs all four stages in sequence: search for imagery, compute NDVI, classify vegetation, and generate the visualization.
  • if __name__ == "__main__" ensures main() only runs when the script is executed directly, not when imported.
  • Save pipeline.py.
  • Create the output directory on your host by running this command:
mkdir -p output
  • Run the full pipeline inside the container by running this command:
podman run -v $(pwd):/app ndvi-pipeline --bbox 12.4 41.8 12.6 42.0 --date 2025-06-01/2025-06-30

What is this command doing?

This runs the ndvi-pipeline container with your project directory mounted at /app. The arguments after the image name are passed directly to pipeline.py through the ENTRYPOINT.

The bounding box 12.4 41.8 12.6 42.0 covers a section of Rome, Italy. You can swap these coordinates for any location you like.

The pipeline takes a minute or two to download the satellite bands from the cloud. You should see progress messages for each of the four stages.

  • Confirm both output files exist by running:
ls output/

You should see ndvi_classified.tif and ndvi_map.png listed.

  • Open output/ndvi_map.png in any image viewer or file browser.

You should see a color-mapped image where green areas represent vegetation and brown or red areas represent bare soil or water.

Pipeline fails with no scenes found?

  • The date range might not have any low-cloud scenes. Try widening it to two or three months.
  • Make sure your bounding box coordinates are in the correct order: west, south, east, north. West must be less than east, and south must be less than north.
  • Check that the bounding box covers land. An ocean-only box will find scenes but produce a mostly-empty NDVI map.

Help me debug my NDVI pipeline run

Seeing a permission denied error on the output directory?

  • Make sure you ran mkdir -p output on your host before running the container.
  • If using SELinux, try adding :Z to the volume mount: -v $(pwd):/app:Z.

Help me fix container volume mount permissions

✔️ Awesome, I've got everything!

Great. Double check you have saved pipeline.py and both output files exist in your output/ directory.

ⓧ I'd like to double check the full code

Here is the complete pipeline.py at this point:

import argparse
import sys

import matplotlib.pyplot as plt
import numpy as np
import rioxarray
from pystac_client import Client

STAC_API_URL = "https://earth-search.aws.element84.com/v1"
COLLECTION = "sentinel-2-l2a"


def search_scene(bbox, date_range, max_cloud_cover=20):
    """Search for the least-cloudy Sentinel-2 scene."""
    client = Client.open(STAC_API_URL)
    search = client.search(
        collections=[COLLECTION],
        bbox=bbox,
        datetime=date_range,
        query={"eo:cloud_cover": {"lt": max_cloud_cover}},
        sortby=[{"field": "properties.eo:cloud_cover", "direction": "asc"}],
        max_items=1,
    )
    items = list(search.items())
    if not items:
        print("No scenes found. Try a wider date range or higher cloud cover limit.")
        sys.exit(1)
    return items[0]


def compute_ndvi(item):
    """Load Red and NIR bands and compute NDVI."""
    red_url = item.assets["red"].href
    nir_url = item.assets["nir"].href

    print(f"Loading Red band: {red_url}")
    red = rioxarray.open_rasterio(red_url, masked=True)

    print(f"Loading NIR band: {nir_url}")
    nir = rioxarray.open_rasterio(nir_url, masked=True)

    nir_matched = nir.rio.reproject_match(red)

    red_float = red.astype("float32")
    nir_float = nir_matched.astype("float32")

    ndvi = (nir_float - red_float) / (nir_float + red_float)
    return ndvi


def classify_ndvi(ndvi):
    """Classify NDVI into vegetation categories."""
    bins = [-1.0, 0.0, 0.2, 0.4, 1.0]
    classified = ndvi.copy()
    classified.values = np.digitize(ndvi.values, bins).astype("int32")
    return classified


def save_geotiff(classified, output_path):
    """Save classified NDVI as a GeoTIFF."""
    classified.rio.to_raster(output_path, dtype="int32")
    print(f"Saved classified GeoTIFF: {output_path}")


def save_png(ndvi, scene_date, output_path):
    """Save a color-mapped PNG visualization."""
    fig, ax = plt.subplots(1, 1, figsize=(10, 8))
    ndvi_plot = ndvi.squeeze()
    im = ax.imshow(
        ndvi_plot.values,
        cmap="RdYlGn",
        vmin=-0.2,
        vmax=0.8,
    )
    ax.set_title(f"NDVI - {scene_date}")
    ax.set_axis_off()
    fig.colorbar(im, ax=ax, label="NDVI", shrink=0.8)
    plt.tight_layout()
    plt.savefig(output_path, dpi=150, bbox_inches="tight")
    plt.close()
    print(f"Saved PNG visualization: {output_path}")


def main():
    parser = argparse.ArgumentParser(description="Sentinel-2 NDVI Pipeline")
    parser.add_argument(
        "--bbox",
        type=float,
        nargs=4,
        required=True,
        metavar=("WEST", "SOUTH", "EAST", "NORTH"),
        help="Bounding box as west south east north",
    )
    parser.add_argument(
        "--date",
        type=str,
        required=True,
        help="Date range as START/END (e.g. 2025-06-01/2025-06-30)",
    )
    parser.add_argument(
        "--output-dir",
        type=str,
        default="/app/output",
        help="Directory for output files",
    )
    args = parser.parse_args()

    print("=== Sentinel-2 NDVI Pipeline ===")
    print(f"Bounding box: {args.bbox}")
    print(f"Date range: {args.date}")

    print("\n[1/4] Searching for satellite imagery...")
    item = search_scene(args.bbox, args.date)
    scene_date = item.properties["datetime"][:10]
    cloud_cover = item.properties["eo:cloud_cover"]
    print(f"Found: {item.id}")
    print(f"  Date: {scene_date}")
    print(f"  Cloud cover: {cloud_cover}%")

    print("\n[2/4] Computing NDVI...")
    ndvi = compute_ndvi(item)

    print("\n[3/4] Classifying vegetation...")
    classified = classify_ndvi(ndvi)
    geotiff_path = f"{args.output_dir}/ndvi_classified.tif"
    save_geotiff(classified, geotiff_path)

    print("\n[4/4] Generating visualization...")
    png_path = f"{args.output_dir}/ndvi_map.png"
    save_png(ndvi, scene_date, png_path)

    print("\n=== Pipeline complete ===")
    print(f"Outputs in: {args.output_dir}")


if __name__ == "__main__":
    main()

Secret mission

Multi-Date Vegetation Change

Your pipeline produces a snapshot of vegetation health for a single date. But real-world analysis often requires comparing two points in time to detect change. Can you extend your pipeline to find scenes from two different periods and produce a difference map showing where vegetation has grown or declined?

Clean Up Your Resources

Clean Up Your Resources

This project runs entirely locally using Podman, so there are no ongoing cloud costs or subscriptions to worry about. Decide whether to keep your resources for future runs or delete them to free up disk space.

Resources you used:

  • Podman image ndvi-pipeline (~400-500 MB on disk).
  • Output files in the output/ directory (ndvi_classified.tif, ndvi_map.png, ndvi_change.png).
  • Project source files (Dockerfile, pipeline.py, requirements.txt).

Keep everything running

No action needed. Choose this if you plan to run the pipeline again with different coordinates or date ranges.

  • The ndvi-pipeline image stays cached locally for instant re-runs. It costs nothing to keep on disk.
  • Your output files remain in the output/ folder, ready for viewing or sharing.
  • To run the pipeline again with a new location, use the same podman run command with different --bbox and --date values.

Pause - I'll come back to this later

Remove output files to free a few MB of disk space while keeping the image for quick re-runs later.

  • Delete the generated output files by running:
  • Keep the ndvi-pipeline image cached so you can re-run the pipeline without rebuilding.
  • Your source files (Dockerfile, pipeline.py, requirements.txt) remain untouched and ready to edit.

Delete - I don't want to use this again

Remove all project resources to start fresh and reclaim disk space.

  • Remove the Podman image by running:
  • Navigate to your Desktop and remove the project folder by running:
  • This removes all source files, the Dockerfile, requirements.txt, pipeline.py, and the output/ directory with all generated maps.

That's a Wrap!

That's a Wrap!

Nice work! You've just built a fully containerized satellite imagery pipeline that fetches real Sentinel-2 data from space and transforms it into a vegetation health map.

You've learned how to:

  • Built a Dockerfile from a geospatial base image and established a container-based development workflow with volume mounts.
  • Used the STAC protocol to programmatically discover satellite imagery and load Cloud-Optimized GeoTIFFs directly from remote URLs.
  • Computed NDVI from spectral bands using array arithmetic and produced a classified vegetation map as both a GeoTIFF and a color-mapped PNG.
  • Secret Mission: Built a multi-date comparison that detects where vegetation has grown or declined over time.

Ready to quiz yourself?