stream_cdf_i.c 136 KB
Newer Older
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1
#if defined (HAVE_CONFIG_H)
2
#include "config.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
3
4
#endif

5
6
#ifdef HAVE_LIBNETCDF

7
//#define TEST_GROUPS 1
Uwe Schulzweida's avatar
Uwe Schulzweida committed
8

Uwe Schulzweida's avatar
Uwe Schulzweida committed
9
#include <ctype.h>
10
#include <limits.h>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
11

Uwe Schulzweida's avatar
Uwe Schulzweida committed
12
#include "dmemory.h"
13
#include "gaussgrid.h"
14
#include "cdi_int.h"
15
#include "cdi_uuid.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
16
17
18
19
#include "stream_cdf.h"
#include "cdf_int.h"
#include "varscan.h"
#include "vlist.h"
20
#include "cdf_util.h"
21
#include "cdf_lazy_grid.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
22
23
24
25
26
27
28


#define  X_AXIS  1
#define  Y_AXIS  2
#define  Z_AXIS  3
#define  T_AXIS  4

29
30
31
#define  POSITIVE_UP    1
#define  POSITIVE_DOWN  2

Uwe Schulzweida's avatar
Uwe Schulzweida committed
32
33
34
35
typedef struct {
  int     ncvarid;
  int     dimtype;
  size_t  len;
36
  char    name[CDI_MAX_NAME];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
37
}
38
ncdim_t;
39
40
41
#define  MAX_COORDVARS  4
#define  MAX_AUXVARS    4

Uwe Schulzweida's avatar
Uwe Schulzweida committed
42
typedef struct {
43
  int      ncid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
44
45
  int      isvar;
  bool     ignore;
46
47
  bool     isx;
  bool     isy;
48
  bool     isc;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
49
50
51
52
53
  bool     islon;
  bool     islat;
  bool     islev;
  bool     istime;
  bool     warn;
54
  bool     calendar;
55
  bool     climatology;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
56
  bool     lformulaterms;
57
  int      tsteptype;
58
  int      param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
59
  int      code;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
60
  int      tabnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
61
62
63
64
65
66
67
68
69
70
71
  int      bounds;
  int      gridID;
  int      zaxisID;
  int      gridtype;
  int      zaxistype;
  int      xdim;
  int      ydim;
  int      zdim;
  int      xvarid;
  int      yvarid;
  int      zvarid;
72
  int      cvarids[MAX_COORDVARS];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
73
  int      tvarid;
74
  int      psvarid;
75
  int      p0varid;
76
  int      ncoordvars;
77
78
79
  int      coordvarids[MAX_COORDVARS];
  int      nauxvars;
  int      auxvarids[MAX_AUXVARS];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
80
81
82
  int      cellarea;
  int      tableID;
  int      truncation;
83
  int      position;
84
85
  bool     defmissval;
  bool     deffillval;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
86
87
88
  int      xtype;
  int      gmapid;
  int      positive;
89
  int      ndims;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
90
91
  int      dimids[8];
  int      dimtype[8];
92
93
94
  int      chunks[8];
  int      chunked;
  int      chunktype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
95
  int      natts;
96
  int      deflate;
97
98
  bool     lunsigned;
  bool     lvalidrange;
Thomas Jahns's avatar
Thomas Jahns committed
99
  int     *atts;
100
101
  size_t   vctsize;
  double  *vct;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
102
  double   missval;
103
  double   fillval;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
104
105
  double   addoffset;
  double   scalefactor;
Deike Kleberg's avatar
Deike Kleberg committed
106
  double   validrange[2];
107
108
109
110
  char     name[CDI_MAX_NAME];
  char     longname[CDI_MAX_NAME];
  char     stdname[CDI_MAX_NAME];
  char     units[CDI_MAX_NAME];
111
  char     extra[CDI_MAX_NAME];
112
  ensinfo_t   *ensdata;    /* Ensemble information */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
113
}
114
115
ncvar_t;

116
117

static
118
void scanTimeString(const char *ptu, int *rdate, int *rtime)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
119
{
120
  int year = 1, month = 1, day = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
121
  int hour = 0, minute = 0, second = 0;
122
  int v1 = 1, v2 = 1, v3 = 1;
123
124
125
126

  *rdate = 0;
  *rtime = 0;

127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
  if ( *ptu )
    {
      v1 = atoi(ptu);
      if ( v1 < 0 ) ptu++;
      while ( isdigit((int) *ptu) ) ptu++;
      if ( *ptu )
        {
          v2 = atoi(++ptu);
          while ( isdigit((int) *ptu) ) ptu++;
          if ( *ptu )
            {
              v3 = atoi(++ptu);
              while ( isdigit((int) *ptu) ) ptu++;
            }
        }
    }
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173

  if ( v3 > 999 && v1 < 32 )
    { year = v3; month = v2; day = v1; }
  else
    { year = v1; month = v2; day = v3; }

  while ( isspace((int) *ptu) ) ptu++;

  if ( *ptu )
    {
      while ( ! isdigit((int) *ptu) ) ptu++;

      hour = atoi(ptu);
      while ( isdigit((int) *ptu) ) ptu++;
      if ( *ptu == ':' )
        {
          ptu++;
          minute = atoi(ptu);
          while ( isdigit((int) *ptu) ) ptu++;
          if ( *ptu == ':' )
            {
              ptu++;
              second = atoi(ptu);
            }
        }
    }

  *rdate = cdiEncodeDate(year, month, day);
  *rtime = cdiEncodeTime(hour, minute, second);
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
174
175
176
static
int scanTimeUnit(const char *unitstr)
{
177
  size_t len = strlen(unitstr);
178
  int timeunit = get_timeunit(len, unitstr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
179
180
181
  if ( timeunit == -1 )
    Message("Unsupported TIMEUNIT: %s!", unitstr);

182
  return timeunit;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
183
184
}

185
186
187
static
void setForecastTime(const char *timestr, taxis_t *taxis)
{
Thomas Jahns's avatar
Thomas Jahns committed
188
189
190
191
192
  size_t len = strlen(timestr);
  if ( len != 0 )
    scanTimeString(timestr, &taxis->fdate, &taxis->ftime);
  else
    taxis->fdate = taxis->ftime = 0;
193
194
195
196
197
}

static
int setBaseTime(const char *timeunits, taxis_t *taxis)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
198
199
200
  int timetype = TAXIS_ABSOLUTE;
  int rdate = -1, rtime = -1;

201
  size_t len = strlen(timeunits);
202
203
  while ( isspace(*timeunits) && len ) { timeunits++; len--; }

Thomas Jahns's avatar
Thomas Jahns committed
204
  char *restrict tu = (char *)Malloc((len+1) * sizeof(char));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
205

Thomas Jahns's avatar
Thomas Jahns committed
206
207
  for ( size_t i = 0; i < len; i++ ) tu[i] = (char)tolower((int)timeunits[i]);
  tu[len] = 0;
208

Thomas Jahns's avatar
Thomas Jahns committed
209
  int timeunit = get_timeunit(len, tu);
210
  if ( timeunit == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
211
    {
212
      Message("Unsupported TIMEUNIT: %s!", timeunits);
213
      return 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
214
215
    }

Thomas Jahns's avatar
Thomas Jahns committed
216
  size_t pos = 0;
217
  while ( pos < len && !isspace(tu[pos]) ) ++pos;
Thomas Jahns's avatar
Thomas Jahns committed
218
  if ( tu[pos] )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
219
    {
Thomas Jahns's avatar
Thomas Jahns committed
220
      while ( isspace(tu[pos]) ) ++pos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
221

Thomas Jahns's avatar
Thomas Jahns committed
222
      if ( str_is_equal(tu+pos, "since") )
223
        timetype = TAXIS_RELATIVE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
224

225
      while ( pos < len && !isspace(tu[pos]) ) ++pos;
Thomas Jahns's avatar
Thomas Jahns committed
226
      if ( tu[pos] )
227
        {
Thomas Jahns's avatar
Thomas Jahns committed
228
          while ( isspace(tu[pos]) ) ++pos;
229
230
231

          if ( timetype == TAXIS_ABSOLUTE )
            {
232
              if ( timeunit == TUNIT_DAY )
233
                {
234
235
236
237
238
                  if ( !str_is_equal(tu+pos, "%y%m%d.%f") )
                    {
                      Message("Unsupported format %s for TIMEUNIT day!", tu+pos);
                      timeunit = -1;
                    }
239
                }
240
              else if ( timeunit == TUNIT_MONTH )
241
                {
242
243
244
245
246
                  if ( !str_is_equal(tu+pos, "%y%m.%f") )
                    {
                      Message("Unsupported format %s for TIMEUNIT month!", tu+pos);
                      timeunit = -1;
                    }
247
248
249
250
                }
            }
          else if ( timetype == TAXIS_RELATIVE )
            {
Thomas Jahns's avatar
Thomas Jahns committed
251
              scanTimeString(tu+pos, &rdate, &rtime);
252

253
254
              taxis->rdate = rdate;
              taxis->rtime = rtime;
255

256
257
258
              if ( CDI_Debug )
                Message("rdate = %d  rtime = %d", rdate, rtime);
            }
259
        }
260
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
261

262
263
  taxis->type = timetype;
  taxis->unit = timeunit;
Deike Kleberg's avatar
Deike Kleberg committed
264

265
  Free(tu);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
266

267
268
  if ( CDI_Debug )
    Message("timetype = %d  unit = %d", timetype, timeunit);
Deike Kleberg's avatar
Deike Kleberg committed
269

270
271
  return 0;
}
272

273
static
274
bool xtypeIsText(int xtype)
275
{
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
  bool isText = ( xtype == NC_CHAR )
#if  defined  (HAVE_NETCDF4)
    || ( xtype == NC_STRING )
#endif
    ;
  return isText;
}

static
bool xtypeIsFloat(nc_type xtype)
{
  bool isFloat = xtype == NC_FLOAT || xtype == NC_DOUBLE;

  return isFloat;
}

static
bool xtypeIsInt(nc_type xtype)
{
  bool isInt = xtype == NC_SHORT || xtype == NC_INT
            || xtype == NC_BYTE
#if  defined  (HAVE_NETCDF4)
            || xtype == NC_USHORT || xtype == NC_UINT
            || xtype == NC_UBYTE
#endif
             ;

  return isInt;
}

static
int cdfInqDatatype(int xtype, bool lunsigned)
{
  int datatype = -1;

#if  defined  (HAVE_NETCDF4)
  if ( xtype == NC_BYTE && lunsigned ) xtype = NC_UBYTE;
#endif

  if      ( xtype == NC_BYTE   )  datatype = CDI_DATATYPE_INT8;
316
  else if ( xtype == NC_CHAR   )  datatype = CDI_DATATYPE_UINT8;
317
318
319
320
321
322
323
324
325
326
327
328
  else if ( xtype == NC_SHORT  )  datatype = CDI_DATATYPE_INT16;
  else if ( xtype == NC_INT    )  datatype = CDI_DATATYPE_INT32;
  else if ( xtype == NC_FLOAT  )  datatype = CDI_DATATYPE_FLT32;
  else if ( xtype == NC_DOUBLE )  datatype = CDI_DATATYPE_FLT64;
#if  defined  (HAVE_NETCDF4)
  else if ( xtype == NC_UBYTE  )  datatype = CDI_DATATYPE_UINT8;
  else if ( xtype == NC_LONG   )  datatype = CDI_DATATYPE_INT32;
  else if ( xtype == NC_USHORT )  datatype = CDI_DATATYPE_UINT16;
  else if ( xtype == NC_UINT   )  datatype = CDI_DATATYPE_UINT32;
  else if ( xtype == NC_INT64  )  datatype = CDI_DATATYPE_FLT64;
  else if ( xtype == NC_UINT64 )  datatype = CDI_DATATYPE_FLT64;
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
329

330
331
332
333
334
335
  return datatype;
}

static
void cdfGetAttInt(int fileID, int ncvarid, const char *attname, size_t attlen, int *attint)
{
336
  *attint = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
337

338
339
  nc_type atttype;
  size_t nc_attlen;
340
341
  cdf_inq_atttype(fileID, ncvarid, attname, &atttype);
  cdf_inq_attlen(fileID, ncvarid, attname, &nc_attlen);
Deike Kleberg's avatar
Deike Kleberg committed
342

343
  if ( xtypeIsFloat(atttype) || xtypeIsInt(atttype) )
344
    {
345
346
      bool lalloc = nc_attlen > attlen;
      int *pintatt = lalloc ? (int *)(Malloc(nc_attlen*sizeof(int))) : attint;
347
      cdf_get_att_int(fileID, ncvarid, attname, pintatt);
348
      if ( lalloc )
349
        {
350
          memcpy(attint, pintatt, attlen*sizeof(int));
351
          Free(pintatt);
352
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
353
354
355
    }
}

356
static
357
void cdfGetAttDouble(int fileID, int ncvarid, char *attname, size_t attlen, double *attdouble)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
358
{
359
  *attdouble = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
360

361
362
  nc_type atttype;
  size_t nc_attlen;
363
364
  cdf_inq_atttype(fileID, ncvarid, attname, &atttype);
  cdf_inq_attlen(fileID, ncvarid, attname, &nc_attlen);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
365

366
  if ( xtypeIsFloat(atttype) || xtypeIsInt(atttype) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
367
    {
368
369
      bool lalloc = nc_attlen > attlen;
      double *pdoubleatt = lalloc ? (double*)Malloc(nc_attlen*sizeof(double)) : attdouble;
370
      cdf_get_att_double(fileID, ncvarid, attname, pdoubleatt);
371
      if ( lalloc )
372
        {
373
          memcpy(attdouble, pdoubleatt, attlen*sizeof(double));
374
375
376
377
          Free(pdoubleatt);
        }
    }
}
378

379
static
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
bool cdfCheckAttText(int fileID, int ncvarid, const char *attname)
{
  bool status = false;
  nc_type atttype;

  int status_nc = nc_inq_atttype(fileID, ncvarid, attname, &atttype);

  if ( status_nc == NC_NOERR
       && (atttype == NC_CHAR
#if  defined  (HAVE_NETCDF4)
           || atttype == NC_STRING
#endif
           ) )
    {
      status = true;
    }

  return status;
}

static
401
void cdfGetAttText(int fileID, int ncvarid, const char *attname, size_t attlen, char *atttext)
402
{
403
404
405
406
407
  nc_type atttype;
  size_t nc_attlen;

  cdf_inq_atttype(fileID, ncvarid, attname, &atttype);
  cdf_inq_attlen(fileID, ncvarid, attname, &nc_attlen);
408

409
  if ( atttype == NC_CHAR )
410
    {
411
412
413
414
      char attbuf[65636];
      if ( nc_attlen < sizeof(attbuf) )
        {
          cdf_get_att_text(fileID, ncvarid, attname, attbuf);
415

416
          if ( nc_attlen > (attlen-1) ) nc_attlen = (attlen-1);
417

418
419
420
421
          attbuf[nc_attlen++] = 0;
          memcpy(atttext, attbuf, nc_attlen);
        }
      else
422
        {
423
          atttext[0] = 0;
424
        }
425
    }
426
427
#if  defined  (HAVE_NETCDF4)
  else if ( atttype == NC_STRING )
428
    {
429
      if ( nc_attlen == 1 )
430
        {
431
432
          char *attbuf = NULL;
          cdf_get_att_string(fileID, ncvarid, attname, &attbuf);
433

434
          size_t ssize = strlen(attbuf) + 1;
435

436
          if ( ssize > attlen ) ssize = attlen;
437
438
439
440
441
442
443
          memcpy(atttext, attbuf, ssize);
          atttext[ssize - 1] = 0;
          Free(attbuf);
        }
      else
        {
          atttext[0] = 0;
444
        }
445
    }
446
447
#endif
}
448

449

450
void cdf_scale_add(size_t size, double *data, double addoffset, double scalefactor)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
451
{
452
453
  bool laddoffset   = IS_NOT_EQUAL(addoffset, 0);
  bool lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
454

455
  if ( laddoffset && lscalefactor )
456
    {
457
458
459
460
461
462
463
464
465
466
467
468
      for (size_t i = 0; i < size; ++i )
        data[i] = data[i] * scalefactor + addoffset;
    }
  else if (lscalefactor)
    {
      for (size_t i = 0; i < size; ++i )
        data[i] *= scalefactor;
    }
  else if (laddoffset)
    {
      for (size_t i = 0; i < size; ++i )
        data[i] += addoffset;
469
470
    }
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
471

472
473
static
void cdfCreateRecords(stream_t *streamptr, int tsID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
474
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
475
  if ( tsID < 0 || (tsID >= streamptr->ntsteps && tsID > 0) ) return;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
476

Uwe Schulzweida's avatar
Uwe Schulzweida committed
477
  if ( streamptr->tsteps[tsID].nallrecs > 0 ) return;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
478

479
480
  int vlistID  = streamptr->vlistID;

481
482
483
484
485
486
487
488
  tsteps_t* sourceTstep = streamptr->tsteps;
  tsteps_t* destTstep = sourceTstep + tsID;

  int nvars = vlistNvars(vlistID);
  int nrecs = vlistNrecs(vlistID);

  if ( nrecs <= 0 ) return;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
489
490
  if ( tsID == 0 )
    {
491
      int nvrecs = nrecs; /* use all records at first timestep */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
492

Uwe Schulzweida's avatar
Uwe Schulzweida committed
493
      streamptr->nrecs += nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
494

495
      destTstep->records    = (record_t *) Malloc((size_t)nrecs*sizeof(record_t));
496
497
498
      destTstep->nrecs      = nrecs;
      destTstep->nallrecs   = nrecs;
      destTstep->recordSize = nrecs;
499
      destTstep->curRecID   = CDI_UNDEFID;
500
      destTstep->recIDs     = (int *) Malloc((size_t)nvrecs*sizeof (int));;
501
      for ( int recID = 0; recID < nvrecs; recID++ ) destTstep->recIDs[recID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
502

503
      record_t *records = destTstep->records;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
504

505
      for ( int varID = 0, recID = 0; varID < nvars; varID++ )
506
        {
507
508
509
          int zaxisID = vlistInqVarZaxis(vlistID, varID);
          int nlev    = zaxisInqSize(zaxisID);
          for ( int levelID = 0; levelID < nlev; levelID++ )
510
511
            {
              recordInitEntry(&records[recID]);
512
513
              records[recID].varID   = (short)varID;
              records[recID].levelID = (short)levelID;
514
515
516
              recID++;
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
517
518
519
    }
  else if ( tsID == 1 )
    {
520
521
      int nvrecs = 0;
      for ( int varID = 0; varID < nvars; varID++ )
522
        {
523
          if ( vlistInqVarTsteptype(vlistID, varID) != TSTEP_CONSTANT )
524
            {
525
              int zaxisID = vlistInqVarZaxis(vlistID, varID);
526
527
528
              nvrecs += zaxisInqSize(zaxisID);
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
529

Uwe Schulzweida's avatar
Uwe Schulzweida committed
530
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
531

532
      destTstep->records    = (record_t *) Malloc((size_t)nrecs*sizeof(record_t));
533
534
535
      destTstep->nrecs      = nvrecs;
      destTstep->nallrecs   = nrecs;
      destTstep->recordSize = nrecs;
536
      destTstep->curRecID   = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
537

538
      memcpy(destTstep->records, sourceTstep->records, (size_t)nrecs*sizeof(record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
539

540
      if ( nvrecs )
541
        {
542
          destTstep->recIDs = (int *) Malloc((size_t)nvrecs * sizeof (int));
543
          for ( int recID = 0, vrecID = 0; recID < nrecs; recID++ )
544
            {
545
              int varID = destTstep->records[recID].varID;
546
              if ( vlistInqVarTsteptype(vlistID, varID) != TSTEP_CONSTANT )
547
                {
548
                  destTstep->recIDs[vrecID++] = recID;
549
550
551
                }
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
552
553
554
    }
  else
    {
555
      if ( streamptr->tsteps[1].records == 0 ) cdfCreateRecords(streamptr, 1);
556

557
      int nvrecs = streamptr->tsteps[1].nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
558

Uwe Schulzweida's avatar
Uwe Schulzweida committed
559
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
560

561
      destTstep->records    = (record_t *) Malloc((size_t)nrecs*sizeof(record_t));
562
563
564
      destTstep->nrecs      = nvrecs;
      destTstep->nallrecs   = nrecs;
      destTstep->recordSize = nrecs;
565
      destTstep->curRecID   = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
566

567
      memcpy(destTstep->records, sourceTstep->records, (size_t)nrecs*sizeof(record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
568

569
      destTstep->recIDs     = (int *) Malloc((size_t)nvrecs * sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
570

571
      memcpy(destTstep->recIDs, streamptr->tsteps[1].recIDs, (size_t)nvrecs*sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
572
573
574
    }
}

575
static
576
int cdf_time_dimid(int fileID, int ndims, int nvars)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
577
{
578
579
  char dimname[80];
  for ( int dimid = 0; dimid < ndims; ++dimid )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
580
    {
581
      dimname[0] = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
582
      cdf_inq_dimname(fileID, dimid, dimname);
583
      if ( str_is_equal(dimname, "time") || str_is_equal(dimname, "Time") ) return dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
584
585
    }

586
  for ( int varid = 0; varid < nvars; ++varid )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
587
    {
588
      nc_type xtype;
589
      int nvdims, nvatts, dimids[9];
590
      cdf_inq_var(fileID, varid, NULL, &xtype, &nvdims, dimids, &nvatts);
591
      if ( nvdims == 1 )
592
        {
593
594
          char sbuf[CDI_MAX_NAME];
          for ( int iatt = 0; iatt < nvatts; ++iatt )
595
            {
596
              sbuf[0] = 0;
597
598
              cdf_inq_attname(fileID, varid, iatt, sbuf);
              if ( strncmp(sbuf, "units", 5) == 0 )
599
                {
600
                  cdfGetAttText(fileID, varid, "units", sizeof(sbuf), sbuf);
601
                  str_tolower(sbuf);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
602

603
                  if ( is_time_units(sbuf) ) return dimids[0];
604
605
606
                }
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
607
608
    }

609
  return CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
610
611
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
612
static
613
void init_ncdims(long ndims, ncdim_t *ncdims)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
614
{
615
  for ( long ncdimid = 0; ncdimid < ndims; ncdimid++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
616
    {
617
618
      ncdims[ncdimid].ncvarid      = CDI_UNDEFID;
      ncdims[ncdimid].dimtype      = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
619
620
621
622
623
624
      ncdims[ncdimid].len          = 0;
      ncdims[ncdimid].name[0]      = 0;
    }
}

static
625
void init_ncvars(long nvars, ncvar_t *ncvars)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
626
{
627
  for ( long ncvarid = 0; ncvarid < nvars; ++ncvarid )
Deike Kleberg's avatar
Deike Kleberg committed
628
    {
629
630
      ncvars[ncvarid].ncid            = CDI_UNDEFID;
      ncvars[ncvarid].isvar           = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
631
      ncvars[ncvarid].ignore          = false;
632
633
      ncvars[ncvarid].isx             = false;
      ncvars[ncvarid].isy             = false;
634
      ncvars[ncvarid].isc             = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
635
636
637
638
639
      ncvars[ncvarid].islon           = false;
      ncvars[ncvarid].islat           = false;
      ncvars[ncvarid].islev           = false;
      ncvars[ncvarid].istime          = false;
      ncvars[ncvarid].warn            = false;
640
      ncvars[ncvarid].calendar        = false;
641
      ncvars[ncvarid].climatology     = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
642
      ncvars[ncvarid].lformulaterms   = false;
643
      ncvars[ncvarid].tsteptype       = TSTEP_CONSTANT;
644
645
      ncvars[ncvarid].param           = CDI_UNDEFID;
      ncvars[ncvarid].code            = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
646
      ncvars[ncvarid].tabnum          = 0;
647
648
649
650
651
652
653
654
655
656
657
658
659
660
      ncvars[ncvarid].bounds          = CDI_UNDEFID;
      ncvars[ncvarid].gridID          = CDI_UNDEFID;
      ncvars[ncvarid].zaxisID         = CDI_UNDEFID;
      ncvars[ncvarid].gridtype        = CDI_UNDEFID;
      ncvars[ncvarid].zaxistype       = CDI_UNDEFID;
      ncvars[ncvarid].xdim            = CDI_UNDEFID;
      ncvars[ncvarid].ydim            = CDI_UNDEFID;
      ncvars[ncvarid].zdim            = CDI_UNDEFID;
      ncvars[ncvarid].xvarid          = CDI_UNDEFID;
      ncvars[ncvarid].yvarid          = CDI_UNDEFID;
      ncvars[ncvarid].zvarid          = CDI_UNDEFID;
      ncvars[ncvarid].tvarid          = CDI_UNDEFID;
      ncvars[ncvarid].psvarid         = CDI_UNDEFID;
      ncvars[ncvarid].p0varid         = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
661
      ncvars[ncvarid].ncoordvars      = 0;
662
      for ( int i = 0; i < MAX_COORDVARS; ++i )
663
664
665
666
        {
          ncvars[ncvarid].coordvarids[i]  = CDI_UNDEFID;
          ncvars[ncvarid].cvarids[i]      = CDI_UNDEFID;
        }
667
668
      ncvars[ncvarid].nauxvars      = 0;
      for ( int i = 0; i < MAX_AUXVARS; ++i )
669
670
671
        ncvars[ncvarid].auxvarids[i]  = CDI_UNDEFID;
      ncvars[ncvarid].cellarea        = CDI_UNDEFID;
      ncvars[ncvarid].tableID         = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
672
673
      ncvars[ncvarid].xtype           = 0;
      ncvars[ncvarid].ndims           = 0;
674
      ncvars[ncvarid].gmapid          = CDI_UNDEFID;
675
676
      ncvars[ncvarid].vctsize         = 0;
      ncvars[ncvarid].vct             = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
677
      ncvars[ncvarid].truncation      = 0;
678
      ncvars[ncvarid].position        = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
679
      ncvars[ncvarid].positive        = 0;
680
      ncvars[ncvarid].chunked         = 0;
681
      ncvars[ncvarid].chunktype       = CDI_UNDEFID;
682
683
      ncvars[ncvarid].defmissval      = false;
      ncvars[ncvarid].deffillval      = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
684
      ncvars[ncvarid].missval         = 0;
685
      ncvars[ncvarid].fillval         = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
686
687
688
689
690
      ncvars[ncvarid].addoffset       = 0;
      ncvars[ncvarid].scalefactor     = 1;
      ncvars[ncvarid].natts           = 0;
      ncvars[ncvarid].atts            = NULL;
      ncvars[ncvarid].deflate         = 0;
691
692
      ncvars[ncvarid].lunsigned       = false;
      ncvars[ncvarid].lvalidrange     = false;
Deike Kleberg's avatar
Deike Kleberg committed
693
694
      ncvars[ncvarid].validrange[0]   = VALIDMISS;
      ncvars[ncvarid].validrange[1]   = VALIDMISS;
695
      ncvars[ncvarid].ensdata         = NULL;
696
697
698
699
700
      memset(ncvars[ncvarid].name, 0, CDI_MAX_NAME);
      memset(ncvars[ncvarid].longname, 0, CDI_MAX_NAME);
      memset(ncvars[ncvarid].stdname, 0, CDI_MAX_NAME);
      memset(ncvars[ncvarid].units, 0, CDI_MAX_NAME);
      memset(ncvars[ncvarid].extra, 0, CDI_MAX_NAME);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
701
702
703
    }
}

704
static
705
void cdf_set_var(ncvar_t *ncvars, int ncvarid, short isvar)
706
{
707
  if ( ncvars[ncvarid].isvar != CDI_UNDEFID &&
708
       ncvars[ncvarid].isvar != isvar   &&
Uwe Schulzweida's avatar
Uwe Schulzweida committed
709
       ncvars[ncvarid].warn  == false )
710
711
712
713
    {
      if ( ! ncvars[ncvarid].ignore )
        Warning("Inconsistent variable definition for %s!", ncvars[ncvarid].name);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
714
      ncvars[ncvarid].warn = true;
715
716
717
718
719
720
721
      isvar = FALSE;
    }

  ncvars[ncvarid].isvar = isvar;
}

static
722
void cdf_set_dim(ncvar_t *ncvars, int ncvarid, int dimid, int dimtype)
723
{
724
  if ( ncvars[ncvarid].dimtype[dimid] != CDI_UNDEFID &&
725
726
727
728
729
730
731
732
733
734
       ncvars[ncvarid].dimtype[dimid] != dimtype )
    {
      Warning("Inconsistent dimension definition for %s! dimid = %d;  type = %d;  newtype = %d",
              ncvars[ncvarid].name, dimid, ncvars[ncvarid].dimtype[dimid], dimtype);
    }

  ncvars[ncvarid].dimtype[dimid] = dimtype;
}

static
735
void scan_hybrid_formulaterms(int ncid, int ncfvarid, int *avarid, int *bvarid, int *psvarid, int *p0varid)
736
{
737
  *avarid  = -1;
738
739
  *bvarid  = -1;
  *psvarid = -1;
740
  *p0varid = -1;
741

742
  char attstring[1024];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
743
744
  cdfGetAttText(ncid, ncfvarid, "formula_terms", sizeof(attstring), attstring);
  char *pstring = attstring;
745

746
  bool lstop = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
747
  for ( int i = 0; i < 4; i++ )
748
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
749
750
751
752
753
754
755
756
757
758
759
760
761
762
      while ( isspace((int) *pstring) ) pstring++;
      if ( *pstring == 0 ) break;
      char *tagname = pstring;
      while ( !isspace((int) *pstring) && *pstring != 0 ) pstring++;
      if ( *pstring == 0 ) lstop = true;
      *pstring++ = 0;

      while ( isspace((int) *pstring) ) pstring++;
      if ( *pstring == 0 ) break;
      char *varname = pstring;
      while ( !isspace((int) *pstring) && *pstring != 0 ) pstring++;
      if ( *pstring == 0 ) lstop = true;
      *pstring++ = 0;

763
      int dimvarid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
764
765
      int status_nc = nc_inq_varid(ncid, varname, &dimvarid);
      if ( status_nc == NC_NOERR )
766
        {
767
768
          if      ( strcmp(tagname, "ap:") == 0 ) *avarid  = dimvarid;
          else if ( strcmp(tagname, "a:")  == 0 ) *avarid  = dimvarid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
769
770
771
          else if ( strcmp(tagname, "b:")  == 0 ) *bvarid  = dimvarid;
          else if ( strcmp(tagname, "ps:") == 0 ) *psvarid = dimvarid;
          else if ( strcmp(tagname, "p0:") == 0 ) *p0varid = dimvarid;
772
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
773
774
775
776
777
778
      else if ( strcmp(tagname, "ps:") != 0 )
        {
          Warning("%s - %s", nc_strerror(status_nc), varname);
        }

      if ( lstop ) break;
779
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
780
781
}

782
static
783
bool isHybridSigmaPressureCoordinate(int ncid, int ncvarid, ncvar_t *ncvars, const ncdim_t *ncdims)
784
{
785
  bool status = false;
786
787
788
789
  ncvar_t *ncvar = &ncvars[ncvarid];

  if ( strcmp(ncvar->stdname, "atmosphere_hybrid_sigma_pressure_coordinate") == 0 )
    {
790
791
      cdiConvention = CDI_CONVENTION_CF;

792
      status = true;
793
      ncvar->zaxistype = ZAXIS_HYBRID;
794
      //int ndims = ncvar->ndims;
795
      int dimid = ncvar->dimids[0];
Thomas Jahns's avatar
Thomas Jahns committed
796
      size_t dimlen = ncdims[dimid].len;
797
      int avarid1 = -1, bvarid1 = -1, psvarid1 = -1, p0varid1 = -1;
798
      int ncfvarid = ncvarid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
799
      if ( ncvars[ncfvarid].lformulaterms )
800
801
        scan_hybrid_formulaterms(ncid, ncfvarid, &avarid1, &bvarid1, &psvarid1, &p0varid1);
      // printf("avarid1, bvarid1, psvarid1, p0varid1 %d %d %d %d\n", avarid1, bvarid1, psvarid1, p0varid1);
802
      if ( avarid1  != -1 ) ncvars[avarid1].isvar = FALSE;
803
804
      if ( bvarid1  != -1 ) ncvars[bvarid1].isvar = FALSE;
      if ( psvarid1 != -1 ) ncvar->psvarid = psvarid1;
805
      if ( p0varid1 != -1 ) ncvar->p0varid = p0varid1;
806

807
      if ( ncvar->bounds != CDI_UNDEFID && ncvars[ncvar->bounds].lformulaterms )
808
809
        {
          ncfvarid = ncvar->bounds;
810
          int avarid2 = -1, bvarid2 = -1, psvarid2 = -1, p0varid2 = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
811
          if ( ncvars[ncfvarid].lformulaterms )
812
813
            scan_hybrid_formulaterms(ncid, ncfvarid, &avarid2, &bvarid2, &psvarid2, &p0varid2);
          // printf("avarid2, bvarid2, psvarid2, p0varid2 %d %d %d %d\n", avarid2, bvarid2, psvarid2, p0varid2);
814
          if ( avarid2 != -1 && bvarid2 != -1 )
815
            {
816
              ncvars[avarid2].isvar = FALSE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
817
              ncvars[bvarid2].isvar = FALSE;
818

819
820
821
822
823
824
              int ndims2 = ncvars[avarid2].ndims;
              int dimid2 = ncvars[avarid2].dimids[0];
              size_t dimlen2 = ncdims[dimid2].len;

              if ( (ndims2 == 2 && dimid == ncvars[avarid2].dimids[0] ) ||
                   (ndims2 == 1 && dimlen == dimlen2-1 ) )
825
                {
826
                  double px = 1;
827
                  if ( p0varid1 != -1 && p0varid1 == p0varid2 )
828
829
                    cdf_get_var_double(ncid, p0varid2, &px);

830
                  double abuf[dimlen*2], bbuf[dimlen*2];
831
                  cdf_get_var_double(ncid, avarid2, abuf);
832
                  cdf_get_var_double(ncid, bvarid2, bbuf);
833

Thomas Jahns's avatar
Thomas Jahns committed
834
                  size_t vctsize = (dimlen+1)*2;
835
836
837
838
839
840
841
842
843
844
845
846
                  double *vct = (double *) Malloc(vctsize*sizeof(double));
                  if ( ndims2 == 2 )
                    {
                      for ( size_t i = 0; i < dimlen; ++i )
                        {
                          vct[i] = abuf[i*2];
                          vct[i+dimlen+1] = bbuf[i*2];
                        }
                      vct[dimlen]     = abuf[dimlen*2-1];
                      vct[dimlen*2+1] = bbuf[dimlen*2-1];
                    }
                  else
847
                    {
848
849
850
851
852
                       for ( size_t i = 0; i < dimlen2; ++i )
                        {
                          vct[i] = abuf[i];
                          vct[i+dimlen+1] = bbuf[i];
                        }
853
854
                    }

855
                  if ( p0varid1 != -1 && IS_NOT_EQUAL(px, 1) )
856
857
                    for ( size_t i = 0; i < dimlen+1; ++i ) vct[i] *= px;

858
859
860
861
862
863
864
865
866
867
                  ncvar->vct = vct;
                  ncvar->vctsize = vctsize;
                }
            }
        }
    }

  return status;
}

868
869
870
871
872
873
874
875
876
877
878
879
880
static
void cdf_set_cdi_attr(int ncid, int ncvarid, int attnum, int cdiID, int varID)
{
  nc_type atttype;
  size_t attlen;
  char attname[CDI_MAX_NAME];

  cdf_inq_attname(ncid, ncvarid, attnum, attname);
  cdf_inq_attlen(ncid, ncvarid, attname, &attlen);
  cdf_inq_atttype(ncid, ncvarid, attname, &atttype);
  if ( xtypeIsInt(atttype) )
    {
      int attint[attlen];
881
      cdfGetAttInt(ncid, ncvarid, attname, attlen, attint);
882
883
      int datatype = (atttype == NC_SHORT)  ? CDI_DATATYPE_INT16 :
                     (atttype == NC_BYTE)   ? CDI_DATATYPE_INT8 :
884
#if  defined  (HAVE_NETCDF4)
885
886
887
                     (atttype == NC_UBYTE)  ? CDI_DATATYPE_UINT8 :
                     (atttype == NC_USHORT) ? CDI_DATATYPE_UINT16 :
                     (atttype == NC_UINT)   ? CDI_DATATYPE_UINT32 :
888
#endif
889
                     CDI_DATATYPE_INT32;
890
      cdiDefAttInt(cdiID, varID, attname, datatype, (int)attlen, attint);
891
892
893
894
    }
  else if ( xtypeIsFloat(atttype) )
    {
      double attflt[attlen];
895
      cdfGetAttDouble(ncid, ncvarid, attname, attlen, attflt);
896
      int datatype = (atttype == NC_FLOAT) ? CDI_DATATYPE_FLT32 : CDI_DATATYPE_FLT64;
897
      cdiDefAttFlt(cdiID, varID, attname, datatype, (int)attlen, attflt);
898
899
900
    }
  else if ( xtypeIsText(atttype) )
    {
901
902
      char attstring[8192];
      cdfGetAttText(ncid, ncvarid, attname, sizeof(attstring), attstring);
903
      cdiDefAttTxt(cdiID, varID, attname, (int)attlen, attstring);
904
905
906
    }
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
907
static
908
void cdf_print_vars(const ncvar_t *ncvars, int nvars, const char *oname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
909
{
910
  char axis[7];
911
  static const char iaxis[] = {'t', 'z', 'y', 'x'};
Uwe Schulzweida's avatar
Uwe Schulzweida committed
912

913
914
  fprintf(stderr, "%s:\n", oname);

915
  for ( int ncvarid = 0; ncvarid < nvars; ncvarid++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
916
    {
917
      int ndim = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
918
      if ( ncvars[ncvarid].isvar )
919
920
        {
          axis[ndim++] = 'v';
921
          axis[ndim++] = ':';
922
          for ( int i = 0; i < ncvars[ncvarid].ndims; i++ )
923
924
925
926
927
928
929
            {/*
              if      ( ncvars[ncvarid].tvarid != -1 ) axis[ndim++] = iaxis[0];
              else if ( ncvars[ncvarid].zvarid != -1 ) axis[ndim++] = iaxis[1];
              else if ( ncvars[ncvarid].yvarid != -1 ) axis[ndim++] = iaxis[2];
              else if ( ncvars[ncvarid].xvarid != -1 ) axis[ndim++] = iaxis[3];
              else
             */
930
931
932
933