Context Navigation

physiq.F @ 197

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

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

Download in other formats:

Original Format