cdf_write.c 41.1 KB
Newer Older
1
2
#ifdef HAVE_CONFIG_H
#include "config.h"
3
4
5
6
7
8
9
10
11
12
13
#endif

#ifdef HAVE_LIBNETCDF

#include "dmemory.h"
#include "cdi.h"
#include "cdi_int.h"
#include "stream_cdf.h"
#include "cdf.h"
#include "cdf_int.h"
#include "vlist.h"
14
#include "vlist_var.h"
15
16
17
18


void cdfDefVarDeflate(int ncid, int ncvarid, int deflate_level)
{
19
#if  defined(HAVE_NETCDF4)
20
21
  int retval;
  /* Set chunking, shuffle, and deflate. */
22
23
  const int shuffle = 1;
  const int deflate = 1;
24
25
26

  if ( deflate_level < 1 || deflate_level > 9 ) deflate_level = 1;

27
  if ( (retval = nc_def_var_deflate(ncid, ncvarid, shuffle, deflate, deflate_level)) )
28
29
30
31
    {
      Error("nc_def_var_deflate failed, status = %d", retval);
    }
#else
32
  static bool lwarn = true;
33
34
  if ( lwarn )
    {
35
      lwarn = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
36
      Warning("Deflate compression failed, NetCDF4 not available!");
37
38
39
40
    }
#endif
}

41
nc_type cdfDefDatatype(int datatype, stream_t *streamptr)
42
{
43
  nc_type xtype = NC_FLOAT;
44

45
  if ( streamptr->filetype == CDI_FILETYPE_NC4 )
46
    {
47
48
49
      if      ( datatype == CDI_DATATYPE_INT8   ) xtype = NC_BYTE;
      else if ( datatype == CDI_DATATYPE_INT16  ) xtype = NC_SHORT;
      else if ( datatype == CDI_DATATYPE_INT32  ) xtype = NC_INT;
50
#ifdef  HAVE_NETCDF4
51
52
53
      else if ( datatype == CDI_DATATYPE_UINT8  ) xtype = NC_UBYTE;
      else if ( datatype == CDI_DATATYPE_UINT16 ) xtype = NC_USHORT;
      else if ( datatype == CDI_DATATYPE_UINT32 ) xtype = NC_UINT;
54
      else if ( datatype == CDI_DATATYPE_CPX32 || datatype == CDI_DATATYPE_CPX64 )
55
        Error("CDI library does not support complex numbers with NetCDF4/HDF5!");
56
#else
57
58
59
      else if ( datatype == CDI_DATATYPE_UINT8  ) xtype = NC_SHORT;
      else if ( datatype == CDI_DATATYPE_UINT16 ) xtype = NC_INT;
      else if ( datatype == CDI_DATATYPE_UINT32 ) xtype = NC_INT;
60
61
      else if ( datatype == CDI_DATATYPE_CPX32 || datatype == CDI_DATATYPE_CPX64 )
        Error("CDI library does not support complex numbers with NetCDF4 classic!");
62
#endif
63
      else if ( datatype == CDI_DATATYPE_FLT64  ) xtype = NC_DOUBLE;
64
      else if ( datatype == CDI_DATATYPE_FLT32  ) xtype = NC_FLOAT;
65
66
67
    }
  else
    {
68
69
70
71
72
73
74
      if      ( datatype == CDI_DATATYPE_INT8   ) xtype = NC_BYTE;
      else if ( datatype == CDI_DATATYPE_INT16  ) xtype = NC_SHORT;
      else if ( datatype == CDI_DATATYPE_INT32  ) xtype = NC_INT;
      else if ( datatype == CDI_DATATYPE_UINT8  ) xtype = NC_SHORT;
      else if ( datatype == CDI_DATATYPE_UINT16 ) xtype = NC_INT;
      else if ( datatype == CDI_DATATYPE_UINT32 ) xtype = NC_INT;
      else if ( datatype == CDI_DATATYPE_FLT64  ) xtype = NC_DOUBLE;
75
76
77
      else if ( datatype == CDI_DATATYPE_FLT32  ) xtype = NC_FLOAT;
      else if ( datatype == CDI_DATATYPE_CPX32 || datatype == CDI_DATATYPE_CPX64 )
        Error("CDI library does not support complex numbers with NetCDF classic!");
78
79
80
81
82
83
84
85
    }

  return xtype;
}

static
void cdfDefVarMissval(stream_t *streamptr, int varID, int dtype, int lcheck)
{
86
  if ( streamptr->vars[varID].defmiss == false )
87
    {
88
89
90
91
      const int vlistID = streamptr->vlistID;
      const int fileID  = streamptr->fileID;
      const int ncvarid = streamptr->vars[varID].ncvarid;
      const double missval = vlistInqVarMissval(vlistID, varID);
92
93
94

      if ( lcheck && streamptr->ncmode == 2 ) cdf_redef(fileID);

95
      nc_type xtype = cdfDefDatatype(dtype, streamptr);
96
97
      if ( xtype == NC_BYTE && missval > 127 && missval < 256 ) xtype = NC_INT;

98
99
      if ( lcheck == 0 ||
           streamptr->ncmode != 2 ||
100
           streamptr->filetype == CDI_FILETYPE_NC ||
Uwe Schulzweida's avatar
Uwe Schulzweida committed
101
102
           streamptr->filetype == CDI_FILETYPE_NC2||
           streamptr->filetype == CDI_FILETYPE_NC5 )
103
        cdf_put_att_double(fileID, ncvarid, "_FillValue", xtype, 1, &missval);
104

105
      cdf_put_att_double(fileID, ncvarid, "missing_value", xtype, 1, &missval);
106
107
108

      if ( lcheck && streamptr->ncmode == 2 ) cdf_enddef(fileID);

109
      streamptr->vars[varID].defmiss = true;
110
111
112
113
    }
}

static
114
void cdfDefInstitut(const stream_t *streamptr)
115
{
116
117
118
  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;
  const int instID  = vlistInqInstitut(vlistID);
119

120
  if ( instID != CDI_UNDEFID )
121
122
123
124
    {
      const char *longname = institutInqLongnamePtr(instID);
      if ( longname )
	{
125
	  const size_t len = strlen(longname);
126
127
128
129
130
131
132
133
134
135
136
	  if ( len > 0 )
	    {
	      if ( streamptr->ncmode == 2 ) cdf_redef(fileID);
	      cdf_put_att_text(fileID, NC_GLOBAL, "institution", len, longname);
	      if ( streamptr->ncmode == 2 ) cdf_enddef(fileID);
	    }
	}
    }
}

static
137
void cdfDefSource(const stream_t *streamptr)
138
{
139
140
141
  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;
  const int modelID = vlistInqModel(vlistID);
142

143
  if ( modelID != CDI_UNDEFID )
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
    {
      const char *longname = modelInqNamePtr(modelID);
      if ( longname )
	{
          size_t len = strlen(longname);
	  if ( len > 0 )
	    {
	      if ( streamptr->ncmode == 2 ) cdf_redef(fileID);
	      cdf_put_att_text(fileID, NC_GLOBAL, "source", len, longname);
	      if ( streamptr->ncmode == 2 ) cdf_enddef(fileID);
	    }
	}
    }
}

