stream_cdf_i.c 131 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;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
48
49
50
51
52
  bool     islon;
  bool     islat;
  bool     islev;
  bool     istime;
  bool     warn;
53
  bool     calendar;
54
  bool     climatology;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
55
  bool     lformulaterms;
56
  int      tsteptype;
57
  int      param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
58
  int      code;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
59
  int      tabnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
60
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;
  int      tvarid;
72
  int      psvarid;
73
  int      p0varid;
74
  int      ncoordvars;
75
76
77
  int      coordvarids[MAX_COORDVARS];
  int      nauxvars;
  int      auxvarids[MAX_AUXVARS];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
78
79
80
  int      cellarea;
  int      tableID;
  int      truncation;
81
  int      position;
82
83
  bool     defmissval;
  bool     deffillval;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
84
85
86
  int      xtype;
  int      gmapid;
  int      positive;
87
  int      ndims;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
88
89
  int      dimids[8];
  int      dimtype[8];
90
91
92
  int      chunks[8];
  int      chunked;
  int      chunktype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
93
  int      natts;
94
  int      deflate;
95
96
  bool     lunsigned;
  bool     lvalidrange;
Thomas Jahns's avatar
Thomas Jahns committed
97
  int     *atts;
98
99
  size_t   vctsize;
  double  *vct;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
100
  double   missval;
101
  double   fillval;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
102
103
  double   addoffset;
  double   scalefactor;
Deike Kleberg's avatar
Deike Kleberg committed
104
  double   validrange[2];
105
106
107
108
  char     name[CDI_MAX_NAME];
  char     longname[CDI_MAX_NAME];
  char     stdname[CDI_MAX_NAME];
  char     units[CDI_MAX_NAME];
109
  char     extra[CDI_MAX_NAME];
110
  ensinfo_t   *ensdata;    /* Ensemble information */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
111
}
112
113
ncvar_t;

114
115

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

  *rdate = 0;
  *rtime = 0;

125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
  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++;
            }
        }
    }
141
142
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

  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
172
173
174
static
int scanTimeUnit(const char *unitstr)
{
175
  size_t len = strlen(unitstr);
176
  int timeunit = get_timeunit(len, unitstr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
177
178
179
  if ( timeunit == -1 )
    Message("Unsupported TIMEUNIT: %s!", unitstr);

180
  return timeunit;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
181
182
}

183
184
185
static
void setForecastTime(const char *timestr, taxis_t *taxis)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
186
187
188
  (*taxis).fdate = 0;
  (*taxis).ftime = 0;

189
  int len = (int) strlen(timestr);
190
191
192
193
194
195
196
197
198
199
200
201
  if ( len == 0 ) return;

  int fdate = 0, ftime = 0;
  scanTimeString(timestr, &fdate, &ftime);

  (*taxis).fdate = fdate;
  (*taxis).ftime = ftime;
}

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

205
  size_t len = strlen(timeunits);
Thomas Jahns's avatar
Thomas Jahns committed
206
  char *restrict tu = (char *)Malloc((len+1) * sizeof(char));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
207

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

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

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

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

Thomas Jahns's avatar
Thomas Jahns committed
227
228
      while ( ! isspace(tu[pos]) && tu[pos] != 0 ) ++pos;
      if ( tu[pos] )
229
        {
Thomas Jahns's avatar
Thomas Jahns committed
230
          while ( isspace(tu[pos]) ) ++pos;
231
232
233

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

249
250
              (*taxis).rdate = rdate;
              (*taxis).rtime = rtime;
251

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

258
259
  (*taxis).type = timetype;
  (*taxis).unit = timeunit;
Deike Kleberg's avatar
Deike Kleberg committed
260

261
  Free(tu);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
262

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

266
267
  return 0;
}
268

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

326
327
328
329
330
331
  return datatype;
}

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

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

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

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

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

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

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

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

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

412
          if ( nc_attlen > (attlen-1) ) nc_attlen = (attlen-1);
413

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

430
          size_t ssize = strlen(attbuf) + 1;
431

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

445

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

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

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

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

475
476
  int vlistID  = streamptr->vlistID;

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
489
      streamptr->nrecs += nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
490

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

499
      record_t *records = destTstep->records;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
500

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
526
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
527

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

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

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
555
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
556

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

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

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

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

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

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

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

605
  return CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
606
607
}

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

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

696
static
697
void cdf_set_var(ncvar_t *ncvars, int ncvarid, short isvar)
698
{
699
  if ( ncvars[ncvarid].isvar != CDI_UNDEFID &&
700
       ncvars[ncvarid].isvar != isvar   &&
Uwe Schulzweida's avatar
Uwe Schulzweida committed
701
       ncvars[ncvarid].warn  == false )
702
703
704
705
    {
      if ( ! ncvars[ncvarid].ignore )
        Warning("Inconsistent variable definition for %s!", ncvars[ncvarid].name);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
706
      ncvars[ncvarid].warn = true;
707
708
709
710
711
712
713
      isvar = FALSE;
    }

  ncvars[ncvarid].isvar = isvar;
}

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

734
  char attstring[1024];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
735
736
  cdfGetAttText(ncid, ncfvarid, "formula_terms", sizeof(attstring), attstring);
  char *pstring = attstring;
737

738
  bool lstop = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
739
  for ( int i = 0; i < 4; i++ )
740
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
741
742
743
744
745
746
747
748
749
750
751
752
753
754
      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;

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

      if ( lstop ) break;
771
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
772
773
}

774
static
775
bool isHybridSigmaPressureCoordinate(int ncid, int ncvarid, ncvar_t *ncvars, const ncdim_t *ncdims)
776
{
777
  bool status = false;
778
779
780
781
  ncvar_t *ncvar = &ncvars[ncvarid];

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

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

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

811
812
813
814
815
816
              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 ) )
817
                {
818
                  double px = 1;
819
                  if ( p0varid1 != -1 && p0varid1 == p0varid2 )
820
821
                    cdf_get_var_double(ncid, p0varid2, &px);

822
                  double abuf[dimlen*2], bbuf[dimlen*2];
823
                  cdf_get_var_double(ncid, avarid2, abuf);
824
                  cdf_get_var_double(ncid, bvarid2, bbuf);
825

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

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

850
851
852
853
854
855
856
857
858
859
                  ncvar->vct = vct;
                  ncvar->vctsize = vctsize;
                }
            }
        }
    }

  return status;
}

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
899
static
900
void cdf_print_vars(const ncvar_t *ncvars, int nvars, const char *oname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
901
{
902
  char axis[7];
903
  static const char iaxis[] = {'t', 'z', 'y', 'x'};
Uwe Schulzweida's avatar
Uwe Schulzweida committed
904

905
906
  fprintf(stderr, "%s:\n", oname);