Context Navigation

physiq.F @ 302

Last change on this file since 302 was 302, checked in by lmdz, 23 years ago
Correction bug dans KE vectorise (FH, JYG) LF
Property svn:eol-style set to `native` Property svn:executable set to ``* Property svn:keywords set to `Author Date Id Revision`
File size: 99.5 KB

Line
1	c $Header$
2	c
3	SUBROUTINE physiq (nlon,nlev,nqmax ,
4	. debut,lafin,rjourvrai,rjour_ecri,gmtime,pdtphys,
5	. paprs,pplay,pphi,pphis,paire,presnivs,clesphy0,
6	. u,v,t,qx,
7	. omega,
8	. d_u, d_v, d_t, d_qx, d_ps)
9	USE ioipsl
10	IMPLICIT none
11	c======================================================================
12	c
13	c Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
14	c
15	c Objet: Moniteur general de la physique du modele
16	cAA Modifications quant aux traceurs :
17	cAA - uniformisation des parametrisations ds phytrac
18	cAA - stockage des moyennes des champs necessaires
19	cAA en mode traceur off-line
20	c======================================================================
21	c modif ( P. Le Van , 12/10/98 )
22	c
23	c Arguments:
24	c
25	c nlon----input-I-nombre de points horizontaux
26	c nlev----input-I-nombre de couches verticales
27	c nqmax---input-I-nombre de traceurs (y compris vapeur d'eau) = 1
28	c debut---input-L-variable logique indiquant le premier passage
29	c lafin---input-L-variable logique indiquant le dernier passage
30	c rjour---input-R-numero du jour de l'experience
31	c gmtime--input-R-temps universel dans la journee (0 a 86400 s)
32	c pdtphys-input-R-pas d'integration pour la physique (seconde)
33	c paprs---input-R-pression pour chaque inter-couche (en Pa)
34	c pplay---input-R-pression pour le mileu de chaque couche (en Pa)
35	c pphi----input-R-geopotentiel de chaque couche (g z) (reference sol)
36	c pphis---input-R-geopotentiel du sol
37	c paire---input-R-aire de chaque maille
38	c presnivs-input_R_pressions approximat. des milieux couches ( en PA)
39	c u-------input-R-vitesse dans la direction X (de O a E) en m/s
40	c v-------input-R-vitesse Y (de S a N) en m/s
41	c t-------input-R-temperature (K)
42	c qx------input-R-humidite specifique (kg/kg) et d'autres traceurs
43	c d_t_dyn-input-R-tendance dynamique pour "t" (K/s)
44	c d_q_dyn-input-R-tendance dynamique pour "q" (kg/kg/s)
45	c omega---input-R-vitesse verticale en Pa/s
46	c
47	c d_u-----output-R-tendance physique de "u" (m/s/s)
48	c d_v-----output-R-tendance physique de "v" (m/s/s)
49	c d_t-----output-R-tendance physique de "t" (K/s)
50	c d_qx----output-R-tendance physique de "qx" (kg/kg/s)
51	c d_ps----output-R-tendance physique de la pression au sol
52	c======================================================================
53	#include "dimensions.h"
54	integer jjmp1
55	parameter (jjmp1=jjm+1-1/jjm)
56	#include "dimphy.h"
57	#include "regdim.h"
58	#include "indicesol.h"
59	#include "dimsoil.h"
60	#include "clesphys.h"
61	#include "control.h"
62	#include "temps.h"
63	c======================================================================
64	LOGICAL check ! Verifier la conservation du modele en eau
65	PARAMETER (check=.FALSE.)
66	LOGICAL ok_stratus ! Ajouter artificiellement les stratus
67	PARAMETER (ok_stratus=.FALSE.)
68	c======================================================================
69	c Parametres lies au coupleur OASIS:
70	#include "oasis.h"
71	INTEGER npas, nexca, itimestep
72	logical rnpb
73	parameter(rnpb=.true.)
74	PARAMETER (npas=1440)
75	PARAMETER (nexca=48)
76	PARAMETER (itimestep=1800)
77	EXTERNAL fromcpl, intocpl, inicma
78	REAL cpl_sst(iim,jjmp1), cpl_sic(iim,jjmp1)
79	REAL cpl_alb_sst(iim,jjmp1), cpl_alb_sic(iim,jjmp1)
80	c======================================================================
81	c ok_ocean indique l'utilisation du modele oceanique "slab ocean",
82	c il faut bien sur s'assurer que le bilan energetique de reference
83	c a la surface de l'ocean est bien present dans le fichier des
84	c conditions aux limites, ainsi que l'indicateur du sol ne contient
85	c pas de glace oceanique (pas de valeur 3).
86	c
87	LOGICAL ok_ocean
88	PARAMETER (ok_ocean=.FALSE.)
89	REAL cyang ! capacite thermique de l'ocean superficiel
90	PARAMETER (cyang=30.0 * 4.228e+06)
91	REAL cbing ! capacite thermique de la glace oceanique
92	PARAMETER (cbing=1.0 * 4.228e+06)
93	REAL cthermiq
94	c======================================================================
95	c Clef controlant l'activation du cycle diurne:
96	ccc LOGICAL cycle_diurne
97	ccc PARAMETER (cycle_diurne=.FALSE.)
98	c======================================================================
99	c Modele thermique du sol, a activer pour le cycle diurne:
100	ccc LOGICAL soil_model
101	ccc PARAMETER (soil_model=.FALSE.)
102	REAL soilcap(klon,nbsrf), soilflux(klon,nbsrf)
103	SAVE soilcap, soilflux
104	c======================================================================
105	c Dans les versions precedentes, l'eau liquide nuageuse utilisee dans
106	c le calcul du rayonnement est celle apres la precipitation des nuages.
107	c Si cette cle new_oliq est activee, ce sera une valeur moyenne entre
108	c la condensation et la precipitation. Cette cle augmente les impacts
109	c radiatifs des nuages.
110	ccc LOGICAL new_oliq
111	ccc PARAMETER (new_oliq=.FALSE.)
112	c======================================================================
113	c Clefs controlant deux parametrisations de l'orographie:
114	cc LOGICAL ok_orodr
115	ccc PARAMETER (ok_orodr=.FALSE.)
116	ccc LOGICAL ok_orolf
117	ccc PARAMETER (ok_orolf=.FALSE.)
118	c======================================================================
119	LOGICAL ok_journe ! sortir le fichier journalier
120	PARAMETER (ok_journe=.FALSE.)
121	c
122	LOGICAL ok_mensuel ! sortir le fichier mensuel
123	PARAMETER (ok_mensuel=.TRUE.)
124	c
125	LOGICAL ok_instan ! sortir le fichier instantane
126	PARAMETER (ok_instan=.FALSE.)
127	c
128	LOGICAL ok_region ! sortir le fichier regional
129	PARAMETER (ok_region=.FALSE.)
130	c======================================================================
131	c
132	INTEGER ivap ! indice de traceurs pour vapeur d'eau
133	PARAMETER (ivap=1)
134	INTEGER iliq ! indice de traceurs pour eau liquide
135	PARAMETER (iliq=2)
136	c
137	INTEGER nvm ! nombre de vegetations
138	PARAMETER (nvm=8)
139	REAL veget(klon,nvm) ! couverture vegetale
140	SAVE veget
141	c
142	c Variables argument:
143	c
144	INTEGER nlon
145	INTEGER nlev
146	INTEGER nqmax
147	REAL rjourvrai, rjour_ecri
148	REAL gmtime
149	REAL pdtphys
150	LOGICAL debut, lafin
151	REAL paprs(klon,klev+1)
152	REAL pplay(klon,klev)
153	REAL pphi(klon,klev)
154	REAL pphis(klon)
155	REAL paire(klon)
156	REAL presnivs(klev)
157	REAL znivsig(klev)
158
159	REAL u(klon,klev)
160	REAL v(klon,klev)
161	REAL t(klon,klev)
162	REAL qx(klon,klev,nqmax)
163
164	REAL t_ancien(klon,klev), q_ancien(klon,klev)
165	SAVE t_ancien, q_ancien
166	LOGICAL ancien_ok
167	SAVE ancien_ok
168
169	REAL d_u_dyn(klon,klev)
170	REAL d_v_dyn(klon,klev)
171	REAL d_t_dyn(klon,klev)
172	REAL d_q_dyn(klon,klev)
173
174	REAL omega(klon,klev)
175
176	REAL d_u(klon,klev)
177	REAL d_v(klon,klev)
178	REAL d_t(klon,klev)
179	REAL d_qx(klon,klev,nqmax)
180	REAL d_ps(klon)
181
182	INTEGER longcles
183	PARAMETER ( longcles = 20 )
184	REAL clesphy0( longcles )
185	c
186	c Variables quasi-arguments
187	c
188	REAL xjour
189	SAVE xjour
190	c
191	c
192	c Variables propres a la physique
193	c
194	REAL dtime
195	SAVE dtime ! pas temporel de la physique
196	c
197	INTEGER radpas
198	SAVE radpas ! frequence d'appel rayonnement
199	c
200	REAL radsol(klon)
201	SAVE radsol ! bilan radiatif au sol
202	c
203	REAL rlat(klon)
204	SAVE rlat ! latitude pour chaque point
205	c
206	REAL rlon(klon)
207	SAVE rlon ! longitude pour chaque point
208	c
209	cc INTEGER iflag_con
210	cc SAVE iflag_con ! indicateur de la convection
211	c
212	INTEGER itap
213	SAVE itap ! compteur pour la physique
214	c
215	REAL co2_ppm
216	SAVE co2_ppm ! concentration du CO2
217	c
218	REAL solaire
219	SAVE solaire ! constante solaire
220	c
221	REAL ftsol(klon,nbsrf)
222	SAVE ftsol ! temperature du sol
223	c
224	REAL ftsoil(klon,nsoilmx,nbsrf)
225	SAVE ftsoil ! temperature dans le sol
226	c
227	REAL deltat(klon)
228	SAVE deltat ! ecart avec la SST de reference
229	c
230	REAL fqsol(klon,nbsrf)
231	SAVE fqsol ! humidite du sol
232	c
233	REAL fsnow(klon,nbsrf)
234	SAVE fsnow ! epaisseur neigeuse
235	c
236	REAL rugmer(klon)
237	SAVE rugmer ! longeur de rugosite sur mer (m)
238	c
239	c Parametres de l'Orographie a l'Echelle Sous-Maille (OESM):
240	c
241	REAL zmea(klon)
242	SAVE zmea ! orographie moyenne
243	c
244	REAL zstd(klon)
245	SAVE zstd ! deviation standard de l'OESM
246	c
247	REAL zsig(klon)
248	SAVE zsig ! pente de l'OESM
249	c
250	REAL zgam(klon)
251	save zgam ! anisotropie de l'OESM
252	c
253	REAL zthe(klon)
254	SAVE zthe ! orientation de l'OESM
255	c
256	REAL zpic(klon)
257	SAVE zpic ! Maximum de l'OESM
258	c
259	REAL zval(klon)
260	SAVE zval ! Minimum de l'OESM
261	c
262	REAL rugoro(klon)
263	SAVE rugoro ! longueur de rugosite de l'OESM
264	c
265	REAL zulow(klon),zvlow(klon),zustr(klon), zvstr(klon)
266	c
267	REAL zuthe(klon),zvthe(klon)
268	SAVE zuthe
269	SAVE zvthe
270	INTEGER igwd,igwdim,idx(klon),itest(klon)
271	c
272	REAL agesno(klon)
273	SAVE agesno ! age de la neige
274	c
275	REAL alb_neig(klon)
276	SAVE alb_neig ! albedo de la neige
277	cKE43
278	c Variables liees a la convection de K. Emanuel (sb):
279	c
280	REAL ema_workcbmf(klon) ! cloud base mass flux
281	SAVE ema_workcbmf
282
283	REAL ema_cbmf(klon) ! cloud base mass flux
284	SAVE ema_cbmf
285
286	REAL ema_pcb(klon) ! cloud base pressure
287	SAVE ema_pcb
288
289	REAL ema_pct(klon) ! cloud top pressure
290	SAVE ema_pct
291
292	REAL bas, top ! cloud base and top levels
293	SAVE bas
294	SAVE top
295
296	REAL Ma(klon,klev) ! undilute upward mass flux
297	SAVE Ma
298	REAL ema_work1(klon, klev), ema_work2(klon, klev)
299	SAVE ema_work1, ema_work2
300	REAL wdn(klon), tdn(klon), qdn(klon)
301	c Variables locales pour la couche limite (al1):
302	c
303	cAl1 REAL pblh(klon) ! Hauteur de couche limite
304	cAl1 SAVE pblh
305	c34EK
306	c
307	c Variables locales:
308	c
309	REAL cdragh(klon) ! drag coefficient pour T and Q
310	REAL cdragm(klon) ! drag coefficient pour vent
311	cAA
312	cAA Pour phytrac
313	cAA
314	REAL ycoefh(klon,klev) ! coef d'echange pour phytrac
315	REAL yu1(klon) ! vents dans la premiere couche U
316	REAL yv1(klon) ! vents dans la premiere couche V
317	LOGICAL offline ! Controle du stockage ds "physique"
318	PARAMETER (offline=.false.)
319	INTEGER physid
320	REAL pfrac_impa(klon,klev)! Produits des coefs lessivage impaction
321	save pfrac_impa
322	REAL pfrac_nucl(klon,klev)! Produits des coefs lessivage nucleation
323	save pfrac_nucl
324	REAL pfrac_1nucl(klon,klev)! Produits des coefs lessi nucl (alpha = 1)
325	save pfrac_1nucl
326	REAL frac_impa(klon,klev) ! fractions d'aerosols lessivees (impaction)
327	REAL frac_nucl(klon,klev) ! idem (nucleation)
328	cAA
329	REAL rain_fall(klon) ! pluie
330	REAL snow_fall(klon) ! neige
331	REAL evap(klon), devap(klon) ! evaporation et sa derivee
332	REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
333	REAL bils(klon) ! bilan de chaleur au sol
334	REAL fder(klon) ! Derive de flux (sensible et latente)
335	REAL ruis(klon) ! ruissellement
336	REAL ve(klon) ! integr. verticale du transport meri. de l'energie
337	REAL vq(klon) ! integr. verticale du transport meri. de l'eau
338	REAL ue(klon) ! integr. verticale du transport zonal de l'energie
339	REAL uq(klon) ! integr. verticale du transport zonal de l'eau
340	c
341	REAL frugs(klon,nbsrf) ! longueur de rugosite
342	REAL zxrugs(klon) ! longueur de rugosite
343	c
344	c Conditions aux limites
345	c
346	INTEGER julien
347	INTEGER idayvrai
348	SAVE idayvrai
349	c
350	INTEGER lmt_pas
351	SAVE lmt_pas ! frequence de mise a jour
352	REAL pctsrf(klon,nbsrf)
353	SAVE pctsrf ! sous-fraction du sol
354	REAL lmt_sst(klon)
355	SAVE lmt_sst ! temperature de la surface ocean
356	REAL lmt_bils(klon)
357	SAVE lmt_bils ! bilan de chaleur au sol
358	REAL lmt_alb(klon)
359	SAVE lmt_alb ! temperature de la surface ocean
360	REAL lmt_rug(klon)
361	SAVE lmt_rug ! longueur de rugosite
362	REAL alb_eau(klon)
363	SAVE alb_eau ! albedo sur l'ocean
364	REAL albsol(klon)
365	SAVE albsol ! albedo du sol total
366	REAL wo(klon,klev)
367	SAVE wo ! ozone
368	c======================================================================
369	c
370	c Declaration des procedures appelees
371	c
372	EXTERNAL angle ! calculer angle zenithal du soleil
373	EXTERNAL alboc ! calculer l'albedo sur ocean
374	EXTERNAL albsno ! calculer albedo sur neige
375	EXTERNAL ajsec ! ajustement sec
376	EXTERNAL clmain ! couche limite
377	EXTERNAL condsurf ! lire les conditions aux limites
378	EXTERNAL conlmd ! convection (schema LMD)
379	cKE43
380	EXTERNAL conema ! convect4.3
381	EXTERNAL fisrtilp ! schema de condensation a grande echelle (pluie)
382	cAA
383	EXTERNAL fisrtilp_tr ! schema de condensation a grande echelle (pluie)
384	c ! stockage des coefficients necessaires au
385	c ! lessivage OFF-LINE et ON-LINE
386	cAA
387	EXTERNAL hgardfou ! verifier les temperatures
388	EXTERNAL hydrol ! hydrologie du sol
389	EXTERNAL nuage ! calculer les proprietes radiatives
390	EXTERNAL o3cm ! initialiser l'ozone
391	EXTERNAL orbite ! calculer l'orbite terrestre
392	EXTERNAL ozonecm ! prescrire l'ozone
393	EXTERNAL phyetat0 ! lire l'etat initial de la physique
394	EXTERNAL phyredem ! ecrire l'etat de redemarrage de la physique
395	EXTERNAL radlwsw ! rayonnements solaire et infrarouge
396	EXTERNAL suphec ! initialiser certaines constantes
397	EXTERNAL transp ! transport total de l'eau et de l'energie
398	EXTERNAL ecribina ! ecrire le fichier binaire global
399	EXTERNAL ecribins ! ecrire le fichier binaire global
400	EXTERNAL ecrirega ! ecrire le fichier binaire regional
401	EXTERNAL ecriregs ! ecrire le fichier binaire regional
402	c
403	c Variables locales
404	c
405	REAL dialiq(klon,klev) ! eau liquide nuageuse
406	REAL diafra(klon,klev) ! fraction nuageuse
407	REAL cldliq(klon,klev) ! eau liquide nuageuse
408	REAL cldfra(klon,klev) ! fraction nuageuse
409	REAL cldtau(klon,klev) ! epaisseur optique
410	REAL cldemi(klon,klev) ! emissivite infrarouge
411	c
412	REAL fluxq(klon,klev) ! flux turbulent d'humidite
413	REAL fluxt(klon,klev) ! flux turbulent de chaleur
414	REAL fluxu(klon,klev) ! flux turbulent de vitesse u
415	REAL fluxv(klon,klev) ! flux turbulent de vitesse v
416	c
417	REAL heat(klon,klev) ! chauffage solaire
418	REAL heat0(klon,klev) ! chauffage solaire ciel clair
419	REAL cool(klon,klev) ! refroidissement infrarouge
420	REAL cool0(klon,klev) ! refroidissement infrarouge ciel clair
421	REAL topsw(klon), toplw(klon), solsw(klon), sollw(klon)
422	REAL topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
423	REAL albpla(klon)
424	c Le rayonnement n'est pas calcule tous les pas, il faut donc
425	c sauvegarder les sorties du rayonnement
426	SAVE heat,cool,albpla,topsw,toplw,solsw,sollw
427	SAVE topsw0,toplw0,solsw0,sollw0, heat0, cool0
428	INTEGER itaprad
429	SAVE itaprad
430	c
431	REAL conv_q(klon,klev) ! convergence de l'humidite (kg/kg/s)
432	REAL conv_t(klon,klev) ! convergence de la temperature(K/s)
433	c
434	REAL cldl(klon),cldm(klon),cldh(klon) !nuages bas, moyen et haut
435	REAL cldt(klon),cldq(klon) !nuage total, eau liquide integree
436	c
437	REAL zx_alb_lic, zx_alb_oce, zx_alb_ter, zx_alb_sic
438	REAL zxtsol(klon), zxqsol(klon), zxsnow(klon)
439	c
440	REAL dist, rmu0(klon), fract(klon)
441	REAL zdtime, zlongi
442	c
443	CHARACTER*2 str2
444	CHARACTER*2 iqn
445	c
446	REAL qcheck
447	REAL z_avant(klon), z_apres(klon), z_factor(klon)
448	LOGICAL zx_ajustq
449	c
450	REAL za, zb
451	REAL zx_t, zx_qs, zdelta, zcor, zfra, zlvdcp, zlsdcp
452	INTEGER i, k, iq, ig, j, nsrf, ll
453	REAL t_coup
454	PARAMETER (t_coup=234.0)
455	c
456	REAL zphi(klon,klev)
457	REAL zx_tmp_x(iim), zx_tmp_yjjmp1
458	REAL zx_relief(iim,jjmp1)
459	REAL zx_aire(iim,jjmp1)
460	cKE43
461	c Variables locales pour la convection de K. Emanuel (sb):
462	c
463	REAL upwd(klon,klev) ! saturated updraft mass flux
464	REAL dnwd(klon,klev) ! saturated downdraft mass flux
465	REAL dnwd0(klon,klev) ! unsaturated downdraft mass flux
466	REAL tvp(klon,klev) ! virtual temp of lifted parcel
467	REAL cape(klon) ! CAPE
468	SAVE cape
469	REAL pbase(klon) ! cloud base pressure
470	SAVE pbase
471	REAL bbase(klon) ! cloud base buoyancy
472	SAVE bbase
473	REAL rflag(klon) ! flag fonctionnement de convect
474	INTEGER iflagctrl(klon) ! flag fonctionnement de convect
475	c -- convect43:
476	INTEGER ntra ! nb traceurs pour convect4.3
477	REAL pori_con(klon) ! pressure at the origin level of lifted parcel
478	REAL plcl_con(klon),dtma_con(klon),dtlcl_con(klon)
479	REAL dtvpdt1(klon,klev), dtvpdq1(klon,klev)
480	REAL dplcldt(klon), dplcldr(klon)
481	c? . condm_con(klon,klev),conda_con(klon,klev),
482	c? . mr_con(klon,klev),ep_con(klon,klev)
483	c? . ,sadiab(klon,klev),wadiab(klon,klev)
484	c --
485	c34EK
486	c
487	c Variables du changement
488	c
489	c con: convection
490	c lsc: condensation a grande echelle (Large-Scale-Condensation)
491	c ajs: ajustement sec
492	c eva: evaporation de l'eau liquide nuageuse
493	c vdf: couche limite (Vertical DiFfusion)
494	REAL d_t_con(klon,klev),d_q_con(klon,klev)
495	REAL d_u_con(klon,klev),d_v_con(klon,klev)
496	REAL d_t_lsc(klon,klev),d_q_lsc(klon,klev),d_ql_lsc(klon,klev)
497	REAL d_t_ajs(klon,klev), d_q_ajs(klon,klev)
498	REAL d_t_eva(klon,klev),d_q_eva(klon,klev)
499	REAL rneb(klon,klev)
500	c
501	REAL pmfu(klon,klev), pmfd(klon,klev)
502	REAL pen_u(klon,klev), pen_d(klon,klev)
503	REAL pde_u(klon,klev), pde_d(klon,klev)
504	INTEGER kcbot(klon), kctop(klon), kdtop(klon)
505	REAL pmflxr(klon,klev+1), pmflxs(klon,klev+1)
506	REAL prfl(klon,klev+1), psfl(klon,klev+1)
507	c
508	INTEGER ibas_con(klon), itop_con(klon)
509	REAL rain_con(klon), rain_lsc(klon)
510	REAL snow_con(klon), snow_lsc(klon)
511	REAL d_ts(klon,nbsrf)
512	c
513	REAL d_u_vdf(klon,klev), d_v_vdf(klon,klev)
514	REAL d_t_vdf(klon,klev), d_q_vdf(klon,klev)
515	c
516	REAL d_u_oro(klon,klev), d_v_oro(klon,klev)
517	REAL d_t_oro(klon,klev)
518	REAL d_u_lif(klon,klev), d_v_lif(klon,klev)
519	REAL d_t_lif(klon,klev)
520
521	REAL ratqs(klon,klev)
522	real zpt_conv(klon,klev)
523
524	c
525	c Variables liees a l'ecriture de la bande histoire physique
526	c
527	INTEGER ecrit_mth
528	SAVE ecrit_mth ! frequence d'ecriture (fichier mensuel)
529	c
530	INTEGER ecrit_day
531	SAVE ecrit_day ! frequence d'ecriture (fichier journalier)
532	c
533	INTEGER ecrit_ins
534	SAVE ecrit_ins ! frequence d'ecriture (fichier instantane)
535	c
536	INTEGER ecrit_reg
537	SAVE ecrit_reg ! frequence d'ecriture
538	c
539	REAL oas_sols(klon), z_sols(iim,jjmp1)
540	SAVE oas_sols
541	REAL oas_nsol(klon), z_nsol(iim,jjmp1)
542	SAVE oas_nsol
543	REAL oas_rain(klon), z_rain(iim,jjmp1)
544	SAVE oas_rain
545	REAL oas_snow(klon), z_snow(iim,jjmp1)
546	SAVE oas_snow
547	REAL oas_evap(klon), z_evap(iim,jjmp1)
548	SAVE oas_evap
549	REAL oas_ruis(klon), z_ruis(iim,jjmp1)
550	SAVE oas_ruis
551	REAL oas_tsol(klon), z_tsol(iim,jjmp1)
552	SAVE oas_tsol
553	REAL oas_fder(klon), z_fder(iim,jjmp1)
554	SAVE oas_fder
555	REAL oas_albe(klon), z_albe(iim,jjmp1)
556	SAVE oas_albe
557	REAL oas_taux(klon), z_taux(iim,jjmp1)
558	SAVE oas_taux
559	REAL oas_tauy(klon), z_tauy(iim,jjmp1)
560	SAVE oas_tauy
561	REAL oas_ruisoce(klon), z_ruisoce(iim,jjmp1)
562	SAVE oas_ruisoce
563	REAL oas_ruisriv(klon), z_ruisriv(iim,jjmp1)
564	SAVE oas_ruisriv
565	c
566	c
567	c Variables locales pour effectuer les appels en serie
568	c
569	REAL t_seri(klon,klev), q_seri(klon,klev)
570	REAL ql_seri(klon,klev)
571	REAL u_seri(klon,klev), v_seri(klon,klev)
572	c
573	REAL tr_seri(klon,klev,nbtr)
574	REAL d_tr(klon,klev,nbtr)
575	REAL source_tr(klon,nbtr)
576
577	REAL zx_rh(klon,klev)
578	REAL dtimeday,dtimecri,dtimexp9,fecri_pas,fecri86400,fecritday
579
580	INTEGER length
581	PARAMETER ( length = 100 )
582	REAL tabcntr0( length )
583	c
584	INTEGER ndex2d(iimjjmp1),ndex3d(iimjjmp1*klev)
585	REAL zx_tmp_fi2d(klon)
586	REAL zx_tmp_2d(iim,jjmp1), zx_tmp_3d(iim,jjmp1,klev)
587	REAL zx_lon(iim,jjmp1), zx_lat(iim,jjmp1)
588	c
589	INTEGER nid_day, nid_mth, nid_ins
590	SAVE nid_day, nid_mth, nid_ins
591	c
592	INTEGER nhori, nvert
593	REAL zsto, zout, zjulian
594
595	character*20 modname
596	character*80 abort_message
597	logical ok_sync
598
599	c
600	c Declaration des constantes et des fonctions thermodynamiques
601	c
602	#include "YOMCST.h"
603	#include "YOETHF.h"
604	#include "FCTTRE.h"
605	c======================================================================
606	modname = 'physiq'
607	ok_sync=.TRUE.
608	IF (nqmax .LT. 2) THEN
609	PRINT*, 'eaux vapeur et liquide sont indispensables'
610	CALL ABORT
611	ENDIF
612	IF (debut) THEN
613	CALL suphec ! initialiser constantes et parametres phys.
614	ENDIF
615	c======================================================================
616	xjour = rjourvrai
617	c
618	c Si c'est le debut, il faut initialiser plusieurs choses
619	c ********
620	c
621	IF (debut) THEN
622	c
623
624	IF (ok_oasis) THEN
625	PRINT*, "Attentions! les parametres suivants sont fixes:"
626	PRINT ,'**********************************************'
627	PRINT*, "npas, nexca, itimestep=", npas, nexca, itimestep
628	PRINT*, "Changer-les manuellement s il le faut"
629	PRINT ,'**********************************************'
630	CALL inicma( npas, nexca, itimestep)
631	ENDIF
632	c
633	IF (ok_ocean) THEN
634	PRINT, '***********************'
635	PRINT*, 'SLAB OCEAN est active, prenez precautions !'
636	PRINT, '***********************'
637	ENDIF
638	c
639	DO k = 2, nvm ! pas de vegetation
640	DO i = 1, klon
641	veget(i,k) = 0.0
642	ENDDO
643	ENDDO
644	DO i = 1, klon
645	veget(i,1) = 1.0 ! il n'y a que du desert
646	ENDDO
647	PRINT*, 'Pas de vegetation; desert partout'
648	c
649	c Initialiser les compteurs:
650	c
651
652	itap = 0
653	itaprad = 0
654	c
655	CALL phyetat0 ("startphy.nc",dtime,co2_ppm,solaire,
656	. rlat,rlon,ftsol,ftsoil,deltat,fqsol,fsnow,
657	. radsol,rugmer,agesno,clesphy0,
658	. zmea,zstd,zsig,zgam,zthe,zpic,zval,rugoro,tabcntr0,
659	. t_ancien, q_ancien, ancien_ok )
660
661	c
662	radpas = NINT( 86400./dtime/nbapp_rad)
663
664	c
665	CALL printflag( tabcntr0,radpas,ok_ocean,ok_oasis ,ok_journe,
666	, ok_instan, ok_region )
667	c
668	IF (ABS(dtime-pdtphys).GT.0.001) THEN
669	PRINT*, 'Pas physique n est pas correcte',dtime,pdtphys
670	abort_message=' See above '
671	call abort_gcm(modname,abort_message,1)
672	ENDIF
673	IF (nlon .NE. klon) THEN
674	PRINT*, 'nlon et klon ne sont pas coherents', nlon, klon
675	abort_message=' See above '
676	call abort_gcm(modname,abort_message,1)
677	ENDIF
678	IF (nlev .NE. klev) THEN
679	PRINT*, 'nlev et klev ne sont pas coherents', nlev, klev
680	abort_message=' See above '
681	call abort_gcm(modname,abort_message,1)
682	ENDIF
683	c
684	IF (dtime*FLOAT(radpas).GT.21600..AND.cycle_diurne) THEN
685	PRINT*, 'Nbre d appels au rayonnement insuffisant'
686	PRINT*, "Au minimum 4 appels par jour si cycle diurne"
687	abort_message=' See above '
688	call abort_gcm(modname,abort_message,1)
689	ENDIF
690	PRINT*, "Clef pour la convection, iflag_con=", iflag_con
691	c
692	cKE43
693	c Initialisation pour la convection de K.E. (sb):
694	IF (iflag_con.GE.3) THEN
695
696	PRINT, "** Convection de Kerry Emanuel 4.3 "
697	PRINT*, "On va utiliser le melange convectif des traceurs qui"
698	PRINT*, "est calcule dans convect4.3"
699	PRINT*, " !!! penser aux logical flags de phytrac"
700
701	DO i = 1, klon
702	ema_cbmf(i) = 0.
703	ema_pcb(i) = 0.
704	ema_pct(i) = 0.
705	ema_workcbmf(i) = 0.
706	ENDDO
707	ENDIF
708	c34EK
709	IF (ok_orodr) THEN
710	DO i=1,klon
711	rugoro(i) = MAX(1.0e-05, zstd(i)*zsig(i)/2.0)
712	ENDDO
713	CALL SUGWD(klon,klev,paprs,pplay)
714	DO i=1,klon
715	zuthe(i)=0.
716	zvthe(i)=0.
717	if(zstd(i).gt.10.)then
718	zuthe(i)=(1.-zgam(i))*cos(zthe(i))
719	zvthe(i)=(1.-zgam(i))*sin(zthe(i))
720	endif
721	ENDDO
722	ENDIF
723	c
724	IF (soil_model) THEN
725	DO nsrf = 1, nbsrf
726	CALL soil(dtime, nsrf, fsnow(1,nsrf),
727	. ftsol(1,nsrf), ftsoil(1,1,nsrf),
728	. soilcap(1,nsrf), soilflux(1,nsrf))
729	ENDDO
730	ENDIF
731	c
732	lmt_pas = NINT(86400./dtime * 1.0) ! tous les jours
733	PRINT*,'La frequence de lecture surface est de ', lmt_pas
734	c
735	ecrit_mth = NINT(86400./dtime *ecritphy) ! tous les ecritphy jours
736	IF (ok_mensuel) THEN
737	PRINT*, 'La frequence de sortie mensuelle est de ', ecrit_mth
738	ENDIF
739	ecrit_day = NINT(86400./dtime *1.0) ! tous les jours
740	IF (ok_journe) THEN
741	PRINT*, 'La frequence de sortie journaliere est de ',ecrit_day
742	ENDIF
743	ccc ecrit_ins = NINT(86400./dtime *0.5) ! 2 fois par jour
744	ecrit_ins = NINT(86400./dtime *0.25) ! tous les jours
745	IF (ok_instan) THEN
746	PRINT*, 'La frequence de sortie instant. est de ', ecrit_ins
747	ENDIF
748	ecrit_reg = NINT(86400./dtime *0.25) ! 4 fois par jour
749	IF (ok_region) THEN
750	PRINT*, 'La frequence de sortie region est de ', ecrit_reg
751	ENDIF
752	c
753	c
754	IF (ok_journe) THEN
755	c
756	CALL ymds2ju(anne_ini, 1, 1, 0.0, zjulian)
757	zjulian = zjulian + day_ini
758	c
759	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
760	DO i = 1, iim
761	zx_lon(i,1) = rlon(i+1)
762	zx_lon(i,jjmp1) = rlon(i+1)
763	ENDDO
764	DO ll=1,klev
765	znivsig(ll)=float(ll)
766	ENDDO
767	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
768	CALL histbeg("histday", iim,zx_lon, jjmp1,zx_lat,
769	. 1,iim,1,jjmp1, 0, zjulian, dtime,
770	. nhori, nid_day)
771	CALL histvert(nid_day, "presnivs", "Vertical levels", "mb",
772	. klev, presnivs, nvert)
773	c call histvert(nid_day, 'sig_s', 'Niveaux sigma','-',
774	c . klev, znivsig, nvert)
775	c
776	zsto = dtime
777	zout = dtime * FLOAT(ecrit_day)
778	c
779	CALL histdef(nid_day, "phis", "Surface geop. height", "-",
780	. iim,jjmp1,nhori, 1,1,1, -99, 32,
781	. "once", zsto,zout)
782	c
783	CALL histdef(nid_day, "aire", "Grid area", "-",
784	. iim,jjmp1,nhori, 1,1,1, -99, 32,
785	. "once", zsto,zout)
786	c
787	c Champs 2D:
788	c
789	CALL histdef(nid_day, "tsol", "Surface Temperature", "K",
790	. iim,jjmp1,nhori, 1,1,1, -99, 32,
791	. "ave(X)", zsto,zout)
792	c
793	CALL histdef(nid_day, "psol", "Surface Pressure", "Pa",
794	. iim,jjmp1,nhori, 1,1,1, -99, 32,
795	. "ave(X)", zsto,zout)
796	c
797	CALL histdef(nid_day, "rain", "Precipitation", "mm/day",
798	. iim,jjmp1,nhori, 1,1,1, -99, 32,
799	. "ave(X)", zsto,zout)
800	c
801	CALL histdef(nid_day, "snow", "Snow fall", "mm/day",
802	. iim,jjmp1,nhori, 1,1,1, -99, 32,
803	. "ave(X)", zsto,zout)
804	c
805	CALL histdef(nid_day, "evap", "Evaporation", "mm/day",
806	. iim,jjmp1,nhori, 1,1,1, -99, 32,
807	. "ave(X)", zsto,zout)
808	c
809	CALL histdef(nid_day, "tops", "Solar rad. at TOA", "W/m2",
810	. iim,jjmp1,nhori, 1,1,1, -99, 32,
811	. "ave(X)", zsto,zout)
812	c
813	CALL histdef(nid_day, "topl", "IR rad. at TOA", "W/m2",
814	. iim,jjmp1,nhori, 1,1,1, -99, 32,
815	. "ave(X)", zsto,zout)
816	c
817	CALL histdef(nid_day, "sols", "Solar rad. at surf.", "W/m2",
818	. iim,jjmp1,nhori, 1,1,1, -99, 32,
819	. "ave(X)", zsto,zout)
820	c
821	CALL histdef(nid_day, "soll", "IR rad. at surface", "W/m2",
822	. iim,jjmp1,nhori, 1,1,1, -99, 32,
823	. "ave(X)", zsto,zout)
824	c
825	CALL histdef(nid_day, "bils", "Surf. total heat flux", "W/m2",
826	. iim,jjmp1,nhori, 1,1,1, -99, 32,
827	. "ave(X)", zsto,zout)
828	c
829	CALL histdef(nid_day, "sens", "Sensible heat flux", "W/m2",
830	. iim,jjmp1,nhori, 1,1,1, -99, 32,
831	. "ave(X)", zsto,zout)
832	c
833	CALL histdef(nid_day, "fder", "Heat flux derivation", "W/m2",
834	. iim,jjmp1,nhori, 1,1,1, -99, 32,
835	. "ave(X)", zsto,zout)
836	c
837	CALL histdef(nid_day, "frtu", "Zonal wind stress", "Pa",
838	. iim,jjmp1,nhori, 1,1,1, -99, 32,
839	. "ave(X)", zsto,zout)
840	c
841	CALL histdef(nid_day, "frtv", "Meridional wind stress", "Pa",
842	. iim,jjmp1,nhori, 1,1,1, -99, 32,
843	. "ave(X)", zsto,zout)
844	c
845	CALL histdef(nid_day, "ruis", "Runoff", "mm/day",
846	. iim,jjmp1,nhori, 1,1,1, -99, 32,
847	. "ave(X)", zsto,zout)
848	c
849	CALL histdef(nid_day, "sicf", "Sea-ice fraction", "-",
850	. iim,jjmp1,nhori, 1,1,1, -99, 32,
851	. "ave(X)", zsto,zout)
852	c
853	CALL histdef(nid_day, "cldl", "Low-level cloudiness", "-",
854	. iim,jjmp1,nhori, 1,1,1, -99, 32,
855	. "ave(X)", zsto,zout)
856	c
857	CALL histdef(nid_day, "cldm", "Mid-level cloudiness", "-",
858	. iim,jjmp1,nhori, 1,1,1, -99, 32,
859	. "ave(X)", zsto,zout)
860	c
861	CALL histdef(nid_day, "cldh", "High-level cloudiness", "-",
862	. iim,jjmp1,nhori, 1,1,1, -99, 32,
863	. "ave(X)", zsto,zout)
864	c
865	CALL histdef(nid_day, "cldt", "Total cloudiness", "-",
866	. iim,jjmp1,nhori, 1,1,1, -99, 32,
867	. "ave(X)", zsto,zout)
868	c
869	CALL histdef(nid_day, "cldq", "Cloud liquid water path", "-",
870	. iim,jjmp1,nhori, 1,1,1, -99, 32,
871	. "ave(X)", zsto,zout)
872	c
873	c Champs 3D:
874	c
875	CALL histdef(nid_day, "temp", "Air temperature", "K",
876	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
877	. "ave(X)", zsto,zout)
878	c
879	CALL histdef(nid_day, "ovap", "Specific humidity", "Kg/Kg",
880	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
881	. "ave(X)", zsto,zout)
882	c
883	CALL histdef(nid_day, "geop", "Geopotential height", "m",
884	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
885	. "ave(X)", zsto,zout)
886	c
887	CALL histdef(nid_day, "vitu", "Zonal wind", "m/s",
888	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
889	. "ave(X)", zsto,zout)
890	c
891	CALL histdef(nid_day, "vitv", "Meridional wind", "m/s",
892	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
893	. "ave(X)", zsto,zout)
894	c
895	CALL histdef(nid_day, "vitw", "Vertical wind", "m/s",
896	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
897	. "ave(X)", zsto,zout)
898	c
899	CALL histdef(nid_day, "pres", "Air pressure", "Pa",
900	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
901	. "ave(X)", zsto,zout)
902	c
903	CALL histend(nid_day)
904	c
905	ndex2d = 0
906	ndex3d = 0
907	c
908	ENDIF ! fin de test sur ok_journe
909	c
910	IF (ok_mensuel) THEN
911	c
912	CALL ymds2ju(anne_ini, 1, 1, 0.0, zjulian)
913	zjulian = zjulian + day_ini
914	c
915	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
916	DO i = 1, iim
917	zx_lon(i,1) = rlon(i+1)
918	zx_lon(i,jjmp1) = rlon(i+1)
919	ENDDO
920	DO ll=1,klev
921	znivsig(ll)=float(ll)
922	ENDDO
923	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
924	CALL histbeg("histmth", iim,zx_lon, jjmp1,zx_lat,
925	. 1,iim,1,jjmp1, 0, zjulian, dtime,
926	. nhori, nid_mth)
927	CALL histvert(nid_mth, "presnivs", "Vertical levels", "mb",
928	. klev, presnivs, nvert)
929	c call histvert(nid_mth, 'sig_s', 'Niveaux sigma','-',
930	c . klev, znivsig, nvert)
931	c
932	zsto = dtime
933	zout = dtime * FLOAT(ecrit_mth)
934	c
935	CALL histdef(nid_mth, "phis", "Surface geop. height", "-",
936	. iim,jjmp1,nhori, 1,1,1, -99, 32,
937	. "once", zsto,zout)
938	c
939	CALL histdef(nid_mth, "aire", "Grid area", "-",
940	. iim,jjmp1,nhori, 1,1,1, -99, 32,
941	. "once", zsto,zout)
942	c
943	c Champs 2D:
944	c
945	CALL histdef(nid_mth, "tsol", "Surface Temperature", "K",
946	. iim,jjmp1,nhori, 1,1,1, -99, 32,
947	. "ave(X)", zsto,zout)
948	c
949	CALL histdef(nid_mth, "psol", "Surface Pressure", "Pa",
950	. iim,jjmp1,nhori, 1,1,1, -99, 32,
951	. "ave(X)", zsto,zout)
952	c
953	CALL histdef(nid_mth, "qsol", "Surface humidity", "mm",
954	. iim,jjmp1,nhori, 1,1,1, -99, 32,
955	. "ave(X)", zsto,zout)
956	c
957	CALL histdef(nid_mth, "rain", "Precipitation", "mm/day",
958	. iim,jjmp1,nhori, 1,1,1, -99, 32,
959	. "ave(X)", zsto,zout)
960	c
961	CALL histdef(nid_mth, "plul", "Large-scale Precip.", "mm/day",
962	. iim,jjmp1,nhori, 1,1,1, -99, 32,
963	. "ave(X)", zsto,zout)
964	c
965	CALL histdef(nid_mth, "pluc", "Convective Precip.", "mm/day",
966	. iim,jjmp1,nhori, 1,1,1, -99, 32,
967	. "ave(X)", zsto,zout)
968	c
969	CALL histdef(nid_mth, "snow", "Snow fall", "mm/day",
970	. iim,jjmp1,nhori, 1,1,1, -99, 32,
971	. "ave(X)", zsto,zout)
972	c
973	CALL histdef(nid_mth, "ages", "Snow age", "day",
974	. iim,jjmp1,nhori, 1,1,1, -99, 32,
975	. "ave(X)", zsto,zout)
976	c
977	CALL histdef(nid_mth, "evap", "Evaporation", "mm/day",
978	. iim,jjmp1,nhori, 1,1,1, -99, 32,
979	. "ave(X)", zsto,zout)
980	c
981	CALL histdef(nid_mth, "tops", "Solar rad. at TOA", "W/m2",
982	. iim,jjmp1,nhori, 1,1,1, -99, 32,
983	. "ave(X)", zsto,zout)
984	c
985	CALL histdef(nid_mth, "topl", "IR rad. at TOA", "W/m2",
986	. iim,jjmp1,nhori, 1,1,1, -99, 32,
987	. "ave(X)", zsto,zout)
988	c
989	CALL histdef(nid_mth, "sols", "Solar rad. at surf.", "W/m2",
990	. iim,jjmp1,nhori, 1,1,1, -99, 32,
991	. "ave(X)", zsto,zout)
992	c
993	CALL histdef(nid_mth, "soll", "IR rad. at surface", "W/m2",
994	. iim,jjmp1,nhori, 1,1,1, -99, 32,
995	. "ave(X)", zsto,zout)
996	c
997	CALL histdef(nid_mth, "tops0", "Solar rad. at TOA", "W/m2",
998	. iim,jjmp1,nhori, 1,1,1, -99, 32,
999	. "ave(X)", zsto,zout)
1000	c
1001	CALL histdef(nid_mth, "topl0", "IR rad. at TOA", "W/m2",
1002	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1003	. "ave(X)", zsto,zout)
1004	c
1005	CALL histdef(nid_mth, "sols0", "Solar rad. at surf.", "W/m2",
1006	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1007	. "ave(X)", zsto,zout)
1008	c
1009	CALL histdef(nid_mth, "soll0", "IR rad. at surface", "W/m2",
1010	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1011	. "ave(X)", zsto,zout)
1012	c
1013	CALL histdef(nid_mth, "bils", "Surf. total heat flux", "W/m2",
1014	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1015	. "ave(X)", zsto,zout)
1016	c
1017	CALL histdef(nid_mth, "sens", "Sensible heat flux", "W/m2",
1018	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1019	. "ave(X)", zsto,zout)
1020	c
1021	CALL histdef(nid_mth, "fder", "Heat flux derivation", "W/m2",
1022	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1023	. "ave(X)", zsto,zout)
1024	c
1025	CALL histdef(nid_mth, "frtu", "Zonal wind stress", "Pa",
1026	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1027	. "ave(X)", zsto,zout)
1028	c
1029	CALL histdef(nid_mth, "frtv", "Meridional wind stress", "Pa",
1030	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1031	. "ave(X)", zsto,zout)
1032	c
1033	CALL histdef(nid_mth, "ruis", "Runoff", "mm/day",
1034	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1035	. "ave(X)", zsto,zout)
1036	c
1037	CALL histdef(nid_mth, "sicf", "Sea-ice fraction", "-",
1038	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1039	. "ave(X)", zsto,zout)
1040	c
1041	CALL histdef(nid_mth, "albs", "Surface albedo", "-",
1042	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1043	. "ave(X)", zsto,zout)
1044	c
1045	CALL histdef(nid_mth, "cdrm", "Momentum drag coef.", "-",
1046	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1047	. "ave(X)", zsto,zout)
1048	c
1049	CALL histdef(nid_mth, "cdrh", "Heat drag coef.", "-",
1050	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1051	. "ave(X)", zsto,zout)
1052	c
1053	CALL histdef(nid_mth, "cldl", "Low-level cloudiness", "-",
1054	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1055	. "ave(X)", zsto,zout)
1056	c
1057	CALL histdef(nid_mth, "cldm", "Mid-level cloudiness", "-",
1058	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1059	. "ave(X)", zsto,zout)
1060	c
1061	CALL histdef(nid_mth, "cldh", "High-level cloudiness", "-",
1062	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1063	. "ave(X)", zsto,zout)
1064	c
1065	CALL histdef(nid_mth, "cldt", "Total cloudiness", "-",
1066	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1067	. "ave(X)", zsto,zout)
1068	c
1069	CALL histdef(nid_mth, "cldq", "Cloud liquid water path", "-",
1070	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1071	. "ave(X)", zsto,zout)
1072	c
1073	CALL histdef(nid_mth, "ue", "Zonal energy transport", "-",
1074	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1075	. "ave(X)", zsto,zout)
1076	c
1077	CALL histdef(nid_mth, "ve", "Merid energy transport", "-",
1078	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1079	. "ave(X)", zsto,zout)
1080	c
1081	CALL histdef(nid_mth, "uq", "Zonal humidity transport", "-",
1082	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1083	. "ave(X)", zsto,zout)
1084	c
1085	CALL histdef(nid_mth, "vq", "Merid humidity transport", "-",
1086	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1087	. "ave(X)", zsto,zout)
1088	cKE43
1089	IF (iflag_con .GE. 3) THEN ! sb
1090	c
1091	CALL histdef(nid_mth, "cape", "Conv avlbl pot ener", "J/Kg",
1092	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1093	. "ave(X)", zsto,zout)
1094	c
1095	CALL histdef(nid_mth, "pbase", "Cld base pressure", "hPa",
1096	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1097	. "ave(X)", zsto,zout)
1098	c
1099	CALL histdef(nid_mth, "ptop", "Cld top pressure", "hPa",
1100	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1101	. "ave(X)", zsto,zout)
1102	c
1103	CALL histdef(nid_mth, "fbase", "Cld base mass flux", "Kg/m2/s",
1104	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1105	. "ave(X)", zsto,zout)
1106	c
1107	c
1108	ENDIF
1109	c34EK
1110	c
1111	c Champs 3D:
1112	c
1113	CALL histdef(nid_mth, "temp", "Air temperature", "K",
1114	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1115	. "ave(X)", zsto,zout)
1116	c
1117	CALL histdef(nid_mth, "ovap", "Specific humidity", "Kg/Kg",
1118	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1119	. "ave(X)", zsto,zout)
1120	c
1121	CALL histdef(nid_mth, "geop", "Geopotential height", "m",
1122	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1123	. "ave(X)", zsto,zout)
1124	c
1125	CALL histdef(nid_mth, "vitu", "Zonal wind", "m/s",
1126	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1127	. "ave(X)", zsto,zout)
1128	c
1129	CALL histdef(nid_mth, "vitv", "Meridional wind", "m/s",
1130	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1131	. "ave(X)", zsto,zout)
1132	c
1133	CALL histdef(nid_mth, "vitw", "Vertical wind", "m/s",
1134	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1135	. "ave(X)", zsto,zout)
1136	c
1137	CALL histdef(nid_mth, "pres", "Air pressure", "Pa",
1138	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1139	. "ave(X)", zsto,zout)
1140	c
1141	CALL histdef(nid_mth, "rneb", "Cloud fraction", "-",
1142	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1143	. "ave(X)", zsto,zout)
1144	c
1145	CALL histdef(nid_mth, "rhum", "Relative humidity", "-",
1146	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1147	. "ave(X)", zsto,zout)
1148	c
1149	CALL histdef(nid_mth, "oliq", "Liquid water content", "kg/kg",
1150	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1151	. "ave(X)", zsto,zout)
1152	c
1153	CALL histdef(nid_mth, "dtdyn", "Dynamics dT", "K/s",
1154	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1155	. "ave(X)", zsto,zout)
1156	c
1157	CALL histdef(nid_mth, "dqdyn", "Dynamics dQ", "Kg/Kg/s",
1158	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1159	. "ave(X)", zsto,zout)
1160	c
1161	CALL histdef(nid_mth, "dtcon", "Convection dT", "K/s",
1162	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1163	. "ave(X)", zsto,zout)
1164	c
1165	CALL histdef(nid_mth, "dqcon", "Convection dQ", "Kg/Kg/s",
1166	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1167	. "ave(X)", zsto,zout)
1168	c
1169	CALL histdef(nid_mth, "dtlsc", "Condensation dT", "K/s",
1170	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1171	. "ave(X)", zsto,zout)
1172	c
1173	CALL histdef(nid_mth, "dqlsc", "Condensation dQ", "Kg/Kg/s",
1174	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1175	. "ave(X)", zsto,zout)
1176	c
1177	CALL histdef(nid_mth, "dtvdf", "Boundary-layer dT", "K/s",
1178	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1179	. "ave(X)", zsto,zout)
1180	c
1181	CALL histdef(nid_mth, "dqvdf", "Boundary-layer dQ", "Kg/Kg/s",
1182	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1183	. "ave(X)", zsto,zout)
1184	c
1185	CALL histdef(nid_mth, "dteva", "Reevaporation dT", "K/s",
1186	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1187	. "ave(X)", zsto,zout)
1188	c
1189	CALL histdef(nid_mth, "dqeva", "Reevaporation dQ", "Kg/Kg/s",
1190	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1191	. "ave(X)", zsto,zout)
1192
1193	CALL histdef(nid_mth, "ptconv", "POINTS CONVECTIFS"," ",
1194	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1195	. "ave(X)", zsto,zout)
1196
1197	CALL histdef(nid_mth, "ratqs", "RATQS"," ",
1198	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1199	. "ave(X)", zsto,zout)
1200
1201	c
1202	CALL histdef(nid_mth, "dtajs", "Dry adjust. dT", "K/s",
1203	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1204	. "ave(X)", zsto,zout)
1205
1206	CALL histdef(nid_mth, "dqajs", "Dry adjust. dQ", "Kg/Kg/s",
1207	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1208	. "ave(X)", zsto,zout)
1209	c
1210	CALL histdef(nid_mth, "dtswr", "SW radiation dT", "K/s",
1211	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1212	. "ave(X)", zsto,zout)
1213	c
1214	CALL histdef(nid_mth, "dtsw0", "SW radiation dT", "K/s",
1215	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1216	. "ave(X)", zsto,zout)
1217	c
1218	CALL histdef(nid_mth, "dtlwr", "LW radiation dT", "K/s",
1219	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1220	. "ave(X)", zsto,zout)
1221	c
1222	CALL histdef(nid_mth, "dtlw0", "LW radiation dT", "K/s",
1223	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1224	. "ave(X)", zsto,zout)
1225	c
1226	CALL histdef(nid_mth, "duvdf", "Boundary-layer dU", "m/s2",
1227	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1228	. "ave(X)", zsto,zout)
1229	c
1230	CALL histdef(nid_mth, "dvvdf", "Boundary-layer dV", "m/s2",
1231	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1232	. "ave(X)", zsto,zout)
1233	c
1234	IF (ok_orodr) THEN
1235	CALL histdef(nid_mth, "duoro", "Orography dU", "m/s2",
1236	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1237	. "ave(X)", zsto,zout)
1238	c
1239	CALL histdef(nid_mth, "dvoro", "Orography dV", "m/s2",
1240	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1241	. "ave(X)", zsto,zout)
1242	c
1243	ENDIF
1244	C
1245	IF (ok_orolf) THEN
1246	CALL histdef(nid_mth, "dulif", "Orography dU", "m/s2",
1247	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1248	. "ave(X)", zsto,zout)
1249	c
1250	CALL histdef(nid_mth, "dvlif", "Orography dV", "m/s2",
1251	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1252	. "ave(X)", zsto,zout)
1253	ENDIF
1254	C
1255	CALL histdef(nid_mth, "ozone", "Ozone concentration", "-",
1256	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1257	. "ave(X)", zsto,zout)
1258	c
1259	if (nqmax.GE.3) THEN
1260	DO iq=1,nqmax-2
1261	IF (iq.LE.99) THEN
1262	WRITE(str2,'(i2.2)') iq
1263	CALL histdef(nid_mth, "trac"//str2, "Tracer No."//str2, "-",
1264	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1265	. "ave(X)", zsto,zout)
1266	ELSE
1267	PRINT*, "Trop de traceurs"
1268	CALL abort
1269	ENDIF
1270	ENDDO
1271	ENDIF
1272	c
1273	cKE43
1274	IF (iflag_con.GE.3) THEN ! (sb)
1275	c
1276	CALL histdef(nid_mth, "upwd", "saturated updraft", "Kg/m2/s",
1277	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1278	. "ave(X)", zsto,zout)
1279	c
1280	CALL histdef(nid_mth, "dnwd", "saturated downdraft","Kg/m2/s",
1281	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1282	. "ave(X)", zsto,zout)
1283	c
1284	CALL histdef(nid_mth, "dnwd0", "unsat. downdraft", "Kg/m2/s",
1285	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1286	. "ave(X)", zsto,zout)
1287	c
1288	CALL histdef(nid_mth,"Ma","undilute adiab updraft","Kg/m2/s",
1289	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1290	. "ave(X)", zsto,zout)
1291	c
1292	c
1293	ENDIF
1294	c34EK
1295	CALL histend(nid_mth)
1296	c
1297	ndex2d = 0
1298	ndex3d = 0
1299	c
1300	ENDIF ! fin de test sur ok_mensuel
1301	c
1302	c
1303	IF (ok_instan) THEN
1304	c
1305	CALL ymds2ju(anne_ini, 1, 1, 0.0, zjulian)
1306	zjulian = zjulian + day_ini
1307	c
1308	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
1309	DO i = 1, iim
1310	zx_lon(i,1) = rlon(i+1)
1311	zx_lon(i,jjmp1) = rlon(i+1)
1312	ENDDO
1313	DO ll=1,klev
1314	znivsig(ll)=float(ll)
1315	ENDDO
1316	CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
1317	CALL histbeg("histins", iim,zx_lon, jjmp1,zx_lat,
1318	. 1,iim,1,jjmp1, 0, zjulian, dtime,
1319	. nhori, nid_ins)
1320	CALL histvert(nid_ins, "presnivs", "Vertical levels", "mb",
1321	. klev, presnivs, nvert)
1322	c call histvert(nid_ins, 'sig_s', 'Niveaux sigma','-',
1323	c . klev, znivsig, nvert)
1324	c
1325	zsto = dtime * ecrit_ins
1326	zout = dtime * ecrit_ins
1327	C
1328	CALL histdef(nid_ins, "phis", "Surface geop. height", "-",
1329	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1330	. "once", zsto,zout)
1331	c
1332	CALL histdef(nid_ins, "aire", "Grid area", "-",
1333	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1334	. "once", zsto,zout)
1335	c
1336	c Champs 2D:
1337	c
1338	CALL histdef(nid_ins, "psol", "Surface Pressure", "Pa",
1339	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1340	. "inst(X)", zsto,zout)
1341	c
1342	CALL histdef(nid_ins, "topl", "OLR", "W/m2",
1343	. iim,jjmp1,nhori, 1,1,1, -99, 32,
1344	. "inst(X)", zsto,zout)
1345	c
1346	c Champs 3D:
1347	c
1348	CALL histdef(nid_ins, "temp", "Temperature", "K",
1349	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1350	. "inst(X)", zsto,zout)
1351	c
1352	CALL histdef(nid_ins, "vitu", "Zonal wind", "m/s",
1353	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1354	. "inst(X)", zsto,zout)
1355	c
1356	CALL histdef(nid_ins, "vitv", "Merid wind", "m/s",
1357	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1358	. "inst(X)", zsto,zout)
1359	c
1360	CALL histdef(nid_ins, "geop", "Geopotential height", "m",
1361	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1362	. "inst(X)", zsto,zout)
1363	c
1364	CALL histdef(nid_ins, "pres", "Air pressure", "Pa",
1365	. iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1366	. "inst(X)", zsto,zout)
1367	c
1368	CALL histend(nid_ins)
1369	c
1370	ndex2d = 0
1371	ndex3d = 0
1372	c
1373	ENDIF
1374	c
1375	c
1376	c
1377	c Prescrire l'ozone dans l'atmosphere
1378	c
1379	c
1380	cc DO i = 1, klon
1381	cc DO k = 1, klev
1382	cc CALL o3cm (paprs(i,k)/100.,paprs(i,k+1)/100., wo(i,k),20)
1383	cc ENDDO
1384	cc ENDDO
1385	c
1386	IF (ok_oasis) THEN
1387	DO i = 1, klon
1388	oas_sols(i) = 0.0
1389	oas_nsol(i) = 0.0
1390	oas_rain(i) = 0.0
1391	oas_snow(i) = 0.0
1392	oas_evap(i) = 0.0
1393	oas_ruis(i) = 0.0
1394	oas_tsol(i) = 0.0
1395	oas_fder(i) = 0.0
1396	oas_albe(i) = 0.0
1397	oas_taux(i) = 0.0
1398	oas_tauy(i) = 0.0
1399	ENDDO
1400	ENDIF
1401	c
1402	ENDIF
1403	c
1404	c ************** Fin de IF ( debut ) *************
1405	c
1406	c
1407	c Mettre a zero des variables de sortie (pour securite)
1408	c
1409	DO i = 1, klon
1410	d_ps(i) = 0.0
1411	ENDDO
1412	DO k = 1, klev
1413	DO i = 1, klon
1414	d_t(i,k) = 0.0
1415	d_u(i,k) = 0.0
1416	d_v(i,k) = 0.0
1417	ENDDO
1418	ENDDO
1419	DO iq = 1, nqmax
1420	DO k = 1, klev
1421	DO i = 1, klon
1422	d_qx(i,k,iq) = 0.0
1423	ENDDO
1424	ENDDO
1425	ENDDO
1426	c
1427	c Ne pas affecter les valeurs entrees de u, v, h, et q
1428	c
1429	DO k = 1, klev
1430	DO i = 1, klon
1431	t_seri(i,k) = t(i,k)
1432	u_seri(i,k) = u(i,k)
1433	v_seri(i,k) = v(i,k)
1434	q_seri(i,k) = qx(i,k,ivap)
1435	ql_seri(i,k) = qx(i,k,iliq)
1436	ENDDO
1437	ENDDO
1438	IF (nqmax.GE.3) THEN
1439	DO iq = 3, nqmax
1440	DO k = 1, klev
1441	DO i = 1, klon
1442	tr_seri(i,k,iq-2) = qx(i,k,iq)
1443	ENDDO
1444	ENDDO
1445	ENDDO
1446	ELSE
1447	DO k = 1, klev
1448	DO i = 1, klon
1449	tr_seri(i,k,1) = 0.0
1450	ENDDO
1451	ENDDO
1452	ENDIF
1453	c
1454	c Diagnostiquer la tendance dynamique
1455	c
1456	IF (ancien_ok) THEN
1457	DO k = 1, klev
1458	DO i = 1, klon
1459	d_t_dyn(i,k) = (t_seri(i,k)-t_ancien(i,k))/dtime
1460	d_q_dyn(i,k) = (q_seri(i,k)-q_ancien(i,k))/dtime
1461	ENDDO
1462	ENDDO
1463	ELSE
1464	DO k = 1, klev
1465	DO i = 1, klon
1466	d_t_dyn(i,k) = 0.0
1467	d_q_dyn(i,k) = 0.0
1468	ENDDO
1469	ENDDO
1470	ancien_ok = .TRUE.
1471	ENDIF
1472	c
1473	c Ajouter le geopotentiel du sol:
1474	c
1475	DO k = 1, klev
1476	DO i = 1, klon
1477	zphi(i,k) = pphi(i,k) + pphis(i)
1478	ENDDO
1479	ENDDO
1480	c
1481	c Verifier les temperatures
1482	c
1483
1484	CALL hgardfou(t_seri,ftsol,'debutphy')
1485	c
1486	c Incrementer le compteur de la physique
1487	c
1488	itap = itap + 1
1489	julien = MOD(NINT(xjour),360)
1490	c
1491	c Mettre en action les conditions aux limites (albedo, sst, etc.).
1492	c Prescrire l'ozone et calculer l'albedo sur l'ocean.
1493	c
1494	IF (MOD(itap-1,lmt_pas) .EQ. 0) THEN
1495	idayvrai = NINT(xjour)
1496	PRINT *,' PHYS cond julien ',julien,idayvrai,ok_limitvrai
1497	CALL condsurf(julien,idayvrai, pctsrf ,
1498	. lmt_sst,lmt_alb,lmt_rug,lmt_bils )
1499	CALL ozonecm( FLOAT(julien), rlat, paprs, wo)
1500	ENDIF
1501	cccccccccc
1502	IF (ok_oasis .AND. MOD(itap-1,nexca).EQ.0) THEN
1503	C
1504	CALL fromcpl(itap,jjmp1*iim,
1505	. cpl_sst,cpl_sic,cpl_alb_sst,cpl_alb_sic)
1506	DO i = 1, iim-1 ! un seul point pour le pole nord
1507	cpl_sst(i,1) = cpl_sst(iim,1)
1508	cpl_sic(i,1) = cpl_sic(iim,1)
1509	cpl_alb_sst(i,1) = cpl_alb_sst(iim,1)
1510	cpl_alb_sic(i,1) = cpl_alb_sic(iim,1)
1511	ENDDO
1512	DO i = 2, iim ! un seul point pour le pole sud
1513	cpl_sst(i,jjmp1) = cpl_sst(1,jjmp1)
1514	cpl_sic(i,jjmp1) = cpl_sic(1,jjmp1)
1515	cpl_alb_sst(i,jjmp1) = cpl_alb_sst(1,jjmp1)
1516	cpl_alb_sic(i,jjmp1) = cpl_alb_sic(1,jjmp1)
1517	ENDDO
1518	c
1519	ig = 1
1520	IF (pctsrf(ig,is_oce).GT.epsfra .OR.
1521	. pctsrf(ig,is_sic).GT.epsfra) THEN
1522	pctsrf(ig,is_oce) = pctsrf(ig,is_oce)
1523	. - (cpl_sic(1,1)-pctsrf(ig,is_sic))
1524	pctsrf(ig,is_sic) = cpl_sic(1,1)
1525	lmt_sst(ig) = cpl_sst(1,1)
1526	ENDIF
1527	DO j = 2, jjm
1528	DO i = 1, iim
1529	ig = ig + 1
1530	IF (pctsrf(ig,is_oce).GT.epsfra .OR.
1531	. pctsrf(ig,is_sic).GT.epsfra) THEN
1532	pctsrf(ig,is_oce) = pctsrf(ig,is_oce)
1533	. - (cpl_sic(i,j)-pctsrf(ig,is_sic))
1534	pctsrf(ig,is_sic) = cpl_sic(i,j)
1535	lmt_sst(ig) = cpl_sst(i,j)
1536	ENDIF
1537	ENDDO
1538	ENDDO
1539	ig = ig + 1
1540	IF (pctsrf(ig,is_oce).GT.epsfra .OR.
1541	. pctsrf(ig,is_sic).GT.epsfra) THEN
1542	pctsrf(ig,is_oce) = pctsrf(ig,is_oce)
1543	. - (cpl_sic(1,jjmp1)-pctsrf(ig,is_sic))
1544	pctsrf(ig,is_sic) = cpl_sic(1,jjmp1)
1545	lmt_sst(ig) = cpl_sst(1,jjmp1)
1546	ENDIF
1547	c
1548	ENDIF ! ok_oasis
1549	cccccccccc
1550	c
1551
1552	IF (ok_ocean) THEN
1553	DO i = 1, klon
1554	ftsol(i,is_oce) = lmt_sst(i) + deltat(i)
1555	ENDDO
1556
1557	ELSE
1558	DO i = 1, klon
1559	ftsol(i,is_oce) = lmt_sst(i)
1560	ENDDO
1561
1562	ENDIF
1563	c
1564	c Re-evaporer l'eau liquide nuageuse
1565	c
1566	DO k = 1, klev ! re-evaporation de l'eau liquide nuageuse
1567	DO i = 1, klon
1568	zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
1569	zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
1570	zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
1571	zb = MAX(0.0,ql_seri(i,k))
1572	za = - MAX(0.0,ql_seri(i,k))
1573	. * (zlvdcp(1.-zdelta)+zlsdcpzdelta)
1574	t_seri(i,k) = t_seri(i,k) + za
1575	q_seri(i,k) = q_seri(i,k) + zb
1576	ql_seri(i,k) = 0.0
1577	d_t_eva(i,k) = za
1578	d_q_eva(i,k) = zb
1579	ENDDO
1580	ENDDO
1581	c
1582	c Appeler la diffusion verticale (programme de couche limite)
1583	c
1584	DO i = 1, klon
1585	frugs(i,is_ter) = SQRT(lmt_rug(i)2+rugoro(i)2)
1586	frugs(i,is_lic) = rugoro(i)
1587	frugs(i,is_oce) = rugmer(i)
1588	frugs(i,is_sic) = 0.001
1589	zxrugs(i) = 0.0
1590	ENDDO
1591	DO nsrf = 1, nbsrf
1592	DO i = 1, klon
1593	frugs(i,nsrf) = MAX(frugs(i,nsrf),0.001)
1594	ENDDO
1595	ENDDO
1596	DO nsrf = 1, nbsrf
1597	DO i = 1, klon
1598	zxrugs(i) = zxrugs(i) + frugs(i,nsrf)*pctsrf(i,nsrf)
1599	ENDDO
1600	ENDDO
1601	c
1602	CALL clmain(dtime,pctsrf,
1603	e t_seri,q_seri,u_seri,v_seri,soil_model,
1604	e ftsol,soilcap,soilflux,paprs,pplay,radsol,
1605	e fsnow,fqsol,
1606	e rlat, frugs,
1607	s d_t_vdf,d_q_vdf,d_u_vdf,d_v_vdf,d_ts,
1608	s fluxt,fluxq,fluxu,fluxv,cdragh,cdragm,rugmer,
1609	s dsens, devap,
1610	s ycoefh,yu1,yv1)
1611	c
1612	DO i = 1, klon
1613	sens(i) = - fluxt(i,1) ! flux de chaleur sensible au sol
1614	evap(i) = - fluxq(i,1) ! flux d'evaporation au sol
1615	fder(i) = dsens(i) + devap(i)
1616	ENDDO
1617	DO k = 1, klev
1618	DO i = 1, klon
1619	t_seri(i,k) = t_seri(i,k) + d_t_vdf(i,k)
1620	q_seri(i,k) = q_seri(i,k) + d_q_vdf(i,k)
1621	u_seri(i,k) = u_seri(i,k) + d_u_vdf(i,k)
1622	v_seri(i,k) = v_seri(i,k) + d_v_vdf(i,k)
1623	ENDDO
1624	ENDDO
1625	c
1626	c Incrementer la temperature du sol
1627	c
1628	DO i = 1, klon
1629	zxtsol(i) = 0.0
1630	IF ( abs( pctsrf(i, is_ter) + pctsrf(i, is_lic) +
1631	$ pctsrf(i, is_oce) + pctsrf(i, is_sic) - 1.) .GT. EPSFRA)
1632	$ THEN
1633	WRITE(,) 'physiq : pb sous surface au point ', i,
1634	$ pctsrf(i, 1 : nbsrf)
1635	ENDIF
1636	ENDDO
1637	DO nsrf = 1, nbsrf
1638	DO i = 1, klon
1639	ftsol(i,nsrf) = ftsol(i,nsrf) + d_ts(i,nsrf)
1640	zxtsol(i) = zxtsol(i) + ftsol(i,nsrf)*pctsrf(i,nsrf)
1641	ENDDO
1642	ENDDO
1643
1644	c
1645	c Si une sous-fraction n'existe pas, elle prend la temp. moyenne
1646	c
1647	DO nsrf = 1, nbsrf
1648	DO i = 1, klon
1649	IF (pctsrf(i,nsrf) .LT. epsfra) ftsol(i,nsrf) = zxtsol(i)
1650	ENDDO
1651	ENDDO
1652
1653	c
1654	c Appeler le modele du sol
1655	c
1656	IF (soil_model) THEN
1657	DO nsrf = 1, nbsrf
1658	CALL soil(dtime, nsrf, fsnow(1,nsrf),
1659	. ftsol(1,nsrf), ftsoil(1,1,nsrf),
1660	. soilcap(1,nsrf), soilflux(1,nsrf))
1661	ENDDO
1662	ENDIF
1663	c
1664	c Calculer la derive du flux infrarouge
1665	c
1666	DO nsrf = 1, nbsrf
1667	DO i = 1, klon
1668	fder(i) = fder(i) - 4.0RSIGMAzxtsol(i)*3
1669	. (ftsol(i,nsrf)-zxtsol(i))
1670	. *pctsrf(i,nsrf)
1671	ENDDO
1672	ENDDO
1673	c
1674	c Appeler la convection (au choix)
1675	c
1676	DO k = 1, klev
1677	DO i = 1, klon
1678	conv_q(i,k) = d_q_dyn(i,k)
1679	. + d_q_vdf(i,k)/dtime
1680	conv_t(i,k) = d_t_dyn(i,k)
1681	. + d_t_vdf(i,k)/dtime
1682	ENDDO
1683	ENDDO
1684	IF (check) THEN
1685	za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
1686	PRINT*, "avantcon=", za
1687	ENDIF
1688	zx_ajustq = .FALSE.
1689	IF (iflag_con.EQ.2) zx_ajustq=.TRUE.
1690	IF (zx_ajustq) THEN
1691	DO i = 1, klon
1692	z_avant(i) = 0.0
1693	ENDDO
1694	DO k = 1, klev
1695	DO i = 1, klon
1696	z_avant(i) = z_avant(i) + (q_seri(i,k)+ql_seri(i,k))
1697	. *(paprs(i,k)-paprs(i,k+1))/RG
1698	ENDDO
1699	ENDDO
1700	ENDIF
1701	IF (iflag_con.EQ.1) THEN
1702	stop'reactiver le call conlmd dans physiq.F'
1703	c CALL conlmd (dtime, paprs, pplay, t_seri, q_seri, conv_q,
1704	c . d_t_con, d_q_con,
1705	c . rain_con, snow_con, ibas_con, itop_con)
1706	ELSE IF (iflag_con.EQ.2) THEN
1707	CALL conflx(dtime, paprs, pplay, t_seri, q_seri,
1708	e conv_t, conv_q, fluxq(1,1), omega,
1709	s d_t_con, d_q_con, rain_con, snow_con,
1710	s pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
1711	s kcbot, kctop, kdtop, pmflxr, pmflxs)
1712	DO i = 1, klon
1713	ibas_con(i) = klev+1 - kcbot(i)
1714	itop_con(i) = klev+1 - kctop(i)
1715	ENDDO
1716	ELSE IF (iflag_con.GE.3) THEN
1717	c nb of tracers for the KE convection:
1718	if (nqmax .GE. 4) then
1719	ntra = nbtr
1720	else
1721	ntra = 1
1722	endif
1723	if (iflag_con.eq.4) then ! vectorise
1724	CALL conemav (dtime,paprs,pplay,t_seri,q_seri,
1725	. u_seri,v_seri,tr_seri,nbtr,
1726	. ema_work1,ema_work2,
1727	. d_t_con,d_q_con,d_u_con,d_v_con,d_tr,
1728	. rain_con, snow_con, ibas_con, itop_con,
1729	. upwd,dnwd,dnwd0,
1730	c . Ma,cape,tvp,(/(nint(rflag(i)),i=1,size(rflag))/),
1731	. Ma,cape,tvp,iflagctrl,
1732	. pbase,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr)
1733
1734	else
1735
1736	CALL conema (dtime,paprs,pplay,t_seri,q_seri,
1737	. u_seri,v_seri,tr_seri,nbtr,
1738	. ema_work1,ema_work2,
1739	. d_t_con,d_q_con,d_u_con,d_v_con,d_tr,
1740	. rain_con, snow_con, ibas_con, itop_con,
1741	. upwd,dnwd,dnwd0,bas,top,
1742	. Ma,cape,tvp,rflag,
1743	. pbase
1744	. ,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr)
1745	endif
1746
1747
1748	DO i = 1, klon
1749	ema_pcb(i) = pbase(i)
1750	ENDDO
1751	DO i = 1, klon
1752	ema_pct(i) = paprs(i,itop_con(i))
1753	ENDDO
1754	DO i = 1, klon
1755	ema_cbmf(i) = ema_workcbmf(i)
1756	ENDDO
1757	ELSE
1758	PRINT*, "iflag_con non-prevu", iflag_con
1759	CALL abort
1760	ENDIF
1761	CALL homogene(paprs, q_seri, d_q_con, u_seri,v_seri,
1762	. d_u_con, d_v_con)
1763	DO k = 1, klev
1764	DO i = 1, klon
1765	t_seri(i,k) = t_seri(i,k) + d_t_con(i,k)
1766	q_seri(i,k) = q_seri(i,k) + d_q_con(i,k)
1767	u_seri(i,k) = u_seri(i,k) + d_u_con(i,k)
1768	v_seri(i,k) = v_seri(i,k) + d_v_con(i,k)
1769	ENDDO
1770	ENDDO
1771	IF (check) THEN
1772	za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
1773	PRINT*, "aprescon=", za
1774	zx_t = 0.0
1775	za = 0.0
1776	DO i = 1, klon
1777	za = za + paire(i)/FLOAT(klon)
1778	zx_t = zx_t + (rain_con(i)+snow_con(i))*paire(i)/FLOAT(klon)
1779	ENDDO
1780	zx_t = zx_t/za*dtime
1781	PRINT*, "Precip=", zx_t
1782	ENDIF
1783	IF (zx_ajustq) THEN
1784	DO i = 1, klon
1785	z_apres(i) = 0.0
1786	ENDDO
1787	DO k = 1, klev
1788	DO i = 1, klon
1789	z_apres(i) = z_apres(i) + (q_seri(i,k)+ql_seri(i,k))
1790	. *(paprs(i,k)-paprs(i,k+1))/RG
1791	ENDDO
1792	ENDDO
1793	DO i = 1, klon
1794	z_factor(i) = (z_avant(i)-(rain_con(i)+snow_con(i))*dtime)
1795	. /z_apres(i)
1796	ENDDO
1797	DO k = 1, klev
1798	DO i = 1, klon
1799	IF (z_factor(i).GT.(1.0+1.0E-08) .OR.
1800	. z_factor(i).LT.(1.0-1.0E-08)) THEN
1801	q_seri(i,k) = q_seri(i,k) * z_factor(i)
1802	ENDIF
1803	ENDDO
1804	ENDDO
1805	ENDIF
1806	zx_ajustq=.FALSE.
1807	c
1808	IF (nqmax.GT.2) THEN !--melange convectif de traceurs
1809	c
1810	IF (iflag_con .NE. 2 .AND. iflag_con .LT. 3 ) THEN
1811	PRINT*, 'Pour l instant, seul conflx fonctionne ',
1812	$ 'avec traceurs', iflag_con
1813	PRINT*,' Mettre iflag_con',
1814	$ ' = 2 ou 4 dans run.def et repasser'
1815	CALL abort
1816	ENDIF
1817	c
1818	ENDIF !--nqmax.GT.2
1819	c
1820	c Appeler l'ajustement sec
1821	c
1822	CALL ajsec(paprs, pplay, t_seri, q_seri, d_t_ajs, d_q_ajs)
1823	DO k = 1, klev
1824	DO i = 1, klon
1825	t_seri(i,k) = t_seri(i,k) + d_t_ajs(i,k)
1826	q_seri(i,k) = q_seri(i,k) + d_q_ajs(i,k)
1827	ENDDO
1828	ENDDO
1829
1830	c RATQS
1831	call calcratqs (
1832	I paprs,pplay,q_seri,d_t_con,d_t_ajs
1833	O ,ratqs,zpt_conv)
1834	c
1835	c Appeler le processus de condensation a grande echelle
1836	c et le processus de precipitation
1837	c
1838	CALL fisrtilp_tr(dtime,paprs,pplay,
1839	. t_seri, q_seri,ratqs,
1840	. d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq,
1841	. rain_lsc, snow_lsc,
1842	. pfrac_impa, pfrac_nucl, pfrac_1nucl,
1843	. frac_impa, frac_nucl,
1844	. prfl, psfl)
1845	DO k = 1, klev
1846	DO i = 1, klon
1847	t_seri(i,k) = t_seri(i,k) + d_t_lsc(i,k)
1848	q_seri(i,k) = q_seri(i,k) + d_q_lsc(i,k)
1849	ql_seri(i,k) = ql_seri(i,k) + d_ql_lsc(i,k)
1850	cldfra(i,k) = rneb(i,k)
1851	IF (.NOT.new_oliq) cldliq(i,k) = ql_seri(i,k)
1852	ENDDO
1853	ENDDO
1854	IF (check) THEN
1855	za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
1856	PRINT*, "apresilp=", za
1857	zx_t = 0.0
1858	za = 0.0
1859	DO i = 1, klon
1860	za = za + paire(i)/FLOAT(klon)
1861	zx_t = zx_t + (rain_lsc(i)+snow_lsc(i))*paire(i)/FLOAT(klon)
1862	ENDDO
1863	zx_t = zx_t/za*dtime
1864	PRINT*, "Precip=", zx_t
1865	ENDIF
1866	c
1867	c Nuages diagnostiques:
1868	c
1869	IF (iflag_con.EQ.2) THEN ! seulement pour Tiedtke
1870	CALL diagcld1(paprs,pplay,
1871	. rain_con,snow_con,ibas_con,itop_con,
1872	. diafra,dialiq)
1873	DO k = 1, klev
1874	DO i = 1, klon
1875	IF (diafra(i,k).GT.cldfra(i,k)) THEN
1876	cldliq(i,k) = dialiq(i,k)
1877	cldfra(i,k) = diafra(i,k)
1878	ENDIF
1879	ENDDO
1880	ENDDO
1881	ENDIF
1882	c
1883	c Nuages stratus artificiels:
1884	c
1885	IF (ok_stratus) THEN
1886	CALL diagcld2(paprs,pplay,t_seri,q_seri, diafra,dialiq)
1887	DO k = 1, klev
1888	DO i = 1, klon
1889	IF (diafra(i,k).GT.cldfra(i,k)) THEN
1890	cldliq(i,k) = dialiq(i,k)
1891	cldfra(i,k) = diafra(i,k)
1892	ENDIF
1893	ENDDO
1894	ENDDO
1895	ENDIF
1896	c
1897	c Precipitation totale
1898	c
1899	DO i = 1, klon
1900	rain_fall(i) = rain_con(i) + rain_lsc(i)
1901	snow_fall(i) = snow_con(i) + snow_lsc(i)
1902	ENDDO
1903	c
1904	c Calculer l'humidite relative pour diagnostique
1905	c
1906	DO k = 1, klev
1907	DO i = 1, klon
1908	zx_t = t_seri(i,k)
1909	IF (thermcep) THEN
1910	zdelta = MAX(0.,SIGN(1.,rtt-zx_t))
1911	zx_qs = r2es * FOEEW(zx_t,zdelta)/pplay(i,k)
1912	zx_qs = MIN(0.5,zx_qs)
1913	zcor = 1./(1.-retv*zx_qs)
1914	zx_qs = zx_qs*zcor
1915	ELSE
1916	IF (zx_t.LT.t_coup) THEN
1917	zx_qs = qsats(zx_t)/pplay(i,k)
1918	ELSE
1919	zx_qs = qsatl(zx_t)/pplay(i,k)
1920	ENDIF
1921	ENDIF
1922	zx_rh(i,k) = q_seri(i,k)/zx_qs
1923	ENDDO
1924	ENDDO
1925	c
1926	c Calculer les parametres optiques des nuages et quelques
1927	c parametres pour diagnostiques:
1928	c
1929	CALL nuage (paprs, pplay,
1930	. t_seri, cldliq, cldfra, cldtau, cldemi,
1931	. cldh, cldl, cldm, cldt, cldq)
1932	c
1933	c Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1934	c
1935	IF (MOD(itaprad,radpas).EQ.0) THEN
1936	CALL orbite(FLOAT(julien),zlongi,dist)
1937	IF (cycle_diurne) THEN
1938	zdtime=dtime*FLOAT(radpas) ! pas de temps du rayonnement (s)
1939	CALL zenang(zlongi,gmtime,zdtime,rlat,rlon,rmu0,fract)
1940	c CALL zenith(zlongi,gmtime,rlat,rlon,rmu0,fract) !va disparaitre
1941	CALL alboc_cd(rmu0,alb_eau)
1942	ELSE
1943	CALL angle(zlongi,rlat,fract,rmu0)
1944	CALL alboc(FLOAT(julien),rlat,alb_eau)
1945	ENDIF
1946	CALL albsno(veget,agesno,alb_neig)
1947	DO i = 1, klon
1948	zx_alb_oce = alb_eau(i)
1949	IF (pctsrf(i,is_oce).GT.epsfra .AND. ftsol(i,is_oce).LT.271.35)
1950	. zx_alb_oce = 0.6 ! pour slab_ocean
1951	zfra = MAX(0.0,MIN(1.0,fsnow(i,is_lic)/(fsnow(i,is_lic)+10.0)))
1952	zx_alb_lic = alb_neig(i)zfra + 0.6(1.0-zfra)
1953	zfra = MAX(0.0,MIN(1.0,fsnow(i,is_ter)/(fsnow(i,is_ter)+10.0)))
1954	zx_alb_ter = alb_neig(i)zfra + lmt_alb(i)(1.0-zfra)
1955	zfra = MAX(0.0,MIN(1.0,fsnow(i,is_sic)/(fsnow(i,is_sic)+10.0)))
1956	zx_alb_sic = alb_neig(i)zfra + 0.6(1.0-zfra)
1957	albsol(i) = zx_alb_oce * pctsrf(i,is_oce)
1958	. + zx_alb_lic * pctsrf(i,is_lic)
1959	. + zx_alb_ter * pctsrf(i,is_ter)
1960	. + zx_alb_sic * pctsrf(i,is_sic)
1961	ENDDO
1962	CALL radlwsw ! nouveau rayonnement (compatible Arpege-IFS)
1963	e (dist, rmu0, fract, co2_ppm, solaire,
1964	e paprs, pplay,zxtsol,albsol, t_seri,q_seri,wo,
1965	e cldfra, cldemi, cldtau,
1966	s heat,heat0,cool,cool0,radsol,albpla,
1967	s topsw,toplw,solsw,sollw,
1968	s topsw0,toplw0,solsw0,sollw0)
1969	itaprad = 0
1970	ENDIF
1971	itaprad = itaprad + 1
1972	c
1973	c Ajouter la tendance des rayonnements (tous les pas)
1974	c
1975	DO k = 1, klev
1976	DO i = 1, klon
1977	t_seri(i,k) = t_seri(i,k)
1978	. + (heat(i,k)-cool(i,k)) * dtime/86400.
1979	ENDDO
1980	ENDDO
1981	c
1982	c Calculer l'hydrologie de la surface
1983	c
1984	CALL hydrol(dtime,pctsrf,rain_fall, snow_fall, evap,
1985	. agesno, ftsol,fqsol,fsnow, ruis)
1986	c
1987	DO i = 1, klon
1988	zxqsol(i) = 0.0
1989	zxsnow(i) = 0.0
1990	ENDDO
1991	DO nsrf = 1, nbsrf
1992	DO i = 1, klon
1993	zxqsol(i) = zxqsol(i) + fqsol(i,nsrf)*pctsrf(i,nsrf)
1994	zxsnow(i) = zxsnow(i) + fsnow(i,nsrf)*pctsrf(i,nsrf)
1995	ENDDO
1996	ENDDO
1997	c
1998	c Si une sous-fraction n'existe pas, elle prend la valeur moyenne
1999	c
2000	DO nsrf = 1, nbsrf
2001	DO i = 1, klon
2002	IF (pctsrf(i,nsrf).LT.epsfra) THEN
2003	fqsol(i,nsrf) = zxqsol(i)
2004	fsnow(i,nsrf) = zxsnow(i)
2005	ENDIF
2006	ENDDO
2007	ENDDO
2008	c
2009	c Calculer le bilan du sol et la derive de temperature (couplage)
2010	c
2011	DO i = 1, klon
2012	bils(i) = radsol(i) - sens(i) - evap(i)*RLVTT
2013	ENDDO
2014	IF (ok_ocean) THEN
2015	DO i = 1, klon
2016	cthermiq = cyang
2017	IF (ftsol(i,is_oce).LT. 271.35) cthermiq = cbing
2018	IF (pctsrf(i,is_oce).GT.epsfra) deltat(i) = deltat(i) +
2019	. (bils(i)-lmt_bils(i))/cthermiq * dtime
2020	IF (deltat(i).GT.15.0 ) deltat(i) = 15.0
2021	IF (deltat(i).LT.-15.0) deltat(i) = -15.0
2022	ENDDO
2023	ENDIF
2024	c
2025	cmoddeblott(jan95)
2026	c Appeler le programme de parametrisation de l'orographie
2027	c a l'echelle sous-maille:
2028	c
2029	IF (ok_orodr) THEN
2030	c
2031	c selection des points pour lesquels le shema est actif:
2032	igwd=0
2033	DO i=1,klon
2034	itest(i)=0
2035	c IF ((zstd(i).gt.10.0)) THEN
2036	IF (((zpic(i)-zmea(i)).GT.100.).AND.(zstd(i).GT.10.0)) THEN
2037	itest(i)=1
2038	igwd=igwd+1
2039	idx(igwd)=i
2040	ENDIF
2041	ENDDO
2042	igwdim=MAX(1,igwd)
2043	c
2044	CALL drag_noro(klon,klev,dtime,paprs,pplay,
2045	e zmea,zstd, zsig, zgam, zthe,zpic,zval,
2046	e igwd,igwdim,idx,itest,
2047	e t_seri, u_seri, v_seri,
2048	s zulow, zvlow, zustr, zvstr,
2049	s d_t_oro, d_u_oro, d_v_oro)
2050	c
2051	c ajout des tendances
2052	DO k = 1, klev
2053	DO i = 1, klon
2054	t_seri(i,k) = t_seri(i,k) + d_t_oro(i,k)
2055	u_seri(i,k) = u_seri(i,k) + d_u_oro(i,k)
2056	v_seri(i,k) = v_seri(i,k) + d_v_oro(i,k)
2057	ENDDO
2058	ENDDO
2059	c
2060	ENDIF ! fin de test sur ok_orodr
2061	c
2062	IF (ok_orolf) THEN
2063	c
2064	c selection des points pour lesquels le shema est actif:
2065	igwd=0
2066	DO i=1,klon
2067	itest(i)=0
2068	IF ((zpic(i)-zmea(i)).GT.100.) THEN
2069	itest(i)=1
2070	igwd=igwd+1
2071	idx(igwd)=i
2072	ENDIF
2073	ENDDO
2074	igwdim=MAX(1,igwd)
2075	c
2076	CALL lift_noro(klon,klev,dtime,paprs,pplay,
2077	e rlat,zmea,zstd, zsig, zgam, zthe,zpic,zval,
2078	e igwd,igwdim,idx,itest,
2079	e t_seri, u_seri, v_seri,
2080	s zulow, zvlow, zustr, zvstr,
2081	s d_t_lif, d_u_lif, d_v_lif)
2082	c
2083	c ajout des tendances
2084	DO k = 1, klev
2085	DO i = 1, klon
2086	t_seri(i,k) = t_seri(i,k) + d_t_lif(i,k)
2087	u_seri(i,k) = u_seri(i,k) + d_u_lif(i,k)
2088	v_seri(i,k) = v_seri(i,k) + d_v_lif(i,k)
2089	ENDDO
2090	ENDDO
2091	c
2092	ENDIF ! fin de test sur ok_orolf
2093	c
2094	cAA
2095	cAA Installation de l'interface online-offline pour traceurs
2096	cAA
2097	c====================================================================
2098	c Calcul des tendances traceurs
2099	c====================================================================
2100	C Pascale : il faut quand meme apeller phytrac car il gere les sorties
2101	cKE43 des traceurs => il faut donc mettre des flags a .false.
2102	IF (iflag_con.GE.3) THEN
2103	c on ajoute les tendances calculees par KE43
2104	DO iq=1, nqmax-2 ! Sandrine a -3 ???
2105	DO k = 1, nlev
2106	DO i = 1, klon
2107	tr_seri(i,k,iq) = tr_seri(i,k,iq) + d_tr(i,k,iq)
2108	ENDDO
2109	ENDDO
2110	WRITE(iqn,'(i2.2)') iq
2111	CALL minmaxqfi(tr_seri(1,1,iq),0.,1.e33,'couche lim iq='//iqn)
2112	ENDDO
2113	CMAF modif pour garder info du nombre de traceurs auxquels
2114	C la physique s'applique
2115	ELSE
2116	CMAF modif pour garder info du nombre de traceurs auxquels
2117	C la physique s'applique
2118	C
2119	call phytrac (rnpb,
2120	I debut,lafin,
2121	I nqmax-2,
2122	I nlon,nlev,dtime,
2123	I t,paprs,pplay,
2124	I pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
2125	I ycoefh,yu1,yv1,ftsol,pctsrf,rlat,
2126	I frac_impa, frac_nucl,
2127	I rlon,presnivs,paire,pphis,
2128	O tr_seri)
2129	ENDIF
2130
2131	IF (offline) THEN
2132
2133	call phystokenc (
2134	I nlon,nlev,pdtphys,rlon,rlat,
2135	I t,pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
2136	I ycoefh,yu1,yv1,ftsol,pctsrf,
2137	I frac_impa, frac_nucl,
2138	I pphis,paire,dtime,itap)
2139
2140
2141	ENDIF
2142
2143	c
2144	c Calculer le transport de l'eau et de l'energie (diagnostique)
2145	c
2146	CALL transp (paprs,zxtsol,
2147	e t_seri, q_seri, u_seri, v_seri, zphi,
2148	s ve, vq, ue, uq)
2149	c
2150	c Accumuler les variables a stocker dans les fichiers histoire:
2151	c
2152	IF (ok_oasis) THEN ! couplage oasis
2153	DO i = 1, klon
2154	oas_sols(i) = oas_sols(i) + solsw(i) / FLOAT(nexca)
2155	oas_nsol(i) = oas_nsol(i) + (bils(i)-solsw(i))/ FLOAT(nexca)
2156	oas_rain(i) = oas_rain(i) + rain_fall(i) / FLOAT(nexca)
2157	oas_snow(i) = oas_snow(i) + snow_fall(i) / FLOAT(nexca)
2158	oas_evap(i) = oas_evap(i) + evap(i) / FLOAT(nexca)
2159	oas_tsol(i) = oas_tsol(i) + zxtsol(i) / FLOAT(nexca)
2160	oas_fder(i) = oas_fder(i) + fder(i) / FLOAT(nexca)
2161	oas_albe(i) = oas_albe(i) + albsol(i) / FLOAT(nexca)
2162	oas_taux(i) = oas_taux(i) + fluxu(i,1) / FLOAT(nexca)
2163	oas_tauy(i) = oas_tauy(i) + fluxv(i,1) / FLOAT(nexca)
2164	oas_ruis(i) = oas_ruis(i) + ruis(i) /FLOAT(nexca)/dtime
2165	ENDDO
2166	ENDIF
2167	c
2168	c
2169	IF (ok_journe) THEN
2170	c
2171	ndex2d = 0
2172	ndex3d = 0
2173	c
2174	c Champs 2D:
2175	c
2176	i = NINT(zout/zsto)
2177	CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
2178	CALL histwrite(nid_day,"phis",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2179	c
2180	i = NINT(zout/zsto)
2181	CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
2182	CALL histwrite(nid_day,"aire",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2183	C
2184	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxtsol,zx_tmp_2d)
2185	CALL histwrite(nid_day,"tsol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2186	c
2187	DO i = 1, klon
2188	zx_tmp_fi2d(i) = paprs(i,1)
2189	ENDDO
2190	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2191	CALL histwrite(nid_day,"psol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2192	c
2193	DO i = 1, klon
2194	zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
2195	ENDDO
2196	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2197	CALL histwrite(nid_day,"rain",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2198	c
2199	CALL gr_fi_ecrit(1, klon,iim,jjmp1, snow_fall,zx_tmp_2d)
2200	CALL histwrite(nid_day,"snow",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2201	c
2202	CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
2203	CALL histwrite(nid_day,"evap",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2204	c
2205	CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw,zx_tmp_2d)
2206	CALL histwrite(nid_day,"tops",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2207	c
2208	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
2209	CALL histwrite(nid_day,"topl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2210	c
2211	CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
2212	CALL histwrite(nid_day,"sols",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2213	c
2214	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
2215	CALL histwrite(nid_day,"soll",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2216	c
2217	CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
2218	CALL histwrite(nid_day,"bils",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2219	c
2220	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
2221	CALL histwrite(nid_day,"sens",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2222	c
2223	CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
2224	CALL histwrite(nid_day,"fder",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2225	c
2226	CALL gr_fi_ecrit(1, klon,iim,jjmp1, ruis,zx_tmp_2d)
2227	CALL histwrite(nid_day,"ruis",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2228	c
2229	DO i = 1, klon
2230	zx_tmp_fi2d(i) = fluxu(i,1)
2231	ENDDO
2232	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2233	CALL histwrite(nid_day,"frtu",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2234	c
2235	DO i = 1, klon
2236	zx_tmp_fi2d(i) = fluxv(i,1)
2237	ENDDO
2238	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2239	CALL histwrite(nid_day,"frtv",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2240	c
2241	DO i = 1, klon
2242	zx_tmp_fi2d(i) = pctsrf(i,is_sic)
2243	ENDDO
2244	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2245	CALL histwrite(nid_day,"sicf",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2246	c
2247	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldl,zx_tmp_2d)
2248	CALL histwrite(nid_day,"cldl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2249	c
2250	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldm,zx_tmp_2d)
2251	CALL histwrite(nid_day,"cldm",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2252	c
2253	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldh,zx_tmp_2d)
2254	CALL histwrite(nid_day,"cldh",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2255	c
2256	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldt,zx_tmp_2d)
2257	CALL histwrite(nid_day,"cldt",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2258	c
2259	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldq,zx_tmp_2d)
2260	CALL histwrite(nid_day,"cldq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2261	c
2262	c Champs 3D:
2263	c
2264	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
2265	CALL histwrite(nid_day,"temp",itap,zx_tmp_3d,
2266	. iimjjmp1klev,ndex3d)
2267	c
2268	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,ivap), zx_tmp_3d)
2269	CALL histwrite(nid_day,"ovap",itap,zx_tmp_3d,
2270	. iimjjmp1klev,ndex3d)
2271	c
2272	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
2273	CALL histwrite(nid_day,"geop",itap,zx_tmp_3d,
2274	. iimjjmp1klev,ndex3d)
2275	c
2276	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
2277	CALL histwrite(nid_day,"vitu",itap,zx_tmp_3d,
2278	. iimjjmp1klev,ndex3d)
2279	c
2280	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
2281	CALL histwrite(nid_day,"vitv",itap,zx_tmp_3d,
2282	. iimjjmp1klev,ndex3d)
2283	c
2284	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, omega, zx_tmp_3d)
2285	CALL histwrite(nid_day,"vitw",itap,zx_tmp_3d,
2286	. iimjjmp1klev,ndex3d)
2287	c
2288	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
2289	CALL histwrite(nid_day,"pres",itap,zx_tmp_3d,
2290	. iimjjmp1klev,ndex3d)
2291	c
2292	if (ok_sync) then
2293	call histsync(nid_day)
2294	endif
2295	ENDIF
2296	C
2297	IF (ok_mensuel) THEN
2298	c
2299	ndex2d = 0
2300	ndex3d = 0
2301	c
2302	c Champs 2D:
2303	c
2304	i = NINT(zout/zsto)
2305	CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
2306	CALL histwrite(nid_mth,"phis",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2307	C
2308	i = NINT(zout/zsto)
2309	CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
2310	CALL histwrite(nid_mth,"aire",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2311
2312	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxtsol,zx_tmp_2d)
2313	CALL histwrite(nid_mth,"tsol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2314	c
2315	DO i = 1, klon
2316	zx_tmp_fi2d(i) = paprs(i,1)
2317	ENDDO
2318	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2319	CALL histwrite(nid_mth,"psol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2320	c
2321	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxqsol,zx_tmp_2d)
2322	CALL histwrite(nid_mth,"qsol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2323	c
2324	DO i = 1, klon
2325	zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
2326	ENDDO
2327	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2328	CALL histwrite(nid_mth,"rain",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2329	c
2330	DO i = 1, klon
2331	zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
2332	ENDDO
2333	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2334	CALL histwrite(nid_mth,"plul",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2335	c
2336	DO i = 1, klon
2337	zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
2338	ENDDO
2339	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2340	CALL histwrite(nid_mth,"pluc",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2341	c
2342	CALL gr_fi_ecrit(1, klon,iim,jjmp1, snow_fall,zx_tmp_2d)
2343	CALL histwrite(nid_mth,"snow",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2344	c
2345	CALL gr_fi_ecrit(1, klon,iim,jjmp1, agesno,zx_tmp_2d)
2346	CALL histwrite(nid_mth,"ages",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2347	c
2348	CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
2349	CALL histwrite(nid_mth,"evap",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2350	c
2351	CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw,zx_tmp_2d)
2352	CALL histwrite(nid_mth,"tops",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2353	c
2354	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
2355	CALL histwrite(nid_mth,"topl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2356	c
2357	CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
2358	CALL histwrite(nid_mth,"sols",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2359	c
2360	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
2361	CALL histwrite(nid_mth,"soll",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2362	c
2363	CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw0,zx_tmp_2d)
2364	CALL histwrite(nid_mth,"tops0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2365	c
2366	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw0,zx_tmp_2d)
2367	CALL histwrite(nid_mth,"topl0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2368	c
2369	CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw0,zx_tmp_2d)
2370	CALL histwrite(nid_mth,"sols0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2371	c
2372	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw0,zx_tmp_2d)
2373	CALL histwrite(nid_mth,"soll0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2374	c
2375	CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
2376	CALL histwrite(nid_mth,"bils",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2377	c
2378	CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
2379	CALL histwrite(nid_mth,"sens",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2380	c
2381	CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
2382	CALL histwrite(nid_mth,"fder",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2383	c
2384	CALL gr_fi_ecrit(1, klon,iim,jjmp1, ruis,zx_tmp_2d)
2385	CALL histwrite(nid_mth,"ruis",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2386	c
2387	DO i = 1, klon
2388	zx_tmp_fi2d(i) = fluxu(i,1)
2389	ENDDO
2390	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2391	CALL histwrite(nid_mth,"frtu",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2392	c
2393	DO i = 1, klon
2394	zx_tmp_fi2d(i) = fluxv(i,1)
2395	ENDDO
2396	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2397	CALL histwrite(nid_mth,"frtv",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2398	c
2399	DO i = 1, klon
2400	zx_tmp_fi2d(i) = pctsrf(i,is_sic)
2401	ENDDO
2402	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2403	CALL histwrite(nid_mth,"sicf",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2404	c
2405	CALL gr_fi_ecrit(1, klon,iim,jjmp1, albsol,zx_tmp_2d)
2406	CALL histwrite(nid_mth,"albs",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2407	c
2408	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragm,zx_tmp_2d)
2409	CALL histwrite(nid_mth,"cdrm",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2410	c
2411	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragh,zx_tmp_2d)
2412	CALL histwrite(nid_mth,"cdrh",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2413	c
2414	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldl,zx_tmp_2d)
2415	CALL histwrite(nid_mth,"cldl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2416	c
2417	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldm,zx_tmp_2d)
2418	CALL histwrite(nid_mth,"cldm",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2419	c
2420	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldh,zx_tmp_2d)
2421	CALL histwrite(nid_mth,"cldh",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2422	c
2423	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldt,zx_tmp_2d)
2424	CALL histwrite(nid_mth,"cldt",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2425	c
2426	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldq,zx_tmp_2d)
2427	CALL histwrite(nid_mth,"cldq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2428	c
2429	CALL gr_fi_ecrit(1, klon,iim,jjmp1, ue,zx_tmp_2d)
2430	CALL histwrite(nid_mth,"ue",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2431	c
2432	CALL gr_fi_ecrit(1, klon,iim,jjmp1, ve,zx_tmp_2d)
2433	CALL histwrite(nid_mth,"ve",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2434	c
2435	CALL gr_fi_ecrit(1, klon,iim,jjmp1, uq,zx_tmp_2d)
2436	CALL histwrite(nid_mth,"uq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2437	c
2438	CALL gr_fi_ecrit(1, klon,iim,jjmp1, vq,zx_tmp_2d)
2439	CALL histwrite(nid_mth,"vq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2440	cKE43
2441	IF (iflag_con .GE. 3) THEN ! sb
2442	c
2443	CALL gr_fi_ecrit(1, klon,iim,jjmp1, cape,zx_tmp_2d)
2444	CALL histwrite(nid_mth,"cape",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2445	c
2446	CALL gr_fi_ecrit(1, klon,iim,jjmp1,pbase,zx_tmp_2d)
2447	CALL histwrite(nid_mth,"pbase",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2448	c
2449	CALL gr_fi_ecrit(1, klon,iim,jjmp1,ema_pct,zx_tmp_2d)
2450	CALL histwrite(nid_mth,"ptop",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2451	c
2452	CALL gr_fi_ecrit(1, klon,iim,jjmp1,ema_cbmf,zx_tmp_2d)
2453	CALL histwrite(nid_mth,"fbase",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2454	c
2455	c
2456	ENDIF
2457	c34EK
2458	c
2459	c Champs 3D:
2460	C
2461	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
2462	CALL histwrite(nid_mth,"temp",itap,zx_tmp_3d,
2463	. iimjjmp1klev,ndex3d)
2464	c
2465	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,ivap), zx_tmp_3d)
2466	CALL histwrite(nid_mth,"ovap",itap,zx_tmp_3d,
2467	. iimjjmp1klev,ndex3d)
2468	c
2469	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
2470	CALL histwrite(nid_mth,"geop",itap,zx_tmp_3d,
2471	. iimjjmp1klev,ndex3d)
2472	c
2473	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
2474	CALL histwrite(nid_mth,"vitu",itap,zx_tmp_3d,
2475	. iimjjmp1klev,ndex3d)
2476	c
2477	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
2478	CALL histwrite(nid_mth,"vitv",itap,zx_tmp_3d,
2479	. iimjjmp1klev,ndex3d)
2480	c
2481	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, omega, zx_tmp_3d)
2482	CALL histwrite(nid_mth,"vitw",itap,zx_tmp_3d,
2483	. iimjjmp1klev,ndex3d)
2484	c
2485	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
2486	CALL histwrite(nid_mth,"pres",itap,zx_tmp_3d,
2487	. iimjjmp1klev,ndex3d)
2488	c
2489	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldfra, zx_tmp_3d)
2490	CALL histwrite(nid_mth,"rneb",itap,zx_tmp_3d,
2491	. iimjjmp1klev,ndex3d)
2492	c
2493	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zx_rh, zx_tmp_3d)
2494	CALL histwrite(nid_mth,"rhum",itap,zx_tmp_3d,
2495	. iimjjmp1klev,ndex3d)
2496	c
2497	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldliq, zx_tmp_3d)
2498	CALL histwrite(nid_mth,"oliq",itap,zx_tmp_3d,
2499	. iimjjmp1klev,ndex3d)
2500	c
2501	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_dyn, zx_tmp_3d)
2502	CALL histwrite(nid_mth,"dtdyn",itap,zx_tmp_3d,
2503	. iimjjmp1klev,ndex3d)
2504	c
2505	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_dyn, zx_tmp_3d)
2506	CALL histwrite(nid_mth,"dqdyn",itap,zx_tmp_3d,
2507	. iimjjmp1klev,ndex3d)
2508	c
2509	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_con, zx_tmp_3d)
2510	CALL histwrite(nid_mth,"dtcon",itap,zx_tmp_3d,
2511	. iimjjmp1klev,ndex3d)
2512	c
2513	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_con, zx_tmp_3d)
2514	CALL histwrite(nid_mth,"dqcon",itap,zx_tmp_3d,
2515	. iimjjmp1klev,ndex3d)
2516	c
2517	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_lsc, zx_tmp_3d)
2518	CALL histwrite(nid_mth,"dtlsc",itap,zx_tmp_3d,
2519	. iimjjmp1klev,ndex3d)
2520	c
2521	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_lsc, zx_tmp_3d)
2522	CALL histwrite(nid_mth,"dqlsc",itap,zx_tmp_3d,
2523	. iimjjmp1klev,ndex3d)
2524	c
2525	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_vdf, zx_tmp_3d)
2526	CALL histwrite(nid_mth,"dtvdf",itap,zx_tmp_3d,
2527	. iimjjmp1klev,ndex3d)
2528	c
2529	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_vdf, zx_tmp_3d)
2530	CALL histwrite(nid_mth,"dqvdf",itap,zx_tmp_3d,
2531	. iimjjmp1klev,ndex3d)
2532	c
2533	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_eva, zx_tmp_3d)
2534	CALL histwrite(nid_mth,"dteva",itap,zx_tmp_3d,
2535	. iimjjmp1klev,ndex3d)
2536	c
2537	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_eva, zx_tmp_3d)
2538	CALL histwrite(nid_mth,"dqeva",itap,zx_tmp_3d,
2539	. iimjjmp1klev,ndex3d)
2540	c
2541	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zpt_conv, zx_tmp_3d)
2542	CALL histwrite(nid_mth,"ptconv",itap,zx_tmp_3d,
2543	. iim(jjmp1)klev,ndex3d)
2544	c
2545	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, ratqs, zx_tmp_3d)
2546	CALL histwrite(nid_mth,"ratqs",itap,zx_tmp_3d,
2547	. iim(jjmp1)klev,ndex3d)
2548	c
2549	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_ajs, zx_tmp_3d)
2550	CALL histwrite(nid_mth,"dtajs",itap,zx_tmp_3d,
2551	. iimjjmp1klev,ndex3d)
2552	c
2553	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_ajs, zx_tmp_3d)
2554	CALL histwrite(nid_mth,"dqajs",itap,zx_tmp_3d,
2555	. iimjjmp1klev,ndex3d)
2556	c
2557	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, heat, zx_tmp_3d)
2558	CALL histwrite(nid_mth,"dtswr",itap,zx_tmp_3d,
2559	. iimjjmp1klev,ndex3d)
2560	c
2561	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, heat0, zx_tmp_3d)
2562	CALL histwrite(nid_mth,"dtsw0",itap,zx_tmp_3d,
2563	. iimjjmp1klev,ndex3d)
2564	c
2565	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cool, zx_tmp_3d)
2566	CALL histwrite(nid_mth,"dtlwr",itap,zx_tmp_3d,
2567	. iimjjmp1klev,ndex3d)
2568	c
2569	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cool0, zx_tmp_3d)
2570	CALL histwrite(nid_mth,"dtlw0",itap,zx_tmp_3d,
2571	. iimjjmp1klev,ndex3d)
2572	c
2573	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_vdf, zx_tmp_3d)
2574	CALL histwrite(nid_mth,"duvdf",itap,zx_tmp_3d,
2575	. iimjjmp1klev,ndex3d)
2576	c
2577	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_vdf, zx_tmp_3d)
2578	CALL histwrite(nid_mth,"dvvdf",itap,zx_tmp_3d,
2579	. iimjjmp1klev,ndex3d)
2580	c
2581	IF (ok_orodr) THEN
2582	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_oro, zx_tmp_3d)
2583	CALL histwrite(nid_mth,"duoro",itap,zx_tmp_3d,
2584	. iimjjmp1klev,ndex3d)
2585	c
2586	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_oro, zx_tmp_3d)
2587	CALL histwrite(nid_mth,"dvoro",itap,zx_tmp_3d,
2588	. iimjjmp1klev,ndex3d)
2589	c
2590	ENDIF
2591	C
2592	IF (ok_orolf) THEN
2593	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_lif, zx_tmp_3d)
2594	CALL histwrite(nid_mth,"dulif",itap,zx_tmp_3d,
2595	. iimjjmp1klev,ndex3d)
2596	c
2597	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_lif, zx_tmp_3d)
2598	CALL histwrite(nid_mth,"dvlif",itap,zx_tmp_3d,
2599	. iimjjmp1klev,ndex3d)
2600	ENDIF
2601	C
2602	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, wo, zx_tmp_3d)
2603	CALL histwrite(nid_mth,"ozone",itap,zx_tmp_3d,
2604	. iimjjmp1klev,ndex3d)
2605	c
2606	IF (nqmax.GE.3) THEN
2607	DO iq=1,nqmax-2
2608	IF (iq.LE.99) THEN
2609	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,iq+2), zx_tmp_3d)
2610	WRITE(str2,'(i2.2)') iq
2611	CALL histwrite(nid_mth,"trac"//str2,itap,zx_tmp_3d,
2612	. iimjjmp1klev,ndex3d)
2613	ELSE
2614	PRINT*, "Trop de traceurs"
2615	CALL abort
2616	ENDIF
2617	ENDDO
2618	ENDIF
2619	cKE43
2620	IF (iflag_con.GE.3) THEN ! (sb)
2621	c
2622	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, upwd, zx_tmp_3d)
2623	CALL histwrite(nid_mth,"upwd",itap,zx_tmp_3d,
2624	. iimjjmp1klev,ndex3d)
2625	c
2626	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, dnwd, zx_tmp_3d)
2627	CALL histwrite(nid_mth,"dnwd",itap,zx_tmp_3d,
2628	. iimjjmp1klev,ndex3d)
2629	c
2630	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, dnwd0, zx_tmp_3d)
2631	CALL histwrite(nid_mth,"dnwd0",itap,zx_tmp_3d,
2632	. iimjjmp1klev,ndex3d)
2633	c
2634	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, Ma, zx_tmp_3d)
2635	CALL histwrite(nid_mth,"Ma",itap,zx_tmp_3d,
2636	. iimjjmp1klev,ndex3d)
2637	c
2638	c
2639	ENDIF
2640	c34EK
2641	c
2642	if (ok_sync) then
2643	call histsync(nid_mth)
2644	endif
2645	ENDIF
2646	c
2647	IF (ok_instan) THEN
2648	c
2649	ndex2d = 0
2650	ndex3d = 0
2651	c
2652	c Champs 2D:
2653	c
2654	i = NINT(zout/zsto)
2655	CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
2656	CALL histwrite(nid_ins,"phis",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2657	c
2658	i = NINT(zout/zsto)
2659	CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
2660	CALL histwrite(nid_ins,"aire",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2661
2662	DO i = 1, klon
2663	zx_tmp_fi2d(i) = paprs(i,1)
2664	ENDDO
2665	CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2666	CALL histwrite(nid_ins,"psol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2667	c
2668	CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
2669	CALL histwrite(nid_ins,"topl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2670	c
2671	c Champs 3D:
2672	c
2673	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
2674	CALL histwrite(nid_ins,"temp",itap,zx_tmp_3d,
2675	. iimjjmp1klev,ndex3d)
2676	c
2677	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
2678	CALL histwrite(nid_ins,"vitu",itap,zx_tmp_3d,
2679	. iimjjmp1klev,ndex3d)
2680	c
2681	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
2682	CALL histwrite(nid_ins,"vitv",itap,zx_tmp_3d,
2683	. iimjjmp1klev,ndex3d)
2684	c
2685	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
2686	CALL histwrite(nid_ins,"geop",itap,zx_tmp_3d,
2687	. iimjjmp1klev,ndex3d)
2688	c
2689	CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
2690	CALL histwrite(nid_ins,"pres",itap,zx_tmp_3d,
2691	. iimjjmp1klev,ndex3d)
2692	c
2693	if (ok_sync) then
2694	call histsync(nid_ins)
2695	endif
2696	ENDIF
2697	c
2698	IF (ok_oasis .AND. mod(itap,nexca).EQ.0) THEN
2699	c
2700	c Je ne traite pas le ruissellement, pour l'instant (qui m'aidera ?)
2701	DO i = 1, klon
2702	oas_ruisoce(i) = 0.0
2703	oas_ruisriv(i) = 0.0
2704	ENDDO
2705	c
2706	ig = 1
2707	DO i = 1, iim
2708	z_sols(i,1) = oas_sols(ig)
2709	z_nsol(i,1) = oas_nsol(ig)
2710	z_rain(i,1) = oas_rain(ig)
2711	z_snow(i,1) = oas_snow(ig)
2712	z_evap(i,1) = oas_evap(ig)
2713	z_ruisoce(i,1) = oas_ruisoce(ig)
2714	z_ruisriv(i,1) = oas_ruisriv(ig)
2715	z_tsol(i,1) = oas_tsol(ig)
2716	z_fder(i,1) = oas_fder(ig)
2717	z_albe(i,1) = oas_albe(ig)
2718	z_taux(i,1) = oas_taux(ig)
2719	z_tauy(i,1) = oas_tauy(ig)
2720	ENDDO
2721	DO j = 2, jjm
2722	DO i = 1, iim
2723	ig = ig + 1
2724	z_sols(i,j) = oas_sols(ig)
2725	z_nsol(i,j) = oas_nsol(ig)
2726	z_rain(i,j) = oas_rain(ig)
2727	z_snow(i,j) = oas_snow(ig)
2728	z_evap(i,j) = oas_evap(ig)
2729	z_ruisoce(i,j) = oas_ruisoce(ig)
2730	z_ruisriv(i,j) = oas_ruisriv(ig)
2731	z_tsol(i,j) = oas_tsol(ig)
2732	z_fder(i,j) = oas_fder(ig)
2733	z_albe(i,j) = oas_albe(ig)
2734	z_taux(i,j) = oas_taux(ig)
2735	z_tauy(i,j) = oas_tauy(ig)
2736	ENDDO
2737	ENDDO
2738	ig = ig + 1
2739	DO i = 1, iim
2740	z_sols(i,jjmp1) = oas_sols(ig)
2741	z_nsol(i,jjmp1) = oas_nsol(ig)
2742	z_rain(i,jjmp1) = oas_rain(ig)
2743	z_snow(i,jjmp1) = oas_snow(ig)
2744	z_evap(i,jjmp1) = oas_evap(ig)
2745	z_ruisoce(i,jjmp1) = oas_ruisoce(ig)
2746	z_ruisriv(i,jjmp1) = oas_ruisriv(ig)
2747	z_tsol(i,jjmp1) = oas_tsol(ig)
2748	z_fder(i,jjmp1) = oas_fder(ig)
2749	z_albe(i,jjmp1) = oas_albe(ig)
2750	z_taux(i,jjmp1) = oas_taux(ig)
2751	z_tauy(i,jjmp1) = oas_tauy(ig)
2752	ENDDO
2753	c
2754	c Passer les champs au coupleur:
2755	c
2756	CALL intocpl(itap,jjmp1*iim,
2757	. z_sols, z_nsol,
2758	. z_rain, z_snow, z_evap,
2759	. z_ruisoce, z_ruisriv,
2760	. z_tsol, z_fder, z_albe,
2761	. z_taux, z_tauy)
2762	DO i = 1, klon
2763	oas_sols(i) = 0.0
2764	oas_nsol(i) = 0.0
2765	oas_rain(i) = 0.0
2766	oas_snow(i) = 0.0
2767	oas_evap(i) = 0.0
2768	oas_ruis(i) = 0.0
2769	oas_tsol(i) = 0.0
2770	oas_fder(i) = 0.0
2771	oas_albe(i) = 0.0
2772	oas_taux(i) = 0.0
2773	oas_tauy(i) = 0.0
2774	ENDDO
2775	ENDIF
2776	c
2777	c Ecrire la bande regionale (binaire grads)
2778	IF (ok_region .AND. mod(itap,ecrit_reg).eq.0) THEN
2779	CALL ecriregs(84,zxtsol)
2780	CALL ecriregs(84,paprs(1,1))
2781	CALL ecriregs(84,topsw)
2782	CALL ecriregs(84,toplw)
2783	CALL ecriregs(84,solsw)
2784	CALL ecriregs(84,sollw)
2785	CALL ecriregs(84,rain_fall)
2786	CALL ecriregs(84,snow_fall)
2787	CALL ecriregs(84,evap)
2788	CALL ecriregs(84,sens)
2789	CALL ecriregs(84,bils)
2790	CALL ecriregs(84,pctsrf(1,is_sic))
2791	CALL ecriregs(84,fluxu(1,1))
2792	CALL ecriregs(84,fluxv(1,1))
2793	CALL ecriregs(84,ue)
2794	CALL ecriregs(84,ve)
2795	CALL ecriregs(84,uq)
2796	CALL ecriregs(84,vq)
2797	c
2798	CALL ecrirega(84,u_seri)
2799	CALL ecrirega(84,v_seri)
2800	CALL ecrirega(84,omega)
2801	CALL ecrirega(84,t_seri)
2802	CALL ecrirega(84,zphi)
2803	CALL ecrirega(84,q_seri)
2804	CALL ecrirega(84,cldfra)
2805	CALL ecrirega(84,cldliq)
2806	CALL ecrirega(84,pplay)
2807
2808
2809	cc CALL ecrirega(84,d_t_dyn)
2810	cc CALL ecrirega(84,d_q_dyn)
2811	cc CALL ecrirega(84,heat)
2812	cc CALL ecrirega(84,cool)
2813	cc CALL ecrirega(84,d_t_con)
2814	cc CALL ecrirega(84,d_q_con)
2815	cc CALL ecrirega(84,d_t_lsc)
2816	cc CALL ecrirega(84,d_q_lsc)
2817	ENDIF
2818	c
2819	c Convertir les incrementations en tendances
2820	c
2821	DO k = 1, klev
2822	DO i = 1, klon
2823	d_u(i,k) = ( u_seri(i,k) - u(i,k) ) / dtime
2824	d_v(i,k) = ( v_seri(i,k) - v(i,k) ) / dtime
2825	d_t(i,k) = ( t_seri(i,k)-t(i,k) ) / dtime
2826	d_qx(i,k,ivap) = ( q_seri(i,k) - qx(i,k,ivap) ) / dtime
2827	d_qx(i,k,iliq) = ( ql_seri(i,k) - qx(i,k,iliq) ) / dtime
2828	ENDDO
2829	ENDDO
2830	c
2831	IF (nqmax.GE.3) THEN
2832	DO iq = 3, nqmax
2833	DO k = 1, klev
2834	DO i = 1, klon
2835	d_qx(i,k,iq) = ( tr_seri(i,k,iq-2) - qx(i,k,iq) ) / dtime
2836	ENDDO
2837	ENDDO
2838	ENDDO
2839	ENDIF
2840	c
2841	c Sauvegarder les valeurs de t et q a la fin de la physique:
2842	c
2843	DO k = 1, klev
2844	DO i = 1, klon
2845	t_ancien(i,k) = t_seri(i,k)
2846	q_ancien(i,k) = q_seri(i,k)
2847	ENDDO
2848	ENDDO
2849	c
2850	c====================================================================
2851	c Si c'est la fin, il faut conserver l'etat de redemarrage
2852	c====================================================================
2853	c
2854	IF (lafin) THEN
2855	ccc IF (ok_oasis) CALL quitcpl
2856	CALL phyredem ("restartphy.nc",dtime,radpas,co2_ppm,solaire,
2857	. rlat,rlon,ftsol,ftsoil,deltat,fqsol,fsnow,
2858	. radsol,rugmer,agesno,
2859	. zmea,zstd,zsig,zgam,zthe,zpic,zval,rugoro,
2860	. t_ancien, q_ancien)
2861	ENDIF
2862
2863	RETURN
2864	END
2865	FUNCTION qcheck(klon,klev,paprs,q,ql,aire)
2866	IMPLICIT none
2867	c
2868	c Calculer et imprimer l'eau totale. A utiliser pour verifier
2869	c la conservation de l'eau
2870	c
2871	#include "YOMCST.h"
2872	INTEGER klon,klev
2873	REAL paprs(klon,klev+1), q(klon,klev), ql(klon,klev)
2874	REAL aire(klon)
2875	REAL qtotal, zx, qcheck
2876	INTEGER i, k
2877	c
2878	zx = 0.0
2879	DO i = 1, klon
2880	zx = zx + aire(i)
2881	ENDDO
2882	qtotal = 0.0
2883	DO k = 1, klev
2884	DO i = 1, klon
2885	qtotal = qtotal + (q(i,k)+ql(i,k)) * aire(i)
2886	. *(paprs(i,k)-paprs(i,k+1))/RG
2887	ENDDO
2888	ENDDO
2889	c
2890	qcheck = qtotal/zx
2891	c
2892	RETURN
2893	END
2894	SUBROUTINE gr_fi_ecrit(nfield,nlon,iim,jjmp1,fi,ecrit)
2895	IMPLICIT none
2896	c
2897	c Tranformer une variable de la grille physique a
2898	c la grille d'ecriture
2899	c
2900	INTEGER nfield,nlon,iim,jjmp1, jjm
2901	REAL fi(nlon,nfield), ecrit(iim*jjmp1,nfield)
2902	c
2903	INTEGER i, j, n, ig
2904	c
2905	jjm = jjmp1 - 1
2906	DO n = 1, nfield
2907	DO i=1,iim
2908	ecrit(i,n) = fi(1,n)
2909	ecrit(i+jjm*iim,n) = fi(nlon,n)
2910	ENDDO
2911	DO ig = 1, nlon - 2
2912	ecrit(iim+ig,n) = fi(1+ig,n)
2913	ENDDO
2914	ENDDO
2915	RETURN
2916	END
2917

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

Download in other formats:

Original Format