Surveillance par séries temporelles de la turbidité#

Mots-clés données utilisées ; wq_annual, données utilisées ; landsat 5, données utilisées ; landsat 7, données utilisées ; landsat 8, données utilisées ; landsat 9, données utilisées ; sentinel-2

Aperçu#

Ce notebook fait partie d’une collection de notebooks d’analyse de la qualité de l’eau qui utilisent des données d’observation de la Terre pour évaluer les conditions des eaux de surface. Ensemble, ces notebooks offrent une vision plus globale de la qualité de l’eau en permettant une analyse à travers de multiples indicateurs.

Ce carnet se concentre sur les matières en suspension totales/turbidité.

Le carnet de jeu de données du Service de surveillance de la qualité de l’eau (WQMS) présente les diverses variables disponibles dans le produit annuel de qualité de l’eau. Cependant, on peut vouloir surveiller ces variables sur des échelles de temps plus courtes, comme nous le faisons avec la turbidité dans ce carnet.

Description#

Ce carnet montre comment surveiller les matières en suspension totales (TSM), un indicateur de turbidité, sur la base d’une scène individuelle et les agréger en médianes spatiales mensuelles.

  1. Charger les données annuelles

  2. Appliquer des algorithmes aux pas de temps individuels

  3. Tracer des séries temporelles

Commencer#

Pour exécuter cette analyse, exécutez toutes les cellules du bloc-notes, en commençant par la cellule « Charger les packages ».

Charger des paquets#

Importez les packages Python utilisés pour l’analyse.

Pour importer le paquet water_quality, veuillez d’abord exécuter git clone -b develop https://github.com/digitalearthafrica/deafrica_water_quality.git dans le terminal et mettre à jour le chemin ci-dessous pour le pointer.

[1]:

%matplotlib inline

import datacube
import numpy as np
import pandas as pd
import xarray as xr
import matplotlib.pyplot as plt
import plotly.graph_objects as go

from deafrica_tools.plotting import display_map
from deafrica_tools.datahandling import load_ard

import sys

sys.path.append("/home/jovyan/deafrica_water_quality/src")
from water_quality.mapping.algorithms import (
    set_wq_algorithms,
    run_wq_algorithms,
    harmonize_wq_variables,
    normalize_wq_variables
)
from water_quality.mapping.pixel_correction import R_correction

Se connecter au datacube#

Connectez-vous au datacube pour que nous puissions accéder aux données de DE Africa. Le paramètre « app » est un nom unique pour l’analyse qui est basé sur le nom du fichier du notebook.

[2]:

dc = datacube.Datacube(app="TSM_monitoring")

Charger les données annuelles#

Nous chargeons d’abord les données du produit wq_annual pour former une référence annuelle. L’emplacement par défaut est le lac Manyara dans le nord de la Tanzanie.

[3]:

query = {
    "product": "wq_annual",
    "x": (35.72, 35.92),
    "y": (-3.90, -3.35),
    "crs": "EPSG:4326",
    "time": ("2020-01-01", "2025-12-31"),
    "output_crs": "EPSG:6933",
    "resolution": (-10, 10),
    "measurements": [
        "tsm",  # Turbidity
        "clear_water",  # Clear water mask
        "water_mask",  # Water presence mask
    ],
}

[4]:

display_map(x=query["x"], y=query["y"])

[4]:
Make this Notebook Trusted to load map: File -> Trust Notebook
[5]:

ds = dc.load(**query)

Appliquer des algorithmes aux pas de temps individuels#

Pour la surveillance à court terme, nous chargeons l’ensemble complet des observations disponibles dans DE Africa à partir des capteurs Thematic Mapper (TM ; Landsat 5 et 7), Operational Land Imager (OLI ; Landsat 8 et 9) et Multispectral Imager (MSI ; Sentinel-2).

Charger des scènes satellites dans la période de temps#

[6]:

# Dictionary to map load_ard bands to those used by water_quality package
measurements = {
    "oli": {
        "blue": "oli02",
        "green": "oli03",
        "red": "oli04",
        "nir": "oli05",
        "swir_1": "oli06",
        "swir_2": "oli07",
        "pixel_quality": "oli_pq",
    },
    "msi": {
        "coastal_aerosol": "msi01",
        "blue": "msi02",
        "green": "msi03",
        "red": "msi04",
        "red_edge_1": "msi05",
        "red_edge_2": "msi06",
        "red_edge_3": "msi07",
        "nir": "msi08",
        "nir_narrow": "msi8a",
        "swir_1": "msi11",
        "swir_2": "msi12",
        "qa": "msi_pq",
    },
    "tm": {
        "blue": "tm01",
        "green": "tm02",
        "red": "tm03",
        "nir": "tm04",
        "swir_1": "tm05",
        "swir_2": "tm07",
        "pixel_quality": "tm_pq",
    },
}

