From c39f9d774517a18c62de72e3994816a5dc9d3d08 Mon Sep 17 00:00:00 2001 From: ninsbl Date: Tue, 6 Feb 2024 11:31:04 +0100 Subject: [PATCH 01/27] first draft --- .../i.sentinel/i.sentinel3.import/Makefile | 7 + .../i.sentinel3.import.html | 137 ++ .../i.sentinel3.import/i.sentinel3.import.py | 1817 +++++++++++++++++ 3 files changed, 1961 insertions(+) create mode 100644 src/imagery/i.sentinel/i.sentinel3.import/Makefile create mode 100644 src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html create mode 100644 src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py diff --git a/src/imagery/i.sentinel/i.sentinel3.import/Makefile b/src/imagery/i.sentinel/i.sentinel3.import/Makefile new file mode 100644 index 0000000000..beb92d7097 --- /dev/null +++ b/src/imagery/i.sentinel/i.sentinel3.import/Makefile @@ -0,0 +1,7 @@ +MODULE_TOPDIR = ../.. + +PGM = i.sentinel3.import + +include $(MODULE_TOPDIR)/include/Make/Script.make + +default: script diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html new file mode 100644 index 0000000000..4dd09781a5 --- /dev/null +++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html @@ -0,0 +1,137 @@ +

DESCRIPTION

+ +The i.sentinel3.import module allows importing Copernicus Sentinel-3 products +downloaded by the i.sentinel.download +module. + +

+By default i.sentinel3.import imports all Sentinel-3 scene files found +in the input directory. The number of scene files can be optionally +reduced by the pattern option or filtered by modification time using the +modified_after and/or modified_before option. In the pattern +option, a regular expression for filtering the file names should be given, e.g. +"0179_076_100_0900_LN2" for importing only specific tiles. +

+The Sentinel-3 products are provided in different formats. Not all of +them are currently directly supported by GDAL. Some, like the e.g. the LST products +consist of several netCDF4 files, where geometry and image information are stored in +different files. Therefore, a GRASS GIS specific import routine has been implemented +for those formats. This routine imports the data pixel-wise after transforming them +from WGS84 into the projection of the current LOCATION. r.in.xyz is then +used to read the coordinates into a GRASS GIS raster map. +Thus, the the user has to set the computational region extent and +resolution appropriately, to import the data of interest in the correct resolution. + +

+By default, ancillary data and quality flag data are imported as well - as artificial +bands. Import of those additional data layers can be deactivated using the ancillary_bands +and quality_bands option. + +

+For each imported band both scene and band specific metadata are written into the map history +(r.support). +In addition, the scene name is stored as source1 and the imported or linked file name as +source2. Also, sensing time is written into the timestamp of the +map. After import, the metadata can be retrieved with r.info -e +as shown below. + + +

NOTES

+ +

+Currently, only import of SL_2_LST__ products from Sentinel-3 SLSTR sensor are supported. +

+By means of the register_file option i.sentinel3.import allows +creating a file which can be used to register imported imagery data +into a space-time raster dateset (STRDS) with +t.register. +See example below. + +

+

Metadata storage

+ +By using the -j flag the band metadata are additionally stored +in JSON format (in the current mapset under cell_misc). + + +

EXAMPLES

+ +

List Sentinel bands

+ +At first, print the list of raster files to be imported with -p flag. For +each file also product content is printed: + +
+i.sentinel3.import -p input=./data/ nprocs=4 modified_before="2021-09-09" \
+product=LST basename=S3_LST
+
+...
+
+ +

Import Sentinel-3 data

+Import all Sentinel-3 data found in data directory and store metadata +as JSON files within the GRASS GIS database directory: +
+i.sentinel3.import -a -c -j input=./data/ nprocs=4 modified_before="2021-09-09" \
+product=LST basename=S3_LST
+
+ + +

Register imported Sentinel-3 data into STRDS

+ +
+i.sentinel3.import -a -c -j input=./data/ register_output=./data/reg.txt nprocs=4 \
+modified_before="2021-09-09" product=LST basename=S3_LST
+
+# register imported data into existing STRDS
+t.register input=Sentinel_3 file=data/reg.txt
+
+ +A register file typically contains two columns: imported raster map +name and timestamp separated by |. In the case of current +development version of GRASS GIS which supports band references concept +(see g.bands module for details) a +register file is extended by a third column containg band reference +information, see the examples below. + +
+# register file produced by stable GRASS GIS 7.8 version
+S3_imp_surface_temperature_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00
+S3_imp_surface_temperature_standard_error_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00
+S3_imp_confidence_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00
+...
+# register file produced by development GRASS GIS 7.9 version
+S3_imp_surface_temperature_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00|S3_surface_temperature
+S3_imp_surface_temperature_standard_error_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00|S3_surface_temperature_standard_error
+S3_imp_confidence_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00|S3_confidence
+
+ +

REQUIREMENTS

+The following Python libraries are required + + +

SEE ALSO

+ + +Overview of i.sentinel toolset + +

+ +i.sentinel.download, +r.in.xyz, +t.register + + +

AUTHOR

