Subversion Repositories pvinversion.ecmwf

Rev

Go to most recent revision | Details | Last modification | View Log | RSS feed

Rev Author Line No. Line
3 michaesp 1
      subroutine wricst(cstnam,datar,aklev,bklev,aklay,bklay,stdate)
2
C------------------------------------------------------------------------
3
 
4
C     Creates the constants file for NetCDF files containing ECMWF
5
C     data. The constants file is compatible with the one created
6
C     for EM data (with subroutine writecst).
7
C
8
C     Input parameters:
9
C
10
C     cstnam    name of constants file
11
C     datar     array contains all required parameters to write file
12
C               datar(1):       number of points along x        
13
C               datar(2):       number of points along y
14
C               datar(3):       maximum latitude of data region (ymax)
15
C               datar(4):       minimum longitude of data region (xmin)
16
C               datar(5):       minimum latitude of data region (ymin)
17
C               datar(6):       maximum longitude of data region (xmax)
18
C               datar(7):       grid increment along x
19
C               datar(8):       grid increment along y
20
C               datar(9):       number of levels        
21
C		datar(10):	data type (forecast or analysis)
22
C		datar(11):	data version
23
C		datar(12):	constants file version
24
C		datar(13):	longitude of pole of coordinate system
25
C		datar(14):	latitude of pole of coordinate system
26
C     aklev     array contains the aklev values
27
C     bklev	array contains the bklev values
28
C     aklay     array contains the aklay values
29
C     bklay     array contains the bklay values
30
C     stdate    array contains date (year,month,day,time,step) of first
31
C               field on file (start-date), dimensionised as stdate(5)
32
C------------------------------------------------------------------------
33
 
34
 
35
      include "netcdf.inc"
36
 
37
      integer   nchar,maxlev
38
 
39
      parameter (nchar=20,maxlev=32)
40
      real	aklev(maxlev),bklev(maxlev)
41
      real      aklay(maxlev),bklay(maxlev)
42
      real	pollat,latmin,latmax
43
      integer   datar(14)
44
      integer	stdate(5)
45
      character*80 cstnam
46
 
47
C     declarations for constants-variables
48
 
49
      integer   nz
50
      integer   dattyp, datver, cstver
51
 
52
C     further declarations
53
 
54
      integer	ierr			! error flag
55
      integer	cdfid			! NetCDF id
56
      integer	xid,yid,zid		! dimension ids
57
      integer	pollonid, pollatid,	! variable ids
58
     >		aklevid, bklevid, aklayid, bklayid,
59
     >		lonminid, lonmaxid, latminid, latmaxid,
60
     >		dellonid, dellatid,
61
     >		startyid, startmid, startdid, starthid, startsid,
62
     >		dattypid, datverid, cstverid
63
 
64
      nz=datar(9)			! number of levels
65
 
66
C     Set data-type and -version, version of cst-file-format
67
 
68
      dattyp=datar(10)
69
      datver=datar(11)
70
      cstver=datar(12)
71
 
72
C     Initially set error to false
73
 
74
      ierr=0
75
 
76
C     Create constants file
77
 
78
      cdfid=nccre(trim(cstnam),NCCLOB,ierr)
79
 
80
C     Define the dimensions
81
 
82
      xid = ncddef (cdfid,'nx',datar(1),ierr)
83
      yid = ncddef (cdfid,'ny',datar(2),ierr)
84
      zid = ncddef (cdfid,'nz',datar(9),ierr)
85
 
86
C     Define integer constants
87
 
88
      pollonid = ncvdef(cdfid,'pollon', NCFLOAT,0,0,ierr)
89
      pollatid = ncvdef(cdfid,'pollat', NCFLOAT,0,0,ierr)
90
 
91
      aklevid = ncvdef (cdfid, 'aklev', NCFLOAT, 1, zid, ierr)
92
      bklevid = ncvdef (cdfid, 'bklev', NCFLOAT, 1, zid, ierr)
93
      aklayid = ncvdef (cdfid, 'aklay', NCFLOAT, 1, zid, ierr)
94
      bklayid = ncvdef (cdfid, 'bklay', NCFLOAT, 1, zid, ierr)
95
 
96
      lonminid = ncvdef (cdfid, 'lonmin', NCFLOAT, 0, 0, ierr)
97
      lonmaxid = ncvdef (cdfid, 'lonmax', NCFLOAT, 0, 0, ierr)
98
      latminid = ncvdef (cdfid, 'latmin', NCFLOAT, 0, 0, ierr)
99
      latmaxid = ncvdef (cdfid, 'latmax', NCFLOAT, 0, 0, ierr)
100
      dellonid = ncvdef (cdfid, 'dellon', NCFLOAT, 0, 0, ierr)
101
      dellatid = ncvdef (cdfid, 'dellat', NCFLOAT, 0, 0, ierr)
102
      startyid = ncvdef (cdfid, 'starty', NCLONG, 0, 0, ierr)
103
      startmid = ncvdef (cdfid, 'startm', NCLONG, 0, 0, ierr)
104
      startdid = ncvdef (cdfid, 'startd', NCLONG, 0, 0, ierr)
105
      starthid = ncvdef (cdfid, 'starth', NCLONG, 0, 0, ierr)
106
      startsid = ncvdef (cdfid, 'starts', NCLONG, 0, 0, ierr)
107
      dattypid = ncvdef (cdfid, 'dattyp', NCLONG, 0, 0, ierr)
108
      datverid = ncvdef (cdfid, 'datver', NCLONG, 0, 0, ierr)
109
      cstverid = ncvdef (cdfid, 'cstver', NCLONG, 0, 0, ierr)
110
 
111
C     Leave define mode
112
 
113
      call ncendf(cdfid,ierr)
114
 
115
C     Store levels
116
      call ncvpt(cdfid, aklevid, 1, nz, aklev, ierr)
117
      call ncvpt(cdfid, bklevid, 1, nz, bklev, ierr)
118
      call ncvpt(cdfid, aklayid, 1, nz, aklay, ierr)
119
      call ncvpt(cdfid, bklayid, 1, nz, bklay, ierr)
120
 
121
C     Store position of pole (trivial for ECMWF data)
122
      call ncvpt1(cdfid, pollonid, 1, real(datar(13))/1000., ierr)
