← Back to branch summary

~javier-lopez/ubuntu/quantal/cdo/sru-1023329

« back to all changes in this revision

Viewing changes to src/IFS2ICON.c

  • Committer: Bazaar Package Importer
  • Author(s): Alastair McKinstry
  • Date: 2010-09-22 15:58:09 UTC
  • mfrom: (1.2.3 upstream)
  • Revision ID: james.westby@ubuntu.com-20100922155809-u1d3vlmlqj02uxjt
Tags: 1.4.6.dfsg.1-1
New upstream release. 

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/IFS2ICON.c
 
1
/*
 
2
  This file is part of CDO. CDO is a collection of Operators to
 
3
  manipulate and analyse Climate model Data.
 
4
 
 
5
  Copyright (C) 2003-2010 Uwe Schulzweida, Uwe.Schulzweida@zmaw.de
 
6
  See COPYING file for copying and redistribution conditions.
 
7
 
 
8
  This program is free software; you can redistribute it and/or modify
 
9
  it under the terms of the GNU General Public License as published by
 
10
  the Free Software Foundation; version 2 of the License.
 
11
 
 
12
  This program is distributed in the hope that it will be useful,
 
13
  but WITHOUT ANY WARRANTY; without even the implied warranty of
 
14
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
15
  GNU General Public License for more details.
 
16
*/
 
17
 
 
18
 
 
19
#include <cdi.h>
 
20
#include "cdo.h"
 
21
#include "cdo_int.h"
 
22
#include "pstream.h"
 
23
#include "namelist.h"
 
24
 
 
25
#define  MAX_PARAM  512
 
26
 
 
27
typedef struct {
 
28
  char *name;
 
29
  int   pos;
 
30
  int   itype;
 
31
} param_t;
 
32
 
 
33
typedef struct {
 
34
  char nml_name[32];
 
35
  int debug_level;
 
36
  int force_remap;
 
37
} nml_control_t;
 
38
 
 
39
 
 
40
static
 
41
void nml_error(int status, const char *nml_name)
 
42
{
 
43
  switch (status)
 
44
    {
 
45
    case 1: 
 
46
      {
 
47
        cdoAbort("Namelist %s not found!", nml_name);
 
48
        break;
 
49
      }
 
50
    default:
 
51
      {
 
52
        cdoAbort("Namelist %s unknown error!", nml_name);
 
53
      }
 
54
    }
 
55
}
 
56
 
 
57
static
 
58
int pos_nml(FILE *fp_nml, const char *nml_name)
 
59
{
 
60
  int status = 1;
 
61
  fpos_t fpos;
 
62
  char line[MAX_LINE_LEN];
 
63
  
 
64
 
 
65
  rewind(fp_nml);
 
66
  
 
67
  fgetpos(fp_nml, &fpos);
 
68
  while ( readline(fp_nml, line, MAX_LINE_LEN) )
 
69
    {
 
70
      if ( line[0] == '#' ) continue;
 
71
      if ( line[0] == '\0' ) continue;
 
72
      if ( line[0] == '&' && strlen(line) > 2 )
 
73
        {
 
74
          strtolower(line+1);
 
75
          if ( strncmp(line+1, nml_name, strlen(nml_name)) == 0 )
 
76
            {
 
77
              status = 0;
 
78
              fsetpos(fp_nml, &fpos);
 
79
              break;
 
80
            }
 
81
        }
 
82
      fgetpos(fp_nml, &fpos);
 
83
    }
 
84
 
 
85
  return status;
 
86
}
 
87
 
 
88
static
 
89
void nml_control_init(nml_control_t *nml_control)
 
90
{
 
91
  strcpy(nml_control->nml_name, "control");
 
92
  nml_control->debug_level = 0;
 
93
  nml_control->force_remap = 0;
 
94
}
 
95
 
 
96
static
 
97
void nml_control_print(nml_control_t nml_control)
 
98
{ 
 
99
  fprintf(stdout, "Namelist: %s\n", nml_control.nml_name);
 
100
  fprintf(stdout, "  debug_level = %d\n", nml_control.debug_level);
 
101
  fprintf(stdout, "  force_remap = %d\n", nml_control.force_remap);
 
102
}
 
103
 
 
104
static
 
105
int nml_control_read(FILE *fp_nml, nml_control_t *nml_control)
 
