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

Thomas Jahns's avatar
Thomas Jahns committed
5
#include <assert.h>
Uwe Schulzweida's avatar
Uwe Schulzweida committed
6
#include <string.h>
7
#include <float.h>  /* FLT_EPSILON */
8
#include <limits.h> /* INT_MAX     */
9

Uwe Schulzweida's avatar
Uwe Schulzweida committed
10
11
#include "dmemory.h"
#include "cdi.h"
12
#include "cdi_cksum.h"
13
#include "cdi_int.h"
14
#include "cdi_uuid.h"
Uwe Schulzweida's avatar
Uwe Schulzweida committed
15
#include "grid.h"
16
#include "gaussgrid.h"
17
#include "resource_handle.h"
18
#include "resource_unpack.h"
19
#include "namespace.h"
20
#include "serialize.h"
21
22
23
24
#include "vlist.h"

#undef  UNDEFID
#define UNDEFID -1
Uwe Schulzweida's avatar
Uwe Schulzweida committed
25

26
27
28
/* the value in the second pair of brackets must match the length of
 * the longest string (including terminating NUL) */
static const char Grids[][17] = {
29
30
31
32
33
34
35
36
37
  /*  0 */  "undefined",
  /*  1 */  "generic",
  /*  2 */  "gaussian",
  /*  3 */  "gaussian reduced",
  /*  4 */  "lonlat",
  /*  5 */  "spectral",
  /*  6 */  "fourier",
  /*  7 */  "gme",
  /*  8 */  "trajectory",
38
  /*  9 */  "unstructured",
39
40
  /* 10 */  "curvilinear",
  /* 11 */  "lcc",
41
42
43
  /* 12 */  "lcc2",
  /* 13 */  "laea",
  /* 14 */  "sinusoidal",
44
  /* 15 */  "projection",
Uwe Schulzweida's avatar
Uwe Schulzweida committed
45
46
};

47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
/* must match table below */
enum xystdname_idx {
  grid_xystdname_grid_latlon,
  grid_xystdname_latlon,
  grid_xystdname_projection,
};
static const char xystdname_tab[][2][24] = {
  [grid_xystdname_grid_latlon] = { "grid_longitude",
                                   "grid_latitude" },
  [grid_xystdname_latlon] = { "longitude",
                              "latitude" },
  [grid_xystdname_projection] = { "projection_x_coordinate",
                                  "projection_y_coordinate" },

};


Uwe Schulzweida's avatar
Uwe Schulzweida committed
64

65
66
67
static int    gridCompareP    ( void * gridptr1, void * gridptr2 );
static void   gridDestroyP    ( void * gridptr );
static void   gridPrintP      ( void * gridptr, FILE * fp );
68
static int    gridGetPackSize ( void * gridptr, void *context);
69
static void   gridPack        ( void * gridptr, void * buff, int size,
70
				int *position, void *context);
Deike Kleberg's avatar
minimal    
Deike Kleberg committed
71
static int    gridTxCode      ( void );
Uwe Schulzweida's avatar
Uwe Schulzweida committed
72

73
static const resOps gridOps = {
74
75
76
77
78
79
  gridCompareP,
  gridDestroyP,
  gridPrintP,
  gridGetPackSize,
  gridPack,
  gridTxCode
80
};
Uwe Schulzweida's avatar
Uwe Schulzweida committed
81

82
static int  GRID_Debug = 0;   /* If set to 1, debugging */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
83

84
85
86
87
grid_t *gridID2Ptr(int gridID)
{
  return (grid_t *)reshGetVal(gridID, &gridOps);
}
88
#define gridID2Ptr(gridID) (grid_t *)reshGetVal(gridID, &gridOps)
89
#define gridMark4Update(gridID) reshSetStatus(gridID, &gridOps, RESH_DESYNC_IN_USE)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
90

91

