source: LMDZ.3.3/branches/rel-LF/libf/phylmd/physiq.F @ 88

Last change on this file since 88 was 88, checked in by lmdzadmin, 25 years ago

Problemes avec les jjm+1

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