stream_gribapi.c 69.6 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"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
15
16
17
18


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


extern int cdiInventoryMode;

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


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

44
  GRIB_CHECK(grib_get_long(gh, "editionNumber", &editionNumber), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
45

46
  if ( editionNumber <= 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
47
    {
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
      GRIB_CHECK(grib_get_long(gh, ">>>dataRepresentation<<<", &lpar), 0);
      gribgridtype = (int) lpar;

      switch (gribgridtype)
	{
	case  GRIB1_GTYPE_LATLON:
	case  GRIB1_GTYPE_LATLON_ROT: { gridtype = GRID_LONLAT;   break; }
	case  GRIB1_GTYPE_LCC:        { gridtype = GRID_LCC;      break; }
	case  GRIB1_GTYPE_GAUSSIAN:   { GRIB_CHECK(grib_get_long(gh, "Ni", &lpar), 0);
					if ( lpar < 0 )
	                                  gridtype = GRID_GAUSSIAN_REDUCED;
                         	        else
				          gridtype = GRID_GAUSSIAN;
          	                        break;
                                      }
	case  GRIB1_GTYPE_SPECTRAL:   { gridtype = GRID_SPECTRAL; break; }
	case  GRIB1_GTYPE_GME:        { gridtype = GRID_GME;      break; }
	}
    }
  else
  */
    {
      GRIB_CHECK(grib_get_long(gh, "gridDefinitionTemplateNumber", &lpar), 0);
      gribgridtype = (int) lpar;

      switch (gribgridtype)
	{
75
76
77
	case  GRIB2_GTYPE_LATLON:     { GRIB_CHECK(grib_get_long(gh, "Ni", &lpar), 0);
	                                if ( lpar < 0 ) break;
                                      }
78
79
80
81
82
83
84
85
86
87
88
	case  GRIB2_GTYPE_LATLON_ROT: { gridtype = GRID_LONLAT;   break; }
	case  GRIB2_GTYPE_LCC:        { gridtype = GRID_LCC;      break; }
	case  GRIB2_GTYPE_GAUSSIAN:   { GRIB_CHECK(grib_get_long(gh, "Ni", &lpar), 0);
					if ( lpar < 0 )
					  gridtype = GRID_GAUSSIAN_REDUCED;
					else
					  gridtype = GRID_GAUSSIAN;
					break;
                                      }
	case  GRIB2_GTYPE_SPECTRAL:   { gridtype = GRID_SPECTRAL; break; }
	case  GRIB2_GTYPE_GME:        { gridtype = GRID_GME;      break; }
89
	case  GRIB2_GTYPE_NUMBER:     { gridtype = GRID_REFERENCE;   break; }
90
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
91
92
93
94
95
96
    }

  return (gridtype);
}
#endif

Uwe Schulzweida's avatar
Uwe Schulzweida committed
97
#if  defined  (HAVE_LIBGRIB_API)
98
static
99
int gribapiGetIsRotated(grib_handle *gh)
100
101
{
  int isRotated = 0;
102
103
104
105
106
107
  int gribgridtype;
  long lpar;
  /*
  long editionNumber;

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

109
110
111
112
113
114
115
116
117
  if ( editionNumber <= 1 )
    {
      GRIB_CHECK(grib_get_long(gh, ">>>dataRepresentation<<<", &lpar), 0);
      gribgridtype = (int) lpar;

      if ( gribgridtype == GRIB1_GTYPE_LATLON_ROT ) isRotated = 1;
    }
  else
  */
118
    {
119
120
121
122
      GRIB_CHECK(grib_get_long(gh, "gridDefinitionTemplateNumber", &lpar), 0);
      gribgridtype = (int) lpar;

      if ( gribgridtype == GRIB2_GTYPE_LATLON_ROT ) isRotated = 1;
123
124
125
126
    }

  return (isRotated);
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
127
#endif
128

Uwe Schulzweida's avatar
Uwe Schulzweida committed
129
130
#if  defined  (HAVE_LIBGRIB_API)
static
131
int gribapiGetZaxisType(long editionNumber, int grib_ltype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
132
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
133
134
135
  int zaxistype = ZAXIS_GENERIC;

  if ( editionNumber <= 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
136
    {
137
      zaxistype = grib1ltypeToZaxisType(grib_ltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
138
139
140
    }
  else
    {
141
      zaxistype = grib2ltypeToZaxisType(grib_ltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
142
143
144
145
146
147
    }

  return (zaxistype);
}
#endif

148
149
150
151
152
153
#if  defined  (HAVE_LIBGRIB_API)
static
int gribapiGetTimeUnits(grib_handle *gh)
{
  int timeunits = -1;
  long lpar;
154
155
  size_t len = 8;
  char stepunits[8];
156
157
  static int lprint = TRUE;

158
159
  GRIB_CHECK(grib_get_string(gh, "stepUnits", stepunits, &len), 0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
160
161
162
163
164
165
166
167
168
169
170
  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;
171
  else if ( lprint )
172
    {
173
      Message("Step units >%s< unsupported!", stepunits);
174
      lprint = FALSE;
175
176
177
178
179
180
181
182
183
184
    }

  return (timeunits);
}
#endif

#if  defined  (HAVE_LIBGRIB_API)
static
int gribapiTimeIsFC(grib_handle *gh)
{
185
  long editionNumber;
186
187
  int isFC = TRUE;

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

190
191
192
193
194
195
196
197
  if ( editionNumber == 2 )
    {
      long sigofrtime;

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

      if ( sigofrtime == 3 ) isFC = FALSE;
    }
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228

  return (isFC);
}
#endif

#if  defined  (HAVE_LIBGRIB_API)
static
int gribapiGetTsteptype(grib_handle *gh)
{
  int tsteptype = 0;
  int timerange;
  long lpar;
  static int lprint = TRUE;

  if ( gribapiTimeIsFC(gh) )
    {
      GRIB_CHECK(grib_get_long(gh, "stepType", &lpar), 0);
      timerange = (int) lpar;

      // printf("timerange %d\n", timerange);

      switch ( timerange )
	{
	case  0:  tsteptype = TSTEP_AVG;    break;
	case  1:  tsteptype = TSTEP_ACCUM;  break;
	case  2:  tsteptype = TSTEP_MIN;    break;
	case  3:  tsteptype = TSTEP_MAX;    break;
	case  4:  tsteptype = TSTEP_DIFF;   break;
	default:
	  if ( lprint )
	    {
229
	      Message("Time range %d unsupported", timerange);
230
231
232
233
234
235
236
237
238
239
	      lprint = FALSE;
	    }
	}
    }

  return (tsteptype);
}
#endif

#if  defined  (HAVE_LIBGRIB_API)
240
static
241
242
243
void gribapiGetValidityDateTime(grib_handle *gh, int *vdate, int *vtime)
{
  long lpar;
244
245
246
247
248
249
250
251
252
  long sigofrtime = 3;
  long editionNumber;

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

  if ( editionNumber == 2 )
    {
      GRIB_CHECK(grib_get_long(gh, "significanceOfReferenceTime", &sigofrtime), 0);
    }
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270

  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
    {
      GRIB_CHECK(grib_get_long(gh, "validityDate", &lpar), 0);
      *vdate = (int) lpar;
      GRIB_CHECK(grib_get_long(gh, "validityTime", &lpar), 0);
      *vtime = (int) lpar*100;
    }
}
#endif

Uwe Schulzweida's avatar
Uwe Schulzweida committed
271
#if  defined  (HAVE_LIBGRIB_API)
272
static
273
void gribapiGetGrid(grib_handle *gh, grid_t *grid)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
274
{
275
  long editionNumber;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
276
  int gridtype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
277
  size_t datasize;
278
279
  long numberOfPoints;
  long lpar;
280
281
282

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

283
  gridtype = gribapiGetGridType(gh);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
284
285
286
287
288
289
290
291
  /*
  if ( streamptr->unreduced && gridtype == GRID_GAUSSIAN_REDUCED )
    {
      gridtype = GRID_GAUSSIAN;
      ISEC2_NumLon = 2*ISEC2_NumLat;
      ISEC4_NumValues = ISEC2_NumLon*ISEC2_NumLat;
    }
  */
292
  memset(grid, 0, sizeof(grid_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
293

Uwe Schulzweida's avatar
Uwe Schulzweida committed
294
295
296
  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
297
298
299
300
301
302
303
304
305
306
307
308
309
  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;

	if ( numberOfPoints != nlon*nlat )
310
	  Error("numberOfPoints (%d) and gridSize (%d) differ!",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
311
		(int)numberOfPoints, nlon*nlat);
312
313
314
315
316
317
318
319
320
321
322
	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
323
	if ( gridtype == GRID_LONLAT )
324
	  GRIB_CHECK(grib_get_double(gh, "jDirectionIncrementInDegrees", &grid->yinc), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
325

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

328
	/* if ( IS_NOT_EQUAL(grid->xfirst, 0) || IS_NOT_EQUAL(grid->xlast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
329
	  {
330
	    if ( grid->xsize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
331
	      {
332
		if ( (grid->xfirst > grid->xlast) && (grid->xfirst >= 180) ) grid->xfirst -= 360;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
333

334
335
336
		if ( editionNumber <= 1 )
		  {
		    /* correct xinc if necessary */
337
		    if ( IS_EQUAL(grid->xfirst, 0) && grid->xlast > 354 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
338
		      {
339
			double xinc = 360. / grid->xsize;
340
		      
341
			if ( fabs(grid->xinc-xinc) > 0.0 )
342
			  {
343
344
			    grid->xinc = xinc;
			    if ( CDI_Debug ) Message("set xinc to %g", grid->xinc);
345
			  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
346
347
348
		      }
		  }
	      }
349
	    grid->xdef   = 2;	    
Uwe Schulzweida's avatar
Uwe Schulzweida committed
350
	  }
351
352
	grid->ydef  = 0;
        /* if ( IS_NOT_EQUAL(grid->yfirst, 0) || IS_NOT_EQUAL(grid->ylast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
353
	  {
354
	    if ( grid->ysize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
355
	      {
356
357
358
		if ( editionNumber <= 1 )
		  {
		  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
359
	      }
360
	    grid->ydef   = 2;	    
Uwe Schulzweida's avatar
Uwe Schulzweida committed
361
362
363
364
365
	  }
	break;
      }
    case GRID_GAUSSIAN_REDUCED:
      {
366
367
368
369
370
371
372
	int nlat, i;
	size_t dummy;
	long *pl;

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

373
	grid->size   = numberOfPoints;
374

375
        grid->rowlon = (int *) malloc(nlat*sizeof(int));
376
377
378
        pl          = (long *) malloc(nlat*sizeof(long));
	dummy       = nlat;
	GRIB_CHECK(grib_get_long_array(gh, "pl", pl, &dummy), 0);
379
	for ( i = 0; i < nlat; ++i ) grid->rowlon[i] = pl[i];
380
381
	free(pl);

382
383
384
385
386
387
388
389
390
	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);
391

392
	if ( IS_EQUAL(grid->xinc, GRIB_MISSING_DOUBLE) ) grid->xinc = 0;
393

394
	/* if ( IS_NOT_EQUAL(grid->xfirst, 0) || IS_NOT_EQUAL(grid->xlast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
395
	  {
396
	    if ( grid->xsize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
397
	      {
398
		if ( (grid->xfirst > grid->xlast) && (grid->xfirst >= 180) ) grid->xfirst -= 360;
399
400
401
402

		if ( editionNumber <= 1 )
		  {
		    /* correct xinc if necessary */
403
		    if ( IS_EQUAL(grid->xfirst, 0) && grid->xlast > 354 )
404
		      {
405
			double xinc = 360. / grid->xsize;
406
		      
407
			if ( fabs(grid->xinc-xinc) > 0.0 )
408
			  {
409
410
			    grid->xinc = xinc;
			    if ( CDI_Debug ) Message("set xinc to %g", grid->xinc);
411
412
413
			  }
		      }
		  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
414
	      }
415
	    grid->xdef   = 2;	    
Uwe Schulzweida's avatar
Uwe Schulzweida committed
416
	  }
417
418
	grid->ydef  = 0;
        /* if ( IS_NOT_EQUAL(grid->yfirst, 0) || IS_NOT_EQUAL(grid->ylast, 0) ) */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
419
	  {
420
	    if ( grid->ysize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
421
	      {
422
423
424
		if ( editionNumber <= 1 )
		  {
		  }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
425
	      }
426
	    grid->ydef   = 2;	    
Uwe Schulzweida's avatar
Uwe Schulzweida committed
427
428
429
	  }
	break;
      }
430
      /*
Uwe Schulzweida's avatar
Uwe Schulzweida committed
431
432
433
    case GRID_LCC:
      {
	if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
434
	  Error("numberOfPoints (%d) and gridSize (%d) differ!",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
435
436
		ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);

437
438
439
	grid->size  = ISEC4_NumValues;
	grid->xsize = ISEC2_NumLon;
	grid->ysize = ISEC2_NumLat;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
440

441
442
443
444
445
446
447
448
449
	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
450

451
452
	grid->xdef   = 0;	    
	grid->ydef   = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
453
454
455

	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
456
      */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
457
458
    case GRID_SPECTRAL:
      {
459
460
	size_t len = 256;
	char typeOfPacking[256];
461
	GRIB_CHECK(grib_get_string(gh, "packingType", typeOfPacking, &len), 0);
462
463
	grid->lcomplex = 0;
	if ( strncmp(typeOfPacking, "spectral_complex", len) == 0 ) grid->lcomplex = 1;
464

465
	grid->size  = datasize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
466
	GRIB_CHECK(grib_get_long(gh, "J", &lpar), 0);
467
	grid->trunc = lpar;
468

Uwe Schulzweida's avatar
Uwe Schulzweida committed
469
470
471
472
	break;
      }
    case GRID_GME:
      {
473
474
475
476
477
	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;
478
479
480

	break;
      }
481
    case GRID_REFERENCE:
482
      {
483
484
485
486
	char reference_link[8192];
	size_t len = sizeof(reference_link);
	reference_link[0] = 0;

487
	grid->size  = numberOfPoints;
488
489
	if ( grib_get_long(gh, "numberOfGridUsed", &lpar) == 0 )
	  {
490
491
	    grid->number   = lpar;
	    if ( grib_get_long(gh, "numberOfGridInReference", &lpar) == 0 ) grid->position = lpar;
492
493
	    if ( grib_get_string(gh, "gridDescriptionFile", reference_link, &len) == 0 )
	      {
494
		if ( strncmp(reference_link, "file://", 7) == 0 )
495
		  grid->reference = strdupx(reference_link);
496
497
	      }
	  }
498

Uwe Schulzweida's avatar
Uwe Schulzweida committed
499
500
501
502
	break;
      }
    case GRID_GENERIC:
      {
503
504
505
506
507
	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;

508
	grid->size  = numberOfPoints;
509

510
	if ( nlon > 0 && nlat > 0 && nlon*nlat == grid->size )
511
	  {
512
513
	    grid->xsize = nlon;
	    grid->ysize = nlat;
514
515
516
	  }
	else
	  {
517
518
	    grid->xsize = 0;
	    grid->ysize = 0;
519
520
	  }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
521
522
523
524
	break;
      }
    default:
      {
525
	Error("Unsupported grid type: %s", gridNamePtr(gridtype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
526
527
528
529
	break;
      }
    }

530
  grid->isRotated = FALSE;
531
  if ( gribapiGetIsRotated(gh) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
532
    {
533
534
535
536
      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);
537
      /* change from south to north pole */
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
      grid->ypole = -grid->ypole;
      grid->xpole =  grid->xpole - 180;
    }

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

#if  defined  (HAVE_LIBGRIB_API)
static
double grib1GetLevel(grib_handle *gh, int leveltype)
{
  double dlevel;

  GRIB_CHECK(grib_get_double(gh, "level", &dlevel), 0);
  if ( leveltype == 100 ) dlevel *= 100;
  if ( dlevel < -2.e9 || dlevel > 2.e9 ) dlevel = 0;

  return (dlevel);
}

static
double grib2GetLevel(grib_handle *gh, int leveltype)
{
  double dlevel;

  GRIB_CHECK(grib_get_double(gh, "level", &dlevel), 0);
  if ( leveltype == 100 ) dlevel *= 100;
  if ( dlevel < -2.e9 || dlevel > 2.e9 ) dlevel = 0;

  return (dlevel);
}

static
void gribapiAddRecord(int streamID, int param, grib_handle *gh,
		      long recsize, off_t position, int prec, int ztype)
{
  long editionNumber;
  int zaxistype;
  int gridID = CDI_UNDEFID, varID;
  int levelID = 0;
  int tsID, recID;
  int level1, level2;
  int numavg;
  int tsteptype;
  int lbounds = 0;
  record_t *record;
  grid_t grid;
  int vlistID;
  stream_t *streamptr;
  int leveltype;
  double dlevel;
  long lpar;
  int status;
  char name[256], longname[256], units[256];
  size_t vlen; 

  streamptr = stream_to_pointer(streamID);

  vlistID = streamInqVlist(streamID);
  tsID    = streamptr->curTsID;
  recID   = recordNewEntry(streamID, tsID);
  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 )
    {
      status = grib_get_long(gh, "indicatorOfTypeOfLevel", &lpar);
      if ( status == 0 )
	{
	  leveltype = (int) lpar;
	  dlevel = grib1GetLevel(gh, leveltype);
	  if ( leveltype == 99 ) leveltype++;
	}
      else 
	{
	  leveltype = 0;
	  dlevel = 0;
	}
    }
  else
    {
      status = grib_get_long(gh, "typeOfFirstFixedSurface", &lpar);
      if ( status == 0 )
	{
	  leveltype = (int) lpar;
	  dlevel = grib2GetLevel(gh, leveltype);
	  if ( leveltype == 99 ) leveltype++;
	}
      else 
	{
	  leveltype = 0;
	  dlevel = 0;
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
639
    }
640

641
642
643
644
645
646
647
648
649
650
651
652
653
  level1 = (int) dlevel;
  level2 = 0;

  // 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
654
655
656

  gridID = varDefGrid(vlistID, grid, 0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
657
  zaxistype = gribapiGetZaxisType(editionNumber, leveltype);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
658

659
  if ( zaxistype == ZAXIS_HYBRID || zaxistype == ZAXIS_HYBRID_HALF )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
660
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
661
      int vctsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
662
      size_t dummy;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
663
      double *vctptr;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
664

Uwe Schulzweida's avatar
Uwe Schulzweida committed
665
666
      GRIB_CHECK(grib_get_long(gh, "NV", &lpar), 0);
      vctsize = lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
667
668
669
670
671
672
673
674
      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);
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
675
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
676

Uwe Schulzweida's avatar
Uwe Schulzweida committed
677
678
679
680
681
  //lbounds = cgribexGetZaxisHasBounds(ISEC1_LevelType);

  if ( prec > 32 ) prec = DATATYPE_PACK32;
  if ( prec <  0 ) prec = DATATYPE_PACK;

682
683
684
685
  name[0] = 0;
  longname[0] = 0;
  units[0] = 0;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
686
687
  vlen = 256;
  GRIB_CHECK(grib_get_string(gh, "shortName", name, &vlen), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
688
689
690
691
692
693
694
695
696
697
698
699
  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
700
  // fprintf(stderr, "param %d name %s %s %s\n", param, name, longname, units); 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
701

Uwe Schulzweida's avatar
Uwe Schulzweida committed
702
  varAddRecord(recID, param, gridID, zaxistype, lbounds, level1, level2,
703
	       prec, &varID, &levelID, 0, numavg, leveltype,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
704
	       name, longname, units);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
705
706
707
708
709
710
711
712

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

  varDefZtype(varID, ztype);

  if ( varInqInst(varID) == CDI_UNDEFID )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
713
714
715
716
717
      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
718
      if ( instID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
719
	instID = institutDef((int)center, (int)subcenter, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
720
721
722
723
724
725
      varDefInst(varID, instID);
    }

  if ( varInqModel(varID) == CDI_UNDEFID )
    {
      int modelID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
726
      long processID;
727
728
729
730
731
732
733
734
      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
735
736
737
738
    }

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
743
744
745
746
      if ( pdis == 255 )
	{
	  int tableID;
	  int tabnum = pcat;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
747

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
750
751
752
	  if ( tableID == CDI_UNDEFID )
	    tableID = tableDef(varInqModel(varID), tabnum, NULL);
	  varDefTable(varID, tableID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
753
754
755
756
757
758
759
	}
    }

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

  if ( CDI_Debug )
760
    Message("varID = %d  param = %d  zaxistype = %d  gridID = %d  levelID = %d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
761
	    varID, param, zaxistype, gridID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
762
763
764
}
#endif

765
int gribapiScanTimestep1(int streamID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
766
767
768
769
770
771
{
#if  defined  (HAVE_LIBGRIB_API)
  off_t recpos = 0;
  unsigned char *gribbuffer = NULL;
  long buffersize = 0;
  int iret = 0, ipunp = 0, iword = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
772
  int rstatus;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
773
774
  int status;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
775
  int rtabnum = 0;
776
777
  int rcode = 0, level1 = 0, level2 = 0;
  int vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
778
  int param = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
779
780
781
782
783
784
785
786
787
788
789
  DateTime datetime, datetime0;
  int tsID;
  int varID;
  size_t readsize;
  int nrecords, nrecs, recID;
  int prec;
  long recsize = 0;
  int warn_time = TRUE;
  int warn_numavg = TRUE;
  int taxisID = -1;
  int rdate = 0, rtime = 0, tunit = 0, fcast = 0;
790
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
791
792
793
794
795
796
797
  int vlistID;
  int ztype;
  long unzipsize;
  compvar2_t compVar, compVar0;
  stream_t *streamptr;
  grib_handle *gh = NULL;
  int leveltype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
798
  int pdis = 0, pcat = 0, pnum = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
799
800
  long editionNumber;
  long lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
801
  int bitsPerValue;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
802
803
804
805
  double dlevel = 0;

  streamptr = stream_to_pointer(streamID);

806
  stream_check_ptr(__func__, streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
807
808
809
810
811
812
813

  streamptr->curTsID = 0;

  tsID  = tstepsNewEntry(streamID);
  taxis = &streamptr->tsteps[tsID].taxis;

  if ( tsID != 0 )
814
    Error("Internal problem! tstepsNewEntry returns %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835

  fileID = streamInqFileID(streamID);

  nrecs = 0;
  while ( TRUE )
    {
      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
836
837
      rstatus = gribRead(fileID, gribbuffer, &readsize);
      if ( rstatus ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
838
839
840
841
842

      ztype = COMPRESS_NONE;
      if ( gribGetZip(recsize, gribbuffer, &unzipsize) > 0 )
	{
	  ztype = COMPRESS_SZIP;
843
	  unzipsize += 100;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
844
845
846
847
848
849
850
851
	  if ( (long) buffersize < unzipsize )
	    {
	      buffersize = unzipsize;
	      gribbuffer = (unsigned char *) realloc(gribbuffer, buffersize);
	    }
	}

      gh = grib_handle_new_from_message(NULL, (void *) gribbuffer, recsize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
852
      GRIB_CHECK(grib_set_double(gh, "missingValue", GRIBAPI_MISSVAL), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
853
854
855
856
857

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

      if ( editionNumber <= 1 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
858
859
	  GRIB_CHECK(grib_get_long(gh, "table2Version", &lpar), 0);
	  rtabnum = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
860
861
	  GRIB_CHECK(grib_get_long(gh, "indicatorOfParameter", &lpar), 0);
	  rcode = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
862

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
865
866
867
868
	  status = grib_get_long(gh, "indicatorOfTypeOfLevel", &lpar);
	  if ( status == 0 )
	    {
	      leveltype = (int) lpar;
869
870
	      dlevel = grib1GetLevel(gh, leveltype);
	      if ( leveltype == 99 ) leveltype++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
871
872
873
874
875
876
877
878
879
	    }
	  else 
	    {
	      leveltype = 0;
	      dlevel = 0;
	    }
	}
      else
	{
880
881
	  size_t len = 256;
	  char typeOfPacking[256];
882

883
	  status = grib_get_string(gh, "packingType", typeOfPacking, &len);
884
885
	  if ( status == 0 )
	    {
886
	      // fprintf(stderr, "packingType %d %s\n", len, typeOfPacking);
887
888
	      if ( strncmp(typeOfPacking, "grid_jpeg", len) == 0 ) ztype = COMPRESS_JPEG;
	    }
889
	  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
890
	  GRIB_CHECK(grib_get_long(gh, "discipline", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
891
892
	  pdis = (int) lpar;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
893
	  GRIB_CHECK(grib_get_long(gh, "parameterCategory", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
894
	  pcat = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
895

Uwe Schulzweida's avatar
Uwe Schulzweida committed
896
	  GRIB_CHECK(grib_get_long(gh, "parameterNumber", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
897
898
	  pnum = (int) lpar;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
899
	  param = cdiEncodeParam(pnum, pcat, pdis);
900

Uwe Schulzweida's avatar
Uwe Schulzweida committed
901
902
903
904
	  status = grib_get_long(gh, "typeOfFirstFixedSurface", &lpar);
	  if ( status == 0 )
	    {
	      leveltype = (int) lpar;
905
906
	      dlevel = grib2GetLevel(gh, leveltype);
	      if ( leveltype == 99 ) leveltype++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
907
908
909
910
911
912
913
914
915
916
917
	    }
	  else 
	    {
	      leveltype = 0;
	      dlevel = 0;
	    }
	}

      level1 = (int) dlevel;
      level2 = 0;

918
919
      gribapiGetValidityDateTime(gh, &vdate, &vtime);
      /*
920
      printf("%d %d %d.%d.%d %d %g\n", vdate, vtime, pnum, pcat, pdis, leveltype, dlevel);
921
      */
922
923
924
925
926
927
928
929
      prec = DATATYPE_PACK;
      status = grib_get_long(gh,"bitsPerValue", &lpar);
      if ( status == 0 )
	{
	  bitsPerValue = (int) lpar;
	  if ( bitsPerValue > 0 && bitsPerValue <= 32 )
	    prec = bitsPerValue;
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
930

Uwe Schulzweida's avatar
Uwe Schulzweida committed
931
932
933
934
      if ( nrecs == 0 )
	{
	  datetime0.date = vdate;
	  datetime0.time = vtime;
935
936
937
938
939
940
	  GRIB_CHECK(grib_get_long(gh, "dataDate", &lpar), 0);
	  rdate = (int) lpar;
	  GRIB_CHECK(grib_get_long(gh, "dataTime", &lpar), 0);
	  rtime = (int) lpar*100;
	  fcast = gribapiTimeIsFC(gh);
	  if ( fcast ) tunit = gribapiGetTimeUnits(gh);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
941
942
943
944
945
	}
      else
	{
	  datetime.date  = vdate;
	  datetime.time  = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
946
947

	  compVar.param  = param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
948
949
          compVar.level1 = level1;
          compVar.level2 = level2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
950
951
	  compVar.ltype  = leveltype;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
952
953
	  for ( recID = 0; recID < nrecs; recID++ )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
954
	      compVar0.param  = streamptr->tsteps[0].records[recID].param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
955
956
957
	      compVar0.level1 = streamptr->tsteps[0].records[recID].ilevel;
	      compVar0.level2 = streamptr->tsteps[0].records[recID].ilevel2;
	      compVar0.ltype  = streamptr->tsteps[0].records[recID].ltype;
958
959
960
961
	      /*
	      printf("var0: %d %d %d %d %d\n", recID, compVar0.param, compVar0.level1, compVar0.level2, compVar0.ltype);
	      printf("var1: %d %d %d %d %d\n", recID, compVar.param, compVar.level1, compVar.level2, compVar.ltype);
	      */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
962
963
964
965
966
967
968
969
970
	      if ( memcmp(&compVar0, &compVar, sizeof(compvar2_t)) == 0 ) break;
	    }

	  if ( cdiInventoryMode == 1 )
	    {
	      if ( recID < nrecs ) break;
	      if ( warn_time )
		if ( memcmp(&datetime, &datetime0, sizeof(DateTime)) != 0 )
		  {
971
972
		    char paramstr[32];
		    cdiParamToString(param, paramstr, sizeof(paramstr));	    
973
		    Warning("Inconsistent verification time (param=%s level=%d)", paramstr, level1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
974
975
976
977
978
979
980
981
982
		    warn_time = FALSE;
		  }
	    }
	  else
	    {
	      if ( memcmp(&datetime, &datetime0, sizeof(DateTime)) != 0 ) break;

	      if ( recID < nrecs )
		{
983
984
		  char paramstr[32];
		  cdiParamToString(param, paramstr, sizeof(paramstr));
985
		  Warning("Param=%s level=%d already exist, skipped!", paramstr, level1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
986
987
988
989
990
991
992
993
994
		  continue;
		}
	    }
	}
      /*
      if ( ISEC1_AvgNum )
	{
	  if (  taxis->numavg && warn_numavg && (taxis->numavg != ISEC1_AvgNum) )
	    {
995
	      Message("Change numavg from %d to %d not allowed!",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
		      taxis->numavg, ISEC1_AvgNum);
	      warn_numavg = FALSE;
	    }
	  else
	    {
	      taxis->numavg = ISEC1_AvgNum;
	    }
	}
      */
      nrecs++;

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1010
      gribapiAddRecord(streamID, param, gh, recsize, recpos, prec, ztype);
1011
1012
1013

      grib_handle_delete(gh);
      gh = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1014
1015
    }

1016
1017
  if ( gh ) grib_handle_delete(gh);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1018
1019
  streamptr->rtsteps = 1;

1020
1021
  if ( nrecs == 0 ) return (CDI_EUFSTRUCT);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
  cdiGenVars(streamID);

  if ( fcast )
    {
      taxisID = taxisCreate(TAXIS_RELATIVE);
      taxis->type  = TAXIS_RELATIVE;
      taxis->rdate = rdate;
      taxis->rtime = rtime;
      taxis->unit  = tunit;
    }
  else
    {
      taxisID = taxisCreate(TAXIS_ABSOLUTE);
      taxis->type  = TAXIS_ABSOLUTE;
    }

  taxis->vdate = datetime0.date;
  taxis->vtime = datetime0.time;

  vlistID = streamInqVlist(streamID);
  vlistDefTaxis(vlistID, taxisID);

  nrecords = streamptr->tsteps[0].nallrecs;
  if ( nrecords < streamptr->tsteps[0].recordSize )
    {
      streamptr->tsteps[0].recordSize = nrecords;
      streamptr->tsteps[0].records =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1049
      (record_t *) realloc(streamptr->tsteps[0].records, nrecords*sizeof(record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
    }

  streamptr->tsteps[0].recIDs = (int *) malloc(nrecords*sizeof(int));
  streamptr->tsteps[0].nrecs = nrecords;
  for ( recID = 0; recID < nrecords; recID++ )
    streamptr->tsteps[0].recIDs[recID] = recID;

  streamptr->record->buffer     = gribbuffer;
  streamptr->record->buffersize = buffersize;

  if ( streamptr->ntsteps == -1 )
    {
      tsID = tstepsNewEntry(streamID);
      if ( tsID != streamptr->rtsteps )
1064
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1065
1066
1067
1068
1069
1070
1071
1072
1073
1074
1075
1076
1077
1078
1079
1080
1081

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

  if ( streamptr->ntsteps == 1 )
    {
      if ( taxis->vdate == 0 && taxis->vtime == 0 )
	{
	  streamptr->ntsteps = 0;
	  for ( varID = 0; varID < streamptr->nvars; varID++ )
	    {
	      vlistDefVarTime(vlistID, varID, TIME_CONSTANT);
	    }
	}
    }
#else
1082
  Error("GRIB_API support not compiled in!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1083
#endif
1084
1085

  return (0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1086
1087
1088
1089
1090
}


int gribapiScanTimestep2(int streamID)
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1091
  int rstatus = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1092
#if  defined  (HAVE_LIBGRIB_API)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1093
  int status;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1094
1095
1096
1097
1098
  off_t recpos = 0;
  unsigned char *gribbuffer = NULL;
  long buffersize = 0;
  int iret = 0, ipunp = 0, iword = 0;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1099
  int rtabnum = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1100
1101
1102
1103
1104
1105
1106
1107
  int rcode = 0, level1 = 0, level2 = 0, vdate = 0, vtime = 0;
  DateTime datetime, datetime0;
  int tsID;
  int varID, gridID;
  size_t readsize;
  int nrecords, nrecs, recID, rindex;
  long recsize = 0;
  int warn_numavg = TRUE;
1108
  int tsteptype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1109
  int taxisID = -1;
1110
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1111
1112
1113
1114
1115
1116
  int vlistID;
  long unzipsize;
  compvar2_t compVar, compVar0;
  stream_t *streamptr;
  grib_handle *gh = NULL;
  int leveltype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1117
1118
  int pdis = 0, pcat = 0, pnum = 0;
  int param = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1119
1120
1121
1122
1123
1124
  long editionNumber;
  long lpar;
  double dlevel = 0;

  streamptr = stream_to_pointer(streamID);

1125
  stream_check_ptr(__func__, streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135
1136
1137

  streamptr->curTsID = 1;

  fileID  = streamInqFileID(streamID);
  vlistID = streamInqVlist(streamID);
  taxisID = vlistInqTaxis(vlistID);

  gribbuffer = (unsigned char *) streamptr->record->buffer;
  buffersize = streamptr->record->buffersize;
                                          
  tsID = streamptr->rtsteps;
  if ( tsID != 1 )
1138
    Error("Internal problem! unexpeceted timestep %d", tsID+1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
1157
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
1170
1171
1172
1173
1174
1175
1176
1177
1178
1179

  taxis = &streamptr->tsteps[tsID].taxis;

  fileSetPos(fileID, streamptr->tsteps[tsID].position, SEEK_SET);

  cdiCreateRecords(streamID, tsID);

  nrecords = streamptr->tsteps[tsID].nallrecs;
  streamptr->tsteps[1].recIDs = (int *) malloc(nrecords*sizeof(int));
  streamptr->tsteps[1].nrecs = 0;
  for ( recID = 0; recID < nrecords; recID++ )
    streamptr->tsteps[1].recIDs[recID] = -1;
      
  for ( recID = 0; recID < nrecords; recID++ )
    {
      varID = streamptr->tsteps[0].records[recID].varID;
      streamptr->tsteps[tsID].records[recID].position = 
	streamptr->tsteps[0].records[recID].position;
      streamptr->tsteps[tsID].records[recID].size     = 
	streamptr->tsteps[0].records[recID].size;
    }

  rindex = 0;
  while ( TRUE )
    {
      if ( rindex > nrecords ) break;

      recsize = gribGetSize(fileID);
      recpos  = fileGetPos(fileID);
      if ( recsize == 0 )
	{
	  streamptr->ntsteps = 2;
	  break;
	}
      if ( recsize > buffersize )
	{
	  buffersize = recsize;
	  gribbuffer = (unsigned char *) realloc(gribbuffer, (size_t)buffersize);
	}

      readsize = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1180
1181
      rstatus = gribRead(fileID, gribbuffer, &readsize);
      if ( rstatus ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1182
1183
1184
1185
1186
1187
1188
1189
1190
1191
1192
1193

      if ( gribGetZip(recsize, gribbuffer, &unzipsize) > 0 )
	{
	  unzipsize += 100; /* need 0 to 1 bytes for rounding of bds */
	  if ( (long) buffersize < unzipsize )
	    {
	      buffersize = unzipsize;
	      gribbuffer = (unsigned char *) realloc(gribbuffer, buffersize);
	    }
	}

      gh = grib_handle_new_from_message(NULL, (void *) gribbuffer, recsize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1194
      GRIB_CHECK(grib_set_double(gh, "missingValue", GRIBAPI_MISSVAL), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1195
1196
1197
1198
1199

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

      if ( editionNumber <= 1 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1200
1201
	  GRIB_CHECK(grib_get_long(gh, "table2Version", &lpar), 0);
	  rtabnum = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1202
1203
	  GRIB_CHECK(grib_get_long(gh, "indicatorOfParameter", &lpar), 0);
	  rcode = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1204

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1207
1208
1209
1210
	  status = grib_get_long(gh, "indicatorOfTypeOfLevel", &lpar);
	  if ( status == 0 )
	    {
	      leveltype = (int) lpar;
1211
1212
	      dlevel = grib1GetLevel(gh, leveltype);
	      if ( leveltype == 99 ) leveltype++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1213
1214
1215
1216
1217
1218
1219
1220
1221
	    }
	  else 
	    {
	      leveltype = 0;
	      dlevel = 0;
	    }
	}
      else
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1222
	  GRIB_CHECK(grib_get_long(gh, "discipline", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1223
	  pdis = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1224

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1225
	  GRIB_CHECK(grib_get_long(gh, "parameterCategory", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1226
	  pcat = (int) lpar;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1227

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1228
	  GRIB_CHECK(grib_get_long(gh, "parameterNumber", &lpar), 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1229
1230
	  pnum = (int) lpar;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1231