1 /*
2 This file is part of CDO. CDO is a collection of Operators to manipulate and analyse Climate model Data.
3
4 Author: Uwe Schulzweida
5
6 */
7
8 /*
9 This module contains the following operators:
10
11 Condc ifthenc If then constant
12 Condc ifnotthenc If not then constant
13 */
14
15 #include <cdi.h>
16
17 #include "process_int.h"
18 #include "cdo_vlist.h"
19 #include "param_conversion.h"
20
21 static void
operator_IFTHENC(size_t n,double mv,const Varray<double> & vIn,const double cVal,Varray<double> & vOut)22 operator_IFTHENC(size_t n, double mv, const Varray<double> &vIn, const double cVal, Varray<double> &vOut)
23 {
24 if (std::isnan(mv))
25 for (size_t i = 0; i < n; i++) vOut[i] = (!dbl_is_equal(vIn[i], mv) && !dbl_is_equal(vIn[i], 0.0)) ? cVal : mv;
26 else
27 for (size_t i = 0; i < n; i++) vOut[i] = (!is_equal(vIn[i], mv) && !is_equal(vIn[i], 0.0)) ? cVal : mv;
28 }
29
30 static void
operator_IFNOTTHENC(size_t n,double mv,const Varray<double> & vIn,const double cVal,Varray<double> & vOut)31 operator_IFNOTTHENC(size_t n, double mv, const Varray<double> &vIn, const double cVal, Varray<double> &vOut)
32 {
33 if (std::isnan(mv))
34 for (size_t i = 0; i < n; i++) vOut[i] = (!dbl_is_equal(vIn[i], mv) && dbl_is_equal(vIn[i], 0.0)) ? cVal : mv;
35 else
36 for (size_t i = 0; i < n; i++) vOut[i] = (!is_equal(vIn[i], mv) && is_equal(vIn[i], 0.0)) ? cVal : mv;
37 }
38
39 void *
Condc(void * process)40 Condc(void *process)
41 {
42 cdo_initialize(process);
43
44 // clang-format off
45 const auto IFTHENC = cdo_operator_add("ifthenc", 0, 0, nullptr);
46 const auto IFNOTTHENC = cdo_operator_add("ifnotthenc", 0, 0, nullptr);
47 // clang-format on
48
49 const auto operatorID = cdo_operator_id();
50
51 operator_input_arg("constant value");
52 const auto rc = parameter_to_double(cdo_operator_argv(0));
53
54 const auto streamID1 = cdo_open_read(0);
55
56 const auto vlistID1 = cdo_stream_inq_vlist(streamID1);
57 const auto vlistID2 = vlistDuplicate(vlistID1);
58
59 const auto taxisID1 = vlistInqTaxis(vlistID1);
60 const auto taxisID2 = taxisDuplicate(taxisID1);
61 vlistDefTaxis(vlistID2, taxisID2);
62
63 VarList varList;
64 varListInit(varList, vlistID1);
65
66 const auto gridsizemax = vlistGridsizeMax(vlistID1);
67 Varray<double> vaIn(gridsizemax), vaOut(gridsizemax);
68
69 const auto streamID2 = cdo_open_write(1);
70 cdo_def_vlist(streamID2, vlistID2);
71
72 int tsID = 0;
73 while (true)
74 {
75 const auto nrecs = cdo_stream_inq_timestep(streamID1, tsID);
76 if (nrecs == 0) break;
77
78 cdo_taxis_copy_timestep(taxisID2, taxisID1);
79 cdo_def_timestep(streamID2, tsID);
80
81 for (int recID = 0; recID < nrecs; recID++)
82 {
83 int varID, levelID;
84 cdo_inq_record(streamID1, &varID, &levelID);
85 size_t nmiss;
86 cdo_read_record(streamID1, vaIn.data(), &nmiss);
87
88 const auto missval = varList[varID].missval;
89 const auto ngp = varList[varID].gridsize;
90
91 if (nmiss > 0) cdo_check_missval(missval, varList[varID].name);
92
93 // clang-format off
94 if (operatorID == IFTHENC) operator_IFTHENC(ngp, missval, vaIn, rc, vaOut);
95 else if (operatorID == IFNOTTHENC) operator_IFNOTTHENC(ngp, missval, vaIn, rc, vaOut);
96 else cdo_abort("Operator not implemented!");
97 // clang-format on
98
99 nmiss = varray_num_mv(ngp, vaOut, missval);
100 cdo_def_record(streamID2, varID, levelID);
101 cdo_write_record(streamID2, vaOut.data(), nmiss);
102 }
103
104 tsID++;
105 }
106
107 cdo_stream_close(streamID2);
108 cdo_stream_close(streamID1);
109
110 cdo_finish();
111
112 return nullptr;
113 }
114