106
{
 
107
  int status = 0;
 
108
 
 
109
  status = pos_nml(fp_nml, nml_control->nml_name);
 
110
 
 
111
  if ( status == 0 )
 
112
    {
 
113
      namelist_t *nml;
 
114
      NML_DEF_INT(debug_level, 1, 0);
 
115
      NML_DEF_INT(force_remap, 1, 0);
 
116
 
 
117
      nml = namelistNew(nml_control->nml_name);
 
118
 
 
119
      NML_ADD_INT(nml, debug_level);
 
120
      NML_ADD_INT(nml, force_remap);
 
121
      
 
122
      namelistRead(fp_nml, nml);
 
123
 
 
124
      if ( (NML_NUM(nml, debug_level)) ) nml_control->debug_level = seldebug_level[0];
 
125
      if ( (NML_NUM(nml, force_remap)) ) nml_control->force_remap = selforce_remap[0];
 
126
 
 
127
      namelistDelete(nml);
 
128
    }
 
129
 
 
130
  return status;
 
131
}
 
132
 
 
133
static
 
134
void read_param(const char *paramfile, param_t *param, int maxparam, int *nparam)
 
135
{
 
136
  static const char *func = "read_param";
 
137
  FILE *fp;
 
138
  namelist_t *nml;
 
139
  int nml_itype, nml_name, nml_pos;
 
140
  int locc, i;
 
141
  int code, new_code, table, pos;
 
142
  int nml_index = 0;
 
143
  int codenum, tabnum, levtype;
 
144
  int numparam = 0;
 
145
  char *datatype = NULL;
 
146
  char *name = NULL, *itype = NULL, *stdname = NULL, longname[256] = "", units[256] = "";
 
147
  char varname[256];
 
148
 
 
149
  fp = fopen(paramfile, "r");
 
150
  if ( fp == NULL ) cdoAbort("Open failed on %s!", paramfile);
 
151
 
 
152
  nml = namelistNew("parameter");
 
153
  nml->dis = 0;
 
154
 
 
155
  nml_name     = namelistAdd(nml, "name",          NML_WORD, 0, &name, 1);
 
156
  nml_pos      = namelistAdd(nml, "pos",           NML_INT,  0, &pos, 1);
 
157
  nml_itype    = namelistAdd(nml, "itype",         NML_WORD, 0, &itype, 1);
 
158
              
 
159
  while ( ! feof(fp) )
 
160
    {
 
161
      namelistReset(nml);
 
162
      namelistRead(fp, nml);
 
163
      
 
164
      locc = FALSE;
 
165
      for ( i = 0; i < nml->size; i++ )
 
166
        {
 
167
          if ( nml->entry[i]->occ ) { locc = TRUE; break; }
 
168
        }
 
169
      
 
170
      if ( locc )
 
171
        {
 
172
          /* namelistPrint(nml); */
 
173
          
 
174
          nml_index++;
 
175
 
 
176
          if ( nml->entry[nml_name]->occ == 0 )
 
177
            {
 
178
              cdoWarning("Parameter %d skipped, variable name not found!", nml_index);
 
179
              continue;
 
180
            }
 
181
 
 
182
          if ( numparam >= maxparam )
 
183
            {
 
184
              cdoWarning("Too many parameter (limit=%d)!", maxparam);
 
185
              break;
 
186
            }
 
187
 
 
188
          param[numparam].name = strdup(name);
 
189
 
 
190
          if ( nml->entry[nml_pos]->occ )
 
191
            param[numparam].pos  = pos;
 
192
          else
 
193
            param[numparam].pos  = 0;
 
194
 
 
195
          if ( nml->entry[nml_itype]->occ )
 
196
            param[numparam].itype = (int) *itype;
 
197
          else
 
198
            param[numparam].itype = 'Q';
 
199
        
 
200
          numparam++;
 
201
          /*
 
202
          for ( varID = 0; varID < nvars; varID++ )
 
203
            {
 
204
              vlistInqVarName(vlistID2, varID, varname);
 
205
              if ( strcmp(varname, name) == 0 ) break;
 
206
            }
 
207
 
 
208
          if ( varID < nvars )
 
209
            {
 
210
              if ( nml->entry[nml_code]->occ     ) vlistDefVarCode(vlistID2, varID, code);
 
211
              if ( nml->entry[nml_new_code]->occ ) vlistDefVarCode(vlistID2, varID, new_code);
 
212
              if ( nml->entry[nml_name]->occ     ) vlistDefVarName(vlistID2, varID, name);
 
213
              if ( nml->entry[nml_new_name]->occ ) vlistDefVarName(vlistID2, varID, new_name);
 
214
              if ( nml->entry[nml_stdname]->occ  ) vlistDefVarStdname(vlistID2, varID, stdname);
 
215
              if ( nml->entry[nml_longname]->occ ) vlistDefVarLongname(vlistID2, varID, longname);
 
216
              if ( nml->entry[nml_units]->occ    ) vlistDefVarUnits(vlistID2, varID, units);
 
217
            }
 
218
          else
 
219
            {
 
220
              if ( cdoVerbose )
 
221
                {
 
222
                  cdoPrint("Variable %s not found!", name);
 
223
                }
 
224
            }
 
225
          */
 
226
        }
 
227
      else
 
228
        break;
 
229
    }
 
230
  
 
231
  namelistDelete(nml);
 
232
 
 
233
  fclose(fp);
 
234
 
 
235
  *nparam = numparam;
 
236
}
 