+ +Stefan Blumentrath, NINA, Norway + + diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py new file mode 100644 index 0000000000..7dfe330a0d --- /dev/null +++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py @@ -0,0 +1,1817 @@ +#!/usr/bin/env python3 + +""" + MODULE: i.sentinel3.import + AUTHOR(S): Stefan Blumentrath + PURPOSE: Import and pre-process Sentinel-3 data from the Copernicus program + COPYRIGHT: (C) 2024 by Norwegian Water and ENergy Directorate, Stefan Blumentrath, + and the GRASS development team + + This program is free software under the GNU General Public + License (>=v2). Read the file COPYING that comes with GRASS + for details. +""" + +# %Module +# % description: Import and pre-process Sentinel-3 data from the Copernicus program +# % keyword: imagery +# % keyword: satellite +# % keyword: Sentinel +# % keyword: import +# % keyword: optical +# % keyword: thermal +# %end + +# %option +# % key: input +# % label: Sentinel-3 input data +# % description: Either a (commaseparated list of) path(s) to Sentinel-3 zip files or a textfile with such paths (one per row) +# % required: yes +# % multiple: yes +# %end + +# %option +# % key: product_type +# % multiple: no +# % options: S3SL1RBT,S3SL2LST +# % answer: S3SL2LST +# % description: Sentinel-3 product type to import (currently, only S3SL1RBT and S3SL2LST are supported) +# % required: yes +# %end + +# %option +# % key: bands +# % multiple: yes +# % answer: all +# % required: yes +# % description: Data bands to import (e.g. LST, default is all available) +# %end + +# %option +# % key: anxillary_bands +# % multiple: yes +# % answer: all +# % required: no +# % description: Anxillary data bands to import (e.g. LST_uncertainty, default is all available) +# %end + +# %option +# % key: flag_bands +# % multiple: yes +# % answer: all +# % required: no +# % description: Quality flag bands to import (e.g. bayes_in, default is all available) +# %end + +# %option +# % key: basename +# % description: Basename used as prefix for map names (default is derived from the input file(s)) +# % required: no +# %end + +# %option G_OPT_F_OUTPUT +# % key: register_output +# % description: Name for output file to use with t.register +# % required: no +# %end + +# %option G_OPT_M_DIR +# % key: metadata +# % description: Name of directory into which Sentinel metadata json dumps are saved +# % required: no +# %end + +# %option G_OPT_M_NPROCS +# %end + +# %option +# % key: memory +# % type: integer +# % required: no +# % multiple: no +# % label: Maximum memory to be used (in MB) +# % description: Cache size for raster rows +# % answer: 300 +# %end + +# %flag +# % key: a +# % description: Apply cloud mask before import (can significantly speed up import) +# % guisection: Settings +# %end + +# %flag +# % key: c +# % description: Import LST in degree celsius (default is kelvin) +# % guisection: Settings +# %end + +# %flag +# % key: d +# % description: Import data with decimals as double precision +# % guisection: Settings +# %end + +# %flag +# % key: e +# % description: Import also elevation from geocoding of stripe +# % guisection: Settings +# %end + +# %flag +# % key: i +# % description: Do not interpolate imported bands +# % guisection: Settings +# %end + +# %flag +# % key: j +# % description: Write metadata json for each band to LOCATION/MAPSET/cell_misc/BAND/description.json +# % guisection: Settings +# %end + +# %flag +# % key: o +# % description: Process oblique view (default is nadir) +# % guisection: Settings +# %end + +# %flag +# % key: p +# % description: Print raster data to be imported and exit +# % guisection: Print +# %end + +# %flag +# % key: r +# % description: Rescale radiance bands to reflectance +# % guisection: Settings +# %end + +# %rules +# % excludes: -p,register_output +# %end + + +# ToDo +# - handle valid range according to: https://docs.unidata.ucar.edu/netcdf-c/current/attribute_conventions.html +# - add evtl. cleaning of external, temporary files +# - Adjust region bounds across tracks (pixel alignment may differ between tracks) +# - Add orphaned pixels +# - parallelize NC conversion + + +import atexit +import json +import os +import re +import sys + +from collections import OrderedDict +from datetime import datetime +from itertools import chain +# from multiprocessing import Pool +from pathlib import Path +from zipfile import ZipFile + +import grass.script as gs +from grass.pygrass.modules import Module, MultiModule, ParallelModuleQueue +from grass.temporal.datetime_math import ( + # adjust_datetime_to_granularity, + create_suffix_from_datetime, + # datetime_to_grass_datetime_string as grass_timestamp, +) +# from grass.temporal.temporal_granularity import check_granularity_string + +S3_SOLAR_FLUX = { + "S3SL1RBT": { + "band": "{}_solar_irradiance_{}", + "nc_file": "{}_quality_{}.nc", + # See: https://github.com/senbox-org/s3tbx/blob/master/s3tbx-rad2refl/src/main/java/org/esa/s3tbx/processor/rad2refl/Rad2ReflConstants.java + "defaults": { + "S1": 1837.39, + "S2": 1525.94, + "S3": 956.17, + "S4": 365.90, + "S5": 248.33, + "S6": 78.33, + }, + }, + "S3SL2LST": None, +} + +S3_SUPPORTED_PRODUCTS = ["S3SL1RBT", "S3SL2LST"] + +S3_FILE_PATTERN = { + "S3OL1ERF": None, + "S3SL1RBT": "S3*SL_1_RBT*.zip", + "S3SL2LST": "S3*SL_2_LST*.zip", +} + +S3_SUN_PARAMTERS = { + "S3OL1ERF": None, + "S3SL1RBT": { + "geocoding": { + "nc_file": "geodetic_tx.nc", + "bands": {"lat": "latitude_tx", "lon": "longitude_tx"}, + }, + "sun_bands": { + "nc_file": "geometry_t{}.nc", + "bands": { + "solar_azimuth": "solar_azimuth_t{}", + "solar_zenith": "solar_zenith_t{}", + }, + }, + }, + "S3SL2LST": None, +} + +S3_GRIDS = { + "S3OL1ERF": None, + "S3SL2LST": {"i": 1000.0}, + "S3SL1RBT": { + "a": 500.0, # Stripe A + "b": 500.0, # Stripe B + "c": 1000.0, # TDI grid + "f": 500.0, # F channel grid + "i": 1000.0, # 1km grid + }, +} + +# Default view is n (nadir), o = oblique +S3_VIEWS = { + "S3OL1ERF": None, + "S3SL1RBT": ["o", "n"], + "S3SL2LST": ["n"], +} + +S3_RADIANCE_ADJUSTMENT = { + "S3_PN_SLSTR_L1_08": { + "o": { + # Nadir + "S1": 0.97, + "S2": 0.98, + "S3": 0.98, + "S5": 1.11, + "S6": 1.13, + }, + "n": { + # Oblique + "S1": 0.94, + "S2": 0.95, + "S3": 0.95, + "S5": 1.04, + "S6": 1.07, + }, + } +} + +S3_BANDS = { + "bands": { + "S3OL1ERF": None, + "S3SL1RBT": { + "F1": { + "geometries": ("f",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "BT", + }, + "F2": { + "geometries": ("i",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "BT", + }, + "S1": { + "geometries": ("a",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "radiance", + "exception": "{}_{}_{}", + }, + "S2": { + "geometries": ("a",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "radiance", + "exception": "{}_{}_{}", + }, + "S3": { + "geometries": ("a",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "radiance", + "exception": "{}_{}_{}", + }, + "S4": { + "geometries": ("a", "b"), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "radiance", + "exception": "{}_{}_{}", + }, + "S5": { + "geometries": ("a", "b"), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "radiance", + "exception": "{}_{}_{}", + }, + "S6": { + "geometries": ("a", "b"), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "radiance", + "exception": "{}_{}_{}", + }, + "S7": { + "geometries": ("i",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "BT", + }, + "S8": { + "geometries": ("i",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "BT", + }, + "S9": { + "geometries": ("i",), + "nc_file": "{}_{}_{}.nc", + "full_name": "{}_{}_{}", + "types": "BT", + }, + }, + "S3SL2LST": { + "LST": { + "geometries": ("i",), + "nc_file": "LST_{}.nc", + "full_name": "{}", + "types": None, + "exception": "exception", + }, + "LST_uncertainty": { + "geometries": ("i",), + "nc_file": "LST_{}.nc", + "full_name": "{}", + "types": None, + "exception": "exception", + }, + }, + }, + "flags": { + "S3OL1ERF": None, + "S3SL1RBT": { + "bayes": { + "geometries": ("a", "b", "f", "i"), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "cloud": { + "geometries": ("a", "b", "f", "i"), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "confidence": { + "geometries": ("a", "b", "f", "i"), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "pointing": { + "geometries": ("a", "b", "f", "i"), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "probability_cloud_dual": { + "geometries": ("i",), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "probability_cloud_single": { + "geometries": ("i",), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + }, + "S3SL2LST": { + "bayes": { + "geometries": ("i",), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "cloud": { + "geometries": ("i",), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "confidence": { + "geometries": ("i",), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "probability_cloud_dual": { + "geometries": ("i",), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + "probability_cloud_single": { + "geometries": ("i",), + "nc_file": "flags_{}.nc", + "full_name": "{}_{}", + "types": None, + }, + }, + }, + "anxillary": { + "S3OL1ERF": None, + "S3SL1RBT": None, + "S3SL2LST": { + "biome": { + "geometries": ("i",), + "nc_file": "LST_ancillary_ds.nc", + "full_name": "{}", + "types": None, + }, + "fraction": { + "geometries": ("i",), + "nc_file": "LST_ancillary_ds.nc", + "full_name": "{}", + "types": None, + }, + "NDVI": { + "geometries": ("i",), + "nc_file": "LST_ancillary_ds.nc", + "full_name": "{}", + "types": None, + }, + "TCWV": { + "geometries": ("i",), + "nc_file": "LST_ancillary_ds.nc", + "full_name": "{}", + "types": None, + }, + }, + }, +} + +# GRASS map precision +DTYPE_TO_GRASS = { + "float32": "FCELL", + "uint16": "CELL", + "uint8": "CELL", + "float64": "DCELL", +} + +OVERWRITE = gs.overwrite() + +GISENV = gs.gisenv() +GISDBASE = GISENV["GISDBASE"] +TMP_FILE = Path(gs.tempfile(create=False)) +TMP_FILE.mkdir(exist_ok=True, parents=True) +TMP_NAME = gs.tempname(12) + +MAPSET = None +SUN_ZENITH_ANGLE = None + + +# From r.import + +# initialize global vars +TMP_REG_NAME = None +TMPLOC = None +SRCGISRC = None +GISDBASE = None +TMP_REG_NAME = None + + +def cleanup(): + """Remove all temporary data""" + # remove temporary maps + if TMP_FILE: + gs.try_remove(TMP_FILE) + # remove temp location + if TMPLOC: + gs.try_rmdir(os.path.join(GISDBASE, TMPLOC)) + if SRCGISRC: + gs.try_remove(SRCGISRC) + if ( + TMP_REG_NAME + and gs.find_file( + name=TMP_REG_NAME, element="vector", mapset=gs.gisenv()["MAPSET"] + )["fullname"] + ): + gs.run_command( + "g.remove", type="vector", name=TMP_REG_NAME, flags="f", quiet=True + ) + + +def create_tmp_location(grass_env=gs.gisenv()): + """Create a temporary location""" + global GISDBASE, TMPLOC, SRCGISRC + GISDBASE = grass_env["GISDBASE"] + TMPLOC = gs.append_node_pid("tmp_s3_import_location") + if not (Path(GISDBASE) / TMPLOC).exists(): + gs.run_command("g.proj", flags="c", location=TMPLOC, epsg=4326) + SRCGISRC, src_env = gs.create_environment(GISDBASE, TMPLOC, "PERMANENT") + return src_env + + +def np_as_scalar(var): + """Return a numpy object as scalar""" + if type(var).__module__ == np.__name__: + if var.size > 1: + return str(var) + return var.item() + return var + + +def write_metadata(json_dict, metadatajson): + """Write extended map metadata to JSON file""" + gs.verbose(_("Writing metadata to maps...")) + with open(metadatajson, "w", encoding="UTF8") as outfile: + json.dump(json_dict, outfile) + + +def write_register_file(filename, register_input): + """Write file to register resulting maps""" + gs.verbose(_("Writing register file <{}>...").format(filename)) + with open(filename, "w", encoding="UTF8") as register_file: + register_file.write("\n".join(chain(*register_input))) + + +def convert_units(np_column, from_u, to_u): + """Converts a numpy column from one unit to + another, convertibility needs to be checked beforehand""" + try: + from cf_units import Unit + except ImportError: + gs.fatal( + _( + "Could not import cf_units. Please install it with:\n" + "'pip install cf_units'!" + ) + ) + + try: + converted_col = Unit(from_u).convert(np_column, Unit(to_u)) + except ValueError: + gs.fatal( + _( + "Warning: Could not convert units from {from_u} to {to_u}.").format( + from_u=from_u, to_u=to_u + ) + ) + converted_col = np_column + + return converted_col + + +def extend_region(region_dict, additional_region_dict): + """Extend region bounds in region_dict to cover also the additional_region_dict""" + region_dict["n"] = max(region_dict["n"], additional_region_dict["n"]) + region_dict["e"] = max(region_dict["e"], additional_region_dict["e"]) + region_dict["s"] = min(region_dict["s"], additional_region_dict["s"]) + region_dict["w"] = min(region_dict["w"], additional_region_dict["w"]) + # region_dict["nsres"] = (region_dict["nsres"] + additional_region_dict["nsres"]) / 2.0 + # region_dict["ewres"] = (region_dict["ewres"] + additional_region_dict["ewres"]) / 2.0 + return region_dict + + +def parse_s3_file_name(file_name): + """Extract info from file name according to naming onvention: + https://sentinels.copernicus.eu/web/sentinel/user-guides/sentinel-3-slstr/naming-convention + Assumes that file name is checked to be a valid / supported Sentinel-3 file name + :param file_name: string representing the file name of a Senintel-3 scene + """ + try: + return { + "mission_id": file_name[0:3], + "instrument": file_name[4:6], + "level": file_name[7], + "product": file_name[9:12], + "start_time": datetime.strptime(file_name[16:31], "%Y%m%dT%H%M%S"), + "end_time": datetime.strptime(file_name[32:47], "%Y%m%dT%H%M%S"), + "ingestion_time": datetime.strptime(file_name[48:63], "%Y%m%dT%H%M%S"), + "duration": file_name[64:68], + "cycle": file_name[69:72], + "relative_orbit": file_name[73:76], + "frame": file_name[77:81], + } + except ValueError: + gs.fatal(_("{} is not a supported Sentinel-3 scene").format(str(file_name))) + + +def group_scenes(s3_files, granularity="1 day"): + """ + Group scenes by information from the file name: + 1. mission ID + 2. product type + 3. temporal granule + 4. duration + 5. cycle + 6. relative orbit + : param s3_filesv: list of pathlib.Path objects with Sentine1-3 files + """ + groups = {} + for sfile in s3_files: + s3_name_dict = parse_s3_file_name(sfile.name) + group_id = "_".join( + [ + s3_name_dict["mission_id"], + options["product_type"][2:], + create_suffix_from_datetime( + s3_name_dict["start_time"], granularity + ).replace("_", ""), + s3_name_dict["duration"], + s3_name_dict["cycle"], + s3_name_dict["relative_orbit"], + ] + ) + if group_id in groups: + groups[group_id]["files"].append(sfile) + if s3_name_dict["start_time"] < groups[group_id]["start_time"]: + groups[group_id]["start_time"] = s3_name_dict["start_time"] + if s3_name_dict["end_time"] > groups[group_id]["end_time"]: + groups[group_id]["end_time"] = s3_name_dict["end_time"] + groups[group_id]["frames"].append(s3_name_dict["frame"]) + else: + groups[group_id] = { + "files": [sfile], + "start_time": s3_name_dict["start_time"], + "end_time": s3_name_dict["end_time"], + "frames": [s3_name_dict["frame"]], + } + return groups + + +def extract_file_info(s3_files, basename=None): + """Extract information from file name according to naming conventions""" + result_dict = {} + for s3_file in s3_files: + file_info = parse_s3_file_name(s3_file.name) + if not result_dict: + result_dict = file_info + result_dict["duration"] = {file_info["duration"]} + result_dict["cycle"] = {file_info["cycle"]} + result_dict["relative_orbit"] = {file_info["relative_orbit"]} + result_dict["frame"] = {file_info["frame"]} + else: + if file_info["mission_id"] != result_dict["mission_id"]: + result_dict["mission_id"] = "S3_" + result_dict["start_time"] = min( + result_dict["start_time"], file_info["start_time"] + ) + result_dict["end_time"] = max( + result_dict["end_time"], file_info["end_time"] + ) + result_dict["duration"].add(file_info["duration"]) + result_dict["cycle"].add(file_info["cycle"]) + result_dict["relative_orbit"].add(file_info["relative_orbit"]) + result_dict["frame"].add(file_info["frame"]) + for key in ["duration", "cycle", "relative_orbit"]: + if len(result_dict[key]) > 1: + gs.warning( + _("Merging {key}s {values}").format( + key=key, values=", ".join(result_dict[key]) + ) + ) + if len(result_dict["mission_id"]) == "3_": + gs.warning( + _("Merging Seninel-3A and Seninel-3B data").format( + key=key, values=", ".join(result_dict[key]) + ) + ) + if not basename: + basename = "_".join( + [ + result_dict["mission_id"], + result_dict["instrument"], + result_dict["level"], + result_dict["product"], + result_dict["start_time"].strftime("%Y%m%d%H%M%S"), + result_dict["end_time"].strftime("%Y%m%d%H%M%S"), + list(result_dict["duration"])[0], + list(result_dict["cycle"])[0], + list(result_dict["relative_orbit"])[0], + ] + ) + return {basename: result_dict} + + +def adjust_region_env(reg, coords): + """Get region bounding box of intersecting area with + coordinates aligned to the current region""" + coords_min = np.min(coords, axis=0) + coords_max = np.max(coords, axis=0) + + diff = coords_max[0:2] / np.array([float(reg["ewres"]), float(reg["nsres"])]) + coords_max_aligned = diff.astype(int) * np.array( + [float(reg["ewres"]), float(reg["nsres"])] + ) + + diff = coords_min[0:2] / np.array([float(reg["ewres"]), float(reg["nsres"])]) + coords_min_aligned = diff.astype(int) * np.array( + [float(reg["ewres"]), float(reg["nsres"])] + ) + + new_bounds = {} + new_bounds["e"], new_bounds["n"] = tuple( + np.min( + np.vstack( + (coords_max_aligned, np.array([reg["e"], reg["n"]]).astype(float)) + ), + axis=0, + ).astype(str) + ) + new_bounds["w"], new_bounds["s"] = tuple( + np.max( + np.vstack( + (coords_min_aligned, np.array([reg["w"], reg["s"]]).astype(float)) + ), + axis=0, + ).astype(str) + ) + + if float(new_bounds["s"]) >= float(new_bounds["n"]) or float( + new_bounds["w"] + ) >= float(new_bounds["e"]): + return None + + compute_env = os.environ.copy() + compute_env["GRASS_region"] = gs.region_env(**new_bounds) + + return compute_env + + +def get_geocoding(zip_file, root, geo_bands_dict, region_bounds=None): + """Get ground control points from NetCDF file""" + member = str(root / geo_bands_dict["nc_file"]) + nc_file_path = zip_file.extract(member, path=TMP_FILE) + with Dataset(nc_file_path) as nc_file_open: + nc_bands = OrderedDict() + print(nc_file_open.resolution) + resolution = nc_file_open.resolution.strip(" ").split(" ")[1:3] + + for band_id, band in geo_bands_dict["bands"].items(): + if band not in nc_file_open.variables: + gs.fatal( + _( + "{s3_file} does not contain a container {container} with band {band}").format( + s3_file=str(Path(nc_file_path).parent.name), + container=geo_bands_dict["nc_file"], + band=", ".join(band), + ) + ) + + # Add band to dict + nc_bands[band_id] = nc_file_open[band][:] + + # Create initial mask + mask = nc_bands["lat"][:].mask + + if region_bounds: + # Mask to region + mask = np.ma.mask_or( + mask, + np.ma.masked_outside( + nc_bands["lon"], + float(region_bounds["ll_w"]), + float(region_bounds["ll_e"]), + ).mask, + ) + mask = np.ma.mask_or( + mask, + np.ma.masked_outside( + nc_bands["lat"], + float(region_bounds["ll_s"]), + float(region_bounds["ll_n"]), + ).mask, + ) + + if mask.all(): + gs.fatal(_("No valid pixels inside computational region")) + + return nc_bands, mask, resolution + + +def setup_import_multi_module( + tmp_ascii, + mapname, + zrange=None, + val_col=None, + data_type=None, + method="mean", + solar_flux=None, + rules=None, + # support_kwargs=None, + # env_=None, +): + """Setup GRASS GIS moduls for importing S3 bands""" + # Basic import module + modules = [ + Module( + "r.in.xyz", + input=str(tmp_ascii), + output=f"{TMP_NAME}_{mapname}", + method=method, + separator=",", + x=1, + y=2, + # Array contains a column z at position 3 (with all 0) + # after coordinate transformation + z=val_col, + flags="i", + type=data_type, + zrange=zrange, + percent=100, + run_=False, + quiet=gs.verbosity() <= 2, + overwrite=OVERWRITE, + # env_=env_, + ) + ] + # Interpolation for missing pixels + nbh_function = "nmedian" if method == "mean" else "nmode" + neighborhoods = ", ".join( + [ + f"{TMP_NAME}_{mapname}{nbh}" + for nbh in [ + "[0,1]", + "[1,1]", + "[1,0]", + "[1,-1]", + "[0,-1]", + "[-1,-1]", + "[-1,0]", + "[-1,1]", + ] + ] + ) + mc_interpolate_template = f"{{mapname}}=if(isnull({TMP_NAME}_{mapname}), {nbh_function}({neighborhoods}), {TMP_NAME}_{mapname})" + + # Add conversion from radiance to reflectance if requested and relevant + if solar_flux: + mc_interpolate = mc_interpolate_template.format( + mapname=f"{TMP_NAME}_{mapname}_rad" + ) + mapname_reflectance = mapname.replace("radiance", "reflectance") + mc_expression = f"eval({mc_interpolate})\n{mapname_reflectance}={TMP_NAME}_{mapname}_rad * {np.pi} / {solar_flux} / cos({SUN_ZENITH_ANGLE})" + mapname = mapname_reflectance + else: + mc_expression = mc_interpolate_template.format(mapname=mapname) + + # Create mapcalc module + modules.append( + Module( + "r.mapcalc", + expression=mc_expression, + quiet=gs.verbosity() <= 2, + overwrite=OVERWRITE, + run_=False, + # env_=env_, + ) + ) + # if support_kwargs: + # modules.extend( + # [ + # Module( + # "r.support", + # **support_kwargs, + # quiet=True, + # run_=False, + # ), + # Module( + # "r.timestamp", + # map=mapname, + # date=start_time, + # quiet=True, + # run_=False, + # ), + # ] + # ) + + # Add categories (for flag datasets) + if rules: + modules.append( + Module( + "r.category", + quiet=True, + map=mapname, + rules="-", + stdin_=rules, + separator=":", + run_=False, + ) + ) + return modules + + +def get_file_metadata(nc_dataset): + """Collect metadata from NetCDF file""" + metadata = { + attr: np_as_scalar(nc_dataset.getncattr(attr)) + for attr in [ + "title", + "creation_time", + "absolute_orbit_number", + "track_offset", + "start_offset", + "institution", + "references", + "resolution", + "source", + "contact", + "comment", + "history", + "processing_baseline", + "product_name", + ] + if attr in nc_dataset.ncattrs() + } + + metadata["start_time"] = parse_timestr(nc_dataset.start_time) + metadata["end_time"] = parse_timestr(nc_dataset.stop_time) + + return metadata + + +def get_band_metadata( + band, nc_variable, fmt, file_metadata=None, basename=None, module_flags=None +): + """Extract band metadata from NetCDF variable""" + metadata = file_metadata.copy() + # Collect metadata + band_attrs = nc_variable.ncattrs() + + # Define variable name + varname_short = band.full_name + datatype = str(nc_variable[:].dtype) + + # Define map name + mapname = f"{basename}_{varname_short}" + metadata["mapname"] = mapname + + band_title = nc_variable.long_name if "long_name" in band_attrs else band.full_name + + # Define unit + unit = nc_variable.units if "units" in band_attrs else None + unit = "degree_celsius" if band.band_id.startswith("LST") and module_flags["c"] else unit + metadata["unit"] = unit + + # Define datatype and import method + if datatype not in DTYPE_TO_GRASS: + gs.fatal( + _("Unsupported datatype {dt} in band {band}").format( + dt=datatype, band=band.full_name + ) + ) + if datatype in ["uint8", "uint16"]: + method = "max" # Unfortunately there is no "mode" in r.in.xyz + fmt += ",%i" + else: + method = "mean" + fmt += ",%.12f" + metadata["method"] = method + metadata["datatype"] = DTYPE_TO_GRASS[datatype] + + # Compile description + for time_reference in ["start_time", "end_time"]: + metadata[time_reference] = metadata[time_reference].isoformat() + description = json.dumps(metadata, separators=["\n", ": "]).lstrip("{").rstrip("}") + + # Get max and min values + if "valid_max" in band_attrs: + min_val = nc_variable.valid_min + max_val = nc_variable.valid_max + if "scale_factor" in band_attrs: + min_val = min_val * nc_variable.scale_factor + max_val = max_val * nc_variable.scale_factor + if "add_offset" in band_attrs: + min_val = min_val + nc_variable.add_offset + max_val = max_val + nc_variable.add_offset + description += f"\n\nvalid_min: {min_val}\nvalid_max: {max_val}" + else: + min_val = nc_variable[:].min() + max_val = nc_variable[:].max() + metadata["description"] = description + + # Define valid input range + if "_FillValue" in band_attrs: + fill_val = nc_variable._FillValue if "_FillValue" in band_attrs else None + if fill_val > 0: + zrange = [ + fill_val - 1 if not min_val else np.max([fill_val - 1, min_val]), + max_val, + ] + else: + zrange = [ + fill_val + 1 if not min_val else np.max([fill_val + 1, min_val]), + max_val, + ] + elif "flag_values" in band_attrs: + zrange = [ + min(nc_variable.flag_values), + max(nc_variable.flag_values), + ] + else: + zrange = None + metadata["zrange"] = zrange + metadata["support_kwargs"] = { + "map": mapname, + "title": f"{band_title} from {metadata['title']}", + "history": None, + "units": unit, + "source1": metadata["product_name"], + "source2": metadata["history"], + "description": description, + "semantic_label": f"S3_{varname_short}", + } + + return metadata, fmt + + +def transform_coordinates(coordinates): + """Tranforms a numy array with coordinates to + projection of the current location""" + + # Create source coordinate reference + s_srs = osr.SpatialReference() + s_srs.ImportFromEPSG(4326) + + # Create target coordinate reference + t_srs = osr.SpatialReference() + t_srs.ImportFromWkt(gs.read_command("g.proj", flags="fw")) + + # Initialize osr transformation + transform = osr.CoordinateTransformation(s_srs, t_srs) + + return ( + coordinates[:, [1, 0]] + if s_srs.IsSame(t_srs) + else np.array(transform.TransformPoints(coordinates))[:, [0, 1]] + ) + + +def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True): + """Write temporary XYZ ascii file""" + # Extract grid coordinates + lon = np.ma.masked_where(mask, nc_bands["lon"][:]).compressed() + lat = np.ma.masked_where(mask, nc_bands["lat"][:]).compressed() + + if project: + # Project coordinates + np_output = transform_coordinates( + np.dstack((lat, lon)).reshape(lat.shape[0], 2) + ) + else: + np_output = np.hstack((lon[:, None], lat[:, None])) + # Fetch, mask and stack requested bands + for band in nc_bands: + if band in ["lat", "lon"]: + continue + add_array = nc_bands[band][:] + if np.ma.is_masked(add_array): + add_array = add_array.filled() + np_output = np.hstack( + (np_output, np.ma.masked_where(mask, add_array).compressed()[:, None]) + ) + print( + band, + np.min(add_array), + np.max(np.ma.masked_where(mask, add_array).compressed()[:, None]), + ) + print( + band, + np.min(np.ma.masked_where(mask, add_array).compressed()[:, None]), + np.max(np.ma.masked_where(mask, add_array).compressed()[:, None]), + ) + # Write to temporary file + if Path(tmp_ascii).exists(): + with open(tmp_ascii, "ab") as tmp_ascii_file: + np.savetxt(tmp_ascii_file, np_output, delimiter=",", fmt=fmt) + else: + np.savetxt(tmp_ascii, np_output, delimiter=",", fmt=fmt) + + stripe_reg = { # + "n": np.ma.max(np_output[:, 1]), # + (nsres / 2.0), + "s": np.ma.min(np_output[:, 1]), # - (nsres / 2.0), + "e": np.ma.max(np_output[:, 0]), # + (ewres / 2.0), + "w": np.ma.min(np_output[:, 0]), # - (ewres / 2.0), + # "nsres": nsres, + # "ewres": ewres, + } + return stripe_reg + + +def import_s3(s3_file, kwargs, s3_product=None): + """Import Sentinel-3 netCDF4 data""" + # Unpack dictionary variables + rmap = list(kwargs["meta_dict"].keys())[0] + region_bounds = kwargs["reg_bounds"] + # current_reg = kwargs["current_reg"] + mod_flags = kwargs["mod_flags"] + json_standard_folder = kwargs["json_standard_folder"] + module_queue = {} + region_dicts = {} + register_strings = {} + + with ZipFile(s3_file) as zip_file: + members = zip_file.namelist() + root = members[0].rsplit(".SEN3", maxsplit=1)[0] + root = Path(f"{root}.SEN3") + # Check if solar flux raster band is required + if mod_flags["r"]: + tmp_ascii_sun = TMP_FILE / f"{rmap}_sun_parameters.txt" + module_queue, sun_region_dict = s3_product.get_sun_parameters( + zip_file, root, rmap, tmp_ascii_sun, region_bounds + ) + region_dicts["sun_parameters"] = sun_region_dict + for stripe, stripe_dict in s3_product.requested_stripe_content.items(): + # Could be parallelized!!! + # Get geocoding (lat/lon, elevation) and initial mask + tmp_ascii = TMP_FILE / f"{rmap}_{stripe}.txt" + nc_bands, mask, resolution = get_geocoding( + zip_file, root, stripe_dict["geocoding"], region_bounds=region_bounds + ) + fmt = ",".join(["%.12f"] * len(nc_bands)) + for nc_file, bands in stripe_dict["containers"].items(): + ( + nc_bands, + module_queue, + register_output, + fmt, + ) = s3_product.process_nc_file( + zip_file, + root, + (nc_file, bands), + nc_bands, + rmap, + module_queue, + mod_flags, + tmp_ascii, + fmt=fmt, + json_standard_folder=json_standard_folder, + ) + register_strings[nc_file] = register_output + region_dicts[stripe] = write_xyz( + tmp_ascii, nc_bands, mask, fmt=fmt, project=True + ) + region_dicts[stripe]["ewres"], region_dicts[stripe]["nsres"] = float( + resolution[0] + ), float(resolution[1]) + + return module_queue, register_strings, region_dicts + + +class S3Product: + """Class to provide information necessary to pre-process Senintel-3 data products + level 1 and 2""" + + def __init__( + self, product_type, view="n", bands=None, flag_bands=None, anxillary_bands=None + ): + self.product_type = product_type + self.available_views = S3_VIEWS[product_type] + self.available_bands = list(S3_BANDS["bands"][product_type].keys()) + self.available_flag_bands = list(S3_BANDS["flags"][product_type].keys()) + self.available_anxillary_bands = ( + list(S3_BANDS["anxillary"][product_type].keys()) + if S3_BANDS["anxillary"][product_type] + else None + ) + self.view = self._check_view(view) + self.bands = self._check_bands(bands) + self.flag_bands = self._check_bands(flag_bands, band_type="flag_bands") + self.anxillary_bands = self._check_bands( + anxillary_bands, band_type="anxillary_bands" + ) + self.grids = S3_GRIDS[product_type] + self.requested_stripe_content = self._collect_requested_stripe_content() + self.file_pattern = S3_FILE_PATTERN[product_type] + + def __str__(self): + return json.dumps(self.__repr__(), indent=2) + + def __repr__(self): + class_dict = self.__dict__.copy() + for band_type in ["bands", "flag_bands", "anxillary_bands"]: + bands_of_type = getattr(self, band_type) + print(bands_of_type) + class_dict[band_type] = ( + { + band: description.__dict__ if description else None + for band, description in bands_of_type.items() + } + if bands_of_type + else None + ) + return json.dumps(class_dict, indent=2) + + def _check_view(self, view): + if view not in self.available_views: + gs.warning(_( + "View {} not available for product type {}").format( + view, self.product_type + ) + ) + return view + + def _check_bands(self, bands, band_type="bands"): + band_types = { + "bands": self.available_bands, + "flag_bands": self.available_flag_bands, + "anxillary_bands": self.available_anxillary_bands, + } + if not band_types[band_type]: + return None + if not bands: + if not bands and band_type == "bands": + gs.warning(_("At least one band band needs to be given")) + return None + if bands == "all": + bands = band_types[band_type] + band_objects = {} + for band in bands: + if band not in band_types[band_type]: + gs.warning(_( + "Band {0} is not available in product_type {1}").format( + band, self.product_type + ) + ) + band_objects[band] = S3Band(self.product_type, band, use_b=False, view="n") + return band_objects + + def _collect_requested_stripe_content(self): + suffixes = {} + for band_type in ["bands", "flag_bands", "anxillary_bands"]: + selected_bands = getattr(self, band_type) + if not selected_bands: + continue + for band_id, band_object in selected_bands.items(): + if band_object.suffix in suffixes: + if "containers" not in suffixes[band_object.suffix]: + suffixes[band_object.suffix]["containers"] = {} + if ( + band_object.nc_file + in suffixes[band_object.suffix]["containers"] + ): + suffixes[band_object.suffix]["containers"][ + band_object.nc_file + ].append(band_id) + else: + suffixes[band_object.suffix]["containers"][ + band_object.nc_file + ] = [band_id] + else: + suffixes[band_object.suffix] = { + "containers": {band_object.nc_file: [band_id]}, + "geocoding": { + "nc_file": f"geodetic_{band_object.suffix}.nc", + "bands": { + "lat": f"latitude_{band_object.suffix}", + "lon": f"longitude_{band_object.suffix}", + "elevation": f"elevation_{band_object.suffix}", + }, + }, + } + return suffixes + + def process_nc_file( + self, + zip_file, + root, + container_dict_items, + nc_bands, + prefix, + module_queue, + module_flags, + tmp_ascii, + fmt=None, + json_standard_folder=None, + ): + """Extract requested bands as numpy arrays from NetCDF file and setup import modules""" + meta_information = {} + member = str(root / container_dict_items[0]) + nc_file_path = zip_file.extract(member, path=TMP_FILE) + with Dataset(nc_file_path) as nc_file_open: + file_metadata = get_file_metadata(nc_file_open) + + for band in container_dict_items[1]: + # Check for band + solar_flux = None + if band in self.bands: + band = self.bands[band] + if band.solar_flux: + solar_flux = band.solar_flux["default"] + elif band in self.flag_bands: + band = self.flag_bands[band] + elif band in self.anxillary_bands: + band = self.anxillary_bands[band] + + if band.full_name not in nc_file_open.variables: + gs.fatal( + _( + "{s3_file} does not contain a container {container} with band {band}").format( + s3_file=str(root), + container=band.nc_file, + band=", ".join(band.full_name), + ) + ) + + nc_variable = nc_file_open[band.full_name] + # Apply radiance adjustment + if band.radiance_adjustment and module_flags["r"]: + add_array = nc_variable[:] * band.radiance_adjustment + else: + add_array = nc_variable[:] + + # Collect band metadata + metadata = get_band_metadata( + band, + nc_variable, + fmt, + file_metadata=file_metadata, + basename=prefix, + module_flags=module_flags, + ) + fmt = metadata[1] + if band.radiance_adjustment and module_flags["r"]: + # Get solar flux for band + metadata[0]["solar_flux"] = band.get_solar_flux( + zip_file, root, + ) + + if band.exception: + if np.ma.is_masked(add_array): + add_array.mask = np.ma.mask_or( + add_array.mask, nc_file_open[band.exception][:] > 0 + ) + else: + add_array = np.ma.masked_array( + add_array, nc_file_open[band.exception][:] > 0 + ) + if np.ma.is_masked(add_array): + add_array = add_array.filled() + + # Rescale temperature variables if requested + if band.full_name.startswith("LST") and module_flags["c"]: + add_array = convert_units(add_array, "K", "degree_celsius") + + # Write metadata json if requested + band_attrs = nc_file_open[band.full_name].ncattrs() + if json_standard_folder: + write_metadata( + { + **{ + a: np_as_scalar(nc_file_open.getncattr(a)) + for a in nc_file_open.ncattrs() + }, + **{"variable": band.full_name}, + **{ + a: np_as_scalar( + nc_file_open[band.full_name].getncattr(a) + ) + for a in nc_file_open[band.full_name].ncattrs() + }, + }, + json_standard_folder.joinpath(metadata[0]["mapname"] + ".json"), + ) + + # Define categories for flag datasets + rules = None + if "flag_masks" in band_attrs: + rules = "\n".join( + [ + ":".join( + [ + str(nc_file_open[band.full_name].flag_masks[idx]), + label, + ] + ) + for idx, label in enumerate( + nc_file_open[band.full_name].flag_meanings.split(" ") + ) + ] + ) + # Setup import modules + module_queue[band.full_name] = setup_import_multi_module( + tmp_ascii, + metadata[0]["mapname"], + zrange=metadata[0]["zrange"], + val_col=len(nc_bands) + 1, + data_type=metadata[0]["datatype"], + method=metadata[0]["method"], + # support_kwargs=metadata[0]["support_kwargs"], + solar_flux=solar_flux, + rules=rules, + # env_=None, + ) + + # Add array with invalid data masked to ordered dict of nc_bands + nc_bands[band.full_name] = add_array + meta_information[metadata[0]["mapname"]] = { + "rules": rules, + "start_time": metadata[0]["start_time"], + "end_time": metadata[0]["end_time"], + "support_kwargs": metadata[0]["support_kwargs"], + } + + return nc_bands, module_queue, meta_information, fmt + + def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds): + """https://github.com/sertit/eoreader/blob/main/eoreader/products/optical/s3_slstr_product.py#L862""" + + # nc_file = sun_azimuth_dict["nc_file"] + # Get values + import_modules = {} + fmt = "%.12f,%.12f" + nc_bands, mask, resolution = get_geocoding( + zip_file, + root, + S3_SUN_PARAMTERS[self.product_type]["geocoding"], + region_bounds=region_bounds, + ) + member = str( + root + / S3_SUN_PARAMTERS[self.product_type]["sun_bands"]["nc_file"].format( + self.view + ) + ) + nc_file_path = zip_file.extract(member, path=TMP_FILE) + sun_parameter_nc = Dataset(nc_file_path) + # sun_metadata = get_file_metadata(sun_parameter_nc) + for band_id, band in S3_SUN_PARAMTERS[self.product_type]["sun_bands"][ + "bands" + ].items(): + fmt += ",%.12f" + sun_parameter_array = sun_parameter_nc[band.format(self.view)] + nc_bands[band_id] = sun_parameter_array[:] * np.pi / 180.0 + map_name = f"{prefix}_{band_id}" + if band_id == "solar_zenith": + global SUN_ZENITH_ANGLE + SUN_ZENITH_ANGLE = map_name + # support_kwargs = { + # "map": map_name, + # "title": f"{sun_parameter_array.long_name} from {root}", + # # "history": nc_file_open.history, + # "units": sun_parameter_array.units, + # "source1": sun_parameter_nc.product_name, + # "source2": None, + # # "description": description, + # "semantic_label": f"S3_{sun_parameter_array.standard_name}", + # } + import_modules[band_id] = setup_import_multi_module( + tmp_ascii, + map_name, + zrange=None, + val_col=len(nc_bands), + data_type="FCELL", + method="mean", + solar_flux=None, + ) + + # Write to temporary file + sun_region_dict = write_xyz(tmp_ascii, nc_bands, mask, fmt=fmt, project=True) + sun_region_dict["ewres"], sun_region_dict["nsres"] = float( + resolution[0] + ), float(resolution[1]) + + return import_modules, sun_region_dict + + +class S3Band: + """Class for properties and methods related to Sentinel-3 bands""" + def __init__( + self, + product_type, + band_id, + use_b=False, + view="n", + radiance_adjustment="S3_PN_SLSTR_L1_08", + ): + self.band_id = band_id + # self.product_type = product_type + # self.use_b = use_b + self.band_type = self._get_band_type(product_type) + # Need to call in this order + # self.view = view + geometry = self._get_geometry(product_type, use_b) + self.resolution = S3_GRIDS[product_type][geometry] + self.suffix = f"{geometry}{view}" + self.full_name = self._get_full_name(product_type) + self.exception = self._get_exception_band(product_type) + self.nc_file = self._get_nc_file(product_type) + self.radiance_adjustment = self._get_radiance_adjustment( + radiance_adjustment, view + ) + self.solar_flux = self._get_solar_flux_dict_for_band(product_type) + + def __str__(self): + return json.dumps(self.__repr__(), indent=2) + + def __repr__(self): + return json.dumps(self.__dict__, indent=2) + + def _get_band_type(self, product_type): + for band_type, s3_bands in S3_BANDS.items(): + if self.band_id in s3_bands[product_type]: + return band_type + gs.fatal(_("Oh no")) + return None + + def _get_geometry(self, product_type, use_b): + band_dict = S3_BANDS[self.band_type][product_type][self.band_id] + return ( + "b" + if use_b and "b" in band_dict["geometries"] + else band_dict["geometries"][0] + ) + + def _get_full_name(self, product_type): + band_dict = S3_BANDS[self.band_type][product_type][self.band_id] + types = band_dict["types"] + if self.band_type != "bands": + return band_dict["full_name"].format(self.band_id, self.suffix) + if types: + return band_dict["full_name"].format(self.band_id, types, self.suffix) + return band_dict["full_name"].format(self.band_id, self.suffix) + + def _get_nc_file(self, product_type): + band_dict = S3_BANDS[self.band_type][product_type][self.band_id] + types = band_dict["types"] + if self.band_type == "flags": + return band_dict["nc_file"].format(self.suffix) + if types: + return band_dict["nc_file"].format(self.band_id, types, self.suffix) + return band_dict["nc_file"].format(self.band_id, self.suffix) + + def _get_radiance_adjustment(self, radiance_adjustment, view): + if radiance_adjustment not in S3_RADIANCE_ADJUSTMENT: + print("Oh no") + if self.band_id in S3_RADIANCE_ADJUSTMENT[radiance_adjustment][view]: + return S3_RADIANCE_ADJUSTMENT[radiance_adjustment][view][self.band_id] + return None + + def _get_exception_band(self, product_type): + if "exception" in S3_BANDS[self.band_type][product_type][self.band_id]: + return S3_BANDS[self.band_type][product_type][self.band_id]["exception"].format( + self.band_id, "exception", self.suffix + ) + return None + + def _get_solar_flux_dict_for_band(self, product_type): + if self.band_id in S3_SOLAR_FLUX[product_type]["defaults"]: + return { + "default": S3_SOLAR_FLUX[product_type]["defaults"][self.band_id], + "band": S3_SOLAR_FLUX[product_type]["band"].format( + self.band_id, self.suffix + ), + "nc_file": S3_SOLAR_FLUX[product_type]["nc_file"].format( + self.band_id, self.suffix + ), + } + return None + + def get_solar_flux(self, zip_file, root_path): + """ + Get solar spectral flux in mW / (m^2 * sr * nm) for band + + Args: + band_object (S3Band): Optical Band with radiance + + Returns: + float: solar Flux + + """ + if not self.solar_flux: + return None + member = str(root_path / self.solar_flux["nc_file"]) + nc_file_path = zip_file.extract(member, path=TMP_FILE) + sf_dataset = Dataset(nc_file_path) + solar_flux = np.nanmean(sf_dataset[self.solar_flux["band"]]) + if np.isnan(solar_flux): + solar_flux = self.solar_flux["default"] + return float(solar_flux) + + def _get_band_geo_points(self, suffix): + """Get ground control point references from suffix""" + if self.suffix.startswith("t"): + return None + return { + "nc_file": f"geodetic_{suffix}.nc", + "lat": f"latitude_{suffix}", + "lon": f"longitude_{suffix}", + "elevation": f"elevation_{suffix}", + } + + def adjust_radiance(self, np_band_array): + """Get radiance adjustment for band object""" + if not self.radiance_adjustment: + return np_band_array * self.radiance_adjustment + return np_band_array + + + +def main(): + """Do the main work""" + pattern = re.compile(S3_FILE_PATTERN[options["product_type"]]) + + # check provided input + s3_files = options["input"].split(",") + if len(s3_files) == 1: + if not re.match(pattern, s3_files[0]): + try: + s3_files = ( + Path(s3_files[0]).read_text(encoding="UTF8").strip().split("\n") + ) + except ValueError: + gs.fatal( + _( + "Input <{}> is neither a supported Sentinel-3 scene nor a text file with scenes" + ).format(options["input"]) + ) + + if not s3_files: + gs.warning("No scenes found to process, please check input.") + sys.exit() + + s3_files = [Path(scene) for scene in s3_files] + for s3_scene in s3_files: + if not s3_scene.exists(): + gs.fatal(_("Input file <{sn}> not found").format(str(sn=s3_scene))) + if re.match(pattern, s3_scene.name): + gs.fatal( + _( + "Input <{sn}> is not a supported Sentinel-3 scene for the requested product type <{pt}>" + ).format(sn=s3_scene, pt=options["product_type"]) + ) + + if flags["f"]: + global DTYPE_TO_GRASS + DTYPE_TO_GRASS["float64"] = "FCELL" + + s3_product = S3Product( + options["product_type"], + view="o" if flags["o"] else "n", + bands=options["bands"], + flag_bands=options["flag_bands"], + anxillary_bands=options["anxillary_bands"], + ) + + meta_info_dict = extract_file_info(s3_files, basename=options["basename"]) + # Group scenes by mission, + # groups = group_scenes(s3_files, granularity="1 day") + if gs.parse_command("g.proj", flags="g")["proj"] == "ll": + gs.fatal(_("Running in lonlat location is not supported")) + + nprocs = int(options["nprocs"]) + + # Get region bounds + region_bounds = gs.parse_command("g.region", flags="ugb", quiet=True) + current_region = gs.parse_command("g.region", flags="ug") + # has_band_ref = float(gs.version()["version"][0:3]) >= 7.9 + + json_standard_folder = None + if flags["j"]: + env = gs.gisenv() + json_standard_folder = Path(env["GISDBASE"]).joinpath( + env["LOCATION_NAME"], env["MAPSET"], "cell_misc" + ) + if not json_standard_folder.exists(): + json_standard_folder.mkdir(parents=True, exist_ok=True) + + # Collect variables for import + import_dict = { + "reg_bounds": dict(region_bounds), + "current_reg": dict(current_region), + "mod_flags": flags, + "json_standard_folder": json_standard_folder, + "meta_dict": meta_info_dict, + } + + if flags["p"]: + print( + "|".join( + [ + "product_file_name", + "nc_file_name", + "nc_file_title", + "nc_file_start_time", + "nc_file_creation_time", + "band", + "band_shape", + "band_title", + "band_standard_name", + "band_long_name", + ] + ) + ) + print("\n".join(chain(*import_result))) + return 0 + + module_queues = [] + register_strings = {} + region_dicts = {} + region_list = [] + + for s3_file in s3_files: + module_list, regs, region_dict = import_s3(s3_file, import_dict, s3_product=s3_product) + module_queues.append(module_list) + register_strings.update(regs) + region_list.append(region_dict) + if not region_dicts: + region_dicts = region_dict.copy() + else: + for stripe_id, region in region_dict.items(): + if stripe_id not in region_dicts: + region_dicts[stripe_id] = region + else: + extend_region(region_dicts[stripe_id], region_dict[stripe_id]) + + stripe_envs = {} + for stripe_id, stripe_region in region_dicts.items(): + stripe_env = os.environ.copy() + stripe_env["GRASS_REGION"] = gs.region_env(**stripe_region) # , env=stripe_env) + stripe_envs[stripe_id] = stripe_env + + queue = ParallelModuleQueue(nprocs) + for solar_parameter in ["solar_azimuth", "solar_zenith"]: # mq.items(): + compute_env = stripe_envs["sun_parameters"] + module_list = [] + for listed_module in module_queues[0][solar_parameter]: + listed_module.env_ = compute_env + listed_module.verbose = True + module_list.append(listed_module) + queue.put(MultiModule(module_list)) + + queue = ParallelModuleQueue(nprocs) + for band_id, modules in module_queues[0].items(): + if band_id in ["solar_azimuth", "solar_zenith"]: + continue + compute_env = stripe_envs[band_id[-2:]] + module_list = [] + for listed_module in modules: + listed_module.env_ = compute_env + listed_module.verbose = True + module_list.append(listed_module) + queue.put(MultiModule(module_list)) + + if options["register_output"]: + # Write t.register file if requested + write_register_file(options["register_output"], import_result) + + return 0 + + +if __name__ == "__main__": + options, flags = gs.parser() + # Lazy imports + try: + from dateutil.parser import isoparse as parse_timestr + except ImportError: + gs.fatal( + _( + "Could not import dateutil. Please install it with:\n" + "'pip install dateutil'!" + ) + ) + try: + from osgeo import osr + except ImportError: + gs.fatal( + _("Could not import gdal. Please install it with:\n" "'pip install GDAL'!") + ) + + try: + from netCDF4 import Dataset + except ImportError: + gs.fatal( + _( + "Could not import netCDF4. Please install it with:\n" + "'pip install netcdf4'!" + ) + ) + + try: + import numpy as np + except ImportError: + gs.fatal( + _( + "Could not import numpy. Please install it with:\n" + "'pip install numpy'!" + ) + ) + + atexit.register(cleanup) + sys.exit(main()) From 9aaf2db295cb1f5acdfbba09af4920cde68277cb Mon Sep 17 00:00:00 2001 From: ninsbl Date: Thu, 8 Feb 2024 15:04:21 +0100 Subject: [PATCH 02/27] fixed some todos --- .../i.sentinel3.import/i.sentinel3.import.py | 389 ++++++++---------- 1 file changed, 175 insertions(+), 214 deletions(-) diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py index 7dfe330a0d..2fb05488e3 100644 --- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py +++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py @@ -154,8 +154,6 @@ # ToDo -# - handle valid range according to: https://docs.unidata.ucar.edu/netcdf-c/current/attribute_conventions.html -# - add evtl. cleaning of external, temporary files # - Adjust region bounds across tracks (pixel alignment may differ between tracks) # - Add orphaned pixels # - parallelize NC conversion @@ -170,18 +168,12 @@ from collections import OrderedDict from datetime import datetime from itertools import chain -# from multiprocessing import Pool from pathlib import Path from zipfile import ZipFile import grass.script as gs from grass.pygrass.modules import Module, MultiModule, ParallelModuleQueue -from grass.temporal.datetime_math import ( - # adjust_datetime_to_granularity, - create_suffix_from_datetime, - # datetime_to_grass_datetime_string as grass_timestamp, -) -# from grass.temporal.temporal_granularity import check_granularity_string + S3_SOLAR_FLUX = { "S3SL1RBT": { @@ -470,7 +462,7 @@ "float32": "FCELL", "uint16": "CELL", "uint8": "CELL", - "float64": "DCELL", + "float64": "FCELL", } OVERWRITE = gs.overwrite() @@ -500,6 +492,16 @@ def cleanup(): # remove temporary maps if TMP_FILE: gs.try_remove(TMP_FILE) + + # Remove external data if mapset uses r.external.out + external = gs.parse_key_val(gs.read_command("r.external.out", flags="p"), sep=": ") + if "directory" in external: + for map_file in Path(external["directory"]).glob( + f"{TMP_NAME}_*{external['extension']}" + ): + if map_file.is_file(): + map_file.unlink() + # remove temp location if TMPLOC: gs.try_rmdir(os.path.join(GISDBASE, TMPLOC)) @@ -515,16 +517,9 @@ def cleanup(): "g.remove", type="vector", name=TMP_REG_NAME, flags="f", quiet=True ) - -def create_tmp_location(grass_env=gs.gisenv()): - """Create a temporary location""" - global GISDBASE, TMPLOC, SRCGISRC - GISDBASE = grass_env["GISDBASE"] - TMPLOC = gs.append_node_pid("tmp_s3_import_location") - if not (Path(GISDBASE) / TMPLOC).exists(): - gs.run_command("g.proj", flags="c", location=TMPLOC, epsg=4326) - SRCGISRC, src_env = gs.create_environment(GISDBASE, TMPLOC, "PERMANENT") - return src_env + gs.run_command( + "g.remove", type="vector", pattern=f"{TMP_NAME}_*", flags="f", quiet=True + ) def np_as_scalar(var): @@ -536,6 +531,32 @@ def np_as_scalar(var): return var +def get_dtype_range(np_datatype): + """Get information to set the valid data range based on dtype + according to: + https://docs.unidata.ucar.edu/netcdf-c/current/attribute_conventions.html#valid_range + """ + if np_datatype == "bool": + return {"min": 0, "max": 1, "reolution": 1} + + dtype = np.dtype(np_datatype) + if not np.issubdtype(dtype, np.floating) and not np.issubdtype(dtype, np.integer): + gs.fatal(_("Unsupported data type {}").format(np_datatype)) + + # Integer dtypes + if "int" in np_datatype: + dt_info = np.iinfo(np_datatype) + return {"min": dt_info.min, "max": dt_info.max, "resolution": 1} + + # Float dtypes + dt_info = np.finfo(np_datatype) + return { + "min": dt_info.min, + "max": dt_info.max, + "resolution": 2.0 * dt_info.resolution, + } + + def write_metadata(json_dict, metadatajson): """Write extended map metadata to JSON file""" gs.verbose(_("Writing metadata to maps...")) @@ -567,11 +588,10 @@ def convert_units(np_column, from_u, to_u): converted_col = Unit(from_u).convert(np_column, Unit(to_u)) except ValueError: gs.fatal( - _( - "Warning: Could not convert units from {from_u} to {to_u}.").format( - from_u=from_u, to_u=to_u - ) + _("Warning: Could not convert units from {from_u} to {to_u}.").format( + from_u=from_u, to_u=to_u ) + ) converted_col = np_column return converted_col @@ -583,8 +603,35 @@ def extend_region(region_dict, additional_region_dict): region_dict["e"] = max(region_dict["e"], additional_region_dict["e"]) region_dict["s"] = min(region_dict["s"], additional_region_dict["s"]) region_dict["w"] = min(region_dict["w"], additional_region_dict["w"]) - # region_dict["nsres"] = (region_dict["nsres"] + additional_region_dict["nsres"]) / 2.0 - # region_dict["ewres"] = (region_dict["ewres"] + additional_region_dict["ewres"]) / 2.0 + return region_dict + + +def adjust_region(region_dict): + """ + Adjust region bounds for r.in.xyz (bounds + half resolution) + aligned to resolution in ew and ns direction + starting from the sw-corner + """ + nsres = int(region_dict["nsres"]) + ewres = int(region_dict["ewres"]) + ns_round_to = 10 ** (len(str(nsres)) - len(str(nsres).rstrip("0"))) + ew_round_to = 10 ** (len(str(ewres)) - len(str(ewres).rstrip("0"))) + region_dict["s"] = ( + np.floor((region_dict["s"] - (nsres / 2.0)) / ns_round_to) * ns_round_to + ) + region_dict["w"] = ( + np.floor((region_dict["w"] - (ewres / 2.0)) / ew_round_to) * ew_round_to + ) + region_dict["n"] = ( + region_dict["s"] + + np.ceil(((region_dict["n"] + (nsres / 2.0) - region_dict["s"]) / nsres)) + * nsres + ) + region_dict["e"] = ( + region_dict["w"] + + np.ceil(((region_dict["e"] + (ewres / 2.0) - region_dict["w"]) / ewres)) + * ewres + ) return region_dict @@ -612,49 +659,6 @@ def parse_s3_file_name(file_name): gs.fatal(_("{} is not a supported Sentinel-3 scene").format(str(file_name))) -def group_scenes(s3_files, granularity="1 day"): - """ - Group scenes by information from the file name: - 1. mission ID - 2. product type - 3. temporal granule - 4. duration - 5. cycle - 6. relative orbit - : param s3_filesv: list of pathlib.Path objects with Sentine1-3 files - """ - groups = {} - for sfile in s3_files: - s3_name_dict = parse_s3_file_name(sfile.name) - group_id = "_".join( - [ - s3_name_dict["mission_id"], - options["product_type"][2:], - create_suffix_from_datetime( - s3_name_dict["start_time"], granularity - ).replace("_", ""), - s3_name_dict["duration"], - s3_name_dict["cycle"], - s3_name_dict["relative_orbit"], - ] - ) - if group_id in groups: - groups[group_id]["files"].append(sfile) - if s3_name_dict["start_time"] < groups[group_id]["start_time"]: - groups[group_id]["start_time"] = s3_name_dict["start_time"] - if s3_name_dict["end_time"] > groups[group_id]["end_time"]: - groups[group_id]["end_time"] = s3_name_dict["end_time"] - groups[group_id]["frames"].append(s3_name_dict["frame"]) - else: - groups[group_id] = { - "files": [sfile], - "start_time": s3_name_dict["start_time"], - "end_time": s3_name_dict["end_time"], - "frames": [s3_name_dict["frame"]], - } - return groups - - def extract_file_info(s3_files, basename=None): """Extract information from file name according to naming conventions""" result_dict = {} @@ -709,70 +713,25 @@ def extract_file_info(s3_files, basename=None): return {basename: result_dict} -def adjust_region_env(reg, coords): - """Get region bounding box of intersecting area with - coordinates aligned to the current region""" - coords_min = np.min(coords, axis=0) - coords_max = np.max(coords, axis=0) - - diff = coords_max[0:2] / np.array([float(reg["ewres"]), float(reg["nsres"])]) - coords_max_aligned = diff.astype(int) * np.array( - [float(reg["ewres"]), float(reg["nsres"])] - ) - - diff = coords_min[0:2] / np.array([float(reg["ewres"]), float(reg["nsres"])]) - coords_min_aligned = diff.astype(int) * np.array( - [float(reg["ewres"]), float(reg["nsres"])] - ) - - new_bounds = {} - new_bounds["e"], new_bounds["n"] = tuple( - np.min( - np.vstack( - (coords_max_aligned, np.array([reg["e"], reg["n"]]).astype(float)) - ), - axis=0, - ).astype(str) - ) - new_bounds["w"], new_bounds["s"] = tuple( - np.max( - np.vstack( - (coords_min_aligned, np.array([reg["w"], reg["s"]]).astype(float)) - ), - axis=0, - ).astype(str) - ) - - if float(new_bounds["s"]) >= float(new_bounds["n"]) or float( - new_bounds["w"] - ) >= float(new_bounds["e"]): - return None - - compute_env = os.environ.copy() - compute_env["GRASS_region"] = gs.region_env(**new_bounds) - - return compute_env - - def get_geocoding(zip_file, root, geo_bands_dict, region_bounds=None): """Get ground control points from NetCDF file""" member = str(root / geo_bands_dict["nc_file"]) nc_file_path = zip_file.extract(member, path=TMP_FILE) with Dataset(nc_file_path) as nc_file_open: nc_bands = OrderedDict() - print(nc_file_open.resolution) resolution = nc_file_open.resolution.strip(" ").split(" ")[1:3] for band_id, band in geo_bands_dict["bands"].items(): if band not in nc_file_open.variables: gs.fatal( _( - "{s3_file} does not contain a container {container} with band {band}").format( - s3_file=str(Path(nc_file_path).parent.name), - container=geo_bands_dict["nc_file"], - band=", ".join(band), - ) + "{s3_file} does not contain a container {container} with band {band}" + ).format( + s3_file=str(Path(nc_file_path).parent.name), + container=geo_bands_dict["nc_file"], + band=", ".join(band), ) + ) # Add band to dict nc_bands[band_id] = nc_file_open[band][:] @@ -966,7 +925,11 @@ def get_band_metadata( # Define unit unit = nc_variable.units if "units" in band_attrs else None - unit = "degree_celsius" if band.band_id.startswith("LST") and module_flags["c"] else unit + unit = ( + "degree_celsius" + if band.band_id.startswith("LST") and module_flags["c"] + else unit + ) metadata["unit"] = unit # Define datatype and import method @@ -990,42 +953,40 @@ def get_band_metadata( metadata[time_reference] = metadata[time_reference].isoformat() description = json.dumps(metadata, separators=["\n", ": "]).lstrip("{").rstrip("}") - # Get max and min values - if "valid_max" in band_attrs: + # Handle offset, scale, and valid data range, see: + # https://docs.unidata.ucar.edu/netcdf-c/current/attribute_conventions.html + zrange = None + if "valid_range" in band_attrs: + min_val, max_val = nc_variable.valid_range + elif "valid_max" in band_attrs: min_val = nc_variable.valid_min max_val = nc_variable.valid_max - if "scale_factor" in band_attrs: - min_val = min_val * nc_variable.scale_factor - max_val = max_val * nc_variable.scale_factor - if "add_offset" in band_attrs: + elif "_FillValue" in band_attrs: + data_range = get_dtype_range(datatype) + fill_val = nc_variable._FillValue + if fill_val > 0: + min_val = data_range["min"] + max_val = fill_val - data_range["resolution"] + elif fill_val < 0: + min_val = fill_val + data_range["resolution"] + max_val = data_range["max"] + else: + min_val, max_val = None, None + else: + min_val, max_val = None, None + + if min_val and max_val: + if "add_offset" in band_attrs and min_val and max_val: min_val = min_val + nc_variable.add_offset max_val = max_val + nc_variable.add_offset + + if "scale_factor" in band_attrs and min_val and max_val: + min_val = min_val * nc_variable.scale_factor + max_val = max_val * nc_variable.scale_factor description += f"\n\nvalid_min: {min_val}\nvalid_max: {max_val}" - else: - min_val = nc_variable[:].min() - max_val = nc_variable[:].max() + zrange = [min_val, max_val] metadata["description"] = description - # Define valid input range - if "_FillValue" in band_attrs: - fill_val = nc_variable._FillValue if "_FillValue" in band_attrs else None - if fill_val > 0: - zrange = [ - fill_val - 1 if not min_val else np.max([fill_val - 1, min_val]), - max_val, - ] - else: - zrange = [ - fill_val + 1 if not min_val else np.max([fill_val + 1, min_val]), - max_val, - ] - elif "flag_values" in band_attrs: - zrange = [ - min(nc_variable.flag_values), - max(nc_variable.flag_values), - ] - else: - zrange = None metadata["zrange"] = zrange metadata["support_kwargs"] = { "map": mapname, @@ -1086,16 +1047,7 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True): np_output = np.hstack( (np_output, np.ma.masked_where(mask, add_array).compressed()[:, None]) ) - print( - band, - np.min(add_array), - np.max(np.ma.masked_where(mask, add_array).compressed()[:, None]), - ) - print( - band, - np.min(np.ma.masked_where(mask, add_array).compressed()[:, None]), - np.max(np.ma.masked_where(mask, add_array).compressed()[:, None]), - ) + # Write to temporary file if Path(tmp_ascii).exists(): with open(tmp_ascii, "ab") as tmp_ascii_file: @@ -1104,12 +1056,10 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True): np.savetxt(tmp_ascii, np_output, delimiter=",", fmt=fmt) stripe_reg = { # - "n": np.ma.max(np_output[:, 1]), # + (nsres / 2.0), - "s": np.ma.min(np_output[:, 1]), # - (nsres / 2.0), - "e": np.ma.max(np_output[:, 0]), # + (ewres / 2.0), - "w": np.ma.min(np_output[:, 0]), # - (ewres / 2.0), - # "nsres": nsres, - # "ewres": ewres, + "n": np.ma.max(np_output[:, 1]), + "s": np.ma.min(np_output[:, 1]), + "e": np.ma.max(np_output[:, 0]), + "w": np.ma.min(np_output[:, 0]), } return stripe_reg @@ -1220,8 +1170,8 @@ def __repr__(self): def _check_view(self, view): if view not in self.available_views: - gs.warning(_( - "View {} not available for product type {}").format( + gs.warning( + _("View {} not available for product type {}").format( view, self.product_type ) ) @@ -1244,8 +1194,8 @@ def _check_bands(self, bands, band_type="bands"): band_objects = {} for band in bands: if band not in band_types[band_type]: - gs.warning(_( - "Band {0} is not available in product_type {1}").format( + gs.warning( + _("Band {0} is not available in product_type {1}").format( band, self.product_type ) ) @@ -1322,12 +1272,13 @@ def process_nc_file( if band.full_name not in nc_file_open.variables: gs.fatal( _( - "{s3_file} does not contain a container {container} with band {band}").format( - s3_file=str(root), - container=band.nc_file, - band=", ".join(band.full_name), - ) + "{s3_file} does not contain a container {container} with band {band}" + ).format( + s3_file=str(root), + container=band.nc_file, + band=", ".join(band.full_name), ) + ) nc_variable = nc_file_open[band.full_name] # Apply radiance adjustment @@ -1349,7 +1300,8 @@ def process_nc_file( if band.radiance_adjustment and module_flags["r"]: # Get solar flux for band metadata[0]["solar_flux"] = band.get_solar_flux( - zip_file, root, + zip_file, + root, ) if band.exception: @@ -1461,16 +1413,6 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds): if band_id == "solar_zenith": global SUN_ZENITH_ANGLE SUN_ZENITH_ANGLE = map_name - # support_kwargs = { - # "map": map_name, - # "title": f"{sun_parameter_array.long_name} from {root}", - # # "history": nc_file_open.history, - # "units": sun_parameter_array.units, - # "source1": sun_parameter_nc.product_name, - # "source2": None, - # # "description": description, - # "semantic_label": f"S3_{sun_parameter_array.standard_name}", - # } import_modules[band_id] = setup_import_multi_module( tmp_ascii, map_name, @@ -1479,7 +1421,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds): data_type="FCELL", method="mean", solar_flux=None, - ) + ) # Write to temporary file sun_region_dict = write_xyz(tmp_ascii, nc_bands, mask, fmt=fmt, project=True) @@ -1492,6 +1434,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds): class S3Band: """Class for properties and methods related to Sentinel-3 bands""" + def __init__( self, product_type, @@ -1501,11 +1444,8 @@ def __init__( radiance_adjustment="S3_PN_SLSTR_L1_08", ): self.band_id = band_id - # self.product_type = product_type - # self.use_b = use_b self.band_type = self._get_band_type(product_type) # Need to call in this order - # self.view = view geometry = self._get_geometry(product_type, use_b) self.resolution = S3_GRIDS[product_type][geometry] self.suffix = f"{geometry}{view}" @@ -1565,9 +1505,9 @@ def _get_radiance_adjustment(self, radiance_adjustment, view): def _get_exception_band(self, product_type): if "exception" in S3_BANDS[self.band_type][product_type][self.band_id]: - return S3_BANDS[self.band_type][product_type][self.band_id]["exception"].format( - self.band_id, "exception", self.suffix - ) + return S3_BANDS[self.band_type][product_type][self.band_id][ + "exception" + ].format(self.band_id, "exception", self.suffix) return None def _get_solar_flux_dict_for_band(self, product_type): @@ -1622,15 +1562,16 @@ def adjust_radiance(self, np_band_array): return np_band_array - def main(): """Do the main work""" - pattern = re.compile(S3_FILE_PATTERN[options["product_type"]]) + pattern = re.compile( + ".*" + S3_FILE_PATTERN[options["product_type"]].replace("*", ".*") + ) # check provided input s3_files = options["input"].split(",") if len(s3_files) == 1: - if not re.match(pattern, s3_files[0]): + if re.match(pattern, s3_files[0]) is None: try: s3_files = ( Path(s3_files[0]).read_text(encoding="UTF8").strip().split("\n") @@ -1650,16 +1591,16 @@ def main(): for s3_scene in s3_files: if not s3_scene.exists(): gs.fatal(_("Input file <{sn}> not found").format(str(sn=s3_scene))) - if re.match(pattern, s3_scene.name): + if not re.match(pattern, s3_scene.name): gs.fatal( _( "Input <{sn}> is not a supported Sentinel-3 scene for the requested product type <{pt}>" ).format(sn=s3_scene, pt=options["product_type"]) ) - if flags["f"]: + if flags["d"]: global DTYPE_TO_GRASS - DTYPE_TO_GRASS["float64"] = "FCELL" + DTYPE_TO_GRASS["float64"] = "DCELL" s3_product = S3Product( options["product_type"], @@ -1673,7 +1614,7 @@ def main(): # Group scenes by mission, # groups = group_scenes(s3_files, granularity="1 day") if gs.parse_command("g.proj", flags="g")["proj"] == "ll": - gs.fatal(_("Running in lonlat location is not supported")) + gs.fatal(_("Running in lonlat location is currently not supported")) nprocs = int(options["nprocs"]) @@ -1718,6 +1659,7 @@ def main(): ) ) print("\n".join(chain(*import_result))) + cleanup() return 0 module_queues = [] @@ -1726,9 +1668,12 @@ def main(): region_list = [] for s3_file in s3_files: - module_list, regs, region_dict = import_s3(s3_file, import_dict, s3_product=s3_product) + gs.verbose(_("Preparing scene {} for import").format(s3_file.name)) + module_list, register_string, region_dict = import_s3( + s3_file, import_dict, s3_product=s3_product + ) module_queues.append(module_list) - register_strings.update(regs) + register_strings.update(register_string) region_list.append(region_dict) if not region_dicts: region_dicts = region_dict.copy() @@ -1742,19 +1687,30 @@ def main(): stripe_envs = {} for stripe_id, stripe_region in region_dicts.items(): stripe_env = os.environ.copy() - stripe_env["GRASS_REGION"] = gs.region_env(**stripe_region) # , env=stripe_env) + stripe_env["GRASS_REGION"] = gs.region_env( + **adjust_region(stripe_region) + ) # , env=stripe_env) stripe_envs[stripe_id] = stripe_env - queue = ParallelModuleQueue(nprocs) - for solar_parameter in ["solar_azimuth", "solar_zenith"]: # mq.items(): compute_env = stripe_envs["sun_parameters"] - module_list = [] - for listed_module in module_queues[0][solar_parameter]: - listed_module.env_ = compute_env - listed_module.verbose = True - module_list.append(listed_module) - queue.put(MultiModule(module_list)) + if flags["r"]: + gs.verbose(_("Importing solar parameter bands")) + queue = ParallelModuleQueue(nprocs) + for solar_parameter in ["solar_azimuth", "solar_zenith"]: + module_list = [] + for solar_module in module_queues[0][solar_parameter]: + solar_module.env_ = compute_env + solar_module.verbose = True + module_list.append(solar_module) + queue.put(MultiModule(module_list)) + queue.wait() + + gs.verbose( + _("Importing scenes {}").format( + "\n" + "\n".join([s3_file.name for s3_file in s3_files]) + ) + ) queue = ParallelModuleQueue(nprocs) for band_id, modules in module_queues[0].items(): if band_id in ["solar_azimuth", "solar_zenith"]: @@ -1765,12 +1721,17 @@ def main(): listed_module.env_ = compute_env listed_module.verbose = True module_list.append(listed_module) - queue.put(MultiModule(module_list)) + queue.put(MultiModule(module_list)) + queue.wait() if options["register_output"]: # Write t.register file if requested - write_register_file(options["register_output"], import_result) - + write_register_file( + options["register_output"], "\n".join(register_strings.values()) + ) + else: + print("\n".join(register_strings.values())) + cleanup() return 0 From 00b4f5b7fdead3b0c198f932eb1fafffd251e8a4 Mon Sep 17 00:00:00 2001 From: ninsbl Date: Fri, 9 Feb 2024 08:41:58 +0100 Subject: [PATCH 03/27] update manual --- .../i.sentinel3.import.html | 76 +++++++++---------- 1 file changed, 38 insertions(+), 38 deletions(-) diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html index 4dd09781a5..ca93f9d3d5 100644 --- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html +++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html @@ -1,31 +1,45 @@

