Add run and finalize methods to marine LETKF task (#2944)

Adds run and finalize methods to marine LETKF task, experiment yaml for
gw-ci in GDASApp, workflow additions, removes bugs found on the way, and
completes the bulk of the work on LETKF task. Conversion of fields to
increments pending.

Partially resolves NOAA-EMC/GDASApp#1091 and NOAA-EMC/GDASApp#1251

Mutual dependency with GDASApp PR
NOAA-EMC/GDASApp#1287 and IC fix file issue
#2944 (comment)

---------

Co-authored-by: Walter Kolczynski - NOAA <Walter.Kolczynski@noaa.gov>
Co-authored-by: David Huber <69919478+DavidHuber-NOAA@users.noreply.github.com>

ci/cases/pr/C48mx500_hybAOWCDA.yaml

@@ -0,0 +1,26 @@
+experiment:
+  system: gfs
+  mode: cycled
+
+arguments:
+  pslot: {{ 'pslot' | getenv }}
+  app: S2S
+  resdetatmos: 48
+  resdetocean: 5.0
+  resensatmos: 48
+  comroot: {{ 'RUNTESTS' | getenv }}/COMROOT
+  expdir: {{ 'RUNTESTS' | getenv }}/EXPDIR
+  icsdir: {{ 'ICSDIR_ROOT' | getenv }}/C48mx500/20240610
+  idate: 2021032412
+  edate: 2021032418
+  nens: 3
+  interval: 0
+  start: warm
+  yaml: {{ HOMEgfs }}/ci/cases/yamls/soca_gfs_defaults_ci.yaml
+
+skip_ci_on_hosts:
+  - wcoss2
+  - orion
+  - hercules
+  - hera
+  - gaea

env/HERA.env

@@ -153,10 +153,10 @@ elif [[ "${step}" = "ocnanalecen" ]]; then
     export NTHREADS_OCNANALECEN=${NTHREADSmax}
     export APRUN_OCNANALECEN="${APRUN_default} --cpus-per-task=${NTHREADS_OCNANALECEN}"
 
-elif [[ "${step}" = "marineanalletkf" ]]; then
+elif [[ "${step}" = "marineanlletkf" ]]; then
 
-    export NTHREADS_MARINEANALLETKF=${NTHREADSmax}
-    export APRUN_MARINEANALLETKF="${APRUN_default} --cpus-per-task=${NTHREADS_MARINEANALLETKF}"
+    export NTHREADS_MARINEANLLETKF=${NTHREADSmax}
+    export APRUN_MARINEANLLETKF=${APRUN_default}
 
 elif [[ "${step}" = "anal" ]] || [[ "${step}" = "analcalc" ]]; then
 

env/ORION.env

@@ -149,10 +149,10 @@ elif [[ "${step}" = "ocnanalecen" ]]; then
     export NTHREADS_OCNANALECEN=${NTHREADSmax}
     export APRUN_OCNANALECEN="${APRUN_default} --cpus-per-task=${NTHREADS_OCNANALECEN}"
 
-elif [[ "${step}" = "marineanalletkf" ]]; then
+elif [[ "${step}" = "marineanlletkf" ]]; then
 
-    export NTHREADS_MARINEANALLETKF=${NTHREADSmax}
-    export APRUN_MARINEANALLETKF="${APRUN_default} --cpus-per-task=${NTHREADS_MARINEANALLETKF}"
+    export NTHREADS_MARINEANLLETKF=${NTHREADSmax}
+    export APRUN_MARINEANLLETKF="${APRUN_default}"
 
 elif [[ "${step}" = "anal" ]] || [[ "${step}" = "analcalc" ]]; then
 

env/WCOSS2.env

@@ -126,10 +126,10 @@ elif [[ "${step}" = "ocnanalecen" ]]; then
     export NTHREADS_OCNANALECEN=${NTHREADSmax}
     export APRUN_OCNANALECEN="${APRUN_default} --cpus-per-task=${NTHREADS_OCNANALECEN}"
 
-elif [[ "${step}" = "marineanalletkf" ]]; then
+elif [[ "${step}" = "marineanlletkf" ]]; then
 
-    export NTHREADS_MARINEANALLETKF=${NTHREADSmax}
-    export APRUN_MARINEANALLETKF="${APRUN_default} --cpus-per-task=${NTHREADS_MARINEANALLETKF}"
+    export NTHREADS_MARINEANLLETKF=${NTHREADSmax}
+    export APRUN_MARINEANLLETKF="${APRUN_default}"
 
 elif [[ "${step}" = "atmanlfv3inc" ]]; then
 

jobs/JGLOBAL_MARINE_ANALYSIS_LETKF

@@ -1,6 +1,13 @@
 #!/bin/bash
 source "${HOMEgfs}/ush/preamble.sh"
-source "${HOMEgfs}/ush/jjob_header.sh" -e "marineanalletkf" -c "base ocnanal marineanalletkf"
+
+export DATAjob="${DATAROOT}/${RUN}marineanalysis.${PDY:-}${cyc}"
+export DATA="${DATAjob}/${jobid}"
+# Create the directory to hold ensemble perturbations
+export DATAens="${DATAjob}/ensdata"
+if [[ ! -d "${DATAens}" ]]; then mkdir -p "${DATAens}"; fi
+
+source "${HOMEgfs}/ush/jjob_header.sh" -e "marineanlletkf" -c "base marineanl marineanlletkf"
 
 ##############################################
 # Set variables used in the script
@@ -11,12 +18,18 @@ GDATE=$(date --utc +%Y%m%d%H -d "${PDY} ${cyc} - ${assim_freq} hours")
 
 gPDY=${GDATE:0:8}
 gcyc=${GDATE:8:2}
+export GDUMP="gdas"
+export GDUMP_ENS="enkf${GDUMP}"
+export OPREFIX="${RUN}.t${cyc}z."
 
-YMD=${gPDY} HH=${gcyc} declare_from_tmpl -rx \
-   COMIN_OCEAN_HISTORY_PREV:COM_OCEAN_HISTORY_TMPL \
-   COMIN_ICE_HISTORY_PREV:COM_ICE_HISTORY_TMPL
+RUN="${GDUMP}" YMD=${gPDY} HH=${gcyc} declare_from_tmpl -rx \
+  COMIN_OCEAN_HISTORY_PREV:COM_OCEAN_HISTORY_TMPL \
+  COMIN_ICE_HISTORY_PREV:COM_ICE_HISTORY_TMPL
 
-YMD=${PDY} HH=${cyc} declare_from_tmpl -rx COMIN_OBS:COM_OBS_TMPL 
+YMD=${PDY} HH=${cyc} declare_from_tmpl -rx \
+  COMIN_OBS:COM_OBS_TMPL \
+  COMOUT_OCEAN_LETKF:COM_OCEAN_LETKF_TMPL \
+  COMOUT_ICE_LETKF:COM_ICE_LETKF_TMPL
 
 ##############################################
 # Begin JOB SPECIFIC work
@@ -25,7 +38,7 @@ YMD=${PDY} HH=${cyc} declare_from_tmpl -rx COMIN_OBS:COM_OBS_TMPL
 ###############################################################
 # Run relevant script
 
-EXSCRIPT=${GDASOCNLETKFPY:-${HOMEgfs}/scripts/exgdas_global_marine_analysis_letkf.py}
+EXSCRIPT=${GDASOCNLETKFPY:-${HOMEgfs}/scripts/exglobal_marine_analysis_letkf.py}
 ${EXSCRIPT}
 status=$?
 [[ ${status} -ne 0 ]] && exit "${status}"

jobs/rocoto/marineanalletkf.sh → jobs/rocoto/marineanlletkf.sh

@@ -8,7 +8,7 @@ source "${HOMEgfs}/ush/preamble.sh"
 status=$?
 [[ ${status} -ne 0 ]] && exit "${status}"
 
-export job="marineanalletkf"
+export job="marineanlletkf"
 export jobid="${job}.$$"
 
 ###############################################################

parm/config/gfs/config.com

@@ -82,12 +82,14 @@ declare -rx COM_OCEAN_HISTORY_TMPL=${COM_BASE}'/model/ocean/history'
 declare -rx COM_OCEAN_RESTART_TMPL=${COM_BASE}'/model/ocean/restart'
 declare -rx COM_OCEAN_INPUT_TMPL=${COM_BASE}'/model/ocean/input'
 declare -rx COM_OCEAN_ANALYSIS_TMPL=${COM_BASE}'/analysis/ocean'
+declare -rx COM_OCEAN_LETKF_TMPL=${COM_BASE}'/analysis/ocean/letkf'
 declare -rx COM_OCEAN_BMATRIX_TMPL=${COM_BASE}'/bmatrix/ocean'
 declare -rx COM_OCEAN_NETCDF_TMPL=${COM_BASE}'/products/ocean/netcdf'
 declare -rx COM_OCEAN_GRIB_TMPL=${COM_BASE}'/products/ocean/grib2'
 declare -rx COM_OCEAN_GRIB_GRID_TMPL=${COM_OCEAN_GRIB_TMPL}'/${GRID}'
 
 declare -rx COM_ICE_ANALYSIS_TMPL=${COM_BASE}'/analysis/ice'
+declare -rx COM_ICE_LETKF_TMPL=${COM_BASE}'/analysis/ice/letkf'
 declare -rx COM_ICE_BMATRIX_TMPL=${COM_BASE}'/bmatrix/ice'
 declare -rx COM_ICE_INPUT_TMPL=${COM_BASE}'/model/ice/input'
 declare -rx COM_ICE_HISTORY_TMPL=${COM_BASE}'/model/ice/history'

parm/config/gfs/config.marineanlletkf

@@ -0,0 +1,20 @@
+#!/bin/bash
+
+########## config.marineanlletkf ##########
+# Ocn Analysis specific
+
+echo "BEGIN: config.marineanlletkf"
+
+# Get task specific resources
+. "${EXPDIR}/config.resources" marineanlletkf
+
+export MARINE_LETKF_EXEC="${EXECgfs}/gdas.x"
+export MARINE_LETKF_YAML_TMPL="${PARMgfs}/gdas/soca/letkf/letkf.yaml.j2"
+export MARINE_LETKF_STAGE_YAML_TMPL="${PARMgfs}/gdas/soca/letkf/letkf_stage.yaml.j2"
+export MARINE_LETKF_SAVE_YAML_TMPL="${PARMgfs}/gdas/soca/letkf/letkf_save.yaml.j2"
+
+export GRIDGEN_EXEC="${EXECgfs}/gdas_soca_gridgen.x"
+export GRIDGEN_YAML="${PARMgfs}/gdas/soca/gridgen/gridgen.yaml"
+export DIST_HALO_SIZE=500000
+
+echo "END: config.marineanlletkf"

parm/config/gfs/config.resources

@@ -601,7 +601,7 @@ case ${step} in
     tasks_per_node=$(( max_tasks_per_node / threads_per_task ))
     ;;
 
-  "marineanalletkf")
+  "marineanlletkf")
     ntasks=16
     case ${OCNRES} in
       "025")

