cdf_write.c 39.9 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
22
23
24
25
26
  int retval;
  /* Set chunking, shuffle, and deflate. */
  int shuffle = 1;
  int deflate = 1;

  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
int cdfDefDatatype(int datatype, stream_t *streamptr)
42
43
44
{
  int xtype = NC_FLOAT;

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
92
93
94
    {
      int vlistID = streamptr->vlistID;
      int fileID  = streamptr->fileID;
      int ncvarid = streamptr->vars[varID].ncvarid;
      double missval = vlistInqVarMissval(vlistID, varID);

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

95
      int xtype = cdfDefDatatype(dtype, streamptr);
96
97
98

      if ( xtype == NC_BYTE && missval > 127 && missval < 256 ) xtype = NC_INT;

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

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

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

110
      streamptr->vars[varID].defmiss = true;
111
112
113
114
115
116
117
118
119
120
    }
}

static
void cdfDefInstitut(stream_t *streamptr)
{
  int vlistID = streamptr->vlistID;
  int fileID  = streamptr->fileID;
  int instID  = vlistInqInstitut(vlistID);

121
  if ( instID != CDI_UNDEFID )
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
    {
      const char *longname = institutInqLongnamePtr(instID);
      if ( longname )
	{
	  size_t len = strlen(longname);
	  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
void cdfDefSource(stream_t *streamptr)
{
  int vlistID = streamptr->vlistID;
  int fileID  = streamptr->fileID;
  int modelID = vlistInqModel(vlistID);

144
  if ( modelID != CDI_UNDEFID )
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
    {
      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)
{
163
  if ( reqSize > *bufSize )
164
165
166
167
168
169
170
    {
      *buf = Realloc(*buf, reqSize);
      *bufSize = reqSize;
    }
  return *buf;
}

171
172

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

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

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

      if ( attlen == 0 ) continue;

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

234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
  Free(attBuf);
}

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

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

  cdfDefSource(streamptr);
  cdfDefInstitut(streamptr);

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

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

253
  cdfDefineAttributes(vlistID, CDI_GLOBAL, fileID, NC_GLOBAL);
254
255
256
257
258
259
260
261
262
263
264
265
266

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

  streamptr->globalatts = 1;
}

static
void cdfDefLocalAtts(stream_t *streamptr)
{
  int vlistID = streamptr->vlistID;
  int fileID  = streamptr->fileID;

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

  streamptr->localatts = 1;

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

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

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

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

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

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

312
313
314
315
316
317
318
319
320
321
322
323
324
325
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;
}

326
327
328
329
static
int cdfDefVar(stream_t *streamptr, int varID)
{
  int ncvarid = -1;
330
  int xid = CDI_UNDEFID, yid = CDI_UNDEFID;
331
332
333
334
335
336
337
338
339
340
341
342
343
344
  size_t xsize = 0, ysize = 0;
  int dims[4];
  size_t chunks[4] = {0,0,0,0};
  int ndims = 0;
  int tablenum;
  int dimorder[3];
  size_t iax = 0;
  char axis[5];

  int fileID  = streamptr->fileID;

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

345
  if ( streamptr->vars[varID].ncvarid != CDI_UNDEFID )
346
347
    return streamptr->vars[varID].ncvarid;

348
349
350
351
352
353
  int vlistID  = streamptr->vlistID;
  int gridID   = vlistInqVarGrid(vlistID, varID);
  int zaxisID  = vlistInqVarZaxis(vlistID, varID);
  int timetype = vlistInqVarTimetype(vlistID, varID);
  int code     = vlistInqVarCode(vlistID, varID);
  int param    = vlistInqVarParam(vlistID, varID);
354
355
356
357
358
359
360
  int pnum, pcat, pdis;
  cdiDecodeParam(param, &pnum, &pcat, &pdis);

  int chunktype = vlistInqVarChunkType(vlistID, varID);

  vlistInqVarDimorder(vlistID, varID, &dimorder);

361
  size_t gridsize = gridInqSize(gridID);
362
  bool lchunk = (gridsize >= 16);
363
  int gridtype  = gridInqType(gridID);
364
  int gridindex = nc_grid_index(streamptr, gridID);
365
366
  if ( gridtype != GRID_TRAJECTORY )
    {
367
      xid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X];
368
      yid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y];
369
370
      if ( xid != CDI_UNDEFID ) cdf_inq_dimlen(fileID, xid, &xsize);
      if ( yid != CDI_UNDEFID ) cdf_inq_dimlen(fileID, yid, &ysize);
371
372
373
374
    }

  int zaxisindex = vlistZaxisIndex(vlistID, zaxisID);
  int zid = streamptr->zaxisID[zaxisindex];
375
  bool zaxis_is_scalar = false;
376
  if ( zid == CDI_UNDEFID ) zaxis_is_scalar = zaxisInqScalar(zaxisID) > 0;
377

378
  if ( dimorder[0] != 3 ) lchunk = false; /* ZYX and ZXY */
379

380
  if ( ((dimorder[0]>0)+(dimorder[1]>0)+(dimorder[2]>0)) < ((xid!=CDI_UNDEFID)+(yid!=CDI_UNDEFID)+(zid!=CDI_UNDEFID)) )
381
382
383
384
385
386
387
    {
      printf("zid=%d  yid=%d  xid=%d\n", zid, yid, xid);
      Error("Internal problem, dimension order missing!");
    }

  int tid = streamptr->basetime.ncdimid;

388
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
389
    {
390
      if ( tid == CDI_UNDEFID ) Error("Internal problem, time undefined!");
391
392
393
394
395
      chunks[ndims] = 1;
      dims[ndims++] = tid;
      axis[iax++] = 'T';
    }
  /*
396
397
398
  if ( zid != CDI_UNDEFID ) axis[iax++] = 'Z';
  if ( zid != CDI_UNDEFID ) chunks[ndims] = 1;
  if ( zid != CDI_UNDEFID ) dims[ndims++] = zid;
399

400
401
  if ( yid != CDI_UNDEFID ) chunks[ndims] = ysize;
  if ( yid != CDI_UNDEFID ) dims[ndims++] = yid;
402

403
404
  if ( xid != CDI_UNDEFID ) chunks[ndims] = xsize;
  if ( xid != CDI_UNDEFID ) dims[ndims++] = xid;
405
  */
406
  size_t chunk_size_max = 65536;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
407
408
409
410
411
412
  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;
    }

