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

#include "dmemory.h"

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


#if defined (HAVE_LIBEXTRA)

21
typedef struct {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
22
  int param;
23
  int level;
24
} extcompvar_t;
25

Uwe Schulzweida's avatar
Uwe Schulzweida committed
26
static
27
int extInqDatatype(int prec, int number)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
28
29
30
{
  int datatype;

31
  if ( number == 2 )
32
    datatype = (prec == EXSE_DOUBLE_PRECISION) ? CDI_DATATYPE_CPX64 : CDI_DATATYPE_CPX32;
33
  else
34
    datatype = (prec == EXSE_DOUBLE_PRECISION) ? CDI_DATATYPE_FLT64 : CDI_DATATYPE_FLT32;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
35

36
  return datatype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
37
38
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
39
static
40
void extDefDatatype(int datatype, int *prec, int *number)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
41
42
{

43
44
45
  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
46

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

49
  *prec = (datatype == CDI_DATATYPE_FLT64 || datatype == CDI_DATATYPE_CPX64) ? EXSE_DOUBLE_PRECISION : EXSE_SINGLE_PRECISION;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
50
51
}

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

63
64
  vlistID = streamptr->vlistID;
  fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
65
66
67
68
69

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

  status = extRead(fileID, extp);
70
  if ( status != 0 ) return 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
71
72
73

  extInqHeader(extp, header);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
74
  icode  = header[1];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
75
76
  ilevel = header[2];

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

79
  if ( *varID == CDI_UNDEFID ) Error("Code %d undefined", icode);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
80
81
82
83

  zaxisID = vlistInqVarZaxis(vlistID, *varID);

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

85
  return 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
86
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
87
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
88

89
void extReadRecord(stream_t *streamptr, double *data, size_t *nmiss)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
90
{
91
92
93
94
95
96
97
  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
98
99
100

  fileSetPos(fileID, recpos, SEEK_SET);

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

105
  int header[4];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
106
107
108
  extInqHeader(extp, header);
  extInqDataDP(extp, data);

109
110
  double missval = vlistInqVarMissval(vlistID, varID);
  int gridID  = vlistInqVarGrid(vlistID, varID);
111
  size_t size = gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
112

Uwe Schulzweida's avatar
Uwe Schulzweida committed
113
  streamptr->numvals += size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
114
115

  *nmiss = 0;
116
117
  if ( vlistInqVarNumber(vlistID, varID) == CDI_REAL )
    {
118
      for ( size_t i = 0; i < size; i++ )
119
120
121
122
123
124
125
126
	if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
	  {
	    data[i] = missval;
	    (*nmiss)++;
	  }
    }
  else
    {
127
      for ( size_t i = 0; i < 2*size; i+=2 )
128
129
130
131
132
133
	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
134
135
136
}


137
void extCopyRecord(stream_t *streamptr2, stream_t *streamptr1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
138
{
139
  streamFCopyRecord(streamptr2, streamptr1, "EXTRA");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
140
141
142
}


143
void extDefRecord(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
144
{
145
  int pdis, pcat, pnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
146
  cdiDecodeParam(streamptr->record->param, &pnum, &pcat, &pdis);
147
148

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

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


162
void extWriteRecord(stream_t *streamptr, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
163
{
164
  int fileID = streamptr->fileID;
165
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
166
167
168
169
170

  extDefDataDP(extp, data);
  extWrite(fileID, extp);
}

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

180
181
182
183
  record->size     = recsize;
  record->position = position;
  record->param    = param;
  record->ilevel   = level;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
184

185
186
  grid_t *grid = (grid_t *)Malloc(sizeof (*grid));
  grid_init(grid);
187
  cdiGridTypeInit(grid, GRID_GENERIC, xysize);
188
189
  grid->x.size = xysize;
  grid->y.size = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
190
  struct addIfNewRes gridAdded = cdiVlistAddGridIfNew(vlistID, grid, 0);
191
192
  int gridID = gridAdded.Id;
  if (!gridAdded.isNew) Free(grid);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
193
194
195
196
  /*
  if ( level == 0 ) leveltype = ZAXIS_SURFACE;
  else              leveltype = ZAXIS_GENERIC;
  */
197
  int leveltype = ZAXIS_GENERIC;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
198

199
  int varID, levelID = 0;
200
  varAddRecord(recID, param, gridID, leveltype, 0, level, 0, 0, 0,
201
               extInqDatatype(prec, number), &varID, &levelID, TSTEP_INSTANT, 0, 0, -1,
202
               NULL, NULL, NULL, NULL, NULL, NULL);
203
204
  record->varID   = (short)varID;
  record->levelID = (short)levelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
205

Uwe Schulzweida's avatar
Uwe Schulzweida committed
206
207
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
208
209

  if ( CDI_Debug )
210
    Message("varID = %d gridID = %d levelID = %d", varID, gridID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
211
212
}

213
static
214
215
void extScanTimestep1(stream_t *streamptr)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
216
  int header[4];
217
  DateTime datetime0 = { LONG_MIN, LONG_MIN };
Uwe Schulzweida's avatar
Uwe Schulzweida committed
218
  int varID;
219
  off_t recpos = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
220
  long recsize;
221
  int recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
222
  extcompvar_t compVar, compVar0;
223
  extrec_t *extp = (extrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
224

Uwe Schulzweida's avatar
Uwe Schulzweida committed
225
  streamptr->curTsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
226

227
228
  int tsID  = tstepsNewEntry(streamptr);
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
229
230

  if ( tsID != 0 )
231
    Error("Internal problem! tstepsNewEntry returns %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
232

233
  int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
234

235
  int nrecs = 0;
236
  while ( true )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
237
238
    {
      recpos = fileGetPos(fileID);
239
      int status = extRead(fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
240
241
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
242
	  streamptr->ntsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
243
244
245
246
247
248
	  break;
	}
      recsize = fileGetPos(fileID) - recpos;

      extInqHeader(extp, header);

249
250
251
252
253
      int vdate   = header[0];
      int vtime   = 0;
      int rcode   = header[1];
      int rlevel  = header[2];
      int rxysize = header[3];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
254

255
      int param = cdiEncodeParam(rcode, 255, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
256

Uwe Schulzweida's avatar
Uwe Schulzweida committed
257
258
259
260
261
262
263
      if ( nrecs == 0 )
	{
	  datetime0.date = vdate;
	  datetime0.time = vtime;
	}
      else
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
264
	  compVar.param = param;
265
          compVar.level = rlevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
266
267
	  for ( recID = 0; recID < nrecs; recID++ )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
268
	      compVar0.param  = streamptr->tsteps[0].records[recID].param;
269
270
	      compVar0.level = streamptr->tsteps[0].records[recID].ilevel;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
271
	      if ( memcmp(&compVar0, &compVar, sizeof(extcompvar_t)) == 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
272
273
	    }
	  if ( recID < nrecs ) break;
274
275
	  DateTime datetime = { .date = vdate, .time = vtime};
	  if ( datetimeCmp(datetime, datetime0) )
276
	    Warning("Inconsistent verification time for code %d level %d", rcode, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
277
278
279
280
281
	}

      nrecs++;

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

284
      extAddRecord(streamptr, param, rlevel, rxysize, (size_t)recsize, recpos, extp->prec, extp->number);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
285
286
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
287
  streamptr->rtsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
288

289
  cdi_generate_vars(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
290

291
  int taxisID = taxisCreate(TAXIS_ABSOLUTE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
292
  taxis->type  = TAXIS_ABSOLUTE;
293
294
  taxis->vdate = (int)datetime0.date;
  taxis->vtime = (int)datetime0.time;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
295

296
  int vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
297
298
  vlistDefTaxis(vlistID, taxisID);

299
  vlist_check_contents(vlistID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
300

301
  int nrecords = streamptr->tsteps[0].nallrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
302
  if ( nrecords < streamptr->tsteps[0].recordSize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
303
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
304
305
      streamptr->tsteps[0].recordSize = nrecords;
      streamptr->tsteps[0].records =
306
        (record_t *) Realloc(streamptr->tsteps[0].records, (size_t)nrecords * sizeof (record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
307
308
    }

309
  streamptr->tsteps[0].recIDs = (int *) Malloc((size_t)nrecords * sizeof (int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
310
  streamptr->tsteps[0].nrecs = nrecords;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
311
  for ( recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
312
    streamptr->tsteps[0].recIDs[recID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
313

Uwe Schulzweida's avatar
Uwe Schulzweida committed
314
  if ( streamptr->ntsteps == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
315
    {
316
      tsID = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
317
      if ( tsID != streamptr->rtsteps )
318
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
319

320
      streamptr->tsteps[tsID-1].next   = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
321
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
322
323
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
324
  if ( streamptr->ntsteps == 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
325
326
327
    {
      if ( taxis->vdate == 0 && taxis->vtime == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
328
329
	  streamptr->ntsteps = 0;
	  for ( varID = 0; varID < streamptr->nvars; varID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
330
	    {
331
	      vlistDefVarTsteptype(vlistID, varID, TSTEP_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
332
333
334
335
336
	    }
	}
    }
}

337
static
338
339
int extScanTimestep2(stream_t *streamptr)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
340
341
342
  int header[4];
  int varID;
  off_t recpos = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
343
  extcompvar_t compVar, compVar0;
344
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
345
346

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

348
349
  int fileID  = streamptr->fileID;
  int vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
350

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

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

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

359
  cdi_create_records(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
360

361
  int nrecords = streamptr->tsteps[0].nallrecs;
362
  streamptr->tsteps[1].recIDs = (int *) Malloc((size_t)nrecords * sizeof (int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
363
  streamptr->tsteps[1].nrecs = 0;
364
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
365
    streamptr->tsteps[1].recIDs[recID] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
366

367
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
368
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
369
      varID = streamptr->tsteps[0].records[recID].varID;
370
      streamptr->tsteps[tsID].records[recID].position =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
371
	streamptr->tsteps[0].records[recID].position;
372
      streamptr->tsteps[tsID].records[recID].size     =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
373
	streamptr->tsteps[0].records[recID].size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
374
375
    }

376
  for ( int rindex = 0; rindex <= nrecords; rindex++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
377
378
    {
      recpos = fileGetPos(fileID);
379
      int status = extRead(fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
380
381
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
382
	  streamptr->ntsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
383
384
	  break;
	}
385
      size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
386
387
388

      extInqHeader(extp, header);

389
390
391
392
      int vdate  = header[0];
      int vtime  = 0;
      int rcode  = header[1];
      int rlevel = header[2];
393
      // rxysize = header[3];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
394

395
      int param = cdiEncodeParam(rcode, 255, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
396

Uwe Schulzweida's avatar
Uwe Schulzweida committed
397
398
399
400
401
402
403
      if ( rindex == 0 )
	{
	  taxis->type  = TAXIS_ABSOLUTE;
	  taxis->vdate = vdate;
	  taxis->vtime = vtime;
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
404
      compVar.param = param;
405
      compVar.level = rlevel;
406
      bool nextstep = false;
407
      int recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
408
409
      for ( recID = 0; recID < nrecords; recID++ )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
410
	  compVar0.param = streamptr->tsteps[tsID].records[recID].param;
411
412
	  compVar0.level = streamptr->tsteps[tsID].records[recID].ilevel;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
413
	  if ( memcmp(&compVar0, &compVar, sizeof(extcompvar_t)) == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
414
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
415
	      if ( streamptr->tsteps[tsID].records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
416
		{
417
		  nextstep = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
418
419
420
		}
	      else
		{
421
		  streamptr->tsteps[tsID].records[recID].used = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
422
		  streamptr->tsteps[tsID].recIDs[rindex] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
423
424
425
426
427
428
		}
	      break;
	    }
	}
      if ( recID == nrecords )
	{
429
	  Warning("Code %d level %d not found at timestep %d", rcode, rlevel, tsID+1);
430
	  return CDI_EUFSTRUCT;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
431
432
433
434
435
	}

      if ( nextstep ) break;

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
438
      streamptr->tsteps[tsID].records[recID].size = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
439

Uwe Schulzweida's avatar
Uwe Schulzweida committed
440
      compVar0.param  = streamptr->tsteps[tsID].records[recID].param;
441
      compVar0.level = streamptr->tsteps[tsID].records[recID].ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
442

Uwe Schulzweida's avatar
Uwe Schulzweida committed
443
      if ( memcmp(&compVar0, &compVar, sizeof(extcompvar_t)) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
444
	{
445
	  Message("tsID = %d recID = %d param = %3d new %3d  level = %3d new %3d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
446
		  tsID, recID,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
447
		  streamptr->tsteps[tsID].records[recID].param, param,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
448
		  streamptr->tsteps[tsID].records[recID].ilevel, rlevel);
449
	  return CDI_EUFSTRUCT;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
450
451
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
452
      streamptr->tsteps[1].records[recID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
453
454
    }

455
456
  int nrecs = 0;
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
457
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
458
      if ( ! streamptr->tsteps[tsID].records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
459
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
460
	  varID = streamptr->tsteps[tsID].records[recID].varID;
461
          vlistDefVarTsteptype(vlistID, varID, TSTEP_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
462
463
464
465
466
467
	}
      else
	{
	  nrecs++;
	}
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
468
  streamptr->tsteps[tsID].nrecs = nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
469

Uwe Schulzweida's avatar
Uwe Schulzweida committed
470
  streamptr->rtsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
471

Uwe Schulzweida's avatar
Uwe Schulzweida committed
472
  if ( streamptr->ntsteps == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
473
    {
474
      tsID = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
475
      if ( tsID != streamptr->rtsteps )
476
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
477

478
      streamptr->tsteps[tsID-1].next   = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
479
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
480
481
    }

482
  return 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
483
484
485
}


486
int extInqContents(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
487
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
488
  streamptr->curTsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
489

490
491
  extScanTimestep1(streamptr);

492
  int status = 0;
493
  if ( streamptr->ntsteps == -1 ) status = extScanTimestep2(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
494

495
  int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
496
497
  fileSetPos(fileID, 0, SEEK_SET);

498
  return status;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
499
500
}

501
static
502
long extScanTimestep(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
503
504
505
506
{
  int header[4];
  off_t recpos = 0;
  int recID;
507
  int nrecs = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
508
  extcompvar_t compVar, compVar0;
509
  void *extp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
510
  /*
Uwe Schulzweida's avatar
Uwe Schulzweida committed
511
512
  if ( CDI_Debug )
    {
513
      Message("streamID = %d", streamptr->self);
514
515
516
      Message("cts = %d", streamptr->curTsID);
      Message("rts = %d", streamptr->rtsteps);
      Message("nts = %d", streamptr->ntsteps);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
517
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
518
  */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
519

Uwe Schulzweida's avatar
Uwe Schulzweida committed
520
521
  int tsID  = streamptr->rtsteps;
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
522

Uwe Schulzweida's avatar
Uwe Schulzweida committed
523
  if ( streamptr->tsteps[tsID].recordSize == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
524
    {
525
      cdi_create_records(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
526

Uwe Schulzweida's avatar
Uwe Schulzweida committed
527
      nrecs = streamptr->tsteps[1].nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
528

Uwe Schulzweida's avatar
Uwe Schulzweida committed
529
      streamptr->tsteps[tsID].nrecs = nrecs;
530
      streamptr->tsteps[tsID].recIDs = (int *) Malloc((size_t)nrecs * sizeof (int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
531
      for ( recID = 0; recID < nrecs; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
532
	streamptr->tsteps[tsID].recIDs[recID] = streamptr->tsteps[1].recIDs[recID];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
533

534
      int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
535

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

538
      for ( int rindex = 0; rindex <= nrecs; rindex++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
539
540
	{
	  recpos = fileGetPos(fileID);
541
	  int status = extRead(fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
542
543
	  if ( status != 0 )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
544
	      streamptr->ntsteps = streamptr->rtsteps + 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
545
546
	      break;
	    }
547
	  size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
548
549
550

	  extInqHeader(extp, header);

551
552
553
554
	  int vdate  = header[0];
	  int vtime  = 0;
	  int rcode  = header[1];
	  int rlevel = header[2];
555
	  // rxysize = header[3];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
556

557
	  int param = cdiEncodeParam(rcode, 255, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
558

559
560
	  // if ( rindex == nrecs ) break; gcc-4.5 internal compiler error
	  if ( rindex == nrecs ) continue;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
561
	  recID = streamptr->tsteps[tsID].recIDs[rindex];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
562
563
564
565
566
567
568

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
570
571
572
	  compVar.param  = param;
          compVar.level  = rlevel;
	  compVar0.param = streamptr->tsteps[tsID].records[recID].param;
573
	  compVar0.level = streamptr->tsteps[tsID].records[recID].ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
574

Uwe Schulzweida's avatar
Uwe Schulzweida committed
575
	  if ( memcmp(&compVar0, &compVar, sizeof(extcompvar_t)) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
576
	    {
577
	      Message("tsID = %d recID = %d param = %3d new %3d  level = %3d new %3d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
578
		      tsID, recID,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
579
		      streamptr->tsteps[tsID].records[recID].param, param,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
580
		      streamptr->tsteps[tsID].records[recID].ilevel, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
581
	      Error("Invalid, unsupported or inconsistent record structure!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
582
583
	    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
584
585
	  streamptr->tsteps[tsID].records[recID].position = recpos;
	  streamptr->tsteps[tsID].records[recID].size = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
586
587

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
591
      streamptr->rtsteps++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
592

Uwe Schulzweida's avatar
Uwe Schulzweida committed
593
      if ( streamptr->ntsteps != streamptr->rtsteps )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
594
	{
595
	  tsID = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
596
	  if ( tsID != streamptr->rtsteps )
597
	    Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
598

599
	  streamptr->tsteps[tsID-1].next   = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
600
	  streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
601
602
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
603
604
      fileSetPos(fileID, streamptr->tsteps[tsID].position, SEEK_SET);
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
605
606
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
607
  if ( nrecs > 0 && nrecs < streamptr->tsteps[tsID].nrecs )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
608
    {
609
      Warning("Incomplete timestep. Stop scanning at timestep %d.", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
610
      streamptr->ntsteps = tsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
611
612
    }

613
  return streamptr->ntsteps;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
614
615
616
}


617
int extInqTimestep(stream_t *streamptr, int tsID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
618
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
619
  if ( tsID == 0 && streamptr->rtsteps == 0 )
620
    Error("Call to cdiInqContents missing!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
621
622

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

625
626
  long ntsteps = CDI_UNDEFID;
  while ( ( tsID + 1 ) > streamptr->rtsteps && ntsteps == CDI_UNDEFID )
627
    ntsteps = extScanTimestep(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
628

629
630
  int nrecs = 0;
  if ( !(tsID >= streamptr->ntsteps && streamptr->ntsteps != CDI_UNDEFID) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
631
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
632
633
      streamptr->curTsID = tsID;
      nrecs = streamptr->tsteps[tsID].nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
634
635
    }

636
  return nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
637
638
639
}


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

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

646
647
  int vlistID  = streamptr->vlistID;
  int fileID   = streamptr->fileID;
648
  /* NOTE: tiles are not supported here! */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
649
  double missval = vlistInqVarMissval(vlistID, varID);
650
  size_t gridsize = gridInqSize(vlistInqVarGrid(vlistID, varID));
651
  int tsid     = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
652

653
  off_t currentfilepos = fileGetPos(fileID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
654

655
  /* NOTE: tiles are not supported here! */
656
657
  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
658
  fileSetPos(fileID, recpos, SEEK_SET);
659
  extRead(fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
660
  int header[4];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
661
662
  extInqHeader(extp, header);
  extInqDataDP(extp, data);
Thomas Jahns's avatar
Thomas Jahns committed
663

Uwe Schulzweida's avatar
Uwe Schulzweida committed
664
665
666
  fileSetPos(fileID, currentfilepos, SEEK_SET);

  *nmiss = 0;
667
668
  if ( vlistInqVarNumber(vlistID, varID) == CDI_REAL )
    {
669
      for ( size_t i = 0; i < gridsize; i++ )
670
671
672
673
674
675
676
677
	if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
	  {
	    data[i] = missval;
	    (*nmiss)++;
	  }
    }
  else
    {
678
      for ( size_t i = 0; i < 2*gridsize; i+=2 )
679
680
681
682
683
684
	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
685
686
687
}


688
void extReadVarDP(stream_t *streamptr, int varID, double *data, size_t *nmiss)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
689
690
691
692
{
  if ( CDI_Debug ) Message("streamID = %d  varID = %d", streamptr->self, varID);

  int vlistID = streamptr->vlistID;
693
  size_t gridsize = gridInqSize(vlistInqVarGrid(vlistID, varID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
694
695
696
697
698
699
700
  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
701
void extWriteVarSliceDP(stream_t *streamptr, int varID, int levID, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
702
{
703
704
  if ( CDI_Debug ) Message("streamID = %d  varID = %d  levID = %d", streamptr->self, varID, levID);

705
706
707
  int vlistID  = streamptr->vlistID;
  int fileID   = streamptr->fileID;
  int tsID     = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
708

709
  int pdis, pcat, pnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
710
  cdiDecodeParam(vlistInqVarParam(vlistID, varID), &pnum, &pcat, &pdis);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
711

Uwe Schulzweida's avatar
Uwe Schulzweida committed
712
  int header[4];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
713
  header[0] = streamptr->tsteps[tsID].taxis.vdate;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
714
  header[1] = pnum;
715
  header[2] = (int)(zaxisInqLevel(vlistInqVarZaxis(vlistID, varID), levID));
716
  int gridID = vlistInqVarGrid(vlistID, varID);
717
  cdi_check_gridsize_int_limit("EXTRA", gridInqSize(gridID));
718
  header[3] = (int) gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
719

Uwe Schulzweida's avatar
Uwe Schulzweida committed
720
  extrec_t *extp = (extrec_t*) streamptr->record->exsep;
721
  extDefDatatype(vlistInqVarDatatype(vlistID, varID), &extp->prec, &extp->number);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
722
  extDefHeader(extp, header);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
723

Uwe Schulzweida's avatar
Uwe Schulzweida committed
724
725
  extDefDataDP(extp, data);
  extWrite(fileID, extp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
726
727
728
}


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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
733
  int vlistID = streamptr->vlistID;
734
  size_t gridsize = gridInqSize(vlistInqVarGrid(vlistID, varID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
735
  size_t nlevs    = (size_t) zaxisInqSize(vlistInqVarZaxis(vlistID, varID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
736

Uwe Schulzweida's avatar
Uwe Schulzweida committed
737
  for ( size_t levID = 0;  levID < nlevs; levID++ )
Thomas Jahns's avatar
Thomas Jahns committed
738
    extWriteVarSliceDP(streamptr, varID, (int)levID, &data[levID*gridsize]);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
739
740
741
}

#endif /* HAVE_LIBEXTRA */
742

743
744
745
746
747
748
749
750
751
/*
 * Local Variables:
 * c-file-style: "Java"
 * c-basic-offset: 2
 * indent-tabs-mode: nil
 * show-trailing-whitespace: t
 * require-trailing-newline: t
 * End:
 */