[7]:

# Map sensors to products
products = {
    "tm": ["ls5_sr", "ls7_sr"],
    "oli": ["ls8_sr", "ls9_sr"],
    "msi": ["s2_l2a"],
}

[8]:

# Sensor end dates
SENSOR_END_DATES = {
    "tm": "2023-12-31",  # Landsat 7 science imaging was suspended in Jan 2024
    "oli": None,  # Landsat 8/9 (ongoing)
    "msi": None,  # Sentinel-2 (ongoing)
}

Nous chargeons les données pour 2023-2025 à 480 m pour les traiter en 20 minutes et dans les limites de mémoire du sandbox standard (16 Go de mémoire). Nous utilisons également des blocs dask pour retarder le chargement en mémoire jusqu’à ce que cela soit nécessaire.

[9]:

# Dictionary to store each sensor dataset
ds_list = {}

x_range = (35.72, 35.92)
y_range = (-3.90, -3.35)  # Lake Manyara
start_time = "2023-01-01"
end_time = "2025-12-31"

for sensor in products.keys():

    print(f"Loading {sensor} products: {products[sensor]}")

    sensor_end = SENSOR_END_DATES[sensor]
    effective_end = min(end_time, sensor_end) if sensor_end else end_time

    if start_time > effective_end:
        print(f"Skipping {sensor} (no valid time range)")
        continue

    print(f"Loading {sensor}")
    print(f"time: {start_time} → {effective_end}")

    ds_list[sensor] = load_ard(
        dc=dc,
        products=products[sensor],
        x=x_range,
        y=y_range,
        time=(start_time, effective_end),
        output_crs="EPSG:6933",
        resolution=(-480, 480),
        verbose=True,
        dask_chunks={},
        measurements=list(measurements[sensor].keys()),
    ).rename(measurements[sensor])

    # Get the QA band name from measurements dict
    qa_band = measurements[sensor].get("qa") or measurements[sensor].get(
        "pixel_quality"
    )

    # The water_quality package assumes reflectances between 0-10,000
    # We only scale Landsat as Sentinel-2 is already scaled coming out of load_ard
    bands_to_scale = [b for b in ds_list[sensor].data_vars if b != qa_band]
    if sensor in ["tm", "oli"]:
        ds_list[sensor][bands_to_scale] = ds_list[sensor][bands_to_scale] * 10000

    print(f"Finished loading {sensor}\n")

Loading tm products: ['ls5_sr', 'ls7_sr']
Loading tm
time: 2023-01-01 → 2023-12-31
Using pixel quality parameters for USGS Collection 2
Finding datasets
    ls5_sr
    ls7_sr
Applying pixel quality/cloud mask
Re-scaling Landsat C2 data
Returning 52 time steps as a dask array
Finished loading tm

Loading oli products: ['ls8_sr', 'ls9_sr']
Loading oli
time: 2023-01-01 → 2025-12-31
Using pixel quality parameters for USGS Collection 2
Finding datasets
    ls8_sr

/opt/venv/lib/python3.12/site-packages/deafrica_tools/datahandling.py:565: FutureWarning: In a future version of xarray the default value for compat will change from compat='no_conflicts' to compat='override'. This is likely to lead to different results when combining overlapping variables with the same name. To opt in to new defaults and get rid of these warnings now use `set_options(use_new_combine_kwarg_defaults=True) or set compat explicitly.
  ds = xr.merge([ds_data, ds_masks])

    ls9_sr
Applying pixel quality/cloud mask
Re-scaling Landsat C2 data
Returning 496 time steps as a dask array
Finished loading oli

Loading msi products: ['s2_l2a']
Loading msi
time: 2023-01-01 → 2025-12-31
Using pixel quality parameters for Sentinel 2
Finding datasets
    s2_l2a

/opt/venv/lib/python3.12/site-packages/deafrica_tools/datahandling.py:565: FutureWarning: In a future version of xarray the default value for compat will change from compat='no_conflicts' to compat='override'. This is likely to lead to different results when combining overlapping variables with the same name. To opt in to new defaults and get rid of these warnings now use `set_options(use_new_combine_kwarg_defaults=True) or set compat explicitly.
  ds = xr.merge([ds_data, ds_masks])

Applying pixel quality/cloud mask
Returning 983 time steps as a dask array
Finished loading msi