413
414
  for ( int id = 0; id < 3; ++id )
    {
415
      if ( dimorder[id] == 3 && zid != CDI_UNDEFID )
416
417
418
419
420
421
        {
          axis[iax++] = 'Z';
          chunks[ndims] = 1;
          dims[ndims] = zid;
          ndims++;
        }
422
      else if ( dimorder[id] == 2 && yid != CDI_UNDEFID )
423
        {
424
425
426
427
          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;
428
429
430
          dims[ndims] = yid;
          ndims++;
        }
431
      else if ( dimorder[id] == 1 && xid != CDI_UNDEFID )
432
        {
433
434
435
          if ( chunktype == CDI_CHUNK_AUTO && yid == CDI_UNDEFID )
            chunks[ndims] = (chunk_size_max > xsize) ? xsize : chunk_size_max;
          else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
436
            chunks[ndims] = (xsize > INT32_MAX) ? INT32_MAX : xsize;
437
438
439
440
441
442
          dims[ndims] = xid;
          ndims++;
        }
    }

  if ( CDI_Debug )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
443
    fprintf(stderr, "chunktype %d  chunks %zu %zu %zu %zu\n", chunktype, chunks[0], chunks[1], chunks[2], chunks[3]);
444

445
446
447
448
449
450
451
452
453
  char varname[CDI_MAX_NAME];
  char name[CDI_MAX_NAME]; name[0] = 0;
  char longname[CDI_MAX_NAME]; longname[0] = 0;
  char stdname[CDI_MAX_NAME]; stdname[0] = 0;
  char units[CDI_MAX_NAME]; units[0] = 0;
  if ( vlistInqVarNamePtr(vlistID, varID) ) vlistInqVarName(vlistID, varID, name);
  vlistInqVarLongname(vlistID, varID, longname);
  vlistInqVarStdname(vlistID, varID, stdname);
  vlistInqVarUnits(vlistID, varID, units);