static inline
void *resizeBuf(void **buf, size_t *bufSize, size_t reqSize)
{
162
  if ( reqSize > *bufSize )
163
164
165
166
167
168
169
    {
      *buf = Realloc(*buf, reqSize);
      *bufSize = reqSize;
    }
  return *buf;
}

170
171

void cdfDefineAttributes(int cdiID, int varID, int fileID, int ncvarID)
172
173
174
175
176
177
178
179
{
  int atttype, attlen;
  size_t len;
  char attname[CDI_MAX_NAME+1];
  void *attBuf = NULL;
  size_t attBufSize = 0;

  int natts;
180
  cdiInqNatts(cdiID, varID, &natts);
181

182
  for ( int iatt = 0; iatt < natts; ++iatt )
183
    {
184
      cdiInqAtt(cdiID, varID, iatt, attname, &atttype, &attlen);
185
186
187

      if ( attlen == 0 ) continue;

188
      if ( atttype == CDI_DATATYPE_TXT )
189
190
191
        {
          size_t attSize = (size_t)attlen*sizeof(char);
          char *atttxt = (char *)resizeBuf(&attBuf, &attBufSize, attSize);
192
          cdiInqAttTxt(cdiID, varID, attname, attlen, atttxt);
193
194
195
          len = (size_t)attlen;
          cdf_put_att_text(fileID, ncvarID, attname, len, atttxt);
        }
196
197
198
      else if ( atttype == CDI_DATATYPE_INT8  || atttype == CDI_DATATYPE_UINT8  ||
                atttype == CDI_DATATYPE_INT16 || atttype == CDI_DATATYPE_UINT16 ||
                atttype == CDI_DATATYPE_INT32 || atttype == CDI_DATATYPE_UINT32 )
199
200
201
        {
          size_t attSize = (size_t)attlen*sizeof(int);
          int *attint = (int *)resizeBuf(&attBuf, &attBufSize, attSize);
202
          cdiInqAttInt(cdiID, varID, attname, attlen, &attint[0]);
203
          len = (size_t)attlen;
204
205
          nc_type xtype = (atttype == CDI_DATATYPE_INT8)  ? NC_BYTE :
                          (atttype == CDI_DATATYPE_INT16) ? NC_SHORT :
206
#if  defined  (HAVE_NETCDF4)
207
208
209
                          (atttype == CDI_DATATYPE_UINT8)  ? NC_UBYTE :
                          (atttype == CDI_DATATYPE_UINT16) ? NC_USHORT :
                          (atttype == CDI_DATATYPE_UINT32) ? NC_UINT :
210
211
212
#endif
                          NC_INT;
          cdf_put_att_int(fileID, ncvarID, attname, xtype, len, attint);
213
        }
214
      else if ( atttype == CDI_DATATYPE_FLT32 || atttype == CDI_DATATYPE_FLT64 )
215
216
217
        {
          size_t attSize = (size_t)attlen * sizeof(double);
          double *attflt = (double *)resizeBuf(&attBuf, &attBufSize, attSize);
218
          cdiInqAttFlt(cdiID, varID, attname, attlen, attflt);
219
          len = (size_t)attlen;
220
          if ( atttype == CDI_DATATYPE_FLT32 )
221
            {
222
223
224
225
226
              float attflt_sp[8];
              float *pattflt_sp = (len > 8) ? (float*) malloc(len*sizeof(float)) : attflt_sp;
              for ( size_t i = 0; i < len; ++i ) pattflt_sp[i] = (float)attflt[i];
              cdf_put_att_float(fileID, ncvarID, attname, NC_FLOAT, len, pattflt_sp);
              if (len > 8) free(pattflt_sp);
227
228
229
230
231
            }
          else
            cdf_put_att_double(fileID, ncvarID, attname, NC_DOUBLE, len, attflt);
        }
    }
232

233
234
235
236
237
238
239
240
  Free(attBuf);
}

static
void cdfDefGlobalAtts(stream_t *streamptr)
{
  if ( streamptr->globalatts ) return;

241
242
  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;
243
244
245
246
247

  cdfDefSource(streamptr);
  cdfDefInstitut(streamptr);

  int natts;
248
  cdiInqNatts(vlistID, CDI_GLOBAL, &natts);
249
250
251

  if ( natts > 0 && streamptr->ncmode == 2 ) cdf_redef(fileID);

252
  cdfDefineAttributes(vlistID, CDI_GLOBAL, fileID, NC_GLOBAL);
253
254
255
256
257
258
259
260
261

  if ( natts > 0 && streamptr->ncmode == 2 ) cdf_enddef(fileID);

  streamptr->globalatts = 1;
}

static
void cdfDefLocalAtts(stream_t *streamptr)
{
262
263
  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;
264
265

  if ( streamptr->localatts ) return;
266
  if ( vlistInqInstitut(vlistID) != CDI_UNDEFID ) return;
267
268
269
270
271
272
273

  streamptr->localatts = 1;

  if ( streamptr->ncmode == 2 ) cdf_redef(fileID);

  for ( int varID = 0; varID < streamptr->nvars; varID++ )
    {
274
      const int instID = vlistInqVarInstitut(vlistID, varID);
275
      if ( instID != CDI_UNDEFID )
276
	{
277
          const int ncvarid = streamptr->vars[varID].ncvarid;
278
279
280
  	  const char *name = institutInqNamePtr(instID);
	  if ( name )
	    {
281
              const size_t len = strlen(name);
282
283
284
285
286
287
288
289
	      cdf_put_att_text(fileID, ncvarid, "institution", len, name);
	    }
	}
      }

  if ( streamptr->ncmode == 2 ) cdf_enddef(fileID);
}

290
291
292
static
void cdf_get_gmapvarname(int gridID, char *gmapvarname)
{
293
  int pgridID = gridID;
294
  char mapping[CDI_MAX_NAME]; mapping[0] = 0;
295
  cdiGridInqKeyStr(pgridID, CDI_KEY_MAPNAME, CDI_MAX_NAME, mapping);
296
297
298
299

  if ( !mapping[0] )
    {
      int projID = gridInqProj(gridID);
300
301
302
      if ( projID != CDI_UNDEFID )
        {
          pgridID = projID;
303
          cdiGridInqKeyStr(pgridID, CDI_KEY_MAPNAME, CDI_MAX_NAME, mapping);
304
        }
305
306
307
    }

  if ( mapping[0] )
308
    cdiGridInqKeyStr(pgridID, CDI_KEY_MAPPING, CDI_MAX_NAME, gmapvarname);
309
310
}

