Echam5ini.c 39 KB
Newer Older
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1
2
3
4
/*
  This file is part of CDO. CDO is a collection of Operators to
  manipulate and analyse Climate model Data.

Uwe Schulzweida's avatar
Uwe Schulzweida committed
5
  Copyright (C) 2003-2007 Uwe Schulzweida, Uwe.Schulzweida@zmaw.de
Uwe Schulzweida's avatar
Uwe Schulzweida committed
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
  See COPYING file for copying and redistribution conditions.

  This program is free software; you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation; version 2 of the License.

  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.
*/

#if  defined  (HAVE_CONFIG_H)
#  include "config.h"
#endif

#include <string.h>

#include "cdi.h"
#include "cdo.h"
#include "cdo_int.h"
#include "pstream.h"


#if  defined  (HAVE_LIBNETCDF)
#  include "netcdf.h"
#endif

static int nvars_ml = 4;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
35
36
static const char strfiletype_ml[]  = "Initial file spectral";
static const char strfiletype_res[] = "Restart history file";
Uwe Schulzweida's avatar
Uwe Schulzweida committed
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53


typedef struct {
  int  gridtype;
  int  zaxistype;
  int  code;
  char *name;
  char *longname;
  char *units;
  int gridID;
  int zaxisID;
  int gridsize;
  int nlev;
  double *ptr;
} VAR;


Uwe Schulzweida's avatar
Uwe Schulzweida committed
54
55
56
57
58
59
60
61
62
63
64
65
66
typedef struct {
  int    naint;
  int    naflt;
  int    natxt;
  char   *aintname[1024];
  int    *aintentry[1024];
  char   *afltname[1024];
  double *afltentry[1024];
  char   *atxtname[1024];
  char   *atxtentry[1024];
} ATTS;


Uwe Schulzweida's avatar
Uwe Schulzweida committed
67
static void iniatts(ATTS *atts)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
68
69
70
71
72
73
74
{
  atts->naint = 0;
  atts->naflt = 0;
  atts->natxt = 0;
}


Uwe Schulzweida's avatar
Uwe Schulzweida committed
75
static void inivar(VAR *var, int gridtype, int zaxistype, int code, const char *name,
Uwe Schulzweida's avatar
Uwe Schulzweida committed
76
77
78
79
80
81
82
	       const char *longname, const char *units)
{
  static char func[] = "inivar_ml";
  
  var->gridtype  = gridtype;
  var->zaxistype = zaxistype;
  var->code      = code;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
83
84
85
86
87
88
  var->name      = NULL;
  var->longname  = NULL;
  var->units     = NULL;
  if ( name )     var->name      = strdup(name);
  if ( longname ) var->longname  = strdup(longname);
  if ( units )    var->units     = strdup(units);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
89
90
}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
91
92

static void inivars_ml(VAR **vars)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
93
94
95
96
97
{
  static char func[] = "inivars_ml";

  *vars = (VAR *) malloc((nvars_ml+1)*sizeof(VAR));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
98
99
100
101
  inivar(&(*vars)[0], GRID_GAUSSIAN, ZAXIS_HYBRID,  133, "Q",   "specific humidity", "kg/kg");
  inivar(&(*vars)[1], GRID_SPECTRAL, ZAXIS_HYBRID,  138, "SVO", "vorticity", "1/s");
  inivar(&(*vars)[2], GRID_SPECTRAL, ZAXIS_HYBRID,  155, "SD",  "divergence", "1/s");
  inivar(&(*vars)[3], GRID_SPECTRAL, ZAXIS_HYBRID,  130, "STP", "temperature", "K");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
102
  /* Don't change the order (lsp must be the last one)! */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
103
  inivar(&(*vars)[4], GRID_SPECTRAL, ZAXIS_SURFACE, 152, "LSP", "log surface pressure", "K");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
}


#if  defined  (HAVE_LIBNETCDF)
static void nce(int istat)
{
  /*
    This routine provides a simple interface to netCDF error message routine.
  */

  if ( istat != NC_NOERR ) cdoAbort(nc_strerror(istat));
}
#endif


Uwe Schulzweida's avatar
Uwe Schulzweida committed
119
static int read_e5ml(const char *filename, VAR **vars)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
120
{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
121
122
123
  static char func[] = "read_e5ml";
  int nvars = 0;
#if  defined  (HAVE_LIBNETCDF)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
124
125
126
127
128
129
130
  int nc_dim_id, nc_var_id;
  size_t dimlen, nvals;
  size_t start[3];
  size_t count[3];
  int nlon, nlat, nlev, nlevp1, nvct, nsp, i, iv;
  int gridIDgp, gridIDsp, zaxisIDml, zaxisIDsfc;
  int gridtype, zaxistype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
131
132
133
  int nc_file_id;
  char filetype[256];
  size_t attlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
134
135
  double *xvals, *yvals, *vct, *levs;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
136
137
138
139
140
141
142
143
144
145
146
  /* open file and check file type */
  nce(nc_open(filename, NC_NOWRITE, &nc_file_id));

  nce(nc_get_att_text(nc_file_id, NC_GLOBAL, "file_type", filetype));
  nce(nc_inq_attlen(nc_file_id, NC_GLOBAL, "file_type", &attlen));
  filetype[attlen] = 0;

  if ( strcmp(filetype, strfiletype_ml) != 0 ) return (0);

  inivars_ml(vars);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
  /* read dimensions */

  nce(nc_inq_dimid(nc_file_id, "lon", &nc_dim_id));
  nce(nc_inq_dimlen(nc_file_id, nc_dim_id, &dimlen));
  nlon = (int) dimlen;

  nce(nc_inq_dimid(nc_file_id, "lat", &nc_dim_id));
  nce(nc_inq_dimlen(nc_file_id, nc_dim_id, &dimlen));
  nlat = (int) dimlen;

  gridIDgp = gridCreate(GRID_GAUSSIAN, nlon*nlat);
  gridDefXsize(gridIDgp, nlon);
  gridDefYsize(gridIDgp, nlat);

  nce(nc_inq_dimid(nc_file_id, "nsp", &nc_dim_id));
  nce(nc_inq_dimlen(nc_file_id, nc_dim_id, &dimlen));
  nsp = (int) dimlen;

  gridIDsp = gridCreate(GRID_SPECTRAL, nsp*2);

  nce(nc_inq_dimid(nc_file_id, "nlev", &nc_dim_id));
  nce(nc_inq_dimlen(nc_file_id, nc_dim_id, &dimlen));
  nlev = (int) dimlen;
  nlevp1 = nlev + 1;
  nvct = nlevp1*2;

  zaxisIDsfc = zaxisCreate(ZAXIS_SURFACE, 1);
  zaxisIDml  = zaxisCreate(ZAXIS_HYBRID, nlev);

  levs = (double *) malloc(nlev*sizeof(double));
  for ( i = 0; i < nlev; i++ ) levs[i] = i+1;
  zaxisDefLevels(zaxisIDml, levs);
  free(levs);

  /* read variables */

  xvals = (double *) malloc(nlon*sizeof(double));
  yvals = (double *) malloc(nlat*sizeof(double));

  nce(nc_inq_varid(nc_file_id, "lon", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, xvals));

  nce(nc_inq_varid(nc_file_id, "lat", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, yvals));

  gridDefXvals(gridIDgp, xvals);
  gridDefYvals(gridIDgp, yvals);

  free(xvals);
  free(yvals);

  vct   = (double *) malloc(nvct*sizeof(double));

  nce(nc_inq_varid(nc_file_id, "vct_a", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, vct));

  nce(nc_inq_varid(nc_file_id, "vct_b", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, vct+nlevp1));

  zaxisDefVct(zaxisIDml, 2*nlevp1, vct);
  free(vct);

  for ( iv = 0; iv < nvars_ml; iv++ )
    {
      nvals = 0;

      gridtype  = (*vars)[iv].gridtype;
      zaxistype = (*vars)[iv].zaxistype;

      if ( gridtype == GRID_GAUSSIAN )
	{
	  (*vars)[iv].gridID = gridIDgp;
	  nvals += nlon*nlat;
	}
      else
	{
	  (*vars)[iv].gridID = gridIDsp;
	  nvals += nsp*2;
	}

      (*vars)[iv].zaxisID   = zaxisIDml;
      (*vars)[iv].gridsize  = nvals;
      (*vars)[iv].nlev      = nlev;

      (*vars)[iv].ptr = (double *) malloc(nlev*nvals*sizeof(double));
      
      for ( i = 0; i < nlev; i++ )
	{
	  if ( gridtype == GRID_GAUSSIAN )
	    {
	      start[0] = 0;     start[1] = i;  start[2] = 0;
	      count[0] = nlat;  count[1] = 1;  count[2] = nlon;     
	    }
	  else
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
242
243
	      start[0] = 0;     start[1] = 0;  start[2] = i;
	      count[0] = nsp;   count[1] = 2;  count[2] = 1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
	    }

	  nce(nc_inq_varid(nc_file_id, (*vars)[iv].name, &nc_var_id));
	  nce(nc_get_vara_double(nc_file_id, nc_var_id, start, count, (*vars)[iv].ptr+i*nvals));
	}
    }

  /* read lsp */

  (*vars)[nvars_ml].gridID    = gridIDsp;
  (*vars)[nvars_ml].zaxisID   = zaxisIDsfc;
  (*vars)[nvars_ml].gridsize  = nsp*2;
  (*vars)[nvars_ml].nlev      = 1;

  start[0] = 0;    start[1] = 0;  start[2] = nlev;
  count[0] = nsp;  count[1] = 2;  count[2] = 1;

  (*vars)[nvars_ml].ptr = (double *) malloc(nsp*2*sizeof(double));

  nce(nc_inq_varid(nc_file_id, "STP", &nc_var_id));
  nce(nc_get_vara_double(nc_file_id, nc_var_id, start, count, (*vars)[nvars_ml].ptr));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