237
 
 
238
 
 
239
void *IFS2ICON(void *argument)
 
240
{
 
241
  static const char *func = "IFS2ICON";
 
242
  int streamID1, streamID2 = CDI_UNDEFID;
 
243
  int nrecs, nvars, newval = -1, tabnum = 0;
 
244
  int tsID1, recID, varID, levelID;
 
245
  int vlistID1, vlistID2;
 
246
  int taxisID1, taxisID2;
 
247
  int nmiss;
 
248
  int gridsize;
 
249
  int index, zaxisID1, zaxisID2, nzaxis, nlevs;
 
250
  int tableID = -1;
 
251
  int tableformat = 0;
 
252
  int zaxistype;
 
253
  int i;
 
254
  int status;
 
255
  int nparam;
 
256
  const char *paramfile;
 
257
  const char *nmlfile;
 
258
  double newlevel = 0;
 
259
  double *levels = NULL;
 
260
  double *array = NULL;
 
261
  param_t param[MAX_PARAM];
 
262
  nml_control_t nml_control;
 
263
  FILE *fp_nml;
 
264
 
 
265
  cdoInitialize(argument);
 
266
 
 
267
  operatorInputArg("namelist file and parameter file");
 
268
  operatorCheckArgc(2);
 
269
 
 
270
  /* read namelists */
 
271
  nmlfile = operatorArgv()[0];
 
272
 
 
273
  fp_nml = fopen(nmlfile, "r");
 
274
  if ( fp_nml == NULL ) cdoAbort("Open failed on %s", nmlfile);
 
275
 
 
276
  nml_control_init(&nml_control);
 
277
  status = nml_control_read(fp_nml, &nml_control);
 
278
  if ( status > 0 ) nml_error(status, nml_control.nml_name);
 
279
 
 
280
  fclose(fp_nml);
 
281
 
 
282
  if ( cdoVerbose )
 
283
    {
 
284
      nml_control_print(nml_control);
 
285
    }
 
286
  
 
287
  paramfile = operatorArgv()[1];
 
288
  read_param(paramfile, param, MAX_PARAM, &nparam);
 
289
 
 
290
  for ( i = 0; i < nparam; ++i )
 
291
    printf("%s %d %c\n", param[i].name, param[i].pos, param[i].itype);
 
292
 
 
293
  streamID1 = streamOpenRead(cdoStreamName(0));
 
294
  if ( streamID1 < 0 ) cdiError(streamID1, "Open failed on %s", cdoStreamName(0));
 
295
 
 
296
  vlistID1 = streamInqVlist(streamID1);
 
297
  vlistID2 = vlistDuplicate(vlistID1);
 
298
  vlistPrint(vlistID2);
 
299
 
 
300
  taxisID1 = vlistInqTaxis(vlistID1);
 
301
  taxisID2 = taxisDuplicate(taxisID1);
 
302
  vlistDefTaxis(vlistID2, taxisID2);
 
303
 
 
304
  streamID2 = streamOpenWrite(cdoStreamName(1), cdoFiletype());
 
305
  if ( streamID2 < 0 ) cdiError(streamID2, "Open failed on %s", cdoStreamName(1));
 
306
 
 
307
  streamDefVlist(streamID2, vlistID2);
 
308
 
 
309
  gridsize = vlistGridsizeMax(vlistID2);
 
310
  array = (double *) malloc(gridsize*sizeof(double));
 
311
 
 
312
  tsID1 = 0;
 
313
  while ( (nrecs = streamInqTimestep(streamID1, tsID1)) )
 
314
    {
 
315
      taxisCopyTimestep(taxisID2, taxisID1);
 
316
 
 
317
      streamDefTimestep(streamID2, tsID1);
 
318
               
 
319
      for ( recID = 0; recID < nrecs; recID++ )
 
320
        {
 
321
          streamInqRecord(streamID1, &varID, &levelID);
 
322
          streamDefRecord(streamID2,  varID,  levelID);
 
323
          
 
324
          streamReadRecord(streamID1, array, &nmiss);
 
325
          streamWriteRecord(streamID2, array, nmiss);
 
326
        }
 
327
      tsID1++;
 
328
    }
 
329
 
 
330
  streamClose(streamID2);
 
331
  streamClose(streamID1);
 
332
 
 
333
  if ( array ) free(array);
 
334
 
 
335
  cdoFinish();
 
336
 
 
337
  return (0);
 
338
}