123
      if (datar(14).gt.0) then
124
        pollat=min(real(datar(14))/1000.,90.)
125
      else
126
        pollat=max(real(datar(14))/1000.,-90.)
127
      endif
128
      call ncvpt1(cdfid, pollatid, 1, pollat, ierr)
129
 
130
C     Store horizontal data borders and grid increments
131
      call ncvpt1(cdfid, lonminid, 1, real(datar(4))/1000., ierr)
132
      call ncvpt1(cdfid, lonmaxid, 1, real(datar(6))/1000., ierr)
133
      latmin=max(real(datar(5))/1000.,-90.)
134
      latmax=min(real(datar(3))/1000.,90.)
135
      call ncvpt1(cdfid, latminid, 1, latmin, ierr)
136
      call ncvpt1(cdfid, latmaxid, 1, latmax, ierr)
137
      call ncvpt1(cdfid, dellonid, 1, real(datar(7))/1000., ierr)
138
      call ncvpt1(cdfid, dellatid, 1, real(datar(8))/1000., ierr)
139
 
140
C     Store date of first field on file (start-date)
141
      call ncvpt1(cdfid, startyid, 1, stdate(1), ierr)
142
      call ncvpt1(cdfid, startmid, 1, stdate(2), ierr)
143
      call ncvpt1(cdfid, startdid, 1, stdate(3), ierr)
144
      call ncvpt1(cdfid, starthid, 1, stdate(4), ierr)
145
      call ncvpt1(cdfid, startsid, 1, stdate(5), ierr)
146
 
147
C     Store datatype and version
148
      call ncvpt1(cdfid, dattypid, 1, dattyp, ierr)
149
      call ncvpt1(cdfid, datverid, 1, datver, ierr)
150
 
151
C     Store version of the constants file format
152
      call ncvpt1(cdfid, cstverid, 1, cstver, ierr)
153
 
154
C     Store strings
155
 
156
      call ncclos(cdfid,ierr)
157
      return
158
 
159
      end
160
      subroutine writelmcst(cdfid,nx,ny,nz,pollon,pollat,lonmin,
161
     &lonmax,latmin,latmax,dellon,dellat,dattyp,datver,cstver,
162
     &psref,tstar,tbeta,pintf,p0top,idate)
163
c     ------------------------------------------------------------------
164
 
165
      implicit none
166
 
167
      integer   cdfid
168
 
169
c     deklarationen der constants-variablen
170
      real       pollon,pollat
171
      real       lonmin,lonmax,latmin,latmax,dellon,dellat
172
      integer    idate(5)
173
      integer    nx,ny,nz
174
      integer    dattyp, datver, cstver
175
      real       psref, tstar, tbeta, pintf, p0top
176
 
177
      include 'netcdf.inc'
178
 
179
* netcdf declaration
180
      integer   iret, k
181
* dimension ids
182
      integer  nxdim, nydim, nzdim
183
* variable ids
184
      integer  startyid, startmid, startdid, starthid
185
* variable shapes, corners and edge lengths
186
      integer dims(1), corner(1), edges(1)
187
 
188
* enter define mode
189
      call ncredf(cdfid, iret)
190
 
191
      startyid = ncvdef (cdfid, 'starty', NCLONG, 0, 0, iret)
192
      startmid = ncvdef (cdfid, 'startm', NCLONG, 0, 0, iret)
193
      startdid = ncvdef (cdfid, 'startd', NCLONG, 0, 0, iret)
194
      starthid = ncvdef (cdfid, 'starth', NCLONG, 0, 0, iret)
195
 
196
* store the rest as global attributes
197
* store nx,ny,nz
198
      call ncapt(cdfid,NCGLOBAL,'nx',NCLONG,1,nx,iret)
199
      call ncapt(cdfid,NCGLOBAL,'ny',NCLONG,1,ny,iret)
200
      call ncapt(cdfid,NCGLOBAL,'nz',NCLONG,1,nz,iret)
201
 
202
* store pollon, pollat
203
      call ncapt(cdfid,NCGLOBAL,'pollon',NCFLOAT,1,pollon,iret)
204
      call ncapt(cdfid,NCGLOBAL,'pollat',NCFLOAT,1,pollat,iret)
205
 
206
* store lonmin, etc
207
      call ncapt(cdfid,NCGLOBAL,'lonmin',NCFLOAT,1,lonmin,iret)
208
      call ncapt(cdfid,NCGLOBAL,'lonmax',NCFLOAT,1,lonmax,iret)
209
      call ncapt(cdfid,NCGLOBAL,'latmin',NCFLOAT,1,latmin,iret)
210
      call ncapt(cdfid,NCGLOBAL,'latmax',NCFLOAT,1,latmax,iret)
211
      call ncapt(cdfid,NCGLOBAL,'dellon',NCFLOAT,1,dellon,iret)
212
      call ncapt(cdfid,NCGLOBAL,'dellat',NCFLOAT,1,dellat,iret)
213
 
214
* store data type and version
215
      call ncapt(cdfid,NCGLOBAL,'dattyp',NCLONG,1,dattyp,iret)
216
      call ncapt(cdfid,NCGLOBAL,'datver',NCLONG,1,datver,iret)
217
      call ncapt(cdfid,NCGLOBAL,'cstver',NCLONG,1,cstver,iret)
218
 
219
* store information of lm model vertical grid
220
      call ncapt(cdfid,NCGLOBAL,'psref',NCFLOAT,1,psref,iret)
221
      call ncapt(cdfid,NCGLOBAL,'tstar',NCFLOAT,1,tstar,iret)
222
      call ncapt(cdfid,NCGLOBAL,'tbeta',NCFLOAT,1,tbeta,iret)
223
      call ncapt(cdfid,NCGLOBAL,'pintf',NCFLOAT,1,pintf,iret)
224
      call ncapt(cdfid,NCGLOBAL,'p0top',NCFLOAT,1,p0top,iret)
225
 
226
* leave define mode
227
      call ncendf(cdfid, iret)
228
 
229
* store starty, etc
230
      corner(1) = 1
231
      edges(1) = 1
