source: LMDZ.3.3/trunk/libf/phylmd/physiq.F @ 519

Last change on this file since 519 was 348, checked in by lmdz, 23 years ago

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