Deike Kleberg's avatar
Deike Kleberg committed
92
void grid_init(grid_t *gridptr)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
93
{
94
  gridptr->self         = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
95
  gridptr->type         = CDI_UNDEFID;
Deike Kleberg's avatar
Deike Kleberg committed
96
  gridptr->proj         = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
97
  gridptr->mask         = NULL;
98
  gridptr->mask_gme     = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
99
100
101
102
103
104
105
  gridptr->xvals        = NULL;
  gridptr->yvals        = NULL;
  gridptr->area         = NULL;
  gridptr->xbounds      = NULL;
  gridptr->ybounds      = NULL;
  gridptr->rowlon       = NULL;
  gridptr->nrowlon      = 0;
106
107
  gridptr->xfirst       = 0.0;
  gridptr->xlast        = 0.0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
108
  gridptr->xinc         = 0.0;
109
110
  gridptr->yfirst       = 0.0;
  gridptr->ylast        = 0.0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
111
  gridptr->yinc         = 0.0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
112
113
114
115
116
117
118
  gridptr->lcc_originLon = 0.0;
  gridptr->lcc_originLat = 0.0;
  gridptr->lcc_lonParY  = 0.0;
  gridptr->lcc_lat1     = 0.0;
  gridptr->lcc_lat2     = 0.0;
  gridptr->lcc_xinc     = 0.0;
  gridptr->lcc_yinc     = 0.0;
119
120
  gridptr->lcc_projflag = 0;
  gridptr->lcc_scanflag = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
121
  gridptr->lcc_defined  = FALSE;
122
123
124
125
126
127
  gridptr->lcc2_lon_0   = 0.0;
  gridptr->lcc2_lat_0   = 0.0;
  gridptr->lcc2_lat_1   = 0.0;
  gridptr->lcc2_lat_2   = 0.0;
  gridptr->lcc2_a       = 0.0;
  gridptr->lcc2_defined = FALSE;
128
129
130
131
  gridptr->laea_lon_0   = 0.0;
  gridptr->laea_lat_0   = 0.0;
  gridptr->laea_a       = 0.0;
  gridptr->laea_defined = FALSE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
132
133
134
135
136
137
  gridptr->trunc        = 0;
  gridptr->nvertex      = 0;
  gridptr->nd           = 0;
  gridptr->ni           = 0;
  gridptr->ni2          = 0;
  gridptr->ni3          = 0;
138
139
140
  gridptr->number       = 0;
  gridptr->position     = 0;
  gridptr->reference    = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
141
142
143
144
  gridptr->prec         = 0;
  gridptr->size         = 0;
  gridptr->xsize        = 0;
  gridptr->ysize        = 0;
145
  gridptr->np           = 0;
146
147
  gridptr->xdef         = 0;
  gridptr->ydef         = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
148
  gridptr->isCyclic     = CDI_UNDEFID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
149
150
151
152
  gridptr->isRotated    = FALSE;
  gridptr->xpole        = 0.0;
  gridptr->ypole        = 0.0;
  gridptr->angle        = 0.0;
153
  gridptr->lcomplex     = 0;
154
  gridptr->hasdims      = TRUE;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
155
156
157
158
159
160
  gridptr->xname[0]     = 0;
  gridptr->yname[0]     = 0;
  gridptr->xlongname[0] = 0;
  gridptr->ylongname[0] = 0;
  gridptr->xunits[0]    = 0;
  gridptr->yunits[0]    = 0;
161
162
  gridptr->xstdname  = NULL;
  gridptr->ystdname  = NULL;
Thomas Jahns's avatar
Thomas Jahns committed
163
  memset(gridptr->uuid, 0, CDI_UUID_SIZE);
Deike Kleberg's avatar
Deike Kleberg committed
164
  gridptr->name         = NULL;
165
  gridptr->vtable       = &cdiGridVtable;
166
  gridptr->extraData    = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
167
168
}

169

170
171
static void
grid_free_components(grid_t *gridptr)
Deike Kleberg's avatar
Deike Kleberg committed
172
{
173
174
175
176
177
178
179
180
  void *p2free[] = { gridptr->mask, gridptr->mask_gme,
                   gridptr->xvals, gridptr->yvals,
                   gridptr->xbounds, gridptr->ybounds,
                   gridptr->rowlon, gridptr->area,
                   gridptr->reference, gridptr->name };
  for (size_t i = 0; i < sizeof (p2free) / sizeof (p2free[0]); ++i)
    if (p2free[i]) Free(p2free[i]);
}
Deike Kleberg's avatar
Deike Kleberg committed
181

182
183
184
void grid_free(grid_t *gridptr)
{
  grid_free_components(gridptr);
Deike Kleberg's avatar
Deike Kleberg committed
185
186
187
  grid_init(gridptr);
}

188
189
static grid_t *
gridNewEntry(cdiResH resH)
190
{
191
  grid_t *gridptr = (grid_t*) Malloc(sizeof(grid_t));
192
193
194
195
196
197
198
199
  grid_init(gridptr);
  if (resH == CDI_UNDEFID)
    gridptr->self = reshPut(gridptr, &gridOps);
  else
    {
      gridptr->self = resH;
      reshReplace(resH, gridptr, &gridOps);
    }
200
  return gridptr;
201
202
}

203
static
204
void gridInit (void)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
205
{
206
207
  static int gridInitialized = 0;
  char *env;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
208

209
  if ( gridInitialized ) return;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
210

Thomas Jahns's avatar
Thomas Jahns committed
211
  gridInitialized = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
212

213
214
  env = getenv("GRID_DEBUG");
  if ( env ) GRID_Debug = atoi(env);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
215
216
}

217
218
static
void grid_copy(grid_t *gridptr2, grid_t *gridptr1)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
219
220
221
222
{
  int gridID2;

  gridID2 = gridptr2->self;
223
  memcpy(gridptr2, gridptr1, sizeof(grid_t));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
224
225
226
  gridptr2->self = gridID2;
}

227
unsigned cdiGridCount(void)
Thomas Jahns's avatar
Thomas Jahns committed
228
{
229
  return reshCountType(&gridOps);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
230
231
}

232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
static inline void
gridSetXname(grid_t *gridptr, const char *xname)
{
  strncpy(gridptr->xname, xname, CDI_MAX_NAME);
  gridptr->xname[CDI_MAX_NAME - 1] = 0;
}

