16 #include "eckit/config/Configuration.h"
18 #include "ioda/ObsDataVector.h"
19 #include "ioda/ObsSpace.h"
20 #include "oops/util/Logger.h"
21 #include "oops/util/missingValues.h"
31 :
FilterBase(obsdb, config, flags, obserr), invars_(config_,
"clw variables") {
47 std::vector<std::vector<bool>> & flagged)
const {
48 oops::Log::trace() <<
"MWCLWCheck postFilter" << std::endl;
50 const oops::Variables observed =
obsdb_.obsvariables();
55 std::vector<float> clw_thresholds =
config_.getFloatVector(
"clw_thresholds");
60 const int clw_option =
config_.getInt(
"clw_option");
64 ASSERT(clw_thresholds.size() == filtervars.
nvars());
66 for (
size_t jv = 0; jv < filtervars.
nvars(); ++jv) {
67 ASSERT(clw_thresholds[jv] !=
missing);
72 ASSERT(clw_option >= 1 && clw_option <=3);
84 std::vector<float> hofx0;
86 std::vector<float> hofx1;
89 auto clw_obs = [&](
size_t jobs) {
return amsua_clw(
90 obs_for_calc[0][jobs], obs_for_calc[1][jobs], sza[0][jobs]); };
92 auto clw_guess = [&](
size_t jobs) {
return amsua_clw(
93 hofx0[jobs], hofx1[jobs], sza[0][jobs]); };
96 for (
size_t jobs = 0; jobs < obs.nlocs(); ++jobs) {
99 clw[0][jobs] = clw_obs(jobs);
102 clw[0][jobs] = clw_guess(jobs);
105 clw_obs_out[0][jobs] = clw_obs(jobs);
106 clw_guess_out[0][jobs] = clw_guess(jobs);
107 if (clw_obs_out[0][jobs] ==
missing || clw_guess_out[0][jobs] ==
missing) {
110 clw[0][jobs] = (clw_obs_out[0][jobs] + clw_guess_out[0][jobs])/2.0;
114 oops::Log::error() <<
"Invalid value for clw_option:" << clw_option << std::endl;
115 ABORT(
"Invalid value for clw_option");
119 for (
size_t jv = 0; jv < filtervars.
nvars(); ++jv) {
122 clw[0][jobs] > clw_thresholds[jv]) flagged[jv][jobs] =
true;
127 clw_obs_out.save(
"Derived");
128 clw_guess_out.save(
"Derived");
134 const float d1 = 0.754;
135 const float d2 = -2.265;
141 tobs1 <= 284.0 && tobs2 <= 284.0 && tobs1 > 0.0 && tobs2 > 0.0) {
142 float cossza = cos(M_PI * sza/180.0);
143 float d0 = 8.240 - (2.622 - 1.846*cossza)*cossza;
144 clw = cossza*(d0 + d1*log(285.0-tobs1)) + d2*log(285.0-tobs2);
145 clw = std::max(0.0f,
clw);
156 os <<
"MWCLWCheck: config = " <<
config_ << std::endl;