Selmulti.cc 49.6 KB
Newer Older
1 2 3 4
/*
  This file is part of CDO. CDO is a collection of Operators to
  manipulate and analyse Climate model Data.

Uwe Schulzweida's avatar
Uwe Schulzweida committed
5
  Copyright (C) 2003-2012 Uwe Schulzweida, <uwe.schulzweida AT mpimet.mpg.de>
6 7 8 9 10 11 12 13 14 15 16 17 18
  See COPYING file for copying and redistribution conditions.

  This program is free software; you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation; version 2 of the License.

  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.
*/

#include <cdi.h>
Oliver Heidmann's avatar
Oliver Heidmann committed
19

20
#include "process_int.h"
21
#include "cdo_zaxis.h"
22
#include "readline.h"
23
#include "cdi_lockedIO.h"
24

25 26
// NOTE: All operators in this module works only on GRIB edition 1 files!

Uwe Schulzweida's avatar
Uwe Schulzweida committed
27
extern "C" void streamGrbChangeParameterIdentification(int code, int ltype, int lev);
28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54

/*
Supported notations:
======================
Selection provided on commandline:
---------------------------------
cdo selmulti,'(33/34;105;10)'
 - or -
cdo selmulti,'(33/34;105;10);(11/6;109;55)'
cdo selmulti,'(33/34;105;10);(11/6;109;40/55)'
cdo selmulti,'(*;105;10);(11/6;109;40/55);(*;105;2)'
cdo selmulti,'{(33/34;105;10);(11/32,8;109;51/52/53/54/55)}'

NOTE: ' .. ' are mandatory !

Selection provided from a text file:
---------------------------------

cdo selmulti,selection_10m_wind.txt

(*A*) Compact general notation, selection file content:

(1; 103; 0)
(33,34; 105; 10)
(11,17; 105; 2)
(71,73,74,75,61,62,65,117,67,122,121,11,131,66,84,111,112; 105; 0)
# If nothing <'sel(' or 'del('>  is specified then
55 56
# the operator -selmulti or -delmulti decides if it will be selection of
extraction or delete # Explicite select or delete is also possible:
57 58 59 60 61 62
#(11; 109; *)
#(*; 105; *)
#del(*; 109; *)
#sel(*; 105; *)
#sel(*; 100; *)

63 64
# BUT simple array arithmetics should be also possible ("*" ~= mulc;  "+' ~=
addc) sel(33,34;105,1000,3000):math(*2;)        # not implemented yet
65 66 67 68 69 70 71
sel(11;105,500,1500,3000):math(+273.15;)  # not implemented yet

(*B*) HIP.X notation (KNMI specific), selection file content:

SELECT, PARAMETER=1, LEVTYPE=103, LEVEL=0
SELECT, PARAMETER=33/34, LEVTYPE=105, LEVEL=10
SELECT, PARAMETER=11/17, LEVTYPE=105, LEVEL=2
72 73 74
SELECT, PARAMETER=71/73/74/75/61/62/65/117/67/122/121/11/131/66/84/111/112,
LEVTYPE=105, LEVEL=0 # Explicite delete is also possible: #DELETE,
PARAMETER=128, LEVTYPE=109, LEVEL=*
75

76 77
# BUT simple array arithmetics should be also possible (SCALE ~= mulc;  OFFSET
~= addc)
78

79 80 81 82
# The following will convert Pressure from Pa into HPa; Temp from Kelvin to
Celsius: SELECT, PARAMETER=1, LEVTYPE= 103, LEVEL=0, SCALE=0.01 SELECT,
PARAMETER=11, LEVTYPE=105, LEVEL=2, OFFSET=273.15 SELECT, PARAMETER=33/34,
LEVTYPE=105, LEVEL=10
83

84 85 86
If SCALE and/or OFFSET are defined, then the data values are scaled as
SCALE*(VALUE-OFFSET). The default value for SCALE is 1.0; the default for OFFSET
is 0.0.
87 88 89

**** changemulti ***********

90 91 92
cdo
changemulti,'{(134;1;*|1;105;*);{(6;1;*|6;105;*)};{(246;*;*|76;*;*)};{(247;*;*|58;*;*)};{(248;*;*|71;*;*)}'
fileIN fileOUT
93 94 95 96


cdo changemulti,'{(134;1;*|1;105;*)}' fileIN fileOUT
# surface pressure has ECMWF code; change it into Hirlam notation ..
97 98 99
grib_set -w indicatorOfParameter=134,indicatorOfTypeOfLevel=1 -s
indicatorOfParameter=1,indicatorOfTypeOfLevel=105 ECMWF_H11_test0.grb
ECMWF_H11_test1.grb
100 101 102

cdo changemulti,'{(6;1;*|6;105;*)}' fileIN fileOUT
# orography has wrong level-type, should be 105
103 104 105
grib_set -w indicatorOfParameter=6,indicatorOfTypeOfLevel=1 -s
indicatorOfParameter=6,indicatorOfTypeOfLevel=105 ECMWF_H11_test1.grb
ECMWF_H11_test2.grb
106 107 108

cdo changemulti,'{(246;*;*|76;*;*)}' fileIN fileOUT
# change code for cloud_water
109 110
grib_set -w indicatorOfParameter=246 -s indicatorOfParameter=76
ECMWF_H11_test2.grb ECMWF_H11_test3.grb
111 112 113

cdo changemulti,'{(247;*;*|58;*;*)}' fileIN fileOUT
# change code for cloud_ice
114 115
grib_set -w indicatorOfParameter=247 -s indicatorOfParameter=58
ECMWF_H11_test3.grb ECMWF_H11_test4.grb
116 117 118

cdo changemulti,'{(248;*;*|71;*;*)}' fileIN fileOUT
# change code for total_cloud_cover
119 120
grib_set -w indicatorOfParameter=248 -s indicatorOfParameter=71
ECMWF_H11_test4.grb ECMWF_H11_test.grb
121 122 123 124


*/

