31
31
32
32
#include " orca-jedi/interpolator/Interpolator.h"
33
33
34
- #include < boost/uuid/uuid.hpp> // uuid class DJL
35
- #include < boost/uuid/uuid_generators.hpp> // generators DJL
36
- #include < boost/uuid/uuid_io.hpp> // streaming operators etc. DJL
37
-
38
- #include " orca-jedi/utilities/IOUtils.h" // DJL
39
-
40
-
41
34
namespace eckit {
42
35
class Configuration ;
43
36
}
@@ -193,7 +186,6 @@ template void Interpolator::executeInterpolation<float>(
193
186
void Interpolator::apply (const oops::Variables& vars, const Increment& inc,
194
187
const std::vector<bool > & mask,
195
188
std::vector<double >& result) const {
196
- // DJL question can the atlas templates help?
197
189
// input is inc output is result
198
190
const size_t nvars = vars.size ();
199
191
@@ -226,7 +218,6 @@ void Interpolator::apply(const oops::Variables& vars, const Increment& inc,
226
218
auto field_view = atlas::array::make_view<double , 2 >(tgt_field);
227
219
atlas::field::MissingValue mv (inc.incrementFields ()[gv_varname]);
228
220
bool has_mv = static_cast <bool >(mv);
229
- oops::Log::debug () << " DJL Interpolator::apply mv " << mv << " has_mv " << has_mv << std::endl;
230
221
for (std::size_t klev=0 ; klev < varSizes[jvar]; ++klev) {
231
222
for (std::size_t iloc=0 ; iloc < nlocs_; iloc++) {
232
223
if (has_mv && mv (field_view (iloc, klev))) {
@@ -254,10 +245,6 @@ void Interpolator::applyAD(const oops::Variables& vars, Increment& inc,
254
245
oops::Log::trace () << " orcamodel::Interpolator::applyAD start "
255
246
<< std::endl;
256
247
257
- oops::Log::debug () << " DJL ** Interpolator::applyAD this needs checking **" << std::endl;
258
-
259
- // ** Trying to do the opposite of apply DJL
260
-
261
248
const size_t nvars = vars.size ();
262
249
263
250
for (size_t j=0 ; j < nvars; ++j) {
@@ -277,65 +264,38 @@ void Interpolator::applyAD(const oops::Variables& vars, Increment& inc,
277
264
inc.geometry ()->variableSizes (vars);
278
265
size_t nvals = 0 ;
279
266
280
- // DJL write debug fields to file
281
- boost::uuids::uuid uuid = boost::uuids::random_generator ()();
282
- std::shared_ptr<const Geometry> geom = inc.geometry ();
283
- writeFieldsToFile (" applyADpre" + boost::uuids::to_string (uuid) +" .nc" , *geom,
284
- inc.validTime (), inc.incrementFields ());
285
-
286
267
for (size_t jvar=0 ; jvar < nvars; ++jvar) nvals += nlocs_ * varSizes[jvar];
287
- // result.resize(nvals);
288
-
289
268
290
269
std::size_t out_idx = 0 ;
291
270
for (size_t jvar=0 ; jvar < nvars; ++jvar) {
292
- oops::Log::debug () << " DJL ** jvar " << jvar << " " << nvars
293
- << " varSizes " << varSizes[jvar]
294
- << std::endl;
295
271
auto gv_varname = vars[jvar].name ();
296
- // atlas::Field tgt_field = atlasObsFuncSpace_.createField<double>(
297
- // atlas::option::name(gv_varname) |
298
- // atlas::option::levels(varSizes[jvar]));
299
- // atlas::Field incField = inc.incrementFields()[jvar];
300
-
301
272
auto tgt_field = atlasObsFuncSpace_.createField <double >(
302
273
atlas::option::name (gv_varname) |
303
274
atlas::option::levels (varSizes[jvar]));
304
275
305
276
// Copying observation array vector to an atlas observation field (tgt_field)
306
- // DJL not sure if the missing value aspect does anything here
307
277
auto field_view = atlas::array::make_view<double , 2 >(tgt_field);
308
- // field_view.assign(0.0);
309
278
atlas::field::MissingValue mv (inc.incrementFields ()[gv_varname]);
310
279
bool has_mv = static_cast <bool >(mv);
311
- oops::Log::debug () << " DJL Interpolator::applyAD mv "
312
- << mv << " has_mv " << has_mv << std::endl;
313
280
314
281
for (std::size_t klev=0 ; klev < varSizes[jvar]; ++klev) {
315
282
for (std::size_t iloc=0 ; iloc < nlocs_; iloc++) {
316
283
if (has_mv && mv (field_view (iloc, klev))) {
317
284
field_view (iloc, klev) = util::missingValue<double >();
318
285
} else {
319
- oops::Log::debug () << " DJL iloc " << iloc << " klev " << klev << " out_idx "
320
- << out_idx << " resultin[out_idx] " << resultin[out_idx] << std::endl;
321
286
field_view (iloc, klev) = resultin[out_idx];
322
287
}
323
288
++out_idx;
324
289
}
325
290
}
326
291
327
- // halo exchange update ghost points DJL
292
+ // halo exchange update ghost points
293
+ std::shared_ptr<const Geometry> geom = inc.geometry ();
328
294
geom->functionSpace ().haloExchange (inc.incrementFields ()[gv_varname]);
329
295
330
296
interpolator_.execute_adjoint (inc.incrementFields ()[gv_varname], tgt_field);
331
297
} // jvar
332
298
333
- // DJL write debug fields to file
334
- // boost::uuids::uuid uuid = boost::uuids::random_generator()();
335
- // std::shared_ptr<const Geometry> geom = inc.geometry();
336
- writeFieldsToFile (" applyAD" + boost::uuids::to_string (uuid) +" .nc" , *geom, inc.validTime (),
337
- inc.incrementFields ());
338
-
339
299
oops::Log::trace () << " orcamodel::Interpolator::applyAD done "
340
300
<< std::endl;
341
301
}
0 commit comments