source: IOIPSL/trunk/src/histcom.f90 @ 75

Last change on this file since 75 was 75, checked in by bellier, 15 years ago

JB: new version of histcom (min-max values) and test

  • Property svn:keywords set to Id
File size: 86.3 KB
Line 
1MODULE histcom
2!-
3!$Id$
4!-
5  USE netcdf
6!-
7  USE stringop,  ONLY : nocomma,cmpblank,findpos,find_str,strlowercase
8  USE mathelp,   ONLY : mathop,moycum,trans_buff,buildop
9  USE fliocom,   ONLY : flio_dom_file,flio_dom_att
10  USE calendar
11  USE errioipsl, ONLY : ipslerr
12!-
13  IMPLICIT NONE
14!-
15  PRIVATE
16  PUBLIC :: histbeg, histdef, histhori, histvert, histend, &
17 &          histwrite, histclo, histsync, ioconf_modname
18!---------------------------------------------------------------------
19!- Some confusing vocabulary in this code !
20!- =========================================
21!-
22!- A REGULAR grid is a grid which is i,j indexes
23!- and thus it is stored in a 2D matrix.
24!- This is opposed to a IRREGULAR grid which is only in a vector
25!- and where we do not know which neighbors we have.
26!- As a consequence we need the bounds for each grid-cell.
27!-
28!- A RECTILINEAR grid is a special case of a regular grid
29!- in which all longitudes for i constant are equal
30!- and all latitudes for j constant.
31!- In other words we do not need the full 2D matrix
32!- to describe the grid, just two vectors.
33!---------------------------------------------------------------------
34  INTERFACE histwrite
35!---------------------------------------------------------------------
36!- The "histwrite" routines will give the data to the I/O system.
37!- It will trigger the operations to be performed,
38!- and the writting to the file if needed
39!-
40!- We test for the work to be done at this time here so that at a
41!- later stage we can call different operation and write subroutine
42!- for the REAL and INTEGER interfaces
43!-
44!- INPUT
45!- pfileid  : The ID of the file on which this variable is to be,
46!-            written. The variable should have been defined in
47!-            this file before.
48!- pvarname : The short name of the variable
49!- pitau    : Current timestep
50!- pdata    : The variable, I mean the real data !
51!- nbindex  : The number of indexes provided. If it is equal to
52!-            the size of the full field as provided in histdef
53!-            then nothing is done.
54!- nindex   : The indices used to expand the variable (pdata)
55!-            onto the full field.
56!---------------------------------------------------------------------
57!- histwrite - we have to prepare different type of fields :
58!-             real and integer, 1,2 or 3D
59    MODULE PROCEDURE histwrite_r1d,histwrite_r2d,histwrite_r3d
60  END INTERFACE
61!-
62  INTERFACE histbeg
63!!  MODULE PROCEDURE histbeg_regular,histbeg_irregular
64    MODULE PROCEDURE histbeg_totreg,histbeg_regular,histbeg_irregular
65  END INTERFACE
66!-
67  INTERFACE histhori
68    MODULE PROCEDURE histhori_regular,histhori_irregular
69  END INTERFACE
70!-
71! Fixed parameter
72!-
73  INTEGER,PARAMETER :: nb_files_max=20, nb_var_max=400, &
74 &                     nb_hax_max=5, nb_zax_max=10, nbopp_max=10
75  REAL,PARAMETER :: missing_val=nf90_fill_real
76!-
77  INTEGER :: bufftmp_max(nb_files_max) = 1
78!-
79! Time variables
80!-
81  INTEGER,SAVE :: itau0(nb_files_max)=0
82  REAL,DIMENSION(nb_files_max),SAVE ::date0,deltat
83!-
84! Counter of elements
85!-
86  INTEGER,SAVE :: nb_files=0
87  INTEGER,DIMENSION(nb_files_max),SAVE :: nb_var=0, nb_tax=0
88!-
89! NETCDF IDs for files and axes
90!-
91  INTEGER,DIMENSION(nb_files_max),SAVE :: ncdf_ids,xid,yid,tid
92  CHARACTER(LEN=500),SAVE :: assc_file=''
93!-
94! General definitions in the NETCDF file
95!-
96  INTEGER,DIMENSION(nb_files_max,2),SAVE :: &
97 &   full_size=0,slab_ori,slab_sz
98!-
99! The horizontal axes
100!-
101  INTEGER,SAVE :: nb_hax(nb_files_max)=0
102  CHARACTER(LEN=25),SAVE :: hax_name(nb_files_max,nb_hax_max,2)
103!-
104! The vertical axes
105!-
106  INTEGER,SAVE :: nb_zax(nb_files_max)=0
107  INTEGER,DIMENSION(nb_files_max,nb_zax_max),SAVE :: &
108  &  zax_size,zax_ids,zax_name_length
109  CHARACTER(LEN=20),SAVE :: zax_name(nb_files_max,nb_zax_max)
110!-
111! Informations on each variable
112!-
113  INTEGER,DIMENSION(nb_files_max,nb_var_max),SAVE :: &
114 &  name_length,nbopp
115  CHARACTER(LEN=20),DIMENSION(nb_files_max,nb_var_max),SAVE :: &
116 &  name,unit_name
117  CHARACTER(LEN=80),DIMENSION(nb_files_max,nb_var_max),SAVE :: &
118 &  title,fullop
119  CHARACTER(LEN=7),SAVE :: topp(nb_files_max,nb_var_max)
120  CHARACTER(LEN=7),SAVE :: sopps(nb_files_max,nb_var_max,nbopp_max)
121  REAL,SAVE :: scal(nb_files_max,nb_var_max,nbopp_max)
122!- Sizes of the associated grid and zommed area
123  INTEGER,DIMENSION(nb_files_max,nb_var_max,3),SAVE :: &
124 &   scsize,zorig,zsize
125!- Sizes for the data as it goes through the various math operations
126  INTEGER,SAVE :: datasz_in(nb_files_max,nb_var_max,3) = -1
127  INTEGER,SAVE :: datasz_max(nb_files_max,nb_var_max) = -1
128!-
129  INTEGER,DIMENSION(nb_files_max,nb_var_max),SAVE :: &
130 &  var_haxid,var_zaxid,var_axid,ncvar_ids
131!-
132  REAL,DIMENSION(nb_files_max,nb_var_max,2),SAVE :: hist_minmax
133  LOGICAL,DIMENSION(nb_files_max,nb_var_max),SAVE :: &
134 &  hist_calc_rng=.FALSE.,hist_wrt_rng=.FALSE.
135!-
136  REAL,DIMENSION(nb_files_max,nb_var_max),SAVE :: &
137 &  freq_opp,freq_wrt
138  INTEGER,DIMENSION(nb_files_max,nb_var_max),SAVE :: &
139 &  last_opp,last_wrt,last_opp_chk,last_wrt_chk,nb_opp,nb_wrt,point
140!-
141! Book keeping for the buffers
142!-
143  INTEGER,SAVE :: buff_pos=0
144  REAL,ALLOCATABLE,SAVE :: buffer(:)
145  LOGICAL,SAVE :: &
146 &  zoom(nb_files_max)=.FALSE., regular(nb_files_max)=.TRUE.
147!-
148! Book keeping of the axes
149!-
150  INTEGER,DIMENSION(nb_files_max,nb_var_max),SAVE :: &
151 &  tdimid,tax_last,tax_name_length
152  CHARACTER(LEN=40),DIMENSION(nb_files_max,nb_var_max),SAVE :: &
153 &  tax_name
154!-
155! A list of functions which require special action
156! (Needs to be updated when functions are added
157!  but they are well located here)
158!-
159  CHARACTER(LEN=120),SAVE :: &
160 &  indchfun = 'scatter, fill, gather, coll', &
161 &  fuchnbout = 'scatter, fill'
162!- Some configurable variables with locks
163  CHARACTER(LEN=80),SAVE :: model_name='An IPSL model'
164  LOGICAL,SAVE :: lock_modname=.FALSE.
165!-
166!===
167CONTAINS
168!===
169!-
170SUBROUTINE histbeg_totreg                     &
171 & (pfilename, pim, plon, pjm, plat,          &
172 &  par_orix, par_szx, par_oriy, par_szy,     &
173 &  pitau0, pdate0, pdeltat, phoriid, pfileid, domain_id)
174!---------------------------------------------------------------------
175!- This is just an interface for histbeg_regular in case when
176!- the user provides plon and plat as vectors. Obviously this can only
177!- be used for very regular grids.
178!-
179!- INPUT
180!-
181!- pfilename : Name of the netcdf file to be created
182!- pim       : Size of arrays in longitude direction
183!- plon      : Coordinates of points in longitude
184!- pjm       : Size of arrays in latitude direction
185!- plat      : Coordinates of points in latitude
186!-
187!- The next 4 arguments allow to define a horizontal zoom
188!- for this file. It is assumed that all variables to come
189!- have the same index space. This can not be assumed for
190!- the z axis and thus we define the zoom in histdef.
191!-
192!- par_orix  : Origin of the slab of data within the X axis (pim)
193!- par_szx   : Size of the slab of data in X
194!- par_oriy  : Origin of the slab of data within the Y axis (pjm)
195!- par_szy   : Size of the slab of data in Y
196!-
197!- pitau0    : time step at which the history tape starts
198!- pdate0    : The Julian date at which the itau was equal to 0
199!- pdeltat   : Time step in seconds. Time step of the counter itau
200!-             used in histwrt for instance
201!-
202!- OUTPUT
203!-
204!- phoriid   : ID of the horizontal grid
205!- pfileid   : ID of the netcdf file
206!-
207!- Optional INPUT arguments
208!-
209!- domain_id       : Domain identifier
210!-
211!- TO DO
212!-
213!- This package should be written in f90
214!- and use the following features :
215!-    - structures for the meta-data of the files and variables
216!-    - memory allocation as needed
217!-    - Pointers
218!-
219!- VERSION
220!-
221!---------------------------------------------------------------------
222  IMPLICIT NONE
223!-
224  CHARACTER(LEN=*) :: pfilename
225  INTEGER,INTENT(IN) :: pim,pjm
226  REAL,DIMENSION(pim),INTENT(IN) :: plon
227  REAL,DIMENSION(pjm),INTENT(IN) :: plat
228  INTEGER,INTENT(IN):: par_orix,par_szx,par_oriy,par_szy
229  INTEGER,INTENT(IN) :: pitau0
230  REAL,INTENT(IN) :: pdate0,pdeltat
231  INTEGER,INTENT(OUT) :: pfileid,phoriid
232  INTEGER,INTENT(IN),OPTIONAL :: domain_id
233!-
234  LOGICAL :: check = .FALSE.
235!-
236  REAL,ALLOCATABLE,DIMENSION(:,:) :: lon_tmp,lat_tmp
237!---------------------------------------------------------------------
238  IF (check) WRITE(*,*) "histbeg_totreg"
239!-
240  ALLOCATE (lon_tmp(pim,pjm),lat_tmp(pim,pjm))
241!-
242  lon_tmp(:,:) = SPREAD(plon(:),2,pjm)
243  lat_tmp(:,:) = SPREAD(plat(:),1,pim)
244!-
245  CALL histbeg_regular &
246 &  (pfilename,pim,lon_tmp,pjm,lat_tmp, &
247 &   par_orix,par_szx,par_oriy,par_szy, &
248 &   pitau0,pdate0,pdeltat,phoriid,pfileid, &
249 &   .TRUE.,domain_id)
250!-
251  DEALLOCATE (lon_tmp, lat_tmp)
252!----------------------------
253END SUBROUTINE histbeg_totreg
254!===
255SUBROUTINE histbeg_regular &
256 & (pfilename, pim, plon, pjm, plat,         &
257 &  par_orix, par_szx, par_oriy, par_szy,    &
258 &  pitau0, pdate0, pdeltat, phoriid, pfileid, &
259 &  opt_rectilinear, domain_id)
260!---------------------------------------------------------------------
261!- This subroutine initializes a netcdf file and returns the ID.
262!- It will set up the geographical space on which the data will be
263!- stored and offers the possibility of seting a zoom.
264!- It also gets the global parameters into the I/O subsystem.
265!-
266!- INPUT
267!-
268!- pfilename : Name of the netcdf file to be created
269!- pim       : Size of arrays in longitude direction
270!- plon      : Coordinates of points in longitude
271!- pjm       : Size of arrays in latitude direction
272!- plat      : Coordinates of points in latitude
273!-
274!- The next 4 arguments allow to define a horizontal zoom
275!- for this file. It is assumed that all variables to come
276!- have the same index space. This can not be assumed for
277!- the z axis and thus we define the zoom in histdef.
278!-
279!- par_orix  : Origin of the slab of data within the X axis (pim)
280!- par_szx   : Size of the slab of data in X
281!- par_oriy  : Origin of the slab of data within the Y axis (pjm)
282!- par_szy   : Size of the slab of data in Y
283!-
284!- pitau0    : time step at which the history tape starts
285!- pdate0    : The Julian date at which the itau was equal to 0
286!- pdeltat   : Time step in seconds. Time step of the counter itau
287!-             used in histwrt for instance
288!-
289!- OUTPUT
290!-
291!- phoriid   : ID of the horizontal grid
292!- pfileid   : ID of the netcdf file
293!-
294!- Optional INPUT arguments
295!-
296!- opt_rectilinear : If true we know the grid is rectilinear
297!- domain_id       : Domain identifier
298!-
299!- TO DO
300!-
301!- This package should be written in F90 and use the following
302!- feature :
303!-   - structures for the meta-data of the files and variables
304!-   - memory allocation as needed
305!-   - Pointers
306!-
307!- VERSION
308!-
309!---------------------------------------------------------------------
310  IMPLICIT NONE
311!-
312  CHARACTER(LEN=*) :: pfilename
313  INTEGER,INTENT(IN) :: pim,pjm
314  REAL,DIMENSION(pim,pjm),INTENT(IN) :: plon,plat
315  INTEGER,INTENT(IN):: par_orix, par_szx, par_oriy, par_szy
316  INTEGER,INTENT(IN) :: pitau0
317  REAL,INTENT(IN) :: pdate0, pdeltat
318  INTEGER,INTENT(OUT) :: pfileid, phoriid
319  LOGICAL,INTENT(IN),OPTIONAL :: opt_rectilinear
320  INTEGER,INTENT(IN),OPTIONAL :: domain_id
321!-
322  INTEGER :: ncid, iret
323  INTEGER :: lengf, lenga
324  CHARACTER(LEN=120) :: file
325  CHARACTER(LEN=30) :: timenow
326  LOGICAL :: rectilinear
327!-
328  LOGICAL :: check = .FALSE.
329!---------------------------------------------------------------------
330  nb_files = nb_files+1
331  pfileid  = nb_files
332!-
333! 1.0 Transfering into the common for future use
334!-
335  IF (check) WRITE(*,*) "histbeg_regular 1.0"
336!-
337  itau0(pfileid) = pitau0
338  date0(pfileid) = pdate0
339  deltat(pfileid) = pdeltat
340!-
341  IF (PRESENT(opt_rectilinear)) THEN
342    rectilinear = opt_rectilinear
343  ELSE
344    rectilinear = .FALSE.
345  ENDIF
346!-
347! 2.0 Initializes all variables for this file
348!-
349  IF (check) WRITE(*,*) "histbeg_regular 2.0"
350!-
351  IF (nb_files > nb_files_max) THEN
352    CALL ipslerr (3,"histbeg", &
353   &  'Table of files too small. You should increase nb_files_max', &
354   &  'in M_HISTCOM.f90 in order to accomodate all these files', ' ')
355  ENDIF
356!-
357  nb_var(pfileid) = 0
358  nb_tax(pfileid) = 0
359  nb_hax(pfileid) = 0
360  nb_zax(pfileid) = 0
361!-
362  slab_ori(pfileid,1:2) = (/ par_orix, par_oriy /)
363  slab_sz(pfileid,1:2)  = (/ par_szx,  par_szy  /)
364!-
365! 3.0 Opening netcdf file and defining dimensions
366!-
367  IF (check) WRITE(*,*) "histbeg_regular 3.0"
368!-
369! Add DOMAIN number and ".nc" suffix in file name if needed
370!-
371  file  = pfilename
372  CALL flio_dom_file (file,domain_id)
373!-
374! Keep track of the name of the files opened
375!-
376  lengf=LEN_TRIM(file)
377  lenga=LEN_TRIM(assc_file)
378  IF (nb_files == 1) THEN
379    assc_file=file(1:lengf)
380  ELSE IF ( (lenga+lengf) < 500) THEN
381    assc_file = assc_file(1:lenga)//' '//file(1:lengf)
382  ELSE IF (     ((lenga+7) < 500) &
383         & .AND.(INDEX(assc_file(1:lenga),'et.al.') < 1) ) THEN
384    assc_file = assc_file(1:lenga)//' et.al.'
385  ELSE
386    CALL ipslerr (2,"histbeg", &
387   & 'The file names do not fit into the associate_file attribute.', &
388   & 'Use shorter names if you wish to keep the information.',' ')
389  ENDIF
390!-
391  iret = NF90_CREATE (file, NF90_CLOBBER, ncid)
392!-
393  IF (rectilinear) THEN
394    iret = NF90_DEF_DIM (ncid, 'lon', par_szx, xid(nb_files))
395    iret = NF90_DEF_DIM (ncid, 'lat', par_szy, yid(nb_files))
396  ELSE
397    iret = NF90_DEF_DIM (ncid, 'x', par_szx, xid(nb_files))
398    iret = NF90_DEF_DIM (ncid, 'y', par_szy, yid(nb_files))
399  ENDIF
400!-
401! 4.0 Declaring the geographical coordinates and other attributes
402!-
403  IF (check) WRITE(*,*) "histbeg_regular 4.0"
404!-
405! 4.3 Global attributes
406!-
407  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'Conventions','GDT 1.3')
408  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'file_name',TRIM(file))
409  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'production',TRIM(model_name))
410  lock_modname = .TRUE.
411  CALL ioget_timestamp (timenow)
412  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'TimeStamp',TRIM(timenow))
413!-
414! Add DOMAIN attributes if needed
415!-
416  CALL flio_dom_att (ncid,domain_id)
417!-
418! 5.0 Saving some important information on this file in the common
419!-
420  IF (check) WRITE(*,*) "histbeg_regular 5.0"
421!-
422  ncdf_ids(pfileid) = ncid
423  full_size(pfileid,1:2) = (/ pim, pjm /)
424!-
425! 6.0 storing the geographical coordinates
426!-
427  IF ( (pim /= par_szx).OR.(pjm /= par_szy) )   zoom(pfileid)=.TRUE.
428  regular(pfileid)=.TRUE.
429!-
430  CALL histhori_regular (pfileid, pim, plon, pjm, plat, &
431 &  ' ' , 'Default grid', phoriid, rectilinear)
432!-----------------------------
433END SUBROUTINE histbeg_regular
434!===
435SUBROUTINE histbeg_irregular &
436 &  (pfilename, pim, plon, plon_bounds, plat, plat_bounds,   &
437 &   pitau0, pdate0, pdeltat, phoriid, pfileid, domain_id)
438!---------------------------------------------------------------------
439!- This subroutine initializes a netcdf file and returns the ID.
440!- This version is for totaly irregular grids. In this case all
441!- all the data comes in as vectors and for the grid we have
442!- the coordinates of the 4 corners.
443!- It also gets the global parameters into the I/O subsystem.
444!-
445!- INPUT
446!-
447!- pfilename   : Name of the netcdf file to be created
448!- pim         : Size of arrays in longitude direction
449!- plon        : Coordinates of points in longitude
450!- plon_bounds : The 2 corners of the grid in longitude
451!- plat        : Coordinates of points in latitude
452!- plat_bounds : The 2 corners of the grid in latitude
453!-
454!- pitau0    : time step at which the history tape starts
455!- pdate0    : The Julian date at which the itau was equal to 0
456!- pdeltat   : Time step in seconds. Time step of the counter itau
457!-             used in histwrt for instance
458!-
459!- OUTPUT
460!-
461!- phoriid   : ID of the horizontal grid
462!- pfileid   : ID of the netcdf file
463!-
464!- Optional INPUT arguments
465!-
466!- domain_id       : Domain identifier
467!-
468!- TO DO
469!-
470!- This package should be written in F90 and use the following
471!- feature :
472!-   - structures for the meta-data of the files and variables
473!-   - memory allocation as needed
474!-   - Pointers
475!-
476!- VERSION
477!-
478!---------------------------------------------------------------------
479  IMPLICIT NONE
480!-
481  CHARACTER(LEN=*) :: pfilename
482  INTEGER,INTENT(IN) :: pim
483  REAL,DIMENSION(pim),INTENT(IN) :: plon,plat
484  REAL,DIMENSION(:,:),INTENT(IN) :: plon_bounds,plat_bounds
485  INTEGER,INTENT(IN) :: pitau0
486  REAL,INTENT(IN) :: pdate0, pdeltat
487  INTEGER,INTENT(OUT) :: pfileid, phoriid
488  INTEGER,INTENT(IN),OPTIONAL :: domain_id
489!-
490  INTEGER :: ncid, iret
491  INTEGER :: lengf, lenga
492  CHARACTER(LEN=120) :: file
493  CHARACTER(LEN=30) :: timenow
494!-
495  LOGICAL :: check = .FALSE.
496!---------------------------------------------------------------------
497  nb_files = nb_files+1
498  pfileid  = nb_files
499!-
500! 1.0 Transfering into the common for future use
501!-
502  IF (check) WRITE(*,*) "histbeg_irregular 1.0"
503!-
504  itau0(pfileid) = pitau0
505  date0(pfileid) = pdate0
506  deltat(pfileid) = pdeltat
507!-
508! 2.0 Initializes all variables for this file
509!-
510  IF (check) WRITE(*,*) "histbeg_irregular 2.0"
511!-
512  IF (nb_files > nb_files_max) THEN
513    CALL ipslerr (3,"histbeg", &
514   &  'Table of files too small. You should increase nb_files_max', &
515   &  'in M_HISTCOM.f90 in order to accomodate all these files', ' ')
516  ENDIF
517!-
518  nb_var(pfileid) = 0
519  nb_tax(pfileid) = 0
520  nb_hax(pfileid) = 0
521  nb_zax(pfileid) = 0
522!-
523  slab_ori(pfileid,1:2) = (/ 1, 1 /)
524  slab_sz(pfileid,1:2)  = (/ pim, 1 /)
525!-
526! 3.0 Opening netcdf file and defining dimensions
527!-
528  IF (check) WRITE(*,*) "histbeg_irregular 3.0"
529!-
530! Add DOMAIN number and ".nc" suffix in file name if needed
531!-
532  file  = pfilename
533  CALL flio_dom_file (file,domain_id)
534!-
535! Keep track of the name of the files opened
536!-
537  lengf=LEN_TRIM(file)
538  lenga=LEN_TRIM(assc_file)
539  IF (nb_files == 1) THEN
540    assc_file=file(1:lengf)
541  ELSE IF ( (lenga+lengf) < 500) THEN
542    assc_file = assc_file(1:lenga)//' '//file(1:lengf)
543  ELSE IF (     ((lenga+7) < 500) &
544         & .AND.(INDEX(assc_file(1:lenga),'et.al.') < 1) ) THEN
545    assc_file = assc_file(1:lenga)//' et.al.'
546  ELSE
547    CALL ipslerr (2,"histbeg", &
548   & 'The file names do not fit into the associate_file attribute.', &
549   & 'Use shorter names if you wish to keep the information.',' ')
550  ENDIF
551!-
552  iret = NF90_CREATE (file, NF90_CLOBBER, ncid)
553!-
554  iret = NF90_DEF_DIM (ncid, 'x', pim, xid(nb_files))
555  yid(nb_files) = 0
556!-
557!- 4.0 Declaring the geographical coordinates and other attributes
558!-
559   IF (check) WRITE(*,*) "histbeg_irregular 4.0"
560!-
561! 4.3 Global attributes
562!-
563  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'Conventions','GDT 1.3')
564  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'file_name',TRIM(file))
565  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'production',TRIM(model_name))
566  lock_modname = .TRUE.
567  CALL ioget_timestamp (timenow)
568  iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'TimeStamp',TRIM(timenow))
569!-
570! Add DOMAIN attributes if needed
571!-
572  CALL flio_dom_att (ncid,domain_id)
573!-
574! 5.0 Saving some important information on this file in the common
575!-
576  IF (check) WRITE(*,*) "histbeg_irregular 5.0"
577!-
578  ncdf_ids(pfileid) = ncid
579  full_size(pfileid,1:2) = (/ pim, 1 /)
580!-
581! 6.0 storing the geographical coordinates
582!-
583  zoom(pfileid)=.FALSE.
584  regular(pfileid)=.FALSE.
585!-
586  CALL histhori_irregular &
587 &  (pfileid, pim, plon, plon_bounds, plat, plat_bounds, &
588 &   ' ' , 'Default grid', phoriid)
589!-------------------------------
590END SUBROUTINE histbeg_irregular
591!===
592SUBROUTINE histhori_regular &
593 &  (pfileid,pim,plon,pjm,plat,phname,phtitle,phid,opt_rectilinear)
594!---------------------------------------------------------------------
595!- This subroutine is made to declare a new horizontale grid.
596!- It has to have the same number of points as
597!- the original and thus in this routine we will only
598!- add two variable (longitude and latitude).
599!- Any variable in the file can thus point to this pair
600!- through an attribute. This routine is very usefull
601!- to allow staggered grids.
602!-
603!- INPUT
604!-
605!- pfileid : The id of the file to which the grid should be added
606!- pim     : Size in the longitude direction
607!- plon    : The longitudes
608!- pjm     : Size in the latitude direction
609!- plat    : The latitudes
610!- phname  : The name of grid
611!- phtitle : The title of the grid
612!-
613!- OUTPUT
614!-
615!- phid    : Id of the created grid
616!-
617!- OPTIONAL
618!-
619!- opt_rectilinear : If true we know the grid is rectilinear.
620!-
621!---------------------------------------------------------------------
622  IMPLICIT NONE
623!-
624  INTEGER,INTENT(IN) :: pfileid, pim, pjm
625  REAL,INTENT(IN),DIMENSION(pim,pjm) :: plon, plat
626  CHARACTER(LEN=*),INTENT(IN) :: phname, phtitle
627  INTEGER,INTENT(OUT) :: phid
628  LOGICAL,INTENT(IN),OPTIONAL :: opt_rectilinear
629!-
630  CHARACTER(LEN=25) :: lon_name, lat_name
631  CHARACTER(LEN=80) :: tmp_title, tmp_name
632  INTEGER :: ndim
633  INTEGER,DIMENSION(2) :: dims(2)
634  INTEGER :: nlonid, nlatid
635  INTEGER :: orix, oriy, par_szx, par_szy
636  INTEGER :: iret, ncid
637  LOGICAL :: rectilinear
638!-
639  LOGICAL :: check = .FALSE.
640!---------------------------------------------------------------------
641!-
642! 1.0 Check that all fits in the buffers
643!-
644  IF (    (pim /= full_size(pfileid,1)) &
645    & .OR.(pjm /= full_size(pfileid,2)) ) THEN
646    CALL ipslerr (3,"histhori", &
647   &  'The new horizontal grid does not have the same size', &
648   &  'as the one provided to histbeg. This is not yet ', &
649   &  'possible in the hist package.')
650  ENDIF
651!-
652  IF (PRESENT(opt_rectilinear)) THEN
653    rectilinear = opt_rectilinear
654  ELSE
655    rectilinear = .FALSE.
656  ENDIF
657!-
658! 1.1 Create all the variables needed
659!-
660  IF (check) WRITE(*,*) "histhori_regular 1.0"
661!-
662  ncid = ncdf_ids(pfileid)
663!-
664  ndim = 2
665  dims(1:2) = (/ xid(pfileid), yid(pfileid) /)
666!-
667  tmp_name = phname
668  IF (rectilinear) THEN
669    IF (nb_hax(pfileid) == 0) THEN
670      lon_name = 'lon'
671      lat_name = 'lat'
672    ELSE
673      lon_name = 'lon_'//TRIM(tmp_name)
674      lat_name = 'lat_'//TRIM(tmp_name)
675    ENDIF
676  ELSE
677    IF (nb_hax(pfileid) == 0) THEN
678      lon_name = 'nav_lon'
679      lat_name = 'nav_lat'
680    ELSE
681      lon_name = 'nav_lon_'//TRIM(tmp_name)
682      lat_name = 'nav_lat_'//TRIM(tmp_name)
683    ENDIF
684  ENDIF
685!-
686! 1.2 Save the informations
687!-
688  phid =  nb_hax(pfileid)+1
689  nb_hax(pfileid) = phid
690!-
691  hax_name(pfileid,phid,1:2) = (/ lon_name, lat_name /)
692  tmp_title = phtitle
693!-
694! 2.0 Longitude
695!-
696  IF (check) WRITE(*,*) "histhori_regular 2.0"
697!-
698  IF (rectilinear) THEN
699    ndim = 1
700    dims(1:1) = (/ xid(pfileid) /)
701  ENDIF
702  iret = NF90_DEF_VAR (ncid,lon_name,NF90_FLOAT,dims(1:ndim),nlonid)
703  iret = NF90_PUT_ATT (ncid,nlonid,'units',"degrees_east")
704  iret = NF90_PUT_ATT (ncid,nlonid,'valid_min', &
705 &                     REAL(MINVAL(plon),KIND=4))
706  iret = NF90_PUT_ATT (ncid,nlonid,'valid_max', &
707 &                     REAL(MAXVAL(plon),KIND=4))
708  iret = NF90_PUT_ATT (ncid,nlonid,'long_name',"Longitude")
709  iret = NF90_PUT_ATT (ncid,nlonid,'nav_model',TRIM(tmp_title))
710!-
711! 3.0 Latitude
712!-
713  IF (check) WRITE(*,*) "histhori_regular 3.0"
714!-
715  IF (rectilinear) THEN
716    ndim = 1
717    dims(1:1) = (/ yid(pfileid) /)
718  ENDIF
719  iret = NF90_DEF_VAR (ncid,lat_name,NF90_FLOAT,dims(1:ndim),nlatid)
720  iret = NF90_PUT_ATT (ncid,nlatid,'units',"degrees_north")
721  iret = NF90_PUT_ATT (ncid,nlatid,'valid_min', &
722 &                     REAL(MINVAL(plat),KIND=4))
723  iret = NF90_PUT_ATT (ncid,nlatid,'valid_max', &
724 &                     REAL(MAXVAL(plat),KIND=4))
725  iret = NF90_PUT_ATT (ncid,nlatid,'long_name',"Latitude")
726  iret = NF90_PUT_ATT (ncid,nlatid,'nav_model',TRIM(tmp_title))
727!-
728  iret = NF90_ENDDEF (ncid)
729!-
730! 4.0 storing the geographical coordinates
731!-
732  IF (check) WRITE(*,*) "histhori_regular 4.0"
733!-
734  orix = slab_ori(pfileid,1)
735  oriy = slab_ori(pfileid,2)
736  par_szx = slab_sz(pfileid,1)
737  par_szy = slab_sz(pfileid,2)
738!-
739! Transfer the longitude
740!-
741  IF (rectilinear) THEN
742    iret = NF90_PUT_VAR (ncid,nlonid,plon(1:par_szx,1))
743  ELSE
744    iret = NF90_PUT_VAR (ncid,nlonid,plon, &
745 &           start=(/orix,oriy/),count=(/par_szx,par_szy/))
746  ENDIF
747!-
748! Transfer the latitude
749!-
750  IF ( rectilinear ) THEN
751    iret = NF90_PUT_VAR (ncid,nlatid,plat(1,1:par_szy))
752  ELSE
753    iret = NF90_PUT_VAR (ncid,nlatid,plat, &
754 &           start=(/orix,oriy/),count=(/par_szx,par_szy/))
755  ENDIF
756!-
757  iret = NF90_REDEF (ncid)
758!------------------------------
759END SUBROUTINE histhori_regular
760!===
761SUBROUTINE histhori_irregular &
762 &  (pfileid, pim, plon, plon_bounds, plat, plat_bounds, &
763 &   phname, phtitle, phid)
764!---------------------------------------------------------------------
765!- This subroutine is made to declare a new horizontale grid.
766!- It has to have the same number of points as
767!- the original and thus in this routine we will only
768!- add two variable (longitude and latitude).
769!- Any variable in the file can thus point to this pair
770!- through an attribute. This routine is very usefull
771!- to allow staggered grids.
772!-
773!- INPUT
774!-
775!- pfileid     : The id of the file to which the grid should be added
776!- pim         : Size in the longitude direction
777!- plon        : The longitudes
778!- plon_bounds : The boundaries of the grid in longitude
779!- plat        : The latitudes
780!- plat_bounds : Boundaries of the grid in latitude
781!- phname      : The name of grid
782!- phtitle     : The title of the grid
783!-
784!- OUTPUT
785!-
786!- phid    : Id of the created grid
787!---------------------------------------------------------------------
788  IMPLICIT NONE
789!-
790  INTEGER,INTENT(IN) :: pfileid, pim
791  REAL,DIMENSION(pim),INTENT(IN) :: plon, plat
792  REAL,DIMENSION(:,:),INTENT(IN) :: plon_bounds, plat_bounds
793  CHARACTER(LEN=*), INTENT(IN) :: phname, phtitle
794  INTEGER,INTENT(OUT) :: phid
795!-
796  CHARACTER(LEN=25) :: lon_name, lat_name
797  CHARACTER(LEN=30) :: lonbound_name, latbound_name
798  CHARACTER(LEN=80) :: tmp_title, tmp_name, longname
799  INTEGER :: ndim, dims(2)
800  INTEGER :: ndimb, dimsb(2)
801  INTEGER :: nbbounds
802  INTEGER :: nlonid, nlatid, nlonidb, nlatidb
803  INTEGER :: iret, ncid, twoid
804  LOGICAL :: transp = .FALSE.
805  REAL, ALLOCATABLE, DIMENSION(:,:) :: bounds_trans
806!-
807  LOGICAL :: check = .FALSE.
808!---------------------------------------------------------------------
809!-
810! 1.0 Check that all fits in the buffers
811!-
812  IF (    (pim /= full_size(pfileid,1)) &
813    & .OR.(full_size(pfileid,2) /= 1) ) THEN
814    CALL ipslerr (3,"histhori", &
815   &  'The new horizontal grid does not have the same size', &
816   &  'as the one provided to histbeg. This is not yet ', &
817   &  'possible in the hist package.')
818  ENDIF
819!-
820! 1.1 Create all the variables needed
821!-
822  IF (check) WRITE(*,*) 'histhori_irregular 1.0'
823!-
824  ncid = ncdf_ids(pfileid)
825!-
826  IF ( SIZE(plon_bounds, DIM=1) == pim ) THEN
827    nbbounds = SIZE(plon_bounds, DIM=2)
828    transp = .TRUE.
829  ELSEIF ( SIZE(plon_bounds, DIM=2) == pim ) THEN
830    nbbounds = SIZE(plon_bounds, DIM=1)
831    transp = .FALSE.
832  ELSE
833    CALL ipslerr (3,"histhori", &
834   &  'The boundary variable does not have any axis corresponding', &
835   &  'to the size of the longitude or latitude variable', &
836   &  '.')
837  ENDIF
838!-
839  IF (.NOT.ALLOCATED(bounds_trans)) THEN
840    ALLOCATE(bounds_trans(nbbounds,pim))
841  ENDIF
842!-
843  iret = NF90_DEF_DIM (ncid, 'nbnd', nbbounds, twoid)
844  ndim = 1
845  dims(1) = xid(pfileid)
846  ndimb = 2
847  dimsb(1:2) = (/ twoid, xid(pfileid) /)
848!-
849  tmp_name = phname
850  IF (nb_hax(pfileid) == 0) THEN
851    lon_name = 'nav_lon'
852    lat_name = 'nav_lat'
853  ELSE
854    lon_name = 'nav_lon_'//TRIM(tmp_name)
855    lat_name = 'nav_lat_'//TRIM(tmp_name)
856  ENDIF
857  lonbound_name = TRIM(lon_name)//'_bounds'
858  latbound_name = TRIM(lat_name)//'_bounds'
859!-
860! 1.2 Save the informations
861!-
862  phid =  nb_hax(pfileid)+1
863  nb_hax(pfileid) = phid
864!-
865  hax_name(pfileid,phid,1:2) = (/ lon_name, lat_name /)
866  tmp_title = phtitle
867!-
868! 2.0 Longitude
869!-
870  IF (check) WRITE(*,*) "histhori_irregular 2.0"
871!-
872  iret = NF90_DEF_VAR (ncid,lon_name,NF90_FLOAT,dims(1:ndim),nlonid)
873  iret = NF90_PUT_ATT (ncid,nlonid,'units',"degrees_east")
874  iret = NF90_PUT_ATT (ncid,nlonid,'valid_min', &
875 &                     REAL(MINVAL(plon),KIND=4))
876  iret = NF90_PUT_ATT (ncid,nlonid,'valid_max', &
877 &                     REAL(MAXVAL(plon),KIND=4))
878  iret = NF90_PUT_ATT (ncid,nlonid,'long_name',"Longitude")
879  iret = NF90_PUT_ATT (ncid,nlonid,'nav_model',TRIM(tmp_title))
880!-
881! 2.1 Longitude bounds
882!-
883  iret = NF90_PUT_ATT (ncid,nlonid,'bounds',TRIM(lonbound_name))
884  iret = NF90_DEF_VAR (ncid,lonbound_name,NF90_FLOAT, &
885 &                     dimsb(1:ndimb),nlonidb)
886  longname = 'Boundaries for coordinate variable '//TRIM(lon_name)
887  iret = NF90_PUT_ATT (ncid,nlonidb,'long_name',TRIM(longname))
888!-
889! 3.0 Latitude
890!-
891  IF (check) WRITE(*,*) "histhori_irregular 3.0"
892!-
893  iret = NF90_DEF_VAR (ncid,lat_name,NF90_FLOAT,dims(1:ndim),nlatid)
894  iret = NF90_PUT_ATT (ncid,nlatid,'units',"degrees_north")
895  iret = NF90_PUT_ATT (ncid,nlatid,'valid_min', &
896 &                     REAL(MINVAL(plat),KIND=4))
897  iret = NF90_PUT_ATT (ncid,nlatid,'valid_max', &
898 &                     REAL(MAXVAL(plat),KIND=4))
899  iret = NF90_PUT_ATT (ncid,nlatid,'long_name',"Latitude")
900  iret = NF90_PUT_ATT (ncid,nlatid,'nav_model',TRIM(tmp_title))
901!-
902! 3.1 Latitude bounds
903!-
904  iret = NF90_PUT_ATT (ncid,nlatid,'bounds',TRIM(latbound_name))
905  iret = NF90_DEF_VAR (ncid,latbound_name,NF90_FLOAT, &
906 &                     dimsb(1:ndimb),nlatidb)
907  longname = 'Boundaries for coordinate variable '//TRIM(lat_name)
908  iret = NF90_PUT_ATT (ncid,nlatidb,'long_name',TRIM(longname))
909!-
910  iret = NF90_ENDDEF (ncid)
911!-
912! 4.0 storing the geographical coordinates
913!-
914  IF (check) WRITE(*,*) "histhori_irregular 4.0"
915!-
916! 4.1 Write the longitude
917!-
918  iret = NF90_PUT_VAR (ncid, nlonid, plon(1:pim))
919!-
920! 4.2 Write the longitude bounds
921!-
922  IF (transp) THEN
923    bounds_trans = TRANSPOSE(plon_bounds)
924  ELSE
925    bounds_trans = plon_bounds
926  ENDIF
927  iret = NF90_PUT_VAR (ncid, nlonidb, bounds_trans(1:nbbounds,1:pim))
928!-
929! 4.3 Write the latitude
930!-
931  iret = NF90_PUT_VAR (ncid, nlatid, plat(1:pim))
932!-
933! 4.4 Write the latitude bounds
934!-
935  IF (transp) THEN
936    bounds_trans = TRANSPOSE(plat_bounds)
937  ELSE
938    bounds_trans = plat_bounds
939  ENDIF
940  iret = NF90_PUT_VAR (ncid,nlatidb,bounds_trans(1:nbbounds,1:pim))
941!-
942  iret = NF90_REDEF (ncid)
943!--------------------------------
944END SUBROUTINE histhori_irregular
945!===
946SUBROUTINE histvert (pfileid, pzaxname, pzaxtitle, &
947 &                   pzaxunit, pzsize, pzvalues, pzaxid, pdirect)
948!---------------------------------------------------------------------
949!- This subroutine defines a vertical axis and returns it s id.
950!- It gives the user the possibility to the user to define many
951!- different vertical axes. For each variable defined with histdef a
952!- vertical axis can be specified with by it s ID.
953!-
954!- INPUT
955!-
956!- pfileid  : ID of the file the variable should be archived in
957!- pzaxname : Name of the vertical axis
958!- pzaxtitle: title of the vertical axis
959!- pzaxunit : Units of the vertical axis
960!- pzsize   : size of the vertical axis
961!- pzvalues : Coordinate values of the vetical axis
962!-
963!- pdirect  : is an optional argument which allows to specify the
964!-            the positive direction of the axis : up or down.
965!- OUTPUT
966!-
967!- pzaxid   : Returns the ID of the axis.
968!-            Note that this is not the netCDF ID !
969!-
970!- VERSION
971!-
972!---------------------------------------------------------------------
973  IMPLICIT NONE
974!-
975  INTEGER,INTENT(IN) :: pfileid,pzsize
976  CHARACTER(LEN=*),INTENT(IN) :: pzaxname,pzaxunit,pzaxtitle
977  REAL,INTENT(IN) :: pzvalues(pzsize)
978  INTEGER, INTENT(OUT) :: pzaxid
979  CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: pdirect
980!-
981  INTEGER :: pos, iv, nb, zdimid, zaxid_tmp
982  CHARACTER(LEN=20) :: str20, tab_str20(nb_zax_max)
983  INTEGER tab_str20_length(nb_zax_max)
984  CHARACTER(LEN=70) :: str70, str71, str72
985  CHARACTER(LEN=80) :: str80
986  CHARACTER(LEN=20) :: direction
987  INTEGER :: iret, leng, ncid
988  LOGICAL :: check = .FALSE.
989!---------------------------------------------------------------------
990!-
991! 1.0 Verifications :
992!    Do we have enough space for an extra axis ?
993!    Is the name already in use ?
994!-
995  IF (check) WRITE(*,*) "histvert : 1.0 Verifications", &
996 &                      pzaxname,'---',pzaxunit,'---',pzaxtitle
997!-
998! - Direction of axis. Can we get if from the user.
999!   If not we put unknown.
1000!-
1001  IF (PRESENT(pdirect)) THEN
1002    direction = TRIM(pdirect)
1003    CALL strlowercase (direction)
1004  ELSE
1005    direction = 'unknown'
1006  ENDIF
1007!-
1008! Check the consistency of the attribute
1009!-
1010  IF (     (direction /= 'unknown') &
1011 &    .AND.(direction /= 'up')      &
1012 &    .AND.(direction /= 'down')   ) THEN
1013    direction = 'unknown'
1014    str80 = 'The specified axis was : '//TRIM(direction)
1015    CALL ipslerr (2,"histvert",&
1016   & "The specified direction for the vertical axis is not possible.",&
1017   & "it is replaced by : unknown", str80)
1018  ENDIF
1019!-
1020  IF ( nb_zax(pfileid)+1 > nb_zax_max) THEN
1021    CALL ipslerr (3,"histvert", &
1022   &  'Table of vertical axes too small. You should increase ',&
1023   &  'nb_zax_max in M_HISTCOM.f90 in order to accomodate all ', &
1024   &  'these variables ')
1025  ENDIF
1026!-
1027  iv = nb_zax(pfileid)
1028  IF ( iv > 1) THEN
1029    str20 = pzaxname
1030    nb = iv-1
1031    tab_str20(1:nb) = zax_name(pfileid,1:nb)
1032    tab_str20_length(1:nb) = zax_name_length(pfileid,1:nb)
1033    CALL find_str(nb, tab_str20, tab_str20_length, str20, pos)
1034  ELSE
1035    pos = 0
1036  ENDIF
1037!-
1038  IF ( pos > 0) THEN
1039    str70 = "Vertical axis already exists"
1040    WRITE(str71,'("Check variable ",a," in file",I3)') str20,pfileid
1041    str72 = "Can also be a wrong file ID in another declaration"
1042    CALL ipslerr (3,"histvert", str70, str71, str72)
1043  ENDIF
1044!-
1045  iv = nb_zax(pfileid)+1
1046!-
1047! 2.0 Add the information to the file
1048!-
1049  IF (check) &
1050 &  WRITE(*,*) "histvert : 2.0 Add the information to the file"
1051!-
1052  ncid = ncdf_ids(pfileid)
1053!-
1054  leng = MIN(LEN_TRIM(pzaxname),20)
1055  iret = NF90_DEF_DIM (ncid,pzaxname(1:leng),pzsize,zaxid_tmp)
1056  iret = NF90_DEF_VAR (ncid,pzaxname(1:leng),NF90_FLOAT, &
1057 &                     zaxid_tmp,zdimid)
1058!-
1059  leng = MIN(LEN_TRIM(pzaxunit),20)
1060  iret = NF90_PUT_ATT (ncid, zdimid, 'units', pzaxunit(1:leng))
1061  iret = NF90_PUT_ATT (ncid, zdimid, 'positive', TRIM(direction))
1062!-
1063  iret = NF90_PUT_ATT (ncid,zdimid,'valid_min', &
1064 &                     REAL(MINVAL(pzvalues(1:pzsize)),KIND=4))
1065  iret = NF90_PUT_ATT (ncid,zdimid,'valid_max', &
1066 &                     REAL(MAXVAL(pzvalues(1:pzsize)),KIND=4))
1067!-
1068  leng = MIN(LEN_TRIM(pzaxname),20)
1069  iret = NF90_PUT_ATT (ncid, zdimid, 'title', pzaxname(1:leng))
1070  leng = MIN(LEN_TRIM(pzaxtitle),80)
1071  iret = NF90_PUT_ATT (ncid, zdimid, 'long_name', pzaxtitle(1:leng))
1072!-
1073  iret = NF90_ENDDEF (ncid)
1074!-
1075  iret = NF90_PUT_VAR (ncid, zdimid, pzvalues(1:pzsize))
1076!-
1077  iret = NF90_REDEF (ncid)
1078!-
1079!- 3.0 add the information to the common
1080!-
1081  IF ( check) &
1082  &  WRITE(*,*) "histvert : 3.0 add the information to the common"
1083!-
1084  nb_zax(pfileid) = iv
1085  zax_size(pfileid, iv) = pzsize
1086  zax_name(pfileid, iv) = pzaxname
1087  zax_name_length(pfileid, iv) = LEN_TRIM(pzaxname)
1088  zax_ids(pfileid, iv) = zaxid_tmp
1089  pzaxid =  iv
1090!----------------------
1091END SUBROUTINE histvert
1092!===
1093SUBROUTINE histdef (pfileid, pvarname, ptitle, punit, &
1094 &                  pxsize, pysize, phoriid, pzsize,  &
1095 &                  par_oriz, par_szz, pzid,          &
1096 &                  pnbbyt, popp, pfreq_opp, pfreq_wrt, var_range)
1097!---------------------------------------------------------------------
1098!- With this subroutine each variable to be archived on the history
1099!- tape should be declared.
1100!-
1101!- It gives the user the choise of operation
1102!- to be performed on the variables, the frequency of this operation
1103!- and finaly the frequency of the archiving.
1104!-
1105!- INPUT
1106!-
1107!- pfileid  : ID of the file the variable should be archived in
1108!- pvarname : Name of the variable, short and easy to remember
1109!- ptitle   : Full name of the variable
1110!- punit    : Units of the variable
1111!-
1112!- The next 3 arguments give the size of that data
1113!- that will be passed to histwrite. The zoom will be
1114!- done there with the horizontal information obtained
1115!- in histbeg and the vertical information to follow.
1116!-
1117!- pxsize   : Size in X direction (size of the data that will be
1118!-            given to histwrite)
1119!- pysize   : Size in Y direction
1120!- phoriid  : ID of the horizontal axis
1121!-
1122!- The next two arguments give the vertical zoom to use.
1123!-
1124!- pzsize   : Size in Z direction (If 1 then no axis is declared
1125!-            for this variable and pzid is not used)
1126!- par_oriz : Off set of the zoom
1127!- par_szz  : Size of the zoom
1128!-
1129!- pzid     : ID of the vertical axis to use. It has to have
1130!-            the size of the zoom.
1131!- pnbbyt   : Number of bytes on which to store in netCDF (Not opp.)
1132!- popp     : Operation to be performed. The following options
1133!-            exist today :
1134!- inst : keeps instantaneous values for writting
1135!- ave  : Computes the average from call between writes
1136!- pfreq_opp: Frequency of this operation (in seconds)
1137!- pfreq_wrt: Frequency at which the variable should be
1138!-            written (in seconds)
1139!- var_range: Range of the variable.
1140!-            If the minimum is greater than the maximum,
1141!-            the values will be calculated.
1142!-
1143!- VERSION
1144!---------------------------------------------------------------------
1145  IMPLICIT NONE
1146!-
1147  INTEGER,INTENT(IN) :: pfileid, pxsize, pysize, pzsize, pzid
1148  INTEGER,INTENT(IN) :: par_oriz, par_szz, pnbbyt, phoriid
1149  CHARACTER(LEN=*),INTENT(IN) :: pvarname, punit, popp
1150  CHARACTER(LEN=*),INTENT(IN) :: ptitle
1151  REAL,INTENT(IN) :: pfreq_opp, pfreq_wrt
1152  REAL,DIMENSION(2),OPTIONAL,INTENT(IN) :: var_range
1153!-
1154  INTEGER :: iv, i, nb
1155  CHARACTER(LEN=70) :: str70, str71, str72
1156  CHARACTER(LEN=20) :: tmp_name
1157  CHARACTER(LEN=20) :: str20, tab_str20(nb_var_max)
1158  INTEGER :: tab_str20_length(nb_var_max)
1159  CHARACTER(LEN=40) :: str40, tab_str40(nb_var_max)
1160  INTEGER :: tab_str40_length(nb_var_max)
1161  CHARACTER(LEN=10) :: str10
1162  CHARACTER(LEN=80) :: tmp_str80
1163  CHARACTER(LEN=7) :: tmp_topp, tmp_sopp(nbopp_max)
1164  CHARACTER(LEN=120) :: ex_topps
1165  REAL :: tmp_scal(nbopp_max), un_an, un_jour, test_fopp, test_fwrt
1166  INTEGER :: pos, buff_sz
1167!-
1168  LOGICAL :: check = .FALSE.
1169!---------------------------------------------------------------------
1170  ex_topps = 'ave, inst, t_min, t_max, t_sum, once, never, l_max, l_min'
1171!-
1172  nb_var(pfileid) = nb_var(pfileid)+1
1173  iv = nb_var(pfileid)
1174!-
1175  IF ( iv > nb_var_max) THEN
1176    CALL ipslerr (3,"histdef", &
1177   &  'Table of variables too small. You should increase nb_var_max',&
1178   &  'in M_HISTCOM.f90 in order to accomodate all these variables', &
1179   &  ' ')
1180  ENDIF
1181!-
1182! 1.0 Transfer informations on the variable to the common
1183!     and verify that it does not already exist
1184!-
1185  IF (check) WRITE(*,*) "histdef : 1.0"
1186!-
1187  IF (iv > 1) THEN
1188    str20 = pvarname
1189    nb = iv-1
1190    tab_str20(1:nb) = name(pfileid,1:nb)
1191    tab_str20_length(1:nb) = name_length(pfileid,1:nb)
1192    CALL find_str (nb, tab_str20, tab_str20_length, str20, pos)
1193  ELSE
1194    pos = 0
1195  ENDIF
1196!-
1197  IF (pos > 0) THEN
1198    str70 = "Variable already exists"
1199    WRITE(str71,'("Check variable  ",a," in file",I3)') str20,pfileid
1200    str72 = "Can also be a wrong file ID in another declaration"
1201    CALL ipslerr (3,"histdef", str70, str71, str72)
1202  ENDIF
1203!-
1204  name(pfileid,iv) = pvarname
1205  name_length(pfileid,iv) = LEN_TRIM(name(pfileid,iv))
1206  title(pfileid,iv) = ptitle
1207  unit_name(pfileid,iv) = punit
1208  tmp_name =  name(pfileid,iv)
1209!-
1210! 1.1 decode the operations
1211!-
1212  fullop(pfileid,iv) = popp
1213  tmp_str80 = popp
1214  CALL buildop &
1215 &  (tmp_str80, ex_topps, tmp_topp, nbopp_max, missing_val, &
1216 &   tmp_sopp, tmp_scal, nbopp(pfileid,iv))
1217!-
1218  topp(pfileid,iv) = tmp_topp
1219  DO i=1,nbopp(pfileid,iv)
1220    sopps(pfileid,iv,i) = tmp_sopp(i)
1221    scal(pfileid,iv,i) = tmp_scal(i)
1222  ENDDO
1223!-
1224! 1.2 If we have an even number of operations
1225!     then we need to add identity
1226!-
1227  IF (2*INT(nbopp(pfileid,iv)/2.0) == nbopp(pfileid,iv)) THEN
1228    nbopp(pfileid,iv) = nbopp(pfileid,iv)+1
1229    sopps(pfileid,iv,nbopp(pfileid,iv)) = 'ident'
1230    scal(pfileid,iv,nbopp(pfileid,iv)) = missing_val
1231  ENDIF
1232!-
1233! 2.0 Put the size of the variable in the common and check
1234!-
1235  IF (check) &
1236 &  WRITE(*,*) "histdef : 2.0", pfileid,iv,nbopp(pfileid,iv), &
1237 &    sopps(pfileid,iv,1:nbopp(pfileid,iv)), &
1238 &    scal(pfileid,iv,1:nbopp(pfileid,iv))
1239!-
1240  scsize(pfileid,iv,1:3) = (/ pxsize, pysize, pzsize /)
1241!-
1242  zorig(pfileid,iv,1:3) = &
1243 &  (/ slab_ori(pfileid,1), slab_ori(pfileid,2), par_oriz /)
1244!-
1245  zsize(pfileid,iv,1:3) = &
1246 &  (/ slab_sz(pfileid,1), slab_sz(pfileid,2), par_szz /)
1247!-
1248! Is the size of the full array the same as that of the coordinates  ?
1249!-
1250  IF (    (pxsize > full_size(pfileid,1)) &
1251 &    .OR.(pysize > full_size(pfileid,2)) ) THEN
1252!-
1253    str70 = "The size of the variable is different "// &
1254 &          "from the one of the coordinates"
1255    WRITE(str71,'("Size of coordinates :", 2I4)') &
1256 &   full_size(pfileid,1), full_size(pfileid,2)
1257    WRITE(str72,'("Size declared for variable ",a," :",2I4)') &
1258 &   TRIM(tmp_name), pxsize, pysize
1259    CALL ipslerr (3,"histdef", str70, str71, str72)
1260  ENDIF
1261!-
1262! Is the size of the zoom smaler than the coordinates ?
1263!-
1264  IF (    (full_size(pfileid,1) < slab_sz(pfileid,1)) &
1265 &    .OR.(full_size(pfileid,2) < slab_sz(pfileid,2)) ) THEN
1266    str70 = &
1267 &   "Size of variable should be greater or equal to those of the zoom"
1268    WRITE(str71,'("Size of XY zoom :", 2I4)') &
1269 &   slab_sz(pfileid,1),slab_sz(pfileid,1)
1270    WRITE(str72,'("Size declared for variable ",a," :",2I4)') &
1271 &   TRIM(tmp_name), pxsize, pysize
1272    CALL ipslerr (3,"histdef", str70, str71, str72)
1273  ENDIF
1274!-
1275! 2.1 We store the horizontal grid information with minimal
1276!     and a fall back onto the default grid
1277!-
1278  IF ( phoriid > 0 .AND. phoriid <= nb_hax(pfileid)) THEN
1279    var_haxid(pfileid,iv) = phoriid
1280  ELSE
1281    var_haxid(pfileid,iv) = 1
1282    CALL ipslerr (2,"histdef", &
1283   &  'We use the default grid for variable as an invalide',&
1284   &  'ID was provided for variable : ', pvarname)
1285  ENDIF
1286!-
1287! 2.2 Check the vertical coordinates if needed
1288!-
1289  IF (par_szz > 1) THEN
1290!-
1291!-- Does the vertical coordinate exist ?
1292!-
1293    IF ( pzid > nb_zax(pfileid)) THEN
1294      WRITE(str70, &
1295 &    '("The vertical coordinate chosen for variable ",a)') &
1296 &     TRIM(tmp_name)
1297      str71 = " Does not exist."
1298      CALL ipslerr (3,"histdef",str70,str71, " ")
1299    ENDIF
1300!-
1301!-- Is the vertical size of the variable equal to that of the axis ?
1302!-
1303    IF (par_szz /= zax_size(pfileid,pzid)) THEN
1304      str20 = zax_name(pfileid,pzid)
1305      str70 = "The size of the zoom does not correspond "// &
1306 &            "to the size of the chosen vertical axis"
1307      WRITE(str71,'("Size of zoom in z :", I4)') par_szz
1308      WRITE(str72,'("Size declared for axis ",a," :",I4)') &
1309 &     TRIM(str20), zax_size(pfileid,pzid)
1310      CALL ipslerr (3,"histdef", str70, str71, str72)
1311    ENDIF
1312!-
1313!-- Is the zoom smaler that the total size of the variable ?
1314!-
1315    IF ( pzsize < par_szz ) THEN
1316      str20 = zax_name(pfileid,pzid)
1317      str70 = "The vertical size of variable "// &
1318 &            "is smaller than that of the zoom."
1319      WRITE(str71,'("Declared vertical size of data :", I5)') pzsize
1320      WRITE(str72,'("Size of zoom for variable ",a," = ",I5)') &
1321 &     TRIM(tmp_name),par_szz
1322      CALL ipslerr (3,"histdef", str70, str71, str72)
1323    ENDIF
1324    var_zaxid(pfileid,iv) = pzid
1325  ELSE
1326    var_zaxid(pfileid,iv) = -99
1327  ENDIF
1328!-
1329! 3.0 Determine the position of the variable in the buffer
1330!     If it is instantaneous output then we do not use the buffer
1331!-
1332  IF (check) WRITE(*,*) "histdef : 3.0"
1333!-
1334! 3.1 We get the size of the arrays histwrite will get and check
1335!     that they fit into the tmp_buffer
1336!-
1337  buff_sz = zsize(pfileid,iv,1)*zsize(pfileid,iv,2)*zsize(pfileid,iv,3)
1338!-
1339! 3.2 move the pointer of the buffer array for operation
1340!     which need bufferisation
1341!-
1342  IF (     (TRIM(tmp_topp) /= "inst") &
1343 &    .AND.(TRIM(tmp_topp) /= "once") &
1344 &    .AND.(TRIM(tmp_topp) /= "never") )THEN
1345    point(pfileid,iv) = buff_pos+1
1346    buff_pos = buff_pos+buff_sz
1347    IF (check) THEN
1348      WRITE(*,*) "histdef : 3.2 bufpos for iv = ",iv, &
1349 &               " pfileid = ",pfileid," is = ",point(pfileid,iv)
1350    ENDIF
1351  ENDIF
1352!-
1353! 4.0 Transfer the frequency of the operations and check
1354!     for validity. We have to pay attention to negative values
1355!     of the frequency which indicate monthly time-steps.
1356!     The strategy is to bring it back to seconds for the tests
1357!-
1358  IF (check) WRITE(*,*) "histdef : 4.0"
1359!-
1360  freq_opp(pfileid,iv) = pfreq_opp
1361  freq_wrt(pfileid,iv) = pfreq_wrt
1362!-
1363  CALL ioget_calendar(un_an, un_jour)
1364  IF ( pfreq_opp < 0) THEN
1365    CALL ioget_calendar(un_an)
1366    test_fopp = pfreq_opp*(-1.)*un_an/12.*un_jour
1367  ELSE
1368    test_fopp = pfreq_opp
1369  ENDIF
1370  IF ( pfreq_wrt < 0) THEN
1371    CALL ioget_calendar(un_an)
1372    test_fwrt = pfreq_wrt*(-1.)*un_an/12.*un_jour
1373  ELSE
1374    test_fwrt = pfreq_wrt
1375  ENDIF
1376!-
1377! 4.1 Frequency of operations and output should be larger than deltat !
1378!-
1379  IF (test_fopp < deltat(pfileid)) THEN
1380    str70 = 'Frequency of operations should be larger than deltat'
1381    WRITE(str71,'("It is not the case for variable ",a," :",F10.4)') &
1382 &   TRIM(tmp_name),pfreq_opp
1383    str72 = "PATCH : frequency set to deltat"
1384!-
1385    CALL ipslerr (2,"histdef", str70, str71, str72)
1386!-
1387    freq_opp(pfileid,iv) = deltat(pfileid)
1388  ENDIF
1389!-
1390  IF (test_fwrt < deltat(pfileid)) THEN
1391    str70 = 'Frequency of output should be larger than deltat'
1392    WRITE(str71,'("It is not the case for variable ",a," :",F10.4)') &
1393 &   TRIM(tmp_name),pfreq_wrt
1394    str72 = "PATCH : frequency set to deltat"
1395!-
1396    CALL ipslerr (2,"histdef", str70, str71, str72)
1397!-
1398    freq_wrt(pfileid,iv) = deltat(pfileid)
1399  ENDIF
1400!-
1401! 4.2 First the existence of the operation is tested and then
1402!     its compaticility with the choice of frequencies
1403!-
1404  IF (TRIM(tmp_topp) == "inst") THEN
1405    IF (test_fopp /= test_fwrt) THEN
1406      str70 = 'For instantaneous output the frequency '// &
1407 &            'of operations and output'
1408      WRITE(str71, &
1409 &     '("should be the same, this was not case for variable ",a)') &
1410 &      TRIM(tmp_name)
1411      str72 = "PATCH : The smalest frequency of both is used"
1412      CALL ipslerr (2,"histdef", str70, str71, str72)
1413      IF ( test_fopp < test_fwrt) THEN
1414        freq_opp(pfileid,iv) = pfreq_opp
1415        freq_wrt(pfileid,iv) = pfreq_opp
1416      ELSE
1417        freq_opp(pfileid,iv) = pfreq_wrt
1418        freq_wrt(pfileid,iv) = pfreq_wrt
1419      ENDIF
1420    ENDIF
1421  ELSE IF (INDEX(ex_topps,TRIM(tmp_topp)) > 0) THEN
1422    IF (test_fopp > test_fwrt) THEN
1423      str70 = 'For averages the frequency of operations '// &
1424&             'should be smaller or equal'
1425      WRITE(str71, &
1426 &     '("to that of output. It is not the case for variable ",a)') &
1427 &     TRIM(tmp_name)
1428      str72 = 'PATCH : The output frequency is used for both'
1429      CALL ipslerr (2,"histdef", str70, str71, str72)
1430      freq_opp(pfileid,iv) = pfreq_wrt
1431    ENDIF
1432  ELSE
1433    WRITE (str70,'("Operation on variable ",a," is unknown")') &
1434&    TRIM(tmp_name)
1435    WRITE (str71, '("operation requested is :",a)') tmp_topp
1436    WRITE (str72, '("File ID :",I3)') pfileid
1437    CALL ipslerr (3,"histdef", str70, str71, str72)
1438  ENDIF
1439!-
1440! 5.0 Initialize other variables of the common
1441!-
1442  IF (check) WRITE(*,*) "histdef : 5.0"
1443!-
1444  hist_wrt_rng(pfileid,iv) = (PRESENT(var_range))
1445  IF (hist_wrt_rng(pfileid,iv)) THEN
1446    hist_calc_rng(pfileid,iv) = (var_range(1) > var_range(2))
1447    IF (hist_calc_rng(pfileid,iv)) THEN
1448      hist_minmax(pfileid,iv,1:2) = &
1449 &      (/ ABS(missing_val), -ABS(missing_val) /)
1450    ELSE
1451      hist_minmax(pfileid,iv,1:2) = var_range(1:2)
1452    ENDIF
1453  ENDIF
1454!-
1455! - freq_opp(pfileid,iv)/2./deltat(pfileid)
1456  last_opp(pfileid,iv) = itau0(pfileid)
1457! - freq_wrt(pfileid,iv)/2./deltat(pfileid)
1458  last_wrt(pfileid,iv) = itau0(pfileid)
1459! - freq_opp(pfileid,iv)/2./deltat(pfileid)
1460  last_opp_chk(pfileid,iv) = itau0(pfileid)
1461! - freq_wrt(pfileid,iv)/2./deltat(pfileid)
1462  last_wrt_chk(pfileid,iv) = itau0(pfileid)
1463  nb_opp(pfileid,iv) = 0
1464  nb_wrt(pfileid,iv) = 0
1465!-
1466! 6.0 Get the time axis for this variable
1467!-
1468  IF (check) WRITE(*,*) "histdef : 6.0"
1469!-
1470  IF ( freq_wrt(pfileid,iv) > 0 ) THEN
1471    WRITE(str10,'(I8.8)') INT(freq_wrt(pfileid,iv))
1472    str40 = TRIM(tmp_topp)//"_"//TRIM(str10)
1473  ELSE
1474    WRITE(str10,'(I2.2,"month")') ABS(INT(freq_wrt(pfileid,iv)))
1475    str40 = TRIM(tmp_topp)//"_"//TRIM(str10)
1476  ENDIF
1477!-
1478  DO i=1,nb_tax(pfileid)
1479    tab_str40(i) = tax_name(pfileid,i)
1480    tab_str40_length(i) = tax_name_length(pfileid,i)
1481  ENDDO
1482!-
1483  CALL find_str (nb_tax(pfileid),tab_str40,tab_str40_length,str40,pos)
1484!-
1485! No time axis for once, l_max, l_min or never operation
1486!-
1487  IF (     (TRIM(tmp_topp) /= 'once')  &
1488 &    .AND.(TRIM(tmp_topp) /= 'never') &
1489 &    .AND.(TRIM(tmp_topp) /= 'l_max') &
1490 &    .AND.(TRIM(tmp_topp) /= 'l_min') ) THEN
1491    IF ( pos < 0) THEN
1492      nb_tax(pfileid) = nb_tax(pfileid)+1
1493      tax_name(pfileid,nb_tax(pfileid)) = str40
1494      tax_name_length(pfileid, nb_tax(pfileid)) = LEN_TRIM(str40)
1495      tax_last(pfileid,nb_tax(pfileid)) = 0
1496      var_axid(pfileid,iv) = nb_tax(pfileid)
1497    ELSE
1498      var_axid(pfileid,iv) = pos
1499    ENDIF
1500  ELSE
1501    IF (check)   WRITE(*,*) "histdef : 7.0 ",TRIM(tmp_topp),'----'
1502    var_axid(pfileid,iv) = -99
1503  ENDIF
1504!-
1505! 7.0 prepare frequence of writing and operation
1506!     for never or once operation
1507!-
1508  IF (    (TRIM(tmp_topp) == 'once')  &
1509 &    .OR.(TRIM(tmp_topp) == 'never') ) THEN
1510    freq_opp(pfileid,iv) = 0.
1511    freq_wrt(pfileid,iv) = 0.
1512  ENDIF
1513!---------------------
1514END SUBROUTINE histdef
1515!===
1516SUBROUTINE histend (pfileid)
1517!---------------------------------------------------------------------
1518!- This subroutine end the decalaration of variables and sets the
1519!- time axes in the netcdf file and puts it into the write mode.
1520!-
1521!- INPUT
1522!-
1523!- pfileid : ID of the file to be worked on
1524!-
1525!- VERSION
1526!-
1527!---------------------------------------------------------------------
1528  IMPLICIT NONE
1529!-
1530  INTEGER, INTENT(IN) :: pfileid
1531!-
1532  INTEGER :: ncid, ncvarid
1533  INTEGER :: iret, ndim, iv, itx, ziv
1534  INTEGER :: itax
1535  INTEGER :: dims(4), dim_cnt
1536  INTEGER :: year, month, day, hours, minutes
1537  REAL :: sec
1538  REAL :: rtime0
1539  CHARACTER(LEN=20) :: tname, tunit
1540  CHARACTER(LEN=30) :: str30
1541  CHARACTER(LEN=80) :: ttitle
1542  CHARACTER(LEN=120) :: assoc
1543  CHARACTER(LEN=70) :: str70
1544  CHARACTER(LEN=3),DIMENSION(12) :: cal =   &
1545 &  (/ 'JAN','FEB','MAR','APR','MAY','JUN', &
1546 &     'JUL','AUG','SEP','OCT','NOV','DEC' /)
1547  CHARACTER(LEN=7) :: tmp_opp
1548!-
1549  LOGICAL :: check = .FALSE.
1550!---------------------------------------------------------------------
1551  ncid = ncdf_ids(pfileid)
1552!-
1553! 1.0 Create the time axes
1554!-
1555  IF (check) WRITE(*,*) "histend : 1.0"
1556!---
1557  iret = NF90_DEF_DIM (ncid,'time_counter',NF90_UNLIMITED,tid(pfileid))
1558!-
1559! 1.1 Define all the time axes needed for this file
1560!-
1561  DO itx=1,nb_tax(pfileid)
1562    dims(1) = tid(pfileid)
1563    IF (nb_tax(pfileid) > 1) THEN
1564      str30 = "t_"//tax_name(pfileid,itx)
1565    ELSE
1566      str30 = "time_counter"
1567    ENDIF
1568    iret = NF90_DEF_VAR (ncid,str30,NF90_FLOAT, &
1569 &                       dims(1),tdimid(pfileid,itx))
1570!---
1571!   To transform the current itau into a real date and take it
1572!   as the origin of the file requires the time counter to change.
1573!   Thus it is an operation the user has to ask for.
1574!   This function should thus only be re-instated
1575!   if there is a ioconf routine to control it.
1576!---
1577!-- rtime0 = itau2date(itau0(pfileid), date0(pfileid), deltat(pfileid))
1578    rtime0 = date0(pfileid)
1579!-
1580    CALL ju2ymds(rtime0, year, month, day, sec)
1581!---
1582!   Catch any error induced by a change in calendar !
1583!---
1584    IF (year < 0) THEN
1585      year = 2000+year
1586    ENDIF
1587!-
1588    hours = INT(sec/(60.*60.))
1589    minutes = INT((sec-hours*60.*60.)/60.)
1590    sec = sec-(hours*60.*60.+minutes*60.)
1591!-
1592    WRITE(str70,7000) year, month, day, hours, minutes, INT(sec)
1593    iret = NF90_PUT_ATT (ncid,tdimid(pfileid,itx),'units',TRIM(str70))
1594!-
1595    CALL ioget_calendar (str30)
1596    iret = NF90_PUT_ATT (ncid,tdimid(pfileid,itx), &
1597 &                       'calendar',TRIM(str30))
1598!-
1599    iret = NF90_PUT_ATT (ncid,tdimid(pfileid,itx),'title','Time')
1600!-
1601    iret = NF90_PUT_ATT (ncid,tdimid(pfileid,itx), &
1602 &                       'long_name','Time axis')
1603!-
1604    WRITE(str70,7001) year, cal(month), day, hours, minutes, INT(sec)
1605    iret = NF90_PUT_ATT (ncid,tdimid(pfileid,itx), &
1606 &                       'time_origin',TRIM(str70))
1607  ENDDO
1608!-
1609! The formats we need
1610!-
16117000 FORMAT('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2)
16127001 FORMAT(' ', I4.4,'-',A3,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2)
1613!-
1614! 2.0 declare the variables
1615!-
1616  IF (check) WRITE(*,*) "histend : 2.0"
1617!-
1618  DO iv=1,nb_var(pfileid)
1619!---
1620    itax = var_axid(pfileid,iv)
1621!---
1622    tname = name(pfileid,iv)
1623    tunit = unit_name(pfileid,iv)
1624    ttitle = title(pfileid,iv)
1625!---
1626    IF ( regular(pfileid) ) THEN
1627      dims(1:2) = (/ xid(pfileid), yid(pfileid) /)
1628      dim_cnt = 2
1629    ELSE
1630      dims(1) = xid(pfileid)
1631      dim_cnt = 1
1632    ENDIF
1633!---
1634    tmp_opp = topp(pfileid,iv)
1635    ziv = var_zaxid(pfileid,iv)
1636!---
1637!   2.1 dimension of field
1638!---
1639    IF ( (TRIM(tmp_opp) /= 'never')) THEN
1640      IF (     (TRIM(tmp_opp) /= 'once')  &
1641     &    .AND.(TRIM(tmp_opp) /= 'l_max') &
1642     &    .AND.(TRIM(tmp_opp) /= 'l_min') ) THEN
1643        IF (ziv == -99) THEN
1644          ndim = dim_cnt+1
1645          dims(dim_cnt+1:dim_cnt+2) = (/ tid(pfileid), 0 /)
1646        ELSE
1647          ndim = dim_cnt+2
1648          dims(dim_cnt+1:dim_cnt+2) = (/ zax_ids(pfileid,ziv), tid(pfileid) /)
1649        ENDIF
1650      ELSE
1651        IF (ziv == -99) THEN
1652          ndim = dim_cnt
1653          dims(dim_cnt+1:dim_cnt+2) = (/ 0, 0 /)
1654        ELSE
1655          ndim = dim_cnt+1
1656          dims(dim_cnt+1:dim_cnt+2) = (/ zax_ids(pfileid,ziv), 0 /)
1657        ENDIF
1658      ENDIF
1659!-
1660      iret = NF90_DEF_VAR (ncid,TRIM(tname),NF90_FLOAT, &
1661 &                         dims(1:ABS(ndim)),ncvarid)
1662!-
1663      ncvar_ids(pfileid,iv) = ncvarid
1664!-
1665      iret = NF90_PUT_ATT (ncid,ncvarid,'units',TRIM(tunit))
1666!-
1667      iret = NF90_PUT_ATT (ncid,ncvarid,'missing_value', &
1668 &                         REAL(missing_val,KIND=4))
1669      IF (hist_wrt_rng(pfileid,iv)) THEN
1670        iret = NF90_PUT_ATT (ncid,ncvarid,'valid_min', &
1671 &                           REAL(hist_minmax(pfileid,iv,1),KIND=4))
1672        iret = NF90_PUT_ATT (ncid,ncvarid,'valid_max', &
1673 &                           REAL(hist_minmax(pfileid,iv,2),KIND=4))
1674      ENDIF
1675!-
1676      iret = NF90_PUT_ATT (ncid,ncvarid,'long_name',TRIM(ttitle))
1677!-
1678      iret = NF90_PUT_ATT (ncid,ncvarid,'short_name',TRIM(tname))
1679!-
1680      iret = NF90_PUT_ATT (ncid,ncvarid,'online_operation', &
1681 &                         TRIM(fullop(pfileid,iv)))
1682!-
1683      SELECT CASE(ndim)
1684      CASE(-3)
1685        str30 = 'ZYX'
1686      CASE(2)
1687        str30 = 'YX'
1688      CASE(3)
1689        str30 = 'TYX'
1690      CASE(4)
1691        str30 = 'TZYX'
1692      CASE DEFAULT
1693        CALL ipslerr (3,"histend", &
1694       &  'less than 2 or more than 4 dimensions are not', &
1695       &  'allowed at this stage',' ')
1696      END SELECT
1697!-
1698      iret = NF90_PUT_ATT (ncid,ncvarid,'axis',TRIM(str30))
1699!-
1700      assoc='nav_lat nav_lon'
1701      ziv = var_zaxid(pfileid, iv)
1702      IF (ziv > 0) THEN
1703        str30 = zax_name(pfileid,ziv)
1704        assoc = TRIM(str30)//' '//TRIM(assoc)
1705      ENDIF
1706!-
1707      IF (itax > 0) THEN
1708        IF ( nb_tax(pfileid) > 1) THEN
1709          str30 = "t_"//tax_name(pfileid,itax)
1710        ELSE
1711          str30 = "time_counter"
1712        ENDIF
1713        assoc = TRIM(str30)//' '//TRIM(assoc)
1714!-
1715        IF (check) THEN
1716          WRITE(*,*) "histend : 2.0.n, freq_opp, freq_wrt", &
1717 &                   freq_opp(pfileid,iv), freq_wrt(pfileid,iv)
1718        ENDIF
1719!-
1720        iret = NF90_PUT_ATT (ncid,ncvarid,'interval_operation', &
1721 &                           REAL(freq_opp(pfileid,iv),KIND=4))
1722        iret = NF90_PUT_ATT (ncid,ncvarid,'interval_write', &
1723 &                           REAL(freq_wrt(pfileid,iv),KIND=4))
1724      ENDIF
1725      iret = NF90_PUT_ATT (ncid,ncvarid,'associate',TRIM(assoc))
1726    ENDIF
1727  ENDDO
1728!-
1729! 3.0 Put the netcdf file into wrte mode
1730!-
1731  IF (check) WRITE(*,*) "histend : 3.0"
1732!-
1733  iret = NF90_ENDDEF (ncid)
1734!-
1735! 4.0 Give some informations to the user
1736!-
1737  IF (check) WRITE(*,*) "histend : 4.0"
1738!-
1739  WRITE(str70,'("All variables have been initialized on file :",I3)') pfileid
1740  CALL ipslerr (1,'histend',str70,'',' ')
1741!---------------------
1742END SUBROUTINE histend
1743!===
1744SUBROUTINE histwrite_r1d (pfileid,pvarname,pitau,pdata,nbindex,nindex)
1745!---------------------------------------------------------------------
1746  IMPLICIT NONE
1747!-
1748  INTEGER,INTENT(IN) :: pfileid, pitau, nbindex, nindex(nbindex)
1749  REAL,DIMENSION(:),INTENT(IN) :: pdata
1750  CHARACTER(LEN=*),INTENT(IN) :: pvarname
1751!-
1752  LOGICAL :: do_oper, do_write, largebuf
1753  INTEGER :: varid, io, nbpt_in, nbpt_out
1754  REAL, ALLOCATABLE,SAVE :: buff_tmp(:)
1755  INTEGER,SAVE :: buff_tmp_sz
1756  CHARACTER(LEN=7) :: tmp_opp
1757!-
1758  LOGICAL :: check = .FALSE.
1759!---------------------------------------------------------------------
1760!-
1761! 1.0 Try to catch errors like specifying the wrong file ID.
1762!     Thanks Marine for showing us what errors users can make !
1763!-
1764  IF ( (pfileid < 1).OR.(pfileid > nb_files) ) THEN
1765    CALL ipslerr (3,"histwrite", &
1766 &    'Illegal file ID in the histwrite of variable',pvarname,' ')
1767  ENDIF
1768!-
1769! 1.1 Find the id of the variable to be written and the real time
1770!-
1771  CALL histvar_seq (pfileid,pvarname,varid)
1772!-
1773! 2.0 do nothing for never operation
1774!-
1775  tmp_opp = topp(pfileid,varid)
1776!-
1777  IF (TRIM(tmp_opp) == "never") THEN
1778    last_opp_chk(pfileid,varid) = -99
1779    last_wrt_chk(pfileid,varid) = -99
1780  ENDIF
1781!-
1782! 3.0 We check if we need to do an operation
1783!-
1784  IF (last_opp_chk(pfileid,varid) == pitau) THEN
1785    CALL ipslerr (3,"histwrite", &
1786 &    'This variable as already been analysed at the present', &
1787 &    'time step',' ')
1788  ENDIF
1789!-
1790  CALL isittime &
1791 &  (pitau, date0(pfileid), deltat(pfileid), freq_opp(pfileid,varid), &
1792 &   last_opp(pfileid,varid), last_opp_chk(pfileid,varid), do_oper)
1793!-
1794! 4.0 We check if we need to write the data
1795!-
1796  IF (last_wrt_chk(pfileid,varid) == pitau) THEN
1797    CALL ipslerr (3,"histwrite", &
1798 &    'This variable as already been written for the present', &
1799 &    'time step',' ')
1800  ENDIF
1801!-
1802  CALL isittime &
1803 &  (pitau, date0(pfileid), deltat(pfileid), freq_wrt(pfileid,varid), &
1804 &   last_wrt(pfileid,varid), last_wrt_chk(pfileid,varid), do_write)
1805!-
1806! 5.0 histwrite called
1807!-
1808  IF (do_oper.OR.do_write) THEN
1809!-
1810!-- 5.1 Get the sizes of the data we will handle
1811!-
1812    IF (datasz_in(pfileid,varid,1) <= 0) THEN
1813!---- There is the risk here that the user has over-sized the array.
1814!---- But how can we catch this ?
1815!---- In the worst case we will do impossible operations
1816!---- on part of the data !
1817      datasz_in(pfileid,varid,1) = SIZE(pdata)
1818      datasz_in(pfileid,varid,2) = -1
1819      datasz_in(pfileid,varid,3) = -1
1820    ENDIF
1821!-
1822!-- 5.2 The maximum size of the data will give the size of the buffer
1823!-
1824    IF (datasz_max(pfileid,varid) <= 0) THEN
1825      largebuf = .FALSE.
1826      DO io=1,nbopp(pfileid,varid)
1827        IF (INDEX(fuchnbout,sopps(pfileid,varid,io)) > 0) THEN
1828          largebuf = .TRUE.
1829        ENDIF
1830      ENDDO
1831      IF (largebuf) THEN
1832        datasz_max(pfileid,varid) = &
1833 &        scsize(pfileid,varid,1) &
1834 &       *scsize(pfileid,varid,2) &
1835 &       *scsize(pfileid,varid,3)
1836      ELSE
1837        datasz_max(pfileid,varid) = &
1838 &        datasz_in(pfileid,varid,1)
1839      ENDIF
1840    ENDIF
1841!-
1842    IF (.NOT.ALLOCATED(buff_tmp)) THEN
1843      IF (check) THEN
1844        WRITE(*,*) &
1845 &       "histwrite_r1d : allocate buff_tmp for buff_sz = ", &
1846 &       datasz_max(pfileid,varid)
1847      ENDIF
1848      ALLOCATE (buff_tmp(datasz_max(pfileid,varid)))
1849      buff_tmp_sz = datasz_max(pfileid,varid)
1850    ELSE IF (datasz_max(pfileid,varid) > buff_tmp_sz) THEN
1851      IF (check) THEN
1852        WRITE(*,*) &
1853 &       "histwrite_r1d : re-allocate buff_tmp for buff_sz = ", &
1854 &       datasz_max(pfileid,varid)
1855      ENDIF
1856      DEALLOCATE (buff_tmp)
1857      ALLOCATE (buff_tmp(datasz_max(pfileid,varid)))
1858      buff_tmp_sz = datasz_max(pfileid,varid)
1859    ENDIF
1860!-
1861!-- We have to do the first operation anyway.
1862!-- Thus we do it here and change the ranke
1863!-- of the data at the same time. This should speed up things.
1864!-
1865    nbpt_in = datasz_in(pfileid,varid,1)
1866    nbpt_out = datasz_max(pfileid,varid)
1867    CALL mathop (sopps(pfileid,varid,1), nbpt_in, pdata, &
1868 &               missing_val, nbindex, nindex, &
1869 &               scal(pfileid,varid,1), nbpt_out, buff_tmp)
1870    CALL histwrite_real (pfileid, varid, pitau, nbpt_out, &
1871 &            buff_tmp, nbindex, nindex, do_oper, do_write)
1872  ENDIF
1873!-
1874! 6.0 Manage time steps
1875!-
1876  IF ((TRIM(tmp_opp) /= "once").AND.(TRIM(tmp_opp) /= "never")) THEN
1877    last_opp_chk(pfileid,varid) = pitau
1878    last_wrt_chk(pfileid,varid) = pitau
1879  ELSE
1880    last_opp_chk(pfileid,varid) = -99
1881    last_wrt_chk(pfileid,varid) = -99
1882  ENDIF
1883!---------------------------
1884END SUBROUTINE histwrite_r1d
1885!===
1886SUBROUTINE histwrite_r2d (pfileid,pvarname,pitau,pdata,nbindex,nindex)
1887!---------------------------------------------------------------------
1888  IMPLICIT NONE
1889!-
1890  INTEGER,INTENT(IN) :: pfileid, pitau, nbindex, nindex(nbindex)
1891  REAL,DIMENSION(:,:),INTENT(IN) :: pdata
1892  CHARACTER(LEN=*),INTENT(IN) :: pvarname
1893!-
1894  LOGICAL :: do_oper, do_write, largebuf
1895  INTEGER :: varid, io, nbpt_in(1:2), nbpt_out
1896  REAL, ALLOCATABLE,SAVE :: buff_tmp(:)
1897  INTEGER,SAVE :: buff_tmp_sz
1898  CHARACTER(LEN=7) :: tmp_opp
1899!-
1900  LOGICAL :: check = .FALSE.
1901!---------------------------------------------------------------------
1902!-
1903! 1.0 Try to catch errors like specifying the wrong file ID.
1904!     Thanks Marine for showing us what errors users can make !
1905!-
1906  IF ( (pfileid < 1).OR.(pfileid > nb_files) ) THEN
1907    CALL ipslerr (3,"histwrite", &
1908 &    'Illegal file ID in the histwrite of variable',pvarname,' ')
1909  ENDIF
1910!-
1911! 1.1 Find the id of the variable to be written and the real time
1912!-
1913  CALL histvar_seq (pfileid,pvarname,varid)
1914!-
1915! 2.0 do nothing for never operation
1916!-
1917  tmp_opp = topp(pfileid,varid)
1918!-
1919  IF (TRIM(tmp_opp) == "never") THEN
1920    last_opp_chk(pfileid,varid) = -99
1921    last_wrt_chk(pfileid,varid) = -99
1922  ENDIF
1923!-
1924! 3.0 We check if we need to do an operation
1925!-
1926  IF (last_opp_chk(pfileid,varid) == pitau) THEN
1927    CALL ipslerr (3,"histwrite", &
1928 &    'This variable as already been analysed at the present', &
1929 &    'time step',' ')
1930  ENDIF
1931!-
1932  CALL isittime &
1933 &  (pitau, date0(pfileid), deltat(pfileid), freq_opp(pfileid,varid), &
1934 &   last_opp(pfileid,varid), last_opp_chk(pfileid,varid), do_oper)
1935!-
1936! 4.0 We check if we need to write the data
1937!-
1938  IF (last_wrt_chk(pfileid,varid) == pitau) THEN
1939    CALL ipslerr (3,"histwrite", &
1940 &    'This variable as already been written for the present', &
1941 &    'time step',' ')
1942  ENDIF
1943!-
1944  CALL isittime &
1945 &  (pitau, date0(pfileid), deltat(pfileid), freq_wrt(pfileid,varid), &
1946 &   last_wrt(pfileid,varid), last_wrt_chk(pfileid,varid), do_write)
1947!-
1948! 5.0 histwrite called
1949!-
1950  IF (do_oper.OR.do_write) THEN
1951!-
1952!-- 5.1 Get the sizes of the data we will handle
1953!-
1954    IF (datasz_in(pfileid,varid,1) <= 0) THEN
1955!---- There is the risk here that the user has over-sized the array.
1956!---- But how can we catch this ?
1957!---- In the worst case we will do impossible operations
1958!---- on part of the data !
1959      datasz_in(pfileid,varid,1) = SIZE(pdata, DIM=1)
1960      datasz_in(pfileid,varid,2) = SIZE(pdata, DIM=2)
1961      datasz_in(pfileid,varid,3) = -1
1962    ENDIF
1963!-
1964!-- 5.2 The maximum size of the data will give the size of the buffer
1965!-
1966    IF (datasz_max(pfileid,varid) <= 0) THEN
1967      largebuf = .FALSE.
1968      DO io=1,nbopp(pfileid,varid)
1969        IF (INDEX(fuchnbout,sopps(pfileid,varid,io)) > 0) THEN
1970          largebuf = .TRUE.
1971        ENDIF
1972      ENDDO
1973      IF (largebuf) THEN
1974        datasz_max(pfileid,varid) = &
1975 &        scsize(pfileid,varid,1) &
1976 &       *scsize(pfileid,varid,2) &
1977 &       *scsize(pfileid,varid,3)
1978      ELSE
1979        datasz_max(pfileid,varid) = &
1980 &        datasz_in(pfileid,varid,1) &
1981 &       *datasz_in(pfileid,varid,2)
1982      ENDIF
1983    ENDIF
1984!-
1985    IF (.NOT.ALLOCATED(buff_tmp)) THEN
1986      IF (check) THEN
1987        WRITE(*,*) &
1988 &       "histwrite_r2d : allocate buff_tmp for buff_sz = ", &
1989 &       datasz_max(pfileid,varid)
1990      ENDIF
1991      ALLOCATE (buff_tmp(datasz_max(pfileid,varid)))
1992      buff_tmp_sz = datasz_max(pfileid,varid)
1993    ELSE IF (datasz_max(pfileid,varid) > buff_tmp_sz) THEN
1994      IF (check) THEN
1995        WRITE(*,*) &
1996 &       "histwrite_r2d : re-allocate buff_tmp for buff_sz = ", &
1997 &       datasz_max(pfileid,varid)
1998      ENDIF
1999      DEALLOCATE (buff_tmp)
2000      ALLOCATE (buff_tmp(datasz_max(pfileid,varid)))
2001      buff_tmp_sz = datasz_max(pfileid,varid)
2002    ENDIF
2003!-
2004!-- We have to do the first operation anyway.
2005!-- Thus we do it here and change the ranke
2006!-- of the data at the same time. This should speed up things.
2007!-
2008    nbpt_in(1:2) = datasz_in(pfileid,varid,1:2)
2009    nbpt_out = datasz_max(pfileid,varid)
2010    CALL mathop (sopps(pfileid,varid,1), nbpt_in, pdata, &
2011 &               missing_val, nbindex, nindex, &
2012 &               scal(pfileid,varid,1), nbpt_out, buff_tmp)
2013    CALL histwrite_real (pfileid, varid, pitau, nbpt_out, &
2014 &            buff_tmp, nbindex, nindex, do_oper, do_write)
2015  ENDIF
2016!-
2017! 6.0 Manage time steps
2018!-
2019  IF ((TRIM(tmp_opp) /= "once").AND.(TRIM(tmp_opp) /= "never")) THEN
2020    last_opp_chk(pfileid,varid) = pitau
2021    last_wrt_chk(pfileid,varid) = pitau
2022  ELSE
2023    last_opp_chk(pfileid,varid) = -99
2024    last_wrt_chk(pfileid,varid) = -99
2025  ENDIF
2026!---------------------------
2027END SUBROUTINE histwrite_r2d
2028!===
2029SUBROUTINE histwrite_r3d (pfileid,pvarname,pitau,pdata,nbindex,nindex)
2030!---------------------------------------------------------------------
2031  IMPLICIT NONE
2032!-
2033  INTEGER,INTENT(IN) :: pfileid, pitau, nbindex, nindex(nbindex)
2034  REAL,DIMENSION(:,:,:),INTENT(IN) :: pdata
2035  CHARACTER(LEN=*),INTENT(IN) :: pvarname
2036!-
2037  LOGICAL :: do_oper, do_write, largebuf
2038  INTEGER :: varid, io, nbpt_in(1:3), nbpt_out
2039  REAL, ALLOCATABLE,SAVE :: buff_tmp(:)
2040  INTEGER,SAVE :: buff_tmp_sz
2041  CHARACTER(LEN=7) :: tmp_opp
2042!-
2043  LOGICAL :: check = .FALSE.
2044!---------------------------------------------------------------------
2045!-
2046! 1.0 Try to catch errors like specifying the wrong file ID.
2047!     Thanks Marine for showing us what errors users can make !
2048!-
2049  IF ( (pfileid < 1).OR.(pfileid > nb_files) ) THEN
2050    CALL ipslerr (3,"histwrite", &
2051 &    'Illegal file ID in the histwrite of variable',pvarname,' ')
2052  ENDIF
2053!-
2054! 1.1 Find the id of the variable to be written and the real time
2055!-
2056  CALL histvar_seq (pfileid,pvarname,varid)
2057!-
2058! 2.0 do nothing for never operation
2059!-
2060  tmp_opp = topp(pfileid,varid)
2061!-
2062  IF (TRIM(tmp_opp) == "never") THEN
2063    last_opp_chk(pfileid,varid) = -99
2064    last_wrt_chk(pfileid,varid) = -99
2065  ENDIF
2066!-
2067! 3.0 We check if we need to do an operation
2068!-
2069  IF (last_opp_chk(pfileid,varid) == pitau) THEN
2070    CALL ipslerr (3,"histwrite", &
2071 &    'This variable as already been analysed at the present', &
2072 &    'time step',' ')
2073  ENDIF
2074!-
2075  CALL isittime &
2076 &  (pitau, date0(pfileid), deltat(pfileid), freq_opp(pfileid,varid), &
2077 &   last_opp(pfileid,varid), last_opp_chk(pfileid,varid), do_oper)
2078!-
2079! 4.0 We check if we need to write the data
2080!-
2081  IF (last_wrt_chk(pfileid,varid) == pitau) THEN
2082    CALL ipslerr (3,"histwrite", &
2083 &    'This variable as already been written for the present', &
2084 &    'time step',' ')
2085  ENDIF
2086!-
2087  CALL isittime &
2088 &  (pitau, date0(pfileid), deltat(pfileid), freq_wrt(pfileid,varid), &
2089 &   last_wrt(pfileid,varid), last_wrt_chk(pfileid,varid), do_write)
2090!-
2091! 5.0 histwrite called
2092!-
2093  IF (do_oper.OR.do_write) THEN
2094!-
2095!-- 5.1 Get the sizes of the data we will handle
2096!-
2097    IF (datasz_in(pfileid,varid,1) <= 0) THEN
2098!---- There is the risk here that the user has over-sized the array.
2099!---- But how can we catch this ?
2100!---- In the worst case we will do impossible operations
2101!---- on part of the data !
2102      datasz_in(pfileid,varid,1) = SIZE(pdata, DIM=1)
2103      datasz_in(pfileid,varid,2) = SIZE(pdata, DIM=2)
2104      datasz_in(pfileid,varid,3) = SIZE(pdata, DIM=3)
2105    ENDIF
2106!-
2107!-- 5.2 The maximum size of the data will give the size of the buffer
2108!-
2109    IF (datasz_max(pfileid,varid) <= 0) THEN
2110      largebuf = .FALSE.
2111      DO io =1,nbopp(pfileid,varid)
2112        IF (INDEX(fuchnbout,sopps(pfileid,varid,io)) > 0) THEN
2113          largebuf = .TRUE.
2114        ENDIF
2115      ENDDO
2116      IF (largebuf) THEN
2117        datasz_max(pfileid,varid) = &
2118 &        scsize(pfileid,varid,1) &
2119 &       *scsize(pfileid,varid,2) &
2120 &       *scsize(pfileid,varid,3)
2121      ELSE
2122        datasz_max(pfileid,varid) = &
2123 &        datasz_in(pfileid,varid,1) &
2124 &       *datasz_in(pfileid,varid,2) &
2125 &       *datasz_in(pfileid,varid,3)
2126      ENDIF
2127    ENDIF
2128!-
2129    IF (.NOT.ALLOCATED(buff_tmp)) THEN
2130      IF (check) THEN
2131        WRITE(*,*) &
2132 &       "histwrite_r1d : allocate buff_tmp for buff_sz = ", &
2133 &       datasz_max(pfileid,varid)
2134      ENDIF
2135      ALLOCATE (buff_tmp(datasz_max(pfileid,varid)))
2136      buff_tmp_sz = datasz_max(pfileid,varid)
2137    ELSE IF (datasz_max(pfileid,varid) > buff_tmp_sz) THEN
2138      IF (check) THEN
2139        WRITE(*,*) &
2140 &       "histwrite_r1d : re-allocate buff_tmp for buff_sz = ", &
2141 &       datasz_max(pfileid,varid)
2142      ENDIF
2143      DEALLOCATE (buff_tmp)
2144      ALLOCATE (buff_tmp(datasz_max(pfileid,varid)))
2145      buff_tmp_sz = datasz_max(pfileid,varid)
2146    ENDIF
2147!-
2148!-- We have to do the first operation anyway.
2149!-- Thus we do it here and change the ranke
2150!-- of the data at the same time. This should speed up things.
2151!-
2152    nbpt_in(1:3) = datasz_in(pfileid,varid,1:3)
2153    nbpt_out = datasz_max(pfileid,varid)
2154    CALL mathop (sopps(pfileid,varid,1), nbpt_in, pdata, &
2155 &               missing_val, nbindex, nindex, &
2156 &               scal(pfileid,varid,1), nbpt_out, buff_tmp)
2157    CALL histwrite_real (pfileid, varid, pitau, nbpt_out, &
2158 &            buff_tmp, nbindex, nindex, do_oper, do_write)
2159  ENDIF
2160!-
2161! 6.0 Manage time steps
2162!-
2163  IF ((TRIM(tmp_opp) /= "once").AND.(TRIM(tmp_opp) /= "never")) THEN
2164    last_opp_chk(pfileid,varid) = pitau
2165    last_wrt_chk(pfileid,varid) = pitau
2166  ELSE
2167    last_opp_chk(pfileid,varid) = -99
2168    last_wrt_chk(pfileid,varid) = -99
2169  ENDIF
2170!---------------------------
2171END SUBROUTINE histwrite_r3d
2172!===
2173SUBROUTINE histwrite_real &
2174 & (pfileid,varid,pitau,nbdpt,buff_tmp,nbindex,nindex,do_oper,do_write)
2175!---------------------------------------------------------------------
2176!- This subroutine is internal and does the calculations and writing
2177!- if needed. At a later stage it should be split into an operation
2178!- and writing subroutines.
2179!---------------------------------------------------------------------
2180  IMPLICIT NONE
2181!-
2182  INTEGER,INTENT(IN) :: pfileid,pitau,varid, &
2183 &                      nbindex,nindex(nbindex),nbdpt
2184  REAL,DIMENSION(:)  :: buff_tmp
2185  LOGICAL,INTENT(IN) :: do_oper,do_write
2186!-
2187  INTEGER :: tsz, ncid, ncvarid
2188  INTEGER :: i, iret, ipt, itax
2189  INTEGER :: io, nbin, nbout
2190  INTEGER,DIMENSION(4) :: corner, edges
2191  INTEGER :: itime
2192!-
2193  REAL :: rtime
2194  CHARACTER(LEN=7) :: tmp_opp
2195  REAL,ALLOCATABLE,SAVE :: buff_tmp2(:)
2196  INTEGER,SAVE          :: buff_tmp2_sz
2197  REAL,ALLOCATABLE,SAVE :: buffer_used(:)
2198  INTEGER,SAVE          :: buffer_sz
2199!-
2200  LOGICAL :: check = .FALSE.
2201!---------------------------------------------------------------------
2202  IF (check) THEN
2203    WRITE(*,*) "histwrite 0.0 :  VAR : ", name(pfileid,varid)
2204    WRITE(*,*) "histwrite 0.0 : nbindex, nindex :", &
2205    &  nbindex,nindex(1:MIN(3,nbindex)),'...',nindex(MAX(1,nbindex-3):nbindex)
2206  ENDIF
2207!-
2208! The sizes which can be encoutered
2209!-
2210  tsz = zsize(pfileid,varid,1)*zsize(pfileid,varid,2)*zsize(pfileid,varid,3)
2211!-
2212! 1.0 We allocate the memory needed to store the data between write
2213!     and the temporary space needed for operations.
2214!     We have to keep precedent buffer if needed
2215!-
2216  IF (.NOT. ALLOCATED(buffer)) THEN
2217    IF (check) WRITE(*,*) "histwrite_real 1.0 allocate buffer ",buff_pos
2218    ALLOCATE(buffer(buff_pos))
2219    buffer_sz = buff_pos
2220    buffer(:)=0.0
2221  ELSE IF (buffer_sz < buff_pos) THEN
2222    IF (check) WRITE(*,*) "histwrite_real 1.0.1 re-allocate buffer for ",buff_pos," instead of ",SIZE(buffer)
2223    IF (SUM(buffer)/=0.0) THEN
2224      IF (check) WRITE (*,*) ' histwrite : buffer has been used. We have to save it before re-allocating it '
2225      ALLOCATE (buffer_used(buffer_sz))
2226      buffer_used(:)=buffer(:)
2227      DEALLOCATE (buffer)
2228      ALLOCATE (buffer(buff_pos))
2229      buffer_sz = buff_pos
2230      buffer(:SIZE(buffer_used))=buffer_used
2231      DEALLOCATE (buffer_used)
2232    ELSE
2233      IF (check) WRITE (*,*) ' histwrite : buffer has not been used. We have just to re-allocating it '
2234      DEALLOCATE (buffer)
2235      ALLOCATE (buffer(buff_pos))
2236      buffer_sz = buff_pos
2237      buffer(:)=0.0
2238    ENDIF
2239  ENDIF
2240!-
2241! The buffers are only deallocated when more space is needed. This
2242! reduces the umber of allocates but increases memory needs.
2243!-
2244  IF (.NOT.ALLOCATED(buff_tmp2)) THEN
2245    IF (check) THEN
2246      WRITE(*,*) "histwrite_real 1.1 allocate buff_tmp2 ",SIZE(buff_tmp)
2247    ENDIF
2248    ALLOCATE (buff_tmp2(datasz_max(pfileid,varid)))
2249    buff_tmp2_sz = datasz_max(pfileid,varid)
2250  ELSE IF ( datasz_max(pfileid,varid) > buff_tmp2_sz) THEN
2251    IF (check) THEN
2252      WRITE(*,*) "histwrite_real 1.2 re-allocate buff_tmp2 : ", &
2253     & SIZE(buff_tmp)," instead of ",SIZE(buff_tmp2)
2254    ENDIF
2255    DEALLOCATE (buff_tmp2)
2256    ALLOCATE (buff_tmp2(datasz_max(pfileid,varid)))
2257    buff_tmp2_sz = datasz_max(pfileid,varid)
2258  ENDIF
2259!-
2260  rtime = pitau * deltat(pfileid)
2261  tmp_opp = topp(pfileid,varid)
2262!-
2263! 3.0 Do the operations or transfer the slab of data into buff_tmp
2264!-
2265  IF (check) WRITE(*,*) "histwrite: 3.0", pfileid
2266!-
2267! 3.1 DO the Operations only if needed
2268!-
2269  IF ( do_oper ) THEN
2270    i = pfileid
2271    nbout = nbdpt
2272!-
2273!-- 3.4 We continue the sequence of operations
2274!--     we started in the interface routine
2275!-
2276    DO io = 2, nbopp(i,varid),2
2277      nbin = nbout
2278      nbout = datasz_max(i,varid)
2279      CALL mathop(sopps(i,varid,io),nbin,buff_tmp,missing_val, &
2280 &      nbindex,nindex,scal(i,varid,io),nbout,buff_tmp2)
2281      IF (check) THEN
2282        WRITE(*,*) &
2283 &       "histwrite: 3.4a nbout : ",nbin,nbout,sopps(i,varid,io)
2284      ENDIF
2285!-
2286      nbin = nbout
2287      nbout = datasz_max(i,varid)
2288      CALL mathop(sopps(i,varid,io+1),nbin,buff_tmp2,missing_val, &
2289 &      nbindex,nindex,scal(i,varid,io+1),nbout,buff_tmp)
2290      IF (check) THEN
2291        WRITE(*,*) &
2292 &       "histwrite: 3.4b nbout : ",nbin,nbout,sopps(i,varid,io+1)
2293      ENDIF
2294    ENDDO
2295!-
2296!   3.5 Zoom into the data
2297!-
2298    IF (check) THEN
2299      WRITE(*,*) &
2300 &     "histwrite: 3.5 size(buff_tmp) : ",SIZE(buff_tmp)
2301      WRITE(*,*) &
2302 &     "histwrite: 3.5 slab in X :",zorig(i,varid,1),zsize(i,varid,1)
2303      WRITE(*,*) &
2304 &     "histwrite: 3.5 slab in Y :",zorig(i,varid,2),zsize(i,varid,2)
2305      WRITE(*,*) &
2306 &     "histwrite: 3.5 slab in Z :",zorig(i,varid,3),zsize(i,varid,3)
2307      WRITE(*,*) &
2308 &     "histwrite: 3.5 slab of input:", &
2309 &     scsize(i,varid,1),scsize(i,varid,2),scsize(i,varid,3)
2310    ENDIF
2311    CALL trans_buff &
2312 &      (zorig(i,varid,1),zsize(i,varid,1), &
2313 &       zorig(i,varid,2),zsize(i,varid,2), &
2314 &       zorig(i,varid,3),zsize(i,varid,3), &
2315 &       scsize(i,varid,1),scsize(i,varid,2),scsize(i,varid,3), &
2316 &       buff_tmp, buff_tmp2_sz,buff_tmp2)
2317!-
2318!-- 4.0 Get the min and max of the field (buff_tmp)
2319!-
2320    IF (check) WRITE(*,*) "histwrite: 4.0 buff_tmp",pfileid,varid, &
2321 &    TRIM(tmp_opp),'----',LEN_TRIM(tmp_opp),nbindex
2322!-
2323    IF (hist_calc_rng(pfileid,varid)) THEN
2324      hist_minmax(pfileid,varid,1) = &
2325 &      MIN(hist_minmax(pfileid,varid,1), &
2326 &      MINVAL(buff_tmp2(1:tsz),MASK=buff_tmp2(1:tsz) /= missing_val))
2327      hist_minmax(pfileid,varid,2) = &
2328 &      MAX(hist_minmax(pfileid,varid,2), &
2329 &      MAXVAL(buff_tmp2(1:tsz),MASK=buff_tmp2(1:tsz) /= missing_val))
2330    ENDIF
2331!-
2332!-- 5.0 Do the operations if needed. In the case of instantaneous
2333!--     output we do not transfer to the buffer.
2334!-
2335    IF (check) WRITE(*,*) "histwrite: 5.0", pfileid, "tsz :", tsz
2336!-
2337    ipt = point(pfileid,varid)
2338!-
2339!   WRITE(*,*) 'OPE ipt, buffer :', pvarname, ipt, varid
2340!-
2341    IF (     (TRIM(tmp_opp) /= "inst") &
2342   &    .AND.(TRIM(tmp_opp) /= "once") ) THEN
2343      CALL moycum(tmp_opp,tsz,buffer(ipt:), &
2344     &       buff_tmp2,nb_opp(pfileid,varid))
2345    ENDIF
2346!-
2347    last_opp(pfileid,varid) = pitau
2348    nb_opp(pfileid,varid) = nb_opp(pfileid,varid)+1
2349!-
2350   ENDIF
2351!-
2352! 6.0 Write to file if needed
2353!-
2354  IF (check) WRITE(*,*) "histwrite: 6.0", pfileid
2355!-
2356  IF ( do_write ) THEN
2357!-
2358    ncvarid = ncvar_ids(pfileid,varid)
2359    ncid = ncdf_ids(pfileid)
2360!-
2361!-- 6.1 Do the operations that are needed before writting
2362!-
2363    IF (check) WRITE(*,*) "histwrite: 6.1", pfileid
2364!-
2365    IF (     (TRIM(tmp_opp) /= "inst") &
2366   &    .AND.(TRIM(tmp_opp) /= "once") ) THEN
2367      rtime = (rtime+last_wrt(pfileid,varid)*deltat(pfileid))/2.0
2368    ENDIF
2369!-
2370!-- 6.2 Add a value to the time axis of this variable if needed
2371!-
2372    IF (     (TRIM(tmp_opp) /= "l_max") &
2373   &    .AND.(TRIM(tmp_opp) /= "l_min") &
2374   &    .AND.(TRIM(tmp_opp) /= "once") ) THEN
2375!-
2376      IF (check) WRITE(*,*) "histwrite: 6.2", pfileid
2377!-
2378      itax = var_axid(pfileid, varid)
2379      itime = nb_wrt(pfileid, varid)+1
2380!-
2381      IF (tax_last(pfileid, itax) < itime) THEN
2382        iret = NF90_PUT_VAR (ncid,tdimid(pfileid,itax),(/ rtime /), &
2383 &                            start=(/ itime /),count=(/ 1 /))
2384        tax_last(pfileid, itax) = itime
2385      ENDIF
2386    ELSE
2387      itime=1
2388    ENDIF
2389!-
2390!-- 6.3 Write the data. Only in the case of instantaneous output
2391!       we do not write the buffer.
2392!-
2393    IF (check) THEN
2394      WRITE(*,*) "histwrite: 6.3",pfileid,ncid,ncvarid,varid,itime
2395    ENDIF
2396!-
2397    IF (scsize(pfileid,varid,3) == 1) THEN
2398      IF (regular(pfileid)) THEN
2399        corner(1:4) = (/ 1, 1, itime, 0 /)
2400        edges(1:4) = (/ zsize(pfileid,varid,1), &
2401 &                      zsize(pfileid,varid,2), &
2402 &                       1, 0 /)
2403      ELSE
2404        corner(1:4) = (/ 1, itime, 0, 0 /)
2405        edges(1:4) = (/ zsize(pfileid,varid,1), 1, 0, 0 /)
2406      ENDIF
2407    ELSE
2408      IF ( regular(pfileid) ) THEN
2409        corner(1:4) = (/ 1, 1, 1, itime /)
2410        edges(1:4) = (/ zsize(pfileid,varid,1), &
2411 &                      zsize(pfileid,varid,2), &
2412 &                      zsize(pfileid,varid,3), 1 /)
2413      ELSE
2414        corner(1:4) = (/ 1, 1, itime, 0 /)
2415        edges(1:4) = (/ zsize(pfileid,varid,1), &
2416 &                      zsize(pfileid,varid,3), 1, 0 /)
2417      ENDIF
2418    ENDIF
2419!-
2420    ipt = point(pfileid,varid)
2421!-
2422    IF (     (TRIM(tmp_opp) /= "inst") &
2423 &      .AND.(TRIM(tmp_opp) /= "once") ) THEN
2424      iret = NF90_PUT_VAR (ncid,ncvarid,buffer(ipt:), &
2425 &                       start=corner(1:4),count=edges(1:4))
2426    ELSE
2427      iret = NF90_PUT_VAR (ncid,ncvarid,buff_tmp2, &
2428 &                       start=corner(1:4),count=edges(1:4))
2429    ENDIF
2430!-
2431    last_wrt(pfileid,varid) = pitau
2432    nb_wrt(pfileid,varid) = nb_wrt(pfileid,varid)+1
2433    nb_opp(pfileid,varid) = 0
2434!---
2435!   After the write the file can be synchronized so that no data is
2436!   lost in case of a crash. This feature gives up on the benefits of
2437!   buffering and should only be used in debuging mode. A flag is
2438!   needed here to switch to this mode.
2439!---
2440!   iret = NF90_SYNC (ncid)
2441!-
2442  ENDIF
2443!----------------------------
2444END SUBROUTINE histwrite_real
2445!===
2446SUBROUTINE histvar_seq (pfid,pvarname,pvid)
2447!---------------------------------------------------------------------
2448!- This subroutine optimized the search for the variable in the table.
2449!- In a first phase it will learn the succession of the variables
2450!- called and then it will use the table to guess what comes next.
2451!- It is the best solution to avoid lengthy searches through array
2452!- vectors.
2453!-
2454!- ARGUMENTS :
2455!-
2456!- pfid  : id of the file on which we work
2457!- pvarname : The name of the variable we are looking for
2458!- pvid     : The var id we found
2459!---------------------------------------------------------------------
2460  IMPLICIT NONE
2461!-
2462  INTEGER,INTENT(in)  :: pfid
2463  CHARACTER(LEN=*),INTENT(IN) :: pvarname
2464  INTEGER,INTENT(out) :: pvid
2465!-
2466  LOGICAL,SAVE :: learning(nb_files_max)=.TRUE.
2467  INTEGER,SAVE :: overlap(nb_files_max) = -1
2468  INTEGER,SAVE :: varseq(nb_files_max, nb_var_max*3)
2469  INTEGER,SAVE :: varseq_len(nb_files_max) = 0
2470  INTEGER,SAVE :: varseq_pos(nb_files_max)
2471  INTEGER,SAVE :: varseq_err(nb_files_max) = 0
2472  INTEGER      :: ib, nb, sp, nx, pos
2473  CHARACTER(LEN=20),DIMENSION(nb_var_max) :: tab_str20
2474  CHARACTER(LEN=20) :: str20
2475  CHARACTER(LEN=70) :: str70
2476  INTEGER      :: tab_str20_length(nb_var_max)
2477!-
2478  LOGICAL :: check = .FALSE.
2479!---------------------------------------------------------------------
2480  nb = nb_var(pfid)
2481!-
2482  IF (check) THEN
2483    WRITE(*,*) 'histvar_seq, start of the subroutine :',learning(pfid)
2484  ENDIF
2485!-
2486  IF (learning(pfid)) THEN
2487!-
2488!-- 1.0 We compute the length over which we are going
2489!--     to check the overlap
2490!-
2491    IF (overlap(pfid) <= 0) THEN
2492      IF (nb_var(pfid) > 6) THEN
2493        overlap(pfid) = nb_var(pfid)/3*2
2494      ELSE
2495        overlap(pfid) = nb_var(pfid)
2496      ENDIF
2497    ENDIF
2498!-
2499!-- 1.1 Find the position of this string
2500!-
2501    str20 = pvarname
2502    tab_str20(1:nb) = name(pfid,1:nb)
2503    tab_str20_length(1:nb) = name_length(pfid,1:nb)
2504!-
2505    CALL find_str (nb, tab_str20, tab_str20_length, str20, pos)
2506!-
2507    IF (pos > 0) THEN
2508      pvid = pos
2509    ELSE
2510      CALL ipslerr (3,"histvar_seq", &
2511 &      'The name of the variable you gave has not been declared', &
2512 &      'You should use subroutine histdef for declaring variable', &
2513 &      TRIM(str20))
2514    ENDIF
2515!-
2516!-- 1.2 If we have not given up we store the position
2517!--     in the sequence of calls
2518!-
2519    IF ( varseq_err(pfid) .GE. 0 ) THEN
2520      sp = varseq_len(pfid)+1
2521      IF (sp <= nb_var_max*3) THEN
2522        varseq(pfid,sp) = pvid
2523        varseq_len(pfid) = sp
2524      ELSE
2525        CALL ipslerr (2,"histvar_seq",&
2526 &       'The learning process has failed and we give up. '// &
2527 &       'Either you sequence is',&
2528 &       'too complex or I am too dumb. '// &
2529 &       'This will only affect the efficiency',&
2530 &       'of your code. Thus if you wish to save time'// &
2531 &       ' contact the IOIPSL team. ')
2532        WRITE(*,*) 'The sequence we have found up to now :'
2533        WRITE(*,*) varseq(pfid,1:sp-1)
2534        varseq_err(pfid) = -1
2535      ENDIF
2536!-
2537!---- 1.3 Check if we have found the right overlap
2538!-
2539      IF (varseq_len(pfid) .GE. overlap(pfid)*2) THEN
2540!-
2541!------ We skip a few variables if needed as they could come
2542!------ from the initialisation of the model.
2543!-
2544        DO ib = 0, sp-overlap(pfid)*2
2545          IF ( learning(pfid) .AND.&
2546            & SUM(ABS(varseq(pfid,ib+1:ib+overlap(pfid)) -&
2547            & varseq(pfid,sp-overlap(pfid)+1:sp))) == 0 ) THEN
2548            learning(pfid) = .FALSE.
2549            varseq_len(pfid) = sp-overlap(pfid)-ib
2550            varseq_pos(pfid) = overlap(pfid)+ib
2551            varseq(pfid,1:varseq_len(pfid)) = &
2552 &            varseq(pfid,ib+1:ib+varseq_len(pfid))
2553          ENDIF
2554        ENDDO
2555      ENDIF
2556    ENDIF
2557  ELSE
2558!-
2559!-- 2.0 Now we know how the calls to histwrite are sequenced
2560!--     and we can get a guess at the var ID
2561!-
2562    nx = varseq_pos(pfid)+1
2563    IF (nx > varseq_len(pfid)) nx = 1
2564!-
2565    pvid = varseq(pfid, nx)
2566!-
2567    IF (    (INDEX(name(pfid,pvid),pvarname) <= 0)         &
2568   &    .OR.(name_length(pfid,pvid) /= len_trim(pvarname)) ) THEN
2569      str20 = pvarname
2570      tab_str20(1:nb) = name(pfid,1:nb)
2571      tab_str20_length(1:nb) = name_length(pfid,1:nb)
2572      CALL find_str (nb,tab_str20,tab_str20_length,str20,pos)
2573      IF (pos > 0) THEN
2574        pvid = pos
2575      ELSE
2576        CALL ipslerr (3,"histvar_seq", &
2577 &  'The name of the variable you gave has not been declared',&
2578 &  'You should use subroutine histdef for declaring variable',str20)
2579      ENDIF
2580      varseq_err(pfid) = varseq_err(pfid)+1
2581    ELSE
2582!-
2583!---- We only keep the new position if we have found the variable
2584!---- this way. This way an out of sequence call to histwrite does
2585!---- not defeat the process.
2586!-
2587      varseq_pos(pfid) = nx
2588    ENDIF
2589!-
2590    IF (varseq_err(pfid) .GE. 10) THEN
2591      WRITE(str70,'("for file ",I3)') pfid
2592      CALL ipslerr (2,"histvar_seq", &
2593 &  'There were 10 errors in the learned sequence of variables',&
2594 &  str70,'This looks like a bug, please report it.')
2595         varseq_err(pfid) = 0
2596    ENDIF
2597  ENDIF
2598!-
2599  IF (check) THEN
2600    WRITE(*,*) &
2601 &   'histvar_seq, end of the subroutine :',TRIM(pvarname),pvid
2602  ENDIF
2603!-------------------------
2604END SUBROUTINE histvar_seq
2605!===
2606SUBROUTINE histsync (file)
2607!---------------------------------------------------------------------
2608!- This subroutine will synchronise all
2609!- (or one if defined) opened files.
2610!-
2611!- VERSION
2612!-
2613!---------------------------------------------------------------------
2614  IMPLICIT NONE
2615!-
2616! file  : optional argument for fileid
2617  INTEGER,INTENT(in),OPTIONAL :: file
2618!-
2619  INTEGER :: ifile,ncid,iret
2620!-
2621  LOGICAL :: file_exists
2622  LOGICAL :: check = .FALSE.
2623!---------------------------------------------------------------------
2624  IF (check) WRITE(*,*) 'Entering loop on files :', nb_files
2625!-
2626! 1.The loop on files to synchronise
2627!-
2628  DO ifile = 1,nb_files
2629!-
2630    IF (PRESENT(file)) THEN
2631      file_exists = (ifile == file)
2632    ELSE
2633      file_exists = .TRUE.
2634    ENDIF
2635!-
2636    IF ( file_exists ) THEN
2637      IF (check) THEN
2638        WRITE(*,*) 'Synchronising specified file number :', file
2639      ENDIF
2640      ncid = ncdf_ids(ifile)
2641      iret = NF90_SYNC (ncid)
2642    ENDIF
2643!-
2644  ENDDO
2645!----------------------
2646END SUBROUTINE histsync
2647!===
2648SUBROUTINE histclo (fid)
2649!---------------------------------------------------------------------
2650!- This subroutine will close all (or one if defined) opened files
2651!-
2652!- VERSION
2653!-
2654!---------------------------------------------------------------------
2655  IMPLICIT NONE
2656!-
2657! fid  : optional argument for fileid
2658  INTEGER,INTENT(in),OPTIONAL :: fid
2659!-
2660  INTEGER :: ifile,ncid,iret,iv,ncvarid
2661  INTEGER :: start_loop,end_loop
2662  CHARACTER(LEN=70) :: str70
2663!-
2664  LOGICAL :: check=.FALSE.
2665!---------------------------------------------------------------------
2666  IF (check) WRITE(*,*) 'Entering loop on files :', nb_files
2667!-
2668  IF (PRESENT(fid)) THEN
2669    start_loop = fid
2670    end_loop = fid
2671  ELSE
2672    start_loop = 1
2673    end_loop = nb_files
2674  ENDIF
2675!-
2676  DO ifile=start_loop,end_loop
2677    IF (check) WRITE(*,*) 'Closing specified file number :', ifile
2678    ncid = ncdf_ids(ifile)
2679    iret = NF90_REDEF (ncid)
2680!-
2681!-- 1. The loop on the number of variables to add
2682!-     some final information
2683!-
2684    IF ( check ) WRITE(*,*) 'Entering loop on vars :', nb_var(ifile)
2685    DO iv = 1,nb_var(ifile)
2686      IF (hist_wrt_rng(ifile,iv)) THEN
2687        IF (check) THEN
2688          WRITE(*,*) 'min value for file :',ifile,' var n. :',iv, &
2689         &           ' is : ',hist_minmax(ifile,iv,1)
2690          WRITE(*,*) 'max value for file :',ifile,' var n. :',iv, &
2691         &           ' is : ',hist_minmax(ifile,iv,2)
2692        ENDIF
2693        IF (hist_calc_rng(ifile,iv)) THEN
2694!-------- 1.1 Put the min and max values on the file
2695          ncvarid = ncvar_ids(ifile,iv)
2696          iret = NF90_PUT_ATT (ncid,ncvarid,'valid_min', &
2697 &                             REAL(hist_minmax(ifile,iv,1),KIND=4))
2698          iret = NF90_PUT_ATT (ncid,ncvarid,'valid_max', &
2699 &                             REAL(hist_minmax(ifile,iv,2),KIND=4))
2700        ENDIF
2701      ENDIF
2702    ENDDO
2703!-
2704!-- 2.0 We list the names of the other files
2705!--     in the associated_file attribute
2706!-
2707    IF (nb_files > 1 ) THEN
2708      iret = NF90_PUT_ATT (ncid,NF90_GLOBAL,'associate_file', &
2709 &                         TRIM(assc_file))
2710    ENDIF
2711    IF ( check ) WRITE(*,*) 'close file :', ncid
2712    iret = NF90_CLOSE (ncid)
2713    IF (iret /= NF90_NOERR) THEN
2714      WRITE(str70,'("This file has been already closed :",I3)') ifile
2715      CALL ipslerr (2,'histclo',str70,'',' ')
2716    ENDIF
2717  ENDDO
2718!---------------------
2719END SUBROUTINE histclo
2720!===
2721SUBROUTINE ioconf_modname (str)
2722!---------------------------------------------------------------------
2723!- This subroutine allows to configure the name
2724!- of the model written into the file
2725!---------------------------------------------------------------------
2726  IMPLICIT NONE
2727!-
2728  CHARACTER(LEN=*),INTENT(IN) :: str
2729!---------------------------------------------------------------------
2730  IF (.NOT.lock_modname) THEN
2731    model_name = str(1:MIN(LEN_TRIM(str),80))
2732    lock_modname = .TRUE.
2733  ELSE
2734    CALL ipslerr (2,"ioconf_modname", &
2735   &  'The model name can only be changed once and only', &
2736   &  'before it is used. It is now set to :',model_name)
2737  ENDIF
2738!----------------------------
2739END SUBROUTINE ioconf_modname
2740!-
2741!===
2742!-
2743!-----------------
2744END MODULE histcom
Note: See TracBrowser for help on using the repository browser.