DESCRIPTION

-The i.sentinel3.import module allows importing Copernicus Sentinel-3 products -downloaded by the i.sentinel.download -module. - -

-By default i.sentinel3.import imports all Sentinel-3 scene files found -in the input directory. The number of scene files can be optionally -reduced by the pattern option or filtered by modification time using the -modified_after and/or modified_before option. In the pattern -option, a regular expression for filtering the file names should be given, e.g. -"0179_076_100_0900_LN2" for importing only specific tiles. +The i.sentinel3.import module allows importing Copernicus Sentinel-3 products.

-The Sentinel-3 products are provided in different formats. Not all of -them are currently directly supported by GDAL. Some, like the e.g. the LST products +The Sentinel-3 products types are provided in different formats. Not all of +them are currently directly supported by GDAL. Some, like the e.g. the SLSTR products consist of several netCDF4 files, where geometry and image information are stored in different files. Therefore, a GRASS GIS specific import routine has been implemented for those formats. This routine imports the data pixel-wise after transforming them from WGS84 into the projection of the current LOCATION. r.in.xyz is then -used to read the coordinates into a GRASS GIS raster map. -Thus, the the user has to set the computational region extent and -resolution appropriately, to import the data of interest in the correct resolution. +used to read the coordinates into a GRASS GIS raster map. Thus, the the user has +to set the computational region extent, to import the data of interest. The resolution +is set from the resolution of the inut bands and pixels are alined to that resolution. +

+Currently, only import of the following product_types is implemented: +

    + +
  • S3SL1RBT: SL_1_RBT__ products from Sentinel-3 SLSTR
    • +
  • S3SL2LST: SL_2_LST__ products from Sentinel-3 SLSTR
    • +
+ + +

+All Sentinel-3 scene files provided in the input option are mosaiced +into one resulting granule with one map for each imported bands. The input +option accepts either one or many path(s) to Sentinel-3 archives (.zip) separated +by comma or a textfile with one Sentinel-3 file per row. The module is designed +to produce temporally non overlaping mosaics along the tracks of the satellite. +Input should therefor not mix data from different tracks or sensing periods. + +

+Users can provide a custom basename, to which band names are appended. +Otherwise, the basename is constructed from the input scenes following the +Sentinel-3 file naming convention. Here, tile id is stripped from the resulting +basename and start of sensing time is set to start of sensing time of the first +input scene and end of sensing time is set to the end of the last input scene.

-By default, ancillary data and quality flag data are imported as well - as artificial -bands. Import of those additional data layers can be deactivated using the ancillary_bands -and quality_bands option. +Anxillary data and quality flag data can be imported as well using the +anxillary_bands and flag_bands option. Here, the user has to specify +the product specific band names or all to imort all available artificial +bands of the given kind.

For each imported band both scene and band specific metadata are written into the map history @@ -38,9 +52,6 @@

DESCRIPTION

NOTES

-

-Currently, only import of SL_2_LST__ products from Sentinel-3 SLSTR sensor are supported. -

By means of the register_file option i.sentinel3.import allows creating a file which can be used to register imported imagery data into a space-time raster dateset (STRDS) with @@ -56,7 +67,7 @@

Metadata storage

EXAMPLES

-

List Sentinel bands

+

Import Sentinel-3 data

Import all Sentinel-3 data found in data directory and store metadata as JSON files within the GRASS GIS database directory: 
-i.sentinel3.import -a -c -j input=./data/ nprocs=4 modified_before="2021-09-09" \
-product=LST basename=S3_LST
+i.sentinel3.import -a -c -j input=./data/ nprocs=4 product_type=S3SL2LST basename=S3_LST
@@ -81,7 +92,7 @@

Register imported Sentinel-3 data into STRDS

 i.sentinel3.import -a -c -j input=./data/ register_output=./data/reg.txt nprocs=4 \
-modified_before="2021-09-09" product=LST basename=S3_LST
+product_type=S3SL2LST basename=S3_LST
 
 # register imported data into existing STRDS
 t.register input=Sentinel_3 file=data/reg.txt
@@ -95,12 +106,6 @@ 
Register imported Sentinel-3 data into STRDS

 information, see the examples below.
 
 
-# register file produced by stable GRASS GIS 7.8 version
-S3_imp_surface_temperature_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00
-S3_imp_surface_temperature_standard_error_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00
-S3_imp_confidence_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00
-...
-# register file produced by development GRASS GIS 7.9 version
 S3_imp_surface_temperature_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00|S3_surface_temperature
 S3_imp_surface_temperature_standard_error_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00|S3_surface_temperature_standard_error
 S3_imp_confidence_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00|S3_confidence
@@ -129,9 +134,4 @@ 
SEE ALSO

 
 
AUTHOR

 
-Stefan Blumentrath, NINA, Norway
-
-
+Stefan Blumentrath, Norway

From 39f75ea6c932203a149e568f1a2e1dcd7b0cd854 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Fri, 23 Feb 2024 14:26:18 +0100
Subject: [PATCH 04/27] test agains reference

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 552 +++++++++---------
 1 file changed, 292 insertions(+), 260 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 2fb05488e3..7aeb7da31a 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -84,21 +84,12 @@
 # %option G_OPT_M_NPROCS
 # %end
 
-# %option
-# % key: memory
-# % type: integer
-# % required: no
-# % multiple: no
-# % label: Maximum memory to be used (in MB)
-# % description: Cache size for raster rows
-# % answer: 300
-# %end
-
-# %flag
-# % key: a
-# % description: Apply cloud mask before import (can significantly speed up import)
-# % guisection: Settings
-# %end
+# to be implemented
+# # %flag
+# # % key: a
+# # % description: Apply cloud mask before import (can significantly speed up import)
+# # % guisection: Settings
+# # %end
 
 # %flag
 # % key: c
@@ -112,21 +103,28 @@
 # % guisection: Settings
 # %end
 
+# to be implemented
+# # %flag
+# # % key: e
+# # % description: Import also elevation from geocoding of stripe
+# # % guisection: Settings
+# # %end
+
 # %flag
-# % key: e
-# % description: Import also elevation from geocoding of stripe
+# % key: n
+# % description: Import data at native resolution of the bands (default is use current region)
 # % guisection: Settings
 # %end
 
 # %flag
-# % key: i
-# % description: Do not interpolate imported bands
+# % key: j
+# % description: Write metadata json for each band to LOCATION/MAPSET/cell_misc/BAND/description.json
 # % guisection: Settings
 # %end
 
 # %flag
-# % key: j
-# % description: Write metadata json for each band to LOCATION/MAPSET/cell_misc/BAND/description.json
+# % key: k
+# % description: Keep original cell values during interpolation (see: r.fill.stats)
 # % guisection: Settings
 # %end
 
@@ -136,11 +134,12 @@
 # % guisection: Settings
 # %end
 
-# %flag
-# % key: p
-# % description: Print raster data to be imported and exit
-# % guisection: Print
-# %end
+# to be implemented
+# # %flag
+# # % key: p
+# # % description: Print raster data to be imported and exit
+# # % guisection: Print
+# # %end
 
 # %flag
 # % key: r
@@ -148,13 +147,16 @@
 # % guisection: Settings
 # %end
 
-# %rules
-# % excludes: -p,register_output
-# %end
+# # %rules
+# # % excludes: -p,register_output
+# # %end
 
 
 # ToDo
-# - Adjust region bounds across tracks (pixel alignment may differ between tracks)
+# - implement printing
+# - implement cloud-masking at import
+# - implement elevation import
+# - do not write to temporary file (feed r.in.xyz from stdin)
 # - Add orphaned pixels
 # - parallelize NC conversion
 
@@ -167,13 +169,17 @@
 
 from collections import OrderedDict
 from datetime import datetime
-from itertools import chain
+# from itertools import chain
 from pathlib import Path
 from zipfile import ZipFile
 
 import grass.script as gs
 from grass.pygrass.modules import Module, MultiModule, ParallelModuleQueue
+from grass.temporal.datetime_math import (
+    datetime_to_grass_datetime_string as grass_timestamp,
+)
 
