Bugfix/sla do not assimilate

src/nemo-feedback/NemoFeedback.cc

@@ -6,8 +6,8 @@
 
 #include "nemo-feedback/NemoFeedback.h"
 
-#include <string>
 #include <string_view>
+#include <string>
 
 #include <algorithm>
 #include <utility>
@@ -46,35 +46,6 @@ namespace nemo_feedback {
 constexpr std::string_view defaultDepthGroup{"MetaData"};
 constexpr std::string_view defaultDepthVariable{"depthBelowWaterSurface"};
 
-/// \brief create whole report variables from a diagnostic flag
-void wholeReportQCFromDiagnosticFlags(const
-    NemoFeedbackDataCreator& creator, const std::string& group, const
-    std::string& name, feedback_io::Data<int32_t>& wholeReportQCData) {
-  if (wholeReportQCData.n_obs() == 0) {
-    wholeReportQCData = feedback_io::Data<int32_t>(creator.indexer(),
-        std::vector<int32_t>(creator.indexer()->n_source_data(), 0));
-  }
-  feedback_io::Data<int32_t> QCData =
-    creator.create(group, name,
-                   ufo::DiagnosticFlag(0), 4, 1);
-  for (size_t iProfile = 0;
-      iProfile < QCData.n_obs(); ++iProfile) {
-    size_t nBadObs = 0;
-    for (size_t iLevel = 0;
-        iLevel < QCData.length(iProfile);
-        ++iLevel) {
-      if (QCData(iProfile, iLevel) == 4) {
-        ++nBadObs;
-      }
-    }
-    if (wholeReportQCData[iProfile] == 0 ||
-        wholeReportQCData[iProfile] == 4) {
-      wholeReportQCData[iProfile] =
-        (nBadObs == QCData.length(iProfile) ? 4 : 1);
-    }
-  }
-}
-
 
 NemoFeedback::NemoFeedback(
     ioda::ObsSpace & obsdb,
@@ -222,10 +193,11 @@ template <typename T>
 void NemoFeedback::write_all_data(feedback_io::Writer<T>& writer,
                                   const NemoFeedbackDataCreator& creator) const
 {
-  feedback_io::Data<int32_t> wholeReportQCData(creator.indexer(),
-      std::vector<int32_t>(creator.indexer()->n_source_data(), 0));
-  feedback_io::Data<int32_t> wholeReportPositionQCData;
-  feedback_io::Data<int32_t> wholeReportTimeQCData;
+  feedback_io::Data<feedback_io::QC::Level> wholeReportQCData(creator.indexer(),
+      std::vector<feedback_io::QC::Level>(creator.indexer()->n_source_data(),
+        feedback_io::QC::Level::None));
+  feedback_io::Data<feedback_io::QC::Level> wholeReportPositionQCData;
+  feedback_io::Data<feedback_io::QC::Level> wholeReportTimeQCData;
   std::vector<ufo::DiagnosticFlag> do_not_assimilate;
   for (const NemoFeedbackVariableParameters& nemoVariableParams :
         parameters_.variables.value()) {
@@ -270,70 +242,78 @@ void NemoFeedback::write_all_data(feedback_io::Writer<T>& writer,
 
     // Quality control rank variables
     if (obsdb_.has("DiagnosticFlags/FinalReject", ufo_name)) {
-
-      feedback_io::Data<int32_t> variableFinalQCData =
+      feedback_io::Data<feedback_io::QC::Level> variableFinalQCData =
         creator.create("DiagnosticFlags/FinalReject", ufo_name,
-                       ufo::DiagnosticFlag(0), 4, 1);
+            ufo::DiagnosticFlag(0), feedback_io::QC::Level::Bad,
+            feedback_io::QC::Level::Good);
 
-      // Add do not assimilate flag if required.
+      // Add do not assimilate flag if required to the final QC information.
       if (obsdb_.has("DiagnosticFlags/DoNotAssimilate", ufo_name)) {
-        feedback_io::Data<int32_t> variableDoNotAssimilateData(
+        feedback_io::Data<feedback_io::QC::Level> variableDoNotAssimilateData(
             creator.create("DiagnosticFlags/DoNotAssimilate", ufo_name,
-              ufo::DiagnosticFlag(0), 128, 0));
+              ufo::DiagnosticFlag(0), feedback_io::QC::Level::DoNotAssimilate,
+              feedback_io::QC::Level::None));
         for (size_t iProfile = 0;
             iProfile < variableDoNotAssimilateData.n_obs(); ++iProfile) {
           for (size_t iLevel = 0;
               iLevel < variableDoNotAssimilateData.length(iProfile);
               ++iLevel) {
-            variableFinalQCData(iProfile, iLevel) +=
-              variableDoNotAssimilateData(iProfile, iLevel);
+            if (variableDoNotAssimilateData(iProfile, iLevel) ==
+                feedback_io::QC::Level::DoNotAssimilate) {
+              variableFinalQCData(iProfile, iLevel) =
+               static_cast<feedback_io::QC::Level>(
+                static_cast<int32_t>(variableFinalQCData(iProfile, iLevel)) +
+                static_cast<int32_t>(feedback_io::QC::Level::DoNotAssimilate));
+            }
           }
         }
       }
 
-      // Per-profile quality rank
-      feedback_io::Data<int32_t> wholeVariableQCData(creator.indexer(),
-          std::vector<int32_t>(creator.indexer()->n_source_data(), 0));
-      for (size_t iProfile = 0;
-          iProfile < variableFinalQCData.n_obs(); ++iProfile) {
-        size_t nBadObs = 0;
-        for (size_t iLevel = 0;
-            iLevel < variableFinalQCData.length(iProfile);
-            ++iLevel) {
-          if (variableFinalQCData(iProfile, iLevel) == 4) {
-            ++nBadObs;
-          }
-        }
-        wholeVariableQCData[iProfile] =
-          (nBadObs == variableFinalQCData.length(iProfile) ? 4 : 1);
-      }
+      // set per-profile quality rank for the variable
+      feedback_io::Data<feedback_io::QC::Level> wholeVariableQCData(
+          creator.indexer(), std::vector<feedback_io::QC::Level>(
+            creator.indexer()->n_source_data(), feedback_io::QC::Level::None));
+      feedback_io::wholeReportFromPerProfile(variableFinalQCData,
+          wholeVariableQCData);
+      // update set per-profile quality rank for the whole report
+      feedback_io::wholeReportFromPerProfile(variableFinalQCData,
+          wholeReportQCData);
 
       writer.write_variable_surf_qc(nemo_name + "_QC",
         wholeVariableQCData);
       writer.write_variable_level_qc(nemo_name + "_LEVEL_QC",
           variableFinalQCData);
 
-      // Whole Observation report QC
+      // Whole Observation Position report QC
       const std::string positionQCGroup("DiagnosticFlags/PositionReject");
       if (obsdb_.has(positionQCGroup, ufo_name)) {
-        wholeReportQCFromDiagnosticFlags(creator, positionQCGroup, ufo_name,
+        feedback_io::Data<feedback_io::QC::Level> QCData = creator.create(
+            positionQCGroup, ufo_name, ufo::DiagnosticFlag(0),
+            feedback_io::QC::Level::Bad, feedback_io::QC::Level::Good);
+        if (wholeReportPositionQCData.n_obs() == 0) {
+          wholeReportPositionQCData = feedback_io::Data<feedback_io::QC::Level>(
+              creator.indexer(), std::vector<feedback_io::QC::Level>(
+                creator.indexer()->n_source_data(),
+                feedback_io::QC::Level::None));
+        }
+        feedback_io::wholeReportFromPerProfile(QCData,
             wholeReportPositionQCData);
       }
 
+      // Whole Observation time report QC
       const std::string timeQCGroup("DiagnosticFlags/TimeReject");
       if (obsdb_.has(timeQCGroup, ufo_name)) {
-        wholeReportQCFromDiagnosticFlags(creator, timeQCGroup, ufo_name,
-            wholeReportTimeQCData);
-      }
-
-      {
-        for (size_t iProfile = 0; iProfile < variableData.n_obs(); ++iProfile) {
-          if (wholeReportQCData[iProfile] == 0 ||
-              wholeReportQCData[iProfile] == 4) {
-            wholeReportQCData[iProfile] =
-              (4 == wholeVariableQCData[iProfile] ? 4 : 1);
-          }
+        feedback_io::Data<feedback_io::QC::Level> QCData =
+          creator.create(timeQCGroup, ufo_name, ufo::DiagnosticFlag(0),
+              feedback_io::QC::Level::Bad, feedback_io::QC::Level::Good);
+        if (wholeReportTimeQCData.n_obs() == 0) {
+          wholeReportTimeQCData = feedback_io::Data<feedback_io::QC::Level>(
+              creator.indexer(), std::vector<feedback_io::QC::Level>(
+                creator.indexer()->n_source_data(),
+                feedback_io::QC::Level::None));
         }
+        feedback_io::wholeReportFromPerProfile(QCData,
+            wholeReportTimeQCData);
       }
     }
 
@@ -356,9 +336,9 @@ void NemoFeedback::write_all_data(feedback_io::Writer<T>& writer,
   const std::string depthVariable = parameters_.depthVariable.value()
     .value_or(static_cast<std::string>(defaultDepthVariable));
   if (obsdb_.has(depthQCGroup, depthVariable)) {
-    feedback_io::Data<int32_t> depthQCData(
-        creator.create(depthQCGroup, depthVariable,
-                       ufo::DiagnosticFlag(0), 4, 1));
+    feedback_io::Data<feedback_io::QC::Level> depthQCData(
+        creator.create(depthQCGroup, depthVariable, ufo::DiagnosticFlag(0),
+          feedback_io::QC::Level::Bad, feedback_io::QC::Level::Good));
     writer.write_variable_level_qc("DEPTH_QC", depthQCData);
   }
   if (wholeReportPositionQCData.n_obs() != 0) {
@@ -424,7 +404,6 @@ std::tuple<feedback_io::Data<std::string>, feedback_io::Data<std::string>>
     stationIdentificationAvailable = true;
   }
 
-  const int32_t buoyIDMissingValueInt = util::missingValue<int32_t>();
   constexpr size_t buoyIDWidth = 8;
   const std::string missingStringFeedback =
     feedback_io::typeToFill::value<std::string>();

src/nemo-feedback/NemoFeedbackDataCreator.cc

@@ -144,14 +144,15 @@ feedback_io::Data<std::string> NemoFeedbackDataCreator::create_from_obsdb(
   return feedback_io::Data<std::string>(indexer_, std::move(data));
 }
 
-feedback_io::Data<int32_t> NemoFeedbackDataCreator::create_from_obsdb(const
-    std::string& obsGroup, const std::string& ufoName, const
-    ufo::DiagnosticFlag TypeInstance, int32_t whenTrue,
-    int32_t whenFalse) const {
+feedback_io::Data<feedback_io::QC::Level>
+ NemoFeedbackDataCreator::create_from_obsdb(
+    const std::string& obsGroup, const std::string& ufoName, const
+    ufo::DiagnosticFlag TypeInstance, feedback_io::QC::Level whenTrue,
+    feedback_io::QC::Level whenFalse) const {
   oops::Log::trace() << NemoFeedbackDataCreator::className()
                      << ":create_from_obsdb ufo::DiagnosticFlag "
                      << obsGroup << "/" << ufoName << std::endl;
-  std::vector<int32_t> data;
+  std::vector<feedback_io::QC::Level> data;
   std::vector<ufo::DiagnosticFlag> flagData(obsdb_.nlocs(), 0);
   obsdb_.get_db(obsGroup, ufoName, flagData);
   for (ufo::DiagnosticFlag flag : flagData) {
@@ -161,7 +162,7 @@ feedback_io::Data<int32_t> NemoFeedbackDataCreator::create_from_obsdb(const
     throw eckit::BadValue(NemoFeedbackDataCreator::className()
         + ":create_from_obsdb ufo::DiagnosticFlag no data.", Here());
 
-  return feedback_io::Data<int32_t>(indexer_, std::move(data));
+  return feedback_io::Data<feedback_io::QC::Level>(indexer_, std::move(data));
 }
 
 template<typename T>
@@ -304,9 +305,10 @@ NemoFeedbackDataCreator::create_altimeter_IDs() const {
         feedback_io::Data<std::string>(indexer_, std::move(station_types)));
 }
 
-feedback_io::Data<int32_t> NemoFeedbackDataCreator::create(const std::string&
-    obsGroup, const std::string& ufoName, const ufo::DiagnosticFlag
-    typeInstance, int32_t whenTrue, int32_t whenFalse) const {
+feedback_io::Data<feedback_io::QC::Level> NemoFeedbackDataCreator::create(
+    const std::string& obsGroup, const std::string& ufoName, const
+    ufo::DiagnosticFlag typeInstance, feedback_io::QC::Level whenTrue,
+    feedback_io::QC::Level whenFalse) const {
   if (obsGroup == "hofx" || obsGroup == "HofX") {
     std::ostringstream err_stream;
     err_stream << NemoFeedbackDataCreator::className()

src/nemo-feedback/NemoFeedbackDataCreator.h

@@ -63,11 +63,11 @@ class NemoFeedbackDataCreator {
     create_altimeter_IDs() const;
 
   /// \brief create DiagnosticFlags feedbackData
-  feedback_io::Data<int32_t> create(const std::string& obsGroup,
+  feedback_io::Data<feedback_io::QC::Level> create(const std::string& obsGroup,
       const std::string& ufoName,
       const ufo::DiagnosticFlag typeInstance,
-      const int32_t whenTrue,
-      const int32_t whenFalse) const;
+      const feedback_io::QC::Level whenTrue,
+      const feedback_io::QC::Level whenFalse) const;
 
   /// \brief create feedbackData
   template<typename T>
@@ -91,11 +91,12 @@ std::shared_ptr<feedback_io::DataIndexer> indexer() const {return indexer_;}
       size_t width, bool leftJustify = false) const;
 
   /// \brief create DiagnosticFlags feedbackData from obsdb
-  feedback_io::Data<int32_t> create_from_obsdb(const std::string& obsGroup,
+  feedback_io::Data<feedback_io::QC::Level> create_from_obsdb(
+      const std::string& obsGroup,
       const std::string& ufoName,
       const ufo::DiagnosticFlag TypeInstance,
-      int32_t whenTrue,
-      int32_t whenFalse) const;
+      feedback_io::QC::Level whenTrue,
+      feedback_io::QC::Level whenFalse) const;
 
   /// \brief create feedbackData from obsdb
   template<typename T>

src/nemo-feedback/feedback_io/Data.cc

@@ -88,7 +88,71 @@ std::vector<T> Data<T>::raw_profile(size_t iProfile) const {
 template class Data<std::string>;
 template class Data<size_t>;
 template class Data<int32_t>;
+template class Data<feedback_io::QC::Level>;
 template class Data<float>;
 template class Data<double>;
+
+/// \brief update the whole report QC information for a profile based on the
+///        current whole variable QC information
+feedback_io::QC::Level wholeReportUpdate(feedback_io::QC::Level current, const
+    size_t length, const size_t nGoodObs, const size_t nBadObs, const size_t
+    nBadDoNotAssimilate, const size_t nDoNotAssimilate) {
+  // If good obs were previously found, the report is good
+  if (current == QC::Level::Good) {
+    return QC::Level::Good;
+  }
+  // Upgrade any reports if they contain any good obs (goodDoNotAssimilate ->
+  // good)
+  if (nGoodObs > 0) {
+    return QC::Level::Good;
+  }
+  // Check currently bad, do-not-assimilate or not-yet-checked reports and
+  // downgrade where necessary
+  if (nBadObs == length) {
+    return QC::Level::Bad;
+  } else if (nBadDoNotAssimilate == length) {
+    // No good observations, all bad-do-not-assimilate
+    return QC::Level::BadDoNotAssimilate;
+  } else if (nDoNotAssimilate == length) {
+    // No good observations, some mix of good and bad do-not-assimilate types
+    //     -> GoodDoNotAssimilate
+    return QC::Level::GoodDoNotAssimilate;
+  } else {
+    // No good observations, some mix of bad and do-not-assimilate -> Bad
+    return QC::Level::Bad;
+  }
+}
+
+void wholeReportFromPerProfile(const Data<feedback_io::QC::Level>& QCData,
+    Data<feedback_io::QC::Level>& wholeReportQCData) {
+  for (size_t iProfile = 0;
+      iProfile < QCData.n_obs(); ++iProfile) {
+    size_t nGoodObs = 0;
+    size_t nBadObs = 0;
+    size_t nBadDoNotAssimilate = 0;
+    size_t nDoNotAssimilate = 0;
+    for (size_t iLevel = 0;
+        iLevel < QCData.length(iProfile);
+        ++iLevel) {
+      if (QCData(iProfile, iLevel) == QC::Level::Bad) {
+        ++nBadObs;
+      } else if (QCData(iProfile, iLevel) ==
+                 QC::Level::BadDoNotAssimilate) {
+        ++nBadDoNotAssimilate;
+      }
+      if (QCData(iProfile, iLevel) == QC::Level::BadDoNotAssimilate ||
+          QCData(iProfile, iLevel) == QC::Level::GoodDoNotAssimilate ||
+          QCData(iProfile, iLevel) == QC::Level::DoNotAssimilate) {
+        ++nDoNotAssimilate;
+      }
+      if (QCData(iProfile, iLevel) == QC::Level::Good) {
+        ++nGoodObs;
+      }
+    }
+    wholeReportQCData[iProfile] = wholeReportUpdate(
+        wholeReportQCData[iProfile], QCData.length(iProfile), nGoodObs,
+        nBadObs, nBadDoNotAssimilate, nDoNotAssimilate);
+  }
+}
 }  // namespace feedback_io
 }  // namespace nemo_feedback

src/nemo-feedback/feedback_io/Data.h

@@ -10,10 +10,12 @@
 
 #include <nemo-feedback/feedback_io/DataIndexer.h>
 
-#include <vector>
 #include <string>
+#include <vector>
 #include <memory>
 
+#include <nemo-feedback/feedback_io/Utils.h>
+
 namespace nemo_feedback {
 namespace feedback_io {
 
@@ -119,5 +121,9 @@ class Data {
   std::shared_ptr<DataIndexer> indexer_;
   std::shared_ptr<std::vector<T>> data_;
 };
+
+/// \brief create whole report data from per-profile QC data
+void wholeReportFromPerProfile(const Data<feedback_io::QC::Level>& QCData,
+    Data<feedback_io::QC::Level>& wholeReportQCData);
 }  // namespace feedback_io
 }  // namespace nemo_feedback

src/nemo-feedback/feedback_io/Utils.cc

@@ -20,6 +20,12 @@
 #include "oops/util/Logger.h"
 #include "oops/util/Duration.h"
 
+#define DOUBLE_FILLVALUE 99999.0
+#define FLOAT_FILLVALUE 99999.0
+#define INT32_FILLVALUE 0
+#define STRING_FILLVALUE "MISSING "
+
+
 namespace nemo_feedback {
 namespace feedback_io {
 
@@ -41,6 +47,9 @@ template <> const float  typeToFill::value<float>() {
 template <> const int32_t  typeToFill::value<int32_t>() {
   return INT32_FILLVALUE;
 }
+template <> const QC::Level  typeToFill::value<QC::Level>() {
+  return QC::Level::None;
+}
 template <> const std::string  typeToFill::value<std::string>() {
   return STRING_FILLVALUE;
 }