Context Navigation

physiq.F @ 395

Last change on this file since 395 was 395, checked in by lmdzadmin, 22 years ago
Changement dans les sorties histoires (rain devient precip pour eviter toute confusion) + unites, commentaires des fichiers histoire LF
Property svn:eol-style set to `native` Property svn:executable set to ``* Property svn:keywords set to `Author Date Id Revision`
File size: 131.4 KB

Line
1	c
2	c $Header$
3	c
4	SUBROUTINE physiq (nlon,nlev,nqmax ,
5	. debut,lafin,rjourvrai,rjour_ecri,gmtime,pdtphys,
6	. paprs,pplay,pphi,pphis,paire,presnivs,clesphy0,
7	. u,v,t,qx,
8	. omega, cufi, cvfi,
9	. d_u, d_v, d_t, d_qx, d_ps)
10	USE ioipsl
11	USE histcom
12	USE writephys
13
14	IMPLICIT none
15	c======================================================================
16	c
17	c Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
18	c
19	c Objet: Moniteur general de la physique du modele
20	cAA Modifications quant aux traceurs :
21	cAA - uniformisation des parametrisations ds phytrac
22	cAA - stockage des moyennes des champs necessaires
23	cAA en mode traceur off-line
24	c======================================================================
25	c modif ( P. Le Van , 12/10/98 )
26	c
27	c Arguments:
28	c
29	c nlon----input-I-nombre de points horizontaux
30	c nlev----input-I-nombre de couches verticales
31	c nqmax---input-I-nombre de traceurs (y compris vapeur d'eau) = 1
32	c debut---input-L-variable logique indiquant le premier passage
33	c lafin---input-L-variable logique indiquant le dernier passage
34	c rjour---input-R-numero du jour de l'experience
35	c gmtime--input-R-temps universel dans la journee (0 a 86400 s)
36	c pdtphys-input-R-pas d'integration pour la physique (seconde)
37	c paprs---input-R-pression pour chaque inter-couche (en Pa)
38	c pplay---input-R-pression pour le mileu de chaque couche (en Pa)
39	c pphi----input-R-geopotentiel de chaque couche (g z) (reference sol)
40	c pphis---input-R-geopotentiel du sol
41	c paire---input-R-aire de chaque maille
42	c presnivs-input_R_pressions approximat. des milieux couches ( en PA)
43	c u-------input-R-vitesse dans la direction X (de O a E) en m/s
44	c v-------input-R-vitesse Y (de S a N) en m/s
45	c t-------input-R-temperature (K)
46	c qx------input-R-humidite specifique (kg/kg) et d'autres traceurs
47	c d_t_dyn-input-R-tendance dynamique pour "t" (K/s)
48	c d_q_dyn-input-R-tendance dynamique pour "q" (kg/kg/s)
49	c omega---input-R-vitesse verticale en Pa/s
50	c cufi----input-R-resolution des mailles en x (m)
51	c cvfi----input-R-resolution des mailles en y (m)
52	c
53	c d_u-----output-R-tendance physique de "u" (m/s/s)
54	c d_v-----output-R-tendance physique de "v" (m/s/s)
55	c d_t-----output-R-tendance physique de "t" (K/s)
56	c d_qx----output-R-tendance physique de "qx" (kg/kg/s)
57	c d_ps----output-R-tendance physique de la pression au sol
58	c======================================================================
59	#include "dimensions.h"
60	integer jjmp1
61	parameter (jjmp1=jjm+1-1/jjm)
62	#include "dimphy.h"
63	#include "regdim.h"
64	#include "indicesol.h"
65	#include "dimsoil.h"
66	#include "clesphys.h"
67	#include "control.h"
68	#include "temps.h"
69	c======================================================================
70	LOGICAL check ! Verifier la conservation du modele en eau
71	PARAMETER (check=.FALSE.)
72	LOGICAL ok_stratus ! Ajouter artificiellement les stratus
73	PARAMETER (ok_stratus=.FALSE.)
74	c======================================================================
75	c Parametres lies au coupleur OASIS:
76	#include "oasis.h"
77	INTEGER,SAVE :: npas, nexca
78	logical rnpb
79	parameter(rnpb=.true.)
80	c PARAMETER (npas=1440)
81	c PARAMETER (nexca=48)
82	EXTERNAL fromcpl, intocpl, inicma
83	c ocean = type de modele ocean a utiliser: force, slab, couple
84	character*6 ocean
85	SAVE ocean
86
87	c parameter (ocean = 'force ')
88	c parameter (ocean = 'couple')
89	logical ok_ocean
90	c======================================================================
91	c Clef controlant l'activation du cycle diurne:
92	ccc LOGICAL cycle_diurne
93	ccc PARAMETER (cycle_diurne=.FALSE.)
94	c======================================================================
95	c Modele thermique du sol, a activer pour le cycle diurne:
96	ccc LOGICAL soil_model
97	ccc PARAMETER (soil_model=.FALSE.)
98	logical ok_veget
99	save ok_veget
100	c parameter (ok_veget = .true.)
101	c parameter (ok_veget = .false.)
102	c======================================================================
103	c Dans les versions precedentes, l'eau liquide nuageuse utilisee dans
104	c le calcul du rayonnement est celle apres la precipitation des nuages.
105	c Si cette cle new_oliq est activee, ce sera une valeur moyenne entre
106	c la condensation et la precipitation. Cette cle augmente les impacts
107	c radiatifs des nuages.
108	ccc LOGICAL new_oliq
109	ccc PARAMETER (new_oliq=.FALSE.)
110	c======================================================================
111	c Clefs controlant deux parametrisations de l'orographie:
112	cc LOGICAL ok_orodr
113	ccc PARAMETER (ok_orodr=.FALSE.)
114	ccc LOGICAL ok_orolf
115	ccc PARAMETER (ok_orolf=.FALSE.)
116	c======================================================================
117	LOGICAL ok_journe ! sortir le fichier journalier
118	save ok_journe
119	c PARAMETER (ok_journe=.true.)
120	c
121	LOGICAL ok_mensuel ! sortir le fichier mensuel
122	save ok_mensuel
123	c PARAMETER (ok_mensuel=.true.)
124	c
125	LOGICAL ok_instan ! sortir le fichier instantane
126	save ok_instan
127	c PARAMETER (ok_instan=.true.)
128	c
129	LOGICAL ok_region ! sortir le fichier regional
130	PARAMETER (ok_region=.FALSE.)
131	c======================================================================
132	c
133	INTEGER ivap ! indice de traceurs pour vapeur d'eau
134	PARAMETER (ivap=1)
135	INTEGER iliq ! indice de traceurs pour eau liquide
136	PARAMETER (iliq=2)
137
138	INTEGER nvm ! nombre de vegetations
139	PARAMETER (nvm=8)
140	REAL veget(klon,nvm) ! couverture vegetale
141	SAVE veget
142
143	c
144	c
145	c Variables argument:
146	c
147	INTEGER nlon
148	INTEGER nlev
149	INTEGER nqmax
150	REAL rjourvrai, rjour_ecri
151	REAL gmtime
152	REAL pdtphys
153	LOGICAL debut, lafin
154	REAL paprs(klon,klev+1)
155	REAL pplay(klon,klev)
156	REAL pphi(klon,klev)
157	REAL pphis(klon)
158	REAL paire(klon)
159	REAL presnivs(klev)
160	REAL znivsig(klev)
161	REAL zsurf(nbsrf)
162	real cufi(klon), cvfi(klon)
163
164	REAL u(klon,klev)
165	REAL v(klon,klev)
166	REAL t(klon,klev)
167	REAL qx(klon,klev,nqmax)
168
169	REAL t_ancien(klon,klev), q_ancien(klon,klev)
170	SAVE t_ancien, q_ancien
171	LOGICAL ancien_ok
172	SAVE ancien_ok
173
174	REAL d_t_dyn(klon,klev)
175	REAL d_q_dyn(klon,klev)
176
177	REAL omega(klon,klev)
178
179	REAL d_u(klon,klev)
180	REAL d_v(klon,klev)
181	REAL d_t(klon,klev)
182	REAL d_qx(klon,klev,nqmax)
183	REAL d_ps(klon)
184
185	INTEGER longcles
186	PARAMETER ( longcles = 20 )
187	REAL clesphy0( longcles )
188	c
189	c Variables quasi-arguments
190	c
191	REAL xjour
192	SAVE xjour
193	c
194	c
195	c Variables propres a la physique
196	c
197	REAL dtime
198	SAVE dtime ! pas temporel de la physique
199	c
200	INTEGER radpas
201	SAVE radpas ! frequence d'appel rayonnement
202	c
203	REAL radsol(klon)
204	SAVE radsol ! bilan radiatif au sol
205	c
206	REAL rlat(klon)
207	SAVE rlat ! latitude pour chaque point
208	c
209	REAL rlon(klon)
210	SAVE rlon ! longitude pour chaque point
211	c
212	cc INTEGER iflag_con
213	cc SAVE iflag_con ! indicateur de la convection
214	c
215	INTEGER itap
216	SAVE itap ! compteur pour la physique
217	c
218	REAL co2_ppm
219	SAVE co2_ppm ! concentration du CO2
220	c
221	REAL solaire
222	SAVE solaire ! constante solaire
223	c
224	REAL ftsol(klon,nbsrf)
225	SAVE ftsol ! temperature du sol
226	c
227	REAL ftsoil(klon,nsoilmx,nbsrf)
228	SAVE ftsoil ! temperature dans le sol
229	c
230	REAL fevap(klon,nbsrf)
231	SAVE fevap ! evaporation
232	REAL fluxlat(klon,nbsrf)
233	SAVE fluxlat
234	c
235	REAL deltat(klon)
236	SAVE deltat ! ecart avec la SST de reference
237	c
238	REAL fqsol(klon,nbsrf)
239	SAVE fqsol ! humidite du sol
240	c
241	REAL fsnow(klon,nbsrf)
242	SAVE fsnow ! epaisseur neigeuse
243	c
244	REAL falbe(klon,nbsrf)
245	SAVE falbe ! albedo par type de surface
246	REAL falblw(klon,nbsrf)
247	SAVE falblw ! albedo par type de surface
248
249	c
250	c
251	c Parametres de l'Orographie a l'Echelle Sous-Maille (OESM):
252	c
253	REAL zmea(klon)
254	SAVE zmea ! orographie moyenne
255	c
256	REAL zstd(klon)
257	SAVE zstd ! deviation standard de l'OESM
258	c
259	REAL zsig(klon)
260	SAVE zsig ! pente de l'OESM
261	c
262	REAL zgam(klon)
263	save zgam ! anisotropie de l'OESM
264	c
265	REAL zthe(klon)
266	SAVE zthe ! orientation de l'OESM
267	c
268	REAL zpic(klon)
269	SAVE zpic ! Maximum de l'OESM
270	c
271	REAL zval(klon)
272	SAVE zval ! Minimum de l'OESM
273	c
274	REAL rugoro(klon)
275	SAVE rugoro ! longueur de rugosite de l'OESM
276	c
277	REAL zulow(klon),zvlow(klon),zustr(klon), zvstr(klon)
278	c
279	REAL zuthe(klon),zvthe(klon)
280	SAVE zuthe
281	SAVE zvthe
282	INTEGER igwd,idx(klon),itest(klon)
283	c
284	REAL agesno(klon,nbsrf)
285	SAVE agesno ! age de la neige
286	c
287	REAL alb_neig(klon)
288	SAVE alb_neig ! albedo de la neige
289	cKE43
290	c Variables liees a la convection de K. Emanuel (sb):
291	c
292	REAL ema_workcbmf(klon) ! cloud base mass flux
293	SAVE ema_workcbmf
294
295	REAL ema_cbmf(klon) ! cloud base mass flux
296	SAVE ema_cbmf
297
298	REAL ema_pcb(klon) ! cloud base pressure
299	SAVE ema_pcb
300
301	REAL ema_pct(klon) ! cloud top pressure
302	SAVE ema_pct
303
304	REAL bas, top ! cloud base and top levels
305	SAVE bas
306	SAVE top
307
308	REAL Ma(klon,klev) ! undilute upward mass flux
309	SAVE Ma
310	REAL ema_work1(klon, klev), ema_work2(klon, klev)
311	SAVE ema_work1, ema_work2
312	REAL wdn(klon), tdn(klon), qdn(klon)
313	c Variables locales pour la couche limite (al1):
314	c
315	cAl1 REAL pblh(klon) ! Hauteur de couche limite
316	cAl1 SAVE pblh
317	c34EK
318	c
319	c Variables locales:
320	c
321	REAL cdragh(klon) ! drag coefficient pour T and Q
322	REAL cdragm(klon) ! drag coefficient pour vent
323	cAA
324	cAA Pour phytrac
325	cAA
326	REAL ycoefh(klon,klev) ! coef d'echange pour phytrac
327	REAL yu1(klon) ! vents dans la premiere couche U
328	REAL yv1(klon) ! vents dans la premiere couche V
329	LOGICAL offline ! Controle du stockage ds "physique"
330	PARAMETER (offline=.false.)
331	INTEGER physid
332	REAL pfrac_impa(klon,klev)! Produits des coefs lessivage impaction
333	save pfrac_impa
334	REAL pfrac_nucl(klon,klev)! Produits des coefs lessivage nucleation
335	save pfrac_nucl
336	REAL pfrac_1nucl(klon,klev)! Produits des coefs lessi nucl (alpha = 1)
337	save pfrac_1nucl
338	REAL frac_impa(klon,klev) ! fractions d'aerosols lessivees (impaction)
339	REAL frac_nucl(klon,klev) ! idem (nucleation)
340	cAA
341	REAL rain_fall(klon) ! pluie
342	REAL snow_fall(klon) ! neige
343	save snow_fall, rain_fall
344	REAL evap(klon), devap(klon) ! evaporation et sa derivee
345	REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
346	REAL dlw(klon) ! derivee infra rouge
347	REAL bils(klon) ! bilan de chaleur au sol
348	REAL fder(klon) ! Derive de flux (sensible et latente)
349	save fder
350	REAL ve(klon) ! integr. verticale du transport meri. de l'energie
351	REAL vq(klon) ! integr. verticale du transport meri. de l'eau
352	REAL ue(klon) ! integr. verticale du transport zonal de l'energie
353	REAL uq(klon) ! integr. verticale du transport zonal de l'eau
354	c
355	REAL frugs(klon,nbsrf) ! longueur de rugosite
356	save frugs
357	REAL zxrugs(klon) ! longueur de rugosite
358	c
359	c Conditions aux limites
360	c
361	INTEGER julien
362	c
363	INTEGER lmt_pas
364	SAVE lmt_pas ! frequence de mise a jour
365	REAL pctsrf(klon,nbsrf)
366	SAVE pctsrf ! sous-fraction du sol
367	REAL albsol(klon)
368	SAVE albsol ! albedo du sol total
369	REAL albsollw(klon)
370	SAVE albsollw ! albedo du sol total
371	REAL albsol1(klon)
372	SAVE albsol1 ! albedo du sol total
373	REAL albsollw1(klon)
374	SAVE albsollw1 ! albedo du sol total
375
376	REAL wo(klon,klev)
377	SAVE wo ! ozone
378	c======================================================================
379	c
380	c Declaration des procedures appelees
381	c
382	EXTERNAL angle ! calculer angle zenithal du soleil
383	EXTERNAL alboc ! calculer l'albedo sur ocean
384	EXTERNAL albsno ! calculer albedo sur neige
385	EXTERNAL ajsec ! ajustement sec
386	EXTERNAL clmain ! couche limite
387	EXTERNAL condsurf ! lire les conditions aux limites
388	EXTERNAL conlmd ! convection (schema LMD)
389	cKE43
390	EXTERNAL conema3 ! convect4.3
391	EXTERNAL fisrtilp ! schema de condensation a grande echelle (pluie)
392	cAA
393	EXTERNAL fisrtilp_tr ! schema de condensation a grande echelle (pluie)
394	c ! stockage des coefficients necessaires au
395	c ! lessivage OFF-LINE et ON-LINE
396	cAA
397	EXTERNAL hgardfou ! verifier les temperatures
398	EXTERNAL nuage ! calculer les proprietes radiatives
399	EXTERNAL o3cm ! initialiser l'ozone
400	EXTERNAL orbite ! calculer l'orbite terrestre
401	EXTERNAL ozonecm ! prescrire l'ozone
402	EXTERNAL phyetat0 ! lire l'etat initial de la physique
403	EXTERNAL phyredem ! ecrire l'etat de redemarrage de la physique
404	EXTERNAL radlwsw ! rayonnements solaire et infrarouge
405	EXTERNAL suphec ! initialiser certaines constantes
406	EXTERNAL transp ! transport total de l'eau et de l'energie
407	EXTERNAL ecribina ! ecrire le fichier binaire global
408	EXTERNAL ecribins ! ecrire le fichier binaire global
409	EXTERNAL ecrirega ! ecrire le fichier binaire regional
410	EXTERNAL ecriregs ! ecrire le fichier binaire regional
411	c
412	c Variables locales
413	c
414	real clwcon(klon,klev),rnebcon(klon,klev)
415	real clwcon0(klon,klev),rnebcon0(klon,klev)
416	save rnebcon, clwcon
417
418	REAL rhcl(klon,klev) ! humiditi relative ciel clair
419	REAL dialiq(klon,klev) ! eau liquide nuageuse
420	REAL diafra(klon,klev) ! fraction nuageuse
421	REAL cldliq(klon,klev) ! eau liquide nuageuse
422	REAL cldfra(klon,klev) ! fraction nuageuse
423	REAL cldtau(klon,klev) ! epaisseur optique
424	REAL cldemi(klon,klev) ! emissivite infrarouge
425	c
426	C§§§ PB
427	REAL fluxq(klon,klev, nbsrf) ! flux turbulent d'humidite
428	REAL fluxt(klon,klev, nbsrf) ! flux turbulent de chaleur
429	REAL fluxu(klon,klev, nbsrf) ! flux turbulent de vitesse u
430	REAL fluxv(klon,klev, nbsrf) ! flux turbulent de vitesse v
431	c
432	REAL zxfluxt(klon, klev)
433	REAL zxfluxq(klon, klev)
434	REAL zxfluxu(klon, klev)
435	REAL zxfluxv(klon, klev)
436	C§§§
437	REAL heat(klon,klev) ! chauffage solaire
438	REAL heat0(klon,klev) ! chauffage solaire ciel clair
439	REAL cool(klon,klev) ! refroidissement infrarouge
440	REAL cool0(klon,klev) ! refroidissement infrarouge ciel clair
441	REAL topsw(klon), toplw(klon), solsw(klon), sollw(klon)
442	real sollwdown(klon) ! downward LW flux at surface
443	REAL topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
444	REAL albpla(klon)
445	c Le rayonnement n'est pas calcule tous les pas, il faut donc
446	c sauvegarder les sorties du rayonnement
447	SAVE heat,cool,albpla,topsw,toplw,solsw,sollw,sollwdown
448	SAVE topsw0,toplw0,solsw0,sollw0, heat0, cool0
449	INTEGER itaprad
450	SAVE itaprad
451	c
452	REAL conv_q(klon,klev) ! convergence de l'humidite (kg/kg/s)
453	REAL conv_t(klon,klev) ! convergence de la temperature(K/s)
454	c
455	REAL cldl(klon),cldm(klon),cldh(klon) !nuages bas, moyen et haut
456	REAL cldt(klon),cldq(klon) !nuage total, eau liquide integree
457	c
458	REAL zxtsol(klon), zxqsol(klon), zxsnow(klon), zxfluxlat(klon)
459	c
460	REAL dist, rmu0(klon), fract(klon)
461	REAL zdtime, zlongi
462	c
463	CHARACTER*2 str2
464	CHARACTER*2 iqn
465	c
466	REAL qcheck
467	REAL z_avant(klon), z_apres(klon), z_factor(klon)
468	LOGICAL zx_ajustq
469	c
470	REAL za, zb
471	REAL zx_t, zx_qs, zdelta, zcor, zfra, zlvdcp, zlsdcp
472	real zqsat(klon,klev)
473	INTEGER i, k, iq, ig, j, nsrf, ll
474	REAL t_coup
475	PARAMETER (t_coup=234.0)
476	c
477	REAL zphi(klon,klev)
478	REAL zx_tmp_x(iim), zx_tmp_yjjmp1
479	REAL zx_relief(iim,jjmp1)
480	REAL zx_aire(iim,jjmp1)
481	cKE43
482	c Variables locales pour la convection de K. Emanuel (sb):
483	c
484	REAL upwd(klon,klev) ! saturated updraft mass flux
485	REAL dnwd(klon,klev) ! saturated downdraft mass flux
486	REAL dnwd0(klon,klev) ! unsaturated downdraft mass flux
487	REAL tvp(klon,klev) ! virtual temp of lifted parcel
488	REAL cape(klon) ! CAPE
489	SAVE cape
490	REAL pbase(klon) ! cloud base pressure
491	SAVE pbase
492	REAL bbase(klon) ! cloud base buoyancy
493	SAVE bbase
494	REAL rflag(klon) ! flag fonctionnement de convect
495	INTEGER iflagctrl(klon) ! flag fonctionnement de convect
496	c -- convect43:
497	INTEGER ntra ! nb traceurs pour convect4.3
498	REAL pori_con(klon) ! pressure at the origin level of lifted parcel
499	REAL plcl_con(klon),dtma_con(klon),dtlcl_con(klon)
500	REAL dtvpdt1(klon,klev), dtvpdq1(klon,klev)
501	REAL dplcldt(klon), dplcldr(klon)
502	c? . condm_con(klon,klev),conda_con(klon,klev),
503	c? . mr_con(klon,klev),ep_con(klon,klev)
504	c? . ,sadiab(klon,klev),wadiab(klon,klev)
505	c --
506	c34EK
507	c
508	c Variables du changement
509	c
510	c con: convection
511	c lsc: condensation a grande echelle (Large-Scale-Condensation)
512	c ajs: ajustement sec
513	c eva: evaporation de l'eau liquide nuageuse
514	c vdf: couche limite (Vertical DiFfusion)
515	REAL d_t_con(klon,klev),d_q_con(klon,klev)
516	REAL d_u_con(klon,klev),d_v_con(klon,klev)
517	REAL d_t_lsc(klon,klev),d_q_lsc(klon,klev),d_ql_lsc(klon,klev)
518	REAL d_t_ajs(klon,klev), d_q_ajs(klon,klev)
519	REAL d_t_eva(klon,klev),d_q_eva(klon,klev)
520	REAL rneb(klon,klev)
521	c
522	REAL pmfu(klon,klev), pmfd(klon,klev)
523	REAL pen_u(klon,klev), pen_d(klon,klev)
524	REAL pde_u(klon,klev), pde_d(klon,klev)
525	INTEGER kcbot(klon), kctop(klon), kdtop(klon)
526	REAL pmflxr(klon,klev+1), pmflxs(klon,klev+1)
527	REAL prfl(klon,klev+1), psfl(klon,klev+1)
528	c
529	INTEGER ibas_con(klon), itop_con(klon)
530	REAL rain_con(klon), rain_lsc(klon)
531	REAL snow_con(klon), snow_lsc(klon)
532	REAL d_ts(klon,nbsrf)
533	c
534	REAL d_u_vdf(klon,klev), d_v_vdf(klon,klev)
535	REAL d_t_vdf(klon,klev), d_q_vdf(klon,klev)
536	c
537	REAL d_u_oro(klon,klev), d_v_oro(klon,klev)
538	REAL d_t_oro(klon,klev)
539	REAL d_u_lif(klon,klev), d_v_lif(klon,klev)
540	REAL d_t_lif(klon,klev)
541
542	REAL ratqs(klon,klev),ratqss(klon,klev),ratqsc(klon,klev)
543	real ratqsbas,ratqshaut
544	save ratqsbas,ratqshaut, ratqs
545	real zpt_conv(klon,klev)
546
547	c Parametres lies au nouveau schema de nuages (SB, PDF)
548	real fact_cldcon
549	real facttemps
550	logical ok_newmicro
551	save ok_newmicro
552	save fact_cldcon,facttemps
553	real facteur
554
555	integer iflag_cldcon
556	save iflag_cldcon
557
558	logical ptconv(klon,klev)
559
560	c
561	c Variables liees a l'ecriture de la bande histoire physique
562	c
563	INTEGER ecrit_mth
564	SAVE ecrit_mth ! frequence d'ecriture (fichier mensuel)
565	c
566	INTEGER ecrit_day
567	SAVE ecrit_day ! frequence d'ecriture (fichier journalier)
568	c
569	INTEGER ecrit_ins
570	SAVE ecrit_ins ! frequence d'ecriture (fichier instantane)
571	c
572	INTEGER ecrit_reg
573	SAVE ecrit_reg ! frequence d'ecriture
574	c
575	integer itau_w ! pas de temps ecriture = itap + itau_phy
576	c
577	c
578	c Variables locales pour effectuer les appels en serie
579	c
580	REAL t_seri(klon,klev), q_seri(klon,klev)
581	REAL ql_seri(klon,klev),qs_seri(klon,klev)
582	REAL u_seri(klon,klev), v_seri(klon,klev)
583	c
584	REAL tr_seri(klon,klev,nbtr)
585	REAL d_tr(klon,klev,nbtr)
586
587	REAL zx_rh(klon,klev)
588
589	INTEGER length
590	PARAMETER ( length = 100 )
591	REAL tabcntr0( length )
592	c
593	INTEGER ndex2d(iimjjmp1),ndex3d(iimjjmp1*klev)
594	REAL zx_tmp_fi2d(klon)
595	REAL zx_tmp_2d(iim,jjmp1), zx_tmp_3d(iim,jjmp1,klev)
596	REAL zx_lon(iim,jjmp1), zx_lat(iim,jjmp1)
597	c
598	INTEGER nid_day, nid_mth, nid_ins
599	SAVE nid_day, nid_mth, nid_ins
600	c
601	INTEGER nhori, nvert
602	REAL zsto, zout
603	real zjulian
604	save zjulian
605
606	character*20 modname
607	character*80 abort_message
608	logical ok_sync
609	real date0
610	integer idayref
611
612	C essai writephys
613	integer fid_day, fid_mth, fid_ins
614	parameter (fid_ins = 1, fid_day = 2, fid_mth = 3)
615	integer prof2d_on, prof3d_on, prof2d_av, prof3d_av
616	parameter (prof2d_on = 1, prof3d_on = 2,
617	. prof2d_av = 3, prof3d_av = 4)
618	character*30 nom_fichier
619	character*10 varname
620	character*40 vartitle
621	character*20 varunits
622	C Variables liees au bilan d'energie et d'enthalpi
623	INTEGER if_ebil ! level for energy conserv. dignostics
624	SAVE if_ebil
625	REAL ztsol(klon)
626	REAL h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot
627	$ , h_qs_tot, qw_tot, ql_tot, qs_tot , ec_tot
628	SAVE h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot
629	$ , h_qs_tot, qw_tot, ql_tot, qs_tot , ec_tot
630	REAL d_h_vcol, d_h_dair, d_qt, d_qw, d_ql, d_qs, d_ec
631	REAL d_h_vcol_phy
632	REAL fs_bound, fq_bound
633	SAVE d_h_vcol_phy
634	REAL zero_v(klon)
635	CHARACTER*15 ztit
636	INTEGER ip_ebil ! PRINT level for energy conserv. diag.
637	SAVE ip_ebil
638	DATA ip_ebil/2/
639	c
640	c Declaration des constantes et des fonctions thermodynamiques
641	c
642	#include "YOMCST.h"
643	#include "YOETHF.h"
644	#include "FCTTRE.h"
645	c======================================================================
646	modname = 'physiq'
647	IF (if_ebil.ge.1) THEN
648	DO i=1,klon
649	zero_v(i)=0.
650	END DO
651	END IF
652	ok_sync=.TRUE.
653	IF (nqmax .LT. 2) THEN
654	PRINT*, 'eaux vapeur et liquide sont indispensables'
655	CALL ABORT
656	ENDIF
657	IF (debut) THEN
658	CALL suphec ! initialiser constantes et parametres phys.
659	ENDIF
660
661
662	c======================================================================
663	xjour = rjourvrai
664	c
665	c Si c'est le debut, il faut initialiser plusieurs choses
666	c ********
667	c
668	IF (debut) THEN
669	C
670	IF (if_ebil.ge.1) d_h_vcol_phy=0.
671	c
672	c appel a la lecture du run.def physique
673	c
674	call conf_phys(ocean, ok_veget, ok_journe, ok_mensuel,
675	. ok_instan, fact_cldcon, facttemps,ok_newmicro,
676	. iflag_cldcon,ratqsbas,ratqshaut, if_ebil)
677
678	DO k = 2, nvm ! pas de vegetation
679	DO i = 1, klon
680	veget(i,k) = 0.0
681	ENDDO
682	ENDDO
683	DO i = 1, klon
684	veget(i,1) = 1.0 ! il n'y a que du desert
685	ENDDO
686	PRINT*, 'Pas de vegetation; desert partout'
687	c
688	c
689	c Initialiser les compteurs:
690	c
691
692	frugs = 0.
693	itap = 0
694	itaprad = 0
695	c
696	CALL phyetat0 ("startphy.nc",dtime,co2_ppm,solaire,
697	. rlat,rlon,pctsrf, ftsol,ftsoil,deltat,fqsol,fsnow,
698	. falbe, fevap, rain_fall,snow_fall,solsw, sollwdown,
699	. dlw,radsol,frugs,agesno,clesphy0,
700	. zmea,zstd,zsig,zgam,zthe,zpic,zval,rugoro,tabcntr0,
701	. t_ancien, q_ancien, ancien_ok, rnebcon, ratqs,clwcon )
702
703	c
704	radpas = NINT( 86400./dtime/nbapp_rad)
705
706	c
707	CALL printflag( tabcntr0,radpas,ok_ocean,ok_oasis ,ok_journe,
708	, ok_instan, ok_region )
709	c
710	IF (ABS(dtime-pdtphys).GT.0.001) THEN
711	PRINT*, 'Pas physique n est pas correcte',dtime,pdtphys
712	abort_message=' See above '
713	call abort_gcm(modname,abort_message,1)
714	ENDIF
715	IF (nlon .NE. klon) THEN
716	PRINT*, 'nlon et klon ne sont pas coherents', nlon, klon
717	abort_message=' See above '
718	call abort_gcm(modname,abort_message,1)
719	ENDIF
720	IF (nlev .NE. klev) THEN
721	PRINT*, 'nlev et klev ne sont pas coherents', nlev, klev
722	abort_message=' See above '
723	call abort_gcm(modname,abort_message,1)
724	ENDIF
725	c
726	IF (dtime*FLOAT(radpas).GT.21600..AND.cycle_diurne) THEN
727	PRINT*, 'Nbre d appels au rayonnement insuffisant'
728	PRINT*, "Au minimum 4 appels par jour si cycle diurne"
729	abort_message=' See above '
730	call abort_gcm(modname,abort_message,1)
731	ENDIF
732	PRINT*, "Clef pour la convection, iflag_con=", iflag_con
733	c
734	cKE43
735	c Initialisation pour la convection de K.E. (sb):
736	IF (iflag_con.GE.3) THEN
737
738	PRINT, "** Convection de Kerry Emanuel 4.3 "
739	PRINT*, "On va utiliser le melange convectif des traceurs qui"
740	PRINT*, "est calcule dans convect4.3"
741	PRINT*, " !!! penser aux logical flags de phytrac"
742
743	DO i = 1, klon
744	ema_cbmf(i) = 0.
745	ema_pcb(i) = 0.
746	ema_pct(i) = 0.
747	ema_workcbmf(i) = 0.
748	ENDDO
749	ENDIF
750	c34EK
751	IF (ok_orodr) THEN
752	DO i=1,klon
753	rugoro(i) = MAX(1.0e-05, zstd(i)*zsig(i)/2.0)
754	ENDDO
755	CALL SUGWD(klon,klev,paprs,pplay)
756	DO i=1,klon
757	zuthe(i)=0.
758	zvthe(i)=0.
759	if(zstd(i).gt.10.)then
760	zuthe(i)=(1.-zgam(i))*cos(zthe(i))
761	zvthe(i)=(1.-zgam(i))*sin(zthe(i))
762	endif
763	ENDDO
764	ENDIF
765	c
766	c
767	lmt_pas = NINT(86400./dtime * 1.0) ! tous les jours
768	PRINT*,'La frequence de lecture surface est de ', lmt_pas
769	c
770	ecrit_mth = NINT(86400./dtime *ecritphy) ! tous les ecritphy jours
771	IF (ok_mensuel) THEN
772	PRINT*, 'La frequence de sortie mensuelle est de ', ecrit_mth
773	ENDIF
774	ecrit_day = NINT(86400./dtime *1.0) ! tous les jours
775	IF (ok_journe) THEN
776	PRINT*, 'La frequence de sortie journaliere est de ',ecrit_day
777	ENDIF
778	ccc ecrit_ins = NINT(86400./dtime *0.5) ! 2 fois par jour
779	ccc ecrit_ins = NINT(86400./dtime *0.25) ! 4 fois par jour
780	ecrit_ins = NINT(86400./dtime/12.) ! toutes les deux heures
781	ecrit_ins = NINT(86400./dtime/48.) ! a chaque pas de temps
782	IF (ok_instan) THEN
783	PRINT*, 'La frequence de sortie instant. est de ', ecrit_ins
784	ENDIF
785	ecrit_reg = NINT(86400./dtime *0.25) ! 4 fois par jour
786	IF (ok_region) THEN
787	PRINT*, 'La frequence de sortie region est de ', ecrit_reg
788	ENDIF
789
790	c
791	c Initialiser le couplage si necessaire
792	c
793	npas = 0
794	nexca = 0
795	if (ocean == 'couple') then
796	npas = itaufin/ iphysiq
797	nexca = 86400 / dtime
798	write(,)' ##### Ocean couple #####'
799	write(,)' Valeurs des pas de temps'
800	write(,)' npas = ', npas
801	write(,)' nexca = ', nexca
802	endif
803	c
804	c
805	c Gestion calendrier
806
807	c
808	IF (ok_journe) THEN
809	c
810	idayref = day_ref
811	CALL ymds2ju(annee_ref, 1, idayref, 0.0, zjulian)
812	c
813	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
814	DO i = 1, iim
815	zx_lon(i,1) = rlon(i+1)
816	zx_lon(i,jjmp1) = rlon(i+1)
817	ENDDO
818	DO ll=1,klev
819	znivsig(ll)=float(ll)
820	ENDDO
821	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
822	write(,)'zx_lon = ',zx_lon(:,1)
823	write(,)'zx_lat = ',zx_lat(1,:)
824	CALL histbeg("histday", iim,zx_lon(:,1), jjmp1,zx_lat(1,:),
825	. 1,iim,1,jjmp1, itau_phy, zjulian, dtime,
826	. nhori, nid_day)
827	write(,)'Journee ', itau_phy, zjulian
828	CALL histvert(nid_day, "presnivs", "Vertical levels", "mb",
829	. klev, presnivs, nvert)
830	c call histvert(nid_day, 'sig_s', 'Niveaux sigma','-',
831	c . klev, znivsig, nvert)
832	c
833	zsto = dtime
834	zout = dtime * FLOAT(ecrit_day)
835	C Essai writephys
836	c nom_fichier = 'histday1'
837	c call writephy_ini(fid_day,nom_fichier,klon,iim,jjmp1,klev,
838	c . rlon,rlat, presnivs,
839	c . zjulian, dtime)
840	c call writephy_def(prof2d_on, fid_day, "once", zsto, zout, 0)
841	c call writephy_def(prof3d_on, fid_day, "once", zsto, zout,
842	c . klev)
843	c call writephy_def(prof2d_av, fid_day, "ave(X)", zsto, zout, 0)
844	c call writephy_def(prof3d_av, fid_day, "ave(X)", zsto, zout,
845	c . klev)
846
847	c
848	CALL histdef(nid_day, "phis", "Surface geop. height", "-",
849	. iim,jjmp1,nhori, 1,1,1, -99, 32,
850	. "once", zsto,zout)
851	c
852	CALL histdef(nid_day, "aire", "Grid area", "-",
853	. iim,jjmp1,nhori, 1,1,1, -99, 32,
854	. "once", zsto,zout)
855	c
856	c Champs 2D:
857	c
858	CALL histdef(nid_day, "tsol", "Surface Temperature", "K",
859	. iim,jjmp1,nhori, 1,1,1, -99, 32,
860	. "ave(X)", zsto,zout)
861	c
862	CALL histdef(nid_day, "tter", "Surface Temperature", "K",
863	. iim,jjmp1,nhori, 1,1,1, -99, 32,
864	. "ave(X)", zsto,zout)
865	c
866	CALL histdef(nid_day, "tlic", "Surface Temperature", "K",
867	. iim,jjmp1,nhori, 1,1,1, -99, 32,
868	. "ave(X)", zsto,zout)
869	c
870	CALL histdef(nid_day, "toce", "Surface Temperature", "K",
871	. iim,jjmp1,nhori, 1,1,1, -99, 32,
872	. "ave(X)", zsto,zout)
873	c
874	CALL histdef(nid_day, "tsic", "Surface Temperature", "K",
875	. iim,jjmp1,nhori, 1,1,1, -99, 32,
876	. "ave(X)", zsto,zout)
877	c
878	CALL histdef(nid_day, "psol", "Surface Pressure", "Pa",
879	. iim,jjmp1,nhori, 1,1,1, -99, 32,
880	. "ave(X)", zsto,zout)
881	c
882	CALL histdef(nid_day, "precip","Precipitation Totale liq+sol"
883	. , "kg/s",
884	. iim,jjmp1,nhori, 1,1,1, -99, 32,
885	. "ave(X)", zsto,zout)
886	c
887	CALL histdef(nid_day, "snow", "Snow fall", "kg/s",
888	. iim,jjmp1,nhori, 1,1,1, -99, 32,
889	. "ave(X)", zsto,zout)
890	c
891	CALL histdef(nid_day, "snow_mass", "Snow Mass", "kg/m2",
892	. iim,jjmp1,nhori, 1,1,1, -99, 32,
893	. "ave(X)", zsto,zout)
894	c
895	CALL histdef(nid_day, "evap", "Evaporation", "kg/s",
896	. iim,jjmp1,nhori, 1,1,1, -99, 32,
897	. "ave(X)", zsto,zout)
898	c
899	CALL histdef(nid_day, "tops", "Solar rad. at TOA", "W/m2",
900	. iim,jjmp1,nhori, 1,1,1, -99, 32,
901	. "ave(X)", zsto,zout)
902	c
903	CALL histdef(nid_day, "topl", "IR rad. at TOA", "W/m2",
904	. iim,jjmp1,nhori, 1,1,1, -99, 32,
905	. "ave(X)", zsto,zout)
906	c
907	CALL histdef(nid_day, "sols", "Net Solar rad. at surf.",
908	. "W/m2",
909	. iim,jjmp1,nhori, 1,1,1, -99, 32,
910	. "ave(X)", zsto,zout)
911	c
912	CALL histdef(nid_day, "soll", "Net IR rad. at surface", "W/m2",
913	. iim,jjmp1,nhori, 1,1,1, -99, 32,
914	. "ave(X)", zsto,zout)
915	c
916	CALL histdef(nid_day, "solldown", "Down. IR rad. at surface",
917	. "W/m2", iim,jjmp1,nhori, 1,1,1, -99, 32,
918	. "ave(X)", zsto,zout)
919	c
920	CALL histdef(nid_day, "bils", "Surf. total heat flux", "W/m2",
921	. iim,jjmp1,nhori, 1,1,1, -99, 32,
922	. "ave(X)", zsto,zout)
923	c
924	CALL histdef(nid_day, "sens", "Sensible heat flux", "W/m2",
925	. iim,jjmp1,nhori, 1,1,1, -99, 32,
926	. "ave(X)", zsto,zout)
927	c
928	CALL histdef(nid_day, "fder", "Heat flux derivation", "W/m2",
929	. iim,jjmp1,nhori, 1,1,1, -99, 32,
930	. "ave(X)", zsto,zout)
931	c
932	CALL histdef(nid_day, "frtu", "Zonal wind stress", "Pa",
933	. iim,jjmp1,nhori, 1,1,1, -99, 32,
934	. "ave(X)", zsto,zout)
935	c
936	CALL histdef(nid_day, "frtv", "Meridional wind stress", "Pa",
937	. iim,jjmp1,nhori, 1,1,1, -99, 32,
938	. "ave(X)", zsto,zout)
939	c
940	C §§§ PB flux pour chauqe sous surface
941	C
942	DO nsrf = 1, nbsrf
943	C
944	call histdef(nid_day, "pourc_"//clnsurf(nsrf),
945	$ "Fraction"//clnsurf(nsrf), "W/m2",
946	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
947	$ "ave(X)", zsto,zout)
948	C
949	call histdef(nid_day, "tsol_"//clnsurf(nsrf),
950	$ "Fraction"//clnsurf(nsrf), "W/m2",
951	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
952	$ "ave(X)", zsto,zout)
953	C
954	call histdef(nid_day, "sens_"//clnsurf(nsrf),
955	$ "Sensible heat flux "//clnsurf(nsrf), "W/m2",
956	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
957	$ "ave(X)", zsto,zout)
958	c
959	call histdef(nid_day, "lat_"//clnsurf(nsrf),
960	$ "Latent heat flux "//clnsurf(nsrf), "W/m2",
961	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
962	$ "ave(X)", zsto,zout)
963	C
964	call histdef(nid_day, "taux_"//clnsurf(nsrf),
965	$ "Zonal wind stress"//clnsurf(nsrf),"Pa",
966	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
967	$ "ave(X)", zsto,zout)
968
969	call histdef(nid_day, "tauy_"//clnsurf(nsrf),
970	$ "Meridional xind stress "//clnsurf(nsrf), "Pa",
971	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
972	$ "ave(X)", zsto,zout)
973	C
974	call histdef(nid_day, "albe_"//clnsurf(nsrf),
975	$ "Albedo surf. "//clnsurf(nsrf), "W/m2",
976	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
977	$ "ave(X)", zsto,zout)
978	C
979	call histdef(nid_day, "rugs_"//clnsurf(nsrf),
980	$ "Latent heat flux "//clnsurf(nsrf), "W/m2",
981	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
982	$ "ave(X)", zsto,zout)
983
984	C§§§
985	END DO
986
987	CALL histdef(nid_day, "sicf", "Sea-ice fraction", "-",
988	. iim,jjmp1,nhori, 1,1,1, -99, 32,
989	. "ave(X)", zsto,zout)
990	c
991	CALL histdef(nid_day, "cldl", "Low-level cloudiness", "-",
992	. iim,jjmp1,nhori, 1,1,1, -99, 32,
993	. "ave(X)", zsto,zout)
994	c
995	CALL histdef(nid_day, "cldm", "Mid-level cloudiness", "-",
996	. iim,jjmp1,nhori, 1,1,1, -99, 32,
997	. "ave(X)", zsto,zout)
998	c
999	CALL histdef(nid_day, "cldh", "High-level cloudiness", "-",
1000	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1001	. "ave(X)", zsto,zout)
1002	c
1003	CALL histdef(nid_day, "cldt", "Total cloudiness", "-",
1004	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1005	. "ave(X)", zsto,zout)
1006	c
1007	CALL histdef(nid_day, "cldq", "Cloud liquid water path", "-",
1008	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1009	. "ave(X)", zsto,zout)
1010	c
1011	c Champs 3D:
1012	c
1013	CALL histdef(nid_day, "temp", "Air temperature", "K",
1014	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1015	. "ave(X)", zsto,zout)
1016	c
1017	CALL histdef(nid_day, "ovap", "Specific humidity", "Kg/Kg",
1018	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1019	. "ave(X)", zsto,zout)
1020	c
1021	CALL histdef(nid_day, "geop", "Geopotential height", "m",
1022	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1023	. "ave(X)", zsto,zout)
1024	c
1025	CALL histdef(nid_day, "vitu", "Zonal wind", "m/s",
1026	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1027	. "ave(X)", zsto,zout)
1028	c
1029	CALL histdef(nid_day, "vitv", "Meridional wind", "m/s",
1030	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1031	. "ave(X)", zsto,zout)
1032	c
1033	CALL histdef(nid_day, "vitw", "Vertical wind", "m/s",
1034	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1035	. "ave(X)", zsto,zout)
1036	c
1037	CALL histdef(nid_day, "pres", "Air pressure", "Pa",
1038	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1039	. "ave(X)", zsto,zout)
1040	c
1041	CALL histend(nid_day)
1042	c
1043	ndex2d = 0
1044	ndex3d = 0
1045	c
1046	ENDIF ! fin de test sur ok_journe
1047	c
1048	IF (ok_mensuel) THEN
1049	c
1050	idayref = day_ref
1051	CALL ymds2ju(annee_ref, 1, idayref, 0.0, zjulian)
1052	c
1053	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
1054	DO i = 1, iim
1055	zx_lon(i,1) = rlon(i+1)
1056	zx_lon(i,jjmp1) = rlon(i+1)
1057	ENDDO
1058	DO ll=1,klev
1059	znivsig(ll)=float(ll)
1060	ENDDO
1061	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
1062	CALL histbeg("histmth.nc", iim,zx_lon(:,1), jjmp1,zx_lat(1,:),
1063	. 1,iim,1,jjmp1, itau_phy, zjulian, dtime,
1064	. nhori, nid_mth)
1065	write(,)'Mensuel ', itau_phy, zjulian
1066	CALL histvert(nid_mth, "presnivs", "Vertical levels", "mb",
1067	. klev, presnivs, nvert)
1068	c call histvert(nid_mth, 'sig_s', 'Niveaux sigma','-',
1069	c . klev, znivsig, nvert)
1070	c
1071	zsto = dtime
1072	zout = dtime * ecrit_mth
1073	c
1074	CALL histdef(nid_mth, "phis", "Surface geop. height", "-",
1075	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1076	. "once", zsto,zout)
1077	c
1078	CALL histdef(nid_mth, "aire", "Grid area", "-",
1079	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1080	. "once", zsto,zout)
1081	c
1082	c Champs 2D:
1083	c
1084	CALL histdef(nid_mth, "tsol", "Surface Temperature", "K",
1085	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1086	. "ave(X)", zsto,zout)
1087	c
1088	CALL histdef(nid_mth, "psol", "Surface Pressure", "Pa",
1089	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1090	. "ave(X)", zsto,zout)
1091	c
1092	CALL histdef(nid_mth, "qsol", "Surface humidity", "mm",
1093	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1094	. "ave(X)", zsto,zout)
1095	c
1096	CALL histdef(nid_mth, "precip", "Precipitation Totale liq+sol",
1097	. "kg/s",
1098	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1099	. "ave(X)", zsto,zout)
1100	c
1101	CALL histdef(nid_mth, "plul", "Large-scale Precip.", "kg/s",
1102	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1103	. "ave(X)", zsto,zout)
1104	c
1105	CALL histdef(nid_mth, "pluc", "Convective Precip.", "kg/s",
1106	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1107	. "ave(X)", zsto,zout)
1108	c
1109	CALL histdef(nid_mth, "snow", "Snow fall", "kg/s",
1110	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1111	. "ave(X)", zsto,zout)
1112	c
1113	CALL histdef(nid_mth, "snow_mass", "Snow Mass", "kg/m2",
1114	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1115	. "ave(X)", zsto,zout)
1116	c
1117	CALL histdef(nid_mth, "evap", "Evaporation", "kg/s",
1118	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1119	. "ave(X)", zsto,zout)
1120	c
1121	CALL histdef(nid_mth, "tops", "Solar rad. at TOA", "W/m2",
1122	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1123	. "ave(X)", zsto,zout)
1124	c
1125	CALL histdef(nid_mth, "topl", "IR rad. at TOA", "W/m2",
1126	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1127	. "ave(X)", zsto,zout)
1128	c
1129	CALL histdef(nid_mth, "sols", "Solar rad. at surf.", "W/m2",
1130	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1131	. "ave(X)", zsto,zout)
1132	c
1133	CALL histdef(nid_mth, "soll", "IR rad. at surface", "W/m2",
1134	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1135	. "ave(X)", zsto,zout)
1136	c
1137	CALL histdef(nid_mth, "solldown", "Down. IR rad. at surface",
1138	. "W/m2", iim,jjmp1,nhori, 1,1,1, -99, 32,
1139	. "ave(X)", zsto,zout)
1140	c
1141	CALL histdef(nid_mth, "tops0", "Solar rad. at TOA", "W/m2",
1142	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1143	. "ave(X)", zsto,zout)
1144	c
1145	CALL histdef(nid_mth, "topl0", "IR rad. at TOA", "W/m2",
1146	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1147	. "ave(X)", zsto,zout)
1148	c
1149	CALL histdef(nid_mth, "sols0", "Solar rad. at surf.", "W/m2",
1150	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1151	. "ave(X)", zsto,zout)
1152	c
1153	CALL histdef(nid_mth, "soll0", "IR rad. at surface", "W/m2",
1154	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1155	. "ave(X)", zsto,zout)
1156	c
1157	CALL histdef(nid_mth, "bils", "Surf. total heat flux", "W/m2",
1158	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1159	. "ave(X)", zsto,zout)
1160	c
1161	CALL histdef(nid_mth, "sens", "Sensible heat flux", "W/m2",
1162	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1163	. "ave(X)", zsto,zout)
1164	c
1165	CALL histdef(nid_mth, "fder", "Heat flux derivation", "W/m2",
1166	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1167	. "ave(X)", zsto,zout)
1168	c
1169	CALL histdef(nid_mth, "frtu", "Zonal wind stress", "Pa",
1170	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1171	. "ave(X)", zsto,zout)
1172	c
1173	CALL histdef(nid_mth, "frtv", "Meridional wind stress", "Pa",
1174	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1175	. "ave(X)", zsto,zout)
1176	c
1177	DO nsrf = 1, nbsrf
1178	C
1179	call histdef(nid_mth, "pourc_"//clnsurf(nsrf),
1180	$ "Fraction "//clnsurf(nsrf), "W/m2",
1181	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1182	$ "ave(X)", zsto,zout)
1183	C
1184	call histdef(nid_mth, "tsol_"//clnsurf(nsrf),
1185	$ "Fraction "//clnsurf(nsrf), "W/m2",
1186	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1187	$ "ave(X)", zsto,zout)
1188	C
1189	call histdef(nid_mth, "sens_"//clnsurf(nsrf),
1190	$ "Sensible heat flux "//clnsurf(nsrf), "W/m2",
1191	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1192	$ "ave(X)", zsto,zout)
1193	c
1194	call histdef(nid_mth, "lat_"//clnsurf(nsrf),
1195	$ "Latent heat flux "//clnsurf(nsrf), "W/m2",
1196	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1197	$ "ave(X)", zsto,zout)
1198	C
1199	call histdef(nid_mth, "taux_"//clnsurf(nsrf),
1200	$ "Zonal wind stress"//clnsurf(nsrf), "Pa",
1201	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1202	$ "ave(X)", zsto,zout)
1203
1204	call histdef(nid_mth, "tauy_"//clnsurf(nsrf),
1205	$ "Meridional xind stress "//clnsurf(nsrf), "Pa",
1206	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1207	$ "ave(X)", zsto,zout)
1208	c
1209	call histdef(nid_mth, "albe_"//clnsurf(nsrf),
1210	$ "Albedo surf. "//clnsurf(nsrf), "W/m2",
1211	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1212	$ "ave(X)", zsto,zout)
1213	c
1214	call histdef(nid_mth, "rugs_"//clnsurf(nsrf),
1215	$ "Latent heat flux "//clnsurf(nsrf), "W/m2",
1216	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1217	$ "ave(X)", zsto,zout)
1218	c
1219	CALL histdef(nid_mth, "ages_"//clnsurf(nsrf), "Snow age","day",
1220	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1221	. "ave(X)", zsto,zout)
1222
1223	END DO
1224	C
1225	CALL histdef(nid_mth, "sicf", "Sea-ice fraction", "-",
1226	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1227	. "ave(X)", zsto,zout)
1228	c
1229	CALL histdef(nid_mth, "albs", "Surface albedo", "-",
1230	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1231	. "ave(X)", zsto,zout)
1232	CALL histdef(nid_mth, "albslw", "Surface albedo LW", "-",
1233	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1234	. "ave(X)", zsto,zout)
1235	c
1236	CALL histdef(nid_mth, "cdrm", "Momentum drag coef.", "-",
1237	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1238	. "ave(X)", zsto,zout)
1239	c
1240	CALL histdef(nid_mth, "cdrh", "Heat drag coef.", "-",
1241	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1242	. "ave(X)", zsto,zout)
1243	c
1244	CALL histdef(nid_mth, "cldl", "Low-level cloudiness", "-",
1245	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1246	. "ave(X)", zsto,zout)
1247	c
1248	CALL histdef(nid_mth, "cldm", "Mid-level cloudiness", "-",
1249	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1250	. "ave(X)", zsto,zout)
1251	c
1252	CALL histdef(nid_mth, "cldh", "High-level cloudiness", "-",
1253	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1254	. "ave(X)", zsto,zout)
1255	c
1256	CALL histdef(nid_mth, "cldt", "Total cloudiness", "-",
1257	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1258	. "ave(X)", zsto,zout)
1259	c
1260	CALL histdef(nid_mth, "cldq", "Cloud liquid water path", "-",
1261	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1262	. "ave(X)", zsto,zout)
1263	c
1264	CALL histdef(nid_mth, "ue", "Zonal energy transport", "-",
1265	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1266	. "ave(X)", zsto,zout)
1267	c
1268	CALL histdef(nid_mth, "ve", "Merid energy transport", "-",
1269	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1270	. "ave(X)", zsto,zout)
1271	c
1272	CALL histdef(nid_mth, "uq", "Zonal humidity transport", "-",
1273	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1274	. "ave(X)", zsto,zout)
1275	c
1276	CALL histdef(nid_mth, "vq", "Merid humidity transport", "-",
1277	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1278	. "ave(X)", zsto,zout)
1279	cKE43
1280	IF (iflag_con .GE. 3) THEN ! sb
1281	c
1282	CALL histdef(nid_mth, "cape", "Conv avlbl pot ener", "J/Kg",
1283	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1284	. "ave(X)", zsto,zout)
1285	c
1286	CALL histdef(nid_mth, "pbase", "Cld base pressure", "hPa",
1287	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1288	. "ave(X)", zsto,zout)
1289	c
1290	CALL histdef(nid_mth, "ptop", "Cld top pressure", "hPa",
1291	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1292	. "ave(X)", zsto,zout)
1293	c
1294	CALL histdef(nid_mth, "fbase", "Cld base mass flux", "Kg/m2/s",
1295	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1296	. "ave(X)", zsto,zout)
1297	c
1298	c
1299	ENDIF
1300	c34EK
1301	c
1302	c Champs 3D:
1303	c
1304	CALL histdef(nid_mth, "temp", "Air temperature", "K",
1305	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1306	. "ave(X)", zsto,zout)
1307	c
1308	CALL histdef(nid_mth, "ovap", "Specific humidity", "Kg/Kg",
1309	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1310	. "ave(X)", zsto,zout)
1311	c
1312	CALL histdef(nid_mth, "geop", "Geopotential height", "m",
1313	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1314	. "ave(X)", zsto,zout)
1315	c
1316	CALL histdef(nid_mth, "vitu", "Zonal wind", "m/s",
1317	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1318	. "ave(X)", zsto,zout)
1319	c
1320	CALL histdef(nid_mth, "vitv", "Meridional wind", "m/s",
1321	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1322	. "ave(X)", zsto,zout)
1323	c
1324	CALL histdef(nid_mth, "vitw", "Vertical wind", "m/s",
1325	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1326	. "ave(X)", zsto,zout)
1327	c
1328	CALL histdef(nid_mth, "pres", "Air pressure", "Pa",
1329	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1330	. "ave(X)", zsto,zout)
1331	c
1332	CALL histdef(nid_mth, "rneb", "Cloud fraction", "-",
1333	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1334	. "ave(X)", zsto,zout)
1335	c
1336	CALL histdef(nid_mth, "rhum", "Relative humidity", "-",
1337	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1338	. "ave(X)", zsto,zout)
1339	c
1340	CALL histdef(nid_mth, "clwcon", "Cloud Liquid water content"
1341	. , "kg/kg",
1342	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1343	. "ave(X)", zsto,zout)
1344	c
1345	CALL histdef(nid_mth, "oliq", "Liquid water content", "kg/kg",
1346	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1347	. "ave(X)", zsto,zout)
1348	c
1349	CALL histdef(nid_mth, "dtdyn", "Dynamics dT", "K/s",
1350	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1351	. "ave(X)", zsto,zout)
1352	c
1353	CALL histdef(nid_mth, "dqdyn", "Dynamics dQ", "Kg/Kg/s",
1354	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1355	. "ave(X)", zsto,zout)
1356	c
1357	CALL histdef(nid_mth, "dtcon", "Convection dT", "K/s",
1358	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1359	. "ave(X)", zsto,zout)
1360	c
1361	CALL histdef(nid_mth, "ducon", "Convection du", "m/s2",
1362	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1363	. "ave(X)", zsto,zout)
1364	c
1365	CALL histdef(nid_mth, "dqcon", "Convection dQ", "Kg/Kg/s",
1366	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1367	. "ave(X)", zsto,zout)
1368	c
1369	CALL histdef(nid_mth, "dtlsc", "Condensation dT", "K/s",
1370	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1371	. "ave(X)", zsto,zout)
1372	c
1373	CALL histdef(nid_mth, "dqlsc", "Condensation dQ", "Kg/Kg/s",
1374	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1375	. "ave(X)", zsto,zout)
1376	c
1377	CALL histdef(nid_mth, "dtvdf", "Boundary-layer dT", "K/s",
1378	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1379	. "ave(X)", zsto,zout)
1380	c
1381	CALL histdef(nid_mth, "dqvdf", "Boundary-layer dQ", "Kg/Kg/s",
1382	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1383	. "ave(X)", zsto,zout)
1384	c
1385	CALL histdef(nid_mth, "dteva", "Reevaporation dT", "K/s",
1386	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1387	. "ave(X)", zsto,zout)
1388	c
1389	CALL histdef(nid_mth, "dqeva", "Reevaporation dQ", "Kg/Kg/s",
1390	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1391	. "ave(X)", zsto,zout)
1392
1393	CALL histdef(nid_mth, "ptconv", "POINTS CONVECTIFS"," ",
1394	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1395	. "ave(X)", zsto,zout)
1396
1397	CALL histdef(nid_mth, "ratqs", "RATQS"," ",
1398	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1399	. "ave(X)", zsto,zout)
1400
1401	c
1402	CALL histdef(nid_mth, "dtajs", "Dry adjust. dT", "K/s",
1403	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1404	. "ave(X)", zsto,zout)
1405
1406	CALL histdef(nid_mth, "dqajs", "Dry adjust. dQ", "Kg/Kg/s",
1407	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1408	. "ave(X)", zsto,zout)
1409	c
1410	CALL histdef(nid_mth, "dtswr", "SW radiation dT", "K/s",
1411	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1412	. "ave(X)", zsto,zout)
1413	c
1414	CALL histdef(nid_mth, "dtsw0", "SW radiation dT", "K/s",
1415	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1416	. "ave(X)", zsto,zout)
1417	c
1418	CALL histdef(nid_mth, "dtlwr", "LW radiation dT", "K/s",
1419	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1420	. "ave(X)", zsto,zout)
1421	c
1422	CALL histdef(nid_mth, "dtlw0", "LW radiation dT", "K/s",
1423	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1424	. "ave(X)", zsto,zout)
1425	c
1426	CALL histdef(nid_mth, "duvdf", "Boundary-layer dU", "m/s2",
1427	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1428	. "ave(X)", zsto,zout)
1429	c
1430	CALL histdef(nid_mth, "dvvdf", "Boundary-layer dV", "m/s2",
1431	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1432	. "ave(X)", zsto,zout)
1433	c
1434	IF (ok_orodr) THEN
1435	CALL histdef(nid_mth, "duoro", "Orography dU", "m/s2",
1436	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1437	. "ave(X)", zsto,zout)
1438	c
1439	CALL histdef(nid_mth, "dvoro", "Orography dV", "m/s2",
1440	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1441	. "ave(X)", zsto,zout)
1442	c
1443	ENDIF
1444	C
1445	IF (ok_orolf) THEN
1446	CALL histdef(nid_mth, "dulif", "Orography dU", "m/s2",
1447	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1448	. "ave(X)", zsto,zout)
1449	c
1450	CALL histdef(nid_mth, "dvlif", "Orography dV", "m/s2",
1451	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1452	. "ave(X)", zsto,zout)
1453	ENDIF
1454	C
1455	CALL histdef(nid_mth, "ozone", "Ozone concentration", "-",
1456	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1457	. "ave(X)", zsto,zout)
1458	c
1459	if (nqmax.GE.3) THEN
1460	DO iq=1,nqmax-2
1461	IF (iq.LE.99) THEN
1462	WRITE(str2,'(i2.2)') iq
1463	CALL histdef(nid_mth, "trac"//str2, "Tracer No."//str2, "-",
1464	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1465	. "ave(X)", zsto,zout)
1466	ELSE
1467	PRINT*, "Trop de traceurs"
1468	CALL abort
1469	ENDIF
1470	ENDDO
1471	ENDIF
1472	c
1473	cKE43
1474	IF (iflag_con.GE.3) THEN ! (sb)
1475	c
1476	CALL histdef(nid_mth, "upwd", "saturated updraft", "Kg/m2/s",
1477	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1478	. "ave(X)", zsto,zout)
1479	c
1480	CALL histdef(nid_mth, "dnwd", "saturated downdraft","Kg/m2/s",
1481	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1482	. "ave(X)", zsto,zout)
1483	c
1484	CALL histdef(nid_mth, "dnwd0", "unsat. downdraft", "Kg/m2/s",
1485	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1486	. "ave(X)", zsto,zout)
1487	c
1488	CALL histdef(nid_mth,"Ma","undilute adiab updraft","Kg/m2/s",
1489	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1490	. "ave(X)", zsto,zout)
1491	c
1492	c
1493	ENDIF
1494	c34EK
1495	CALL histend(nid_mth)
1496	c
1497	ndex2d = 0
1498	ndex3d = 0
1499	c
1500	ENDIF ! fin de test sur ok_mensuel
1501	c
1502	c
1503	IF (ok_instan) THEN
1504	c
1505	idayref = day_ref
1506	CALL ymds2ju(annee_ref, 1, idayref, 0.0, zjulian)
1507	c
1508	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
1509	DO i = 1, iim
1510	zx_lon(i,1) = rlon(i+1)
1511	zx_lon(i,jjmp1) = rlon(i+1)
1512	ENDDO
1513	DO ll=1,klev
1514	znivsig(ll)=float(ll)
1515	ENDDO
1516	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
1517	CALL histbeg("histins", iim,zx_lon(:,1), jjmp1,zx_lat(1,:),
1518	. 1,iim,1,jjmp1, itau_phy, zjulian, dtime,
1519	. nhori, nid_ins)
1520	write(,)'Inst ', itau_phy, zjulian
1521	CALL histvert(nid_ins, "presnivs", "Vertical levels", "mb",
1522	. klev, presnivs, nvert)
1523	c call histvert(nid_ins, 'sig_s', 'Niveaux sigma','-',
1524	c . klev, znivsig, nvert)
1525	c
1526	c
1527	zsto = dtime * ecrit_ins
1528	zout = dtime * ecrit_ins
1529	C
1530	CALL histdef(nid_ins, "phis", "Surface geop. height", "-",
1531	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1532	. "once", zsto,zout)
1533	c
1534	CALL histdef(nid_ins, "aire", "Grid area", "-",
1535	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1536	. "once", zsto,zout)
1537	c
1538	c Champs 2D:
1539	c
1540	CALL histdef(nid_ins, "tsol", "Surface Temperature", "K",
1541	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1542	. "inst(X)", zsto,zout)
1543	c
1544	CALL histdef(nid_ins, "psol", "Surface Pressure", "Pa",
1545	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1546	. "inst(X)", zsto,zout)
1547	c
1548	CALL histdef(nid_ins, "plul", "Large-scale Precip.", "mm/day",
1549	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1550	. "inst(X)", zsto,zout)
1551	c
1552	CALL histdef(nid_ins, "pluc", "Convective Precip.", "mm/day",
1553	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1554	. "inst(X)", zsto,zout)
1555
1556	CALL histdef(nid_ins, "qsol", "Surface humidity", "mm",
1557	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1558	. "inst(X)", zsto,zout)
1559	c
1560	CALL histdef(nid_ins, "cdrm", "Momentum drag coef.", "-",
1561	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1562	. "inst(X)", zsto,zout)
1563	c
1564	CALL histdef(nid_ins, "cdrh", "Heat drag coef.", "-",
1565	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1566	. "inst(X)", zsto,zout)
1567	c
1568	CALL histdef(nid_ins, "precip", "Precipitation Totale liq+sol",
1569	. "kg/s",
1570	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1571	. "inst(X)", zsto,zout)
1572	c
1573	CALL histdef(nid_ins, "snow", "Snow fall", "kg/s",
1574	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1575	. "inst(X)", zsto,zout)
1576	c
1577	CALL histdef(nid_ins, "snow_mass", "Snow Mass", "kg/m2",
1578	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1579	. "inst(X)", zsto,zout)
1580	c
1581	CALL histdef(nid_ins, "topl", "OLR", "W/m2",
1582	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1583	. "inst(X)", zsto,zout)
1584	c
1585	CALL histdef(nid_ins, "evap", "Evaporation", "kg/s",
1586	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1587	. "inst(X)", zsto,zout)
1588	c
1589	CALL histdef(nid_ins, "sols", "Solar rad. at surf.", "W/m2",
1590	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1591	. "inst(X)", zsto,zout)
1592	c
1593	CALL histdef(nid_ins, "soll", "IR rad. at surface", "W/m2",
1594	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1595	. "inst(X)", zsto,zout)
1596	c
1597	CALL histdef(nid_ins, "solldown", "Down. IR rad. at surface",
1598	. "W/m2", iim,jjmp1,nhori, 1,1,1, -99, 32,
1599	. "inst(X)", zsto,zout)
1600	c
1601	CALL histdef(nid_ins, "bils", "Surf. total heat flux", "W/m2",
1602	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1603	. "inst(X)", zsto,zout)
1604	c
1605	CALL histdef(nid_ins, "sens", "Sensible heat flux", "W/m2",
1606	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1607	. "inst(X)", zsto,zout)
1608	c
1609	CALL histdef(nid_ins, "fder", "Heat flux derivation", "W/m2",
1610	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1611	. "inst(X)", zsto,zout)
1612	c
1613	CALL histdef(nid_ins, "dtsvdfo", "Boundary-layer dTs(o)", "K/s",
1614	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1615	. "inst(X)", zsto,zout)
1616	c
1617	CALL histdef(nid_ins, "dtsvdft", "Boundary-layer dTs(t)", "K/s",
1618	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1619	. "inst(X)", zsto,zout)
1620	c
1621	CALL histdef(nid_ins, "dtsvdfg", "Boundary-layer dTs(g)", "K/s",
1622	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1623	. "inst(X)", zsto,zout)
1624	c
1625	CALL histdef(nid_ins, "dtsvdfi", "Boundary-layer dTs(g)", "K/s",
1626	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1627	. "inst(X)", zsto,zout)
1628
1629	DO nsrf = 1, nbsrf
1630	C
1631	call histdef(nid_ins, "pourc_"//clnsurf(nsrf),
1632	$ "Fraction"//clnsurf(nsrf), "W/m2",
1633	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1634	$ "inst(X)", zsto,zout)
1635
1636	call histdef(nid_ins, "sens_"//clnsurf(nsrf),
1637	$ "Sensible heat flux "//clnsurf(nsrf), "W/m2",
1638	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1639	$ "inst(X)", zsto,zout)
1640	c
1641	call histdef(nid_ins, "tsol_"//clnsurf(nsrf),
1642	$ "Surface Temperature"//clnsurf(nsrf), "W/m2",
1643	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1644	$ "inst(X)", zsto,zout)
1645	c
1646	call histdef(nid_ins, "lat_"//clnsurf(nsrf),
1647	$ "Latent heat flux "//clnsurf(nsrf), "W/m2",
1648	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1649	$ "inst(X)", zsto,zout)
1650	C
1651	call histdef(nid_ins, "taux_"//clnsurf(nsrf),
1652	$ "Zonal wind stress"//clnsurf(nsrf),"Pa",
1653	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1654	$ "inst(X)", zsto,zout)
1655
1656	call histdef(nid_ins, "tauy_"//clnsurf(nsrf),
1657	$ "Meridional xind stress "//clnsurf(nsrf), "Pa",
1658	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1659	$ "inst(X)", zsto,zout)
1660	c
1661	call histdef(nid_ins, "albe_"//clnsurf(nsrf),
1662	$ "Albedo "//clnsurf(nsrf), "-",
1663	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1664	$ "inst(X)", zsto,zout)
1665	c
1666	call histdef(nid_ins, "rugs_"//clnsurf(nsrf),
1667	$ "rugosite "//clnsurf(nsrf), "-",
1668	$ iim,jjmp1,nhori, 1,1,1, -99, 32,
1669	$ "inst(X)", zsto,zout)
1670	C§§§
1671	END DO
1672	CALL histdef(nid_ins, "rugs", "rugosity", "-",
1673	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1674	. "inst(X)", zsto,zout)
1675
1676	c
1677	CALL histdef(nid_ins, "albs", "Surface albedo", "-",
1678	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1679	. "inst(X)", zsto,zout)
1680	CALL histdef(nid_ins, "albslw", "Surface albedo LW", "-",
1681	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1682	. "inst(X)", zsto,zout)
1683	c
1684	c
1685	c Champs 3D:
1686	c
1687	CALL histdef(nid_ins, "temp", "Temperature", "K",
1688	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1689	. "inst(X)", zsto,zout)
1690	c
1691	CALL histdef(nid_ins, "vitu", "Zonal wind", "m/s",
1692	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1693	. "inst(X)", zsto,zout)
1694	c
1695	CALL histdef(nid_ins, "vitv", "Merid wind", "m/s",
1696	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1697	. "inst(X)", zsto,zout)
1698	c
1699	CALL histdef(nid_ins, "geop", "Geopotential height", "m",
1700	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1701	. "inst(X)", zsto,zout)
1702	c
1703	CALL histdef(nid_ins, "pres", "Air pressure", "Pa",
1704	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1705	. "inst(X)", zsto,zout)
1706	c
1707	CALL histdef(nid_ins, "dtvdf", "Boundary-layer dT", "K/s",
1708	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1709	. "inst(X)", zsto,zout)
1710	c
1711	CALL histdef(nid_ins, "dqvdf", "Boundary-layer dQ", "Kg/Kg/s",
1712	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1713	. "inst(X)", zsto,zout)
1714	c
1715
1716	CALL histend(nid_ins)
1717	c
1718	ndex2d = 0
1719	ndex3d = 0
1720	c
1721	ENDIF
1722
1723	c$$$PB Positionner date0 pour initialisation de ORCHIDEE
1724	date0 = zjulian
1725	C date0 = day_ini
1726	WRITE(,) 'physiq date0 : ',date0
1727	c
1728	c
1729	c
1730	c Prescrire l'ozone dans l'atmosphere
1731	c
1732	c
1733	cc DO i = 1, klon
1734	cc DO k = 1, klev
1735	cc CALL o3cm (paprs(i,k)/100.,paprs(i,k+1)/100., wo(i,k),20)
1736	cc ENDDO
1737	cc ENDDO
1738	c
1739	c
1740	ENDIF
1741	c
1742	c ************** Fin de IF ( debut ) *************
1743	c
1744	c
1745	c Mettre a zero des variables de sortie (pour securite)
1746	c
1747	DO i = 1, klon
1748	d_ps(i) = 0.0
1749	ENDDO
1750	DO k = 1, klev
1751	DO i = 1, klon
1752	d_t(i,k) = 0.0
1753	d_u(i,k) = 0.0
1754	d_v(i,k) = 0.0
1755	ENDDO
1756	ENDDO
1757	DO iq = 1, nqmax
1758	DO k = 1, klev
1759	DO i = 1, klon
1760	d_qx(i,k,iq) = 0.0
1761	ENDDO
1762	ENDDO
1763	ENDDO
1764	c
1765	c Ne pas affecter les valeurs entrees de u, v, h, et q
1766	c
1767	DO k = 1, klev
1768	DO i = 1, klon
1769	t_seri(i,k) = t(i,k)
1770	u_seri(i,k) = u(i,k)
1771	v_seri(i,k) = v(i,k)
1772	q_seri(i,k) = qx(i,k,ivap)
1773	ql_seri(i,k) = qx(i,k,iliq)
1774	qs_seri(i,k) = 0.
1775	ENDDO
1776	ENDDO
1777	IF (nqmax.GE.3) THEN
1778	DO iq = 3, nqmax
1779	DO k = 1, klev
1780	DO i = 1, klon
1781	tr_seri(i,k,iq-2) = qx(i,k,iq)
1782	ENDDO
1783	ENDDO
1784	ENDDO
1785	ELSE
1786	DO k = 1, klev
1787	DO i = 1, klon
1788	tr_seri(i,k,1) = 0.0
1789	ENDDO
1790	ENDDO
1791	ENDIF
1792	C
1793	IF (if_ebil.ge.1) THEN
1794	DO i = 1, klon
1795	ztsol(i) = 0.
1796	ENDDO
1797	DO nsrf = 1, nbsrf
1798	DO i = 1, klon
1799	ztsol(i) = ztsol(i) + ftsol(i,nsrf)*pctsrf(i,nsrf)
1800	ENDDO
1801	ENDDO
1802	ztit='after dynamic'
1803	CALL diagetpq(paire,ztit,ip_ebil,1,1,dtime
1804	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
1805	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1806	C Comme les tendances de la physique sont ajoute dans la dynamique,
1807	C on devrait avoir que la variation d'entalpie par la dynamique
1808	C est egale a la variation de la physique au pas de temps precedent.
1809	C Donc la somme de ces 2 variations devrait etre nulle.
1810	call diagphy(paire,ztit,ip_ebil
1811	e , zero_v, zero_v, zero_v, zero_v, zero_v
1812	e , zero_v, zero_v, zero_v, ztsol
1813	e , d_h_vcol+d_h_vcol_phy, d_qt, 0.
1814	s , fs_bound, fq_bound )
1815	END IF
1816
1817	c Diagnostiquer la tendance dynamique
1818	c
1819	IF (ancien_ok) THEN
1820	DO k = 1, klev
1821	DO i = 1, klon
1822	d_t_dyn(i,k) = (t_seri(i,k)-t_ancien(i,k))/dtime
1823	d_q_dyn(i,k) = (q_seri(i,k)-q_ancien(i,k))/dtime
1824	ENDDO
1825	ENDDO
1826	ELSE
1827	DO k = 1, klev
1828	DO i = 1, klon
1829	d_t_dyn(i,k) = 0.0
1830	d_q_dyn(i,k) = 0.0
1831	ENDDO
1832	ENDDO
1833	ancien_ok = .TRUE.
1834	ENDIF
1835	c
1836	c Ajouter le geopotentiel du sol:
1837	c
1838	DO k = 1, klev
1839	DO i = 1, klon
1840	zphi(i,k) = pphi(i,k) + pphis(i)
1841	ENDDO
1842	ENDDO
1843	c
1844	c Verifier les temperatures
1845	c
1846	CALL hgardfou(t_seri,ftsol,'debutphy')
1847	c
1848	c Incrementer le compteur de la physique
1849	c
1850	itap = itap + 1
1851	julien = MOD(NINT(xjour),360)
1852	c
1853	c Mettre en action les conditions aux limites (albedo, sst, etc.).
1854	c Prescrire l'ozone et calculer l'albedo sur l'ocean.
1855	c
1856	IF (MOD(itap-1,lmt_pas) .EQ. 0) THEN
1857	PRINT *,' PHYS cond julien ',julien
1858	CALL ozonecm( FLOAT(julien), rlat, paprs, wo)
1859	ENDIF
1860	c
1861	c Re-evaporer l'eau liquide nuageuse
1862	c
1863	DO k = 1, klev ! re-evaporation de l'eau liquide nuageuse
1864	DO i = 1, klon
1865	zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
1866	c zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
1867	zlsdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
1868	zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
1869	zb = MAX(0.0,ql_seri(i,k))
1870	za = - MAX(0.0,ql_seri(i,k))
1871	. * (zlvdcp(1.-zdelta)+zlsdcpzdelta)
1872	t_seri(i,k) = t_seri(i,k) + za
1873	q_seri(i,k) = q_seri(i,k) + zb
1874	ql_seri(i,k) = 0.0
1875	d_t_eva(i,k) = za
1876	d_q_eva(i,k) = zb
1877	ENDDO
1878	ENDDO
1879	c
1880	IF (if_ebil.ge.2) THEN
1881	ztit='after reevap'
1882	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
1883	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
1884	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1885	call diagphy(paire,ztit,ip_ebil
1886	e , zero_v, zero_v, zero_v, zero_v, zero_v
1887	e , zero_v, zero_v, zero_v, ztsol
1888	e , d_h_vcol, d_qt, d_ec
1889	s , fs_bound, fq_bound )
1890	C
1891	END IF
1892	C
1893	c
1894	c Appeler la diffusion verticale (programme de couche limite)
1895	c
1896	DO i = 1, klon
1897	c if (.not. ok_veget) then
1898	c frugs(i,is_ter) = SQRT(frugs(i,is_ter)2+rugoro(i)2)
1899	c endif
1900	c frugs(i,is_lic) = rugoro(i)
1901	c frugs(i,is_oce) = rugmer(i)
1902	c frugs(i,is_sic) = 0.001
1903	zxrugs(i) = 0.0
1904	ENDDO
1905	DO nsrf = 1, nbsrf
1906	DO i = 1, klon
1907	frugs(i,nsrf) = MAX(frugs(i,nsrf),0.001)
1908	cccc frugs(i,nsrf) = MAX(frugs(i,nsrf),0.000015)
1909	ENDDO
1910	ENDDO
1911	DO nsrf = 1, nbsrf
1912	DO i = 1, klon
1913	zxrugs(i) = zxrugs(i) + frugs(i,nsrf)*pctsrf(i,nsrf)
1914	ENDDO
1915	ENDDO
1916	c
1917	C calculs necessaires au calcul de l'albedo dans l'interface
1918	c
1919	CALL orbite(FLOAT(julien),zlongi,dist)
1920	IF (cycle_diurne) THEN
1921	zdtime=dtime*FLOAT(radpas) ! pas de temps du rayonnement (s)
1922	CALL zenang(zlongi,gmtime,zdtime,rlat,rlon,rmu0,fract)
1923	ELSE
1924	rmu0 = -999.999
1925	ENDIF
1926
1927	fder = dlw
1928
1929	CALL clmain(dtime,itap,date0,pctsrf,
1930	e t_seri,q_seri,u_seri,v_seri,
1931	e julien, rmu0,
1932	e ok_veget, ocean, npas, nexca, ftsol,
1933	$ soil_model,ftsoil,
1934	$ paprs,pplay,radsol, fsnow,fqsol,fevap,falbe,falblw,
1935	$ fluxlat,
1936	e rain_fall, snow_fall, solsw, sollw, sollwdown, fder,
1937	e rlon, rlat, cufi, cvfi, frugs,
1938	e debut, lafin, agesno,rugoro ,
1939	s d_t_vdf,d_q_vdf,d_u_vdf,d_v_vdf,d_ts,
1940	s fluxt,fluxq,fluxu,fluxv,cdragh,cdragm,
1941	s dsens, devap,
1942	s ycoefh,yu1,yv1)
1943
1944	c
1945	C§§§ PB
1946	C§§§ Incrementation des flux
1947	C§§
1948	zxfluxt=0.
1949	zxfluxq=0.
1950	zxfluxu=0.
1951	zxfluxv=0.
1952	DO nsrf = 1, nbsrf
1953	DO k = 1, klev
1954	DO i = 1, klon
1955	zxfluxt(i,k) = zxfluxt(i,k) +
1956	$ fluxt(i,k,nsrf) * pctsrf( i, nsrf)
1957	zxfluxq(i,k) = zxfluxq(i,k) +
1958	$ fluxq(i,k,nsrf) * pctsrf( i, nsrf)
1959	zxfluxu(i,k) = zxfluxu(i,k) +
1960	$ fluxu(i,k,nsrf) * pctsrf( i, nsrf)
1961	zxfluxv(i,k) = zxfluxv(i,k) +
1962	$ fluxv(i,k,nsrf) * pctsrf( i, nsrf)
1963	END DO
1964	END DO
1965	END DO
1966	DO i = 1, klon
1967	sens(i) = - zxfluxt(i,1) ! flux de chaleur sensible au sol
1968	c evap(i) = - fluxq(i,1) ! flux d'evaporation au sol
1969	evap(i) = - zxfluxq(i,1) ! flux d'evaporation au sol
1970	fder(i) = dlw(i) + dsens(i) + devap(i)
1971	ENDDO
1972
1973
1974	DO k = 1, klev
1975	DO i = 1, klon
1976	t_seri(i,k) = t_seri(i,k) + d_t_vdf(i,k)
1977	q_seri(i,k) = q_seri(i,k) + d_q_vdf(i,k)
1978	u_seri(i,k) = u_seri(i,k) + d_u_vdf(i,k)
1979	v_seri(i,k) = v_seri(i,k) + d_v_vdf(i,k)
1980	ENDDO
1981	ENDDO
1982	c
1983	IF (if_ebil.ge.2) THEN
1984	ztit='after clmain'
1985	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
1986	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
1987	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1988	call diagphy(paire,ztit,ip_ebil
1989	e , zero_v, zero_v, zero_v, zero_v, sens
1990	e , evap , zero_v, zero_v, ztsol
1991	e , d_h_vcol, d_qt, d_ec
1992	s , fs_bound, fq_bound )
1993	END IF
1994	C
1995	c
1996	c Incrementer la temperature du sol
1997	c
1998	DO i = 1, klon
1999	zxtsol(i) = 0.0
2000	zxfluxlat(i) = 0.0
2001	IF ( abs( pctsrf(i, is_ter) + pctsrf(i, is_lic) +
2002	$ pctsrf(i, is_oce) + pctsrf(i, is_sic) - 1.) .GT. EPSFRA)
2003	$ THEN
2004	WRITE(,) 'physiq : pb sous surface au point ', i,
2005	$ pctsrf(i, 1 : nbsrf)
2006	ENDIF
2007	ENDDO
2008	DO nsrf = 1, nbsrf
2009	DO i = 1, klon
2010	c$$$ IF (pctsrf(i,nsrf) .GE. EPSFRA) THEN
2011	ftsol(i,nsrf) = ftsol(i,nsrf) + d_ts(i,nsrf)
2012	zxtsol(i) = zxtsol(i) + ftsol(i,nsrf)*pctsrf(i,nsrf)
2013	zxfluxlat(i) = zxfluxlat(i) + fluxlat(i,nsrf)*pctsrf(i,nsrf)
2014	c$$$ ENDIF
2015	ENDDO
2016	ENDDO
2017
2018	c
2019	c Si une sous-fraction n'existe pas, elle prend la temp. moyenne
2020	c
2021	DO nsrf = 1, nbsrf
2022	DO i = 1, klon
2023	IF (pctsrf(i,nsrf) .LT. epsfra) ftsol(i,nsrf) = zxtsol(i)
2024	ENDDO
2025	ENDDO
2026
2027	c
2028	c Calculer la derive du flux infrarouge
2029	c
2030	c$$$ DO nsrf = 1, nbsrf
2031	DO i = 1, klon
2032	c$$$ IF (pctsrf(i,nsrf) .GE. EPSFRA) THEN
2033	dlw(i) = - 4.0RSIGMAzxtsol(i)**3
2034	c$$$ . *(ftsol(i,nsrf)-zxtsol(i))
2035	c$$$ . *pctsrf(i,nsrf)
2036	c$$$ ENDIF
2037	c$$$ ENDDO
2038	ENDDO
2039	c
2040	c Appeler la convection (au choix)
2041	c
2042	DO k = 1, klev
2043	DO i = 1, klon
2044	conv_q(i,k) = d_q_dyn(i,k)
2045	. + d_q_vdf(i,k)/dtime
2046	conv_t(i,k) = d_t_dyn(i,k)
2047	. + d_t_vdf(i,k)/dtime
2048	ENDDO
2049	ENDDO
2050	IF (check) THEN
2051	za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
2052	PRINT*, "avantcon=", za
2053	ENDIF
2054	zx_ajustq = .FALSE.
2055	IF (iflag_con.EQ.2) zx_ajustq=.TRUE.
2056	IF (zx_ajustq) THEN
2057	DO i = 1, klon
2058	z_avant(i) = 0.0
2059	ENDDO
2060	DO k = 1, klev
2061	DO i = 1, klon
2062	z_avant(i) = z_avant(i) + (q_seri(i,k)+ql_seri(i,k))
2063	. *(paprs(i,k)-paprs(i,k+1))/RG
2064	ENDDO
2065	ENDDO
2066	ENDIF
2067	IF (iflag_con.EQ.1) THEN
2068	stop'reactiver le call conlmd dans physiq.F'
2069	c CALL conlmd (dtime, paprs, pplay, t_seri, q_seri, conv_q,
2070	c . d_t_con, d_q_con,
2071	c . rain_con, snow_con, ibas_con, itop_con)
2072	ELSE IF (iflag_con.EQ.2) THEN
2073	CALL conflx(dtime, paprs, pplay, t_seri, q_seri,
2074	e conv_t, conv_q, zxfluxq(1,1), omega,
2075	s d_t_con, d_q_con, rain_con, snow_con,
2076	s pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
2077	s kcbot, kctop, kdtop, pmflxr, pmflxs)
2078	WHERE (rain_con < 0.) rain_con = 0.
2079	WHERE (snow_con < 0.) snow_con = 0.
2080	DO i = 1, klon
2081	ibas_con(i) = klev+1 - kcbot(i)
2082	itop_con(i) = klev+1 - kctop(i)
2083	ENDDO
2084	ELSE IF (iflag_con.GE.3) THEN
2085	c nb of tracers for the KE convection:
2086	if (nqmax .GE. 4) then
2087	ntra = nbtr
2088	else
2089	ntra = 1
2090	endif
2091	if (iflag_con.eq.4) then ! vectorise
2092	CALL conemav (dtime,paprs,pplay,t_seri,q_seri,
2093	. u_seri,v_seri,tr_seri,nbtr,
2094	. ema_work1,ema_work2,
2095	. d_t_con,d_q_con,d_u_con,d_v_con,d_tr,
2096	. rain_con, snow_con, ibas_con, itop_con,
2097	. upwd,dnwd,dnwd0,
2098	c . Ma,cape,tvp,(/(nint(rflag(i)),i=1,size(rflag))/),
2099	. Ma,cape,tvp,iflagctrl,
2100	. pbase,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr)
2101
2102	else
2103
2104	c print*,'Avant conema OUI'
2105	CALL conema3 (dtime,
2106	. paprs,pplay,t_seri,q_seri,
2107	. u_seri,v_seri,tr_seri,nbtr,
2108	. ema_work1,ema_work2,
2109	. d_t_con,d_q_con,d_u_con,d_v_con,d_tr,
2110	. rain_con, snow_con, ibas_con, itop_con,
2111	. upwd,dnwd,dnwd0,bas,top,
2112	. Ma,cape,tvp,rflag,
2113	. pbase
2114	. ,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr
2115	. ,clwcon0)
2116	print*,'Apres conema3 '
2117
2118	c Calculer l'humidite relative pour diagnostique
2119	c
2120	DO k = 1, klev
2121	DO i = 1, klon
2122	zx_t = t_seri(i,k)
2123	IF (thermcep) THEN
2124	zdelta = MAX(0.,SIGN(1.,rtt-zx_t))
2125	zx_qs = r2es * FOEEW(zx_t,zdelta)/pplay(i,k)
2126	zx_qs = MIN(0.5,zx_qs)
2127	zcor = 1./(1.-retv*zx_qs)
2128	zx_qs = zx_qs*zcor
2129	ELSE
2130	IF (zx_t.LT.t_coup) THEN
2131	zx_qs = qsats(zx_t)/pplay(i,k)
2132	ELSE
2133	zx_qs = qsatl(zx_t)/pplay(i,k)
2134	ENDIF
2135	ENDIF
2136	zqsat(i,k)=zx_qs
2137	ENDDO
2138	ENDDO
2139
2140	c calcul des propriétés des nuages convectifs
2141	clwcon0(:,:)=fact_cldcon*clwcon0(:,:)
2142	call clouds_gno
2143	s (klon,klev,q_seri,zqsat,clwcon0,ptconv,ratqsc,rnebcon0)
2144
2145	endif
2146	DO i = 1, klon
2147	ema_pcb(i) = pbase(i)
2148	ENDDO
2149	DO i = 1, klon
2150	ema_pct(i) = paprs(i,itop_con(i))
2151	ENDDO
2152	DO i = 1, klon
2153	ema_cbmf(i) = ema_workcbmf(i)
2154	ENDDO
2155	ELSE
2156	PRINT*, "iflag_con non-prevu", iflag_con
2157	CALL abort
2158	ENDIF
2159
2160	c CALL homogene(paprs, q_seri, d_q_con, u_seri,v_seri,
2161	c . d_u_con, d_v_con)
2162	DO k = 1, klev
2163	DO i = 1, klon
2164	t_seri(i,k) = t_seri(i,k) + d_t_con(i,k)
2165	q_seri(i,k) = q_seri(i,k) + d_q_con(i,k)
2166	u_seri(i,k) = u_seri(i,k) + d_u_con(i,k)
2167	v_seri(i,k) = v_seri(i,k) + d_v_con(i,k)
2168	ENDDO
2169	ENDDO
2170	c
2171	IF (if_ebil.ge.2) THEN
2172	ztit='after convect'
2173	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
2174	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
2175	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
2176	call diagphy(paire,ztit,ip_ebil
2177	e , zero_v, zero_v, zero_v, zero_v, zero_v
2178	e , zero_v, rain_con, snow_con, ztsol
2179	e , d_h_vcol, d_qt, d_ec
2180	s , fs_bound, fq_bound )
2181	END IF
2182	C
2183	IF (check) THEN
2184	za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
2185	PRINT*, "aprescon=", za
2186	zx_t = 0.0
2187	za = 0.0
2188	DO i = 1, klon
2189	za = za + paire(i)/FLOAT(klon)
2190	zx_t = zx_t + (rain_con(i)+snow_con(i))*paire(i)/FLOAT(klon)
2191	ENDDO
2192	zx_t = zx_t/za*dtime
2193	PRINT*, "Precip=", zx_t
2194	ENDIF
2195	IF (zx_ajustq) THEN
2196	DO i = 1, klon
2197	z_apres(i) = 0.0
2198	ENDDO
2199	DO k = 1, klev
2200	DO i = 1, klon
2201	z_apres(i) = z_apres(i) + (q_seri(i,k)+ql_seri(i,k))
2202	. *(paprs(i,k)-paprs(i,k+1))/RG
2203	ENDDO
2204	ENDDO
2205	DO i = 1, klon
2206	z_factor(i) = (z_avant(i)-(rain_con(i)+snow_con(i))*dtime)
2207	. /z_apres(i)
2208	ENDDO
2209	DO k = 1, klev
2210	DO i = 1, klon
2211	IF (z_factor(i).GT.(1.0+1.0E-08) .OR.
2212	. z_factor(i).LT.(1.0-1.0E-08)) THEN
2213	q_seri(i,k) = q_seri(i,k) * z_factor(i)
2214	ENDIF
2215	ENDDO
2216	ENDDO
2217	ENDIF
2218	zx_ajustq=.FALSE.
2219	c
2220	IF (nqmax.GT.2) THEN !--melange convectif de traceurs
2221	c
2222	IF (iflag_con .NE. 2 .AND. debut) THEN
2223	PRINT*, 'Pour l instant, seul conflx fonctionne ',
2224	$ 'avec traceurs', iflag_con
2225	PRINT*,' Mettre iflag_con',
2226	$ ' = 2 dans run.def et repasser'
2227	c CALL abort
2228	ENDIF
2229	c
2230	ENDIF !--nqmax.GT.2
2231	c
2232	c Appeler l'ajustement sec
2233	c
2234	CALL ajsec(paprs, pplay, t_seri, q_seri, d_t_ajs, d_q_ajs)
2235	DO k = 1, klev
2236	DO i = 1, klon
2237	t_seri(i,k) = t_seri(i,k) + d_t_ajs(i,k)
2238	q_seri(i,k) = q_seri(i,k) + d_q_ajs(i,k)
2239	ENDDO
2240	ENDDO
2241	c
2242	IF (if_ebil.ge.2) THEN
2243	ztit='after dry_adjust'
2244	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
2245	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
2246	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
2247	END IF
2248
2249
2250	c-------------------------------------------------------------------------
2251	c Caclul des ratqs
2252	c-------------------------------------------------------------------------
2253
2254	c print*,'calcul des ratqs'
2255	c ratqs convectifs a l'ancienne en fonction de q(z=0)-q / q
2256	c ----------------
2257	c on ecrase le tableau ratqsc calcule par clouds_gno
2258	if (iflag_cldcon.eq.1) then
2259	do k=1,klev
2260	do i=1,klon
2261	if(ptconv(i,k)) then
2262	ratqsc(i,k)=ratqsbas
2263	s +fact_cldcon*(q_seri(i,1)-q_seri(i,k))/q_seri(i,k)
2264	else
2265	ratqsc(i,k)=0.
2266	endif
2267	enddo
2268	enddo
2269	endif
2270
2271	c ratqs stables
2272	c -------------
2273	do k=1,klev
2274	ratqss(:,k)=ratqsbas+(ratqshaut-ratqsbas)*
2275	s min((paprs(:,1)-pplay(:,k))/(paprs(:,1)-30000.),1.)
2276	enddo
2277
2278
2279	c ratqs final
2280	c -----------
2281	if (iflag_cldcon.eq.1 .or.iflag_cldcon.eq.2) then
2282	c les ratqs sont une conbinaison de ratqss et ratqsc
2283	c ratqs final
2284	c 1e4 (en gros 3 heures), en dur pour le moment, est le temps de
2285	c relaxation des ratqs
2286	c facttemps=exp(-pdtphys/1.e4)
2287	facteur=exp(-pdtphys*facttemps)
2288	ratqs(:,:)=max(ratqs(:,:)*facteur,ratqss(:,:))
2289	ratqs(:,:)=max(ratqs(:,:),ratqsc(:,:))
2290	c print*,'calcul des ratqs fini'
2291	else
2292	c on ne prend que le ratqs stable pour fisrtilp
2293	ratqs(:,:)=ratqss(:,:)
2294	endif
2295
2296
2297	c
2298	c Appeler le processus de condensation a grande echelle
2299	c et le processus de precipitation
2300	c-------------------------------------------------------------------------
2301	CALL fisrtilp(dtime,paprs,pplay,
2302	. t_seri, q_seri,ptconv,ratqs,
2303	. d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq,
2304	. rain_lsc, snow_lsc,
2305	. pfrac_impa, pfrac_nucl, pfrac_1nucl,
2306	. frac_impa, frac_nucl,
2307	. prfl, psfl, rhcl)
2308
2309	WHERE (rain_lsc < 0) rain_lsc = 0.
2310	WHERE (snow_lsc < 0) snow_lsc = 0.
2311	DO k = 1, klev
2312	DO i = 1, klon
2313	t_seri(i,k) = t_seri(i,k) + d_t_lsc(i,k)
2314	q_seri(i,k) = q_seri(i,k) + d_q_lsc(i,k)
2315	ql_seri(i,k) = ql_seri(i,k) + d_ql_lsc(i,k)
2316	cldfra(i,k) = rneb(i,k)
2317	IF (.NOT.new_oliq) cldliq(i,k) = ql_seri(i,k)
2318	ENDDO
2319	ENDDO
2320	IF (check) THEN
2321	za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
2322	PRINT*, "apresilp=", za
2323	zx_t = 0.0
2324	za = 0.0
2325	DO i = 1, klon
2326	za = za + paire(i)/FLOAT(klon)
2327	zx_t = zx_t + (rain_lsc(i)+snow_lsc(i))*paire(i)/FLOAT(klon)
2328	ENDDO
2329	zx_t = zx_t/za*dtime
2330	PRINT*, "Precip=", zx_t
2331	ENDIF
2332	c
2333	IF (if_ebil.ge.2) THEN
2334	ztit='after fisrt'
2335	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
2336	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
2337	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
2338	call diagphy(paire,ztit,ip_ebil
2339	e , zero_v, zero_v, zero_v, zero_v, zero_v
2340	e , zero_v, rain_lsc, snow_lsc, ztsol
2341	e , d_h_vcol, d_qt, d_ec
2342	s , fs_bound, fq_bound )
2343	END IF
2344	c
2345	c-------------------------------------------------------------------
2346	c PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
2347	c-------------------------------------------------------------------
2348
2349	c 1. NUAGES CONVECTIFS
2350	c
2351	IF (iflag_cldcon.eq.-1) THEN ! seulement pour Tiedtke
2352
2353	c Nuages diagnostiques pour Tiedtke
2354	CALL diagcld1(paprs,pplay,
2355	. rain_con,snow_con,ibas_con,itop_con,
2356	. diafra,dialiq)
2357	DO k = 1, klev
2358	DO i = 1, klon
2359	IF (diafra(i,k).GT.cldfra(i,k)) THEN
2360	cldliq(i,k) = dialiq(i,k)
2361	cldfra(i,k) = diafra(i,k)
2362	ENDIF
2363	ENDDO
2364	ENDDO
2365
2366	ELSE IF (iflag_cldcon.eq.3) THEN
2367	c On prend pour les nuages convectifs le max du calcul de la
2368	c convection et du calcul du pas de temps précédent diminué d'un facteur
2369	c facttemps
2370	c facttemps=pdtphys/1.e4
2371	facteur = pdtphys *facttemps
2372	do k=1,klev
2373	do i=1,klon
2374	rnebcon(i,k)=rnebcon(i,k)*facteur
2375	if (rnebcon0(i,k)clwcon0(i,k).gt.rnebcon(i,k)clwcon(i,k))
2376	s then
2377	rnebcon(i,k)=rnebcon0(i,k)
2378	clwcon(i,k)=clwcon0(i,k)
2379	endif
2380	enddo
2381	enddo
2382
2383	c On prend la somme des fractions nuageuses et des contenus en eau
2384	cldfra(:,:)=min(max(cldfra(:,:),rnebcon(:,:)),1.)
2385	cldliq(:,:)=cldliq(:,:)+rnebcon(:,:)*clwcon(:,:)
2386
2387
2388	ENDIF
2389	c
2390	c 2. NUAGES STARTIFORMES
2391	c
2392	IF (ok_stratus) THEN
2393	CALL diagcld2(paprs,pplay,t_seri,q_seri, diafra,dialiq)
2394	DO k = 1, klev
2395	DO i = 1, klon
2396	IF (diafra(i,k).GT.cldfra(i,k)) THEN
2397	cldliq(i,k) = dialiq(i,k)
2398	cldfra(i,k) = diafra(i,k)
2399	ENDIF
2400	ENDDO
2401	ENDDO
2402	ENDIF
2403	c
2404	c Precipitation totale
2405	c
2406	DO i = 1, klon
2407	rain_fall(i) = rain_con(i) + rain_lsc(i)
2408	snow_fall(i) = snow_con(i) + snow_lsc(i)
2409	ENDDO
2410	c
2411	IF (if_ebil.ge.2) THEN
2412	ztit="after diagcld"
2413	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
2414	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
2415	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
2416	END IF
2417	c
2418	c Calculer l'humidite relative pour diagnostique
2419	c
2420	DO k = 1, klev
2421	DO i = 1, klon
2422	zx_t = t_seri(i,k)
2423	IF (thermcep) THEN
2424	zdelta = MAX(0.,SIGN(1.,rtt-zx_t))
2425	zx_qs = r2es * FOEEW(zx_t,zdelta)/pplay(i,k)
2426	zx_qs = MIN(0.5,zx_qs)
2427	zcor = 1./(1.-retv*zx_qs)
2428	zx_qs = zx_qs*zcor
2429	ELSE
2430	IF (zx_t.LT.t_coup) THEN
2431	zx_qs = qsats(zx_t)/pplay(i,k)
2432	ELSE
2433	zx_qs = qsatl(zx_t)/pplay(i,k)
2434	ENDIF
2435	ENDIF
2436	zx_rh(i,k) = q_seri(i,k)/zx_qs
2437	zqsat(i,k)=zx_qs
2438	ENDDO
2439	ENDDO
2440	c
2441	c Calculer les parametres optiques des nuages et quelques
2442	c parametres pour diagnostiques:
2443	c
2444	if (ok_newmicro) then
2445	CALL newmicro (paprs, pplay,ok_newmicro,
2446	. t_seri, cldliq, cldfra, cldtau, cldemi,
2447	. cldh, cldl, cldm, cldt, cldq)
2448	else
2449	CALL nuage (paprs, pplay,
2450	. t_seri, cldliq, cldfra, cldtau, cldemi,
2451	. cldh, cldl, cldm, cldt, cldq)
2452	endif
2453	c
2454	c Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
2455	c
2456	IF (MOD(itaprad,radpas).EQ.0) THEN
2457	DO i = 1, klon
2458	albsol(i) = falbe(i,is_oce) * pctsrf(i,is_oce)
2459	. + falbe(i,is_lic) * pctsrf(i,is_lic)
2460	. + falbe(i,is_ter) * pctsrf(i,is_ter)
2461	. + falbe(i,is_sic) * pctsrf(i,is_sic)
2462	albsollw(i) = falblw(i,is_oce) * pctsrf(i,is_oce)
2463	. + falblw(i,is_lic) * pctsrf(i,is_lic)
2464	. + falblw(i,is_ter) * pctsrf(i,is_ter)
2465	. + falblw(i,is_sic) * pctsrf(i,is_sic)
2466	ENDDO
2467	! if (debut) then
2468	! albsol1 = albsol
2469	! albsollw1 = albsollw
2470	! endif
2471	! albsol = albsol1
2472	! albsollw = albsollw1
2473	CALL radlwsw ! nouveau rayonnement (compatible Arpege-IFS)
2474	e (dist, rmu0, fract, co2_ppm, solaire,
2475	e paprs, pplay,zxtsol,albsol, albsollw, t_seri,q_seri,
2476	e wo,
2477	e cldfra, cldemi, cldtau,
2478	s heat,heat0,cool,cool0,radsol,albpla,
2479	s topsw,toplw,solsw,sollw,
2480	s sollwdown,
2481	s topsw0,toplw0,solsw0,sollw0)
2482	itaprad = 0
2483	ENDIF
2484	itaprad = itaprad + 1
2485	c
2486	c Ajouter la tendance des rayonnements (tous les pas)
2487	c
2488	DO k = 1, klev
2489	DO i = 1, klon
2490	t_seri(i,k) = t_seri(i,k)
2491	. + (heat(i,k)-cool(i,k)) * dtime/86400.
2492	ENDDO
2493	ENDDO
2494	c
2495	IF (if_ebil.ge.2) THEN
2496	ztit='after rad'
2497	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
2498	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
2499	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
2500	call diagphy(paire,ztit,ip_ebil
2501	e , topsw, toplw, solsw, sollw, zero_v
2502	e , zero_v, zero_v, zero_v, ztsol
2503	e , d_h_vcol, d_qt, d_ec
2504	s , fs_bound, fq_bound )
2505	END IF
2506	c
2507	c
2508	c Calculer l'hydrologie de la surface
2509	c
2510	c CALL hydrol(dtime,pctsrf,rain_fall, snow_fall, zxevap,
2511	c . agesno, ftsol,fqsol,fsnow, ruis)
2512	c
2513	DO i = 1, klon
2514	zxqsol(i) = 0.0
2515	zxsnow(i) = 0.0
2516	ENDDO
2517	DO nsrf = 1, nbsrf
2518	DO i = 1, klon
2519	zxqsol(i) = zxqsol(i) + fqsol(i,nsrf)*pctsrf(i,nsrf)
2520	zxsnow(i) = zxsnow(i) + fsnow(i,nsrf)*pctsrf(i,nsrf)
2521	ENDDO
2522	ENDDO
2523	c
2524	c Si une sous-fraction n'existe pas, elle prend la valeur moyenne
2525	c
2526	c$$$ DO nsrf = 1, nbsrf
2527	c$$$ DO i = 1, klon
2528	c$$$ IF (pctsrf(i,nsrf).LT.epsfra) THEN
2529	c$$$ fqsol(i,nsrf) = zxqsol(i)
2530	c$$$ fsnow(i,nsrf) = zxsnow(i)
2531	c$$$ ENDIF
2532	c$$$ ENDDO
2533	c$$$ ENDDO
2534	c
2535	c Calculer le bilan du sol et la derive de temperature (couplage)
2536	c
2537	DO i = 1, klon
2538	c bils(i) = radsol(i) - sens(i) - evap(i)*RLVTT
2539	c a la demande de JLD
2540	bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
2541	ENDDO
2542	c
2543	cmoddeblott(jan95)
2544	c Appeler le programme de parametrisation de l'orographie
2545	c a l'echelle sous-maille:
2546	c
2547	IF (ok_orodr) THEN
2548	c
2549	c selection des points pour lesquels le shema est actif:
2550	igwd=0
2551	DO i=1,klon
2552	itest(i)=0
2553	c IF ((zstd(i).gt.10.0)) THEN
2554	IF (((zpic(i)-zmea(i)).GT.100.).AND.(zstd(i).GT.10.0)) THEN
2555	itest(i)=1
2556	igwd=igwd+1
2557	idx(igwd)=i
2558	ENDIF
2559	ENDDO
2560	c igwdim=MAX(1,igwd)
2561	c
2562	CALL drag_noro(klon,klev,dtime,paprs,pplay,
2563	e zmea,zstd, zsig, zgam, zthe,zpic,zval,
2564	e igwd,idx,itest,
2565	e t_seri, u_seri, v_seri,
2566	s zulow, zvlow, zustr, zvstr,
2567	s d_t_oro, d_u_oro, d_v_oro)
2568	c
2569	c ajout des tendances
2570	DO k = 1, klev
2571	DO i = 1, klon
2572	t_seri(i,k) = t_seri(i,k) + d_t_oro(i,k)
2573	u_seri(i,k) = u_seri(i,k) + d_u_oro(i,k)
2574	v_seri(i,k) = v_seri(i,k) + d_v_oro(i,k)
2575	ENDDO
2576	ENDDO
2577	c
2578	ENDIF ! fin de test sur ok_orodr
2579	c
2580	IF (ok_orolf) THEN
2581	c
2582	c selection des points pour lesquels le shema est actif:
2583	igwd=0
2584	DO i=1,klon
2585	itest(i)=0
2586	IF ((zpic(i)-zmea(i)).GT.100.) THEN
2587	itest(i)=1
2588	igwd=igwd+1
2589	idx(igwd)=i
2590	ENDIF
2591	ENDDO
2592	c igwdim=MAX(1,igwd)
2593	c
2594	CALL lift_noro(klon,klev,dtime,paprs,pplay,
2595	e rlat,zmea,zstd,zpic,
2596	e itest,
2597	e t_seri, u_seri, v_seri,
2598	s zulow, zvlow, zustr, zvstr,
2599	s d_t_lif, d_u_lif, d_v_lif)
2600	c
2601	c ajout des tendances
2602	DO k = 1, klev
2603	DO i = 1, klon
2604	t_seri(i,k) = t_seri(i,k) + d_t_lif(i,k)
2605	u_seri(i,k) = u_seri(i,k) + d_u_lif(i,k)
2606	v_seri(i,k) = v_seri(i,k) + d_v_lif(i,k)
2607	ENDDO
2608	ENDDO
2609	c
2610	ENDIF ! fin de test sur ok_orolf
2611	c
2612	IF (if_ebil.ge.2) THEN
2613	ztit='after orography'
2614	CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
2615	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
2616	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
2617	END IF
2618	c
2619	c
2620	cAA
2621	cAA Installation de l'interface online-offline pour traceurs
2622	cAA
2623	c====================================================================
2624	c Calcul des tendances traceurs
2625	c====================================================================
2626	C Pascale : il faut quand meme apeller phytrac car il gere les sorties
2627	cKE43 des traceurs => il faut donc mettre des flags a .false.
2628	IF (iflag_con.GE.3) THEN
2629	c on ajoute les tendances calculees par KE43
2630	c$$$ OM on onhibe la convection sur les traceurs
2631	DO iq=1, nqmax-2 ! Sandrine a -3 ???
2632	c$$$ OM on inhibe la convection sur les traceur
2633	c$$$ DO k = 1, nlev
2634	c$$$ DO i = 1, klon
2635	c$$$ tr_seri(i,k,iq) = tr_seri(i,k,iq) + d_tr(i,k,iq)
2636	c$$$ ENDDO
2637	c$$$ ENDDO
2638	WRITE(iqn,'(i2.2)') iq
2639	CALL minmaxqfi(tr_seri(1,1,iq),0.,1.e33,'couche lim iq='//iqn)
2640	ENDDO
2641	CMAF modif pour garder info du nombre de traceurs auxquels
2642	C la physique s'applique
2643	ELSE
2644	CMAF modif pour garder info du nombre de traceurs auxquels
2645	C la physique s'applique
2646	C
2647	call phytrac (rnpb,
2648	I debut,lafin,
2649	I nqmax-2,
2650	I nlon,nlev,dtime,
2651	I t,paprs,pplay,
2652	I pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
2653	I ycoefh,yu1,yv1,ftsol,pctsrf,rlat,
2654	I frac_impa, frac_nucl,
2655	I rlon,presnivs,paire,pphis,
2656	O tr_seri)
2657	ENDIF
2658
2659	IF (offline) THEN
2660
2661	call phystokenc (
2662	I nlon,nlev,pdtphys,rlon,rlat,
2663	I t,pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
2664	I ycoefh,yu1,yv1,ftsol,pctsrf,
2665	I frac_impa, frac_nucl,
2666	I pphis,paire,dtime,itap)
2667
2668
2669	ENDIF
2670
2671	c
2672	c Calculer le transport de l'eau et de l'energie (diagnostique)
2673	c
2674	CALL transp (paprs,zxtsol,
2675	e t_seri, q_seri, u_seri, v_seri, zphi,
2676	s ve, vq, ue, uq)
2677	c
2678	c Accumuler les variables a stocker dans les fichiers histoire:
2679	c
2680	c
2681	c
2682	IF (if_ebil.ge.1) THEN
2683	ztit='after physic'
2684	CALL diagetpq(paire,ztit,ip_ebil,1,1,dtime
2685	e , t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,paprs,pplay
2686	s , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
2687	C Comme les tendances de la physique sont ajoute dans la dynamique,
2688	C on devrait avoir que la variation d'entalpie par la dynamique
2689	C est egale a la variation de la physique au pas de temps precedent.
2690	C Donc la somme de ces 2 variations devrait etre nulle.
2691	call diagphy(paire,ztit,ip_ebil
2692	e , topsw, toplw, solsw, sollw, sens
2693	e , evap, rain_fall, snow_fall, ztsol
2694	e , d_h_vcol, d_qt, d_ec
2695	s , fs_bound, fq_bound )
2696	C
2697	d_h_vcol_phy=d_h_vcol
2698	C
2699	END IF
2700	C
2701	IF (ok_journe) THEN
2702	c
2703	ndex2d = 0
2704	ndex3d = 0
2705	c
2706	c Champs 2D:
2707	c
2708	zsto = dtime
2709	zout = dtime * FLOAT(ecrit_day)
2710	itau_w = itau_phy + itap
2711
2712	i = NINT(zout/zsto)
2713	CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
2714	CALL histwrite(nid_day,"phis",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2715	varname = 'phis'
2716	vartitle= 'Surface geop. height'
2717	varunits= '-'
2718	c call writephy(fid_day,prof2d_on,varname,pphis,vartitle,
2719	c . varunits)
2720	c
2721	i = NINT(zout/zsto)
2722	CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
2723	CALL histwrite(nid_day,"aire",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2724	varname = 'aire'
2725	vartitle= 'Grid area'
2726	varunits= '-'
2727	c call writephy(fid_day,prof2d_on,varname,paire,vartitle,
2728	c . varunits)
2729	C
2730	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxtsol,zx_tmp_2d)
2731	CALL histwrite(nid_day,"tsol",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2732	c call writephy(fid_day,prof2d_av,'tsol',zxtsol,
2733	c . 'Surface Temperature','K')
2734	c
2735	C
2736	zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, is_ter)
2737	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d ,zx_tmp_2d)
2738	CALL histwrite(nid_day,"tter",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2739	c call writephy(fid_day,prof2d_av,'tter',ftsol(1 : klon, is_ter),
2740	c . 'Surface Temperature','K')
2741	C
2742	zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, is_lic)
2743	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2744	CALL histwrite(nid_day,"tlic",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2745	c call writephy(fid_day,prof2d_av,'tlic',ftsol(1 : klon, is_lic),
2746	c . 'Surface Temperature','K')
2747	C
2748	zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, is_oce)
2749	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2750	CALL histwrite(nid_day,"toce",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2751	c call writephy(fid_day,prof2d_av,'toce',ftsol(1 : klon, is_oce),
2752	c . 'Surface Temperature','K')
2753	C
2754	zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, is_sic)
2755	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2756	CALL histwrite(nid_day,"tsic",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2757	c call writephy(fid_day,prof2d_av,'tsic',ftsol(1 : klon, is_sic),
2758	c . 'Surface Temperature','K')
2759	C
2760	DO i = 1, klon
2761	zx_tmp_fi2d(i) = paprs(i,1)
2762	ENDDO
2763	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2764	CALL histwrite(nid_day,"psol",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2765	c Essai writephys
2766	varname = 'psol'
2767	vartitle= 'pression au sol'
2768	varunits= 'hPa'
2769	c call writephy(fid_day,prof2d_av,varname,zx_tmp_fi2d,vartitle,
2770	c . varunits)
2771	c
2772	DO i = 1, klon
2773	zx_tmp_fi2d(i) = (rain_fall(i) + snow_fall(i))
2774	ENDDO
2775	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2776	CALL histwrite(nid_day,"precip",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2777	c call writephy(fid_day,prof2d_av,'rain',zx_tmp_fi2d,
2778	c . 'Precipitation','mm/day')
2779
2780
2781	c
2782	CALL gr_fi_ecrit(1, klon,iim,jjmp1, snow_fall,zx_tmp_2d)
2783	CALL histwrite(nid_day,"snow",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2784	c call writephy(fid_day,prof2d_av,'snow',snow_fall,
2785	c . 'Snow','mm/day')
2786	c
2787	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxsnow,zx_tmp_2d)
2788	CALL histwrite(nid_day,"snow_mass",itau_w,zx_tmp_2d,iim*jjmp1,
2789	. ndex2d)
2790	c call writephy(fid_day,prof2d_av,'snow_mass',zxsnow,
2791	c . 'Snow cover','mm')
2792	c
2793	CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
2794	CALL histwrite(nid_day,"evap",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2795	c call writephy(fid_day,prof2d_av,'evap',evap,
2796	c . 'Evaporation','mm/day')
2797	c
2798	CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw,zx_tmp_2d)
2799	CALL histwrite(nid_day,"tops",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2800	c call writephy(fid_day,prof2d_av,'tops',topsw,
2801	c . 'Solar rad. at TOA','W/m2')
2802	c
2803	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
2804	CALL histwrite(nid_day,"topl",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2805	c call writephy(fid_day,prof2d_av,'topl',toplw,
2806	c . 'IR rad. at TOA','W/m2')
2807	c
2808	CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
2809	CALL histwrite(nid_day,"sols",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2810	c call writephy(fid_day,prof2d_av,'sols',solsw,
2811	c . 'Solar rad. at surf.','W/m2')
2812	c
2813	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
2814	CALL histwrite(nid_day,"soll",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2815	c call writephy(fid_day,prof2d_av,'soll',sollw,
2816	c . 'IR rad. at surface','W/m2')
2817	c
2818	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollwdown,zx_tmp_2d)
2819	CALL histwrite(nid_day,"solldown",itau_w,zx_tmp_2d,iim*jjmp1,
2820	. ndex2d)
2821	c call writephy(fid_day,prof2d_av,'solldown',sollwdown,
2822	c . 'Down. IR rad. at surface','W/m2')
2823	c
2824	CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
2825	CALL histwrite(nid_day,"bils",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2826	c call writephy(fid_day,prof2d_av,'bils',bils,
2827	c . 'Surf. total heat flux','W/m2')
2828	c
2829	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
2830	CALL histwrite(nid_day,"sens",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2831	c call writephy(fid_day,prof2d_av,'sens',sens,
2832	c . 'Sensible heat flux','W/m2')
2833	c
2834	CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
2835	CALL histwrite(nid_day,"fder",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2836	c call writephy(fid_day,prof2d_av,'fder',fder,
2837	c . 'Heat flux derivation','W/m2')
2838	c
2839	c
2840	DO nsrf = 1, nbsrf
2841	C§§§
2842	zx_tmp_fi2d(1 : klon) = pctsrf( 1 : klon, nsrf)
2843	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2844	CALL histwrite(nid_day,"pourc_"//clnsurf(nsrf),itau_w,
2845	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2846	c call writephy(fid_day,prof2d_av,'pourc_'//clnsurf(nsrf),
2847	c . pctsrf( 1 : klon, nsrf),
2848	c . 'Fraction'//clnsurf(nsrf),'-')
2849	C
2850	zx_tmp_fi2d(1 : klon) = ftsol( 1 : klon, nsrf)
2851	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2852	CALL histwrite(nid_day,"tsol_"//clnsurf(nsrf),itau_w,
2853	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2854	c call writephy(fid_day,prof2d_av,'tsol_'//clnsurf(nsrf),
2855	c . ftsol( 1 : klon, nsrf),
2856	c . 'Surf. Temp'//clnsurf(nsrf),'K')
2857	C
2858	zx_tmp_fi2d(1 : klon) = fluxt( 1 : klon, 1, nsrf)
2859	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2860	CALL histwrite(nid_day,"sens_"//clnsurf(nsrf),itau_w,
2861	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2862	c call writephy(fid_day,prof2d_av,'sens_'//clnsurf(nsrf),
2863	c . fluxt( 1 : klon, 1, nsrf),
2864	c . 'Sensible heat flux '//clnsurf(nsrf),'W/m2')
2865	C
2866	zx_tmp_fi2d(1 : klon) = fluxlat( 1 : klon, nsrf)
2867	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2868	CALL histwrite(nid_day,"lat_"//clnsurf(nsrf),itau_w,
2869	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2870	c call writephy(fid_day,prof2d_av,'lat_'//clnsurf(nsrf),
2871	c . fluxlat( 1 : klon, nsrf),
2872	c . 'Latent heat flux '//clnsurf(nsrf),'W/m2')
2873	C
2874	zx_tmp_fi2d(1 : klon) = fluxu( 1 : klon, 1, nsrf)
2875	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2876	CALL histwrite(nid_day,"taux_"//clnsurf(nsrf),itau_w,
2877	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2878	c call writephy(fid_day,prof2d_av,'taux_'//clnsurf(nsrf),
2879	c . fluxu( 1 : klon, 1, nsrf),
2880	c . 'Zonal wind stress '//clnsurf(nsrf),'Pa')
2881	C
2882	zx_tmp_fi2d(1 : klon) = fluxv( 1 : klon, 1, nsrf)
2883	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2884	CALL histwrite(nid_day,"tauy_"//clnsurf(nsrf),itau_w,
2885	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2886	c call writephy(fid_day,prof2d_av,'tauy_'//clnsurf(nsrf),
2887	c . fluxv( 1 : klon, 1, nsrf),
2888	c . 'Meridional wind stress '//clnsurf(nsrf),'Pa')
2889	C
2890	zx_tmp_fi2d(1 : klon) = falbe( 1 : klon, nsrf)
2891	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2892	CALL histwrite(nid_day,"albe_"//clnsurf(nsrf),itau_w,
2893	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2894	c call writephy(fid_day,prof2d_av,'albe_'//clnsurf(nsrf),
2895	c . falbe( 1 : klon, nsrf),
2896	c . 'Albedo surf. SW'//clnsurf(nsrf),'-')
2897	c call writephy(fid_day,prof2d_av,'alblw_'//clnsurf(nsrf),
2898	c . falblw( 1 : klon, nsrf),
2899	c . 'Albedo surf. LW'//clnsurf(nsrf),'-')
2900	C
2901	zx_tmp_fi2d(1 : klon) = frugs( 1 : klon, nsrf)
2902	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
2903	CALL histwrite(nid_day,"rugs_"//clnsurf(nsrf),itau_w,
2904	$ zx_tmp_2d,iim*jjmp1,ndex2d)
2905	c call writephy(fid_day,prof2d_av,'rugs_'//clnsurf(nsrf),
2906	c . frugs( 1 : klon, nsrf),
2907	c . 'Rugosity '//clnsurf(nsrf),' - ')
2908	C
2909	END DO
2910	C
2911	c$$$ DO i = 1, klon
2912	c$$$ zx_tmp_fi2d(i) = pctsrf(i,is_sic)
2913	c$$$ ENDDO
2914	c$$$ CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2915	c$$$ CALL histwrite(nid_day,"sicf",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2916	c
2917	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldl,zx_tmp_2d)
2918	CALL histwrite(nid_day,"cldl",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2919	c call writephy(fid_day,prof2d_av,'cldl',cldl,
2920	c . 'Low-level cloudiness','-')
2921	c
2922	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldm,zx_tmp_2d)
2923	CALL histwrite(nid_day,"cldm",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2924	c call writephy(fid_day,prof2d_av,'cldm',cldm,
2925	c . 'Mid-level cloudiness','-')
2926	c
2927	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldh,zx_tmp_2d)
2928	CALL histwrite(nid_day,"cldh",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2929	c call writephy(fid_day,prof2d_av,'cldh',cldh,
2930	c . 'High-level cloudiness','-')
2931	c
2932	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldt,zx_tmp_2d)
2933	CALL histwrite(nid_day,"cldt",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2934	c call writephy(fid_day,prof2d_av,'cldt',cldt,
2935	c . 'Total cloudiness','-')
2936	c
2937	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldq,zx_tmp_2d)
2938	CALL histwrite(nid_day,"cldq",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
2939	c call writephy(fid_day,prof2d_av,'cldq',cldq,
2940	c . 'Cloud liquid water path','-')
2941	c
2942	c Champs 3D:
2943	c
2944	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
2945	CALL histwrite(nid_day,"temp",itau_w,zx_tmp_3d,
2946	. iimjjmp1klev,ndex3d)
2947	c Essai writephys
2948	varname = 'temp'
2949	vartitle= 'temperature 3D'
2950	varunits= 'K'
2951	c call writephy(fid_day,prof3d_av,varname,t_seri,vartitle,varunits)
2952	c
2953	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,ivap), zx_tmp_3d)
2954	CALL histwrite(nid_day,"ovap",itau_w,zx_tmp_3d,
2955	. iimjjmp1klev,ndex3d)
2956	c call writephy(fid_day,prof3d_av,'ovap',qx(1,1,ivap),
2957	c . 'Specific humidity','Kg/Kg')
2958	c
2959	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
2960	CALL histwrite(nid_day,"geop",itau_w,zx_tmp_3d,
2961	. iimjjmp1klev,ndex3d)
2962	c call writephy(fid_day,prof3d_av,'geop',zphi,
2963	c . 'Geopotential height','m')
2964	c
2965	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
2966	CALL histwrite(nid_day,"vitu",itau_w,zx_tmp_3d,
2967	. iimjjmp1klev,ndex3d)
2968	c call writephy(fid_day,prof3d_av,'vitu',u_seri,
2969	c . 'Zonal wind','m/s')
2970	c
2971	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
2972	CALL histwrite(nid_day,"vitv",itau_w,zx_tmp_3d,
2973	. iimjjmp1klev,ndex3d)
2974	c call writephy(fid_day,prof3d_av,'vitv',v_seri,
2975	c . 'Meridional wind','m/s')
2976	c
2977	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, omega, zx_tmp_3d)
2978	CALL histwrite(nid_day,"vitw",itau_w,zx_tmp_3d,
2979	. iimjjmp1klev,ndex3d)
2980	c call writephy(fid_day,prof3d_av,'vitw',omega,
2981	c . 'Vertical wind','m/s')
2982	c
2983	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
2984	CALL histwrite(nid_day,"pres",itau_w,zx_tmp_3d,
2985	. iimjjmp1klev,ndex3d)
2986	c call writephy(fid_day,prof3d_av,'pres',pplay,
2987	c . 'Air pressure','Pa')
2988
2989	c
2990	if (ok_sync) then
2991	c call writephy_sync(fid_day)
2992	call histsync(nid_day)
2993	endif
2994	ENDIF
2995	C
2996	IF (ok_mensuel) THEN
2997	c
2998	ndex2d = 0
2999	ndex3d = 0
3000	c
3001	c Champs 2D:
3002	c
3003	zsto = dtime
3004	zout = dtime * ecrit_mth
3005	itau_w = itau_phy + itap
3006
3007	i = NINT(zout/zsto)
3008	CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
3009	CALL histwrite(nid_mth,"phis",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3010	C
3011	i = NINT(zout/zsto)
3012	CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
3013	CALL histwrite(nid_mth,"aire",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3014
3015	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxtsol,zx_tmp_2d)
3016	CALL histwrite(nid_mth,"tsol",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3017	c
3018	DO i = 1, klon
3019	zx_tmp_fi2d(i) = paprs(i,1)
3020	ENDDO
3021	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3022	CALL histwrite(nid_mth,"psol",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3023	c
3024	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxqsol,zx_tmp_2d)
3025	CALL histwrite(nid_mth,"qsol",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3026	c
3027	DO i = 1, klon
3028	zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
3029	ENDDO
3030	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3031	CALL histwrite(nid_mth,"precip",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3032	c
3033	DO i = 1, klon
3034	zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
3035	ENDDO
3036	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3037	CALL histwrite(nid_mth,"plul",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3038	c
3039	DO i = 1, klon
3040	zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
3041	ENDDO
3042	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3043	CALL histwrite(nid_mth,"pluc",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3044	c
3045	CALL gr_fi_ecrit(1, klon,iim,jjmp1, snow_fall,zx_tmp_2d)
3046	CALL histwrite(nid_mth,"snow",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3047	c
3048	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxsnow,zx_tmp_2d)
3049	CALL histwrite(nid_mth,"snow_mass",itau_w,zx_tmp_2d,iim*jjmp1,
3050	. ndex2d)
3051	c
3052	CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
3053	CALL histwrite(nid_mth,"evap",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3054	c
3055	CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw,zx_tmp_2d)
3056	CALL histwrite(nid_mth,"tops",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3057	c
3058	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
3059	CALL histwrite(nid_mth,"topl",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3060	c
3061	CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
3062	CALL histwrite(nid_mth,"sols",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3063	c
3064	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
3065	CALL histwrite(nid_mth,"soll",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3066	c
3067	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollwdown,zx_tmp_2d)
3068	CALL histwrite(nid_mth,"solldown",itau_w,zx_tmp_2d,iim*jjmp1,
3069	. ndex2d)
3070	c
3071	CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw0,zx_tmp_2d)
3072	CALL histwrite(nid_mth,"tops0",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3073	c
3074	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw0,zx_tmp_2d)
3075	CALL histwrite(nid_mth,"topl0",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3076	c
3077	CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw0,zx_tmp_2d)
3078	CALL histwrite(nid_mth,"sols0",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3079	c
3080	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw0,zx_tmp_2d)
3081	CALL histwrite(nid_mth,"soll0",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3082	c
3083	CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
3084	CALL histwrite(nid_mth,"bils",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3085	c
3086	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
3087	CALL histwrite(nid_mth,"sens",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3088	c
3089	CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
3090	CALL histwrite(nid_mth,"fder",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3091	c
3092	c
3093	c DO i = 1, klon
3094	c zx_tmp_fi2d(i) = fluxu(i,1)
3095	c ENDDO
3096	c CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3097	c CALL histwrite(nid_mth,"frtu",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3098	c
3099	c DO i = 1, klon
3100	c zx_tmp_fi2d(i) = fluxv(i,1)
3101	c ENDDO
3102	c CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3103	c CALL histwrite(nid_mth,"frtv",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3104	c
3105	DO nsrf = 1, nbsrf
3106	C§§§
3107	zx_tmp_fi2d(1 : klon) = pctsrf( 1 : klon, nsrf)
3108	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3109	CALL histwrite(nid_mth,"pourc_"//clnsurf(nsrf),itau_w,
3110	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3111	C
3112	zx_tmp_fi2d(1 : klon) = ftsol( 1 : klon, nsrf)
3113	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3114	CALL histwrite(nid_mth,"tsol_"//clnsurf(nsrf),itau_w,
3115	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3116	C
3117	zx_tmp_fi2d(1 : klon) = fluxt( 1 : klon, 1, nsrf)
3118	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3119	CALL histwrite(nid_mth,"sens_"//clnsurf(nsrf),itau_w,
3120	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3121	C
3122	zx_tmp_fi2d(1 : klon) = fluxlat( 1 : klon, nsrf)
3123	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3124	CALL histwrite(nid_mth,"lat_"//clnsurf(nsrf),itau_w,
3125	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3126	C
3127	zx_tmp_fi2d(1 : klon) = fluxu( 1 : klon, 1, nsrf)
3128	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3129	CALL histwrite(nid_mth,"taux_"//clnsurf(nsrf),itau_w,
3130	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3131	C
3132	zx_tmp_fi2d(1 : klon) = fluxv( 1 : klon, 1, nsrf)
3133	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3134	CALL histwrite(nid_mth,"tauy_"//clnsurf(nsrf),itau_w,
3135	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3136	C
3137	zx_tmp_fi2d(1 : klon) = falbe( 1 : klon, nsrf)
3138	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3139	CALL histwrite(nid_mth,"albe_"//clnsurf(nsrf),itau_w,
3140	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3141	C
3142	zx_tmp_fi2d(1 : klon) = frugs( 1 : klon, nsrf)
3143	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3144	CALL histwrite(nid_mth,"rugs_"//clnsurf(nsrf),itau_w,
3145	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3146	c
3147	zx_tmp_fi2d(1 : klon) = agesno( 1 : klon, nsrf)
3148	CALL gr_fi_ecrit(1, klon,iim,jjmp1, agesno,zx_tmp_2d)
3149	CALL histwrite(nid_mth,"ages_"//clnsurf(nsrf),itau_w
3150	$ ,zx_tmp_2d,iim*jjmp1,ndex2d)
3151
3152	END DO
3153	c$$$ DO i = 1, klon
3154	c$$$ zx_tmp_fi2d(i) = pctsrf(i,is_sic)
3155	c$$$ ENDDO
3156	c$$$ CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3157	c$$$ CALL histwrite(nid_mth,"sicf",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3158	c
3159	CALL gr_fi_ecrit(1, klon,iim,jjmp1, albsol,zx_tmp_2d)
3160	CALL histwrite(nid_mth,"albs",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3161	CALL gr_fi_ecrit(1, klon,iim,jjmp1, albsollw,zx_tmp_2d)
3162	CALL histwrite(nid_mth,"albslw",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3163	c
3164	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragm,zx_tmp_2d)
3165	CALL histwrite(nid_mth,"cdrm",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3166	c
3167	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragh,zx_tmp_2d)
3168	CALL histwrite(nid_mth,"cdrh",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3169	c
3170	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldl,zx_tmp_2d)
3171	CALL histwrite(nid_mth,"cldl",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3172	c
3173	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldm,zx_tmp_2d)
3174	CALL histwrite(nid_mth,"cldm",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3175	c
3176	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldh,zx_tmp_2d)
3177	CALL histwrite(nid_mth,"cldh",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3178	c
3179	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldt,zx_tmp_2d)
3180	CALL histwrite(nid_mth,"cldt",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3181	c
3182	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldq,zx_tmp_2d)
3183	CALL histwrite(nid_mth,"cldq",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3184	c
3185	CALL gr_fi_ecrit(1, klon,iim,jjmp1, ue,zx_tmp_2d)
3186	CALL histwrite(nid_mth,"ue",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3187	c
3188	CALL gr_fi_ecrit(1, klon,iim,jjmp1, ve,zx_tmp_2d)
3189	CALL histwrite(nid_mth,"ve",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3190	c
3191	CALL gr_fi_ecrit(1, klon,iim,jjmp1, uq,zx_tmp_2d)
3192	CALL histwrite(nid_mth,"uq",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3193	c
3194	CALL gr_fi_ecrit(1, klon,iim,jjmp1, vq,zx_tmp_2d)
3195	CALL histwrite(nid_mth,"vq",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3196	cKE43
3197	IF (iflag_con .GE. 3) THEN ! sb
3198	c
3199	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cape,zx_tmp_2d)
3200	CALL histwrite(nid_mth,"cape",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3201	c
3202	CALL gr_fi_ecrit(1, klon,iim,jjmp1,pbase,zx_tmp_2d)
3203	CALL histwrite(nid_mth,"pbase",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3204	c
3205	CALL gr_fi_ecrit(1, klon,iim,jjmp1,ema_pct,zx_tmp_2d)
3206	CALL histwrite(nid_mth,"ptop",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3207	c
3208	CALL gr_fi_ecrit(1, klon,iim,jjmp1,ema_cbmf,zx_tmp_2d)
3209	CALL histwrite(nid_mth,"fbase",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3210	c
3211	c
3212	ENDIF
3213	c34EK
3214	c
3215	c Champs 3D:
3216	C
3217	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
3218	CALL histwrite(nid_mth,"temp",itau_w,zx_tmp_3d,
3219	. iimjjmp1klev,ndex3d)
3220	c
3221	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,ivap), zx_tmp_3d)
3222	CALL histwrite(nid_mth,"ovap",itau_w,zx_tmp_3d,
3223	. iimjjmp1klev,ndex3d)
3224	c
3225	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
3226	CALL histwrite(nid_mth,"geop",itau_w,zx_tmp_3d,
3227	. iimjjmp1klev,ndex3d)
3228	c
3229	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
3230	CALL histwrite(nid_mth,"vitu",itau_w,zx_tmp_3d,
3231	. iimjjmp1klev,ndex3d)
3232	c
3233	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
3234	CALL histwrite(nid_mth,"vitv",itau_w,zx_tmp_3d,
3235	. iimjjmp1klev,ndex3d)
3236	c
3237	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, omega, zx_tmp_3d)
3238	CALL histwrite(nid_mth,"vitw",itau_w,zx_tmp_3d,
3239	. iimjjmp1klev,ndex3d)
3240	c
3241	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
3242	CALL histwrite(nid_mth,"pres",itau_w,zx_tmp_3d,
3243	. iimjjmp1klev,ndex3d)
3244	c
3245	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldfra, zx_tmp_3d)
3246	CALL histwrite(nid_mth,"rneb",itau_w,zx_tmp_3d,
3247	. iimjjmp1klev,ndex3d)
3248	c
3249	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zx_rh, zx_tmp_3d)
3250	CALL histwrite(nid_mth,"rhum",itau_w,zx_tmp_3d,
3251	. iimjjmp1klev,ndex3d)
3252	c
3253	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldliq, zx_tmp_3d)
3254	CALL histwrite(nid_mth,"oliq",itau_w,zx_tmp_3d,
3255	. iimjjmp1klev,ndex3d)
3256	c
3257	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, clwcon0, zx_tmp_3d)
3258	CALL histwrite(nid_mth,"clwcon",itau_w,zx_tmp_3d,
3259	. iimjjmp1klev,ndex3d)
3260	c
3261	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_dyn, zx_tmp_3d)
3262	CALL histwrite(nid_mth,"dtdyn",itau_w,zx_tmp_3d,
3263	. iimjjmp1klev,ndex3d)
3264	c
3265	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_dyn, zx_tmp_3d)
3266	CALL histwrite(nid_mth,"dqdyn",itau_w,zx_tmp_3d,
3267	. iimjjmp1klev,ndex3d)
3268	c
3269	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_con, zx_tmp_3d)
3270	CALL histwrite(nid_mth,"dtcon",itau_w,zx_tmp_3d,
3271	. iimjjmp1klev,ndex3d)
3272	c
3273	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_con, zx_tmp_3d)
3274	CALL histwrite(nid_mth,"dqcon",itau_w,zx_tmp_3d,
3275	. iimjjmp1klev,ndex3d)
3276	c
3277	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_lsc, zx_tmp_3d)
3278	CALL histwrite(nid_mth,"dtlsc",itau_w,zx_tmp_3d,
3279	. iimjjmp1klev,ndex3d)
3280	c
3281	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_lsc, zx_tmp_3d)
3282	CALL histwrite(nid_mth,"dqlsc",itau_w,zx_tmp_3d,
3283	. iimjjmp1klev,ndex3d)
3284	c
3285	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_vdf, zx_tmp_3d)
3286	CALL histwrite(nid_mth,"dtvdf",itau_w,zx_tmp_3d,
3287	. iimjjmp1klev,ndex3d)
3288	c
3289	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_vdf, zx_tmp_3d)
3290	CALL histwrite(nid_mth,"dqvdf",itau_w,zx_tmp_3d,
3291	. iimjjmp1klev,ndex3d)
3292	c
3293	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_eva, zx_tmp_3d)
3294	CALL histwrite(nid_mth,"dteva",itau_w,zx_tmp_3d,
3295	. iimjjmp1klev,ndex3d)
3296	c
3297	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_eva, zx_tmp_3d)
3298	CALL histwrite(nid_mth,"dqeva",itau_w,zx_tmp_3d,
3299	. iimjjmp1klev,ndex3d)
3300	c
3301	zpt_conv = 0.
3302	where (ptconv) zpt_conv = 1.
3303	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zpt_conv, zx_tmp_3d)
3304	CALL histwrite(nid_mth,"ptconv",itau_w,zx_tmp_3d,
3305	. iim(jjmp1)klev,ndex3d)
3306	c
3307	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, ratqs, zx_tmp_3d)
3308	CALL histwrite(nid_mth,"ratqs",itau_w,zx_tmp_3d,
3309	. iim(jjmp1)klev,ndex3d)
3310	c
3311	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_ajs, zx_tmp_3d)
3312	CALL histwrite(nid_mth,"dtajs",itau_w,zx_tmp_3d,
3313	. iimjjmp1klev,ndex3d)
3314	c
3315	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_ajs, zx_tmp_3d)
3316	CALL histwrite(nid_mth,"dqajs",itau_w,zx_tmp_3d,
3317	. iimjjmp1klev,ndex3d)
3318	c
3319	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, heat, zx_tmp_3d)
3320	CALL histwrite(nid_mth,"dtswr",itau_w,zx_tmp_3d,
3321	. iimjjmp1klev,ndex3d)
3322	c
3323	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, heat0, zx_tmp_3d)
3324	CALL histwrite(nid_mth,"dtsw0",itau_w,zx_tmp_3d,
3325	. iimjjmp1klev,ndex3d)
3326	c
3327	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cool, zx_tmp_3d)
3328	CALL histwrite(nid_mth,"dtlwr",itau_w,zx_tmp_3d,
3329	. iimjjmp1klev,ndex3d)
3330	c
3331	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cool0, zx_tmp_3d)
3332	CALL histwrite(nid_mth,"dtlw0",itau_w,zx_tmp_3d,
3333	. iimjjmp1klev,ndex3d)
3334	c
3335	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_vdf, zx_tmp_3d)
3336	CALL histwrite(nid_mth,"duvdf",itau_w,zx_tmp_3d,
3337	. iimjjmp1klev,ndex3d)
3338	c
3339	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_vdf, zx_tmp_3d)
3340	CALL histwrite(nid_mth,"dvvdf",itau_w,zx_tmp_3d,
3341	. iimjjmp1klev,ndex3d)
3342	c
3343	IF (ok_orodr) THEN
3344	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_oro, zx_tmp_3d)
3345	CALL histwrite(nid_mth,"duoro",itau_w,zx_tmp_3d,
3346	. iimjjmp1klev,ndex3d)
3347	c
3348	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_oro, zx_tmp_3d)
3349	CALL histwrite(nid_mth,"dvoro",itau_w,zx_tmp_3d,
3350	. iimjjmp1klev,ndex3d)
3351	c
3352	ENDIF
3353	C
3354	IF (ok_orolf) THEN
3355	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_lif, zx_tmp_3d)
3356	CALL histwrite(nid_mth,"dulif",itau_w,zx_tmp_3d,
3357	. iimjjmp1klev,ndex3d)
3358	c
3359	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_lif, zx_tmp_3d)
3360	CALL histwrite(nid_mth,"dvlif",itau_w,zx_tmp_3d,
3361	. iimjjmp1klev,ndex3d)
3362	ENDIF
3363	C
3364	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, wo, zx_tmp_3d)
3365	CALL histwrite(nid_mth,"ozone",itau_w,zx_tmp_3d,
3366	. iimjjmp1klev,ndex3d)
3367	c
3368	IF (nqmax.GE.3) THEN
3369	DO iq=1,nqmax-2
3370	IF (iq.LE.99) THEN
3371	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,iq+2), zx_tmp_3d)
3372	WRITE(str2,'(i2.2)') iq
3373	CALL histwrite(nid_mth,"trac"//str2,itau_w,zx_tmp_3d,
3374	. iimjjmp1klev,ndex3d)
3375	ELSE
3376	PRINT*, "Trop de traceurs"
3377	CALL abort
3378	ENDIF
3379	ENDDO
3380	ENDIF
3381	cKE43
3382	IF (iflag_con.GE.3) THEN ! (sb)
3383	c
3384	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, upwd, zx_tmp_3d)
3385	CALL histwrite(nid_mth,"upwd",itau_w,zx_tmp_3d,
3386	. iimjjmp1klev,ndex3d)
3387	c
3388	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, dnwd, zx_tmp_3d)
3389	CALL histwrite(nid_mth,"dnwd",itau_w,zx_tmp_3d,
3390	. iimjjmp1klev,ndex3d)
3391	c
3392	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, dnwd0, zx_tmp_3d)
3393	CALL histwrite(nid_mth,"dnwd0",itau_w,zx_tmp_3d,
3394	. iimjjmp1klev,ndex3d)
3395	c
3396	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, Ma, zx_tmp_3d)
3397	CALL histwrite(nid_mth,"Ma",itau_w,zx_tmp_3d,
3398	. iimjjmp1klev,ndex3d)
3399	c
3400	c
3401	ENDIF
3402	c34EK
3403	c
3404	if (ok_sync) then
3405	call histsync(nid_mth)
3406	endif
3407	ENDIF
3408	c
3409	IF (ok_instan) THEN
3410	c
3411	ndex2d = 0
3412	ndex3d = 0
3413	c
3414	c Champs 2D:
3415	c
3416	zsto = dtime * ecrit_ins
3417	zout = dtime * ecrit_ins
3418	itau_w = itau_phy + itap
3419
3420	i = NINT(zout/zsto)
3421	CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
3422	CALL histwrite(nid_ins,"phis",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3423	c
3424	i = NINT(zout/zsto)
3425	CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
3426	CALL histwrite(nid_ins,"aire",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3427
3428	DO i = 1, klon
3429	zx_tmp_fi2d(i) = paprs(i,1)
3430	ENDDO
3431	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3432	CALL histwrite(nid_ins,"psol",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3433	c
3434	DO i = 1, klon
3435	zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
3436	ENDDO
3437	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3438	CALL histwrite(nid_ins,"precip",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3439	c
3440	DO i = 1, klon
3441	zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
3442	ENDDO
3443	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3444	CALL histwrite(nid_ins,"plul",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3445	c
3446	DO i = 1, klon
3447	zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
3448	ENDDO
3449	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
3450	CALL histwrite(nid_ins,"pluc",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3451
3452	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxtsol,zx_tmp_2d)
3453	CALL histwrite(nid_ins,"tsol",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3454	c
3455	CALL gr_fi_ecrit(1, klon,iim,jjmp1, snow_fall,zx_tmp_2d)
3456	CALL histwrite(nid_ins,"snow",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3457
3458	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragm,zx_tmp_2d)
3459	CALL histwrite(nid_ins,"cdrm",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3460	c
3461	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragh,zx_tmp_2d)
3462	CALL histwrite(nid_ins,"cdrh",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3463	c
3464	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
3465	CALL histwrite(nid_ins,"topl",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3466	c
3467	CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
3468	CALL histwrite(nid_ins,"evap",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3469	c
3470	CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
3471	CALL histwrite(nid_ins,"sols",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3472	c
3473	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
3474	CALL histwrite(nid_ins,"soll",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3475	c
3476	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollwdown,zx_tmp_2d)
3477	CALL histwrite(nid_ins,"solldown",itau_w,zx_tmp_2d,iim*jjmp1,
3478	. ndex2d)
3479	c
3480	CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
3481	CALL histwrite(nid_ins,"bils",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3482	c
3483	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
3484	CALL histwrite(nid_ins,"sens",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3485	c
3486	CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
3487	CALL histwrite(nid_ins,"fder",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3488	c
3489	CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_oce),zx_tmp_2d)
3490	CALL histwrite(nid_ins,"dtsvdfo",itau_w,zx_tmp_2d,iim*jjmp1,
3491	. ndex2d)
3492	c
3493	CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_ter),zx_tmp_2d)
3494	CALL histwrite(nid_ins,"dtsvdft",itau_w,zx_tmp_2d,iim*jjmp1,
3495	. ndex2d)
3496	c
3497	CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_lic),zx_tmp_2d)
3498	CALL histwrite(nid_ins,"dtsvdfg",itau_w,zx_tmp_2d,iim*jjmp1,
3499	. ndex2d)
3500	c
3501	CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_sic),zx_tmp_2d)
3502	CALL histwrite(nid_ins,"dtsvdfi",itau_w,zx_tmp_2d,iim*jjmp1,
3503	. ndex2d)
3504
3505	DO nsrf = 1, nbsrf
3506	C§§§
3507	zx_tmp_fi2d(1 : klon) = pctsrf( 1 : klon, nsrf)
3508	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3509	CALL histwrite(nid_ins,"pourc_"//clnsurf(nsrf),itau_w,
3510	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3511	C
3512	zx_tmp_fi2d(1 : klon) = fluxt( 1 : klon, 1, nsrf)
3513	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3514	CALL histwrite(nid_ins,"sens_"//clnsurf(nsrf),itau_w,
3515	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3516	C
3517	zx_tmp_fi2d(1 : klon) = fluxlat( 1 : klon, nsrf)
3518	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3519	CALL histwrite(nid_ins,"lat_"//clnsurf(nsrf),itau_w,
3520	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3521	C
3522	zx_tmp_fi2d(1 : klon) = ftsol( 1 : klon, nsrf)
3523	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3524	CALL histwrite(nid_ins,"tsol_"//clnsurf(nsrf),itau_w,
3525	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3526	C
3527	zx_tmp_fi2d(1 : klon) = fluxu( 1 : klon, 1, nsrf)
3528	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3529	CALL histwrite(nid_ins,"taux_"//clnsurf(nsrf),itau_w,
3530	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3531	C
3532	zx_tmp_fi2d(1 : klon) = fluxv( 1 : klon, 1, nsrf)
3533	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3534	CALL histwrite(nid_ins,"tauy_"//clnsurf(nsrf),itau_w,
3535	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3536	C
3537	zx_tmp_fi2d(1 : klon) = frugs( 1 : klon, nsrf)
3538	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3539	CALL histwrite(nid_ins,"rugs_"//clnsurf(nsrf),itau_w,
3540	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3541	C
3542	zx_tmp_fi2d(1 : klon) = falbe( 1 : klon, nsrf)
3543	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d , zx_tmp_2d)
3544	CALL histwrite(nid_ins,"albe_"//clnsurf(nsrf),itau_w,
3545	$ zx_tmp_2d,iim*jjmp1,ndex2d)
3546	C
3547	END DO
3548	CALL gr_fi_ecrit(1, klon,iim,jjmp1, albsol,zx_tmp_2d)
3549	CALL histwrite(nid_ins,"albs",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3550	CALL gr_fi_ecrit(1, klon,iim,jjmp1, albsollw,zx_tmp_2d)
3551	CALL histwrite(nid_ins,"albslw",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3552	c
3553	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxsnow,zx_tmp_2d)
3554	CALL histwrite(nid_ins,"snow_mass",itau_w,zx_tmp_2d,iim*jjmp1,
3555	. ndex2d)
3556	c
3557	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxrugs,zx_tmp_2d)
3558	CALL histwrite(nid_ins,"rugs",itau_w,zx_tmp_2d,iim*jjmp1,ndex2d)
3559	c
3560	c Champs 3D:
3561	c
3562	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
3563	CALL histwrite(nid_ins,"temp",itau_w,zx_tmp_3d,
3564	. iimjjmp1klev,ndex3d)
3565	c
3566	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
3567	CALL histwrite(nid_ins,"vitu",itau_w,zx_tmp_3d,
3568	. iimjjmp1klev,ndex3d)
3569	c
3570	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
3571	CALL histwrite(nid_ins,"vitv",itau_w,zx_tmp_3d,
3572	. iimjjmp1klev,ndex3d)
3573	c
3574	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
3575	CALL histwrite(nid_ins,"geop",itau_w,zx_tmp_3d,
3576	. iimjjmp1klev,ndex3d)
3577	c
3578	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
3579	CALL histwrite(nid_ins,"pres",itau_w,zx_tmp_3d,
3580	. iimjjmp1klev,ndex3d)
3581	c
3582	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_vdf, zx_tmp_3d)
3583	CALL histwrite(nid_ins,"dtvdf",itau_w,zx_tmp_3d,
3584	. iimjjmp1klev,ndex3d)
3585	c
3586	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_vdf, zx_tmp_3d)
3587	CALL histwrite(nid_ins,"dqvdf",itau_w,zx_tmp_3d,
3588	. iimjjmp1klev,ndex3d)
3589
3590	c
3591	if (ok_sync) then
3592	call histsync(nid_ins)
3593	endif
3594	ENDIF
3595	c
3596	c
3597	c Ecrire la bande regionale (binaire grads)
3598	IF (ok_region .AND. mod(itap,ecrit_reg).eq.0) THEN
3599	CALL ecriregs(84,zxtsol)
3600	CALL ecriregs(84,paprs(1,1))
3601	CALL ecriregs(84,topsw)
3602	CALL ecriregs(84,toplw)
3603	CALL ecriregs(84,solsw)
3604	CALL ecriregs(84,sollw)
3605	CALL ecriregs(84,rain_fall)
3606	CALL ecriregs(84,snow_fall)
3607	CALL ecriregs(84,evap)
3608	CALL ecriregs(84,sens)
3609	CALL ecriregs(84,bils)
3610	CALL ecriregs(84,pctsrf(1,is_sic))
3611	CALL ecriregs(84,zxfluxu(1,1))
3612	CALL ecriregs(84,zxfluxv(1,1))
3613	CALL ecriregs(84,ue)
3614	CALL ecriregs(84,ve)
3615	CALL ecriregs(84,uq)
3616	CALL ecriregs(84,vq)
3617	c
3618	CALL ecrirega(84,u_seri)
3619	CALL ecrirega(84,v_seri)
3620	CALL ecrirega(84,omega)
3621	CALL ecrirega(84,t_seri)
3622	CALL ecrirega(84,zphi)
3623	CALL ecrirega(84,q_seri)
3624	CALL ecrirega(84,cldfra)
3625	CALL ecrirega(84,cldliq)
3626	CALL ecrirega(84,pplay)
3627
3628
3629	cc CALL ecrirega(84,d_t_dyn)
3630	cc CALL ecrirega(84,d_q_dyn)
3631	cc CALL ecrirega(84,heat)
3632	cc CALL ecrirega(84,cool)
3633	cc CALL ecrirega(84,d_t_con)
3634	cc CALL ecrirega(84,d_q_con)
3635	cc CALL ecrirega(84,d_t_lsc)
3636	cc CALL ecrirega(84,d_q_lsc)
3637	ENDIF
3638	c
3639	c Convertir les incrementations en tendances
3640	c
3641	DO k = 1, klev
3642	DO i = 1, klon
3643	d_u(i,k) = ( u_seri(i,k) - u(i,k) ) / dtime
3644	d_v(i,k) = ( v_seri(i,k) - v(i,k) ) / dtime
3645	d_t(i,k) = ( t_seri(i,k)-t(i,k) ) / dtime
3646	d_qx(i,k,ivap) = ( q_seri(i,k) - qx(i,k,ivap) ) / dtime
3647	d_qx(i,k,iliq) = ( ql_seri(i,k) - qx(i,k,iliq) ) / dtime
3648	ENDDO
3649	ENDDO
3650	c
3651	IF (nqmax.GE.3) THEN
3652	DO iq = 3, nqmax
3653	DO k = 1, klev
3654	DO i = 1, klon
3655	d_qx(i,k,iq) = ( tr_seri(i,k,iq-2) - qx(i,k,iq) ) / dtime
3656	ENDDO
3657	ENDDO
3658	ENDDO
3659	ENDIF
3660	c
3661	c Sauvegarder les valeurs de t et q a la fin de la physique:
3662	c
3663	DO k = 1, klev
3664	DO i = 1, klon
3665	t_ancien(i,k) = t_seri(i,k)
3666	q_ancien(i,k) = q_seri(i,k)
3667	ENDDO
3668	ENDDO
3669	c
3670	c====================================================================
3671	c Si c'est la fin, il faut conserver l'etat de redemarrage
3672	c====================================================================
3673	c
3674	IF (lafin) THEN
3675	itau_phy = itau_phy + itap
3676	ccc IF (ok_oasis) CALL quitcpl
3677	CALL phyredem ("restartphy.nc",dtime,radpas,co2_ppm,solaire,
3678	. rlat, rlon, pctsrf, ftsol, ftsoil, deltat, fqsol, fsnow,
3679	. falbe, fevap, rain_fall, snow_fall,
3680	. solsw, sollwdown,dlw,
3681	. radsol,frugs,agesno,
3682	. zmea,zstd,zsig,zgam,zthe,zpic,zval,rugoro,
3683	. t_ancien, q_ancien, rnebcon, ratqs, clwcon)
3684	ENDIF
3685
3686	RETURN
3687	END
3688	FUNCTION qcheck(klon,klev,paprs,q,ql,aire)
3689	IMPLICIT none
3690	c
3691	c Calculer et imprimer l'eau totale. A utiliser pour verifier
3692	c la conservation de l'eau
3693	c
3694	#include "YOMCST.h"
3695	INTEGER klon,klev
3696	REAL paprs(klon,klev+1), q(klon,klev), ql(klon,klev)
3697	REAL aire(klon)
3698	REAL qtotal, zx, qcheck
3699	INTEGER i, k
3700	c
3701	zx = 0.0
3702	DO i = 1, klon
3703	zx = zx + aire(i)
3704	ENDDO
3705	qtotal = 0.0
3706	DO k = 1, klev
3707	DO i = 1, klon
3708	qtotal = qtotal + (q(i,k)+ql(i,k)) * aire(i)
3709	. *(paprs(i,k)-paprs(i,k+1))/RG
3710	ENDDO
3711	ENDDO
3712	c
3713	qcheck = qtotal/zx
3714	c
3715	RETURN
3716	END
3717	SUBROUTINE gr_fi_ecrit(nfield,nlon,iim,jjmp1,fi,ecrit)
3718	IMPLICIT none
3719	c
3720	c Tranformer une variable de la grille physique a
3721	c la grille d'ecriture
3722	c
3723	INTEGER nfield,nlon,iim,jjmp1, jjm
3724	REAL fi(nlon,nfield), ecrit(iim*jjmp1,nfield)
3725	c
3726	INTEGER i, n, ig
3727	c
3728	jjm = jjmp1 - 1
3729	DO n = 1, nfield
3730	DO i=1,iim
3731	ecrit(i,n) = fi(1,n)
3732	ecrit(i+jjm*iim,n) = fi(nlon,n)
3733	ENDDO
3734	DO ig = 1, nlon - 2
3735	ecrit(iim+ig,n) = fi(1+ig,n)
3736	ENDDO
3737	ENDDO
3738	RETURN
3739	END
3740

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

Download in other formats:

Original Format