Uwe Schulzweida's avatar
Uwe Schulzweida committed
125
struct TUPLEREC
126
{
127
  std::vector<int> codeLST;
128
  int ncodes;
129

130
  std::vector<int> levelTypeLST;
131
  int nlevelTypes;
132

133
  std::vector<int> levelLST;
134 135 136
  int nlevels;
  int sel_or_del_or_change;  // sel_or_del_or_change:  0:  operator decides,
                             // 1:select , 2:delete, 3:change
Uwe Schulzweida's avatar
Uwe Schulzweida committed
137
  int simpleMath;            // 1:  simple array arithmetics ( *,+), 0: do nothing
138 139 140
  float scale;
  float offset;

141
  int changedCode;  // used only changemulti mode
142 143
  int changedLevelType;
  int changedLevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
144
};
145

146
int checkListContainsInt(int value, std::vector<int> &list, int arraylen);
147 148 149 150 151 152 153

#define MAX_TUPLES 1000
static int NUMTUPLES = 0;
static TUPLEREC *SelTUPLEREC[MAX_TUPLES];

TUPLEREC *TUPLERECNew();
void push_backSelTuple(TUPLEREC *tp);
154
TUPLEREC *getSelTuple(int index);
155 156 157 158 159 160 161

void printSelectionTuples();
int getNumberOfSelectionTuples();
int getNumberOfDeleteSelectionTuples();

int multiSelectionParser(const char *filenameOrString);