232
      call ncvpt(cdfid, startyid, corner, edges, idate(1), iret)
233
      call ncvpt(cdfid, startmid, corner, edges, idate(2), iret)
234
      call ncvpt(cdfid, startdid, corner, edges, idate(3), iret)
235
      call ncvpt(cdfid, starthid, corner, edges, idate(4), iret)
236
 
237
      end
238
      subroutine globcst(cdfnam,datar,aklev,bklev,aklay,bklay,stdate)
239
C------------------------------------------------------------------------
240
C+
241
C NAME:
242
C     subroutine globcst
243
C
244
C PURPOSE:
245
C     instead of writing a constants-file (*_cst), the information
246
C     is added to the netCDF file as global variables
247
C     the data format is compatible with the one requested by
248
C     the IVE ETH/MIT version, contact author about details
249
C
250
C CATEGORY:
251
C     model,netCDF
252
C
253
C CALLING SEQUENCE:
254
C     subroutine globcst(cdfnam,datar,aklev,bklev,aklay,bklay,stdate)
255
C
256
C INPUTS:
257
C     cdfnam    name of netCDF file
258
C               The file needs to exist, otherwise an ERROR occurs,
259
C               i.e. nothing is done
260
C     datar     array contains all required parameters to write file
261
C               datar(1):       number of points along x
262
C               datar(2):       number of points along y
263
C               datar(3):       maximum latitude of data region (ymax)
264
C               datar(4):       minimum longitude of data region (xmin)
265
C               datar(5):       minimum latitude of data region (ymin)
266
C               datar(6):       maximum longitude of data region (xmax)
267
C               datar(7):       grid increment along x
268
C               datar(8):       grid increment along y
269
C               datar(9):       number of levels
270
C               datar(10):      data type (forecast or analysis)
271
C               datar(11):      data version
272
C               datar(12):      constants file version
273
C               datar(13):      longitude of pole of coordinate system
274
C               datar(14):      latitude of pole of coordinate system
275
C     aklev     array contains the aklev values
276
C     bklev     array contains the bklev values
277
C     aklay     array contains the aklay values
278
C     bklay     array contains the bklay values
279
C     stdate    array contains date (year,month,day,time,step) of first
280
C               field on file (start-date), dimensionised as stdate(5)
281
C     list    the griblist-ASCII-file
282
C     varno   the GRIB code number
283
C
284
C OUTPUTS:
285
C     Adds cdf-information to EXISTING netCDF-file
286
C
287
C MODIFICATION HISTORY:
288
C
289
C     June  93    Christoph Schaer (ETHZ) created
290
C     Nov   93    Heini Wernli (ETHZ) wricst
291
C     Nov   98    David N. Bresch (MIT) wricst to globcst
292
C-
293
 
294
C     Sun include statement.
295
      include "netcdf.inc"
296
 
297
      integer   nchar,maxlev
298
 
299
      parameter (nchar=20,maxlev=32)
300
      real      aklev(maxlev),bklev(maxlev)
301
      real      aklay(maxlev),bklay(maxlev)
302
      integer   datar(14)
303
      integer   stdate(5)
304
      character*80 cdfnam
305
 
306
C     declarations for constants-variables
307
 
308
      integer   nz
309
      integer   dattyp, datver, cstver
310
 
311
C     further declarations
312
 
313
      integer   ierr                    ! error flag
314
      integer   cdfid                   ! NetCDF id
315
      integer   xid,yid,zid             ! dimension ids
316
      integer   pollonid, pollatid,     ! variable ids
317
     >          aklevid, bklevid, aklayid, bklayid,
318
     >          lonminid, lonmaxid, latminid, latmaxid,
319
     >          dellonid, dellatid,
320
     >          startyid, startmid, startdid, starthid, startsid,
321
     >          dattypid, datverid, cstverid
322
 
323
      nz=datar(9)                       ! number of levels
324
 
325
C     Set data-type and -version, version of cst-file-format
326
 
327
      dattyp=datar(10)
328
      datver=datar(11)
329
      cstver=datar(12)
330
 
331
C     Initially set error to false
332
 
333
      ierr=0
334
 
335
C     open the netCDF-file:
336
 
337
      call cdfwopn(cdfnam,cdfid,ierr)
338
      if (ierr.ne.0) then
339
         print*,'ERROR opening netCDF-file ',cdfnam
340
         return
341
      endif
342
 
343
C     Put file into define mode
344
      call ncredf(cdfid,ierr)
345
      if (ierr.ne.0) then
346
         print*,'ERROR switching to netCDF redefine mode'
347
         return
348
      endif
349
 
350
C     Define the dimensions
351
 
352
      xid = ncddef (cdfid,'nx',datar(1),ierr)
353
      yid = ncddef (cdfid,'ny',datar(2),ierr)
354
      zid = ncddef (cdfid,'nz',datar(9),ierr)
355
 
356
C     Define integer constants
357
 
358
      pollonid = ncvdef(cdfid,'pollon', NCFLOAT,0,0,ierr)
359
      pollatid = ncvdef(cdfid,'pollat', NCFLOAT,0,0,ierr)
360
 
361
      aklevid = ncvdef (cdfid, 'aklev', NCFLOAT, 1, zid, ierr)
362
      bklevid = ncvdef (cdfid, 'bklev', NCFLOAT, 1, zid, ierr)
363
      aklayid = ncvdef (cdfid, 'aklay', NCFLOAT, 1, zid, ierr)
364
      bklayid = ncvdef (cdfid, 'bklay', NCFLOAT, 1, zid, ierr)
365
 
366
      lonminid = ncvdef (cdfid, 'lonmin', NCFLOAT, 0, 0, ierr)
367
      lonmaxid = ncvdef (cdfid, 'lonmax', NCFLOAT, 0, 0, ierr)
368
      latminid = ncvdef (cdfid, 'latmin', NCFLOAT, 0, 0, ierr)
369
      latmaxid = ncvdef (cdfid, 'latmax', NCFLOAT, 0, 0, ierr)
370
      dellonid = ncvdef (cdfid, 'dellon', NCFLOAT, 0, 0, ierr)
371
      dellatid = ncvdef (cdfid, 'dellat', NCFLOAT, 0, 0, ierr)
