stream_ieg.c 35.2 KB
Newer Older
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1
2
3
4
#if defined (HAVE_CONFIG_H)
#  include "config.h"
#endif

5
#include <limits.h>
Thomas Jahns's avatar
Thomas Jahns committed
6
#include <stdbool.h>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
7
#include <stdio.h>
8
#include <stdlib.h>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
9
10
11
12
13
14
15
16
17
#include <string.h>
#include <float.h>
#include <math.h>

#include "dmemory.h"

#include "error.h"
#include "file.h"
#include "cdi.h"
18
#include "cdi_int.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
19
20
21
#include "varscan.h"
#include "datetime.h"
#include "ieg.h"
22
#include "stream_ieg.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
23
#include "vlist.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
24
#include "exse.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
25
26
27
28
29
30
31
32


#undef  UNDEFID
#define UNDEFID  CDI_UNDEFID


#if defined (HAVE_LIBIEG)

33
34

typedef struct {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
35
  int param;
36
  int level;
37
} iegcompvar_t;
38
39


40
41
static
int iegInqDatatype(int prec)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
42
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
43
  return (prec == EXSE_DOUBLE_PRECISION) ? DATATYPE_FLT64 : DATATYPE_FLT32;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
44
45
}

46
47
static
int iegDefDatatype(int datatype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
48
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
49
  if ( datatype == DATATYPE_CPX32 || datatype == DATATYPE_CPX64 )
50
    Error("CDI/IEG library does not support complex numbers!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
51

52
53
  if ( datatype != DATATYPE_FLT32 && datatype != DATATYPE_FLT64 )
    datatype = DATATYPE_FLT32;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
54

Uwe Schulzweida's avatar
Uwe Schulzweida committed
55
  return (datatype == DATATYPE_FLT64) ? EXSE_DOUBLE_PRECISION : EXSE_SINGLE_PRECISION;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
56
57
}

58
/* not used
59
int iegInqRecord(stream_t *streamptr, int *varID, int *levelID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
60
61
62
{
  int status;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
63
  int icode, ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
64
65
  int zaxisID = -1;
  int vlistID;
66
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
67

Uwe Schulzweida's avatar
Uwe Schulzweida committed
68
69
  vlistID = streamptr->vlistID;
  fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
70
71
72
73
74

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

  status = iegRead(fileID, iegp);
75
  if ( status != 0 ) return 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
76

Uwe Schulzweida's avatar
Uwe Schulzweida committed
77
  icode  = IEG_P_Parameter(iegp->ipdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
78
79
80
81
82
  if ( IEG_P_LevelType(iegp->ipdb) == IEG_LTYPE_HYBRID_LAYER )
    ilevel = IEG_P_Level1(iegp->ipdb);
  else
    ilevel = IEG_P_Level2(iegp->ipdb);

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

85
  if ( *varID == UNDEFID ) Error("Code %d undefined", icode);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
86
87
88
89

  zaxisID = vlistInqVarZaxis(vlistID, *varID);

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

91
  return 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
92
}
93
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
94

95
void iegReadRecord(stream_t *streamptr, double *data, int *nmiss)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
96
{
97
98
99
100
101
102
103
  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
104
105
106

  fileSetPos(fileID, recpos, SEEK_SET);

107
  void *iegp = streamptr->record->exsep;
108
  int status = iegRead(fileID, iegp);
109
110
  if ( status != 0 )
    Error("Could not read IEG record!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
111
112
113

  iegInqDataDP(iegp, data);

114
115
116
  double missval = vlistInqVarMissval(vlistID, varID);
  int gridID  = vlistInqVarGrid(vlistID, varID);
  int size    = gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
117

Uwe Schulzweida's avatar
Uwe Schulzweida committed
118
  streamptr->numvals += size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
119
120

  *nmiss = 0;
121
  for ( int i = 0; i < size; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
122
123
124
125
126
127
128
    if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
      {
	data[i] = missval;
	(*nmiss)++;
      }
}

129
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
int iegGetZaxisType(int iegleveltype)
{
  int leveltype = 0;

  switch ( iegleveltype )
    {
    case IEG_LTYPE_SURFACE:
      {
	leveltype = ZAXIS_SURFACE;
	break;
      }
    case IEG_LTYPE_99:
    case IEG_LTYPE_ISOBARIC:
      {
	leveltype = ZAXIS_PRESSURE;
	break;
      }
    case IEG_LTYPE_HEIGHT:
      {
	leveltype = ZAXIS_HEIGHT;
	break;
      }
    case IEG_LTYPE_ALTITUDE:
      {
	leveltype = ZAXIS_ALTITUDE;
	break;
      }
    case IEG_LTYPE_HYBRID:
    case IEG_LTYPE_HYBRID_LAYER:
      {
	leveltype = ZAXIS_HYBRID;
	break;
      }
    case IEG_LTYPE_LANDDEPTH:
    case IEG_LTYPE_LANDDEPTH_LAYER:
      {
	leveltype = ZAXIS_DEPTH_BELOW_LAND;
	break;
      }
    case IEG_LTYPE_SEADEPTH:
      {
	leveltype = ZAXIS_DEPTH_BELOW_SEA;
	break;
      }
    default:
      {
	leveltype = ZAXIS_GENERIC;
	break;
      }
    }

181
  return leveltype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
182
183
184
}


Thomas Jahns's avatar
Thomas Jahns committed
185
static void iegDefTime(int *pdb, int date, int time, int taxisID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
186
187
{
  int timetype = -1;
188
  if ( taxisID != -1 ) timetype = taxisInqType(taxisID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
189
190
191

  if ( timetype == TAXIS_ABSOLUTE || timetype == TAXIS_RELATIVE )
    {
192
      int year, month, day, hour, minute, second;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
193
194
      cdiDecodeDate(date, &year, &month, &day);
      cdiDecodeTime(time, &hour, &minute, &second);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211

      IEG_P_Year(pdb)     = year;
      IEG_P_Month(pdb)    = month;
      IEG_P_Day(pdb)      = day;
      IEG_P_Hour(pdb)     = hour;
      IEG_P_Minute(pdb)   = minute;

      pdb[15] = 1;
      pdb[16] = 0;
      pdb[17] = 0;
      pdb[18] = 10;
      pdb[36] = 1;
    }

  pdb[5] = 128;
}

Thomas Jahns's avatar
Thomas Jahns committed
212
213
214
215
216
/* find smallest power of 10 in [1000,10000000] that upon
 * multiplication results in fractional part close to zero for all
 * arguments */
static double
calc_resfac(double xfirst, double xlast, double xinc, double yfirst, double ylast, double yinc)
217
{
Thomas Jahns's avatar
Thomas Jahns committed
218
219
220
221
222
223
224
225
226
227
  double resfac = 1000.0;
  enum {
    nPwrOf10 = 5,
    nMultTests = 6,
  };
  static const double scaleFactors[nPwrOf10]
    = { 1000, 10000, 100000, 1000000, 10000000 };
  double vals[nMultTests] = { xfirst, xlast, xinc, yfirst, ylast, yinc };

  for (size_t j = 0; j < nPwrOf10; ++j )
228
    {
Thomas Jahns's avatar
Thomas Jahns committed
229
230
231
      double scaleBy = scaleFactors[j];
      bool fractionalScale = false;
      for (size_t i = 0; i < nMultTests; ++i )
232
        {
Thomas Jahns's avatar
Thomas Jahns committed
233
234
          fractionalScale = fractionalScale
            || fabs(vals[i]*scaleBy - round(vals[i]*scaleBy)) > FLT_EPSILON;
235
        }
Thomas Jahns's avatar
Thomas Jahns committed
236
      if ( !fractionalScale )
237
        {
Thomas Jahns's avatar
Thomas Jahns committed
238
          resfac = scaleBy;
239
240
241
242
          break;
        }
    }

243
  return resfac;
244
245
}

246
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
247
248
void iegDefGrid(int *gdb, int gridID)
{
249
  int gridtype = gridInqType(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
250
251
252

  if ( gridtype == GRID_GENERIC )
    {
253
254
      int xsize = gridInqXsize(gridID);
      int ysize = gridInqYsize(gridID);
255
256

      if ( (ysize == 32  || ysize == 48 || ysize == 64 ||
257
	    ysize == 96  || ysize == 160) &&
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
	   (xsize == 2*ysize || xsize == 1) )
	{
	  gridtype = GRID_GAUSSIAN;
	  gridChangeType(gridID, gridtype);
	}
      else if ( (xsize == 1 && ysize == 1) || (xsize == 0 && ysize == 0) )
	{
	  gridtype = GRID_LONLAT;
	  gridChangeType(gridID, gridtype);
	}
      else if ( gridInqXvals(gridID, NULL) && gridInqYvals(gridID, NULL) )
	{
	  gridtype = GRID_LONLAT;
	  gridChangeType(gridID, gridtype);
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
273
274
275
276
277
278
    }
  else if ( gridtype == GRID_CURVILINEAR )
    {
      gridtype = GRID_LONLAT;
    }

279
  if ( gridtype == GRID_LONLAT || gridtype == GRID_GAUSSIAN )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
280
    {
281
282
      double xfirst = 0, xlast = 0, xinc = 0;
      double yfirst = 0, ylast = 0, yinc = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
283

284
285
      int nlon = gridInqXsize(gridID),
        nlat = gridInqYsize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
286

287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
      if ( nlon == 0 )
	{
	  nlon = 1;
	}
      else
	{
	  xfirst = gridInqXval(gridID,      0);
	  xlast  = gridInqXval(gridID, nlon-1);
	  xinc   = gridInqXinc(gridID);
	}

      if ( nlat == 0 )
	{
	  nlat = 1;
	}
      else
	{
	  yfirst = gridInqYval(gridID,      0);
	  ylast  = gridInqYval(gridID, nlat-1);
	  yinc   = gridInqYinc(gridID);
	}

309
310
311
      if ( gridtype == GRID_GAUSSIAN )
	IEG_G_GridType(gdb) = 4;
      else if ( gridtype == GRID_LONLAT && gridIsRotated(gridID) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
312
313
314
315
	IEG_G_GridType(gdb) = 10;
      else
	IEG_G_GridType(gdb) = 0;

Thomas Jahns's avatar
Thomas Jahns committed
316
317
      double resfac = calc_resfac(xfirst, xlast, xinc, yfirst, ylast, yinc);
      int iresfac = (int)resfac;
318
319
320
321
      if ( iresfac == 1000 ) iresfac = 0;

      IEG_G_ResFac(gdb)   = iresfac;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
322
323
      IEG_G_NumLon(gdb)   = nlon;
      IEG_G_NumLat(gdb)   = nlat;
324
325
326
327
328
      IEG_G_FirstLat(gdb) = (int)lround(yfirst*resfac);
      IEG_G_LastLat(gdb)  = (int)lround(ylast*resfac);
      IEG_G_FirstLon(gdb) = (int)lround(xfirst*resfac);
      IEG_G_LastLon(gdb)  = (int)lround(xlast*resfac);
      IEG_G_LonIncr(gdb)  = (int)lround(xinc*resfac);
329
      if ( fabs(xinc*resfac - IEG_G_LonIncr(gdb)) > FLT_EPSILON )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
330
331
	IEG_G_LonIncr(gdb) = 0;

332
333
334
335
      if ( gridtype == GRID_GAUSSIAN )
	IEG_G_LatIncr(gdb) = nlat/2;
      else
	{
336
	  IEG_G_LatIncr(gdb) = (int)lround(yinc*resfac);
337
	  if ( fabs(yinc*resfac - IEG_G_LatIncr(gdb)) > FLT_EPSILON )
338
	    IEG_G_LatIncr(gdb) = 0;
339
340

	  if ( IEG_G_LatIncr(gdb) < 0 ) IEG_G_LatIncr(gdb) = -IEG_G_LatIncr(gdb);
341
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
342

343
      if ( IEG_G_NumLon(gdb) > 1 && IEG_G_NumLat(gdb) == 1 )
344
345
	if ( IEG_G_LonIncr(gdb) != 0 && IEG_G_LatIncr(gdb) == 0 ) IEG_G_LatIncr(gdb) = IEG_G_LonIncr(gdb);

346
      if ( IEG_G_NumLon(gdb) == 1 && IEG_G_NumLat(gdb) > 1 )
347
348
349
350
351
352
353
	if ( IEG_G_LonIncr(gdb) == 0 && IEG_G_LatIncr(gdb) != 0 ) IEG_G_LonIncr(gdb) = IEG_G_LatIncr(gdb);

      if ( IEG_G_LatIncr(gdb) == 0 || IEG_G_LonIncr(gdb) == 0 )
	IEG_G_ResFlag(gdb) = 0;
      else
	IEG_G_ResFlag(gdb) = 128;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
354
355
      if ( gridIsRotated(gridID) )
	{
356
357
	  IEG_G_LatSP(gdb) = - (int)lround(gridInqYpole(gridID) * resfac);
	  IEG_G_LonSP(gdb) =   (int)lround((gridInqXpole(gridID) + 180) * resfac);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
358
359
360
361
362
363
364
365
366
	  IEG_G_Size(gdb)  = 42;
	}
      else
	{
	  IEG_G_Size(gdb)  = 32;
	}
    }
  else
    {
367
      Error("Unsupported grid type: %s", gridNamePtr(gridtype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
368
369
370
371
372
    }

  IEG_G_ScanFlag(gdb) = 64;
}

373
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
374
375
376
void iegDefLevel(int *pdb, int *gdb, double *vct, int zaxisID, int levelID)
{
  double level;
377
  int ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
378

379
  int leveltype = zaxisInqType(zaxisID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
380
381
382

  if ( leveltype == ZAXIS_GENERIC )
    {
383
      Message("Changed zaxis type from %s to %s",
384
385
	      zaxisNamePtr(leveltype),
	      zaxisNamePtr(ZAXIS_PRESSURE));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
386
387
      leveltype = ZAXIS_PRESSURE;
      zaxisChangeType(zaxisID, leveltype);
388
      zaxisDefUnits(zaxisID, "Pa");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
389
390
391
392
393
394
395
396
397
398
    }

  /*  IEG_G_NumVCP(gdb) = 0; */

  switch (leveltype)
    {
    case ZAXIS_SURFACE:
      {
	IEG_P_LevelType(pdb) = IEG_LTYPE_SURFACE;
	IEG_P_Level1(pdb)    = 0;
399
	IEG_P_Level2(pdb)    = (int)(zaxisInqLevel(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
400
401
402
403
404
405
406
	break;
      }
    case ZAXIS_HYBRID:
      {
	if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) )
	  {
	    IEG_P_LevelType(pdb) = IEG_LTYPE_HYBRID_LAYER;
407
408
	    IEG_P_Level1(pdb)    = (int)(zaxisInqLbound(zaxisID, levelID));
	    IEG_P_Level2(pdb)    = (int)(zaxisInqUbound(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
409
410
411
412
413
	  }
	else
	  {
	    IEG_P_LevelType(pdb) = IEG_LTYPE_HYBRID;
	    IEG_P_Level1(pdb)    = 0;
414
	    IEG_P_Level2(pdb)    = (int)(zaxisInqLevel(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
415
416
	  }

417
	int vctsize = zaxisInqVctSize(zaxisID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
418
419
	if ( vctsize > 100 )
	  {
420
            static bool vct_warning = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
421
422
423
	    /*	    IEG_G_NumVCP(gdb) = 0; */
	    if ( vct_warning )
	      {
424
		Warning("VCT size of %d is too large (maximum is 100). Set to 0!", vctsize);
425
		vct_warning = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
426
427
428
429
430
	      }
	  }
	else
	  {
	    IEG_G_Size(gdb) += (vctsize*4);
431
432
	    memcpy(vct, zaxisInqVctPtr(zaxisID), (size_t)vctsize/2*sizeof(double));
	    memcpy(vct+50, zaxisInqVctPtr(zaxisID)+vctsize/2, (size_t)vctsize/2*sizeof(double));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
433
434
435
436
437
438
439
440
441
442
	  }
	break;
      }
    case ZAXIS_PRESSURE:
      {
	double dum;
	char units[128];

	level = zaxisInqLevel(zaxisID, levelID);
	if ( level < 0 )
443
	  Warning("pressure level of %f Pa is below 0.", level);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
444
445

	zaxisInqUnits(zaxisID, units);
446
	if ( memcmp(units, "hPa", 3) == 0 || memcmp(units, "mb",2 ) == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
	  level = level*100;

	ilevel = (int) level;
	if ( level < 32768 && (level < 100 || modf(level/100, &dum) > 0) )
	  {
	    IEG_P_LevelType(pdb) = IEG_LTYPE_99;
	    IEG_P_Level1(pdb)    = 0;
	    IEG_P_Level2(pdb)    = ilevel;
	  }
	else
	  {
	    IEG_P_LevelType(pdb) = IEG_LTYPE_ISOBARIC;
	    IEG_P_Level1(pdb)    = 0;
	    IEG_P_Level2(pdb)    = ilevel/100;
	  }
	break;
      }
    case ZAXIS_HEIGHT:
      {
	level = zaxisInqLevel(zaxisID, levelID);

	ilevel = (int) level;
	IEG_P_LevelType(pdb) = IEG_LTYPE_HEIGHT;
	IEG_P_Level1(pdb)    = 0;
	IEG_P_Level2(pdb)    = ilevel;

	break;
      }
    case ZAXIS_ALTITUDE:
      {
	level = zaxisInqLevel(zaxisID, levelID);

	ilevel = (int) level;
	IEG_P_LevelType(pdb) = IEG_LTYPE_ALTITUDE;
	IEG_P_Level1(pdb)    = 0;
	IEG_P_Level2(pdb)    = ilevel;

	break;
      }
    case ZAXIS_DEPTH_BELOW_LAND:
      {
	if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) )
	  {
	    IEG_P_LevelType(pdb) = IEG_LTYPE_LANDDEPTH_LAYER;
491
492
	    IEG_P_Level1(pdb)    = (int)(zaxisInqLbound(zaxisID, levelID));
	    IEG_P_Level2(pdb)    = (int)(zaxisInqUbound(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
	  }
	else
	  {
	    level = zaxisInqLevel(zaxisID, levelID);

	    ilevel = (int) level;
	    IEG_P_LevelType(pdb) = IEG_LTYPE_LANDDEPTH;
	    IEG_P_Level1(pdb)    = 0;
	    IEG_P_Level2(pdb)    = ilevel;
	  }

	break;
      }
    case ZAXIS_DEPTH_BELOW_SEA:
      {
	level = zaxisInqLevel(zaxisID, levelID);

	ilevel = (int) level;
	IEG_P_LevelType(pdb) = IEG_LTYPE_SEADEPTH;
	IEG_P_Level1(pdb)    = 0;
	IEG_P_Level2(pdb)    = ilevel;

	break;
      }
    case ZAXIS_ISENTROPIC:
      {
	level = zaxisInqLevel(zaxisID, levelID);

	ilevel = (int) level;
	IEG_P_LevelType(pdb) = 113;
	IEG_P_Level1(pdb)    = 0;
	IEG_P_Level2(pdb)    = ilevel;

	break;
      }
    default:
      {
530
	Error("Unsupported zaxis type: %s", zaxisNamePtr(leveltype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
531
532
533
534
535
536
	break;
      }
    }
}


537
void iegCopyRecord(stream_t *streamptr2, stream_t *streamptr1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
538
{
539
  streamFCopyRecord(streamptr2, streamptr1, "IEG");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
540
541
542
}


543
void iegDefRecord(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
544
{
545
546
  Record *record = streamptr->record;
  iegrec_t *iegp = (iegrec_t*) record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
547

548
549
  int vlistID = streamptr->vlistID;
  int byteorder = streamptr->byteorder;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
550

551
552
  int varID   = record->varID;
  int levelID = record->levelID;
553
  int tsID    = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
554

555
556
  int gridID  = vlistInqVarGrid(vlistID, varID);
  int zaxisID = vlistInqVarZaxis(vlistID, varID);
557

558
  iegInitMem(iegp);
559
  for ( int i = 0; i < 37; i++ ) iegp->ipdb[i] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
560
561
562

  iegp->byteswap = getByteswap(byteorder);

563
  int param = vlistInqVarParam(vlistID, varID);
564
  int pdis, pcat, pnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
565
  cdiDecodeParam(param, &pnum, &pcat, &pdis);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
566
567
  IEG_P_Parameter(iegp->ipdb) = pnum;
  if ( pdis == 255 ) IEG_P_CodeTable(iegp->ipdb) = pcat;
568
569
  int date = streamptr->tsteps[tsID].taxis.vdate;
  int time = streamptr->tsteps[tsID].taxis.vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
570
571
572
573
574

  iegDefTime(iegp->ipdb, date, time, vlistInqTaxis(vlistID));
  iegDefGrid(iegp->igdb, gridID);
  iegDefLevel(iegp->ipdb, iegp->igdb, iegp->vct, zaxisID, levelID);

575
  int datatype = record->prec;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
576
577
578
579
  iegp->dprec = iegDefDatatype(datatype);
}


580
void iegWriteRecord(stream_t *streamptr, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
581
{
582
583
  Record *record = streamptr->record;
  iegrec_t *iegp = (iegrec_t*) record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
584

585
  int fileID = streamptr->fileID;
586
  int gridID = record->gridID;
587

588
  int gridsize = gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
589

590
591
  double refval = data[0];
  for ( int i = 1; i < gridsize; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
592
593
594
595
596
597
598
599
    if ( data[i] < refval ) refval = data[i];

  iegp->refval = refval;

  iegDefDataDP(iegp, data);
  iegWrite(fileID, iegp);
}

600
static
601
void iegAddRecord(stream_t *streamptr, int param, int *pdb, int *gdb, double *vct,
602
		  size_t recsize, off_t position, int prec)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
603
{
604
605
606
607
  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
608

609
  int level1, level2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
610
611
612
613
614
615
616
617
618
  if ( IEG_P_LevelType(pdb) == IEG_LTYPE_HYBRID_LAYER )
    {
      level1 = IEG_P_Level1(pdb);
      level2 = IEG_P_Level2(pdb);
    }
  else
    {
      level1 = IEG_P_Level2(pdb);
      level2 = 0;
619
      if ( IEG_P_LevelType(pdb) == 100 ) level1 *= 100;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
620
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
621

622
623
624
625
626
627
  record->size     = recsize;
  record->position = position;
  record->param    = param;
  record->ilevel   = level1;
  record->ilevel2  = level2;
  record->ltype    = IEG_P_LevelType(pdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
628

629
630
631
  int gridtype =
   ( IEG_G_GridType(gdb) == 0 || IEG_G_GridType(gdb) == 10 ) ? GRID_LONLAT :
    ( IEG_G_GridType(gdb) == 4 ) ? GRID_GAUSSIAN : GRID_GENERIC;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
632

633
634
  grid_t *grid = (grid_t *)Malloc(sizeof (*grid));
  grid_init(grid);
635
  cdiGridTypeInit(grid, gridtype, IEG_G_NumLon(gdb)*IEG_G_NumLat(gdb));
636
637
638
639
  int xsize = IEG_G_NumLon(gdb);
  int ysize = IEG_G_NumLat(gdb);
  grid->x.size = xsize;
  grid->y.size = ysize;
640
641
642
  grid->x.inc  = 0;
  grid->y.inc  = 0;
  grid->x.flag = 0;
643
644
645
646
647

  int iresfac = IEG_G_ResFac(gdb);
  if ( iresfac == 0 ) iresfac = 1000;
  double resfac = 1./(double) iresfac;

648
649
  /* if ( IEG_G_FirstLon != 0 || IEG_G_LastLon != 0 ) */
  {
650
    if ( xsize > 1 )
651
652
      {
	if ( IEG_G_ResFlag(gdb) && IEG_G_LonIncr(gdb) > 0 )
653
	  grid->x.inc = IEG_G_LonIncr(gdb) * resfac;
654
	else
655
	  grid->x.inc = (IEG_G_LastLon(gdb) - IEG_G_FirstLon(gdb)) * resfac / (xsize - 1);
656
657
658
659

	/* correct xinc if necessary */
	if ( IEG_G_FirstLon(gdb) == 0 && IEG_G_LastLon(gdb) > 354000 )
	  {
660
	    double xinc = 360. / xsize;
661
662
            /* FIXME: why not use grid->x.inc != xinc as condition? */
	    if ( fabs(grid->x.inc-xinc) > 0.0 )
663
	      {
664
665
		grid->x.inc = xinc;
		if ( CDI_Debug ) Message("set xinc to %g", grid->x.inc);
666
667
668
	      }
	  }
      }
669
670
671
    grid->x.first = IEG_G_FirstLon(gdb) * resfac;
    grid->x.last  = IEG_G_LastLon(gdb)  * resfac;
    grid->x.flag  = 2;
672
  }
673
  grid->y.flag = 0;
674
675
  /* if ( IEG_G_FirstLat != 0 || IEG_G_LastLat != 0 ) */
  {
676
    if ( ysize > 1 )
677
678
      {
	if ( IEG_G_ResFlag(gdb) && IEG_G_LatIncr(gdb) > 0 )
679
	  grid->y.inc = IEG_G_LatIncr(gdb) * resfac;
680
	else
681
	  grid->y.inc = (IEG_G_LastLat(gdb) - IEG_G_FirstLat(gdb)) * resfac / (ysize - 1);
682
      }
683
684
685
    grid->y.first = IEG_G_FirstLat(gdb) * resfac;
    grid->y.last  = IEG_G_LastLat(gdb)  * resfac;
    grid->y.flag  = 2;
686
  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
687

688
  grid->isRotated = FALSE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
689
690
  if ( IEG_G_GridType(gdb) == 10 )
    {
691
      grid->isRotated = TRUE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
692
693
694
      grid->rll.ypole = - IEG_G_LatSP(gdb) * resfac;
      grid->rll.xpole =   IEG_G_LonSP(gdb) * resfac - 180;
      grid->rll.angle = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
695
696
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
697
  struct addIfNewRes gridAdded = cdiVlistAddGridIfNew(vlistID, grid, 0);
698
699
  int gridID = gridAdded.Id;
  if (!gridAdded.isNew) Free(grid);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
700

701
  int leveltype = iegGetZaxisType(IEG_P_LevelType(pdb));
702

Uwe Schulzweida's avatar
Uwe Schulzweida committed
703
704
705
  if ( leveltype == ZAXIS_HYBRID )
    {
      double tmpvct[100];
706
      size_t vctsize = (size_t)IEG_G_NumVCP(gdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
707

708
709
      for (size_t i = 0; i < vctsize/2; i++ ) tmpvct[i] = vct[i];
      for (size_t i = 0; i < vctsize/2; i++ ) tmpvct[i+vctsize/2] = vct[i+50];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
710
711
712
713

      varDefVCT(vctsize, tmpvct);
    }

714
  int lbounds = IEG_P_LevelType(pdb) == IEG_LTYPE_HYBRID_LAYER ? 1 : 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
715

716
  int datatype = iegInqDatatype(prec);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
717

718
  int varID;
719
  int levelID = 0;
720
  varAddRecord(recID, param, gridID, leveltype, lbounds, level1, level2, 0, 0,
721
722
	       datatype, &varID, &levelID, TSTEP_INSTANT, 0, 0, -1,
               NULL, NULL, NULL, NULL, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
723

724
725
  record->varID   = (short)varID;
  record->levelID = (short)levelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
726

Uwe Schulzweida's avatar
Uwe Schulzweida committed
727
728
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
729
730

  if ( CDI_Debug )
731
    Message("varID = %d gridID = %d levelID = %d", varID, gridID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
732
733
}

Thomas Jahns's avatar
Thomas Jahns committed
734
#if 0
735
static
736
void iegCmpRecord(stream_t *streamptr, int tsID, int recID, off_t position, int param,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
737
738
739
740
		  int level, int xsize, int ysize)
{
  int varID = 0;
  int levelID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
741

742
  record_t *record  = &streamptr->tsteps[tsID].records[recID];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
743

Uwe Schulzweida's avatar
Uwe Schulzweida committed
744
  if ( param != (*record).param || level != (*record).ilevel )
745
    Error("inconsistent timestep");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
746
747
748
749
750
751

  (*record).position = position;
  /*
  varID   = (*record).varID;
  levelID = (*record).levelID;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
752
  streamptr->vars[varID].level[levelID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
753

Uwe Schulzweida's avatar
Uwe Schulzweida committed
754
755
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
756
757
  */
  if ( CDI_Debug )
758
    Message("varID = %d levelID = %d", varID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
759
}
Thomas Jahns's avatar
Thomas Jahns committed
760
#endif
761

762
763
static
void iegDateTime(int *pdb, int *date, int *time)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
764
{
765
766
767
  int ryear   = IEG_P_Year(pdb);
  int rmonth  = IEG_P_Month(pdb);
  int rday    = IEG_P_Day(pdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
768

769
770
  int rhour   = IEG_P_Hour(pdb);
  int rminute = IEG_P_Minute(pdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
771
772
773

  if ( rminute == -1 ) rminute = 0;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
774
775
  *date = cdiEncodeDate(ryear, rmonth, rday);
  *time = cdiEncodeTime(rhour, rminute, 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
776
777
}

778
static
779
780
void iegScanTimestep1(stream_t *streamptr)
{
781
  DateTime datetime0 = { LONG_MIN, LONG_MIN };
Uwe Schulzweida's avatar
Uwe Schulzweida committed
782
  off_t recpos;
783
  iegcompvar_t compVar, compVar0;
784
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
785

Uwe Schulzweida's avatar
Uwe Schulzweida committed
786
  streamptr->curTsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
787

788
  int tsID = tstepsNewEntry(streamptr);
789
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
790
791

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

794
  int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
795

796
  int nrecs = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
797
798
799
  while ( TRUE )
    {
      recpos = fileGetPos(fileID);
800
      int status = iegRead(fileID, iegp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
801
802
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
803
	  streamptr->ntsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
804
805
	  break;
	}
806
      size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
807

808
809
810
811
      int prec   = iegp->dprec;
      int rcode  = IEG_P_Parameter(iegp->ipdb);
      int tabnum = IEG_P_CodeTable(iegp->ipdb);
      int param  = cdiEncodeParam(rcode, tabnum, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
812

813
      int rlevel = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
814
815
816
817
818
      if ( IEG_P_LevelType(iegp->ipdb) == IEG_LTYPE_HYBRID_LAYER )
	rlevel = IEG_P_Level1(iegp->ipdb);
      else
	rlevel = IEG_P_Level2(iegp->ipdb);

819
820
      if ( IEG_P_LevelType(iegp->ipdb) == 100 ) rlevel *= 100;

821
      int vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
822
823
824
825
826
827
828
829
830
      iegDateTime(iegp->ipdb, &vdate, &vtime);

      if ( nrecs == 0 )
	{
	  datetime0.date = vdate;
	  datetime0.time = vtime;
	}
      else
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
831
	  compVar.param = param;
832
          compVar.level = rlevel;
833
          int recID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
834
835
	  for ( recID = 0; recID < nrecs; recID++ )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
836
	      compVar0.param = streamptr->tsteps[0].records[recID].param;
837
838
	      compVar0.level = streamptr->tsteps[0].records[recID].ilevel;

839
	      if ( memcmp(&compVar0, &compVar, sizeof(iegcompvar_t)) == 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
840
841
	    }
	  if ( recID < nrecs ) break;
842
843
	  DateTime datetime = { .date = vdate, .time = vtime};
	  if ( datetimeCmp(datetime, datetime0) )
844
	    Warning("Inconsistent verification time for param %d level %d", param, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
845
846
847
848
849
	}

      nrecs++;

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

852
      iegAddRecord(streamptr, param, iegp->ipdb, iegp->igdb, iegp->vct, recsize, recpos, prec);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
853
854
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
855
  streamptr->rtsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
856

857
  cdi_generate_vars(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
858

859
  int taxisID = taxisCreate(TAXIS_ABSOLUTE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
860
  taxis->type  = TAXIS_ABSOLUTE;
861
862
  taxis->vdate = (int)datetime0.date;
  taxis->vtime = (int)datetime0.time;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
863

864
  int vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
865
866
  vlistDefTaxis(vlistID, taxisID);

867
  vlist_check_contents(vlistID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
868

869
  int nrecords = streamptr->tsteps[0].nallrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
870
  if ( nrecords < streamptr->tsteps[0].recordSize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
871
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
872
873
      streamptr->tsteps[0].recordSize = nrecords;
      streamptr->tsteps[0].records =
874
	(record_t *) Realloc(streamptr->tsteps[0].records,
875
                             (size_t)nrecords * sizeof (record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
876
877
    }

878
  streamptr->tsteps[0].recIDs = (int *) Malloc((size_t)nrecords * sizeof (int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
879
  streamptr->tsteps[0].nrecs = nrecords;
880
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
881
    streamptr->tsteps[0].recIDs[recID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
882

Uwe Schulzweida's avatar
Uwe Schulzweida committed
883
  if ( streamptr->ntsteps == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
884
    {
885
      tsID = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
886
      if ( tsID != streamptr->rtsteps )
887
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
888

Uwe Schulzweida's avatar
Uwe Schulzweida committed
889
890
      streamptr->tsteps[tsID-1].next   = TRUE;
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
891
892
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
893
  if ( streamptr->ntsteps == 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
894
895
896
    {
      if ( taxis->vdate == 0 && taxis->vtime == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
897
	  streamptr->ntsteps = 0;
898
	  for ( int varID = 0; varID < streamptr->nvars; varID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
899
	    {
900
	      vlistDefVarTsteptype(vlistID, varID, TSTEP_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
901
902
903
904
905
	    }
	}
    }
}

906
static
907
int iegScanTimestep2(stream_t *streamptr)
908
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
909
  off_t recpos = 0;
910
  iegcompvar_t compVar, compVar0;
911
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
912
913

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

915
916
  int vlistID = streamptr->vlistID;
  int fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
917

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

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

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

926
  cdi_create_records(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
927

928
  int nrecords = streamptr->tsteps[0].nallrecs;
929
  streamptr->tsteps[1].recIDs = (int *) Malloc((size_t)nrecords * sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
930
  streamptr->tsteps[1].nrecs = 0;
931
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
932
    streamptr->tsteps[1].recIDs[recID] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
933

934
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
935
    {
936
      streamptr->tsteps[tsID].records[recID].position =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
937
	streamptr->tsteps[0].records[recID].position;
938
      streamptr->tsteps[tsID].records[recID].size     =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
939
	streamptr->tsteps[0].records[recID].size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
940
941
    }

942
  for ( int rindex = 0; rindex <= nrecords; rindex++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
943
944
    {
      recpos = fileGetPos(fileID);
945
      int status = iegRead(fileID, iegp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
946
947
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
948
	  streamptr->ntsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
949
950
	  break;
	}
951
      size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
952

953
954
955
      int rcode  = IEG_P_Parameter(iegp->ipdb);
      int tabnum = IEG_P_CodeTable(iegp->ipdb);
      int param  = cdiEncodeParam(rcode, tabnum, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
956

957
      int rlevel = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
958
959
960
961
962
      if ( IEG_P_LevelType(iegp->ipdb) == IEG_LTYPE_HYBRID_LAYER )
	rlevel = IEG_P_Level1(iegp->ipdb);
      else
	rlevel = IEG_P_Level2(iegp->ipdb);

963
964
      if ( IEG_P_LevelType(iegp->ipdb) == 100 ) rlevel *= 100;

965
      int vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
966
967
968
969
970
971
972
973
974
      iegDateTime(iegp->ipdb, &vdate, &vtime);

      if ( rindex == 0 )
	{
	  taxis->type  = TAXIS_ABSOLUTE;
	  taxis->vdate = vdate;
	  taxis->vtime = vtime;
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
975
      compVar.param = param;
976
      compVar.level = rlevel;
977
978
      bool nextstep = false;
      int recID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
979
980
      for ( recID = 0; recID < nrecords; recID++ )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
981
	  compVar0.param = streamptr->tsteps[tsID].records[recID].param;
982
983
	  compVar0.level = streamptr->tsteps[tsID].records[recID].ilevel;

984
	  if ( memcmp(&compVar0, &compVar, sizeof(iegcompvar_t)) == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
985
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
986
	      if ( streamptr->tsteps[tsID].records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
987
		{
988
		  nextstep = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
989
990
991
		}
	      else
		{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
992
993
		  streamptr->tsteps[tsID].records[recID].used = TRUE;
		  streamptr->tsteps[tsID].recIDs[rindex] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
994
995
996
997
998
999
		}
	      break;
	    }
	}
      if ( recID == nrecords )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1000
1001
	  char paramstr[32];
	  cdiParamToString(param, paramstr, sizeof(paramstr));
1002
	  Warning("param %s level %d not defined at timestep 1", paramstr, rlevel);
1003
	  return CDI_EUFSTRUCT;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1004
1005
1006
1007
1008
	}

      if ( nextstep ) break;

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1013
      compVar0.param = streamptr->tsteps[tsID].records[recID].param;
1014
      compVar0.level = streamptr->tsteps[tsID].records[recID].ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1015

1016
      if ( memcmp(&compVar0</