Context Navigation

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

Last change on this file since 1011 was 1011, checked in by bellier, 14 years ago
New handling for time axes
Property svn:keywords set to `Id`
File size: 77.7 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats: