stream_cgribex.c 58.4 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>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
6
7
8
9
#include <stdio.h>

#include "dmemory.h"
#include "cdi.h"
10
#include "cdi_int.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
11
12
13
14
#include "file.h"
#include "varscan.h"
#include "datetime.h"
#include "vlist.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
15
#include "stream_grb.h"
16
#include "stream_cgribex.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
17

Uwe Schulzweida's avatar
Uwe Schulzweida committed
18
19
20
#if  defined  (HAVE_LIBCGRIBEX)
#  include "cgribex.h"
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
21
22

typedef struct {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
23
  int param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
24
25
26
  int level1;
  int level2;
  int ltype;
27
  int tsteptype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
28
} compvar_t;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
29
30


Uwe Schulzweida's avatar
Uwe Schulzweida committed
31
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
32
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
33
int cgribexGetGridType(int *isec2)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
34
{
35
  int gridtype = GRID_GENERIC;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
36
37
38

  switch (ISEC2_GridType)
    {
39
    case  GRIB1_GTYPE_LATLON:     { gridtype = GRID_LONLAT;     break; }
40
41
    case  GRIB1_GTYPE_LATLON_ROT: { gridtype = GRID_PROJECTION; break; }
    case  GRIB1_GTYPE_LCC:        { gridtype = GRID_LCC;        break; }
42
43
44
45
46
47
    case  GRIB1_GTYPE_GAUSSIAN:   { if ( ISEC2_Reduced )
	                              gridtype = GRID_GAUSSIAN_REDUCED;
                         	    else
				      gridtype = GRID_GAUSSIAN;
          	                    break;
                                  }
48
49
    case  GRIB1_GTYPE_SPECTRAL:   { gridtype = GRID_SPECTRAL;   break; }
    case  GRIB1_GTYPE_GME:        { gridtype = GRID_GME;        break; }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
50
51
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
52
  return gridtype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
53
54
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
55
static
56
bool cgribexGetIsRotated(int *isec2)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
57
{
58
  return (ISEC2_GridType == GRIB1_GTYPE_LATLON_ROT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
59
60
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
61
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
62
int cgribexGetZaxisHasBounds(int grb_ltype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
63
64
65
{
  int lbounds = 0;

66
  switch (grb_ltype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
67
    {
68
    case GRIB1_LTYPE_SIGMA_LAYER:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
69
70
    case GRIB1_LTYPE_HYBRID_LAYER:
    case GRIB1_LTYPE_LANDDEPTH_LAYER:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
71
72
73
74
75
76
      {
	lbounds = 1;
	break;
      }
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
77
  return lbounds;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
78
79
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
80
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
81
int cgribexGetTimeUnit(int *isec1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
82
{
83
  int timeunit = TUNIT_HOUR;
84
  static bool lprint = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
85
86
87

  switch ( ISEC1_TimeUnit )
    {
88
89
90
91
92
93
94
95
    case ISEC1_TABLE4_MINUTE:    timeunit = TUNIT_MINUTE;    break;
    case ISEC1_TABLE4_QUARTER:   timeunit = TUNIT_QUARTER;   break;
    case ISEC1_TABLE4_30MINUTES: timeunit = TUNIT_30MINUTES; break;
    case ISEC1_TABLE4_HOUR:      timeunit = TUNIT_HOUR;      break;
    case ISEC1_TABLE4_3HOURS:    timeunit = TUNIT_3HOURS;    break;
    case ISEC1_TABLE4_6HOURS:    timeunit = TUNIT_6HOURS;    break;
    case ISEC1_TABLE4_12HOURS:   timeunit = TUNIT_12HOURS;   break;
    case ISEC1_TABLE4_DAY:       timeunit = TUNIT_DAY;       break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
96
97
98
    default:
      if ( lprint )
	{
99
	  Message("GRIB time unit %d unsupported!", ISEC1_TimeUnit);
100
	  lprint = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
101
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
102
      break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
103
104
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
105
  return timeunit;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
106
107
}

108
static
109
bool cgribexTimeIsFC(int *isec1)
110
{
111
  bool isFC = (ISEC1_TimeRange == 10 && ISEC1_TimePeriod1 == 0 && ISEC1_TimePeriod2 == 0) ? false : true;
112

Uwe Schulzweida's avatar
Uwe Schulzweida committed
113
  return isFC;
114
115
116
117
118
}

static
int cgribexGetTsteptype(int timerange)
{
119
  int tsteptype = TSTEP_INSTANT;
120
  static bool lprint = true;
121
122
123
124
125
126
127
128
129
130
131
132
133

  switch ( timerange )
    {
    case  0:  tsteptype = TSTEP_INSTANT;  break;
    case  1:  tsteptype = TSTEP_INSTANT2; break;
    case  2:  tsteptype = TSTEP_RANGE;    break;
    case  3:  tsteptype = TSTEP_AVG;      break;
    case  4:  tsteptype = TSTEP_ACCUM;    break;
    case  5:  tsteptype = TSTEP_DIFF;     break;
    case 10:  tsteptype = TSTEP_INSTANT3; break;
    default:
      if ( lprint )
	{
134
	  Message("Time range indicator %d unsupported, set to 0!", timerange);
135
	  lprint = false;
136
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
137
      break;
138
139
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
140
  return tsteptype;
141
142
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
143
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
144
void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, int iret)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
145
{
146
  bool compyinc = true;
147
  int gridtype = cgribexGetGridType(isec2);
148
  int projtype = (gridtype == GRID_PROJECTION && cgribexGetIsRotated(isec2)) ? CDI_PROJ_RLL : CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
149

150
  if ( streamptr->unreduced && gridtype == GRID_GAUSSIAN_REDUCED && iret != -801 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
151
    {
152
153
154
      int ilat, nlon = 0;
      for ( ilat = 0; ilat < ISEC2_NumLat; ++ilat )
        if ( ISEC2_RowLon(ilat) > nlon ) nlon = ISEC2_RowLon(ilat);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
155
      gridtype = GRID_GAUSSIAN;
156
157
      ISEC2_NumLon = nlon;
      ISEC4_NumValues = nlon*ISEC2_NumLat;
158
      compyinc = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
159
160
    }

161
  grid_init(grid);
162
  cdiGridTypeInit(grid, gridtype, 0);
163

164
  if ( gridtype == GRID_LONLAT || gridtype == GRID_GAUSSIAN || projtype == CDI_PROJ_RLL )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
165
    {
166
167
168
169
170
171
172
173
174
175
      if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
        Error("numberOfPoints (%d) and gridSize (%d) differ!", ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);
      grid->size  = ISEC4_NumValues;
      grid->x.size = ISEC2_NumLon;
      grid->y.size = ISEC2_NumLat;
      if ( gridtype == GRID_GAUSSIAN ) grid->np = ISEC2_NumPar;
      grid->x.inc  = 0;
      grid->y.inc  = 0;
      grid->x.flag = 0;
      /* if ( ISEC2_FirstLon != 0 || ISEC2_LastLon != 0 ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
176
      {
177
178
179
        if ( grid->x.size > 1 )
          {
            bool recompinc = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
180

181
182
183
184
185
            if ( ISEC2_LastLon < ISEC2_FirstLon && ISEC2_LastLon < 0 ) ISEC2_LastLon += 360000;

            if ( ISEC2_ResFlag && ISEC2_LonIncr > 0 )
              {
                if ( abs(ISEC2_LastLon - (ISEC2_FirstLon+ISEC2_LonIncr*(grid->x.size-1))) <= 2 )
186
                  {
187
188
                    recompinc = false;
                    grid->x.inc = ISEC2_LonIncr * 0.001;
189
                  }
190
              }
191

192
193
            /* recompute xinc if necessary */
            if ( recompinc ) grid->x.inc = (ISEC2_LastLon - ISEC2_FirstLon) * 0.001 / (grid->x.size-1);
194

195
196
197
198
199
200
201
202
203
204
205
206
207
208
            /* correct xinc if necessary */
            if ( ISEC2_FirstLon == 0 && ISEC2_LastLon > 354000 && ISEC2_LastLon < 360000 )
              {
                double xinc = 360. / grid->x.size;
                if ( fabs(grid->x.inc-xinc) > 0.0 )
                  {
                    grid->x.inc = xinc;
                    if ( CDI_Debug ) Message("set xinc to %g", grid->x.inc);
                  }
              }
          }
        grid->x.first = ISEC2_FirstLon * 0.001;
        grid->x.last  = ISEC2_LastLon  * 0.001;
        grid->x.flag  = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
209
      }
210
211
      grid->y.flag = 0;
      /* if ( ISEC2_FirstLat != 0 || ISEC2_LastLat != 0 ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
212
      {
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
        if ( grid->y.size > 1 && compyinc )
          {
            bool recompinc = true;
            if ( ISEC2_ResFlag && ISEC2_LatIncr > 0 )
              {
                if ( abs(ISEC2_LastLat - (ISEC2_FirstLat+ISEC2_LatIncr*(grid->y.size-1))) <= 2 )
                  {
                    recompinc = false;
                    grid->y.inc = ISEC2_LatIncr * 0.001;
                  }
              }

            /* recompute yinc if necessary */
            if ( recompinc ) grid->y.inc = (ISEC2_LastLat - ISEC2_FirstLat) * 0.001 / (grid->y.size - 1);
          }
        grid->y.first = ISEC2_FirstLat * 0.001;
        grid->y.last  = ISEC2_LastLat  * 0.001;
        grid->y.flag  = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
231
      }
232
233
234
235
236
237
238
239
240
241
242
243
    }
  else if ( gridtype == GRID_GAUSSIAN_REDUCED )
    {
      grid->np      = ISEC2_NumPar;
      grid->size    = ISEC4_NumValues;
      grid->rowlon  = ISEC2_RowLonPtr;
      grid->nrowlon = ISEC2_NumLat;
      grid->y.size  = ISEC2_NumLat;
      grid->x.inc   = 0;
      grid->y.inc   = 0;
      grid->x.flag  = 0;
      /* if ( ISEC2_FirstLon != 0 || ISEC2_LastLon != 0 ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
244
      {
245
246
247
248
249
250
251
252
253
254
255
256
        if ( grid->x.size > 1 )
          {
            if ( ISEC2_LastLon < ISEC2_FirstLon && ISEC2_LastLon < 0 ) ISEC2_LastLon += 360000;

            if ( ISEC2_ResFlag && ISEC2_LonIncr > 0 )
              grid->x.inc = ISEC2_LonIncr * 0.001;
            else
              grid->x.inc = (ISEC2_LastLon - ISEC2_FirstLon) * 0.001 / (grid->x.size - 1);
          }
        grid->x.first = ISEC2_FirstLon * 0.001;
        grid->x.last  = ISEC2_LastLon  * 0.001;
        grid->x.flag  = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
257
      }
258
259
      grid->y.flag = 0;
      /* if ( ISEC2_FirstLat != 0 || ISEC2_LastLat != 0 ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
260
      {
261
262
263
264
265
266
267
268
269
270
        if ( grid->y.size > 1 )
          {
            if ( ISEC2_ResFlag && ISEC2_LatIncr > 0 )
              grid->y.inc = ISEC2_LatIncr * 0.001;
            else
              grid->y.inc = (ISEC2_LastLat - ISEC2_FirstLat) * 0.001 / (grid->y.size - 1);
          }
        grid->y.first = ISEC2_FirstLat * 0.001;
        grid->y.last  = ISEC2_LastLat  * 0.001;
        grid->y.flag  = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
271
272
      }
    }
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
  else if ( gridtype == GRID_LCC )
    {
      if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
        Error("numberOfPoints (%d) and gridSize (%d) differ!",
              ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);

      grid->size  = ISEC4_NumValues;
      grid->x.size = ISEC2_NumLon;
      grid->y.size = ISEC2_NumLat;

      grid->lcc.xinc      = ISEC2_Lambert_dx;
      grid->lcc.yinc      = ISEC2_Lambert_dy;
      grid->lcc.originLon = ISEC2_FirstLon * 0.001;
      grid->lcc.originLat = ISEC2_FirstLat * 0.001;
      grid->lcc.lonParY   = ISEC2_Lambert_Lov * 0.001;
      grid->lcc.lat1      = ISEC2_Lambert_LatS1 * 0.001;
      grid->lcc.lat2      = ISEC2_Lambert_LatS2 * 0.001;
      grid->lcc.projflag  = ISEC2_Lambert_ProjFlag;
      grid->lcc.scanflag  = ISEC2_ScanFlag;

      grid->x.flag = 0;
      grid->y.flag = 0;
    }
  else if ( gridtype == GRID_SPECTRAL )
    {
      grid->size  = ISEC4_NumValues;
      grid->trunc = ISEC2_PentaJ;
      if ( ISEC2_RepMode == 2 )
        grid->lcomplex = 1;
      else
        grid->lcomplex = 0;
    }
  else if ( gridtype == GRID_GME )
    {
      grid->size  = ISEC4_NumValues;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
308
309
310
311
      grid->gme.nd  = ISEC2_GME_ND;
      grid->gme.ni  = ISEC2_GME_NI;
      grid->gme.ni2 = ISEC2_GME_NI2;
      grid->gme.ni3 = ISEC2_GME_NI3;
312
313
314
315
316
317
318
319
320
321
322
    }
  else if ( gridtype == GRID_GENERIC )
    {
      grid->size  = ISEC4_NumValues;
      grid->x.size = 0;
      grid->y.size = 0;
    }
  else
    {
      Error("Unsupported grid type: %s", gridNamePtr(gridtype));
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
323

324
325
  grid->type = gridtype;
  grid->projtype = projtype;
326
327
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
328
329
330
331
332
333
334
335
336
337
static
void cgribexGetLevel(int *isec1, int *leveltype, int *level1, int *level2)
{
  *leveltype = ISEC1_LevelType;
  *level1 = ISEC1_Level1;
  *level2 = ISEC1_Level2;
  if ( *leveltype == GRIB1_LTYPE_ISOBARIC ) *level1 *= 100;
  else if ( *leveltype == GRIB1_LTYPE_99 || *leveltype == GRIB1_LTYPE_ISOBARIC_PA ) *leveltype = GRIB1_LTYPE_ISOBARIC;
}

338
static
339
void cgribexAddRecord(stream_t * streamptr, int param, int *isec1, int *isec2, double *fsec2, double *fsec3,
340
		      int *isec4, size_t recsize, off_t position, int datatype, int comptype, int lmv, int iret)
341
{
342
  int varID;
343
344
  int levelID = 0;

345
346
347
  int vlistID = streamptr->vlistID;
  int tsID    = streamptr->curTsID;
  int recID   = recordNewEntry(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
348
  record_t *record = &streamptr->tsteps[tsID].records[recID];
349

350
351
  int tsteptype = cgribexGetTsteptype(ISEC1_TimeRange);
  int numavg    = ISEC1_AvgNum;
352

Uwe Schulzweida's avatar
Uwe Schulzweida committed
353
354
  int leveltype, level1, level2;
  cgribexGetLevel(isec1, &leveltype, &level1, &level2);
355

Uwe Schulzweida's avatar
Uwe Schulzweida committed
356
  /* fprintf(stderr, "param %d %d %d %d\n", param, level1, level2, leveltype); */
357

358
  record->size      = recsize;
359
360
361
362
  record->position  = position;
  record->param     = param;
  record->ilevel    = level1;
  record->ilevel2   = level2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
363
  record->ltype     = leveltype;
Thomas Jahns's avatar
Thomas Jahns committed
364
  record->tsteptype = (short)tsteptype;
365

Uwe Schulzweida's avatar
Uwe Schulzweida committed
366
  grid_t *gridptr = (grid_t*) Malloc(sizeof(*gridptr));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
367
  cgribexGetGrid(streamptr, isec2, isec4, gridptr, iret);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
368

Uwe Schulzweida's avatar
Uwe Schulzweida committed
369
  struct addIfNewRes gridAdded = cdiVlistAddGridIfNew(vlistID, gridptr, 0);
370
  int gridID = gridAdded.Id;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
371
372
373
374
375
376
377
378
379
  if ( gridAdded.isNew )
    {
      if ( gridptr->nrowlon )
        {
          size_t nrowlon = (size_t) gridptr->nrowlon;
          int *rowlon = gridptr->rowlon;
          gridptr->rowlon = (int*) Malloc(nrowlon * sizeof(int));
          memcpy(gridptr->rowlon, rowlon, nrowlon * sizeof(int));
        }
380
381
382
383
384
      else if ( gridptr->projtype == CDI_PROJ_RLL )
        {
          double xpole =   ISEC2_LonSP*0.001 - 180;
          double ypole = - ISEC2_LatSP*0.001;
          double angle = - FSEC2_RotAngle;
385
          gridDefParamRLL(gridID, xpole, ypole, angle);
386
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
387
388
389
    }
  else
    Free(gridptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
390

Uwe Schulzweida's avatar
Uwe Schulzweida committed
391
  int zaxistype = grib1ltypeToZaxisType(leveltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
392

393
  if ( zaxistype == ZAXIS_HYBRID || zaxistype == ZAXIS_HYBRID_HALF )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
394
    {
395
      size_t vctsize = (size_t)ISEC2_NumVCP;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
396
397
398
399
400
      double *vctptr = &fsec2[10];

      varDefVCT(vctsize, vctptr);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
401
  int lbounds = cgribexGetZaxisHasBounds(leveltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
402

403
404
  if ( datatype > 32 ) datatype = CDI_DATATYPE_PACK32;
  if ( datatype <  0 ) datatype = CDI_DATATYPE_PACK;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
405

406
  varAddRecord(recID, param, gridID, zaxistype, lbounds, level1, level2, 0, 0,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
407
	       datatype, &varID, &levelID, tsteptype, numavg, leveltype, -1,
408
               NULL, NULL, NULL, NULL, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
409

Thomas Jahns's avatar
Thomas Jahns committed
410
411
  record->varID   = (short)varID;
  record->levelID = (short)levelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
412

413
  varDefCompType(varID, comptype);
414

415
416
417
418
419
420
421
  if ( ISEC1_LocalFLag )
    {
      if      ( ISEC1_CenterID == 78  && isec1[36] == 253 ) // DWD local extension
        varDefEnsembleInfo(varID, isec1[54], isec1[53], isec1[52]);
      else if ( ISEC1_CenterID == 252 && isec1[36] ==   1 ) // MPIM local extension
        varDefEnsembleInfo(varID, isec1[38], isec1[39], isec1[37]);
    }
422

423
424
  if ( lmv ) varDefMissval(varID, FSEC3_MissVal);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
425
  if ( varInqInst(varID) == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
426
    {
427
428
429
      int center    = ISEC1_CenterID;
      int subcenter = ISEC1_SubCenterID;
      int instID    = institutInq(center, subcenter, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
430
      if ( instID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
431
432
433
434
	instID = institutDef(center, subcenter, NULL, NULL);
      varDefInst(varID, instID);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
435
  if ( varInqModel(varID) == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
436
    {
437
      int modelID = modelInq(varInqInst(varID), ISEC1_ModelID, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
438
      if ( modelID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
439
440
441
442
	modelID = modelDef(varInqInst(varID), ISEC1_ModelID, NULL);
      varDefModel(varID, modelID);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
443
  if ( varInqTable(varID) == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
444
    {
445
      int tableID = tableInq(varInqModel(varID), ISEC1_CodeTable, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
446
      if ( tableID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
447
448
449
450
	tableID = tableDef(varInqModel(varID), ISEC1_CodeTable, NULL);
      varDefTable(varID, tableID);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
451
452
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
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
static
void MCH_get_undef(int *isec1, double *undef_pds, double *undef_eps)
{
  /* 2010-01-13: Oliver Fuhrer */
  if ( ISEC1_CenterID == 215 ) {
    if (isec1[34] != 0 && isec1[34] != 255) {
      if (isec1[34] & 2) {
        if (isec1[34] & 1) {
          *undef_pds = -0.99*pow(10.0,-isec1[35]);
        } else {
          *undef_pds = +0.99*pow(10.0,-isec1[35]);
        }
        *undef_eps = pow(10.0,-isec1[35]-1);
      } else {
        if (isec1[34] & 1) {
          *undef_pds = -0.99*pow(10.0,+isec1[35]);
        } else {
          *undef_pds = +0.99*pow(10.0,+isec1[35]);
        }
        *undef_eps = pow(10.0,isec1[35]-1);
      }
    }
  }
}

static
void cgribexDecodeHeader(int *isec0, int *isec1, int *isec2, double *fsec2,
482
			 int *isec3, double *fsec3, int *isec4, double *fsec4,
483
			 int *gribbuffer, int recsize, int *lmv, int *iret)
484
{
485
  int ipunp = 0, iword = 0;
486

487
488
  memset(isec1, 0, 256*sizeof(int));

489
  gribExDP(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, fsec4,
490
	   ipunp, (int *) gribbuffer, recsize, &iword, "J", iret);
491
492
493
494
495
496
497
498
499
500
501
502

  *lmv = 0;

  if ( ISEC1_CenterID == 215 && (isec1[34] != 0 && isec1[34] != 255) )
    {
      double undef_pds, undef_eps;

      MCH_get_undef(isec1, &undef_pds, &undef_eps);
      FSEC3_MissVal = undef_pds;
      *lmv = 1;
    }
}
503
504

static
505
compvar_t cgribexVarSet(int param, int level1, int level2, int leveltype, int trange)
506
507
{
  compvar_t compVar;
508
  int tsteptype = cgribexGetTsteptype(trange);
509

510
511
512
513
514
  compVar.param     = param;
  compVar.level1    = level1;
  compVar.level2    = level2;
  compVar.ltype     = leveltype;
  compVar.tsteptype = tsteptype;
515

Uwe Schulzweida's avatar
Uwe Schulzweida committed
516
  return compVar;
517
518
}

Thomas Jahns's avatar
Thomas Jahns committed
519
520
static inline int
cgribexVarCompare(compvar_t compVar, record_t record, int flag)
521
{
522
523
524
525
  bool vinst = compVar.tsteptype == TSTEP_INSTANT || compVar.tsteptype == TSTEP_INSTANT2 || compVar.tsteptype == TSTEP_INSTANT3;
  bool rinst = record.tsteptype == TSTEP_INSTANT || record.tsteptype == TSTEP_INSTANT2 || record.tsteptype == TSTEP_INSTANT3;
  int tstepDiff = (!((flag == 0) & (vinst && rinst)))
                & (compVar.tsteptype != record.tsteptype);
Thomas Jahns's avatar
Thomas Jahns committed
526
527
528
529
530
  int rstatus = (compVar.param != record.param)
    |           (compVar.level1 != record.ilevel)
    |           (compVar.level2 != record.ilevel2)
    |           (compVar.ltype != record.ltype)
    |           tstepDiff;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
531
  return rstatus;
532
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
533
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
534

535
536
537
#define gribWarning(text, nrecs, timestep, paramstr, level1, level2) \
            Warning("Record %2d (id=%s lev1=%d lev2=%d) timestep %d: %s", nrecs, paramstr, level1, level2, timestep, text)

538
#if  defined  (HAVE_LIBCGRIBEX)
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562

static inline void
cgribexScanTsFixNtsteps(stream_t *streamptr, off_t recpos)
{
  if ( streamptr->ntsteps == -1 )
    {
      int tsID = tstepsNewEntry(streamptr);
      if ( tsID != streamptr->rtsteps )
	Error("Internal error. tsID = %d", tsID);

      streamptr->tsteps[tsID-1].next   = TRUE;
      streamptr->tsteps[tsID].position = recpos;
    }
}

static inline void
cgribexScanTsConstAdjust(stream_t *streamptr, taxis_t *taxis)
{
  int vlistID = streamptr->vlistID;
  if ( streamptr->ntsteps == 1 )
    {
      if ( taxis->vdate == 0 && taxis->vtime == 0 )
	{
	  streamptr->ntsteps = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
563
564
	  for ( int varID = 0; varID < streamptr->nvars; varID++ )
            vlistDefVarTsteptype(vlistID, varID, TSTEP_CONSTANT);
565
566
567
568
569
	}
    }
}


Uwe Schulzweida's avatar
Uwe Schulzweida committed
570
int cgribexScanTimestep1(stream_t *streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
571
572
{
  double fsec2[512], fsec3[2], *fsec4 = NULL;
573
  int lmv = 0, iret = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
574
  off_t recpos = 0;
575
  void *gribbuffer = NULL;
576
  size_t buffersize = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
577
  int rstatus;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
578
  int param = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
579
  int leveltype = 0, level1 = 0, level2 = 0, vdate = 0, vtime = 0;
580
  DateTime datetime, datetime0 = { LONG_MIN, LONG_MIN };
Uwe Schulzweida's avatar
Uwe Schulzweida committed
581
  size_t readsize;
582
  unsigned nrecords, recID;
583
  int nrecs_scanned = 0;
584
  int datatype;
585
  size_t recsize = 0;
586
587
  bool warn_time = true;
  bool warn_numavg = true;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
588
  int taxisID = -1;
589
590
  int rdate = 0, rtime = 0, tunit = 0;
  bool fcast = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
591
  int vlistID;
592
  int comptype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
593
  size_t unzipsize;
594
  char paramstr[32];
595
  int nskip = cdiSkipRecords;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
596
597

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

599
600
601
602
603
  int *isec0 = streamptr->record->sec0;
  int *isec1 = streamptr->record->sec1;
  int *isec2 = streamptr->record->sec2;
  int *isec3 = streamptr->record->sec3;
  int *isec4 = streamptr->record->sec4;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
604

605
606
  int tsID  = tstepsNewEntry(streamptr);
  taxis_t *taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
607
608

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

611
  int fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
612

613
614
615
  while ( nskip-- > 0 )
    {
      recsize = gribGetSize(fileID);
616
      if ( recsize == 0 )
617
	Error("Skipping of %d records failed!", cdiSkipRecords);
618

619
      recpos  = fileGetPos(fileID);
620
      fileSetPos(fileID, (off_t)recsize, SEEK_CUR);
621
622
    }

623
  unsigned nrecs = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
624
625
626
627
628
629
630
  while ( TRUE )
    {
      recsize = gribGetSize(fileID);
      recpos  = fileGetPos(fileID);

      if ( recsize == 0 )
	{
631
	  if ( nrecs == 0 )
632
	    Error("No GRIB records found!");
633

Uwe Schulzweida's avatar
Uwe Schulzweida committed
634
	  streamptr->ntsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
635
636
	  break;
	}
637
      if ( recsize > buffersize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
638
	{
639
	  buffersize = recsize;
640
	  gribbuffer = Realloc(gribbuffer, buffersize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
641
642
	}

643
      readsize = recsize;
644
      rstatus = gribRead(fileID, (unsigned char *)gribbuffer, &readsize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
645
      if ( rstatus ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
646

647
      comptype = CDI_COMPRESS_NONE;
648
      if ( gribGetZip(recsize, (unsigned char *)gribbuffer, &unzipsize) > 0 )
649
	{
650
	  comptype = CDI_COMPRESS_SZIP;
651
	  unzipsize += 100; /* need 0 to 1 bytes for rounding of bds */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
652
	  if ( buffersize < unzipsize )
653
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
654
	      buffersize = unzipsize;
655
	      gribbuffer = Realloc(gribbuffer, buffersize);
656
657
658
	    }
	}

659
      nrecs_scanned++;
660
      cgribexDecodeHeader(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, fsec4,
661
			  (int *) gribbuffer, (int)recsize, &lmv, &iret);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
662

Uwe Schulzweida's avatar
Uwe Schulzweida committed
663
      param = cdiEncodeParam(ISEC1_Parameter, ISEC1_CodeTable, 255);
664
665
      cdiParamToString(param, paramstr, sizeof(paramstr));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
666
      cgribexGetLevel(isec1, &leveltype, &level1, &level2);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
667

Uwe Schulzweida's avatar
Uwe Schulzweida committed
668
669
      gribDateTime(isec1, &vdate, &vtime);

670
      if ( ISEC4_NumBits > 0 && ISEC4_NumBits <= 32 )
671
	datatype = ISEC4_NumBits;
672
      else
673
        datatype = CDI_DATATYPE_PACK;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
674
675
676
677
678
679
680

      if ( nrecs == 0 )
	{
	  datetime0.date = vdate;
	  datetime0.time = vtime;
	  rdate = gribRefDate(isec1);
	  rtime = gribRefTime(isec1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
681
	  tunit = cgribexGetTimeUnit(isec1);
682
	  fcast = cgribexTimeIsFC(isec1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
683
684
685
	}
      else
	{
686
687
	  datetime.date = vdate;
	  datetime.time = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
688
	  compvar_t compVar = cgribexVarSet(param, level1, level2, leveltype, ISEC1_TimeRange);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
689
690
	  for ( recID = 0; recID < nrecs; recID++ )
	    {
691
	      if ( cgribexVarCompare(compVar, streamptr->tsteps[0].records[recID], 0) == 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
692
	    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
693
694
695
696
697

	  if ( cdiInventoryMode == 1 )
	    {
	      if ( recID < nrecs ) break;
	      if ( warn_time )
698
		if ( datetimeCmp(datetime, datetime0) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
699
		  {
700
                    gribWarning("Inconsistent verification time!", nrecs_scanned, tsID+1, paramstr, level1, level2);
701
		    warn_time = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
702
703
704
705
		  }
	    }
	  else
	    {
706
	      if ( datetimeCmp(datetime, datetime0) != 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
707
708
709

	      if ( recID < nrecs )
		{
710
		  gribWarning("Parameter already exist, skipped!", nrecs_scanned, tsID+1, paramstr, level1, level2);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
711
712
713
		  continue;
		}
	    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
714
715
716
717
718
719
	}

      if ( ISEC1_AvgNum )
	{
	  if (  taxis->numavg && warn_numavg && (taxis->numavg != ISEC1_AvgNum) )
	    {
720
	      Warning("Changing numavg from %d to %d not supported!", taxis->numavg, ISEC1_AvgNum);
721
	      warn_numavg = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
722
723
724
725
726
727
728
729
730
731
	    }
	  else
	    {
	      taxis->numavg = ISEC1_AvgNum;
	    }
	}

      nrecs++;

      if ( CDI_Debug )
732
	Message("Read record %2d (id=%s lev1=%d lev2=%d) %8d %6d", nrecs_scanned, paramstr, level1, level2, vdate, vtime);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
733

734
      cgribexAddRecord(streamptr, param, isec1, isec2, fsec2, fsec3,
735
		       isec4, recsize, recpos, datatype, comptype, lmv, iret);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
736
737
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
738
  streamptr->rtsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
739

Uwe Schulzweida's avatar
Uwe Schulzweida committed
740
  if ( nrecs == 0 ) return CDI_EUFSTRUCT;
741

742
  cdi_generate_vars(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
743
744
745

  if ( fcast )
    {
746
      taxisID = taxisCreate(TAXIS_RELATIVE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
747
748
749
750
751
752
753
      taxis->type  = TAXIS_RELATIVE;
      taxis->rdate = rdate;
      taxis->rtime = rtime;
      taxis->unit  = tunit;
    }
  else
    {
754
      taxisID = taxisCreate(TAXIS_ABSOLUTE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
755
      taxis->type  = TAXIS_ABSOLUTE;
756
      taxis->unit  = tunit;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
757
758
    }

759
760
  taxis->vdate = (int)datetime0.date;
  taxis->vtime = (int)datetime0.time;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
761

Uwe Schulzweida's avatar
Uwe Schulzweida committed
762
  vlistID = streamptr->vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
763
764
  vlistDefTaxis(vlistID, taxisID);

765
766
  nrecords = (unsigned)streamptr->tsteps[0].nallrecs;
  if ( nrecords < (unsigned)streamptr->tsteps[0].recordSize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
767
    {
768
      streamptr->tsteps[0].recordSize = (int)nrecords;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
769
      streamptr->tsteps[0].records =
770
      (record_t *) Realloc(streamptr->tsteps[0].records, nrecords*sizeof(record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
771
772
    }

773
  streamptr->tsteps[0].recIDs = (int *) Malloc(nrecords*sizeof(int));
774
  streamptr->tsteps[0].nrecs = (int)nrecords;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
775
  for ( recID = 0; recID < nrecords; recID++ )
776
    streamptr->tsteps[0].recIDs[recID] = (int)recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
777

Uwe Schulzweida's avatar
Uwe Schulzweida committed
778
  streamptr->record->buffer     = gribbuffer;
779
  streamptr->record->buffersize = buffersize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
780

781
782
  cgribexScanTsFixNtsteps(streamptr, recpos);
  cgribexScanTsConstAdjust(streamptr, taxis);
783

Uwe Schulzweida's avatar
Uwe Schulzweida committed
784
  return 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
785
786
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
787

788
int cgribexScanTimestep2(stream_t * streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
789
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
790
  int rstatus = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
791
  double fsec2[512], fsec3[2], *fsec4 = NULL;
792
  int lmv = 0, iret = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
793
  off_t recpos = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
794
  int param = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
795
  int leveltype = 0, level1 = 0, level2 = 0, vdate = 0, vtime = 0;
796
  DateTime datetime, datetime0 = { LONG_MIN, LONG_MIN };
Uwe Schulzweida's avatar
Uwe Schulzweida committed
797
798
  int varID, gridID;
  size_t readsize;
799
  int nrecs, recID;
800
  size_t recsize = 0;
801
  bool warn_numavg = true;
802
  int tsteptype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
803
  size_t unzipsize;
804
  char paramstr[32];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
805

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

808
809
810
811
812
  int *isec0 = streamptr->record->sec0;
  int *isec1 = streamptr->record->sec1;
  int *isec2 = streamptr->record->sec2;
  int *isec3 = streamptr->record->sec3;
  int *isec4 = streamptr->record->sec4;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
813

814
815
816
  int fileID  = streamptr->fileID;
  int vlistID = streamptr->vlistID;
  int taxisID = vlistInqTaxis(vlistID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
817

818
819
  void *gribbuffer = streamptr->record->buffer;
  size_t buffersize = streamptr->record->buffersize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
820

821
  int tsID = streamptr->rtsteps;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
822
  if ( tsID != 1 )
Thomas Jahns's avatar
Thomas Jahns committed
823
    Error("Internal problem! unexpected timestep %d", tsID+1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
824

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

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

829
  cdi_create_records(streamptr, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
830

831
  int nrecords = streamptr->tsteps[tsID].nallrecs;
832
  if ( nrecords ) streamptr->tsteps[1].recIDs = (int *) Malloc((size_t)nrecords * sizeof(int));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
833
  streamptr->tsteps[1].nrecs = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
834
  for ( recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
835
    streamptr->tsteps[1].recIDs[recID] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
836

Uwe Schulzweida's avatar
Uwe Schulzweida committed
837
838
  for ( recID = 0; recID < nrecords; recID++ )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
839
      varID = streamptr->tsteps[0].records[recID].varID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
840
841
      streamptr->tsteps[tsID].records[recID].position =	streamptr->tsteps[0].records[recID].position;
      streamptr->tsteps[tsID].records[recID].size     =	streamptr->tsteps[0].records[recID].size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
842
843
    }

844
845
  int nrecs_scanned = nrecords;
  int rindex = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
846
  while ( TRUE )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
847
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
848
849
      if ( rindex > nrecords ) break;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
850
851
852
853
      recsize = gribGetSize(fileID);
      recpos  = fileGetPos(fileID);
      if ( recsize == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
854
	  streamptr->ntsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
855
856
	  break;
	}
857
      if ( recsize > buffersize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
858
	{
859
	  buffersize = recsize;
860
	  gribbuffer = Realloc(gribbuffer, buffersize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
861
862
	}

863
      readsize = recsize;
864
      rstatus = gribRead(fileID, (unsigned char *)gribbuffer, &readsize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
865
      if ( rstatus ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
866

867
      if ( gribGetZip(recsize, (unsigned char *)gribbuffer, &unzipsize) > 0 )
868
869
	{
	  unzipsize += 100; /* need 0 to 1 bytes for rounding of bds */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
870
	  if ( buffersize < unzipsize )
871
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
872
	      buffersize = unzipsize;
873
	      gribbuffer = Realloc(gribbuffer, buffersize);
874
875
876
	    }
	}

877
      cgribexDecodeHeader(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, fsec4,
878
			  (int *) gribbuffer, (int)recsize, &lmv, &iret);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
879

880
881
      nrecs_scanned++;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
882
      param = cdiEncodeParam(ISEC1_Parameter, ISEC1_CodeTable, 255);
883
884
      cdiParamToString(param, paramstr, sizeof(paramstr));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
885
      cgribexGetLevel(isec1, &leveltype, &level1, &level2);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900

      gribDateTime(isec1, &vdate, &vtime);

      if ( rindex == 0 )
	{
	  if ( taxisInqType(taxisID) == TAXIS_RELATIVE )
	    {
	      taxis->type  = TAXIS_RELATIVE;
	      taxis->rdate = gribRefDate(isec1);
	      taxis->rtime = gribRefTime(isec1);
	    }
	  else
	    {
	      taxis->type  = TAXIS_ABSOLUTE;
	    }
901
	  taxis->unit  = cgribexGetTimeUnit(isec1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
902
903
	  taxis->vdate = vdate;
	  taxis->vtime = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
904
905
906

	  datetime0.date = vdate;
	  datetime0.time = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
907
908
	}

909
910
      tsteptype = cgribexGetTsteptype(ISEC1_TimeRange);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
911
912
913
      if ( ISEC1_AvgNum )
	{
	  if (  taxis->numavg && warn_numavg &&
914
        	(taxis->numavg != ISEC1_AvgNum) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
915
916
	    {
	  /*
917
	      Warning("Changing numavg from %d to %d not supported!", taxis->numavg, ISEC1_AvgNum);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
918
	  */
919
	      warn_numavg = false;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
920
921
922
923
924
925
926
	    }
	  else
	    {
	      taxis->numavg = ISEC1_AvgNum;
	    }
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
927
928
      datetime.date  = vdate;
      datetime.time  = vtime;
929

Uwe Schulzweida's avatar
Uwe Schulzweida committed
930
      compvar_t compVar = cgribexVarSet(param, level1, level2, leveltype, ISEC1_TimeRange);
931

Uwe Schulzweida's avatar
Uwe Schulzweida committed
932
933
      for ( recID = 0; recID < nrecords; recID++ )
	{
934
	  if ( cgribexVarCompare(compVar, streamptr->tsteps[tsID].records[recID], 0)