311
312
313
314
315
316
317
318
319
320
321
322
323
324
static
int nc_grid_index(stream_t *streamptr, int gridID)
{
  int index = 0;
  int vlistID = streamptr->vlistID;
  int ngrids = vlistNgrids(vlistID);
  for ( index = 0; index < ngrids; ++index )
    if ( streamptr->ncgrid[index].gridID == gridID ) break;

  assert(index < ngrids);

  return index;
}

325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
static
void cdfDefVarCompression(const stream_t *streamptr, int ncvarid, bool lchunk)
{
  if ( streamptr->comptype == CDI_COMPRESS_ZIP )
    {
      if ( lchunk && (streamptr->filetype == CDI_FILETYPE_NC4 || streamptr->filetype == CDI_FILETYPE_NC4C) )
        {
          cdfDefVarDeflate(streamptr->fileID, ncvarid, streamptr->complevel);
        }
      else
        {
          if ( lchunk )
            {
              static bool lwarn = true;
              if ( lwarn )
                {
                  lwarn = false;
                  Warning("Deflate compression is only available for NetCDF4!");
                }
            }
        }
    }
}

349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
static
void cdfDefVarPacking(const stream_t *streamptr, int ncvarid, nc_type xtype, int vlistID, int varID)
{
  /*  if ( xtype == NC_BYTE || xtype == NC_SHORT || xtype == NC_INT ) */
    {
      const double addoffset   = vlistInqVarAddoffset(vlistID, varID);
      const double scalefactor = vlistInqVarScalefactor(vlistID, varID);
      const bool laddoffset   = IS_NOT_EQUAL(addoffset, 0);
      const bool lscalefactor = IS_NOT_EQUAL(scalefactor, 1);

      if ( laddoffset || lscalefactor )
        {
          nc_type astype = (xtype == NC_FLOAT) ? NC_FLOAT : NC_DOUBLE;
          if ( (astype == NC_DOUBLE) &&
               IS_EQUAL(addoffset,   (double) ((float) addoffset)) &&
               IS_EQUAL(scalefactor, (double) ((float) scalefactor)) )
            {
              astype = NC_FLOAT;
            }

          cdf_put_att_double(streamptr->fileID, ncvarid, "add_offset",   astype, 1, &addoffset);
          cdf_put_att_double(streamptr->fileID, ncvarid, "scale_factor", astype, 1, &scalefactor);
        }
    }
}

375
376
377
378
379
380
381
382
383
384
385
static
void cdfAppendCoordinates(int fileID, int ncvarid, char coordinates[CDI_MAX_NAME])
{
  if (ncvarid != CDI_UNDEFID)
    {
      size_t len = strlen(coordinates);
      if (len) coordinates[len++] = ' ';
      cdf_inq_varname(fileID, ncvarid, coordinates+len);
    }
}

386
static
387
void cdfDefineCoordinates(const stream_t *streamptr, int ncvarid, int nczvarID, int gridtype, int gridID, int gridindex, int xid, int yid, size_t gridsize, char axis[5], size_t iax)
388
{
389
390
  const int fileID  = streamptr->fileID;

391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
  char coordinates[CDI_MAX_NAME]; coordinates[0] = 0;

  if (nczvarID != CDI_UNDEFID) cdfAppendCoordinates(fileID, nczvarID, coordinates);

  if ( gridtype != GRID_GENERIC && gridtype != GRID_LONLAT && gridtype != GRID_GAUSSIAN &&
       gridtype != GRID_PROJECTION && gridtype != GRID_CURVILINEAR && gridtype != GRID_CHARXY )
    {
      const size_t len = strlen(gridNamePtr(gridtype));
      if ( len > 0 ) cdf_put_att_text(fileID, ncvarid, "CDI_grid_type", len, gridNamePtr(gridtype));
    }

  char gmapvarname[CDI_MAX_NAME]; gmapvarname[0] = 0;
  cdf_get_gmapvarname(gridID, gmapvarname);
  if ( gmapvarname[0] ) cdf_put_att_text(fileID, ncvarid, "grid_mapping", strlen(gmapvarname), gmapvarname);

  if ( gridtype == GRID_TRAJECTORY )
    {
      cdf_put_att_text(fileID, ncvarid, "coordinates", 9, "tlon tlat" );
    }
  else if ( gridtype == GRID_LONLAT && xid == CDI_UNDEFID && yid == CDI_UNDEFID && gridsize == 1 )
    {
      const int ncxvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_X];
      const int ncyvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_Y];
      cdfAppendCoordinates(fileID, ncyvarID, coordinates);
      cdfAppendCoordinates(fileID, ncxvarID, coordinates);
    }
  else if ( gridtype == GRID_UNSTRUCTURED || gridtype == GRID_CURVILINEAR )
    {
      char cellarea[CDI_MAX_NAME] = "area: ";
      const int ncxvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_X];
      const int ncyvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_Y];
      const int ncavarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_A];
      // CMOR order: coordinates = "lat lon"
      if ( cdiCoordinatesLonLat )
        {
          cdfAppendCoordinates(fileID, ncxvarID, coordinates);
          cdfAppendCoordinates(fileID, ncyvarID, coordinates);
        }
      else
        {
          cdfAppendCoordinates(fileID, ncyvarID, coordinates);
          cdfAppendCoordinates(fileID, ncxvarID, coordinates);
        }

      if ( ncavarID != CDI_UNDEFID )
        {
          size_t len = strlen(cellarea);
          cdf_inq_varname(fileID, ncavarID, cellarea+len);
          len = strlen(cellarea);
          cdf_put_att_text(fileID, ncvarid, "cell_measures", len, cellarea);
        }

      if ( gridtype == GRID_UNSTRUCTURED )
        {
          int position = gridInqPosition(gridID);
          if ( position > 0 )
            cdf_put_att_int(fileID, ncvarid, "number_of_grid_in_reference", NC_INT, 1, &position);
        }
    }
  else if ( gridtype == GRID_SPECTRAL || gridtype == GRID_FOURIER )
    {
      axis[iax++] = '-';
      axis[iax++] = '-';
      cdf_put_att_text(fileID, ncvarid, "axis", iax, axis);
      int gridTruncation = gridInqTrunc(gridID);
      cdf_put_att_int(fileID, ncvarid, "truncation", NC_INT, 1, &gridTruncation);
    }
  else if ( gridtype == GRID_CHARXY )
    {
      if ( gridInqXIsc(gridID) )
        {
          const int ncxvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_X];
          cdfAppendCoordinates(fileID, ncxvarID, coordinates);
        }
      else if ( gridInqYIsc(gridID) )
        {
          const int ncyvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_Y];
          cdfAppendCoordinates(fileID, ncyvarID, coordinates);
        }
    }

  size_t len = strlen(coordinates);
  if ( len ) cdf_put_att_text(fileID, ncvarid, "coordinates", len, coordinates);
}