266
267
268
269
270
271
  /*close input file */
  nce(nc_close(nc_file_id));

  nvars = nvars_ml + 1;

#else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
272
  cdoAbort("netCDF support not compiled in!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
273
274
275
#endif

  return (nvars);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
276
}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
277
278


Uwe Schulzweida's avatar
Uwe Schulzweida committed
279
280
281
282
static void write_e5ml(const char *filename, VAR *vars, int nvars)
{
  static char func[] = "write_e5ml";
#if  defined  (HAVE_LIBNETCDF)
283
284
  int nc_var_id;
  size_t nvals;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
285
286
287
288
289
  size_t start[3], count[3];
  int dimidsp[9];
  int varid;
  int ilev;
  int lon, lat;
290
  int nlon, nlat, nlev, nlevp1, nvct, nsp, n2, i, nvclev;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
291
  int lat_dimid, lon_dimid, nlev_dimid, nlevp1_dimid, nsp_dimid, nvclev_dimid, n2_dimid;
292
  int gridIDgp = -1, gridIDsp, zaxisIDml = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
293
294
295
  int gridtype, zaxistype;
  int nc_file_id;
  int nc_stpid, lspid;
296
  double *xvals, *yvals;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
297
298
299
300
301
302
303
304
305
  const double *vct;

  /* create file */
  nce(nc_create(filename, NC_CLOBBER, &nc_file_id));

  nce(nc_put_att_text(nc_file_id, NC_GLOBAL, "file_type", strlen(strfiletype_ml), strfiletype_ml));

  n2 = 2;

306
  lon = 0; lat = 0; nsp = 0; nlev = 0; nlevp1 = 0; nvclev = 0;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
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
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
  for ( varid = 0; varid < nvars; ++varid )
    {
      gridtype  = vars[varid].gridtype;
      zaxistype = vars[varid].zaxistype;

      if ( gridtype == GRID_GAUSSIAN && lat == 0 )
	{
	  gridIDgp = vars[varid].gridID;
	  lon = gridInqXsize(vars[varid].gridID);
	  lat = gridInqYsize(vars[varid].gridID);
	}
      else if ( gridtype == GRID_SPECTRAL && nsp == 0 )
	{
	  gridIDsp = vars[varid].gridID;
	  nsp = gridInqSize(vars[varid].gridID);
	  nsp = nsp/2;
	}

      if ( zaxistype == ZAXIS_HYBRID && nlev == 0 )
	{
	  zaxisIDml = vars[varid].zaxisID;
	  nlev = zaxisInqSize(vars[varid].zaxisID);
	  nlevp1 = nlev + 1;
	  nvclev = nlev + 1;
	}
    }

  if ( lat == 0 ) cdoAbort("Gaussian grid not found!");
  if ( nsp == 0 ) cdoAbort("Spectral data not found!");
  if ( nlev == 0 ) cdoAbort("Hybrid level not found!");

  nlon = lon;
  nlat = lat;
  
  nce(nc_def_dim(nc_file_id, "lat", lat, &lat_dimid));

  nce(nc_def_dim(nc_file_id, "lon", lon, &lon_dimid));

  nce(nc_def_dim(nc_file_id, "nlev", nlev, &nlev_dimid));
  nce(nc_def_dim(nc_file_id, "nlevp1", nlevp1, &nlevp1_dimid));

  nce(nc_def_dim(nc_file_id, "nsp", nsp, &nsp_dimid));

  nce(nc_def_dim(nc_file_id, "nvclev", nvclev, &nvclev_dimid));

  nce(nc_def_dim(nc_file_id, "n2", n2, &n2_dimid));

  nce(nc_enddef(nc_file_id));

  /* define gaussian grid */

  xvals = (double *) malloc(nlon*sizeof(double));
  yvals = (double *) malloc(nlat*sizeof(double));

  gridInqXvals(gridIDgp, xvals);
  gridInqYvals(gridIDgp, yvals);

  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "lat", NC_DOUBLE, 1, &lat_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, yvals));
   
  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "lon", NC_DOUBLE, 1, &lon_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, xvals));

  free(xvals);
  free(yvals);

  /* define model level */

  nvct = nvclev*2;

  /* vct   = (double *) malloc(nvct*sizeof(double)); */

  vct = zaxisInqVctPtr(zaxisIDml);

  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "vct_a", NC_DOUBLE, 1, &nvclev_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, vct));

  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "vct_b", NC_DOUBLE, 1, &nvclev_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, vct+nlevp1));

  /* free(vct); */

  lspid = -1;
  nc_stpid = -1;

  for ( varid = 0; varid < nvars; varid++ )
    {
      nvals = 0;

      gridtype  = vars[varid].gridtype;
      zaxistype = vars[varid].zaxistype;

      ilev = zaxisInqSize(vars[varid].zaxisID);

      if ( ilev == 1 )
	{
	  lspid = varid;
412
413
	  if ( gridtype != GRID_SPECTRAL )
	    cdoAbort("%s has wrong gridtype!", vars[varid].name);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
	  continue;
	}

      if ( nlev != ilev )
	cdoAbort("Unexpected number of level %d!", ilev);

      if ( gridtype == GRID_GAUSSIAN )
	{
	  nvals = nlon*nlat;

	  dimidsp[0] = lat_dimid;
	  dimidsp[1] = nlev_dimid;
	  dimidsp[2] = lon_dimid;
	}
      else if ( gridtype == GRID_SPECTRAL )
	{
	  nvals = nsp*2;

	  dimidsp[0] = nsp_dimid;
	  dimidsp[1] = n2_dimid;

	  if ( strcmp(vars[varid].name, "STP") == 0 || strcmp(vars[varid].name, "T") == 0 )
	    dimidsp[2] = nlevp1_dimid;
	  else
	    dimidsp[2] = nlev_dimid;
	}
      else
	cdoAbort("Unsupported grid!");

      nce(nc_redef(nc_file_id));
      nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 3, dimidsp, &nc_var_id));
      /*
	nce(nc_put_att_text(nc_file_id, nc_var_id, "long_name", strlen(vars[varid].longname), vars[varid].longname));*/
      nce(nc_enddef(nc_file_id));

      if (  dimidsp[2] == nlevp1_dimid ) nc_stpid = nc_var_id;

      for ( i = 0; i < nlev; i++ )
	{
	  if ( gridtype == GRID_GAUSSIAN )
	    {
	      start[0] = 0;     start[1] = i;  start[2] = 0;
	      count[0] = nlat;  count[1] = 1;  count[2] = nlon;     
	    }
	  else
	    {
	      start[0] = 0;     start[1] = 0;  start[2] = i;
	      count[0] = nsp;   count[1] = 2;  count[2] = 1;
	    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
463
	  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
464
465
466
467
468
469
470
471
472
473
474
	  nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	}
    }

  if ( lspid == -1 ) cdoAbort("LSP not found!");
  if ( nc_stpid == -1 ) cdoAbort("STP not found!");

  /* write lsp */
  start[0] = 0;    start[1] = 0;  start[2] = nlev;
  count[0] = nsp;  count[1] = 2;  count[2] = 1;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