372
      startyid = ncvdef (cdfid, 'starty', NCLONG, 0, 0, ierr)
373
      startmid = ncvdef (cdfid, 'startm', NCLONG, 0, 0, ierr)
374
      startdid = ncvdef (cdfid, 'startd', NCLONG, 0, 0, ierr)
375
      starthid = ncvdef (cdfid, 'starth', NCLONG, 0, 0, ierr)
376
      startsid = ncvdef (cdfid, 'starts', NCLONG, 0, 0, ierr)
377
      dattypid = ncvdef (cdfid, 'dattyp', NCLONG, 0, 0, ierr)
378
      datverid = ncvdef (cdfid, 'datver', NCLONG, 0, 0, ierr)
379
      cstverid = ncvdef (cdfid, 'cstver', NCLONG, 0, 0, ierr)
380
 
381
C     Leave define mode
382
 
383
      call ncendf(cdfid,ierr)
384
      if (ierr.ne.0) then
385
         print*,'ERROR exiting define mode'
386
         return
387
      endif
388
 
389
C     Store levels
390
      call ncvpt(cdfid, aklevid, 1, nz, aklev, ierr)
391
      call ncvpt(cdfid, bklevid, 1, nz, bklev, ierr)
392
      call ncvpt(cdfid, aklayid, 1, nz, aklay, ierr)
393
      call ncvpt(cdfid, bklayid, 1, nz, bklay, ierr)
394
 
395
C     Store position of pole (trivial for ECMWF data)
396
      call ncvpt1(cdfid, pollonid, 1, real(datar(13))/1000., ierr)
397
      call ncvpt1(cdfid, pollatid, 1, real(datar(14))/1000., ierr)
398
 
399
C     Store horizontal data borders and grid increments
400
      call ncvpt1(cdfid, lonminid, 1, real(datar(4))/1000., ierr)
401
      call ncvpt1(cdfid, lonmaxid, 1, real(datar(6))/1000., ierr)
402
      call ncvpt1(cdfid, latminid, 1, real(datar(5))/1000., ierr)
403
      call ncvpt1(cdfid, latmaxid, 1, real(datar(3))/1000., ierr)
404
      call ncvpt1(cdfid, dellonid, 1, real(datar(7))/1000., ierr)
405
      call ncvpt1(cdfid, dellatid, 1, real(datar(8))/1000., ierr)
406
 
407
C     Store date of first field on file (start-date)
408
      call ncvpt1(cdfid, startyid, 1, stdate(1), ierr)
409
      call ncvpt1(cdfid, startmid, 1, stdate(2), ierr)
410
      call ncvpt1(cdfid, startdid, 1, stdate(3), ierr)
411
      call ncvpt1(cdfid, starthid, 1, stdate(4), ierr)
412
      call ncvpt1(cdfid, startsid, 1, stdate(5), ierr)
413
 
414
C     Store datatype and version
415
      call ncvpt1(cdfid, dattypid, 1, dattyp, ierr)
416
      call ncvpt1(cdfid, datverid, 1, datver, ierr)
417
 
418
C     Store version of the constants file format
419
      call ncvpt1(cdfid, cstverid, 1, cstver, ierr)
420
 
421
      if (ierr.ne.0) then
422
         print*,'ERROR adding cst-date as global variables'
423
         return
424
      endif
425
 
426
C     Store strings
427
 
428
      call ncclos(cdfid,ierr)
429
      if (ierr.ne.0) then
430
         print*,'ERROR closing netCDF file'
431
      endif
432
 
433
      return
434
      end
435
 
436
      subroutine getsdat(cdfid,varnam,time,ix,iy,iz,sx,sy,sz,dat,error)
437
c-----------------------------------------------------------------------
438
c     Purpose:
439
c        This routine is called to read the data within a selected
440
c	 domain of a variable from an IVE-NetCDF file.
441
c        Prior to calling this routine, the file must be opened with
442
c        a call to opncdf (for extension) or crecdf (for creation) or
443
c        readcdf (for readonly).
444
c     Arguments:
445
c        cdfid   int   input   file-identifier
446
c                              (must be obtained by calling routine
447
c                              opncdf,readcdf  or crecdf)
448
c        varnam  char  input   the user-supplied variable name
449
c        time    real  input   the user-supplied time-level of the
450
c                              data to be read from the file (the time-
451
c                              levels stored in the file can be obtained
452
c                              with a call to gettimes).
453
c        ix/y/z  int   input   indices of lower left corner of selected
454
c			       data volume.
455
c	 sx/y/z  int   input   size of selected data volume
456
c        dat     real  output  data-array with dimensions (sx,sy,sz).
457
c        error   int   output  indicates possible errors found in this
458
c                              routine.
459
c                              error = 0   no errors detected.
460
c                              error = 1   the variable is not present on
461
c                                          the file.
462
c                              error = 2   the value of 'time' is not
463
c                                          known.to the file.
464
c			       error = 6,7,8   data volume too large
465
c                              error =10   another error.
466
c     History:
467
c       June  93    Christoph Schaer (ETHZ)  Created getdat
468
c	Nov   93    Heini Wernli (ETHZ)	     Created getsdat
469
c-----------------------------------------------------------------------
470
 
471
      include "netcdf.inc"
472
 
473
C     Declaration of local variables
474
      character*(*) varnam
475
      character*(20) chars
476
      integer cdfid
477
 
478
      integer     ix,iy,iz,sx,sy,sz
479
      real        dat(sx,sy,sz)
480
      real        misdat,varmin(4),varmax(4),stag(4)
481
      real        time, timeval
482
 
483
      integer     corner(4),edgeln(4),didtim,vardim(4),ndims
484
      integer     error, ierr
485
      integer     ntime
486
      integer     idtime,idvar,iflag
487
      integer     i
488
 
489
      call ncpopt(NCVERBOS)
490
 
491
c     access the variable
492
      call getdef (cdfid, trim(varnam), ndims, misdat,
493
     &                           vardim, varmin, varmax, stag, ierr)
494
      if (ierr.ne.0) then
495
        print *,'*ERROR* in getdef in getdat'
496
        error=1
497
        return
498
      endif
499
      idvar=ncvid(cdfid,trim(varnam),ierr)
