stream_cgribex.c 53.8 KB
Newer Older
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
#if defined (HAVE_CONFIG_H)
#  include "config.h"
#endif

#include <stdio.h>
#include <float.h>  /* FLT_EPSILON */

#include "dmemory.h"
#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
15
#include "stream_grb.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
16

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
21
extern int cdiInventoryMode;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
22
23

typedef struct {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
24
  int param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
25
26
27
  int level1;
  int level2;
  int ltype;
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
36
37
38
{
  int gridtype = 0;

  switch (ISEC2_GridType)
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
39
40
    case  GRIB1_GTYPE_LATLON:
    case  GRIB1_GTYPE_LATLON_ROT:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
41
42
43
44
      {
	gridtype = GRID_LONLAT;
	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
45
    case  GRIB1_GTYPE_LCC:
46
      {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
47
	gridtype = GRID_LCC;
48
49
	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
50
    case  GRIB1_GTYPE_GAUSSIAN:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
51
52
53
54
55
56
57
      {
	if ( ISEC2_Reduced )
	  gridtype = GRID_GAUSSIAN_REDUCED;
	else
	  gridtype = GRID_GAUSSIAN;
	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
58
    case  GRIB1_GTYPE_SPECTRAL:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
59
60
61
62
      {
	gridtype = GRID_SPECTRAL;
	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
63
    case  GRIB1_GTYPE_GME:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
64
65
66
67
68
69
70
71
72
73
74
75
76
      {
	gridtype = GRID_GME;
	break;
      }
    default:
      {
	gridtype = GRID_GENERIC;
	break;
      }
    }

  return (gridtype);
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
77
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
78

Uwe Schulzweida's avatar
Uwe Schulzweida committed
79
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
80
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
81
int cgribexGetIsRotated(int *isec2)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
82
83
84
85
86
87
88
89
90
91
{
  int isRotated = 0;

  if ( ISEC2_GridType == 10 )
    {
      isRotated = 1;
    }

  return (isRotated);
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
92
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
93

Uwe Schulzweida's avatar
Uwe Schulzweida committed
94
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
95
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
96
int cgribexGetZaxisHasBounds(int grb_ltype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
97
98
99
{
  int lbounds = 0;

100
  switch (grb_ltype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
101
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
102
103
    case GRIB1_LTYPE_HYBRID_LAYER:
    case GRIB1_LTYPE_LANDDEPTH_LAYER:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
104
105
106
107
108
109
110
111
      {
	lbounds = 1;
	break;
      }
    }

  return (lbounds);
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
112
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
113

Uwe Schulzweida's avatar
Uwe Schulzweida committed
114
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
115
static
Uwe Schulzweida's avatar
Uwe Schulzweida committed
116
int cgribexGetTimeUnit(int *isec1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
117
118
119
120
121
122
{
  int timeunit = -1;
  static int lprint = TRUE;

  switch ( ISEC1_TimeUnit )
    {
123
124
125
    case ISEC1_TABLE4_MINUTE:  timeunit = TUNIT_MINUTE;  break;
    case ISEC1_TABLE4_QUARTER: timeunit = TUNIT_QUARTER; break;
    case ISEC1_TABLE4_HOUR:    timeunit = TUNIT_HOUR;    break;
126
127
128
    case ISEC1_TABLE4_3HOURS:  timeunit = TUNIT_3HOURS;  break;
    case ISEC1_TABLE4_6HOURS:  timeunit = TUNIT_6HOURS;  break;
    case ISEC1_TABLE4_12HOURS: timeunit = TUNIT_12HOURS; break;
129
    case ISEC1_TABLE4_DAY:     timeunit = TUNIT_DAY;     break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
130
131
132
    default:
      if ( lprint )
	{
133
	  Message("Time unit %d unsupported", ISEC1_TimeUnit);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
134
135
136
137
138
139
	  lprint = FALSE;
	}
    }

  return (timeunit);
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
140
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
141

Uwe Schulzweida's avatar
Uwe Schulzweida committed
142
#if  defined  (HAVE_LIBCGRIBEX)
143
144
145
146
147
148
149
150
151
152
static
int cgribexTimeIsFC(int *isec1)
{
  int isFC = TRUE;

  if ( ISEC1_TimeRange == 10 && ISEC1_TimePeriod1 == 0 && ISEC1_TimePeriod2 == 0 )
    isFC = FALSE;

  return (isFC);
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
153
#endif
154

Uwe Schulzweida's avatar
Uwe Schulzweida committed
155
#if  defined  (HAVE_LIBCGRIBEX)
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
static
int cgribexGetTsteptype(int timerange)
{
  int tsteptype = 0;
  static int lprint = TRUE;

  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 )
	{
174
	  Message("Time range %d unsupported", timerange);
175
176
177
178
179
180
	  lprint = FALSE;
	}
    }

  return (tsteptype);
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
181
#endif
182

Uwe Schulzweida's avatar
Uwe Schulzweida committed
183
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
184
static
185
186
void cgribexAddRecord(int streamID, int param, int *isec1, int *isec2, double *fsec2, double *fsec3,
		      int *isec4, long recsize, off_t position, int prec, int ztype, int lmv)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
187
188
{
  int gridtype;
189
  int zaxistype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
190
  int gridID = CDI_UNDEFID, varID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
191
192
193
194
  int levelID = 0;
  int tsID, recID;
  int level1, level2;
  int numavg;
195
  int tsteptype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
196
  int lbounds = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
197
  record_t *record;
198
  grid_t grid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
199
  int vlistID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
200
  stream_t *streamptr;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
201
202

  streamptr = stream_to_pointer(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
203
204

  vlistID = streamInqVlist(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
205
  tsID    = streamptr->curTsID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
206
  recID   = recordNewEntry(streamID, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
207
  record  = &streamptr->tsteps[tsID].records[recID];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
208

209
210
  tsteptype = cgribexGetTsteptype(ISEC1_TimeRange);
  numavg    = ISEC1_AvgNum;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
211
212
213
214

  level1  = ISEC1_Level1;
  level2  = ISEC1_Level2;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
215
  /* fprintf(stderr, "param %d %d %d %d\n", param, level1, level2, ISEC1_LevelType); */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
216
217
218

  (*record).size     = recsize;
  (*record).position = position;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
219
  (*record).param    = param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
220
221
222
223
  (*record).ilevel   = level1;
  (*record).ilevel2  = level2;
  (*record).ltype    = ISEC1_LevelType;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
224
  gridtype = cgribexGetGridType(isec2);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
225

Uwe Schulzweida's avatar
Uwe Schulzweida committed
226
  if ( streamptr->unreduced && gridtype == GRID_GAUSSIAN_REDUCED )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
227
228
229
230
231
232
    {
      gridtype = GRID_GAUSSIAN;
      ISEC2_NumLon = 2*ISEC2_NumLat;
      ISEC4_NumValues = ISEC2_NumLon*ISEC2_NumLat;
    }

233
  memset(&grid, 0, sizeof(grid_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
234
235
236
237
238
239
  switch (gridtype)
    {
    case GRID_LONLAT:
    case GRID_GAUSSIAN:
      {
	if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
240
	  Error("numberOfPoints (%d) and gridSize (%d) differ!",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
241
242
243
244
245
246
247
248
249
250
251
		ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);
	grid.size  = ISEC4_NumValues;
	grid.xsize = ISEC2_NumLon;
	grid.ysize = ISEC2_NumLat;
	grid.xinc  = 0;
	grid.yinc  = 0;
	grid.xdef  = 0;
	/* if ( ISEC2_FirstLon != 0 || ISEC2_LastLon != 0 ) */
	  {
	    if ( grid.xsize > 1 )
	      {
252
253
254
		if ( ISEC2_ResFlag && ISEC2_LonIncr > 0 )
		  if ( ISEC2_LonIncr*(grid.xsize-1) > 360000 ) ISEC2_LonIncr = 0;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
255
256
257
258
259
260
261
262
263
264
265
266
267
		if ( ISEC2_ResFlag && ISEC2_LonIncr > 0 )
		  grid.xinc = ISEC2_LonIncr * 0.001;
		else
		  grid.xinc = (ISEC2_LastLon - ISEC2_FirstLon) * 0.001 / (grid.xsize - 1);

		/* correct xinc if necessary */
		if ( ISEC2_FirstLon == 0 && ISEC2_LastLon > 354000 )
		  {
		    double xinc = 360. / grid.xsize;

		    if ( fabs(grid.xinc-xinc) > 0.0 )
		      {
			grid.xinc = xinc;
268
			if ( CDI_Debug ) Message("set xinc to %g", grid.xinc);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
269
270
271
272
273
274
275
276
277
278
279
280
		      }
		  }
	      }
	    grid.xfirst = ISEC2_FirstLon * 0.001;
	    grid.xlast  = ISEC2_LastLon  * 0.001;
	    grid.xdef   = 2;	    
	  }
	grid.ydef  = 0;
	/* if ( ISEC2_FirstLat != 0 || ISEC2_LastLat != 0 ) */
	  {
	    if ( grid.ysize > 1 )
	      {
281
282
283
		if ( ISEC2_ResFlag && ISEC2_LatIncr > 0 )
		  if ( ISEC2_LatIncr*(grid.ysize-1) > 180000 ) ISEC2_LatIncr = 0;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
		if ( ISEC2_ResFlag && ISEC2_LatIncr > 0 )
		  grid.yinc = ISEC2_LatIncr * 0.001;
		else
		  grid.yinc = (ISEC2_LastLat - ISEC2_FirstLat) * 0.001 / (grid.ysize - 1);
	      }
	    grid.yfirst = ISEC2_FirstLat * 0.001;
	    grid.ylast  = ISEC2_LastLat  * 0.001;
	    grid.ydef   = 2;	    
	  }
	break;
      }
    case GRID_GAUSSIAN_REDUCED:
      {
	grid.size   = ISEC4_NumValues;
        grid.rowlon = ISEC2_RowLonPtr;
	grid.ysize  = ISEC2_NumLat;
	grid.xinc  = 0;
	grid.yinc  = 0;
	grid.xdef  = 0;
	/* if ( ISEC2_FirstLon != 0 || ISEC2_LastLon != 0 ) */
	  {
	    if ( grid.xsize > 1 )
	      {
		if ( ISEC2_ResFlag && ISEC2_LonIncr > 0 )
		  grid.xinc = ISEC2_LonIncr * 0.001;
		else
		  grid.xinc = (ISEC2_LastLon - ISEC2_FirstLon) * 0.001 / (grid.xsize - 1);
	      }
	    grid.xfirst = ISEC2_FirstLon * 0.001;
	    grid.xlast  = ISEC2_LastLon  * 0.001;
	    grid.xdef   = 2;	    
	  }
	grid.ydef  = 0;
	/* if ( ISEC2_FirstLat != 0 || ISEC2_LastLat != 0 ) */
	  {
	    if ( grid.ysize > 1 )
	      {
		if ( ISEC2_ResFlag && ISEC2_LatIncr > 0 )
		  grid.yinc = ISEC2_LatIncr * 0.001;
		else
		  grid.yinc = (ISEC2_LastLat - ISEC2_FirstLat) * 0.001 / (grid.ysize - 1);
	      }
	    grid.yfirst = ISEC2_FirstLat * 0.001;
	    grid.ylast  = ISEC2_LastLat  * 0.001;
	    grid.ydef   = 2;	    
	  }
	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
332
    case GRID_LCC:
333
334
      {
	if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
335
	  Error("numberOfPoints (%d) and gridSize (%d) differ!",
336
		ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
337

338
339
340
341
	grid.size  = ISEC4_NumValues;
	grid.xsize = ISEC2_NumLon;
	grid.ysize = ISEC2_NumLat;

342
343
	grid.lcc_xinc      = ISEC2_Lambert_dx;
	grid.lcc_yinc      = ISEC2_Lambert_dy;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
344
345
	grid.lcc_originLon = ISEC2_FirstLon * 0.001;
	grid.lcc_originLat = ISEC2_FirstLat * 0.001;
346
347
348
349
350
	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;
351

352
	grid.xdef   = 0;	    
353
354
355
356
	grid.ydef   = 0;

	break;
      }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
357
358
359
360
    case GRID_SPECTRAL:
      {
	grid.size  = ISEC4_NumValues;
	grid.trunc = ISEC2_PentaJ;
361
362
363
364
365
	if ( ISEC2_RepMode == 2 )
	  grid.lcomplex = 1;
	else
	  grid.lcomplex = 0;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
	break;
      }
    case GRID_GME:
      {
	grid.size  = ISEC4_NumValues;
	grid.nd    = ISEC2_GME_ND;
	grid.ni    = ISEC2_GME_NI;
	grid.ni2   = ISEC2_GME_NI2;
	grid.ni3   = ISEC2_GME_NI3;
	break;
      }
    case GRID_GENERIC:
      {
	grid.size  = ISEC4_NumValues;
	grid.xsize = 0;
	grid.ysize = 0;
	break;
      }
    default:
      {
386
	Error("Unsupported grid type: %s", gridNamePtr(gridtype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
387
388
389
390
391
	break;
      }
    }

  grid.isRotated = FALSE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
392
  if ( cgribexGetIsRotated(isec2) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
393
394
395
396
397
398
399
400
401
402
403
404
405
    {
      grid.isRotated = TRUE;
      grid.ypole     = - ISEC2_LatSP * 0.001;
      grid.xpole     =   ISEC2_LonSP * 0.001 - 180;
      grid.angle     = 0;
    }

  grid.xvals = NULL;
  grid.yvals = NULL;
  grid.type  = gridtype;

  gridID = varDefGrid(vlistID, grid, 0);

406
  zaxistype = grib1ltypeToZaxisType(ISEC1_LevelType);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
407

408
  if ( zaxistype == ZAXIS_HYBRID || zaxistype == ZAXIS_HYBRID_HALF )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
409
410
411
412
413
414
415
    {
      int vctsize = ISEC2_NumVCP;
      double *vctptr = &fsec2[10];

      varDefVCT(vctsize, vctptr);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
416
  lbounds = cgribexGetZaxisHasBounds(ISEC1_LevelType);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
417

418
419
  if ( prec > 32 ) prec = DATATYPE_PACK32;
  if ( prec <  0 ) prec = DATATYPE_PACK;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
420

Uwe Schulzweida's avatar
Uwe Schulzweida committed
421
  varAddRecord(recID, param, gridID, zaxistype, lbounds, level1, level2,
422
	       prec, &varID, &levelID, tsteptype, numavg, ISEC1_LevelType, NULL, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
423
424
425
426

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

427
  varDefZtype(varID, ztype);
428

429
430
  if ( lmv ) varDefMissval(varID, FSEC3_MissVal);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
431
  if ( varInqInst(varID) == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
432
433
434
435
436
    {
      int center, subcenter, instID;
      center    = ISEC1_CenterID;
      subcenter = ISEC1_SubCenterID;
      instID    = institutInq(center, subcenter, NULL, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
437
      if ( instID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
438
439
440
441
	instID = institutDef(center, subcenter, NULL, NULL);
      varDefInst(varID, instID);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
442
  if ( varInqModel(varID) == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
443
444
445
    {
      int modelID;
      modelID = modelInq(varInqInst(varID), ISEC1_ModelID, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
446
      if ( modelID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
447
448
449
450
	modelID = modelDef(varInqInst(varID), ISEC1_ModelID, NULL);
      varDefModel(varID, modelID);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
451
  if ( varInqTable(varID) == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
452
453
    {
      int tableID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
454

Uwe Schulzweida's avatar
Uwe Schulzweida committed
455
      tableID = tableInq(varInqModel(varID), ISEC1_CodeTable, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
456

Uwe Schulzweida's avatar
Uwe Schulzweida committed
457
      if ( tableID == CDI_UNDEFID )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
458
459
460
461
	tableID = tableDef(varInqModel(varID), ISEC1_CodeTable, NULL);
      varDefTable(varID, tableID);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
462
463
  streamptr->tsteps[tsID].nallrecs++;
  streamptr->nrecs++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
464
465

  if ( CDI_Debug )
466
    Message("varID = %d  param = %d  zaxistype = %d  gridID = %d  levelID = %d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
467
	    varID, param, zaxistype, gridID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
468
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
469
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
470

Uwe Schulzweida's avatar
Uwe Schulzweida committed
471
#if  defined  (HAVE_LIBCGRIBEX)
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
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);
      }
    }
  }
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
496
#endif
497

Uwe Schulzweida's avatar
Uwe Schulzweida committed
498
#if  defined  (HAVE_LIBCGRIBEX)
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
static
void cgribexDecodeHeader(int *isec0, int *isec1, int *isec2, double *fsec2,
			 int *isec3, double *fsec3, int *isec4, double *fsec4, 
			 int *gribbuffer, int recsize, int *lmv)
{
  int iret = 0, ipunp = 0, iword = 0;

  gribExDP(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, fsec4,
	   ipunp, (int *) gribbuffer, recsize, &iword, "J", &iret);

  *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;
    }
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
520
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
521

522
int cgribexScanTimestep1(int streamID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
523
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
524
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
525
526
  int *isec0, *isec1, *isec2, *isec3, *isec4;
  double fsec2[512], fsec3[2], *fsec4 = NULL;
527
  int lmv = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
528
529
530
  off_t recpos = 0;
  unsigned char *gribbuffer = NULL;
  long buffersize = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
531
  int rstatus;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
532
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
533
534
  int param = 0;
  int level1 = 0, level2 = 0, vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
535
536
537
538
539
540
541
  DateTime datetime, datetime0;
  int tsID;
  int varID;
  size_t readsize;
  int nrecords, nrecs, recID;
  int prec;
  long recsize = 0;
542
  int warn_time = TRUE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
543
544
545
  int warn_numavg = TRUE;
  int taxisID = -1;
  int rdate = 0, rtime = 0, tunit = 0, fcast = 0;
546
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
547
  int vlistID;
548
  int ztype;
549
  long unzipsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
550
551
  compvar_t compVar, compVar0;
  stream_t *streamptr;
552
553
  extern int cdiSkipRecords;
  int nskip = cdiSkipRecords;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
554
555

  streamptr = stream_to_pointer(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
556

557
  stream_check_ptr(__func__, streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
558

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
561
562
563
564
565
  isec0 = streamptr->record->sec0;
  isec1 = streamptr->record->sec1;
  isec2 = streamptr->record->sec2;
  isec3 = streamptr->record->sec3;
  isec4 = streamptr->record->sec4;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
566
567

  tsID  = tstepsNewEntry(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
568
  taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
569
570

  if ( tsID != 0 )
571
    Error("Internal problem! tstepsNewEntry returns %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
572
573
574

  fileID = streamInqFileID(streamID);

575
576
577
  while ( nskip-- > 0 )
    {
      recsize = gribGetSize(fileID);
578
      if ( recsize == 0 )
579
	Error("Skipping of %d records failed!", cdiSkipRecords);
580

581
582
583
584
      recpos  = fileGetPos(fileID);
      fileSetPos(fileID, recsize, SEEK_CUR);
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
585
586
587
588
589
590
591
592
  nrecs = 0;
  while ( TRUE )
    {
      recsize = gribGetSize(fileID);
      recpos  = fileGetPos(fileID);

      if ( recsize == 0 )
	{
593
	  if ( nrecs == 0 )
594
	    Error("No GRIB records found!");
595

Uwe Schulzweida's avatar
Uwe Schulzweida committed
596
	  streamptr->ntsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
597
598
599
600
601
602
603
604
605
	  break;
	}
      if ( recsize > buffersize )
	{
	  buffersize = recsize;
	  gribbuffer = (unsigned char *) realloc(gribbuffer, buffersize);
	}

      readsize = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
606
607
      rstatus = gribRead(fileID, gribbuffer, &readsize);
      if ( rstatus ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
608

609
      ztype = COMPRESS_NONE;
610
611
      if ( gribGetZip(recsize, gribbuffer, &unzipsize) > 0 )
	{
612
	  ztype = COMPRESS_SZIP;
613
614
615
616
617
618
619
620
	  unzipsize += 100; /* need 0 to 1 bytes for rounding of bds */
	  if ( (long) buffersize < unzipsize )
	    {
	      buffersize = unzipsize;
	      gribbuffer = (unsigned char *) realloc(gribbuffer, buffersize);
	    }
	}

621
622
      cgribexDecodeHeader(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, fsec4,
			  (int *) gribbuffer, recsize, &lmv);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
623

Uwe Schulzweida's avatar
Uwe Schulzweida committed
624
      param = cdiEncodeParam(ISEC1_Parameter, ISEC1_CodeTable, 255);
625
626
      if ( ISEC1_LevelType == 100 ) ISEC1_Level1 *= 100;
      if ( ISEC1_LevelType ==  99 ) ISEC1_LevelType = 100;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
627
628
629
      level1   = ISEC1_Level1;
      level2   = ISEC1_Level2;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
630
631
      gribDateTime(isec1, &vdate, &vtime);

632
633
634
635
      if ( ISEC4_NumBits > 0 && ISEC4_NumBits <= 32 )
	prec = ISEC4_NumBits;
      else
        prec = DATATYPE_PACK;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
636
637
638
639
640
641
642

      if ( nrecs == 0 )
	{
	  datetime0.date = vdate;
	  datetime0.time = vtime;
	  rdate = gribRefDate(isec1);
	  rtime = gribRefTime(isec1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
643
	  tunit = cgribexGetTimeUnit(isec1);
644
	  fcast = cgribexTimeIsFC(isec1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
645
646
647
648
649
	}
      else
	{
	  datetime.date  = vdate;
	  datetime.time  = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
650
	  compVar.param  = param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
651
652
653
654
655
          compVar.level1 = level1;
          compVar.level2 = level2;
          compVar.ltype  = ISEC1_LevelType;
	  for ( recID = 0; recID < nrecs; recID++ )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
656
	      compVar0.param  = streamptr->tsteps[0].records[recID].param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
657
658
659
	      compVar0.level1 = streamptr->tsteps[0].records[recID].ilevel;
	      compVar0.level2 = streamptr->tsteps[0].records[recID].ilevel2;
	      compVar0.ltype  = streamptr->tsteps[0].records[recID].ltype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
660

Uwe Schulzweida's avatar
Uwe Schulzweida committed
661
	      if ( memcmp(&compVar0, &compVar, sizeof(compvar_t)) == 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
662
	    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
663
664
665
666
667
668
669

	  if ( cdiInventoryMode == 1 )
	    {
	      if ( recID < nrecs ) break;
	      if ( warn_time )
		if ( memcmp(&datetime, &datetime0, sizeof(DateTime)) != 0 )
		  {
670
671
		    char paramstr[32];
		    cdiParamToString(param, paramstr, sizeof(paramstr));
672
		    Warning("Inconsistent verification time (param=%s level=%d)", paramstr, level1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
673
674
675
676
677
678
679
680
681
		    warn_time = FALSE;
		  }
	    }
	  else
	    {
	      if ( memcmp(&datetime, &datetime0, sizeof(DateTime)) != 0 ) break;

	      if ( recID < nrecs )
		{
682
683
		  char paramstr[32];
		  cdiParamToString(param, paramstr, sizeof(paramstr));
684
		  Warning("Param=%s level=%d already exist, skipped!", paramstr, level1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
685
686
687
		  continue;
		}
	    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
688
689
690
691
692
693
	}

      if ( ISEC1_AvgNum )
	{
	  if (  taxis->numavg && warn_numavg && (taxis->numavg != ISEC1_AvgNum) )
	    {
694
	      Warning("Changing numavg from %d to %d not supported!",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
695
696
697
698
699
700
701
702
703
704
705
706
		      taxis->numavg, ISEC1_AvgNum);
	      warn_numavg = FALSE;
	    }
	  else
	    {
	      taxis->numavg = ISEC1_AvgNum;
	    }
	}

      nrecs++;

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

709
710
      cgribexAddRecord(streamID, param, isec1, isec2, fsec2, fsec3,
		       isec4, recsize, recpos, prec, ztype, lmv);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
711
712
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
713
  streamptr->rtsteps = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
714

715
716
  if ( nrecs == 0 ) return (CDI_EUFSTRUCT);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
717
718
719
720
  cdiGenVars(streamID);

  if ( fcast )
    {
721
      taxisID = taxisCreate(TAXIS_RELATIVE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
722
723
724
725
726
727
728
      taxis->type  = TAXIS_RELATIVE;
      taxis->rdate = rdate;
      taxis->rtime = rtime;
      taxis->unit  = tunit;
    }
  else
    {
729
      taxisID = taxisCreate(TAXIS_ABSOLUTE);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
730
      taxis->type  = TAXIS_ABSOLUTE;
731
      taxis->unit  = tunit;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
732
733
734
735
736
737
738
739
    }

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
740
741
  nrecords = streamptr->tsteps[0].nallrecs;
  if ( nrecords < streamptr->tsteps[0].recordSize )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
742
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
743
744
      streamptr->tsteps[0].recordSize = nrecords;
      streamptr->tsteps[0].records =
Uwe Schulzweida's avatar
Uwe Schulzweida committed
745
      (record_t *) realloc(streamptr->tsteps[0].records, nrecords*sizeof(record_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
746
747
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
748
749
  streamptr->tsteps[0].recIDs = (int *) malloc(nrecords*sizeof(int));
  streamptr->tsteps[0].nrecs = nrecords;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
750
  for ( recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
751
    streamptr->tsteps[0].recIDs[recID] = recID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
752

Uwe Schulzweida's avatar
Uwe Schulzweida committed
753
754
  streamptr->record->buffer     = gribbuffer;
  streamptr->record->buffersize = buffersize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
755

Uwe Schulzweida's avatar
Uwe Schulzweida committed
756
  if ( streamptr->ntsteps == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
757
758
    {
      tsID = tstepsNewEntry(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
759
      if ( tsID != streamptr->rtsteps )
760
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
761

Uwe Schulzweida's avatar
Uwe Schulzweida committed
762
763
      streamptr->tsteps[tsID-1].next   = TRUE;
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
764
765
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
766
  if ( streamptr->ntsteps == 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
767
768
769
    {
      if ( taxis->vdate == 0 && taxis->vtime == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
770
771
	  streamptr->ntsteps = 0;
	  for ( varID = 0; varID < streamptr->nvars; varID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
772
773
774
775
776
	    {
	      vlistDefVarTime(vlistID, varID, TIME_CONSTANT);
	    }
	}
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
777
778
779
#else
  Error("CGRIBEX support not compiled in!");
#endif
780
781

  return (0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
782
783
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
784
785

int cgribexScanTimestep2(int streamID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
786
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
787
788
  int rstatus = 0;
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
789
790
  int *isec0, *isec1, *isec2, *isec3, *isec4;
  double fsec2[512], fsec3[2], *fsec4 = NULL;
791
  int lmv = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
792
793
794
795
  off_t recpos = 0;
  unsigned char *gribbuffer = NULL;
  long buffersize = 0;
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
796
797
  int param = 0;
  int level1 = 0, level2 = 0, vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
798
  DateTime datetime, datetime0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
799
800
801
802
803
804
  int tsID;
  int varID, gridID;
  size_t readsize;
  int nrecords, nrecs, recID, rindex;
  long recsize = 0;
  int warn_numavg = TRUE;
805
  int tsteptype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
806
  int taxisID = -1;
807
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
808
  int vlistID;
809
  long unzipsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
810
811
  compvar_t compVar, compVar0;
  stream_t *streamptr;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
812

Uwe Schulzweida's avatar
Uwe Schulzweida committed
813
  streamptr = stream_to_pointer(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
814

815
  stream_check_ptr(__func__, streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
816

Uwe Schulzweida's avatar
Uwe Schulzweida committed
817
818
819
820
821
822
823
  streamptr->curTsID = 1;

  isec0 = streamptr->record->sec0;
  isec1 = streamptr->record->sec1;
  isec2 = streamptr->record->sec2;
  isec3 = streamptr->record->sec3;
  isec4 = streamptr->record->sec4;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
824
825
826
827
828

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
829
830
  gribbuffer = (unsigned char *) streamptr->record->buffer;
  buffersize = streamptr->record->buffersize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
831
                                          
Uwe Schulzweida's avatar
Uwe Schulzweida committed
832
  tsID = streamptr->rtsteps;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
833
  if ( tsID != 1 )
834
    Error("Internal problem! unexpeceted timestep %d", tsID+1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
835

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
838
  fileSetPos(fileID, streamptr->tsteps[tsID].position, SEEK_SET);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
839
840
841

  cdiCreateRecords(streamID, tsID);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
842
843
844
  nrecords = streamptr->tsteps[tsID].nallrecs;
  streamptr->tsteps[1].recIDs = (int *) malloc(nrecords*sizeof(int));
  streamptr->tsteps[1].nrecs = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
845
  for ( recID = 0; recID < nrecords; recID++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
846
    streamptr->tsteps[1].recIDs[recID] = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
847
848
849
      
  for ( recID = 0; recID < nrecords; recID++ )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
850
851
852
853
854
      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;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
855
856
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
857
858
  rindex = 0;
  while ( TRUE )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
859
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
860
861
      if ( rindex > nrecords ) break;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
862
863
864
865
      recsize = gribGetSize(fileID);
      recpos  = fileGetPos(fileID);
      if ( recsize == 0 )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
866
	  streamptr->ntsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
867
868
869
870
871
	  break;
	}
      if ( recsize > buffersize )
	{
	  buffersize = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
872
	  gribbuffer = (unsigned char *) realloc(gribbuffer, (size_t)buffersize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
873
874
875
	}

      readsize = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
876
877
      rstatus = gribRead(fileID, gribbuffer, &readsize);
      if ( rstatus ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
878

879
880
881
882
883
884
885
886
887
888
      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);
	    }
	}

889
890
      cgribexDecodeHeader(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, fsec4,
			  (int *) gribbuffer, recsize, &lmv);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
891

Uwe Schulzweida's avatar
Uwe Schulzweida committed
892
      param = cdiEncodeParam(ISEC1_Parameter, ISEC1_CodeTable, 255);
893
894
      if ( ISEC1_LevelType == 100 ) ISEC1_Level1 *= 100;
      if ( ISEC1_LevelType ==  99 ) ISEC1_LevelType = 100;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
      level1    = ISEC1_Level1;
      level2    = ISEC1_Level2;

      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;
	    }
912
	  taxis->unit  = cgribexGetTimeUnit(isec1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
913
914
	  taxis->vdate = vdate;
	  taxis->vtime = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
915
916
917

	  datetime0.date = vdate;
	  datetime0.time = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
918
919
	}

920
921
      tsteptype = cgribexGetTsteptype(ISEC1_TimeRange);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
922
923
924
925
926
927
      if ( ISEC1_AvgNum )
	{
	  if (  taxis->numavg && warn_numavg &&
		(taxis->numavg != ISEC1_AvgNum) )
	    {
	  /*
928
	      Warning("Changing numavg from %d to %d not supported!",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
929
930
931
932
933
934
935
936
937
938
		      taxis->numavg, ISEC1_AvgNum);
	  */
	      warn_numavg = FALSE;
	    }
	  else
	    {
	      taxis->numavg = ISEC1_AvgNum;
	    }
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
939
940
      datetime.date  = vdate;
      datetime.time  = vtime;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
941
      compVar.param  = param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
942
943
944
945
946
      compVar.level1 = level1;
      compVar.level2 = level2;
      compVar.ltype  = ISEC1_LevelType;
      for ( recID = 0; recID < nrecords; recID++ )
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
947
	  compVar0.param  = streamptr->tsteps[tsID].records[recID].param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
948
949
950
	  compVar0.level1 = streamptr->tsteps[tsID].records[recID].ilevel;
	  compVar0.level2 = streamptr->tsteps[tsID].records[recID].ilevel2;
	  compVar0.ltype  = streamptr->tsteps[tsID].records[recID].ltype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
951

Uwe Schulzweida's avatar
Uwe Schulzweida committed
952
	  if ( memcmp(&compVar0, &compVar, sizeof(compvar_t)) == 0 ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
953
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
954

Uwe Schulzweida's avatar
Uwe Schulzweida committed
955
956
      if ( recID == nrecords )
	{
957
958
	  char paramstr[32];
	  cdiParamToString(param, paramstr, sizeof(paramstr));
959
	  Warning("Param=%s level=%d not defined at timestep 1!", paramstr, level1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
960
961
962
	  return (CDI_EUFSTRUCT);
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
      if ( cdiInventoryMode == 1 )
	{
	  if ( streamptr->tsteps[tsID].records[recID].used )
	    {
	      break;
	    }
	  else
	    {
	      streamptr->tsteps[tsID].records[recID].used = TRUE;
	      streamptr->tsteps[tsID].recIDs[rindex] = recID;
	    }    
	}
      else
	{
	  if ( streamptr->tsteps[tsID].records[recID].used )
	    {
979
980
981
	      char paramstr[32];
	      cdiParamToString(param, paramstr, sizeof(paramstr));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
982
983
	      if ( memcmp(&datetime, &datetime0, sizeof(DateTime)) != 0 ) break;

984
	      Warning("Param=%s level=%d already exist, skipped!", paramstr, level1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
985
986
987
988
989
990
991
992
	      continue;
	    }
	  else
	    {
	      streamptr->tsteps[tsID].records[recID].used = TRUE;
	      streamptr->tsteps[tsID].recIDs[rindex] = recID;
	    }    
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
993
994

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

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
999
      compVar0.param  = streamptr->tsteps[tsID].records[recID].param;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1000
1001
1002
      compVar0.level1 = streamptr->tsteps[tsID].records[recID].ilevel;
      compVar0.level2 = streamptr->tsteps[tsID].records[recID].ilevel2;
      compVar0.ltype  = streamptr->tsteps[tsID].records[recID].ltype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1003

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1004
      if ( memcmp(&compVar0, &compVar, sizeof(compvar_t)) != 0 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1005
	{
1006
	  Message("tsID = %d recID = %d param = %3d new %3d  level = %3d new %3d",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1007
		  tsID, recID,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1008
		  streamptr->tsteps[tsID].records[recID].param, param,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1009
		  streamptr->tsteps[tsID].records[recID].ilevel, level1);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1010
1011
1012
	  return (CDI_EUFSTRUCT);
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1013
1014
      streamptr->tsteps[1].records[recID].position = recpos;
      varID = streamptr->tsteps[tsID].records[recID].varID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1015
1016
1017
      gridID = vlistInqVarGrid(vlistID, varID);
      if ( gridInqSize(gridID) == 1 && gridInqType(gridID) == GRID_LONLAT )
	{
1018
1019
	  if ( IS_NOT_EQUAL(gridInqXval(gridID, 0),ISEC2_FirstLon*0.001) ||
	       IS_NOT_EQUAL(gridInqYval(gridID, 0),ISEC2_FirstLat*0.001) )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1020
1021
	    gridChangeType(gridID, GRID_TRAJECTORY);
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1022

1023
1024
1025
      if ( tsteptype != vlistInqVarTsteptype(vlistID, varID) )
	vlistDefVarTsteptype(vlistID, varID, tsteptype);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1026
      rindex++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1027
1028
1029
1030
1031
    }

  nrecs = 0;
  for ( recID = 0; recID < nrecords; recID++ )
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1032
      if ( ! streamptr->tsteps[tsID].records[recID].used )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1033
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1034
	  varID = streamptr->tsteps[tsID].records[recID].varID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1035
1036
1037
1038
1039
1040
1041
	  vlistDefVarTime(vlistID, varID, TIME_CONSTANT);
	}
      else
	{
	  nrecs++;
	}
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1042
  streamptr->tsteps[tsID].nrecs = nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1043

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1044
  streamptr->rtsteps = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1045

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1046
  if ( streamptr->ntsteps == -1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1047
1048
    {
      tsID = tstepsNewEntry(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1049
      if ( tsID != streamptr->rtsteps )
1050
	Error("Internal error. tsID = %d", tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1051

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1052
1053
      streamptr->tsteps[tsID-1].next   = TRUE;
      streamptr->tsteps[tsID].position = recpos;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1054
1055
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1056
1057
  streamptr->record->buffer     = gribbuffer;
  streamptr->record->buffersize = buffersize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1058
#endif
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1059

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1060
  return (rstatus);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1061
1062
1063
}


Uwe Schulzweida's avatar
Uwe Schulzweida committed
1064
int cgribexScanTimestep(int streamID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1065
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1066
1067
  int rstatus = 0;
#if  defined  (HAVE_LIBCGRIBEX)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1068
1069
  int *isec0, *isec1, *isec2, *isec3, *isec4;
  double fsec2[512], fsec3[2], *fsec4 = NULL;
1070
  int lmv = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1071
1072
1073
1074
  long recsize = 0;
  off_t recpos = 0;
  unsigned char *gribbuffer;
  long buffersize = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1075
  int fileID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1076
1077
  int param = 0;
  int level1 = 0, level2 = 0, vdate = 0, vtime = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1078
1079
1080
1081
1082
1083
  DateTime datetime, datetime0;
  int tsID;
  int vrecID, recID;
  int warn_numavg = TRUE;
  size_t readsize;
  int taxisID = -1;
1084
  taxis_t *taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1085
1086
1087
1088
  int vlistID;
  int rindex, nrecs = 0;
  long unzipsize;
  compvar_t compVar, compVar0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1089
  stream_t *streamptr;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1090
1091

  streamptr = stream_to_pointer(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1092

1093
  stream_check_ptr(__func__, streamptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1094

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1095
  vlistID = streamInqVlist(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1096

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1097
1098
  if ( CDI_Debug )
    {
1099
1100
1101
1102
      Message("streamID = %d", streamID);
      Message("cts = %d", streamptr->curTsID);
      Message("rts = %d", streamptr->rtsteps);
      Message("nts = %d", streamptr->ntsteps);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1103
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1104

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1105
1106
1107
1108
1109
  isec0 = streamptr->record->sec0;
  isec1 = streamptr->record->sec1;
  isec2 = streamptr->record->sec2;
  isec3 = streamptr->record->sec3;
  isec4 = streamptr->record->sec4;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1110

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1111
1112
  tsID  = streamptr->rtsteps;
  taxis = &streamptr->tsteps[tsID].taxis;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1113

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1114
1115
1116
1117
  if ( streamptr->tsteps[tsID].recordSize == 0 )
    {
      gribbuffer = (unsigned char *) streamptr->record->buffer;
      buffersize = streamptr->record->buffersize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1118

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1119
      cdiCreateRecords(streamID, tsID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1120

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1121
      nrecs = streamptr->tsteps[1].nrecs;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1122

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1123
1124
1125
1126
      streamptr->tsteps[tsID].nrecs = nrecs;
      streamptr->tsteps[tsID].recIDs = (int *) malloc(nrecs*sizeof(int));
      for ( recID = 0; recID < nrecs; recID++ )
	streamptr->tsteps[tsID].recIDs[recID] = streamptr->tsteps[1].recIDs[recID];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1127

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1128
      fileID = streamInqFileID(streamID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1129

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

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1132
1133
      rindex = 0;
      while ( TRUE )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1134
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1135
	  if ( rindex > nrecs ) break;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1136

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1137
1138
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
	  recsize = gribGetSize(fileID);
	  recpos  = fileGetPos(fileID);
	  if ( recsize == 0 )
	    {
	      streamptr->ntsteps = streamptr->rtsteps + 1;
	      break;
	    }
	  if ( recsize > buffersize )
	    {
	      buffersize = recsize;
	      gribbuffer = (unsigned char *) realloc(gribbuffer, buffersize);
	    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1149

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1150
	  readsize = recsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1151
1152
	  rstatus = gribRead(fileID, gribbuffer, &readsize);
	  if ( rstatus )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1153
	    {
1154
	      Error("Inconsistent timestep %d (GRIB record %d/%d)!\n", tsID+1, rindex+1,