stream_cdf_i.c 139 KB
Newer Older
1
#ifdef HAVE_CONFIG_H
2
#include "config.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
3
4
#endif

5
6
#ifdef HAVE_LIBNETCDF

Uwe Schulzweida's avatar
Uwe Schulzweida committed
7
#include <ctype.h>
8
#include <limits.h>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
9

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


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

26
27
28
#define  POSITIVE_UP    1
#define  POSITIVE_DOWN  2

Uwe Schulzweida's avatar
Uwe Schulzweida committed
29
typedef struct {
30
  int     dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
31
32
33
  int     ncvarid;
  int     dimtype;
  size_t  len;
34
  char    name[CDI_MAX_NAME];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
35
}
36
ncdim_t;
37

38
39
40
#define  MAX_COORDVARS  4
#define  MAX_AUXVARS    4

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

117
118

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
125
  if ( *ptu ) sscanf(ptu, "%d-%d-%d %d:%d:%d", &v1, &v2, &v3, &hour, &minute, &second);
126
127
128
129
130
131
132
133
134
135

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
136
137
138
static
int scanTimeUnit(const char *unitstr)
{
139
140
141
  const size_t len = strlen(unitstr);
  const int timeunit = get_timeunit(len, unitstr);
  if ( timeunit == -1 ) Message("Unsupported TIMEUNIT: %s!", unitstr);
142
  return timeunit;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
143
144
}

145
146
147
static
void setForecastTime(const char *timestr, taxis_t *taxis)
{
148
  const size_t len = strlen(timestr);
Thomas Jahns's avatar
Thomas Jahns committed
149
150
151
152
  if ( len != 0 )
    scanTimeString(timestr, &taxis->fdate, &taxis->ftime);
  else
    taxis->fdate = taxis->ftime = 0;
153
154
155
156
157
}