476
static
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
int cdfDefineDimsAndChunks(const stream_t *streamptr, int varID, int xid, int yid, int zid, size_t gridsize, const int dimorder[3], int dims[4], bool lchunk, size_t chunks[4], char axis[5], size_t *piax)
{
  const int fileID  = streamptr->fileID;
  const int vlistID = streamptr->vlistID;

  size_t iax = *piax;
  int ndims = 0;

  for (int i = 0; i < 4; ++i) chunks[i] = 0;

  size_t xsize = 0, ysize = 0;
  if ( xid != CDI_UNDEFID ) cdf_inq_dimlen(fileID, xid, &xsize);
  if ( yid != CDI_UNDEFID ) cdf_inq_dimlen(fileID, yid, &ysize);

  const int timetype = vlistInqVarTimetype(vlistID, varID);
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
    {
      const int tid = streamptr->basetime.ncdimid;
      if (tid == CDI_UNDEFID) Error("Internal problem, time undefined!");
      axis[iax++] = 'T';
      chunks[ndims] = 1;
      dims[ndims] = tid;
      ndims++;
    }

  int chunktype = vlistInqVarChunkType(vlistID, varID);
  size_t chunk_size_max = 65536;
  if ( chunktype != CDI_CHUNK_LINES && gridsize > INT32_MAX )
    {
      if ( CDI_Debug ) fprintf(stderr, "gridsize > %d, changed chunktype to CDI_CHUNK_LINES!\n", INT32_MAX);
      chunktype = CDI_CHUNK_LINES;
    }

  for ( int id = 0; id < 3; ++id )
    {
      if ( dimorder[id] == 3 && zid != CDI_UNDEFID )
        {
          axis[iax++] = 'Z';
          chunks[ndims] = 1;
          dims[ndims] = zid;
          ndims++;
        }
      else if ( dimorder[id] == 2 && yid != CDI_UNDEFID )
        {
          if ( chunktype == CDI_CHUNK_AUTO )
            chunks[ndims] = (chunk_size_max > gridsize) ? ysize : chunk_size_max/xsize;
          else
            chunks[ndims] = (chunktype == CDI_CHUNK_LINES) ? 1 : ysize;
          dims[ndims] = yid;
          ndims++;
        }
      else if ( dimorder[id] == 1 && xid != CDI_UNDEFID )
        {
          if ( chunktype == CDI_CHUNK_AUTO && yid == CDI_UNDEFID )
            chunks[ndims] = (chunk_size_max > xsize) ? xsize : chunk_size_max;
          else
            chunks[ndims] = (xsize > INT32_MAX) ? INT32_MAX : xsize;
          dims[ndims] = xid;
          ndims++;
        }
    }

  if ( CDI_Debug )
    fprintf(stderr, "lchunk %d chunktype %d  chunks %zu %zu %zu %zu\n", lchunk, chunktype, chunks[0], chunks[1], chunks[2], chunks[3]);

  iax = *piax;
  return ndims;
}

static
void cdfDefineAttrLeveltype(int fileID, int ncvarid, int zaxisID, int zaxistype)
548
549
550
551
552
553
554
555
556
557
558
559
560
{
  if ( zaxistype == ZAXIS_CLOUD_BASE          ||
       zaxistype == ZAXIS_CLOUD_TOP           ||
       zaxistype == ZAXIS_ISOTHERM_ZERO       ||
       zaxistype == ZAXIS_TOA                 ||
       zaxistype == ZAXIS_SEA_BOTTOM          ||
       zaxistype == ZAXIS_LAKE_BOTTOM         ||
       zaxistype == ZAXIS_SEDIMENT_BOTTOM     ||
       zaxistype == ZAXIS_SEDIMENT_BOTTOM_TA  ||
       zaxistype == ZAXIS_SEDIMENT_BOTTOM_TW  ||
       zaxistype == ZAXIS_MIX_LAYER           ||
       zaxistype == ZAXIS_ATMOSPHERE )
    {
561
      char varname[CDI_MAX_NAME];
562
563
564
565
566
      zaxisInqName(zaxisID, varname);
      cdf_put_att_text(fileID, ncvarid, "level_type", strlen(varname), varname);
    }
}

567
568
569
570
571
572
573
574
575
576
577
578
static
void cdfDefineAttrEnsemble(int fileID, int ncvarid, int vlistID, int varID)
{
  int perturbationNumber, numberOfForecastsInEnsemble, typeOfEnsembleForecast;
  if ( cdiInqKeyInt(vlistID, varID, CDI_KEY_PERTURBATIONNUMBER, &perturbationNumber) == 0 )
    cdf_put_att_int(fileID, ncvarid, "realization", NC_INT, 1, &perturbationNumber);
  if ( cdiInqKeyInt(vlistID, varID, CDI_KEY_NUMBEROFFORECASTSINENSEMBLE, &numberOfForecastsInEnsemble) == 0 )
    cdf_put_att_int(fileID, ncvarid, "ensemble_members", NC_INT, 1, &numberOfForecastsInEnsemble);
  if ( cdiInqKeyInt(vlistID, varID, CDI_KEY_TYPEOFENSEMBLEFORECAST, &typeOfEnsembleForecast) == 0 )
    cdf_put_att_int(fileID, ncvarid, "forecast_init_type", NC_INT, 1, &typeOfEnsembleForecast);
}

579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
static
void cdfCheckVarname(int fileID, char name[CDI_MAX_NAME])
{
  if ( name[0] )
    {
      int ncvarid;
      char varname[CDI_MAX_NAME];
      sprintf(varname, "%s", name);
      char *varname2 = varname+strlen(varname);
      unsigned iz = 0;

      do
        {
          if ( iz ) sprintf(varname2, "_%u", iz+1);

          if ( nc_inq_varid(fileID, varname, &ncvarid) != NC_NOERR ) break;

          ++iz;
        }
      while ( iz <= 99 );

      if (iz > 99) Error("Variable name %s already exsist!", name);

      if ( strcmp(name, varname) != 0 )
        Warning("Changed %s entry of variable name '%s' to '%s'!", (iz==1)?"double":"multiple", name, varname);

      strcpy(name, varname);
    }
}

static
void cdfGenVarname(int fileID, char name[CDI_MAX_NAME], int pnum, int pcat, int *pdis, int *pcode)
{
  char varname[CDI_MAX_NAME];

  int code = *pcode;
  if ( code < 0 ) code = -code;
  if ( pnum < 0 ) pnum = -pnum;

  if ( *pdis == 255 )
    sprintf(varname, "var%d", code);
  else
    sprintf(varname, "param%d.%d.%d", pnum, pcat, *pdis);

  char *varname2 = varname+strlen(varname);
  int ncvarid;
  unsigned iz = 0;

  do
    {
      if ( iz ) sprintf(varname2, "_%u", iz+1);

      if ( nc_inq_varid(fileID, varname, &ncvarid) != NC_NOERR ) break;

      ++iz;
    }
  while ( iz <= 99 );

  if (iz > 99) Error("Variable name %s already exsist!", name);

  strcpy(name, varname);
  *pcode = 0;
  *pdis = 255;
}

