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

#include <stdio.h>

Uwe Schulzweida's avatar
Uwe Schulzweida committed
7
#include "dmemory.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
8
9
10
11
12
13
#include "cdi.h"
#include "stream_int.h"
#include "file.h"
#include "varscan.h"
#include "datetime.h"
#include "vlist.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
14
#include "stream_grb.h"
15
#include "calendar.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
16
17
18
19


#if  defined  (HAVE_LIBGRIB_API)
#  include "cgribex.h"      /* gribGetSize, gribRead, gribGetZip */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
20
#  include "gribapi.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
21
22
23
24
25
26
27
#  include "grib_api.h"
#endif


extern int cdiInventoryMode;

typedef struct {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
28
  int param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
29
30
31
  int level1;
  int level2;
  int ltype;
32
} compvar2_t;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
33
34
35
36


#if  defined  (HAVE_LIBGRIB_API)
static
37
int gribapiGetGridType(grib_handle *gh)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
38
{
39
  int gridtype = GRID_GENERIC;
40
  int gribgridtype = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
41
42
  long lpar;

43
    {
44
45
46
47
      int status;
      status = grib_get_long(gh, "gridDefinitionTemplateNumber", &lpar);

      if ( status ==  0 ) gribgridtype = (int) lpar;
48
49
50

      switch (gribgridtype)
	{
51
	case  GRIB2_GTYPE_LATLON:     { GRIB_CHECK(grib_get_long(gh, "Ni", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
52
	                                if ( lpar == (long) GRIB_MISSING_LONG ) break;
53
                                      }
54
55
	case  GRIB2_GTYPE_LATLON_ROT: { gridtype = GRID_LONLAT;    break; }
	case  GRIB2_GTYPE_LCC:        { gridtype = GRID_LCC;       break; }
56
	case  GRIB2_GTYPE_GAUSSIAN:   { GRIB_CHECK(grib_get_long(gh, "Ni", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
57
	                                if ( lpar == (long) GRIB_MISSING_LONG )
58
59
60
61
62
					  gridtype = GRID_GAUSSIAN_REDUCED;
					else
					  gridtype = GRID_GAUSSIAN;
					break;
                                      }
63
64
65
	case  GRIB2_GTYPE_SPECTRAL:   { gridtype = GRID_SPECTRAL;  break; }
	case  GRIB2_GTYPE_GME:        { gridtype = GRID_GME;       break; }
	case  GRIB2_GTYPE_NUMBER:     { gridtype = GRID_REFERENCE; break; }
66
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
67
68
69
70
71
    }

  return (gridtype);
}

72
static
73
int gribapiGetIsRotated(grib_handle *gh)
74
75
{
  int isRotated = 0;
76
  int gribgridtype = -1;
77
  long lpar;
78
  int status;
79

80
  status = grib_get_long(gh, "gridDefinitionTemplateNumber", &lpar);
81

82
83
84
  if ( status ==  0 ) gribgridtype = (int) lpar;

  if ( gribgridtype == GRIB2_GTYPE_LATLON_ROT ) isRotated = 1;
85
86
87
88

  return (isRotated);
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
89
static
90
int gribapiGetZaxisType(long editionNumber, int grib_ltype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
91
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
92
93
94
  int zaxistype = ZAXIS_GENERIC;

  if ( editionNumber <= 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
95
    {
96
      zaxistype = grib1ltypeToZaxisType(grib_ltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
97
98
99
    }
  else
    {
100
      zaxistype = grib2ltypeToZaxisType(grib_ltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
101
102
103
104
105
    }

  return (zaxistype);
}

106
static
107
int getTimeunits(long unitsOfTime)
108
109
{
  int timeunits = -1;
110
111
112
113
114
115
116
117
118
119
120
121
122

  switch (unitsOfTime)
    {
    case 13:  timeunits = TUNIT_SECOND;  break;
    case  0:  timeunits = TUNIT_MINUTE;  break;
    case  1:  timeunits = TUNIT_HOUR;    break;
    case 10:  timeunits = TUNIT_3HOURS;  break;
    case 11:  timeunits = TUNIT_6HOURS;  break;
    case 12:  timeunits = TUNIT_12HOURS; break;
    case  2:  timeunits = TUNIT_DAY;     break;
    default:  timeunits = TUNIT_HOUR;    break;
    }

123
124
125
126
127
128
129
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
  return (timeunits);
}

static
double timeunit_factor(int tu1, int tu2)
{
  double factor = 1;

  if ( tu2 == TUNIT_HOUR )
    {
      switch (tu1)
        {
        case TUNIT_SECOND:  factor = 3600;   break;
        case TUNIT_MINUTE:  factor = 60;     break;
        case TUNIT_HOUR:    factor = 1;      break;
        case TUNIT_3HOURS:  factor = 1./3;   break;
        case TUNIT_6HOURS:  factor = 1./6;   break;
        case TUNIT_12HOURS: factor = 1./12;  break;
        case TUNIT_DAY:     factor = 1./24;  break;
        }
    }

  return (factor);
}

static
int gribapiGetEndStep(grib_handle *gh, int startStep, int timeunits)
{
  int endStep = startStep;
  int timeunits2;
  long unitsOfTime;
  long lpar;

  GRIB_CHECK(grib_get_long(gh, "stepUnits", &unitsOfTime), 0);

  timeunits2 = getTimeunits(unitsOfTime);

  GRIB_CHECK(grib_get_long(gh, "endStep", &lpar), 0);

  endStep = (int)  ((lpar * timeunit_factor(timeunits, timeunits2)) + 0.5);

  return (endStep);
}

static
int gribapiGetTimeUnits(grib_handle *gh)
{
  int timeunits = -1;
  long unitsOfTime;
  // size_t len = 8;
  //char stepunits[8];
  //static int lprint = TRUE;

  GRIB_CHECK(grib_get_long(gh, "indicatorOfUnitOfTimeRange", &unitsOfTime), 0);

  timeunits = getTimeunits(unitsOfTime);

180
  /*
181
182
  GRIB_CHECK(grib_get_string(gh, "stepUnits", stepunits, &len), 0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
183
184
185
186
187
188
189
190
191
192
193
  len--;

  if      ( memcmp(stepunits, "s",   len) == 0 ) timeunits = TUNIT_SECOND;
  else if ( memcmp(stepunits, "m",   len) == 0 ) timeunits = TUNIT_MINUTE;
  else if ( memcmp(stepunits, "h",   len) == 0 ) timeunits = TUNIT_HOUR;
  else if ( memcmp(stepunits, "3h",  len) == 0 ) timeunits = TUNIT_3HOURS;
  else if ( memcmp(stepunits, "6h",  len) == 0 ) timeunits = TUNIT_6HOURS;
  else if ( memcmp(stepunits, "12h", len) == 0 ) timeunits = TUNIT_12HOURS;
  else if ( memcmp(stepunits, "D",   len) == 0 ) timeunits = TUNIT_DAY;
  else if ( memcmp(stepunits, "M",   len) == 0 ) timeunits = TUNIT_MONTH;
  else if ( memcmp(stepunits, "Y",   len) == 0 ) timeunits = TUNIT_YEAR;
194
  else if ( lprint )
195
    {
196
      Message("Step units >%s< unsupported!", stepunits);
197
      lprint = FALSE;
198
    }
199
  */
200
201
202
203
204
205
206

  return (timeunits);
}

static
int gribapiTimeIsFC(grib_handle *gh)
{
207
  long editionNumber;
208
209
  int isFC = TRUE;

210
  GRIB_CHECK(grib_get_long(gh, "editionNumber", &editionNumber), 0);
211

212
  if ( editionNumber > 1 )
213
214
215
216
217
218
219
    {
      long sigofrtime;

      GRIB_CHECK(grib_get_long(gh, "significanceOfReferenceTime", &sigofrtime), 0);

      if ( sigofrtime == 3 ) isFC = FALSE;
    }
220
221
222
223
224
225
226

  return (isFC);
}

static
int gribapiGetTsteptype(grib_handle *gh)
{
227
  int tsteptype = TSTEP_INSTANT;
228
229
230
231
  static int lprint = TRUE;

  if ( gribapiTimeIsFC(gh) )
    {
232
      int status;
233
234
      size_t len = 256;
      char stepType[256];
235

236
237
      status = grib_get_string(gh, "stepType", stepType, &len);
      if ( status == 0 && len > 1 && len < 256 )
238
	{
239
	  if      ( strncmp("instant", stepType, len) == 0 ) tsteptype = TSTEP_INSTANT;
240
	  else if ( strncmp("avg",     stepType, len) == 0 ) tsteptype = TSTEP_AVG;
241
	  else if ( strncmp("accum",   stepType, len) == 0 ) tsteptype = TSTEP_ACCUM;
242
243
244
245
246
247
248
	  else if ( strncmp("max",     stepType, len) == 0 ) tsteptype = TSTEP_MAX;
	  else if ( strncmp("min",     stepType, len) == 0 ) tsteptype = TSTEP_MIN;
	  else if ( strncmp("diff",    stepType, len) == 0 ) tsteptype = TSTEP_DIFF;
	  else if ( strncmp("rms",     stepType, len) == 0 ) tsteptype = TSTEP_RMS;
	  else if ( strncmp("sd",      stepType, len) == 0 ) tsteptype = TSTEP_SD;
	  else if ( strncmp("cov",     stepType, len) == 0 ) tsteptype = TSTEP_COV;
	  else if ( strncmp("ratio",   stepType, len) == 0 ) tsteptype = TSTEP_RATIO;
249
	  else if ( lprint )
250
	    {
251
252
	      Message("stepType %s unsupported, set to instant!", stepType);
	      lprint = FALSE;
253
	    }
254
255

	  // printf("stepType: %s %ld %d\n", stepType, len, tsteptype);
256
257
258
259
260
261
	}
    }

  return (tsteptype);
}

262
static
263
int gribapiGetValidityDateTime(grib_handle *gh, int *vdate, int *vtime)
264
{
265
  int rdate, rtime;
266
267
268
  int timeUnits, startStep, endStep;
  int tstepRange = 0;
  int range;
269
  long lpar;
270
271
272
273
274
  long sigofrtime = 3;
  long editionNumber;

  GRIB_CHECK(grib_get_long(gh, "editionNumber", &editionNumber), 0);

275
  if ( editionNumber > 1 )
276
277
278
    {
      GRIB_CHECK(grib_get_long(gh, "significanceOfReferenceTime", &sigofrtime), 0);
    }
279
280
281
282
283
284
285
286
287
288

  if ( sigofrtime == 3 )
    {
      GRIB_CHECK(grib_get_long(gh, "dataDate", &lpar), 0);
      *vdate = (int) lpar;
      GRIB_CHECK(grib_get_long(gh, "dataTime", &lpar), 0);
      *vtime = (int) lpar*100;
    }
  else
    {
289
290
291
      GRIB_CHECK(grib_get_long(gh, "dataDate", &lpar), 0);
      rdate = (int) lpar;
      GRIB_CHECK(grib_get_long(gh, "dataTime", &lpar), 0);
292
      rtime = (int) lpar*100;
293
      GRIB_CHECK(grib_get_long(gh, "forecastTime", &lpar), 0);
294
      startStep = (int) lpar;
295
      timeUnits = gribapiGetTimeUnits(gh);
296
      endStep = gribapiGetEndStep(gh, startStep, timeUnits);
297
298
299
300
301
302
303
304
305

      range = endStep - startStep;

      if ( range > 0 )
	{
	  if ( startStep == 0 ) tstepRange = -1;
	  else                  tstepRange =  1;
	}

306
307
308
309
      {
	static int lprint = TRUE;
	extern int grib_calendar;
	int ryear, rmonth, rday, rhour, rminute, rsecond;
310
	int time_period = endStep;
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
	int julday, secofday, addsec;

	cdiDecodeDate(rdate, &ryear, &rmonth, &rday);
	cdiDecodeTime(rtime, &rhour, &rminute, &rsecond);

	encode_caldaysec(grib_calendar, ryear, rmonth, rday, rhour, rminute, rsecond, &julday, &secofday);

	addsec = 0;
	switch ( timeUnits )
	  {
	  case TUNIT_SECOND:  addsec =         time_period; break;
	  case TUNIT_MINUTE:  addsec =    60 * time_period; break;
	  case TUNIT_HOUR:    addsec =  3600 * time_period; break;
	  case TUNIT_3HOURS:  addsec = 10800 * time_period; break;
	  case TUNIT_6HOURS:  addsec = 21600 * time_period; break;
	  case TUNIT_12HOURS: addsec = 43200 * time_period; break;
	  case TUNIT_DAY:     addsec = 86400 * time_period; break;
	  default:
	    if ( lprint )
	      {
	        Warning("Time unit %d unsupported", timeUnits);
		lprint = FALSE;
	      }
	  }

	julday_add_seconds(addsec, &julday, &secofday);

	decode_caldaysec(grib_calendar, julday, secofday, &ryear, &rmonth, &rday, &rhour, &rminute, &rsecond);
	/*
	  printf("new %d/%d/%d %d:%d\n", ryear, rmonth, rday, rhour, rminute);
	*/
	*vdate = cdiEncodeDate(ryear, rmonth, rday);
	*vtime = cdiEncodeTime(rhour, rminute, rsecond);
      }
345
    }
346
347

  return (tstepRange);
348
349
}

350
static
351
void gribapiGetGrid(grib_handle *gh, grid_t *grid)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
352
{
353
  long editionNumber;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
354
  int gridtype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
355
  size_t datasize;
356
357
  long numberOfPoints;
  long lpar;
358
359
360

  GRIB_CHECK(grib_get_long(gh, "editionNumber", &editionNumber), 0);

361
  gridtype = gribapiGetGridType(gh);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
362
363
364
365
366
367
368
369
  /*
  if ( streamptr->unreduced && gridtype == GRID_GAUSSIAN_REDUCED )
    {
      gridtype = GRID_GAUSSIAN;
      ISEC2_NumLon = 2*ISEC2_NumLat;
      ISEC4_NumValues = ISEC2_NumLon*ISEC2_NumLat;
    }
  */
370
  memset(grid, 0, sizeof(grid_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
371

Uwe Schulzweida's avatar
Uwe Schulzweida committed
372
373
374
  GRIB_CHECK(grib_get_size(gh, "values", &datasize), 0);
  GRIB_CHECK(grib_get_long(gh, "numberOfPoints", &numberOfPoints), 0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
375
376
377
378
379
380
381
382
383
384
385
386
  switch (gridtype)
    {
    case GRID_LONLAT:
    case GRID_GAUSSIAN:
      {
	int nlon, nlat;

	GRIB_CHECK(grib_get_long(gh, "Ni", &lpar), 0);
	nlon = lpar;
	GRIB_CHECK(grib_get_long(gh, "Nj", &lpar), 0);
	nlat = lpar;

387
388
389
390
391
392
	if ( gridtype == GRID_GAUSSIAN )
          {
            GRIB_CHECK(grib_get_long(gh, "numberOfParallelsBetweenAPoleAndTheEquator", &lpar), 0);
            grid->np = lpar;
          }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
393
	if ( numberOfPoints != nlon*nlat )
394
	  Error("numberOfPoints (%d) and gridSize (%d) differ!", (int)numberOfPoints, nlon*nlat);
395

396
397
398
399
400
401
402
403
404
405
406
	grid->size  = numberOfPoints;
	grid->xsize = nlon;
	grid->ysize = nlat;
	grid->xinc  = 0;
	grid->yinc  = 0;
	grid->xdef  = 0;
	GRIB_CHECK(grib_get_double(gh, "longitudeOfFirstGridPointInDegrees", &grid->xfirst), 0);
	GRIB_CHECK(grib_get_double(gh, "longitudeOfLastGridPointInDegrees",  &grid->xlast), 0);
	GRIB_CHECK(grib_get_double(gh, "latitudeOfFirstGridPointInDegrees",  &grid->yfirst), 0);
	GRIB_CHECK(grib_get_double(gh, "latitudeOfLastGridPointInDegrees",   &grid->ylast), 0);
	GRIB_CHECK(grib_get_double(gh, "iDirectionIncrementInDegrees", &grid->xinc), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
407
	if ( gridtype == GRID_LONLAT )
408
	  GRIB_CHECK(grib_get_double(gh, "jDirectionIncrementInDegrees", &grid->yinc), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
409

410
	if ( IS_EQUAL(grid->xinc, GRIB_MISSING_DOUBLE) ) grid->xinc = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
411

412
	/* if ( IS_NOT_EQUAL(grid->xfirst, 0) || IS_NOT_EQUAL(grid->xlast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
413
	  {
414
	    if ( grid->xsize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
415
	      {
416
		if ( (grid->xfirst >= grid->xlast) && (grid->xfirst >= 180) ) grid->xfirst -= 360;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
417

418
419
420
		if ( editionNumber <= 1 )
		  {
		    /* correct xinc if necessary */
421
		    if ( IS_EQUAL(grid->xfirst, 0) && grid->xlast > 354 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
422
		      {
423
			double xinc = 360. / grid->xsize;
424

425
			if ( fabs(grid->xinc-xinc) > 0.0 )
426
			  {
427
428
			    grid->xinc = xinc;
			    if ( CDI_Debug ) Message("set xinc to %g", grid->xinc);
429
			  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
430
431
432
		      }
		  }
	      }
433
	    grid->xdef = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
434
	  }
435
	grid->ydef = 0;
436
        /* if ( IS_NOT_EQUAL(grid->yfirst, 0) || IS_NOT_EQUAL(grid->ylast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
437
	  {
438
	    if ( grid->ysize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
439
	      {
440
441
442
		if ( editionNumber <= 1 )
		  {
		  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
443
	      }
444
	    grid->ydef = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
445
446
447
448
449
	  }
	break;
      }
    case GRID_GAUSSIAN_REDUCED:
      {
450
451
452
453
	int nlat, i;
	size_t dummy;
	long *pl;

454
455
456
        GRIB_CHECK(grib_get_long(gh, "numberOfParallelsBetweenAPoleAndTheEquator", &lpar), 0);
        grid->np = lpar;

457
458
459
	GRIB_CHECK(grib_get_long(gh, "Nj", &lpar), 0);
	nlat = lpar;

460
	grid->size   = numberOfPoints;
461

462
        grid->rowlon = (int *) malloc(nlat*sizeof(int));
463
464
465
        pl          = (long *) malloc(nlat*sizeof(long));
	dummy       = nlat;
	GRIB_CHECK(grib_get_long_array(gh, "pl", pl, &dummy), 0);
466
	for ( i = 0; i < nlat; ++i ) grid->rowlon[i] = pl[i];
467
468
	free(pl);

469
470
471
472
473
474
475
476
477
	grid->ysize  = nlat;
	grid->xinc   = 0;
	grid->yinc   = 0;
	grid->xdef   = 0;
	GRIB_CHECK(grib_get_double(gh, "longitudeOfFirstGridPointInDegrees", &grid->xfirst), 0);
	GRIB_CHECK(grib_get_double(gh, "longitudeOfLastGridPointInDegrees",  &grid->xlast), 0);
	GRIB_CHECK(grib_get_double(gh, "latitudeOfFirstGridPointInDegrees",  &grid->yfirst), 0);
	GRIB_CHECK(grib_get_double(gh, "latitudeOfLastGridPointInDegrees",   &grid->ylast), 0);
	GRIB_CHECK(grib_get_double(gh, "iDirectionIncrementInDegrees", &grid->xinc), 0);
478

479
	if ( IS_EQUAL(grid->xinc, GRIB_MISSING_DOUBLE) ) grid->xinc = 0;
480

481
	/* if ( IS_NOT_EQUAL(grid->xfirst, 0) || IS_NOT_EQUAL(grid->xlast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
482
	  {
483
	    if ( grid->xsize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
484
	      {
485
		if ( (grid->xfirst > grid->xlast) && (grid->xfirst >= 180) ) grid->xfirst -= 360;
486
487
488
489

		if ( editionNumber <= 1 )
		  {
		    /* correct xinc if necessary */
490
		    if ( IS_EQUAL(grid->xfirst, 0) && grid->xlast > 354 )
491
		      {
492
			double xinc = 360. / grid->xsize;
493

494
			if ( fabs(grid->xinc-xinc) > 0.0 )
495
			  {
496
497
			    grid->xinc = xinc;
			    if ( CDI_Debug ) Message("set xinc to %g", grid->xinc);
498
499
500
			  }
		      }
		  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
501
	      }
502
	    grid->xdef = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
503
	  }
504
505
	grid->ydef  = 0;
        /* if ( IS_NOT_EQUAL(grid->yfirst, 0) || IS_NOT_EQUAL(grid->ylast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
506
	  {
507
	    if ( grid->ysize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
508
	      {
509
510
511
		if ( editionNumber <= 1 )
		  {
		  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
512
	      }
513
	    grid->ydef = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
514
515
516
	  }
	break;
      }
517
      /*
Uwe Schulzweida's avatar
Uwe Schulzweida committed
518
519
520
    case GRID_LCC:
      {
	if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
521
	  Error("numberOfPoints (%d) and gridSize (%d) differ!",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
522
523
		ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);

524
525
526
	grid->size  = ISEC4_NumValues;
	grid->xsize = ISEC2_NumLon;
	grid->ysize = ISEC2_NumLat;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
527

528
529
530
531
532
533
534
535
536
	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;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
537

538
	grid->xdef   = 0;
539
	grid->ydef   = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
540
541
542

	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
543
      */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
544
545
    case GRID_SPECTRAL:
      {
546
547
	size_t len = 256;
	char typeOfPacking[256];
548
	GRIB_CHECK(grib_get_string(gh, "packingType", typeOfPacking, &len), 0);
549
550
	grid->lcomplex = 0;
	if ( strncmp(typeOfPacking, "spectral_complex", len) == 0 ) grid->lcomplex = 1;
551

552
	grid->size  = datasize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
553
	GRIB_CHECK(grib_get_long(gh, "J", &lpar), 0);
554
	grid->trunc = lpar;
555

Uwe Schulzweida's avatar
Uwe Schulzweida committed
556
557
558
559
	break;
      }
    case GRID_GME:
      {
560
561
562
563
564
	grid->size  = numberOfPoints;
	if ( grib_get_long(gh, "nd", &lpar) == 0 ) grid->nd  = lpar;
	if ( grib_get_long(gh, "Ni", &lpar) == 0 ) grid->ni  = lpar;
	if ( grib_get_long(gh, "n2", &lpar) == 0 ) grid->ni2 = lpar;
	if ( grib_get_long(gh, "n3", &lpar) == 0 ) grid->ni3 = lpar;
565
566
567

	break;
      }
568
    case GRID_REFERENCE:
569
      {
570
        char uuid[17];
571
572
573
574
	char reference_link[8192];
	size_t len = sizeof(reference_link);
	reference_link[0] = 0;

575
	grid->size  = numberOfPoints;
576
577
	if ( grib_get_long(gh, "numberOfGridUsed", &lpar) == 0 )
	  {
578
579
	    grid->number   = lpar;
	    if ( grib_get_long(gh, "numberOfGridInReference", &lpar) == 0 ) grid->position = lpar;
580
581
	    if ( grib_get_string(gh, "gridDescriptionFile", reference_link, &len) == 0 )
	      {
582
		if ( strncmp(reference_link, "file://", 7) == 0 )
583
		  grid->reference = strdupx(reference_link);
584
	      }
585
586
587
588
589
            len = (size_t) 16;
            if ( grib_get_string(gh, "uuidOfHGrid", uuid, &len) == 0)
              {
                strncpy(grid->uuid, uuid, 16);
              }
590
	  }
591

Uwe Schulzweida's avatar
Uwe Schulzweida committed
592
593
594
595
	break;
      }
    case GRID_GENERIC:
      {
596
597
598
599
600
	int nlon = 0, nlat = 0;

	if ( grib_get_long(gh, "Ni", &lpar) == 0 ) nlon = lpar;
	if ( grib_get_long(gh, "Nj", &lpar) == 0 ) nlat = lpar;

601
	grid->size  = numberOfPoints;
602

603
	if ( nlon > 0 && nlat > 0 && nlon*nlat == grid->size )
604
	  {
605
606
	    grid->xsize = nlon;
	    grid->ysize = nlat;
607
608
609
	  }
	else
	  {
610
611
	    grid->xsize = 0;
	    grid->ysize = 0;
612
613
	  }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
614
615
616
617
	break;
      }
    default:
      {
618
	Error("Unsupported grid type: %s", gridNamePtr(gridtype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
619
620
621
622
	break;
      }
    }

623
  grid->isRotated = FALSE;
624
  if ( gribapiGetIsRotated(gh) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
625
    {
626
627
628
629
      grid->isRotated = TRUE;
      GRIB_CHECK(grib_get_double(gh, "latitudeOfSouthernPoleInDegrees",  &grid->ypole), 0);
      GRIB_CHECK(grib_get_double(gh, "longitudeOfSouthernPoleInDegrees", &grid->xpole), 0);
      GRIB_CHECK(grib_get_double(gh, "angleOfRotation", &grid->angle), 0);
630
      /* change from south to north pole */
631
632
633
634
635
636
637
638
639
640
      grid->ypole = -grid->ypole;
      grid->xpole =  grid->xpole - 180;
    }

  grid->xvals = NULL;
  grid->yvals = NULL;
  grid->type  = gridtype;
}

static
641
void grib1GetLevel(grib_handle *gh, int *leveltype, int *lbounds, int *level1, int *level2)
642
{
643
644
645
  int status;
  long lpar;
  double dlevel;
646

647
648
649
650
  *leveltype = 0;
  *lbounds = 0;
  *level1  = 0;
  *level2  = 0;
651

652
653
  status = grib_get_long(gh, "indicatorOfTypeOfLevel", &lpar);
  if ( status == 0 )
654
    {
655
      *leveltype = (int) lpar;
656

657
      switch (*leveltype)
658
659
660
661
	{
	case GRIB1_LTYPE_SIGMA_LAYER:
	case GRIB1_LTYPE_HYBRID_LAYER:
	case GRIB1_LTYPE_LANDDEPTH_LAYER:
662
	  { *lbounds = 1; break; }
663
	}
664

665
666
667
668
669
670
      if ( *lbounds == 0 )
	{
	  GRIB_CHECK(grib_get_double(gh, "level", &dlevel), 0);
	  if ( *leveltype == 100 ) dlevel *= 100;
	  if ( dlevel < -2.e9 || dlevel > 2.e9 ) dlevel = 0;
	  if ( *leveltype == 99 ) *leveltype = 100;
671

672
673
674
675
676
677
678
679
680
681
	  *level1 = (int) dlevel;
	  *level2 = 0;
	}
      else
	{
	  GRIB_CHECK(grib_get_long(gh, "topLevel", &lpar), 0);
	  *level1 = lpar;
	  GRIB_CHECK(grib_get_long(gh, "bottomLevel", &lpar), 0);
	  *level2 = lpar;
	}
682
    }
683
684
685
}

static
686
void grib2GetLevel(grib_handle *gh, int *leveltype, int *lbounds, int *level1, int *level2)
687
{
688
689
690
691
  int status;
  int leveltype2 = -1;
  long lpar;
  long factor;
692
693
  double dlevel;

694
695
696
697
698
699
700
  *leveltype = 0;
  *lbounds = 0;
  *level1  = 0;
  *level2  = 0;

  status = grib_get_long(gh, "typeOfFirstFixedSurface", &lpar);
  if ( status == 0 )
701
    {
702
      *leveltype = (int) lpar;
703

704
705
      status = grib_get_long(gh, "typeOfSecondFixedSurface", &lpar);
      if ( status == 0 ) leveltype2 = lpar;
706

707
      if ( *leveltype == leveltype2 && *leveltype != 255 ) *lbounds = 1;
708

709
      if ( *lbounds == 0 )
710
	{
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
	  if ( *leveltype == GRIB2_LTYPE_LANDDEPTH )
	    {
	      GRIB_CHECK(grib_get_long(gh, "scaleFactorOfFirstFixedSurface", &factor), 0);
	      GRIB_CHECK(grib_get_double(gh, "scaledValueOfFirstFixedSurface", &dlevel), 0);
	      if      ( factor == 0 ) dlevel *= 100;   //  m to cm
	      else if ( factor == 1 ) dlevel *=  10;   // dm to cm
	      else if ( factor == 3 ) dlevel *=   0.1; // mm to cm
	    }
	  else
	    {
	      GRIB_CHECK(grib_get_double(gh, "level", &dlevel), 0);
	      if ( *leveltype == GRIB2_LTYPE_ISOBARIC ) dlevel *= 100;
	      if ( dlevel < -2.e9 || dlevel > 2.e9 ) dlevel = 0;
	      if ( *leveltype == 99 ) *leveltype = 100;
	    }

727
	  *level1 = (int) dlevel;
728
	  *level2 = 0;
729
730
731
	}
      else
	{
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
	  if ( *leveltype == GRIB2_LTYPE_LANDDEPTH )
	    {
	      GRIB_CHECK(grib_get_long(gh, "scaleFactorOfFirstFixedSurface", &factor), 0);
	      GRIB_CHECK(grib_get_double(gh, "scaledValueOfFirstFixedSurface", &dlevel), 0);
	      if      ( factor == 0 ) dlevel *= 100;   //  m to cm
	      else if ( factor == 1 ) dlevel *=  10;   // dm to cm
	      else if ( factor == 3 ) dlevel *=   0.1; // mm to cm
	      *level1 = (int) dlevel;
	      GRIB_CHECK(grib_get_long(gh, "scaleFactorOfSecondFixedSurface", &factor), 0);
	      GRIB_CHECK(grib_get_double(gh, "scaledValueOfSecondFixedSurface", &dlevel), 0);
	      if      ( factor == 0 ) dlevel *= 100;   //  m to cm
	      else if ( factor == 1 ) dlevel *=  10;   // dm to cm
	      else if ( factor == 3 ) dlevel *=   0.1; // mm to cm
	      *level2 = (int) dlevel;
	    }
	  else
	    {
	      GRIB_CHECK(grib_get_long(gh, "topLevel", &lpar), 0);
	      *level1 = lpar;
	      GRIB_CHECK(grib_get_long(gh, "bottomLevel", &lpar), 0);
	      *level2 = lpar;
	    }
754
755
	}
    }
756
757
758
}

static
759
void gribapiAddRecord(stream_t * streamptr, int param, grib_handle *gh,
760
		      long recsize, off_t position, int datatype, int comptype)
761
762
763
764
765
766
{
  long editionNumber;
  int zaxistype;
  int gridID = CDI_UNDEFID, varID;
  int levelID = 0;
  int tsID, recID;
767
  int level1 = 0, level2 = 0;
768
769
770
771
772
773
774
775
776
777
  int numavg;
  int tsteptype;
  int lbounds = 0;
  record_t *record;
  grid_t grid;
  int vlistID;
  int leveltype;
  long lpar;
  int status;
  char name[256], longname[256], units[256];
778
  size_t vlen;
779
  long ens_index = 0, ens_count = 0, ens_forecast_type = 0;
780

Uwe Schulzweida's avatar
Uwe Schulzweida committed
781
  vlistID = streamptr->vlistID;
782
  tsID    = streamptr->curTsID;
783
  recID   = recordNewEntry(streamptr, tsID);
784
785
786
787
788
789
790
791
792
  record  = &streamptr->tsteps[tsID].records[recID];

  tsteptype = gribapiGetTsteptype(gh);
  // numavg  = ISEC1_AvgNum;
  numavg  = 0;

  GRIB_CHECK(grib_get_long(gh, "editionNumber", &editionNumber), 0);

  if ( editionNumber <= 1 )
793
    grib1GetLevel(gh, &leveltype, &lbounds, &level1, &level2);
794
  else
795
    grib2GetLevel(gh, &leveltype, &lbounds, &level1, &level2);
796

797
798
799
800
801
802
803
804
805
806
  // fprintf(stderr, "param %d %d %d %d\n", param, level1, level2, leveltype);

  (*record).size     = recsize;
  (*record).position = position;
  (*record).param    = param;
  (*record).ilevel   = level1;
  (*record).ilevel2  = level2;
  (*record).ltype    = leveltype;

  gribapiGetGrid(gh, &grid);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
807
808
809

  gridID = varDefGrid(vlistID, grid, 0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
810
  zaxistype = gribapiGetZaxisType(editionNumber, leveltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
811

812
  switch (zaxistype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
813
    {
814
815
816
817
818
819
    case ZAXIS_HYBRID:
    case ZAXIS_HYBRID_HALF:
      {
        int vctsize;
        size_t dummy;
        double *vctptr;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
820

821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
        GRIB_CHECK(grib_get_long(gh, "NV", &lpar), 0);
        vctsize = lpar;
        if ( vctsize > 0 )
          {
            vctptr = (double *) malloc(vctsize*sizeof(double));
            dummy = vctsize;
            GRIB_CHECK(grib_get_double_array(gh, "pv", vctptr, &dummy), 0);
            varDefVCT(vctsize, vctptr);
            free(vctptr);
          }
        break;
      }
    case ZAXIS_REFERENCE:
      {
        size_t len;
        char uuid[17];
        long nlev, nvgrid;

        GRIB_CHECK(grib_get_long(gh, "NV", &lpar), 0);
840
        if ( lpar != 3 )
841
842
843
844
845
846
847
848
849
850
851
852
          {
            fprintf(stderr, "Warning ...\n");
          }
        GRIB_CHECK(grib_get_long(gh, "nlev", &nlev), 0);
        GRIB_CHECK(grib_get_long(gh, "numberOfVGridUsed", &nvgrid), 0);

        len = (size_t) 16;
        uuid[16] = 0;
        GRIB_CHECK(grib_get_string(gh, "uuidOfVGrid", uuid, &len), 0);
        varDefZAxisReference((int) nlev, (int) nvgrid, uuid);
        break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
853
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
854

855
  // if ( datatype > 32 ) datatype = DATATYPE_PACK32;
856
  if ( datatype <  0 ) datatype = DATATYPE_PACK;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
857

858
859
860
861
  name[0] = 0;
  longname[0] = 0;
  units[0] = 0;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
862
863
  vlen = 256;
  GRIB_CHECK(grib_get_string(gh, "shortName", name, &vlen), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
864
865
866
867
868
869
870
871
872
873
874
875
  if      ( vlen == 8 && memcmp(name, "unknown", vlen) == 0 ) name[0] = 0;
  else if ( vlen == 2 && memcmp(name, "~", vlen)       == 0 ) name[0] = 0;

  if ( name[0] != 0 )
    {
      vlen = 256;
      GRIB_CHECK(grib_get_string(gh, "name", longname, &vlen), 0);
      if ( vlen == 8 && memcmp(longname, "unknown", vlen) == 0 ) longname[0] = 0;
      vlen = 256;
      GRIB_CHECK(grib_get_string(gh, "units", units, &vlen), 0);
      if ( vlen == 8 && memcmp(units, "unknown", vlen) == 0 ) units[0] = 0;
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
876
  // fprintf(stderr, "param %d name %s %s %s\n", param, name, longname, units); 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
877

Uwe Schulzweida's avatar
Uwe Schulzweida committed
878
  varAddRecord(recID, param, gridID, zaxistype, lbounds, level1, level2,
879
	       datatype, &varID, &levelID, tsteptype, numavg, leveltype,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
880
	       name, longname, units);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
881
882
883
884

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

885
  varDefCompType(varID, comptype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
886

887
888
  /*
    Get the ensemble Info from the grib-2 Tables and update the intermediate datastructure.
889
    Further update to the "vlist" is handled in the same way as for GRIB-1 by "cdi_generate_vars"
890
  */
891
892
893
894
895
896
897
898
899
  {
    int status;
    status = grib_get_long(gh, "typeOfEnsembleForecast", &ens_forecast_type );
    if ( status == 0 )
      {
        GRIB_CHECK(grib_get_long(gh, "numberOfForecastsInEnsemble", &ens_count ), 0);
        GRIB_CHECK(grib_get_long(gh, "perturbationNumber", &ens_index ), 0);
      }
  }
900
901
902
903
  if ( ens_index > 0 )
    varDefEnsembleInfo(varID, (int)ens_index, (int)ens_count, (int)ens_forecast_type);


Uwe Schulzweida's avatar
Uwe Schulzweida committed
904
905
  if ( varInqInst(varID) == CDI_UNDEFID )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
906
907
908
909
910
      long center, subcenter;
      int instID;
      GRIB_CHECK(grib_get_long(gh, "centre", &center), 0);
      GRIB_CHECK(grib_get_long(gh, "subCentre", &subcenter), 0);
      instID    = institutInq((int)center, (int)subcenter, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
911
      if ( instID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
912
	instID = institutDef((int)center, (int)subcenter, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
913
914
915
916
917
918
      varDefInst(varID, instID);
    }

  if ( varInqModel(varID) == CDI_UNDEFID )
    {
      int modelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
919
      long processID;
920
921
922
923
924
925
926
927
      status = grib_get_long(gh, "generatingProcessIdentifier", &processID);
      if ( status == 0 )
	{
	  modelID = modelInq(varInqInst(varID), processID, NULL);
	  if ( modelID == CDI_UNDEFID )
	    modelID = modelDef(varInqInst(varID), processID, NULL);
	  varDefModel(varID, modelID);
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
928
929
930
931
    }

  if ( varInqTable(varID) == CDI_UNDEFID )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
932
      int pdis, pcat, pnum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
933

Uwe Schulzweida's avatar
Uwe Schulzweida committed
934
      cdiDecodeParam(param, &pnum, &pcat, &pdis);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
935

Uwe Schulzweida's avatar
Uwe Schulzweida committed
936
937
938
939
      if ( pdis == 255 )
	{
	  int tableID;
	  int tabnum = pcat;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
940

Uwe Schulzweida's avatar
Uwe Schulzweida committed
941
	  tableID = tableInq(varInqModel(varID), tabnum, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
942

Uwe Schulzweida's avatar
Uwe Schulzweida committed
943
944
945
	  if ( tableID == CDI_UNDEFID )
	    tableID = tableDef(varInqModel(varID), tabnum, NULL);
	  varDefTable(varID, tableID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
946
947
948
949
950
951
952
	}
    }

  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;

  if ( CDI_Debug )
953
    Message("varID = %d  param = %d  zaxistype = %d  gridID = %d  levelID = %d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
954
	    varID, param, zaxistype, gridID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
955
956
957
}
#endif

958
int gribapiScanTimestep1(stream_t * streamptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
959
960
961
962
963
{
#if  defined  (HAVE_LIBGRIB_API)
  off_t recpos = 0;
  unsigned char *gribbuffer = NULL;
  long buffersize = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
964
  int rstatus;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
965
966
  int status;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
967
  int rtabnum = 0;
968
969
  int rcode = 0, level1 = 0, level2 = 0;
  int vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
970
  int param = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
971
972
973
974
975
  DateTime datetime, datetime0;
  int tsID;
  int varID;
  size_t readsize;
  int nrecords, nrecs, recID;
976
  int datatype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
977
978
  long recsize = 0;
  int warn_time = TRUE;
979
  // int warn_numavg = TRUE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
980
981
  int taxisID = -1;
  int rdate = 0, rtime = 0, tunit = 0, fcast = 0;
982
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
983
  int vlistID;
984
  int comptype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
985
986
987
988
  long unzipsize;
  compvar2_t compVar, compVar0;
  grib_handle *gh = NULL;
  int leveltype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
989
  int pdis = 0, pcat = 0, pnum = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
990
991
  long editionNumber;
  long lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
992
  int bitsPerValue;
993
  int lieee = FALSE;
994
  int lbounds;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
995
996
997

  streamptr->curTsID = 0;

998
  tsID  = tstepsNewEntry(streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
999
1000
1001
  taxis = &streamptr->tsteps[tsID].taxis;

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1004
  fileID = streamptr->fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1005
1006
1007
1008

  nrecs = 0;
  while ( TRUE )
    {
1009
1010
      level1 = 0;
      level2 = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
      recsize = gribGetSize(fileID);
      recpos  = fileGetPos(fileID);

      if ( recsize == 0 )
	{
	  streamptr->ntsteps = 1;
	  break;
	}
      if ( recsize > buffersize )
	{
	  buffersize = recsize;
	  gribbuffer = (unsigned char *) realloc(gribbuffer, buffersize);
	}

      readsize = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1026
1027
      rstatus = gribRead(fileID, gribbuffer, &readsize);
      if ( rstatus ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1028

1029
1030
      lieee = FALSE;

1031
      comptype = COMPRESS_NONE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1032
1033
      if ( gribGetZip(recsize, gribbuffer, &unzipsize) > 0 )
	{
1034
	  comptype = COMPRESS_SZIP;
1035
	  unzipsize += 100;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1036
1037
1038
1039
1040
1041
1042
1043
	  if ( (long) buffersize < unzipsize )
	    {
	      buffersize = unzipsize;
	      gribbuffer = (unsigned char *) realloc(gribbuffer, buffersize);
	    }
	}

      gh = grib_handle_new_from_message(NULL, (void *) gribbuffer, recsize);
1044
      GRIB_CHECK(grib_set_double(gh, "missingValue", cdiDefaultMissval), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1045
1046
1047
1048
1049

      GRIB_CHECK(grib_get_long(gh, "editionNumber", &editionNumber), 0);

      if ( editionNumber <= 1 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1050
1051
	  GRIB_CHECK(grib_get_long(gh, "table2Version", &lpar), 0);
	  rtabnum = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1052
1053
	  GRIB_CHECK(grib_get_long(gh, "indicatorOfParameter", &lpar), 0);
	  rcode = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1054

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1055
	  param = cdiEncodeParam(rcode, rtabnum, 255);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1056