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 #include <cdi.h>
9 
10 #include "process_int.h"
11 #include "cdo_default_values.h"  // Namespace CdoDefault
12 
13 void *
Tocomplex(void * process)14 Tocomplex(void *process)
15 {
16   cdo_initialize(process);
17 
18   const auto RETOCOMPLEX = cdo_operator_add("retocomplex", 0, 0, nullptr);
19   const auto IMTOCOMPLEX = cdo_operator_add("imtocomplex", 0, 0, nullptr);
20 
21   const auto operatorID = cdo_operator_id();
22 
23   operator_check_argc(0);
24 
25   const auto streamID1 = cdo_open_read(0);
26 
27   const auto vlistID1 = cdo_stream_inq_vlist(streamID1);
28   const auto vlistID2 = vlistDuplicate(vlistID1);
29 
30   const auto nvars = vlistNvars(vlistID2);
31   for (int varID = 0; varID < nvars; ++varID)
32     {
33       auto datatype = vlistInqVarDatatype(vlistID2, varID);
34       datatype = (datatype == CDI_DATATYPE_FLT64) ? CDI_DATATYPE_CPX64 : CDI_DATATYPE_CPX32;
35       vlistDefVarDatatype(vlistID2, varID, datatype);
36     }
37 
38   const auto taxisID1 = vlistInqTaxis(vlistID1);
39   const auto taxisID2 = taxisDuplicate(taxisID1);
40   vlistDefTaxis(vlistID2, taxisID2);
41 
42   // if (CdoDefault::FileType != CDI_FILETYPE_EXT) cdo_abort("Complex numbers need EXTRA format; used CDO option -f ext!");
43   const auto streamID2 = cdo_open_write(1);
44   cdo_def_vlist(streamID2, vlistID2);
45 
46   const auto gridsizemax = vlistGridsizeMax(vlistID1);
47   Varray<double> array1(gridsizemax), array2(2 * gridsizemax);
48 
49   VarList varList1;
50   varListInit(varList1, vlistID1);
51 
52   int tsID = 0;
53   int tsID2 = 0;
54   while (true)
55     {
56       const auto nrecs = cdo_stream_inq_timestep(streamID1, tsID);
57       if (nrecs == 0) break;
58 
59       cdo_taxis_copy_timestep(taxisID2, taxisID1);
60       cdo_def_timestep(streamID2, tsID2++);
61 
62       for (int recID = 0; recID < nrecs; recID++)
63         {
64           size_t nmiss;
65           int varID, levelID;
66           cdo_inq_record(streamID1, &varID, &levelID);
67           cdo_def_record(streamID2, varID, levelID);
68 
69           cdo_read_record(streamID1, array1.data(), &nmiss);
70 
71           const auto gridsize = varList1[varID].gridsize;
72           if (operatorID == RETOCOMPLEX)
73             {
74               for (size_t i = 0; i < gridsize; ++i)
75                 {
76                   array2[2 * i] = array1[i];
77                   array2[2 * i + 1] = 0;
78                 }
79             }
80           else if (operatorID == IMTOCOMPLEX)
81             {
82               for (size_t i = 0; i < gridsize; ++i)
83                 {
84                   array2[2 * i] = 0;
85                   array2[2 * i + 1] = array1[i];
86                 }
87             }
88 
89           cdo_write_record(streamID2, array2.data(), nmiss);
90         }
91 
92       tsID++;
93     }
94 
95   cdo_stream_close(streamID2);
96   cdo_stream_close(streamID1);
97 
98   vlistDestroy(vlistID2);
99 
100   cdo_finish();
101 
102   return nullptr;
103 }
104