Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add reader for CAMEL emissivity datasets. #3023

Open
wants to merge 6 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
147 changes: 147 additions & 0 deletions satpy/etc/readers/camel_l3_nc.yaml
Original file line number Diff line number Diff line change
@@ -0,0 +1,147 @@
reader:
name: camel_l3_nc
short_name: CAMEL L3
long_name: CAMEL emissivity level 3 data in netCDF4 format.
description: >
Reader for the CAMEL emissivity product, produced from various L2/L3
datasets on a monthly basis. More details
`here <https://lpdaac.usgs.gov/products/cam5k30emv002/>`_.
status: Nominal
supports_fsspec: false
sensors: [combined]
reader: !!python/name:satpy.readers.yaml_reader.FileYAMLReader

file_types:
camel_emis_file:
file_reader: !!python/name:satpy.readers.camel_l3_nc.CAMELL3NCFileHandler
file_patterns:
- 'CAM5K30EM_emis_{start_period:%Y%m}_V{version:3s}.nc'


datasets:
# QA products
aster_ndvi:
name: aster_ndvi
file_key: aster_ndvi
file_type: [ camel_emis_file ]
aster_qflag:
name: aster_qflag
file_key: aster_qflag
file_type: [ camel_emis_file ]
bfemis_qflag:
name: bfemis_qflag
file_key: bfemis_qflag
file_type: [ camel_emis_file ]
camel_qflag:
name: camel_qflag
file_key: camel_qflag
file_type: [ camel_emis_file ]
snow_fraction:
name: snow_fraction
file_key: snow_fraction
file_type: [ camel_emis_file ]

# Emissivity bands
camel_emis_b1:
name: camel_emis_b1
file_key: camel_emis
band_id: 0
file_type: [ camel_emis_file ]
wavelength: 3.6
resolution: 0.05

camel_emis_b2:
name: camel_emis_b2
file_key: camel_emis
band_id: 1
file_type: [ camel_emis_file ]
wavelength: 4.3
resolution: 0.05

camel_emis_b3:
name: camel_emis_b3
file_key: camel_emis
band_id: 2
file_type: [ camel_emis_file ]
wavelength: 5.0
resolution: 0.05

camel_emis_b4:
name: camel_emis_b4
file_key: camel_emis
band_id: 3
file_type: [ camel_emis_file ]
wavelength: 5.8
resolution: 0.05

camel_emis_b5:
name: camel_emis_b5
file_key: camel_emis
band_id: 4
file_type: [ camel_emis_file ]
wavelength: 7.6
resolution: 0.05

camel_emis_b6:
name: camel_emis_b6
file_key: camel_emis
band_id: 5
file_type: [ camel_emis_file ]
wavelength: 8.3
resolution: 0.05

camel_emis_b7:
name: camel_emis_b7
file_key: camel_emis
band_id: 6
file_type: [ camel_emis_file ]
wavelength: 8.6
resolution: 0.05

camel_emis_b8:
name: camel_emis_b8
file_key: camel_emis
band_id: 7
file_type: [ camel_emis_file ]
wavelength: 9.1
resolution: 0.05

camel_emis_b9:
name: camel_emis_b9
file_key: camel_emis
band_id: 8
file_type: [ camel_emis_file ]
wavelength: 10.6
resolution: 0.05

camel_emis_b10:
name: camel_emis_b10
file_key: camel_emis
band_id: 9
file_type: [ camel_emis_file ]
wavelength: 10.8
resolution: 0.05

camel_emis_b11:
name: camel_emis_b11
file_key: camel_emis
band_id: 10
file_type: [ camel_emis_file ]
wavelength: 11.3
resolution: 0.05

camel_emis_b12:
name: camel_emis_b12
file_key: camel_emis
band_id: 11
file_type: [ camel_emis_file ]
wavelength: 12.1
resolution: 0.05

camel_emis_b13:
name: camel_emis_b13
file_key: camel_emis
band_id: 12
file_type: [ camel_emis_file ]
wavelength: 14.3
resolution: 0.05
106 changes: 106 additions & 0 deletions satpy/readers/camel_l3_nc.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,106 @@
#!/usr/bin/env python
# -*- coding: utf-8 -*-
# Copyright (c) 2024 Satpy developers
#
# This file is part of satpy.
#
# satpy is free software: you can redistribute it and/or modify it under the
simonrp84 marked this conversation as resolved.
Show resolved Hide resolved
# terms of the GNU General Public License as published by the Free Software
# Foundation, either version 3 of the License, or (at your option) any later
# version.
#
# satpy is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License along with
# satpy. If not, see <http://www.gnu.org/licenses/>.

