stream_ext.c 16.4 KB
Newer Older
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1
2
#ifdef HAVE_CONFIG_H
#include "config.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
3
4
5
6
7
8
9
#endif

#include "dmemory.h"

#include "error.h"
#include "file.h"
#include "cdi.h"
10
#include "cdi_int.h"
11
#include "stream_scan.h"
12
#include "stream_ext.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
13
14
15
#include "varscan.h"
#include "datetime.h"
#include "extra.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
16
#include "exse.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
17
18


Uwe Schulzweida's avatar
Uwe Schulzweida committed
19
#ifdef HAVE_LIBEXTRA
Uwe Schulzweida's avatar
Uwe Schulzweida committed
20

Uwe Schulzweida's avatar
Uwe Schulzweida committed
21
static
22
int extInqDatatype(int prec, int number)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
23
24
25
{
  int datatype;

26
  if ( number == 2 )
27
    datatype = (prec == EXSE_DOUBLE_PRECISION) ? CDI_DATATYPE_CPX64 : CDI_DATATYPE_CPX32;
28
  else
29
    datatype = (prec == EXSE_DOUBLE_PRECISION) ? CDI_DATATYPE_FLT64 : CDI_DATATYPE_FLT32;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
30

31
  return datatype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
32
33
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
34
static
35
void extDefDatatype(int datatype, int *prec, int *number)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
36
37
{

38
39
40
  if ( datatype != CDI_DATATYPE_FLT32 && datatype != CDI_DATATYPE_FLT64 &&
       datatype != CDI_DATATYPE_CPX32 && datatype != CDI_DATATYPE_CPX64 )
    datatype = CDI_DATATYPE_FLT32;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
41

42
  *number = (datatype == CDI_DATATYPE_CPX32 || datatype == CDI_DATATYPE_CPX64) ? 2 : 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
43

44
  *prec = (datatype == CDI_DATATYPE_FLT64 || datatype == CDI_DATATYPE_CPX64) ? EXSE_DOUBLE_PRECISION : EXSE_SINGLE_PRECISION;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
45
46
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
47
/* not used
48
int extInqRecord(stream_t *streamptr, int *varID, int *levelID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
49
50
51
{
  int status;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
52
  int icode, ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
53
54
55
  int zaxisID = -1;
  int header[4];
  int vlistID;
56
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
57

58
59
  vlistID = streamptr->vlistID;
  fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
60
61
62
63
64

  *varID   = -1;
  *levelID = -1;

  status = extRead(fileID, extp);
65
  if ( status != 0 ) return 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
66
67
68

  extInqHeader(extp, header);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
69
  icode  = header[1];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
70
71
  ilevel = header[2];

Uwe Schulzweida's avatar
Uwe Schulzweida committed
72
  *varID = vlistInqVarID(vlistID, icode);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
73

74
  if ( *varID == CDI_UNDEFID ) Error("Code %d undefined", icode);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
75
76
77
78

  zaxisID = vlistInqVarZaxis(vlistID, *varID);

  *levelID = zaxisInqLevelID(zaxisID, (double) ilevel);
79

80
  return 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
81
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
82
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
83

84
void extReadRecord(stream_t *streamptr, double *data, size_t *nmiss)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
85
{
86
87
88
89
90
91
92
  int vlistID = streamptr->vlistID;
  int fileID  = streamptr->fileID;
  int tsID    = streamptr->curTsID;
  int vrecID  = streamptr->tsteps[tsID].curRecID;
  int recID   = streamptr->tsteps[tsID].recIDs[vrecID];
  int varID   = streamptr->tsteps[tsID].records[recID].varID;
  off_t recpos = streamptr->tsteps[tsID].records[recID].position;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
93
94
95

  fileSetPos(fileID, recpos, SEEK_SET);

96
  void *extp = streamptr->record->exsep;
97
98
  int status = extRead(fileID, extp);
  if ( status != 0 ) Error("Failed to read EXTRA record");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
99

100
  int header[4];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
101
102
103
  extInqHeader(extp, header);
  extInqDataDP(extp, data);

104
105
  double missval = vlistInqVarMissval(vlistID, varID);
  int gridID  = vlistInqVarGrid(vlistID, varID);
106
  size_t size = gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
107

Uwe Schulzweida's avatar
Uwe Schulzweida committed
108
  streamptr->numvals += size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
109
110

  *nmiss = 0;
111
112
  if ( vlistInqVarNumber(vlistID, varID) == CDI_REAL )
    {
113
      for ( size_t i = 0; i < size; i++ )
114
115
116
117
118
119
120
121
	if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
	  {
	    data[i] = missval;
	    (*nmiss)++;
	  }
    }
  else
    {
122
      for ( size_t i = 0; i < 2*size; i+=2 )
123
124
125
126
127
128
	if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
	  {
	    data[i] = missval;
	    (*nmiss)++;
	  }
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
129
130
131
}


132
void extCopyRecord(stream_t *streamptr2, stream_t *streamptr1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
133
{
134
  streamFCopyRecord(streamptr2, streamptr1, "EXTRA");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
135
136
137
}


138
void extDefRecord(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
139
{
140
141
  Record *record = streamptr->record;

142
  int pdis, pcat, pnum;
143
  cdiDecodeParam(record->param, &pnum, &pcat, &pdis);
144
145

  int header[4];
146
  header[0] = record->date;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
147
  header[1] = pnum;
148
149
  header[2] = record->level;
  int gridID = record->gridID;
150
  cdi_check_gridsize_int_limit("EXTRA", gridInqSize(gridID));
151
  header[3] = (int) gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
152

153
154
  extrec_t *extp = (extrec_t*) record->exsep;
  extDefDatatype(record->prec, &extp->prec, &extp->number);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
155
156
157
158
  extDefHeader(extp, header);
}


159
void extWriteRecord(stream_t *streamptr, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
160
{
161
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
162
  extDefDataDP(extp, data);
163
  extWrite(streamptr->fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
164
165
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
166
static
167
void extAddRecord(stream_t *streamptr, int param, int level, size_t xysize,
168
		  size_t recsize, off_t position, int prec, int number)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
169
{
170
171
172
  const int vlistID = streamptr->vlistID;
  const int tsID    = streamptr->curTsID;
  const int recID   = recordNewEntry(streamptr, tsID);
173
  record_t *record = &streamptr->tsteps[tsID].records[recID];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
174

175
176
177
178
  record->size     = recsize;
  record->position = position;
  record->param    = param;
  record->ilevel   = level;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
179

180
181
  grid_t *grid = (grid_t *)Malloc(sizeof (*grid));
  grid_init(grid);
182
  cdiGridTypeInit(grid, GRID_GENERIC, xysize);
183
184
  grid->x.size = xysize;
  grid->y.size = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
185
  struct addIfNewRes gridAdded = cdiVlistAddGridIfNew(vlistID, grid, 0);
186
  const int gridID = gridAdded.Id;
187
  if (!gridAdded.isNew) Free(grid);
188

189
190
  const int leveltype = ZAXIS_GENERIC;
  const int datatype = extInqDatatype(prec, number);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
191

192
  int varID, levelID = 0;
193
  varAddRecord(recID, param, gridID, leveltype, 0, level, 0, 0, 0,
194
               datatype, &varID, &levelID, TSTEP_INSTANT, 0, 0, -1,
195
               NULL, NULL, NULL, NULL);
196

197
198
  record->varID   = (short)varID;
  record->levelID = (short)levelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
199

Uwe Schulzweida's avatar
Uwe Schulzweida committed
200
201
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
202
203

  if ( CDI_Debug )
204
    Message("varID = %d gridID = %d levelID = %d", varID, gridID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
205
206
}

207
static
208
209
void extScanTimestep1(stream_t *streamptr)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
210
  int header[4];
211
  DateTime datetime0 = { LONG_MIN, LONG_MIN };
212
  off_t recpos = 0;
213
  extrec_t *extp = (extrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
214

Uwe Schulzweida's avatar
Uwe Schulzweida committed
215
  streamptr->curTsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
216

217
  int tsID  = tstepsNewEntry(streamptr);
218
  if ( tsID != 0 ) Error("Internal problem! tstepsNewEntry returns %d", tsID);
219
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
220

221
  const int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
222

223
  int nrecs = 0;
224
  while ( true )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
225
226
    {
      recpos = fileGetPos(fileID);
227
      if ( extRead(fileID, extp) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
228
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
229
	  streamptr->ntsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
230
231
232
	  break;
	}

233
      const size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
234

235
      extInqHeader(extp, header);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
236

237
238
239
240
241
242
      const int vdate   = header[0];
      const int vtime   = 0;
      const int rcode   = header[1];
      const int rlevel  = header[2];
      const int rxysize = header[3];
      const int param   = cdiEncodeParam(rcode, 255, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
243
      DateTime datetime = { .date = vdate, .time = vtime};
Uwe Schulzweida's avatar
Uwe Schulzweida committed
244

Uwe Schulzweida's avatar
Uwe Schulzweida committed
245
246
      if ( nrecs == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
247
248
249
	  datetime0 = datetime;
          taxis->vdate = vdate;
          taxis->vtime = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
250
251
252
	}
      else
	{
253
254
255
256
          record_t *records = streamptr->tsteps[tsID].records;
	  for ( int recID = 0; recID < nrecs; recID++ )
            if ( param == records[recID].param && rlevel == records[recID].ilevel )
              goto tstepScanLoopFinished;
257

258
	  if ( datetimeDiffer(datetime, datetime0) )
259
	    Warning("Inconsistent verification time for code %d level %d", rcode, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
260
261
262
263
264
	}

      nrecs++;

      if ( CDI_Debug )
265
	Message("%4d%8d%4d%8d%8d%6d", nrecs, (int)recpos, rcode, rlevel, vdate, vtime);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
266

267
      extAddRecord(streamptr, param, rlevel, rxysize, recsize, recpos, extp->prec, extp->number);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
268
269
    }

270
  tstepScanLoopFinished:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
271
  streamptr->rtsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
272

273
  cdi_generate_vars(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
274

275
  const int taxisID = taxisCreate(TAXIS_ABSOLUTE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
276
  taxis->type  = TAXIS_ABSOLUTE;
277
278
  taxis->rdate = taxis->vdate;
  taxis->rtime = taxis->vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
279

280
  const int vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
281
282
  vlistDefTaxis(vlistID, taxisID);

283
  vlist_check_contents(vlistID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
284

285
  streamScanResizeRecords1(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
286

287
  streamScanTsFixNtsteps(streamptr, recpos);
288
  streamScanTimeConstAdjust(streamptr, taxis);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
289
290
}

291
static
292
293
int extScanTimestep2(stream_t *streamptr)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
294
295
  int header[4];
  off_t recpos = 0;
296
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
297
298

  streamptr->curTsID = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
299

300
301
  const int fileID  = streamptr->fileID;
  const int vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
302

303
  int tsID = streamptr->rtsteps;
304
  if ( tsID != 1 ) Error("Internal problem! unexpected timestep %d", tsID+1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
305

306
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
307

Uwe Schulzweida's avatar
Uwe Schulzweida committed
308
  fileSetPos(fileID, streamptr->tsteps[tsID].position, SEEK_SET);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
309

310
  cdi_create_records(streamptr, tsID);
311
  record_t *records = streamptr->tsteps[tsID].records;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
312

313
  const int nrecords = streamScanInitRecords2(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
314

315
  for ( int rindex = 0; rindex <= nrecords; rindex++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
316
317
    {
      recpos = fileGetPos(fileID);
318
      if ( extRead(fileID, extp) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
319
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
320
	  streamptr->ntsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
321
322
323
	  break;
	}

324
      const size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
325

326
      extInqHeader(extp, header);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
327

328
329
330
331
332
      const int vdate  = header[0];
      const int vtime  = 0;
      const int rcode  = header[1];
      const int rlevel = header[2];
      const int param  = cdiEncodeParam(rcode, 255, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
333

Uwe Schulzweida's avatar
Uwe Schulzweida committed
334
335
336
337
338
339
340
      if ( rindex == 0 )
	{
	  taxis->type  = TAXIS_ABSOLUTE;
	  taxis->vdate = vdate;
	  taxis->vtime = vtime;
	}

341
      bool nextstep = false;
342
      int recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
343
344
      for ( recID = 0; recID < nrecords; recID++ )
	{
345
	  if ( param == records[recID].param && rlevel == records[recID].ilevel )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
346
	    {
347
	      if ( records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
348
		{
349
		  nextstep = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
350
351
352
		}
	      else
		{
353
		  records[recID].used = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
354
		  streamptr->tsteps[tsID].recIDs[rindex] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
355
356
357
358
359
360
		}
	      break;
	    }
	}
      if ( recID == nrecords )
	{
361
	  Warning("Code %d level %d not found at timestep %d", rcode, rlevel, tsID+1);
362
	  return CDI_EUFSTRUCT;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
363
364
365
366
367
	}

      if ( nextstep ) break;

      if ( CDI_Debug )
368
	Message("%4d%8d%4d%8d%8d%6d", rindex+1, (int)recpos, rcode, rlevel, vdate, vtime);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
369

370
      if ( param != records[recID].param || rlevel != records[recID].ilevel )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
371
	{
372
	  Message("tsID = %d recID = %d param = %3d new %3d  level = %3d new %3d",
373
		  tsID, recID, records[recID].param, param, records[recID].ilevel, rlevel);
374
	  return CDI_EUFSTRUCT;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
375
376
	}

377
378
      records[recID].position = recpos;
      records[recID].size = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
379
380
    }

381
382
  int nrecs = 0;
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
383
    {
384
      if ( ! records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
385
	{
386
          vlistDefVarTimetype(vlistID, records[recID].varID, TIME_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
387
388
389
390
391
392
	}
      else
	{
	  nrecs++;
	}
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
393
  streamptr->tsteps[tsID].nrecs = nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
394

Uwe Schulzweida's avatar
Uwe Schulzweida committed
395
  streamptr->rtsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
396

397
  streamScanTsFixNtsteps(streamptr, recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
398

399
  return 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
400
401
402
}


403
int extInqContents(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
404
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
405
  streamptr->curTsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
406

407
408
  extScanTimestep1(streamptr);

409
  const int status = (streamptr->ntsteps == -1) ? extScanTimestep2(streamptr) : 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
410

411
  fileSetPos(streamptr->fileID, 0, SEEK_SET);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
412

413
  return status;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
414
415
}

416
static
417
long extScanTimestep(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
418
419
420
{
  int header[4];
  off_t recpos = 0;
421
  int nrecs = 0;
422
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
423

424
  int tsID = streamptr->rtsteps;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
425
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
426

Uwe Schulzweida's avatar
Uwe Schulzweida committed
427
  if ( streamptr->tsteps[tsID].recordSize == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
428
    {
429
      cdi_create_records(streamptr, tsID);
430
      record_t *records = streamptr->tsteps[tsID].records;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
431

432
      nrecs = streamScanInitRecords(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
433

434
      const int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
435

Uwe Schulzweida's avatar
Uwe Schulzweida committed
436
      fileSetPos(fileID, streamptr->tsteps[tsID].position, SEEK_SET);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
437

438
      for ( int rindex = 0; rindex <= nrecs; rindex++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
439
440
	{
	  recpos = fileGetPos(fileID);
441
	  if ( extRead(fileID, extp) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
442
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
443
	      streamptr->ntsteps = streamptr->rtsteps + 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
444
445
446
	      break;
	    }

447
	  const size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
448

449
	  extInqHeader(extp, header);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
450

451
452
453
454
455
	  const int vdate  = header[0];
	  const int vtime  = 0;
	  const int rcode  = header[1];
	  const int rlevel = header[2];
	  const int param  = cdiEncodeParam(rcode, 255, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
456

457
458
	  // if ( rindex == nrecs ) break; gcc-4.5 internal compiler error
	  if ( rindex == nrecs ) continue;
459
	  const int recID = streamptr->tsteps[tsID].recIDs[rindex];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
460
461
462
463
464
465
466

	  if ( rindex == 0 )
	    {
	      taxis->type  = TAXIS_ABSOLUTE;
	      taxis->vdate = vdate;
	      taxis->vtime = vtime;
	    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
467

468
          if ( param != records[recID].param || rlevel != records[recID].ilevel )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
469
	    {
470
	      Message("tsID = %d recID = %d param = %3d new %3d  level = %3d new %3d",
471
		      tsID, recID, records[recID].param, param, records[recID].ilevel, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
472
	      Error("Invalid, unsupported or inconsistent record structure!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
473
474
	    }

475
476
	  records[recID].position = recpos;
	  records[recID].size = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
477
478

	  if ( CDI_Debug )
479
	    Message("%4d%8d%4d%8d%8d%6d", rindex, (int)recpos, rcode, rlevel, vdate, vtime);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
480
481
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
482
      streamptr->rtsteps++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
483

Uwe Schulzweida's avatar
Uwe Schulzweida committed
484
      if ( streamptr->ntsteps != streamptr->rtsteps )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
485
	{
486
	  tsID = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
487
	  if ( tsID != streamptr->rtsteps )
488
	    Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
489

490
	  streamptr->tsteps[tsID-1].next   = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
491
	  streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
492
493
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
494
495
      fileSetPos(fileID, streamptr->tsteps[tsID].position, SEEK_SET);
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
496
497
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
498
  if ( nrecs > 0 && nrecs < streamptr->tsteps[tsID].nrecs )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
499
    {
500
      Warning("Incomplete timestep. Stop scanning at timestep %d.", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
501
      streamptr->ntsteps = tsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
502
503
    }

504
  return streamptr->ntsteps;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
505
506
507
}


508
int extInqTimestep(stream_t *streamptr, int tsID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
509
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
510
  if ( tsID == 0 && streamptr->rtsteps == 0 )
511
    Error("Call to cdiInqContents missing!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
512
513

  if ( CDI_Debug )
514
    Message("tsID = %d rtsteps = %d", tsID, streamptr->rtsteps);
515

516
517
  long ntsteps = CDI_UNDEFID;
  while ( ( tsID + 1 ) > streamptr->rtsteps && ntsteps == CDI_UNDEFID )
518
    ntsteps = extScanTimestep(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
519

520
521
  int nrecs = 0;
  if ( !(tsID >= streamptr->ntsteps && streamptr->ntsteps != CDI_UNDEFID) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
522
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
523
524
      streamptr->curTsID = tsID;
      nrecs = streamptr->tsteps[tsID].nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
525
526
    }

527
  return nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
528
529
530
}


531
void extReadVarSliceDP(stream_t *streamptr, int varID, int levID, double *data, size_t *nmiss)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
532
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
533
534
  if ( CDI_Debug ) Message("streamID = %d  varID = %d  levID = %d", streamptr->self, varID, levID);

535
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
536

537
538
  int vlistID  = streamptr->vlistID;
  int fileID   = streamptr->fileID;
539
  /* NOTE: tiles are not supported here! */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
540
  double missval = vlistInqVarMissval(vlistID, varID);
541
  size_t gridsize = gridInqSize(vlistInqVarGrid(vlistID, varID));
542
  int tsid = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
543

544
  off_t currentfilepos = fileGetPos(fileID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
545

546
  /* NOTE: tiles are not supported here! */
547
548
  int recID = streamptr->vars[varID].recordTable[0].recordID[levID];
  off_t recpos = streamptr->tsteps[tsid].records[recID].position;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
549
  fileSetPos(fileID, recpos, SEEK_SET);
550
  extRead(fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
551
  int header[4];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
552
553
  extInqHeader(extp, header);
  extInqDataDP(extp, data);
Thomas Jahns's avatar
Thomas Jahns committed
554

Uwe Schulzweida's avatar
Uwe Schulzweida committed
555
556
557
  fileSetPos(fileID, currentfilepos, SEEK_SET);

  *nmiss = 0;
558
559
  if ( vlistInqVarNumber(vlistID, varID) == CDI_REAL )
    {
560
      for ( size_t i = 0; i < gridsize; i++ )
561
562
563
564
565
566
567
568
	if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
	  {
	    data[i] = missval;
	    (*nmiss)++;
	  }
    }
  else
    {
569
      for ( size_t i = 0; i < 2*gridsize; i+=2 )
570
571
572
573
574
575
	if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
	  {
	    data[i] = missval;
	    (*nmiss)++;
	  }
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
576
577
578
}


579
void extReadVarDP(stream_t *streamptr, int varID, double *data, size_t *nmiss)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
580
581
582
583
{
  if ( CDI_Debug ) Message("streamID = %d  varID = %d", streamptr->self, varID);

  int vlistID = streamptr->vlistID;
584
  size_t gridsize = gridInqSize(vlistInqVarGrid(vlistID, varID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
585
586
587
588
589
590
591
  size_t nlevs    = (size_t) streamptr->vars[varID].recordTable[0].nlevs;

  for ( size_t levID = 0; levID < nlevs; levID++)
    extReadVarSliceDP(streamptr, varID, (int)levID, &data[levID*gridsize], nmiss);
}


Uwe Schulzweida's avatar
Uwe Schulzweida committed
592
void extWriteVarSliceDP(stream_t *streamptr, int varID, int levID, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
593
{
594
595
  if ( CDI_Debug ) Message("streamID = %d  varID = %d  levID = %d", streamptr->self, varID, levID);

596
597
598
  int vlistID  = streamptr->vlistID;
  int fileID   = streamptr->fileID;
  int tsID     = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
599

600
  int pdis, pcat, pnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
601
  cdiDecodeParam(vlistInqVarParam(vlistID, varID), &pnum, &pcat, &pdis);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
602

Uwe Schulzweida's avatar
Uwe Schulzweida committed
603
  int header[4];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
604
  header[0] = streamptr->tsteps[tsID].taxis.vdate;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
605
  header[1] = pnum;
606
  header[2] = (int)(zaxisInqLevel(vlistInqVarZaxis(vlistID, varID), levID));
607
  int gridID = vlistInqVarGrid(vlistID, varID);
608
  cdi_check_gridsize_int_limit("EXTRA", gridInqSize(gridID));
609
  header[3] = (int) gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
610

Uwe Schulzweida's avatar
Uwe Schulzweida committed
611
  extrec_t *extp = (extrec_t*) streamptr->record->exsep;
612
  extDefDatatype(vlistInqVarDatatype(vlistID, varID), &extp->prec, &extp->number);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
613
  extDefHeader(extp, header);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
614

Uwe Schulzweida's avatar
Uwe Schulzweida committed
615
616
  extDefDataDP(extp, data);
  extWrite(fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
617
618
619
}


Uwe Schulzweida's avatar
Uwe Schulzweida committed
620
void extWriteVarDP(stream_t *streamptr, int varID, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
621
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
622
  if ( CDI_Debug ) Message("streamID = %d  varID = %d", streamptr->self, varID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
623

Uwe Schulzweida's avatar
Uwe Schulzweida committed
624
  int vlistID = streamptr->vlistID;
625
  size_t gridsize = gridInqSize(vlistInqVarGrid(vlistID, varID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
626
  size_t nlevs    = (size_t) zaxisInqSize(vlistInqVarZaxis(vlistID, varID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
627

Uwe Schulzweida's avatar
Uwe Schulzweida committed
628
  for ( size_t levID = 0;  levID < nlevs; levID++ )
Thomas Jahns's avatar
Thomas Jahns committed
629
    extWriteVarSliceDP(streamptr, varID, (int)levID, &data[levID*gridsize]);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
630
631
632
}

#endif /* HAVE_LIBEXTRA */
633

634
635
636
637
638
639
640
641
642
/*
 * Local Variables:
 * c-file-style: "Java"
 * c-basic-offset: 2
 * indent-tabs-mode: nil
 * show-trailing-whitespace: t
 * require-trailing-newline: t
 * End:
 */