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

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

Rajout de la derivee des flux dans le fichier restart physique + menage
LF

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