stream_ieg.c 34.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
250
251
  int projID = gridInqProj(gridID);
  if ( projID != CDI_UNDEFID && gridInqProjType(projID) == CDI_PROJ_RLL ) gridID = projID;

252
  int gridtype = gridInqType(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
253

254
255
256
257
258
259
  int projtype = CDI_UNDEFID;
  if ( gridtype == GRID_PROJECTION && gridInqProjType(gridID) == CDI_PROJ_RLL ) projtype = CDI_PROJ_RLL;

  int xsize = gridInqXsize(gridID);
  int ysize = gridInqYsize(gridID);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
260
261
  if ( gridtype == GRID_GENERIC )
    {
262
      if ( (ysize == 32  || ysize == 48 || ysize == 64 ||
263
	    ysize == 96  || ysize == 160) &&
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
	   (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
279
280
281
282
283
284
    }
  else if ( gridtype == GRID_CURVILINEAR )
    {
      gridtype = GRID_LONLAT;
    }

285
286
287
  bool lrotated = ((gridtype == GRID_LONLAT && gridIsRotated(gridID)) || projtype == CDI_PROJ_RLL);

  if ( gridtype == GRID_LONLAT || gridtype == GRID_GAUSSIAN || projtype == CDI_PROJ_RLL )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
288
    {
289
290
      double xfirst = 0, xlast = 0, xinc = 0;
      double yfirst = 0, ylast = 0, yinc = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
291

292
      if ( xsize == 0 ) xsize = 1;
293
294
      else
	{
295
296
	  xfirst = gridInqXval(gridID,       0);
	  xlast  = gridInqXval(gridID, xsize-1);
297
298
299
	  xinc   = gridInqXinc(gridID);
	}

300
      if ( ysize == 0 ) ysize = 1;
301
302
      else
	{
303
304
	  yfirst = gridInqYval(gridID,       0);
	  ylast  = gridInqYval(gridID, ysize-1);
305
306
307
	  yinc   = gridInqYinc(gridID);
	}

308
309
      if ( gridtype == GRID_GAUSSIAN )
	IEG_G_GridType(gdb) = 4;
310
      else if ( lrotated )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
311
312
313
314
	IEG_G_GridType(gdb) = 10;
      else
	IEG_G_GridType(gdb) = 0;

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

      IEG_G_ResFac(gdb)   = iresfac;

321
322
      IEG_G_NumLon(gdb)   = xsize;
      IEG_G_NumLat(gdb)   = ysize;
323
324
325
326
327
      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);
328
      if ( fabs(xinc*resfac - IEG_G_LonIncr(gdb)) > FLT_EPSILON )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
329
330
	IEG_G_LonIncr(gdb) = 0;

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

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

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

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

348
      IEG_G_ResFlag(gdb) = (IEG_G_LatIncr(gdb) == 0 || IEG_G_LonIncr(gdb) == 0) ? 0 : 128;
349

350
      if ( lrotated )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
351
	{
352
353
          double xpole = 0, ypole = 0, angle = 0;
          if ( projtype == CDI_PROJ_RLL )
354
            gridInqParamRLL(gridID, &xpole, &ypole, &angle);
355
356
357
358
359
360
361
362
363
          else
            {
              xpole = gridInqXpole(gridID);
              ypole = gridInqYpole(gridID);
              angle = gridInqAngle(gridID);
            }

	  IEG_G_LatSP(gdb) = - (int)lround(ypole * resfac);
	  IEG_G_LonSP(gdb) =   (int)lround((xpole + 180) * resfac);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
364
365
366
367
368
369
370
371
372
	  IEG_G_Size(gdb)  = 42;
	}
      else
	{
	  IEG_G_Size(gdb)  = 32;
	}
    }
  else
    {
373
      Error("Unsupported grid type: %s", gridNamePtr(gridtype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
374
375
376
377
378
    }

  IEG_G_ScanFlag(gdb) = 64;
}

379
380
381
382
383
384
385
386
static
void pdbDefLevel(int *pdb, int leveltype, int level1, int level2)
{
  IEG_P_LevelType(pdb) = leveltype;
  IEG_P_Level1(pdb)    = level1;
  IEG_P_Level2(pdb)    = level2;
}

387
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
388
389
390
void iegDefLevel(int *pdb, int *gdb, double *vct, int zaxisID, int levelID)
{
  double level;
391
  int ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
392

393
  int leveltype = zaxisInqType(zaxisID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
394
395
396

  if ( leveltype == ZAXIS_GENERIC )
    {
397
      Message("Changed zaxis type from %s to %s",
398
	      zaxisNamePtr(leveltype), zaxisNamePtr(ZAXIS_PRESSURE));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
399
400
      leveltype = ZAXIS_PRESSURE;
      zaxisChangeType(zaxisID, leveltype);
401
      zaxisDefUnits(zaxisID, "Pa");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
402
403
404
405
406
407
408
409
    }

  /*  IEG_G_NumVCP(gdb) = 0; */

  switch (leveltype)
    {
    case ZAXIS_SURFACE:
      {
410
        pdbDefLevel(pdb, IEG_LTYPE_SURFACE, 0, (int)zaxisInqLevel(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
411
412
413
414
415
	break;
      }
    case ZAXIS_HYBRID:
      {
	if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) )
416
417
          pdbDefLevel(pdb, IEG_LTYPE_HYBRID_LAYER, (int)zaxisInqLbound(zaxisID, levelID),
                      (int)zaxisInqUbound(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
418
	else
419
          pdbDefLevel(pdb, IEG_LTYPE_HYBRID, 0, (int)zaxisInqLevel(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
420

421
	int vctsize = zaxisInqVctSize(zaxisID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
422
423
	if ( vctsize > 100 )
	  {
424
            static bool vct_warning = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
425
426
427
	    /*	    IEG_G_NumVCP(gdb) = 0; */
	    if ( vct_warning )
	      {
428
		Warning("VCT size of %d is too large (maximum is 100). Set to 0!", vctsize);
429
		vct_warning = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
430
431
432
433
434
	      }
	  }
	else
	  {
	    IEG_G_Size(gdb) += (vctsize*4);
435
436
	    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
437
438
439
440
441
442
443
444
445
	  }
	break;
      }
    case ZAXIS_PRESSURE:
      {
	double dum;
	char units[128];

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

	zaxisInqUnits(zaxisID, units);
449
	if ( memcmp(units, "hPa", 3) == 0 || memcmp(units, "mb",2 ) == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
450
451
452
453
	  level = level*100;

	ilevel = (int) level;
	if ( level < 32768 && (level < 100 || modf(level/100, &dum) > 0) )
454
          pdbDefLevel(pdb, IEG_LTYPE_99, 0, ilevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
455
	else
456
457
458
          pdbDefLevel(pdb, IEG_LTYPE_ISOBARIC, 0, ilevel/100);

        break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
459
460
461
462
      }
    case ZAXIS_HEIGHT:
      {
	level = zaxisInqLevel(zaxisID, levelID);
463
        pdbDefLevel(pdb, IEG_LTYPE_HEIGHT, 0, (int)level);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
464
465
466
467
468
	break;
      }
    case ZAXIS_ALTITUDE:
      {
	level = zaxisInqLevel(zaxisID, levelID);
469
        pdbDefLevel(pdb, IEG_LTYPE_ALTITUDE, 0, (int)level);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
470
471
472
473
474
	break;
      }
    case ZAXIS_DEPTH_BELOW_LAND:
      {
	if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) )
475
          pdbDefLevel(pdb, IEG_LTYPE_LANDDEPTH_LAYER, (int)zaxisInqLbound(zaxisID, levelID), (int)zaxisInqUbound(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
476
	else
477
          pdbDefLevel(pdb, IEG_LTYPE_LANDDEPTH, 0, (int)zaxisInqLevel(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
478
479
480
481
482
483

	break;
      }
    case ZAXIS_DEPTH_BELOW_SEA:
      {
	level = zaxisInqLevel(zaxisID, levelID);
484
        pdbDefLevel(pdb, IEG_LTYPE_SEADEPTH, 0, (int)level);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
485
486
487
488
489
	break;
      }
    case ZAXIS_ISENTROPIC:
      {
	level = zaxisInqLevel(zaxisID, levelID);
490
        pdbDefLevel(pdb, 113, 0, (int)level);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
491
492
493
494
	break;
      }
    default:
      {
495
	Error("Unsupported zaxis type: %s", zaxisNamePtr(leveltype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
496
497
498
499
500
501
	break;
      }
    }
}


502
void iegCopyRecord(stream_t *streamptr2, stream_t *streamptr1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
503
{
504
  streamFCopyRecord(streamptr2, streamptr1, "IEG");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
505
506
507
}


508
void iegDefRecord(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
509
{
510
511
  Record *record = streamptr->record;
  iegrec_t *iegp = (iegrec_t*) record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
512

513
514
  int vlistID = streamptr->vlistID;
  int byteorder = streamptr->byteorder;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
515

516
517
  int varID   = record->varID;
  int levelID = record->levelID;
518
  int tsID    = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
519

520
521
  int gridID  = vlistInqVarGrid(vlistID, varID);
  int zaxisID = vlistInqVarZaxis(vlistID, varID);
522

523
  iegInitMem(iegp);
524
  for ( int i = 0; i < 37; i++ ) iegp->ipdb[i] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
525
526
527

  iegp->byteswap = getByteswap(byteorder);

528
  int param = vlistInqVarParam(vlistID, varID);
529
  int pdis, pcat, pnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
530
  cdiDecodeParam(param, &pnum, &pcat, &pdis);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
531
532
  IEG_P_Parameter(iegp->ipdb) = pnum;
  if ( pdis == 255 ) IEG_P_CodeTable(iegp->ipdb) = pcat;
533

534
535
  int date = streamptr->tsteps[tsID].taxis.vdate;
  int time = streamptr->tsteps[tsID].taxis.vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
536
537
538
539
  iegDefTime(iegp->ipdb, date, time, vlistInqTaxis(vlistID));
  iegDefGrid(iegp->igdb, gridID);
  iegDefLevel(iegp->ipdb, iegp->igdb, iegp->vct, zaxisID, levelID);

540
  iegp->dprec = iegDefDatatype(record->prec);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
541
542
543
}


544
void iegWriteRecord(stream_t *streamptr, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
545
{
546
547
  Record *record = streamptr->record;
  iegrec_t *iegp = (iegrec_t*) record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
548

549
  int fileID = streamptr->fileID;
550
  int gridsize = gridInqSize(record->gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
551

552
553
  double refval = data[0];
  for ( int i = 1; i < gridsize; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
554
555
556
557
558
559
560
561
    if ( data[i] < refval ) refval = data[i];

  iegp->refval = refval;

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

562
static
563
void iegAddRecord(stream_t *streamptr, int param, int *pdb, int *gdb, double *vct,
564
		  size_t recsize, off_t position, int prec)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
565
{
566
567
568
569
  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
570

571
  int level1, level2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
572
573
574
575
576
577
578
579
580
  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;
581
      if ( IEG_P_LevelType(pdb) == 100 ) level1 *= 100;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
582
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
583

584
585
586
587
588
589
  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
590

591
592
593
594
595
596
597
598
#ifdef TEST_PROJECTION
  int gridtype = (IEG_G_GridType(gdb) == 0) ? GRID_LONLAT :
                 (IEG_G_GridType(gdb) == 10) ? GRID_PROJECTION :
                 (IEG_G_GridType(gdb) == 4) ? GRID_GAUSSIAN : GRID_GENERIC;
#else
  int gridtype = (IEG_G_GridType(gdb) == 0 || IEG_G_GridType(gdb) == 10) ? GRID_LONLAT :
                 (IEG_G_GridType(gdb) == 4) ? GRID_GAUSSIAN : GRID_GENERIC;
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
599

600
601
  grid_t *grid = (grid_t *)Malloc(sizeof (*grid));
  grid_init(grid);
602
  cdiGridTypeInit(grid, gridtype, IEG_G_NumLon(gdb)*IEG_G_NumLat(gdb));
603
604
605
606
  int xsize = IEG_G_NumLon(gdb);
  int ysize = IEG_G_NumLat(gdb);
  grid->x.size = xsize;
  grid->y.size = ysize;
607
608
609
  grid->x.inc  = 0;
  grid->y.inc  = 0;
  grid->x.flag = 0;
610
611
612
613
614

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

615
616
  /* if ( IEG_G_FirstLon != 0 || IEG_G_LastLon != 0 ) */
  {
617
    if ( xsize > 1 )
618
619
      {
	if ( IEG_G_ResFlag(gdb) && IEG_G_LonIncr(gdb) > 0 )
620
	  grid->x.inc = IEG_G_LonIncr(gdb) * resfac;
621
	else
622
	  grid->x.inc = (IEG_G_LastLon(gdb) - IEG_G_FirstLon(gdb)) * resfac / (xsize - 1);
623
624
625
626

	/* correct xinc if necessary */
	if ( IEG_G_FirstLon(gdb) == 0 && IEG_G_LastLon(gdb) > 354000 )
	  {
627
	    double xinc = 360. / xsize;
628
629
            /* FIXME: why not use grid->x.inc != xinc as condition? */
	    if ( fabs(grid->x.inc-xinc) > 0.0 )
630
	      {
631
632
		grid->x.inc = xinc;
		if ( CDI_Debug ) Message("set xinc to %g", grid->x.inc);
633
634
635
	      }
	  }
      }
636
637
638
    grid->x.first = IEG_G_FirstLon(gdb) * resfac;
    grid->x.last  = IEG_G_LastLon(gdb)  * resfac;
    grid->x.flag  = 2;
639
  }
640
  grid->y.flag = 0;
641
642
  /* if ( IEG_G_FirstLat != 0 || IEG_G_LastLat != 0 ) */
  {
643
    if ( ysize > 1 )
644
645
      {
	if ( IEG_G_ResFlag(gdb) && IEG_G_LatIncr(gdb) > 0 )
646
	  grid->y.inc = IEG_G_LatIncr(gdb) * resfac;
647
	else
648
	  grid->y.inc = (IEG_G_LastLat(gdb) - IEG_G_FirstLat(gdb)) * resfac / (ysize - 1);
649
      }
650
651
652
    grid->y.first = IEG_G_FirstLat(gdb) * resfac;
    grid->y.last  = IEG_G_LastLat(gdb)  * resfac;
    grid->y.flag  = 2;
653
  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
654

655
  grid->isRotated = FALSE;
656
  double xpole = 0, ypole = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
657
658
  if ( IEG_G_GridType(gdb) == 10 )
    {
659
660
661
662
663
664
665
666
667
      xpole =   IEG_G_LonSP(gdb) * resfac - 180;
      ypole = - IEG_G_LatSP(gdb) * resfac;
      if ( gridtype == GRID_LONLAT )
        {
          grid->isRotated = TRUE;
          grid->rll.xpole = xpole;
          grid->rll.ypole = ypole;
          grid->rll.angle = 0;
        }
668
669
      else
        grid->projtype = CDI_PROJ_RLL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
670
671
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
672
  struct addIfNewRes gridAdded = cdiVlistAddGridIfNew(vlistID, grid, 0);
673
  int gridID = gridAdded.Id;
674
675
  if ( !gridAdded.isNew ) Free(grid);
  else if ( gridtype == GRID_PROJECTION ) gridDefProjParamRLL(gridID, xpole, ypole, 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
676

677
  int leveltype = iegGetZaxisType(IEG_P_LevelType(pdb));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
678
679
680
  if ( leveltype == ZAXIS_HYBRID )
    {
      double tmpvct[100];
681
      size_t vctsize = (size_t)IEG_G_NumVCP(gdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
682

683
      for ( size_t i = 0; i < vctsize/2; i++ ) tmpvct[i] = vct[i];
684
      for ( size_t i = 0; i < vctsize/2; i++ ) tmpvct[i+vctsize/2] = vct[i+50];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
685
686
687
688

      varDefVCT(vctsize, tmpvct);
    }

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

691
  int datatype = iegInqDatatype(prec);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
692

693
  int varID;
694
  int levelID = 0;
695
  varAddRecord(recID, param, gridID, leveltype, lbounds, level1, level2, 0, 0,
696
697
	       datatype, &varID, &levelID, TSTEP_INSTANT, 0, 0, -1,
               NULL, NULL, NULL, NULL, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
698

699
700
  record->varID   = (short)varID;
  record->levelID = (short)levelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
701

Uwe Schulzweida's avatar
Uwe Schulzweida committed
702
703
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
704
705

  if ( CDI_Debug )
706
    Message("varID = %d gridID = %d levelID = %d", varID, gridID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
707
708
}

Thomas Jahns's avatar
Thomas Jahns committed
709
#if 0
710
static
711
void iegCmpRecord(stream_t *streamptr, int tsID, int recID, off_t position, int param,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
712
713
714
715
		  int level, int xsize, int ysize)
{
  int varID = 0;
  int levelID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
716

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
719
  if ( param != (*record).param || level != (*record).ilevel )
720
    Error("inconsistent timestep");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
721
722
723
724
725
726

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
729
730
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
731
732
  */
  if ( CDI_Debug )
733
    Message("varID = %d levelID = %d", varID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
734
}
Thomas Jahns's avatar
Thomas Jahns committed
735
#endif
736

737
738
static
void iegDateTime(int *pdb, int *date, int *time)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
739
{
740
741
742
  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
743

744
745
  int rhour   = IEG_P_Hour(pdb);
  int rminute = IEG_P_Minute(pdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
746
747
748

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
749
750
  *date = cdiEncodeDate(ryear, rmonth, rday);
  *time = cdiEncodeTime(rhour, rminute, 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
751
752
}

753
static
754
755
void iegScanTimestep1(stream_t *streamptr)
{
756
  DateTime datetime0 = { LONG_MIN, LONG_MIN };
Uwe Schulzweida's avatar
Uwe Schulzweida committed
757
  off_t recpos;
758
  iegcompvar_t compVar, compVar0;
759
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
760

Uwe Schulzweida's avatar
Uwe Schulzweida committed
761
  streamptr->curTsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
762

763
  int tsID = tstepsNewEntry(streamptr);
764
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
765
766

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

769
  int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
770

771
  int nrecs = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
772
773
774
  while ( TRUE )
    {
      recpos = fileGetPos(fileID);
775
      int status = iegRead(fileID, iegp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
776
777
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
778
	  streamptr->ntsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
779
780
	  break;
	}
781
      size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
782

783
784
785
786
      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
787

788
      int rlevel = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
789
790
791
792
793
      if ( IEG_P_LevelType(iegp->ipdb) == IEG_LTYPE_HYBRID_LAYER )
	rlevel = IEG_P_Level1(iegp->ipdb);
      else
	rlevel = IEG_P_Level2(iegp->ipdb);

794
795
      if ( IEG_P_LevelType(iegp->ipdb) == 100 ) rlevel *= 100;

796
      int vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
797
798
799
800
801
802
803
804
805
      iegDateTime(iegp->ipdb, &vdate, &vtime);

      if ( nrecs == 0 )
	{
	  datetime0.date = vdate;
	  datetime0.time = vtime;
	}
      else
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
806
	  compVar.param = param;
807
          compVar.level = rlevel;
808
          int recID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
809
810
	  for ( recID = 0; recID < nrecs; recID++ )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
811
	      compVar0.param = streamptr->tsteps[0].records[recID].param;
812
813
	      compVar0.level = streamptr->tsteps[0].records[recID].ilevel;

814
	      if ( memcmp(&compVar0, &compVar, sizeof(iegcompvar_t)) == 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
815
816
	    }
	  if ( recID < nrecs ) break;
817
818
	  DateTime datetime = { .date = vdate, .time = vtime};
	  if ( datetimeCmp(datetime, datetime0) )
819
	    Warning("Inconsistent verification time for param %d level %d", param, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
820
821
822
823
824
	}

      nrecs++;

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

827
      iegAddRecord(streamptr, param, iegp->ipdb, iegp->igdb, iegp->vct, recsize, recpos, prec);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
828
829
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
830
  streamptr->rtsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
831

832
  cdi_generate_vars(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
833

834
  int taxisID = taxisCreate(TAXIS_ABSOLUTE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
835
  taxis->type  = TAXIS_ABSOLUTE;
836
837
  taxis->vdate = (int)datetime0.date;
  taxis->vtime = (int)datetime0.time;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
838

839
  int vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
840
841
  vlistDefTaxis(vlistID, taxisID);

842
  vlist_check_contents(vlistID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
843

844
  int nrecords = streamptr->tsteps[0].nallrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
845
  if ( nrecords < streamptr->tsteps[0].recordSize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
846
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
847
848
      streamptr->tsteps[0].recordSize = nrecords;
      streamptr->tsteps[0].records =
849
	(record_t *) Realloc(streamptr->tsteps[0].records,
850
                             (size_t)nrecords * sizeof (record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
851
852
    }

853
  streamptr->tsteps[0].recIDs = (int *) Malloc((size_t)nrecords * sizeof (int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
854
  streamptr->tsteps[0].nrecs = nrecords;
855
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
856
    streamptr->tsteps[0].recIDs[recID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
857

Uwe Schulzweida's avatar
Uwe Schulzweida committed
858
  if ( streamptr->ntsteps == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
859
    {
860
      tsID = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
861
      if ( tsID != streamptr->rtsteps )
862
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
863

Uwe Schulzweida's avatar
Uwe Schulzweida committed
864
865
      streamptr->tsteps[tsID-1].next   = TRUE;
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
866
867
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
868
  if ( streamptr->ntsteps == 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
869
870
871
    {
      if ( taxis->vdate == 0 && taxis->vtime == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
872
	  streamptr->ntsteps = 0;
873
	  for ( int varID = 0; varID < streamptr->nvars; varID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
874
	    {
875
	      vlistDefVarTsteptype(vlistID, varID, TSTEP_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
876
877
878
879
880
	    }
	}
    }
}

881
static
882
int iegScanTimestep2(stream_t *streamptr)
883
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
884
  off_t recpos = 0;
885
  iegcompvar_t compVar, compVar0;
886
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
887
888

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

890
891
  int vlistID = streamptr->vlistID;
  int fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
892

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

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

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

901
  cdi_create_records(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
902

903
  int nrecords = streamptr->tsteps[0].nallrecs;
904
  streamptr->tsteps[1].recIDs = (int *) Malloc((size_t)nrecords * sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
905
  streamptr->tsteps[1].nrecs = 0;
906
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
907
    streamptr->tsteps[1].recIDs[recID] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
908

909
  for ( int recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
910
    {
911
      streamptr->tsteps[tsID].records[recID].position =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
912
	streamptr->tsteps[0].records[recID].position;
913
      streamptr->tsteps[tsID].records[recID].size     =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
914
	streamptr->tsteps[0].records[recID].size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
915
916
    }

917
  for ( int rindex = 0; rindex <= nrecords; rindex++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
918
919
    {
      recpos = fileGetPos(fileID);
920
      int status = iegRead(fileID, iegp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
921
922
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
923
	  streamptr->ntsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
924
925
	  break;
	}
926
      size_t recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
927

928
929
930
      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
931

932
      int rlevel = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
933
934
935
936
937
      if ( IEG_P_LevelType(iegp->ipdb) == IEG_LTYPE_HYBRID_LAYER )
	rlevel = IEG_P_Level1(iegp->ipdb);
      else
	rlevel = IEG_P_Level2(iegp->ipdb);

938
939
      if ( IEG_P_LevelType(iegp->ipdb) == 100 ) rlevel *= 100;

940
      int vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
941
942
943
944
945
946
947
948
949
      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
950
      compVar.param = param;
951
      compVar.level = rlevel;
952
953
      bool nextstep = false;
      int recID = 0;