500
      if (ierr.ne.0) then
501
        print *,'*ERROR* in ncvid in getsdat'
502
        error=1
503
        return
504
      endif
505
 
506
C     Get times-array
507
      didtim=ncdid(cdfid,'time',ierr)
508
      if (ierr.ne.0) then
509
        print *,'*ERROR* didtim in getsdat'
510
        error=10
511
        return
512
      endif
513
      call ncdinq(cdfid,didtim,chars,ntime,ierr)
514
      if (ierr.ne.0) then
515
        print *,'*ERROR* in ncdinq in getsdat'
516
        error=10
517
        return
518
      endif
519
      idtime=ncvid(cdfid,'time',ierr)
520
      if (ierr.ne.0) then
521
        print *,'*ERROR* in ncvid for time in getsdat'
522
        error=10
523
        return
524
      endif
525
c     find appropriate time-index
526
      iflag=0
527
      do i=1,ntime
528
        call ncvgt1(cdfid,idtime,i,timeval,ierr)
529
        if (ierr.ne.0) print *,'*ERROR* in ncvgt1 in getsdat'
530
        if (time.eq.timeval) iflag=i
531
      enddo
532
      if (iflag.eq.0) then
533
        error=2
534
        print *,'Error: Unknown time in getsdat'
535
        print *,time,timeval
536
        return
537
      endif
538
 
539
C     Define data volume to be written (index space)
540
      corner(1)=ix
541
      corner(2)=iy
542
      corner(3)=iz
543
      corner(4)=iflag
544
      edgeln(1)=sx
545
      edgeln(2)=sy
546
      edgeln(3)=sz
547
      edgeln(4)=1
548
 
549
C     Check if data volume is within data domain
550
 
551
      if (ix+sx-1.gt.vardim(1)) then
552
        error=7
553
        print *,'Error: data volume too large in x-direction'
554
        print *,ix,sx,vardim(1)
555
        return
556
      endif
557
      if (iy+sy-1.gt.vardim(2)) then
558
        error=8
559
        print *,'Error: data volume too large in y-direction'
560
        return
561
      endif
562
      if (iz+sz-1.gt.vardim(3)) then
563
        error=9
564
        print *,'Error: data volume too large in z-direction'
565
        return
566
      endif
567
 
568
C     Read data from NetCDF file
569
 
570
      call ncvgt(cdfid,idvar,corner,edgeln,dat,error)
571
      if (error.ne.0) then
572
        print *, 'corner ',corner(1),corner(2),corner(3)
573
        print *, 'edgeln ',edgeln(1),edgeln(2),edgeln(3)
574
        print *, '*ERROR* in ncvgt in getsdat'
575
        error=10
576
      endif
577
      end
578
      subroutine getlevs(cstid,nlev,aklev,bklev,aklay,bklay,error)
579
c-----------------------------------------------------------------------
580
c     Purpose:
581
c     	This routine is called to get the level arrays aklev and
582
c	bklev from a NetCDF constants file.
583
c     Arguments:
584
c	cstid     int	input   identifier for NetCDF constants file
585
c	nlev	  int	input	number of levels
586
c	aklev     real	output  array contains all aklev values
587
c       bklev     real  output  array contains all bklev values
588
c	aklay	  real  output	array contains all aklay values
589
c	bklay	  real	output	array contains all bklay values
590
c	error	  int	output	error flag
591
c				error = 0   no errors detected
592
c				error = 1   error detected
593
c     History:
594
c	Aug. 93	  Heini Wernli		Created.
595
c-----------------------------------------------------------------------
596
 
597
      integer   error
598
 
599
      integer   cstid
600
      integer   ncdid,ncvid		! NetCDF functions
601
      integer   didz,idak,idbk,idaky,idbky
602
      integer   nlev
603
      real      aklev(nlev),bklev(nlev),aklay(nlev),bklay(nlev)
604
      character*(20) dimnam
605
      integer   i
606
 
607
      didz	=ncdid(cstid,'nz',error)
608
      if (error.ne.0) goto 920
609
      idak	=ncvid(cstid,'aklev',error)
610
      if (error.ne.0) goto 920
611
      idbk	=ncvid(cstid,'bklev',error)
612
      if (error.ne.0) goto 920
613
      idaky     =ncvid(cstid,'aklay',error)
614
      if (error.ne.0) goto 920
615
      idbky     =ncvid(cstid,'bklay',error)
616
      if (error.ne.0) goto 920
617
 
618
      call ncdinq(cstid,didz,dimnam,nlev,error)	! read number of levels
619
      if (error.ne.0) goto 920
620
 
621
      do 10 i=1,nlev
622
        call ncvgt1(cstid,idak,i,aklev(i),error)      ! get aklev
623
        call ncvgt1(cstid,idbk,i,bklev(i),error)      ! get bklev
624
        call ncvgt1(cstid,idaky,i,aklay(i),error)      ! get aklay
625
        call ncvgt1(cstid,idbky,i,bklay(i),error)      ! get bklay
626
        if (error.ne.0) goto 920
627
   10 continue
628
 
629
      return
630
 
631
c     Error exits.
632
  920 write(*,*)'*ERROR*: An error occured in subroutine getlevs'
633
      return
634
 
635
      end
636
      subroutine getntim(cdfid,ntimes,ierr)
637
C------------------------------------------------------------------------
638
C     Purpose:
639
C        Get number of times on the specified NetCDF file
640
C     Arguments:
641
C        cdfid  int  input   identifier for NetCDF file
642
C        ntimes int  output  number of times on the file
643
C        error  int  output  errorflag
644
C     History:
645
C        Heini Wernli, ETHZ
646
C------------------------------------------------------------------------
647
 
648
      include "netcdf.inc"
649
 
650
      integer   ierr
651
      integer didtim,ntimes
652
 
653
      integer   cdfid,idtime
654
      integer   ncopts
655
      character*(20) dimnam
656
 
657
c     Get current value of error options, and make sure netCDF-errors do
658
c     not abort execution
659
      call ncgopt (ncopts)
660
      call ncpopt(NCVERBOS)
661
 
662
      didtim=ncdid(cdfid,'time',ierr)   ! inquire id for time dimension