static
int setBaseTime(const char *timeunits, taxis_t *taxis)
{
158
  int taxistype = TAXIS_ABSOLUTE;
159
160
  int64_t rdate = -1;
  int rtime = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
161

162
  size_t len = strlen(timeunits);
163
164
  while ( isspace(*timeunits) && len ) { timeunits++; len--; }

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

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

Thomas Jahns's avatar
Thomas Jahns committed
170
  int timeunit = get_timeunit(len, tu);
171
  if ( timeunit == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
172
    {
173
      Message("Unsupported TIMEUNIT: %s!", timeunits);
174
      return 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
175
176
    }

Thomas Jahns's avatar
Thomas Jahns committed
177
  size_t pos = 0;
178
  while ( pos < len && !isspace(tu[pos]) ) ++pos;
Thomas Jahns's avatar
Thomas Jahns committed
179
  if ( tu[pos] )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
180
    {
Thomas Jahns's avatar
Thomas Jahns committed
181
      while ( isspace(tu[pos]) ) ++pos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
182

183
      if ( strStartsWith(tu+pos, "since") )
184
        taxistype = TAXIS_RELATIVE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
185

186
      while ( pos < len && !isspace(tu[pos]) ) ++pos;
Thomas Jahns's avatar
Thomas Jahns committed
187
      if ( tu[pos] )
188
        {
Thomas Jahns's avatar
Thomas Jahns committed
189
          while ( isspace(tu[pos]) ) ++pos;
190

191
          if ( taxistype == TAXIS_ABSOLUTE )
192
            {
193
              if ( timeunit == TUNIT_DAY )
194
                {
195
                  if ( !strStartsWith(tu+pos, "%y%m%d.%f") )
196
197
198
199
                    {
                      Message("Unsupported format %s for TIMEUNIT day!", tu+pos);
                      timeunit = -1;
                    }
200
                }
201
              else if ( timeunit == TUNIT_MONTH )
202
                {
203
                  if ( !strStartsWith(tu+pos, "%y%m.%f") )
204
205
206
207
                    {
                      Message("Unsupported format %s for TIMEUNIT month!", tu+pos);
                      timeunit = -1;
                    }
208
209
                }
            }
210
          else if ( taxistype == TAXIS_RELATIVE )
211
            {
Thomas Jahns's avatar
Thomas Jahns committed
212
              scanTimeString(tu+pos, &rdate, &rtime);
213

214
215
              taxis->rdate = rdate;
              taxis->rtime = rtime;
216

217
              if ( CDI_Debug )
218
                Message("rdate = %lld  rtime = %d", rdate, rtime);
219
            }
220
        }
221
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
222

223
  taxis->type = taxistype;
224
  taxis->unit = timeunit;
Deike Kleberg's avatar
Deike Kleberg committed
225

226
  Free(tu);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
227

228
  if ( CDI_Debug )
229
    Message("taxistype = %d  unit = %d", taxistype, timeunit);
Deike Kleberg's avatar
Deike Kleberg committed
230

231
232
  return 0;
}
233

234
static
235
bool xtypeIsText(int xtype)
236
{
237
  const bool isText = (xtype == NC_CHAR)
238
#ifdef HAVE_NETCDF4
Uwe Schulzweida's avatar
Uwe Schulzweida committed
239
    || (xtype == NC_STRING)
240
241
242
243
244
245
246
247
#endif
    ;
  return isText;
}

static
bool xtypeIsFloat(nc_type xtype)
{
248
  const bool isFloat = xtype == NC_FLOAT || xtype == NC_DOUBLE;
249
250
251
252
253
254
  return isFloat;
}

static
bool xtypeIsInt(nc_type xtype)
{
255
  const bool isInt = xtype == NC_SHORT || xtype == NC_INT
256
            || xtype == NC_BYTE
257
#ifdef HAVE_NETCDF4
258
259
260
261
262
263
264
265
266
267
268
269
270
            || xtype == NC_USHORT || xtype == NC_UINT
            || xtype == NC_UBYTE
#endif
             ;

  return isInt;
}

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

271
#ifdef HAVE_NETCDF4
272
273
274
275
  if ( xtype == NC_BYTE && lunsigned ) xtype = NC_UBYTE;
#endif

  if      ( xtype == NC_BYTE   )  datatype = CDI_DATATYPE_INT8;
276
  else if ( xtype == NC_CHAR   )  datatype = CDI_DATATYPE_UINT8;
277
278
279
280
  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;
281
#ifdef HAVE_NETCDF4
282
283
284
285
286
287
288
  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
289

290
291
292
293
294
295
  return datatype;
}

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

298
299
  nc_type atttype;
  size_t nc_attlen;
300
301
  cdf_inq_atttype(fileID, ncvarid, attname, &atttype);
  cdf_inq_attlen(fileID, ncvarid, attname, &nc_attlen);
Deike Kleberg's avatar
Deike Kleberg committed
302

303
  if ( xtypeIsFloat(atttype) || xtypeIsInt(atttype) )
304
    {
305
      const bool lalloc = nc_attlen > attlen;
306
      int *pintatt = lalloc ? (int *)(Malloc(nc_attlen*sizeof(int))) : attint;
307
      cdf_get_att_int(fileID, ncvarid, attname, pintatt);
308
      if ( lalloc )
309
        {
310
          memcpy(attint, pintatt, attlen*sizeof(int));
311
          Free(pintatt);
312
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
313
314
315
    }
}

316
static
317
void cdfGetAttDouble(int fileID, int ncvarid, char *attname, size_t attlen, double *attdouble)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
318
{
319
  *attdouble = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
320

321
322
  nc_type atttype;
  size_t nc_attlen;
323
324
  cdf_inq_atttype(fileID, ncvarid, attname, &atttype);
  cdf_inq_attlen(fileID, ncvarid, attname, &nc_attlen);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
325

326
  if ( xtypeIsFloat(atttype) || xtypeIsInt(atttype) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
327
    {
328
      const bool lalloc = nc_attlen > attlen;
329
      double *pdoubleatt = lalloc ? (double*)Malloc(nc_attlen*sizeof(double)) : attdouble;
330
      cdf_get_att_double(fileID, ncvarid, attname, pdoubleatt);
331
      if ( lalloc )
332
        {
333
          memcpy(attdouble, pdoubleatt, attlen*sizeof(double));
334
335
336
337
          Free(pdoubleatt);
        }
    }
}
338

339
static
340
341
342
343
344
bool cdfCheckAttText(int fileID, int ncvarid, const char *attname)
{
  bool status = false;
  nc_type atttype;

345
346
  const int status_nc = nc_inq_atttype(fileID, ncvarid, attname, &atttype);
  if ( status_nc == NC_NOERR && (atttype == NC_CHAR
347
#ifdef HAVE_NETCDF4
348
349
350
351
352
353
354
355
356
357
358
           || atttype == NC_STRING
#endif
           ) )
    {
      status = true;
    }

  return status;
}

static
359
void cdfGetAttText(int fileID, int ncvarid, const char *attname, size_t attlen, char *atttext)
360
{
361
362
  atttext[0] = 0;

363
364
365
366
  nc_type atttype;
  size_t nc_attlen;
  cdf_inq_atttype(fileID, ncvarid, attname, &atttype);
  cdf_inq_attlen(fileID, ncvarid, attname, &nc_attlen);
367

368
  if ( atttype == NC_CHAR )
369
    {
370
371
372
373
      char attbuf[65636];
      if ( nc_attlen < sizeof(attbuf) )
        {
          cdf_get_att_text(fileID, ncvarid, attname, attbuf);
374

375
          if ( nc_attlen > (attlen-1) ) nc_attlen = (attlen-1);
376

377
378
379
          attbuf[nc_attlen++] = 0;
          memcpy(atttext, attbuf, nc_attlen);
        }
380
    }
381
#ifdef HAVE_NETCDF4
382
  else if ( atttype == NC_STRING )
383
    {
384
      if ( nc_attlen == 1 )
385
        {
386
387
          char *attbuf = NULL;
          cdf_get_att_string(fileID, ncvarid, attname, &attbuf);
388

389
          size_t ssize = strlen(attbuf) + 1;
390
          if ( ssize > attlen ) ssize = attlen;
391
392
393
394
          memcpy(atttext, attbuf, ssize);
          atttext[ssize - 1] = 0;
          Free(attbuf);
        }
395
    }
396
397
#endif
}
398

399

400
void cdf_scale_add(size_t size, double *data, double addoffset, double scalefactor)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
401
{
402
403
  const bool laddoffset   = IS_NOT_EQUAL(addoffset, 0);
  const bool lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
404

405
  if ( laddoffset && lscalefactor )
406
    {
407
      for (size_t i = 0; i < size; ++i) data[i] = data[i] * scalefactor + addoffset;
408
409
410
    }
  else if (lscalefactor)
    {
411
      for (size_t i = 0; i < size; ++i) data[i] *= scalefactor;
412
413
414
    }
  else if (laddoffset)
    {
415
      for (size_t i = 0; i < size; ++i) data[i] += addoffset;
416
417
    }
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
418

419
420
static
void cdfCreateRecords(stream_t *streamptr, int tsID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
421
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
422
  if ( tsID < 0 || (tsID >= streamptr->ntsteps && tsID > 0) ) return;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
423

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

426
427
  int vlistID  = streamptr->vlistID;

428
429
430
  tsteps_t* sourceTstep = streamptr->tsteps;
  tsteps_t* destTstep = sourceTstep + tsID;

431
432
  const int nvars = vlistNvars(vlistID);
  const int nrecs = vlistNrecs(vlistID);
433
434
  if ( nrecs <= 0 ) return;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
435
436
  if ( tsID == 0 )
    {
437
      const int nvrecs = nrecs; /* use all records at first timestep */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
438

Uwe Schulzweida's avatar
Uwe Schulzweida committed
439
      streamptr->nrecs += nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
440

441
      destTstep->records    = (record_t *) Malloc((size_t)nrecs*sizeof(record_t));
442
443
444
      destTstep->nrecs      = nrecs;
      destTstep->nallrecs   = nrecs;
      destTstep->recordSize = nrecs;
445
      destTstep->curRecID   = CDI_UNDEFID;
446
      destTstep->recIDs     = (int *) Malloc((size_t)nvrecs*sizeof (int));;
447
      for ( int recID = 0; recID < nvrecs; recID++ ) destTstep->recIDs[recID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
448

449
      record_t *records = destTstep->records;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
450

451
      for ( int varID = 0, recID = 0; varID < nvars; varID++ )
452
        {
453
454
          const int zaxisID = vlistInqVarZaxis(vlistID, varID);
          const int nlev    = zaxisInqSize(zaxisID);
455
          for ( int levelID = 0; levelID < nlev; levelID++ )
456
457
            {
              recordInitEntry(&records[recID]);
458
459
              records[recID].varID   = (short)varID;
              records[recID].levelID = (short)levelID;
460
461
462
              recID++;
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
463
464
465
    }
  else if ( tsID == 1 )
    {
466
467
      int nvrecs = 0;
      for ( int varID = 0; varID < nvars; varID++ )
468
        {
469
          if ( vlistInqVarTimetype(vlistID, varID) != TIME_CONSTANT )
470
            {
471
              int zaxisID = vlistInqVarZaxis(vlistID, varID);
472
473
474
              nvrecs += zaxisInqSize(zaxisID);
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
475

Uwe Schulzweida's avatar
Uwe Schulzweida committed
476
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
477

478
      destTstep->records    = (record_t *) Malloc((size_t)nrecs*sizeof(record_t));
479
480
481
      destTstep->nrecs      = nvrecs;
      destTstep->nallrecs   = nrecs;
      destTstep->recordSize = nrecs;
482
      destTstep->curRecID   = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
483

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

486
      if ( nvrecs )
487
        {
488
          destTstep->recIDs = (int *) Malloc((size_t)nvrecs * sizeof (int));
489
          for ( int recID = 0, vrecID = 0; recID < nrecs; recID++ )
490
            {
491
              int varID = destTstep->records[recID].varID;
492
              if ( vlistInqVarTimetype(vlistID, varID) != TIME_CONSTANT )
493
                {
494
                  destTstep->recIDs[vrecID++] = recID;
495
496
497
                }
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
498
499
500
    }
  else
    {
501
      if ( streamptr->tsteps[1].records == 0 ) cdfCreateRecords(streamptr, 1);
502

503
      const int nvrecs = streamptr->tsteps[1].nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
504

Uwe Schulzweida's avatar
Uwe Schulzweida committed
505
      streamptr->nrecs += nvrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
506

507
      destTstep->records    = (record_t *) Malloc((size_t)nrecs*sizeof(record_t));
508
509
510
      destTstep->nrecs      = nvrecs;
      destTstep->nallrecs   = nrecs;
      destTstep->recordSize = nrecs;
511
      destTstep->curRecID   = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
512

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

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

517
      memcpy(destTstep->recIDs, streamptr->tsteps[1].recIDs, (size_t)nvrecs*sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
518
519
520
    }
}

521
static
522
int cdf_time_dimid(int fileID, int ndims, int nvars, ncdim_t *ncdims)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
523
{
524
525
  char dimname[CDI_MAX_NAME];

526
  for ( int dimid = 0; dimid < ndims; ++dimid )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
527
    {
528
      dimname[0] = 0;
529
      cdf_inq_dimname(fileID, ncdims[dimid].dimid, dimname);
530
      if ( strIsEqual(dimname, "time") || strIsEqual(dimname, "Time") ) return dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
531
532
    }

533
  for ( int varid = 0; varid < nvars; ++varid )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
534
    {
535
      nc_type xtype;
536
      int nvdims, nvatts, dimids[9];
537
      cdf_inq_var(fileID, varid, NULL, &xtype, &nvdims, dimids, &nvatts);
538
539
540
541
542
543
544
545
      for ( int i = 0; i < nvdims; ++i )
        for ( int dimid = 0; dimid < ndims; ++dimid )
          if ( ncdims[dimid].dimid == dimids[i] )
            {
              dimids[i] = dimid;
              break;
            }

546
      if ( nvdims == 1 )
547
        {
548
549
          char sbuf[CDI_MAX_NAME];
          for ( int iatt = 0; iatt < nvatts; ++iatt )
550
            {
551
              sbuf[0] = 0;
552
553
              cdf_inq_attname(fileID, varid, iatt, sbuf);
              if ( strncmp(sbuf, "units", 5) == 0 )
554
                {
555
                  cdfGetAttText(fileID, varid, "units", sizeof(sbuf), sbuf);
556
                  strToLower(sbuf);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
557

558
                  if ( is_time_units(sbuf) ) return dimids[0];
559
560
561
                }
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
562
563
    }

564
  return CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
565
566
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
567
static
568
void init_ncdims(long ndims, ncdim_t *ncdims)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
569
{
570
  for ( long ncdimid = 0; ncdimid < ndims; ncdimid++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
571
    {
572
      ncdims[ncdimid].dimid        = CDI_UNDEFID;
573
574
      ncdims[ncdimid].ncvarid      = CDI_UNDEFID;
      ncdims[ncdimid].dimtype      = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
575
576
577
578
579
580
      ncdims[ncdimid].len          = 0;
      ncdims[ncdimid].name[0]      = 0;
    }
}

static
581
void init_ncvars(long nvars, ncvar_t *ncvars)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
582
{
583
  for ( long ncvarid = 0; ncvarid < nvars; ++ncvarid )
Deike Kleberg's avatar
Deike Kleberg committed
584
    {
585
586
      ncvars[ncvarid].ncid            = CDI_UNDEFID;
      ncvars[ncvarid].isvar           = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
587
      ncvars[ncvarid].ignore          = false;
588
589
      ncvars[ncvarid].isx             = false;
      ncvars[ncvarid].isy             = false;
590
      ncvars[ncvarid].isc             = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
591
592
593
594
595
      ncvars[ncvarid].islon           = false;
      ncvars[ncvarid].islat           = false;
      ncvars[ncvarid].islev           = false;
      ncvars[ncvarid].istime          = false;
      ncvars[ncvarid].warn            = false;
596
      ncvars[ncvarid].calendar        = false;
597
      ncvars[ncvarid].climatology     = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
598
      ncvars[ncvarid].lformulaterms   = false;
599
      ncvars[ncvarid].timetype        = TIME_CONSTANT;
600
601
      ncvars[ncvarid].param           = CDI_UNDEFID;
      ncvars[ncvarid].code            = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
602
      ncvars[ncvarid].tabnum          = 0;
603
604
605
606
607
608
609
610
611
612
613
614
615
616
      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
617
      ncvars[ncvarid].ncoordvars      = 0;
618
      for ( int i = 0; i < MAX_COORDVARS; ++i )
619
620
621
622
        {
          ncvars[ncvarid].coordvarids[i]  = CDI_UNDEFID;
          ncvars[ncvarid].cvarids[i]      = CDI_UNDEFID;
        }
623
624
      ncvars[ncvarid].nauxvars      = 0;
      for ( int i = 0; i < MAX_AUXVARS; ++i )
625
626
627
        ncvars[ncvarid].auxvarids[i]  = CDI_UNDEFID;
      ncvars[ncvarid].cellarea        = CDI_UNDEFID;
      ncvars[ncvarid].tableID         = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
628
629
      ncvars[ncvarid].xtype           = 0;
      ncvars[ncvarid].ndims           = 0;
630
      ncvars[ncvarid].gmapid          = CDI_UNDEFID;
631
632
      ncvars[ncvarid].vctsize         = 0;
      ncvars[ncvarid].vct             = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
633
      ncvars[ncvarid].truncation      = 0;
634
      ncvars[ncvarid].position        = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
635
      ncvars[ncvarid].positive        = 0;
636
      ncvars[ncvarid].chunked         = 0;
637
      ncvars[ncvarid].chunktype       = CDI_UNDEFID;
638
639
      ncvars[ncvarid].defmissval      = false;
      ncvars[ncvarid].deffillval      = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
640
      ncvars[ncvarid].missval         = 0;
641
      ncvars[ncvarid].fillval         = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
642
643
644
645
646
      ncvars[ncvarid].addoffset       = 0;
      ncvars[ncvarid].scalefactor     = 1;
      ncvars[ncvarid].natts           = 0;
      ncvars[ncvarid].atts            = NULL;
      ncvars[ncvarid].deflate         = 0;
647
648
      ncvars[ncvarid].lunsigned       = false;
      ncvars[ncvarid].lvalidrange     = false;
Deike Kleberg's avatar
Deike Kleberg committed
649
650
      ncvars[ncvarid].validrange[0]   = VALIDMISS;
      ncvars[ncvarid].validrange[1]   = VALIDMISS;
651
652
653
      ncvars[ncvarid].typeOfEnsembleForecast      = -1;
      ncvars[ncvarid].numberOfForecastsInEnsemble = -1;
      ncvars[ncvarid].perturbationNumber          = -1;
654
655
656
657
658
      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
659
660
661
    }
}

662
static
663
void cdf_set_var(ncvar_t *ncvars, int ncvarid, short isvar)
664
{
665
  if ( ncvars[ncvarid].isvar != CDI_UNDEFID &&
666
       ncvars[ncvarid].isvar != isvar   &&
Uwe Schulzweida's avatar
Uwe Schulzweida committed
667
       ncvars[ncvarid].warn  == false )
668
669
670
671
    {
      if ( ! ncvars[ncvarid].ignore )
        Warning("Inconsistent variable definition for %s!", ncvars[ncvarid].name);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
672
      ncvars[ncvarid].warn = true;
673
674
675
676
677
678
679
      isvar = FALSE;
    }

  ncvars[ncvarid].isvar = isvar;
}

static
680
void cdf_set_dim(ncvar_t *ncvars, int ncvarid, int dimid, int dimtype)
681
{
682
  if ( ncvars[ncvarid].dimtype[dimid] != CDI_UNDEFID &&
683
684
685
686
687
688
689
690
691
692
       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
693
void scan_hybrid_formulaterms(int ncid, int ncfvarid, int *avarid, int *bvarid, int *psvarid, int *p0varid)
694
{
695
  *avarid  = -1;
696
697
  *bvarid  = -1;
  *psvarid = -1;
698
  *p0varid = -1;
699

700
  char attstring[1024];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
701
702
  cdfGetAttText(ncid, ncfvarid, "formula_terms", sizeof(attstring), attstring);
  char *pstring = attstring;
703

704
  bool lstop = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
705
  for ( int i = 0; i < 4; i++ )
706
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
707
708
709
710
711
712
713
714
715
716
717
718
719
720
      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;

721
      int dimvarid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
722
723
      int status_nc = nc_inq_varid(ncid, varname, &dimvarid);
      if ( status_nc == NC_NOERR )
724
        {
725
726
727
728
729
          if      ( strIsEqual(tagname, "ap:") ) *avarid  = dimvarid;
          else if ( strIsEqual(tagname, "a:")  ) *avarid  = dimvarid;
          else if ( strIsEqual(tagname, "b:")  ) *bvarid  = dimvarid;
          else if ( strIsEqual(tagname, "ps:") ) *psvarid = dimvarid;
          else if ( strIsEqual(tagname, "p0:") ) *p0varid = dimvarid;
730
        }
731
      else if ( !strIsEqual(tagname, "ps:") )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
732
733
734
735
736
        {
          Warning("%s - %s", nc_strerror(status_nc), varname);
        }

      if ( lstop ) break;
737
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
738
739
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
static
void readVCT(int ncid, int ndims2, size_t dimlen, size_t dimlen2, int avarid2, int bvarid2, double *vct)
{
  double *abuf = (double*) Malloc(dimlen*2*sizeof(double));
  double *bbuf = (double*) Malloc(dimlen*2*sizeof(double));
  cdf_get_var_double(ncid, avarid2, abuf);
  cdf_get_var_double(ncid, bvarid2, bbuf);

  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
    {
      for ( size_t i = 0; i < dimlen2; ++i )
        {
          vct[i] = abuf[i];
          vct[i+dimlen+1] = bbuf[i];
        }
    }

  Free(abuf);
  Free(bbuf);
}

771
static
772
bool isHybridSigmaPressureCoordinate(int ncid, int ncvarid, ncvar_t *ncvars, const ncdim_t *ncdims)
773
{
774
  bool status = false;
775
776
  ncvar_t *ncvar = &ncvars[ncvarid];

777
  if ( strIsEqual(ncvar->stdname, "atmosphere_hybrid_sigma_pressure_coordinate") )
778
    {
779
      CDI_convention = CDI_CONVENTION_CF;
780

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

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

808
809
810
              const int ndims2 = ncvars[avarid2].ndims;
              const int dimid2 = ncvars[avarid2].dimids[0];
              const size_t dimlen2 = ncdims[dimid2].len;
811
812
813

              if ( (ndims2 == 2 && dimid == ncvars[avarid2].dimids[0] ) ||
                   (ndims2 == 1 && dimlen == dimlen2-1 ) )
814
                {
815
                  double px = 1;
816
                  if ( p0varid1 != -1 && p0varid1 == p0varid2 )
817
818