Initial commit

parm/atm/jcb-base.yaml.j2

@@ -9,12 +9,12 @@ app_path_observation_chronicle: "{{PARMgfs}}/gdas/jcb-gdas/observation_chronicle
 
 # Places where we deviate from the generic file name of a yaml
 # ------------------------------------------------------------
-#final_increment_file: final_increment_gaussian
-final_increment_file: atmosphere_final_increment_cubed_sphere
-output_ensemble_increments_file: atmosphere_output_ensemble_increments_cubed_sphere
+final_increment_file: atmosphere_final_increment_fms_nonrestart
+output_ensemble_increments_file: atmosphere_output_ensemble_increments_fms_nonrestart
 model_file: atmosphere_model_pseudo
 initial_condition_file: atmosphere_background  # Initial conditions for 4D apps is background
 background_error_file: "{{BERROR_YAML}}"
+posterior_output_file: atmosphere_posterior_output_cubed_sphere
 
 # Assimilation standard things (not prepended with model choice)
 # ----------------------------
@@ -24,11 +24,15 @@ bound_to_include: begin
 minimizer: DRPCG
 final_diagnostics_departures: anlmob
 final_prints_frequency: PT3H
+analysis_variables: &anvars [ua,va,t,ps,sphum,ice_wat,liq_wat,o3mr]
 number_of_outer_loops: 2
-analysis_variables: [ua,va,t,ps,sphum,ice_wat,liq_wat,o3mr]
 
 # Model things
 # ------------
+
+# Analysis variables
+atmosphere_analysis_variables: *anvars
+
 # Geometry
 atmosphere_layout_x: {{ layout_x | default(1, true) }}
 atmosphere_layout_y: {{ layout_y | default(1, true) }}
@@ -62,8 +66,6 @@ atmosphere_forecast_timestep: "{{ BKG_TSTEP }}"
 # Write final increment on Guassian grid in variational
 atmosphere_final_increment_prefix: "./anl/atminc."
 
-atmosphere_variational_history_prefix: "{{GPREFIX}}"
-
 # Observation things
 # ------------------
 observations: all_observations
@@ -125,7 +127,7 @@ inflation_mult: 1.0
 
 # Driver
 driver_update_obs_config_with_geometry_info: true
-driver_save_posterior_mean: false
+driver_save_posterior_mean: true
 driver_save_posterior_ensemble: false
 driver_save_prior_mean: false
 driver_save_posterior_mean_increment: false

utils/fv3jedi/CMakeLists.txt

@@ -3,3 +3,9 @@ ecbuild_add_executable( TARGET fv3jedi_fv3inc.x
                         SOURCES fv3jedi_fv3inc.cc fv3jedi_fv3inc.h)
 target_compile_features( fv3jedi_fv3inc.x PUBLIC cxx_std_17)
 target_link_libraries( fv3jedi_fv3inc.x PUBLIC NetCDF::NetCDF_CXX oops fv3jedi)
+
+# Ensemble increment recentering
+ecbuild_add_executable( TARGET fv3jedi_ecen.x
+                        SOURCES fv3jedi_ecen.cc fv3jedi_ecen.h)
+target_compile_features( fv3jedi_ecen.x PUBLIC cxx_std_17)
+target_link_libraries( fv3jedi_ecen.x PUBLIC NetCDF::NetCDF_CXX oops fv3jedi)

utils/fv3jedi/fv3jedi_ecen.cc

@@ -0,0 +1,8 @@
+#include "fv3jedi_ecen.h"
+#include "oops/runs/Run.h"
+
+int main(int argc,  char ** argv) {
+  oops::Run run(argc, argv);
+  gdasapp::ecen ecen;
+  return run.execute(ecen);
+}

utils/fv3jedi/fv3jedi_ecen.h