454

455
  int tableID  = vlistInqVarTable(vlistID, varID);
456
  if ( !name[0] ) tableInqEntry(tableID, code, -1, name, longname, units);
457
  if ( name[0] )
458
459
    {
      sprintf(varname, "%s", name);
460
461
      char *varname2 = varname+strlen(varname);
      unsigned iz = 0;
462

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

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

469
          ++iz;
470
        }
471
472
473
      while ( iz <= 99 );

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

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

478
      strcpy(name, varname);
479
480
481
482
483
484
485
486
487
488
489
490
    }
  else
    {
      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);
491
      unsigned iz = 0;
492

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

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

499
          ++iz;
500
        }
501
502
503
      while ( iz <= 99 );

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

505
      strcpy(name, varname);
506
507
508
509
510
511
512
      code = 0;
      pdis = 255;
    }

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

  int dtype = vlistInqVarDatatype(vlistID, varID);
513
  int xtype = cdfDefDatatype(dtype, streamptr);
514
515
516
517

  cdf_def_var(fileID, name, (nc_type) xtype, ndims, dims, &ncvarid);

#if  defined  (HAVE_NETCDF4)
518
  if ( lchunk && (streamptr->filetype == CDI_FILETYPE_NC4 || streamptr->filetype == CDI_FILETYPE_NC4C) )
519
    cdf_def_var_chunking(fileID, ncvarid, NC_CHUNKED, chunks);
520
521
#endif

522
  if ( streamptr->comptype == CDI_COMPRESS_ZIP )
523
    {
524
      if ( lchunk && (streamptr->filetype == CDI_FILETYPE_NC4 || streamptr->filetype == CDI_FILETYPE_NC4C) )
525
526
527
528
529
530
531
        {
          cdfDefVarDeflate(fileID, ncvarid, streamptr->complevel);
        }
      else
        {
          if ( lchunk )
            {
532
              static bool lwarn = true;
533
534
              if ( lwarn )
                {
535
                  lwarn = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
536
                  Warning("Deflate compression is only available for NetCDF4!");
537
538
539
540
541
                }
            }
        }
    }

542
543
544
  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);
545
546
547
548
549
550
551
552
553
554
555

  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);
    }

556
  if ( tableID != CDI_UNDEFID )
557
558
559
560
561
562
    {
      tablenum = tableInqNum(tableID);
      if ( tablenum > 0 )
        cdf_put_att_int(fileID, ncvarid, "table", NC_INT, 1, &tablenum);
    }

563
  char coordinates[CDI_MAX_NAME]; coordinates[0] = 0;
564

565
  if ( zaxis_is_scalar || zaxisInqType(zaxisID) == ZAXIS_CHAR )
566
567
568
569
570
571
572
573
574
575
    {
      int nczvarID = streamptr->nczvarID[zaxisindex];
      if ( nczvarID != CDI_UNDEFID )
        {
          size_t len = strlen(coordinates);
          if ( len ) coordinates[len++] = ' ';
          cdf_inq_varname(fileID, nczvarID, coordinates+len);
        }
    }

576
  if ( gridtype != GRID_GENERIC && gridtype != GRID_LONLAT && gridtype != GRID_GAUSSIAN &&
577
       gridtype != GRID_PROJECTION && gridtype != GRID_CURVILINEAR && gridtype != GRID_CHARXY )
