source: IOIPSL/trunk/src/flincom.f90 @ 386

Last change on this file since 386 was 386, checked in by bellier, 16 years ago

Added CeCILL License information

  • Property svn:keywords set to Id
File size: 57.0 KB
Line 
1MODULE flincom
2!-
3!$Id$
4!-
5! This software is governed by the CeCILL license
6! See IOIPSL/IOIPSL_License_CeCILL.txt
7!---------------------------------------------------------------------
8  USE netcdf
9!-
10  USE calendar,  ONLY : ju2ymds, ymds2ju, ioconf_calendar
11  USE errioipsl, ONLY : histerr
12  USE stringop,  ONLY : strlowercase
13!-
14  IMPLICIT NONE
15!-
16  PRIVATE
17  PUBLIC :: flinput, flincre, flinget, flinclo, &
18            flinopen, flininfo, flininspect, flinquery_var
19!-
20  INTERFACE flinopen
21!---------------------------------------------------------------------
22!- The "flinopen" routines will open an input file
23!-
24!- INPUT
25!-
26!- filename  : Name of the netCDF file to be opened
27!-
28!- iideb     : index i for zoom     !
29!- iilen     : length  of  zoom     !   for
30!- jjdeb     : index j for zoom     !   zoom
31!- jjlen     : length  of  zoom     !
32!-
33!- do_test   : A flag that enables the testing of the content
34!-             of the file against the input from the model
35!-
36!- INPUT if do_test=TRUE      OUTPUT else
37!-
38!- iim       : size in the x direction in the file (longitude)
39!- jjm       : size in the y direction
40!- llm       : number of levels
41!-             (llm = 0 means no axis to be expected)
42!- lon       : array of (iilen,jjlen) (zoom), or (iim,jjm) (no zoom),
43!-             that contains the longitude of each point
44!- lat       : same for latitude
45!- lev       : An array of llm for the latitude
46!-
47!- WARNING :
48!- In the case of do_test=FALSE it is for the user to check
49!- that the dimensions of lon lat and lev are correct when passed to
50!- flinopen. This can be done after the call when iim and jjm have
51!- been retrieved from the netCDF file. In F90 this problem will
52!- be solved with an internal assign
53!- IF iim, jjm, llm or ttm are parameters in the calling program and
54!- you use the option do_test=FALSE it will create a segmentation fault
55!-
56!-  OUTPUT
57!-
58!-  ttm       : size of time axis
59!-  itaus     : Time steps within this file
60!-  date0     : Julian date at which itau = 0
61!-  dt        : length of the time steps of the data
62!-  fid       : returned file ID which is later used to read the data
63!---------------------------------------------------------------------
64    MODULE PROCEDURE flinopen_zoom2d, flinopen_nozoom
65  END INTERFACE
66!-
67  INTERFACE flinput
68!---------------------------------------------------------------------
69!- The "flinput" routines will put a variable
70!- on the netCDF file created by flincre.
71!- If the sizes of the axis do not match the one of the IDs
72!- then a new axis is created.
73!- That is we loose the possibility of writting hyperslabs of data.
74!-
75!- Again here if iim = jjm = llm = ttm = 0
76!- then a global attribute is added to the file.
77!-
78!- INPUT
79!-
80!- fid      : Identification of the file in which we will write
81!- varname  : Name of variable to be written
82!- iim      : size in x of variable
83!- nlonid   : ID of x axis which could fit for this axis
84!- jjm      : size in y of variable
85!- nlatid   : ID of y axis which could fit for this axis
86!- llm      : size in z of variable
87!- zdimid   : ID of z axis which could fit for this axis
88!- ttm      : size in t of variable
89!- tdimid   : ID of t axis which could fit for this axis
90!-
91!- OUTPUT
92!-
93!- NONE
94!---------------------------------------------------------------------
95    MODULE PROCEDURE flinput_r4d, flinput_r3d, flinput_r2d, &
96                     flinput_r1d, flinput_scal
97  END INTERFACE
98!-
99  INTERFACE flinget
100    MODULE PROCEDURE flinget_r4d, flinget_r3d, flinget_r2d, &
101                     flinget_r1d, flinget_scal, &
102                     flinget_r4d_zoom2d, flinget_r3d_zoom2d, &
103                     flinget_r2d_zoom2d
104  END INTERFACE
105!-
106! This is the data we keep on each file we open
107!-
108  INTEGER, PARAMETER :: nbfile_max = 200
109  INTEGER, SAVE :: nbfiles = 0
110  INTEGER, SAVE :: ncids(nbfile_max), ncnbd(nbfile_max), &
111                   ncfunli(nbfile_max), ncnba(nbfile_max)
112  INTEGER, SAVE :: ncnbva(nbfile_max), ncdims(nbfile_max,4)
113  LOGICAL, SAVE :: ncfileopen(nbfile_max)=.FALSE.
114!-
115  INTEGER, SAVE :: cind_vid, cind_fid, cind_len
116  INTEGER,DIMENSION(:),ALLOCATABLE,SAVE :: cindex
117!-
118  INTEGER,DIMENSION(4) :: w_sta, w_len, w_dim
119!-
120CONTAINS
121!-
122!===
123!-
124SUBROUTINE flincre &
125  (filename, iim1, jjm1, lon1, lat1, llm1, lev1, ttm1, itaus, &
126   time0, dt, fid_out, nlonid1, nlatid1, zdimid1, tdimid1)
127!---------------------------------------------------------------------
128!- This is a "low level" subroutine for opening netCDF files wich
129!- contain the major coordinate system of the model.
130!- Other coordinates needed for other variables
131!- will be added as they are needed.
132!-
133!- INPUT
134!-
135!- filename    : Name of the file to be created
136!- iim1, jjm1  : Horizontal size of the grid
137!-               which will be stored in the file
138!- lon1, lat1  : Horizontal grids
139!- llm1        : Size of the vertical grid
140!- lev1        : Vertical grid
141!- ttm1        : Size of time axis
142!- itaus       : time steps on the time axis
143!- time0       : Time in julian days at which itau = 0
144!- dt          : time step in seconds between itaus
145!-               (one step of itau)
146!-
147!- OUTPUT
148!-
149!- fid         : File identification
150!- nlonid1     : Identification of longitudinal axis
151!- nlatid1     : Identification of latitudinal axis
152!- zdimid1     : ID of vertical axis
153!- tdimid1     : ID of time axis
154!---------------------------------------------------------------------
155  IMPLICIT NONE
156!-
157! ARGUMENTS
158!-
159  CHARACTER(LEN=*) :: filename
160  INTEGER :: iim1, jjm1, llm1, ttm1
161  REAL :: lon1(iim1,jjm1)
162  REAL :: lat1(iim1,jjm1)
163  REAL :: lev1(llm1)
164  INTEGER :: itaus(ttm1)
165  REAL :: time0
166  REAL :: dt
167  INTEGER :: fid_out, zdimid1, nlonid1, nlatid1, tdimid1
168!-
169! LOCAL
170!-
171  INTEGER :: iret, lll, fid
172  INTEGER :: lonid, latid, levid, timeid
173  INTEGER :: year, month, day
174  REAL :: sec
175  CHARACTER(LEN=250):: name
176!-
177  LOGICAL :: check = .FALSE.
178!---------------------------------------------------------------------
179  lll = LEN_TRIM(filename)
180  IF (filename(lll-2:lll) /= '.nc') THEN
181    name=filename(1:lll)//'.nc'
182  ELSE
183    name=filename(1:lll)
184  ENDIF
185!-
186  iret = NF90_CREATE (name, NF90_CLOBBER, fid)
187!-
188  iret = NF90_DEF_DIM (fid, 'x',     iim1, nlonid1)
189  iret = NF90_DEF_DIM (fid, 'y',     jjm1, nlatid1)
190  iret = NF90_DEF_DIM (fid, 'lev',   llm1, zdimid1)
191  iret = NF90_DEF_DIM (fid, 'tstep', ttm1, tdimid1)
192!-
193! Vertical axis
194!-
195  IF (check) WRITE(*,*) 'flincre Vertical axis'
196!-
197  iret = NF90_DEF_VAR (fid, 'lev', NF90_FLOAT, zdimid1, levid)
198  iret = NF90_PUT_ATT (fid, levid, 'units',     '-')
199  iret = NF90_PUT_ATT (fid, levid, 'title',     'levels')
200  iret = NF90_PUT_ATT (fid, levid, 'long_name', 'Sigma Levels')
201!-
202! Time axis
203!-
204  IF (check) WRITE(*,*) 'flincre time axis'
205!-
206  iret = NF90_DEF_VAR (fid, 'tstep', NF90_FLOAT, tdimid1, timeid)
207  iret = NF90_PUT_ATT (fid, timeid, 'units',     '-')
208  iret = NF90_PUT_ATT (fid, timeid, 'title',     'time')
209  iret = NF90_PUT_ATT (fid, timeid, 'long_name', 'time steps')
210!-
211! The longitude
212!-
213  IF (check) WRITE(*,*) 'flincre Longitude axis'
214!-
215  iret = NF90_DEF_VAR (fid, "nav_lon", NF90_FLOAT, &
216                       (/ nlonid1, nlatid1 /), lonid)
217  iret = NF90_PUT_ATT (fid, lonid, 'units', "degrees_east")
218  iret = NF90_PUT_ATT (fid, lonid, 'title', "Longitude")
219  iret = NF90_PUT_ATT (fid, lonid, 'nav_model', &
220                       "Lambert projection of PROMES")
221  iret = NF90_PUT_ATT (fid, lonid, 'valid_min', &
222                       REAL(MINVAL(lon1),KIND=4))
223  iret = NF90_PUT_ATT (fid, lonid, 'valid_max', &
224                       REAL(MAXVAL(lon1),KIND=4))
225!-
226! The Latitude
227!-
228  IF (check) WRITE(*,*) 'flincre Latitude axis'
229!-
230  iret = NF90_DEF_VAR (fid, "nav_lat", NF90_FLOAT, &
231                     (/ nlonid1, nlatid1 /), latid)
232  iret = NF90_PUT_ATT (fid, latid, 'units', "degrees_north")
233  iret = NF90_PUT_ATT (fid, latid, 'title', "Latitude")
234  iret = NF90_PUT_ATT (fid, latid, 'nav_model', &
235                       "Lambert projection of PROMES")
236  iret = NF90_PUT_ATT (fid, latid, 'valid_min', &
237                       REAL(MINVAL(lat1),KIND=4))
238  iret = NF90_PUT_ATT (fid, latid, 'valid_max', &
239                       REAL(MAXVAL(lat1),KIND=4))
240!-
241! The time coordinates
242!-
243  iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'delta_tstep_sec', &
244                       REAL(dt,KIND=4))
245!-
246  CALL ju2ymds (time0, year, month, day, sec)
247!-
248  iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'year0',  REAL(year,KIND=4))
249  iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'month0', REAL(month,KIND=4))
250  iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'day0',   REAL(day,KIND=4))
251  iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'sec0',   REAL(sec,KIND=4))
252!-
253  iret = NF90_ENDDEF (fid)
254!-
255  IF (check) WRITE(*,*) 'flincre Variable'
256!-
257  iret = NF90_PUT_VAR (fid, levid, lev1(1:llm1))
258!-
259  IF (check) WRITE(*,*) 'flincre Time Variable'
260!-
261  iret = NF90_PUT_VAR (fid, timeid, REAL(itaus(1:ttm1)))
262!-
263  IF (check) WRITE(*,*) 'flincre Longitude'
264!-
265  iret = NF90_PUT_VAR (fid, lonid, lon1(1:iim1,1:jjm1))
266!-
267  IF (check) WRITE(*,*) 'flincre Latitude'
268!-
269  iret = NF90_PUT_VAR (fid, latid, lat1(1:iim1,1:jjm1))
270!-
271! Keep all this information
272!-
273  nbfiles = nbfiles+1
274!-
275  IF (nbfiles > nbfile_max) THEN
276    CALL histerr (3,'flincre', &
277     'Too many files. Please increase nbfil_max', &
278     'in program flincom.F90.',' ')
279  ENDIF
280!-
281  ncids(nbfiles) = fid
282  ncnbd(nbfiles) = 4
283!-
284  ncdims(nbfiles,1:4) = (/ iim1, jjm1, llm1, ttm1 /)
285!-
286  ncfunli(nbfiles) = -1
287  ncnba(nbfiles)   =  4
288  ncnbva(nbfiles)  =  0
289  ncfileopen(nbfiles) = .TRUE.
290!-
291  fid_out = nbfiles
292!---------------------
293END SUBROUTINE flincre
294!-
295!===
296!-
297SUBROUTINE flinopen_zoom2d &
298  (filename, iideb, iilen, jjdeb, jjlen, do_test, &
299  iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out)
300!---------------------------------------------------------------------
301  IMPLICIT NONE
302!-
303! ARGUMENTS
304!-
305  CHARACTER(LEN=*) :: filename
306  LOGICAL :: do_test
307  INTEGER :: iim, jjm, llm, ttm, iideb, iilen, jjdeb, jjlen
308  REAL :: lon(iilen,jjlen), lat(iilen,jjlen), lev(llm)
309  INTEGER :: itaus(ttm)
310  REAL :: date0, dt
311  INTEGER :: fid_out
312!-
313  LOGICAL :: check = .FALSE.
314!---------------------------------------------------------------------
315  IF (check) WRITE (*,*) ' iideb, iilen, jjdeb, jjlen, iim, jjm ', &
316                           iideb, iilen, jjdeb, jjlen, iim, jjm
317  IF (check) WRITE (*,*) ' lon ', lon(1,1), lon(iilen,jjlen)
318  IF (check) WRITE (*,*) ' lat ', lat(1,1), lat(iilen,jjlen)
319!-
320  CALL flinopen_work &
321    (filename, iideb, iilen, jjdeb, jjlen, do_test, &
322     iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out)
323!-----------------------------
324END SUBROUTINE flinopen_zoom2d
325!-
326!===
327!-
328SUBROUTINE flinopen_nozoom &
329  (filename, do_test, iim, jjm, llm, lon, lat, lev, ttm, &
330   itaus, date0, dt, fid_out)
331!---------------------------------------------------------------------
332  IMPLICIT NONE
333!-
334! ARGUMENTS
335!-
336  CHARACTER(LEN=*) :: filename
337  LOGICAL :: do_test
338  INTEGER :: iim, jjm, llm, ttm
339  REAL :: lon(iim,jjm), lat(iim,jjm), lev(llm)
340  INTEGER :: itaus(ttm)
341  REAL :: date0, dt
342  INTEGER :: fid_out
343!---------------------------------------------------------------------
344  CALL flinopen_work &
345    (filename, 1, iim, 1, jjm, do_test, &
346     iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out)
347!-------------------------
348END SUBROUTINE flinopen_nozoom
349!-
350!===
351!-
352SUBROUTINE flinopen_work &
353  (filename, iideb, iilen, jjdeb, jjlen, do_test, &
354   iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out)
355!---------------------------------------------------------------------
356  IMPLICIT NONE
357!-
358! ARGUMENTS
359!-
360  CHARACTER(LEN=*) :: filename
361  LOGICAL :: do_test
362  INTEGER :: iim, jjm, llm, ttm, iideb, iilen, jjdeb, jjlen
363  REAL :: lon(iilen,jjlen), lat(iilen,jjlen), lev(llm)
364  INTEGER :: itaus(ttm)
365  REAL :: date0, dt
366  INTEGER :: fid_out
367!-
368! LOCAL
369!-
370  REAL, PARAMETER :: eps = 1.e-4
371!-
372  INTEGER :: iret, vid, fid, nbdim, i, iilast, jjlast
373  INTEGER :: gdtt_id, old_id, iv, gdtmaf_id
374  CHARACTER(LEN=250) :: name
375  CHARACTER(LEN=80) :: units, calendar
376  INTEGER :: tmp_iim, tmp_jjm, tmp_llm, tmp_ttm
377  REAL :: x_first, x_last
378  INTEGER :: year, month, day
379  REAL :: r_year, r_month, r_day
380  INTEGER :: year0, month0, day0, hours0, minutes0, seci
381  REAL :: sec, sec0
382  CHARACTER :: strc
383!-
384  REAL,DIMENSION(:),ALLOCATABLE :: vec_tmp
385!-
386  LOGICAL :: open_file
387  LOGICAL :: check = .FALSE.
388!---------------------------------------------------------------------
389  iilast = iideb+iilen-1
390  jjlast = jjdeb+jjlen-1
391  IF (check) WRITE (*,*) &
392    ' flinopen_work zoom 2D information '// &
393    ' iideb, iilen, iilast, jjdeb, jjlen, jjlast ', &
394      iideb, iilen, iilast, jjdeb, jjlen, jjlast
395!-
396! 1.0 get all infos on the file
397!-
398! Either the fid_out has not been initialized (0 or very large)
399! then we have to open anyway. Else we only need to open the file
400! if it has not been opened before.
401!-
402  IF ( (fid_out < 1).OR.(fid_out > nbfile_max) ) THEN
403    open_file = .TRUE.
404  ELSE IF (.NOT.ncfileopen(fid_out)) THEN
405    open_file = .TRUE.
406  ELSE
407    open_file = .FALSE.
408  ENDIF
409!-
410  IF (open_file) THEN
411    CALL flininfo (filename,tmp_iim,tmp_jjm,tmp_llm,tmp_ttm,fid_out)
412  ELSE
413!-- The user has already opened the file
414!-- and we trust that he knows the dimensions
415    tmp_iim = iim
416    tmp_jjm = jjm
417    tmp_llm = llm
418    tmp_ttm = ttm
419  ENDIF
420!-
421  IF (check) &
422    WRITE(*,*) 'OUT OF flininfo :',tmp_iim,tmp_jjm,tmp_llm,tmp_ttm
423!-
424  fid = ncids(fid_out)
425!-
426! 2.0 get the sizes and names of the different coordinates
427!     and do a first set of verification.
428!-
429! 2.2 We test the axis if we have to.
430!-
431  IF (check) &
432    WRITE(*,*) 'flininfo 2.2 We test if we have to test : ',do_test
433!-
434  IF (do_test) THEN
435    IF      (iim /= tmp_iim) THEN
436      CALL histerr (3,'flinopen', &
437        'file '//filename//' does not have the ', &
438        'required dimension in x direction (longitude)',' ')
439    ELSE IF (jjm /= tmp_jjm) THEN
440      CALL histerr (3,'flinopen', &
441        'file '//filename//' does not have the ', &
442        'required dimension in y direction (latitude)',' ')
443    ELSE IF ( llm /= tmp_llm .AND. llm > 0 ) THEN
444      CALL histerr (3,'flinopen', &
445        'file '//filename//' does not have the ', &
446        'required dimension in the vertical',' ')
447    ENDIF
448  ELSE
449!---
450!-- 2.3 Else the sizes of the axes are returned to the user
451!---
452    IF (check) WRITE(*,*) 'flinopen 2.3 Else sizes are returned'
453!---
454    iim = tmp_iim
455    jjm = tmp_jjm
456    llm = tmp_llm
457  ENDIF
458!-
459  ttm = tmp_ttm
460!-
461! 3.0 Check if we are realy talking about the same coodinate system
462!     if not then we get the lon, lat and lev variables from the file
463!-
464  IF (check) WRITE(*,*) 'flinopen 3.0 we are realy talking'
465!-
466  IF (do_test) THEN
467!---
468    CALL flinfindcood (fid_out, 'lon', vid, nbdim)
469    iret = NF90_GET_VAR (fid, vid, x_first, start=(/ iideb, jjdeb /))
470    iret = NF90_GET_VAR (fid, vid, x_last, start=(/ iilast, jjlast /))
471!---
472    IF (check) &
473      WRITE(*,*) 'from file lon first and last, modulo 360. ', &
474        x_first, x_last, MODULO(x_first,360.), MODULO(x_last,360.)
475    IF (check) &
476      WRITE(*,*) 'from model lon first and last, modulo 360. ', &
477        lon(1,1),lon(iilen,jjlen), &
478        MODULO(lon(1,1),360.), MODULO(lon(iilen,jjlen),360.)
479    IF (    (ABS( MODULO(x_first,360.) &
480                 -MODULO(lon(1,1),360.)) > eps) &
481        .OR.(ABS( MODULO(x_last,360.) &
482                 -MODULO(lon(iilen ,jjlen),360.)) > eps ) ) THEN
483      CALL histerr (3,'flinopen', &
484        'file '//filename//' and the model do not', &
485        'share the same longitude coordinate', &
486        'Obtained by comparing the first and last values ')
487    ENDIF
488!---
489    CALL flinfindcood (fid_out, 'lat', vid, nbdim)
490    iret = NF90_GET_VAR (fid, vid, x_first, start=(/ iideb, jjdeb /))
491    iret = NF90_GET_VAR (fid, vid, x_last, start=(/ iilast, jjlast /))
492!---
493    IF (check) WRITE(*,*) &
494      'from file lat first and last ',x_first,x_last
495    IF (check) WRITE(*,*) &
496      'from model lat first and last ',lat(1,1),lat(iilen,jjlen)
497!---
498    IF (    (ABS(x_first-lat(1,1)) > eps) &
499        .OR.(ABS(x_last-lat(iilen,jjlen)) > eps) ) THEN
500      CALL histerr (3,'flinopen', &
501        'file '//filename//' and the model do not', &
502        'share the same latitude coordinate', &
503        'Obtained by comparing the first and last values ')
504    ENDIF
505!---
506    IF (llm > 0) THEN
507      CALL flinfindcood (fid_out, 'lev', vid, nbdim)
508      iret = NF90_GET_VAR (fid, vid, x_first, start=(/ 1 /))
509      iret = NF90_GET_VAR (fid, vid, x_last, start=(/ llm /))
510!-----
511      IF (check) WRITE(*,*) &
512        'from file lev first and last ',x_first ,x_last
513      IF (check) WRITE(*,*) &
514        'from model lev first and last ',lev(1),lev(llm)
515!-----
516      IF (    (ABS(x_first-lev(1)) > eps) &
517          .OR.(ABS(x_last-lev(llm)) > eps) ) THEN
518        CALL histerr (3,'flinopen', &
519          'file '//filename//' and the model do not', &
520          'share the same vertical coordinate', &
521          'Obtained by comparing the first and last values')
522      ENDIF
523    ENDIF
524!---
525  ELSE
526!---
527!-- 4.0 extracting the coordinates if we do not check
528!---
529    IF (check) WRITE(*,*) 'flinopen 4.0 extracting the coordinates'
530!---
531    CALL flinfindcood (fid_out, 'lon', vid, nbdim)
532    IF (nbdim == 2) THEN
533      iret = NF90_GET_VAR (fid, vid, lon, &
534               start=(/ iideb, jjdeb /), count=(/ iilen, jjlen /))
535    ELSE
536      ALLOCATE(vec_tmp(iilen))
537      iret = NF90_GET_VAR (fid, vid, vec_tmp, &
538               start=(/ iideb /), count=(/ iilen /))
539      DO i=1,jjlen
540        lon(:,i) = vec_tmp(:)
541      ENDDO
542      DEALLOCATE(vec_tmp)
543    ENDIF
544!---
545    CALL flinfindcood (fid_out, 'lat', vid, nbdim)
546    IF (nbdim == 2) THEN
547      iret = NF90_GET_VAR (fid, vid, lat, &
548               start=(/ iideb, jjdeb /), count=(/ iilen, jjlen /))
549    ELSE
550      ALLOCATE(vec_tmp(jjlen))
551      iret = NF90_GET_VAR (fid, vid, vec_tmp, &
552               start=(/ jjdeb /), count=(/ jjlen /))
553      DO i=1,iilen
554        lat(i,:) = vec_tmp(:)
555      ENDDO
556      DEALLOCATE(vec_tmp)
557    ENDIF
558!---
559    IF (llm > 0) THEN
560      CALL flinfindcood (fid_out, 'lev', vid, nbdim)
561      IF (nbdim == 1) THEN
562        iret = NF90_GET_VAR (fid, vid, lev, &
563                 start=(/ 1 /), count=(/ llm /))
564      ELSE
565        CALL histerr (3,'flinopen', &
566          'Can not handle vertical coordinates that have more',&
567          'than 1 dimension',' ')
568      ENDIF
569    ENDIF
570  ENDIF
571!-
572! 5.0 Get all the details for the time if possible needed
573!-
574  IF (check) WRITE(*,*) 'flinopen 5.0 Get time'
575!-
576  IF (ttm > 0) THEN
577!---
578!-- 5.1 Find the time axis. Prefered method is the 'timestep since'
579!---
580    gdtmaf_id = -1
581    gdtt_id = -1
582    old_id = -1
583    DO iv=1,ncnbva(fid_out)
584      name=''
585      iret = NF90_INQUIRE_VARIABLE (fid, iv, name=name)
586      units=''
587      iret = NF90_GET_ATT (fid, iv, 'units', units)
588      IF (INDEX(units,'seconds since') > 0) gdtmaf_id = iv
589      IF (INDEX(units,'timesteps since') > 0) gdtt_id = iv
590      IF (INDEX(name, 'tstep') > 0) old_id = iv
591    ENDDO
592!---
593    IF (gdtt_id > 0) THEN
594      vid = gdtt_id
595    ELSE IF (gdtmaf_id > 0) THEN
596      vid = gdtmaf_id
597    ELSE IF (old_id > 0) THEN
598      vid = old_id
599    ELSE
600      CALL histerr (3, 'flinopen', 'No time axis found',' ',' ')
601    ENDIF
602!---
603    ALLOCATE(vec_tmp(ttm))
604    iret = NF90_GET_VAR (fid,vid,vec_tmp,start=(/ 1 /),count=(/ ttm /))
605    itaus(1:ttm) = NINT(vec_tmp(1:ttm))
606    DEALLOCATE(vec_tmp)
607!---
608    IF (check) WRITE(*,*) 'flinopen 5.1 Times ',itaus
609!---
610!-- Getting all the details for the time axis
611!---
612!-- Find the calendar
613    iret = NF90_GET_ATT (fid,gdtmaf_id,'calendar',calendar)
614    IF (iret == NF90_NOERR) THEN
615      CALL ioconf_calendar(calendar)
616    ENDIF
617!--
618    units = ''
619    iret = NF90_GET_ATT (fid,vid,'units',units)
620    IF (gdtt_id > 0) THEN
621      units = units(INDEX(units,'since')+6:LEN_TRIM(units))
622      READ (units,'(I4.4,5(a,I2.2))') &
623        year0, strc, month0, strc, day0, &
624               strc, hours0, strc, minutes0, strc, seci
625      sec0 = hours0*3600. + minutes0*60. + seci
626      CALL ymds2ju (year0, month0, day0, sec0, date0)
627      IF (check) &
628        WRITE(*,*) 'flinopen 5.1 gdtt_id year0 ... date0 ', &
629                   year0, month0, day0, sec0, date0
630!-----
631      iret = NF90_GET_ATT (fid, gdtt_id, 'tstep_sec', dt)
632    ELSE IF (gdtmaf_id > 0) THEN
633      units = units(INDEX(units,'since')+6:LEN_TRIM(units))
634      READ (units,'(I4.4,5(a,I2.2))') &
635        year0, strc, month0, strc, day0, &
636               strc, hours0, strc, minutes0, strc, seci
637      sec0 = hours0*3600. + minutes0*60. + seci
638      CALL ymds2ju (year0, month0, day0, sec0, date0)
639!-----
640      IF (check) &
641        WRITE(*,*) 'flinopen 5.1 gdtmaf_id year0 ... date0 ', &
642                   year0, month0, day0, sec0, date0
643    ELSE IF (old_id > 0) THEN
644      iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'delta_tstep_sec', dt)
645      iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'day0', r_day)
646      iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'sec0', sec)
647      iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'year0', r_year)
648      iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'month0', r_month)
649!-----
650      day = INT(r_day)
651      month = INT(r_month)
652      year = INT(r_year)
653!-----
654      CALL ymds2ju (year, month, day, sec, date0)
655    ENDIF
656  ENDIF
657!-
658  IF (check) WRITE(*,*) 'flinopen 6.0 File opened', date0, dt
659!---------------------------
660END SUBROUTINE flinopen_work
661!-
662!===
663!-
664SUBROUTINE flininfo (filename, iim, jjm, llm, ttm, fid_out)
665!---------------------------------------------------------------------
666!- This subroutine allows to get some information.
667!- It is usualy done within flinopen but the user may want to call
668!- it before in order to allocate the space needed to extract the
669!- data from the file.
670!---------------------------------------------------------------------
671  IMPLICIT NONE
672!-
673! ARGUMENTS
674!-
675  CHARACTER(LEN=*) :: filename
676  INTEGER :: iim, jjm, llm, ttm, fid_out
677!-
678! LOCAL
679!-
680  INTEGER :: iret, fid, ndims, nvars, nb_atts, id_unlim
681  INTEGER :: iv, lll
682  INTEGER :: xid, yid, zid, tid
683  CHARACTER(LEN=80) :: name
684  CHARACTER(LEN=30) :: axname
685!-
686  LOGICAL :: check = .FALSE.
687!---------------------------------------------------------------------
688  lll = LEN_TRIM(filename)
689  IF (filename(lll-2:lll) /= '.nc') THEN
690    name = filename(1:lll)//'.nc'
691  ELSE
692    name = filename(1:lll)
693  ENDIF
694!-
695  iret = NF90_OPEN (name, NF90_NOWRITE, fid)
696  IF (iret /= NF90_NOERR) THEN
697    CALL histerr(3, 'flininfo','Could not open file :',TRIM(name),' ')
698  ENDIF
699!-
700  iret = NF90_INQUIRE (fid, nDimensions=ndims, nVariables=nvars, &
701                      nAttributes=nb_atts, unlimitedDimId=id_unlim)
702!-
703  xid = -1; iim = 0;
704  yid = -1; jjm = 0;
705  zid = -1; llm = 0;
706  tid = -1; ttm = 0;
707!-
708  DO iv=1,ndims
709!---
710    iret = NF90_INQUIRE_DIMENSION (fid, iv, name=axname, len=lll)
711    CALL strlowercase (axname)
712    axname = ADJUSTL(axname)
713!---
714    IF (check) WRITE(*,*) &
715      'flininfo - getting axname',iv,axname,lll
716!---
717    IF      (    (INDEX(axname,'x') == 1) &
718             .OR.(INDEX(axname,'lon') == 1) ) THEN
719      xid = iv; iim = lll;
720    ELSE IF (    (INDEX(axname,'y') == 1) &
721             .OR.(INDEX(axname,'lat') == 1) ) THEN
722      yid = iv; jjm = lll;
723    ELSE IF (    (INDEX(axname,'lev') == 1) &
724             .OR.(INDEX(axname,'plev') == 1) &
725             .OR.(INDEX(axname,'z') == 1) &
726             .OR.(INDEX(axname,'depth') == 1) ) THEN
727      zid = iv; llm = lll;
728    ELSE IF (    (INDEX(axname,'tstep') == 1) &
729             .OR.(INDEX(axname,'time_counter') == 1) ) THEN
730!---- For the time we certainly need to allow for other names
731      tid = iv; ttm = lll;
732    ELSE IF (ndims == 1) THEN
733!---- Nothing was found and ndims=1 then we have a vector of data
734      xid = 1; iim = lll;
735    ENDIF
736!---
737  ENDDO
738!-
739! Keep all this information
740!-
741  nbfiles = nbfiles+1
742!-
743  IF (nbfiles > nbfile_max) THEN
744    CALL histerr (3,'flininfo', &
745      'Too many files. Please increase nbfil_max', &
746      'in program flincom.F90.',' ')
747  ENDIF
748!-
749  ncids(nbfiles) = fid
750  ncnbd(nbfiles) = ndims
751!-
752  ncdims(nbfiles,1:4) = (/ iim, jjm, llm, ttm /)
753!-
754  ncfunli(nbfiles) = id_unlim
755  ncnba(nbfiles)   = nb_atts
756  ncnbva(nbfiles)  = nvars
757  ncfileopen(nbfiles) = .TRUE.
758!-
759  fid_out = nbfiles
760!----------------------
761END SUBROUTINE flininfo
762!-
763!===
764!-
765SUBROUTINE flinput_r1d &
766  (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var)
767!---------------------------------------------------------------------
768  IMPLICIT NONE
769!-
770  INTEGER :: fid_in
771  CHARACTER(LEN=*) :: varname
772  INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid
773  REAL :: var(:)
774!-
775  INTEGER :: fid, ncvarid, ndim, iret
776  LOGICAL :: check = .FALSE.
777!---------------------------------------------------------------------
778  IF (check) WRITE(*,*) &
779     "flinput_r1d : SIZE(var) = ",SIZE(var)
780!-
781  CALL flinput_mat &
782    (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, &
783     fid,ncvarid,ndim)
784!-
785  iret = NF90_PUT_VAR (fid, ncvarid, var, &
786           start=w_sta(1:ndim), count=w_len(1:ndim))
787!-------------------------
788END SUBROUTINE flinput_r1d
789!-
790!===
791!-
792SUBROUTINE flinput_r2d &
793  (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var)
794!---------------------------------------------------------------------
795  IMPLICIT NONE
796!-
797  INTEGER :: fid_in
798  CHARACTER(LEN=*) :: varname
799  INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid
800  REAL :: var(:,:)
801!-
802  INTEGER :: fid, ncvarid, ndim, iret
803  LOGICAL :: check = .FALSE.
804!---------------------------------------------------------------------
805  IF (check) WRITE(*,*) &
806     "flinput_r2d : SIZE(var) = ",SIZE(var)
807!-
808  CALL flinput_mat &
809    (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, &
810     fid,ncvarid,ndim)
811!-
812  iret = NF90_PUT_VAR (fid, ncvarid, var, &
813           start=w_sta(1:ndim), count=w_len(1:ndim))
814!-------------------------
815END SUBROUTINE flinput_r2d
816!-
817!===
818!-
819SUBROUTINE flinput_r3d &
820  (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var)
821!---------------------------------------------------------------------
822  IMPLICIT NONE
823!-
824  INTEGER :: fid_in
825  CHARACTER(LEN=*) :: varname
826  INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid
827  REAL :: var(:,:,:)
828!-
829  INTEGER :: fid, ncvarid, ndim, iret
830  LOGICAL :: check = .FALSE.
831!---------------------------------------------------------------------
832  IF (check) WRITE(*,*) &
833     "flinput_r3d : SIZE(var) = ",SIZE(var)
834!-
835  CALL flinput_mat &
836    (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, &
837     fid,ncvarid,ndim)
838!-
839  iret = NF90_PUT_VAR (fid, ncvarid, var, &
840           start=w_sta(1:ndim), count=w_len(1:ndim))
841!-------------------------
842END SUBROUTINE flinput_r3d
843!-
844!===
845!-
846SUBROUTINE flinput_r4d &
847  (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var)
848!---------------------------------------------------------------------
849  IMPLICIT NONE
850!-
851  INTEGER :: fid_in
852  CHARACTER(LEN=*) :: varname
853  INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid
854  REAL :: var(:,:,:,:)
855!-
856  INTEGER :: fid, ncvarid, ndim, iret
857  LOGICAL :: check = .FALSE.
858!---------------------------------------------------------------------
859  IF (check) WRITE(*,*) &
860     "flinput_r4d : SIZE(var) = ",SIZE(var)
861!-
862  CALL flinput_mat &
863    (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, &
864     fid,ncvarid,ndim)
865!-
866  iret = NF90_PUT_VAR (fid, ncvarid, var, &
867           start=w_sta(1:ndim), count=w_len(1:ndim))
868!-------------------------
869END SUBROUTINE flinput_r4d
870!-
871!===
872!-
873SUBROUTINE flinput_mat &
874  (fid_in,varname,iim,nlonid,jjm,nlatid, &
875                  llm,zdimid,ttm,tdimid,fid,ncvarid,ndim)
876!---------------------------------------------------------------------
877  IMPLICIT NONE
878!-
879  INTEGER :: fid_in
880  CHARACTER(LEN=*) :: varname
881  INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid
882  INTEGER :: fid, ncvarid, ndim
883!-
884! LOCAL
885!-
886  INTEGER :: iret
887!---------------------------------------------------------------------
888  fid = ncids(fid_in)
889!-
890  w_sta(1:4) = (/      1,      1,  1,  1 /)
891  w_len(1:2) = (/    iim,    jjm /)
892  w_dim(1:2) = (/ nlonid, nlatid /)
893!-
894  IF ( (llm > 0).AND.(ttm > 0) ) THEN
895    ndim = 4
896    w_len(3:4) = (/    llm,    ttm /)
897    w_dim(3:4) = (/ zdimid, tdimid /)
898  ELSE IF (llm > 0) THEN
899    ndim = 3
900    w_dim(3) = zdimid
901    w_len(3) = llm
902  ELSE IF (ttm > 0) THEN
903    ndim = 3
904    w_dim(3) = tdimid
905    w_len(3) = ttm
906  ELSE
907    ndim = 2
908  ENDIF
909!-
910  iret = NF90_REDEF   (fid)
911  iret = NF90_DEF_VAR (fid,varname,NF90_FLOAT,w_dim(1:ndim),ncvarid)
912  iret = NF90_PUT_ATT (fid,ncvarid,'short_name',TRIM(varname))
913  iret = NF90_ENDDEF  (fid)
914!--------------------------
915END  SUBROUTINE flinput_mat
916!-
917!===
918!-
919SUBROUTINE flinput_scal &
920  (fid_in, varname, iim, nlonid, jjm, nlatid, &
921                    llm, zdimid, ttm, tdimid, var)
922!---------------------------------------------------------------------
923  IMPLICIT NONE
924!-
925  INTEGER :: fid_in
926  CHARACTER(LEN=*) :: varname
927  INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid
928  REAL :: var
929!-
930! LOCAL
931!-
932  INTEGER :: fid, iret
933!---------------------------------------------------------------------
934  fid = ncids(fid_in)
935!-
936  iret = NF90_REDEF   (fid)
937  iret = NF90_PUT_ATT (fid, NF90_GLOBAL, varname, REAL(var,KIND=4))
938  iret = NF90_ENDDEF  (fid)
939!---------------------------
940END  SUBROUTINE flinput_scal
941!-
942!===
943!-
944SUBROUTINE flinget_r1d &
945  (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var)
946!---------------------------------------------------------------------
947  IMPLICIT NONE
948!-
949  INTEGER :: fid_in
950  CHARACTER(LEN=*) :: varname
951  INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin
952  REAL :: var(:)
953!-
954  INTEGER :: jl, ji
955  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp
956  LOGICAL :: check = .FALSE.
957!---------------------------------------------------------------------
958  IF (.NOT.ALLOCATED(buff_tmp)) THEN
959    IF (check) WRITE(*,*) &
960      "flinget_r1d : allocate buff_tmp for buff_sz = ",SIZE(var)
961    ALLOCATE (buff_tmp(SIZE(var)))
962  ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN
963    IF (check) WRITE(*,*) &
964      "flinget_r1d : re-allocate buff_tmp for buff_sz = ",SIZE(var)
965    DEALLOCATE (buff_tmp)
966    ALLOCATE (buff_tmp(SIZE(var)))
967  ENDIF
968!-
969  CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, &
970                    itau_dep,itau_fin,1,iim,1,jjm,buff_tmp)
971!-
972  jl=0
973  DO ji=1,SIZE(var,1)
974    jl=jl+1
975    var(ji) = buff_tmp(jl)
976  ENDDO
977!-------------------------
978END SUBROUTINE flinget_r1d
979!-
980!===
981!-
982SUBROUTINE flinget_r2d &
983  (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var)
984!---------------------------------------------------------------------
985  IMPLICIT NONE
986!-
987  INTEGER :: fid_in
988  CHARACTER(LEN=*) :: varname
989  INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin
990  REAL :: var(:,:)
991!-
992  INTEGER :: jl, jj, ji
993  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp
994  LOGICAL :: check = .FALSE.
995!---------------------------------------------------------------------
996  IF (.NOT.ALLOCATED(buff_tmp)) THEN
997    IF (check) WRITE(*,*) &
998      "flinget_r2d : allocate buff_tmp for buff_sz = ",SIZE(var)
999    ALLOCATE (buff_tmp(SIZE(var)))
1000  ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN
1001    IF (check) WRITE(*,*) &
1002      "flinget_r2d : re-allocate buff_tmp for buff_sz = ",SIZE(var)
1003    DEALLOCATE (buff_tmp)
1004    ALLOCATE (buff_tmp(SIZE(var)))
1005  ENDIF
1006!-
1007  CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, &
1008                    itau_dep,itau_fin,1,iim,1,jjm,buff_tmp)
1009!-
1010  jl=0
1011  DO jj=1,SIZE(var,2)
1012    DO ji=1,SIZE(var,1)
1013      jl=jl+1
1014      var(ji,jj) = buff_tmp(jl)
1015    ENDDO
1016  ENDDO
1017!-------------------------
1018END SUBROUTINE flinget_r2d
1019!-
1020!===
1021!-
1022SUBROUTINE flinget_r2d_zoom2d &
1023  (fid_in,varname,iim,jjm,llm,ttm, &
1024   itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var)
1025!---------------------------------------------------------------------
1026  IMPLICIT NONE
1027!-
1028  INTEGER :: fid_in
1029  CHARACTER(LEN=*) :: varname
1030  INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen
1031  REAL :: var(:,:)
1032!-
1033  INTEGER :: jl, jj, ji
1034  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp
1035  LOGICAL :: check = .FALSE.
1036!---------------------------------------------------------------------
1037  IF (.NOT.ALLOCATED(buff_tmp)) THEN
1038    IF (check) WRITE(*,*) &
1039      "flinget_r2d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var)
1040    ALLOCATE (buff_tmp(SIZE(var)))
1041  ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN
1042    IF (check) WRITE(*,*) &
1043      "flinget_r2d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var)
1044    DEALLOCATE (buff_tmp)
1045    ALLOCATE (buff_tmp(SIZE(var)))
1046  ENDIF
1047!-
1048  CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, &
1049                    itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp)
1050!-
1051  jl=0
1052  DO jj=1,SIZE(var,2)
1053    DO ji=1,SIZE(var,1)
1054      jl=jl+1
1055      var(ji,jj) = buff_tmp(jl)
1056    ENDDO
1057  ENDDO
1058!--------------------------------
1059END SUBROUTINE flinget_r2d_zoom2d
1060!-
1061!===
1062!-
1063SUBROUTINE flinget_r3d &
1064  (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var)
1065!---------------------------------------------------------------------
1066  IMPLICIT NONE
1067!-
1068  INTEGER :: fid_in
1069  CHARACTER(LEN=*) :: varname
1070  INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin
1071  REAL :: var(:,:,:)
1072!-
1073  INTEGER :: jl, jk, jj, ji
1074  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp
1075  LOGICAL :: check = .FALSE.
1076!---------------------------------------------------------------------
1077  IF (.NOT.ALLOCATED(buff_tmp)) THEN
1078    IF (check) WRITE(*,*) &
1079      "flinget_r3d : allocate buff_tmp for buff_sz = ",SIZE(var)
1080    ALLOCATE (buff_tmp(SIZE(var)))
1081  ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN
1082    IF (check) WRITE(*,*) &
1083      "flinget_r3d : re-allocate buff_tmp for buff_sz = ",SIZE(var)
1084    DEALLOCATE (buff_tmp)
1085    ALLOCATE (buff_tmp(SIZE(var)))
1086  ENDIF
1087!-
1088  CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, &
1089                    itau_dep,itau_fin,1,iim,1,jjm,buff_tmp)
1090!-
1091  jl=0
1092  DO jk=1,SIZE(var,3)
1093    DO jj=1,SIZE(var,2)
1094      DO ji=1,SIZE(var,1)
1095        jl=jl+1
1096        var(ji,jj,jk) = buff_tmp(jl)
1097      ENDDO
1098    ENDDO
1099  ENDDO
1100!-------------------------
1101END SUBROUTINE flinget_r3d
1102!-
1103!===
1104!-
1105SUBROUTINE flinget_r3d_zoom2d &
1106  (fid_in,varname,iim,jjm,llm,ttm, &
1107   itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var)
1108!---------------------------------------------------------------------
1109  IMPLICIT NONE
1110!-
1111  INTEGER :: fid_in
1112  CHARACTER(LEN=*) :: varname
1113  INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen
1114  REAL :: var(:,:,:)
1115!-
1116  INTEGER :: jl, jk, jj, ji
1117  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp
1118  LOGICAL :: check = .FALSE.
1119!---------------------------------------------------------------------
1120  IF (.NOT.ALLOCATED(buff_tmp)) THEN
1121    IF (check) WRITE(*,*) &
1122      "flinget_r3d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var)
1123    ALLOCATE (buff_tmp(SIZE(var)))
1124  ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN
1125    IF (check) WRITE(*,*) &
1126      "flinget_r3d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var)
1127    DEALLOCATE (buff_tmp)
1128    ALLOCATE (buff_tmp(SIZE(var)))
1129  ENDIF
1130!-
1131  CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, &
1132                    itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp)
1133!-
1134  jl=0
1135  DO jk=1,SIZE(var,3)
1136    DO jj=1,SIZE(var,2)
1137      DO ji=1,SIZE(var,1)
1138        jl=jl+1
1139        var(ji,jj,jk) = buff_tmp(jl)
1140      ENDDO
1141    ENDDO
1142  ENDDO
1143!--------------------------------
1144END SUBROUTINE flinget_r3d_zoom2d
1145!-
1146!===
1147!-
1148SUBROUTINE flinget_r4d &
1149  (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var)
1150!---------------------------------------------------------------------
1151  IMPLICIT NONE
1152!-
1153  INTEGER :: fid_in
1154  CHARACTER(LEN=*) :: varname
1155  INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin
1156  REAL :: var(:,:,:,:)
1157!-
1158  INTEGER :: jl, jk, jj, ji, jm
1159  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp
1160  LOGICAL :: check = .FALSE.
1161!---------------------------------------------------------------------
1162  IF (.NOT.ALLOCATED(buff_tmp)) THEN
1163    IF (check) WRITE(*,*) &
1164      "flinget_r4d : allocate buff_tmp for buff_sz = ",SIZE(var)
1165    ALLOCATE (buff_tmp(SIZE(var)))
1166  ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN
1167    IF (check) WRITE(*,*) &
1168      "flinget_r4d : re-allocate buff_tmp for buff_sz = ",SIZE(var)
1169    DEALLOCATE (buff_tmp)
1170    ALLOCATE (buff_tmp(SIZE(var)))
1171  ENDIF
1172!-
1173  CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, &
1174                    itau_dep,itau_fin,1,iim,1,jjm,buff_tmp)
1175!-
1176  jl=0
1177  DO jm=1,SIZE(var,4)
1178    DO jk=1,SIZE(var,3)
1179      DO jj=1,SIZE(var,2)
1180        DO ji=1,SIZE(var,1)
1181          jl=jl+1
1182          var(ji,jj,jk,jm) = buff_tmp(jl)
1183        ENDDO
1184      ENDDO
1185    ENDDO
1186  ENDDO
1187!-------------------------
1188END SUBROUTINE flinget_r4d
1189!-
1190!===
1191!-
1192SUBROUTINE flinget_r4d_zoom2d &
1193  (fid_in,varname,iim,jjm,llm,ttm, &
1194   itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var)
1195!---------------------------------------------------------------------
1196  IMPLICIT NONE
1197!-
1198  INTEGER :: fid_in
1199  CHARACTER(LEN=*) :: varname
1200  INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen
1201  REAL :: var(:,:,:,:)
1202!-
1203  INTEGER :: jl, jk, jj, ji, jm
1204  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp
1205  LOGICAL :: check = .FALSE.
1206!---------------------------------------------------------------------
1207  IF (.NOT.ALLOCATED(buff_tmp)) THEN
1208    IF (check) WRITE(*,*) &
1209      "flinget_r4d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var)
1210    ALLOCATE (buff_tmp(SIZE(var)))
1211  ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN
1212    IF (check) WRITE(*,*) &
1213      "flinget_r4d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var)
1214    DEALLOCATE (buff_tmp)
1215    ALLOCATE (buff_tmp(SIZE(var)))
1216  ENDIF
1217!-
1218  CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, &
1219                    itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp)
1220!-
1221  jl=0
1222  DO jm=1,SIZE(var,4)
1223    DO jk=1,SIZE(var,3)
1224      DO jj=1,SIZE(var,2)
1225        DO ji=1,SIZE(var,1)
1226          jl=jl+1
1227          var(ji,jj,jk,jm) = buff_tmp(jl)
1228        ENDDO
1229      ENDDO
1230    ENDDO
1231  ENDDO
1232!--------------------------------
1233END SUBROUTINE flinget_r4d_zoom2d
1234!-
1235!===
1236!-
1237SUBROUTINE flinget_mat &
1238  (fid_in, varname, iim, jjm, llm, ttm, itau_dep, &
1239   itau_fin, iideb, iilen, jjdeb, jjlen, var)
1240!---------------------------------------------------------------------
1241!- This subroutine will read the variable named varname from
1242!- the file previously opened by flinopen and identified by fid
1243!-
1244!- It is checked that the dimensions of the variable to be read
1245!- correspond to what the user requested when he specified
1246!- iim, jjm and llm. The only exception which is allowed is
1247!- for compressed data where the horizontal grid is not expected
1248!- to be iim x jjm.
1249!-
1250!- If variable is of size zero a global attribute is read.
1251!- This global attribute will be of type real
1252!-
1253!- INPUT
1254!-
1255!- fid      : File ID returned by flinopen
1256!- varname  : Name of the variable to be read from the file
1257!- iim      : | These three variables give the size of the variables
1258!- jjm      : | to be read. It will be verified that the variables
1259!- llm      : | fits in there.
1260!- ttm      : |
1261!- itau_dep : Time step at which we will start to read
1262!- itau_fin : Time step until which we are going to read
1263!-            For the moment this is done on indexes
1264!-            but it should be in the physical space.
1265!-            If there is no time-axis in the file then use a
1266!-            itau_fin < itau_dep, this will tell flinget not to
1267!-            expect a time-axis in the file.
1268!- iideb    : index i for zoom
1269!- iilen    : length of zoom
1270!- jjdeb    : index j for zoom
1271!- jjlen    : length of zoom
1272!-
1273!- OUTPUT
1274!-
1275!- var      : array that will contain the data
1276!---------------------------------------------------------------------
1277  IMPLICIT NONE
1278!-
1279! ARGUMENTS
1280!-
1281  INTEGER :: fid_in
1282  CHARACTER(LEN=*) :: varname
1283  INTEGER :: iim, jjm, llm, ttm
1284  INTEGER :: itau_dep, itau_fin, iideb, iilen, jjdeb, jjlen
1285  REAL :: var(:)
1286!-
1287! LOCAL
1288!-
1289  INTEGER :: iret, fid
1290  INTEGER :: vid, cvid, clen
1291  CHARACTER(LEN=70) :: str1
1292  CHARACTER(LEN=250) :: att_n, tmp_n
1293  CHARACTER(LEN=5) :: axs_l
1294  INTEGER :: tmp_i
1295  REAL,SAVE :: mis_v=0.
1296  REAL :: tmp_r
1297  INTEGER :: ndims, x_typ, nb_atts
1298  INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: dimids
1299  INTEGER :: i, nvars, i2d, cnd
1300  REAL,DIMENSION(:),ALLOCATABLE,SAVE :: var_tmp
1301  LOGICAL :: uncompress = .FALSE.
1302  LOGICAL :: check = .FALSE.
1303!---------------------------------------------------------------------
1304  fid = ncids(fid_in)
1305!-
1306  IF (check) THEN
1307    WRITE(*,*) &
1308    'flinget_mat : fid_in, fid, varname :', fid_in, fid, TRIM(varname)
1309    WRITE(*,*) &
1310    'flinget_mat : iim, jjm, llm, ttm, itau_dep, itau_fin :', &
1311    iim, jjm, llm, ttm, itau_dep, itau_fin
1312    WRITE(*,*) &
1313    'flinget_mat : iideb, iilen, jjdeb, jjlen :', &
1314    iideb, iilen, jjdeb, jjlen
1315  ENDIF
1316!-
1317  uncompress = .FALSE.
1318!-
1319! 1.0 We get first all the details on this variable from the file
1320!-
1321  nvars = ncnbva(fid_in)
1322!-
1323  vid = -1
1324  iret = NF90_INQ_VARID (fid, varname, vid)
1325!-
1326  IF (vid < 0 .OR. iret /= NF90_NOERR) THEN
1327    CALL histerr (3,'flinget', &
1328      'Variable '//TRIM(varname)//' not found in file',' ',' ')
1329  ENDIF
1330!-
1331  iret = NF90_INQUIRE_VARIABLE (fid, vid, &
1332           ndims=ndims, dimids=dimids, nAtts=nb_atts)
1333  IF (check) THEN
1334    WRITE(*,*) &
1335    'flinget_mat : fid, vid :', fid, vid
1336    WRITE(*,*) &
1337    'flinget_mat : ndims, dimids(1:ndims), nb_atts :', &
1338    ndims, dimids(1:ndims), nb_atts
1339  ENDIF
1340!-
1341  w_dim(:) = 0
1342  DO i=1,ndims
1343    iret  = NF90_INQUIRE_DIMENSION (fid, dimids(i), len=w_dim(i))
1344  ENDDO
1345  IF (check) WRITE(*,*) &
1346    'flinget_mat : w_dim :', w_dim(1:ndims)
1347!-
1348  mis_v = 0.0; axs_l = ' ';
1349!-
1350  IF (nb_atts > 0) THEN
1351    IF (check) THEN
1352      WRITE(*,*) 'flinget_mat : attributes for variable :'
1353    ENDIF
1354  ENDIF
1355  DO i=1,nb_atts
1356    iret = NF90_INQ_ATTNAME (fid, vid, i, att_n)
1357    iret = NF90_INQUIRE_ATTRIBUTE (fid, vid, att_n, xtype=x_typ)
1358    CALL strlowercase (att_n)
1359    IF      (    (x_typ == NF90_INT).OR.(x_typ == NF90_SHORT) &
1360             .OR.(x_typ == NF90_BYTE) ) THEN
1361      iret = NF90_GET_ATT (fid, vid, att_n, tmp_i)
1362      IF (check) THEN
1363        WRITE(*,*) '   ',TRIM(att_n),' : ',tmp_i
1364      ENDIF
1365    ELSE IF ( (x_typ == NF90_FLOAT).OR.(x_typ == NF90_DOUBLE) ) THEN
1366      iret = NF90_GET_ATT (fid, vid, att_n, tmp_r)
1367      IF (check) THEN
1368        WRITE(*,*) '   ',TRIM(att_n),' : ',tmp_r
1369      ENDIF
1370      IF (index(att_n,'missing_value') > 0) THEN
1371        mis_v = tmp_r
1372      ENDIF
1373    ELSE
1374      tmp_n = ''
1375      iret = NF90_GET_ATT (fid, vid, att_n, tmp_n)
1376      IF (check) THEN
1377        WRITE(*,*) '   ',TRIM(att_n),' : ',TRIM(tmp_n)
1378      ENDIF
1379      IF (index(att_n,'axis') > 0) THEN
1380        axs_l = tmp_n
1381      ENDIF
1382    ENDIF
1383  ENDDO
1384!?
1385!!!!!!!!!! We will need a verification on the type of the variable
1386!?
1387!-
1388! 2.0 The dimensions are analysed to determine what is to be read
1389!-
1390! 2.1 the longitudes
1391!-
1392  IF ( w_dim(1) /= iim .OR. w_dim(2) /= jjm) THEN
1393!---
1394!-- There is a possibility that we have to deal with a compressed axis !
1395!---
1396    iret = NF90_INQUIRE_DIMENSION (fid, dimids(1), &
1397             name=tmp_n, len=clen)
1398    iret = NF90_INQ_VARID (fid, tmp_n, cvid)
1399!---
1400    IF (check) WRITE(*,*) &
1401      'Dimname, iret , NF90_NOERR : ',TRIM(tmp_n),iret,NF90_NOERR
1402!---
1403!-- If we have an axis which has the same name
1404!-- as the dimension we can see if it is compressed
1405!---
1406!-- TODO TODO for zoom2d
1407!---
1408    IF (iret == NF90_NOERR) THEN
1409      iret = NF90_GET_ATT (fid, cvid, 'compress', str1)
1410!-----
1411      IF (iret == NF90_NOERR) THEN
1412        iret = NF90_INQUIRE_VARIABLE (fid,cvid,xtype=x_typ,ndims=cnd)
1413!-------
1414        IF ( cnd /= 1 .AND. x_typ /= NF90_INT) THEN
1415          CALL histerr (3,'flinget', &
1416            'Variable '//TRIM(tmp_n)//' can not be a compressed axis', &
1417            'Either it has too many dimensions'// &
1418            ' or it is not of type integer', ' ')
1419        ELSE
1420!---------
1421!-------- Let us see if we already have that index table
1422!---------
1423          IF (    (cind_len /= clen).OR.(cind_vid /= cvid) &
1424              .OR.(cind_fid /= fid) ) THEN
1425            IF (ALLOCATED(cindex))   DEALLOCATE(cindex)
1426            ALLOCATE(cindex(clen))
1427            cind_len = clen
1428            cind_vid = cvid
1429            cind_fid = fid
1430            iret = NF90_GET_VAR (fid, cvid, cindex)
1431          ENDIF
1432!---------
1433!-------- In any case we need to set the slab of data to be read
1434!---------
1435          uncompress = .TRUE.
1436          w_sta(1) = 1
1437          w_len(1) = clen
1438          i2d = 1
1439        ENDIF
1440      ELSE
1441        str1 = 'The horizontal dimensions of '//varname
1442        CALL histerr (3,'flinget',str1, &
1443          'is not compressed and does not'// &
1444          ' correspond to the requested size',' ')
1445      ENDIF
1446    ELSE
1447      IF (w_dim(1) /= iim) THEN
1448        str1 = 'The longitude dimension of '//varname
1449        CALL histerr (3,'flinget',str1, &
1450          'in the file is not equal to the dimension', &
1451          'that should be read')
1452      ENDIF
1453      IF (w_dim(2) /= jjm) THEN
1454        str1 = 'The latitude dimension of '//varname
1455        CALL histerr (3,'flinget',str1, &
1456          'in the file is not equal to the dimension', &
1457          'that should be read')
1458      ENDIF
1459    ENDIF
1460  ELSE
1461    w_sta(1:2) = (/ iideb, jjdeb /)
1462    w_len(1:2) = (/ iilen, jjlen /)
1463    i2d = 2
1464  ENDIF
1465!-
1466! 2.3 Now the difficult part, the 3rd dimension which can be
1467! time or levels.
1468!-
1469! Priority is given to the time axis if only three axes are present.
1470!-
1471  IF (ndims > i2d) THEN
1472!---
1473!-- 2.3.1 We have a vertical axis
1474!---
1475    IF (llm == 1 .AND. ndims == i2d+2 .OR. llm == w_dim(i2d+1)) THEN
1476!-----
1477      IF (w_dim(i2d+1) /= llm) THEN
1478        CALL histerr (3,'flinget', &
1479          'The vertical dimension of '//varname, &
1480          'in the file is not equal to the dimension', &
1481          'that should be read')
1482      ELSE
1483        w_sta(i2d+1) = 1
1484        IF (llm > 0) THEN
1485          w_len(i2d+1) = llm
1486        ELSE
1487          w_len(i2d+1) = w_sta(i2d+1)
1488        ENDIF
1489      ENDIF
1490!-----
1491      IF ((itau_fin-itau_dep) >= 0) THEN
1492        IF      (ndims /= i2d+2) THEN
1493          CALL histerr (3,'flinget', &
1494            'You attempt to read a time slab', &
1495            'but there is no time axis on this variable', varname)
1496        ELSE IF ((itau_fin - itau_dep) <= w_dim(i2d+2)) THEN
1497          w_sta(i2d+2) = itau_dep
1498          w_len(i2d+2) = itau_fin-itau_dep+1
1499        ELSE
1500          CALL histerr (3,'flinget', &
1501            'The time step you try to read is not', &
1502            'in the file (1)', varname)
1503        ENDIF
1504      ELSE IF (ndims == i2d+2 .AND. w_dim(i2d+2) > 1) THEN
1505        CALL histerr (3,'flinget', &
1506          'There is a time axis in the file but no', &
1507          'time step give in the call', varname)
1508      ELSE
1509        w_sta(i2d+2) = 1
1510        w_len(i2d+2) = 1
1511      ENDIF
1512    ELSE
1513!-----
1514!---- 2.3.2 We do not have any vertical axis
1515!-----
1516      IF (ndims == i2d+2) THEN
1517        CALL histerr (3,'flinget', &
1518          'The file contains 4 dimensions', &
1519          'but only 3 are requestes for variable ', varname)
1520      ENDIF
1521      IF ((itau_fin-itau_dep) >= 0) THEN
1522        IF (ndims == i2d+1) THEN
1523          IF ((itau_fin-itau_dep) < w_dim(i2d+1) ) THEN
1524            w_sta(i2d+1) = itau_dep
1525            w_len(i2d+1) = itau_fin-itau_dep+1
1526          ELSE
1527            CALL histerr (3,'flinget', &
1528              'The time step you try to read is not', &
1529              'in the file (2)', varname)
1530          ENDIF
1531        ELSE
1532          CALL histerr (3,'flinget', &
1533            'From your input you sould have 3 dimensions', &
1534            'in the file but there are 4', varname)
1535        ENDIF
1536      ELSE
1537        IF (ndims == i2d+1 .AND. w_dim(i2d+1) > 1) THEN
1538          CALL histerr (3,'flinget', &
1539            'There is a time axis in the file but no', &
1540            'time step given in the call', varname)
1541        ELSE
1542          w_sta(i2d+1) = 1
1543          w_len(i2d+1) = 1
1544        ENDIF
1545      ENDIF
1546    ENDIF
1547  ELSE
1548!---
1549!-- 2.3.3 We do not have any vertical axis
1550!---
1551    w_sta(i2d+1:i2d+2) = (/ 0, 0 /)
1552    w_len(i2d+1:i2d+2) = (/ 0, 0 /)
1553  ENDIF
1554!-
1555! 3.0 Reading the data
1556!-
1557  IF (check) WRITE(*,*) &
1558    'flinget_mat 3.0 : ', uncompress, w_sta, w_len
1559!---
1560  IF (uncompress) THEN
1561!---
1562    IF (ALLOCATED(var_tmp)) THEN
1563      IF (SIZE(var_tmp) < clen) THEN
1564        DEALLOCATE(var_tmp)
1565        ALLOCATE(var_tmp(clen))
1566      ENDIF
1567    ELSE
1568      ALLOCATE(var_tmp(clen))
1569    ENDIF
1570!---
1571    iret = NF90_GET_VAR (fid, vid, var_tmp, &
1572             start=w_sta(:), count=w_len(:))
1573!---
1574    var(:) = mis_v
1575    var(cindex(:)) = var_tmp(:)
1576!---
1577  ELSE
1578    iret = NF90_GET_VAR (fid, vid, var, &
1579             start=w_sta(:), count=w_len(:))
1580  ENDIF
1581!-
1582  IF (check) WRITE(*,*) 'flinget_mat 3.1 : ',NF90_STRERROR (iret)
1583!--------------------------
1584END  SUBROUTINE flinget_mat
1585!-
1586!===
1587!-
1588SUBROUTINE flinget_scal &
1589  (fid_in, varname, iim, jjm, llm, ttm, itau_dep, itau_fin, var)
1590!---------------------------------------------------------------------
1591!- This subroutine will read the variable named varname from
1592!- the file previously opened by flinopen and identified by fid
1593!-
1594!- If variable is of size zero a global attribute is read. This
1595!- global attribute will be of type real
1596!-
1597!- INPUT
1598!-
1599!- fid      : File ID returned by flinopen
1600!- varname  : Name of the variable to be read from the file
1601!- iim      : | These three variables give the size of the variables
1602!- jjm      : | to be read. It will be verified that the variables
1603!- llm      : | fits in there.
1604!- ttm      : |
1605!- itau_dep : Time step at which we will start to read
1606!- itau_fin : Time step until which we are going to read
1607!-           For the moment this is done on indeces but it should be
1608!-           in the physical space
1609!-           If there is no time-axis in the file then use a
1610!-           itau_fin < itau_dep, this will tell flinget not to
1611!-           expect a time-axis in the file.
1612!-
1613!- OUTPUT
1614!-
1615!- var      : scalar that will contain the data
1616!---------------------------------------------------------------------
1617  IMPLICIT NONE
1618!-
1619! ARGUMENTS
1620!-
1621  INTEGER :: fid_in
1622  CHARACTER(LEN=*) :: varname
1623  INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin
1624  REAL :: var
1625!-
1626! LOCAL
1627!-
1628  INTEGER :: iret, fid
1629!-
1630  LOGICAL :: check = .FALSE.
1631!---------------------------------------------------------------------
1632  IF (check) THEN
1633    WRITE (*,*) 'flinget_scal in file with id ',fid_in
1634  ENDIF
1635!-
1636  fid = ncids(fid_in)
1637!-
1638! 1.0 Reading a global attribute
1639!-
1640  iret = NF90_GET_ATT (fid, NF90_GLOBAL, varname, var)
1641!---------------------------
1642END  SUBROUTINE flinget_scal
1643!-
1644!===
1645!-
1646SUBROUTINE flinfindcood (fid_in, axtype, vid, ndim)
1647!---------------------------------------------------------------------
1648!- This subroutine explores the file in order to find
1649!- the coordinate according to a number of rules
1650!---------------------------------------------------------------------
1651  IMPLICIT NONE
1652!-
1653! ARGUMENTS
1654!-
1655  INTEGER :: fid_in, vid, ndim
1656  CHARACTER(LEN=3) :: axtype
1657!-
1658! LOCAL
1659!-
1660  INTEGER :: iv, iret, dimnb
1661  CHARACTER(LEN=40) :: dimname, dimuni1, dimuni2, dimuni3
1662  CHARACTER(LEN=30) :: str1
1663  LOGICAL :: found_rule = .FALSE.
1664!---------------------------------------------------------------------
1665  vid = -1
1666!-
1667! Make sure all strings are invalid
1668!-
1669  dimname = '?-?'
1670  dimuni1 = '?-?'
1671  dimuni2 = '?-?'
1672  dimuni3 = '?-?'
1673!-
1674! First rule : we look for the correct units
1675! lon : east
1676! lat : north
1677! We make an exact check as it would be too easy to mistake
1678! some units by just comparing the substrings.
1679!-
1680  SELECTCASE(axtype)
1681  CASE ('lon')
1682    dimuni1 = 'degree_e'
1683    dimuni2 = 'degrees_e'
1684    found_rule = .TRUE.
1685  CASE('lat')
1686    dimuni1 = 'degree_n'
1687    dimuni2 = 'degrees_n'
1688    found_rule = .TRUE.
1689  CASE('lev')
1690    dimuni1 = 'm'
1691    dimuni2 = 'km'
1692    dimuni3 = 'hpa'
1693    found_rule = .TRUE.
1694  CASE DEFAULT
1695    found_rule = .FALSE.
1696  END SELECT
1697!-
1698  IF (found_rule) THEN
1699    iv = 0
1700    DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) )
1701      iv = iv+1
1702      str1 = ''
1703      iret = NF90_GET_ATT (ncids(fid_in), iv, 'units', str1)
1704      IF (iret == NF90_NOERR) THEN
1705        CALL strlowercase (str1)
1706        IF (    (INDEX(str1, TRIM(dimuni1)) == 1) &
1707            .OR.(INDEX(str1, TRIM(dimuni2)) == 1) &
1708            .OR.(INDEX(str1, TRIM(dimuni3)) == 1) ) THEN
1709          vid = iv
1710          iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, ndims=ndim)
1711        ENDIF
1712      ENDIF
1713    ENDDO
1714  ENDIF
1715!-
1716! Second rule : we find specific names :
1717! lon : nav_lon
1718! lat : nav_lat
1719! Here we can check if we find the substring as the
1720! names are more specific.
1721!-
1722  SELECTCASE(axtype)
1723  CASE ('lon')
1724    dimname = 'nav_lon lon longitude'
1725    found_rule = .TRUE.
1726  CASE('lat')
1727    dimname = 'nav_lat lat latitude'
1728    found_rule = .TRUE.
1729  CASE('lev')
1730    dimname = 'plev level depth deptht'
1731    found_rule = .TRUE.
1732  CASE DEFAULT
1733    found_rule = .FALSE.
1734  END SELECT
1735!-
1736  IF (found_rule) THEN
1737    iv = 0
1738    DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) )
1739      iv = iv+1
1740      str1=''
1741      iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, &
1742                                    name=str1, ndims=ndim)
1743      IF (INDEX(dimname,TRIM(str1)) >= 1) THEN
1744        vid = iv
1745      ENDIF
1746    ENDDO
1747  ENDIF
1748!-
1749! Third rule : we find a variable with the same name as the dimension
1750! lon = 1
1751! lat = 2
1752! lev = 3
1753!-
1754  IF (vid < 0) THEN
1755    SELECTCASE(axtype)
1756    CASE ('lon')
1757      dimnb = 1
1758      found_rule = .TRUE.
1759    CASE('lat')
1760      dimnb = 2
1761      found_rule = .TRUE.
1762    CASE('lev')
1763      dimnb = 3
1764      found_rule = .TRUE.
1765    CASE DEFAULT
1766      found_rule = .FALSE.
1767    END SELECT
1768!---
1769    IF (found_rule) THEN
1770      iret = NF90_INQUIRE_DIMENSION (ncids(fid_in), dimnb, name=dimname)
1771      iv = 0
1772      DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) )
1773        iv = iv+1
1774        str1=''
1775        iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, &
1776                                      name=str1, ndims=ndim)
1777        IF (INDEX(dimname,TRIM(str1)) == 1) THEN
1778          vid = iv
1779        ENDIF
1780      ENDDO
1781    ENDIF
1782  ENDIF
1783!-
1784! Stop the program if no coordinate was found
1785!-
1786  IF (vid < 0) THEN
1787    CALL histerr (3,'flinfindcood', &
1788           'No coordinate axis was found in the file', &
1789           'The data in this file can not be used', axtype)
1790  ENDIF
1791!--------------------------
1792END SUBROUTINE flinfindcood
1793!-
1794!===
1795!-
1796SUBROUTINE flinclo (fid_in)
1797!---------------------------------------------------------------------
1798  IMPLICIT NONE
1799!-
1800  INTEGER :: fid_in
1801!-
1802  INTEGER :: iret
1803!---------------------------------------------------------------------
1804  iret = NF90_CLOSE (ncids(fid_in))
1805  ncfileopen(fid_in) = .FALSE.
1806!---------------------
1807END SUBROUTINE flinclo
1808!-
1809!===
1810!-
1811SUBROUTINE flinquery_var(fid_in, varname, exists)
1812!---------------------------------------------------------------------
1813!- Queries the existance of a variable in the file.
1814!---------------------------------------------------------------------
1815  IMPLICIT NONE
1816!-
1817  INTEGER :: fid_in
1818  CHARACTER(LEN=*) varname
1819  LOGICAL :: exists
1820!-
1821  INTEGER :: iret, fid, vid
1822!---------------------------------------------------------------------
1823  fid = ncids(fid_in)
1824  vid = -1
1825  iret = NF90_INQ_VARID (fid, varname, vid)
1826!-
1827  exists = ( (vid >= 0).AND.(iret == NF90_NOERR) )
1828!---------------------------
1829END SUBROUTINE flinquery_var
1830!-
1831!===
1832!-
1833SUBROUTINE flininspect (fid_in)
1834!---------------------------------------------------------------------
1835  IMPLICIT NONE
1836!-
1837! fid : File id to inspect
1838!-
1839  INTEGER :: fid_in
1840!-
1841!- LOCAL
1842!-
1843  INTEGER :: iim, jjm, llm, ttm, fid_out
1844  INTEGER :: iret, fid, ndims, nvars, nb_atts, id_unlim
1845  INTEGER :: iv, in, lll
1846  INTEGER :: xid, yid, zid, tid
1847  INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: idimid
1848  CHARACTER(LEN=80) :: name
1849  CHARACTER(LEN=30) :: axname
1850!---------------------------------------------------------------------
1851  fid = ncids(fid_in)
1852!-
1853  iret = NF90_INQUIRE (fid, nDimensions=ndims, nVariables=nvars, &
1854                       nAttributes=nb_atts, unlimitedDimId=id_unlim)
1855!-
1856  WRITE (*,*) 'IOIPSL ID                   : ',fid_in
1857  WRITE (*,*) 'NetCDF ID                   : ',fid
1858  WRITE (*,*) 'Number of dimensions        : ',ndims
1859  WRITE (*,*) 'Number of variables         : ',nvars
1860  WRITE (*,*) 'Number of global attributes : ',nb_atts
1861  WRITE (*,*) 'ID unlimited                : ',id_unlim
1862!-
1863  xid = -1; iim = 0;
1864  yid = -1; jjm = 0;
1865  zid = -1; llm = 0;
1866  tid = -1; ttm = 0;
1867!-
1868  DO iv=1,ndims
1869!---
1870    iret = NF90_INQUIRE_DIMENSION (fid, iv, name=axname, len=lll)
1871    CALL strlowercase (axname)
1872    axname = ADJUSTL(axname)
1873!---
1874    WRITE (*,*) 'Dimension number : ',iv
1875    WRITE (*,*) 'Dimension name   : ',TRIM(axname)
1876!---
1877    IF      (    (INDEX(axname,'x') == 1) &
1878             .OR.(INDEX(axname,'lon') == 1)) THEN
1879      xid = iv; iim = lll;
1880      WRITE (*,*) 'Dimension X size   : ',iim
1881    ELSE IF (    (INDEX(axname,'y') == 1) &
1882             .OR.(INDEX(axname,'lat') == 1)) THEN
1883      yid = iv; jjm = lll;
1884      WRITE (*,*) 'Dimension Y size   : ',jjm
1885    ELSE IF (    (INDEX(axname,'lev') == 1) &
1886             .OR.(INDEX(axname,'plev') == 1) &
1887             .OR.(INDEX(axname,'z') == 1) &
1888             .OR.(INDEX(axname,'depth') == 1)) THEN
1889      zid = iv; llm = lll;
1890      WRITE (*,*) 'Dimension Z size   : ',llm
1891    ELSE IF (    (INDEX(axname,'tstep') == 1) &
1892             .OR.(INDEX(axname,'time_counter') == 1)) THEN
1893!---- For the time we certainly need to allow for other names
1894      tid = iv; ttm = lll;
1895    ELSE IF (ndims == 1) THEN
1896!---- Nothing was found and ndims=1 then we have a vector of data
1897      xid = 1; iim = lll;
1898    ENDIF
1899!---
1900  ENDDO
1901!-
1902! Keep all this information
1903!-
1904  nbfiles = nbfiles+1
1905!-
1906  IF (nbfiles > nbfile_max) THEN
1907    CALL histerr(3,'flininspect', &
1908      'Too many files. Please increase nbfil_max', &
1909      'in program flincom.F90.',' ')
1910  ENDIF
1911!-
1912  ncids(nbfiles) = fid
1913  ncnbd(nbfiles) = ndims
1914!-
1915  ncdims(nbfiles,1:4) = (/ iim, jjm, llm, ttm /)
1916!-
1917  ncfunli(nbfiles) = id_unlim
1918  ncnba(nbfiles)   = nb_atts
1919  ncnbva(nbfiles)  = nvars
1920  ncfileopen(nbfiles) = .TRUE.
1921!-
1922  fid_out = nbfiles
1923!-
1924  DO in=1,nvars
1925    iret = NF90_INQUIRE_VARIABLE (fid, in, &
1926             name=name, ndims=ndims, dimids=idimid, nAtts=nb_atts)
1927    WRITE (*,*) 'Variable number  ------------ > ', in
1928    WRITE (*,*) 'Variable name        : ', TRIM(name)
1929    WRITE (*,*) 'Number of dimensions : ', ndims
1930    WRITE (*,*) 'Dimensions ID''s     : ', idimid(1:ndims)
1931    WRITE (*,*) 'Number of attributes : ', nb_atts
1932  ENDDO
1933!-------------------------
1934END SUBROUTINE flininspect
1935!-
1936!===
1937!-
1938END MODULE flincom
Note: See TracBrowser for help on using the repository browser.