Context Navigation

physiq.F @ 15

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

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

Download in other formats:

Original Format