578
579
580
    {
      size_t len = strlen(gridNamePtr(gridtype));
      if ( len > 0 )
581
        cdf_put_att_text(fileID, ncvarid, "CDI_grid_type", len, gridNamePtr(gridtype));
582
583
    }

584
  char gmapvarname[CDI_MAX_NAME]; gmapvarname[0] = 0;
585

586
  cdf_get_gmapvarname(gridID, gmapvarname);
587
588
589
590

  if ( gmapvarname[0] ) cdf_put_att_text(fileID, ncvarid, "grid_mapping", strlen(gmapvarname), gmapvarname);

  if ( gridtype == GRID_TRAJECTORY )
591
    {
592
      cdf_put_att_text(fileID, ncvarid, "coordinates", 9, "tlon tlat" );
593
    }
594
  else if ( gridtype == GRID_LONLAT && xid == CDI_UNDEFID && yid == CDI_UNDEFID && gridsize == 1 )
595
    {
596
      int ncxvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_X];
597
      int ncyvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_Y];
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
      if ( ncyvarID != CDI_UNDEFID )
        {
          size_t len = strlen(coordinates);
          if ( len ) coordinates[len++] = ' ';
          cdf_inq_varname(fileID, ncyvarID, coordinates+len);
        }
      if ( ncxvarID != CDI_UNDEFID )
        {
          size_t len = strlen(coordinates);
          if ( len ) coordinates[len++] = ' ';
          cdf_inq_varname(fileID, ncxvarID, coordinates+len);
        }
    }
  else if ( gridtype == GRID_UNSTRUCTURED || gridtype == GRID_CURVILINEAR )
    {
      char cellarea[CDI_MAX_NAME] = "area: ";
614
      int ncxvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_X];
615
      int ncyvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_Y];
616
      int ncavarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_A];
617
618
      // CMOR order: coordinates = "lat lon"
      if ( cdiCoordinatesLonLat )
619
        {
620
621
622
623
624
625
626
627
628
629
630
631
          if ( ncxvarID != CDI_UNDEFID )
            {
              size_t len = strlen(coordinates);
              if ( len ) coordinates[len++] = ' ';
              cdf_inq_varname(fileID, ncxvarID, coordinates+len);
            }
          if ( ncyvarID != CDI_UNDEFID )
            {
              size_t len = strlen(coordinates);
              if ( len ) coordinates[len++] = ' ';
              cdf_inq_varname(fileID, ncyvarID, coordinates+len);
            }
632
        }
633
      else
634
        {
635
636
637
638
639
640
641
642
643
644
645
646
          if ( ncyvarID != CDI_UNDEFID )
            {
              size_t len = strlen(coordinates);
              if ( len ) coordinates[len++] = ' ';
              cdf_inq_varname(fileID, ncyvarID, coordinates+len);
            }
          if ( ncxvarID != CDI_UNDEFID )
            {
              size_t len = strlen(coordinates);
              if ( len ) coordinates[len++] = ' ';
              cdf_inq_varname(fileID, ncxvarID, coordinates+len);
            }
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
        }

      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 )
    {
      int gridTruncation = gridInqTrunc(gridID);
      axis[iax++] = '-';
      axis[iax++] = '-';
      cdf_put_att_text(fileID, ncvarid, "axis", iax, axis);
      cdf_put_att_int(fileID, ncvarid, "truncation", NC_INT, 1, &gridTruncation);
    }
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
  else if ( gridtype == GRID_CHARXY )
    {
      if ( gridInqXIsc(gridID) )
        {
          int ncxvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_X];
          size_t len = strlen(coordinates);
          if ( len ) coordinates[len++] = ' ';
          cdf_inq_varname(fileID, ncxvarID, coordinates+len);
        }
      else if ( gridInqYIsc(gridID) )
        {
          int ncyvarID = streamptr->ncgrid[gridindex].ncIDs[CDF_VARID_Y];
          size_t len = strlen(coordinates);
          if ( len ) coordinates[len++] = ' ';
          cdf_inq_varname(fileID, ncyvarID, coordinates+len);
        }
    }
