Refactor ObsADT::Change georef to include only QCed observations

src/mains/RunCRTM.h

@@ -76,6 +76,11 @@ template <typename OBS> class RunCRTM : public oops::Application {
       ObsDiags_ diag(obsdb[jj], hop.locations(), diagvars);
 
       hop.simulateObs(gval, hofx, ybias, qcflags, bias, diag);
+
+      const double zz = hofx.rms();
+      const double xx = conf[jj].getDouble("rms ref");
+      const double tol = conf[jj].getDouble("tolerance");
+//      BOOST_CHECK_CLOSE(xx, zz, tol);
     }
 
     return 0;

src/mains/RunPreProcess.h

@@ -57,8 +57,10 @@ class RunPreProcess : public oops::Application {
                          util::TimeWindow(params.timeWindow.value()),
                          this->getComm());
 
-    ioda::ObsDataVector<float> obserr_(obsdb, obsdb.obsvariables(), "ObsError");
-    ioda::ObsDataVector<int> qcflags_(obsdb, obsdb.obsvariables());
+    std::shared_ptr<ioda::ObsDataVector<float>> obserr_ =
+      std::make_shared<ioda::ObsDataVector<float>>(obsdb, obsdb.obsvariables(), "ObsError");
+    std::shared_ptr<ioda::ObsDataVector<int>> qcflags_ =
+      std::make_shared<ioda::ObsDataVector<int>>(obsdb, obsdb.obsvariables());
 
     ufo::ObsFilters filters(obsdb,
                             params.filtersParams.toConfiguration(),

src/ufo/CMakeLists.txt

@@ -195,7 +195,7 @@ if(geos-aero_FOUND)
 endif()
 
 if(gsw_FOUND)
-    target_link_libraries(ufo PUBLIC gsw ${gsw_LIBRARIES})
+    target_link_libraries(ufo PUBLIC gsw)
     target_compile_definitions(ufo PUBLIC GSW_FOUND)
 endif()
 

src/ufo/ObsFilter.cc

@@ -21,7 +21,7 @@ namespace ufo {
 
 ObsFilter::ObsFilter(ioda::ObsSpace & os,
                      const ObsFilterParametersBase & parameters,
-                     ioda::ObsDataVector<int> & flags, ioda::ObsDataVector<float> & obserr)
+                     ObsDataPtr_<int> flags, ObsDataPtr_<float> obserr)
   : ofilt_(), filterName_("ObsFilter::"+parameters.filter.value().value())
 {
   oops::Log::trace() << "ObsFilter::ObsFilter starting" << std::endl;

src/ufo/ObsFilter.h

@@ -30,6 +30,8 @@ class ObsDiagnostics;
 /// \brief A filter processing observations.
 class ObsFilter : public util::Printable,
                   private util::ObjectCounter<ObsFilter > {
+  template <typename DATA> using ObsDataPtr_ = std::shared_ptr<ioda::ObsDataVector<DATA> >;
+
  public:
   static const std::string classname() {return "ufo::ObsFilter";}
 
@@ -39,13 +41,13 @@ class ObsFilter : public util::Printable,
   ///   Space containing the observations to process.
   /// \param params
   ///   The filter's configuration parameters.
-  /// \param flags
+  /// \param qcflags
   ///   Quality control flags. They may be modified by the filter.
-  /// \param obserr
+  /// \param obserrors
   ///   Estimates of the standard deviations of observation errors. They may be modified by the
   ///   filter.
   ObsFilter(ioda::ObsSpace &obsspace, const ObsFilterParametersBase &params,
-            ioda::ObsDataVector<int> & flags, ioda::ObsDataVector<float> & obserr);
+            ObsDataPtr_<int> qcflags, ObsDataPtr_<float> obserrors);
   ObsFilter(const ObsFilter &) = delete;
   ObsFilter(ObsFilter &&) = default;
   ObsFilter& operator=(const ObsFilter &) = delete;

src/ufo/ObsFilterBase.cc

@@ -33,7 +33,7 @@ FilterFactory::FilterFactory(const std::string & name) {
 
 std::unique_ptr<ObsFilterBase>
 FilterFactory::create(ioda::ObsSpace & os, const ObsFilterParametersBase & params,
-                      ioda::ObsDataVector<int> & flags, ioda::ObsDataVector<float> & obserr) {
+                      ObsDataPtr_<int> flags, ObsDataPtr_<float> obserr) {
   oops::Log::trace() << "FilterFactory::create starting" << std::endl;
   const std::string &id = params.filter.value().value();
   typename std::map<std::string, FilterFactory*>::iterator

src/ufo/ObsFilterBase.h

@@ -59,15 +59,15 @@ enum class FilterStage {AUTO,
 /// In the former case, they should provide a constructor with the following signature:
 ///
 ///    ObsFilter(ObsSpace_ &, const eckit::Configuration &,
-///              ioda::ObsDataVector<int> &, ioda::ObsDataVector<float> &);
+///              ObsDataPtr_<int>, ObsDataPtr_<float>);
 ///
 /// In the latter case, the implementer should first define a subclass of ObsFilterParametersBase
 /// holding the settings of the filter in question. The implementation of the ObsFilter interface
 /// should then typedef `Parameters_` to the name of that subclass and provide a constructor with
 /// the following signature:
 ///
 ///    ObsFilter(ObsSpace_ &, const Parameters_ &,
-///              ioda::ObsDataVector<int> &, ioda::ObsDataVector<float> &);
+///              ObsDataPtr_<int>, ObsDataPtr_<float>);
 // TODO(someone): check whether we can just keep Parameters here and remove the
 //                Configuration option.
 class ObsFilterBase : public util::Printable,
@@ -136,6 +136,8 @@ class ObsFilterParametersWrapper : public oops::Parameters {
 
 /// ObsFilter factory
 class FilterFactory {
+  template <typename DATA> using ObsDataPtr_ = std::shared_ptr<ioda::ObsDataVector<DATA> >;
+
  public:
   /// \brief Create and return a new observation filter.
   ///
@@ -144,8 +146,8 @@ class FilterFactory {
   /// otherwise an exception will be thrown.
   static std::unique_ptr<ObsFilterBase> create(ioda::ObsSpace &,
                                                const ObsFilterParametersBase & params,
-                                               ioda::ObsDataVector<int> & flags,
-                                               ioda::ObsDataVector<float> & obserr);
+                                               ObsDataPtr_<int> flags = ObsDataPtr_<int>(),
+                                               ObsDataPtr_<float> obserr = ObsDataPtr_<float>());
 
   /// \brief Create and return an instance of the subclass of ObsFilterParametersBase
   /// storing parameters of observation filters of the specified type.
@@ -167,8 +169,7 @@ class FilterFactory {
  private:
   virtual std::unique_ptr<ObsFilterBase> make(ioda::ObsSpace &,
                                               const ObsFilterParametersBase &,
-                                              ioda::ObsDataVector<int> &,
-                                              ioda::ObsDataVector<float> &) = 0;
+                                              ObsDataPtr_<int>, ObsDataPtr_<float>) = 0;
 
   virtual std::unique_ptr<ObsFilterParametersBase> makeParameters() const = 0;
 
@@ -187,10 +188,12 @@ class FilterMaker : public FilterFactory {
   typedef oops::TParameters_IfAvailableElseFallbackType_t<T, GenericObsFilterParameters>
           Parameters_;
 
+  template <typename DATA> using ObsDataPtr_ = std::shared_ptr<ioda::ObsDataVector<DATA> >;
+
   std::unique_ptr<ObsFilterBase> make(ioda::ObsSpace & os,
                                       const ObsFilterParametersBase & params,
-                                      ioda::ObsDataVector<int> & flags,
-                                      ioda::ObsDataVector<float> & obserr) override {
+                                      ObsDataPtr_<int> flags,
+                                      ObsDataPtr_<float> obserr) override {
         const auto &stronglyTypedParams = dynamic_cast<const Parameters_&>(params);
         return std::make_unique<T>(os,
                oops::parametersOrConfiguration<oops::HasParameters_<T>::value>(stronglyTypedParams),

src/ufo/ObsFilters.cc

@@ -22,11 +22,10 @@ namespace ufo {
 
 ObsFilters::ObsFilters(ioda::ObsSpace & os,
                        const eckit::Configuration & config,
-                       ioda::ObsDataVector<int> & flags,
-                       ioda::ObsDataVector<float> & obserr,
+                       ObsDataPtr_<int> qcflags, ObsDataPtr_<float> obserr,
                        const int iteration)
-  : obsspace_(os), qcmanager_(nullptr), autoFilters_(), preFilters_(), priorFilters_(),
-    postFilters_(), geovars_(), diagvars_(), qcflags_(flags), obserr_(obserr),
+  : obsspace_(os), autoFilters_(), preFilters_(), priorFilters_(), postFilters_(),
+    geovars_(), diagvars_(), qcflags_(qcflags), obserrfilter_(obserr),
     iteration_(iteration) {
   oops::Log::trace() << "ObsFilters::ObsFilters starting" << std::endl;
 
@@ -52,21 +51,22 @@ ObsFilters::ObsFilters(ioda::ObsSpace & os,
 
   // If at least one filter has been configured, prepend the QC manager and
   // append the Final Check to the list of filters.
-  atLeastOneFilterConfigured_ =
-      atLeastOneAutoFilterConfigured || atLeastOnePrePriorPostFilterConfigured;
+  const bool atLeastOneFilterConfigured =
+    atLeastOneAutoFilterConfigured || atLeastOnePrePriorPostFilterConfigured;
 
-  // Prepare QC handling and statistics if any filters are present
-  if (atLeastOneFilterConfigured_) {
+// Prepare QC handling and statistics if any filters are present
+  if (atLeastOneFilterConfigured) {
     eckit::LocalConfiguration conf;
     conf.set("filter", "QCmanager");
     ObsFilterParametersWrapper filterParams;
     filterParams.validateAndDeserialize(conf);
-    qcmanager_ = std::make_unique<ObsFilter>(os, filterParams.filterParameters,
-                                             qcflags_, obserr_);
+    if (atLeastOneAutoFilterConfigured)
+      autoFilters_.emplace_back(os, filterParams.filterParameters, qcflags_, obserrfilter_);
+    else
+      preFilters_.emplace_back(os, filterParams.filterParameters, qcflags_, obserrfilter_);
   }
 
-  // Create the filters, only at 0-th iteration, or at iterations specified in
-  // "apply at iterations"
+// Create the filters, only at 0-th iteration, or at iterations specified in "apply at iterations"
   if (atLeastOneAutoFilterConfigured)
     appendToFiltersList(autoFiltersParams, autoFilters_);
   if (atLeastOnePrePriorPostFilterConfigured) {
@@ -75,39 +75,35 @@ ObsFilters::ObsFilters(ioda::ObsSpace & os,
     appendToFiltersList(postFiltersParams, postFilters_);
   }
 
-  // Create the final filter run at the end of the pipeline
-  if (atLeastOneFilterConfigured_) {
+// Create the final filter run at the end of the pipeline
+  if (atLeastOneFilterConfigured) {
     eckit::LocalConfiguration conf;
     conf.set("filter", "Final Check");
     ObsFilterParametersWrapper filterParams;
     filterParams.validateAndDeserialize(conf);
     if (atLeastOneAutoFilterConfigured)
-      autoFilters_.emplace_back(os, filterParams.filterParameters, qcflags_,
-                                obserr_);
+      autoFilters_.emplace_back(os, filterParams.filterParameters, qcflags_, obserrfilter_);
     else
-      postFilters_.emplace_back(os, filterParams.filterParameters, qcflags_,
-                                obserr_);
+      postFilters_.emplace_back(os, filterParams.filterParameters, qcflags_, obserrfilter_);
   }
 
   oops::Log::trace() << "ObsFilters::ObsFilters done" << std::endl;
 }
 
 // -----------------------------------------------------------------------------
 
-void ObsFilters::appendToFiltersList(const FilterParams_& filtersParams,
-                                     std::vector<ObsFilter>& filters) {
-  for (const ObsFilterParametersWrapper& filterParams : filtersParams) {
+void ObsFilters::appendToFiltersList(const FilterParams_ & filtersParams,
+                                     std::vector<ObsFilter> & filters) {
+  for (const ObsFilterParametersWrapper & filterParams : filtersParams) {
     // Only create filters for the 0-th iteration by default
     bool apply = (iteration_ == 0);
     // If "apply at iterations" is set, check if this is the right iteration
     if (filterParams.applyAtIterations.value() != boost::none) {
-      const std::set<int> iters =
-          oops::parseIntSet(*filterParams.applyAtIterations.value());
+      const std::set<int> iters = oops::parseIntSet(*filterParams.applyAtIterations.value());
       apply = oops::contains(iters, iteration_);
     }
     if (apply) {
-      filters.emplace_back(obsspace_, filterParams.filterParameters,
-                           qcflags_, obserr_);
+      filters.emplace_back(obsspace_, filterParams.filterParameters, qcflags_, obserrfilter_);
       geovars_ += filters.back().requiredVars();
       diagvars_ += filters.back().requiredHdiagnostics();
     }
@@ -117,78 +113,71 @@ void ObsFilters::appendToFiltersList(const FilterParams_& filtersParams,
 // -----------------------------------------------------------------------------
 
 void ObsFilters::preProcess() {
-  if (atLeastOneFilterConfigured_) {
-    if (autoFilters_.size() > 0) {
-      qcmanager_->preProcess();  // No operation, placeholder
-      for (ObsFilter& filter : autoFilters_) {
-        filter.checkFilterData(FilterStage::AUTO);
-        filter.preProcess();
-      }
-    } else {
-      qcmanager_->preProcess();  // No operation, placeholder
-      for (ObsFilter& filter : preFilters_) {
-        filter.checkFilterData(FilterStage::PRE);
-        filter.preProcess();
-      }
+  if (autoFilters_.size() > 0) {
+    for (ObsFilter & filter : autoFilters_) {
+      filter.checkFilterData(FilterStage::AUTO);
+      filter.preProcess();
+    }
+  } else {
+    for (ObsFilter & filter : preFilters_) {
+      filter.checkFilterData(FilterStage::PRE);
+      filter.preProcess();
     }
   }
 }
 
 // -----------------------------------------------------------------------------
 
-void ObsFilters::priorFilter(const GeoVaLs& gv) {
-  if (atLeastOneFilterConfigured_) {
-    qcmanager_->priorFilter(gv);  // No operation, placeholder
-    if (autoFilters_.size() > 0) {
-      for (ObsFilter& filter : autoFilters_) {
-        filter.checkFilterData(FilterStage::AUTO);
-        filter.priorFilter(gv);
-      }
-    } else {
-      qcmanager_->priorFilter(gv);  // No operation, placeholder
-      for (ObsFilter& filter : priorFilters_) {
-        filter.checkFilterData(FilterStage::PRIOR);
-        filter.priorFilter(gv);
-      }
+void ObsFilters::priorFilter(const GeoVaLs & gv) {
+  if (autoFilters_.size() > 0) {
+    for (ObsFilter & filter : autoFilters_) {
+      filter.checkFilterData(FilterStage::AUTO);
+      filter.priorFilter(gv);
+    }
+  } else {
+    for (ObsFilter & filter : priorFilters_) {
+      filter.checkFilterData(FilterStage::PRIOR);
+      filter.priorFilter(gv);
     }
   }
 }
 
 // -----------------------------------------------------------------------------
 
-void ObsFilters::postFilter(const GeoVaLs& gv, const ioda::ObsVector& hofx,
-                            const ioda::ObsVector& bias,
-                            const ObsDiagnostics& diags) {
-  if (atLeastOneFilterConfigured_) {
-    if (autoFilters_.size() > 0) {
-      qcmanager_->postFilter(gv, hofx, bias, diags);
-      for (ObsFilter& filter : autoFilters_) {
-        filter.checkFilterData(FilterStage::AUTO);
-        filter.postFilter(gv, hofx, bias, diags);
-      }
-    } else {
-      qcmanager_->postFilter(gv, hofx, bias, diags);
-      for (ObsFilter& filter : postFilters_) {
-        filter.checkFilterData(FilterStage::POST);
-        filter.postFilter(gv, hofx, bias, diags);
-      }
+void ObsFilters::postFilter(const GeoVaLs & gv,
+                            const ioda::ObsVector & hofx,
+                            const ioda::ObsVector & bias,
+                            const ObsDiagnostics & diags) {
+  if (autoFilters_.size() > 0) {
+    for (ObsFilter & filter : autoFilters_) {
+      filter.checkFilterData(FilterStage::AUTO);
+      filter.postFilter(gv, hofx, bias, diags);
+    }
+  } else {
+    for (ObsFilter & filter : postFilters_) {
+      filter.checkFilterData(FilterStage::POST);
+      filter.postFilter(gv, hofx, bias, diags);
     }
   }
 }
 
 // -----------------------------------------------------------------------------
 
-void ObsFilters::print(std::ostream& os) const {
+void ObsFilters::print(std::ostream & os) const {
   if (autoFilters_.size() > 0) {
     os << "ObsFilters: " << autoFilters_.size() << " elements:" << std::endl;
-    for (const ObsFilter& filter : autoFilters_) os << filter << std::endl;
+    for (const ObsFilter & filter : autoFilters_)
+      os << filter << std::endl;
   } else {
     os << "Pre filters: " << preFilters_.size() << " elements:" << std::endl;
-    for (const ObsFilter& preFilter : preFilters_) os << preFilter << std::endl;
+    for (const ObsFilter & preFilter : preFilters_)
+      os << preFilter << std::endl;
     os << "Prior filters: " << priorFilters_.size() << " elements:" << std::endl;
-    for (const ObsFilter& priorFilter : priorFilters_) os << priorFilter << std::endl;
+    for (const ObsFilter & priorFilter : priorFilters_)
+      os << priorFilter << std::endl;
     os << "Post filters: " << postFilters_.size() << " elements:" << std::endl;
-    for (const ObsFilter& postFilter : postFilters_) os << postFilter << std::endl;
+    for (const ObsFilter & postFilter : postFilters_)
+      os << postFilter << std::endl;
   }
 }