162 163
void *
Selmulti(void *process)
164
{
165 166 167 168 169 170
  int varID, levelID;
  int nlevs, code, zaxisID;
  int ltype = 0;
  int varID2, levelID2;
  int sellevel, selcode, selltype;
  bool lcopy = false;
171
  size_t nmiss;
172
  int simpleMath = 0;  // 1:  simple array arithmetics ( *,+), 0: do nothing
173
  float scale = 1.0;
174 175 176 177
  float offset = 0.0;  // If SCALE and/or OFFSET are defined, then the data
                       // values are scaled as SCALE*(VALUE-OFFSET). The default
                       // value for SCALE is 1.0; the default for OFFSET is 0.0.
  double missval;
178

179
  cdoInitialize(process);
180

Uwe Schulzweida's avatar
Uwe Schulzweida committed
181
  // clang-format off
182 183 184
  const auto SELMULTI    = cdoOperatorAdd("selmulti",    0, 0, "filename/string with selection specification ");
  const auto DELMULTI    = cdoOperatorAdd("delmulti",    0, 0, "filename/string with selection specification ");
  const auto CHANGEMULTI = cdoOperatorAdd("changemulti", 0, 0, "filename/string with selection specification ");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
185
  // clang-format on
186

187
  const auto operatorID = cdoOperatorID();
188 189 190

  operatorInputArg(cdoOperatorEnter(operatorID));

191
  const char *filenameOrString;
192

193
  // operatorCheckArgc(1);
194
  filenameOrString = cdoOperatorArgv(0).c_str();
Oliver Heidmann's avatar
Oliver Heidmann committed
195
  if (cdoDebugExt)
196
    {
197
      printf("Given operator arguments (nr=%d): \n", operatorArgc());
198
      for (int i = 0; i < operatorArgc(); i++) printf("%s", cdoOperatorArgv(i).c_str());
199
      printf("\n");
200
    }
201
  if (!multiSelectionParser(filenameOrString))
202
    cdoWarning("Error processing file with selection description!\n%s", filenameOrString);
203

204
  if (operatorID == SELMULTI)
205
    if (getNumberOfSelectionTuples() == 0)
206 207
      cdoAbort("Error! You must provide at lease ONE selection tuple!\n"
               "Notations: 'SELECT,  .. or sel(/;;) or (/;;)'\nCheck the file: %s",
208
               filenameOrString);
209

210
  if (operatorID == DELMULTI)
211
    if (getNumberOfDeleteSelectionTuples() == 0)
212 213
      cdoAbort("Error! You must provide at lease ONE selection tuple!\n"
               "Notations: 'DELETE,  .. or del(/;;) or (/;;)'\nCheck the file: %s",
214
               filenameOrString);
215

216
  if (operatorID == CHANGEMULTI)
217
    if (getNumberOfSelectionTuples() == 0)
218 219
      cdoAbort("Error! You must provide at lease ONE selection tuple!\n"
               "Notations: 'CHANGE,  .. or (/;;|;;;)'\nCheck the file: %s",
220
               filenameOrString);
221

222
  const auto streamID1 = cdoOpenRead(0);
223

224
  const auto vlistID1 = cdoStreamInqVlist(streamID1);
225

226
  vlistClearFlag(vlistID1);
227
  auto nvars = vlistNvars(vlistID1);
228

Uwe Schulzweida's avatar
Uwe Schulzweida committed
229
  Debug(cdoDebugExt, " Total number of variables: %d", nvars);
230

231
  for (varID = 0; varID < nvars; varID++)
232
    {
233
      code = vlistInqVarCode(vlistID1, varID);
234 235
      zaxisID = vlistInqVarZaxis(vlistID1, varID);
      ltype = zaxis2ltype(zaxisID);
236
      nlevs = zaxisInqSize(zaxisID);
237

238
      for (levelID = 0; levelID < nlevs; levelID++)
239
        {
240
          const double level = zaxisInqLevel(zaxisID, levelID);
241

242
          if (operatorID == DELMULTI)
243
            vlistDefFlag(vlistID1, varID, levelID, true);  // set initially, override bellow if in selection
244 245
          if (operatorID == CHANGEMULTI)
            {
246
              vlistDefFlag(vlistID1, varID, levelID, true);  // change operation copies all fields
247 248
              continue;
            }
249

250
          for (int ii = 0; ii < NUMTUPLES; ii++)
251
            {
252
              TUPLEREC *tuplerec = getSelTuple(ii);
Oliver Heidmann's avatar
Oliver Heidmann committed
253
              // if ( cdoDebugExt ) cdoPrint(" Processing: (code %d,
254 255 256 257
              // ltype %d, level %d);  nvars=%d, varID=%d", code, ltype,
              // (int)level, nvars, varID);
              // Note: When the list is Empty then function
              // checkListContainsInt() also returns true !
Uwe Schulzweida's avatar
Uwe Schulzweida committed
258 259 260
              selcode = checkListContainsInt(code, tuplerec->codeLST, tuplerec->ncodes);
              selltype = checkListContainsInt(ltype, tuplerec->levelTypeLST, tuplerec->nlevelTypes);
              sellevel = checkListContainsInt((int) level, tuplerec->levelLST, tuplerec->nlevels);
261
              if (selcode && selltype && sellevel)
262
                {
263
                  if (operatorID == SELMULTI)
264
                    {
265 266
                      switch (tuplerec->sel_or_del_or_change)
                        {
267
                        case 0:  // operator decides ...
268
                          vlistDefFlag(vlistID1, varID, levelID, true);
Oliver Heidmann's avatar
Oliver Heidmann committed
269
                          if (cdoDebugExt)
270 271
                            {
                              if (!tuplerec->simpleMath)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
272 273
                                cdoPrint(" Selecting : (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", code, ltype,
                                         (int) (level), varID, levelID);
274
                              else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
275
                                cdoPrint(" Selecting : (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]; SCALE=%f; "
Uwe Schulzweida's avatar
Uwe Schulzweida committed
276 277
                                         "OFFSET=%f",
                                         code, ltype, (int) (level), varID, levelID, tuplerec->scale, tuplerec->offset);
278 279
                            }
                          break;
280
                        case 1:
281
                          vlistDefFlag(vlistID1, varID, levelID, true);
Oliver Heidmann's avatar
Oliver Heidmann committed
282
                          if (cdoDebugExt)
283 284
                            {
                              if (!tuplerec->simpleMath)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
285 286
                                cdoPrint(" Selecting : (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", code, ltype,
                                         (int) (level), varID, levelID);
287
                              else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
288
                                cdoPrint(" Selecting : (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]; SCALE=%f; "
Uwe Schulzweida's avatar
Uwe Schulzweida committed
289 290
                                         "OFFSET=%f",
                                         code, ltype, (int) (level), varID, levelID, tuplerec->scale, tuplerec->offset);
291 292
                            }
                          break;
293
                        case 2:
294
                          vlistDefFlag(vlistID1, varID, levelID, false);
Oliver Heidmann's avatar
Oliver Heidmann committed
295
                          if (cdoDebugExt)
296 297
                            {
                              if (!tuplerec->simpleMath)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
298 299
                                cdoPrint(" Selecting for removal: (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", code,
                                         ltype, (int) (level), varID, levelID);
300
                              else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
301 302
                                cdoPrint(" Selecting for removal: (code %3i, ltype %3i, level %3i)   "
                                         "[varID(%d),levelID(%d)]; SCALE=%f; OFFSET=%f",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
303
                                         code, ltype, (int) (level), varID, levelID, tuplerec->scale, tuplerec->offset);
304 305 306
                            }
                          break;
                        }
307
                    }
308
                  else if (operatorID == DELMULTI)
309
                    {
310 311
                      switch (tuplerec->sel_or_del_or_change)
                        {
312
                        case 0:  // operator decides ...
313
                          vlistDefFlag(vlistID1, varID, levelID, false);
Oliver Heidmann's avatar
Oliver Heidmann committed
314
                          if (cdoDebugExt)
315 316
                            {
                              if (!tuplerec->simpleMath)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
317 318
                                cdoPrint(" Selecting for removal: (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", code,
                                         ltype, (int) (level), varID, levelID);
319
                              else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
320 321
                                cdoPrint(" Selecting for removal: (code %3i, ltype %3i, level %3i)   "
                                         "[varID(%d),levelID(%d)]; SCALE=%f; OFFSET=%f",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
322
                                         code, ltype, (int) (level), varID, levelID, tuplerec->scale, tuplerec->offset);
323 324
                            }
                          break;
325
                        case 1:
326
                          vlistDefFlag(vlistID1, varID, levelID, true);
Oliver Heidmann's avatar
Oliver Heidmann committed
327
                          if (cdoDebugExt)
328 329
                            {
                              if (!tuplerec->simpleMath)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
330 331
                                cdoPrint(" Selecting : (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", code, ltype,
                                         (int) (level), varID, levelID);
332
                              else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
333
                                cdoPrint(" Selecting : (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]; SCALE=%f; "
Uwe Schulzweida's avatar
Uwe Schulzweida committed
334 335
                                         "OFFSET=%f",
                                         code, ltype, (int) (level), varID, levelID, tuplerec->scale, tuplerec->offset);
336 337
                            }
                          break;
338
                        case 2:
339
                          vlistDefFlag(vlistID1, varID, levelID, false);
Oliver Heidmann's avatar
Oliver Heidmann committed
340
                          if (cdoDebugExt)
341 342
                            {
                              if (!tuplerec->simpleMath)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
343 344
                                cdoPrint(" Selecting for removal: (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", code,
                                         ltype, (int) (level), varID, levelID);
345
                              else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
346 347
                                cdoPrint(" Selecting for removal: (code %3i, ltype %3i, level %3i)   "
                                         "[varID(%d),levelID(%d)]; SCALE=%f; OFFSET=%f",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
348
                                         code, ltype, (int) (level), varID, levelID, tuplerec->scale, tuplerec->offset);
349 350 351
                            }
                          break;
                        }
352
                    }
353
                  break;
354
                }
355 356 357
            }  // end for ( .. NUMTUPLES
        }      // end for ( levelID
    }          // end for ( varID
358

Uwe Schulzweida's avatar
Uwe Schulzweida committed
359
  Debug(cdoDebugExt, " Writing the selected fields ...");
360

361
  auto vlistID2 = vlistCreate();
362
  cdoVlistCopyFlag(vlistID2, vlistID1);
363

364
  nvars = vlistNvars(vlistID2);
365 366 367
  for (varID = 0; varID < nvars; ++varID)
    if (vlistInqVarTimetype(vlistID2, varID) != TIME_CONSTANT) break;
  if (varID == nvars) vlistDefNtsteps(vlistID2, 0);
368

369 370
  const auto taxisID1 = vlistInqTaxis(vlistID1);
  const auto taxisID2 = taxisDuplicate(taxisID1);
371
  vlistDefTaxis(vlistID2, taxisID2);
372

373
  const auto streamID2 = cdoOpenWrite(1);
374
  cdoDefVlist(streamID2, vlistID2);
375

376
  auto gridsize = vlistGridsizeMax(vlistID1);
377
  if (vlistNumber(vlistID1) != CDI_REAL) gridsize *= 2;
378
  Varray<double> array(gridsize);
379

380
  int tsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
381
  while (true)
382
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
383 384 385
      const auto nrecs = cdoStreamInqTimestep(streamID1, tsID);
      if (nrecs == 0) break;

386
      taxisCopyTimestep(taxisID2, taxisID1);
387
      cdoDefTimestep(streamID2, tsID);
388

389
      for (int recID = 0; recID < nrecs; recID++)
390
        {
391
          cdoInqRecord(streamID1, &varID, &levelID);
392
          missval = vlistInqVarMissval(vlistID1, varID);
393

394
          if (vlistInqFlag(vlistID1, varID, levelID) == true)
395
            {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
396
              simpleMath = 0;  // 1:  simple array arithmetics ( *,+), 0: do nothing
397 398
              scale = 1.0;
              offset = 0.0;
399
              code = vlistInqVarCode(vlistID1, varID);
400
              zaxisID = vlistInqVarZaxis(vlistID1, varID);
401
              const auto level = zaxisInqLevel(zaxisID, levelID);
402
              ltype = zaxis2ltype(zaxisID);
403
              for (int ii = 0; ii < NUMTUPLES; ii++)
404
                {
405
                  TUPLEREC *tuplerec = getSelTuple(ii);
406 407
                  // Note: When the list is Empty then function
                  // checkListContainsInt() also returns true !
Uwe Schulzweida's avatar
Uwe Schulzweida committed
408 409 410
                  selcode = checkListContainsInt(code, tuplerec->codeLST, tuplerec->ncodes);
                  selltype = checkListContainsInt(ltype, tuplerec->levelTypeLST, tuplerec->nlevelTypes);
                  sellevel = checkListContainsInt((int) level, tuplerec->levelLST, tuplerec->nlevels);
411
                  lcopy = true;
412
                  if (selcode && selltype && sellevel)
413
                    {
414
                      if (operatorID == CHANGEMULTI)
415
                        {
Oliver Heidmann's avatar
Oliver Heidmann committed
416
                          if (cdoDebugExt)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
417
                            cdoPrint(" Processing: (code %d, ltype %d, level %d);  nvars=%d, varID=%d => (selcode %d, selltype "
Uwe Schulzweida's avatar
Uwe Schulzweida committed
418
                                     "%d, sellevel %d) => change (%d,%d,%d)",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
419 420
                                     code, ltype, (int) level, nvars, varID, selcode, selltype, sellevel, tuplerec->changedCode,
                                     tuplerec->changedLevelType, tuplerec->changedLevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
421
                          if ((tuplerec->changedCode == -1) && (tuplerec->changedLevelType == -1) && (tuplerec->changedLevel == -1))
422 423
                            cdoPrint(" WARNING: Cannot CHANGE identification!");
                          else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
424 425
                            streamGrbChangeParameterIdentification(tuplerec->changedCode, tuplerec->changedLevelType,
                                                                   tuplerec->changedLevel);
426 427 428 429 430 431
                          // Calling PROXY function
                          // streamGrbChangeParameterIdentification() which
                          // results in later calling func.
                          // gribapiChangeParameterIdentification(); see
                          // stream_gribapi.c The change happens during the last
                          // step of writing a grib-record into a file.
432
                        }
433
                      else
434
                        {
Oliver Heidmann's avatar
Oliver Heidmann committed
435
                          if (cdoDebugExt)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
436
                            cdoPrint(" Processing: (code %d, ltype %d, level %d);  nvars=%d, varID=%d => (selcode %d, "
437
                                     "selltype %d, sellevel %d)",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
438
                                     code, ltype, (int) level, nvars, varID, selcode, selltype, sellevel);
439
                        }
440
                      simpleMath = tuplerec->simpleMath;  // 1:  simple array arithmetics ( *,+),
441
                                                          // 0: do nothing
442
                      if (simpleMath)
443
                        {
444 445 446
                          scale = tuplerec->scale;
                          offset = tuplerec->offset;
                          lcopy = false;
447
                        }
448
                      break;  // get out of this for loop
449
                    }
450
                }  // end of for (ii=0; ii<NUMTUPLES ..
451

452
              varID2 = vlistFindVar(vlistID2, varID);
453
              levelID2 = vlistFindLevel(vlistID2, varID, levelID);
454

455
              // tijdelijk PATCH M.K.
456 457 458
              if ((varID2 == -1) || (levelID2 == -1))
                {
                  cdoPrint(" Warning: Missing varID or levelID with (code %3i, "
Uwe Schulzweida's avatar
Uwe Schulzweida committed
459
                           "ltype %3i, level %3i)   [varID(%d),levelID(%d)] .. #2[varID(%d),levelID(%d)]",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
460
                           code, ltype, (int) (level), varID, levelID, varID2, levelID2);
461 462
                  continue;
                }
463
              cdoDefRecord(streamID2, varID2, levelID2);
464

465
              if (lcopy)
466
                {
Oliver Heidmann's avatar
Oliver Heidmann committed
467
                  if (cdoDebugExt)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
468 469
                    cdoPrint(" Copying record [%4d] with (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", recID, code,
                             ltype, (int) (level), varID, levelID);
470
                  cdoCopyRecord(streamID2, streamID1);
471
                }
472
              else
473
                {
474
                  cdoReadRecord(streamID1, array.data(), &nmiss);
475

476
                  if (!simpleMath)
477
                    {
Oliver Heidmann's avatar
Oliver Heidmann committed
478
                      if (cdoDebugExt)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
479 480
                        cdoPrint(" Writing record [%4d] with (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]", recID,
                                 code, ltype, (int) (level), varID, levelID);
481
                      cdoWriteRecord(streamID2, array.data(), nmiss);
482
                    }
483
                  else  // 1:  simple array arithmetics ( *,+)
484
                    {
Oliver Heidmann's avatar
Oliver Heidmann committed
485
                      if (cdoDebugExt)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
486
                        cdoPrint(" Writing record [%4d] with (code %3i, ltype %3i, level %3i)   [varID(%d),levelID(%d)]; "
487
                                 "SCALE=%f; OFFSET=%f",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
488
                                 recID, code, ltype, (int) (level), varID, levelID, scale, offset);
489 490
                      for (size_t li = 0; li < gridsize; ++li)
                        if (!DBL_IS_EQUAL(array[li], missval))
491
                          {
492 493 494
                            array[li] = scale * (array[li] - offset);
                            // If SCALE and/or OFFSET are defined, then the data
                            // values are scaled as SCALE*(VALUE-OFFSET).
495
                          }
496
                      cdoWriteRecord(streamID2, array.data(), nmiss);
497
                    }
498
                }  // end of else ( lcopy )
499
            }      // end if ( vlistInqFlag(vlistID1, varID, levelID) == true )
500
        }          // end for ( recID ..
501

502
      tsID++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
503
    }  // end while
504

505 506
  cdoStreamClose(streamID1);
  cdoStreamClose(streamID2);
507

508
  vlistDestroy(vlistID2);
509

510
  cdoFinish();
511

Oliver Heidmann's avatar
Oliver Heidmann committed
512
  cdoDebugExt = 0;
513

514
  return nullptr;
515 516
}