"""Reader for CAMEL Level 3 emissivity files in netCDF4 format.

For more information about the data, see: <https://lpdaac.usgs.gov/products/cam5k30emv002/>.

NOTE: This reader only supports the global 0.05 degree grid data.
"""


import datetime as dt
import logging

import xarray as xr
from pyresample import geometry

from satpy.readers.file_handlers import BaseFileHandler

logger = logging.getLogger(__name__)

# Area extent for the CAMEL product (global)
GLOB_AREA_EXT = [-180, -90, 180, 90]


class CAMELL3NCFileHandler(BaseFileHandler):
"""File handler for Himawari L2 NOAA enterprise data in netCDF format."""
simonrp84 marked this conversation as resolved.
Show resolved Hide resolved

def __init__(self, filename, filename_info, filetype_info):
"""Initialize the reader."""
super().__init__(filename, filename_info, filetype_info)
self.nc = xr.open_dataset(self.filename,
decode_cf=True,
mask_and_scale=True,
chunks={"xc": "auto", "yc": "auto"})

if self.nc.attrs["geospatial_lon_resolution"] != "0.05 degree grid ":
raise ValueError("Only 0.05 degree grid data is supported.")
if self.nc.attrs["geospatial_lat_resolution"] != "0.05 degree grid ":
raise ValueError("Only 0.05 degree grid data is supported.")
simonrp84 marked this conversation as resolved.
Show resolved Hide resolved

self.nlines = self.nc.sizes["latitude"]
self.ncols = self.nc.sizes["longitude"]
self.area = None


@property
def start_time(self):
"""Start timestamp of the dataset."""
date_str = self.nc.attrs["time_coverage_start"]
return dt.datetime.strptime(date_str, "%Y-%m-%d %H:%M:%SZ")

@property
def end_time(self):
"""End timestamp of the dataset."""
date_str = self.nc.attrs["time_coverage_end"]
return dt.datetime.strptime(date_str, "%Y-%m-%d %H:%M:%SZ")


def get_dataset(self, key, info):
"""Load a dataset."""
var = info["file_key"]
logger.debug("Reading in get_dataset %s.", var)
variable = self.nc[var]

# For the emissivity there are multiple bands, so we need to select the correct one
if var == "camel_emis":
if info["band_id"] >= variable.shape[2]:
raise ValueError("Band id requested is larger than dataset.")
variable = variable[:, :, info["band_id"]]

# Rename the latitude and longitude dimensions to x and y
variable = variable.rename({"latitude": "y", "longitude": "x"})

variable.attrs.update(key.to_dict())
return variable


def get_area_def(self, dsid):
"""Get the area definition, a global lat/lon area for this type of dataset."""
proj_param = "EPSG:4326"

area = geometry.AreaDefinition("gridded_camel",
"A global gridded area",
"longlat",
proj_param,
self.ncols,
self.nlines,
GLOB_AREA_EXT)
self.area = area
return area
133 changes: 133 additions & 0 deletions satpy/tests/reader_tests/test_camel_l3_nc.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,133 @@
"""Tests for the CAMEL L3 netCDF reader."""

import datetime as dt

import numpy as np
import pytest
import xarray as xr

from satpy.readers.camel_l3_nc import CAMELL3NCFileHandler
from satpy.tests.utils import make_dataid

rng = np.random.default_rng()
ndvi_data = rng.integers(0, 1000, (3600, 7200), dtype=np.int16)
emis_data = rng.integers(0, 1000, (3600, 7200, 5), dtype=np.int16)

lon_data = np.arange(-180, 180, 0.05)
lat_data = np.arange(-90, 90, 0.05)

start_time = dt.datetime(2023, 8, 1, 0, 0, 0)
end_time = dt.datetime(2023, 9, 1, 0, 0, 0)

fill_val = -999
scale_val = 0.001

dimensions = {"longitude": 7200, "latitude": 3600, "spectra": 13}

