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 Pardup pardup Duplicate parameters
12 Pardup parmul Multiply parameters
13 */
14
15 #include <cdi.h>
16
17 #include "cdo_vlist.h"
18 #include "process_int.h"
19 #include "param_conversion.h"
20
21 void *
Pardup(void * process)22 Pardup(void *process)
23 {
24 cdo_initialize(process);
25
26 const auto PARDUP = cdo_operator_add("pardup", 0, 0, nullptr);
27 const auto PARMUL = cdo_operator_add("parmul", 0, 0, nullptr);
28
29 const auto operatorID = cdo_operator_id();
30
31 int nmul = 0;
32 if (operatorID == PARDUP)
33 {
34 nmul = 2;
35 }
36 else if (operatorID == PARMUL)
37 {
38 operator_input_arg("number of multiply");
39 nmul = parameter_to_int(cdo_operator_argv(0));
40 }
41 else
42 cdo_abort("operator not implemented!");
43
44 const auto streamID1 = cdo_open_read(0);
45
46 const auto vlistID1 = cdo_stream_inq_vlist(streamID1);
47 const auto vlistID2 = vlistDuplicate(vlistID1);
48
49 const auto taxisID1 = vlistInqTaxis(vlistID1);
50 const auto taxisID2 = taxisDuplicate(taxisID1);
51 vlistDefTaxis(vlistID2, taxisID2);
52
53 VarList varList1;
54 varListInit(varList1, vlistID1);
55
56 const auto nvars = vlistNvars(vlistID1);
57
58 const auto maxrecs = vlistNrecs(vlistID1);
59 std::vector<RecordInfo> recList(maxrecs);
60
61 const auto gridsizemax = vlistGridsizeMax(vlistID1);
62 Varray<double> array(gridsizemax);
63 Varray2D<double> vardata(nvars);
64 std::vector<std::vector<size_t>> varnmiss(nvars);
65
66 for (int varID = 0; varID < nvars; varID++)
67 {
68 const auto gridsize = varList1[varID].gridsize;
69 const auto nlevels = varList1[varID].nlevels;
70 vardata[varID].resize(gridsize * nlevels);
71 varnmiss[varID].resize(nlevels);
72 }
73
74 for (int i = 1; i < nmul; i++)
75 {
76 vlistCat(vlistID2, vlistID1);
77 for (int varID = 0; varID < nvars; varID++)
78 vlistDefVarParam(vlistID2, varID + nvars * i, cdiEncodeParam(-(varID + nvars * i + 1), 255, 255));
79 }
80
81 const auto streamID2 = cdo_open_write(1);
82 cdo_def_vlist(streamID2, vlistID2);
83
84 int tsID = 0;
85 while (true)
86 {
87 const auto nrecs = cdo_stream_inq_timestep(streamID1, tsID);
88 if (nrecs == 0) break;
89
90 cdo_taxis_copy_timestep(taxisID2, taxisID1);
91 cdo_def_timestep(streamID2, tsID);
92
93 for (int recID = 0; recID < nrecs; recID++)
94 {
95 int varID, levelID;
96 cdo_inq_record(streamID1, &varID, &levelID);
97
98 recList[recID].varID = varID;
99 recList[recID].levelID = levelID;
100 recList[recID].lconst = (varList1[varID].timetype == TIME_CONSTANT);
101
102 const auto gridsize = varList1[varID].gridsize;
103 const auto offset = gridsize * levelID;
104 auto single = &vardata[varID][offset];
105
106 size_t nmiss;
107 cdo_read_record(streamID1, single, &nmiss);
108 varnmiss[varID][levelID] = nmiss;
109 }
110
111 for (int i = 0; i < nmul; i++)
112 for (int recID = 0; recID < nrecs; recID++)
113 {
114 const auto varID = recList[recID].varID;
115 const auto levelID = recList[recID].levelID;
116 const auto varID2 = varID + i * nvars;
117
118 const auto gridsize = varList1[varID].gridsize;
119 const auto offset = gridsize * levelID;
120 auto single = &vardata[varID][offset];
121 const auto nmiss = varnmiss[varID][levelID];
122
123 array_copy(gridsize, single, array.data());
124 cdo_def_record(streamID2, varID2, levelID);
125 cdo_write_record(streamID2, array.data(), nmiss);
126 }
127
128 tsID++;
129 }
130
131 cdo_stream_close(streamID2);
132 cdo_stream_close(streamID1);
133
134 cdo_finish();
135
136 return nullptr;
137 }
138