689
690
691
692
693
694
695
696

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

  /*  if ( xtype == NC_BYTE || xtype == NC_SHORT || xtype == NC_INT ) */
    {
      int astype = NC_DOUBLE;

697
698
699
700
      double addoffset   = vlistInqVarAddoffset(vlistID, varID);
      double scalefactor = vlistInqVarScalefactor(vlistID, varID);
      bool laddoffset   = IS_NOT_EQUAL(addoffset, 0);
      bool lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716

      if ( laddoffset || lscalefactor )
        {
          if ( IS_EQUAL(addoffset,   (double) ((float) addoffset)) &&
               IS_EQUAL(scalefactor, (double) ((float) scalefactor)) )
            {
              astype = NC_FLOAT;
            }

          if ( xtype == (int) NC_FLOAT ) astype = NC_FLOAT;

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

717
  if ( dtype == CDI_DATATYPE_UINT8 && xtype == NC_BYTE )
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
    {
      int validrange[2] = {0, 255};
      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);

  if ( zid == -1 )
    {
      if ( zaxisInqType(zaxisID) == ZAXIS_CLOUD_BASE          ||
           zaxisInqType(zaxisID) == ZAXIS_CLOUD_TOP           ||
           zaxisInqType(zaxisID) == ZAXIS_ISOTHERM_ZERO       ||
           zaxisInqType(zaxisID) == ZAXIS_TOA                 ||
           zaxisInqType(zaxisID) == ZAXIS_SEA_BOTTOM          ||
           zaxisInqType(zaxisID) == ZAXIS_LAKE_BOTTOM         ||
           zaxisInqType(zaxisID) == ZAXIS_SEDIMENT_BOTTOM     ||
           zaxisInqType(zaxisID) == ZAXIS_SEDIMENT_BOTTOM_TA  ||
           zaxisInqType(zaxisID) == ZAXIS_SEDIMENT_BOTTOM_TW  ||
           zaxisInqType(zaxisID) == ZAXIS_MIX_LAYER           ||
           zaxisInqType(zaxisID) == ZAXIS_ATMOSPHERE )
        {
          zaxisInqName(zaxisID, varname);
          cdf_put_att_text(fileID, ncvarid, "level_type", strlen(varname), varname);
        }
    }

748
749
750
751
752
753
754
  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);
755
756

  /* Attributes */
757
  cdfDefineAttributes(vlistID, varID, fileID, ncvarid);
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793

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

  return ncvarid;
}


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

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

      int nvars =  streamptr->nvars;
      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
  int gridindex = nc_grid_index(streamptr, gridID);
795
  int lonID = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X];
796
  int latID = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y];
797
798
799
800
801
802

  double xlon = gridInqXval(gridID, 0);
  double xlat = gridInqYval(gridID, 0);
  int tsID = streamptr->curTsID;
  size_t index = (size_t)tsID;

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

static
809
void cdf_write_var_data(int fileID, int vlistID, int varID, int ncvarid, int dtype, size_t nvals, size_t xsize, size_t ysize,
810
                        bool swapxy, size_t *start, size_t *count, int memtype, const void *data, size_t nmiss)