static inline void
gridSetXlongname(grid_t *gridptr, const char *xlongname)
{
  strncpy(gridptr->xlongname, xlongname, CDI_MAX_NAME);
  gridptr->xlongname[CDI_MAX_NAME - 1] = 0;
}

static inline void
gridSetXunits(grid_t *gridptr, const char *xunits)
{
  strncpy(gridptr->xunits, xunits, CDI_MAX_NAME);
  gridptr->xunits[CDI_MAX_NAME - 1] = 0;
}

static inline void
gridSetYname(grid_t *gridptr, const char *yname)
{
  strncpy(gridptr->yname, yname, CDI_MAX_NAME);
  gridptr->yname[CDI_MAX_NAME - 1] = 0;
}

static inline void
gridSetYlongname(grid_t *gridptr, const char *ylongname)
{
  strncpy(gridptr->ylongname, ylongname, CDI_MAX_NAME);
  gridptr->ylongname[CDI_MAX_NAME - 1] = 0;
}

static inline void
gridSetYunits(grid_t *gridptr, const char *yunits)
{
  strncpy(gridptr->yunits, yunits, CDI_MAX_NAME);
  gridptr->yunits[CDI_MAX_NAME - 1] = 0;
}

274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
void
cdiGridTypeInit(grid_t *gridptr, int gridtype, int size)
{

  gridptr->type = gridtype;
  gridptr->size = size;

  switch (gridtype)
    {
    case GRID_CURVILINEAR:
      gridptr->nvertex = 4;
      /* Fall through */
    case GRID_LONLAT:
    case GRID_GAUSSIAN:
    case GRID_GAUSSIAN_REDUCED:
    case GRID_TRAJECTORY:
      {
        if ( gridtype == GRID_TRAJECTORY )
          {
            gridSetXname(gridptr, "tlon");
            gridSetYname(gridptr, "tlat");
          }
        else
          {
            gridSetXname(gridptr, "lon");
            gridSetYname(gridptr, "lat");
          }
        gridSetXlongname(gridptr, "longitude");
        gridSetYlongname(gridptr, "latitude");

        /*
        if ( gridtype == GRID_CURVILINEAR )
          {
            gridptr->xstdname = xystdname_tab[grid_xystdname_grid_latlon][0];
            gridptr->ystdname = xystdname_tab[grid_xystdname_grid_latlon][1];
            gridDefXunits(gridID, "degrees");
            gridDefYunits(gridID, "degrees");
          }
        else
        */
          {
            gridptr->xstdname = xystdname_tab[grid_xystdname_latlon][0];
            gridptr->ystdname = xystdname_tab[grid_xystdname_latlon][1];
            gridSetXunits(gridptr, "degrees_east");
            gridSetYunits(gridptr, "degrees_north");
          }

        break;
      }
    case GRID_UNSTRUCTURED:
      gridptr->xsize = size;
      /* Fall through */
    case GRID_GME:
      {
        gridSetXname(gridptr, "lon");
        gridSetYname(gridptr, "lat");
        gridptr->xstdname = xystdname_tab[grid_xystdname_latlon][0];
        gridptr->ystdname = xystdname_tab[grid_xystdname_latlon][1];
        gridSetXunits(gridptr, "degrees_east");
        gridSetYunits(gridptr, "degrees_north");
        break;
      }
    case GRID_GENERIC:
      {

        /* gridptr->xsize = size; */
        gridSetXname(gridptr, "x");
        gridSetYname(gridptr, "y");
        /*
        strcpy(gridptr->xstdname, "grid_longitude");
        strcpy(gridptr->ystdname, "grid_latitude");
        gridptr->xstdname = xystdname_tab[grid_xystdname_grid_latlon][0];
        gridptr->ystdname = xystdname_tab[grid_xystdname_grid_latlon][1];
        gridDefXunits(gridID, "degrees");
        gridDefYunits(gridID, "degrees");
        */
        break;
      }
    case GRID_LCC2:
    case GRID_SINUSOIDAL:
    case GRID_LAEA:
      {
        gridSetXname(gridptr, "x");
        gridSetYname(gridptr, "y");
        gridptr->xstdname = xystdname_tab[grid_xystdname_projection][0];
        gridptr->ystdname = xystdname_tab[grid_xystdname_projection][1];
        gridSetXunits(gridptr, "m");
        gridSetYunits(gridptr, "m");
        break;
      }
    }

}


