Context Navigation

physiq.F @ 53

Last change on this file since 53 was 53, checked in by lmdzadmin, 24 years ago
Correction sur le calendrier. LF
Property svn:eol-style set to `native` Property svn:executable set to ``* Property svn:keywords set to `Author Date Id Revision`
File size: 91.5 KB

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

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

Download in other formats:

Original Format