/opt/venv/lib/python3.12/site-packages/deafrica_tools/datahandling.py:565: FutureWarning: In a future version of xarray the default value for compat will change from compat='no_conflicts' to compat='override'. This is likely to lead to different results when combining overlapping variables with the same name. To opt in to new defaults and get rid of these warnings now use `set_options(use_new_combine_kwarg_defaults=True) or set compat explicitly.
  ds = xr.merge([ds_data, ds_masks])

Appliquer l’algorithme TSM#

Nous fournissons des fonctions optionnelles pour calibrer les valeurs avec des données in-situ, dont la dérivation est mise en évidence dans le Guide de l’utilisateur de la qualité de l’eau de DE Africa. L’utilisateur peut souhaiter modifier ces fonctions.

Remarque : Les cellules de traitement TSM ci-dessous prennent \(\sim\) 20 minutes avec une résolution de 480 m pour 2023-2025 sur le lac Manyara dans le Sandbox (2 cœurs, 16 Go de mémoire), qui a une surface de 470 :math:textrm{km}^2` <https://en.wikipedia.org/wiki/Lake_Manyara>`__.

[10]:

def calibrate_tsm(data):
    a = 0
    b = 0.5824
    c = -2.0559
    return a * data**2 + b * data + c

[11]:

# No suffix for individual time steps
suffix = ""

# Set group of algorithms for TSM
_, ALGORITHMS_TSM = set_wq_algorithms(suffix)

Appliquer les étapes de prétraitement et l’algorithme TSM#

Les étapes sont appliquées séquentiellement en utilisant les fonctions du paquet water_quality.

Nous alignons d’abord le masque d’eau claire du produit annuel de qualité de l’eau (où l’indice des algues flottantes est inférieur à un seuil de 0,05) pour correspondre aux dimensions spatio-temporelles des observations. Nous appliquons ensuite une correction de Rayleigh des pixels sombres atmosphériques. Enfin, nous exécutons les algorithmes TSM et appliquons des étapes d’harmonisation, de normalisation et de calibration.

Plus de détails sur les étapes sont disponibles dans le Guide de l’utilisateur de la qualité de l’eau de DE Africa.

[12]:

%%time

ds_dict_tsm = {}

for sensor in products.keys():
    print(sensor)

    print("Aligning clear water mask...")
    water_mask_aligned = ds.clear_water.reindex_like(ds_list[sensor], method="nearest")

    print("Applying Rayleigh correction...")
    ds_corrected = R_correction(
        ds_list[sensor], water_mask_aligned, instrument=sensor, drop=True
    )

    print("Applying clear water mask...")
    ds_masked = ds_corrected.where(water_mask_aligned)

    print("Running TSM algorithms...")
    tsm_ds = run_wq_algorithms(
        instrument_data={sensor: ds_masked}, algorithms_group=ALGORITHMS_TSM
    )

    print("Applying harmonisation...")
    tsm_ds = harmonize_wq_variables(tsm_ds)

    print("Applying normalisation...")
    tsm_ds = normalize_wq_variables(tsm_ds)

    print("Applying calibration...")
    for var in tsm_ds.data_vars:
        tsm_ds[var] = calibrate_tsm(tsm_ds[var])

    ds_dict_tsm[sensor] = tsm_ds

    print(f"Finished {sensor}\n")

tm
Aligning clear water mask...
Applying Rayleigh correction...
Applying clear water mask...
Running TSM algorithms...
Applying harmonisation...

/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: overflow encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in sqrt
  return self.func(*new_argspec)
Variable oli01 expected but not found in the dataset - (non-fatal error) for instrument oli

Applying normalisation...
Applying calibration...
Finished tm

oli
Aligning clear water mask...
Applying Rayleigh correction...
Applying clear water mask...
Running TSM algorithms...
Applying harmonisation...

/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/dask/_task_spec.py:768: RuntimeWarning: invalid value encountered in power
  return self.func(*new_argspec)
