Ingestion of synoptic observations

src/CSET/cset_workflow/.gitignore

@@ -2,6 +2,7 @@
 site/
 opt/
 demo_pointstat/
+*metdb*
 restricted*
 
 # User workflow configuration.

src/CSET/cset_workflow/app/fetch_fcst/bin/fetch_data.py

@@ -15,6 +15,7 @@
 """Retrieve the files from the filesystem for the current cycle point."""
 
 import abc
+import ast
 import glob
 import itertools
 import logging
@@ -177,6 +178,28 @@ def _get_needed_environment_variables() -> dict:
     return variables
 
 
+def _get_needed_environment_variables_obs() -> dict:
+    """Load the needed variables from the environment."""
+    variables = {
+        "subtype": os.environ["OBS_SUBTYPE"],
+        "data_time": datetime.fromisoformat(os.environ["CYLC_TASK_CYCLE_POINT"]),
+        "forecast_length": isodate.parse_duration(os.environ["ANALYSIS_LENGTH"]),
+        "obs_fields": ast.literal_eval(os.environ["SURFACE_SYNOP_FIELDS"]),
+        "model_identifier": "OBS",
+        "wmo_nmbrs": ast.literal_eval(os.environ.get("WMO_BLOCK_STTN_NMBRS"))
+        if len(os.environ.get("WMO_BLOCK_STTN_NMBRS")) > 0
+        else None,
+        "subarea_extent": ast.literal_eval(os.environ.get("SUBAREA_EXTENT"))
+        if len(os.environ.get("SUBAREA_EXTENT")) > 0
+        else None,
+        "obs_interval": isodate.parse_duration(os.environ["SURFACE_SYNOP_INTERVAL"]),
+        "obs_offset": isodate.parse_duration(os.environ["SURFACE_SYNOP_OFFSET"]),
+        "rose_datac": os.environ["ROSE_DATAC"],
+    }
+    logging.debug("Environment variables loaded: %s", variables)
+    return variables
+
+
 def _template_file_path(
     raw_path: str,
     date_type: Literal["validity", "initiation"],
@@ -270,3 +293,72 @@ def fetch_data(file_retriever: FileRetrieverABC):
     # exiting the with block.
     if not files_found:
         raise FileNotFoundError("No files found for model!")
+
+
+def fetch_obs(obs_retriever: FileRetrieverABC):
+    """Fetch the observations corresponding to a model run.
+
+    The following environment variables need to be set:
+    * ANALYSIS_OFFSET
+    * ANALYSIS_LENGTH
+    * CYLC_TASK_CYCLE_POINT
+    * DATA_PATH
+    * DATA_PERIOD
+    * DATE_TYPE
+    * MODEL_IDENTIFIER
+    * ROSE_DATAC
+
+    Parameters
+    ----------
+    obs_retriever: ObsRetriever
+        ObsRetriever implementation to use. Defaults to FilesystemFileRetriever.
+
+    Raises
+    ------
+    FileNotFound:
+        If no observations are available.
+    """
+    v = _get_needed_environment_variables_obs()
+
+    # Prepare output directory.
+    cycle_obs_dir = f"{v['rose_datac']}/data/OBS"
+    os.makedirs(cycle_obs_dir, exist_ok=True)
+    logging.debug("Output directory: %s", cycle_obs_dir)
+
+    # We will get just one file for now, but follow the templating
+    # syntax for the model for consistency.
+    obs_base_path = (
+        v["subtype"]
+        + "_"
+        + "%Y%m%dT%H%MZ_dt_"
+        + str(int(v["forecast_length"].total_seconds() // 3600)).zfill(3)
+        + ".nc"
+    )
+    paths = _template_file_path(
+        obs_base_path,
+        "initiation",
+        v["data_time"],
+        v["forecast_length"],
+        timedelta(seconds=0),
+        v["obs_interval"],
+    )
+    logging.info("Retrieving paths:\n%s", "\n".join(paths))
+
+    # Use obs retriever to transfer data with multiple threads.
+    # We shouldn't need to iterate as we do for the forecast data
+    # because these files will be smaller.
+    try:
+        obs_retriever.get_file(
+            paths[0],
+            v["subtype"],
+            v["obs_fields"],
+            v["data_time"],
+            v["obs_offset"],
+            v["forecast_length"],
+            v["obs_interval"],
+            cycle_obs_dir,
+            wmo_nmbrs=v["wmo_nmbrs"],
+            subarea_extent=v["subarea_extent"],
+        )
+    except Exception as exc:
+        raise ValueError("No observations available.") from exc

src/CSET/cset_workflow/app/fetch_obs/rose-app.conf

@@ -0,0 +1,3 @@
+[command]
+default=echo "Please set ROSE_APP_COMMAND_KEY to your database."; false
+metdb=app_env_wrapper fetch-obs-metdb.py

src/CSET/cset_workflow/flow.cylc

@@ -42,6 +42,7 @@ final cycle point = {{CSET_TRIAL_END_DATE}}
             R1/{{date}} = """
             setup_complete[^] =>
             FETCH_DATA:succeed-all =>
+            FETCH_OBS:succeed-all =>
             fetch_complete =>
             PROCESS:{% if LOGLEVEL == "DEBUG" %}succeed-all{% else %}finish-all{% endif %} =>
             cycle_complete
@@ -52,6 +53,7 @@ final cycle point = {{CSET_TRIAL_END_DATE}}
         {{CSET_TRIAL_CYCLE_PERIOD}} = """
         setup_complete[^] =>
         FETCH_DATA:succeed-all =>
+        FETCH_OBS:succeed-all =>
         fetch_complete =>
         PROCESS:{% if LOGLEVEL == "DEBUG" %}succeed-all{% else %}finish-all{% endif %} =>
         cycle_complete
@@ -109,13 +111,16 @@ final cycle point = {{CSET_TRIAL_END_DATE}}
     execution retry delays = PT1M, PT15M, PT1H
     [[[directives]]]
         --mem=10000
+    [[[environment]]]
+        PLOTTING_PROJECTION = {{PLOTTING_PROJECTION | default("")}}
 
     [[PROCESS_CASE_AGGREGATION]]
     script = rose task-run -v --app-key=run_cset_recipe
     execution time limit = PT3H
         [[[directives]]]
         --mem=100000
         [[[environment]]]
+        PLOTTING_PROJECTION = {{PLOTTING_PROJECTION | default("")}}
         # As this process is used for case aggregation we hard code to True.
         DO_CASE_AGGREGATION = True
 
@@ -125,6 +130,12 @@ final cycle point = {{CSET_TRIAL_END_DATE}}
         [[[environment]]]
         ANALYSIS_LENGTH = {{ANALYSIS_LENGTH}}
 
+    [[FETCH_OBS]]
+    script = rose task-run -v --app-key=fetch_obs
+    execution time limit = PT1H
+        [[[environment]]]
+        ANALYSIS_LENGTH = {{ANALYSIS_LENGTH}}
+
     [[METPLUS]]
         [[[environment]]]
         {% if METPLUS_GRID_STAT|default(False) %}
@@ -150,6 +161,9 @@ final cycle point = {{CSET_TRIAL_END_DATE}}
     [[cycle_complete]]
     inherit = DUMMY_TASK
 
+    [[dummy_fetch_obs]]
+    inherit = DUMMY_TASK, FETCH_OBS
+
     [[dummy_process]]
     inherit = DUMMY_TASK, PROCESS
 
@@ -185,6 +199,21 @@ final cycle point = {{CSET_TRIAL_END_DATE}}
         ANALYSIS_OFFSET = {{model["analysis_offset"]}}
     {% endfor %}
 
+    {% if SURFACE_SYNOP_OBS %}
+    [[fetch_obs]]
+    # Fetch observations from the MetDB.
+    inherit = FETCH_OBS
+        [[[environment]]]
+        MODEL_IDENTIFIER = "OBS"
+        ROSE_APP_COMMAND_KEY = "metdb"
+        SURFACE_SYNOP_FIELDS = {{SURFACE_SYNOP_FIELDS}}
+        SURFACE_SYNOP_INTERVAL = {{SURFACE_SYNOP_INTERVAL}}
+        SURFACE_SYNOP_OFFSET = {{SURFACE_SYNOP_OFFSET}}
+        WMO_BLOCK_STTN_NMBRS = {{WMO_BLOCK_STTN_NMBRS | default("")}}
+        SUBAREA_EXTENT = {{SUBAREA_EXTENT | default("")}}
+        OBS_SUBTYPE = "LNDSYB"
+    {% endif %}
+
     [[housekeeping]]
     # Housekeep input data files.
         [[[environment]]]

src/CSET/cset_workflow/includes/spatial_surface_model_obs_diference.cylc

@@ -0,0 +1,25 @@
+{% if SURFACE_SYNOP_OBS and SURFACE_SYNOP_DIFFS %}
+{% for model in models %}
+{% for obs_field in (SURFACE_SYNOP_FIELDS | list) %}
+[runtime]
+    [[obs_model_difference_scatterplot_m{{model["id"]}}_{{obs_field}}]]
+    inherit = PROCESS
+    execution time limit = PT6H
+        [[[environment]]]
+        CSET_RECIPE_NAME = "generic_model_obs_difference_scatterplot.yaml"
+        CSET_ADDOPTS = """
+            --OBSVARNAME={{obs_field}}
+            --VARNAME={{(SURFACE_SYNOP_MODEL_FIELDS | list)[loop.index0]}}
+            --BASE_MODEL='{{model["id"]}}'
+            --OBSERVATIONS='OBS'
+            --MODEL_NAME='{{model["name"]}}'
+            --DIFFERENCE=True
+            --SUBAREA_TYPE={{SUBAREA_TYPE | default("")}}
+            --SUBAREA_EXTENT={{SUBAREA_TYPE | default("")}}
+            --PLOTTING_PROJECTION={{PLOTTING_PROJECTION | default("")}}
+        """
+        MODEL_IDENTIFIERS = 'OBS {{model["id"]}}'
+
+{% endfor %}
+{% endfor %}
+{% endif %}

src/CSET/cset_workflow/includes/spatial_surface_obs_field.cylc

@@ -0,0 +1,18 @@
+{% if SURFACE_SYNOP_OBS %}
+{% for obs_field in SURFACE_SYNOP_FIELDS %}
+[runtime]
+    [[generic_obs_scatterplot_{{obs_field}}]]
+    inherit = PROCESS
+    execution time limit = PT1H
+        [[[environment]]]
+        CSET_RECIPE_NAME = "generic_obs_scatterplot.yaml"
+        CSET_ADDOPTS = """
+            --OBSVARNAME={{obs_field}}
+            --SUBAREA_TYPE={{SUBAREA_TYPE | default("")}}
+            --SUBAREA_EXTENT={{SUBAREA_TYPE | default("")}}
+            --PLOTTING_PROJECTION={{PLOTTING_PROJECTION | default("")}}
+        """
+        MODEL_IDENTIFIERS = "OBS"
+
+{% endfor %}
+{% endif %}

src/CSET/cset_workflow/meta/observations/rose-meta.conf

@@ -0,0 +1,105 @@
+################################################################################
+# Observations
+################################################################################
+
+[Observations]
+ns=Observations
+sort-key=sec-d
+title=Observational Comparisons
+
+################################################################################
+# Spatial (2D map) fields
+################################################################################
+
+[Observations/Types]
+ns=Observations/Types
+sort-key=sec-d1
+title=Types of observations
+
+# Surface synoptic fields.
+[template variables=SURFACE_SYNOP_OBS]
+ns=Observations/Types
+description=Compare with surface synoptic observations
+type=python_boolean
+compulsory=true
+sort-key=0synop1a
+trigger=template variables=SURFACE_SYNOP_FIELDS: True;
+        template variables=USE_WMO_STATION_NUMBERS: True;
+        template variables=SURFACE_SYNOP_INTERVAL: True;
+        template variables=SURFACE_SYNOP_OFFSET: True;
+
+################################################################
+
+[Observations/Types/Synoptic]
+ns=Observations/Types/Synoptic
+sort-key=sec-d2
+title=Synoptic observations
+
+[template variables=SURFACE_SYNOP_FIELDS]
+ns=Observations/Types/Synoptic
+title=Surface synoptic observations
+description=Synoptic observations
+help=Name of field in the style of the MetDB.
+compulsory=true
+type=python_list
+sort-key=0synop1b
+
+[template variables=SURFACE_SYNOP_DIFFS]
+ns=Observations/Types/Synoptic
+description=Calculate the corresponding differences.
+type=python_boolean
+compulsory=true
+sort-key=0synop1c
+trigger=template variables=SURFACE_SYNOP_MODEL_FIELDS: True;
+
+[template variables=SURFACE_SYNOP_MODEL_FIELDS]
+ns=Observations/Types/Synoptic
+title=Corresponding model variables
+description=Synoptic observations
+help=Name of field in the model: the order must match that of the synoptic observations
+compulsory=true
+type=python_list
+sort-key=0synop1d
+
+[template variables=USE_WMO_STATION_NUMBERS]
+ns=Observations/Types/Synoptic
+description=Use WMO numbers
+help=Select to set stations by WMO station or block number.
+compulsory=true
+type=python_boolean
+sort-key=0synop1e
+trigger=template variables=WMO_BLOCK_STTN_NMBRS: True;
+
+[template variables=WMO_BLOCK_STTN_NMBRS]
+ns=Observations/Types/Synoptic
+title=WMO Block or station numbers.
+description=WMO numbers
+help=WMO block numbers (2 digits) or station numbers (5 digits).
+compulsory=true
+type=python_list
+sort-key=0synop1f
+
+[template variables=SURFACE_SYNOP_INTERVAL]
+ns=Observations/Types/Synoptic
+title=Interval between synoptic observations
+description=Synoptic observations
+help=The time interval between rerievals of synoptic observations.
+compulsory=true
+type=quoted
+sort-key=0synop1g
+
+[template variables=SURFACE_SYNOP_OFFSET]
+ns=Observations/Types/Synoptic
+title=Offset for synoptic observations
+description=Offset for observations
+help=The offset of the first retrieved observations relative to the initial time.
+compulsory=true
+type=quoted
+sort-key=0synop1h
+
+################################################################
+
+[Observations/Types/Placeholder]
+ns=Observations/Types/Placeholder
+sort-key=sec-e2
+title=Observations from an undefined source.

src/CSET/cset_workflow/meta/rose-meta.conf.jinja2

@@ -158,6 +158,17 @@ type=quoted
 compulsory=true
 sort-key=setup-h-out2
 
+[template variables=PLOTTING_PROJECTION]
+ns=Setup
+title=Projection for the plots
+description=Geographical projection of plots
+help=This is is the projection used for the final plots, specified as a string.
+    If no value is given, native latitudes and longitudes are used.
+compulsory=true
+values="", "NP_Stereo"
+value-titles="", "North polar stereographic"
+sort-key=setup-h-out3a
+
 [template variables=PLOT_RESOLUTION]
 ns=Setup
 title=Plot resolution