475
  nce(nc_put_vara_double(nc_file_id, nc_stpid, start, count, vars[lspid].ptr));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
476
477
478
479
480
481
482
483
484
485

  /*close input file */
  nce(nc_close(nc_file_id));

#else
  cdoAbort("netCDF support not compiled in!");
#endif
}


Uwe Schulzweida's avatar
Uwe Schulzweida committed
486
#if  defined  (HAVE_LIBNETCDF)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
487
static void read_gg3d(int nc_file_id, const char *name, VAR *var, int gridID, int zaxisID)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
{
  static char func[] = "read_gg3d";
  int nlev, nlat, nlon, gridsize, i;
  int gridtype, zaxistype;
  int nc_var_id;
  size_t start[3], count[3];

  gridtype = gridInqType(gridID);
  zaxistype = zaxisInqType(zaxisID);

  inivar(var, gridtype, zaxistype,  0, name, NULL, NULL);

  nlon = gridInqXsize(gridID);
  nlat = gridInqYsize(gridID);
  nlev = zaxisInqSize(zaxisID);

  gridsize = nlon*nlat;

  var->gridID    = gridID;
  var->zaxisID   = zaxisID;
  var->gridsize  = gridsize;
  var->nlev      = nlev;

  var->ptr = (double *) malloc(nlev*gridsize*sizeof(double));
  
  for ( i = 0; i < nlev; i++ )
    {
      start[0] = 0;     start[1] = i;  start[2] = 0;
      count[0] = nlat;  count[1] = 1;  count[2] = nlon;     

      nce(nc_inq_varid(nc_file_id, name, &nc_var_id));
      nce(nc_get_vara_double(nc_file_id, nc_var_id, start, count, var->ptr+i*gridsize));
    }
}
#endif

#if  defined  (HAVE_LIBNETCDF)
static void read_fc4d(int nc_file_id, const char *name, VAR *var, int gridID, int zaxisID, int nhgl, int nmp1)
{
  static char func[] = "read_fc4d";
  int nlev, nfc, i;
  int gridtype, zaxistype;
  int nc_var_id;
  size_t start[4], count[4];

  gridtype = gridInqType(gridID);
  zaxistype = zaxisInqType(zaxisID);

  inivar(var, gridtype, zaxistype,  0, name, NULL, NULL);

  nfc  = gridInqSize(gridID);
  nlev = zaxisInqSize(zaxisID);

  if ( nfc != nhgl*nmp1*2 ) cdoAbort("%s: inconsistent dimension length!", func);

  var->gridID    = gridID;
  var->zaxisID   = zaxisID;
  var->gridsize  = nfc;
  var->nlev      = nlev;

  var->ptr = (double *) malloc(nlev*nfc*sizeof(double));
  
  for ( i = 0; i < nlev; i++ )
    {
      start[0] = 0;     start[1] = 0;    start[2] = 0;  start[3] = i;
      count[0] = nhgl;  count[1] = nmp1; count[2] = 2;  count[3] = 1;     

      nce(nc_inq_varid(nc_file_id, name, &nc_var_id));
      nce(nc_get_vara_double(nc_file_id, nc_var_id, start, count, var->ptr+i*nfc));
    }
}
#endif

Uwe Schulzweida's avatar
Uwe Schulzweida committed
561