+S3_SUPPORTED_PRODUCTS = ["S3SL1RBT", "S3SL2LST"]
 
 S3_SOLAR_FLUX = {
     "S3SL1RBT": {
@@ -192,8 +198,6 @@
     "S3SL2LST": None,
 }
 
-S3_SUPPORTED_PRODUCTS = ["S3SL1RBT", "S3SL2LST"]
-
 S3_FILE_PATTERN = {
     "S3OL1ERF": None,
     "S3SL1RBT": "S3*SL_1_RBT*.zip",
@@ -239,7 +243,7 @@
 
 S3_RADIANCE_ADJUSTMENT = {
     "S3_PN_SLSTR_L1_08": {
-        "o": {
+        "n": {
             # Nadir
             "S1": 0.97,
             "S2": 0.98,
@@ -247,7 +251,7 @@
             "S5": 1.11,
             "S6": 1.13,
         },
-        "n": {
+        "o": {
             # Oblique
             "S1": 0.94,
             "S2": 0.95,
@@ -468,25 +472,13 @@
 OVERWRITE = gs.overwrite()
 
 GISENV = gs.gisenv()
-GISDBASE = GISENV["GISDBASE"]
 TMP_FILE = Path(gs.tempfile(create=False))
 TMP_FILE.mkdir(exist_ok=True, parents=True)
 TMP_NAME = gs.tempname(12)
 
-MAPSET = None
 SUN_ZENITH_ANGLE = None
 
 
-# From r.import
-
-# initialize global vars
-TMP_REG_NAME = None
-TMPLOC = None
-SRCGISRC = None
-GISDBASE = None
-TMP_REG_NAME = None
-
-
 def cleanup():
     """Remove all temporary data"""
     # remove temporary maps
@@ -502,23 +494,8 @@ def cleanup():
             if map_file.is_file():
                 map_file.unlink()
 
-    # remove temp location
-    if TMPLOC:
-        gs.try_rmdir(os.path.join(GISDBASE, TMPLOC))
-    if SRCGISRC:
-        gs.try_remove(SRCGISRC)
-    if (
-        TMP_REG_NAME
-        and gs.find_file(
-            name=TMP_REG_NAME, element="vector", mapset=gs.gisenv()["MAPSET"]
-        )["fullname"]
-    ):
-        gs.run_command(
-            "g.remove", type="vector", name=TMP_REG_NAME, flags="f", quiet=True
-        )
-
     gs.run_command(
-        "g.remove", type="vector", pattern=f"{TMP_NAME}_*", flags="f", quiet=True
+        "g.remove", type="raster", pattern=f"{TMP_NAME}_*", flags="f", quiet=True
     )
 
 
@@ -557,20 +534,37 @@ def get_dtype_range(np_datatype):
     }
 
 
+def consolidat_metadata_dicts(metadata_dicts):
+    """Consolidate a list of metadata dictionaries for all
+    input scenes into unified metadata for the resulting map"""
+    map_metadata_dicts = {}
+    for map_dict in metadata_dicts:
+        for map_name, meta_dict in map_dict.items():
+            if map_name not in map_metadata_dicts:
+                map_metadata_dicts[map_name] = {}
+            for meta_key, meta_value in meta_dict.items():
+                if meta_key not in map_metadata_dicts[map_name]:
+                    map_metadata_dicts[map_name][meta_key] = meta_value
+                else:
+                    if (
+                        isinstance(map_metadata_dicts[map_name][meta_key], list)
+                        and meta_value not in map_metadata_dicts[map_name][meta_key]
+                    ):
+                        map_metadata_dicts[map_name][meta_key].append(meta_value)
+                    elif meta_value != map_metadata_dicts[map_name][meta_key]:
+                        map_metadata_dicts[map_name][meta_key] = [
+                            map_metadata_dicts[map_name][meta_key],
+                            meta_value,
+                        ]
+    return map_metadata_dicts
+
+
 def write_metadata(json_dict, metadatajson):
     """Write extended map metadata to JSON file"""
-    gs.verbose(_("Writing metadata to maps..."))
     with open(metadatajson, "w", encoding="UTF8") as outfile:
         json.dump(json_dict, outfile)
 
 
-def write_register_file(filename, register_input):
-    """Write file to register resulting maps"""
-    gs.verbose(_("Writing register file <{}>...").format(filename))
-    with open(filename, "w", encoding="UTF8") as register_file:
-        register_file.write("\n".join(chain(*register_input)))
-
-
 def convert_units(np_column, from_u, to_u):
     """Converts a numpy column from one unit to
     another, convertibility needs to be checked beforehand"""
@@ -654,6 +648,10 @@ def parse_s3_file_name(file_name):
             "cycle": file_name[69:72],
             "relative_orbit": file_name[73:76],
             "frame": file_name[77:81],
+            "producer": file_name[82:85],
+            "product_class": file_name[86:87],
+            "timeliness": file_name[88:90],
+            "baseline_collection": file_name[91:94],
         }
     except ValueError:
         gs.fatal(_("{} is not a supported Sentinel-3 scene").format(str(file_name)))
@@ -670,6 +668,8 @@ def extract_file_info(s3_files, basename=None):
             result_dict["cycle"] = {file_info["cycle"]}
             result_dict["relative_orbit"] = {file_info["relative_orbit"]}
             result_dict["frame"] = {file_info["frame"]}
+            result_dict["start_time"] = file_info["start_time"]
+            result_dict["end_time"] = file_info["end_time"]
         else:
             if file_info["mission_id"] != result_dict["mission_id"]:
                 result_dict["mission_id"] = "S3_"
@@ -683,19 +683,15 @@ def extract_file_info(s3_files, basename=None):
             result_dict["cycle"].add(file_info["cycle"])
             result_dict["relative_orbit"].add(file_info["relative_orbit"])
             result_dict["frame"].add(file_info["frame"])
-    for key in ["duration", "cycle", "relative_orbit"]:
+    for key in ["cycle", "relative_orbit"]:
         if len(result_dict[key]) > 1:
             gs.warning(
                 _("Merging {key}s {values}").format(
                     key=key, values=", ".join(result_dict[key])
                 )
             )
-    if len(result_dict["mission_id"]) == "3_":
-        gs.warning(
-            _("Merging Seninel-3A and Seninel-3B data").format(
-                key=key, values=", ".join(result_dict[key])
-            )
-        )
+    if result_dict["mission_id"] == "3_":
+        gs.warning(_("Merging Seninel-3A and Seninel-3B data"))
     if not basename:
         basename = "_".join(
             [
@@ -710,7 +706,10 @@ def extract_file_info(s3_files, basename=None):
                 list(result_dict["relative_orbit"])[0],
             ]
         )
-    return {basename: result_dict}
+    return (
+        basename,
+        result_dict,
+    )
 
 
 def get_geocoding(zip_file, root, geo_bands_dict, region_bounds=None):
@@ -767,14 +766,14 @@ def get_geocoding(zip_file, root, geo_bands_dict, region_bounds=None):
 def setup_import_multi_module(
     tmp_ascii,
     mapname,
+    distance=None,
+    fill_flags=False,
     zrange=None,
     val_col=None,
     data_type=None,
     method="mean",
     solar_flux=None,
     rules=None,
-    # support_kwargs=None,
-    # env_=None,
 ):
     """Setup GRASS GIS moduls for importing S3 bands"""
     # Basic import module
@@ -795,70 +794,44 @@ def setup_import_multi_module(
             zrange=zrange,
             percent=100,
             run_=False,
-            quiet=gs.verbosity() <= 2,
+            verbose=False,
+            quiet=True,
             overwrite=OVERWRITE,
-            # env_=env_,
         )
     ]
-    # Interpolation for missing pixels
-    nbh_function = "nmedian" if method == "mean" else "nmode"
-    neighborhoods = ", ".join(
-        [
-            f"{TMP_NAME}_{mapname}{nbh}"
-            for nbh in [
-                "[0,1]",
-                "[1,1]",
-                "[1,0]",
-                "[1,-1]",
-                "[0,-1]",
-                "[-1,-1]",
-                "[-1,0]",
-                "[-1,1]",
-            ]
-        ]
-    )
-    mc_interpolate_template = f"{{mapname}}=if(isnull({TMP_NAME}_{mapname}), {nbh_function}({neighborhoods}), {TMP_NAME}_{mapname})"
 
+    # Interpolation of missing / empty pixels
+    interp_mod = Module(
+        "r.fill.stats",
+        flags=fill_flags,
+        input=f"{TMP_NAME}_{mapname}",
+        mode="wmean" if method == "mean" else "mode",
+        cells=3,
+        distance=distance,
+        power=2.0,
+        run_=False,
+        quiet=True,
+        overwrite=OVERWRITE,
+    )
     # Add conversion from radiance to reflectance if requested and relevant
     if solar_flux:
-        mc_interpolate = mc_interpolate_template.format(
-            mapname=f"{TMP_NAME}_{mapname}_rad"
-        )
+        interp_mod.outputs.output = f"{TMP_NAME}_{mapname}_rad"
+        modules.append(interp_mod)
         mapname_reflectance = mapname.replace("radiance", "reflectance")
-        mc_expression = f"eval({mc_interpolate})\n{mapname_reflectance}={TMP_NAME}_{mapname}_rad * {np.pi} / {solar_flux} / cos({SUN_ZENITH_ANGLE})"
+        # Create mapcalc module
+        modules.append(
+            Module(
+                "r.mapcalc",
+                expression=f"{mapname_reflectance}={TMP_NAME}_{mapname}_rad * ({np.pi} / {solar_flux} / cos({SUN_ZENITH_ANGLE}))",
+                quiet=True,
+                overwrite=OVERWRITE,
+                run_=False,
+            )
+        )
         mapname = mapname_reflectance
     else:
-        mc_expression = mc_interpolate_template.format(mapname=mapname)
-
-    # Create mapcalc module
-    modules.append(
-        Module(
-            "r.mapcalc",
-            expression=mc_expression,
-            quiet=gs.verbosity() <= 2,
-            overwrite=OVERWRITE,
-            run_=False,
-            # env_=env_,
-        )
-    )
-    # if support_kwargs:
-    #     modules.extend(
-    #         [
-    #             Module(
-    #                 "r.support",
-    #                 **support_kwargs,
-    #                 quiet=True,
-    #                 run_=False,
-    #             ),
-    #             Module(
-    #                 "r.timestamp",
-    #                 map=mapname,
-    #                 date=start_time,
-    #                 quiet=True,
-    #                 run_=False,
-    #             ),
-    #         ]
-    #     )
+        interp_mod.outputs.output = mapname
+        modules.append(interp_mod)
 
     # Add categories (for flag datasets)
     if rules:
@@ -901,6 +874,7 @@ def get_file_metadata(nc_dataset):
 
     metadata["start_time"] = parse_timestr(nc_dataset.start_time)
     metadata["end_time"] = parse_timestr(nc_dataset.stop_time)
+    metadata["history"] = metadata["history"].strip()
 
     return metadata
 
@@ -921,7 +895,7 @@ def get_band_metadata(
     mapname = f"{basename}_{varname_short}"
     metadata["mapname"] = mapname
 
-    band_title = nc_variable.long_name if "long_name" in band_attrs else band.full_name
+    # band_title = nc_variable.long_name if "long_name" in band_attrs else band.full_name
 
     # Define unit
     unit = nc_variable.units if "units" in band_attrs else None
@@ -940,22 +914,23 @@ def get_band_metadata(
             )
         )
     if datatype in ["uint8", "uint16"]:
-        method = "max"  # Unfortunately there is no "mode" in r.in.xyz
+        metadata["method"] = "max"  # Unfortunately there is no "mode" in r.in.xyz
         fmt += ",%i"
     else:
-        method = "mean"
+        metadata["method"] = "mean"
         fmt += ",%.12f"
-    metadata["method"] = method
     metadata["datatype"] = DTYPE_TO_GRASS[datatype]
 
     # Compile description
     for time_reference in ["start_time", "end_time"]:
         metadata[time_reference] = metadata[time_reference].isoformat()
-    description = json.dumps(metadata, separators=["\n", ": "]).lstrip("{").rstrip("}")
+    metadata["description"] = (
+        json.dumps(metadata, separators=["\n", ": "]).lstrip("{").rstrip("}")
+    )
 
     # Handle offset, scale, and valid data range, see:
     # https://docs.unidata.ucar.edu/netcdf-c/current/attribute_conventions.html
-    zrange = None
+    metadata["zrange"] = None
     if "valid_range" in band_attrs:
         min_val, max_val = nc_variable.valid_range
     elif "valid_max" in band_attrs:
@@ -983,21 +958,10 @@ def get_band_metadata(
         if "scale_factor" in band_attrs and min_val and max_val:
             min_val = min_val * nc_variable.scale_factor
             max_val = max_val * nc_variable.scale_factor
-        description += f"\n\nvalid_min: {min_val}\nvalid_max: {max_val}"
-        zrange = [min_val, max_val]
-    metadata["description"] = description
-
-    metadata["zrange"] = zrange
-    metadata["support_kwargs"] = {
-        "map": mapname,
-        "title": f"{band_title} from {metadata['title']}",
-        "history": None,
-        "units": unit,
-        "source1": metadata["product_name"],
-        "source2": metadata["history"],
-        "description": description,
-        "semantic_label": f"S3_{varname_short}",
-    }
+        metadata["description"] += f"\n\nvalid_min: {min_val}\nvalid_max: {max_val}"
+        metadata["zrange"] = [min_val, max_val]
+
+    metadata["semantic_label"] = f"S3_{varname_short}"
 
     return metadata, fmt
 
@@ -1067,14 +1031,12 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
 def import_s3(s3_file, kwargs, s3_product=None):
     """Import Sentinel-3 netCDF4 data"""
     # Unpack dictionary variables
-    rmap = list(kwargs["meta_dict"].keys())[0]
+    rmap = kwargs["meta_dict"][0]
     region_bounds = kwargs["reg_bounds"]
-    # current_reg = kwargs["current_reg"]
     mod_flags = kwargs["mod_flags"]
-    json_standard_folder = kwargs["json_standard_folder"]
     module_queue = {}
     region_dicts = {}
-    register_strings = {}
+    register_strings = []
 
     with ZipFile(s3_file) as zip_file:
         members = zip_file.namelist()
@@ -1111,9 +1073,8 @@ def import_s3(s3_file, kwargs, s3_product=None):
                     mod_flags,
                     tmp_ascii,
                     fmt=fmt,
-                    json_standard_folder=json_standard_folder,
                 )
-                register_strings[nc_file] = register_output
+                register_strings.append(register_output)
             region_dicts[stripe] = write_xyz(
                 tmp_ascii, nc_bands, mask, fmt=fmt, project=True
             )
@@ -1157,7 +1118,6 @@ def __repr__(self):
         class_dict = self.__dict__.copy()
         for band_type in ["bands", "flag_bands", "anxillary_bands"]:
             bands_of_type = getattr(self, band_type)
-            print(bands_of_type)
             class_dict[band_type] = (
                 {
                     band: description.__dict__ if description else None
@@ -1248,7 +1208,6 @@ def process_nc_file(
         module_flags,
         tmp_ascii,
         fmt=None,
-        json_standard_folder=None,
     ):
         """Extract requested bands as numpy arrays from NetCDF file and setup import modules"""
         meta_information = {}
@@ -1262,8 +1221,11 @@ def process_nc_file(
                 solar_flux = None
                 if band in self.bands:
                     band = self.bands[band]
-                    if band.solar_flux:
-                        solar_flux = band.solar_flux["default"]
+                    if band.solar_flux and module_flags["r"]:
+                        solar_flux = band.get_solar_flux(
+                            zip_file,
+                            root,
+                        )
                 elif band in self.flag_bands:
                     band = self.flag_bands[band]
                 elif band in self.anxillary_bands:
@@ -1297,12 +1259,9 @@ def process_nc_file(
                     module_flags=module_flags,
                 )
                 fmt = metadata[1]
-                if band.radiance_adjustment and module_flags["r"]:
-                    # Get solar flux for band
-                    metadata[0]["solar_flux"] = band.get_solar_flux(
-                        zip_file,
-                        root,
-                    )
+
+                # Get solar flux for band
+                metadata[0]["solar_flux"] = solar_flux
 
                 if band.exception:
                     if np.ma.is_masked(add_array):
@@ -1322,23 +1281,6 @@ def process_nc_file(
 
                 # Write metadata json if requested
                 band_attrs = nc_file_open[band.full_name].ncattrs()
-                if json_standard_folder:
-                    write_metadata(
-                        {
-                            **{
-                                a: np_as_scalar(nc_file_open.getncattr(a))
-                                for a in nc_file_open.ncattrs()
-                            },
-                            **{"variable": band.full_name},
-                            **{
-                                a: np_as_scalar(
-                                    nc_file_open[band.full_name].getncattr(a)
-                                )
-                                for a in nc_file_open[band.full_name].ncattrs()
-                            },
-                        },
-                        json_standard_folder.joinpath(metadata[0]["mapname"] + ".json"),
-                    )
 
                 # Define categories for flag datasets
                 rules = None
@@ -1357,26 +1299,39 @@ def process_nc_file(
                         ]
                     )
                 # Setup import modules
+                fill_flags = "m"
+                if not module_flags["d"] and metadata[0]["datatype"] != "CELL":
+                    fill_flags += "s"
+                if module_flags["k"]:
+                    fill_flags += "k"
                 module_queue[band.full_name] = setup_import_multi_module(
                     tmp_ascii,
                     metadata[0]["mapname"],
+                    distance=2.0 * band.resolution,
+                    fill_flags=fill_flags,
                     zrange=metadata[0]["zrange"],
                     val_col=len(nc_bands) + 1,
                     data_type=metadata[0]["datatype"],
                     method=metadata[0]["method"],
-                    # support_kwargs=metadata[0]["support_kwargs"],
                     solar_flux=solar_flux,
                     rules=rules,
-                    # env_=None,
                 )
 
                 # Add array with invalid data masked to ordered dict of nc_bands
                 nc_bands[band.full_name] = add_array
+
                 meta_information[metadata[0]["mapname"]] = {
-                    "rules": rules,
-                    "start_time": metadata[0]["start_time"],
-                    "end_time": metadata[0]["end_time"],
-                    "support_kwargs": metadata[0]["support_kwargs"],
+                    "semantic_label": metadata[0]["semantic_label"],
+                    "unit": metadata[0]["unit"],
+                    **{
+                        a: np_as_scalar(nc_file_open.getncattr(a))
+                        for a in nc_file_open.ncattrs()
+                    },
+                    **{"variable": band.full_name},
+                    **{
+                        a: np_as_scalar(nc_file_open[band.full_name].getncattr(a))
+                        for a in nc_file_open[band.full_name].ncattrs()
+                    },
                 }
 
             return nc_bands, module_queue, meta_information, fmt
@@ -1408,7 +1363,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds):
         ].items():
             fmt += ",%.12f"
             sun_parameter_array = sun_parameter_nc[band.format(self.view)]
-            nc_bands[band_id] = sun_parameter_array[:] * np.pi / 180.0
+            nc_bands[band_id] = sun_parameter_array[:]
             map_name = f"{prefix}_{band_id}"
             if band_id == "solar_zenith":
                 global SUN_ZENITH_ANGLE
@@ -1416,6 +1371,8 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds):
             import_modules[band_id] = setup_import_multi_module(
                 tmp_ascii,
                 map_name,
+                distance=3,
+                fill_flags="ks",
                 zrange=None,
                 val_col=len(nc_bands),
                 data_type="FCELL",
@@ -1498,7 +1455,7 @@ def _get_nc_file(self, product_type):
 
     def _get_radiance_adjustment(self, radiance_adjustment, view):
         if radiance_adjustment not in S3_RADIANCE_ADJUSTMENT:
-            print("Oh no")
+            gs.fatal("Missing information on radiance adjustment")
         if self.band_id in S3_RADIANCE_ADJUSTMENT[radiance_adjustment][view]:
             return S3_RADIANCE_ADJUSTMENT[radiance_adjustment][view][self.band_id]
         return None
@@ -1526,12 +1483,8 @@ def _get_solar_flux_dict_for_band(self, product_type):
     def get_solar_flux(self, zip_file, root_path):
         """
         Get solar spectral flux in mW / (m^2 * sr * nm) for band
-
-        Args:
-            band_object (S3Band): Optical Band with radiance
-
-        Returns:
-            float: solar Flux
+        :returns: solar Flux
+        :type: float
 
         """
         if not self.solar_flux:
@@ -1590,7 +1543,7 @@ def main():
     s3_files = [Path(scene) for scene in s3_files]
     for s3_scene in s3_files:
         if not s3_scene.exists():
-            gs.fatal(_("Input file <{sn}> not found").format(str(sn=s3_scene)))
+            gs.fatal(_("Input file <{sn}> not found").format(sn=str(s3_scene)))
         if not re.match(pattern, s3_scene.name):
             gs.fatal(
                 _(
@@ -1611,8 +1564,7 @@ def main():
     )
 
     meta_info_dict = extract_file_info(s3_files, basename=options["basename"])
-    # Group scenes by mission,
-    # groups = group_scenes(s3_files, granularity="1 day")
+
     if gs.parse_command("g.proj", flags="g")["proj"] == "ll":
         gs.fatal(_("Running in lonlat location is currently not supported"))
 
@@ -1621,49 +1573,38 @@ def main():
     # Get region bounds
     region_bounds = gs.parse_command("g.region", flags="ugb", quiet=True)
     current_region = gs.parse_command("g.region", flags="ug")
-    # has_band_ref = float(gs.version()["version"][0:3]) >= 7.9
-
-    json_standard_folder = None
-    if flags["j"]:
-        env = gs.gisenv()
-        json_standard_folder = Path(env["GISDBASE"]).joinpath(
-            env["LOCATION_NAME"], env["MAPSET"], "cell_misc"
-        )
-        if not json_standard_folder.exists():
-            json_standard_folder.mkdir(parents=True, exist_ok=True)
 
     # Collect variables for import
     import_dict = {
         "reg_bounds": dict(region_bounds),
         "current_reg": dict(current_region),
         "mod_flags": flags,
-        "json_standard_folder": json_standard_folder,
         "meta_dict": meta_info_dict,
     }
 
-    if flags["p"]:
-        print(
-            "|".join(
-                [
-                    "product_file_name",
-                    "nc_file_name",
-                    "nc_file_title",
-                    "nc_file_start_time",
-                    "nc_file_creation_time",
-                    "band",
-                    "band_shape",
-                    "band_title",
-                    "band_standard_name",
-                    "band_long_name",
-                ]
-            )
-        )
-        print("\n".join(chain(*import_result)))
-        cleanup()
-        return 0
+    # if flags["p"]:
+    #     print(
+    #         "|".join(
+    #             [
+    #                 "product_file_name",
+    #                 "nc_file_name",
+    #                 "nc_file_title",
+    #                 "nc_file_start_time",
+    #                 "nc_file_creation_time",
+    #                 "band",
+    #                 "band_shape",
+    #                 "band_title",
+    #                 "band_standard_name",
+    #                 "band_long_name",
+    #             ]
+    #         )
+    #     )
+    #     print("\n".join(chain(*import_result)))
+    #     cleanup()
+    #     return 0
 
     module_queues = []
-    register_strings = {}
+    register_strings = []
     region_dicts = {}
     region_list = []
 
@@ -1673,7 +1614,7 @@ def main():
             s3_file, import_dict, s3_product=s3_product
         )
         module_queues.append(module_list)
-        register_strings.update(register_string)
+        register_strings.extend(register_string)
         region_list.append(region_dict)
         if not region_dicts:
             region_dicts = region_dict.copy()
@@ -1686,22 +1627,20 @@ def main():
 
     stripe_envs = {}
     for stripe_id, stripe_region in region_dicts.items():
-        stripe_env = os.environ.copy()
-        stripe_env["GRASS_REGION"] = gs.region_env(
-            **adjust_region(stripe_region)
-        )  # , env=stripe_env)
-        stripe_envs[stripe_id] = stripe_env
-
-        compute_env = stripe_envs["sun_parameters"]
+        if flags["n"] or stripe_id == "sun_parameters":
+            stripe_env = os.environ.copy()
+            stripe_env["GRASS_REGION"] = gs.region_env(**adjust_region(stripe_region))
+            stripe_envs[stripe_id] = stripe_env
 
     if flags["r"]:
         gs.verbose(_("Importing solar parameter bands"))
         queue = ParallelModuleQueue(nprocs)
+        compute_env = stripe_envs["sun_parameters"]
         for solar_parameter in ["solar_azimuth", "solar_zenith"]:
             module_list = []
             for solar_module in module_queues[0][solar_parameter]:
                 solar_module.env_ = compute_env
-                solar_module.verbose = True
+                # solar_module.verbose = True
                 module_list.append(solar_module)
             queue.put(MultiModule(module_list))
         queue.wait()
@@ -1715,23 +1654,116 @@ def main():
     for band_id, modules in module_queues[0].items():
         if band_id in ["solar_azimuth", "solar_zenith"]:
             continue
-        compute_env = stripe_envs[band_id[-2:]]
-        module_list = []
-        for listed_module in modules:
-            listed_module.env_ = compute_env
-            listed_module.verbose = True
-            module_list.append(listed_module)
-        queue.put(MultiModule(module_list))
+        if flags["n"]:
+            compute_env = stripe_envs[band_id[-2:]]
+            module_list = []
+            for listed_module in modules:
+                listed_module.env_ = compute_env
+                # listed_module.verbose = True
+                module_list.append(listed_module)
+            queue.put(MultiModule(module_list))
+        else:
+            queue.put(MultiModule(modules))
     queue.wait()
 
+    # Update map support information
+    t_register_strings = []
+    queue = ParallelModuleQueue(nprocs)
+    for mapname, metadata in consolidat_metadata_dicts(register_strings).items():
+        mapname = (
+            mapname if not flags["r"] else mapname.replace("radiance", "reflectance")
+        )
+        metadata["semantic_label"] = (
+            metadata["semantic_label"]
+            if not flags["r"]
+            else metadata["semantic_label"].replace("radiance", "reflectance")
+        )
+        # Write raster history
+        gs.raster_history(mapname, overwrite=True)
+
+        # Write extended metadata if requested
+        if flags["j"]:
+            gs.verbose(_("Writing metadata to maps..."))
+            json_standard_folder = Path(GISENV["GISDBASE"]).joinpath(
+                GISENV["LOCATION_NAME"], GISENV["MAPSET"], "cell_misc", mapname
+            )
+            if not json_standard_folder.exists():
+                json_standard_folder.mkdir(parents=True, exist_ok=True)
+
+            write_metadata(metadata, str(json_standard_folder / "description.json"))
+
+        description = (
+            json.dumps(metadata, separators=["\n", ": "]).lstrip("{").rstrip("}")
+        )
+        support_kwargs = {
+            "map": mapname,
+            "title": ",".join(metadata["title"])
+            if isinstance(metadata["title"], list)
+            else metadata["title"],
+            "history": ",".join(metadata["history"])
+            if isinstance(metadata["history"], list)
+            else metadata["history"],
+            "units": metadata["unit"],
+            "source1": ",".join(metadata["product_name"])
+            if isinstance(metadata["product_name"], list)
+            else metadata["product_name"],
+            "source2": ",".join(metadata["processing_baseline"])
+            if isinstance(metadata["processing_baseline"], list)
+            else metadata["processing_baseline"],
+            "description": description,
+            "semantic_label": metadata["semantic_label"],
+        }
+
+        queue.put(
+            MultiModule(
+                [
+                    Module(
+                        "r.support",
+                        **support_kwargs,
+                        quiet=True,
+                        run_=False,
+                    ),
+                    Module(
+                        "r.timestamp",
+                        map=mapname,
+                        date="/".join(
+                            [
+                                grass_timestamp(meta_info_dict[1]["start_time"]),
+                                grass_timestamp(meta_info_dict[1]["end_time"]),
+                            ]
+                        ),
+                        quiet=True,
+                        run_=False,
+                    ),
+                ]
+            )
+        )
+        t_register_strings.append(
+            "|".join(
+                [
+                    f"{mapname}@{GISENV['MAPSET']}",
+                    meta_info_dict[1]["start_time"].isoformat(
+                        sep=" ", timespec="seconds"
+                    ),
+                    meta_info_dict[1]["end_time"].isoformat(
+                        sep=" ", timespec="seconds"
+                    ),
+                    metadata["semantic_label"],
+                ]
+            )
+        )
+    queue.wait()
+
+    # Write t.register file if requested
     if options["register_output"]:
-        # Write t.register file if requested
-        write_register_file(
-            options["register_output"], "\n".join(register_strings.values())
+        gs.verbose(
+            _("Writing register file <{}>...").format(options["register_output"])
+        )
+        Path(options["register_output"]).write_text(
+            "\n".join(t_register_strings) + "\n", encoding="UTF8"
         )
     else:
-        print("\n".join(register_strings.values()))
-    cleanup()
+        print("\n".join(t_register_strings))
     return 0
 
 

From 1f102361ea4c1cf7d2036a8708f803f7b3376b7c Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Fri, 23 Feb 2024 15:11:36 +0100
Subject: [PATCH 05/27] extend basename

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 33 +++++++++++--------
 1 file changed, 19 insertions(+), 14 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 7aeb7da31a..9c972c70f8 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -169,6 +169,7 @@
 
 from collections import OrderedDict
 from datetime import datetime
+
 # from itertools import chain
 from pathlib import Path
 from zipfile import ZipFile
@@ -660,16 +661,23 @@ def parse_s3_file_name(file_name):
 def extract_file_info(s3_files, basename=None):
     """Extract information from file name according to naming conventions"""
     result_dict = {}
+    product_track_ids = [
+        "cycle",
+        "relative_orbit",
+        "producer",
+        "product_class",
+        "timeliness",
+        "baseline_collection",
+    ]
     for s3_file in s3_files:
         file_info = parse_s3_file_name(s3_file.name)
         if not result_dict:
             result_dict = file_info
-            result_dict["duration"] = {file_info["duration"]}
-            result_dict["cycle"] = {file_info["cycle"]}
-            result_dict["relative_orbit"] = {file_info["relative_orbit"]}
-            result_dict["frame"] = {file_info["frame"]}
             result_dict["start_time"] = file_info["start_time"]
             result_dict["end_time"] = file_info["end_time"]
+            result_dict["frame"] = {file_info["frame"]}
+            for pid in product_track_ids:
+                result_dict[pid] = {file_info[pid]}
         else:
             if file_info["mission_id"] != result_dict["mission_id"]:
                 result_dict["mission_id"] = "S3_"
@@ -679,15 +687,14 @@ def extract_file_info(s3_files, basename=None):
             result_dict["end_time"] = max(
                 result_dict["end_time"], file_info["end_time"]
             )
-            result_dict["duration"].add(file_info["duration"])
-            result_dict["cycle"].add(file_info["cycle"])
-            result_dict["relative_orbit"].add(file_info["relative_orbit"])
+            for pid in product_track_ids:
+                result_dict[pid].add(file_info[pid])
             result_dict["frame"].add(file_info["frame"])
-    for key in ["cycle", "relative_orbit"]:
-        if len(result_dict[key]) > 1:
+    for pid in product_track_ids:
+        if len(result_dict[pid]) > 1:
             gs.warning(
-                _("Merging {key}s {values}").format(
-                    key=key, values=", ".join(result_dict[key])
+                _("Merging {key} {values}").format(
+                    key=pid, values=", ".join(result_dict[pid])
                 )
             )
     if result_dict["mission_id"] == "3_":
@@ -701,9 +708,7 @@ def extract_file_info(s3_files, basename=None):
                 result_dict["product"],
                 result_dict["start_time"].strftime("%Y%m%d%H%M%S"),
                 result_dict["end_time"].strftime("%Y%m%d%H%M%S"),
-                list(result_dict["duration"])[0],
-                list(result_dict["cycle"])[0],
-                list(result_dict["relative_orbit"])[0],
+                *[list(result_dict[pid])[0] for pid in product_track_ids],
             ]
         )
     return (

From a385bc962a981b80e2973393813ec7c196b5aa93 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Fri, 23 Feb 2024 15:25:31 +0100
Subject: [PATCH 06/27] update to new version

---
 .../i.sentinel3.import.html                   | 51 ++++++++++++-------
 1 file changed, 33 insertions(+), 18 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html
index ca93f9d3d5..ca08095a58 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.html
@@ -1,6 +1,9 @@
 
DESCRIPTION

 
-The i.sentinel3.import module allows importing Copernicus Sentinel-3 products.
+The i.sentinel3.import module allows importing
+
+Copernicus Sentinel-3 products.
+
 

 The Sentinel-3 products types are provided in different formats. Not all of
 them are currently directly supported by GDAL. Some, like the e.g. the SLSTR products
@@ -22,11 +25,13 @@ 
DESCRIPTION

 
 

 All Sentinel-3 scene files provided in the input option are mosaiced
-into one resulting granule with one map for each imported bands. The input
+into one resulting granule with one map for each imported band. The input
 option accepts either one or many path(s) to Sentinel-3 archives (.zip) separated
-by comma or a textfile with  one Sentinel-3 file per row. The module is designed
-to produce temporally non overlaping mosaics along the tracks of the satellite.
-Input should therefor not mix data from different tracks or sensing periods.
+by comma or a text-file with one Sentinel-3 file per row. The module is designed
+to produce temporally non-overlaping mosaics along the tracks of the satellite.
+Input should therefore not mix data from different tracks, sensing periods,
+or missions (Sentine-3A / Sentinel-3B). The module gives a warning if data from
+different tracks is used as input.
 
 

 Users can provide a custom basename, to which band names are appended.
@@ -49,9 +54,17 @@ 
DESCRIPTION

 map. After import, the metadata can be retrieved with r.info -e
 as shown below.
 
-
+

+i.sentinel3.import imports the selected bands with the resolution and extent of the current
+region (except for the solar bands), unless the n-flag is used to import bands with their
+native resolution (usually 500m or 1000m).
+    
 
NOTES

 
+i.sentinel3.import does ONLY support projected coordinate systems and assumes map units to
+be meter.
+
+

 By means of the register_file option i.sentinel3.import allows
 creating a file which can be used to register imported imagery data
 into a space-time raster dateset (STRDS) with
@@ -81,29 +94,30 @@ 
EXAMPLES

 -->
 
 
Import Sentinel-3 data

-Import all Sentinel-3 data found in data directory and store metadata
-as JSON files within the GRASS GIS database directory:
+Import all Sentinel-3B RBT data of a given start date and track found in data directory, store metadata
+as JSON files within the GRASS GIS database directory, keep pixel values during interpolation (no smoothing),
+import bands at native resolution and rescale radiance bands to reflectance:
 
-i.sentinel3.import -a -c -j input=./data/ nprocs=4 product_type=S3SL2LST basename=S3_LST
+i.sentinel3.import -k -n -r -j input="$(find data/ -iname "S3B*RBT*____20240129*_089_094_*.zip" | tr '\n' ',' | sed 's/,$/\n/')" \
+    product_type=S3SL1RBT register_output=data/reg.txt nprocs=8 --v --o
+
+# register imported data into existing STRDS
+t.register input=Sentinel_3 file=data/reg.txt
 

 
 
 
Register imported Sentinel-3 data into STRDS

 
 
-i.sentinel3.import -a -c -j input=./data/ register_output=./data/reg.txt nprocs=4 \
-product_type=S3SL2LST basename=S3_LST
+i.sentinel3.import -k -n -c -j  input="$(find data/ -iname "S3B*LST*____20240129*_089_094_*.zip" | tr '\n' ',' | sed 's/,$/\n/')" \
+product_type=S3SL2LST register_output=data/reg.txt nprocs=8 --v --o
 
 # register imported data into existing STRDS
-t.register input=Sentinel_3 file=data/reg.txt
+t.register input=Sentinel_3_LST file=data/reg.txt
 

 
-A register file typically contains two columns: imported raster map
-name and timestamp separated by |. In the case of current
-development version of GRASS GIS which supports band references concept
-(see g.bands module for details) a
-register file is extended by a third column containg band reference
-information, see the examples below.
+A register file typically contains the following columns: imported raster map
+name, start and end timestamp as well as semantic label, all separated by |, see the examples below:
 
 
 S3_imp_surface_temperature_20210907T205026|2021-09-07T20:50:26.285902+00:00|2021-09-07T20:53:25.977846+00:00|S3_surface_temperature
@@ -129,6 +143,7 @@ 
SEE ALSO

 
 i.sentinel.download,
 r.in.xyz,
+r.fill.stats,
 t.register
 
 

From cf7d30e288ec520969d3dce181124c0f33660b46 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Sat, 24 Feb 2024 17:08:19 +0100
Subject: [PATCH 07/27] zrange and LST

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 67 +++++++++----------
 1 file changed, 32 insertions(+), 35 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 9c972c70f8..740b3d2be5 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -343,14 +343,14 @@
         "S3SL2LST": {
             "LST": {
                 "geometries": ("i",),
-                "nc_file": "LST_{}.nc",
+                "nc_file": "{}_{}.nc",
                 "full_name": "{}",
                 "types": None,
                 "exception": "exception",
             },
             "LST_uncertainty": {
                 "geometries": ("i",),
-                "nc_file": "LST_{}.nc",
+                "nc_file": "LST_in.nc",
                 "full_name": "{}",
                 "types": None,
                 "exception": "exception",
@@ -936,36 +936,32 @@ def get_band_metadata(
     # Handle offset, scale, and valid data range, see:
     # https://docs.unidata.ucar.edu/netcdf-c/current/attribute_conventions.html
     metadata["zrange"] = None
+    min_val, max_val = None, None
     if "valid_range" in band_attrs:
         min_val, max_val = nc_variable.valid_range
-    elif "valid_max" in band_attrs:
+    if "valid_min" in band_attrs:
         min_val = nc_variable.valid_min
+    if "valid_max" in band_attrs:
         max_val = nc_variable.valid_max
-    elif "_FillValue" in band_attrs:
-        data_range = get_dtype_range(datatype)
-        fill_val = nc_variable._FillValue
+    data_range = get_dtype_range(datatype)
+    if "_FillValue" in band_attrs:
+        fill_val = np.float64(nc_variable._FillValue)
         if fill_val > 0:
-            min_val = data_range["min"]
-            max_val = fill_val - data_range["resolution"]
-        elif fill_val < 0:
-            min_val = fill_val + data_range["resolution"]
-            max_val = data_range["max"]
-        else:
-            min_val, max_val = None, None
+            min_val = data_range["min"] if min_val is None else min_val
+            max_val = fill_val - 1 if max_val is None else max_val  # data_range["resolution"]
+        elif fill_val <= 0:
+            min_val = fill_val + 1 if min_val is None else min_val  # data_range["resolution"]
+            max_val = data_range["max"] if max_val is None else max_val
     else:
-        min_val, max_val = None, None
-
-    if min_val and max_val:
-        if "add_offset" in band_attrs and min_val and max_val:
-            min_val = min_val + nc_variable.add_offset
-            max_val = max_val + nc_variable.add_offset
+        min_val = data_range["min"] + 1 if min_val is None else min_val  # data_range["resolution"]
+        max_val = data_range["max"] if max_val is None else max_val
 
-        if "scale_factor" in band_attrs and min_val and max_val:
-            min_val = min_val * nc_variable.scale_factor
-            max_val = max_val * nc_variable.scale_factor
-        metadata["description"] += f"\n\nvalid_min: {min_val}\nvalid_max: {max_val}"
-        metadata["zrange"] = [min_val, max_val]
+    if "flag_masks" in band_attrs:
+        min_val = min(nc_variable.flag_masks)
+        max_val = max(nc_variable.flag_masks)
 
+    metadata["description"] += f"\n\nvalid_min: {min_val}\nvalid_max: {max_val}"
+    metadata["zrange"] = [min_val, max_val]
     metadata["semantic_label"] = f"S3_{varname_short}"
 
     return metadata, fmt
@@ -1220,7 +1216,6 @@ def process_nc_file(
         nc_file_path = zip_file.extract(member, path=TMP_FILE)
         with Dataset(nc_file_path) as nc_file_open:
             file_metadata = get_file_metadata(nc_file_open)
-
             for band in container_dict_items[1]:
                 # Check for band
                 solar_flux = None
@@ -1235,7 +1230,6 @@ def process_nc_file(
                     band = self.flag_bands[band]
                 elif band in self.anxillary_bands:
                     band = self.anxillary_bands[band]
-
                 if band.full_name not in nc_file_open.variables:
                     gs.fatal(
                         _(
@@ -1271,19 +1265,19 @@ def process_nc_file(
                 if band.exception:
                     if np.ma.is_masked(add_array):
                         add_array.mask = np.ma.mask_or(
-                            add_array.mask, nc_file_open[band.exception][:] > 0
+                            add_array.mask, nc_file_open[band.exception][:] >= 3
                         )
                     else:
                         add_array = np.ma.masked_array(
-                            add_array, nc_file_open[band.exception][:] > 0
+                            add_array, nc_file_open[band.exception][:] >= 3
                         )
                 if np.ma.is_masked(add_array):
                     add_array = add_array.filled()
 
                 # Rescale temperature variables if requested
-                if band.full_name.startswith("LST") and module_flags["c"]:
+                if band.full_name == "LST" and module_flags["c"]:
                     add_array = convert_units(add_array, "K", "degree_celsius")
-
+ 
                 # Write metadata json if requested
                 band_attrs = nc_file_open[band.full_name].ncattrs()
 
@@ -1378,7 +1372,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds):
                 map_name,
                 distance=3,
                 fill_flags="ks",
-                zrange=None,
+                zrange=[-32767,32768],
                 val_col=len(nc_bands),
                 data_type="FCELL",
                 method="mean",
@@ -1473,7 +1467,7 @@ def _get_exception_band(self, product_type):
         return None
 
     def _get_solar_flux_dict_for_band(self, product_type):
-        if self.band_id in S3_SOLAR_FLUX[product_type]["defaults"]:
+        if S3_SOLAR_FLUX[product_type] and self.band_id in S3_SOLAR_FLUX[product_type]["defaults"]:
             return {
                 "default": S3_SOLAR_FLUX[product_type]["defaults"][self.band_id],
                 "band": S3_SOLAR_FLUX[product_type]["band"].format(
@@ -1660,11 +1654,14 @@ def main():
         if band_id in ["solar_azimuth", "solar_zenith"]:
             continue
         if flags["n"]:
-            compute_env = stripe_envs[band_id[-2:]]
+            # S3SL2LST product has only "in" stripe and band names have no suffix
+            if band_id[-2:] not in stripe_envs:
+                compute_env = stripe_envs["in"]
+            else:
+                compute_env = stripe_envs[band_id[-2:]]
             module_list = []
             for listed_module in modules:
                 listed_module.env_ = compute_env
-                # listed_module.verbose = True
                 module_list.append(listed_module)
             queue.put(MultiModule(module_list))
         else:
@@ -1674,6 +1671,7 @@ def main():
     # Update map support information
     t_register_strings = []
     queue = ParallelModuleQueue(nprocs)
+    gs.verbose(_("Writing metadata to maps..."))
     for mapname, metadata in consolidat_metadata_dicts(register_strings).items():
         mapname = (
             mapname if not flags["r"] else mapname.replace("radiance", "reflectance")
@@ -1688,7 +1686,6 @@ def main():
 
         # Write extended metadata if requested
         if flags["j"]:
-            gs.verbose(_("Writing metadata to maps..."))
             json_standard_folder = Path(GISENV["GISDBASE"]).joinpath(
                 GISENV["LOCATION_NAME"], GISENV["MAPSET"], "cell_misc", mapname
             )

From e558826be5244f4486a1d7a2e4765e29472d4f10 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Sat, 24 Feb 2024 23:39:57 +0100
Subject: [PATCH 08/27] fix typo

---
 .../i.sentinel/i.sentinel3.import/i.sentinel3.import.py     | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 740b3d2be5..50955a6e93 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -535,7 +535,7 @@ def get_dtype_range(np_datatype):
     }
 
 
-def consolidat_metadata_dicts(metadata_dicts):
+def consolidate_metadata_dicts(metadata_dicts):
     """Consolidate a list of metadata dictionaries for all
     input scenes into unified metadata for the resulting map"""
     map_metadata_dicts = {}
@@ -945,7 +945,7 @@ def get_band_metadata(
         max_val = nc_variable.valid_max
     data_range = get_dtype_range(datatype)
     if "_FillValue" in band_attrs:
-        fill_val = np.float64(nc_variable._FillValue)
+        fill_val = nc_variable._FillValue
         if fill_val > 0:
             min_val = data_range["min"] if min_val is None else min_val
             max_val = fill_val - 1 if max_val is None else max_val  # data_range["resolution"]
@@ -1672,7 +1672,7 @@ def main():
     t_register_strings = []
     queue = ParallelModuleQueue(nprocs)
     gs.verbose(_("Writing metadata to maps..."))
-    for mapname, metadata in consolidat_metadata_dicts(register_strings).items():
+    for mapname, metadata in consolidate_metadata_dicts(register_strings).items():
         mapname = (
             mapname if not flags["r"] else mapname.replace("radiance", "reflectance")
         )

From e6cf555c86f0e570de0890bab095b00b9a4a7965 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Mon, 4 Mar 2024 23:30:57 +0100
Subject: [PATCH 09/27] import intersection

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 23 +++++++++++++++++--
 1 file changed, 21 insertions(+), 2 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 50955a6e93..bbee928104 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -630,6 +630,23 @@ def adjust_region(region_dict):
     return region_dict
 
 
+def intersect_region(region_a, region_b):
+    """
+    Adjust region bounds for r.in.xyz (bounds + half resolution)
+    aligned to resolution in ew and ns direction
+    starting from the sw-corner
+    """
+    relevant_keys = ["n", "s", "e", "w", "nsres", "ewres"]
+    region_dict = {"nsres": region_a["nsres"], "ewres": region_a["ewres"]}
+    region_a = {region_key: float(region_a[region_key]) for region_key in relevant_keys}
+    region_b = {region_key: float(region_b[region_key]) for region_key in relevant_keys}
+    region_dict["s"] = region_a["s"] if region_a["s"] >= region_b["s"] else region_a["s"] + np.ceil((region_a["s"] - region_b["s"]) / region_a["nsres"]) * region_a["nsres"]
+    region_dict["w"] = region_a["w"] if region_a["w"] >= region_b["w"] else region_a["w"] + np.ceil((region_a["w"] - region_b["w"]) / region_a["ewres"]) * region_a["ewres"]
+    region_dict["n"] = region_a["n"] if region_a["n"] <= region_b["n"] else region_a["n"] - np.floor((region_a["n"] - region_b["n"]) / region_a["nsres"]) * region_a["nsres"]
+    region_dict["e"] = region_a["e"] if region_a["e"] <= region_b["e"] else region_a["e"] - np.floor((region_a["e"] - region_b["e"]) / region_a["ewres"]) * region_a["ewres"]
+    return region_dict
+
+
 def parse_s3_file_name(file_name):
     """Extract info from file name according to naming onvention:
     https://sentinels.copernicus.eu/web/sentinel/user-guides/sentinel-3-slstr/naming-convention
@@ -1626,10 +1643,12 @@ def main():
 
     stripe_envs = {}
     for stripe_id, stripe_region in region_dicts.items():
+        stripe_env = os.environ.copy()
         if flags["n"] or stripe_id == "sun_parameters":
-            stripe_env = os.environ.copy()
             stripe_env["GRASS_REGION"] = gs.region_env(**adjust_region(stripe_region))
-            stripe_envs[stripe_id] = stripe_env
+        else:
+            stripe_env["GRASS_REGION"] = gs.region_env(**intersect_region(dict(current_region), stripe_region))
+        stripe_envs[stripe_id] = stripe_env
 
     if flags["r"]:
         gs.verbose(_("Importing solar parameter bands"))

From 90027dcc387c696f531e07012b3290594e47d982 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Tue, 12 Mar 2024 13:03:28 +0100
Subject: [PATCH 10/27] fix intersection, add max sun angle

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 215 +++++++++++-------
 1 file changed, 136 insertions(+), 79 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index bbee928104..6bd0cef8d6 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -81,6 +81,13 @@
 # % required: no
 # %end
 
+# %option
+# % key: maximum_solar_angle
+# % type: double
+# % description: Import only pixels where solar angle is lower or equal to the given maximum
+# % required: no
+# %end
+
 # %option G_OPT_M_NPROCS
 # %end
 
@@ -630,20 +637,28 @@ def adjust_region(region_dict):
     return region_dict
 
 
-def intersect_region(region_a, region_b):
+def intersect_region(region_current, region_stripe, align_current=True):
     """
-    Adjust region bounds for r.in.xyz (bounds + half resolution)
+    Adjust region bounds for r.in.xyz to the intersection
+    of the current region (bounds + half resolution)
     aligned to resolution in ew and ns direction
     starting from the sw-corner
     """
     relevant_keys = ["n", "s", "e", "w", "nsres", "ewres"]
-    region_dict = {"nsres": region_a["nsres"], "ewres": region_a["ewres"]}
-    region_a = {region_key: float(region_a[region_key]) for region_key in relevant_keys}
-    region_b = {region_key: float(region_b[region_key]) for region_key in relevant_keys}
-    region_dict["s"] = region_a["s"] if region_a["s"] >= region_b["s"] else region_a["s"] + np.ceil((region_a["s"] - region_b["s"]) / region_a["nsres"]) * region_a["nsres"]
-    region_dict["w"] = region_a["w"] if region_a["w"] >= region_b["w"] else region_a["w"] + np.ceil((region_a["w"] - region_b["w"]) / region_a["ewres"]) * region_a["ewres"]
-    region_dict["n"] = region_a["n"] if region_a["n"] <= region_b["n"] else region_a["n"] - np.floor((region_a["n"] - region_b["n"]) / region_a["nsres"]) * region_a["nsres"]
-    region_dict["e"] = region_a["e"] if region_a["e"] <= region_b["e"] else region_a["e"] - np.floor((region_a["e"] - region_b["e"]) / region_a["ewres"]) * region_a["ewres"]
+    if align_current:
+        region_dict = {"nsres": region_current["nsres"], "ewres": region_current["ewres"]}
+    else:
+        region_dict = {"nsres": region_stripe["nsres"], "ewres": region_stripe["ewres"]}
+    region_current = {region_key: float(region_current[region_key]) for region_key in relevant_keys}
+    region_stripe = {region_key: float(region_stripe[region_key]) for region_key in relevant_keys}
+    stripe_s = region_stripe["s"] - region_stripe["nsres"] / 2.0
+    stripe_n = region_stripe["n"] + region_stripe["nsres"] / 2.0
+    stripe_w = region_stripe["w"] - region_stripe["ewres"] / 2.0
+    stripe_e = region_stripe["e"] + region_stripe["ewres"] / 2.0
+    region_dict["s"] = region_current["s"] if region_current["s"] >= stripe_s else region_current["s"] + np.ceil((region_current["s"] - stripe_s) / float(region_dict["nsres"])) * float(region_dict["nsres"])
+    region_dict["w"] = region_current["w"] if region_current["w"] >= stripe_w else region_current["w"] + np.ceil((region_current["w"] - stripe_w) / float(region_dict["ewres"])) * float(region_dict["ewres"])
+    region_dict["n"] = region_current["n"] if region_current["n"] <= stripe_n else region_current["n"] - np.floor((region_current["n"] - stripe_n) / float(region_dict["nsres"])) * float(region_dict["nsres"])
+    region_dict["e"] = region_current["e"] if region_current["e"] <= stripe_e else region_current["e"] - np.floor((region_current["e"] - stripe_e) / float(region_dict["ewres"])) * float(region_dict["ewres"])
     return region_dict
 
 
@@ -734,7 +749,7 @@ def extract_file_info(s3_files, basename=None):
     )
 
 
-def get_geocoding(zip_file, root, geo_bands_dict, region_bounds=None):
+def get_geocoding(zip_file, root, geo_bands_dict, sun_mask=None, region_bounds=None):
     """Get ground control points from NetCDF file"""
     member = str(root / geo_bands_dict["nc_file"])
     nc_file_path = zip_file.extract(member, path=TMP_FILE)
@@ -758,7 +773,10 @@ def get_geocoding(zip_file, root, geo_bands_dict, region_bounds=None):
             nc_bands[band_id] = nc_file_open[band][:]
 
         # Create initial mask
-        mask = nc_bands["lat"][:].mask
+        if sun_mask is not None:
+            mask = np.logical_or(sun_mask, nc_bands["lat"][:].mask)
+        else:
+            mask = nc_bands["lat"][:].mask
 
         if region_bounds:
             # Mask to region
@@ -780,7 +798,8 @@ def get_geocoding(zip_file, root, geo_bands_dict, region_bounds=None):
             )
 
             if mask.all():
-                gs.fatal(_("No valid pixels inside computational region"))
+                gs.warning(_("No valid pixels inside computational region found in {}").format(str(zip_file.filename)))
+                return None, None, None
 
     return nc_bands, mask, resolution
 
@@ -1061,19 +1080,27 @@ def import_s3(s3_file, kwargs, s3_product=None):
         root = members[0].rsplit(".SEN3", maxsplit=1)[0]
         root = Path(f"{root}.SEN3")
         # Check if solar flux raster band is required
-        if mod_flags["r"]:
+        sun_region_bounds = region_bounds.copy()
+        if mod_flags["r"] or kwargs["maximum_solar_angle"]:
             tmp_ascii_sun = TMP_FILE / f"{rmap}_sun_parameters.txt"
             module_queue, sun_region_dict = s3_product.get_sun_parameters(
-                zip_file, root, rmap, tmp_ascii_sun, region_bounds
+                zip_file, root, rmap, tmp_ascii_sun, region_bounds, maximum_solar_angle=kwargs["maximum_solar_angle"]
             )
+            if module_queue is None:
+                return None, None, None
             region_dicts["sun_parameters"] = sun_region_dict
+
+            if kwargs["maximum_solar_angle"]:
+                sun_region_bounds = gs.parse_command("g.region", flags="ugb", quiet=True, **intersect_region(dict(kwargs["current_reg"]), sun_region_dict, align_current=False))
         for stripe, stripe_dict in s3_product.requested_stripe_content.items():
             # Could be parallelized!!!
             # Get geocoding (lat/lon, elevation) and initial mask
             tmp_ascii = TMP_FILE / f"{rmap}_{stripe}.txt"
             nc_bands, mask, resolution = get_geocoding(
-                zip_file, root, stripe_dict["geocoding"], region_bounds=region_bounds
+                zip_file, root, stripe_dict["geocoding"], region_bounds=sun_region_bounds
             )
+            if nc_bands is None:
+                return None, None, None
             fmt = ",".join(["%.12f"] * len(nc_bands))
             for nc_file, bands in stripe_dict["containers"].items():
                 (
@@ -1093,12 +1120,16 @@ def import_s3(s3_file, kwargs, s3_product=None):
                     fmt=fmt,
                 )
                 register_strings.append(register_output)
-            region_dicts[stripe] = write_xyz(
+            stripe_region_dict = write_xyz(
                 tmp_ascii, nc_bands, mask, fmt=fmt, project=True
             )
-            region_dicts[stripe]["ewres"], region_dicts[stripe]["nsres"] = float(
-                resolution[0]
-            ), float(resolution[1])
+            stripe_region_dict["ewres"] = float(resolution[0])
+            stripe_region_dict["nsres"] = float(resolution[1])
+            sun_region_dict_intersect = intersect_region(dict(kwargs["current_reg"]), sun_region_dict, align_current=False)
+            if kwargs["maximum_solar_angle"]:
+                region_dicts[stripe] = intersect_region(sun_region_dict_intersect, stripe_region_dict, align_current=False)
+            else:
+                region_dicts[stripe] = stripe_region_dict
 
     return module_queue, register_strings, region_dicts
 
@@ -1352,19 +1383,13 @@ def process_nc_file(
 
             return nc_bands, module_queue, meta_information, fmt
 
-    def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds):
+    def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, maximum_solar_angle=None):
         """https://github.com/sertit/eoreader/blob/main/eoreader/products/optical/s3_slstr_product.py#L862"""
 
         # nc_file = sun_azimuth_dict["nc_file"]
         # Get values
         import_modules = {}
         fmt = "%.12f,%.12f"
-        nc_bands, mask, resolution = get_geocoding(
-            zip_file,
-            root,
-            S3_SUN_PARAMTERS[self.product_type]["geocoding"],
-            region_bounds=region_bounds,
-        )
         member = str(
             root
             / S3_SUN_PARAMTERS[self.product_type]["sun_bands"]["nc_file"].format(
@@ -1372,35 +1397,48 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds):
             )
         )
         nc_file_path = zip_file.extract(member, path=TMP_FILE)
-        sun_parameter_nc = Dataset(nc_file_path)
-        # sun_metadata = get_file_metadata(sun_parameter_nc)
-        for band_id, band in S3_SUN_PARAMTERS[self.product_type]["sun_bands"][
-            "bands"
-        ].items():
-            fmt += ",%.12f"
-            sun_parameter_array = sun_parameter_nc[band.format(self.view)]
-            nc_bands[band_id] = sun_parameter_array[:]
-            map_name = f"{prefix}_{band_id}"
-            if band_id == "solar_zenith":
-                global SUN_ZENITH_ANGLE
-                SUN_ZENITH_ANGLE = map_name
-            import_modules[band_id] = setup_import_multi_module(
-                tmp_ascii,
-                map_name,
-                distance=3,
-                fill_flags="ks",
-                zrange=[-32767,32768],
-                val_col=len(nc_bands),
-                data_type="FCELL",
-                method="mean",
-                solar_flux=None,
+        # sun_region_bounds = region_bounds.copy()
+        with Dataset(nc_file_path) as sun_parameter_nc:
+            if maximum_solar_angle:
+                sun_mask = sun_parameter_nc[S3_SUN_PARAMTERS[self.product_type]["sun_bands"]["bands"]["solar_zenith"].format(self.view)][:] >= float(maximum_solar_angle)
+            nc_bands, mask, resolution = get_geocoding(
+                zip_file,
+                root,
+                S3_SUN_PARAMTERS[self.product_type]["geocoding"],
+                sun_mask=sun_mask if maximum_solar_angle and np.any(sun_mask) else None,
+                region_bounds=region_bounds,
             )
+            if nc_bands is None:
+                return None, None
+
+            # sun_metadata = get_file_metadata(sun_parameter_nc)
+            for band_id, band in S3_SUN_PARAMTERS[self.product_type]["sun_bands"][
+                "bands"
+            ].items():
+                fmt += ",%.12f"
+                sun_parameter_array = sun_parameter_nc[band.format(self.view)]
+                nc_bands[band_id] = sun_parameter_array[:]
+                map_name = f"{prefix}_{band_id}"
+                if band_id == "solar_zenith":
+                    global SUN_ZENITH_ANGLE
+                    SUN_ZENITH_ANGLE = map_name
+                import_modules[band_id] = setup_import_multi_module(
+                    tmp_ascii,
+                    map_name,
+                    distance=3,
+                    fill_flags="ks",
+                    zrange=[-32767,32768],
+                    val_col=len(nc_bands),
+                    data_type="FCELL",
+                    method="mean",
+                    solar_flux=None,
+                )
 
-        # Write to temporary file
-        sun_region_dict = write_xyz(tmp_ascii, nc_bands, mask, fmt=fmt, project=True)
-        sun_region_dict["ewres"], sun_region_dict["nsres"] = float(
-            resolution[0]
-        ), float(resolution[1])
+            # Write to temporary file
+            sun_region_dict = write_xyz(tmp_ascii, nc_bands, mask, fmt=fmt, project=True)
+            sun_region_dict["ewres"], sun_region_dict["nsres"] = float(
+                resolution[0]
+            ), float(resolution[1])
 
         return import_modules, sun_region_dict
 
@@ -1531,6 +1569,19 @@ def adjust_radiance(self, np_band_array):
         return np_band_array
 
 
+def get_solar_angle_bounds(region_bounds, sun_region_dict):
+    """"""
+    solar_bounds = gs.parse_command("g.region", flags="ugl", **sun_region_dict)
+    print(solar_bounds)
+    solar_bounds["ll_n"] = max(float(solar_bounds["nw_lat"]), float(solar_bounds["ne_lat"]))
+    solar_bounds["ll_s"] = min(float(solar_bounds["sw_lat"]), float(solar_bounds["se_lat"]))
+    if solar_bounds["ll_n"] < float(region_bounds["ll_n"]):
+        region_bounds["ll_n"] = solar_bounds["ll_n"]
+    if float(solar_bounds["ll_s"]) > float(region_bounds["ll_s"]):
+        region_bounds["ll_s"] = solar_bounds["ll_s"]
+    return region_bounds
+
+
 def main():
     """Do the main work"""
     pattern = re.compile(
@@ -1596,6 +1647,7 @@ def main():
         "current_reg": dict(current_region),
         "mod_flags": flags,
         "meta_dict": meta_info_dict,
+        "maximum_solar_angle": float(options["maximum_solar_angle"]) if options["maximum_solar_angle"] else None,
     }
 
     # if flags["p"]:
@@ -1629,25 +1681,28 @@ def main():
         module_list, register_string, region_dict = import_s3(
             s3_file, import_dict, s3_product=s3_product
         )
-        module_queues.append(module_list)
-        register_strings.extend(register_string)
-        region_list.append(region_dict)
-        if not region_dicts:
-            region_dicts = region_dict.copy()
-        else:
-            for stripe_id, region in region_dict.items():
-                if stripe_id not in region_dicts:
-                    region_dicts[stripe_id] = region
-                else:
-                    extend_region(region_dicts[stripe_id], region_dict[stripe_id])
-
+        if module_list is not None:
+            module_queues.append(module_list)
+            register_strings.extend(register_string)
+            region_list.append(region_dict)
+            if not region_dicts:
+                region_dicts = region_dict.copy()
+            else:
+                for stripe_id, region in region_dict.items():
+                    if stripe_id not in region_dicts:
+                        region_dicts[stripe_id] = region
+                    else:
+                        extend_region(region_dicts[stripe_id], region_dict[stripe_id])
+    if not module_queues:
+        gs.warning(_("Nothing to import with the given input"))
+        sys.exit(0)
     stripe_envs = {}
     for stripe_id, stripe_region in region_dicts.items():
         stripe_env = os.environ.copy()
-        if flags["n"] or stripe_id == "sun_parameters":
+        if flags["n"]:
             stripe_env["GRASS_REGION"] = gs.region_env(**adjust_region(stripe_region))
         else:
-            stripe_env["GRASS_REGION"] = gs.region_env(**intersect_region(dict(current_region), stripe_region))
+            stripe_env["GRASS_REGION"] = gs.region_env(**intersect_region(dict(current_region), stripe_region, align_current=stripe_id != "sun_parameters"))
         stripe_envs[stripe_id] = stripe_env
 
     if flags["r"]:
@@ -1662,6 +1717,11 @@ def main():
                 module_list.append(solar_module)
             queue.put(MultiModule(module_list))
         queue.wait()
+        # Here one could zoom to non-null pixels in solar angle map
+        # g.region zoom=zolar_zenith 
+        # if maximum_solar_angle is given in consistently limit 
+        # imported pixels
+
 
     gs.verbose(
         _("Importing scenes {}").format(
@@ -1672,19 +1732,16 @@ def main():
     for band_id, modules in module_queues[0].items():
         if band_id in ["solar_azimuth", "solar_zenith"]:
             continue
-        if flags["n"]:
-            # S3SL2LST product has only "in" stripe and band names have no suffix
-            if band_id[-2:] not in stripe_envs:
-                compute_env = stripe_envs["in"]
-            else:
-                compute_env = stripe_envs[band_id[-2:]]
-            module_list = []
-            for listed_module in modules:
-                listed_module.env_ = compute_env
-                module_list.append(listed_module)
-            queue.put(MultiModule(module_list))
+        # S3SL2LST product has only "in" stripe and band names have no suffix
+        if band_id[-2:] not in stripe_envs:
+            compute_env = stripe_envs["in"]
         else:
-            queue.put(MultiModule(modules))
+            compute_env = stripe_envs[band_id[-2:]]
+        module_list = []
+        for listed_module in modules:
+            listed_module.env_ = compute_env
+            module_list.append(listed_module)
+        queue.put(MultiModule(module_list))
     queue.wait()
 
     # Update map support information

From d0e8e8082062c65ee508bb2b58882b519d894738 Mon Sep 17 00:00:00 2001
From: Stefan Blumentrath 
Date: Thu, 11 Apr 2024 12:41:33 +0200
Subject: [PATCH 11/27] i.sentinel3.import: do not import all anxillary bands

---
 .../i.sentinel/i.sentinel3.import/i.sentinel3.import.py     | 6 ++----
 1 file changed, 2 insertions(+), 4 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 6bd0cef8d6..f33129da95 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -50,17 +50,15 @@
 # %option
 # % key: anxillary_bands
 # % multiple: yes
-# % answer: all
 # % required: no
-# % description: Anxillary data bands to import (e.g. LST_uncertainty, default is all available)
+# % description: Anxillary data bands to import (e.g. LST_uncertainty, default is None, use "all" to import all available)
 # %end
 
 # %option
 # % key: flag_bands
 # % multiple: yes
-# % answer: all
 # % required: no
-# % description: Quality flag bands to import (e.g. bayes_in, default is all available)
+# % description: Quality flag bands to import (e.g. bayes_in, default is None, use "all" to import all available)
 # %end
 
 # %option

From 99ed69a2bac0fc1300759b6909f474ea6d847b29 Mon Sep 17 00:00:00 2001
From: Stefan Blumentrath 
Date: Thu, 25 Apr 2024 14:10:03 +0200
Subject: [PATCH 12/27] Update i.sentinel3.import.py

---
 .../i.sentinel3.import/i.sentinel3.import.py     | 16 ++++++++++------
 1 file changed, 10 insertions(+), 6 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index f33129da95..876d4bc02c 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -1742,14 +1742,19 @@ def main():
         queue.put(MultiModule(module_list))
     queue.wait()
 
-    # Update map support information
+    # Update map support information and filter out empty maps
+    # Empty maps may occure because of a lack of valid data within the comutational region
     t_register_strings = []
     queue = ParallelModuleQueue(nprocs)
-    gs.verbose(_("Writing metadata to maps..."))
+    gs.verbose(_("Writing metadata to maps and filtering empty maps..."))
     for mapname, metadata in consolidate_metadata_dicts(register_strings).items():
         mapname = (
             mapname if not flags["r"] else mapname.replace("radiance", "reflectance")
         )
+        if not gs.raster_info(mapname)["max"]:
+            gs.warning(("Raster map {} is empty. Removing...").format(mapname))
+            gs.run_command("g.remove", type="raster", name=mapname, flags="f", quiet=True)
+            continue
         metadata["semantic_label"] = (
             metadata["semantic_label"]
             if not flags["r"]
@@ -1832,14 +1837,13 @@ def main():
 
     # Write t.register file if requested
     if options["register_output"]:
+        t_register_strings = "\n".join(t_register_strings) if t_register_strings else ""
         gs.verbose(
             _("Writing register file <{}>...").format(options["register_output"])
         )
-        Path(options["register_output"]).write_text(
-            "\n".join(t_register_strings) + "\n", encoding="UTF8"
-        )
+        Path(options["register_output"]).write_text(t_register_strings + "\n", encoding="UTF8")
     else:
-        print("\n".join(t_register_strings))
+        print(t_register_strings)
     return 0
 
 

From 42644a5180c62e6f79432f434733ee06577285cb Mon Sep 17 00:00:00 2001
From: Stefan Blumentrath 
Date: Fri, 26 Apr 2024 13:41:14 +0200
Subject: [PATCH 13/27] Update i.sentinel3.import.py

---
 src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 876d4bc02c..2a50116bff 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -1836,8 +1836,8 @@ def main():
     queue.wait()
 
     # Write t.register file if requested
+    t_register_strings = "\n".join(t_register_strings) if t_register_strings else ""
     if options["register_output"]:
-        t_register_strings = "\n".join(t_register_strings) if t_register_strings else ""
         gs.verbose(
             _("Writing register file <{}>...").format(options["register_output"])
         )

From 5dc24a8fa84b9aae5fd6f03d71d06b145a604d6f Mon Sep 17 00:00:00 2001
From: Stefan Blumentrath 
Date: Fri, 26 Apr 2024 13:47:44 +0200
Subject: [PATCH 14/27] typo

---
 src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 2a50116bff..d08b5b5172 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -4,7 +4,7 @@
  MODULE:      i.sentinel3.import
  AUTHOR(S):   Stefan Blumentrath
  PURPOSE:     Import and pre-process Sentinel-3 data from the Copernicus program
- COPYRIGHT:   (C) 2024 by Norwegian Water and ENergy Directorate, Stefan Blumentrath,
+ COPYRIGHT:   (C) 2024 by Norwegian Water and Energy Directorate, Stefan Blumentrath,
               and the GRASS development team
 
               This program is free software under the GNU General Public

From ca65f2902328f0be530b7efa01a553b3728daab1 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Wed, 15 May 2024 11:05:56 +0200
Subject: [PATCH 15/27] more fixes for metadata

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 50 ++++++++++++++-----
 1 file changed, 38 insertions(+), 12 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index d08b5b5172..d55413255a 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -919,7 +919,7 @@ def get_file_metadata(nc_dataset):
 
 
 def get_band_metadata(
-    band, nc_variable, fmt, file_metadata=None, basename=None, module_flags=None
+    band_tuple, nc_variable, fmt, file_metadata=None, basename=None, to_celsius=False
 ):
     """Extract band metadata from NetCDF variable"""
     metadata = file_metadata.copy()
@@ -927,20 +927,20 @@ def get_band_metadata(
     band_attrs = nc_variable.ncattrs()
 
     # Define variable name
-    varname_short = band.full_name
+    varname_short = band_tuple[1]
     datatype = str(nc_variable[:].dtype)
 
     # Define map name
     mapname = f"{basename}_{varname_short}"
     metadata["mapname"] = mapname
 
-    # band_title = nc_variable.long_name if "long_name" in band_attrs else band.full_name
+    # band_title = nc_variable.long_name if "long_name" in band_attrs else band_tuple[1]
 
     # Define unit
     unit = nc_variable.units if "units" in band_attrs else None
     unit = (
         "degree_celsius"
-        if band.band_id.startswith("LST") and module_flags["c"]
+        if band_tuple[0].startswith("LST") and to_celsius
         else unit
     )
     metadata["unit"] = unit
@@ -949,7 +949,7 @@ def get_band_metadata(
     if datatype not in DTYPE_TO_GRASS:
         gs.fatal(
             _("Unsupported datatype {dt} in band {band}").format(
-                dt=datatype, band=band.full_name
+                dt=datatype, band=band_tuple[1]
             )
         )
     if datatype in ["uint8", "uint16"]:
@@ -1081,9 +1081,10 @@ def import_s3(s3_file, kwargs, s3_product=None):
         sun_region_bounds = region_bounds.copy()
         if mod_flags["r"] or kwargs["maximum_solar_angle"]:
             tmp_ascii_sun = TMP_FILE / f"{rmap}_sun_parameters.txt"
-            module_queue, sun_region_dict = s3_product.get_sun_parameters(
+            module_queue, register_output, sun_region_dict = s3_product.get_sun_parameters(
                 zip_file, root, rmap, tmp_ascii_sun, region_bounds, maximum_solar_angle=kwargs["maximum_solar_angle"]
             )
+            register_strings.append(register_output)
             if module_queue is None:
                 return None, None, None
             region_dicts["sun_parameters"] = sun_region_dict
@@ -1296,12 +1297,12 @@ def process_nc_file(
 
                 # Collect band metadata
                 metadata = get_band_metadata(
-                    band,
+                    (band.band_id, band.full_name),
                     nc_variable,
                     fmt,
                     file_metadata=file_metadata,
                     basename=prefix,
-                    module_flags=module_flags,
+                    to_celsius=module_flags["c"],
                 )
                 fmt = metadata[1]
 
@@ -1387,6 +1388,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
         # nc_file = sun_azimuth_dict["nc_file"]
         # Get values
         import_modules = {}
+        meta_information = {}
         fmt = "%.12f,%.12f"
         member = str(
             root
@@ -1397,6 +1399,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
         nc_file_path = zip_file.extract(member, path=TMP_FILE)
         # sun_region_bounds = region_bounds.copy()
         with Dataset(nc_file_path) as sun_parameter_nc:
+            file_metadata = get_file_metadata(sun_parameter_nc)
             if maximum_solar_angle:
                 sun_mask = sun_parameter_nc[S3_SUN_PARAMTERS[self.product_type]["sun_bands"]["bands"]["solar_zenith"].format(self.view)][:] >= float(maximum_solar_angle)
             nc_bands, mask, resolution = get_geocoding(
@@ -1407,14 +1410,24 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
                 region_bounds=region_bounds,
             )
             if nc_bands is None:
-                return None, None
+                return None, None, None
 
             # sun_metadata = get_file_metadata(sun_parameter_nc)
             for band_id, band in S3_SUN_PARAMTERS[self.product_type]["sun_bands"][
                 "bands"
             ].items():
+                band = band.format(self.view)
+                sun_parameter_array = sun_parameter_nc[band]
+
+                metadata = get_band_metadata(
+                    (band_id, band),
+                    sun_parameter_array,
+                    fmt,
+                    file_metadata=file_metadata,
+                    basename=prefix,
+                )
+
                 fmt += ",%.12f"
-                sun_parameter_array = sun_parameter_nc[band.format(self.view)]
                 nc_bands[band_id] = sun_parameter_array[:]
                 map_name = f"{prefix}_{band_id}"
                 if band_id == "solar_zenith":
@@ -1431,6 +1444,20 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
                     method="mean",
                     solar_flux=None,
                 )
+                meta_information[map_name] = {
+                    "semantic_label": f"S3_{band_id}",
+                    "unit": metadata[0]["unit"],
+                    **{
+                        a: np_as_scalar(sun_parameter_nc.getncattr(a))
+                        for a in sun_parameter_nc.ncattrs()
+                    },
+                    **{"variable": band},
+                    **{
+                        a: np_as_scalar(sun_parameter_nc[band].getncattr(a))
+                        for a in sun_parameter_nc[band].ncattrs()
+                    },
+                }
+
 
             # Write to temporary file
             sun_region_dict = write_xyz(tmp_ascii, nc_bands, mask, fmt=fmt, project=True)
@@ -1438,7 +1465,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
                 resolution[0]
             ), float(resolution[1])
 
-        return import_modules, sun_region_dict
+        return import_modules, meta_information, sun_region_dict
 
 
 class S3Band:
@@ -1570,7 +1597,6 @@ def adjust_radiance(self, np_band_array):
 def get_solar_angle_bounds(region_bounds, sun_region_dict):
     """"""
     solar_bounds = gs.parse_command("g.region", flags="ugl", **sun_region_dict)
-    print(solar_bounds)
     solar_bounds["ll_n"] = max(float(solar_bounds["nw_lat"]), float(solar_bounds["ne_lat"]))
     solar_bounds["ll_s"] = min(float(solar_bounds["sw_lat"]), float(solar_bounds["se_lat"]))
     if solar_bounds["ll_n"] < float(region_bounds["ll_n"]):

From 7a64aa11d03d6ae309bda80db420978651053723 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Fri, 26 Jul 2024 14:23:26 +0200
Subject: [PATCH 16/27] handle meridian wrapping

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 167 +++++++++++++-----
 1 file changed, 126 insertions(+), 41 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index d55413255a..577a5ad87b 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -643,20 +643,49 @@ def intersect_region(region_current, region_stripe, align_current=True):
     starting from the sw-corner
     """
     relevant_keys = ["n", "s", "e", "w", "nsres", "ewres"]
-    if align_current:
-        region_dict = {"nsres": region_current["nsres"], "ewres": region_current["ewres"]}
-    else:
-        region_dict = {"nsres": region_stripe["nsres"], "ewres": region_stripe["ewres"]}
-    region_current = {region_key: float(region_current[region_key]) for region_key in relevant_keys}
-    region_stripe = {region_key: float(region_stripe[region_key]) for region_key in relevant_keys}
+    region_current = {
+        region_key: float(region_current[region_key]) for region_key in relevant_keys
+    }
+    region_stripe = {
+        region_key: float(region_stripe[region_key]) for region_key in relevant_keys
+    }
+    region_dict = region_current.copy()
+    if not align_current:
+        region_dict.update(
+            {"nsres": region_stripe["nsres"], "ewres": region_stripe["ewres"]}
+        )
     stripe_s = region_stripe["s"] - region_stripe["nsres"] / 2.0
     stripe_n = region_stripe["n"] + region_stripe["nsres"] / 2.0
     stripe_w = region_stripe["w"] - region_stripe["ewres"] / 2.0
     stripe_e = region_stripe["e"] + region_stripe["ewres"] / 2.0
-    region_dict["s"] = region_current["s"] if region_current["s"] >= stripe_s else region_current["s"] + np.ceil((region_current["s"] - stripe_s) / float(region_dict["nsres"])) * float(region_dict["nsres"])
-    region_dict["w"] = region_current["w"] if region_current["w"] >= stripe_w else region_current["w"] + np.ceil((region_current["w"] - stripe_w) / float(region_dict["ewres"])) * float(region_dict["ewres"])
-    region_dict["n"] = region_current["n"] if region_current["n"] <= stripe_n else region_current["n"] - np.floor((region_current["n"] - stripe_n) / float(region_dict["nsres"])) * float(region_dict["nsres"])
-    region_dict["e"] = region_current["e"] if region_current["e"] <= stripe_e else region_current["e"] - np.floor((region_current["e"] - stripe_e) / float(region_dict["ewres"])) * float(region_dict["ewres"])
+    region_dict["s"] = max(
+        region_current["s"],
+        region_current["s"]
+        + np.ceil((region_current["s"] - stripe_s) / float(region_dict["nsres"]))
+        * float(region_dict["nsres"]),
+    )
+    region_dict["n"] = min(
+        region_current["n"],
+        region_current["n"]
+        - np.floor((region_current["n"] - stripe_n) / float(region_dict["nsres"]))
+        * float(region_dict["nsres"]),
+    )
+    if stripe_e > stripe_w:
+        new_w = max(
+            region_current["w"],
+            region_current["w"]
+            + np.ceil((region_current["w"] - stripe_w) / float(region_dict["ewres"]))
+            * float(region_dict["ewres"]),
+        )
+        new_e = min(
+            region_current["e"],
+            region_current["e"]
+            - np.floor((region_current["e"] - stripe_e) / float(region_dict["ewres"]))
+            * float(region_dict["ewres"]),
+        )
+        if new_w < new_e:
+            region_dict["w"] = new_w
+            region_dict["e"] = new_e
     return region_dict
 
 
@@ -796,7 +825,11 @@ def get_geocoding(zip_file, root, geo_bands_dict, sun_mask=None, region_bounds=N
             )
 
             if mask.all():
-                gs.warning(_("No valid pixels inside computational region found in {}").format(str(zip_file.filename)))
+                gs.warning(
+                    _("No valid pixels inside computational region found in {}").format(
+                        str(zip_file.filename)
+                    )
+                )
                 return None, None, None
 
     return nc_bands, mask, resolution
@@ -938,11 +971,7 @@ def get_band_metadata(
 
     # Define unit
     unit = nc_variable.units if "units" in band_attrs else None
-    unit = (
-        "degree_celsius"
-        if band_tuple[0].startswith("LST") and to_celsius
-        else unit
-    )
+    unit = "degree_celsius" if band_tuple[0].startswith("LST") and to_celsius else unit
     metadata["unit"] = unit
 
     # Define datatype and import method
@@ -982,12 +1011,18 @@ def get_band_metadata(
         fill_val = nc_variable._FillValue
         if fill_val > 0:
             min_val = data_range["min"] if min_val is None else min_val
-            max_val = fill_val - 1 if max_val is None else max_val  # data_range["resolution"]
+            max_val = (
+                fill_val - 1 if max_val is None else max_val
+            )  # data_range["resolution"]
         elif fill_val <= 0:
-            min_val = fill_val + 1 if min_val is None else min_val  # data_range["resolution"]
+            min_val = (
+                fill_val + 1 if min_val is None else min_val
+            )  # data_range["resolution"]
             max_val = data_range["max"] if max_val is None else max_val
     else:
-        min_val = data_range["min"] + 1 if min_val is None else min_val  # data_range["resolution"]
+        min_val = (
+            data_range["min"] + 1 if min_val is None else min_val
+        )  # data_range["resolution"]
         max_val = data_range["max"] if max_val is None else max_val
 
     if "flag_masks" in band_attrs:
@@ -1081,8 +1116,17 @@ def import_s3(s3_file, kwargs, s3_product=None):
         sun_region_bounds = region_bounds.copy()
         if mod_flags["r"] or kwargs["maximum_solar_angle"]:
             tmp_ascii_sun = TMP_FILE / f"{rmap}_sun_parameters.txt"
-            module_queue, register_output, sun_region_dict = s3_product.get_sun_parameters(
-                zip_file, root, rmap, tmp_ascii_sun, region_bounds, maximum_solar_angle=kwargs["maximum_solar_angle"]
+            (
+                module_queue,
+                register_output,
+                sun_region_dict,
+            ) = s3_product.get_sun_parameters(
+                zip_file,
+                root,
+                rmap,
+                tmp_ascii_sun,
+                region_bounds,
+                maximum_solar_angle=kwargs["maximum_solar_angle"],
             )
             register_strings.append(register_output)
             if module_queue is None:
@@ -1090,13 +1134,25 @@ def import_s3(s3_file, kwargs, s3_product=None):
             region_dicts["sun_parameters"] = sun_region_dict
 
             if kwargs["maximum_solar_angle"]:
-                sun_region_bounds = gs.parse_command("g.region", flags="ugb", quiet=True, **intersect_region(dict(kwargs["current_reg"]), sun_region_dict, align_current=False))
+                sun_region_bounds = gs.parse_command(
+                    "g.region",
+                    flags="ugb",
+                    quiet=True,
+                    **intersect_region(
+                        dict(kwargs["current_reg"]),
+                        sun_region_dict,
+                        align_current=False,
+                    ),
+                )
         for stripe, stripe_dict in s3_product.requested_stripe_content.items():
             # Could be parallelized!!!
             # Get geocoding (lat/lon, elevation) and initial mask
             tmp_ascii = TMP_FILE / f"{rmap}_{stripe}.txt"
             nc_bands, mask, resolution = get_geocoding(
-                zip_file, root, stripe_dict["geocoding"], region_bounds=sun_region_bounds
+                zip_file,
+                root,
+                stripe_dict["geocoding"],
+                region_bounds=sun_region_bounds,
             )
             if nc_bands is None:
                 return None, None, None
@@ -1124,9 +1180,13 @@ def import_s3(s3_file, kwargs, s3_product=None):
             )
             stripe_region_dict["ewres"] = float(resolution[0])
             stripe_region_dict["nsres"] = float(resolution[1])
-            sun_region_dict_intersect = intersect_region(dict(kwargs["current_reg"]), sun_region_dict, align_current=False)
+            sun_region_dict_intersect = intersect_region(
+                dict(kwargs["current_reg"]), sun_region_dict, align_current=False
+            )
             if kwargs["maximum_solar_angle"]:
-                region_dicts[stripe] = intersect_region(sun_region_dict_intersect, stripe_region_dict, align_current=False)
+                region_dicts[stripe] = intersect_region(
+                    sun_region_dict_intersect, stripe_region_dict, align_current=False
+                )
             else:
                 region_dicts[stripe] = stripe_region_dict
 
@@ -1324,7 +1384,7 @@ def process_nc_file(
                 # Rescale temperature variables if requested
                 if band.full_name == "LST" and module_flags["c"]:
                     add_array = convert_units(add_array, "K", "degree_celsius")
- 
+
                 # Write metadata json if requested
                 band_attrs = nc_file_open[band.full_name].ncattrs()
 
@@ -1382,7 +1442,9 @@ def process_nc_file(
 
             return nc_bands, module_queue, meta_information, fmt
 
-    def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, maximum_solar_angle=None):
+    def get_sun_parameters(
+        self, zip_file, root, prefix, tmp_ascii, region_bounds, maximum_solar_angle=None
+    ):
         """https://github.com/sertit/eoreader/blob/main/eoreader/products/optical/s3_slstr_product.py#L862"""
 
         # nc_file = sun_azimuth_dict["nc_file"]
@@ -1401,7 +1463,11 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
         with Dataset(nc_file_path) as sun_parameter_nc:
             file_metadata = get_file_metadata(sun_parameter_nc)
             if maximum_solar_angle:
-                sun_mask = sun_parameter_nc[S3_SUN_PARAMTERS[self.product_type]["sun_bands"]["bands"]["solar_zenith"].format(self.view)][:] >= float(maximum_solar_angle)
+                sun_mask = sun_parameter_nc[
+                    S3_SUN_PARAMTERS[self.product_type]["sun_bands"]["bands"][
+                        "solar_zenith"
+                    ].format(self.view)
+                ][:] >= float(maximum_solar_angle)
             nc_bands, mask, resolution = get_geocoding(
                 zip_file,
                 root,
@@ -1438,7 +1504,7 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
                     map_name,
                     distance=3,
                     fill_flags="ks",
-                    zrange=[-32767,32768],
+                    zrange=[-32767, 32768],
                     val_col=len(nc_bands),
                     data_type="FCELL",
                     method="mean",
@@ -1458,9 +1524,10 @@ def get_sun_parameters(self, zip_file, root, prefix, tmp_ascii, region_bounds, m
                     },
                 }
 
-
             # Write to temporary file
-            sun_region_dict = write_xyz(tmp_ascii, nc_bands, mask, fmt=fmt, project=True)
+            sun_region_dict = write_xyz(
+                tmp_ascii, nc_bands, mask, fmt=fmt, project=True
+            )
             sun_region_dict["ewres"], sun_region_dict["nsres"] = float(
                 resolution[0]
             ), float(resolution[1])
@@ -1547,7 +1614,10 @@ def _get_exception_band(self, product_type):
         return None
 
     def _get_solar_flux_dict_for_band(self, product_type):
-        if S3_SOLAR_FLUX[product_type] and self.band_id in S3_SOLAR_FLUX[product_type]["defaults"]:
+        if (
+            S3_SOLAR_FLUX[product_type]
+            and self.band_id in S3_SOLAR_FLUX[product_type]["defaults"]
+        ):
             return {
                 "default": S3_SOLAR_FLUX[product_type]["defaults"][self.band_id],
                 "band": S3_SOLAR_FLUX[product_type]["band"].format(
@@ -1597,8 +1667,12 @@ def adjust_radiance(self, np_band_array):
 def get_solar_angle_bounds(region_bounds, sun_region_dict):
     """"""
     solar_bounds = gs.parse_command("g.region", flags="ugl", **sun_region_dict)
-    solar_bounds["ll_n"] = max(float(solar_bounds["nw_lat"]), float(solar_bounds["ne_lat"]))
-    solar_bounds["ll_s"] = min(float(solar_bounds["sw_lat"]), float(solar_bounds["se_lat"]))
+    solar_bounds["ll_n"] = max(
+        float(solar_bounds["nw_lat"]), float(solar_bounds["ne_lat"])
+    )
+    solar_bounds["ll_s"] = min(
+        float(solar_bounds["sw_lat"]), float(solar_bounds["se_lat"])
+    )
     if solar_bounds["ll_n"] < float(region_bounds["ll_n"]):
         region_bounds["ll_n"] = solar_bounds["ll_n"]
     if float(solar_bounds["ll_s"]) > float(region_bounds["ll_s"]):
@@ -1671,7 +1745,9 @@ def main():
         "current_reg": dict(current_region),
         "mod_flags": flags,
         "meta_dict": meta_info_dict,
-        "maximum_solar_angle": float(options["maximum_solar_angle"]) if options["maximum_solar_angle"] else None,
+        "maximum_solar_angle": float(options["maximum_solar_angle"])
+        if options["maximum_solar_angle"]
+        else None,
     }
 
     # if flags["p"]:
@@ -1726,7 +1802,13 @@ def main():
         if flags["n"]:
             stripe_env["GRASS_REGION"] = gs.region_env(**adjust_region(stripe_region))
         else:
-            stripe_env["GRASS_REGION"] = gs.region_env(**intersect_region(dict(current_region), stripe_region, align_current=stripe_id != "sun_parameters"))
+            stripe_env["GRASS_REGION"] = gs.region_env(
+                **intersect_region(
+                    dict(current_region),
+                    stripe_region,
+                    align_current=stripe_id != "sun_parameters",
+                )
+            )
         stripe_envs[stripe_id] = stripe_env
 
     if flags["r"]:
@@ -1742,11 +1824,10 @@ def main():
             queue.put(MultiModule(module_list))
         queue.wait()
         # Here one could zoom to non-null pixels in solar angle map
-        # g.region zoom=zolar_zenith 
-        # if maximum_solar_angle is given in consistently limit 
+        # g.region zoom=zolar_zenith
+        # if maximum_solar_angle is given in consistently limit
         # imported pixels
 
-
     gs.verbose(
         _("Importing scenes {}").format(
             "\n" + "\n".join([s3_file.name for s3_file in s3_files])
@@ -1779,7 +1860,9 @@ def main():
         )
         if not gs.raster_info(mapname)["max"]:
             gs.warning(("Raster map {} is empty. Removing...").format(mapname))
-            gs.run_command("g.remove", type="raster", name=mapname, flags="f", quiet=True)
+            gs.run_command(
+                "g.remove", type="raster", name=mapname, flags="f", quiet=True
+            )
             continue
         metadata["semantic_label"] = (
             metadata["semantic_label"]
@@ -1867,7 +1950,9 @@ def main():
         gs.verbose(
             _("Writing register file <{}>...").format(options["register_output"])
         )
-        Path(options["register_output"]).write_text(t_register_strings + "\n", encoding="UTF8")
+        Path(options["register_output"]).write_text(
+            t_register_strings + "\n", encoding="UTF8"
+        )
     else:
         print(t_register_strings)
     return 0

From f70b32bcda5a0e52c9f722aad3300b53904b36d7 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Tue, 18 Mar 2025 10:27:28 +0100
Subject: [PATCH 17/27] fix region intersection

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 68 +++++++++----------
 1 file changed, 32 insertions(+), 36 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 577a5ad87b..93748992ba 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -635,13 +635,16 @@ def adjust_region(region_dict):
     return region_dict
 
 
-def intersect_region(region_current, region_stripe, align_current=True):
-    """
+def intersect_region(region_current, region_stripe, align_current=True) -> dict:
+    """Intersect stripe region with current region aligned to resolution and extended for r.in.xyz.
+
     Adjust region bounds for r.in.xyz to the intersection
     of the current region (bounds + half resolution)
     aligned to resolution in ew and ns direction
     starting from the sw-corner
     """
+
+    # Remove potential, irrelevant / conflicting dict-keys and convert to float
     relevant_keys = ["n", "s", "e", "w", "nsres", "ewres"]
     region_current = {
         region_key: float(region_current[region_key]) for region_key in relevant_keys
@@ -650,49 +653,40 @@ def intersect_region(region_current, region_stripe, align_current=True):
         region_key: float(region_stripe[region_key]) for region_key in relevant_keys
     }
     region_dict = region_current.copy()
+    # Choose pixel alignment with stripe or region
     if not align_current:
         region_dict.update(
             {"nsres": region_stripe["nsres"], "ewres": region_stripe["ewres"]}
         )
-    stripe_s = region_stripe["s"] - region_stripe["nsres"] / 2.0
-    stripe_n = region_stripe["n"] + region_stripe["nsres"] / 2.0
-    stripe_w = region_stripe["w"] - region_stripe["ewres"] / 2.0
-    stripe_e = region_stripe["e"] + region_stripe["ewres"] / 2.0
-    region_dict["s"] = max(
-        region_current["s"],
-        region_current["s"]
-        + np.ceil((region_current["s"] - stripe_s) / float(region_dict["nsres"]))
-        * float(region_dict["nsres"]),
-    )
-    region_dict["n"] = min(
-        region_current["n"],
-        region_current["n"]
-        - np.floor((region_current["n"] - stripe_n) / float(region_dict["nsres"]))
-        * float(region_dict["nsres"]),
-    )
-    if stripe_e > stripe_w:
-        new_w = max(
-            region_current["w"],
-            region_current["w"]
-            + np.ceil((region_current["w"] - stripe_w) / float(region_dict["ewres"]))
-            * float(region_dict["ewres"]),
-        )
-        new_e = min(
-            region_current["e"],
-            region_current["e"]
-            - np.floor((region_current["e"] - stripe_e) / float(region_dict["ewres"]))
-            * float(region_dict["ewres"]),
-        )
-        if new_w < new_e:
-            region_dict["w"] = new_w
-            region_dict["e"] = new_e
+    for direction in "nsew":
+        resolution = region_dict["nsres"] if direction in "ns" else region_dict["ewres"]
+        if direction in "ne":
+            # Get the intersection of the stripe and the current region
+            region_dict[direction] = min(
+                region_current[direction], region_stripe[direction]
+            )
+            # Align region to required resolution and extend by one pixel for r.in.xyz
+            region_dict[direction] = (
+                np.ceil(region_dict[direction] / resolution) * resolution + resolution
+            )
+        else:
+            region_dict[direction] = max(
+                region_current[direction], region_stripe[direction]
+            )
+            region_dict[direction] = (
+                np.floor(region_dict[direction] / resolution) * resolution - resolution
+            )
+
     return region_dict
 
 
-def parse_s3_file_name(file_name):
-    """Extract info from file name according to naming onvention:
+def parse_s3_file_name(file_name: str) -> dict:
+    """Extract info from file name according to ESA naming onvention.
+
+    Naming convention is documented here:
     https://sentinels.copernicus.eu/web/sentinel/user-guides/sentinel-3-slstr/naming-convention
     Assumes that file name is checked to be a valid / supported Sentinel-3 file name
+
     :param file_name: string representing the file name of a Senintel-3 scene
     """
     try:
@@ -1972,6 +1966,8 @@ def main():
         )
     try:
         from osgeo import osr
+
+        osr.UseExceptions()
     except ImportError:
         gs.fatal(
             _("Could not import gdal. Please install it with:\n" "'pip install GDAL'!")

From 0cfc36572877dad8d2eb4083bd0488ca0ad1b045 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Tue, 18 Mar 2025 11:23:58 +0100
Subject: [PATCH 18/27] some ruff fixes

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 161 ++++++++++--------
 1 file changed, 91 insertions(+), 70 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 93748992ba..b7f1a0327c 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -485,8 +485,8 @@
 SUN_ZENITH_ANGLE = None
 
 
-def cleanup():
-    """Remove all temporary data"""
+def cleanup() -> None:
+    """Remove all temporary data."""
     # remove temporary maps
     if TMP_FILE:
         gs.try_remove(TMP_FILE)
@@ -506,7 +506,7 @@ def cleanup():
 
 
 def np_as_scalar(var):
-    """Return a numpy object as scalar"""
+    """Return a numpy object as scalar."""
     if type(var).__module__ == np.__name__:
         if var.size > 1:
             return str(var)
@@ -514,9 +514,10 @@ def np_as_scalar(var):
     return var
 
 
-def get_dtype_range(np_datatype):
-    """Get information to set the valid data range based on dtype
-    according to:
+def get_dtype_range(np_datatype: str) -> dict:
+    """Get information to set the valid data range based on dtype.
+
+    dtype range according to:
     https://docs.unidata.ucar.edu/netcdf-c/current/attribute_conventions.html#valid_range
     """
     if np_datatype == "bool":
@@ -541,8 +542,11 @@ def get_dtype_range(np_datatype):
 
 
 def consolidate_metadata_dicts(metadata_dicts):
-    """Consolidate a list of metadata dictionaries for all
-    input scenes into unified metadata for the resulting map"""
+    """Consolidate a list of metadata dictionaries for all input scenes.
+
+    Consolidate a list of metadata dictionaries for all
+    input scenes into unified metadata for the resulting map
+    """
     map_metadata_dicts = {}
     for map_dict in metadata_dicts:
         for map_name, meta_dict in map_dict.items():
@@ -565,15 +569,18 @@ def consolidate_metadata_dicts(metadata_dicts):
     return map_metadata_dicts
 
 
-def write_metadata(json_dict, metadatajson):
-    """Write extended map metadata to JSON file"""
+def write_metadata(json_dict: dict, metadatajson: str) -> None:
+    """Write extended map metadata to JSON file."""
     with open(metadatajson, "w", encoding="UTF8") as outfile:
         json.dump(json_dict, outfile)
 
 
-def convert_units(np_column, from_u, to_u):
-    """Converts a numpy column from one unit to
-    another, convertibility needs to be checked beforehand"""
+def convert_units(np_column, from_u: str, to_u: str):
+    """Convert units of a numpy array.
+
+    Converts a numpy column from one unit to
+    another, convertibility needs to be checked beforehand.
+    """
     try:
         from cf_units import Unit
     except ImportError:
@@ -597,8 +604,8 @@ def convert_units(np_column, from_u, to_u):
     return converted_col
 
 
-def extend_region(region_dict, additional_region_dict):
-    """Extend region bounds in region_dict to cover also the additional_region_dict"""
+def extend_region(region_dict: dict, additional_region_dict: dict) -> dict:
+    """Extend region bounds in region_dict to cover also the additional_region_dict."""
     region_dict["n"] = max(region_dict["n"], additional_region_dict["n"])
     region_dict["e"] = max(region_dict["e"], additional_region_dict["e"])
     region_dict["s"] = min(region_dict["s"], additional_region_dict["s"])
@@ -606,8 +613,9 @@ def extend_region(region_dict, additional_region_dict):
     return region_dict
 
 
-def adjust_region(region_dict):
-    """
+def adjust_region(region_dict: dict) -> dict:
+    """Adjust region bounds for r.in.xyz.
+
     Adjust region bounds for r.in.xyz (bounds + half resolution)
     aligned to resolution in ew and ns direction
     starting from the sw-corner
@@ -635,7 +643,9 @@ def adjust_region(region_dict):
     return region_dict
 
 
-def intersect_region(region_current, region_stripe, align_current=True) -> dict:
+def intersect_region(
+    region_current: dict, region_stripe: dict, *, align_current: bool = True
+) -> dict:
     """Intersect stripe region with current region aligned to resolution and extended for r.in.xyz.
 
     Adjust region bounds for r.in.xyz to the intersection
@@ -711,8 +721,8 @@ def parse_s3_file_name(file_name: str) -> dict:
         gs.fatal(_("{} is not a supported Sentinel-3 scene").format(str(file_name)))
 
 
-def extract_file_info(s3_files, basename=None):
-    """Extract information from file name according to naming conventions"""
+def extract_file_info(s3_files: list, basename: str = None) -> tuple(str, dict):
+    """Extract information from file name according to naming conventions."""
     result_dict = {}
     product_track_ids = [
         "cycle",
@@ -771,7 +781,7 @@ def extract_file_info(s3_files, basename=None):
 
 
 def get_geocoding(zip_file, root, geo_bands_dict, sun_mask=None, region_bounds=None):
-    """Get ground control points from NetCDF file"""
+    """Get ground control points from NetCDF file."""
     member = str(root / geo_bands_dict["nc_file"])
     nc_file_path = zip_file.extract(member, path=TMP_FILE)
     with Dataset(nc_file_path) as nc_file_open:
@@ -841,7 +851,7 @@ def setup_import_multi_module(
     solar_flux=None,
     rules=None,
 ):
-    """Setup GRASS GIS moduls for importing S3 bands"""
+    """Setup GRASS GIS moduls for importing S3 bands."""
     # Basic import module
     modules = [
         Module(
@@ -915,8 +925,8 @@ def setup_import_multi_module(
     return modules
 
 
-def get_file_metadata(nc_dataset):
-    """Collect metadata from NetCDF file"""
+def get_file_metadata(nc_dataset) -> dict:
+    """Collect metadata from NetCDF file."""
     metadata = {
         attr: np_as_scalar(nc_dataset.getncattr(attr))
         for attr in [
@@ -948,7 +958,7 @@ def get_file_metadata(nc_dataset):
 def get_band_metadata(
     band_tuple, nc_variable, fmt, file_metadata=None, basename=None, to_celsius=False
 ):
-    """Extract band metadata from NetCDF variable"""
+    """Extract band metadata from NetCDF variable."""
     metadata = file_metadata.copy()
     # Collect metadata
     band_attrs = nc_variable.ncattrs()
@@ -975,7 +985,7 @@ def get_band_metadata(
                 dt=datatype, band=band_tuple[1]
             )
         )
-    if datatype in ["uint8", "uint16"]:
+    if datatype in {"uint8", "uint16"}:
         metadata["method"] = "max"  # Unfortunately there is no "mode" in r.in.xyz
         fmt += ",%i"
     else:
@@ -1031,8 +1041,11 @@ def get_band_metadata(
 
 
 def transform_coordinates(coordinates):
-    """Tranforms a numy array with coordinates to
-    projection of the current location"""
+    """Tranform coordinates to projection of the current location.
+
+    Tranforms a numpy array with coordinates to the
+    projection of the current location
+    """
 
     # Create source coordinate reference
     s_srs = osr.SpatialReference()
@@ -1053,7 +1066,7 @@ def transform_coordinates(coordinates):
 
 
 def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
-    """Write temporary XYZ ascii file"""
+    """Write temporary XYZ ascii file."""
     # Extract grid coordinates
     lon = np.ma.masked_where(mask, nc_bands["lon"][:]).compressed()
     lat = np.ma.masked_where(mask, nc_bands["lat"][:]).compressed()
@@ -1067,7 +1080,7 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
         np_output = np.hstack((lon[:, None], lat[:, None]))
     # Fetch, mask and stack requested bands
     for band in nc_bands:
-        if band in ["lat", "lon"]:
+        if band in {"lat", "lon"}:
             continue
         add_array = nc_bands[band][:]
         if np.ma.is_masked(add_array):
@@ -1083,7 +1096,7 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
     else:
         np.savetxt(tmp_ascii, np_output, delimiter=",", fmt=fmt)
 
-    stripe_reg = {  #
+    stripe_reg = {
         "n": np.ma.max(np_output[:, 1]),
         "s": np.ma.min(np_output[:, 1]),
         "e": np.ma.max(np_output[:, 0]),
@@ -1092,8 +1105,8 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
     return stripe_reg
 
 
-def import_s3(s3_file, kwargs, s3_product=None):
-    """Import Sentinel-3 netCDF4 data"""
+def import_s3(s3_file: str, kwargs: dict, s3_product: str = None):
+    """Import Sentinel-3 netCDF4 data."""
     # Unpack dictionary variables
     rmap = kwargs["meta_dict"][0]
     region_bounds = kwargs["reg_bounds"]
@@ -1188,12 +1201,16 @@ def import_s3(s3_file, kwargs, s3_product=None):
 
 
 class S3Product:
-    """Class to provide information necessary to pre-process Senintel-3 data products
-    level 1 and 2"""
+    """Class for Senintel-3 data products.
+
+    The class provides information necessary to pre-process Senintel-3 data products
+    level 1 and 2.
+    """
 
     def __init__(
         self, product_type, view="n", bands=None, flag_bands=None, anxillary_bands=None
     ):
+        """Initialize a Sentinel-3 product."""
         self.product_type = product_type
         self.available_views = S3_VIEWS[product_type]
         self.available_bands = list(S3_BANDS["bands"][product_type].keys())
@@ -1214,9 +1231,11 @@ def __init__(
         self.file_pattern = S3_FILE_PATTERN[product_type]
 
     def __str__(self):
+        """Return Sentinel-3 product information as JSON."""
         return json.dumps(self.__repr__(), indent=2)
 
     def __repr__(self):
+        """Return Sentinel-3 product information as JSON."""
         class_dict = self.__dict__.copy()
         for band_type in ["bands", "flag_bands", "anxillary_bands"]:
             bands_of_type = getattr(self, band_type)
@@ -1231,6 +1250,7 @@ def __repr__(self):
         return json.dumps(class_dict, indent=2)
 
     def _check_view(self, view):
+        """Check requested Satellite view to process."""
         if view not in self.available_views:
             gs.warning(
                 _("View {} not available for product type {}").format(
@@ -1240,6 +1260,7 @@ def _check_view(self, view):
         return view
 
     def _check_bands(self, bands, band_type="bands"):
+        """Check requested Satellite bands to process."""
         band_types = {
             "bands": self.available_bands,
             "flag_bands": self.available_flag_bands,
@@ -1265,6 +1286,7 @@ def _check_bands(self, bands, band_type="bands"):
         return band_objects
 
     def _collect_requested_stripe_content(self):
+        """Collect requested stripe content."""
         suffixes = {}
         for band_type in ["bands", "flag_bands", "anxillary_bands"]:
             selected_bands = getattr(self, band_type)
@@ -1311,7 +1333,7 @@ def process_nc_file(
         tmp_ascii,
         fmt=None,
     ):
-        """Extract requested bands as numpy arrays from NetCDF file and setup import modules"""
+        """Extract requested bands as numpy arrays from NetCDF file and setup import modules."""
         meta_information = {}
         member = str(root / container_dict_items[0])
         nc_file_path = zip_file.extract(member, path=TMP_FILE)
@@ -1364,13 +1386,16 @@ def process_nc_file(
                 metadata[0]["solar_flux"] = solar_flux
 
                 if band.exception:
+                    maximum_exception = 3
                     if np.ma.is_masked(add_array):
                         add_array.mask = np.ma.mask_or(
-                            add_array.mask, nc_file_open[band.exception][:] >= 3
+                            add_array.mask,
+                            nc_file_open[band.exception][:] >= maximum_exception,
                         )
                     else:
                         add_array = np.ma.masked_array(
-                            add_array, nc_file_open[band.exception][:] >= 3
+                            add_array,
+                            nc_file_open[band.exception][:] >= maximum_exception,
                         )
                 if np.ma.is_masked(add_array):
                     add_array = add_array.filled()
@@ -1439,9 +1464,11 @@ def process_nc_file(
     def get_sun_parameters(
         self, zip_file, root, prefix, tmp_ascii, region_bounds, maximum_solar_angle=None
     ):
-        """https://github.com/sertit/eoreader/blob/main/eoreader/products/optical/s3_slstr_product.py#L862"""
+        """Extract sun parameters from S3 SLSTR product.
 
-        # nc_file = sun_azimuth_dict["nc_file"]
+        Oriented on:
+        https://github.com/sertit/eoreader/blob/main/eoreader/products/optical/s3_slstr_product.py#L862
+        """
         # Get values
         import_modules = {}
         meta_information = {}
@@ -1453,7 +1480,6 @@ def get_sun_parameters(
             )
         )
         nc_file_path = zip_file.extract(member, path=TMP_FILE)
-        # sun_region_bounds = region_bounds.copy()
         with Dataset(nc_file_path) as sun_parameter_nc:
             file_metadata = get_file_metadata(sun_parameter_nc)
             if maximum_solar_angle:
@@ -1472,7 +1498,6 @@ def get_sun_parameters(
             if nc_bands is None:
                 return None, None, None
 
-            # sun_metadata = get_file_metadata(sun_parameter_nc)
             for band_id, band in S3_SUN_PARAMTERS[self.product_type]["sun_bands"][
                 "bands"
             ].items():
@@ -1530,7 +1555,7 @@ def get_sun_parameters(
 
 
 class S3Band:
-    """Class for properties and methods related to Sentinel-3 bands"""
+    """Class for properties and methods related to Sentinel-3 bands."""
 
     def __init__(
         self,
@@ -1623,9 +1648,9 @@ def _get_solar_flux_dict_for_band(self, product_type):
             }
         return None
 
-    def get_solar_flux(self, zip_file, root_path):
-        """
-        Get solar spectral flux in mW / (m^2 * sr * nm) for band
+    def get_solar_flux(self, zip_file: ZipFile, root_path: Path) -> float:
+        """Get solar spectral flux in mW / (m^2 * sr * nm) for band.
+
         :returns: solar Flux
         :type: float
 
@@ -1640,8 +1665,8 @@ def get_solar_flux(self, zip_file, root_path):
             solar_flux = self.solar_flux["default"]
         return float(solar_flux)
 
-    def _get_band_geo_points(self, suffix):
-        """Get ground control point references from suffix"""
+    def _get_band_geo_points(self, suffix: str) -> dict:
+        """Get ground control point references from suffix."""
         if self.suffix.startswith("t"):
             return None
         return {
@@ -1652,14 +1677,14 @@ def _get_band_geo_points(self, suffix):
         }
 
     def adjust_radiance(self, np_band_array):
-        """Get radiance adjustment for band object"""
+        """Get radiance adjustment for band object."""
         if not self.radiance_adjustment:
             return np_band_array * self.radiance_adjustment
         return np_band_array
 
 
-def get_solar_angle_bounds(region_bounds, sun_region_dict):
-    """"""
+def get_solar_angle_bounds(region_bounds: dict, sun_region_dict: dict) -> dict:
+    """Get latlon coordinates of bounds in solar angle bands."""
     solar_bounds = gs.parse_command("g.region", flags="ugl", **sun_region_dict)
     solar_bounds["ll_n"] = max(
         float(solar_bounds["nw_lat"]), float(solar_bounds["ne_lat"])
@@ -1674,26 +1699,23 @@ def get_solar_angle_bounds(region_bounds, sun_region_dict):
     return region_bounds
 
 
-def main():
-    """Do the main work"""
+def main() -> None:
+    """Do the main work."""
     pattern = re.compile(
         ".*" + S3_FILE_PATTERN[options["product_type"]].replace("*", ".*")
     )
 
     # check provided input
     s3_files = options["input"].split(",")
-    if len(s3_files) == 1:
-        if re.match(pattern, s3_files[0]) is None:
-            try:
-                s3_files = (
-                    Path(s3_files[0]).read_text(encoding="UTF8").strip().split("\n")
-                )
-            except ValueError:
-                gs.fatal(
-                    _(
-                        "Input <{}> is neither a supported Sentinel-3 scene nor a text file with scenes"
-                    ).format(options["input"])
-                )
+    if len(s3_files) == 1 and re.match(pattern, s3_files[0]) is None:
+        try:
+            s3_files = Path(s3_files[0]).read_text(encoding="UTF8").strip().split("\n")
+        except ValueError:
+            gs.fatal(
+                _(
+                    "Input <{}> is neither a supported Sentinel-3 scene nor a text file with scenes"
+                ).format(options["input"])
+            )
 
     if not s3_files:
         gs.warning("No scenes found to process, please check input.")
@@ -1786,7 +1808,7 @@ def main():
                     if stripe_id not in region_dicts:
                         region_dicts[stripe_id] = region
                     else:
-                        extend_region(region_dicts[stripe_id], region_dict[stripe_id])
+                        extend_region(region_dicts[stripe_id], region)
     if not module_queues:
         gs.warning(_("Nothing to import with the given input"))
         sys.exit(0)
@@ -1809,11 +1831,10 @@ def main():
         gs.verbose(_("Importing solar parameter bands"))
         queue = ParallelModuleQueue(nprocs)
         compute_env = stripe_envs["sun_parameters"]
-        for solar_parameter in ["solar_azimuth", "solar_zenith"]:
+        for solar_parameter in ("solar_azimuth", "solar_zenith"):
             module_list = []
             for solar_module in module_queues[0][solar_parameter]:
                 solar_module.env_ = compute_env
-                # solar_module.verbose = True
                 module_list.append(solar_module)
             queue.put(MultiModule(module_list))
         queue.wait()
@@ -1829,7 +1850,7 @@ def main():
     )
     queue = ParallelModuleQueue(nprocs)
     for band_id, modules in module_queues[0].items():
-        if band_id in ["solar_azimuth", "solar_zenith"]:
+        if band_id in {"solar_azimuth", "solar_zenith"}:
             continue
         # S3SL2LST product has only "in" stripe and band names have no suffix
         if band_id[-2:] not in stripe_envs:
@@ -1970,7 +1991,7 @@ def main():
         osr.UseExceptions()
     except ImportError:
         gs.fatal(
-            _("Could not import gdal. Please install it with:\n" "'pip install GDAL'!")
+            _("Could not import gdal. Please install it with:\n'pip install GDAL'!")
         )
 
     try:

From fd10386b3ba9936a355d4d3d20ae6badf238cd6f Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Tue, 18 Mar 2025 12:51:18 +0100
Subject: [PATCH 19/27] some ruff fixes

---
 src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index b7f1a0327c..952b8043f3 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -721,7 +721,7 @@ def parse_s3_file_name(file_name: str) -> dict:
         gs.fatal(_("{} is not a supported Sentinel-3 scene").format(str(file_name)))
 
 
-def extract_file_info(s3_files: list, basename: str = None) -> tuple(str, dict):
+def extract_file_info(s3_files: list, basename: str = None):
     """Extract information from file name according to naming conventions."""
     result_dict = {}
     product_track_ids = [

From ea1b0bc330e7f8b36c0716d2d54aefc3043f5e9b Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Mon, 24 Mar 2025 08:34:34 +0100
Subject: [PATCH 20/27] update copyright

---
 .../i.sentinel/i.sentinel3.import/i.sentinel3.import.py       | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 952b8043f3..9e50d3904d 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -4,8 +4,8 @@
  MODULE:      i.sentinel3.import
  AUTHOR(S):   Stefan Blumentrath
  PURPOSE:     Import and pre-process Sentinel-3 data from the Copernicus program
- COPYRIGHT:   (C) 2024 by Norwegian Water and Energy Directorate, Stefan Blumentrath,
-              and the GRASS development team
+ COPYRIGHT:   (C) 2024-202-20255 BY Norwegian Water and Energy Directorate,
+              Stefan Blumentrath, and the GRASS development team
 
               This program is free software under the GNU General Public
               License (>=v2). Read the file COPYING that comes with GRASS

From 3d7bea2db3321b23917b98d7ef72337852554222 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Thu, 10 Apr 2025 11:29:38 +0200
Subject: [PATCH 21/27] handle out of bounds + ruff

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 502 ++++++++++--------
 1 file changed, 285 insertions(+), 217 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 9e50d3904d..8a6ace92e3 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -1,15 +1,14 @@
 #!/usr/bin/env python3
 
-"""
- MODULE:      i.sentinel3.import
- AUTHOR(S):   Stefan Blumentrath
- PURPOSE:     Import and pre-process Sentinel-3 data from the Copernicus program
- COPYRIGHT:   (C) 2024-202-20255 BY Norwegian Water and Energy Directorate,
-              Stefan Blumentrath, and the GRASS development team
-
-              This program is free software under the GNU General Public
-              License (>=v2). Read the file COPYING that comes with GRASS
-              for details.
+"""MODULE:      i.sentinel3.import
+AUTHOR(S):   Stefan Blumentrath
+PURPOSE:     Import and pre-process Sentinel-3 data from the Copernicus program
+COPYRIGHT:   (C) 2024-202-20255 BY Norwegian Water and Energy Directorate,
+Stefan Blumentrath, and the GRASS development team
+
+This program is free software under the GNU General Public
+License (>=v2). Read the file COPYING that comes with GRASS
+for details.
 """
 
 # %Module
@@ -139,13 +138,6 @@
 # % guisection: Settings
 # %end
 
-# to be implemented
-# # %flag
-# # % key: p
-# # % description: Print raster data to be imported and exit
-# # % guisection: Print
-# # %end
-
 # %flag
 # % key: r
 # % description: Rescale radiance bands to reflectance
@@ -157,25 +149,22 @@
 # # %end
 
 
-# ToDo
+# TODO
 # - implement printing
 # - implement cloud-masking at import
 # - implement elevation import
 # - do not write to temporary file (feed r.in.xyz from stdin)
 # - Add orphaned pixels
-# - parallelize NC conversion
 
+from __future__ import annotations
 
 import atexit
 import json
 import os
 import re
 import sys
-
 from collections import OrderedDict
 from datetime import datetime
-
-# from itertools import chain
 from pathlib import Path
 from zipfile import ZipFile
 
@@ -265,7 +254,7 @@
             "S5": 1.04,
             "S6": 1.07,
         },
-    }
+    },
 }
 
 S3_BANDS = {
@@ -495,17 +484,21 @@ def cleanup() -> None:
     external = gs.parse_key_val(gs.read_command("r.external.out", flags="p"), sep=": ")
     if "directory" in external:
         for map_file in Path(external["directory"]).glob(
-            f"{TMP_NAME}_*{external['extension']}"
+            f"{TMP_NAME}_*{external['extension']}",
         ):
             if map_file.is_file():
                 map_file.unlink()
 
     gs.run_command(
-        "g.remove", type="raster", pattern=f"{TMP_NAME}_*", flags="f", quiet=True
+        "g.remove",
+        type="raster",
+        pattern=f"{TMP_NAME}_*",
+        flags="f",
+        quiet=True,
     )
 
 
-def np_as_scalar(var):
+def np_as_scalar(var: np.dtype) -> str | int | float:
     """Return a numpy object as scalar."""
     if type(var).__module__ == np.__name__:
         if var.size > 1:
@@ -541,7 +534,7 @@ def get_dtype_range(np_datatype: str) -> dict:
     }
 
 
-def consolidate_metadata_dicts(metadata_dicts):
+def consolidate_metadata_dicts(metadata_dicts: dict) -> dict:
     """Consolidate a list of metadata dictionaries for all input scenes.
 
     Consolidate a list of metadata dictionaries for all
@@ -555,27 +548,26 @@ def consolidate_metadata_dicts(metadata_dicts):
             for meta_key, meta_value in meta_dict.items():
                 if meta_key not in map_metadata_dicts[map_name]:
                     map_metadata_dicts[map_name][meta_key] = meta_value
-                else:
-                    if (
-                        isinstance(map_metadata_dicts[map_name][meta_key], list)
-                        and meta_value not in map_metadata_dicts[map_name][meta_key]
-                    ):
-                        map_metadata_dicts[map_name][meta_key].append(meta_value)
-                    elif meta_value != map_metadata_dicts[map_name][meta_key]:
-                        map_metadata_dicts[map_name][meta_key] = [
-                            map_metadata_dicts[map_name][meta_key],
-                            meta_value,
-                        ]
+                elif (
+                    isinstance(map_metadata_dicts[map_name][meta_key], list)
+                    and meta_value not in map_metadata_dicts[map_name][meta_key]
+                ):
+                    map_metadata_dicts[map_name][meta_key].append(meta_value)
+                elif meta_value != map_metadata_dicts[map_name][meta_key]:
+                    map_metadata_dicts[map_name][meta_key] = [
+                        map_metadata_dicts[map_name][meta_key],
+                        meta_value,
+                    ]
     return map_metadata_dicts
 
 
 def write_metadata(json_dict: dict, metadatajson: str) -> None:
     """Write extended map metadata to JSON file."""
-    with open(metadatajson, "w", encoding="UTF8") as outfile:
+    with Path(metadatajson).open("w", encoding="UTF8") as outfile:
         json.dump(json_dict, outfile)
 
 
-def convert_units(np_column, from_u: str, to_u: str):
+def convert_units(np_column: np.array, from_u: str, to_u: str) -> np.array:
     """Convert units of a numpy array.
 
     Converts a numpy column from one unit to
@@ -587,8 +579,8 @@ def convert_units(np_column, from_u: str, to_u: str):
         gs.fatal(
             _(
                 "Could not import cf_units. Please install it with:\n"
-                "'pip install cf_units'!"
-            )
+                "'pip install cf_units'!",
+            ),
         )
 
     try:
@@ -596,8 +588,9 @@ def convert_units(np_column, from_u: str, to_u: str):
     except ValueError:
         gs.fatal(
             _("Warning: Could not convert units from {from_u} to {to_u}.").format(
-                from_u=from_u, to_u=to_u
-            )
+                from_u=from_u,
+                to_u=to_u,
+            ),
         )
         converted_col = np_column
 
@@ -632,19 +625,20 @@ def adjust_region(region_dict: dict) -> dict:
     )
     region_dict["n"] = (
         region_dict["s"]
-        + np.ceil(((region_dict["n"] + (nsres / 2.0) - region_dict["s"]) / nsres))
-        * nsres
+        + np.ceil((region_dict["n"] + (nsres / 2.0) - region_dict["s"]) / nsres) * nsres
     )
     region_dict["e"] = (
         region_dict["w"]
-        + np.ceil(((region_dict["e"] + (ewres / 2.0) - region_dict["w"]) / ewres))
-        * ewres
+        + np.ceil((region_dict["e"] + (ewres / 2.0) - region_dict["w"]) / ewres) * ewres
     )
     return region_dict
 
 
 def intersect_region(
-    region_current: dict, region_stripe: dict, *, align_current: bool = True
+    region_current: dict,
+    region_stripe: dict,
+    *,
+    align_current: bool = True,
 ) -> dict:
     """Intersect stripe region with current region aligned to resolution and extended for r.in.xyz.
 
@@ -653,7 +647,6 @@ def intersect_region(
     aligned to resolution in ew and ns direction
     starting from the sw-corner
     """
-
     # Remove potential, irrelevant / conflicting dict-keys and convert to float
     relevant_keys = ["n", "s", "e", "w", "nsres", "ewres"]
     region_current = {
@@ -666,14 +659,15 @@ def intersect_region(
     # Choose pixel alignment with stripe or region
     if not align_current:
         region_dict.update(
-            {"nsres": region_stripe["nsres"], "ewres": region_stripe["ewres"]}
+            {"nsres": region_stripe["nsres"], "ewres": region_stripe["ewres"]},
         )
     for direction in "nsew":
         resolution = region_dict["nsres"] if direction in "ns" else region_dict["ewres"]
         if direction in "ne":
             # Get the intersection of the stripe and the current region
             region_dict[direction] = min(
-                region_current[direction], region_stripe[direction]
+                region_current[direction],
+                region_stripe[direction],
             )
             # Align region to required resolution and extend by one pixel for r.in.xyz
             region_dict[direction] = (
@@ -681,12 +675,12 @@ def intersect_region(
             )
         else:
             region_dict[direction] = max(
-                region_current[direction], region_stripe[direction]
+                region_current[direction],
+                region_stripe[direction],
             )
             region_dict[direction] = (
                 np.floor(region_dict[direction] / resolution) * resolution - resolution
             )
-
     return region_dict
 
 
@@ -721,7 +715,7 @@ def parse_s3_file_name(file_name: str) -> dict:
         gs.fatal(_("{} is not a supported Sentinel-3 scene").format(str(file_name)))
 
 
-def extract_file_info(s3_files: list, basename: str = None):
+def extract_file_info(s3_files: list, basename: str = None) -> tuple[str, dict]:
     """Extract information from file name according to naming conventions."""
     result_dict = {}
     product_track_ids = [
@@ -745,10 +739,12 @@ def extract_file_info(s3_files: list, basename: str = None):
             if file_info["mission_id"] != result_dict["mission_id"]:
                 result_dict["mission_id"] = "S3_"
             result_dict["start_time"] = min(
-                result_dict["start_time"], file_info["start_time"]
+                result_dict["start_time"],
+                file_info["start_time"],
             )
             result_dict["end_time"] = max(
-                result_dict["end_time"], file_info["end_time"]
+                result_dict["end_time"],
+                file_info["end_time"],
             )
             for pid in product_track_ids:
                 result_dict[pid].add(file_info[pid])
@@ -757,8 +753,9 @@ def extract_file_info(s3_files: list, basename: str = None):
         if len(result_dict[pid]) > 1:
             gs.warning(
                 _("Merging {key} {values}").format(
-                    key=pid, values=", ".join(result_dict[pid])
-                )
+                    key=pid,
+                    values=", ".join(result_dict[pid]),
+                ),
             )
     if result_dict["mission_id"] == "3_":
         gs.warning(_("Merging Seninel-3A and Seninel-3B data"))
@@ -772,7 +769,7 @@ def extract_file_info(s3_files: list, basename: str = None):
                 result_dict["start_time"].strftime("%Y%m%d%H%M%S"),
                 result_dict["end_time"].strftime("%Y%m%d%H%M%S"),
                 *[list(result_dict[pid])[0] for pid in product_track_ids],
-            ]
+            ],
         )
     return (
         basename,
@@ -780,7 +777,13 @@ def extract_file_info(s3_files: list, basename: str = None):
     )
 
 
-def get_geocoding(zip_file, root, geo_bands_dict, sun_mask=None, region_bounds=None):
+def get_geocoding(
+    zip_file: ZipFile,
+    root: Path,
+    geo_bands_dict: dict,
+    sun_mask: np.array = None,
+    region_bounds: dict = None,
+) -> tuple[dict, np.array, list[int]]:
     """Get ground control points from NetCDF file."""
     member = str(root / geo_bands_dict["nc_file"])
     nc_file_path = zip_file.extract(member, path=TMP_FILE)
@@ -792,12 +795,12 @@ def get_geocoding(zip_file, root, geo_bands_dict, sun_mask=None, region_bounds=N
             if band not in nc_file_open.variables:
                 gs.fatal(
                     _(
-                        "{s3_file} does not contain a container {container} with band {band}"
+                        "{s3_file} does not contain a container {container} with band {band}",
                     ).format(
                         s3_file=str(Path(nc_file_path).parent.name),
                         container=geo_bands_dict["nc_file"],
                         band=", ".join(band),
-                    )
+                    ),
                 )
 
             # Add band to dict
@@ -831,27 +834,26 @@ def get_geocoding(zip_file, root, geo_bands_dict, sun_mask=None, region_bounds=N
             if mask.all():
                 gs.warning(
                     _("No valid pixels inside computational region found in {}").format(
-                        str(zip_file.filename)
-                    )
+                        str(zip_file.filename),
+                    ),
                 )
                 return None, None, None
-
     return nc_bands, mask, resolution
 
 
 def setup_import_multi_module(
-    tmp_ascii,
-    mapname,
-    distance=None,
-    fill_flags=False,
-    zrange=None,
-    val_col=None,
-    data_type=None,
-    method="mean",
-    solar_flux=None,
-    rules=None,
-):
-    """Setup GRASS GIS moduls for importing S3 bands."""
+    tmp_ascii: Path,
+    mapname: str,
+    distance: float = None,
+    fill_flags: str | bool = False,
+    zrange: tuple[float, float] = None,
+    val_col: int = None,
+    data_type: str = None,
+    method: str = "mean",
+    solar_flux: float | None = None,
+    rules: str | None = None,
+) -> list[Module]:
+    """Set up GRASS GIS moduls for importing S3 bands."""
     # Basic import module
     modules = [
         Module(
@@ -873,7 +875,7 @@ def setup_import_multi_module(
             verbose=False,
             quiet=True,
             overwrite=OVERWRITE,
-        )
+        ),
     ]
 
     # Interpolation of missing / empty pixels
@@ -902,7 +904,7 @@ def setup_import_multi_module(
                 quiet=True,
                 overwrite=OVERWRITE,
                 run_=False,
-            )
+            ),
         )
         mapname = mapname_reflectance
     else:
@@ -920,12 +922,12 @@ def setup_import_multi_module(
                 stdin_=rules,
                 separator=":",
                 run_=False,
-            )
+            ),
         )
     return modules
 
 
-def get_file_metadata(nc_dataset) -> dict:
+def get_file_metadata(nc_dataset: Dataset) -> dict:
     """Collect metadata from NetCDF file."""
     metadata = {
         attr: np_as_scalar(nc_dataset.getncattr(attr))
@@ -956,8 +958,13 @@ def get_file_metadata(nc_dataset) -> dict:
 
 
 def get_band_metadata(
-    band_tuple, nc_variable, fmt, file_metadata=None, basename=None, to_celsius=False
-):
+    band_tuple: tuple,
+    nc_variable: np.array,
+    fmt: str,
+    file_metadata: dict = None,
+    basename: str = None,
+    to_celsius: bool = False,
+) -> tuple[dict, str]:
     """Extract band metadata from NetCDF variable."""
     metadata = file_metadata.copy()
     # Collect metadata
@@ -982,8 +989,9 @@ def get_band_metadata(
     if datatype not in DTYPE_TO_GRASS:
         gs.fatal(
             _("Unsupported datatype {dt} in band {band}").format(
-                dt=datatype, band=band_tuple[1]
-            )
+                dt=datatype,
+                band=band_tuple[1],
+            ),
         )
     if datatype in {"uint8", "uint16"}:
         metadata["method"] = "max"  # Unfortunately there is no "mode" in r.in.xyz
@@ -1040,13 +1048,12 @@ def get_band_metadata(
     return metadata, fmt
 
 
-def transform_coordinates(coordinates):
+def transform_coordinates(coordinates: np.array) -> np.array:
     """Tranform coordinates to projection of the current location.
 
     Tranforms a numpy array with coordinates to the
     projection of the current location
     """
-
     # Create source coordinate reference
     s_srs = osr.SpatialReference()
     s_srs.ImportFromEPSG(4326)
@@ -1065,7 +1072,13 @@ def transform_coordinates(coordinates):
     )
 
 
-def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
+def write_xyz(
+    tmp_ascii: str,
+    nc_bands: np.array,
+    mask: np.array,
+    fmt: str = None,
+    project: bool = True,
+) -> None:
     """Write temporary XYZ ascii file."""
     # Extract grid coordinates
     lon = np.ma.masked_where(mask, nc_bands["lon"][:]).compressed()
@@ -1074,7 +1087,7 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
     if project:
         # Project coordinates
         np_output = transform_coordinates(
-            np.dstack((lat, lon)).reshape(lat.shape[0], 2)
+            np.dstack((lat, lon)).reshape(lat.shape[0], 2),
         )
     else:
         np_output = np.hstack((lon[:, None], lat[:, None]))
@@ -1086,26 +1099,30 @@ def write_xyz(tmp_ascii, nc_bands, mask, fmt=None, project=True):
         if np.ma.is_masked(add_array):
             add_array = add_array.filled()
         np_output = np.hstack(
-            (np_output, np.ma.masked_where(mask, add_array).compressed()[:, None])
+            (np_output, np.ma.masked_where(mask, add_array).compressed()[:, None]),
         )
 
     # Write to temporary file
-    if Path(tmp_ascii).exists():
-        with open(tmp_ascii, "ab") as tmp_ascii_file:
+    tmp_ascii_path = Path(tmp_ascii)
+    if tmp_ascii_path.exists():
+        with tmp_ascii_path.open("ab", encoding="UTF8") as tmp_ascii_file:
             np.savetxt(tmp_ascii_file, np_output, delimiter=",", fmt=fmt)
     else:
         np.savetxt(tmp_ascii, np_output, delimiter=",", fmt=fmt)
 
-    stripe_reg = {
+    return {
         "n": np.ma.max(np_output[:, 1]),
         "s": np.ma.min(np_output[:, 1]),
         "e": np.ma.max(np_output[:, 0]),
         "w": np.ma.min(np_output[:, 0]),
     }
-    return stripe_reg
 
 
-def import_s3(s3_file: str, kwargs: dict, s3_product: str = None):
+def import_s3(
+    s3_file: str,
+    kwargs: dict,
+    s3_product: str = None,
+) -> tuple[dict, list[str], dict]:
     """Import Sentinel-3 netCDF4 data."""
     # Unpack dictionary variables
     rmap = kwargs["meta_dict"][0]
@@ -1151,6 +1168,7 @@ def import_s3(s3_file: str, kwargs: dict, s3_product: str = None):
                         align_current=False,
                     ),
                 )
+
         for stripe, stripe_dict in s3_product.requested_stripe_content.items():
             # Could be parallelized!!!
             # Get geocoding (lat/lon, elevation) and initial mask
@@ -1183,20 +1201,29 @@ def import_s3(s3_file: str, kwargs: dict, s3_product: str = None):
                 )
                 register_strings.append(register_output)
             stripe_region_dict = write_xyz(
-                tmp_ascii, nc_bands, mask, fmt=fmt, project=True
+                tmp_ascii,
+                nc_bands,
+                mask,
+                fmt=fmt,
+                project=True,
             )
             stripe_region_dict["ewres"] = float(resolution[0])
             stripe_region_dict["nsres"] = float(resolution[1])
-            sun_region_dict_intersect = intersect_region(
-                dict(kwargs["current_reg"]), sun_region_dict, align_current=False
-            )
+
             if kwargs["maximum_solar_angle"]:
+                sun_region_dict_intersect = intersect_region(
+                    dict(kwargs["current_reg"]),
+                    sun_region_dict,
+                    align_current=False,
+                )
+
                 region_dicts[stripe] = intersect_region(
-                    sun_region_dict_intersect, stripe_region_dict, align_current=False
+                    sun_region_dict_intersect,
+                    stripe_region_dict,
+                    align_current=False,
                 )
             else:
                 region_dicts[stripe] = stripe_region_dict
-
     return module_queue, register_strings, region_dicts
 
 
@@ -1208,8 +1235,13 @@ class S3Product:
     """
 
     def __init__(
-        self, product_type, view="n", bands=None, flag_bands=None, anxillary_bands=None
-    ):
+        self,
+        product_type: str,
+        view: str = "n",
+        bands: list[str] | str = None,
+        flag_bands: list[str] | str = None,
+        anxillary_bands: list[str] | str = None,
+    ) -> None:
         """Initialize a Sentinel-3 product."""
         self.product_type = product_type
         self.available_views = S3_VIEWS[product_type]
@@ -1224,17 +1256,18 @@ def __init__(
         self.bands = self._check_bands(bands)
         self.flag_bands = self._check_bands(flag_bands, band_type="flag_bands")
         self.anxillary_bands = self._check_bands(
-            anxillary_bands, band_type="anxillary_bands"
+            anxillary_bands,
+            band_type="anxillary_bands",
         )
         self.grids = S3_GRIDS[product_type]
         self.requested_stripe_content = self._collect_requested_stripe_content()
         self.file_pattern = S3_FILE_PATTERN[product_type]
 
-    def __str__(self):
+    def __str__(self) -> str:
         """Return Sentinel-3 product information as JSON."""
         return json.dumps(self.__repr__(), indent=2)
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         """Return Sentinel-3 product information as JSON."""
         class_dict = self.__dict__.copy()
         for band_type in ["bands", "flag_bands", "anxillary_bands"]:
@@ -1249,17 +1282,22 @@ def __repr__(self):
             )
         return json.dumps(class_dict, indent=2)
 
-    def _check_view(self, view):
+    def _check_view(self, view: str) -> str:
         """Check requested Satellite view to process."""
         if view not in self.available_views:
             gs.warning(
                 _("View {} not available for product type {}").format(
-                    view, self.product_type
-                )
+                    view,
+                    self.product_type,
+                ),
             )
         return view
 
-    def _check_bands(self, bands, band_type="bands"):
+    def _check_bands(
+        self,
+        bands: list[str] | str,
+        band_type: str = "bands",
+    ) -> list[S3Band]:
         """Check requested Satellite bands to process."""
         band_types = {
             "bands": self.available_bands,
@@ -1279,13 +1317,14 @@ def _check_bands(self, bands, band_type="bands"):
             if band not in band_types[band_type]:
                 gs.warning(
                     _("Band {0} is not available in product_type {1}").format(
-                        band, self.product_type
-                    )
+                        band,
+                        self.product_type,
+                    ),
                 )
             band_objects[band] = S3Band(self.product_type, band, use_b=False, view="n")
         return band_objects
 
-    def _collect_requested_stripe_content(self):
+    def _collect_requested_stripe_content(self) -> dict:
         """Collect requested stripe content."""
         suffixes = {}
         for band_type in ["bands", "flag_bands", "anxillary_bands"]:
@@ -1323,16 +1362,16 @@ def _collect_requested_stripe_content(self):
 
     def process_nc_file(
         self,
-        zip_file,
-        root,
-        container_dict_items,
-        nc_bands,
-        prefix,
-        module_queue,
-        module_flags,
-        tmp_ascii,
-        fmt=None,
-    ):
+        zip_file: ZipFile,
+        root: Path,
+        container_dict_items: tuple,
+        nc_bands: list,
+        prefix: str,
+        module_queue: dict,
+        module_flags: dict,
+        tmp_ascii: str,
+        fmt: str = None,
+    ) -> tuple[list, dict, dict, str]:
         """Extract requested bands as numpy arrays from NetCDF file and setup import modules."""
         meta_information = {}
         member = str(root / container_dict_items[0])
@@ -1356,12 +1395,12 @@ def process_nc_file(
                 if band.full_name not in nc_file_open.variables:
                     gs.fatal(
                         _(
-                            "{s3_file} does not contain a container {container} with band {band}"
+                            "{s3_file} does not contain a container {container} with band {band}",
                         ).format(
                             s3_file=str(root),
                             container=band.nc_file,
                             band=", ".join(band.full_name),
-                        )
+                        ),
                     )
 
                 nc_variable = nc_file_open[band.full_name]
@@ -1416,12 +1455,12 @@ def process_nc_file(
                                 [
                                     str(nc_file_open[band.full_name].flag_masks[idx]),
                                     label,
-                                ]
+                                ],
                             )
                             for idx, label in enumerate(
-                                nc_file_open[band.full_name].flag_meanings.split(" ")
+                                nc_file_open[band.full_name].flag_meanings.split(" "),
                             )
-                        ]
+                        ],
                     )
                 # Setup import modules
                 fill_flags = "m"
@@ -1452,7 +1491,7 @@ def process_nc_file(
                         a: np_as_scalar(nc_file_open.getncattr(a))
                         for a in nc_file_open.ncattrs()
                     },
-                    **{"variable": band.full_name},
+                    "variable": band.full_name,
                     **{
                         a: np_as_scalar(nc_file_open[band.full_name].getncattr(a))
                         for a in nc_file_open[band.full_name].ncattrs()
@@ -1462,8 +1501,14 @@ def process_nc_file(
             return nc_bands, module_queue, meta_information, fmt
 
     def get_sun_parameters(
-        self, zip_file, root, prefix, tmp_ascii, region_bounds, maximum_solar_angle=None
-    ):
+        self,
+        zip_file: ZipFile,
+        root: Path,
+        prefix: str,
+        tmp_ascii: str,
+        region_bounds: dict,
+        maximum_solar_angle: float = None,
+    ) -> tuple[dict, dict, dict]:
         """Extract sun parameters from S3 SLSTR product.
 
         Oriented on:
@@ -1476,8 +1521,8 @@ def get_sun_parameters(
         member = str(
             root
             / S3_SUN_PARAMTERS[self.product_type]["sun_bands"]["nc_file"].format(
-                self.view
-            )
+                self.view,
+            ),
         )
         nc_file_path = zip_file.extract(member, path=TMP_FILE)
         with Dataset(nc_file_path) as sun_parameter_nc:
@@ -1536,7 +1581,7 @@ def get_sun_parameters(
                         a: np_as_scalar(sun_parameter_nc.getncattr(a))
                         for a in sun_parameter_nc.ncattrs()
                     },
-                    **{"variable": band},
+                    "variable": band,
                     **{
                         a: np_as_scalar(sun_parameter_nc[band].getncattr(a))
                         for a in sun_parameter_nc[band].ncattrs()
@@ -1545,11 +1590,18 @@ def get_sun_parameters(
 
             # Write to temporary file
             sun_region_dict = write_xyz(
-                tmp_ascii, nc_bands, mask, fmt=fmt, project=True
+                tmp_ascii,
+                nc_bands,
+                mask,
+                fmt=fmt,
+                project=True,
+            )
+            sun_region_dict["ewres"], sun_region_dict["nsres"] = (
+                float(
+                    resolution[0],
+                ),
+                float(resolution[1]),
             )
-            sun_region_dict["ewres"], sun_region_dict["nsres"] = float(
-                resolution[0]
-            ), float(resolution[1])
 
         return import_modules, meta_information, sun_region_dict
 
@@ -1559,12 +1611,13 @@ class S3Band:
 
     def __init__(
         self,
-        product_type,
-        band_id,
-        use_b=False,
-        view="n",
-        radiance_adjustment="S3_PN_SLSTR_L1_08",
-    ):
+        product_type: str,
+        band_id: str,
+        use_b: bool = False,
+        view: str = "n",
+        radiance_adjustment: str = "S3_PN_SLSTR_L1_08",
+    ) -> None:
+        """Initialize a Sentinel-3 band."""
         self.band_id = band_id
         self.band_type = self._get_band_type(product_type)
         # Need to call in this order
@@ -1575,24 +1628,25 @@ def __init__(
         self.exception = self._get_exception_band(product_type)
         self.nc_file = self._get_nc_file(product_type)
         self.radiance_adjustment = self._get_radiance_adjustment(
-            radiance_adjustment, view
+            radiance_adjustment,
+            view,
         )
         self.solar_flux = self._get_solar_flux_dict_for_band(product_type)
 
-    def __str__(self):
+    def __str__(self) -> str:
         return json.dumps(self.__repr__(), indent=2)
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         return json.dumps(self.__dict__, indent=2)
 
-    def _get_band_type(self, product_type):
+    def _get_band_type(self, product_type: str) -> str | None:
         for band_type, s3_bands in S3_BANDS.items():
             if self.band_id in s3_bands[product_type]:
                 return band_type
         gs.fatal(_("Oh no"))
         return None
 
-    def _get_geometry(self, product_type, use_b):
+    def _get_geometry(self, product_type: str, use_b: bool) -> str:
         band_dict = S3_BANDS[self.band_type][product_type][self.band_id]
         return (
             "b"
@@ -1600,7 +1654,7 @@ def _get_geometry(self, product_type, use_b):
             else band_dict["geometries"][0]
         )
 
-    def _get_full_name(self, product_type):
+    def _get_full_name(self, product_type: str) -> str:
         band_dict = S3_BANDS[self.band_type][product_type][self.band_id]
         types = band_dict["types"]
         if self.band_type != "bands":
@@ -1609,7 +1663,7 @@ def _get_full_name(self, product_type):
             return band_dict["full_name"].format(self.band_id, types, self.suffix)
         return band_dict["full_name"].format(self.band_id, self.suffix)
 
-    def _get_nc_file(self, product_type):
+    def _get_nc_file(self, product_type: str) -> str:
         band_dict = S3_BANDS[self.band_type][product_type][self.band_id]
         types = band_dict["types"]
         if self.band_type == "flags":
@@ -1618,21 +1672,25 @@ def _get_nc_file(self, product_type):
             return band_dict["nc_file"].format(self.band_id, types, self.suffix)
         return band_dict["nc_file"].format(self.band_id, self.suffix)
 
-    def _get_radiance_adjustment(self, radiance_adjustment, view):
+    def _get_radiance_adjustment(
+        self,
+        radiance_adjustment: str,
+        view: str,
+    ) -> float | None:
         if radiance_adjustment not in S3_RADIANCE_ADJUSTMENT:
             gs.fatal("Missing information on radiance adjustment")
         if self.band_id in S3_RADIANCE_ADJUSTMENT[radiance_adjustment][view]:
             return S3_RADIANCE_ADJUSTMENT[radiance_adjustment][view][self.band_id]
         return None
 
-    def _get_exception_band(self, product_type):
+    def _get_exception_band(self, product_type: str) -> dict | None:
         if "exception" in S3_BANDS[self.band_type][product_type][self.band_id]:
             return S3_BANDS[self.band_type][product_type][self.band_id][
                 "exception"
             ].format(self.band_id, "exception", self.suffix)
         return None
 
-    def _get_solar_flux_dict_for_band(self, product_type):
+    def _get_solar_flux_dict_for_band(self, product_type: str) -> dict | None:
         if (
             S3_SOLAR_FLUX[product_type]
             and self.band_id in S3_SOLAR_FLUX[product_type]["defaults"]
@@ -1640,10 +1698,12 @@ def _get_solar_flux_dict_for_band(self, product_type):
             return {
                 "default": S3_SOLAR_FLUX[product_type]["defaults"][self.band_id],
                 "band": S3_SOLAR_FLUX[product_type]["band"].format(
-                    self.band_id, self.suffix
+                    self.band_id,
+                    self.suffix,
                 ),
                 "nc_file": S3_SOLAR_FLUX[product_type]["nc_file"].format(
-                    self.band_id, self.suffix
+                    self.band_id,
+                    self.suffix,
                 ),
             }
         return None
@@ -1676,7 +1736,7 @@ def _get_band_geo_points(self, suffix: str) -> dict:
             "elevation": f"elevation_{suffix}",
         }
 
-    def adjust_radiance(self, np_band_array):
+    def adjust_radiance(self, np_band_array: np.array) -> np.array:
         """Get radiance adjustment for band object."""
         if not self.radiance_adjustment:
             return np_band_array * self.radiance_adjustment
@@ -1687,10 +1747,12 @@ def get_solar_angle_bounds(region_bounds: dict, sun_region_dict: dict) -> dict:
     """Get latlon coordinates of bounds in solar angle bands."""
     solar_bounds = gs.parse_command("g.region", flags="ugl", **sun_region_dict)
     solar_bounds["ll_n"] = max(
-        float(solar_bounds["nw_lat"]), float(solar_bounds["ne_lat"])
+        float(solar_bounds["nw_lat"]),
+        float(solar_bounds["ne_lat"]),
     )
     solar_bounds["ll_s"] = min(
-        float(solar_bounds["sw_lat"]), float(solar_bounds["se_lat"])
+        float(solar_bounds["sw_lat"]),
+        float(solar_bounds["se_lat"]),
     )
     if solar_bounds["ll_n"] < float(region_bounds["ll_n"]):
         region_bounds["ll_n"] = solar_bounds["ll_n"]
@@ -1702,7 +1764,7 @@ def get_solar_angle_bounds(region_bounds: dict, sun_region_dict: dict) -> dict:
 def main() -> None:
     """Do the main work."""
     pattern = re.compile(
-        ".*" + S3_FILE_PATTERN[options["product_type"]].replace("*", ".*")
+        ".*" + S3_FILE_PATTERN[options["product_type"]].replace("*", ".*"),
     )
 
     # check provided input
@@ -1713,8 +1775,8 @@ def main() -> None:
         except ValueError:
             gs.fatal(
                 _(
-                    "Input <{}> is neither a supported Sentinel-3 scene nor a text file with scenes"
-                ).format(options["input"])
+                    "Input <{}> is neither a supported Sentinel-3 scene nor a text file with scenes",
+                ).format(options["input"]),
             )
 
     if not s3_files:
@@ -1728,14 +1790,17 @@ def main() -> None:
         if not re.match(pattern, s3_scene.name):
             gs.fatal(
                 _(
-                    "Input <{sn}> is not a supported Sentinel-3 scene for the requested product type <{pt}>"
-                ).format(sn=s3_scene, pt=options["product_type"])
+                    "Input <{sn}> is not a supported Sentinel-3 scene for the requested product type <{pt}>",
+                ).format(sn=s3_scene, pt=options["product_type"]),
             )
 
     if flags["d"]:
         global DTYPE_TO_GRASS
         DTYPE_TO_GRASS["float64"] = "DCELL"
 
+    solar_bands = S3_SUN_PARAMTERS[options["product_type"]]
+    solar_bands = solar_bands.get("sun_bands").get("bands") if solar_bands else {}
+
     s3_product = S3Product(
         options["product_type"],
         view="o" if flags["o"] else "n",
@@ -1766,27 +1831,6 @@ def main() -> None:
         else None,
     }
 
-    # if flags["p"]:
-    #     print(
-    #         "|".join(
-    #             [
-    #                 "product_file_name",
-    #                 "nc_file_name",
-    #                 "nc_file_title",
-    #                 "nc_file_start_time",
-    #                 "nc_file_creation_time",
-    #                 "band",
-    #                 "band_shape",
-    #                 "band_title",
-    #                 "band_standard_name",
-    #                 "band_long_name",
-    #             ]
-    #         )
-    #     )
-    #     print("\n".join(chain(*import_result)))
-    #     cleanup()
-    #     return 0
-
     module_queues = []
     register_strings = []
     region_dicts = {}
@@ -1795,7 +1839,9 @@ def main() -> None:
     for s3_file in s3_files:
         gs.verbose(_("Preparing scene {} for import").format(s3_file.name))
         module_list, register_string, region_dict = import_s3(
-            s3_file, import_dict, s3_product=s3_product
+            s3_file,
+            import_dict,
+            s3_product=s3_product,
         )
         if module_list is not None:
             module_queues.append(module_list)
@@ -1810,20 +1856,33 @@ def main() -> None:
                     else:
                         extend_region(region_dicts[stripe_id], region)
     if not module_queues:
-        gs.warning(_("Nothing to import with the given input"))
+        gs.warning(_("Nothing to import with the given input."))
         sys.exit(0)
     stripe_envs = {}
     for stripe_id, stripe_region in region_dicts.items():
         stripe_env = os.environ.copy()
+        if (
+            stripe_region["s"] >= stripe_region["n"]
+            or stripe_region["w"] >= stripe_region["e"]
+        ):
+            gs.warning(
+                _(
+                    "No valid data found in data stripe {}.\n"
+                    "Nothing to import with the given input."
+                ).format(stripe_id),
+            )
+            sys.exit(0)
+
         if flags["n"]:
             stripe_env["GRASS_REGION"] = gs.region_env(**adjust_region(stripe_region))
         else:
+            stripe_region_dict = intersect_region(
+                dict(current_region),
+                stripe_region,
+                align_current=stripe_id != "sun_parameters",
+            )
             stripe_env["GRASS_REGION"] = gs.region_env(
-                **intersect_region(
-                    dict(current_region),
-                    stripe_region,
-                    align_current=stripe_id != "sun_parameters",
-                )
+                **stripe_region_dict,
             )
         stripe_envs[stripe_id] = stripe_env
 
@@ -1831,7 +1890,7 @@ def main() -> None:
         gs.verbose(_("Importing solar parameter bands"))
         queue = ParallelModuleQueue(nprocs)
         compute_env = stripe_envs["sun_parameters"]
-        for solar_parameter in ("solar_azimuth", "solar_zenith"):
+        for solar_parameter in solar_bands:
             module_list = []
             for solar_module in module_queues[0][solar_parameter]:
                 solar_module.env_ = compute_env
@@ -1845,12 +1904,12 @@ def main() -> None:
 
     gs.verbose(
         _("Importing scenes {}").format(
-            "\n" + "\n".join([s3_file.name for s3_file in s3_files])
-        )
+            "\n" + "\n".join([s3_file.name for s3_file in s3_files]),
+        ),
     )
     queue = ParallelModuleQueue(nprocs)
     for band_id, modules in module_queues[0].items():
-        if band_id in {"solar_azimuth", "solar_zenith"}:
+        if band_id in solar_bands:
             continue
         # S3SL2LST product has only "in" stripe and band names have no suffix
         if band_id[-2:] not in stripe_envs:
@@ -1874,9 +1933,13 @@ def main() -> None:
             mapname if not flags["r"] else mapname.replace("radiance", "reflectance")
         )
         if not gs.raster_info(mapname)["max"]:
-            gs.warning(("Raster map {} is empty. Removing...").format(mapname))
+            gs.warning(f"Raster map {mapname} is empty. Removing...")
             gs.run_command(
-                "g.remove", type="raster", name=mapname, flags="f", quiet=True
+                "g.remove",
+                type="raster",
+                name=mapname,
+                flags="f",
+                quiet=True,
             )
             continue
         metadata["semantic_label"] = (
@@ -1890,7 +1953,10 @@ def main() -> None:
         # Write extended metadata if requested
         if flags["j"]:
             json_standard_folder = Path(GISENV["GISDBASE"]).joinpath(
-                GISENV["LOCATION_NAME"], GISENV["MAPSET"], "cell_misc", mapname
+                GISENV["LOCATION_NAME"],
+                GISENV["MAPSET"],
+                "cell_misc",
+                mapname,
             )
             if not json_standard_folder.exists():
                 json_standard_folder.mkdir(parents=True, exist_ok=True)
@@ -1935,27 +2001,29 @@ def main() -> None:
                             [
                                 grass_timestamp(meta_info_dict[1]["start_time"]),
                                 grass_timestamp(meta_info_dict[1]["end_time"]),
-                            ]
+                            ],
                         ),
                         quiet=True,
                         run_=False,
                     ),
-                ]
-            )
+                ],
+            ),
         )
         t_register_strings.append(
             "|".join(
                 [
                     f"{mapname}@{GISENV['MAPSET']}",
                     meta_info_dict[1]["start_time"].isoformat(
-                        sep=" ", timespec="seconds"
+                        sep=" ",
+                        timespec="seconds",
                     ),
                     meta_info_dict[1]["end_time"].isoformat(
-                        sep=" ", timespec="seconds"
+                        sep=" ",
+                        timespec="seconds",
                     ),
                     metadata["semantic_label"],
-                ]
-            )
+                ],
+            ),
         )
     queue.wait()
 
@@ -1963,10 +2031,11 @@ def main() -> None:
     t_register_strings = "\n".join(t_register_strings) if t_register_strings else ""
     if options["register_output"]:
         gs.verbose(
-            _("Writing register file <{}>...").format(options["register_output"])
+            _("Writing register file <{}>...").format(options["register_output"]),
         )
         Path(options["register_output"]).write_text(
-            t_register_strings + "\n", encoding="UTF8"
+            t_register_strings + "\n",
+            encoding="UTF8",
         )
     else:
         print(t_register_strings)
@@ -1982,8 +2051,8 @@ def main() -> None:
         gs.fatal(
             _(
                 "Could not import dateutil. Please install it with:\n"
-                "'pip install dateutil'!"
-            )
+                "'pip install dateutil'!",
+            ),
         )
     try:
         from osgeo import osr
@@ -1991,7 +2060,7 @@ def main() -> None:
         osr.UseExceptions()
     except ImportError:
         gs.fatal(
-            _("Could not import gdal. Please install it with:\n'pip install GDAL'!")
+            _("Could not import gdal. Please install it with:\n'pip install GDAL'!"),
         )
 
     try:
@@ -2000,8 +2069,8 @@ def main() -> None:
         gs.fatal(
             _(
                 "Could not import netCDF4. Please install it with:\n"
-                "'pip install netcdf4'!"
-            )
+                "'pip install netcdf4'!",
+            ),
         )
 
     try:
@@ -2009,9 +2078,8 @@ def main() -> None:
     except ImportError:
         gs.fatal(
             _(
-                "Could not import numpy. Please install it with:\n"
-                "'pip install numpy'!"
-            )
+                "Could not import numpy. Please install it with:\n'pip install numpy'!",
+            ),
         )
 
     atexit.register(cleanup)

From 1f8325bd91fb7bd88e9915d0f01c095328d33c10 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Fri, 11 Apr 2025 09:42:45 +0200
Subject: [PATCH 22/27] fix open() regression

---
 src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 8a6ace92e3..c90a1cdf72 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -1105,7 +1105,7 @@ def write_xyz(
     # Write to temporary file
     tmp_ascii_path = Path(tmp_ascii)
     if tmp_ascii_path.exists():
-        with tmp_ascii_path.open("ab", encoding="UTF8") as tmp_ascii_file:
+        with tmp_ascii_path.open("ab") as tmp_ascii_file:
             np.savetxt(tmp_ascii_file, np_output, delimiter=",", fmt=fmt)
     else:
         np.savetxt(tmp_ascii, np_output, delimiter=",", fmt=fmt)

From e293d940d72bafbb600bc043011c2d6a863c0154 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Fri, 25 Apr 2025 14:54:47 +0200
Subject: [PATCH 23/27] fix corner cases

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 36 +++++++++++++------
 1 file changed, 25 insertions(+), 11 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index c90a1cdf72..8824068b19 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -552,12 +552,19 @@ def consolidate_metadata_dicts(metadata_dicts: dict) -> dict:
                     isinstance(map_metadata_dicts[map_name][meta_key], list)
                     and meta_value not in map_metadata_dicts[map_name][meta_key]
                 ):
-                    map_metadata_dicts[map_name][meta_key].append(meta_value)
+                    if isinstance(meta_value, (list, set)):
+                        map_metadata_dicts[map_name][meta_key].extend(meta_value)
+                    else:
+                        map_metadata_dicts[map_name][meta_key].append(meta_value)
                 elif meta_value != map_metadata_dicts[map_name][meta_key]:
-                    map_metadata_dicts[map_name][meta_key] = [
-                        map_metadata_dicts[map_name][meta_key],
-                        meta_value,
-                    ]
+                    if isinstance(map_metadata_dicts[map_name][meta_key], (list, set)):
+                        if meta_value not in map_metadata_dicts[map_name][meta_key]:
+                            map_metadata_dicts[map_name][meta_key].append(meta_value)
+                    else:
+                        map_metadata_dicts[map_name][meta_key] = [
+                            map_metadata_dicts[map_name][meta_key],
+                            meta_value,
+                        ]
     return map_metadata_dicts
 
 
@@ -949,7 +956,6 @@ def get_file_metadata(nc_dataset: Dataset) -> dict:
         ]
         if attr in nc_dataset.ncattrs()
     }
-
     metadata["start_time"] = parse_timestr(nc_dataset.start_time)
     metadata["end_time"] = parse_timestr(nc_dataset.stop_time)
     metadata["history"] = metadata["history"].strip()
@@ -1158,15 +1164,23 @@ def import_s3(
             region_dicts["sun_parameters"] = sun_region_dict
 
             if kwargs["maximum_solar_angle"]:
+                sun_region_dict = intersect_region(
+                    dict(kwargs["current_reg"]), sun_region_dict, align_current=False
+                )
+                print(sun_region_dict)
+                if sun_region_dict["e"] <= sun_region_dict["w"]:
+                    # East is wrapped around meridian
+                    if sun_region_dict["e"] <= float(kwargs["current_reg"]["w"]):
+                        sun_region_dict["e"] = float(kwargs["current_reg"]["e"])
+                    # West is wrapped around meridian
+                    if sun_region_dict["w"] >= float(kwargs["current_reg"]["e"]):
+                        sun_region_dict["w"] = float(kwargs["current_reg"]["w"])
+
                 sun_region_bounds = gs.parse_command(
                     "g.region",
                     flags="ugb",
                     quiet=True,
-                    **intersect_region(
-                        dict(kwargs["current_reg"]),
-                        sun_region_dict,
-                        align_current=False,
-                    ),
+                    **sun_region_dict,
                 )
 
         for stripe, stripe_dict in s3_product.requested_stripe_content.items():

From 0d0073f9db129e9a53e10380c012acf035535a99 Mon Sep 17 00:00:00 2001
From: Stefan Blumentrath 
Date: Tue, 29 Apr 2025 10:51:03 +0200
Subject: [PATCH 24/27] remove left-over print

---
 src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 8824068b19..e7b7b6cc0d 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -1167,7 +1167,7 @@ def import_s3(
                 sun_region_dict = intersect_region(
                     dict(kwargs["current_reg"]), sun_region_dict, align_current=False
                 )
-                print(sun_region_dict)
+
                 if sun_region_dict["e"] <= sun_region_dict["w"]:
                     # East is wrapped around meridian
                     if sun_region_dict["e"] <= float(kwargs["current_reg"]["w"]):

From 45554c32132a0791122caa83700062ec5748b710 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Tue, 13 May 2025 10:16:08 +0200
Subject: [PATCH 25/27] check region validity

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 28 +++++++++++--------
 1 file changed, 17 insertions(+), 11 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index e7b7b6cc0d..396741c432 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -1775,6 +1775,21 @@ def get_solar_angle_bounds(region_bounds: dict, sun_region_dict: dict) -> dict:
     return region_bounds
 
 
+def check_region_validity(stripe_id: str, stripe_region: dict) -> None:
+    """Check if the computational region is valid and exit if not."""
+    if (
+        stripe_region["s"] >= stripe_region["n"]
+        or stripe_region["w"] >= stripe_region["e"]
+    ):
+        gs.warning(
+            _(
+                "No valid data found in data stripe {}.\n"
+                "Nothing to import with the given input."
+            ).format(stripe_id),
+        )
+        sys.exit(0)
+
+
 def main() -> None:
     """Do the main work."""
     pattern = re.compile(
@@ -1875,17 +1890,7 @@ def main() -> None:
     stripe_envs = {}
     for stripe_id, stripe_region in region_dicts.items():
         stripe_env = os.environ.copy()
-        if (
-            stripe_region["s"] >= stripe_region["n"]
-            or stripe_region["w"] >= stripe_region["e"]
-        ):
-            gs.warning(
-                _(
-                    "No valid data found in data stripe {}.\n"
-                    "Nothing to import with the given input."
-                ).format(stripe_id),
-            )
-            sys.exit(0)
+        check_region_validity(stripe_id, stripe_region)
 
         if flags["n"]:
             stripe_env["GRASS_REGION"] = gs.region_env(**adjust_region(stripe_region))
@@ -1895,6 +1900,7 @@ def main() -> None:
                 stripe_region,
                 align_current=stripe_id != "sun_parameters",
             )
+            check_region_validity(stripe_id, stripe_region_dict)
             stripe_env["GRASS_REGION"] = gs.region_env(
                 **stripe_region_dict,
             )

From 4bd53b2d478e6d6fdc12c5237bf9f24836721ae4 Mon Sep 17 00:00:00 2001
From: ninsbl 
Date: Wed, 25 Jun 2025 11:05:46 +0200
Subject: [PATCH 26/27] avoid metadata key errors

---
 .../i.sentinel3.import/i.sentinel3.import.py  | 42 ++++++++++---------
 1 file changed, 22 insertions(+), 20 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 396741c432..79df3c9a91 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -722,7 +722,7 @@ def parse_s3_file_name(file_name: str) -> dict:
         gs.fatal(_("{} is not a supported Sentinel-3 scene").format(str(file_name)))
 
 
-def extract_file_info(s3_files: list, basename: str = None) -> tuple[str, dict]:
+def extract_file_info(s3_files: list, basename: str | None = None) -> tuple[str, dict]:
     """Extract information from file name according to naming conventions."""
     result_dict = {}
     product_track_ids = [
@@ -775,7 +775,7 @@ def extract_file_info(s3_files: list, basename: str = None) -> tuple[str, dict]:
                 result_dict["product"],
                 result_dict["start_time"].strftime("%Y%m%d%H%M%S"),
                 result_dict["end_time"].strftime("%Y%m%d%H%M%S"),
-                *[list(result_dict[pid])[0] for pid in product_track_ids],
+                *[next(iter(result_dict[pid])) for pid in product_track_ids],
             ],
         )
     return (
@@ -984,8 +984,6 @@ def get_band_metadata(
     mapname = f"{basename}_{varname_short}"
     metadata["mapname"] = mapname
 
-    # band_title = nc_variable.long_name if "long_name" in band_attrs else band_tuple[1]
-
     # Define unit
     unit = nc_variable.units if "units" in band_attrs else None
     unit = "degree_celsius" if band_tuple[0].startswith("LST") and to_celsius else unit
@@ -1165,7 +1163,9 @@ def import_s3(
 
             if kwargs["maximum_solar_angle"]:
                 sun_region_dict = intersect_region(
-                    dict(kwargs["current_reg"]), sun_region_dict, align_current=False
+                    dict(kwargs["current_reg"]),
+                    sun_region_dict,
+                    align_current=False,
                 )
 
                 if sun_region_dict["e"] <= sun_region_dict["w"]:
@@ -1648,9 +1648,11 @@ def __init__(
         self.solar_flux = self._get_solar_flux_dict_for_band(product_type)
 
     def __str__(self) -> str:
+        """Return a JSON representation of the S3Band instance."""
         return json.dumps(self.__repr__(), indent=2)
 
     def __repr__(self) -> str:
+        """Return a JSON representation of the S3Band instance."""
         return json.dumps(self.__dict__, indent=2)
 
     def _get_band_type(self, product_type: str) -> str | None:
@@ -1784,12 +1786,19 @@ def check_region_validity(stripe_id: str, stripe_region: dict) -> None:
         gs.warning(
             _(
                 "No valid data found in data stripe {}.\n"
-                "Nothing to import with the given input."
+                "Nothing to import with the given input.",
             ).format(stripe_id),
         )
         sys.exit(0)
 
 
+def _get_grass_metadata(grass_md: dict, md_field: str) -> str | None:
+    md_field = grass_md.get(md_field)
+    if isinstance(md_field, list):
+        return ",".join(md_field)
+    return md_field
+
+
 def main() -> None:
     """Do the main work."""
     pattern = re.compile(
@@ -1986,23 +1995,16 @@ def main() -> None:
         description = (
             json.dumps(metadata, separators=["\n", ": "]).lstrip("{").rstrip("}")
         )
+
         support_kwargs = {
             "map": mapname,
-            "title": ",".join(metadata["title"])
-            if isinstance(metadata["title"], list)
-            else metadata["title"],
-            "history": ",".join(metadata["history"])
-            if isinstance(metadata["history"], list)
-            else metadata["history"],
-            "units": metadata["unit"],
-            "source1": ",".join(metadata["product_name"])
-            if isinstance(metadata["product_name"], list)
-            else metadata["product_name"],
-            "source2": ",".join(metadata["processing_baseline"])
-            if isinstance(metadata["processing_baseline"], list)
-            else metadata["processing_baseline"],
+            "title": _get_grass_metadata(metadata, "title"),
+            "history": _get_grass_metadata(metadata, "history"),
+            "units": metadata.get("unit"),
+            "source1": _get_grass_metadata(metadata, "product_name"),
+            "source2": _get_grass_metadata(metadata, "processing_baseline"),
             "description": description,
-            "semantic_label": metadata["semantic_label"],
+            "semantic_label": metadata.get("semantic_label"),
         }
 
         queue.put(

From 4b7bf9809d93bef6f434ea235abc4127ef6390b4 Mon Sep 17 00:00:00 2001
From: Stefan Blumentrath 
Date: Wed, 10 Dec 2025 15:27:59 +0100
Subject: [PATCH 27/27] Apply suggestions from code review

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>
---
 .../i.sentinel3.import/i.sentinel3.import.py  | 29 ++++++++++---------
 1 file changed, 15 insertions(+), 14 deletions(-)

diff --git a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
index 79df3c9a91..276302c45b 100644
--- a/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
+++ b/src/imagery/i.sentinel/i.sentinel3.import/i.sentinel3.import.py
@@ -173,6 +173,7 @@
 from grass.temporal.datetime_math import (
     datetime_to_grass_datetime_string as grass_timestamp,
 )
+from typing import Optional
 
 S3_SUPPORTED_PRODUCTS = ["S3SL1RBT", "S3SL2LST"]
 
@@ -789,7 +790,7 @@ def get_geocoding(
     root: Path,
     geo_bands_dict: dict,
     sun_mask: np.array = None,
-    region_bounds: dict = None,
+    region_bounds: dict | None = None,
 ) -> tuple[dict, np.array, list[int]]:
     """Get ground control points from NetCDF file."""
     member = str(root / geo_bands_dict["nc_file"])
@@ -851,11 +852,11 @@ def get_geocoding(
 def setup_import_multi_module(
     tmp_ascii: Path,
     mapname: str,
-    distance: float = None,
+    distance: float | None = None,
     fill_flags: str | bool = False,
-    zrange: tuple[float, float] = None,
-    val_col: int = None,
-    data_type: str = None,
+    zrange: tuple[float, float] | None = None,
+    val_col: int | None = None,
+    data_type: str | None = None,
     method: str = "mean",
     solar_flux: float | None = None,
     rules: str | None = None,
@@ -967,8 +968,8 @@ def get_band_metadata(
     band_tuple: tuple,
     nc_variable: np.array,
     fmt: str,
-    file_metadata: dict = None,
-    basename: str = None,
+    file_metadata: dict | None = None,
+    basename: str | None = None,
     to_celsius: bool = False,
 ) -> tuple[dict, str]:
     """Extract band metadata from NetCDF variable."""
@@ -1080,7 +1081,7 @@ def write_xyz(
     tmp_ascii: str,
     nc_bands: np.array,
     mask: np.array,
-    fmt: str = None,
+    fmt: str | None = None,
     project: bool = True,
 ) -> None:
     """Write temporary XYZ ascii file."""
@@ -1125,7 +1126,7 @@ def write_xyz(
 def import_s3(
     s3_file: str,
     kwargs: dict,
-    s3_product: str = None,
+    s3_product: str | None = None,
 ) -> tuple[dict, list[str], dict]:
     """Import Sentinel-3 netCDF4 data."""
     # Unpack dictionary variables
@@ -1252,9 +1253,9 @@ def __init__(
         self,
         product_type: str,
         view: str = "n",
-        bands: list[str] | str = None,
-        flag_bands: list[str] | str = None,
-        anxillary_bands: list[str] | str = None,
+        bands: list[str] | str | None = None,
+        flag_bands: list[str] | str | None = None,
+        anxillary_bands: list[str] | str | None = None,
     ) -> None:
         """Initialize a Sentinel-3 product."""
         self.product_type = product_type
@@ -1384,7 +1385,7 @@ def process_nc_file(
         module_queue: dict,
         module_flags: dict,
         tmp_ascii: str,
-        fmt: str = None,
+        fmt: str | None = None,
     ) -> tuple[list, dict, dict, str]:
         """Extract requested bands as numpy arrays from NetCDF file and setup import modules."""
         meta_information = {}
@@ -1521,7 +1522,7 @@ def get_sun_parameters(
         prefix: str,
         tmp_ascii: str,
         region_bounds: dict,
-        maximum_solar_angle: float = None,
+        maximum_solar_angle: float | None = None,
     ) -> tuple[dict, dict, dict]:
         """Extract sun parameters from S3 SLSTR product.