644
645
646
static
int cdfDefVar(stream_t *streamptr, int varID)
{
647
648
  if (streamptr->vars[varID].ncvarid != CDI_UNDEFID)
    return streamptr->vars[varID].ncvarid;
649

650
  const int fileID  = streamptr->fileID;
651
  if (CDI_Debug) Message("streamID = %d, fileID = %d, varID = %d", streamptr->self, fileID, varID);
652

653
654
655
  const int vlistID = streamptr->vlistID;
  const int param = vlistInqVarParam(vlistID, varID);
  int code = vlistInqVarCode(vlistID, varID);
656
657
658
  int pnum, pcat, pdis;
  cdiDecodeParam(param, &pnum, &pcat, &pdis);

659
  const int gridID = vlistInqVarGrid(vlistID, varID);
660
661
662
  const size_t gridsize = gridInqSize(gridID);
  const int gridtype = gridInqType(gridID);
  const int gridindex = nc_grid_index(streamptr, gridID);
663
664
  const int xid = (gridtype != GRID_TRAJECTORY) ? streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X] : CDI_UNDEFID;
  const int yid = (gridtype != GRID_TRAJECTORY) ? streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y] : CDI_UNDEFID;
665

666
  const int zaxisID = vlistInqVarZaxis(vlistID, varID);
667
  const int zaxistype = zaxisInqType(zaxisID);
668
669
  const int zaxisindex = vlistZaxisIndex(vlistID, zaxisID);
  const int zid = streamptr->zaxisID[zaxisindex];
670

671
672
673
  int dimorder[3]; // ZYX and ZXY
  vlistInqVarDimorder(vlistID, varID, &dimorder);
  const bool lchunk = (dimorder[0] == 3) ? (gridsize >= 16) : false;
674

675
  if ( ((dimorder[0]>0)+(dimorder[1]>0)+(dimorder[2]>0)) < ((xid!=CDI_UNDEFID)+(yid!=CDI_UNDEFID)+(zid!=CDI_UNDEFID)) )
676
677
678
679
680
    {
      printf("zid=%d  yid=%d  xid=%d\n", zid, yid, xid);
      Error("Internal problem, dimension order missing!");
    }

681
682
683
684
685
  size_t iax = 0;
  char axis[5];
  int dims[4];
  size_t chunks[4];
  int ndims = cdfDefineDimsAndChunks(streamptr, varID, xid, yid, zid, gridsize, dimorder, dims, lchunk, chunks, axis, &iax);
686

687
  char name[CDI_MAX_NAME]; name[0] = 0;
688
689
  if ( vlistInqVarNamePtr(vlistID, varID) ) vlistInqVarName(vlistID, varID, name);

690
691
692
693
694
695
  char longname[CDI_MAX_NAME]; longname[0] = 0;
  char stdname[CDI_MAX_NAME]; stdname[0] = 0;
  char units[CDI_MAX_NAME]; units[0] = 0;
  vlistInqVarLongname(vlistID, varID, longname);
  vlistInqVarStdname(vlistID, varID, stdname);
  vlistInqVarUnits(vlistID, varID, units);
696

697
  const int tableID  = vlistInqVarTable(vlistID, varID);
698
699
700
  if (!name[0]) tableInqEntry(tableID, code, -1, name, longname, units);
  if (name[0]) cdfCheckVarname(fileID, name);
  else         cdfGenVarname(fileID, name, pnum, pcat, &pdis, &code);
701

702
703
  const int dtype = vlistInqVarDatatype(vlistID, varID);
  const nc_type xtype = cdfDefDatatype(dtype, streamptr);
704

705
  int ncvarid = -1;
706
  cdf_def_var(fileID, name, xtype, ndims, dims, &ncvarid);
707
708

#if  defined  (HAVE_NETCDF4)
709
  if ( lchunk && (streamptr->filetype == CDI_FILETYPE_NC4 || streamptr->filetype == CDI_FILETYPE_NC4C) )
710
    cdf_def_var_chunking(fileID, ncvarid, NC_CHUNKED, chunks);
711
712
#endif

713
  cdfDefVarCompression(streamptr, ncvarid, lchunk);
714

715
716
717
  if ( *stdname )  cdf_put_att_text(fileID, ncvarid, "standard_name", strlen(stdname), stdname);
  if ( *longname ) cdf_put_att_text(fileID, ncvarid, "long_name", strlen(longname), longname);
  if ( *units )    cdf_put_att_text(fileID, ncvarid, "units", strlen(units), units);
718
719
720
721
722
723
724
725
726
727
728

  if ( code > 0 && pdis == 255 )
    cdf_put_att_int(fileID, ncvarid, "code", NC_INT, 1, &code);

  if ( pdis != 255 )
    {
      char paramstr[32];
      cdiParamToString(param, paramstr, sizeof(paramstr));
      cdf_put_att_text(fileID, ncvarid, "param", strlen(paramstr), paramstr);
    }

729
  if ( tableID != CDI_UNDEFID )
730
    {
731
732
      int tablenum = tableInqNum(tableID);
      if (tablenum > 0) cdf_put_att_int(fileID, ncvarid, "table", NC_INT, 1, &tablenum);
733
734
    }

735
  const bool zaxisIsScalar = (zid == CDI_UNDEFID) ? (zaxisInqScalar(zaxisID) > 0) : false;
736
  const int nczvarID = (zaxisIsScalar || zaxistype == ZAXIS_CHAR) ? streamptr->nczvarID[zaxisindex] : CDI_UNDEFID;
737

738
  cdfDefineCoordinates(streamptr, ncvarid, nczvarID, gridtype, gridID, gridindex, xid, yid, gridsize, axis, iax);
739

740
  cdfDefVarPacking(streamptr, ncvarid, xtype, vlistID, varID);
741

742
  if ( dtype == CDI_DATATYPE_UINT8 && xtype == NC_BYTE )
743
    {
744
      const int validrange[2] = {0, 255};
745
746
747
748
749
750
751
752
753
      cdf_put_att_int(fileID, ncvarid, "valid_range", NC_SHORT, 2, validrange);
      cdf_put_att_text(fileID, ncvarid, "_Unsigned", 4, "true");
    }

  streamptr->vars[varID].ncvarid = ncvarid;

  if ( vlistInqVarMissvalUsed(vlistID, varID) )
    cdfDefVarMissval(streamptr, varID, vlistInqVarDatatype(vlistID, varID), 0);

754
  if (zid == CDI_UNDEFID) cdfDefineAttrLeveltype(fileID, ncvarid, zaxisID, zaxistype);
755

756
  cdfDefineAttrEnsemble(fileID, ncvarid, vlistID, varID);
757

758
  // Attributes
759
  cdfDefineAttributes(vlistID, varID, fileID, ncvarid);
760
761
762
763
764
765
766
767
768
769
770
771

  return ncvarid;
}