811
812
813
814
815
816
817
818
{
  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;

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

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

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

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

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

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

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

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

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

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

953
  if ( swapxy ) // implemented only for cdf_write_var_slice()
954
    {
955
      size_t gridsize = xsize*ysize;
956
957
      if ( memtype == MEMTYPE_FLOAT )
        {
958
          sdata_sp = (float *) Malloc(gridsize*sizeof(float));
959
960
961
962
963
964
965
          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
        {
966
          sdata_dp = (double *) Malloc(gridsize*sizeof (double));
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
          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);
}


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

990
991
992
993
  size_t xsize = 0, ysize = 0;
  size_t size;
  size_t start[5];
  size_t count[5];
994
  bool swapxy = false;
Thomas Jahns's avatar
Thomas Jahns committed
995
  size_t ndims = 0;
996
997
998
999
1000
1001
1002
1003
1004
1005
1006

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

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

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

  int ncvarid = cdfDefVar(streamptr, varID);

1007
1008
1009
  int gridID   = vlistInqVarGrid(vlistID, varID);
  int zaxisID  = vlistInqVarZaxis(vlistID, varID);
  int timetype = vlistInqVarTimetype(vlistID, varID);
1010

1011
  int xid = CDI_UNDEFID, yid = CDI_UNDEFID;
1012
1013
1014
1015
1016
1017
  if ( gridInqType(gridID) == GRID_TRAJECTORY )
    {
      cdfWriteGridTraj(streamptr, gridID);
    }
  else
    {
1018
      int gridindex = nc_grid_index(streamptr, gridID);
1019
      xid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X];
1020
      yid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y];
1021
1022
1023
1024
1025
    }

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

1026
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
1027
1028
1029
1030
1031
    {
      start[ndims] = (size_t)ntsteps - 1;
      count[ndims] = 1;
      ndims++;
    }
1032

1033
  if ( zid != CDI_UNDEFID )
1034
1035
1036
1037
1038
    {
      start[ndims] = 0;
      count[ndims] = (size_t)zaxisInqSize(zaxisID);
      ndims++;
    }
1039

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

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

  if ( CDI_Debug )
Thomas Jahns's avatar
Thomas Jahns committed
1059
    for (size_t idim = 0; idim < ndims; idim++)
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069
1070
1071
      Message("dim = %d  start = %d  count = %d", idim, start[idim], count[idim]);

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

  int dtype = vlistInqVarDatatype(vlistID, varID);

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

1072
  size_t nvals = gridInqSize(gridID) * (size_t)(zaxisInqSize(zaxisID));
1073
1074
1075
1076
1077
1078

  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,
1079
                         const int rect[][2], const void *data, size_t nmiss)
1080
{
1081
1082
  if ( streamptr->accessmode == 0 ) cdfEndDef(streamptr);

1083
  int xid = CDI_UNDEFID, yid = CDI_UNDEFID;
1084
1085
1086
  size_t xsize = 0, ysize = 0;
  size_t start[5];
  size_t count[5];
1087
  bool swapxy = false;
Thomas Jahns's avatar
Thomas Jahns committed
1088
  size_t ndims = 0;
1089
1090
1091
1092
1093
1094
1095
1096
1097
  int streamID = streamptr->self;

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

  int vlistID = streamInqVlist(streamID);
  int fileID  = streamInqFileID(streamID);

  long ntsteps = streamptr->ntsteps;
1098
  if ( CDI_Debug ) Message("ntsteps = %ld", ntsteps);
1099
1100
1101

  int ncvarid = cdfDefVar(streamptr, varID);

1102
1103
1104
  int gridID   = vlistInqVarGrid(vlistID, varID);
  int zaxisID  = vlistInqVarZaxis(vlistID, varID);
  int timetype = vlistInqVarTimetype(vlistID, varID);
1105
1106
1107
1108
1109
1110
1111

  if ( gridInqType(gridID) == GRID_TRAJECTORY )
    {
      cdfWriteGridTraj(streamptr, gridID);
    }
  else
    {
1112
      int gridindex = nc_grid_index(streamptr, gridID);
1113
      xid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X];
1114
      yid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y];
1115
1116
1117
1118
1119
    }

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

1120
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
1121
1122
1123
1124
1125
    {
      start[ndims] = (size_t)ntsteps - 1;
      count[ndims] = 1;
      ndims++;
    }
1126
  if ( zid != CDI_UNDEFID )
1127
1128
1129
1130
1131
1132
1133
1134
    {
      int size = zaxisInqSize(zaxisID);
      xassert(rect[2][0] >= 0 && rect[2][0] <= rect[2][1]
              && rect[2][1] <= size);
      start[ndims] = (size_t)rect[2][0];
      count[ndims] = (size_t)rect[2][1] - (size_t)rect[2][0] + 1;
      ndims++;
    }