parm/stage/ocean_ens_perturbations.yaml.j2

@@ -9,5 +9,5 @@ ocean_ens_perturbation:
         {% for mem in range(first_mem + 1, last_mem + 1) %}
         {% set imem = mem - first_mem %}
         {% set COMOUT_OCEAN_ANALYSIS_MEM = COMOUT_OCEAN_ANALYSIS_MEM_list[imem] %}
-        - ["{{ ICSDIR }}/{{ COMOUT_OCEAN_ANALYSIS_MEM | relpath(ROTDIR) }}/{{ m_prefix }}.mom6_perturbation.nc", "{{ COMOUT_OCEAN_ANALYSIS_MEM }}/mom6_increment.nc"]
+        - ["{{ ICSDIR }}/{{ COMOUT_OCEAN_ANALYSIS_MEM | relpath(ROTDIR) }}/{{ m_prefix }}.mom6_perturbation.nc", "{{ COMOUT_OCEAN_ANALYSIS_MEM }}/{{ RUN }}.t{{ current_cycle_HH }}z.ocninc.nc"]
         {% endfor %} # mem loop

ush/forecast_postdet.sh

@@ -466,12 +466,13 @@ MOM6_postdet() {
     fi
 
     # GEFS perturbations
-    # TODO if [[ $RUN} == "gefs" ]] block maybe be needed
-    #     to ensure it does not interfere with the GFS when ensemble is updated in the GFS
-    if (( MEMBER > 0 )) && [[ "${ODA_INCUPD:-False}" == "True" ]]; then
-      ${NCP} "${COMIN_OCEAN_ANALYSIS}/mom6_increment.nc" "${DATA}/INPUT/mom6_increment.nc" \
-      || ( echo "FATAL ERROR: Unable to copy ensemble MOM6 increment, ABORT!"; exit 1 )
-    fi
+    if [[ "${RUN}" == "gefs" ]]; then 
+    #     to ensure it does not interfere with the GFS
+      if (( MEMBER > 0 )) && [[ "${ODA_INCUPD:-False}" == "True" ]]; then
+        ${NCP} "${COMIN_OCEAN_ANALYSIS}/${RUN}.t${cyc}z.ocninc.nc" "${DATA}/INPUT/mom6_increment.nc" \
+        || ( echo "FATAL ERROR: Unable to copy ensemble MOM6 increment, ABORT!"; exit 1 )
+      fi
+    fi # if [[ "${RUN}" == "gefs" ]]; then 
   fi  # if [[ "${RERUN}" == "NO" ]]; then
 
   # Link output files

ush/python/pygfs/task/marine_bmat.py

@@ -318,7 +318,7 @@ def finalize(self: Task) -> None:
         FileHandler({'copy': diagb_list}).sync()
 
         # Copy the ensemble perturbation diagnostics to the ROTDIR
-        if self.task_config.DOHYBVAR == "YES" or self.task_config.NMEM_ENS > 3:
+        if self.task_config.DOHYBVAR == "YES" or self.task_config.NMEM_ENS > 2:
             window_middle_iso = self.task_config.MARINE_WINDOW_MIDDLE.strftime('%Y-%m-%dT%H:%M:%SZ')
             weight_list = []
             src = os.path.join(self.task_config.DATA, f"ocn.ens_weights.incr.{window_middle_iso}.nc")

ush/python/pygfs/task/marine_letkf.py

@@ -1,11 +1,13 @@
 #!/usr/bin/env python3
 
 import f90nml
+import pygfs.utils.marine_da_utils as mdau
 from logging import getLogger
 import os
 from pygfs.task.analysis import Analysis
 from typing import Dict
 from wxflow import (AttrDict,
+                    Executable,
                     FileHandler,
                     logit,
                     parse_j2yaml,
@@ -41,6 +43,8 @@ def __init__(self, config: Dict) -> None:
                             'soca',
                             'localensembleda',
                             _letkf_yaml_file]
+        # compute the relative path from self.task_config.DATA to self.task_config.DATAenspert
+        _enspert_relpath = os.path.relpath(self.task_config.DATAens, self.task_config.DATA)
 
         self.task_config.WINDOW_MIDDLE = self.task_config.current_cycle
         self.task_config.WINDOW_BEGIN = self.task_config.current_cycle - _half_assim_freq
@@ -49,6 +53,7 @@ def __init__(self, config: Dict) -> None:
         self.task_config.mom_input_nml_tmpl = os.path.join(self.task_config.DATA, 'mom_input.nml.tmpl')
         self.task_config.mom_input_nml = os.path.join(self.task_config.DATA, 'mom_input.nml')
         self.task_config.obs_dir = os.path.join(self.task_config.DATA, 'obs')
+        self.task_config.ENSPERT_RELPATH = _enspert_relpath
 
     @logit(logger)
     def initialize(self):
@@ -64,26 +69,50 @@ def initialize(self):
         logger.info("initialize")
 
         # make directories and stage ensemble background files
-        ensbkgconf = AttrDict()
-        keys = ['previous_cycle', 'current_cycle', 'DATA', 'NMEM_ENS',
-                'PARMgfs', 'ROTDIR', 'COM_OCEAN_HISTORY_TMPL', 'COM_ICE_HISTORY_TMPL']
-        for key in keys:
-            ensbkgconf[key] = self.task_config[key]
-        ensbkgconf.RUN = 'enkfgdas'
-        soca_ens_bkg_stage_list = parse_j2yaml(self.task_config.SOCA_ENS_BKG_STAGE_YAML_TMPL, ensbkgconf)
-        FileHandler(soca_ens_bkg_stage_list).sync()
         soca_fix_stage_list = parse_j2yaml(self.task_config.SOCA_FIX_YAML_TMPL, self.task_config)
         FileHandler(soca_fix_stage_list).sync()
-        letkf_stage_list = parse_j2yaml(self.task_config.MARINE_LETKF_STAGE_YAML_TMPL, self.task_config)
+        stageconf = AttrDict()
+        keys = ['current_cycle',
+                'previous_cycle',
+                'COM_ICE_LETKF_TMPL',
+                'COM_OCEAN_LETKF_TMPL',
+                'COM_ICE_HISTORY_TMPL',
+                'COM_OCEAN_HISTORY_TMPL',
+                'COMIN_OCEAN_HISTORY_PREV',
+                'COMIN_ICE_HISTORY_PREV',
+                'COMOUT_ICE_LETKF',
+                'COMOUT_OCEAN_LETKF',
+                'DATA',
+                'ENSPERT_RELPATH',
+                'GDUMP_ENS',
+                'NMEM_ENS',
+                'OPREFIX',
+                'PARMgfs',
+                'ROTDIR',
+                'RUN',
+                'WINDOW_BEGIN',
+                'WINDOW_MIDDLE']
+        for key in keys:
+            stageconf[key] = self.task_config[key]
+
+        # stage ensemble background files
+        soca_ens_bkg_stage_list = parse_j2yaml(self.task_config.MARINE_ENSDA_STAGE_BKG_YAML_TMPL, stageconf)
+        FileHandler(soca_ens_bkg_stage_list).sync()
+
+        # stage letkf-specific files
+        letkf_stage_list = parse_j2yaml(self.task_config.MARINE_LETKF_STAGE_YAML_TMPL, stageconf)
         FileHandler(letkf_stage_list).sync()
 
-        obs_list = parse_j2yaml(self.task_config.OBS_YAML, self.task_config)
+        obs_list = parse_j2yaml(self.task_config.MARINE_OBS_LIST_YAML, self.task_config)
 
         # get the list of observations
         obs_files = []
         for ob in obs_list['observers']:
             obs_name = ob['obs space']['name'].lower()
-            obs_filename = f"{self.task_config.RUN}.t{self.task_config.cyc}z.{obs_name}.{to_YMDH(self.task_config.current_cycle)}.nc"
+            # TODO(AFE) - this should be removed when the obs config yamls are jinjafied
+            if 'distribution' not in ob['obs space']:
+                ob['obs space']['distribution'] = {'name': 'Halo', 'halo size': self.task_config['DIST_HALO_SIZE']}
+            obs_filename = f"{self.task_config.RUN}.t{self.task_config.cyc}z.{obs_name}.{to_YMDH(self.task_config.current_cycle)}.nc4"
             obs_files.append((obs_filename, ob))
 
         obs_files_to_copy = []
@@ -102,12 +131,7 @@ def initialize(self):
         FileHandler({'copy': obs_files_to_copy}).sync()
 
         # make the letkf.yaml
-        letkfconf = AttrDict()
-        keys = ['WINDOW_BEGIN', 'WINDOW_MIDDLE', 'RUN', 'gcyc', 'NMEM_ENS']
-        for key in keys:
-            letkfconf[key] = self.task_config[key]
-        letkfconf.RUN = 'enkfgdas'
-        letkf_yaml = parse_j2yaml(self.task_config.MARINE_LETKF_YAML_TMPL, letkfconf)
+        letkf_yaml = parse_j2yaml(self.task_config.MARINE_LETKF_YAML_TMPL, stageconf)
         letkf_yaml.observations.observers = obs_to_use
         letkf_yaml.save(self.task_config.letkf_yaml_file)
 
@@ -133,6 +157,18 @@ def run(self):
 
         logger.info("run")
 
+        exec_cmd_gridgen = Executable(self.task_config.APRUN_MARINEANLLETKF)
+        exec_cmd_gridgen.add_default_arg(self.task_config.GRIDGEN_EXEC)
+        exec_cmd_gridgen.add_default_arg(self.task_config.GRIDGEN_YAML)
+
+        mdau.run(exec_cmd_gridgen)
+
+        exec_cmd_letkf = Executable(self.task_config.APRUN_MARINEANLLETKF)
+        for letkf_exec_arg in self.task_config.letkf_exec_args:
+            exec_cmd_letkf.add_default_arg(letkf_exec_arg)
+
+        mdau.run(exec_cmd_letkf)
+
     @logit(logger)
     def finalize(self):
         """Method finalize for ocean and sea ice LETKF task
@@ -145,3 +181,11 @@ def finalize(self):
         """
 
         logger.info("finalize")
+
+        letkfsaveconf = AttrDict()
+        keys = ['current_cycle', 'DATA', 'NMEM_ENS', 'WINDOW_BEGIN', 'GDUMP_ENS',
+                'PARMgfs', 'ROTDIR', 'COM_OCEAN_LETKF_TMPL', 'COM_ICE_LETKF_TMPL']
+        for key in keys:
+            letkfsaveconf[key] = self.task_config[key]
+        letkf_save_list = parse_j2yaml(self.task_config.MARINE_LETKF_SAVE_YAML_TMPL, letkfsaveconf)
+        FileHandler(letkf_save_list).sync()

workflow/applications/gfs_cycled.py

@@ -46,7 +46,7 @@ def _get_app_configs(self):
         if self.do_jediocnvar:
             configs += ['prepoceanobs', 'marineanlinit', 'marinebmat', 'marineanlvar']
             if self.do_hybvar:
-                configs += ['ocnanalecen']
+                configs += ['marineanlletkf', 'ocnanalecen']
             configs += ['marineanlchkpt', 'marineanlfinal']
             if self.do_vrfy_oceanda:
                 configs += ['ocnanalvrfy']
@@ -148,7 +148,7 @@ def get_task_names(self):
         if self.do_jediocnvar:
             gdas_gfs_common_tasks_before_fcst += ['prepoceanobs', 'marineanlinit', 'marinebmat', 'marineanlvar']
             if self.do_hybvar:
-                gdas_gfs_common_tasks_before_fcst += ['ocnanalecen']
+                gdas_gfs_common_tasks_before_fcst += ['marineanlletkf', 'ocnanalecen']
             gdas_gfs_common_tasks_before_fcst += ['marineanlchkpt', 'marineanlfinal']
             if self.do_vrfy_oceanda:
                 gdas_gfs_common_tasks_before_fcst += ['ocnanalvrfy']