void cdfEndDef(stream_t *streamptr)
{
  cdfDefGlobalAtts(streamptr);
  cdfDefLocalAtts(streamptr);

  if ( streamptr->accessmode == 0 )
    {
772
      const int fileID = streamptr->fileID;
773
774
      if ( streamptr->ncmode == 2 ) cdf_redef(fileID);

775
      const int nvars = streamptr->nvars;
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
      for ( int varID = 0; varID < nvars; varID++ )
	cdfDefVar(streamptr, varID);

      if ( streamptr->ncmode == 2 )
        {
          if ( CDI_netcdf_hdr_pad == 0UL )
            cdf_enddef(fileID);
          else
            cdf__enddef(fileID, CDI_netcdf_hdr_pad);
        }

      streamptr->accessmode = 1;
    }
}

static
void cdfWriteGridTraj(stream_t *streamptr, int gridID)
{
794
795
796
  const int gridindex = nc_grid_index(streamptr, gridID);
  const int lonID = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X];
  const int latID = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y];
797
798
  const int tsID = streamptr->curTsID;
  size_t index = (size_t)tsID;
799
800
801
802

  double xlon = gridInqXval(gridID, 0);
  double xlat = gridInqYval(gridID, 0);

803
804
  cdf_put_var1_double(streamptr->fileID, lonID, &index, &xlon);
  cdf_put_var1_double(streamptr->fileID, latID, &index, &xlat);
805
806
807
}

static
808
void cdf_write_var_data(int fileID, int vlistID, int varID, int ncvarid, int dtype, size_t nvals, size_t xsize, size_t ysize,
809
                        bool swapxy, size_t *start, size_t *count, int memtype, const void *data, size_t nmiss)
810
811
812
813
814
815
816
817
{
  const double *pdata_dp = (const double *) data;
  double *mdata_dp = NULL;
  double *sdata_dp = NULL;
  const float *pdata_sp = (const float *) data;
  float *mdata_sp = NULL;
  float *sdata_sp = NULL;

818
  /*  if ( dtype == CDI_DATATYPE_INT8 || dtype == CDI_DATATYPE_INT16 || dtype == CDI_DATATYPE_INT32 ) */
819
    {
820
821
822
823
824
      const double missval      = vlistInqVarMissval(vlistID, varID);
      const double addoffset    = vlistInqVarAddoffset(vlistID, varID);
      const double scalefactor  = vlistInqVarScalefactor(vlistID, varID);
      const bool laddoffset     = IS_NOT_EQUAL(addoffset, 0);
      const bool lscalefactor   = IS_NOT_EQUAL(scalefactor, 1);
825
826
827
828
829

      if ( laddoffset || lscalefactor )
        {
          if ( memtype == MEMTYPE_FLOAT )
            {
830
831
              mdata_sp = (float *) Malloc(nvals*sizeof(float));
              memcpy(mdata_sp, pdata_sp, nvals*sizeof(float));
832
833
834
835
              pdata_sp = mdata_sp;

              if ( nmiss > 0 )
                {
836
                  for ( size_t i = 0; i < nvals; i++ )
837
838
839
840
841
842
843
844
845
846
847
848
                    {
                      double temp = mdata_sp[i];
                      if ( !DBL_IS_EQUAL(temp, missval) )
                        {
                          if ( laddoffset )   temp -= addoffset;
                          if ( lscalefactor ) temp /= scalefactor;
                          mdata_sp[i] = (float)temp;
                        }
                    }
                }
              else
                {
849
                  for ( size_t i = 0; i < nvals; i++ )
850
851
852
853
854
855
856
857
858
859
                    {
                      double temp = mdata_sp[i];
                      if ( laddoffset )   temp -= addoffset;
                      if ( lscalefactor ) temp /= scalefactor;
                      mdata_sp[i] = (float)temp;
                    }
                }
            }
          else
            {
860
861
              mdata_dp = (double *) Malloc(nvals*sizeof(double));
              memcpy(mdata_dp, pdata_dp, nvals*sizeof(double));
862
863
864
865
              pdata_dp = mdata_dp;

              if ( nmiss > 0 )
                {
866
                  for ( size_t i = 0; i < nvals; i++ )
867
868
869
870
871
872
873
874
875
876
                    {
                      if ( !DBL_IS_EQUAL(mdata_dp[i], missval) )
                        {
                          if ( laddoffset )   mdata_dp[i] -= addoffset;
                          if ( lscalefactor ) mdata_dp[i] /= scalefactor;
                        }
                    }
                }
              else
                {
877
                  for ( size_t i = 0; i < nvals; i++ )
878
879
880
881
882
883
884
885
                    {
                      if ( laddoffset )   mdata_dp[i] -= addoffset;
                      if ( lscalefactor ) mdata_dp[i] /= scalefactor;
                    }
                }
            }
        }

886
887
      if ( dtype == CDI_DATATYPE_UINT8 || dtype == CDI_DATATYPE_INT8 ||
           dtype == CDI_DATATYPE_INT16 || dtype == CDI_DATATYPE_INT32 )
888
889
890
891
892
        {
          if ( memtype == MEMTYPE_FLOAT )
            {
              if ( mdata_sp == NULL )
                {
893
894
                  mdata_sp = (float *) Malloc(nvals*sizeof(float));
                  memcpy(mdata_sp, pdata_sp, nvals*sizeof(float));
895
896
897
                  pdata_sp = mdata_sp;
                }

898
              for ( size_t i = 0; i < nvals; i++ ) mdata_sp[i] = roundf(mdata_sp[i]);
899

900
              if ( dtype == CDI_DATATYPE_UINT8 )
901
902
903
904
905
                {
                  nc_type xtype;
                  cdf_inq_vartype(fileID, ncvarid, &xtype);
                  if ( xtype == NC_BYTE )
                    {
906
                      for ( size_t i = 0; i < nvals; ++i )
907
908
909
910
911
912
913
914
                        if ( mdata_sp[i] > 127 ) mdata_sp[i] -= 256;
                    }
                }
            }
          else
            {
              if ( mdata_dp == NULL )
                {
915
916
                  mdata_dp = (double *) Malloc(nvals*sizeof(double));
                  memcpy(mdata_dp, pdata_dp, nvals*sizeof(double));
917
918
919
                  pdata_dp = mdata_dp;
                }

920
              for ( size_t i = 0; i < nvals; i++ ) mdata_dp[i] = round(mdata_dp[i]);
921

922
              if ( dtype == CDI_DATATYPE_UINT8 )
923
924
925
926
927
                {
                  nc_type xtype;
                  cdf_inq_vartype(fileID, ncvarid, &xtype);
                  if ( xtype == NC_BYTE )
                    {
928
                      for ( size_t i = 0; i < nvals; ++i )
929
930
931
932
933
934
935
936
                        if ( mdata_dp[i] > 127 ) mdata_dp[i] -= 256;
                    }
                }
            }
        }

      if ( CDF_Debug && memtype != MEMTYPE_FLOAT )
        {
937
938
          double fmin =  1.0e200;
          double fmax = -1.0e200;
939
          for ( size_t i = 0; i < nvals; ++i )
940
941
942
943
944
945
946
            {
              if ( !DBL_IS_EQUAL(pdata_dp[i], missval) )
                {
                  if ( pdata_dp[i] < fmin ) fmin = pdata_dp[i];
                  if ( pdata_dp[i] > fmax ) fmax = pdata_dp[i];
                }
            }
947
          Message("nvals = %zu, nmiss = %d, missval = %g, minval = %g, maxval = %g",
948
949
950
951
                  nvals, nmiss, missval, fmin, fmax);
        }
    }

952
  if ( swapxy ) // implemented only for cdf_write_var_slice()
953
    {
954
      size_t gridsize = xsize*ysize;
955
956
      if ( memtype == MEMTYPE_FLOAT )
        {
957
          sdata_sp = (float *) Malloc(gridsize*sizeof(float));
958
959
960
961
962
963
964
          for ( size_t j = 0; j < ysize; ++j )
            for ( size_t i = 0; i < xsize; ++i )
              sdata_sp[i*ysize+j] = pdata_sp[j*xsize+i];
          pdata_sp = sdata_sp;
        }
      else
        {
965
          sdata_dp = (double *) Malloc(gridsize*sizeof (double));
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
          for ( size_t j = 0; j < ysize; ++j )
            for ( size_t i = 0; i < xsize; ++i )
              sdata_dp[i*ysize+j] = pdata_dp[j*xsize+i];
          pdata_dp = sdata_dp;
        }
    }

  if ( memtype == MEMTYPE_FLOAT )
    cdf_put_vara_float(fileID, ncvarid, start, count, pdata_sp);
  else
    cdf_put_vara_double(fileID, ncvarid, start, count, pdata_dp);

  if ( mdata_dp ) Free(mdata_dp);
  if ( sdata_dp ) Free(sdata_dp);
  if ( mdata_sp ) Free(mdata_sp);
  if ( sdata_sp ) Free(sdata_sp);
}