exp_ext = (-180.0, -90.0, 180.0, 90.0)

global_attrs = {"time_coverage_start": start_time.strftime("%Y-%m-%d %H:%M:%SZ"),
"time_coverage_end": end_time.strftime("%Y-%m-%d %H:%M:%SZ"),
"geospatial_lon_resolution": "0.05 degree grid ",
"geospatial_lat_resolution": "0.05 degree grid ",
}

bad_attrs1 = global_attrs.copy()
bad_attrs1["geospatial_lon_resolution"] = "0.1 degree grid "
bad_attrs2 = global_attrs.copy()
bad_attrs2["geospatial_lat_resolution"] = "0.1 degree grid "


def _make_ds(the_attrs, tmp_factory):
"""Make a dataset for use in tests."""
fname = f'{tmp_factory.mktemp("data")}/CAM5K30EM_emis_202308_V003.nc'
ds = xr.Dataset({"aster_ndvi": (["Rows", "Columns"], ndvi_data),
"camel_emis": (["latitude", "longitude", "spectra"], emis_data)},
coords={"latitude": (["Rows"], lat_data),
"longitude": (["Columns"], lon_data)},
attrs=the_attrs)
ds.to_netcdf(fname)
return fname


def camel_l3_filehandler(fname):
"""Instantiate a Filehandler."""
fileinfo = {"start_period": "202308",
"version": "003"}
filetype = None
fh = CAMELL3NCFileHandler(fname, fileinfo, filetype)
return fh


@pytest.fixture(scope="session")
def camel_filename(tmp_path_factory):
"""Create a fake camel l3 file."""
return _make_ds(global_attrs, tmp_path_factory)


@pytest.fixture(scope="session")
def camel_filename_bad1(tmp_path_factory):
"""Create a fake camel l3 file."""
return _make_ds(bad_attrs1, tmp_path_factory)


@pytest.fixture(scope="session")
def camel_filename_bad2(tmp_path_factory):
"""Create a fake camel l3 file."""
return _make_ds(bad_attrs2, tmp_path_factory)


def test_startend(camel_filename):
"""Test start and end times are set correctly."""
fh = camel_l3_filehandler(camel_filename)
assert fh.start_time == start_time
assert fh.end_time == end_time


def test_camel_l3_area_def(camel_filename, caplog):
"""Test reader handles area definition correctly."""
ps = "+proj=longlat +datum=WGS84 +no_defs +type=crs"

# Check case where input data is correct size.
fh = camel_l3_filehandler(camel_filename)
ndvi_id = make_dataid(name="aster_ndvi")
area_def = fh.get_area_def(ndvi_id)
assert area_def.width == dimensions["longitude"]
assert area_def.height == dimensions["latitude"]
assert np.allclose(area_def.area_extent, exp_ext)

assert area_def.proj4_string == ps


def test_bad_longitude(camel_filename_bad1):
"""Check case where longitude grid is not correct."""
with pytest.raises(ValueError, match="Only 0.05 degree grid data is supported."):
camel_l3_filehandler(camel_filename_bad1)


def test_bad_latitude(camel_filename_bad2):
"""Check case where latitude grid is not correct."""
with pytest.raises(ValueError, match="Only 0.05 degree grid data is supported."):
camel_l3_filehandler(camel_filename_bad2)


def test_load_ndvi_data(camel_filename):
"""Test that data is loaded successfully."""
fh = camel_l3_filehandler(camel_filename)
ndvi_id = make_dataid(name="aster_ndvi")
ndvi = fh.get_dataset(ndvi_id, {"file_key": "aster_ndvi"})
assert np.allclose(ndvi.data, ndvi_data)


def test_load_emis_data(camel_filename):
"""Test that data is loaded successfully."""
fh = camel_l3_filehandler(camel_filename)
emis_id = make_dataid(name="camel_emis")

# This is correct data
emis = fh.get_dataset(emis_id, {"file_key": "camel_emis", "band_id": 2})
assert np.allclose(emis.data, emis_data[:, :, 2])

# This will fail as we are requesting a band too high data
with pytest.raises(ValueError, match="Band id requested is larger than dataset."):
fh.get_dataset(emis_id, {"file_key": "camel_emis", "band_id": 12})
Loading