From b82cda9dcf48b6faaf5e0d7a9cd3db2253179724 Mon Sep 17 00:00:00 2001 From: Simon Proud Date: Thu, 25 Apr 2024 08:53:09 +0200 Subject: [PATCH] Initial commit for H-SAF products in NC format --- satpy/etc/composites/hsaf.yaml | 8 +- satpy/etc/enhancements/generic.yaml | 13 +++ satpy/etc/enhancements/hsaf.yaml | 13 +++ satpy/etc/readers/hsaf_nc.yaml | 34 +++++++ satpy/readers/hsaf_nc.py | 152 ++++++++++++++++++++++++++++ 5 files changed, 219 insertions(+), 1 deletion(-) create mode 100644 satpy/etc/enhancements/hsaf.yaml create mode 100644 satpy/etc/readers/hsaf_nc.yaml create mode 100644 satpy/readers/hsaf_nc.py diff --git a/satpy/etc/composites/hsaf.yaml b/satpy/etc/composites/hsaf.yaml index ac1b850656..76f62d39ad 100644 --- a/satpy/etc/composites/hsaf.yaml +++ b/satpy/etc/composites/hsaf.yaml @@ -1,4 +1,4 @@ -sensor_name: hsaf +sensor_name: visir/hsaf composites: instantaneous_rainrate_3: @@ -21,3 +21,9 @@ composites: prerequisites: - name: h05B standard_name: accum_rainrate_5b + + instantaneous_rainfall_rate: + compositor: !!python/name:satpy.composites.GenericCompositor + prerequisites: + - name: h63 + standard_name: instantaneous_rainfall_rate \ No newline at end of file diff --git a/satpy/etc/enhancements/generic.yaml b/satpy/etc/enhancements/generic.yaml index b3ec45501c..cbc3ffb19f 100644 --- a/satpy/etc/enhancements/generic.yaml +++ b/satpy/etc/enhancements/generic.yaml @@ -1253,3 +1253,16 @@ enhancements: stretch: crude min_stretch: [0,0,0] max_stretch: [1,1,1] + + instantaneous_rainfall_rate: + standard_name: instantaneous_rainfall_rate + operations: + - name: colorize + method: !!python/name:satpy.enhancements.colorize + kwargs: + palettes: + - {colors: 'spectral', + reverse: true, + min_value: 0.0, + max_value: 65.0, + } \ No newline at end of file diff --git a/satpy/etc/enhancements/hsaf.yaml b/satpy/etc/enhancements/hsaf.yaml new file mode 100644 index 0000000000..bc4dc2713d --- /dev/null +++ b/satpy/etc/enhancements/hsaf.yaml @@ -0,0 +1,13 @@ +enhancements: + instantaneous_rainfall_rate: + standard_name: instantaneous_rainfall_rate + operations: + - name: colorize + method: !!python/name:satpy.enhancements.colorize + kwargs: + palettes: + - {colors: 'spectral', + reverse: true, + min_value: 0.0, + max_value: 65.0, + } \ No newline at end of file diff --git a/satpy/etc/readers/hsaf_nc.yaml b/satpy/etc/readers/hsaf_nc.yaml new file mode 100644 index 0000000000..fda3ea5e13 --- /dev/null +++ b/satpy/etc/readers/hsaf_nc.yaml @@ -0,0 +1,34 @@ +reader: + name: hsaf_nc + short_name: Hydrology SAF + long_name: Hydrology SAF products in NetCDF format + description: Reader for Hydrology SAF products + status: Beta, only h63 currently supported + supports_fsspec: false + reader: !!python/name:satpy.readers.yaml_reader.GEOFlippableFileYAMLReader + sensors: [hsaf] + +file_types: + hsafnc_63: + file_reader: !!python/name:satpy.readers.hsaf_nc.HSAFNCFileHandler + file_patterns: ['h63_{start_time:%Y%m%d_%H%M}_{region:s}.nc'] + hsafnc_90: + file_reader: !!python/name:satpy.readers.hsaf_nc.HSAFNCFileHandler + file_patterns: ['h90_{start_time:%Y%m%d_%H%M}_{acc_time:2s}_{region:s}.nc'] +datasets: + h63: + name: h63 + file_key: rr + sensor: hsaf + resolution: 3000 + standard_name: instantaneous_rainfall_rate + units: mm/h + file_type: hsafnc_63 + h90: + name: h90 + file_key: acc_rr + sensor: hsaf + resolution: 3000 + standard_name: accumulated_precipitation_estimate + units: mm + file_type: hsafnc_90 diff --git a/satpy/readers/hsaf_nc.py b/satpy/readers/hsaf_nc.py new file mode 100644 index 0000000000..f41b44172b --- /dev/null +++ b/satpy/readers/hsaf_nc.py @@ -0,0 +1,152 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +# Copyright (c) 2024. +# +# This file is part of satpy. +# +# satpy is free software: you can redistribute it and/or modify it under the +# 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 . +"""A reader for NetCDF Hydrology SAF products. In this beta version, only H63 is supported.""" +import logging +import os +from contextlib import suppress +from datetime import timedelta +from pyresample import geometry + +import dask.array as da +import numpy as np +import xarray as xr +from satpy.readers._geos_area import get_area_definition, get_area_extent + +from satpy.readers.file_handlers import BaseFileHandler +from satpy.readers.utils import unzip_file +from satpy.utils import get_chunk_size_limit + +logger = logging.getLogger(__name__) + +CHUNK_SIZE = get_chunk_size_limit() + +PLATFORM_NAMES = {"MSG1": "Meteosat-8", + "MSG2": "Meteosat-9", + "MSG3": "Meteosat-10", + "MSG4": "Meteosat-11", + "GOES16": "GOES-16", + "GOES17": "GOES-17", + } + +class HSAFNCFileHandler(BaseFileHandler): + """NWCSAF PPS&MSG NetCDF reader.""" + + def __init__(self, filename, filename_info, filetype_info): + """Init method.""" + super(HSAFNCFileHandler, self).__init__(filename, filename_info, + filetype_info) + + self._unzipped = unzip_file(self.filename) + if self._unzipped: + self.filename = self._unzipped + + self.cache = {} + self.nc = xr.open_dataset(self.filename, + decode_cf=True, + mask_and_scale=False, + chunks=CHUNK_SIZE) + + if 'xc' in self.nc.dims: + self.nc = self.nc.rename({"xc": "x", "yc": "y"}) + elif 'nx' in self.nc.dims: + self.nc = self.nc.rename({"nx": "x", "ny": "y"}) + + try: + kwrgs = {"sat_id": self.nc.attrs["satellite_identifier"]} + except KeyError: + kwrgs = {"sat_id": None} + + self.set_platform_and_sensor(**kwrgs) + + def __del__(self): + """Delete the instance.""" + if self._unzipped: + try: + os.remove(self._unzipped) + except OSError: + pass + + def set_platform_and_sensor(self, **kwargs): + """Set some metadata: platform_name, sensors, and pps (identifying PPS or Geo).""" + self.platform_name = PLATFORM_NAMES.get(kwargs["sat_id"], 'N/A') + + self.sensor = "seviri" + + def get_dataset(self, dsid, info): + """Load a dataset.""" + dsid_name = info["file_key"] + if dsid_name in self.cache: + logger.debug("Get the data set from cache: %s.", dsid_name) + return self.cache[dsid_name] + + logger.debug("Reading %s.", dsid_name) + variable = self.nc[dsid_name] + + # Data is transposed in file, fix it here + variable.data = variable.data.T + + variable.attrs["start_time"] = self.start_time + variable.attrs["end_time"] = self.end_time + + # Fill value is not defined as an attribute, manually specify + variable.data = da.where(variable.data > 0, variable.data, np.nan) + + return variable + + def get_area_def(self, dsid): + """Get the area definition of the band.""" + + val_dict = {} + for keyval in self.nc.attrs['cgms_projection'].split(): + try: + key, val = keyval.split('=') + val_dict[key] = val + except: + pass + + pdict = {} + pdict['scandir'] = 'N2S' + pdict["ssp_lon"] = np.float32(self.nc.attrs["sub-satellite_longitude"][:-1]) + pdict['a'] = float(val_dict['+r_eq'])*1000 + pdict['b'] = float(val_dict['+r_pol'])*1000 + pdict['h'] = float(val_dict['+h'])*1000 - pdict['a'] + pdict['loff'] = float(val_dict['+loff']) + pdict['coff'] = float(val_dict['+coff']) + pdict['lfac'] = -float(val_dict['+lfac']) + pdict['cfac'] = -float(val_dict['+cfac']) + pdict['ncols'] = self.nc.x.size + pdict['nlines'] = self.nc.y.size + pdict["a_name"] = "seviri_geos_fds" + pdict["a_desc"] = "SEVIRI full disk area at native resolution" + pdict["p_id"] = "seviri_fixed_grid" + + area_extent = get_area_extent(pdict) + fg_area_def = get_area_definition(pdict, area_extent) + return fg_area_def + + @property + def start_time(self): + """Return the start time of the object.""" + return self.filename_info["start_time"] + + @property + def end_time(self): + """Return the end time of the object. + + The product does not provide the end time, so the start time is used.""" + return self.filename_info["start_time"] \ No newline at end of file