984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
1001
1002
1003
1004
static
void cdfGetXYZid(stream_t *streamptr, int gridID, int zaxisID, int *xid, int *yid, int *zid)
{
  *xid = CDI_UNDEFID;
  *yid = CDI_UNDEFID;

  if (gridInqType(gridID) == GRID_TRAJECTORY)
    {
      cdfWriteGridTraj(streamptr, gridID);
    }
  else
    {
      const int gridindex = nc_grid_index(streamptr, gridID);
      *xid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X];
      *yid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y];
    }

  const int vlistID = streamptr->vlistID;
  const int zaxisindex = vlistZaxisIndex(vlistID, zaxisID);
  *zid = streamptr->zaxisID[zaxisindex];
}
1005

1006
1007
static
void cdfDefineStartAndCount(stream_t *streamptr, int varID, int xid, int yid, int zid, size_t start[5], size_t count[5], size_t *xsize, size_t *ysize)
1008
{
Thomas Jahns's avatar
Thomas Jahns committed
1009
  size_t ndims = 0;
1010
1011
  *xsize = 0;
  *ysize = 0;
1012

1013
1014
  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;
1015

1016
  const long ntsteps = streamptr->ntsteps;
1017
1018
  if ( CDI_Debug ) Message("ntsteps = %ld", ntsteps);

1019
  const int timetype = vlistInqVarTimetype(vlistID, varID);
1020

1021
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
1022
1023
1024
1025
1026
    {
      start[ndims] = (size_t)ntsteps - 1;
      count[ndims] = 1;
      ndims++;
    }
1027

1028
  if ( zid != CDI_UNDEFID )
1029
    {
1030
      const int zaxisID = vlistInqVarZaxis(vlistID, varID);
1031
1032
1033
1034
      start[ndims] = 0;
      count[ndims] = (size_t)zaxisInqSize(zaxisID);
      ndims++;
    }
1035

1036
  if ( yid != CDI_UNDEFID )
1037
1038
    {
      start[ndims] = 0;
1039
      size_t size;
1040
1041
1042
1043
1044
      cdf_inq_dimlen(fileID, yid, &size);
      /*      count[ndims] = gridInqYsize(gridID); */
      count[ndims] = size;
      ndims++;
    }
1045

1046
  if ( xid != CDI_UNDEFID )
1047
1048
    {
      start[ndims] = 0;
1049
      size_t size;
1050
1051
1052
1053
1054
1055
1056
      cdf_inq_dimlen(fileID, xid, &size);
      /*      count[ndims] = gridInqXsize(gridID); */
      count[ndims] = size;
      ndims++;
    }

  if ( CDI_Debug )
Thomas Jahns's avatar
Thomas Jahns committed
1057
    for (size_t idim = 0; idim < ndims; idim++)
1058
      Message("dim = %d  start = %d  count = %d", idim, start[idim], count[idim]);
1059
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069
1070
1071
1072
1073
1074
1075
1076
1077
1078
1079
1080
1081
1082
1083
}

void cdf_write_var(stream_t *streamptr, int varID, int memtype, const void *data, size_t nmiss)
{
  if ( streamptr->accessmode == 0 ) cdfEndDef(streamptr);

  if ( CDI_Debug ) Message("streamID = %d  varID = %d", streamptr->self, varID);

  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;

  const long ntsteps = streamptr->ntsteps;
  if ( CDI_Debug ) Message("ntsteps = %ld", ntsteps);

  const int ncvarid = cdfDefVar(streamptr, varID);

  const int gridID   = vlistInqVarGrid(vlistID, varID);
  const int zaxisID  = vlistInqVarZaxis(vlistID, varID);

  int xid, yid, zid;
  cdfGetXYZid(streamptr, gridID, zaxisID, &xid, &yid, &zid);

  size_t xsize, ysize;
  size_t start[5], count[5];
  cdfDefineStartAndCount(streamptr, varID, xid, yid, zid, start, count, &xsize, &ysize);
1084
1085
1086
1087
1088
1089
1090

  if ( streamptr->ncmode == 1 )
    {
      cdf_enddef(fileID);
      streamptr->ncmode = 2;
    }

1091
  const int dtype = vlistInqVarDatatype(vlistID, varID);
1092
1093
1094

  if ( nmiss > 0 ) cdfDefVarMissval(streamptr, varID, dtype, 1);

1095
  const size_t nvals = gridInqSize(gridID) * (size_t)(zaxisInqSize(zaxisID));
1096

1097
  bool swapxy = false;
1098
1099
1100
1101
1102
  cdf_write_var_data(fileID, vlistID, varID, ncvarid, dtype, nvals, xsize, ysize, swapxy, start, count, memtype, data, nmiss);
}


void cdf_write_var_chunk(stream_t *streamptr, int varID, int memtype,
1103
                         const int rect[][2], const void *data, size_t nmiss)
