Skip to content

Run matched filters in PRISMA product

Author: Gonzalo Mateo-GarcĂ­a

See example of using PRISMA reader in georeader.

The Window-Wide Matched Filter is proposed in: > J. Roger, L. Guanter, J. Gorroño, and I. Irakulis-Loitxate, “Exploiting the entire near-infrared spectral range to improve the detection of methane plumes with high-resolution imaging spectrometers,” Atmospheric Measurement Techniques, vol. 17, no. 4, pp. 1333–1346, Feb. 2024, doi: 10.5194/amt-17-1333-2024.

Download PRISMA image from UNEP IMEO Cloud Storage

This requires the adlfs package:

pip install adlfs

import os

SAS_TOKEN = ""
AZURE_STORAGE_ACCOUNT=""
CONTAINER = ""

folder_dest = "../tests/data/PRISMA"
tile = "PRS_L1_STD_OFFL_20250316071846_20250316071850_0001"
path_raster_dst = os.path.join(folder_dest, f"{tile}.he5")

NEED_DOWNLOAD = not os.path.exists(path_raster_dst)
print(f"File '{path_raster_dst}' exists: {not NEED_DOWNLOAD}")

if NEED_DOWNLOAD:
    import adlfs
    fs = adlfs.AzureBlobFileSystem(account_name=AZURE_STORAGE_ACCOUNT, sas_token=SAS_TOKEN, container=CONTAINER)

if NEED_DOWNLOAD:
    os.makedirs(folder_dest, exist_ok=True)

    path_raster = f"az://{CONTAINER}/PRISMA/raster/{tile}.he5"
    if not fs.exists(path_raster):
        raise FileNotFoundError(f"File {path_raster} does not exist")

    print(f"Downloading {path_raster}...")
    fs.get(path_raster, path_raster_dst)
else:
    print(f"Skipping download, '{path_raster_dst}' already exists")

!tree {folder_dest}

PRISMA
└── PRS_L1_STD_OFFL_20250316071846_20250316071850_0001.he5

0 directories, 1 file

Load image

from georeader.readers import prisma

rst = prisma.PRISMA(path_raster_dst)
rst

/home/gonzalo/git/marshsi/.venv/lib/python3.10/site-packages/georeader/reflectance.py:9: UserWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html. The pkg_resources package is slated for removal as early as 2025-11-30. Refrain from using this package or pin to Setuptools<81.
  import pkg_resources



        File: PRISMA/PRS_L1_STD_OFFL_20250316071846_20250316071850_0001.he5
        Bounds: (53.92854309082031, 38.101470947265625, 54.339393615722656, 38.42490768432617)
        Time: 2025-03-16 07:18:46.598999+00:00
        VNIR Range: (np.float32(406.9934), np.float32(977.3654)) 63 bands
        SWIR Range: (np.float32(943.3579), np.float32(2497.1155)) 171 bands

Run WW-MF retrieval on PRISMA product

Run wide matched filter retrieval based on the work of Roger et al. 2024.

from marshsi.prismaenmap import retrieval_upv_prisma_enmap

cube_crs = retrieval_upv_prisma_enmap.MF_sunglint_combo_prisma(rst)
cube_crs

INFO:marshsi.prismaenmap.retrieval_upv_prisma_enmap:Running WMF in wavelength range: 978.68798828125 - 2462.7236328125
INFO:marshsi.prismaenmap.retrieval_upv_prisma_enmap:Running MF2100 in wavelength range: 2102.565673828125 - 2462.7236328125


 
         Transform: | 30.00, 0.00, 230898.88|
| 0.00,-30.00, 4257208.59|
| 0.00, 0.00, 1.00|
         Shape: (7, 1221, 1223)
         Resolution: (30.0, 30.0)
         Bounds: (230898.877003546, 4220578.5899621155, 267588.877003546, 4257208.5899621155)
         CRS: EPSG:32640
         fill_value_default: -1
        

rgb_local = rst.load_rgb(as_reflectance=True, raw=False)
rgb_local

 
         Transform: | 30.00, 0.00, 230898.88|
| 0.00,-30.00, 4257208.59|
| 0.00, 0.00, 1.00|
         Shape: (3, 1221, 1223)
         Resolution: (30.0, 30.0)
         Bounds: (230898.877003546, 4220578.5899621155, 267588.877003546, 4257208.5899621155)
         CRS: EPSG:32640
         fill_value_default: -1

Plot results

from georeader import plot
import matplotlib.pyplot as plt

fig, ax = plt.subplots(1,3, figsize=(15,5), sharey=True)
plot.show(rgb_local, ax=ax[0], mask=True, title="RGB")

mf = cube_crs.isel({"band": 1})
plot.show(mf, ax=ax[1], vmin=0, vmax=.3, title= "MF $\Delta X$CH$_4$ [ppm]", mask=True, add_colorbar_next_to=True)
wmf = cube_crs.isel({"band": 2})
plot.show(wmf, ax=ax[2], vmin=0, vmax=.3, title= "WWMF $\Delta X$CH$_4$ [ppm]", mask=True,add_colorbar_next_to=True)

WARNING:matplotlib.image:Clipping input data to the valid range for imshow with RGB data ([0..1] for floats or [0..255] for integers). Got range [-1.0..1.0].


<Axes: title={'center': 'WWMF $\\Delta X$CH$_4$ [ppm]'}>
No description has been provided for this image

The marshsi package and tutorials are released under a Creative Commons non-commercial Share-Alike licence. For using the codes in comercial pipelines written consent by the authors must be provided.

If you find this work useful please cite:

@Article{roger_2024,
    AUTHOR = {Roger, J. and Guanter, L. and Gorroño, J. and Irakulis-Loitxate, I.},
    TITLE = {Exploiting the entire near-infrared spectral range to improve the detection of methane plumes with high-resolution imaging spectrometers},
    JOURNAL = {Atmospheric Measurement Techniques},
    VOLUME = {17},
    YEAR = {2024},
    NUMBER = {4},
    PAGES = {1333--1346},
    URL = {https://amt.copernicus.org/articles/17/1333/2024/},
    DOI = {10.5194/amt-17-1333-2024}
}