369
// used also in CDO
370
void gridGenXvals(int xsize, double xfirst, double xlast, double xinc, double *restrict xvals)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
371
{
372
  if ( (! (fabs(xinc) > 0)) && xsize > 1 )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
373
374
    {
      if ( xfirst >= xlast )
Thomas Jahns's avatar
Thomas Jahns committed
375
376
377
378
        {
          while ( xfirst >= xlast ) xlast += 360;
          xinc = (xlast-xfirst)/(xsize);
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
379
      else
Thomas Jahns's avatar
Thomas Jahns committed
380
381
382
        {
          xinc = (xlast-xfirst)/(xsize-1);
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
383
384
    }

385
  for ( int i = 0; i < xsize; ++i )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
386
387
388
    xvals[i] = xfirst + i*xinc;
}

389
static
390
void calc_gaussgrid(double *restrict yvals, int ysize, double yfirst, double ylast)
391
{
392
  double *restrict yw = (double *) Malloc((size_t)ysize * sizeof(double));
393
  gaussaw(yvals, yw, (size_t)ysize);
394
  Free(yw);
395
  for (int i = 0; i < ysize; i++ )
396
397
398
399
    yvals[i] = asin(yvals[i])/M_PI*180.0;

  if ( yfirst < ylast && yfirst > -90.0 && ylast < 90.0 )
    {
400
401
      int yhsize = ysize/2;
      for (int i = 0; i < yhsize; i++ )
Thomas Jahns's avatar
Thomas Jahns committed
402
        {
403
          double ytmp = yvals[i];
Thomas Jahns's avatar
Thomas Jahns committed
404
405
406
          yvals[i] = yvals[ysize-i-1];
          yvals[ysize-i-1] = ytmp;
        }
407
408
409
    }
}

410
// used also in CDO
Thomas Jahns's avatar
Thomas Jahns committed
411
void gridGenYvals(int gridtype, int ysize, double yfirst, double ylast, double yinc, double *restrict yvals)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
412
{
413
  const double deleps = 0.002;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
414
415
416
417

  if ( gridtype == GRID_GAUSSIAN || gridtype == GRID_GAUSSIAN_REDUCED )
    {
      if ( ysize > 2 )
Deike Kleberg's avatar
Deike Kleberg committed
418
419
420
421
422
423
	{
	  calc_gaussgrid(yvals, ysize, yfirst, ylast);

	  if ( ! (IS_EQUAL(yfirst, 0) && IS_EQUAL(ylast, 0)) )
	    if ( fabs(yvals[0] - yfirst) > deleps || fabs(yvals[ysize-1] - ylast) > deleps )
	      {
424
		double *restrict ytmp = NULL;
Deike Kleberg's avatar
Deike Kleberg committed
425
		int nstart, lfound = 0;
426
		int ny = (int) (180./fabs(ylast-yfirst)/(ysize-1) + 0.5);
Deike Kleberg's avatar
Deike Kleberg committed
427
428
429
		ny -= ny%2;
		if ( ny > ysize && ny < 4096 )
		  {
430
		    ytmp = (double *) Malloc((size_t)ny * sizeof (double));
Deike Kleberg's avatar
Deike Kleberg committed
431
		    calc_gaussgrid(ytmp, ny, yfirst, ylast);
432
433
434
435
436
437
                    {
                      int i;
                      for ( i = 0; i < (ny-ysize); i++ )
                        if ( fabs(ytmp[i] - yfirst) < deleps ) break;
                      nstart = i;
                    }
Deike Kleberg's avatar
Deike Kleberg committed
438

439
440
		    lfound = (nstart+ysize-1) < ny
                      && fabs(ytmp[nstart+ysize-1] - ylast) < deleps;
441
442
443
444
                    if ( lfound )
                      {
                        for (int i = 0; i < ysize; i++) yvals[i] = ytmp[i+nstart];
                      }
Deike Kleberg's avatar
Deike Kleberg committed
445
446
		  }

447
		if ( !lfound )
Deike Kleberg's avatar
Deike Kleberg committed
448
449
		  {
		    Warning("Cannot calculate gaussian latitudes for lat1 = %g latn = %g!", yfirst, ylast);
450
		    for (int i = 0; i < ysize; i++ ) yvals[i] = 0;
Deike Kleberg's avatar
Deike Kleberg committed
451
452
453
454
		    yvals[0] = yfirst;
		    yvals[ysize-1] = ylast;
		  }

455
		if ( ytmp ) Free(ytmp);
Deike Kleberg's avatar
Deike Kleberg committed
456
457
	      }
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
458
      else
Thomas Jahns's avatar
Thomas Jahns committed
459
460
461
462
        {
          yvals[0] = yfirst;
          yvals[ysize-1] = ylast;
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
463
464
465
466
    }
  /*     else if ( gridtype == GRID_LONLAT || gridtype == GRID_GENERIC ) */
  else
    {
467
      if ( (! (fabs(yinc) > 0)) && ysize > 1 )
Thomas Jahns's avatar
Thomas Jahns committed
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
        {
          if ( IS_EQUAL(yfirst, ylast) && IS_NOT_EQUAL(yfirst, 0) ) ylast *= -1;

          if ( yfirst > ylast )
            yinc = (yfirst-ylast)/(ysize-1);
          else if ( yfirst < ylast )
            yinc = (ylast-yfirst)/(ysize-1);
          else
            {
              if ( ysize%2 != 0 )
                {
                  yinc = 180.0/(ysize-1);
                  yfirst = -90;
                }
              else
                {
                  yinc = 180.0/ysize;
                  yfirst = -90 + yinc/2;
                }
            }
        }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
489
490
491

      if ( yfirst > ylast && yinc > 0 ) yinc = -yinc;

492
      for (int i = 0; i < ysize; i++ )
Thomas Jahns's avatar
Thomas Jahns committed
493
        yvals[i] = yfirst + i*yinc;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
494
495
496
    }
  /*
    else
497
    Error("unable to calculate values for %s grid!", gridNamePtr(gridtype));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
498
499
500
501
  */
}

/*
502
@Function  gridCreate
503
@Title     Create a horizontal Grid
Uwe Schulzweida's avatar
Uwe Schulzweida committed
504

505
@Prototype int gridCreate(int gridtype, int size)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
506
507
@Parameter
    @Item  gridtype  The type of the grid, one of the set of predefined CDI grid types.
508
                     The valid CDI grid types are @func{GRID_GENERIC}, @func{GRID_GAUSSIAN},
Uwe Schulzweida's avatar
Uwe Schulzweida committed
509
                     @func{GRID_LONLAT}, @func{GRID_LCC}, @func{GRID_SPECTRAL},
510
                     @func{GRID_GME}, @func{GRID_CURVILINEAR} and @func{GRID_UNSTRUCTURED} and.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
511
512
513
    @Item  size      Number of gridpoints.

@Description
514
The function @func{gridCreate} creates a horizontal Grid.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
515
516

@Result
517
@func{gridCreate} returns an identifier to the Grid.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
518
519

@Example
520
Here is an example using @func{gridCreate} to create a regular lon/lat Grid:
Uwe Schulzweida's avatar
Uwe Schulzweida committed
521
522

@Source
Uwe Schulzweida's avatar
Uwe Schulzweida committed
523
524
#include "cdi.h"
   ...
525
526
#define  nlon  12
#define  nlat   6
Uwe Schulzweida's avatar
Uwe Schulzweida committed
527
   ...
528
529
double lons[nlon] = {0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330};
double lats[nlat] = {-75, -45, -15, 15, 45, 75};
Uwe Schulzweida's avatar
Uwe Schulzweida committed
530
531
int gridID;
   ...
532
533
534
gridID = gridCreate(GRID_LONLAT, nlon*nlat);
gridDefXsize(gridID, nlon);
gridDefYsize(gridID, nlat);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
535
536
537
538
gridDefXvals(gridID, lons);
gridDefYvals(gridID, lats);
   ...
@EndSource
Uwe Schulzweida's avatar
Uwe Schulzweida committed
539
540
@EndFunction
*/
541
int gridCreate(int gridtype, int size)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
542
{
543
  if ( CDI_Debug ) Message("gridtype=%s  size=%d", gridNamePtr(gridtype), size);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
544

545
  if ( size < 0 || size > INT_MAX ) Error("Grid size (%d) out of bounds (0 - %d)!", size, INT_MAX);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
546

547
  gridInit();
548

Uwe Schulzweida's avatar
Uwe Schulzweida committed
549
  grid_t *gridptr = gridNewEntry(CDI_UNDEFID);
550
  if ( ! gridptr ) Error("No memory");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
551

Uwe Schulzweida's avatar
Uwe Schulzweida committed
552
  int gridID = gridptr->self;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
553

554
  if ( CDI_Debug ) Message("gridID: %d", gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
555

556
  cdiGridTypeInit(gridptr, gridtype, size);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
557
558
559
560

  return (gridID);
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
561
562
static
void gridDestroyKernel( grid_t * gridptr )
563
564
565
{
  int id;

Deike Kleberg's avatar
Deike Kleberg committed
566
  xassert ( gridptr );
567
568
569

  id = gridptr->self;

570
  grid_free_components(gridptr);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
571
  Free( gridptr );
572

573
  reshRemove ( id, &gridOps );
574
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
575

576
577
578
579
580
581
/*
@Function  gridDestroy
@Title     Destroy a horizontal Grid

@Prototype void gridDestroy(int gridID)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
582
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
583
584
585
586
587

@EndFunction
*/
void gridDestroy(int gridID)
{
588
  grid_t *gridptr = gridID2Ptr(gridID);
589
  gridptr->vtable->destroy(gridptr);
590
591
}

592
void gridDestroyP ( void * gridptr )
593
{
594
  ((grid_t *)gridptr)->vtable->destroy((grid_t *)gridptr);
595
596
597
}


598
const char *gridNamePtr(int gridtype)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
599
{
600
  int size = (int) (sizeof(Grids)/sizeof(Grids[0]));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
601

602
  const char *name = gridtype >= 0 && gridtype < size ? Grids[gridtype] : Grids[GRID_GENERIC];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618

  return (name);
}


void gridName(int gridtype, char *gridname)
{
  strcpy(gridname, gridNamePtr(gridtype));
}

/*
@Function  gridDefXname
@Title     Define the name of a X-axis

@Prototype void gridDefXname(int gridID, const char *name)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
619
620
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
    @Item  name     Name of the X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
621
622

@Description
623
The function @func{gridDefXname} defines the name of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
624
625
626

@EndFunction
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
627
void gridDefXname(int gridID, const char *xname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
628
{
629
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
630

631
632
  if ( xname )
    {
633
      gridSetXname(gridptr, xname);
634
      gridMark4Update(gridID);
635
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
636
637
}

638

Uwe Schulzweida's avatar
Uwe Schulzweida committed
639
640
641
642
643
644
/*
@Function  gridDefXlongname
@Title     Define the longname of a X-axis

@Prototype void gridDefXlongname(int gridID, const char *longname)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
645
646
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
    @Item  longname Longname of the X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
647
648

@Description
649
The function @func{gridDefXlongname} defines the longname of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
650
651
652

@EndFunction
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
653
void gridDefXlongname(int gridID, const char *xlongname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
654
{
655
  grid_t *gridptr = gridID2Ptr(gridID);
656
657
  if ( xlongname )
    {
658
      gridSetXlongname(gridptr, xlongname);
659
      gridMark4Update(gridID);
660
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
661
662
663
664
665
666
667
668
}

/*
@Function  gridDefXunits
@Title     Define the units of a X-axis

@Prototype void gridDefXunits(int gridID, const char *units)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
669
670
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
    @Item  units    Units of the X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
671
672

@Description
673
The function @func{gridDefXunits} defines the units of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
674
675
676

@EndFunction
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
677
void gridDefXunits(int gridID, const char *xunits)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
678
{
679
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
680

681
682
  if ( xunits )
    {
683
      gridSetXunits(gridptr, xunits);
684
      gridMark4Update(gridID);
685
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
686
687
688
689
690
691
692
693
}

/*
@Function  gridDefYname
@Title     Define the name of a Y-axis

@Prototype void gridDefYname(int gridID, const char *name)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
694
695
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
    @Item  name     Name of the Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
696
697

@Description
698
The function @func{gridDefYname} defines the name of a Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
699
700
701

@EndFunction
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
702
void gridDefYname(int gridID, const char *yname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
703
{
704
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
705

706
707
  if ( yname )
    {
708
      gridSetYname(gridptr, yname);
709
      gridMark4Update(gridID);
710
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
711
712
713
714
715
716
717
718
}

/*
@Function  gridDefYlongname
@Title     Define the longname of a Y-axis

@Prototype void gridDefYlongname(int gridID, const char *longname)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
719
720
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
    @Item  longname Longname of the Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
721
722

@Description
723
The function @func{gridDefYlongname} defines the longname of a Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
724
725
726

@EndFunction
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
727
void gridDefYlongname(int gridID, const char *ylongname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
728
{
729
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
730

731
732
  if ( ylongname )
    {
733
      gridSetYlongname(gridptr, ylongname);
734
      gridMark4Update(gridID);
735
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
736
737
738
739
740
741
742
743
}

/*
@Function  gridDefYunits
@Title     Define the units of a Y-axis

@Prototype void gridDefYunits(int gridID, const char *units)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
744
745
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
    @Item  units    Units of the Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
746
747

@Description
748
The function @func{gridDefYunits} defines the units of a Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
749
750
751

@EndFunction
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
752
void gridDefYunits(int gridID, const char *yunits)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
753
{
754
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
755

756
757
  if ( yunits )
    {
758
      gridSetYunits(gridptr, yunits);
759
      gridMark4Update(gridID);
760
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
761
762
763
764
765
766
}

/*
@Function  gridInqXname
@Title     Get the name of a X-axis

Deike Kleberg's avatar
Deike Kleberg committed
767
@Prototype void gridInqXname(int gridID, char *name)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
768
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
769
770
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
    @Item  name     Name of the X-axis. The caller must allocate space for the
Deike Kleberg's avatar
Deike Kleberg committed
771
                    returned string. The maximum possible length, in characters, of
772
                    the string is given by the predefined constant @func{CDI_MAX_NAME}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
773
774

@Description
775
The function @func{gridInqXname} returns the name of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
776
777

@Result
778
@func{gridInqXname} returns the name of the X-axis to the parameter name.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
779
780
781
782
783

@EndFunction
*/
void gridInqXname(int gridID, char *xname)
{
784
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
785
786
787
788
789
790
791
792

  strcpy(xname, gridptr->xname);
}

/*
@Function  gridInqXlongname
@Title     Get the longname of a X-axis

Deike Kleberg's avatar
Deike Kleberg committed
793
@Prototype void gridInqXlongname(int gridID, char *longname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
794
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
795
796
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
    @Item  longname Longname of the X-axis. The caller must allocate space for the
Deike Kleberg's avatar
Deike Kleberg committed
797
                    returned string. The maximum possible length, in characters, of
798
                    the string is given by the predefined constant @func{CDI_MAX_NAME}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
799
800

@Description
801
The function @func{gridInqXlongname} returns the longname of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
802
803

@Result
804
@func{gridInqXlongname} returns the longname of the X-axis to the parameter longname.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
805
806
807
808
809

@EndFunction
*/
void gridInqXlongname(int gridID, char *xlongname)
{
810
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
811
812
813
814
815
816
817
818

  strcpy(xlongname, gridptr->xlongname);
}

/*
@Function  gridInqXunits
@Title     Get the units of a X-axis

Deike Kleberg's avatar
Deike Kleberg committed
819
@Prototype void gridInqXunits(int gridID, char *units)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
820
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
821
822
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
    @Item  units    Units of the X-axis. The caller must allocate space for the
Deike Kleberg's avatar
Deike Kleberg committed
823
                    returned string. The maximum possible length, in characters, of
824
                    the string is given by the predefined constant @func{CDI_MAX_NAME}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
825
826

@Description
827
The function @func{gridInqXunits} returns the units of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
828
829

@Result
830
@func{gridInqXunits} returns the units of the X-axis to the parameter units.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
831
832
833
834
835

@EndFunction
*/
void gridInqXunits(int gridID, char *xunits)
{
836
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
837
838
839
840
841
842
843

  strcpy(xunits, gridptr->xunits);
}


void gridInqXstdname(int gridID, char *xstdname)
{
844
  grid_t *gridptr = gridID2Ptr(gridID);
845
846
847
848
  if ( gridptr->xstdname )
    strcpy(xstdname, gridptr->xstdname);
  else
    xstdname[0] = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
849
850
851
852
853
854
}

/*
@Function  gridInqYname
@Title     Get the name of a Y-axis

Deike Kleberg's avatar
Deike Kleberg committed
855
@Prototype void gridInqYname(int gridID, char *name)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
856
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
857
858
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
    @Item  name     Name of the Y-axis. The caller must allocate space for the
Deike Kleberg's avatar
Deike Kleberg committed
859
                    returned string. The maximum possible length, in characters, of
860
                    the string is given by the predefined constant @func{CDI_MAX_NAME}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
861
862

@Description
863
The function @func{gridInqYname} returns the name of a Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
864
865

@Result
866
@func{gridInqYname} returns the name of the Y-axis to the parameter name.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
867
868
869
870
871

@EndFunction
*/
void gridInqYname(int gridID, char *yname)
{
872
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
873
874
875
876
877
878
879
880

  strcpy(yname, gridptr->yname);
}

/*
@Function  gridInqYlongname
@Title     Get the longname of a Y-axis

Deike Kleberg's avatar
Deike Kleberg committed
881
@Prototype void gridInqXlongname(int gridID, char *longname)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
882
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
883
884
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
    @Item  longname Longname of the Y-axis. The caller must allocate space for the
Deike Kleberg's avatar
Deike Kleberg committed
885
                    returned string. The maximum possible length, in characters, of
886
                    the string is given by the predefined constant @func{CDI_MAX_NAME}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
887
888

@Description
889
The function @func{gridInqYlongname} returns the longname of a Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
890
891

@Result
892
@func{gridInqYlongname} returns the longname of the Y-axis to the parameter longname.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
893
894
895
896
897

@EndFunction
*/
void gridInqYlongname(int gridID, char *ylongname)
{
898
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
899
900
901
902
903
904
905
906

  strcpy(ylongname, gridptr->ylongname);
}

/*
@Function  gridInqYunits
@Title     Get the units of a Y-axis

Deike Kleberg's avatar
Deike Kleberg committed
907
@Prototype void gridInqYunits(int gridID, char *units)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
908
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
909
910
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
    @Item  units    Units of the Y-axis. The caller must allocate space for the
Deike Kleberg's avatar
Deike Kleberg committed
911
                    returned string. The maximum possible length, in characters, of
912
                    the string is given by the predefined constant @func{CDI_MAX_NAME}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
913
914

@Description
915
The function @func{gridInqYunits} returns the units of a Y-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
916
917

@Result
918
@func{gridInqYunits} returns the units of the Y-axis to the parameter units.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
919
920
921
922
923

@EndFunction
*/
void gridInqYunits(int gridID, char *yunits)
{
924
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
925
926
927
928
929
930

  strcpy(yunits, gridptr->yunits);
}

void gridInqYstdname(int gridID, char *ystdname)
{
931
  grid_t *gridptr = gridID2Ptr(gridID);
932
933
934
935
  if ( gridptr->ystdname )
    strcpy(ystdname, gridptr->ystdname);
  else
    ystdname[0] = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
936
937
938
939
940
941
942
943
}

/*
@Function  gridInqType
@Title     Get the type of a Grid

@Prototype int gridInqType(int gridID)
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
944
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
945
946

@Description
947
The function @func{gridInqType} returns the type of a Grid.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
948
949

@Result
950
@func{gridInqType} returns the type of the grid,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
951
one of the set of predefined CDI grid types.
952
The valid CDI grid types are @func{GRID_GENERIC}, @func{GRID_GAUSSIAN},
Uwe Schulzweida's avatar
Uwe Schulzweida committed
953
@func{GRID_LONLAT}, @func{GRID_LCC}, @func{GRID_SPECTRAL}, @func{GRID_GME},
954
@func{GRID_CURVILINEAR} and @func{GRID_UNSTRUCTURED}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
955
956
957
958
959

@EndFunction
*/
int gridInqType(int gridID)
{
960
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
961
962
963
964
965
966
967
968
969
970
971

  return (gridptr->type);
}


/*
@Function  gridInqSize
@Title     Get the size of a Grid

@Prototype int gridInqSize(int gridID)
@Parameter
Uwe Schulzweida's avatar
Uwe Schulzweida committed
972
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate} or @fref{vlistInqVarGrid}.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
973
974

@Description
975
The function @func{gridInqSize} returns the size of a Grid.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
976
977

@Result
978
@func{gridInqSize} returns the number of grid points of a Grid.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
979
980
981
982
983

@EndFunction
*/
int gridInqSize(int gridID)
{
984
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
985

Uwe Schulzweida's avatar
Uwe Schulzweida committed
986
  int size = gridptr->size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
987
988
989
990
991
992
993
994
995

  if ( ! size )
    {
      int xsize, ysize;

      xsize = gridptr->xsize;
      ysize = gridptr->ysize;

      if ( ysize )
996
        size = xsize * ysize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
997
      else
Thomas Jahns's avatar
Thomas Jahns committed
998
        size = xsize;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
999

Thomas Jahns's avatar
Thomas Jahns committed
1000
      gridptr->size = size;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1001
1002
1003
1004
1005
    }

  return (size);
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1006
1007
static
int nsp2trunc(int nsp)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1008
1009
1010
1011
1012
1013
1014
{
  /*  nsp = (trunc+1)*(trunc+1)              */
  /*      => trunc^2 + 3*trunc - (x-2) = 0   */
  /*                                         */
  /*  with:  y^2 + p*y + q = 0               */
  /*         y = -p/2 +- sqrt((p/2)^2 - q)   */
  /*         p = 3 and q = - (x-2)           */
1015
  int trunc = (int) (sqrt(nsp*4 + 1.) - 3) / 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1016
1017
1018
1019
1020
1021
  return (trunc);
}


int gridInqTrunc(int gridID)
{
1022
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1023
1024
1025
1026

  if ( gridptr->trunc == 0 )
    {
      if ( gridptr->type == GRID_SPECTRAL )
Thomas Jahns's avatar
Thomas Jahns committed
1027
        gridptr->trunc = nsp2trunc(gridptr->size);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1028
1029
      /*
      else if      ( gridptr->type == GRID_GAUSSIAN )
Thomas Jahns's avatar
Thomas Jahns committed
1030
        gridptr->trunc = nlat2trunc(gridptr->ysize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1031
1032
1033
1034
1035
1036
      */
    }

  return (gridptr->trunc);
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1037

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1038
1039
void gridDefTrunc(int gridID, int trunc)
{
1040
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1041

1042
1043
  if (gridptr->trunc != trunc)
    {
1044
      gridMark4Update(gridID);
1045
1046
      gridptr->trunc = trunc;
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1047
1048
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1049
1050
1051
1052
1053
1054
/*
@Function  gridDefXsize
@Title     Define the number of values of a X-axis

@Prototype void gridDefXsize(int gridID, int xsize)
@Parameter
Deike Kleberg's avatar
Deike Kleberg committed
1055
1056
    @Item  gridID   Grid ID, from a previous call to @fref{gridCreate}.
    @Item  xsize    Number of values of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1057
1058

@Description
1059
The function @func{gridDefXsize} defines the number of values of a X-axis.
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1060
1061
1062

@EndFunction
*/
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1063
1064
void gridDefXsize(int gridID, int xsize)
{
1065
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1066

1067
1068
1069
  int gridSize = gridInqSize(gridID);
  if ( xsize > gridSize )
    Error("xsize %d is greater then gridsize %d", xsize, gridSize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1070

1071
1072
  int gridType = gridInqType(gridID);
  if ( gridType == GRID_UNSTRUCTURED && xsize != gridSize )
1073
    Error("xsize %d must be equal to gridsize %d for gridtype: UNSTRUCTURED", xsize, gridSize);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1074

1075
1076
  if (gridptr->xsize != xsize)
    {
1077
      gridMark4Update(gridID);
1078
1079
      gridptr->xsize = xsize;
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1080

1081
  if ( gridType != GRID_UNSTRUCTURED )
1082
    {
1083
1084
      long axisproduct = gridptr->xsize*gridptr->ysize;
      if ( axisproduct > 0 && axisproduct != gridSize )
1085
        Error("Inconsistent grid declaration! (xsize=%d ysize=%d gridsize=%d)",
1086
              gridptr->xsize, gridptr->ysize, gridSize);
1087
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1088
1089
1090
}

/*
Thomas Jahns's avatar
Thomas Jahns committed
1091
1092
@Function
@Title
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1093

Thomas Jahns's avatar
Thomas Jahns committed
1094
@Prototype
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1095
1096
1097
1098
1099
1100
1101
@Parameter
    @Item  Grid identifier

@EndFunction
*/
void gridDefPrec(int gridID, int prec)
{
1102
  grid_t *gridptr = gridID2Ptr(gridID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1103

1104
1105
  if (gridptr->prec != prec)
    {
1106
      gridMark4Update(gridID);