1104
{
1105
1106
  if ( streamptr->accessmode == 0 ) cdfEndDef(streamptr);

1107
  size_t xsize = 0, ysize = 0;
1108
  size_t start[5],  count[5];
Thomas Jahns's avatar
Thomas Jahns committed
1109
  size_t ndims = 0;
1110
  const int streamID = streamptr->self;
1111
1112
1113
1114

  if ( CDI_Debug )
    Message("streamID = %d  varID = %d", streamID, varID);

1115
1116
  const int vlistID = streamInqVlist(streamID);
  const int fileID  = streamInqFileID(streamID);
1117

1118
  const long ntsteps = streamptr->ntsteps;
1119
  if ( CDI_Debug ) Message("ntsteps = %ld", ntsteps);
1120

1121
  const int ncvarid = cdfDefVar(streamptr, varID);
1122

1123
1124
1125
  const int gridID   = vlistInqVarGrid(vlistID, varID);
  const int zaxisID  = vlistInqVarZaxis(vlistID, varID);
  const int timetype = vlistInqVarTimetype(vlistID, varID);
1126

1127
1128
  int xid, yid, zid;
  cdfGetXYZid(streamptr, gridID, zaxisID, &xid, &yid, &zid);
1129

1130
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
1131
1132
1133
1134
1135
    {
      start[ndims] = (size_t)ntsteps - 1;
      count[ndims] = 1;
      ndims++;
    }
1136
  if ( zid != CDI_UNDEFID )
1137
1138
    {
      int size = zaxisInqSize(zaxisID);
1139
      xassert(rect[2][0] >= 0 && rect[2][0] <= rect[2][1] && rect[2][1] <= size);
1140
1141
1142
1143
      start[ndims] = (size_t)rect[2][0];
      count[ndims] = (size_t)rect[2][1] - (size_t)rect[2][0] + 1;
      ndims++;
    }
1144
  if ( yid != CDI_UNDEFID )
1145
1146
1147
    {
      size_t size;
      cdf_inq_dimlen(fileID, yid, &size);
1148
      xassert(rect[1][0] >= 0 && rect[1][0] <= rect[1][1] && (size_t)rect[1][1] <= size);
1149
1150
1151
1152
      start[ndims] = (size_t)rect[1][0];
      count[ndims] = (size_t)rect[1][1] - (size_t)rect[1][0] + 1;
      ndims++;
    }
1153
  if ( xid != CDI_UNDEFID )
1154
1155
1156
    {
      size_t size;
      cdf_inq_dimlen(fileID, xid, &size);
1157
      xassert(rect[0][0] >= 0 && rect[0][0] <= rect[0][1] && (size_t)rect[0][1] <= size);
1158
1159
1160
1161
1162
1163
      start[ndims] = (size_t)rect[0][0];
      count[ndims] = (size_t)rect[0][1] - (size_t)rect[0][0] + 1;
      ndims++;
    }

  if ( CDI_Debug )
Thomas Jahns's avatar
Thomas Jahns committed
1164
    for (size_t idim = 0; idim < ndims; idim++)
1165
1166
1167
1168
1169
1170
1171
1172
      Message("dim = %d  start = %d  count = %d", idim, start[idim], count[idim]);

  if ( streamptr->ncmode == 1 )
    {
      cdf_enddef(fileID);
      streamptr->ncmode = 2;
    }

1173
  const int dtype = vlistInqVarDatatype(vlistID, varID);
1174
1175
1176

  if ( nmiss > 0 ) cdfDefVarMissval(streamptr, varID, dtype, 1);

1177
  const size_t nvals = gridInqSize(gridID) * (size_t)(zaxisInqSize(zaxisID));
1178

1179
  bool swapxy = false;
1180
1181
1182
1183
  cdf_write_var_data(fileID, vlistID, varID, ncvarid, dtype, nvals,
                     xsize, ysize, swapxy, start, count, memtype, data, nmiss);
}

1184
1185
static
void cdfDefineStartAndCountSlice(stream_t *streamptr, int varID, int levelID, int dimorder[3], int xid, int yid, int zid, size_t start[5], size_t count[5], size_t *xsize, size_t *ysize)
1186
{
1187
1188
1189
  size_t ndims = 0;
  *xsize = 0;
  *ysize = 0;
1190

1191
1192
  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;
1193

1194
  const long ntsteps = streamptr->ntsteps;
1195
1196
  if ( CDI_Debug ) Message("ntsteps = %ld", ntsteps);

1197
  const int timetype = vlistInqVarTimetype(vlistID, varID);
1198

1199
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
1200
1201
1202
1203
1204
1205
1206
1207
    {
      start[ndims] = (size_t)ntsteps - 1;
      count[ndims] = 1;
      ndims++;
    }

  for ( int id = 0; id < 3; ++id )
    {
1208
      if ( dimorder[id] == 3 && zid != CDI_UNDEFID )
1209
1210
1211
1212
1213
        {
          start[ndims] = (size_t)levelID;
          count[ndims] = 1;
          ndims++;
        }
1214
      else if ( dimorder[id] == 2 && yid != CDI_UNDEFID )
1215
1216
        {
          start[ndims] = 0;
1217
1218
          cdf_inq_dimlen(fileID, yid, ysize);
          count[ndims] = *ysize;
1219
1220
          ndims++;
        }
1221
      else if ( dimorder[id] == 1 && xid != CDI_UNDEFID )
1222
1223
        {
          start[ndims] = 0;
1224
1225
          cdf_inq_dimlen(fileID, xid, xsize);
          count[ndims] = *xsize;
1226
1227
1228
1229
1230
1231
1232
          ndims++;
        }
    }

  if ( CDI_Debug )
    for (size_t idim = 0; idim < ndims; idim++)
      Message("dim = %d  start = %d  count = %d", idim, start[idim], count[idim]);
1233
1234
1235
1236
1237
1238
1239
1240
1241
1242
1243
1244
1245
1246
1247
1248
1249
1250
1251
1252
1253
1254
1255
1256
1257
1258
1259
1260
1261
}

void cdf_write_var_slice(stream_t *streamptr, int varID, int levelID, int memtype, const void *data, size_t nmiss)
{
  if ( streamptr->accessmode == 0 ) cdfEndDef(streamptr);

  if ( CDI_Debug ) Message("streamID = %d  varID = %d", streamptr->self, varID);

  const int vlistID = streamptr->vlistID;
  const int fileID  = streamptr->fileID;

  const long ntsteps = streamptr->ntsteps;
  if ( CDI_Debug ) Message("ntsteps = %ld", ntsteps);

  const int ncvarid = cdfDefVar(streamptr, varID);

  const int gridID   = vlistInqVarGrid(vlistID, varID);
  const int zaxisID  = vlistInqVarZaxis(vlistID, varID);

  int xid, yid, zid;
  cdfGetXYZid(streamptr, gridID, zaxisID, &xid, &yid, &zid);

  int dimorder[3];
  vlistInqVarDimorder(vlistID, <