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

Last change on this file since 80 was 80, checked in by lmdzadmin, 24 years ago

Rajout de sorties de champs instantanes pour le debugging. LF

  • 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
554c======================================================================
555      xjour = rjourvrai
556c
557c Si c'est le debut, il faut initialiser plusieurs choses
558c          ********
559c
560       IF (debut) THEN
561c
562 
563         IF (ok_oasis) THEN
564            PRINT*, "Attentions! les parametres suivants sont fixes:"
565            PRINT *,'***********************************************'
566            PRINT*, "npas, nexca, itimestep=", npas, nexca, itimestep
567            PRINT*, "Changer-les manuellement s il le faut"
568            PRINT *,'***********************************************'
569            CALL inicma( npas, nexca, itimestep)
570         ENDIF
571c
572         IF (ok_ocean) THEN
573            PRINT*, '************************'
574            PRINT*, 'SLAB OCEAN est active, prenez precautions !'
575            PRINT*, '************************'
576         ENDIF
577c
578         DO k = 2, nvm          ! pas de vegetation
579            DO i = 1, klon
580               veget(i,k) = 0.0
581            ENDDO
582         ENDDO
583         DO i = 1, klon
584            veget(i,1) = 1.0    ! il n'y a que du desert
585         ENDDO
586         PRINT*, 'Pas de vegetation; desert partout'
587c
588c Initialiser les compteurs:
589c
590
591         itap    = 0
592         itaprad = 0
593c
594         CALL phyetat0 ("startphy.nc",dtime,co2_ppm,solaire,
595     .         rlat,rlon,ftsol,ftsoil,deltat,fqsol,fsnow,
596     .        radsol,rugmer,agesno,clesphy0,
597     .       zmea,zstd,zsig,zgam,zthe,zpic,zval,rugoro,tabcntr0,
598     .       t_ancien, q_ancien, ancien_ok )
599
600c
601         radpas = NINT( 86400./dtime/nbapp_rad)
602
603c
604         CALL printflag( tabcntr0,radpas,ok_ocean,ok_oasis ,ok_journe,
605     ,                    ok_instan, ok_region )
606c
607         IF (ABS(dtime-pdtphys).GT.0.001) THEN
608            PRINT*, 'Pas physique n est pas correcte',dtime,pdtphys
609            abort_message=' See above '
610            call abort_gcm(modname,abort_message,1)
611         ENDIF
612         IF (nlon .NE. klon) THEN
613            PRINT*, 'nlon et klon ne sont pas coherents', nlon, klon
614            abort_message=' See above '
615            call abort_gcm(modname,abort_message,1)
616         ENDIF
617         IF (nlev .NE. klev) THEN
618            PRINT*, 'nlev et klev ne sont pas coherents', nlev, klev
619            abort_message=' See above '
620            call abort_gcm(modname,abort_message,1)
621         ENDIF
622c
623         IF (dtime*FLOAT(radpas).GT.21600..AND.cycle_diurne) THEN
624           PRINT*, 'Nbre d appels au rayonnement insuffisant'
625           PRINT*, "Au minimum 4 appels par jour si cycle diurne"
626           abort_message=' See above '
627           call abort_gcm(modname,abort_message,1)
628         ENDIF
629         PRINT*, "Clef pour la convection, iflag_con=", iflag_con
630c
631         IF (ok_orodr) THEN
632         DO i=1,klon
633         rugoro(i) = MAX(1.0e-05, zstd(i)*zsig(i)/2.0)
634         ENDDO
635         CALL SUGWD(klon,klev,paprs,pplay)
636         DO i=1,klon
637         zuthe(i)=0.
638         zvthe(i)=0.
639         if(zstd(i).gt.10.)then
640           zuthe(i)=(1.-zgam(i))*cos(zthe(i))
641           zvthe(i)=(1.-zgam(i))*sin(zthe(i))
642         endif
643         ENDDO
644         ENDIF
645c
646         IF (soil_model) THEN
647            DO nsrf = 1, nbsrf
648            CALL soil(dtime, nsrf, fsnow(1,nsrf),
649     .                ftsol(1,nsrf), ftsoil(1,1,nsrf),
650     .                soilcap(1,nsrf), soilflux(1,nsrf))
651            ENDDO
652         ENDIF
653c
654         lmt_pas = NINT(86400./dtime * 1.0)   ! tous les jours
655         PRINT*,'La frequence de lecture surface est de ', lmt_pas
656c
657         ecrit_mth = NINT(86400./dtime *ecritphy)  ! tous les ecritphy jours
658         IF (ok_mensuel) THEN
659         PRINT*, 'La frequence de sortie mensuelle est de ', ecrit_mth
660         ENDIF
661         ecrit_day = NINT(86400./dtime *1.0)  ! tous les jours
662         IF (ok_journe) THEN
663         PRINT*, 'La frequence de sortie journaliere est de ',ecrit_day
664         ENDIF
665ccc         ecrit_ins = NINT(86400./dtime *0.5)  ! 2 fois par jour
666ccc         ecrit_ins = NINT(86400./dtime *0.25)  ! 4 fois par jour
667         ecrit_ins = NINT(86400./dtime/48.)  ! a chaque pas de temps
668         IF (ok_instan) THEN
669         PRINT*, 'La frequence de sortie instant. est de ', ecrit_ins
670         ENDIF
671         ecrit_reg = NINT(86400./dtime *0.25)  ! 4 fois par jour
672         IF (ok_region) THEN
673         PRINT*, 'La frequence de sortie region est de ', ecrit_reg
674         ENDIF
675c
676c
677      IF (ok_journe) THEN
678c
679         CALL ymds2ju(anne_ini, 1, 1, 0.0, zjulian)
680         zjulian = zjulian + day_ini
681c
682         CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
683         DO i = 1, iim
684            zx_lon(i,1) = rlon(i+1)
685            zx_lon(i,jjmp1) = rlon(i+1)
686         ENDDO
687         DO ll=1,klev
688            znivsig(ll)=float(ll)
689         ENDDO
690         CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
691         CALL histbeg("histday", iim,zx_lon, jjmp1,zx_lat,
692     .                 1,iim,1,jjmp1, 0, zjulian, dtime,
693     .                 nhori, nid_day)
694c         CALL histvert(nid_day, "presnivs", "Vertical levels", "mb",
695c     .                 klev, presnivs, nvert)
696         call histvert(nid_day, 'sig_s', 'Niveaux sigma','-',
697     .              klev, znivsig, nvert)
698c
699         zsto = dtime
700         zout = dtime * FLOAT(ecrit_day)
701c
702         CALL histdef(nid_day, "phis", "Surface geop. height", "-",
703     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
704     .                "once", zsto,zout)
705c
706         CALL histdef(nid_day, "aire", "Grid area", "-",
707     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
708     .                "once", zsto,zout)
709c
710c Champs 2D:
711c
712         CALL histdef(nid_day, "tsol", "Surface Temperature", "K",
713     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
714     .                "ave(X)", zsto,zout)
715c
716         CALL histdef(nid_day, "psol", "Surface Pressure", "Pa",
717     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
718     .                "ave(X)", zsto,zout)
719c
720         CALL histdef(nid_day, "rain", "Precipitation", "mm/day",
721     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
722     .                "ave(X)", zsto,zout)
723c
724         CALL histdef(nid_day, "snow", "Snow fall", "mm/day",
725     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
726     .                "ave(X)", zsto,zout)
727c
728         CALL histdef(nid_day, "evap", "Evaporation", "mm/day",
729     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
730     .                "ave(X)", zsto,zout)
731c
732         CALL histdef(nid_day, "tops", "Solar rad. at TOA", "W/m2",
733     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
734     .                "ave(X)", zsto,zout)
735c
736         CALL histdef(nid_day, "topl", "IR rad. at TOA", "W/m2",
737     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
738     .                "ave(X)", zsto,zout)
739c
740         CALL histdef(nid_day, "sols", "Solar rad. at surf.", "W/m2",
741     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
742     .                "ave(X)", zsto,zout)
743c
744         CALL histdef(nid_day, "soll", "IR rad. at surface", "W/m2",
745     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
746     .                "ave(X)", zsto,zout)
747c
748         CALL histdef(nid_day, "bils", "Surf. total heat flux", "W/m2",
749     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
750     .                "ave(X)", zsto,zout)
751c
752         CALL histdef(nid_day, "sens", "Sensible heat flux", "W/m2",
753     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
754     .                "ave(X)", zsto,zout)
755c
756         CALL histdef(nid_day, "fder", "Heat flux derivation", "W/m2",
757     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
758     .                "ave(X)", zsto,zout)
759c
760         CALL histdef(nid_day, "frtu", "Zonal wind stress", "Pa",
761     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
762     .                "ave(X)", zsto,zout)
763c
764         CALL histdef(nid_day, "frtv", "Meridional wind stress", "Pa",
765     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
766     .                "ave(X)", zsto,zout)
767c
768         CALL histdef(nid_day, "ruis", "Runoff", "mm/day",
769     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
770     .                "ave(X)", zsto,zout)
771c
772         CALL histdef(nid_day, "sicf", "Sea-ice fraction", "-",
773     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
774     .                "ave(X)", zsto,zout)
775c
776         CALL histdef(nid_day, "cldl", "Low-level cloudiness", "-",
777     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
778     .                "ave(X)", zsto,zout)
779c
780         CALL histdef(nid_day, "cldm", "Mid-level cloudiness", "-",
781     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
782     .                "ave(X)", zsto,zout)
783c
784         CALL histdef(nid_day, "cldh", "High-level cloudiness", "-",
785     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
786     .                "ave(X)", zsto,zout)
787c
788         CALL histdef(nid_day, "cldt", "Total cloudiness", "-",
789     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
790     .                "ave(X)", zsto,zout)
791c
792         CALL histdef(nid_day, "cldq", "Cloud liquid water path", "-",
793     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
794     .                "ave(X)", zsto,zout)
795c
796c Champs 3D:
797c
798         CALL histdef(nid_day, "temp", "Air temperature", "K",
799     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
800     .                "ave(X)", zsto,zout)
801c
802         CALL histdef(nid_day, "ovap", "Specific humidity", "Kg/Kg",
803     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
804     .                "ave(X)", zsto,zout)
805c
806         CALL histdef(nid_day, "geop", "Geopotential height", "m",
807     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
808     .                "ave(X)", zsto,zout)
809c
810         CALL histdef(nid_day, "vitu", "Zonal wind", "m/s",
811     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
812     .                "ave(X)", zsto,zout)
813c
814         CALL histdef(nid_day, "vitv", "Meridional wind", "m/s",
815     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
816     .                "ave(X)", zsto,zout)
817c
818         CALL histdef(nid_day, "vitw", "Vertical wind", "m/s",
819     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
820     .                "ave(X)", zsto,zout)
821c
822         CALL histdef(nid_day, "pres", "Air pressure", "Pa",
823     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
824     .                "ave(X)", zsto,zout)
825c
826         CALL histend(nid_day)
827c
828         ndex2d = 0
829         ndex3d = 0
830c
831      ENDIF ! fin de test sur ok_journe
832c
833      IF (ok_mensuel) THEN
834c
835         CALL ymds2ju(anne_ini, 1, 1, 0.0, zjulian)
836         zjulian = zjulian + day_ini
837c
838         CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
839         DO i = 1, iim
840            zx_lon(i,1) = rlon(i+1)
841            zx_lon(i,jjmp1) = rlon(i+1)
842         ENDDO
843         DO ll=1,klev
844            znivsig(ll)=float(ll)
845         ENDDO
846         CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
847         CALL histbeg("histmth", iim,zx_lon, jjmp1,zx_lat,
848     .                 1,iim,1,jjmp1, 0, zjulian, dtime,
849     .                 nhori, nid_mth)
850c         CALL histvert(nid_mth, "presnivs", "Vertical levels", "mb",
851c     .                 klev, presnivs, nvert)
852         call histvert(nid_mth, 'sig_s', 'Niveaux sigma','-',
853     .              klev, znivsig, nvert)
854c
855         zsto = dtime
856         zout = dtime * ecrit_mth
857c
858         CALL histdef(nid_mth, "phis", "Surface geop. height", "-",
859     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
860     .                "once",  zsto,zout)
861c
862         CALL histdef(nid_mth, "aire", "Grid area", "-",
863     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
864     .                "once",  zsto,zout)
865c
866c Champs 2D:
867c
868         CALL histdef(nid_mth, "tsol", "Surface Temperature", "K",
869     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
870     .                "ave(X)", zsto,zout)
871c
872         CALL histdef(nid_mth, "psol", "Surface Pressure", "Pa",
873     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
874     .                "ave(X)", zsto,zout)
875c
876         CALL histdef(nid_mth, "qsol", "Surface humidity", "mm",
877     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
878     .                "ave(X)", zsto,zout)
879c
880         CALL histdef(nid_mth, "rain", "Precipitation", "mm/day",
881     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
882     .                "ave(X)", zsto,zout)
883c
884         CALL histdef(nid_mth, "plul", "Large-scale Precip.", "mm/day",
885     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
886     .                "ave(X)", zsto,zout)
887c
888         CALL histdef(nid_mth, "pluc", "Convective Precip.", "mm/day",
889     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
890     .                "ave(X)", zsto,zout)
891c
892         CALL histdef(nid_mth, "snow", "Snow fall", "mm/day",
893     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
894     .                "ave(X)", zsto,zout)
895c
896         CALL histdef(nid_mth, "ages", "Snow age", "day",
897     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
898     .                "ave(X)", zsto,zout)
899c
900         CALL histdef(nid_mth, "evap", "Evaporation", "mm/day",
901     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
902     .                "ave(X)", zsto,zout)
903c
904         CALL histdef(nid_mth, "tops", "Solar rad. at TOA", "W/m2",
905     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
906     .                "ave(X)", zsto,zout)
907c
908         CALL histdef(nid_mth, "topl", "IR rad. at TOA", "W/m2",
909     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
910     .                "ave(X)", zsto,zout)
911c
912         CALL histdef(nid_mth, "sols", "Solar rad. at surf.", "W/m2",
913     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
914     .                "ave(X)", zsto,zout)
915c
916         CALL histdef(nid_mth, "soll", "IR rad. at surface", "W/m2",
917     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
918     .                "ave(X)", zsto,zout)
919c
920         CALL histdef(nid_mth, "tops0", "Solar rad. at TOA", "W/m2",
921     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
922     .                "ave(X)", zsto,zout)
923c
924         CALL histdef(nid_mth, "topl0", "IR rad. at TOA", "W/m2",
925     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
926     .                "ave(X)", zsto,zout)
927c
928         CALL histdef(nid_mth, "sols0", "Solar rad. at surf.", "W/m2",
929     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
930     .                "ave(X)", zsto,zout)
931c
932         CALL histdef(nid_mth, "soll0", "IR rad. at surface", "W/m2",
933     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
934     .                "ave(X)", zsto,zout)
935c
936         CALL histdef(nid_mth, "bils", "Surf. total heat flux", "W/m2",
937     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
938     .                "ave(X)", zsto,zout)
939c
940         CALL histdef(nid_mth, "sens", "Sensible heat flux", "W/m2",
941     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
942     .                "ave(X)", zsto,zout)
943c
944         CALL histdef(nid_mth, "fder", "Heat flux derivation", "W/m2",
945     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
946     .                "ave(X)", zsto,zout)
947c
948         CALL histdef(nid_mth, "frtu", "Zonal wind stress", "Pa",
949     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
950     .                "ave(X)", zsto,zout)
951c
952         CALL histdef(nid_mth, "frtv", "Meridional wind stress", "Pa",
953     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
954     .                "ave(X)", zsto,zout)
955c
956         CALL histdef(nid_mth, "ruis", "Runoff", "mm/day",
957     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
958     .                "ave(X)", zsto,zout)
959c
960         CALL histdef(nid_mth, "sicf", "Sea-ice fraction", "-",
961     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
962     .                "ave(X)", zsto,zout)
963c
964         CALL histdef(nid_mth, "albs", "Surface albedo", "-",
965     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
966     .                "ave(X)", zsto,zout)
967c
968         CALL histdef(nid_mth, "cdrm", "Momentum drag coef.", "-",
969     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
970     .                "ave(X)", zsto,zout)
971c
972         CALL histdef(nid_mth, "cdrh", "Heat drag coef.", "-",
973     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
974     .                "ave(X)", zsto,zout)
975c
976         CALL histdef(nid_mth, "cldl", "Low-level cloudiness", "-",
977     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
978     .                "ave(X)", zsto,zout)
979c
980         CALL histdef(nid_mth, "cldm", "Mid-level cloudiness", "-",
981     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
982     .                "ave(X)", zsto,zout)
983c
984         CALL histdef(nid_mth, "cldh", "High-level cloudiness", "-",
985     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
986     .                "ave(X)", zsto,zout)
987c
988         CALL histdef(nid_mth, "cldt", "Total cloudiness", "-",
989     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
990     .                "ave(X)", zsto,zout)
991c
992         CALL histdef(nid_mth, "cldq", "Cloud liquid water path", "-",
993     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
994     .                "ave(X)", zsto,zout)
995c
996         CALL histdef(nid_mth, "ue", "Zonal energy transport", "-",
997     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
998     .                "ave(X)", zsto,zout)
999c
1000         CALL histdef(nid_mth, "ve", "Merid energy transport", "-",
1001     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1002     .                "ave(X)", zsto,zout)
1003c
1004         CALL histdef(nid_mth, "uq", "Zonal humidity transport", "-",
1005     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1006     .                "ave(X)", zsto,zout)
1007c
1008         CALL histdef(nid_mth, "vq", "Merid humidity transport", "-",
1009     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1010     .                "ave(X)", zsto,zout)
1011c
1012c Champs 3D:
1013c
1014         CALL histdef(nid_mth, "temp", "Air temperature", "K",
1015     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1016     .                "ave(X)", zsto,zout)
1017c
1018         CALL histdef(nid_mth, "ovap", "Specific humidity", "Kg/Kg",
1019     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1020     .                "ave(X)", zsto,zout)
1021c
1022         CALL histdef(nid_mth, "geop", "Geopotential height", "m",
1023     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1024     .                "ave(X)", zsto,zout)
1025c
1026         CALL histdef(nid_mth, "vitu", "Zonal wind", "m/s",
1027     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1028     .                "ave(X)", zsto,zout)
1029c
1030         CALL histdef(nid_mth, "vitv", "Meridional wind", "m/s",
1031     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1032     .                "ave(X)", zsto,zout)
1033c
1034         CALL histdef(nid_mth, "vitw", "Vertical wind", "m/s",
1035     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1036     .                "ave(X)", zsto,zout)
1037c
1038         CALL histdef(nid_mth, "pres", "Air pressure", "Pa",
1039     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1040     .                "ave(X)", zsto,zout)
1041c
1042         CALL histdef(nid_mth, "rneb", "Cloud fraction", "-",
1043     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1044     .                "ave(X)", zsto,zout)
1045c
1046         CALL histdef(nid_mth, "rhum", "Relative humidity", "-",
1047     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1048     .                "ave(X)", zsto,zout)
1049c
1050         CALL histdef(nid_mth, "oliq", "Liquid water content", "kg/kg",
1051     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1052     .                "ave(X)", zsto,zout)
1053c
1054         CALL histdef(nid_mth, "dtdyn", "Dynamics dT", "K/s",
1055     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1056     .                "ave(X)", zsto,zout)
1057c
1058         CALL histdef(nid_mth, "dqdyn", "Dynamics dQ", "Kg/Kg/s",
1059     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1060     .                "ave(X)", zsto,zout)
1061c
1062         CALL histdef(nid_mth, "dtcon", "Convection dT", "K/s",
1063     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1064     .                "ave(X)", zsto,zout)
1065c
1066         CALL histdef(nid_mth, "dqcon", "Convection dQ", "Kg/Kg/s",
1067     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1068     .                "ave(X)", zsto,zout)
1069c
1070         CALL histdef(nid_mth, "dtlsc", "Condensation dT", "K/s",
1071     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1072     .                "ave(X)", zsto,zout)
1073c
1074         CALL histdef(nid_mth, "dqlsc", "Condensation dQ", "Kg/Kg/s",
1075     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1076     .                "ave(X)", zsto,zout)
1077c
1078         CALL histdef(nid_mth, "dtvdf", "Boundary-layer dT", "K/s",
1079     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1080     .                "ave(X)", zsto,zout)
1081c
1082         CALL histdef(nid_mth, "dqvdf", "Boundary-layer dQ", "Kg/Kg/s",
1083     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1084     .                "ave(X)", zsto,zout)
1085c
1086         CALL histdef(nid_mth, "dteva", "Reevaporation dT", "K/s",
1087     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1088     .                "ave(X)", zsto,zout)
1089c
1090         CALL histdef(nid_mth, "dqeva", "Reevaporation dQ", "Kg/Kg/s",
1091     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1092     .                "ave(X)", zsto,zout)
1093c
1094         CALL histdef(nid_mth, "dtajs", "Dry adjust. dT", "K/s",
1095     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1096     .                "ave(X)", zsto,zout)
1097
1098         CALL histdef(nid_mth, "dqajs", "Dry adjust. dQ", "Kg/Kg/s",
1099     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1100     .                "ave(X)", zsto,zout)
1101c
1102         CALL histdef(nid_mth, "dtswr", "SW radiation dT", "K/s",
1103     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1104     .                "ave(X)", zsto,zout)
1105c
1106         CALL histdef(nid_mth, "dtsw0", "SW radiation dT", "K/s",
1107     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1108     .                "ave(X)", zsto,zout)
1109c
1110         CALL histdef(nid_mth, "dtlwr", "LW radiation dT", "K/s",
1111     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1112     .                "ave(X)", zsto,zout)
1113c
1114         CALL histdef(nid_mth, "dtlw0", "LW radiation dT", "K/s",
1115     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1116     .                "ave(X)", zsto,zout)
1117c
1118         CALL histdef(nid_mth, "duvdf", "Boundary-layer dU", "m/s2",
1119     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1120     .                "ave(X)", zsto,zout)
1121c
1122         CALL histdef(nid_mth, "dvvdf", "Boundary-layer dV", "m/s2",
1123     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1124     .                "ave(X)", zsto,zout)
1125c
1126         IF (ok_orodr) THEN
1127         CALL histdef(nid_mth, "duoro", "Orography dU", "m/s2",
1128     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1129     .                "ave(X)", zsto,zout)
1130c
1131         CALL histdef(nid_mth, "dvoro", "Orography dV", "m/s2",
1132     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1133     .                "ave(X)", zsto,zout)
1134c
1135         ENDIF
1136C
1137         IF (ok_orolf) THEN
1138         CALL histdef(nid_mth, "dulif", "Orography dU", "m/s2",
1139     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1140     .                "ave(X)", zsto,zout)
1141c
1142         CALL histdef(nid_mth, "dvlif", "Orography dV", "m/s2",
1143     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1144     .                "ave(X)", zsto,zout)
1145         ENDIF
1146C
1147         CALL histdef(nid_mth, "ozone", "Ozone concentration", "-",
1148     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1149     .                "ave(X)", zsto,zout)
1150c
1151         if (nqmax.GE.3) THEN
1152         DO iq=1,nqmax-2
1153         IF (iq.LE.99) THEN
1154         WRITE(str2,'(i2.2)') iq
1155         CALL histdef(nid_mth, "trac"//str2, "Tracer No."//str2, "-",
1156     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1157     .                "ave(X)", zsto,zout)
1158         ELSE
1159         PRINT*, "Trop de traceurs"
1160         CALL abort
1161         ENDIF
1162         ENDDO
1163         ENDIF
1164c
1165         CALL histend(nid_mth)
1166c
1167         ndex2d = 0
1168         ndex3d = 0
1169c
1170      ENDIF ! fin de test sur ok_mensuel
1171c
1172c
1173      IF (ok_instan) THEN
1174c
1175         CALL ymds2ju(anne_ini, 1, 1, 0.0, zjulian)
1176         zjulian = zjulian + day_ini
1177c
1178         CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon,zx_lon)
1179         DO i = 1, iim
1180            zx_lon(i,1) = rlon(i+1)
1181            zx_lon(i,jjmp1) = rlon(i+1)
1182         ENDDO
1183         DO ll=1,klev
1184            znivsig(ll)=float(ll)
1185         ENDDO
1186         CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat,zx_lat)
1187         CALL histbeg("histins", iim,zx_lon, jjmp1,zx_lat,
1188     .                 1,iim,1,jjmp1, 0, zjulian, dtime,
1189     .                 nhori, nid_ins)
1190c         CALL histvert(nid_ins, "presnivs", "Vertical levels", "mb",
1191c     .                 klev, presnivs, nvert)
1192         call histvert(nid_ins, 'sig_s', 'Niveaux sigma','-',
1193     .              klev, znivsig, nvert)
1194c
1195c
1196         zsto = dtime * ecrit_ins
1197         zout = dtime * ecrit_ins
1198C
1199         CALL histdef(nid_ins, "phis", "Surface geop. height", "-",
1200     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1201     .                "once", zsto,zout)
1202c
1203         CALL histdef(nid_ins, "aire", "Grid area", "-",
1204     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1205     .                "once", zsto,zout)
1206c
1207c Champs 2D:
1208c
1209         CALL histdef(nid_ins, "psol", "Surface Pressure", "Pa",
1210     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1211     .                "inst(X)", zsto,zout)
1212c
1213         CALL histdef(nid_ins, "topl", "OLR", "W/m2",
1214     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1215     .                "inst(X)", zsto,zout)
1216c
1217         CALL histdef(nid_ins, "evap", "Evaporation", "mm/day",
1218     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1219     .                "inst(X)", zsto,zout)
1220c
1221         CALL histdef(nid_ins, "sols", "Solar rad. at surf.", "W/m2",
1222     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1223     .                "inst(X)", zsto,zout)
1224c
1225         CALL histdef(nid_ins, "soll", "IR rad. at surface", "W/m2",
1226     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1227     .                "inst(X)", zsto,zout)
1228c
1229         CALL histdef(nid_ins, "bils", "Surf. total heat flux", "W/m2",
1230     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1231     .                "inst(X)", zsto,zout)
1232c
1233         CALL histdef(nid_ins, "sens", "Sensible heat flux", "W/m2",
1234     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1235     .                "inst(X)", zsto,zout)
1236c
1237         CALL histdef(nid_ins, "fder", "Heat flux derivation", "W/m2",
1238     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1239     .                "inst(X)", zsto,zout)
1240c
1241      CALL histdef(nid_ins, "dtsvdfo", "Boundary-layer dTs(o)", "K/s",
1242     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1243     .                "inst(X)", zsto,zout)
1244c
1245      CALL histdef(nid_ins, "dtsvdft", "Boundary-layer dTs(t)", "K/s",
1246     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1247     .                "inst(X)", zsto,zout)
1248c
1249      CALL histdef(nid_ins, "dtsvdfg", "Boundary-layer dTs(g)", "K/s",
1250     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1251     .                "inst(X)", zsto,zout)
1252c
1253      CALL histdef(nid_ins, "dtsvdfi", "Boundary-layer dTs(g)", "K/s",
1254     .                iim,jjmp1,nhori, 1,1,1, -99, 32,
1255     .                "inst(X)", zsto,zout)
1256c
1257c Champs 3D:
1258c
1259         CALL histdef(nid_ins, "temp", "Temperature", "K",
1260     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1261     .                "inst(X)", zsto,zout)
1262c
1263         CALL histdef(nid_ins, "vitu", "Zonal wind", "m/s",
1264     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1265     .                "inst(X)", zsto,zout)
1266c
1267         CALL histdef(nid_ins, "vitv", "Merid wind", "m/s",
1268     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1269     .                "inst(X)", zsto,zout)
1270c
1271         CALL histdef(nid_ins, "geop", "Geopotential height", "m",
1272     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1273     .                "inst(X)", zsto,zout)
1274c
1275         CALL histdef(nid_ins, "pres", "Air pressure", "Pa",
1276     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1277     .                "inst(X)", zsto,zout)
1278c
1279         CALL histdef(nid_ins, "dtvdf", "Boundary-layer dT", "K/s",
1280     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1281     .                "inst(X)", zsto,zout)
1282c
1283         CALL histdef(nid_ins, "dqvdf", "Boundary-layer dQ", "Kg/Kg/s",
1284     .                iim,jjmp1,nhori, klev,1,klev,nvert, 32,
1285     .                "inst(X)", zsto,zout)
1286c
1287
1288         CALL histend(nid_ins)
1289c
1290         ndex2d = 0
1291         ndex3d = 0
1292c
1293      ENDIF
1294c
1295c
1296c
1297c Prescrire l'ozone dans l'atmosphere
1298c
1299c
1300cc         DO i = 1, klon
1301cc         DO k = 1, klev
1302cc            CALL o3cm (paprs(i,k)/100.,paprs(i,k+1)/100., wo(i,k),20)
1303cc         ENDDO
1304cc         ENDDO
1305c
1306         IF (ok_oasis) THEN
1307         DO i = 1, klon
1308           oas_sols(i) = 0.0
1309           oas_nsol(i) = 0.0
1310           oas_rain(i) = 0.0
1311           oas_snow(i) = 0.0
1312           oas_evap(i) = 0.0
1313           oas_ruis(i) = 0.0
1314           oas_tsol(i) = 0.0
1315           oas_fder(i) = 0.0
1316           oas_albe(i) = 0.0
1317           oas_taux(i) = 0.0
1318           oas_tauy(i) = 0.0
1319         ENDDO
1320         ENDIF
1321c
1322      ENDIF
1323c
1324c   ****************     Fin  de   IF ( debut  )   ***************
1325c
1326c
1327c Mettre a zero des variables de sortie (pour securite)
1328c
1329      DO i = 1, klon
1330         d_ps(i) = 0.0
1331      ENDDO
1332      DO k = 1, klev
1333      DO i = 1, klon
1334         d_t(i,k) = 0.0
1335         d_u(i,k) = 0.0
1336         d_v(i,k) = 0.0
1337      ENDDO
1338      ENDDO
1339      DO iq = 1, nqmax
1340      DO k = 1, klev
1341      DO i = 1, klon
1342         d_qx(i,k,iq) = 0.0
1343      ENDDO
1344      ENDDO
1345      ENDDO
1346c
1347c Ne pas affecter les valeurs entrees de u, v, h, et q
1348c
1349      DO k = 1, klev
1350      DO i = 1, klon
1351         t_seri(i,k)  = t(i,k)
1352         u_seri(i,k)  = u(i,k)
1353         v_seri(i,k)  = v(i,k)
1354         q_seri(i,k)  = qx(i,k,ivap)
1355         ql_seri(i,k) = qx(i,k,iliq)
1356      ENDDO
1357      ENDDO
1358      IF (nqmax.GE.3) THEN
1359      DO iq = 3, nqmax
1360      DO  k = 1, klev
1361      DO  i = 1, klon
1362         tr_seri(i,k,iq-2) = qx(i,k,iq)
1363      ENDDO
1364      ENDDO
1365      ENDDO
1366      ELSE
1367      DO k = 1, klev
1368      DO i = 1, klon
1369         tr_seri(i,k,1) = 0.0
1370      ENDDO
1371      ENDDO
1372      ENDIF
1373c
1374c Diagnostiquer la tendance dynamique
1375c
1376      IF (ancien_ok) THEN
1377         DO k = 1, klev
1378         DO i = 1, klon
1379            d_t_dyn(i,k) = (t_seri(i,k)-t_ancien(i,k))/dtime
1380            d_q_dyn(i,k) = (q_seri(i,k)-q_ancien(i,k))/dtime
1381         ENDDO
1382         ENDDO
1383      ELSE
1384         DO k = 1, klev
1385         DO i = 1, klon
1386            d_t_dyn(i,k) = 0.0
1387            d_q_dyn(i,k) = 0.0
1388         ENDDO
1389         ENDDO
1390         ancien_ok = .TRUE.
1391      ENDIF
1392c
1393c Ajouter le geopotentiel du sol:
1394c
1395      DO k = 1, klev
1396      DO i = 1, klon
1397         zphi(i,k) = pphi(i,k) + pphis(i)
1398      ENDDO
1399      ENDDO
1400c
1401c Verifier les temperatures
1402c
1403
1404      CALL hgardfou(t_seri,ftsol,'debutphy')
1405c
1406c Incrementer le compteur de la physique
1407c
1408      itap   = itap + 1
1409      julien = MOD(NINT(xjour),360)
1410c
1411c Mettre en action les conditions aux limites (albedo, sst, etc.).
1412c Prescrire l'ozone et calculer l'albedo sur l'ocean.
1413c
1414      IF (MOD(itap-1,lmt_pas) .EQ. 0) THEN
1415         idayvrai = NINT(xjour)
1416         PRINT *,' PHYS cond  julien ',julien,idayvrai
1417         CALL condsurf(julien,idayvrai, pctsrf ,
1418     .                  lmt_sst,lmt_alb,lmt_rug,lmt_bils  )
1419         CALL ozonecm( FLOAT(julien), rlat, paprs, wo)
1420      ENDIF
1421cccccccccc
1422      IF (ok_oasis .AND. MOD(itap-1,nexca).EQ.0) THEN
1423C
1424         CALL fromcpl(itap,jjmp1*iim,
1425     .        cpl_sst,cpl_sic,cpl_alb_sst,cpl_alb_sic)
1426         DO i = 1, iim-1 ! un seul point pour le pole nord
1427            cpl_sst(i,1) = cpl_sst(iim,1)
1428            cpl_sic(i,1) = cpl_sic(iim,1)
1429            cpl_alb_sst(i,1) = cpl_alb_sst(iim,1)
1430            cpl_alb_sic(i,1) = cpl_alb_sic(iim,1)
1431         ENDDO
1432         DO i = 2, iim ! un seul point pour le pole sud
1433            cpl_sst(i,jjmp1) = cpl_sst(1,jjmp1)
1434            cpl_sic(i,jjmp1) = cpl_sic(1,jjmp1)
1435            cpl_alb_sst(i,jjmp1) = cpl_alb_sst(1,jjmp1)
1436            cpl_alb_sic(i,jjmp1) = cpl_alb_sic(1,jjmp1)
1437         ENDDO
1438c
1439         ig = 1
1440         IF (pctsrf(ig,is_oce).GT.epsfra .OR.
1441     .       pctsrf(ig,is_sic).GT.epsfra) THEN
1442            pctsrf(ig,is_oce) = pctsrf(ig,is_oce)
1443     .                        - (cpl_sic(1,1)-pctsrf(ig,is_sic))
1444            pctsrf(ig,is_sic) = cpl_sic(1,1)
1445            lmt_sst(ig) = cpl_sst(1,1)
1446         ENDIF
1447         DO j = 2, jjm
1448         DO i = 1, iim
1449         ig = ig + 1
1450         IF (pctsrf(ig,is_oce).GT.epsfra .OR.
1451     .       pctsrf(ig,is_sic).GT.epsfra) THEN
1452           pctsrf(ig,is_oce) = pctsrf(ig,is_oce)
1453     .                       - (cpl_sic(i,j)-pctsrf(ig,is_sic))
1454           pctsrf(ig,is_sic) = cpl_sic(i,j)
1455           lmt_sst(ig) = cpl_sst(i,j)
1456         ENDIF
1457         ENDDO
1458         ENDDO
1459         ig = ig + 1
1460         IF (pctsrf(ig,is_oce).GT.epsfra .OR.
1461     .       pctsrf(ig,is_sic).GT.epsfra) THEN
1462            pctsrf(ig,is_oce) = pctsrf(ig,is_oce)
1463     .                        - (cpl_sic(1,jjmp1)-pctsrf(ig,is_sic))
1464            pctsrf(ig,is_sic) = cpl_sic(1,jjmp1)
1465            lmt_sst(ig) = cpl_sst(1,jjmp1)
1466         ENDIF
1467c
1468      ENDIF ! ok_oasis
1469cccccccccc
1470c
1471 
1472      IF (ok_ocean) THEN
1473         DO i = 1, klon
1474            ftsol(i,is_oce) = lmt_sst(i) + deltat(i)
1475         ENDDO
1476
1477      ELSE
1478         DO i = 1, klon
1479            ftsol(i,is_oce) = lmt_sst(i)
1480         ENDDO
1481
1482      ENDIF
1483c
1484c Re-evaporer l'eau liquide nuageuse
1485c
1486      DO k = 1, klev  ! re-evaporation de l'eau liquide nuageuse
1487      DO i = 1, klon
1488         zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
1489         zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
1490         zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
1491         zb = MAX(0.0,ql_seri(i,k))
1492         za = - MAX(0.0,ql_seri(i,k))
1493     .                  * (zlvdcp*(1.-zdelta)+zlsdcp*zdelta)
1494         t_seri(i,k) = t_seri(i,k) + za
1495         q_seri(i,k) = q_seri(i,k) + zb
1496         ql_seri(i,k) = 0.0
1497         d_t_eva(i,k) = za
1498         d_q_eva(i,k) = zb
1499      ENDDO
1500      ENDDO
1501c
1502c Appeler la diffusion verticale (programme de couche limite)
1503c
1504      DO i = 1, klon
1505         frugs(i,is_ter) = SQRT(lmt_rug(i)**2+rugoro(i)**2)
1506         frugs(i,is_lic) = rugoro(i)
1507         frugs(i,is_oce) = rugmer(i)
1508         frugs(i,is_sic) = 0.001
1509         zxrugs(i) = 0.0
1510      ENDDO
1511      DO nsrf = 1, nbsrf
1512      DO i = 1, klon
1513         frugs(i,nsrf) = MAX(frugs(i,nsrf),0.001)
1514      ENDDO
1515      ENDDO
1516      DO nsrf = 1, nbsrf
1517      DO i = 1, klon
1518         zxrugs(i) = zxrugs(i) + frugs(i,nsrf)*pctsrf(i,nsrf)
1519      ENDDO
1520      ENDDO
1521c
1522      CALL clmain(dtime,pctsrf,
1523     e            t_seri,q_seri,u_seri,v_seri,soil_model,
1524     e            ftsol,soilcap,soilflux,paprs,pplay,radsol,
1525     e            fsnow,fqsol,
1526     e            rlat, frugs,
1527     s            d_t_vdf,d_q_vdf,d_u_vdf,d_v_vdf,d_ts,
1528     s            fluxt,fluxq,fluxu,fluxv,cdragh,cdragm,rugmer,
1529     s            dsens, devap,
1530     s            ycoefh,yu1,yv1)
1531c
1532      DO i = 1, klon
1533         sens(i) = - fluxt(i,1) ! flux de chaleur sensible au sol
1534         evap(i) = - fluxq(i,1) ! flux d'evaporation au sol
1535         fder(i) = dsens(i) + devap(i)
1536      ENDDO
1537
1538      DO k = 1, klev
1539      DO i = 1, klon
1540         t_seri(i,k) = t_seri(i,k) + d_t_vdf(i,k)
1541         q_seri(i,k) = q_seri(i,k) + d_q_vdf(i,k)
1542         u_seri(i,k) = u_seri(i,k) + d_u_vdf(i,k)
1543         v_seri(i,k) = v_seri(i,k) + d_v_vdf(i,k)
1544      ENDDO
1545      ENDDO
1546c
1547c Incrementer la temperature du sol
1548c
1549      DO i = 1, klon
1550         zxtsol(i) = 0.0
1551      ENDDO
1552      DO nsrf = 1, nbsrf
1553      DO i = 1, klon
1554         ftsol(i,nsrf) = ftsol(i,nsrf) + d_ts(i,nsrf)
1555         zxtsol(i) = zxtsol(i) + ftsol(i,nsrf)*pctsrf(i,nsrf)
1556      ENDDO
1557      ENDDO
1558
1559c
1560c Si une sous-fraction n'existe pas, elle prend la temp. moyenne
1561c
1562      DO nsrf = 1, nbsrf
1563      DO i = 1, klon
1564         IF (pctsrf(i,nsrf).LT.epsfra) ftsol(i,nsrf) = zxtsol(i)
1565      ENDDO
1566      ENDDO
1567
1568c
1569c Appeler le modele du sol
1570c
1571      IF (soil_model) THEN
1572         DO nsrf = 1, nbsrf
1573         CALL soil(dtime, nsrf, fsnow(1,nsrf),
1574     .             ftsol(1,nsrf), ftsoil(1,1,nsrf),
1575     .             soilcap(1,nsrf), soilflux(1,nsrf))
1576         ENDDO
1577      ENDIF
1578c
1579c Calculer la derive du flux infrarouge
1580c
1581      DO nsrf = 1, nbsrf
1582      DO i = 1, klon
1583         fder(i) = fder(i) - 4.0*RSIGMA*zxtsol(i)**3 *
1584     .                       (ftsol(i,nsrf)-zxtsol(i))
1585     .                      *pctsrf(i,nsrf)
1586      ENDDO
1587      ENDDO
1588c
1589c Appeler la convection (au choix)
1590c
1591      DO k = 1, klev
1592      DO i = 1, klon
1593         conv_q(i,k) = d_q_dyn(i,k)
1594     .               + d_q_vdf(i,k)/dtime
1595         conv_t(i,k) = d_t_dyn(i,k)
1596     .               + d_t_vdf(i,k)/dtime
1597      ENDDO
1598      ENDDO
1599      IF (check) THEN
1600         za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
1601         PRINT*, "avantcon=", za
1602      ENDIF
1603      zx_ajustq = .FALSE.
1604      IF (iflag_con.EQ.2) zx_ajustq=.TRUE.
1605      IF (zx_ajustq) THEN
1606         DO i = 1, klon
1607            z_avant(i) = 0.0
1608         ENDDO
1609         DO k = 1, klev
1610         DO i = 1, klon
1611            z_avant(i) = z_avant(i) + (q_seri(i,k)+ql_seri(i,k))
1612     .                        *(paprs(i,k)-paprs(i,k+1))/RG
1613         ENDDO
1614         ENDDO
1615      ENDIF
1616      IF (iflag_con.EQ.1) THEN
1617          stop'reactiver le call conlmd dans physiq.F'
1618c     CALL conlmd (dtime, paprs, pplay, t_seri, q_seri, conv_q,
1619c    .             d_t_con, d_q_con,
1620c    .             rain_con, snow_con, ibas_con, itop_con)
1621      ELSE IF (iflag_con.EQ.2) THEN
1622      CALL conflx(dtime, paprs, pplay, t_seri, q_seri,
1623     e            conv_t, conv_q, fluxq(1,1), omega,
1624     s            d_t_con, d_q_con, rain_con, snow_con,
1625     s            pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
1626     s            kcbot, kctop, kdtop, pmflxr, pmflxs)
1627      DO i = 1, klon
1628         ibas_con(i) = klev+1 - kcbot(i)
1629         itop_con(i) = klev+1 - kctop(i)
1630      ENDDO
1631      ELSE IF (iflag_con.EQ.3) THEN
1632          stop'reactiver le call conlmd dans physiq.F'
1633c     CALL conccm (dtime,paprs,pplay,t_seri,q_seri,conv_q,
1634c    s             d_t_con, d_q_con,
1635c    s             rain_con, snow_con, ibas_con, itop_con)
1636      ELSE
1637      PRINT*, "iflag_con non-prevu", iflag_con
1638      CALL abort
1639      ENDIF
1640      CALL homogene(paprs, q_seri, d_q_con, u_seri,v_seri,
1641     .              d_u_con, d_v_con)
1642      DO k = 1, klev
1643      DO i = 1, klon
1644         t_seri(i,k) = t_seri(i,k) + d_t_con(i,k)
1645         q_seri(i,k) = q_seri(i,k) + d_q_con(i,k)
1646         u_seri(i,k) = u_seri(i,k) + d_u_con(i,k)
1647         v_seri(i,k) = v_seri(i,k) + d_v_con(i,k)
1648      ENDDO
1649      ENDDO
1650      IF (check) THEN
1651         za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
1652         PRINT*, "aprescon=", za
1653         zx_t = 0.0
1654         za = 0.0
1655         DO i = 1, klon
1656            za = za + paire(i)/FLOAT(klon)
1657            zx_t = zx_t + (rain_con(i)+snow_con(i))*paire(i)/FLOAT(klon)
1658        ENDDO
1659         zx_t = zx_t/za*dtime
1660         PRINT*, "Precip=", zx_t
1661      ENDIF
1662      IF (zx_ajustq) THEN
1663         DO i = 1, klon
1664            z_apres(i) = 0.0
1665         ENDDO
1666         DO k = 1, klev
1667         DO i = 1, klon
1668            z_apres(i) = z_apres(i) + (q_seri(i,k)+ql_seri(i,k))
1669     .                        *(paprs(i,k)-paprs(i,k+1))/RG
1670         ENDDO
1671         ENDDO
1672         DO i = 1, klon
1673         z_factor(i) = (z_avant(i)-(rain_con(i)+snow_con(i))*dtime)
1674     .                /z_apres(i)
1675         ENDDO
1676         DO k = 1, klev
1677         DO i = 1, klon
1678         IF (z_factor(i).GT.(1.0+1.0E-08) .OR.
1679     .       z_factor(i).LT.(1.0-1.0E-08)) THEN
1680               q_seri(i,k) = q_seri(i,k) * z_factor(i)
1681         ENDIF
1682         ENDDO
1683         ENDDO
1684      ENDIF
1685      zx_ajustq=.FALSE.
1686c
1687      IF (nqmax.GT.2) THEN !--melange convectif de traceurs
1688c
1689      IF (iflag_con.NE.2) THEN
1690         PRINT*, "Pour l instant, seul conflx fonctionne avec traceurs"
1691         PRINT*,' Mettre iflag_con = 2  dans  run.def et repasser  !'
1692         CALL abort
1693      ENDIF
1694c
1695      ENDIF !--nqmax.GT.2
1696c
1697c Appeler l'ajustement sec
1698c
1699      CALL ajsec(paprs, pplay, t_seri, q_seri, d_t_ajs, d_q_ajs)
1700      DO k = 1, klev
1701      DO i = 1, klon
1702         t_seri(i,k) = t_seri(i,k) + d_t_ajs(i,k)
1703         q_seri(i,k) = q_seri(i,k) + d_q_ajs(i,k)
1704      ENDDO
1705      ENDDO
1706c
1707c Appeler le processus de condensation a grande echelle
1708c et le processus de precipitation
1709c
1710      CALL fisrtilp_tr(dtime,paprs,pplay,
1711     .           t_seri, q_seri,
1712     .           d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq,
1713     .           rain_lsc, snow_lsc,
1714     .           pfrac_impa, pfrac_nucl, pfrac_1nucl,
1715     .           frac_impa, frac_nucl,
1716     .           prfl, psfl)
1717      DO k = 1, klev
1718      DO i = 1, klon
1719         t_seri(i,k) = t_seri(i,k) + d_t_lsc(i,k)
1720         q_seri(i,k) = q_seri(i,k) + d_q_lsc(i,k)
1721         ql_seri(i,k) = ql_seri(i,k) + d_ql_lsc(i,k)
1722         cldfra(i,k) = rneb(i,k)
1723         IF (.NOT.new_oliq) cldliq(i,k) = ql_seri(i,k)
1724      ENDDO
1725      ENDDO
1726      IF (check) THEN
1727         za = qcheck(klon,klev,paprs,q_seri,ql_seri,paire)
1728         PRINT*, "apresilp=", za
1729         zx_t = 0.0
1730         za = 0.0
1731         DO i = 1, klon
1732            za = za + paire(i)/FLOAT(klon)
1733            zx_t = zx_t + (rain_lsc(i)+snow_lsc(i))*paire(i)/FLOAT(klon)
1734        ENDDO
1735         zx_t = zx_t/za*dtime
1736         PRINT*, "Precip=", zx_t
1737      ENDIF
1738c
1739c Nuages diagnostiques:
1740c
1741      IF (iflag_con.EQ.2) THEN ! seulement pour Tiedtke
1742      CALL diagcld1(paprs,pplay,
1743     .             rain_con,snow_con,ibas_con,itop_con,
1744     .             diafra,dialiq)
1745      DO k = 1, klev
1746      DO i = 1, klon
1747      IF (diafra(i,k).GT.cldfra(i,k)) THEN
1748         cldliq(i,k) = dialiq(i,k)
1749         cldfra(i,k) = diafra(i,k)
1750      ENDIF
1751      ENDDO
1752      ENDDO
1753      ENDIF
1754c
1755c Nuages stratus artificiels:
1756c
1757      IF (ok_stratus) THEN
1758      CALL diagcld2(paprs,pplay,t_seri,q_seri, diafra,dialiq)
1759      DO k = 1, klev
1760      DO i = 1, klon
1761      IF (diafra(i,k).GT.cldfra(i,k)) THEN
1762         cldliq(i,k) = dialiq(i,k)
1763         cldfra(i,k) = diafra(i,k)
1764      ENDIF
1765      ENDDO
1766      ENDDO
1767      ENDIF
1768c
1769c Precipitation totale
1770c
1771      DO i = 1, klon
1772         rain_fall(i) = rain_con(i) + rain_lsc(i)
1773         snow_fall(i) = snow_con(i) + snow_lsc(i)
1774      ENDDO
1775c
1776c Calculer l'humidite relative pour diagnostique
1777c
1778      DO k = 1, klev
1779      DO i = 1, klon
1780         zx_t = t_seri(i,k)
1781         IF (thermcep) THEN
1782            zdelta = MAX(0.,SIGN(1.,rtt-zx_t))
1783            zx_qs  = r2es * FOEEW(zx_t,zdelta)/pplay(i,k)
1784            zx_qs  = MIN(0.5,zx_qs)
1785            zcor   = 1./(1.-retv*zx_qs)
1786            zx_qs  = zx_qs*zcor
1787         ELSE
1788           IF (zx_t.LT.t_coup) THEN
1789              zx_qs = qsats(zx_t)/pplay(i,k)
1790           ELSE
1791              zx_qs = qsatl(zx_t)/pplay(i,k)
1792           ENDIF
1793         ENDIF
1794         zx_rh(i,k) = q_seri(i,k)/zx_qs
1795      ENDDO
1796      ENDDO
1797c
1798c Calculer les parametres optiques des nuages et quelques
1799c parametres pour diagnostiques:
1800c
1801      CALL nuage (paprs, pplay,
1802     .            t_seri, cldliq, cldfra, cldtau, cldemi,
1803     .            cldh, cldl, cldm, cldt, cldq)
1804c
1805c Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1806c
1807      IF (MOD(itaprad,radpas).EQ.0) THEN
1808      CALL orbite(FLOAT(julien),zlongi,dist)
1809      IF (cycle_diurne) THEN
1810        zdtime=dtime*FLOAT(radpas) ! pas de temps du rayonnement (s)
1811        CALL zenang(zlongi,gmtime,zdtime,rlat,rlon,rmu0,fract)
1812c        CALL zenith(zlongi,gmtime,rlat,rlon,rmu0,fract) !va disparaitre
1813        CALL alboc_cd(rmu0,alb_eau)
1814      ELSE
1815        CALL angle(zlongi,rlat,fract,rmu0)
1816        CALL alboc(FLOAT(julien),rlat,alb_eau)
1817      ENDIF
1818      CALL albsno(veget,agesno,alb_neig)
1819      DO i = 1, klon
1820         zx_alb_oce = alb_eau(i)
1821         IF (pctsrf(i,is_oce).GT.epsfra .AND. ftsol(i,is_oce).LT.271.35)
1822     .   zx_alb_oce = 0.6 ! pour slab_ocean
1823         zfra = MAX(0.0,MIN(1.0,fsnow(i,is_lic)/(fsnow(i,is_lic)+10.0)))
1824         zx_alb_lic = alb_neig(i)*zfra + 0.6*(1.0-zfra)
1825         zfra = MAX(0.0,MIN(1.0,fsnow(i,is_ter)/(fsnow(i,is_ter)+10.0)))
1826         zx_alb_ter = alb_neig(i)*zfra + lmt_alb(i)*(1.0-zfra)
1827         zfra = MAX(0.0,MIN(1.0,fsnow(i,is_sic)/(fsnow(i,is_sic)+10.0)))
1828         zx_alb_sic = alb_neig(i)*zfra + 0.6*(1.0-zfra)
1829         albsol(i) = zx_alb_oce * pctsrf(i,is_oce)
1830     .             + zx_alb_lic * pctsrf(i,is_lic)
1831     .             + zx_alb_ter * pctsrf(i,is_ter)
1832     .             + zx_alb_sic * pctsrf(i,is_sic)
1833      ENDDO
1834      CALL radlwsw ! nouveau rayonnement (compatible Arpege-IFS)
1835     e            (dist, rmu0, fract, co2_ppm, solaire,
1836     e             paprs, pplay,zxtsol,albsol, t_seri,q_seri,wo,
1837     e             cldfra, cldemi, cldtau,
1838     s             heat,heat0,cool,cool0,radsol,albpla,
1839     s             topsw,toplw,solsw,sollw,
1840     s             topsw0,toplw0,solsw0,sollw0)
1841      itaprad = 0
1842      ENDIF
1843      itaprad = itaprad + 1
1844c
1845c Ajouter la tendance des rayonnements (tous les pas)
1846c
1847      DO k = 1, klev
1848      DO i = 1, klon
1849         t_seri(i,k) = t_seri(i,k)
1850     .               + (heat(i,k)-cool(i,k)) * dtime/86400.
1851      ENDDO
1852      ENDDO
1853c
1854c Calculer l'hydrologie de la surface
1855c
1856      CALL hydrol(dtime,pctsrf,rain_fall, snow_fall, evap,
1857     .            agesno, ftsol,fqsol,fsnow, ruis)
1858c
1859      DO i = 1, klon
1860         zxqsol(i) = 0.0
1861         zxsnow(i) = 0.0
1862      ENDDO
1863      DO nsrf = 1, nbsrf
1864      DO i = 1, klon
1865         zxqsol(i) = zxqsol(i) + fqsol(i,nsrf)*pctsrf(i,nsrf)
1866         zxsnow(i) = zxsnow(i) + fsnow(i,nsrf)*pctsrf(i,nsrf)
1867      ENDDO
1868      ENDDO
1869c
1870c Si une sous-fraction n'existe pas, elle prend la valeur moyenne
1871c
1872      DO nsrf = 1, nbsrf
1873      DO i = 1, klon
1874         IF (pctsrf(i,nsrf).LT.epsfra) THEN
1875            fqsol(i,nsrf) = zxqsol(i)
1876            fsnow(i,nsrf) = zxsnow(i)
1877         ENDIF
1878      ENDDO
1879      ENDDO
1880c
1881c Calculer le bilan du sol et la derive de temperature (couplage)
1882c
1883      DO i = 1, klon
1884         bils(i) = radsol(i) - sens(i) - evap(i)*RLVTT
1885      ENDDO
1886      IF (ok_ocean) THEN
1887         DO i = 1, klon
1888            cthermiq = cyang
1889            IF (ftsol(i,is_oce).LT. 271.35) cthermiq = cbing
1890            IF (pctsrf(i,is_oce).GT.epsfra) deltat(i) = deltat(i) +
1891     .                          (bils(i)-lmt_bils(i))/cthermiq * dtime
1892            IF (deltat(i).GT.15.0 ) deltat(i) = 15.0
1893            IF (deltat(i).LT.-15.0) deltat(i) = -15.0
1894         ENDDO
1895      ENDIF
1896c
1897cmoddeblott(jan95)
1898c Appeler le programme de parametrisation de l'orographie
1899c a l'echelle sous-maille:
1900c
1901      IF (ok_orodr) THEN
1902c
1903c  selection des points pour lesquels le shema est actif:
1904        igwd=0
1905        DO i=1,klon
1906        itest(i)=0
1907c        IF ((zstd(i).gt.10.0)) THEN
1908        IF (((zpic(i)-zmea(i)).GT.100.).AND.(zstd(i).GT.10.0)) THEN
1909          itest(i)=1
1910          igwd=igwd+1
1911          idx(igwd)=i
1912        ENDIF
1913        ENDDO
1914        igwdim=MAX(1,igwd)
1915c
1916        CALL drag_noro(klon,klev,dtime,paprs,pplay,
1917     e                   zmea,zstd, zsig, zgam, zthe,zpic,zval,
1918     e                   igwd,igwdim,idx,itest,
1919     e                   t_seri, u_seri, v_seri,
1920     s                   zulow, zvlow, zustr, zvstr,
1921     s                   d_t_oro, d_u_oro, d_v_oro)
1922c
1923c  ajout des tendances
1924        DO k = 1, klev
1925        DO i = 1, klon
1926           t_seri(i,k) = t_seri(i,k) + d_t_oro(i,k)
1927           u_seri(i,k) = u_seri(i,k) + d_u_oro(i,k)
1928           v_seri(i,k) = v_seri(i,k) + d_v_oro(i,k)
1929        ENDDO
1930        ENDDO
1931c
1932      ENDIF ! fin de test sur ok_orodr
1933c
1934      IF (ok_orolf) THEN
1935c
1936c  selection des points pour lesquels le shema est actif:
1937        igwd=0
1938        DO i=1,klon
1939        itest(i)=0
1940        IF ((zpic(i)-zmea(i)).GT.100.) THEN
1941          itest(i)=1
1942          igwd=igwd+1
1943          idx(igwd)=i
1944        ENDIF
1945        ENDDO
1946        igwdim=MAX(1,igwd)
1947c
1948        CALL lift_noro(klon,klev,dtime,paprs,pplay,
1949     e                   rlat,zmea,zstd, zsig, zgam, zthe,zpic,zval,
1950     e                   igwd,igwdim,idx,itest,
1951     e                   t_seri, u_seri, v_seri,
1952     s                   zulow, zvlow, zustr, zvstr,
1953     s                   d_t_lif, d_u_lif, d_v_lif)
1954c
1955c  ajout des tendances
1956        DO k = 1, klev
1957        DO i = 1, klon
1958           t_seri(i,k) = t_seri(i,k) + d_t_lif(i,k)
1959           u_seri(i,k) = u_seri(i,k) + d_u_lif(i,k)
1960           v_seri(i,k) = v_seri(i,k) + d_v_lif(i,k)
1961        ENDDO
1962        ENDDO
1963c
1964      ENDIF ! fin de test sur ok_orolf
1965c
1966cAA
1967cAA Installation de l'interface online-offline pour traceurs
1968cAA
1969c====================================================================
1970c   Calcul  des tendances traceurs
1971c====================================================================
1972CMAF modif pour garder info du nombre de traceurs auxquels
1973C la physique s'applique
1974C
1975      call phytrac (rnpb,
1976     I                   debut,
1977     I                   nqmax-2,
1978     I                   nlon,nlev,dtime,
1979     I                   t,paprs,pplay,
1980     I                   pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
1981     I                   ycoefh,yu1,yv1,ftsol,pctsrf,rlat,
1982     I                   frac_impa, frac_nucl,
1983     I                   rlon,presnivs,paire,pphis,
1984     O                   tr_seri)
1985
1986      IF (offline) THEN
1987
1988         call phystokenc (
1989     I                   nlon,nlev,pdtphys,rlon,rlat,
1990     I                   pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
1991     I                   ycoefh,yu1,yv1,ftsol,pctsrf,
1992     I                   frac_impa, frac_nucl,
1993     I                   pphis,paire,dtime,itap,
1994     O                   physid)
1995
1996      ENDIF
1997
1998c
1999c Calculer le transport de l'eau et de l'energie (diagnostique)
2000c
2001      CALL transp (paprs,zxtsol,
2002     e                   t_seri, q_seri, u_seri, v_seri, zphi,
2003     s                   ve, vq, ue, uq)
2004c
2005c Accumuler les variables a stocker dans les fichiers histoire:
2006c
2007      IF (ok_oasis) THEN ! couplage oasis
2008      DO i = 1, klon
2009        oas_sols(i) = oas_sols(i) + solsw(i)          / FLOAT(nexca)
2010        oas_nsol(i) = oas_nsol(i) + (bils(i)-solsw(i))/ FLOAT(nexca)
2011        oas_rain(i) = oas_rain(i) + rain_fall(i)      / FLOAT(nexca)
2012        oas_snow(i) = oas_snow(i) + snow_fall(i)      / FLOAT(nexca)
2013        oas_evap(i) = oas_evap(i) + evap(i)           / FLOAT(nexca)
2014        oas_tsol(i) = oas_tsol(i) + zxtsol(i)         / FLOAT(nexca)
2015        oas_fder(i) = oas_fder(i) + fder(i)           / FLOAT(nexca)
2016        oas_albe(i) = oas_albe(i) + albsol(i)         / FLOAT(nexca)
2017        oas_taux(i) = oas_taux(i) + fluxu(i,1)        / FLOAT(nexca)
2018        oas_tauy(i) = oas_tauy(i) + fluxv(i,1)        / FLOAT(nexca)
2019        oas_ruis(i) = oas_ruis(i) + ruis(i)       /FLOAT(nexca)/dtime
2020      ENDDO
2021      ENDIF
2022c
2023c
2024      IF (ok_journe) THEN
2025c
2026      ndex2d = 0
2027      ndex3d = 0
2028c
2029c Champs 2D:
2030c
2031         i = NINT(zout/zsto)
2032         CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
2033         CALL histwrite(nid_day,"phis",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2034c
2035         i = NINT(zout/zsto)
2036         CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
2037         CALL histwrite(nid_day,"aire",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2038C
2039      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxtsol,zx_tmp_2d)
2040      CALL histwrite(nid_day,"tsol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2041c
2042      DO i = 1, klon
2043         zx_tmp_fi2d(i) = paprs(i,1)
2044      ENDDO
2045      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2046      CALL histwrite(nid_day,"psol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2047c
2048      DO i = 1, klon
2049         zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
2050      ENDDO
2051      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2052      CALL histwrite(nid_day,"rain",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2053c
2054      CALL gr_fi_ecrit(1, klon,iim,jjmp1, snow_fall,zx_tmp_2d)
2055      CALL histwrite(nid_day,"snow",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2056c
2057      CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
2058      CALL histwrite(nid_day,"evap",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2059c
2060      CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw,zx_tmp_2d)
2061      CALL histwrite(nid_day,"tops",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2062c
2063      CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
2064      CALL histwrite(nid_day,"topl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2065c
2066      CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
2067      CALL histwrite(nid_day,"sols",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2068c
2069      CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
2070      CALL histwrite(nid_day,"soll",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2071c
2072      CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
2073      CALL histwrite(nid_day,"bils",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2074c
2075      CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
2076      CALL histwrite(nid_day,"sens",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2077c
2078      CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
2079      CALL histwrite(nid_day,"fder",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2080c
2081      CALL gr_fi_ecrit(1, klon,iim,jjmp1, ruis,zx_tmp_2d)
2082      CALL histwrite(nid_day,"ruis",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2083c
2084      DO i = 1, klon
2085         zx_tmp_fi2d(i) = fluxu(i,1)
2086      ENDDO
2087      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2088      CALL histwrite(nid_day,"frtu",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2089c
2090      DO i = 1, klon
2091         zx_tmp_fi2d(i) = fluxv(i,1)
2092      ENDDO
2093      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2094      CALL histwrite(nid_day,"frtv",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2095c
2096      DO i = 1, klon
2097         zx_tmp_fi2d(i) = pctsrf(i,is_sic)
2098      ENDDO
2099      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2100      CALL histwrite(nid_day,"sicf",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2101c
2102      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldl,zx_tmp_2d)
2103      CALL histwrite(nid_day,"cldl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2104c
2105      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldm,zx_tmp_2d)
2106      CALL histwrite(nid_day,"cldm",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2107c
2108      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldh,zx_tmp_2d)
2109      CALL histwrite(nid_day,"cldh",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2110c
2111      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldt,zx_tmp_2d)
2112      CALL histwrite(nid_day,"cldt",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2113c
2114      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldq,zx_tmp_2d)
2115      CALL histwrite(nid_day,"cldq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2116c
2117c Champs 3D:
2118c
2119      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
2120      CALL histwrite(nid_day,"temp",itap,zx_tmp_3d,
2121     .                                   iim*jjmp1*klev,ndex3d)
2122c
2123      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,ivap), zx_tmp_3d)
2124      CALL histwrite(nid_day,"ovap",itap,zx_tmp_3d,
2125     .                                   iim*jjmp1*klev,ndex3d)
2126c
2127      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
2128      CALL histwrite(nid_day,"geop",itap,zx_tmp_3d,
2129     .                                   iim*jjmp1*klev,ndex3d)
2130c
2131      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
2132      CALL histwrite(nid_day,"vitu",itap,zx_tmp_3d,
2133     .                                   iim*jjmp1*klev,ndex3d)
2134c
2135      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
2136      CALL histwrite(nid_day,"vitv",itap,zx_tmp_3d,
2137     .                                   iim*jjmp1*klev,ndex3d)
2138c
2139      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, omega, zx_tmp_3d)
2140      CALL histwrite(nid_day,"vitw",itap,zx_tmp_3d,
2141     .                                   iim*jjmp1*klev,ndex3d)
2142c
2143      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
2144      CALL histwrite(nid_day,"pres",itap,zx_tmp_3d,
2145     .                                   iim*jjmp1*klev,ndex3d)
2146c
2147      if (ok_sync) then
2148        call histsync(nid_day)
2149      endif
2150      ENDIF
2151C
2152      IF (ok_mensuel) THEN
2153c
2154      ndex2d = 0
2155      ndex3d = 0
2156c
2157c Champs 2D:
2158c
2159         i = NINT(zout/zsto)
2160         CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
2161         CALL histwrite(nid_mth,"phis",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2162C
2163         i = NINT(zout/zsto)
2164         CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
2165         CALL histwrite(nid_mth,"aire",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2166
2167      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxtsol,zx_tmp_2d)
2168      CALL histwrite(nid_mth,"tsol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2169c
2170      DO i = 1, klon
2171         zx_tmp_fi2d(i) = paprs(i,1)
2172      ENDDO
2173      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2174      CALL histwrite(nid_mth,"psol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2175c
2176      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zxqsol,zx_tmp_2d)
2177      CALL histwrite(nid_mth,"qsol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2178c
2179      DO i = 1, klon
2180         zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
2181      ENDDO
2182      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2183      CALL histwrite(nid_mth,"rain",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2184c
2185      DO i = 1, klon
2186         zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
2187      ENDDO
2188      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2189      CALL histwrite(nid_mth,"plul",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2190c
2191      DO i = 1, klon
2192         zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
2193      ENDDO
2194      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2195      CALL histwrite(nid_mth,"pluc",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2196c
2197      CALL gr_fi_ecrit(1, klon,iim,jjmp1, snow_fall,zx_tmp_2d)
2198      CALL histwrite(nid_mth,"snow",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2199c
2200      CALL gr_fi_ecrit(1, klon,iim,jjmp1, agesno,zx_tmp_2d)
2201      CALL histwrite(nid_mth,"ages",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2202c
2203      CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
2204      CALL histwrite(nid_mth,"evap",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2205c
2206      CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw,zx_tmp_2d)
2207      CALL histwrite(nid_mth,"tops",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2208c
2209      CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
2210      CALL histwrite(nid_mth,"topl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2211c
2212      CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
2213      CALL histwrite(nid_mth,"sols",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2214c
2215      CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
2216      CALL histwrite(nid_mth,"soll",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2217c
2218      CALL gr_fi_ecrit(1, klon,iim,jjmp1, topsw0,zx_tmp_2d)
2219      CALL histwrite(nid_mth,"tops0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2220c
2221      CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw0,zx_tmp_2d)
2222      CALL histwrite(nid_mth,"topl0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2223c
2224      CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw0,zx_tmp_2d)
2225      CALL histwrite(nid_mth,"sols0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2226c
2227      CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw0,zx_tmp_2d)
2228      CALL histwrite(nid_mth,"soll0",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2229c
2230      CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
2231      CALL histwrite(nid_mth,"bils",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2232c
2233      CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
2234      CALL histwrite(nid_mth,"sens",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2235c
2236      CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
2237      CALL histwrite(nid_mth,"fder",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2238c
2239      CALL gr_fi_ecrit(1, klon,iim,jjmp1, ruis,zx_tmp_2d)
2240      CALL histwrite(nid_mth,"ruis",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2241c
2242      DO i = 1, klon
2243         zx_tmp_fi2d(i) = fluxu(i,1)
2244      ENDDO
2245      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2246      CALL histwrite(nid_mth,"frtu",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2247c
2248      DO i = 1, klon
2249         zx_tmp_fi2d(i) = fluxv(i,1)
2250      ENDDO
2251      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2252      CALL histwrite(nid_mth,"frtv",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2253c
2254      DO i = 1, klon
2255         zx_tmp_fi2d(i) = pctsrf(i,is_sic)
2256      ENDDO
2257      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2258      CALL histwrite(nid_mth,"sicf",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2259c
2260      CALL gr_fi_ecrit(1, klon,iim,jjmp1, albsol,zx_tmp_2d)
2261      CALL histwrite(nid_mth,"albs",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2262c
2263      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragm,zx_tmp_2d)
2264      CALL histwrite(nid_mth,"cdrm",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2265c
2266      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cdragh,zx_tmp_2d)
2267      CALL histwrite(nid_mth,"cdrh",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2268c
2269      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldl,zx_tmp_2d)
2270      CALL histwrite(nid_mth,"cldl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2271c
2272      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldm,zx_tmp_2d)
2273      CALL histwrite(nid_mth,"cldm",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2274c
2275      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldh,zx_tmp_2d)
2276      CALL histwrite(nid_mth,"cldh",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2277c
2278      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldt,zx_tmp_2d)
2279      CALL histwrite(nid_mth,"cldt",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2280c
2281      CALL gr_fi_ecrit(1, klon,iim,jjmp1, cldq,zx_tmp_2d)
2282      CALL histwrite(nid_mth,"cldq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2283c
2284      CALL gr_fi_ecrit(1, klon,iim,jjmp1, ue,zx_tmp_2d)
2285      CALL histwrite(nid_mth,"ue",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2286c
2287      CALL gr_fi_ecrit(1, klon,iim,jjmp1, ve,zx_tmp_2d)
2288      CALL histwrite(nid_mth,"ve",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2289c
2290      CALL gr_fi_ecrit(1, klon,iim,jjmp1, uq,zx_tmp_2d)
2291      CALL histwrite(nid_mth,"uq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2292c
2293      CALL gr_fi_ecrit(1, klon,iim,jjmp1, vq,zx_tmp_2d)
2294      CALL histwrite(nid_mth,"vq",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2295c
2296c Champs 3D:
2297C
2298      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
2299      CALL histwrite(nid_mth,"temp",itap,zx_tmp_3d,
2300     .                                   iim*jjmp1*klev,ndex3d)
2301c
2302      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,ivap), zx_tmp_3d)
2303      CALL histwrite(nid_mth,"ovap",itap,zx_tmp_3d,
2304     .                                   iim*jjmp1*klev,ndex3d)
2305c
2306      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zphi, zx_tmp_3d)
2307      CALL histwrite(nid_mth,"geop",itap,zx_tmp_3d,
2308     .                                   iim*jjmp1*klev,ndex3d)
2309c
2310      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
2311      CALL histwrite(nid_mth,"vitu",itap,zx_tmp_3d,
2312     .                                   iim*jjmp1*klev,ndex3d)
2313c
2314      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
2315      CALL histwrite(nid_mth,"vitv",itap,zx_tmp_3d,
2316     .                                   iim*jjmp1*klev,ndex3d)
2317c
2318      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, omega, zx_tmp_3d)
2319      CALL histwrite(nid_mth,"vitw",itap,zx_tmp_3d,
2320     .                                   iim*jjmp1*klev,ndex3d)
2321c
2322      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, pplay, zx_tmp_3d)
2323      CALL histwrite(nid_mth,"pres",itap,zx_tmp_3d,
2324     .                                   iim*jjmp1*klev,ndex3d)
2325c
2326      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldfra, zx_tmp_3d)
2327      CALL histwrite(nid_mth,"rneb",itap,zx_tmp_3d,
2328     .                                   iim*jjmp1*klev,ndex3d)
2329c
2330      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, zx_rh, zx_tmp_3d)
2331      CALL histwrite(nid_mth,"rhum",itap,zx_tmp_3d,
2332     .                                   iim*jjmp1*klev,ndex3d)
2333c
2334      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cldliq, zx_tmp_3d)
2335      CALL histwrite(nid_mth,"oliq",itap,zx_tmp_3d,
2336     .                                   iim*jjmp1*klev,ndex3d)
2337c
2338      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_dyn, zx_tmp_3d)
2339      CALL histwrite(nid_mth,"dtdyn",itap,zx_tmp_3d,
2340     .                                   iim*jjmp1*klev,ndex3d)
2341c
2342      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_dyn, zx_tmp_3d)
2343      CALL histwrite(nid_mth,"dqdyn",itap,zx_tmp_3d,
2344     .                                   iim*jjmp1*klev,ndex3d)
2345c
2346      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_con, zx_tmp_3d)
2347      CALL histwrite(nid_mth,"dtcon",itap,zx_tmp_3d,
2348     .                                   iim*jjmp1*klev,ndex3d)
2349c
2350      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_con, zx_tmp_3d)
2351      CALL histwrite(nid_mth,"dqcon",itap,zx_tmp_3d,
2352     .                                   iim*jjmp1*klev,ndex3d)
2353c
2354      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_lsc, zx_tmp_3d)
2355      CALL histwrite(nid_mth,"dtlsc",itap,zx_tmp_3d,
2356     .                                   iim*jjmp1*klev,ndex3d)
2357c
2358      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_lsc, zx_tmp_3d)
2359      CALL histwrite(nid_mth,"dqlsc",itap,zx_tmp_3d,
2360     .                                   iim*jjmp1*klev,ndex3d)
2361c
2362      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_vdf, zx_tmp_3d)
2363      CALL histwrite(nid_mth,"dtvdf",itap,zx_tmp_3d,
2364     .                                   iim*jjmp1*klev,ndex3d)
2365c
2366      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_vdf, zx_tmp_3d)
2367      CALL histwrite(nid_mth,"dqvdf",itap,zx_tmp_3d,
2368     .                                   iim*jjmp1*klev,ndex3d)
2369c
2370      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_eva, zx_tmp_3d)
2371      CALL histwrite(nid_mth,"dteva",itap,zx_tmp_3d,
2372     .                                   iim*jjmp1*klev,ndex3d)
2373c
2374      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_eva, zx_tmp_3d)
2375      CALL histwrite(nid_mth,"dqeva",itap,zx_tmp_3d,
2376     .                                   iim*jjmp1*klev,ndex3d)
2377c
2378      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_ajs, zx_tmp_3d)
2379      CALL histwrite(nid_mth,"dtajs",itap,zx_tmp_3d,
2380     .                                   iim*jjmp1*klev,ndex3d)
2381c
2382      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_ajs, zx_tmp_3d)
2383      CALL histwrite(nid_mth,"dqajs",itap,zx_tmp_3d,
2384     .                                   iim*jjmp1*klev,ndex3d)
2385c
2386      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, heat, zx_tmp_3d)
2387      CALL histwrite(nid_mth,"dtswr",itap,zx_tmp_3d,
2388     .                                   iim*jjmp1*klev,ndex3d)
2389c
2390      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, heat0, zx_tmp_3d)
2391      CALL histwrite(nid_mth,"dtsw0",itap,zx_tmp_3d,
2392     .                                   iim*jjmp1*klev,ndex3d)
2393c
2394      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cool, zx_tmp_3d)
2395      CALL histwrite(nid_mth,"dtlwr",itap,zx_tmp_3d,
2396     .                                   iim*jjmp1*klev,ndex3d)
2397c
2398      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, cool0, zx_tmp_3d)
2399      CALL histwrite(nid_mth,"dtlw0",itap,zx_tmp_3d,
2400     .                                   iim*jjmp1*klev,ndex3d)
2401c
2402      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_vdf, zx_tmp_3d)
2403      CALL histwrite(nid_mth,"duvdf",itap,zx_tmp_3d,
2404     .                                   iim*jjmp1*klev,ndex3d)
2405c
2406      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_vdf, zx_tmp_3d)
2407      CALL histwrite(nid_mth,"dvvdf",itap,zx_tmp_3d,
2408     .                                   iim*jjmp1*klev,ndex3d)
2409c
2410      IF (ok_orodr) THEN
2411      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_oro, zx_tmp_3d)
2412      CALL histwrite(nid_mth,"duoro",itap,zx_tmp_3d,
2413     .                                   iim*jjmp1*klev,ndex3d)
2414c
2415      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_oro, zx_tmp_3d)
2416      CALL histwrite(nid_mth,"dvoro",itap,zx_tmp_3d,
2417     .                                   iim*jjmp1*klev,ndex3d)
2418c
2419      ENDIF
2420C
2421      IF (ok_orolf) THEN
2422      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_u_lif, zx_tmp_3d)
2423      CALL histwrite(nid_mth,"dulif",itap,zx_tmp_3d,
2424     .                                   iim*jjmp1*klev,ndex3d)
2425c
2426      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_v_lif, zx_tmp_3d)
2427      CALL histwrite(nid_mth,"dvlif",itap,zx_tmp_3d,
2428     .                                   iim*jjmp1*klev,ndex3d)
2429      ENDIF
2430C
2431      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, wo, zx_tmp_3d)
2432      CALL histwrite(nid_mth,"ozone",itap,zx_tmp_3d,
2433     .                                   iim*jjmp1*klev,ndex3d)
2434c
2435      IF (nqmax.GE.3) THEN
2436      DO iq=1,nqmax-2
2437      IF (iq.LE.99) THEN
2438         CALL gr_fi_ecrit(klev,klon,iim,jjmp1, qx(1,1,iq+2), zx_tmp_3d)
2439         WRITE(str2,'(i2.2)') iq
2440         CALL histwrite(nid_mth,"trac"//str2,itap,zx_tmp_3d,
2441     .                                   iim*jjmp1*klev,ndex3d)
2442      ELSE
2443         PRINT*, "Trop de traceurs"
2444         CALL abort
2445      ENDIF
2446      ENDDO
2447      ENDIF
2448c
2449      if (ok_sync) then
2450        call histsync(nid_mth)
2451      endif
2452      ENDIF
2453c
2454      IF (ok_instan) THEN
2455c
2456      ndex2d = 0
2457      ndex3d = 0
2458c
2459c Champs 2D:
2460c
2461         i = NINT(zout/zsto)
2462         CALL gr_fi_ecrit(1,klon,iim,jjmp1,pphis,zx_tmp_2d)
2463         CALL histwrite(nid_ins,"phis",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2464c
2465         i = NINT(zout/zsto)
2466         CALL gr_fi_ecrit(1,klon,iim,jjmp1,paire,zx_tmp_2d)
2467         CALL histwrite(nid_ins,"aire",i,zx_tmp_2d,iim*jjmp1,ndex2d)
2468
2469      DO i = 1, klon
2470         zx_tmp_fi2d(i) = paprs(i,1)
2471      ENDDO
2472      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zx_tmp_fi2d,zx_tmp_2d)
2473      CALL histwrite(nid_ins,"psol",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2474c
2475      CALL gr_fi_ecrit(1, klon,iim,jjmp1, toplw,zx_tmp_2d)
2476      CALL histwrite(nid_ins,"topl",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2477c
2478      CALL gr_fi_ecrit(1, klon,iim,jjmp1, evap,zx_tmp_2d)
2479      CALL histwrite(nid_ins,"evap",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2480c
2481      CALL gr_fi_ecrit(1, klon,iim,jjmp1, solsw,zx_tmp_2d)
2482      CALL histwrite(nid_ins,"sols",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2483c
2484      CALL gr_fi_ecrit(1, klon,iim,jjmp1, sollw,zx_tmp_2d)
2485      CALL histwrite(nid_ins,"soll",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2486c
2487      CALL gr_fi_ecrit(1, klon,iim,jjmp1, bils,zx_tmp_2d)
2488      CALL histwrite(nid_ins,"bils",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2489c
2490      CALL gr_fi_ecrit(1, klon,iim,jjmp1, sens,zx_tmp_2d)
2491      CALL histwrite(nid_ins,"sens",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2492c
2493      CALL gr_fi_ecrit(1, klon,iim,jjmp1, fder,zx_tmp_2d)
2494      CALL histwrite(nid_ins,"fder",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2495c
2496      CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_oce),zx_tmp_2d)
2497      CALL histwrite(nid_ins,"dtsvdfo",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2498c
2499      CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_ter),zx_tmp_2d)
2500      CALL histwrite(nid_ins,"dtsvdft",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2501c
2502      CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_lic),zx_tmp_2d)
2503      CALL histwrite(nid_ins,"dtsvdfg",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2504c
2505      CALL gr_fi_ecrit(1, klon,iim,jjmp1, d_ts(1,is_sic),zx_tmp_2d)
2506      CALL histwrite(nid_ins,"dtsvdfi",itap,zx_tmp_2d,iim*jjmp1,ndex2d)
2507
2508c
2509c Champs 3D:
2510c
2511      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, t_seri, zx_tmp_3d)
2512      CALL histwrite(nid_ins,"temp",itap,zx_tmp_3d,
2513     .                                   iim*jjmp1*klev,ndex3d)
2514c
2515      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, u_seri, zx_tmp_3d)
2516      CALL histwrite(nid_ins,"vitu",itap,zx_tmp_3d,
2517     .                                   iim*jjmp1*klev,ndex3d)
2518c
2519      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, v_seri, zx_tmp_3d)
2520      CALL histwrite(nid_ins,"vitv",itap,zx_tmp_3d,
2521     .                                   iim*jjmp1*klev,ndex3d)
2522c
2523      CALL gr_fi_ecrit(klev,klon,iim,jjm+1, zphi, zx_tmp_3d)
2524      CALL histwrite(nid_ins,"geop",itap,zx_tmp_3d,
2525     .                                   iim*jjmp1*klev,ndex3d)
2526c
2527      CALL gr_fi_ecrit(klev,klon,iim,jjm+1, pplay, zx_tmp_3d)
2528      CALL histwrite(nid_ins,"pres",itap,zx_tmp_3d,
2529     .                                   iim*jjmp1*klev,ndex3d)
2530c
2531      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_t_vdf, zx_tmp_3d)
2532      CALL histwrite(nid_ins,"dtvdf",itap,zx_tmp_3d,
2533     .                                   iim*jjmp1*klev,ndex3d)
2534c
2535      CALL gr_fi_ecrit(klev,klon,iim,jjmp1, d_q_vdf, zx_tmp_3d)
2536      CALL histwrite(nid_ins,"dqvdf",itap,zx_tmp_3d,
2537     .                                   iim*jjmp1*klev,ndex3d)
2538
2539c
2540      if (ok_sync) then
2541        call histsync(nid_ins)
2542      endif
2543      ENDIF
2544c
2545      IF (ok_oasis .AND. mod(itap,nexca).EQ.0) THEN
2546c
2547c Je ne traite pas le ruissellement, pour l'instant (qui m'aidera ?)
2548         DO i = 1, klon
2549            oas_ruisoce(i) = 0.0
2550            oas_ruisriv(i) = 0.0
2551         ENDDO
2552c
2553         ig = 1
2554         DO i = 1, iim
2555            z_sols(i,1) = oas_sols(ig)
2556            z_nsol(i,1) = oas_nsol(ig)
2557            z_rain(i,1) = oas_rain(ig)
2558            z_snow(i,1) = oas_snow(ig)
2559            z_evap(i,1) = oas_evap(ig)
2560            z_ruisoce(i,1) = oas_ruisoce(ig)
2561            z_ruisriv(i,1) = oas_ruisriv(ig)
2562            z_tsol(i,1) = oas_tsol(ig)
2563            z_fder(i,1) = oas_fder(ig)
2564            z_albe(i,1) = oas_albe(ig)
2565            z_taux(i,1) = oas_taux(ig)
2566            z_tauy(i,1) = oas_tauy(ig)
2567         ENDDO
2568         DO j = 2, jjm
2569         DO i = 1, iim
2570            ig = ig + 1
2571            z_sols(i,j) = oas_sols(ig)
2572            z_nsol(i,j) = oas_nsol(ig)
2573            z_rain(i,j) = oas_rain(ig)
2574            z_snow(i,j) = oas_snow(ig)
2575            z_evap(i,j) = oas_evap(ig)
2576            z_ruisoce(i,j) = oas_ruisoce(ig)
2577            z_ruisriv(i,j) = oas_ruisriv(ig)
2578            z_tsol(i,j) = oas_tsol(ig)
2579            z_fder(i,j) = oas_fder(ig)
2580            z_albe(i,j) = oas_albe(ig)
2581            z_taux(i,j) = oas_taux(ig)
2582            z_tauy(i,j) = oas_tauy(ig)
2583         ENDDO
2584         ENDDO
2585         ig = ig + 1
2586         DO i = 1, iim
2587            z_sols(i,jjmp1)    = oas_sols(ig)
2588            z_nsol(i,jjmp1)    = oas_nsol(ig)
2589            z_rain(i,jjmp1)    = oas_rain(ig)
2590            z_snow(i,jjmp1)    = oas_snow(ig)
2591            z_evap(i,jjmp1)    = oas_evap(ig)
2592            z_ruisoce(i,jjmp1) = oas_ruisoce(ig)
2593            z_ruisriv(i,jjmp1) = oas_ruisriv(ig)
2594            z_tsol(i,jjmp1)    = oas_tsol(ig)
2595            z_fder(i,jjmp1)    = oas_fder(ig)
2596            z_albe(i,jjmp1)    = oas_albe(ig)
2597            z_taux(i,jjmp1)    = oas_taux(ig)
2598            z_tauy(i,jjmp1)    = oas_tauy(ig)
2599         ENDDO
2600c
2601c Passer les champs au coupleur:
2602c
2603         CALL intocpl(itap,jjmp1*iim,
2604     .                   z_sols, z_nsol,
2605     .                   z_rain, z_snow, z_evap,
2606     .                   z_ruisoce, z_ruisriv,
2607     .                   z_tsol, z_fder, z_albe,
2608     .                   z_taux, z_tauy)
2609         DO i = 1, klon
2610           oas_sols(i) = 0.0
2611           oas_nsol(i) = 0.0
2612           oas_rain(i) = 0.0
2613           oas_snow(i) = 0.0
2614           oas_evap(i) = 0.0
2615           oas_ruis(i) = 0.0
2616           oas_tsol(i) = 0.0
2617           oas_fder(i) = 0.0
2618           oas_albe(i) = 0.0
2619           oas_taux(i) = 0.0
2620           oas_tauy(i) = 0.0
2621         ENDDO
2622      ENDIF
2623c
2624c Ecrire la bande regionale (binaire grads)
2625      IF (ok_region .AND. mod(itap,ecrit_reg).eq.0) THEN
2626         CALL ecriregs(84,zxtsol)
2627         CALL ecriregs(84,paprs(1,1))
2628         CALL ecriregs(84,topsw)
2629         CALL ecriregs(84,toplw)
2630         CALL ecriregs(84,solsw)
2631         CALL ecriregs(84,sollw)
2632         CALL ecriregs(84,rain_fall)
2633         CALL ecriregs(84,snow_fall)
2634         CALL ecriregs(84,evap)
2635         CALL ecriregs(84,sens)
2636         CALL ecriregs(84,bils)
2637         CALL ecriregs(84,pctsrf(1,is_sic))
2638         CALL ecriregs(84,fluxu(1,1))
2639         CALL ecriregs(84,fluxv(1,1))
2640         CALL ecriregs(84,ue)
2641         CALL ecriregs(84,ve)
2642         CALL ecriregs(84,uq)
2643         CALL ecriregs(84,vq)
2644c
2645         CALL ecrirega(84,u_seri)
2646         CALL ecrirega(84,v_seri)
2647         CALL ecrirega(84,omega)
2648         CALL ecrirega(84,t_seri)
2649         CALL ecrirega(84,zphi)
2650         CALL ecrirega(84,q_seri)
2651         CALL ecrirega(84,cldfra)
2652         CALL ecrirega(84,cldliq)
2653         CALL ecrirega(84,pplay)
2654
2655
2656cc         CALL ecrirega(84,d_t_dyn)
2657cc         CALL ecrirega(84,d_q_dyn)
2658cc         CALL ecrirega(84,heat)
2659cc         CALL ecrirega(84,cool)
2660cc         CALL ecrirega(84,d_t_con)
2661cc         CALL ecrirega(84,d_q_con)
2662cc         CALL ecrirega(84,d_t_lsc)
2663cc         CALL ecrirega(84,d_q_lsc)
2664      ENDIF
2665c
2666c Convertir les incrementations en tendances
2667c
2668      DO k = 1, klev
2669      DO i = 1, klon
2670         d_u(i,k) = ( u_seri(i,k) - u(i,k) ) / dtime
2671         d_v(i,k) = ( v_seri(i,k) - v(i,k) ) / dtime
2672         d_t(i,k) = ( t_seri(i,k)-t(i,k) ) / dtime
2673         d_qx(i,k,ivap) = ( q_seri(i,k) - qx(i,k,ivap) ) / dtime
2674         d_qx(i,k,iliq) = ( ql_seri(i,k) - qx(i,k,iliq) ) / dtime
2675      ENDDO
2676      ENDDO
2677c
2678      IF (nqmax.GE.3) THEN
2679      DO iq = 3, nqmax
2680      DO  k = 1, klev
2681      DO  i = 1, klon
2682         d_qx(i,k,iq) = ( tr_seri(i,k,iq-2) - qx(i,k,iq) ) / dtime
2683      ENDDO
2684      ENDDO
2685      ENDDO
2686      ENDIF
2687c
2688c Sauvegarder les valeurs de t et q a la fin de la physique:
2689c
2690      DO k = 1, klev
2691      DO i = 1, klon
2692         t_ancien(i,k) = t_seri(i,k)
2693         q_ancien(i,k) = q_seri(i,k)
2694      ENDDO
2695      ENDDO
2696c
2697c====================================================================
2698c Si c'est la fin, il faut conserver l'etat de redemarrage
2699c====================================================================
2700c
2701      IF (lafin) THEN
2702ccc         IF (ok_oasis) CALL quitcpl
2703         CALL phyredem ("restartphy.nc",dtime,radpas,co2_ppm,solaire,
2704     .        rlat,rlon,ftsol,ftsoil,deltat,fqsol,fsnow,
2705     .        radsol,rugmer,agesno,
2706     .        zmea,zstd,zsig,zgam,zthe,zpic,zval,rugoro,
2707     .        t_ancien, q_ancien)
2708      ENDIF
2709
2710      RETURN
2711      END
2712      FUNCTION qcheck(klon,klev,paprs,q,ql,aire)
2713      IMPLICIT none
2714c
2715c Calculer et imprimer l'eau totale. A utiliser pour verifier
2716c la conservation de l'eau
2717c
2718#include "YOMCST.h"
2719      INTEGER klon,klev
2720      REAL paprs(klon,klev+1), q(klon,klev), ql(klon,klev)
2721      REAL aire(klon)
2722      REAL qtotal, zx, qcheck
2723      INTEGER i, k
2724c
2725      zx = 0.0
2726      DO i = 1, klon
2727         zx = zx + aire(i)
2728      ENDDO
2729      qtotal = 0.0
2730      DO k = 1, klev
2731      DO i = 1, klon
2732         qtotal = qtotal + (q(i,k)+ql(i,k)) * aire(i)
2733     .                     *(paprs(i,k)-paprs(i,k+1))/RG
2734      ENDDO
2735      ENDDO
2736c
2737      qcheck = qtotal/zx
2738c
2739      RETURN
2740      END
2741      SUBROUTINE gr_fi_ecrit(nfield,nlon,iim,jjmp1,fi,ecrit)
2742      IMPLICIT none
2743c
2744c Tranformer une variable de la grille physique a
2745c la grille d'ecriture
2746c
2747      INTEGER nfield,nlon,iim,jjmp1, jjm
2748      REAL fi(nlon,nfield), ecrit(iim*jjmp1,nfield)
2749c
2750      INTEGER i, j, n, ig
2751c
2752      jjm = jjmp1 - 1
2753      DO n = 1, nfield
2754         DO i=1,iim
2755            ecrit(i,n) = fi(1,n)
2756            ecrit(i+jjm*iim,n) = fi(nlon,n)
2757         ENDDO
2758         DO ig = 1, nlon - 2
2759           ecrit(iim+ig,n) = fi(1+ig,n)
2760         ENDDO
2761      ENDDO
2762      RETURN
2763      END
2764
Note: See TracBrowser for help on using the repository browser.