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

Last change on this file since 171 was 171, checked in by lmdzadmin, 23 years ago

Changements necessaires pour l'appel a orchidee
LF

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