1135
  if ( yid != CDI_UNDEFID )
1136
1137
1138
1139
1140
1141
1142
1143
1144
    {
      size_t size;
      cdf_inq_dimlen(fileID, yid, &size);
      xassert(rect[1][0] >= 0 && rect[1][0] <= rect[1][1]
              && (size_t)rect[1][1] <= size);
      start[ndims] = (size_t)rect[1][0];
      count[ndims] = (size_t)rect[1][1] - (size_t)rect[1][0] + 1;
      ndims++;
    }
1145
  if ( xid != CDI_UNDEFID )
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
    {
      size_t size;
      cdf_inq_dimlen(fileID, xid, &size);
      xassert(rect[0][0] >= 0 && rect[0][0] <= rect[0][1]
              && (size_t)rect[0][1] <= size);
      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
1157
    for (size_t idim = 0; idim < ndims; idim++)
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
      Message("dim = %d  start = %d  count = %d", idim, start[idim], count[idim]);

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

  int dtype = vlistInqVarDatatype(vlistID, varID);

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

1170
  size_t nvals = gridInqSize(gridID) * (size_t)(zaxisInqSize(zaxisID));
1171
1172
1173
1174
1175
1176

  cdf_write_var_data(fileID, vlistID, varID, ncvarid, dtype, nvals,
                     xsize, ysize, swapxy, start, count, memtype, data, nmiss);
}


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

1181
1182
1183
1184
  size_t xsize = 0, ysize = 0;
  size_t start[5];
  size_t count[5];
  int dimorder[3];
1185
  int xid = CDI_UNDEFID, yid = CDI_UNDEFID;
1186
1187
1188
1189
1190
1191
1192
1193
1194
1195
1196

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

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

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

  int ncvarid = cdfDefVar(streamptr, varID);

1197
1198
1199
  int gridID   = vlistInqVarGrid(vlistID, varID);
  int zaxisID  = vlistInqVarZaxis(vlistID, varID);
  int timetype = vlistInqVarTimetype(vlistID, varID);
1200
1201
1202
1203
1204
1205
1206
1207
  vlistInqVarDimorder(vlistID, varID, &dimorder);

  if ( gridInqType(gridID) == GRID_TRAJECTORY )
    {
      cdfWriteGridTraj(streamptr, gridID);
    }
  else
    {
1208
      int gridindex = nc_grid_index(streamptr, gridID);
1209
      xid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_X];
1210
      yid = streamptr->ncgrid[gridindex].ncIDs[CDF_DIMID_Y];
1211
1212
1213
1214
1215
    }

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

1216
  bool swapxy = (dimorder[2] == 2 || dimorder[0] == 1) && xid != CDI_UNDEFID && yid != CDI_UNDEFID;
1217
1218

  size_t ndims = 0;
1219
  if ( vlistHasTime(vlistID) && timetype != TIME_CONSTANT )
1220
1221
1222
1223
1224
1225
1226
1227
    {
      start[ndims] = (size_t)ntsteps - 1;
      count[ndims] = 1;
      ndims++;
    }

  for ( int id = 0; id < 3; ++id )
    {
1228
      if ( dimorder[id] == 3 && zid != CDI_UNDEFID )
1229
1230
1231
1232
1233
        {
          start[ndims] = (size_t)levelID;
          count[ndims] = 1;
          ndims++;
        }
1234
      else if ( dimorder[id] == 2 && yid != CDI_UNDEFID )
1235
1236
1237
1238
1239
1240
        {
          start[ndims] = 0;
          cdf_inq_dimlen(fileID, yid, &ysize);
          count[ndims] = ysize;
          ndims++;
        }
1241
      else if ( dimorder[id] == 1 && xid != CDI_UNDEFID )
1242
1243
1244