Uwe Schulzweida's avatar
Uwe Schulzweida committed
562
static int read_e5res(const char *filename, VAR **vars, ATTS *atts)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
563
564
565
566
{
  static char func[] = "read_e5res";
  int nvars = 0;
#if  defined  (HAVE_LIBNETCDF)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
567
  int nc_var_id;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
568
  int varid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
569
570
571
572
  size_t nvals;
  size_t start[3], count[3];
  int nlon, nlat, nlev, nvct, nfc, i;
  int gridIDgp, gridIDfc, gridIDhgl, zaxisIDml, zaxisIDmlh, zaxisIDsfc, zaxisIDbsfc, zaxisIDn2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
573
574
575
576
577
578
579
580
581
582
583
  int nc_file_id;
  char filetype[256];
  double *xvals, *yvals, *vct, *levs;
  int ncvarid;
  int ndims, ngatts, unlimdimid;
  int nvdims, nvatts;
  int dimidsp[9];
  int max_vars = 4096;
  char name[256];
  int lon_dimid, lat_dimid, nhgl_dimid, nlevp1_dimid, spc_dimid, nvclev_dimid;
  int complex_dimid, nmp1_dimid, belowsurface_dimid, lev_dimid, ilev_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
584
  int /* surface_dimid, height2m_dimid, height10m_dimid,*/ n2_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
585
586
  int lon, lat, nhgl, nlevp1, spc, nvclev;
  int complex, nmp1, belowsurface, lev, ilev;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
587
588
589
590
  int /* surface, height2m, height10m,*/ n2;
  int iatt;
  int attint;
  double attflt;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
591
  size_t attlen, dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
592
593
594
  nc_type xtype;
  char attname[256];
  const int attstringlen = 8192; char attstring[8192];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
595
596
597
598
599
600
601
602
603
604

  /* open file and check file type */
  nce(nc_open(filename, NC_NOWRITE, &nc_file_id));

  nce(nc_get_att_text(nc_file_id, NC_GLOBAL, "file_type", filetype));
  nce(nc_inq_attlen(nc_file_id, NC_GLOBAL, "file_type", &attlen));
  filetype[attlen] = 0;

  if ( strncmp(filetype, strfiletype_res, strlen(strfiletype_res)) != 0 ) return (0);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
605
  /* printf("%s\n", filetype); */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
606
607

  nce(nc_inq_dimid(nc_file_id, "lon", &lon_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
608
609
  nce(nc_inq_dimlen(nc_file_id, lon_dimid, &dimlen));
  lon = dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
610
611

  nce(nc_inq_dimid(nc_file_id, "lat", &lat_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
612
613
  nce(nc_inq_dimlen(nc_file_id, lat_dimid, &dimlen));
  lat = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
614
615

  nce(nc_inq_dimid(nc_file_id, "nhgl", &nhgl_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
616
617
  nce(nc_inq_dimlen(nc_file_id, nhgl_dimid, &dimlen));
  nhgl = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
618
619

  nce(nc_inq_dimid(nc_file_id, "nlevp1", &nlevp1_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
620
621
  nce(nc_inq_dimlen(nc_file_id, nlevp1_dimid, &dimlen));
  nlevp1 = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
622
623

  nce(nc_inq_dimid(nc_file_id, "spc", &spc_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
624
625
  nce(nc_inq_dimlen(nc_file_id, spc_dimid, &dimlen));
  spc = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
626
627

  nce(nc_inq_dimid(nc_file_id, "nvclev", &nvclev_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
628
629
  nce(nc_inq_dimlen(nc_file_id, nvclev_dimid, &dimlen));
  nvclev = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
630
631

  nce(nc_inq_dimid(nc_file_id, "complex", &complex_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
632
633
  nce(nc_inq_dimlen(nc_file_id, complex_dimid, &dimlen));
  complex = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
634
635

  nce(nc_inq_dimid(nc_file_id, "nmp1", &nmp1_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
636
637
  nce(nc_inq_dimlen(nc_file_id, nmp1_dimid, &dimlen));
  nmp1 = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
638
639

  nce(nc_inq_dimid(nc_file_id, "belowsurface", &belowsurface_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
640
641
  nce(nc_inq_dimlen(nc_file_id, belowsurface_dimid, &dimlen));
  belowsurface = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
642
643

  nce(nc_inq_dimid(nc_file_id, "lev", &lev_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
644
645
  nce(nc_inq_dimlen(nc_file_id, lev_dimid, &dimlen));
  lev = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
646
647

  nce(nc_inq_dimid(nc_file_id, "ilev", &ilev_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
648
649
  nce(nc_inq_dimlen(nc_file_id, ilev_dimid, &dimlen));
  ilev = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
650

Uwe Schulzweida's avatar
Uwe Schulzweida committed
651
  /*
Uwe Schulzweida's avatar
Uwe Schulzweida committed
652
653
  nce(nc_inq_dimid(nc_file_id, "surface", &surface_dimid));
  nce(nc_inq_dimlen(nc_file_id, surface_dimid, &surface));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
654

Uwe Schulzweida's avatar
Uwe Schulzweida committed
655
656
657
658
659
660
661
  nce(nc_inq_dimid(nc_file_id, "height2m", &height2m_dimid));
  nce(nc_inq_dimlen(nc_file_id, height2m_dimid, &height2m));

  nce(nc_inq_dimid(nc_file_id, "height10m", &height10m_dimid));
  nce(nc_inq_dimlen(nc_file_id, height10m_dimid, &height10m));
  */
  nce(nc_inq_dimid(nc_file_id, "n2", &n2_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
662
663
  nce(nc_inq_dimlen(nc_file_id, n2_dimid, &dimlen));
  n2 = (int) dimlen;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688

  /* define gaussian grid */

  nlon = lon;
  nlat = lat;

  gridIDgp = gridCreate(GRID_GAUSSIAN, nlon*nlat);
  gridDefXsize(gridIDgp, nlon);
  gridDefYsize(gridIDgp, nlat);

  xvals = (double *) malloc(nlon*sizeof(double));
  yvals = (double *) malloc(nlat*sizeof(double));

  nce(nc_inq_varid(nc_file_id, "lon", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, xvals));

  nce(nc_inq_varid(nc_file_id, "lat", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, yvals));

  gridDefXvals(gridIDgp, xvals);
  gridDefYvals(gridIDgp, yvals);

  free(xvals);
  free(yvals);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
689
690
691
692
693
694
695
696
697
698
699
  /* define fourier grid */

  nfc = nhgl*nmp1*complex;

  gridIDfc = gridCreate(GRID_FOURIER, nfc);

  /* define nhgl grid */

  gridIDhgl = gridCreate(GRID_GENERIC, nhgl);


Uwe Schulzweida's avatar
Uwe Schulzweida committed
700
701
702
703
  /* define surface level */

  zaxisIDsfc = zaxisCreate(ZAXIS_SURFACE, 1);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
704
  /* define below surface level */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
705

Uwe Schulzweida's avatar
Uwe Schulzweida committed
706
707
  nlev = belowsurface;
  zaxisIDbsfc = zaxisCreate(ZAXIS_DEPTH_BELOW_LAND, nlev);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
708

Uwe Schulzweida's avatar
Uwe Schulzweida committed
709
710
711
712
713
714
715
716
717
718
  levs = (double *) malloc(nlev*sizeof(double));
  for ( i = 0; i < nlev; i++ ) levs[i] = 0;
  zaxisDefLevels(zaxisIDbsfc, levs);
  free(levs);


  /* define n2 level */

  nlev = n2;
  zaxisIDn2 = zaxisCreate(ZAXIS_GENERIC, nlev);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
719
720

  levs = (double *) malloc(nlev*sizeof(double));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
721
722
  for ( i = 0; i < nlev; i++ ) levs[i] = 0;
  zaxisDefLevels(zaxisIDn2, levs);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
723
724
  free(levs);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
725
726
727
728
729
  /* define model level */

  nlev = lev;
  nvct = nvclev*2;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
730
731
732
733
734
735
736
737
  vct   = (double *) malloc(nvct*sizeof(double));

  nce(nc_inq_varid(nc_file_id, "vct_a", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, vct));

  nce(nc_inq_varid(nc_file_id, "vct_b", &nc_var_id));
  nce(nc_get_var_double(nc_file_id, nc_var_id, vct+nlevp1));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
738
739
740
741
742
743
744
745
746
  /* ZAXIS_HYBRID */ 

  zaxisIDml  = zaxisCreate(ZAXIS_HYBRID, nlev);

  levs = (double *) malloc(nlev*sizeof(double));
  for ( i = 0; i < nlev; i++ ) levs[i] = i+1;
  zaxisDefLevels(zaxisIDml, levs);
  free(levs);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
747
  zaxisDefVct(zaxisIDml, 2*nlevp1, vct);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
748
749
750
751
752
753
754
755
756
757
758
759

  /* ZAXIS_HYBRID_HALF */ 

  zaxisIDmlh  = zaxisCreate(ZAXIS_HYBRID_HALF, nlevp1);

  levs = (double *) malloc(nlevp1*sizeof(double));
  for ( i = 0; i < nlevp1; i++ ) levs[i] = i;
  zaxisDefLevels(zaxisIDmlh, levs);
  free(levs);

  zaxisDefVct(zaxisIDmlh, 2*nlevp1, vct);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
760
761
762
763
764
  free(vct);


  nce(nc_inq(nc_file_id, &ndims, &nvars, &ngatts, &unlimdimid));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
  /* read global attributtes*/
  for ( iatt = 0; iatt < ngatts; iatt++ )
    {
      nce(nc_inq_attname(nc_file_id, NC_GLOBAL, iatt, attname));
      nce(nc_inq_atttype(nc_file_id, NC_GLOBAL, attname, &xtype));
      nce(nc_inq_attlen(nc_file_id, NC_GLOBAL, attname, &attlen));

      if ( xtype == NC_CHAR )
	{
	  if ( (int)attlen > attstringlen )
	    {
	      fprintf(stderr, "Attribute %s too large, skipped!\n", attname);
	      continue;
	    }
	  nce(nc_get_att_text(nc_file_id, NC_GLOBAL, attname, attstring));
	  attstring[attlen] = 0;

	  atts->atxtname[atts->natxt]  = strdup(attname);
	  atts->atxtentry[atts->natxt] = strdup(attstring);
	  atts->natxt++;
	}
      else if ( xtype == NC_INT )
	{
	  if ( attlen > 1 )
	    {
	      fprintf(stderr, "Attribute %s too large, skipped!\n", attname);
	      continue;
	    }
	  nce(nc_get_att_int(nc_file_id, NC_GLOBAL, attname, &attint));
	}
      else if ( xtype == NC_DOUBLE )
	{
	  if ( attlen > 1 )
	    {
	      fprintf(stderr, "Attribute %s too large, skipped!\n", attname);
	      continue;
	    }
	  nce(nc_get_att_double(nc_file_id, NC_GLOBAL, attname, &attflt));
	}
      else
	{
	  fprintf(stderr, "Unsupported attribute type for %s\n", attname);
	}
    }

Uwe Schulzweida's avatar
Uwe Schulzweida committed
810
811
812
813
814
815
816
  *vars = (VAR *) malloc(max_vars*sizeof(VAR));

  varid = 0;
  for ( ncvarid = 0; ncvarid < nvars; ncvarid++ )
    {      
      nce(nc_inq_var(nc_file_id, ncvarid, name, &xtype, &nvdims, dimidsp, &nvatts));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
817
818
      if ( varid >= max_vars ) cdoAbort("Too many variables (max = %d)!", max_vars);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
819
820
821
822
823
      if ( nvdims == 4 )
	{
	  if ( dimidsp[0] == nhgl_dimid    && dimidsp[1] == nmp1_dimid &&
	       dimidsp[2] == complex_dimid && dimidsp[3] == lev_dimid )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
824
825
	      read_fc4d(nc_file_id, name, &(*vars)[varid], gridIDfc, zaxisIDml, nhgl, nmp1);
	      varid++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
826
827
828
829
	    }
	  else if ( dimidsp[0] == nhgl_dimid    && dimidsp[1] == nmp1_dimid &&
	            dimidsp[2] == complex_dimid && dimidsp[3] == ilev_dimid )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
830
831
	      read_fc4d(nc_file_id, name, &(*vars)[varid], gridIDfc, zaxisIDmlh, nhgl, nmp1);
	      varid++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
832
833
834
	    }
	  else
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
835
	      fprintf(stderr, "4D structure of %s unsupported!\n", name);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
836
837
838
839
840
841
842
	    }
	}
      else if ( nvdims == 3 )
	{
	  if ( dimidsp[0] == lat_dimid && dimidsp[1] == lev_dimid &&
	       dimidsp[2] == lon_dimid)
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
843
	      read_gg3d(nc_file_id, name, &(*vars)[varid], gridIDgp, zaxisIDml);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
844
845
846
847
	      varid++;
	    }
	  else if ( dimidsp[0] == lat_dimid && dimidsp[1] == ilev_dimid &&
		    dimidsp[2] == lon_dimid)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
848
849
850
	    {  
	      read_gg3d(nc_file_id, name, &(*vars)[varid], gridIDgp, zaxisIDmlh);
	      varid++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
851
852
853
854
	    }
	  else if ( dimidsp[0] == lat_dimid && dimidsp[1] == belowsurface_dimid &&
		    dimidsp[2] == lon_dimid)
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
855
856
	      read_gg3d(nc_file_id, name, &(*vars)[varid], gridIDgp, zaxisIDbsfc);
	      varid++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
857
858
859
860
	    }
	  else if ( dimidsp[0] == lat_dimid && dimidsp[1] == n2_dimid &&
		    dimidsp[2] == lon_dimid)
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
861
862
	      read_gg3d(nc_file_id, name, &(*vars)[varid], gridIDgp, zaxisIDn2);
	      varid++;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
863
864
865
	    }
	  else
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
866
	      fprintf(stderr, "3D structure of %s unsupported!\n", name);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
	    }
	}
      else if ( nvdims == 2 )
	{
	  if ( dimidsp[0] == lat_dimid && dimidsp[1] == lon_dimid)
	    {
	      nvals = nlon*nlat;
  
	      inivar(&(*vars)[varid], GRID_GAUSSIAN, ZAXIS_SURFACE,  0, name, NULL, NULL);

	      (*vars)[varid].gridID    = gridIDgp;
	      (*vars)[varid].zaxisID   = zaxisIDsfc;
	      (*vars)[varid].gridsize  = nvals;
	      (*vars)[varid].nlev      = 1;

	      (*vars)[varid].ptr = (double *) malloc(nvals*sizeof(double));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
	      nce(nc_inq_varid(nc_file_id, name, &nc_var_id));
	      nce(nc_get_var_double(nc_file_id, nc_var_id, (*vars)[varid].ptr));

	      varid++;
	    }
	  else if ( dimidsp[0] == nhgl_dimid && dimidsp[1] == lev_dimid)
	    {
	      nvals = nhgl;
  
	      inivar(&(*vars)[varid], GRID_GENERIC, ZAXIS_HYBRID,  0, name, NULL, NULL);

	      (*vars)[varid].gridID    = gridIDhgl;
	      (*vars)[varid].zaxisID   = zaxisIDml;
	      (*vars)[varid].gridsize  = nvals;
	      (*vars)[varid].nlev      = nlev;

	      (*vars)[varid].ptr = (double *) malloc(nvals*nlev*sizeof(double));

Uwe Schulzweida's avatar
Uwe Schulzweida committed
902
903
	      for ( i = 0; i < nlev; i++ )
		{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
904
905
906
		  start[0] = 0;      start[1] = i;;
		  count[0] = nvals;  count[1] = 1;;     

Uwe Schulzweida's avatar
Uwe Schulzweida committed
907
		  nce(nc_inq_varid(nc_file_id, name, &nc_var_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
908
		  nce(nc_get_vara_double(nc_file_id, nc_var_id, start, count, (*vars)[varid].ptr+i*nvals));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
909
910
911
912
913
914
		}

	      varid++;
	    }
	  else
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
915
	      fprintf(stderr, "2D structure of %s unsupported!\n", name);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
	    }
	}
      else if ( nvdims == 1 )
	{
	  if ( dimidsp[0] == lat_dimid && strcmp(name, "lat") == 0 )
	    {
	    }
	  else if ( dimidsp[0] == lon_dimid && strcmp(name, "lon") == 0 )
	    {
	    }
	  else if ( dimidsp[0] == nvclev_dimid && strcmp(name, "vct_a") == 0 )
	    {
	    }
	  else if ( dimidsp[0] == nvclev_dimid && strcmp(name, "vct_b") == 0 )
	    {
	    }
	  else
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
934
	      fprintf(stderr, "1D structure of %s unsupported!\n", name);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
935
936
937
938
	    }
	}
      else
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
939
	  fprintf(stderr, "structure of %s unsupported!\n", name);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
940
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
941
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
942
943
944
945
946
947
948

  /*close input file */
  nce(nc_close(nc_file_id));

  nvars = varid;

#else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
949
  cdoAbort("netCDF support not compiled in!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
950
951
952
953
954
955
#endif

  return (nvars);
}


Uwe Schulzweida's avatar
Uwe Schulzweida committed
956
static void write_e5res(const char *filename, VAR *vars, int nvars)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
957
958
959
{
  static char func[] = "write_e5res";
#if  defined  (HAVE_LIBNETCDF)
Uwe Schulzweida's avatar
Uwe Schulzweida committed
960
  int nc_var_id;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
961
  int varid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
962
  size_t nvals;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
963
  size_t start[4], count[4];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
964
965
  int nlon, nlat, nlev, nvct, nfc, i;
  int gridIDgp = -1, zaxisIDml = -1;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
966
  int nc_file_id;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
967
968
  double *xvals, *yvals;
  const double *vct;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
969
970
971
  int dimidsp[9];
  int lon_dimid, lat_dimid, nhgl_dimid, nlevp1_dimid, spc_dimid, nvclev_dimid;
  int complex_dimid, nmp1_dimid, belowsurface_dimid, lev_dimid, ilev_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
972
973
974
975
976
  int /* surface_dimid, height2m_dimid, height10m_dimid,*/ n2_dimid;
  int lon = 0, lat = 0, nhgl = 0, nlevp1 = 0, spc = 0, nvclev = 0;
  int complex, nmp1, belowsurface, lev = 0, ilev = 0;
  int /* surface, height2m, height10m,*/ n2;
  int gridtype, zaxistype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
977

Uwe Schulzweida's avatar
Uwe Schulzweida committed
978
979
  /* create file */
  nce(nc_create(filename, NC_CLOBBER, &nc_file_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
980

Uwe Schulzweida's avatar
Uwe Schulzweida committed
981
  nce(nc_put_att_text(nc_file_id, NC_GLOBAL, "file_type", strlen(strfiletype_res), strfiletype_res));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
982

Uwe Schulzweida's avatar
Uwe Schulzweida committed
983
984
  n2 = 2;
  complex = 2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
985

Uwe Schulzweida's avatar
Uwe Schulzweida committed
986
987
988
989
990
  lon = 0; lat = 0; nfc = 0; lev = 0; belowsurface = 0;
  for ( varid = 0; varid < nvars; ++varid )
    {
      gridtype  = vars[varid].gridtype;
      zaxistype = vars[varid].zaxistype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
991

Uwe Schulzweida's avatar
Uwe Schulzweida committed
992
993
994
995
996
997
998
999
1000
1001
1002
      if ( gridtype == GRID_GAUSSIAN && lat == 0 )
	{
	  gridIDgp = vars[varid].gridID;
	  lon = gridInqXsize(vars[varid].gridID);
	  lat = gridInqYsize(vars[varid].gridID);
	  nhgl = lat/2;
	}
      else if ( gridtype == GRID_FOURIER && nfc == 0 )
	{
	  nfc = gridInqSize(vars[varid].gridID);
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1003

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
      if ( zaxistype == ZAXIS_HYBRID && lev == 0 )
	{
	  zaxisIDml = vars[varid].zaxisID;
	  lev = zaxisInqSize(vars[varid].zaxisID);
	  nlevp1 = lev + 1;
	  ilev = lev + 1;
	  nvclev = ilev;
	}
      else if ( zaxistype == ZAXIS_DEPTH_BELOW_LAND && belowsurface == 0 )
	{
	  belowsurface = zaxisInqSize(vars[varid].zaxisID);
	}
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1017

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1018
1019
1020
  if ( lat == 0 ) cdoAbort("Gaussian grid not found!");
  if ( nfc == 0 ) cdoAbort("Fourier data not found!");
  if ( lev == 0 ) cdoAbort("Hybrid level not found!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1021

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1022
1023
  nmp1 = (nfc/nhgl)/2;
  spc  = nmp1*(nmp1+1)/2;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1024

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1025
1026
1027
1028
1029
  nlon = lon;
  nlat = lat;
  nlev = lev;
  
  nce(nc_def_dim(nc_file_id, "lon", lon, &lon_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1030

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1031
  nce(nc_def_dim(nc_file_id, "lat", lat, &lat_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1032

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1033
  nce(nc_def_dim(nc_file_id, "nhgl", nhgl, &nhgl_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1034

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1035
  nce(nc_def_dim(nc_file_id, "nlevp1", nlevp1, &nlevp1_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1036

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1037
  nce(nc_def_dim(nc_file_id, "spc", spc, &spc_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1038

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1039
  nce(nc_def_dim(nc_file_id, "nvclev", nvclev, &nvclev_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1040

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1041
  nce(nc_def_dim(nc_file_id, "complex", complex, &complex_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1042

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1043
  nce(nc_def_dim(nc_file_id, "nmp1", nmp1, &nmp1_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1044

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1045
  nce(nc_def_dim(nc_file_id, "belowsurface", belowsurface, &belowsurface_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1046

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1047
  nce(nc_def_dim(nc_file_id, "lev", lev, &lev_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1048

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1049
  nce(nc_def_dim(nc_file_id, "ilev", ilev, &ilev_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1050

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1051
1052
  /*
  nce(nc_def_dim(nc_file_id, "surface", surface, &surface_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1053

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1054
  nce(nc_def_dim(nc_file_id, "height2m", height2m, &height2m_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1055

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1056
1057
1058
  nce(nc_def_dim(nc_file_id, "height10m", height10m, &height10m_dimid));
  */
  nce(nc_def_dim(nc_file_id, "n2", n2, &n2_dimid));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1059

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1060
  nce(nc_enddef(nc_file_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1061
1062


Uwe Schulzweida's avatar
Uwe Schulzweida committed
1063
1064
1065
1066
  for ( varid = 0; varid < nvars; varid++ )
    {
      gridtype  = vars[varid].gridtype;
      zaxistype = vars[varid].zaxistype;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1067

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1068
1069
1070
      if ( gridtype == GRID_FOURIER && zaxistype == ZAXIS_HYBRID )
	{
	  nvals = nfc;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1071

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1072
1073
1074
1075
	  dimidsp[0] = nhgl_dimid;
	  dimidsp[1] = nmp1_dimid;
	  dimidsp[2] = complex_dimid;
	  dimidsp[3] = lev_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1076

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1077
1078
1079
	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 4, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1080

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1081
1082
1083
1084
	  for ( i = 0; i < nlev; i++ )
	    {
	      start[0] = 0;     start[1] = 0;    start[2] = 0;  start[3] = i;
	      count[0] = nhgl;  count[1] = nmp1; count[2] = 2;  count[3] = 1; 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1085

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1086
1087
1088
1089
1090
1091
	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
	}
      else if ( gridtype == GRID_FOURIER && zaxistype == ZAXIS_HYBRID_HALF )
	{
	  nvals = nfc;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1092

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1093
1094
1095
1096
	  dimidsp[0] = nhgl_dimid;
	  dimidsp[1] = nmp1_dimid;
	  dimidsp[2] = complex_dimid;
	  dimidsp[3] = ilev_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1097

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1098
1099
1100
	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 4, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1101

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1102
1103
1104
1105
	  for ( i = 0; i < ilev; i++ )
	    {
	      start[0] = 0;     start[1] = 0;    start[2] = 0;  start[3] = i;
	      count[0] = nhgl;  count[1] = nmp1; count[2] = 2;  count[3] = 1; 
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1106

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1107
1108
1109
1110
1111
1112
	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
	}
      else if ( gridtype == GRID_GAUSSIAN && zaxistype == ZAXIS_HYBRID )
	{
	  nvals = lat*lon;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1113

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1114
1115
1116
	  dimidsp[0] = lat_dimid;
	  dimidsp[1] = lev_dimid;
	  dimidsp[2] = lon_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1117

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1118
1119
1120
	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 3, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1121

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1122
1123
1124
1125
	  for ( i = 0; i < nlev; i++ )
	    {
	      start[0] = 0;    start[1] = i;  start[2] = 0;
	      count[0] = lat;  count[1] = 1;  count[2] = lon;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1126

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1127
1128
1129
1130
1131
1132
	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
	}
      else if ( gridtype == GRID_GAUSSIAN && zaxistype == ZAXIS_HYBRID_HALF )
	{
	  nvals = lat*lon;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1133

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1134
1135
1136
	  dimidsp[0] = lat_dimid;
	  dimidsp[1] = ilev_dimid;
	  dimidsp[2] = lon_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1137

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1138
1139
1140
1141
1142
	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 3, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));

	  for ( i = 0; i < ilev; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1143
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1144
1145
1146
1147
1148
	      start[0] = 0;    start[1] = i;  start[2] = 0;
	      count[0] = lat;  count[1] = 1;  count[2] = lon;

	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1149
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1150
      else if ( gridtype == GRID_GAUSSIAN && zaxistype == ZAXIS_DEPTH_BELOW_LAND )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1151
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1152
	  nvals = lat*lon;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1153

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1154
1155
1156
	  dimidsp[0] = lat_dimid;
	  dimidsp[1] = belowsurface_dimid;
	  dimidsp[2] = lon_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1157

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1158
1159
1160
	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 3, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1161

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1162
	  for ( i = 0; i < belowsurface; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1163
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1164
1165
1166
1167
1168
	      start[0] = 0;    start[1] = i;  start[2] = 0;
	      count[0] = lat;  count[1] = 1;  count[2] = lon;

	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1169
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1170
      else if ( gridtype == GRID_GAUSSIAN && zaxistype == ZAXIS_GENERIC )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1171
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1172
1173
1174
1175
1176
1177
1178
1179
1180
1181
1182
	  nvals = lat*lon;

	  dimidsp[0] = lat_dimid;
	  dimidsp[1] = n2_dimid;
	  dimidsp[2] = lon_dimid;

	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 3, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));

	  for ( i = 0; i < n2; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1183
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1184
1185
	      start[0] = 0;    start[1] = i;  start[2] = 0;
	      count[0] = lat;  count[1] = 1;  count[2] = lon;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1186

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1187
1188
1189
1190
1191
1192
	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
	}
      else if ( gridtype == GRID_GAUSSIAN && zaxistype == ZAXIS_HYBRID_HALF )
	{
	  nvals = lat*lon;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1193

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1194
1195
1196
	  dimidsp[0] = lat_dimid;
	  dimidsp[1] = ilev_dimid;
	  dimidsp[2] = lon_dimid;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1197

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1198
1199
1200
	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 3, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1201

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1202
	  for ( i = 0; i < ilev; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1203
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1204
1205
1206
1207
1208
	      start[0] = 0;    start[1] = i;  start[2] = 0;
	      count[0] = lat;  count[1] = 1;  count[2] = lon;

	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1209
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1210
      else if ( gridtype == GRID_GAUSSIAN && zaxistype == ZAXIS_SURFACE )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1211
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1212
1213
1214
1215
1216
1217
1218
1219
1220
1221
1222
1223
1224
1225
1226
1227
1228
1229
1230
1231
	  dimidsp[0] = lat_dimid;
	  dimidsp[1] = lon_dimid;

	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 2, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));
	  nce(nc_put_var_double(nc_file_id, nc_var_id, vars[varid].ptr));
	}
      else if ( gridtype == GRID_GENERIC && zaxistype == ZAXIS_HYBRID )
	{
	  nvals = nhgl;

	  dimidsp[0] = nhgl_dimid;
	  dimidsp[1] = lev_dimid;

	  nce(nc_redef(nc_file_id));
	  nce(nc_def_var(nc_file_id, vars[varid].name, NC_DOUBLE, 2, dimidsp, &nc_var_id));
	  nce(nc_enddef(nc_file_id));

	  for ( i = 0; i < nlev; i++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1232
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1233
1234
1235
1236
1237
	      start[0] = 0;      start[1] = i;
	      count[0] = nvals;  count[1] = 1;

	      nce(nc_put_vara_double(nc_file_id, nc_var_id, start, count, vars[varid].ptr+i*nvals));
	    }  
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1238
1239
1240
	}
      else
	{
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1241
	  fprintf(stderr, "structure of %s unsupported!\n", vars[varid].name);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1242
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1243
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1244
1245
1246
1247
1248
1249
1250
1251
1252
1253
1254
1255
1256
1257
1258
1259
1260
1261
1262
1263
1264
1265
1266
1267
1268
1269
1270
1271
1272
1273
1274
1275
1276
1277
1278
1279
1280
1281
1282
1283
1284
1285
1286
1287
1288
1289

  /* define gaussian grid */

  nlon = lon;
  nlat = lat;

  xvals = (double *) malloc(nlon*sizeof(double));
  yvals = (double *) malloc(nlat*sizeof(double));

  gridInqXvals(gridIDgp, xvals);
  gridInqYvals(gridIDgp, yvals);


  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "lat", NC_DOUBLE, 1, &lat_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, yvals));
   
  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "lon", NC_DOUBLE, 1, &lon_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, xvals));

  free(xvals);
  free(yvals);

  /* define model level */

  nlev = lev;
  nvct = nvclev*2;

  /* vct   = (double *) malloc(nvct*sizeof(double)); */

  vct = zaxisInqVctPtr(zaxisIDml);

  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "vct_a", NC_DOUBLE, 1, &nvclev_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, vct));

  nce(nc_redef(nc_file_id));
  nce(nc_def_var(nc_file_id, "vct_b", NC_DOUBLE, 1, &nvclev_dimid, &nc_var_id));
  nce(nc_enddef(nc_file_id));
  nce(nc_put_var_double(nc_file_id, nc_var_id, vct+nlevp1));

  /* free(vct); */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1290
1291
1292
1293
1294

  /*close input file */
  nce(nc_close(nc_file_id));

#else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1295
  cdoAbort("netCDF support not compiled in!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1296
1297
1298
1299
1300
1301
1302
1303
#endif
}


void *Echam5ini(void *argument)
{
  static char func[] = "Echam5ini";
  int operatorID;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1304
  int operfunc;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1305
1306
  int READ_E5ML,  READ_E5RES;
  int WRITE_E5ML, WRITE_E5RES;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1307
1308
1309
1310
1311
1312
1313
1314
1315
1316
1317
  int streamID1, streamID2 = CDI_UNDEFID;
  int nrecs = 0;
  int recID, varID, levelID;
  int vlistID1, vlistID2;
  int nvars = 0;
  int iv, nlev;
  int gridsize, nmiss;
  int taxisID, tsID;

  cdoInitialize(argument);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1318
1319
1320
1321
  READ_E5ML   = cdoOperatorAdd("read_e5ml",   func_read,  0, NULL);
  READ_E5RES  = cdoOperatorAdd("read_e5res",  func_read,  0, NULL);
  WRITE_E5ML  = cdoOperatorAdd("write_e5ml",  func_write, 0, NULL);
  WRITE_E5RES = cdoOperatorAdd("write_e5res", func_write, 0, NULL);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1322
1323

  operatorID = cdoOperatorID();
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1324
  operfunc = cdoOperatorFunc(operatorID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1325

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1326
  if ( operfunc == func_read )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1327
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1328
      VAR *vars = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1329
1330
1331
1332
      ATTS atts;
      int iatt;

      iniatts(&atts);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1333

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1334
1335
1336
      if ( operatorID == READ_E5ML )
	nvars = read_e5ml(cdoStreamName(0), &vars);
      else if ( operatorID == READ_E5RES )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1337
	nvars = read_e5res(cdoStreamName(0), &vars, &atts);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1338
      else
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1339
	cdoAbort("Operator not implemented!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1340

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1341
      if ( nvars == 0 ) cdoAbort("Unsupported file type!");
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1342
1343
      
      vlistID2 = vlistCreate();
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1344
      vlistDefNtsteps(vlistID2, 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1345
1346

      for ( iv = 0; iv < nvars; iv++ )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1347
1348
1349
	{/*
	  fprintf(stderr, "%d %s %d %d %d %d\n", iv, vars[iv].name, vars[iv].gridID, vars[iv].zaxisID, gridInqSize(vars[iv].gridID), zaxisInqSize(vars[iv].zaxisID));
	 */
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1350
	  varID = vlistDefVar(vlistID2, vars[iv].gridID, vars[iv].zaxisID, TIME_CONSTANT);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1351
1352
1353
1354
1355
	  if ( vars[iv].code > 0 ) vlistDefVarCode(vlistID2, varID, vars[iv].code);
	  if ( vars[iv].name )     vlistDefVarName(vlistID2, varID, vars[iv].name);
	  if ( vars[iv].longname ) vlistDefVarLongname(vlistID2, varID, vars[iv].longname);
          if ( vars[iv].units )    vlistDefVarUnits(vlistID2, varID, vars[iv].units);
	  vlistDefVarDatatype(vlistID2, varID, DATATYPE_FLT64);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1356
1357
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1358
1359
1360
1361
1362
1363
      for ( iatt = 0; iatt < atts.natxt; ++iatt )
	{
	  /* printf("%s: %s\n", atts.atxtname[iatt], atts.atxtentry[iatt]); */
	  vlistDefAttribute(vlistID2, atts.atxtname[iatt], atts.atxtentry[iatt]);
	}

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1364
1365
1366
      taxisID = taxisCreate(TAXIS_ABSOLUTE);
      vlistDefTaxis(vlistID2, taxisID);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1367
1368
1369
      if ( cdoDefaultFileType == CDI_UNDEFID )
	cdoDefaultFileType = FILETYPE_NC;

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1370
1371
1372
1373
1374
1375
1376
1377
1378
1379
1380
1381
1382
1383
1384
      streamID2 = streamOpenWrite(cdoStreamName(1), cdoFiletype());
      if ( streamID2 < 0 ) cdiError(streamID2, "Open failed on %s", cdoStreamName(1));

      streamDefVlist(streamID2, vlistID2);

      tsID = 0;
      streamDefTimestep(streamID2, tsID);

      for ( varID = 0; varID < nvars; varID++ )
	{
	  gridsize = vars[varID].gridsize;
	  nlev     = vars[varID].nlev;

	  for ( levelID = 0; levelID < nlev; levelID++ )
	    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1385
	      streamDefRecord(streamID2, varID, levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1386
1387
1388
1389
1390
1391
1392
1393
	      streamWriteRecord(streamID2, vars[varID].ptr+levelID*gridsize, 0);
	    }
	}

      streamClose(streamID2);

      vlistDestroy(vlistID2);
    }
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1394
  else if ( operfunc == func_write )
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1395
    {
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1396
      VAR *vars = NULL;
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1397
1398
      int gridID, zaxisID, gridtype, zaxistype, gridsize, nlev;
      char name[256];
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1399
1400
1401
1402
1403
1404

      streamID1 = streamOpenRead(cdoStreamName(0));
      if ( streamID1 < 0 ) cdiError(streamID1, "Open failed on %s", cdoStreamName(0));

      vlistID1 = streamInqVlist(streamID1);

Uwe Schulzweida's avatar
Uwe Schulzweida committed
1405
1406
1407
1408
1409
1410
1411
1412
1413
1414
1415
1416
1417
1418
1419
1420
1421
1422
1423
1424
1425
1426
1427
1428
1429
1430
1431
1432
      nvars = vlistNvars(vlistID1);

      vars = (VAR *) malloc(nvars*sizeof(VAR));

      for ( varID = 0; varID < nvars; ++varID )
	{
	  vlistInqVarName(vlistID1, varID, name);

	  gridID = vlistInqVarGrid(vlistID1, varID);
	  zaxisID = vlistInqVarZaxis(vlistID1, varID);

	  gridtype = gridInqType(gridID);
	  zaxistype = zaxisInqType(zaxisID);

	  gridsize = gridInqSize(gridID);
	  nlev     = zaxisInqSize(zaxisID);

	  inivar(&vars[varID], gridtype, zaxistype,  0, name, NULL, NULL);
	  
	  vars[varID].gridID = gridID;
	  vars[varID].zaxisID   = zaxisID;
	  vars[varID].gridsize  = gridsize;
	  vars[varID].nlev      = nlev;

	  vars[varID].ptr = (double *) malloc(nlev*gridsize*sizeof(double));
	}

      nrecs = streamInqTimestep(streamID1, 0);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1433
1434
1435
1436

      for ( recID = 0; recID < nrecs; recID++ )
	{
	  streamInqRecord(streamID1, &varID, &levelID);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1437
1438
1439
1440
1441

	  gridID = vlistInqVarGrid(vlistID1, varID);
	  gridsize = gridInqSize(gridID);

	  streamReadRecord(streamID1, vars[varID].ptr+levelID*gridsize, &nmiss);
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1442
	}
Uwe Schulzweida's avatar
Uwe Schulzweida committed
1443
1444

      streamClose(streamID1);