663
      if (ierr.ne.0) goto 900
664
      idtime=ncvid(cdfid,'time',ierr)   ! inquire id for time array
665
      if (ierr.ne.0) goto 900
666
      call ncdinq(cdfid,didtim,dimnam,ntimes,ierr)      ! inquire # of times
667
      if (ierr.ne.0) goto 900
668
 
669
c     normal exit
670
      call ncpopt (ncopts)
671
      return
672
 
673
c     error exit
674
 900  ntimes=1
675
      call ncpopt (ncopts)
676
      end
677
      subroutine getstart(cdfid,idate,ierr)
678
C------------------------------------------------------------------------
679
C     Purpose:
680
C	Get start date for fields on specified NetCDF file
681
C     Arguments:
682
C	cdfid	int	input	identifier for NetCDF file
683
C	idate	int	output	array contains date (year,month,day,time,step)
684
C				dimensioned as idate(5)
685
C	ierr	int	output	error flag
686
C------------------------------------------------------------------------
687
 
688
      include "netcdf.inc"
689
 
690
c     variable declarations
691
      integer   ierr
692
      integer   idate(5)
693
      integer   cdfid,ncopts,idvar,nvars
694
      integer   ndims,ngatts,recdim,i,vartyp,nvatts,vardim(4)
695
      character*20 vnam(100)
696
 
697
c     Get current value of error options, and make sure NetCDF-errors do
698
c     not abort execution
699
      call ncgopt (ncopts)
700
      call ncpopt (NCVERBOS)
701
 
702
      idvar=ncvid(cdfid,'starty',ierr)
703
      if (ierr.ne.0) goto 930
704
      call ncvgt1(cdfid,idvar,1,idate(1),ierr)
705
      if (ierr.ne.0) goto 920
706
 
707
      idvar=ncvid(cdfid,'startm',ierr)
708
      if (ierr.ne.0) goto 920
709
      call ncvgt1(cdfid,idvar,1,idate(2),ierr)
710
      if (ierr.ne.0) goto 920
711
 
712
      idvar=ncvid(cdfid,'startd',ierr)
713
      if (ierr.ne.0) goto920
714
      call ncvgt1(cdfid,idvar,1,idate(3),ierr)
715
      if (ierr.ne.0) goto 920
716
 
717
      idvar=ncvid(cdfid,'starth',ierr)
718
      if (ierr.ne.0) goto 920
719
      call ncvgt1(cdfid,idvar,1,idate(4),ierr)
720
      if (ierr.ne.0) goto 920
721
 
722
C     Starts is not defined on all files
723
C     Only ask for it if it exists
724
C     Inquire number of dimensions, variables and attributes
725
 
726
      idate(5)=0
727
      call ncinq(cdfid,ndims,nvars,ngatts,recdim,ierr)
728
      do i=1,nvars
729
        call ncvinq(cdfid,i,vnam(i),vartyp,ndims,vardim,nvatts,ierr)
730
        if (vnam(i).eq.'starts') then
731
          idvar=ncvid(cdfid,'starts',ierr)
732
          call ncvgt1(cdfid,idvar,1,idate(5),ierr)
733
          if (ierr.ne.0) goto 920
734
        endif
735
      enddo
736
 
737
c     normal exit
738
      call ncpopt (ncopts)
739
      return
740
 
741
c     error exit
742
 920  continue
743
      write (6, *) 'ERROR: An error occurred while attempting to ',
744
     &             'read the starting-time in subroutine putstart.'
745
 930  continue
746
      call ncpopt (ncopts)
747
 
748
      end
749
      subroutine putstart(cdfid,idate,ierr)
750
C----------------------------------------------------------------------
751
C     Purpose:
752
C        Puts the 'starting-time' on the specified NetCDF file.
753
C     Arguments:
754
C        cdfid   int     input   identifier for NetCDF file
755
C        idate   int     input   array contains date (year,month,day,time,step)
756
C                                dimensioned as idate(5)
757
C        ierr    int     output  error flag
758
C------------------------------------------------------------------------
759
 
760
      include "netcdf.inc"
761
 
762
c     variable declarations
763
      integer   ierr,idate(5),startid(5),cdfid,ncopts,i
764
 
765
c     Get current value of error options, and make sure NetCDF-errors do
766
c     not abort execution
767
      call ncgopt (ncopts)
768
      call ncpopt (NCVERBOS)
769
 
770
c     define variables
771
      call ncredf(cdfid,ierr)
772
      if (ierr.ne.0) goto 920
773
      startid(1) = ncvdef (cdfid, 'starty', NCLONG, 0, 0, ierr)
774
      if (ierr.ne.0) goto 920
775
      startid(2) = ncvdef (cdfid, 'startm', NCLONG, 0, 0, ierr)
776
      if (ierr.ne.0) goto 920
777
      startid(3) = ncvdef (cdfid, 'startd', NCLONG, 0, 0, ierr)
778
      if (ierr.ne.0) goto 920
779
      startid(4) = ncvdef (cdfid, 'starth', NCLONG, 0, 0, ierr)
780
      if (ierr.ne.0) goto 920
781
      startid(5) = ncvdef (cdfid, 'starts', NCLONG, 0, 0, ierr)
782
      if (ierr.ne.0) goto 920
783
      call ncendf(cdfid, ierr)
784
      if (ierr.ne.0) goto 920
785
 
786
c     store variables
787
      do i=1,5
788
        call ncvpt1(cdfid,startid(i),1,idate(i),ierr)
789
        if (ierr.ne.0) goto 920
790
      enddo
791
 
792
c     synchronyse output to disk, revert to previous error-mode, and exit
793
      call ncsnc (cdfid,ierr)
794
      call ncpopt (ncopts)
795
      return
796
 
797
c     error exit
798
 920  write (6, *) 'ERROR: An error occurred while attempting to ',
799
     &             'write the starting-time in subroutine putstart.'
800
      call ncpopt (ncopts)
801
      call ncclos (cdfid, ierr)
802
 
803
      end
804
      subroutine getgrid(cdfid,dx,dy,ierr)