517 518
TUPLEREC *
TUPLERECNew()
519 520 521 522 523 524 525 526
{
  TUPLEREC *tpl = (TUPLEREC *) malloc(sizeof(TUPLEREC));

  tpl->ncodes = 0;

  tpl->nlevelTypes = 0;

  tpl->nlevels = 0;
527 528
  tpl->sel_or_del_or_change = 0;  // sel_or_del_or_change:  0:  operator
                                  // decides, 1:select , 2:delete, 3:change
529

530
  tpl->simpleMath = 0;  // 1:  simple array arithmetics ( *,+), 0: do nothing
531 532 533
  tpl->scale = 1.0;
  tpl->offset = 0.0;

534
  tpl->changedCode = -1;  // used only changemulti mode
535 536 537
  tpl->changedLevelType = -1;
  tpl->changedLevel = -1;

538
  return tpl;
539 540
}

541 542
void
push_backSelTuple(TUPLEREC *tp)
543
{
544 545 546 547 548
  if (NUMTUPLES < MAX_TUPLES)
    {
      SelTUPLEREC[NUMTUPLES] = tp;
      NUMTUPLES++;
    }
549 550
}

551 552
TUPLEREC *
getSelTuple(int index)
553
{
554 555 556 557
  if (index < NUMTUPLES)
    {
      return SelTUPLEREC[index];
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
558
  return nullptr;
559 560
}

561
int
562
checkListContainsInt(int value, std::vector<int> &list, int arraylen)
563 564
// Note: When the list is Empty it also returns true !
{
565 566 567 568
  int i;

  if (arraylen == 0)
    {
Oliver Heidmann's avatar
Oliver Heidmann committed
569
      // if ( cdoDebugExt ) cdoPrint(" value=%d: found. (list empty)");
570 571 572 573
      return 1;
    }
  for (i = 0; i < arraylen; i++)
    {
574
      if (list[i] == -1)
575 576
        {  // this is for '*' selection; can be any code, any level or any
           // level-type
Oliver Heidmann's avatar
Oliver Heidmann committed
577
          // if ( cdoDebugExt ) cdoPrint(" value=%d: found.");
578 579
          return 1;
        }
580

581
      if (list[i] == value)
582
        {
Oliver Heidmann's avatar
Oliver Heidmann committed
583
          // if ( cdoDebugExt ) cdoPrint(" value=%d: found.");
584 585 586
          return 1;
        }
    }
Oliver Heidmann's avatar
Oliver Heidmann committed
587
  // if ( cdoDebugExt ) cdoPrint(" value=%d: NOT found.");
588 589
  return 0;
}
590 591 592

#define MAX_LINE_LEN 65536

593 594
static char *
removeSpaces(char *pline)
595
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
596
  if (pline == nullptr) return nullptr;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
597
  while (isspace((int) *pline)) pline++;
598
  return pline;
599 600
}

