stream_ieg.c 36.8 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_fcommon.h"
23
#include "stream_ieg.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
24
25
26
27
28
29
30
31
32
33
34
#include "vlist.h"


#undef  UNDEFID
#define UNDEFID  CDI_UNDEFID

#define SINGLE_PRECISION  4
#define DOUBLE_PRECISION  8

#if defined (HAVE_LIBIEG)

35
36

typedef struct {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
37
  int param;
38
  int level;
39
} IEGCOMPVAR;
40
41


Thomas Jahns's avatar
Thomas Jahns committed
42
static int iegInqDatatype(int prec)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
43
44
45
{
  int datatype;

46
47
  if ( prec == DOUBLE_PRECISION ) datatype = DATATYPE_FLT64;
  else                            datatype = DATATYPE_FLT32;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
48
49
50
51
52

  return (datatype);
}


Thomas Jahns's avatar
Thomas Jahns committed
53
static int iegDefDatatype(int datatype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
54
55
56
{
  int prec;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
57
  if ( datatype == DATATYPE_CPX32 || datatype == DATATYPE_CPX64 )
58
    Error("CDI/IEG library does not support complex numbers!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
59

60
61
  if ( datatype != DATATYPE_FLT32 && datatype != DATATYPE_FLT64 )
    datatype = DATATYPE_FLT32;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
62

63
  if ( datatype == DATATYPE_FLT64 ) prec = DOUBLE_PRECISION;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
64
65
66
67
68
  else                              prec = SINGLE_PRECISION;

  return (prec);
}

69
/* not used
70
int iegInqRecord(stream_t *streamptr, int *varID, int *levelID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
71
72
73
{
  int status;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
74
  int icode, ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
75
76
  int zaxisID = -1;
  int vlistID;
77
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
78

Uwe Schulzweida's avatar
Uwe Schulzweida committed
79
80
  vlistID = streamptr->vlistID;
  fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
81
82
83
84
85
86
87

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

  status = iegRead(fileID, iegp);
  if ( status != 0 ) return (0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
88
  icode  = IEG_P_Parameter(iegp->ipdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
89
90
91
92
93
  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
94
  *varID = vlistInqVarID(vlistID, icode);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
95

96
  if ( *varID == UNDEFID ) Error("Code %d undefined", icode);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
97
98
99
100

  zaxisID = vlistInqVarZaxis(vlistID, *varID);

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
102
103
  return (1);
}
104
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
105

106
void iegReadRecord(stream_t *streamptr, double *data, int *nmiss)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
107
108
109
110
111
112
113
114
{
  int vlistID, fileID;
  int status;
  int recID, vrecID, tsID;
  off_t recpos;
  int varID, gridID;
  int i, size;
  double missval;
115
  void *iegp = streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
116

Uwe Schulzweida's avatar
Uwe Schulzweida committed
117
118
  vlistID = streamptr->vlistID;
  fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
119
120
121
122
123
  tsID    = streamptr->curTsID;
  vrecID  = streamptr->tsteps[tsID].curRecID;
  recID   = streamptr->tsteps[tsID].recIDs[vrecID];
  recpos  = streamptr->tsteps[tsID].records[recID].position;
  varID   = streamptr->tsteps[tsID].records[recID].varID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
124
125
126
127

  fileSetPos(fileID, recpos, SEEK_SET);

  status = iegRead(fileID, iegp);
128
129
  if ( status != 0 )
    Error("Could not read IEG record!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
130
131
132

  iegInqDataDP(iegp, data);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
133
  missval = vlistInqVarMissval(vlistID, varID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
134
135
  gridID  = vlistInqVarGrid(vlistID, varID);
  size    = gridInqSize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
136

Uwe Schulzweida's avatar
Uwe Schulzweida committed
137
  streamptr->numvals += size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
138
139
140
141
142
143
144
145
146
147

  *nmiss = 0;
  for ( i = 0; i < size; i++ )
    if ( DBL_IS_EQUAL(data[i], missval) || DBL_IS_EQUAL(data[i], (float)missval) )
      {
	data[i] = missval;
	(*nmiss)++;
      }
}

148
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
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
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
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;
      }
    }

  return (leveltype);
}


Thomas Jahns's avatar
Thomas Jahns committed
204
static void iegDefTime(int *pdb, int date, int time, int taxisID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
205
{
206
  int year, month, day, hour, minute, second;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
207
208
  int timetype = -1;

209
  if ( taxisID != -1 ) timetype = taxisInqType(taxisID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
210
211
212

  if ( timetype == TAXIS_ABSOLUTE || timetype == TAXIS_RELATIVE )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
213
214
      cdiDecodeDate(date, &year, &month, &day);
      cdiDecodeTime(time, &hour, &minute, &second);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231

      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
232
233
234
235
236
/* 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)
237
{
Thomas Jahns's avatar
Thomas Jahns committed
238
239
240
241
242
243
244
245
246
247
  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 )
248
    {
Thomas Jahns's avatar
Thomas Jahns committed
249
250
251
      double scaleBy = scaleFactors[j];
      bool fractionalScale = false;
      for (size_t i = 0; i < nMultTests; ++i )
252
        {
Thomas Jahns's avatar
Thomas Jahns committed
253
254
          fractionalScale = fractionalScale
            || fabs(vals[i]*scaleBy - round(vals[i]*scaleBy)) > FLT_EPSILON;
255
        }
Thomas Jahns's avatar
Thomas Jahns committed
256
      if ( !fractionalScale )
257
        {
Thomas Jahns's avatar
Thomas Jahns committed
258
          resfac = scaleBy;
259
260
261
262
          break;
        }
    }

Thomas Jahns's avatar
Thomas Jahns committed
263
  return (resfac);
264
265
}

266
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
267
268
void iegDefGrid(int *gdb, int gridID)
{
269
  int gridtype = gridInqType(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
270
271
272

  if ( gridtype == GRID_GENERIC )
    {
273
274
275
276
277
278
      int xsize, ysize;

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

      if ( (ysize == 32  || ysize == 48 || ysize == 64 ||
279
	    ysize == 96  || ysize == 160) &&
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
	   (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
295
296
297
298
299
300
    }
  else if ( gridtype == GRID_CURVILINEAR )
    {
      gridtype = GRID_LONLAT;
    }

301
  if ( gridtype == GRID_LONLAT || gridtype == GRID_GAUSSIAN )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
302
    {
303
304
      double xfirst = 0, xlast = 0, xinc = 0;
      double yfirst = 0, ylast = 0, yinc = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
305

306
307
      int nlon = gridInqXsize(gridID),
        nlat = gridInqYsize(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
308

309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
      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);
	}

331
332
333
      if ( gridtype == GRID_GAUSSIAN )
	IEG_G_GridType(gdb) = 4;
      else if ( gridtype == GRID_LONLAT && gridIsRotated(gridID) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
334
335
336
337
	IEG_G_GridType(gdb) = 10;
      else
	IEG_G_GridType(gdb) = 0;

Thomas Jahns's avatar
Thomas Jahns committed
338
339
      double resfac = calc_resfac(xfirst, xlast, xinc, yfirst, ylast, yinc);
      int iresfac = (int)resfac;
340
341
342
343
      if ( iresfac == 1000 ) iresfac = 0;

      IEG_G_ResFac(gdb)   = iresfac;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
344
345
      IEG_G_NumLon(gdb)   = nlon;
      IEG_G_NumLat(gdb)   = nlat;
346
347
348
349
350
      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);
351
      if ( fabs(xinc*resfac - IEG_G_LonIncr(gdb)) > FLT_EPSILON )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
352
353
	IEG_G_LonIncr(gdb) = 0;

354
355
356
357
      if ( gridtype == GRID_GAUSSIAN )
	IEG_G_LatIncr(gdb) = nlat/2;
      else
	{
358
	  IEG_G_LatIncr(gdb) = (int)lround(yinc*resfac);
359
	  if ( fabs(yinc*resfac - IEG_G_LatIncr(gdb)) > FLT_EPSILON )
360
	    IEG_G_LatIncr(gdb) = 0;
361
362

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

365
      if ( IEG_G_NumLon(gdb) > 1 && IEG_G_NumLat(gdb) == 1 )
366
367
	if ( IEG_G_LonIncr(gdb) != 0 && IEG_G_LatIncr(gdb) == 0 ) IEG_G_LatIncr(gdb) = IEG_G_LonIncr(gdb);

368
      if ( IEG_G_NumLon(gdb) == 1 && IEG_G_NumLat(gdb) > 1 )
369
370
371
372
373
374
375
	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
376
377
      if ( gridIsRotated(gridID) )
	{
378
379
	  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
380
381
382
383
384
385
386
387
388
	  IEG_G_Size(gdb)  = 42;
	}
      else
	{
	  IEG_G_Size(gdb)  = 32;
	}
    }
  else
    {
389
      Error("Unsupported grid type: %s", gridNamePtr(gridtype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
390
391
392
393
394
    }

  IEG_G_ScanFlag(gdb) = 64;
}

395
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
396
397
398
399
400
401
402
403
404
405
void iegDefLevel(int *pdb, int *gdb, double *vct, int zaxisID, int levelID)
{
  double level;
  int ilevel, leveltype;
  static int vct_warning = 1;

  leveltype = zaxisInqType(zaxisID);

  if ( leveltype == ZAXIS_GENERIC )
    {
406
      Message("Changed zaxis type from %s to %s",
407
408
	      zaxisNamePtr(leveltype),
	      zaxisNamePtr(ZAXIS_PRESSURE));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
409
410
      leveltype = ZAXIS_PRESSURE;
      zaxisChangeType(zaxisID, leveltype);
411
      zaxisDefUnits(zaxisID, "Pa");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
412
413
414
415
416
417
418
419
420
421
    }

  /*  IEG_G_NumVCP(gdb) = 0; */

  switch (leveltype)
    {
    case ZAXIS_SURFACE:
      {
	IEG_P_LevelType(pdb) = IEG_LTYPE_SURFACE;
	IEG_P_Level1(pdb)    = 0;
422
	IEG_P_Level2(pdb)    = (int)(zaxisInqLevel(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
423
424
425
426
427
428
429
430
431
	break;
      }
    case ZAXIS_HYBRID:
      {
	int vctsize;

	if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) )
	  {
	    IEG_P_LevelType(pdb) = IEG_LTYPE_HYBRID_LAYER;
432
433
	    IEG_P_Level1(pdb)    = (int)(zaxisInqLbound(zaxisID, levelID));
	    IEG_P_Level2(pdb)    = (int)(zaxisInqUbound(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
434
435
436
437
438
	  }
	else
	  {
	    IEG_P_LevelType(pdb) = IEG_LTYPE_HYBRID;
	    IEG_P_Level1(pdb)    = 0;
439
	    IEG_P_Level2(pdb)    = (int)(zaxisInqLevel(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
440
441
442
443
444
445
446
447
	  }

	vctsize = zaxisInqVctSize(zaxisID);
	if ( vctsize > 100 )
	  {
	    /*	    IEG_G_NumVCP(gdb) = 0; */
	    if ( vct_warning )
	      {
448
		Warning("VCT size of %d is too large (maximum is 100). Set to 0!", vctsize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
449
450
451
452
453
454
		vct_warning = 0;
	      }
	  }
	else
	  {
	    IEG_G_Size(gdb) += (vctsize*4);
455
456
	    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
457
458
459
460
461
462
463
464
465
466
	  }
	break;
      }
    case ZAXIS_PRESSURE:
      {
	double dum;
	char units[128];

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

	zaxisInqUnits(zaxisID, units);
470
	if ( memcmp(units, "hPa", 3) == 0 || memcmp(units, "mb",2 ) == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
	  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;
515
516
	    IEG_P_Level1(pdb)    = (int)(zaxisInqLbound(zaxisID, levelID));
	    IEG_P_Level2(pdb)    = (int)(zaxisInqUbound(zaxisID, levelID));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
	  }
	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:
      {
554
	Error("Unsupported zaxis type: %s", zaxisNamePtr(leveltype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
555
556
557
558
559
560
	break;
      }
    }
}


561
void iegCopyRecord(stream_t *streamptr2, stream_t *streamptr1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
562
{
563
  streamFCopyRecord(streamptr2, streamptr1, "IEG");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
564
565
566
}


567
void iegDefRecord(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
568
{
569
  int vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
570
571
572
573
  int gridID;
  int date, time;
  int datatype;
  int i;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
574
  int param, pdis, pcat, pnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
575
576
  int varID, levelID, tsID, zaxisID;
  int byteorder;
577
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
578

Uwe Schulzweida's avatar
Uwe Schulzweida committed
579
  vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
580
  byteorder = streamptr->byteorder;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
581

Uwe Schulzweida's avatar
Uwe Schulzweida committed
582
583
584
  varID   = streamptr->record->varID;
  levelID = streamptr->record->levelID;
  tsID    = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
585
586
587

  gridID  = vlistInqVarGrid(vlistID, varID);
  zaxisID = vlistInqVarZaxis(vlistID, varID);
588

589
  iegInitMem(iegp);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
590
591
592
593
  for ( i = 0; i < 37; i++ ) iegp->ipdb[i] = -1;

  iegp->byteswap = getByteswap(byteorder);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
594
  param =  vlistInqVarParam(vlistID, varID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
595
  cdiDecodeParam(param, &pnum, &pcat, &pdis);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
596
597
  IEG_P_Parameter(iegp->ipdb) = pnum;
  if ( pdis == 255 ) IEG_P_CodeTable(iegp->ipdb) = pcat;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
598
599
  date     = streamptr->tsteps[tsID].taxis.vdate;
  time     = streamptr->tsteps[tsID].taxis.vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
600
601
602
603
604

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
605
  datatype = streamptr->record->prec;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
606
607
608
609
610

  iegp->dprec = iegDefDatatype(datatype);
}


611
void iegWriteRecord(stream_t *streamptr, const double *data)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
612
613
614
615
{
  int fileID;
  int i, gridsize, gridID;
  double refval;
616
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
617

Uwe Schulzweida's avatar
Uwe Schulzweida committed
618
  fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
619
  gridID = streamptr->record->gridID;
620

Uwe Schulzweida's avatar
Uwe Schulzweida committed
621
622
623
624
625
626
627
628
629
630
631
632
633
  gridsize = gridInqSize(gridID);

  refval = data[0];
  for ( i = 1; i < gridsize; i++ )
    if ( data[i] < refval ) refval = data[i];

  iegp->refval = refval;

  iegDefDataDP(iegp, data);

  iegWrite(fileID, iegp);
}

634
static
635
void iegAddRecord(stream_t *streamptr, int param, int *pdb, int *gdb, double *vct,
636
		  size_t recsize, off_t position, int prec)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
637
638
{
  int levelID = 0;
639
640
641
642
  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
643

644
  int level1, level2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
645
646
647
648
649
650
651
652
653
  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;
654
      if ( IEG_P_LevelType(pdb) == 100 ) level1 *= 100;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
655
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
656

657
658
659
660
661
662
  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
663

664
665
666
  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
667

668
669
  grid_t *grid = (grid_t *)Malloc(sizeof (*grid));
  grid_init(grid);
670
  cdiGridTypeInit(grid, gridtype, IEG_G_NumLon(gdb)*IEG_G_NumLat(gdb));
671
672
673
674
675
  grid->xsize = IEG_G_NumLon(gdb);
  grid->ysize = IEG_G_NumLat(gdb);
  grid->xinc  = 0;
  grid->yinc  = 0;
  grid->xdef  = 0;
676
677
678
679
680

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

681
682
  /* if ( IEG_G_FirstLon != 0 || IEG_G_LastLon != 0 ) */
  {
683
    if ( grid->xsize > 1 )
684
685
      {
	if ( IEG_G_ResFlag(gdb) && IEG_G_LonIncr(gdb) > 0 )
686
	  grid->xinc = IEG_G_LonIncr(gdb) * resfac;
687
	else
688
	  grid->xinc = (IEG_G_LastLon(gdb) - IEG_G_FirstLon(gdb)) * resfac / (grid->xsize - 1);
689
690
691
692

	/* correct xinc if necessary */
	if ( IEG_G_FirstLon(gdb) == 0 && IEG_G_LastLon(gdb) > 354000 )
	  {
693
694
695
	    double xinc = 360. / grid->xsize;
            /* FIXME: why not use grid->xinc != xinc as condition? */
	    if ( fabs(grid->xinc-xinc) > 0.0 )
696
	      {
697
698
		grid->xinc = xinc;
		if ( CDI_Debug ) Message("set xinc to %g", grid->xinc);
699
700
701
	      }
	  }
      }
702
703
704
    grid->xfirst = IEG_G_FirstLon(gdb) * resfac;
    grid->xlast  = IEG_G_LastLon(gdb)  * resfac;
    grid->xdef   = 2;
705
  }
706
  grid->ydef  = 0;
707
708
  /* if ( IEG_G_FirstLat != 0 || IEG_G_LastLat != 0 ) */
  {
709
    if ( grid->ysize > 1 )
710
711
      {
	if ( IEG_G_ResFlag(gdb) && IEG_G_LatIncr(gdb) > 0 )
712
	  grid->yinc = IEG_G_LatIncr(gdb) * resfac;
713
	else
714
	  grid->yinc = (IEG_G_LastLat(gdb) - IEG_G_FirstLat(gdb)) * resfac / (grid->ysize - 1);
715
      }
716
717
718
    grid->yfirst = IEG_G_FirstLat(gdb) * resfac;
    grid->ylast  = IEG_G_LastLat(gdb)  * resfac;
    grid->ydef   = 2;
719
720
  }
  /*
721
722
723
724
725
726
727
728
  grid->xfirst= IEG_G_FirstLon(gdb) * resfac;
  grid->xlast = IEG_G_LastLon(gdb) * resfac;
  grid->xinc  = IEG_G_LonIncr(gdb) * resfac;
  grid->xdef  = 2;
  grid->yfirst= IEG_G_FirstLat(gdb) * resfac;
  grid->ylast = IEG_G_LastLat(gdb) * resfac;
  grid->yinc  = IEG_G_LatIncr(gdb) * resfac;
  grid->ydef  = 2;
729
  */
730
731
  grid->xvals = NULL;
  grid->yvals = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
732

733
  grid->isRotated = FALSE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
734
735
  if ( IEG_G_GridType(gdb) == 10 )
    {
736
737
738
739
      grid->isRotated = TRUE;
      grid->ypole     = - IEG_G_LatSP(gdb) * resfac;
      grid->xpole     =   IEG_G_LonSP(gdb) * resfac - 180;
      grid->angle     = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
740
741
    }

742
743
744
  struct addIffNewRes gridAdded = cdiVlistAddGridIfNew(vlistID, grid, 0);
  int gridID = gridAdded.Id;
  if (!gridAdded.isNew) Free(grid);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
745

746
  int leveltype = iegGetZaxisType(IEG_P_LevelType(pdb));
747

Uwe Schulzweida's avatar
Uwe Schulzweida committed
748
749
750
  if ( leveltype == ZAXIS_HYBRID )
    {
      double tmpvct[100];
751
      size_t vctsize = (size_t)IEG_G_NumVCP(gdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
752

753
754
      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
755
756
757
758

      varDefVCT(vctsize, tmpvct);
    }

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

761
  int datatype = iegInqDatatype(prec);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
762

763
  int varID;
764
  varAddRecord(recID, param, gridID, leveltype, lbounds, level1, level2, 0, 0,
765
766
	       datatype, &varID, &levelID, TSTEP_INSTANT, 0, 0, -1,
               NULL, NULL, NULL, NULL, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
767

768
769
  record->varID   = (short)varID;
  record->levelID = (short)levelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
770

Uwe Schulzweida's avatar
Uwe Schulzweida committed
771
772
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
773
774

  if ( CDI_Debug )
775
    Message("varID = %d gridID = %d levelID = %d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
776
777
778
	    varID, gridID, levelID);
}

Thomas Jahns's avatar
Thomas Jahns committed
779
#if 0
780
static
781
void iegCmpRecord(stream_t *streamptr, int tsID, int recID, off_t position, int param,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
782
783
784
785
		  int level, int xsize, int ysize)
{
  int varID = 0;
  int levelID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
786
  record_t *record;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
787
788

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
790
  if ( param != (*record).param || level != (*record).ilevel )
791
    Error("inconsistent timestep");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
792
793
794
795
796
797

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
800
801
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
802
803
  */
  if ( CDI_Debug )
804
    Message("varID = %d levelID = %d", varID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
805
}
Thomas Jahns's avatar
Thomas Jahns committed
806
#endif
807

Thomas Jahns's avatar
Thomas Jahns committed
808
static void iegDateTime(int *pdb, int *date, int *time)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
809
810
811
812
813
814
815
816
817
818
819
820
821
{
  int ryear, rmonth, rday, rhour, rminute;

  ryear   = IEG_P_Year(pdb);

  rmonth  = IEG_P_Month(pdb);
  rday    = IEG_P_Day(pdb);

  rhour   = IEG_P_Hour(pdb);
  rminute = IEG_P_Minute(pdb);

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
822
823
  *date = cdiEncodeDate(ryear, rmonth, rday);
  *time = cdiEncodeTime(rhour, rminute, 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
824
825
}

826
static
827
828
void iegScanTimestep1(stream_t *streamptr)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
829
830
831
  int prec = 0;
  int status;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
832
833
  int tabnum;
  int param = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
834
  int rcode = 0, rlevel = 0, vdate = 0, vtime = 0;
835
  DateTime datetime0 = { LONG_MIN, LONG_MIN };
Uwe Schulzweida's avatar
Uwe Schulzweida committed
836
837
  int tsID;
  int varID;
838
  size_t recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
839
840
841
  off_t recpos;
  int nrecords, nrecs, recID;
  int taxisID = -1;
842
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
843
  int vlistID;
844
  IEGCOMPVAR compVar, compVar0;
845
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
846

Uwe Schulzweida's avatar
Uwe Schulzweida committed
847
  streamptr->curTsID = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
848

849
  tsID  = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
850
  taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
851
852

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
855
  fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
856
857
858
859
860
861
862
863

  nrecs = 0;
  while ( TRUE )
    {
      recpos = fileGetPos(fileID);
      status = iegRead(fileID, iegp);
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
864
	  streamptr->ntsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
865
866
	  break;
	}
867
      recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
868
869
870

      prec   = iegp->dprec;
      rcode  = IEG_P_Parameter(iegp->ipdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
871
      tabnum = IEG_P_CodeTable(iegp->ipdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
872
      param  = cdiEncodeParam(rcode, tabnum, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
873

Uwe Schulzweida's avatar
Uwe Schulzweida committed
874
875
876
877
878
      if ( IEG_P_LevelType(iegp->ipdb) == IEG_LTYPE_HYBRID_LAYER )
	rlevel = IEG_P_Level1(iegp->ipdb);
      else
	rlevel = IEG_P_Level2(iegp->ipdb);

879
880
      if ( IEG_P_LevelType(iegp->ipdb) == 100 ) rlevel *= 100;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
881
882
883
884
885
886
887
888
889
      iegDateTime(iegp->ipdb, &vdate, &vtime);

      if ( nrecs == 0 )
	{
	  datetime0.date = vdate;
	  datetime0.time = vtime;
	}
      else
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
890
	  compVar.param = param;
891
          compVar.level = rlevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
892
893
	  for ( recID = 0; recID < nrecs; recID++ )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
894
	      compVar0.param = streamptr->tsteps[0].records[recID].param;
895
896
897
	      compVar0.level = streamptr->tsteps[0].records[recID].ilevel;

	      if ( memcmp(&compVar0, &compVar, sizeof(IEGCOMPVAR)) == 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
898
899
	    }
	  if ( recID < nrecs ) break;
900
901
	  DateTime datetime = { .date = vdate, .time = vtime};
	  if ( datetimeCmp(datetime, datetime0) )
902
	    Warning("Inconsistent verification time for param %d level %d", param, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
903
904
905
906
907
	}

      nrecs++;

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

910
      iegAddRecord(streamptr, param, iegp->ipdb, iegp->igdb, iegp->vct, recsize, recpos, prec);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
911
912
    }

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

915
  cdi_generate_vars(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
916

917
  taxisID = taxisCreate(TAXIS_ABSOLUTE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
918
  taxis->type  = TAXIS_ABSOLUTE;
919
920
  taxis->vdate = (int)datetime0.date;
  taxis->vtime = (int)datetime0.time;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
921

Uwe Schulzweida's avatar
Uwe Schulzweida committed
922
  vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
923
924
  vlistDefTaxis(vlistID, taxisID);

925
  vlist_check_contents(vlistID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
926

Uwe Schulzweida's avatar
Uwe Schulzweida committed
927
928
  nrecords = streamptr->tsteps[0].nallrecs;
  if ( nrecords < streamptr->tsteps[0].recordSize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
929
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
930
931
      streamptr->tsteps[0].recordSize = nrecords;
      streamptr->tsteps[0].records =
932
	(record_t *) Realloc(streamptr->tsteps[0].records,
933
                             (size_t)nrecords * sizeof (record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
934
935
    }

936
  streamptr->tsteps[0].recIDs = (int *) Malloc((size_t)nrecords * sizeof (int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
937
  streamptr->tsteps[0].nrecs = nrecords;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
938
  for ( recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
939
    streamptr->tsteps[0].recIDs[recID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
940

Uwe Schulzweida's avatar
Uwe Schulzweida committed
941
  if ( streamptr->ntsteps == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
942
    {
943
      tsID = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
944
      if ( tsID != streamptr->rtsteps )
945
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
946

Uwe Schulzweida's avatar
Uwe Schulzweida committed
947
948
      streamptr->tsteps[tsID-1].next   = TRUE;
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
949
950
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
951
  if ( streamptr->ntsteps == 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
952
953
954
    {
      if ( taxis->vdate == 0 && taxis->vtime == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
955
956
	  streamptr->ntsteps = 0;
	  for ( varID = 0; varID < streamptr->nvars; varID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
957
	    {
958
	      vlistDefVarTsteptype(vlistID, varID, TSTEP_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
959
960
961
962
963
	    }
	}
    }
}

964
static
965
int iegScanTimestep2(stream_t *streamptr)
966
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
967
968
  int status;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
969
970
  int tabnum;
  int param = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
971
972
973
  int rcode = 0, rlevel = 0, vdate = 0, vtime = 0;
  int tsID;
  int varID;
974
  size_t recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
975
976
977
  off_t recpos = 0;
  int nrecords, nrecs, recID, rindex;
  int nextstep;
978
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
979
  int vlistID;
980
  IEGCOMPVAR compVar, compVar0;
981
  iegrec_t *iegp = (iegrec_t*) streamptr->record->exsep;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
982
983

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
985
986
  vlistID = streamptr->vlistID;
  fileID  = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
987

Uwe Schulzweida's avatar
Uwe Schulzweida committed
988
  tsID = streamptr->rtsteps;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
989
  if ( tsID != 1 )
990
    Error("Internal problem! unexpected timestep %d", tsID+1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
991

Uwe Schulzweida's avatar
Uwe Schulzweida committed
992
  taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
993

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

996
  cdi_create_records(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
997

Uwe Schulzweida's avatar
Uwe Schulzweida committed
998
  nrecords = streamptr->tsteps[0].nallrecs;
999
  streamptr->tsteps[1].recIDs = (int *) Malloc((size_t)nrecords * sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1000
  streamptr->tsteps[1].nrecs = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1001
  for ( recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1002
    streamptr->tsteps[1].recIDs[recID] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1003
1004
1005

  for ( recID = 0; recID < nrecords; recID++ )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1006
      varID = streamptr->tsteps[0].records[recID].varID;
1007
      streamptr->tsteps[tsID].records[recID].position =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1008
	streamptr->tsteps[0].records[recID].position;
1009
      streamptr->tsteps[tsID].records[recID].size     =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1010
	streamptr->tsteps[0].records[recID].size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1011
1012
1013
1014
1015
1016
1017
1018
    }

  for ( rindex = 0; rindex <= nrecords; rindex++ )
    {
      recpos = fileGetPos(fileID);
      status = iegRead(fileID, iegp);
      if ( status != 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1019
	  streamptr->ntsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1020
1021
	  break;
	}
1022
      recsize = (size_t)(fileGetPos(fileID) - recpos);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1023
1024

      rcode  = IEG_P_Parameter(iegp->ipdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1025
      tabnum = IEG_P_CodeTable(iegp->ipdb);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1026
      param  = cdiEncodeParam(rcode, tabnum, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1027

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1028
1029
1030
1031
1032
      if ( IEG_P_LevelType(iegp->ipdb) == IEG_LTYPE_HYBRID_LAYER )
	rlevel = IEG_P_Level1(iegp->ipdb);
      else
	rlevel = IEG_P_Level2(iegp->ipdb);

1033
1034
      if ( IEG_P_LevelType(iegp->ipdb) == 100 ) rlevel *= 100;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1035
1036
1037
1038
1039
1040
1041
1042
1043
      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
1044
      compVar.param = param;
1045
      compVar.level = rlevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1046
1047
1048
      nextstep = FALSE;
      for ( recID = 0; recID < nrecords; recID++ )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1049
	  compVar0.param = streamptr->tsteps[tsID].records[recID].param;
1050
1051
1052
	  compVar0.level = streamptr->tsteps[tsID].records[recID].ilevel;

	  if ( memcmp(&compVar0, &compVar, sizeof(IEGCOMPVAR)) == 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1053
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1054
	      if ( streamptr->tsteps[tsID].records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1055
1056
1057
1058
1059
		{
		  nextstep = TRUE;
		}
	      else
		{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1060
1061
		  streamptr->tsteps[tsID].records[recID].used = TRUE;
		  streamptr->tsteps[tsID].recIDs[rindex] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1062
1063
1064
1065
1066
1067
		}
	      break;
	    }
	}
      if ( recID == nrecords )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1068
1069
	  char paramstr[32];
	  cdiParamToString(param, paramstr, sizeof(paramstr));
1070
	  Warning("param %s level %d not defined at timestep 1", paramstr, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1071
1072
1073
1074
1075
1076
	  return (CDI_EUFSTRUCT);
	}

      if ( nextstep ) break;

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1081
      compVar0.param = streamptr->tsteps[tsID].records[recID].param;
1082
      compVar0.level = streamptr->tsteps[tsID].records[recID].ilevel;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1083

1084
      if ( memcmp(&compVar0, &compVar, sizeof(IEGCOMPVAR)) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1085
	{
1086
	  Message("tsID = %d recID = %d param = %3d new %3d  level = %3d new %3d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1087
		  tsID, recID,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1088
		  streamptr->tsteps[tsID].records[recID].param, param,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1089
		  streamptr->tsteps[tsID].records[recID].ilevel, rlevel);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1090
1091
1092
	  return (CDI_EUFSTRUCT);
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1093
      streamptr->tsteps[1].records[recID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1094
1095
1096
1097
1098
    }

  nrecs = 0;
  for ( recID = 0; recID < nrecords; recID++ )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1099
      if ( ! streamptr->tsteps[tsID].records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1100
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1101
	  varID = streamptr->tsteps[tsID].records[recID].varID;
1102
          vlistDefVarTsteptype(vlistID, varID, TSTEP_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1103
1104
1105
1106
1107
1108
	}
      else
	{