805
C------------------------------------------------------------------------
806
C     Purpose:
807
C       Get grid increments for fields on specified NetCDF file
808
C     Arguments:
809
C       cdfid   int     input   identifier for NetCDF file
810
C	dx	real	output	grid increment along latitude
811
C	dy	real	output	grid increment along longitude
812
C       ierr    int     output  error flag
813
C------------------------------------------------------------------------
814
 
815
      integer   ierr
816
 
817
      integer   cdfid
818
      integer   ncvid
819
 
820
      integer   idilon,idilat
821
      real	dx,dy
822
 
823
      idilon    =ncvid(cdfid,'dellon',ierr)
824
      if (ierr.ne.0) return
825
      idilat    =ncvid(cdfid,'dellat',ierr)
826
      if (ierr.ne.0) return
827
 
828
      call ncvgt1(cdfid,idilon,1,dx,ierr)
829
      if (ierr.ne.0) return
830
      call ncvgt1(cdfid,idilat,1,dy,ierr)
831
      if (ierr.ne.0) return
832
 
833
      end
834
      subroutine getdattyp(cdfid,typ,ierr)
835
C------------------------------------------------------------------------
836
C     Purpose:
837
C       Get data type for specified NetCDF file
838
C     Arguments:
839
C       cdfid   int     input   identifier for NetCDF file
840
C       typ     int     output  data type: 1 (52) for pressure (theta) coord
841
C       ierr    int     output  error flag
842
C------------------------------------------------------------------------
843
 
844
      integer   ierr
845
 
846
      integer   cdfid
847
      integer   ncvid
848
 
849
      integer   idtyp,typ
850
 
851
      idtyp    =ncvid(cdfid,'dattyp',ierr)
852
      if (ierr.ne.0) return
853
 
854
      call ncvgt1(cdfid,idtyp,1,typ,ierr)
855
      if (ierr.ne.0) return
856
 
857
      end
858
      subroutine getpole(cdfid,pollon,pollat,ierr)
859
C------------------------------------------------------------------------
860
C     Purpose:
861
C       Get physical coordinates of pole of coordinate system
862
C     Arguments:
863
C       cdfid   int     input   identifier for NetCDF file
864
C	pollon	real	output	longitude of pole
865
C	pollat	real	output	latitude of pole
866
C       ierr    int     output  error flag
867
C------------------------------------------------------------------------
868
 
869
      integer   ierr
870
 
871
      integer   cdfid
872
      integer   ncvid
873
 
874
      integer   idplon,idplat
875
      real      pollon,pollat
876
 
877
      idplon    =ncvid(cdfid,'pollon',ierr)
878
      if (ierr.ne.0) return
879
      idplat    =ncvid(cdfid,'pollat',ierr)
880
      if (ierr.ne.0) return
881
 
882
      call ncvgt1(cdfid,idplon,1,pollon,ierr)
883
      if (ierr.ne.0) return
884
      call ncvgt1(cdfid,idplat,1,pollat,ierr)
885
      if (ierr.ne.0) return
886
 
887
      end
888
      subroutine getmc2grid(cdfid,polx,poly,delx,shem,phi0,lam0,ierr)
889
C------------------------------------------------------------------------
890
C     Purpose:
891
C       Get physical coordinates of pole of coordinate system
892
C     Arguments:
893
C       cdfid   int     input   identifier for NetCDF file
894
C       ierr    int     output  error flag
895
C------------------------------------------------------------------------
896
 
897
      integer   ierr
898
 
899
      integer   cdfid
900
      integer   ncvid
901
 
902
      integer   idpolx,idpoly,iddelx,idshem,idphi0,idlam0
903
      real      polx,poly,delx,shem,phi0,lam0
904
 
905
      idpolx    =ncvid(cdfid,'polx',ierr)
906
      if (ierr.ne.0) return
907
      idpoly    =ncvid(cdfid,'poly',ierr)
908
      if (ierr.ne.0) return
909
      iddelx	=ncvid(cdfid,'delx',ierr)
910
      if (ierr.ne.0) return
911
      idshem	=ncvid(cdfid,'shem',ierr)
912
      if (ierr.ne.0) return
913
      idphi0	=ncvid(cdfid,'phi0',ierr)
914
      if (ierr.ne.0) return
915
      idlam0	=ncvid(cdfid,'lam0',ierr)
916
      if (ierr.ne.0) return
917
 
918
      call ncvgt1(cdfid,idpolx,1,polx,ierr)
919
      if (ierr.ne.0) return
920
      call ncvgt1(cdfid,idpoly,1,poly,ierr)
921
      if (ierr.ne.0) return
922
      call ncvgt1(cdfid,iddelx,1,delx,ierr)
923
      if (ierr.ne.0) return
924
      call ncvgt1(cdfid,idshem,1,shem,ierr)
925
      if (ierr.ne.0) return
926
      call ncvgt1(cdfid,idphi0,1,phi0,ierr)
927
      if (ierr.ne.0) return
928
      call ncvgt1(cdfid,idlam0,1,lam0,ierr)
929
      if (ierr.ne.0) return
930
 
931
      end
932
 
933
 
934
      subroutine getcfn(cdfid,cfn,ierr)
935
C------------------------------------------------------------------------
936
C     Purpose:
937
C       Get name of constants file
938
C     Arguments:
939
C       cdfid   int     input   identifier for NetCDF file
940
C       cfn     char    output  name of constants file
941
C       ierr    int     output  error flag
942
C------------------------------------------------------------------------
943
 
944
      include 'netcdf.inc'
945
 
946
      integer   ierr
947
      integer   cdfid,lenstr
948
      character*80 cfn
949
 
950
      lenstr=80
951
      call ncagtc(cdfid,NCGLOBAL,"constants_file_name",cfn,lenstr,ierr)
952
      if (ierr.ne.0) write(*,*)'error in SR getcfn'
953
 
954
      end
955
 
956
 
957
      subroutine getdim (cdfid, varnam, nx, ny, nz, error)