601 602
static char *
skipSeparator(char *pline)
603
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
604
  if (pline == nullptr) return nullptr;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
605
  while (isspace((int) *pline)) pline++;
606
  if (*pline == '=' || *pline == ':' || *pline == '/' || *pline == ',') pline++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
607
  while (isspace((int) *pline)) pline++;
608
  return pline;
609 610
}

611 612
static char *
goToNextSeparator(char *pline)
613
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
614
  if (pline == nullptr) return nullptr;
615 616 617
  int separatorFound = 0;
  while (isspace((int) *pline) || !separatorFound)
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
618
      if (*pline == '\0') return nullptr;
619
      pline++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
620
      if (*pline == '|') return nullptr;
621 622 623 624 625 626 627 628 629
      if (*pline == '=' || *pline == ':' || *pline == '/' || *pline == ',')
        separatorFound = 1;
      else if (*pline == ';')
        {
          pline++;
          pline = removeSpaces(pline);
          return (pline);
        }
    }
630
  if (separatorFound) pline++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
631
  Debug(cdoDebugExt >= 100, "goToNextSeparator():  pline= ('%s') ", pline);
632
  // while ( isspace((int) *pline) ) pline++;
633
  pline = removeSpaces(pline);
634
  return pline;
635 636
}

637 638
static char *
strContains(char *str, const char *substr)
639
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
640 641
  if (str == nullptr) return nullptr;
  if (substr == nullptr) return nullptr;