@@ -0,0 +1,118 @@
+#pragma once
+
+#include <iostream>
+#include <memory>
+#include <string>
+#include <vector>
+
+#include "eckit/config/LocalConfiguration.h"
+
+#include "fv3jedi/Geometry/Geometry.h"
+#include "fv3jedi/Increment/Increment.h"
+#include "fv3jedi/State/State.h"
+
+#include "oops/mpi/mpi.h"
+#include "oops/runs/Application.h"
+#include "oops/util/ConfigFunctions.h"
+#include "oops/util/DateTime.h"
+#include "oops/util/Duration.h"
+#include "oops/util/Logger.h"
+
+namespace gdasapp {
+
+  // Main application class
+  class ecen : public oops::Application {
+   public:
+    explicit ecen(const eckit::mpi::Comm & comm = oops::mpi::world())
+      : Application(comm) {}
+    static const std::string classname() {return "gdasapp::ecen";}
+
+    int execute(const eckit::Configuration & fullConfig) const {
+      // Get analysis parameters
+      std::vector<std::string> fcstHours;
+      std::string windowBeginStr;
+      fullConfig.get("forecast hours", fcstHours);
+      fullConfig.get("window begin", windowBeginStr);
+      const util::DateTime windowBegin(windowBeginStr);
+      const oops::Variables stateVars(fullConfig, "state variables");
+      const oops::Variables incrVars(fullConfig, "increment variables");
+
+      // Get geometry configurations
+      const eckit::LocalConfiguration bkgGeomConfig(fullConfig, "background geometry");
+      const eckit::LocalConfiguration incrGeomConfig(fullConfig, "increment geometry");
+      const eckit::LocalConfiguration anlEnsMeanGeomConfig(fullConfig, \
+                                                           "ensemble mean analysis geometry");
+      const eckit::LocalConfiguration incrCorGeomConfig(fullConfig, \
+                                                        "correction increment geometry");
+
+      // Setup geometries
+      const fv3jedi::Geometry incrGeom(incrGeomConfig, this->getComm());
+      const fv3jedi::Geometry bkgGeom(bkgGeomConfig, this->getComm());
+      const fv3jedi::Geometry anlEnsMeanGeom(anlEnsMeanGeomConfig, this->getComm());
+      const fv3jedi::Geometry incrCorGeom(incrCorGeomConfig, this->getComm());
+
+      // Get additions configuration
+      int nhrs = fcstHours.size();
+      std::vector<eckit::LocalConfiguration> additionsConfig;
+      eckit::LocalConfiguration additionsFromTemplateConfig(fullConfig, \
+                                                            "additions from template");
+      eckit::LocalConfiguration templateConfig(additionsFromTemplateConfig, "template");
+      std::string pattern;
+      additionsFromTemplateConfig.get("pattern", pattern);
+
+      for ( int ihrs = 0; ihrs < nhrs; ihrs++ ) {
+        eckit::LocalConfiguration thisAdditionsConfig(templateConfig);
+        util::seekAndReplace(thisAdditionsConfig, pattern, fcstHours[ihrs]);
+        additionsConfig.push_back(thisAdditionsConfig);
+      }
+
+      // Loops through forecast hours
+      for ( int ihrs = 0; ihrs < nhrs; ihrs++ ) {
+        // Get DateTime for forecast hour
+        const util::Duration fcstHour(3600*std::stoi(fcstHours[ihrs]) - 3*3600);
+        util::DateTime currentCycle = windowBegin + fcstHour;
+
+        // Get elements of individual additions configurations
+        const eckit::LocalConfiguration bkgConfig(additionsConfig[ihrs], "background");
+        const eckit::LocalConfiguration incrConfig(additionsConfig[ihrs], "increment");
+        const eckit::LocalConfiguration anlEnsMeanConfig(additionsConfig[ihrs], \
+                                                         "ensemble mean analysis");
+        const eckit::LocalConfiguration incrCorConfig(additionsConfig[ihrs], \
+                                                      "correction increment");
+
+        // Read background
+        fv3jedi::State xxBkg(bkgGeom, stateVars, currentCycle);
+        xxBkg.read(bkgConfig);
+
+        // Read increment
+        fv3jedi::Increment dx(incrGeom, incrVars, currentCycle);
+        dx.read(incrConfig);
+
+        // Read ensemble mean analysis
+        fv3jedi::State xxAnlEnsMean(anlEnsMeanGeom, incrVars, currentCycle);
+        xxAnlEnsMean.read(anlEnsMeanConfig);
+
+        // Compute deterministic analysis
+        fv3jedi::State xxAnlDet(bkgGeom, xxBkg);
+        xxAnlDet += dx;
+
+        // Interpolate full resolution analysis to ensemble resolution and then change variables
+        fv3jedi::State xxAnlDetEnsRes(incrCorGeom, fv3jedi::State(incrVars, xxAnlDet));
+
+        // Compute correction increment
+        fv3jedi::Increment dxCor(incrCorGeom, incrVars, xxBkg.validTime());
+        dxCor.diff(xxAnlEnsMean, xxAnlDetEnsRes);
+
+        // Write correction increment
+        dxCor.write(incrCorConfig);
+      }
+
+      return 0;
+    }
+
+   private:
+    std::string appname() const {
+      return "gdasapp::ecen";
+    }
+  };
+}  // namespace gdasapp