958
c-------------------------------------------------------------------------
959
c     Purpose:
960
c        This routine is called to get the dimensions of
961
c        a variable from an IVE-NetCDF file for use with the IVE plotting
962
c        package. Prior to calling this routine, the file must be opened
963
c        with a call to opncdf.
964
c     Arguments:
965
c        cdfid   int   input   file-identifier
966
c                              (can be obtained by calling routine
967
c                              opncdf)
968
c        varnam  char  input   the user-supplied variable name.
969
c                              (can be obtained by calling routine
970
c                              opncdf)
971
c        nx      int   output  the zonal dimension of the variable.
972
c        ny      int   output  the meridional dimension of the variable.
973
c        nz      int   output  the vertical dimension of the variable.
974
c
975
c        error   int   output  indicates possible errors found in this
976
c                              routine.
977
c                              error = 0   no errors detected.
978
c                              error = 1   the variable is not on the file.
979
c                              error =10   other errors.
980
c     History:
981
c        March 2000    Heini Wernli (ETHZ)     Created.
982
c-----------------------------------------------------------------------
983
 
984
      include  "netcdf.inc"
985
 
986
c     Argument declarations.
987
      character *(*) varnam
988
      integer        vardim(4), ndim, error, cdfid
989
      integer        nx,ny,nz
990
 
991
c     Local variable declarations.
992
      character *(20) dimnam(MAXNCDIM),vnam
993
      integer         id,i,k
994
      integer         ndims,nvars,ngatts,recdim,dimsiz(MAXNCDIM)
995
      integer         vartyp,nvatts, ncopts
996
 
997
c     Get current value of error options.
998
      call ncgopt (ncopts)
999
 
1000
c     make sure NetCDF-errors do not abort execution
1001
      call ncpopt(NCVERBOS)
1002
 
1003
c     Initially set error to indicate no errors.
1004
      error = 0
1005
 
1006
c     inquire for number of dimensions
1007
      call ncinq(cdfid,ndims,nvars,ngatts,recdim,error)
1008
      if (error.eq.1) goto 920
1009
 
1010
c     read dimension-table
1011
      do i=1,ndims
1012
        call ncdinq(cdfid,i,dimnam(i),dimsiz(i),error)
1013
        if (error.gt.0) goto 920
1014
      enddo
1015
 
1016
c     get id of the variable
1017
      id=ncvid(cdfid,varnam,error)
1018
      if (error.eq.1) goto 910
1019
 
1020
c     inquire about variable
1021
      call ncvinq(cdfid,id,vnam,vartyp,ndim,vardim,nvatts,error)
1022
      if (vartyp.ne.NCFLOAT) error=1
1023
      if (error.gt.0) goto 920
1024
 
1025
c     get dimensions from dimension-table
1026
      do k=1,ndim
1027
        vardim(k)=dimsiz(vardim(k))
1028
      enddo
1029
 
1030
      nx=vardim(1)
1031
      ny=vardim(2)
1032
      nz=vardim(3)
1033
 
1034
c     normal exit
1035
      call ncpopt (ncopts)
1036
      return
1037
 
1038
c     Error exits.
1039
 910  write (6, *) '*ERROR*: The selected variable could not be found ',
1040
     &             'in the file by getdim.'
1041
      call ncpopt (ncopts)
1042
      call ncclos (cdfid, error)
1043
      return
1044
 
1045
 920  write (6, *) '*ERROR*: An error occurred while attempting to ',
1046
     &             'read the data file in subroutine getcdf.'
1047
      call ncpopt (ncopts)
1048
      call ncclos (cdfid, error)
1049
      end
1050
 
1051
 
1052
      subroutine new_gettra(cdfid,varnam,ix,ntimes,array,ierr)
1053
C------------------------------------------------------------------------
1054
C
1055
C     Reads the time-evolution for one grid-point of the variable
1056
C     indicated by varnam.
1057
C
1058
C     cdfid     int     input   identifier for NetCDF file
1059
C     varnam    char    input   name of variable
1060
C     ix        int     input   index for trajectory to read
1061
C     ntimes    int     input   number of time-indices to read
1062
C     array     real    output  array contains the readed values
1063
C     ierr      int     output  error flag
1064
C------------------------------------------------------------------------
1065
 
1066
C     Declaration of attributes
1067
 
1068
      integer   cdfid
1069
      character*(*) varnam
1070
      integer   ix
1071
      integer   ntimes
1072
      real      array(ntimes)
1073
 
1074
C     Declaration of local variables
1075
 
1076
      integer   corner(4),edgeln(4)
1077
      integer   idvar,ierr
1078
      integer   ncvid
1079
      integer	strend
1080
 
1081
      corner(1)=ix
1082
      corner(2)=1
1083
      corner(3)=1
1084
      corner(4)=1
1085
      edgeln(1)=1
1086
      edgeln(2)=1
1087
      edgeln(3)=1
1088
      edgeln(4)=ntimes
1089
 
1090
      idvar =ncvid(cdfid,varnam(1:strend(varnam)),ierr)
1091
      call ncvgt(cdfid,idvar,corner,edgeln,array,ierr)
1092
      if (ierr.ne.0) goto 991
1093
 
1094
      return
1095
  991 stop 'Variable not found on NetCDF file in SR new_gettra'
1096
      end
1097
      subroutine getvars(cdfid,nvars,vnam,ierr)
1098
C------------------------------------------------------------------------
1099
 
1100
C     Opens the NetCDF file 'filnam' and returns its identifier cdfid.
1101
 
1102
C     filnam    char    input   name of NetCDF file to open
1103
C     nvars     int     output  number of variables on file
1104
C     vnam	char	output  array with variable names
1105
C     ierr      int     output  error flag
1106
C------------------------------------------------------------------------
1107
 
1108
      include 'netcdf.inc'
1109
 
1110
      integer   cdfid,ierr,nvars
1111
      character*(*) vnam(*)
1112
 
1113
      integer	ndims,ngatts,recdim,i,vartyp,nvatts,vardim(4)
1114
 
1115
      call ncpopt(NCVERBOS)
1116
 
1117
C     Inquire number of dimensions, variables and attributes
1118
 
1119
      call ncinq(cdfid,ndims,nvars,ngatts,recdim,ierr)
1120
 
1121
C     Inquire variable names from NetCDF file
1122
 
1123
      do i=1,nvars
1124
        call ncvinq(cdfid,i,vnam(i),vartyp,ndims,vardim,nvatts,ierr)
1125
      enddo
1126
 
1127
      return
1128
      end