Initial attempt to implement an increment interpolator unit test.

Author: frld

.github/workflows/ci.yml

@@ -37,7 +37,7 @@ jobs:
         with:
           path: ci/saber
           repository: JCSDA-internal/saber
-          token: ${{ secrets.GH_PAT }}
+          token: ${{ secrets.GHCR_PAT }}
 
       - name: checkout ioda
         uses: actions/checkout@v3

src/orca-jedi/increment/Increment.cc

@@ -37,6 +37,11 @@
 #include "atlas/mesh.h"
 #include "atlas-orca/grid/OrcaGrid.h"
 
+#include <boost/uuid/uuid.hpp>            // uuid class DJL
+#include <boost/uuid/uuid_generators.hpp> // generators DJL
+#include <boost/uuid/uuid_io.hpp>         // streaming operators etc. DJL
+
+
 namespace orcamodel {
 
 // -----------------------------------------------------------------------------
@@ -264,10 +269,17 @@ void Increment::setval(const double & val) {
           if (!has_mv || (has_mv && !mv(field_view(j, k)))) {
             field_view(j, k) = val;
           }
+        } else {
+          field_view(j, k) = -99;            /// DJL
         }
       }
     }
   }
+  
+  // DJL write debug fields to file
+    boost::uuids::uuid uuid = boost::uuids::random_generator()();    
+    writeFieldsToFile("incsetval"+ boost::uuids::to_string(uuid) +".nc", *geom_, validTime(), incrementFields_);
+
 
   oops::Log::trace() << "Increment(ORCA)::setval done" << std::endl;
 }

src/orca-jedi/interpolator/Interpolator.cc

@@ -31,6 +31,13 @@
 
 #include "orca-jedi/interpolator/Interpolator.h"
 
+#include <boost/uuid/uuid.hpp>            // uuid class DJL
+#include <boost/uuid/uuid_generators.hpp> // generators DJL
+#include <boost/uuid/uuid_io.hpp>         // streaming operators etc. DJL
+
+#include "orca-jedi/utilities/IOUtils.h" // DJL
+
+
 namespace eckit {
 class Configuration;
 }
@@ -178,6 +185,11 @@ template void Interpolator::executeInterpolation<float>(
     const std::vector<bool> & mask,
     std::vector<double>::iterator& iter) const;
 
+/// \brief Interpolate from model space to observation space
+/// \param vars Oops variables
+/// \param inc Increment object (input)
+/// \param mask Mask vector in observation space
+/// \param result Result (output) vector in observation space
 void Interpolator::apply(const oops::Variables& vars, const Increment& inc,
            const std::vector<bool> & mask,
            std::vector<double>& result) const {
@@ -230,6 +242,11 @@ void Interpolator::apply(const oops::Variables& vars, const Increment& inc,
   }
 }
 
+/// \brief Interpolate from observation space to model space 
+/// \param vars Oops variables
+/// \param inc Increment object (output)
+/// \param mask Mask (observation space) vector
+/// \param resultin Values (observation space) vector (input)
 void Interpolator::applyAD(const oops::Variables& vars, Increment& inc,
              const std::vector<bool> & mask,
              const std::vector<double> & resultin) const