/opt/venv/lib/python3.12/site-packages/rasterio/warp.py:385: NotGeoreferencedWarning: Dataset has no geotransform, gcps, or rpcs. The identity matrix will be returned.
  dest = _reproject(

Applying normalisation...
Applying calibration...
Finished oli

msi
Aligning clear water mask...
Applying Rayleigh correction...
Applying clear water mask...
Running TSM algorithms...
Applying harmonisation...
Applying normalisation...
Applying calibration...
Finished msi

CPU times: user 3min, sys: 11.5 s, total: 3min 11s
Wall time: 17min 46s

[13]:

# Concatenate values across sensors
tsm_ds = xr.concat(ds_dict_tsm.values(), dim="time").sortby("time")

[14]:

# Stack values across the variables dimension
tsm_da = tsm_ds.to_stacked_array(
    new_dim="tsm_measures",  # New dimension for stacking
    sample_dims=["time", "y", "x"],  # Keep spatial + temporal dims
    variable_dim="tsm_wq_vars",  # Stores the names of the algorithms
    name="tsm",
)

[15]:

# Take the median across the variables dimension
tsm_final = tsm_da.median(dim="tsm_measures")

[16]:

# Calculate the monthly mean image
tsm_monthly = tsm_final.resample(time="1M").mean()

[17]:

# Take the median over area
tsm_month_median = tsm_monthly.median(dim=("x", "y"))

Tracer des séries temporelles#

Ensuite, nous traçons une série temporelle interactive du TSM mensuel et du CHLA. Nous superposons également les résultats annuels sous forme de graphique en escalier pour une surveillance à plus long terme.

[18]:

# Calibrate annual data to field data
ds["tsm"] = calibrate_tsm(ds.tsm)

[19]:

# Define series for plotting monthly results
x = tsm_month_median.time.values
y = tsm_month_median.values

# Calculate quantile values for TSM
Q = np.array([0, 0.21, 0.79, 1])
quantile_vals = tsm_month_median.quantile(Q, skipna=True).values
q1, q2 = quantile_vals[1], quantile_vals[2]

# Classify points within quantiles in a vectorised way
status = xr.where(y < q1, "Low", xr.where(y < q2, "Normal", "High"))

[20]:

# Original yearly data
year_x = pd.to_datetime(ds.tsm.time.values).to_period("Y").to_timestamp()
year_y = ds.tsm.median(dim=("x", "y")).values

# Build monthly timeline covering full span
monthly_index = pd.date_range(
    start=year_x[0], end=year_x[-1] + pd.offsets.YearEnd(0), freq="ME"  # Month end
)

# Map each month to its year value
year_series = pd.Series(year_y, index=year_x)

# Forward-fill within each year
year_monthly = year_series.reindex(monthly_index, method="ffill")

# Final arrays
year_x_monthly = year_monthly.index.values
year_y_monthly = year_monthly.values

[21]:

fig = go.Figure()

# Low band
fig.add_hrect(
    y0=quantile_vals[0],
    y1=quantile_vals[1],
    fillcolor="blue",
    opacity=0.1,
    line_width=0,
)

# Normal band
fig.add_hrect(
    y0=quantile_vals[1],
    y1=quantile_vals[2],
    fillcolor="green",
    opacity=0.1,
    line_width=0,
)

# High band
fig.add_hrect(
    y0=quantile_vals[2], y1=quantile_vals[3], fillcolor="red", opacity=0.1, line_width=0
)

# Yearly step line
fig.add_trace(
    go.Scatter(
        x=year_x_monthly,
        y=year_y_monthly,
        mode="lines",
        line_shape="hv",
        line=dict(color="steelblue", width=4, dash="solid"),
        name="Yearly spatial median",
        hovertemplate="<b>Year</b>: %{x|%Y}<br>" "<b>TSM</b>: %{y:.2f}<extra></extra>",
    )
)

# Monthly points
fig.add_trace(
    go.Scatter(
        x=x,
        y=y,
        mode="lines+markers",
        name="Monthly median TSM",
        line=dict(color="grey", dash="dash"),
        marker=dict(color="black", size=8),
        customdata=status,
        hovertemplate="<b>Date</b>: %{x}<br>"
        "<b>TSM</b>: %{y:.2f}<br>"
        "<b>Status</b>: %{customdata}<extra></extra>",
    )
)

fig.update_layout(
    title="Median TSM",
    yaxis_title="TSM (relative units)",
    template="plotly_white",
    hovermode="x unified",
    dragmode="zoom",
    xaxis=dict(title="Time", type="date"),
)

fig.show()
../../../_images/sandbox_notebooks_Real_world_examples_Water_quality_monthly_monitoring_turbidity_34_0.png

Informations Complémentaires#

Licence : Le code de ce carnet est sous licence Apache, version 2.0 <https://www.apache.org/licenses/LICENSE-2.0>. Les données de Digital Earth Africa sont sous licence Creative Commons par attribution 4.0 <https://creativecommons.org/licenses/by/4.0/>.

Contact : Si vous avez besoin d’aide, veuillez poster une question sur le canal Slack Open Data Cube <http://slack.opendatacube.org/>`__ ou sur le GIS Stack Exchange en utilisant la balise open-data-cube (vous pouvez consulter les questions posées précédemment ici). Si vous souhaitez signaler un problème avec ce bloc-notes, vous pouvez en déposer un sur Github.

Version de Datacube compatible :

[22]:

print(datacube.__version__)

1.9.13

Dernier test :

[23]:

from datetime import datetime

datetime.today().strftime("%Y-%m-%d")

[23]:

'2026-03-24'