642 643 644 645 646 647

  str = removeSpaces(str);

  size_t lensub = strlen(substr);
  size_t lenstr = strlen(str);

648 649
  if (lensub > lenstr)
    {
Oliver Heidmann's avatar
Oliver Heidmann committed
650
      if (cdoDebugExt >= 100)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
651
        cdoPrint("strContains():  substr('%s') NOT found in str('%s');  lensub(%zu)>lenstr(%zu) ", substr, str, lensub, lenstr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
652
      return nullptr;
653
    }
654
  char *rv = strstr(str, substr);
655 656
  if (rv)
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
657
      Debug(cdoDebugExt >= 100, "strContains():  substr('%s') FOUND in str('%s')", substr, str);
658 659 660 661
      return (rv + lensub);  // points after subStr ..
    }
  else
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
662
      Debug(cdoDebugExt >= 100, "strContains():  substr('%s') NOT found in str('%s')", substr, str);
663 664
      return rv;
    }
665 666
}

667 668
static char *
findParamEnd(char *str)
669 670
{
  char *ptr = str;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
671
  char *ptrEnding = nullptr;
672

Uwe Schulzweida's avatar
Uwe Schulzweida committed
673
  if (str == nullptr) return nullptr;
674
  // supported endings are: ", " or ";"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
675 676 677
  if (ptrEnding == nullptr) ptrEnding = strContains(ptr, ", ");  // HIP notation
  if (ptrEnding == nullptr) ptrEnding = strContains(ptr, ";");   // compact notation
  if (ptrEnding != nullptr)
678 679
    {
      ptrEnding = removeSpaces(ptrEnding);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
680
      Debug(cdoDebugExt >= 100, " ptrEnding='%s'", ptrEnding);
681 682
      return ptrEnding;
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
683
  Debug(cdoDebugExt >= 100, " ptrEnding=end-of-string");
684
  size_t lenstr = strlen(str);
685
  ptrEnding = str + lenstr;
686 687 688 689
  ptrEnding = removeSpaces(ptrEnding);
  return ptrEnding;
}

690 691
static char *
findTupleEnd(char *str)
692 693 694 695
{
  char *ptr = str;
  char *ptrEnding;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
696
  if (str == nullptr) return nullptr;
697
  // supported endings are: ")" or ");"
698
  ptrEnding = strContains(ptr, ")");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
699 700
  if (ptrEnding == nullptr) ptrEnding = strContains(ptr, ");");
  if (ptrEnding != nullptr)
701 702
    {
      ptrEnding = removeSpaces(ptrEnding);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
703
      Debug(cdoDebugExt >= 100, " findTupleEnd='%s'", ptrEnding);
704 705
      return ptrEnding;
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
706
  Debug(cdoDebugExt >= 100, " findTupleEnd=end-of-string");
707
  size_t lenstr = strlen(str);
708
  ptrEnding = str + lenstr;
709 710 711 712
  ptrEnding = removeSpaces(ptrEnding);
  return ptrEnding;
}

713 714
static char *
readlineForParsing(FILE *gfp, char *strToParsePtr, char *line)
715
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
716
  if (gfp != nullptr)  // file is open => we parse text from a file
717
    {