add StdDev variable change (#60)

src/umdsst/CMakeLists.txt

@@ -31,5 +31,6 @@ add_subdirectory(Fields)
 add_subdirectory(Geometry)
 add_subdirectory(GetValues)
 add_subdirectory(Increment)
+add_subdirectory(LinearVariableChange)
 add_subdirectory(ModelAux)
 add_subdirectory(State)

src/umdsst/Increment/Increment.cc

@@ -175,6 +175,16 @@ namespace umdsst {
       fd(i, 0) *= fd_rhs(i, 0);
   }
 
+// ----------------------------------------------------------------------------
+
+  void Increment::schur_product_with_inv(const Increment &rhs ) {
+    auto fd = make_view<double, 2>(atlasFieldSet_->field(0));
+    auto fd_rhs = make_view<double, 2>(rhs.atlasFieldSet()->field(0));
+
+    const int size = geom_->atlasFunctionSpace()->size();
+    for (int i = 0; i < size; i++)
+      fd(i, 0) *= 1.0 / fd_rhs(i, 0);
+  }
 // ----------------------------------------------------------------------------
 
   void Increment::zero() {

src/umdsst/Increment/Increment.h

@@ -59,6 +59,7 @@ namespace umdsst {
     void ones();
     void random();
     void schur_product_with(const Increment &);
+    void schur_product_with_inv(const Increment &);
     void zero();
     void zero(const util::DateTime &);
 

src/umdsst/LinearVariableChange/CMakeLists.txt

@@ -0,0 +1,4 @@
+umdsst_target_sources(
+    StdDev.cc
+    StdDev.h
+)

src/umdsst/LinearVariableChange/StdDev.cc

@@ -0,0 +1,74 @@
+/*
+ * (C) Copyright 2021-2021 UCAR
+ *
+ * This software is licensed under the terms of the Apache Licence Version 2.0
+ * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+ */
+
+#include "umdsst/LinearVariableChange/StdDev.h"
+#include "umdsst/State/State.h"
+#include "umdsst/Geometry/Geometry.h"
+#include "umdsst/Traits.h"
+#include "umdsst/Increment/Increment.h"
+
+#include "oops/interface/LinearVariableChange.h"
+#include "oops/util/abor1_cpp.h"
+
+namespace umdsst {
+
+// ----------------------------------------------------------------------------
+
+StdDev::StdDev(const State &bkg, const State &traj, const Geometry &geom,
+               const eckit::Configuration &conf) {
+  stddev_.reset(new Increment(geom, bkg.variables(), bkg.validTime()));
+  stddev_->ones();
+
+  if ( conf.has("fixed") ) {
+    // a single global fixed value
+    double val;
+    conf.get("fixed", val);
+    *stddev_ *= val;
+  } else {
+    util::abor1_cpp("StdDev::StdDev() no standard deviation "
+                    "method specified", __FILE__, __LINE__);
+  }
+
+  // TODO(someone) place the real horizontally varying background error
+  // standard deviation in stddev_
+}
+
+// ----------------------------------------------------------------------------
+
+void StdDev::multiply(const Increment &dxin, Increment &dxout) const {
+  dxout = dxin;
+  dxout.schur_product_with(*stddev_);
+}
+
+// ----------------------------------------------------------------------------
+
+void StdDev::multiplyAD(const Increment &dxin, Increment &dxout) const {
+  multiply(dxin, dxout);
+}
+
+// ----------------------------------------------------------------------------
+
+void StdDev::multiplyInverse(const Increment &dxin, Increment &dxout) const {
+  dxout = dxin;
+  dxout.schur_product_with_inv(*stddev_);
+}
+
+// ----------------------------------------------------------------------------
+
+void StdDev::multiplyInverseAD(const Increment &dxin, Increment &dxout) const {
+  multiplyInverse(dxin, dxout);
+}
+
+// ----------------------------------------------------------------------------
+
+oops::LinearVariableChangeMaker<Traits,
+                                oops::LinearVariableChange<Traits, StdDev> >
+  makerLinearVariableChangeStdDev_("umdsstStdDev");
+
+// ----------------------------------------------------------------------------
+
+}  // namespace umdsst

src/umdsst/LinearVariableChange/StdDev.h

@@ -0,0 +1,47 @@
+/*
+ * (C) Copyright 2021-2021 UCAR
+ *
+ * This software is licensed under the terms of the Apache Licence Version 2.0
+ * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+ */
+
+#ifndef UMDSST_LINEARVARIABLECHANGE_STDDEV_H_
+#define UMDSST_LINEARVARIABLECHANGE_STDDEV_H_
+
+#include <memory>
+#include <string>
+
+#include "oops/util/Printable.h"
+
+// Forward Declaration
+namespace eckit {
+  class Configuration;
+}
+namespace umdsst {
+  class Geometry;
+  class Increment;
+  class State;
+}
+
+namespace umdsst {
+
+class StdDev: public util::Printable {
+ public:
+  static const std::string classname() {return "umdsst:StdDev";}
+
+  explicit StdDev(const State &, const State &, const Geometry &,
+                  const eckit::Configuration &);
+
+  void multiply(const Increment &, Increment &) const;
+  void multiplyInverse(const Increment &, Increment &) const;
+  void multiplyAD(const Increment &, Increment &) const;
+  void multiplyInverseAD(const Increment &, Increment &) const;
+
+ private:
+  void print(std::ostream &) const override {}
+  std::unique_ptr<Increment> stddev_;
+};
+
+}  // namespace umdsst
+
+#endif  // UMDSST_LINEARVARIABLECHANGE_STDDEV_H_

test/CMakeLists.txt

@@ -5,6 +5,7 @@ list( APPEND umdsst_test_input
   testinput/hofx_nomodel.yml
   testinput/increment.yml
   testinput/lineargetvalues.yml
+  testinput/linearvarchange_stddev.yml
   testinput/modelaux.yml
   testinput/state.yml
   testinput/dirac.yml
@@ -154,6 +155,13 @@ endfunction()
      MPI     ${MPI_PES}
      LIBS    umdsst )
 
+   ecbuild_add_test(
+    TARGET  test_umdsst_linearvarchange_stddev
+    SOURCES executables/TestLinearVariableChange.cc
+    ARGS    testinput/linearvarchange_stddev.yml
+    MPI     ${MPI_PES}
+    LIBS    umdsst )
+
    ecbuild_add_test(
      TARGET  test_umdsst_errorcovariance
      SOURCES executables/TestErrorCovariance.cc

test/executables/TestLinearVariableChange.cc

@@ -0,0 +1,17 @@
+/*
+ * (C) Copyright 2021-2021 UCAR.
+ *
+ * This software is licensed under the terms of the Apache Licence Version 2.0
+ * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+ */
+
+#include "umdsst/Traits.h"
+#include "oops/runs/Run.h"
+#include "test/interface/LinearVariableChange.h"
+
+int main(int argc,  char ** argv) {
+  oops::Run run(argc, argv);
+  test::LinearVariableChange<umdsst::Traits> tests;
+  return run.execute(tests);
+}
+

test/testinput/dirac.yml

@@ -25,6 +25,11 @@ background error:
     prefix: bump_sst
     mpicom: 2
     strategy: specific_univariate
+  variable changes:
+  - variable change: umdsstStdDev
+    input variables: *vars
+    output variables: *vars
+    fixed: 1.0
 
 dirac:
   date: 2010-01-01T12:00:00Z

test/testinput/linearvarchange_stddev.yml

@@ -0,0 +1,19 @@
+geometry:
+  grid: S360x180
+  domain:
+    type: global
+    west: -180
+  landmask:
+    filename: Data/landmask.nc
+
+background:
+  state variables: &vars [sea_surface_temperature]
+  date: 2018-04-15T00:00:00Z
+
+linear variable change tests:
+- variable change: umdsstStdDev
+  tolerance inverse: 1e-12
+  test inverse: 1
+  input variables: *vars
+  output variables: *vars
+  fixed: 2.0

test/testinput/var.yml

@@ -30,6 +30,11 @@ cost function:
       prefix: bump_sst
       mpicom: 2
       strategy: specific_univariate
+    variable changes:
+    - variable change: umdsstStdDev
+      input variables: *vars
+      output variables: *vars
+      fixed: 1.0
 
   observations:
   - obs error: