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);
Thomas Jahns's avatar
Thomas Jahns committed
202
  char *restrict tu = (char *)Malloc((len+1) * sizeof(char));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
203

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

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

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

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

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

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

251
252
              taxis->rdate = rdate;
              taxis->rtime = rtime;
253

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

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

263
  Free(tu);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
264

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

268
269
  return 0;
}
270

271
static
272
bool xtypeIsText(int xtype)
273
{
274
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
  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;
314
  else if ( xtype == NC_CHAR   )  datatype = CDI_DATATYPE_UINT8;
315
316
317
318
319
320
321
322
323
324
325
326
  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
327

328
329
330
331
332
333
  return datatype;
}

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

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

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

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

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

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

377
static
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
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
399
void cdfGetAttText(int fileID, int ncvarid, const char *attname, size_t attlen, char *atttext)
400
{
401
402
403
404
405
  nc_type atttype;
  size_t nc_attlen;

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

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

414
          if ( nc_attlen > (attlen-1) ) nc_attlen = (attlen-1);
415

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

432
          size_t ssize = strlen(attbuf) + 1;
433

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

447

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

453
  if ( laddoffset && lscalefactor )
454
    {
455
456
457
458
459
460
461
462
463
464
465
466
      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;
467
468
    }
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
469

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

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

477
478
  int vlistID  = streamptr->vlistID;

479
480
481
482
483
484
485
486
  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
487
488
  if ( tsID == 0 )
    {
489
      int nvrecs = nrecs; /* use all records at first timestep */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
490

Uwe Schulzweida's avatar
Uwe Schulzweida committed
491
      streamptr->nrecs += nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
492

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

501
      record_t *records = destTstep->records;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
502

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
528
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
529

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

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

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
557
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
558

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

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

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

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

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

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

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

607
  return CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
608
609
}

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

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

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

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

  ncvars[ncvarid].isvar = isvar;
}

static
720
void cdf_set_dim(ncvar_t *ncvars, int ncvarid, int dimid, int dimtype)
721
{
722
  if ( ncvars[ncvarid].dimtype[dimid] != CDI_UNDEFID &&
723
724
725
726
727
728
729
730
731
732
       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
733
void scan_hybrid_formulaterms(int ncid, int ncfvarid, int *avarid, int *bvarid, int *psvarid, int *p0varid)
734
{
735
  *avarid  = -1;
736
737
  *bvarid  = -1;
  *psvarid = -1;
738
  *p0varid = -1;
739

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

744
  bool lstop = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
745
  for ( int i = 0; i < 4; i++ )
746
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
747
748
749
750
751
752
753
754
755
756
757
758
759
760
      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;

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

      if ( lstop ) break;
777
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
778
779
}

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

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

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

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

817
818
819
820
821
822
              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 ) )
823
                {
824
                  double px = 1;
825
                  if ( p0varid1 != -1 && p0varid1 == p0varid2 )
826
827
                    cdf_get_var_double(ncid, p0varid2, &px);

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

Thomas Jahns's avatar
Thomas Jahns committed
832
                  size_t vctsize = (dimlen+1)*2;
833
834
835
836
837
838
839
840
841
842
843
844
                  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
845
                    {
846
847
848
849
850
                       for ( size_t i = 0; i < dimlen2; ++i )
                        {
                          vct[i] = abuf[i];
                          vct[i+dimlen+1] = bbuf[i];
                        }
851
852
                    }

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

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

  return status;
}

866
867
868
869
870
871
872
873
874
875
876
877
878
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];
879
      cdfGetAttInt(ncid, ncvarid, attname, attlen, attint);
880
881
      int datatype = (atttype == NC_SHORT)  ? CDI_DATATYPE_INT16 :
                     (atttype == NC_BYTE)   ? CDI_DATATYPE_INT8 :
882
#if  defined  (HAVE_NETCDF4)
883
884
885
                     (atttype == NC_UBYTE)  ? CDI_DATATYPE_UINT8 :
                     (atttype == NC_USHORT) ? CDI_DATATYPE_UINT16 :
                     (atttype == NC_UINT)   ? CDI_DATATYPE_UINT32 :
886
#endif
887
                     CDI_DATATYPE_INT32;
888
      cdiDefAttInt(cdiID, varID, attname, datatype, (int)attlen, attint);
889
890
891
892
    }
  else if ( xtypeIsFloat(atttype) )
    {
      double attflt[attlen];
893
      cdfGetAttDouble(ncid, ncvarid, attname, attlen, attflt);
894
      int datatype = (atttype == NC_FLOAT) ? CDI_DATATYPE_FLT32 : CDI_DATATYPE_FLT64;
895
      cdiDefAttFlt(cdiID, varID, attname, datatype, (int)attlen, attflt);
896
897
898
    }
  else if ( xtypeIsText(atttype) )
    {
899
900
      char attstring[8192];
      cdfGetAttText(ncid, ncvarid, attname, sizeof(attstring), attstring);
901
      cdiDefAttTxt(cdiID, varID, attname, (int)attlen, attstring);
902