@@ -242,7 +259,7 @@ void Interpolator::applyAD(const oops::Variables& vars, Increment& inc,
 
   oops::Log::debug() << "DJL ** Interpolator::applyAD this needs checking **" << std::endl;
 
-// ** Not sure what I'm doing yet - Just trying to do the opposite of applyAD
+// ** Trying to do the opposite of apply DJL
 
   const size_t nvars = vars.size();
 
@@ -263,13 +280,15 @@ void Interpolator::applyAD(const oops::Variables& vars, Increment& inc,
     inc.geometry()->variableSizes(vars);
   size_t nvals = 0;
 
-//    boost::uuids::uuid uuid = boost::uuids::random_generator()();    
-//    std::shared_ptr<const Geometry> geom = inc.geometry();
-//    writeGenFieldsToFile("applyADpre"+ boost::uuids::to_string(uuid) +".nc", *geom, inc.validTime(), inc.incrementFields());
+// DJL write debug fields to file
+    boost::uuids::uuid uuid = boost::uuids::random_generator()();    
+    std::shared_ptr<const Geometry> geom = inc.geometry();
+    writeFieldsToFile("applyADpre"+ boost::uuids::to_string(uuid) +".nc", *geom, inc.validTime(), inc.incrementFields());
 
   for (size_t jvar=0; jvar < nvars; ++jvar) nvals += nlocs_ * varSizes[jvar];
 //    result.resize(nvals);
 
+
   std::size_t out_idx = 0;
   for (size_t jvar=0; jvar < nvars; ++jvar) {
     oops::Log::debug() << "DJL ** jvar " << jvar << " " << nvars
@@ -301,21 +320,19 @@ void Interpolator::applyAD(const oops::Variables& vars, Increment& inc,
       }
     }
 
+// halo exchange update ghost points DJL
+    geom->functionSpace().haloExchange(inc.incrementFields()[gv_varname]);
 
 //      atlas::field::MissingValue mv(inc.incrementFields()[gv_varname]);
 //      bool has_mv = static_cast<bool>(mv);
     interpolator_.execute_adjoint(inc.incrementFields()[gv_varname], tgt_field);
 
   }    // jvar
 
-  // DJL write to file
-
-//  boost::uuids::uuid uuid = boost::uuids::random_generator()();    
-//  std::shared_ptr<const Geometry> geom = inc.geometry();
-//    std::ostringstream out;
-//    out << std::setfill('0') << std::setw(6) << uuid;
-//  writeGenFieldsToFile("applyAD"+ boost::uuids::to_string(uuid) +".nc", *geom, inc.validTime(), inc.incrementFields());
-//    fileCounter_++;
+// DJL write debug fields to file
+//    boost::uuids::uuid uuid = boost::uuids::random_generator()();    
+//    std::shared_ptr<const Geometry> geom = inc.geometry();
+    writeFieldsToFile("applyAD"+ boost::uuids::to_string(uuid) +".nc", *geom, inc.validTime(), inc.incrementFields());
 
   oops::Log::trace() << "orcamodel::Interpolator::applyAD done "
                      << std::endl;

src/tests/orca-jedi/test_interpolator.cc

@@ -34,6 +34,7 @@ struct InterpTestSettingsFixture {
  public:
   eckit::LocalConfiguration geometry_config;
   eckit::LocalConfiguration state_config;
+  eckit::LocalConfiguration increment_config;
   eckit::LocalConfiguration interpolator_config;
   size_t nlocs, nlevs;
   std::vector<double> lons;
@@ -87,6 +88,9 @@ CASE("test  interpolator") {
     settings_map["ORCA2_T"].state_config.set("nemo field file", "../Data/orca2_t_nemo.nc");
     settings_map["ORCA2_T"].state_config.set("variance field file", "../Data/orca2_t_bkg_var.nc");
 
+    settings_map["ORCA2_T"].increment_config.set("date", "2021-06-30T00:00:00Z");
+    settings_map["ORCA2_T"].increment_config.set("output path", "../testoutput/orca2_t_increment_interptest.nc");
+
     settings_map["ORCA2_T"].lons = std::vector<double>{0, 120, 270};
     settings_map["ORCA2_T"].lats = std::vector<double>{88, 0, 30};
 
@@ -137,6 +141,9 @@ CASE("test  interpolator") {
     settings_map["AMM1"].state_config.set("nemo field file", "../Data/amm1_nemo.nc");
     settings_map["AMM1"].state_config.set("variance field file", "../Data/amm1_nemo.nc");
 
+    settings_map["AMM1"].increment_config.set("date", "2021-06-30T00:00:00Z");
+    settings_map["AMM1"].increment_config.set("output path", "../testoutput/amm1_t_increment_interptest.nc");
+
     settings_map["AMM1"].lons = std::vector<double>{-17.5, -6.78, -16.1};
     settings_map["AMM1"].lats = std::vector<double>{58.16, 58.91, 63.55};
 
@@ -209,6 +216,40 @@ CASE("test  interpolator") {
         EXPECT(std::abs(vals[i] - settings.vol_values[i]) < ATOL);
       }
     }
+
+    
+    SECTION("test " + key + " interpolator.apply/applyAD with increment") {
+
+      OrcaIncrementParameters incrementParams;
+      incrementParams.validateAndDeserialize(settings.increment_config);
+
+      Increment increment(geometry, settings.surf_vars, incrementParams.date);
+      increment.ones();  
+
+      // two variables at n locations
+      std::vector<double> vals(2*settings.nlocs);
+      std::vector<bool> mask(settings.nlocs, true);
+
+      // increment -> observation space
+      interpolator.apply(settings.surf_vars, increment, mask, vals);
+
+      for (size_t i=0; i < settings.surf_values.size(); ++i) {
+        std::cout << "vals[" << i << "] " << std::setprecision(12) << vals[i]
+                  << " kgo_values[" << i << "] " << settings.surf_values[i] << std::endl;
+      }
+      //for (size_t i=0; i < settings.surf_values.size(); ++i) {
+      //  EXPECT(std::abs(vals[i] - settings.surf_values[i]) < ATOL);
+      //}
+      
+      Increment incrementout(geometry, settings.surf_vars, incrementParams.date);
+      // observation space -> increment
+      interpolator.applyAD(settings.surf_vars, incrementout, mask, vals);
+      
+      incrementout.write(incrementParams);
+            
+    }
+    
+    
   }
 }