Changeset 344


Ignore:
Timestamp:
Mar 6, 2002, 3:58:31 PM (23 years ago)
Author:
lmdz
Message:

Inclusion des modifs de D. Hauglustaine pour la version 1 de INCA
LF

Location:
LMDZ.3.3/trunk/libf
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • LMDZ.3.3/trunk/libf/dyn3d/calfis.F

    r343 r344  
    1       SUBROUTINE calfis(nq, lafin, rdayvrai,rday_ecri, heure,
    2      $            pucov,pvcov,pteta,pq,pmasse,pps,pp,ppk,pphis,pphi,
    3      $            pducov,pdvcov,pdteta,pdq,pw, clesphy0,
    4      $            pdufi,pdvfi,pdhfi,pdqfi,pdpsfi                   )
     1      SUBROUTINE calfis(nq,
     2     $                  lafin,
     3     $                  rdayvrai,
     4     $                  rday_ecri,
     5     $                  heure,
     6     $                  pucov,
     7     $                  pvcov,
     8     $                  pteta,
     9     $                  pq,
     10     $                  pmasse,
     11     $                  pps,
     12     $                  pp,
     13     $                  ppk,
     14     $                  pphis,
     15     $                  pphi,
     16     $                  pducov,
     17     $                  pdvcov,
     18     $                  pdteta,
     19     $                  pdq,
     20     $                  pw,
     21#ifdef INCA_CH4
     22     $                  flxw,   
     23#endif
     24     $                  clesphy0,
     25     $                  pdufi,
     26     $                  pdvfi,
     27     $                  pdhfi,
     28     $                  pdqfi,
     29     $                  pdpsfi)
    530c
    631c    Auteur :  P. Le Van, F. Hourdin
    … …  
    132157      REAL zsinbis(iim),zcosbis(iim),z1bis(iim)
    133158      REAL unskap, pksurcp
    134 c
     159 
     160#ifdef INCA_CH4
     161      REAL flxw(iip1,jjp1,llm)
     162      REAL flxwfi(ngridmx,llm)
     163#endif
    135164     
    136165      EXTERNAL gr_dyn_fi,gr_fi_dyn
    … …  
    421450      ENDDO
    422451
     452#ifdef INCA_CH4
     453      CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,flxw,flxwfi)
     454#endif
     455
    423456c-----------------------------------------------------------------------
    424457c   Appel de la physique:
    … …  
    426459
    427460
    428       CALL physiq (ngridmx,llm,nq,debut,lafin,
    429      ,   rdayvrai,rday_ecri,heure,dtphys,zplev,zplay,zphi,zphis,airefi,
    430      ,             presnivs,clesphy0,  zufi, zvfi,ztfi, zqfi, 
    431 ccc     ,             pcvgu, pcvgv, pcvgt, pcvgq,
    432      ,             pvervel,
    433 C - sorties
    434      s             zdufi, zdvfi, zdtfi, zdqfi,zdpsrf              )
     461      CALL physiq (ngridmx,
     462     .             llm,
     463     .             nq,
     464     .             debut,
     465     .             lafin,
     466     .             rdayvrai,
     467     .             rday_ecri,
     468     .             heure,
     469     .             dtphys,
     470     .             zplev,
     471     .             zplay,
     472     .             zphi,
     473     .             zphis,
     474     .             airefi,
     475     .             presnivs,
     476     .             clesphy0,
     477     .             zufi,
     478     .             zvfi,
     479     .             ztfi,
     480     .             zqfi, 
     481     .             pvervel,
     482#ifdef INCA_CH4
     483     .             flxwfi,
     484#endif
     485     .             zdufi,
     486     .             zdvfi,
     487     .             zdtfi,
     488     .             zdqfi,
     489     .             zdpsrf)
    435490
    436491500   CONTINUE

  • LMDZ.3.3/trunk/libf/dyn3d/gcm.F

    r195 r344  
    22      PROGRAM gcm
    33
     4#ifdef INCA
     5      USE transport_controls, ONLY : adv_flg, mmt_adj
     6#endif 
    47      USE IOIPSL
    58
    … …  
    127130      REAL alpha(ip1jmp1,llm),beta(ip1jmp1,llm)
    128131
     132#ifdef INCA_CH4
     133      REAL :: flxw(ip1jmp1,llm)
     134#endif
     135
    129136      LOGICAL offline  ! Controle du stockage ds "fluxmass"
    130137      PARAMETER (offline=.true.)
    … …  
    409416cccc      de masse pour  Van-Leer dans la routine  tracvl  .
    410417c
    411             CALL vanleer(numvanle,iapp_tracvl,nqmx,q,pbaru,pbarv,
    412      *                      p, masse, dq,  iadv(1), teta, pk     )
     418      CALL vanleer(numvanle,
     419     .             iapp_tracvl,
     420     .             nqmx,
     421     .             q,
     422     .             pbaru,
     423     .             pbarv,
     424     .             p,
     425     .             masse,
     426     .             dq, 
     427     .             iadv(1),
     428     .             teta,
     429#ifdef INCA_CH4
     430     .             flxw,
     431     .             pk,
     432     .             mmt_adj,
     433     .             adv_flg)
     434#else
     435     .             pk)
     436#endif
    413437c
    414438c                   ...  Modif  F.Codron  ....
    … …  
    470494           ENDIF
    471495c
    472         CALL calfis( nqmx, lafin ,rdayvrai,rday_ecri,time  ,
    473      $                 ucov,vcov,teta,q,masse,ps,p,pk,phis,phi ,
    474      $     du,dv,dteta,dq,w,clesphy0, dufi,dvfi,dhfi,dqfi,dpfi  )
     496        CALL calfis( nqmx,
     497     $               lafin,
     498     $               rdayvrai,
     499     $               rday_ecri,
     500     $               time,
     501     $               ucov,
     502     $               vcov,
     503     $               teta,
     504     $               q,
     505     $               masse,
     506     $               ps,
     507     $               p,
     508     $               pk,
     509     $               phis,
     510     $               phi,
     511     $               du,
     512     $               dv,
     513     $               dteta,
     514     $               dq,
     515     $               w,
     516#ifdef INCA_CH4
     517     $               flxw,
     518#endif
     519     $               clesphy0,
     520     $               dufi,
     521     $               dvfi,
     522     $               dhfi,
     523     $               dqfi,
     524     $               dpfi)
    475525
    476526c      ajout des tendances physiques:

  • LMDZ.3.3/trunk/libf/dyn3d/tracvl.F

    r78 r344  
    1       SUBROUTINE tracvl(numvanle,iapp_tracvl,nq,pbaru,pbarv ,
    2      *                    p, masse , q, iapptrac, iadv1, teta, pk  )
     1      SUBROUTINE tracvl(numvanle,
     2     *                  iapp_tracvl,
     3     *                  nq,
     4     *                  pbaru,
     5     *                  pbarv,
     6     *                  p,
     7     *                  masse,
     8     *                  q,
     9     *                  iapptrac,
     10     *                  iadv1,
     11     *                  teta,
     12#ifdef INCA_CH4
     13     *                  flxw,
     14     *                  pk,
     15     *                  mmt_adj,
     16     *                  adv_flg)
     17#else
     18     *                  pk)
     19#endif
    320c
    421c     Auteur :  F. Hourdin
    … …  
    2037c
    2138      INTEGER numvanle, nq, iapp_tracvl, iapptrac, iadv1
    22 
    2339      REAL pbaru(ip1jmp1,llm),pbarv(ip1jm,llm)
    2440      REAL q(ip1jmp1,llm,nq),masse(ip1jmp1,llm)
    2541      REAL p( ip1jmp1,llmp1 ),teta(ip1jmp1,llm)
    2642      REAL pk(ip1jmp1,llm)
    27 
     43#ifdef INCA_CH4
     44      INTEGER, PARAMETER :: ntra   = 1
     45      INTEGER, PARAMETER :: nprath = 1
     46      INTEGER            :: adv_flg(nq)
     47      REAL               :: mmt_adj(ip1jmp1,llm,nprath)
     48      REAL, SAVE         :: qpente(ip1jmp1,llm,10,nprath)
     49      REAL               :: flxw(ip1jmp1,llm)
     50#endif
    2851c     ....  var. locales  .....
    2952c
    … …  
    82105        CALL groupe( massem, pbaruc,pbarvc, pbarug,pbarvg,wg )
    83106
     107#ifdef INCA_CH4
     108      ! ... Flux de masse diaganostiques traceurs
     109      flxw = wg / FLOAT(iapp_tracvl)
     110#endif
    84111
    85112c  test sur l'eventuelle creation de valeurs negatives de la masse
    … …  
    115142         ENDIF
    116143
    117          DO iq = numvan, nq
    118           CALL vlsplt( q(1,1,iq), 2. ,massem,wg,pbarug,pbarvg,dtvr )
    119          ENDDO
     144#ifdef INCA_CH4
     145      do iq = 2, 10
     146        qpente(:,:,iq,1)=qpente(:,:,iq,1)*mmt_adj(:,:,1)
     147      enddo
     148#endif
    120149
    121          iadvtr=0
     150      DO iq = numvan, 2
     151#ifdef INCA
     152      IF (adv_flg(iq) == 0) CYCLE
     153#endif
     154      CALL vlsplt(q(1,1,iq),2.,massem,wg,pbarug,pbarvg,dtvr)
     155      ENDDO
     156
     157#ifdef INCA_CH4
     158!     CALL prather(q(1,1,3),wg,massem,pbarug,pbarvg,ntra,qpente(1,1,1,1))
     159#endif
     160
     161      DO  iq =3,nq
     162#ifdef INCA
     163      IF (adv_flg(iq) == 0) CYCLE
     164#endif
     165      CALL vlsplt(q(1,1,iq),2.,massem,wg,pbarug,pbarvg,dtvr)
     166      ENDDO
     167
     168      iadvtr=0
    122169
    123170c   on reinitialise a zero les flux de masse cumules.

  • LMDZ.3.3/trunk/libf/dyn3d/vanleer.F

    r78 r344  
    1       SUBROUTINE vanleer(numvanle,iapp_tracvl,nq,q,pbaru,pbarv ,
    2      *                     p ,masse, dq ,  iadv1, teta, pk      )
     1      SUBROUTINE vanleer(numvanle,
     2     *                   iapp_tracvl,
     3     *                   nq,
     4     *                   q,
     5     *                   pbaru,
     6     *                   pbarv,
     7     *                   p,
     8     *                   masse,
     9     *                   dq,
     10     *                   iadv1,
     11     *                   teta,
     12#ifdef INCA_CH4
     13     *                   flxw,
     14     *                   pk,
     15     *                   mmt_adj,
     16     *                   adv_flg)
     17#else
     18     *                   pk)
     19#endif
    320c
    421      IMPLICIT NONE
    … …  
    2441      REAL p( ip1jmp1,llmp1),q( ip1jmp1,llm,nq),dq( ip1jmp1,llm,nq )
    2542      REAL teta( ip1jmp1,llm),pk( ip1jmp1,llm)
     43#ifdef INCA_CH4
     44      INTEGER, PARAMETER :: nprath=1
     45      INTEGER            :: adv_flg(nq)
     46      REAL               :: mmt_adj(iip1,jjp1,llm,nprath)
     47      REAL               :: flxw(ip1jmp1,llm)
     48#endif
    2649c  ..................................................................
    2750c
    … …  
    6487c   advection
    6588
    66       CALL tracvl( numvanle,iapp_tracvl,nq,pbaru,pbarv,p , masse,q  ,
    67      *                      iapptrac, iadv1, teta ,pk              )
     89      CALL tracvl( numvanle,
     90     .             iapp_tracvl,
     91     .             nq,
     92     .             pbaru,
     93     .             pbarv,
     94     .             p,
     95     .             masse,
     96     .             q,
     97     .             iapptrac,
     98     .             iadv1,
     99     .             teta,
     100#ifdef INCA_CH4
     101     .             flxw,
     102     .             pk,
     103     .             mmt_adj,
     104     .             adv_flg)
     105#else
     106     .             pk)
     107#endif
    68108
    69109      IF( numvanle.EQ.1 ) THEN

  • LMDZ.3.3/trunk/libf/phylmd/physiq.F

    r303 r344  
    11c $Header$
    22c
    3       SUBROUTINE physiq (nlon,nlev,nqmax  ,
    4      .            debut,lafin,rjourvrai,rjour_ecri,gmtime,pdtphys,
    5      .            paprs,pplay,pphi,pphis,paire,presnivs,clesphy0,
    6      .            u,v,t,qx,
    7      .            omega,
    8      .            d_u, d_v, d_t, d_qx, d_ps)
     3      SUBROUTINE physiq (nlon,
     4     .                   nlev,
     5     .                   nqmax,
     6     .                   debut,
     7     .                   lafin,
     8     .                   rjourvrai,
     9     .                   rjour_ecri,
     10     .                   gmtime,
     11     .                   pdtphys,
     12     .                   paprs,
     13     .                   pplay,
     14     .                   pphi,
     15     .                   pphis,
     16     .                   paire,
     17     .                   presnivs,
     18     .                   clesphy0,
     19     .                   u,
     20     .                   v,
     21     .                   t,
     22     .                   qx,
     23     .                   omega,
     24#ifdef INCA_CH4
     25     .                   flxmass_w,
     26#endif
     27     .                   d_u,
     28     .                   d_v,
     29     .                   d_t,
     30     .                   d_qx,
     31     .                   d_ps)
    932      USE ioipsl
     33#ifdef INCA
     34      USE chemshut
     35      USE obs_pos
     36#endif
    1037      IMPLICIT none
    1138c======================================================================
    … …  
    7198      INTEGER npas, nexca, itimestep
    7299      logical rnpb
     100#ifdef INCA
     101      parameter(rnpb=.false.)
     102#else
    73103      parameter(rnpb=.true.)
     104#endif
    74105      PARAMETER (npas=1440)
    75106      PARAMETER (nexca=48)
    … …  
    174205      REAL omega(klon,klev)
    175206
     207      REAL flxmass_w(klon,klev)
     208
    176209      REAL d_u(klon,klev)
    177210      REAL d_v(klon,klev)
    … …  
    326359      REAL frac_impa(klon,klev) ! fractions d'aerosols lessivees (impaction)
    327360      REAL frac_nucl(klon,klev) ! idem (nucleation)
     361
     362#ifdef INCA
     363      REAL          :: calday
     364#endif
     365
    328366cAA
    329367      REAL rain_fall(klon) ! pluie
    … …  
    403441c Variables locales
    404442c
     443      REAL rhcl(klon,klev)    ! humiditi relative ciel clair
    405444      REAL dialiq(klon,klev)  ! eau liquide nuageuse
    406445      REAL diafra(klon,klev)  ! fraction nuageuse
    … …  
    14011440         ENDIF
    14021441c
     1442
     1443#ifdef INCA
     1444           calday = FLOAT(MOD(NINT(xjour),360)) + gmtime
     1445           print *, 'initial time ', xjour, calday
     1446#ifdef INCAINFO
     1447           PRINT *, 'Appel CHEMINI ...'
     1448#endif
     1449           CALL chemini( rpi,
     1450     $                   rg,
     1451     $                   ra,
     1452     $                   paire,
     1453     $                   rlat,
     1454     $                   rlon,
     1455     $                   calday,
     1456     $                   tracnam,
     1457     $                   natsnam,
     1458     $                   mxoutflds,
     1459     $                   outinst,
     1460     $                   outtimav,
     1461     $                   klon,
     1462     $                   nqmax)
     1463#ifdef INCAINFO
     1464           PRINT *, 'OK.'
     1465#endif
     1466#endif
     1467
    14031468      ENDIF
    14041469c
    … …  
    18481913     .           pfrac_impa, pfrac_nucl, pfrac_1nucl,
    18491914     .           frac_impa, frac_nucl,
    1850      .           prfl, psfl)
     1915     .           prfl, psfl, rhcl)
    18511916      DO k = 1, klev
    18521917      DO i = 1, klon
    … …  
    19291994      ENDDO
    19301995      ENDDO
     1996
     1997#ifdef INCA
     1998           calday = FLOAT(julien) + gmtime
     1999
     2000#ifdef INCA_AER
     2001      call AEROSOL_METEO_CALC(calday,pdtphys,pplay,paprs,t,pmflxr,pmflxs,
     2002     &   prfl,psfl,pctsrf,paire)
     2003#endif
     2004
     2005#ifdef INCAINFO
     2006           PRINT *, 'Appel CHEMHOOK_BEGIN ...'
     2007#endif
     2008           CALL chemhook_begin (calday,
     2009     $                          pctsrf(1,3),
     2010     $                          rlat,
     2011     $                          rlon,
     2012     $                          paire,
     2013     $                          paprs,
     2014     $                          pplay,
     2015     $                          ycoefh,
     2016     $                          pphi,
     2017     $                          t_seri,
     2018     $                          u,
     2019     $                          v,
     2020     $                          wo,
     2021     $                          q_seri,
     2022     $                          zxtsol,
     2023     $                          zxsnow,
     2024     $                          solsw,
     2025     $                          albsol,
     2026     $                          rain_fall,
     2027     $                          snow_fall,
     2028     $                          itop_con,
     2029     $                          ibas_con,
     2030     $                          cldfra,
     2031     $                          iim,
     2032     $                          jjm,
     2033     $                          tr_seri)
     2034#ifdef INCAINFO
     2035           PRINT *, 'OK.'
     2036#endif
     2037#endif
     2038
    19312039c
    19322040c Calculer les parametres optiques des nuages et quelques
    … …  
    21242232C
    21252233      call phytrac (rnpb,
     2234     I                   itap, julien, gmtime,
    21262235     I                   debut,lafin,
    21272236     I                   nqmax-2,
    21282237     I                   nlon,nlev,dtime,
    2129      I                   t,paprs,pplay,
     2238     I                   u,v,t,paprs,pplay,
    21302239     I                   pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
    21312240     I                   ycoefh,yu1,yv1,ftsol,pctsrf,rlat,
    21322241     I                   frac_impa, frac_nucl,
    2133      I                   rlon,presnivs,paire,pphis,
     2242     I                   rlon,paire,pphis,pphi,
     2243     I                   albsol,
     2244     I                   qx(1,1,1), rhcl,
     2245     I                   cldfra, rneb, diafra, cldliq, itop_con,
     2246     I                   ibas_con,
     2247     I                   pmflxr,pmflxs,prfl,psfl,flxmass_w,
    21342248     O                   tr_seri)
    21352249      ENDIF
    … …  
    28222936cc         CALL ecrirega(84,d_q_lsc)
    28232937      ENDIF
     2938
     2939#ifdef INCA
     2940#ifdef INCAINFO
     2941           PRINT *, 'Appel CHEMHOOK_END ...'
     2942#endif
     2943           CALL chemhook_end (calday,
     2944     $                        dtime,
     2945     $                        pplay,
     2946     $                        t_seri,
     2947     $                        tr_seri,
     2948     $                        nbtr,
     2949     $                        paprs,
     2950     $                        anne_ini,
     2951     $                        day_ini,
     2952     $                        xjour)
     2953#ifdef INCAINFO
     2954           PRINT *, 'OK.'
     2955#endif
     2956#endif
     2957
    28242958c
    28252959c Convertir les incrementations en tendances

  • LMDZ.3.3/trunk/libf/phylmd/phytrac.F

    r209 r344  
    33c
    44      SUBROUTINE phytrac (rnpb,
    5      I                   debutphy,lafin,
    6      I                   nqmax,
    7      I                   nlon,nlev,pdtphys,
    8      I                   t_seri,paprs,pplay,
    9      I                   pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
    10      I                   coefh,yu1,yv1,ftsol,pctsrf,xlat,
    11      I                   frac_impa, frac_nucl,
    12      I                   xlon,presnivs,paire,pphis,
    13      O                   tr_seri)
     5     I                    nstep,
     6     I                    julien,
     7     I                    gmtime,
     8     I                    debutphy,
     9     I                    lafin,
     10     I                    nqmax,
     11     I                    nlon,
     12     I                    nlev,
     13     I                    pdtphys,
     14     I                    u,
     15     I                    v,
     16     I                    t_seri,
     17     I                    paprs,
     18     I                    pplay,
     19     I                    pmfu,
     20     I                    pmfd,
     21     I                    pen_u,
     22     I                    pde_u,
     23     I                    pen_d,
     24     I                    pde_d,
     25     I                    coefh,
     26     I                    yu1,
     27     I                    yv1,
     28     I                    ftsol,
     29     I                    pctsrf,
     30     I                    xlat,
     31     I                    frac_impa,
     32     I                    frac_nucl,
     33     I                    xlon,
     34     I                    paire,
     35     I                    pphis,
     36     I                    pphi,
     37     I                    albsol,
     38     I                    sh,
     39     I                    rh,
     40     I                    cldfra,
     41     I                    rneb,
     42     I                    diafra,
     43     I                    cldliq,
     44     I                    itop_con,
     45     I                    ibas_con,
     46     I                    pmflxr,
     47     I                    pmflxs,
     48     I                    prfl,
     49     I                    psfl,
     50     I                    flxmass_w,
     51     O                    tr_seri)
    1452      USE ioipsl
     53
     54#ifdef INCA
     55      USE sflx
     56      USE chem_tracnm
     57      USE species_names
     58      USE chem_mods
     59      USE pht_tables, ONLY : jrates
     60      USE transport_controls, ONLY : conv_flg, pbl_flg
     61      USE airplane_src, ONLY : ptrop
     62      USE lightning, ONLY : prod_light
     63#ifdef INCA_AER
     64      USE AEROSOL_DIAG
     65#endif
     66#endif
     67
    1568      IMPLICIT none
    1669c======================================================================
    … …  
    2881#include "YOMCST.h"
    2982#include "dimensions.h"
     83#include "paramet.h"
    3084#include "dimphy.h"
    3185#include "indicesol.h"
     86#include "comvert.h"
    3287#include "temps.h"
    3388#include "control.h"
    … …  
    4297c   -------
    4398c
    44       integer nlon  ! nombre de points horizontaux
    45       integer nlev  ! nombre de couches verticales
    46       integer nqmax ! nombre de traceurs auxquels on applique la physique
     99      integer nlon   ! nombre de points horizontaux
     100      integer nlev   ! nombre de couches verticales
     101      integer nqmax  ! nombre de traceurs auxquels on applique la physique
     102      integer nstep  ! appel physique
     103      integer julien !jour julien
     104      integer itop_con(nlon)
     105      integer ibas_con(nlon)
     106      real gmtime
    47107      real pdtphys  ! pas d'integration pour la physique (seconde)
    48108      real t_seri(nlon,nlev) ! temperature
     109      real u(klon,klev)
     110      real v(klon,klev)
     111      real sh(nlon,nlev)     ! humidite specifique
     112      real rh(nlon,nlev)     ! humidite relative
     113      real cldliq(nlon,nlev) ! eau liquide nuageuse
     114      real cldfra(nlon,nlev) ! fraction nuageuse (tous les nuages)
     115      real diafra(nlon,nlev) ! fraction nuageuse (convection ou stratus artificiels)
     116      real rneb(nlon,nlev)   ! fraction nuageuse (grande echelle)
     117      real albsol(nlon)      ! albedo surface
    49118      real tr_seri(nlon,nlev,nbtr) ! traceur 
    50119      real paprs(nlon,nlev+1)  ! pression pour chaque inter-couche (en Pa)
     120      real ps(nlon)  ! pression surface
    51121      real pplay(nlon,nlev)  ! pression pour le mileu de chaque couche (en Pa)
    52       real presnivs(klev) ! pressions approximat. des milieux couches ( en PA)
     122      real pphi(nlon,klev) ! geopotentiel
     123      real znivsig(klev) ! niveaux sigma
    53124      real paire(klon)
    54125      real pphis(klon)
     126      REAL pmflxr(klon,klev+1), pmflxs(klon,klev+1)   !--lessivage convection
     127      REAL prfl(klon,klev+1),   psfl(klon,klev+1)     !--lessivage large-scale
     128      REAL flxmass_w(klon,klev)
    55129      logical debutphy       ! le flag de l'initialisation de la physique
    56130      logical lafin          ! le flag de la fin de la physique
    … …  
    131205      SAVE scavtr
    132206cAA
    133       CHARACTER*1 itn
     207      CHARACTER*2 itn
    134208C maf ioipsl
    135209      CHARACTER*2 str2
    136       INTEGER nhori, nvert
     210      INTEGER nhori, nvert, nverta, nvertb, nverts
    137211      REAL zsto, zout, zjulian
    138212      INTEGER nid_tra
    … …  
    207281      data first,couchelimite,convection,lessivage,sorties
    208282     s     /.true.,.true.,.true.,.true.,.true./
     283
     284#ifdef INCA
     285      INTEGER           :: ncsec
     286      INTEGER           :: prt_flag_ts(nbtr)=(/1,1,1
     287#ifdef INCA_CH4
     288     .                                              ,0,0,1,1,1,1,1,
     289     .                                         0,1,0,0,0,0,0,1,0,0,
     290     .                                         0,1,1,1,1,0,1,1,1,0,
     291     .                                         1,1,1,1,1,1,1,1,1,1,
     292     .                                         1,0,0
     293#endif
     294#ifdef INCA_AER
     295     .                                              ,1,1,1,1,1,1
     296#endif
     297     .                                         /)
     298      REAL, PARAMETER   :: dry_mass = 28.966
     299      REAL, POINTER     :: hbuf(:)
     300      REAL, ALLOCATABLE :: obuf(:)
     301      REAL              :: calday
     302      REAL              :: pdel(klon,klev)
     303      REAL              :: dummy(klon,klev) = 0.
     304#endif
     305
    209306c
    210307c======================================================================
    211308c         ecrit_tra = NINT(86400./pdtphys *1.0)  ! tous les jours
    212309         modname='phytrac'
     310         ecrit_tra = NINT(86400./pdtphys *ecritphy)   
     311!DH      print*,'dans phytrac ',pdtphys,ecritphy,ecrit_tra
     312         ps(:) = paprs(:,1)
    213313
    214314c        print*,'DANS PHYTRAC debutphy=',debutphy
    … …  
    228328
    229329         inirnpb=rnpb
    230          PRINT*, 'La frequence de sortie traceurs est  ', ecrit_tra
    231          itra=0
    232          itap=0
     330!DH      PRINT*, 'La frequence de sortie traceurs est  ', ecrit_tra
     331
    233332C         
    234333         CALL ymds2ju(anne_ini, 1, 1, 0.0, zjulian)
    235          zjulian = zjulian + day_ini
    236 c
     334         itra = NINT(FLOAT(day_ini)*86400./pdtphys)
     335         itap = NINT(FLOAT(day_ini)*86400./pdtphys)
     336!        itra=0
     337!        itap=0
     338!        zjulian = zjulian + day_ini
     339 
    237340         CALL gr_fi_ecrit(1,klon,iim,jjm+1,xlon,zx_lon)
    238341         DO i = 1, iim
    … …  
    240343            zx_lon(i,jjm+1) = xlon(i+1)
    241344         ENDDO
    242 c         DO ll=1,klev
    243 c            znivsig(ll)=float(ll)
    244 c         ENDDO
     345         DO ll=1,klev
     346            znivsig(ll)=float(ll)
     347         ENDDO
    245348         CALL gr_fi_ecrit(1,klon,iim,jjm+1,xlat,zx_lat)
    246349         CALL histbeg("histrac", iim,zx_lon, jjm+1,zx_lat,
    247      .                 1,iim,1,jjm+1, 0, zjulian, pdtphys,
     350     .                 1,iim,1,jjm+1, itra, zjulian, pdtphys,
    248351     .                 nhori, nid_tra)
    249          CALL histvert(nid_tra, "presnivs", "Vertical levels", "mb",
     352
     353         call histvert(nid_tra, "presnivs", "presnivs", "mb",
    250354     .                 klev, presnivs, nvert)
     355#ifdef INCA
     356!        call histvert(nid_tra, "ap", "Hybrid A parameter", "-",
     357!    .                 klev+1, ap, nverta)
     358!        call histvert(nid_tra, "bp", "Hybrid B parameter", "-",
     359!    .                 klev+1, bp, nvertb)
     360#endif
     361
    251362         zsto = pdtphys
    252363         zout = pdtphys * FLOAT(ecrit_tra)
    … …  
    260371     .                "once",  zsto,zout)
    261372
    262         goto 666
    263          CALL histdef(nid_tra, "pyu1", "Vent niv 1", "-",
    264      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    265      .                "inst(X)",  zsto,zout)
    266 
    267          CALL histdef(nid_tra, "pyv1", "Vent niv 1", "-",
    268      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    269      .                "inst(X)",  zsto,zout)
    270          CALL histdef(nid_tra, "psrf1", "nature sol", "-",
    271      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    272      .                "inst(X)",  zsto,zout)
    273          CALL histdef(nid_tra, "psrf2", "nature sol", "-",
    274      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    275      .                "inst(X)",  zsto,zout)
    276          CALL histdef(nid_tra, "psrf3", "nature sol", "-",
    277      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    278      .                "inst(X)",  zsto,zout)
    279          CALL histdef(nid_tra, "psrf4", "nature sol", "-",
    280      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    281      .                "inst(X)",  zsto,zout)
    282          CALL histdef(nid_tra, "ftsol1", "temper sol", "-",
    283      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    284      .                "inst(X)",  zsto,zout)
    285          CALL histdef(nid_tra, "ftsol2", "temper sol", "-",
    286      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    287      .                "inst(X)",  zsto,zout)
    288          CALL histdef(nid_tra, "ftsol3", "temper sol", "-",
    289      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    290      .                "inst",  zsto,zout)
    291          CALL histdef(nid_tra, "ftsol4", "temper sol", "-",
    292      .                iim,jjm+1,nhori, 1,1,1, -99, 32,
    293      .                "inst(X)",  zsto,zout)
    294          CALL histdef(nid_tra, "pplay", "flux u mont","-",
    295      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    296      .                "inst(X)", zsto,zout)
    297          CALL histdef(nid_tra, "t", "flux u mont","-",
    298      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    299      .                "inst(X)", zsto,zout)
    300          CALL histdef(nid_tra, "mfu", "flux u mont","-",
    301      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    302      .                "ave(X)", zsto,zout)
    303          CALL histdef(nid_tra, "mfd", "flux u decen","-",
    304      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    305      .                "ave(X)", zsto,zout)
    306          CALL histdef(nid_tra, "en_u", "flux u mont","-",
    307      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    308      .                "ave(X)", zsto,zout)
    309          CALL histdef(nid_tra, "en_d", "flux u mont","-",
    310      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    311      .                "ave(X)", zsto,zout)
    312          CALL histdef(nid_tra, "de_u", "flux u mont","-",
    313      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    314      .                "ave(X)", zsto,zout)
    315          CALL histdef(nid_tra, "de_d", "flux u mont","-",
    316      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    317      .                "ave(X)", zsto,zout)
    318          CALL histdef(nid_tra, "coefh", "turbulent coef","-",
    319      .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    320      .                "ave(X)", zsto,zout)
    321 
    322 666     continue
    323 c
     373#ifdef INCA
     374         CALL histdef(nid_tra, "ps", "Surface pressure", "Pa",
     375     .                iim,jjm+1,nhori, 1,1,1,-99, 32,
     376     .                "ave(X)", zsto,zout)
     377
     378         CALL histdef(nid_tra, "ptrop", "Tropopause pressure", "Pa",
     379     .                iim,jjm+1,nhori, 1,1,1,-99, 32,
     380     .                "ave(X)", zsto,zout)
     381
     382         CALL histdef(nid_tra, "temp", "Air temperature", "K",
     383     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     384     .                "ave(X)", zsto,zout)
     385
     386         CALL histdef(nid_tra, "u", "zonal wind component", "m/s",
     387     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     388     .                "ave(X)", zsto,zout)
     389
     390         CALL histdef(nid_tra, "v", "zonal wind component", "m/s",
     391     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     392     .                "ave(X)", zsto,zout)
     393
     394         CALL histdef(nid_tra, "h2o", "Specific Humidity", "MMR",
     395     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     396     .                "ave(X)", zsto,zout)
     397
     398         CALL histdef(nid_tra, "pmid", "Pressure", "Pa",
     399     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     400     .                "ave(X)", zsto,zout)
     401
     402         CALL histdef(nid_tra, "pdel", "Delta Pressure", "Pa",
     403     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     404     .                "ave(X)", zsto,zout)
     405#ifdef INCA_CH4
     406#ifdef INCAINFO
     407         DO it=1, phtcnt
     408         WRITE(str2,'(i2.2)') it
     409         CALL histdef(nid_tra, "j"//str2,"j"//str2, "CM-3 S-1",
     410     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     411     .                "ave(X)", zsto,zout)
     412         ENDDO
     413         DO it=1, hetcnt
     414         WRITE(str2,'(i2.2)') it
     415         CALL histdef(nid_tra, "w"//str2,"w"//str2, "S-1",
     416     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     417     .                "ave(X)", zsto,zout)
     418         ENDDO
     419
     420         DO it=1, extcnt
     421         WRITE(str2,'(i2.2)') it
     422        CALL histdef(nid_tra, "ext"//str2,"ext"//str2, "CM-3 S-1",
     423     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     424     .                "ave(X)", zsto,zout)
     425         ENDDO
     426
     427         DO it=1, nfs
     428         WRITE(str2,'(i2.2)') it
     429         CALL histdef(nid_tra, "INV"//str2, "INV"//str2, "CM-3",
     430     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     431     .                "ave(X)", zsto,zout)
     432         ENDDO
     433
     434#else
     435         CALL histdef(nid_tra, "jO3","jO3", "CM-3 S-1",
     436     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     437     .                "ave(X)", zsto,zout)
     438         CALL histdef(nid_tra, "jNO2","jNO2", "CM-3 S-1",
     439     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     440     .                "ave(X)", zsto,zout)
     441         CALL histdef(nid_tra, "jH2O2","jH2O2", "CM-3 S-1",
     442     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     443     .                "ave(X)", zsto,zout)
     444         CALL histdef(nid_tra, "wHNO3","wHNO3", "S-1",
     445     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     446     .                "ave(X)", zsto,zout)
     447         CALL histdef(nid_tra, "kN2O5", "kN2O5","CM-3 S-1",
     448     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     449     .                "ave(X)", zsto,zout)
     450         CALL histdef(nid_tra, "LghtNO","LghtNO", "CM-3 S-1",
     451     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     452     .                "ave(X)", zsto,zout)
     453#endif
     454
     455         DO it=1, grpcnt
     456         CALL histdef(nid_tra, grpsym(it), grpsym(it), "VMR",
     457     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     458     .                "ave(X)", zsto,zout)
     459         ENDDO
     460#endif
     461#endif
     462
     463#ifdef INCA_AER
     464C   3d FIELDS
     465
     466        CALL histdef(nid_tra, "scavcoef_st","scavcoef_st", "S-1",
     467     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     468     .                "ave(X)", zsto,zout)
     469        CALL histdef(nid_tra, "scavcoef_cv","scavcoef_cv", "S-1",
     470     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     471     .                "ave(X)", zsto,zout)
     472#endif
     473 
    324474         DO it=1,nqmax
    325475         IF (it.LE.99) THEN
    326476         WRITE(str2,'(i2.2)') it
    327 C champ 3D
     477C champ 2D
     478
     479#ifdef INCA
     480         IF ( prt_flag_ts(it) == 0 ) CYCLE
     481
     482         CALL histdef(nid_tra, "Emi_"//solsym(it), "Emi_"//solsym(it),
     483     .           "kg/m2/s", iim,jjm+1,nhori, 1,1,1, -99, 32,
     484     .           "ave(X)", zsto,zout)
     485         CALL histdef(nid_tra, "Dep_"//solsym(it), "Dep_"//solsym(it),
     486     .           "cm/s", iim,jjm+1,nhori, 1,1,1, -99, 32,
     487     .           "ave(X)", zsto,zout)
     488#ifdef INCA_AER
     489        IF (solsym(it) .eq. 'CIDUSTM'.or.solsym(it) .eq. 'CIN'
     490     .  .or.solsym(it) .eq. 'CSSSM'  .or.solsym(it) .eq. 'CSN'
     491     .  .or.solsym(it) .eq. 'ASSSM'  .or.solsym(it) .eq. 'ASN' ) THEN
     492         CALL histdef(nid_tra, "Sed_"//solsym(it), "Sed_"//solsym(it),
     493     .           "kg/m2", iim,jjm+1,nhori, 1,1,1, -99, 32,
     494     .           "ave(X)", zsto,zout)
     495        endif
     496        IF (solsym(it) .eq. 'CIDUSTM') THEN
     497         CALL histdef(nid_tra, "OD_"//solsym(it), "OD_"//solsym(it),
     498     .           "opt. depth", iim,jjm+1,nhori, 1,1,1, -99, 32,
     499     .           "ave(X)", zsto,zout)
     500        endif
     501
     502#endif
     503         CALL histdef(nid_tra, solsym(it), solsym(it), "VMR",
     504     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     505     .                "ave(X)", zsto,zout)
     506#else
    328507         CALL histdef(nid_tra, "tr"//str2, "Tracer No."//str2, "U/kga",
    329508     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
    … …  
    334513     .                "ave(X)", zsto,zout)
    335514         ENDIF
     515#endif
    336516         ELSE
    337517         PRINT*, "Trop de traceurs"
    … …  
    339519         ENDIF
    340520         ENDDO
     521
     522#ifdef INCA_CH4
     523         CALL histdef(nid_tra, "O3_column", "O3_column",
     524     .           "DU", iim,jjm+1,nhori, 1,1,1, -99, 32,
     525     .           "ave(X)", zsto,zout)
     526         CALL histdef(nid_tra, "CO_column", "CO_column",
     527     .           "10^18 CM-2", iim,jjm+1,nhori, 1,1,1, -99, 32,
     528     .           "ave(X)", zsto,zout)
     529         CALL histdef(nid_tra, "CH4_column", "CH4_column",
     530     .           "10^18 CM-2", iim,jjm+1,nhori, 1,1,1, -99, 32,
     531     .           "ave(X)", zsto,zout)
     532         CALL histdef(nid_tra, "NO2_column", "NO2_column",
     533     .           "10^15 CM-2", iim,jjm+1,nhori, 1,1,1, -99, 32,
     534     .           "ave(X)", zsto,zout)
     535         CALL histdef(nid_tra, "O3_ste", "O3_ste",
     536     .           "CM-2 S-1", iim,jjm+1,nhori, 1,1,1, -99, 32,
     537     .           "ave(X)", zsto,zout)
     538         CALL histdef(nid_tra, "O3_prod", "O3_prod", "CM-3 S-1",
     539     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     540     .                "ave(X)", zsto,zout)
     541         CALL histdef(nid_tra, "O3_loss", "O3_loss", "CM-3 S-1",
     542     .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     543     .                "ave(X)", zsto,zout)
     544
     545!        Special variables for daytime averaging
     546!        CALL histdef(nid_tra, "day_cnt", "day_cnt", "-",
     547!    .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     548!    .                "t_sum(X)", zsto,zout)
     549!        CALL histdef(nid_tra, "NO_day", "NO_day", "VMR",
     550!    .                iim,jjm+1,nhori, klev,1,klev,nvert, 32,
     551!    .                "t_sum(X)", zsto,zout)
     552
     553#endif
     554
    341555         CALL histend(nid_tra)
    342556         ndex(1) = 0
    … …  
    367581         read(99,*) (trs(i,1),i=1,klon)
    368582999      close(99)
    369          print*, 'apres starttrac'
     583!DH      print*, 'apres starttrac'
    370584
    371585c Initialisation de la fraction d'aerosols lessivee
    … …  
    402616         inirnpb=.false.
    403617      endif
    404       if(nqmax.gt.2) aerosol(3)=.true.
    405 
    406 
    407 c  abder
    408         goto 777
    409             do i=1,nlon
    410                pftsol1(i) = ftsol(i,1)
    411                pftsol2(i) = ftsol(i,2)
    412                pftsol3(i) = ftsol(i,3)
    413                pftsol4(i) = ftsol(i,4)
    414 
    415                ppsrf1(i) = pctsrf(i,1)
    416                ppsrf2(i) = pctsrf(i,2)
    417                ppsrf3(i) = pctsrf(i,3)
    418                ppsrf4(i) = pctsrf(i,4)
    419 
    420             enddo
    421          ndex(1)=0
    422          itap=itap+1
    423          CALL gr_fi_ecrit(1,klon,iim,jjm+1,yu1,zx_tmp_2d)
    424          CALL histwrite(nid_tra,"pyu1",itap,zx_tmp_2d,
    425      s                                  iim*(jjm+1),ndex)
    426          
    427          CALL gr_fi_ecrit(1,klon,iim,jjm+1,yv1,zx_tmp_2d)
    428          CALL histwrite(nid_tra,"pyv1",itap,zx_tmp_2d,
    429      s                                  iim*(jjm+1),ndex)
    430 
    431          CALL gr_fi_ecrit(1,klon,iim,jjm+1,pftsol1,zx_tmp_2d)
    432          CALL histwrite(nid_tra,"ftsol1",itap,zx_tmp_2d,
    433      s                                       iim*(jjm+1),ndex)
    434 
    435          CALL gr_fi_ecrit(1,klon,iim,jjm+1,pftsol2,zx_tmp_2d)
    436          CALL histwrite(nid_tra,"ftsol2",itap,zx_tmp_2d,
    437      s                                       iim*(jjm+1),ndex)
    438 
    439          CALL gr_fi_ecrit(1,klon,iim,jjm+1,pftsol3,zx_tmp_2d)
    440          CALL histwrite(nid_tra,"ftsol3",itap,zx_tmp_2d,
    441      s                                      iim*(jjm+1),ndex)
    442 
    443          CALL gr_fi_ecrit(1,klon,iim,jjm+1,pftsol4,zx_tmp_2d)
    444          CALL histwrite(nid_tra,"ftsol4",itap,zx_tmp_2d,
    445      s                                      iim*(jjm+1),ndex)
    446 
    447          CALL gr_fi_ecrit(1,klon,iim,jjm+1,ppsrf1,zx_tmp_2d)
    448          CALL histwrite(nid_tra,"psrf1",itap,zx_tmp_2d,
    449      s                                     iim*(jjm+1),ndex)
    450 
    451          CALL gr_fi_ecrit(1,klon,iim,jjm+1,ppsrf2,zx_tmp_2d)
    452          CALL histwrite(nid_tra,"psrf2",itap,zx_tmp_2d,
    453      s                                     iim*(jjm+1),ndex)
    454 
    455          CALL gr_fi_ecrit(1,klon,iim,jjm+1,ppsrf3,zx_tmp_2d)
    456          CALL histwrite(nid_tra,"psrf3",itap,zx_tmp_2d,
    457      s                                     iim*(jjm+1),ndex)
    458 
    459          CALL gr_fi_ecrit(1,klon,iim,jjm+1,ppsrf4,zx_tmp_2d)
    460          CALL histwrite(nid_tra,"psrf4",itap,zx_tmp_2d,
    461      s                                     iim*(jjm+1),ndex)
    462 777     continue
     618
     619#ifdef INCA
     620!======================================================================
     621!     Chimie
     622!======================================================================
     623
     624        calday = FLOAT(julien) + gmtime
     625        ncsec  = NINT (86400.*gmtime)
     626
     627        DO k = 1, nlev
     628        pdel(:,k) = paprs(:,k) - paprs (:,k+1)
     629        END DO
     630
     631#ifdef INCAINFO
     632        PRINT *, 'CHEMMAIN @ ', calday, ' ... '
     633        DO it = 1, nbtr
     634        PRINT *, solsym(it), MINVAL(tr_seri(:,:,it)),
     635     $                       MAXVAL(tr_seri(:,:,it))
     636      END DO
     637#endif
     638
     639
     640#ifdef INCA_AER
     641        CALL aerosolmain (tr_seri,
     642     $                 pdtphys,
     643     $                 pplay,
     644     $                 prfl,
     645     $                 pmflxr,
     646     $                 psfl,
     647     $                 pmflxs,
     648     $                 pmfu,
     649     $                 itop_con,
     650     $                 ibas_con,
     651     $                 pphi,
     652     $                 paire,
     653     $                 nstep)
     654#endif
     655
     656        CALL chemmain (tr_seri,
     657     $                 nas,
     658     $                 nstep,
     659     $                 calday,
     660     $                 julien,
     661     $                 ncsec,
     662     $                 1,
     663     $                 pdtphys,
     664     $                 paprs(1,1),
     665     $                 pplay,
     666     $                 pdel,
     667     $                 pctsrf(1,3),
     668     $                 ftsol,
     669     $                 albsol,
     670     $                 pphi,
     671     $                 pphis,
     672     $                 cldfra,
     673     $                 rneb,
     674     $                 diafra,
     675     $                 itop_con,
     676     $                 cldliq,
     677     $                 prfl,
     678     $                 pmflxr,
     679     $                 psfl,
     680     $                 pmflxs,
     681     $                 pmfu,
     682     $                 flxmass_w,
     683     $                 t_seri,
     684     $                 sh,
     685     $                 rh,
     686     $                 .false.,
     687     $                 hbuf,
     688     $                 obuf,
     689     $                 iip1,
     690     $                 jjp1)
     691#ifdef INCAINFO
     692      PRINT *, 'OK.'
     693      DO it = 1, nbtr
     694      PRINT *, solsym(it), MINVAL(tr_seri(:,:,it)),
     695     $                     MAXVAL(tr_seri(:,:,it))
     696      END DO
     697#endif
     698#endif
     699C
     700
    463701c======================================================================
    464702c   Calcul de l'effet de la convection
    465703c======================================================================
    466         print*,'Avant convection'
     704!DH   print*,'Avant convection'
    467705      do it=1,nqmax
    468706         WRITE(itn,'(i1)') it
    … …  
    472710      if (convection) then
    473711
    474       print*,'Pas de temps dans phytrac : ',pdtphys
     712!DH   print*,'Pas de temps dans phytrac : ',pdtphys
    475713      DO it=1, nqmax
     714#ifdef INCA
     715      IF ( conv_flg(it) == 0 ) CYCLE
     716#endif
    476717      CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d,
    477718     .            pplay, paprs, tr_seri(1,1,it), d_tr_cv)
    … …  
    481722      ENDDO
    482723      ENDDO
    483 c      WRITE(itn,'(i1)') it
    484 c      CALL minmaxqfi(tr_seri(1,1,it),0.,1.e33,'convection it='//itn)
    485       ENDDO
    486 c      print*,'apres nflxtr'
    487 
     724c     WRITE(itn,'(i2)') it
     725#ifdef INCA
     726      CALL minmaxqfi(tr_seri(1,1,it),0.,1.e33,'convection it = '//solsym(it))
     727#else
     728      CALL minmaxqfi(tr_seri(1,1,it),0.,1.e33,'convection it = '//itn)
     729#endif
     730      ENDDO
     731c     print*,'apres nflxtr'
    488732
    489733      endif ! convection
    … …  
    513757C maf modif pour tenir compte du cas rnpb + traceur
    514758      DO it=1, nqmax
    515 c     print *,'it',it,clsol(it)
     759#ifdef INCA
     760      IF ( pbl_flg(it) == 0 ) CYCLE
     761#endif
     762C     print *,'it',it,clsol(it)
    516763      if (clsol(it)) then  ! couche limite avec quantite dans le sol calculee
    517764          CALL cltracrn(it, pdtphys, yu1, yv1,
    … …  
    539786C         CALL minmaxqfi(tr_seri(1,1,it),0.,1.e33,'cltracrn it='//itn)
    540787      else ! couche limite avec flux prescrit
     788
     789#ifdef INCA
     790        DO k =  1, klon
     791          source(k) = eflux(k,it)-dflux(k,it)
     792        END DO
     793#else
    541794Cmaf provisoire source / traceur a creer
    542795        DO i=1, klon
    543796          source(i) = 0.0 ! pas de source, pour l'instant
    544797        ENDDO
    545 C
     798#endif
    546799          CALL cltrac(pdtphys, coefh,t_seri,
    547800     s               tr_seri(1,1,it), source,
    … …  
    574827c si radio=true mais pour l'instant radiornpb propre au cas rnpb
    575828      if(rnpb) then
     829        print *, 'decroissance radiactive activee'
    576830        call radiornpb (tr_seri,pdtphys,tautr,d_tr_dec)
    577831C
    … …  
    589843c
    590844      endif ! rnpb decroissance  radioactive
    591 C
     845
    592846c======================================================================
    593847c   Calcul de l'effet de la precipitation
    594848c======================================================================
    595849
    596       print*,'LESSIVAGE =',lessivage
     850!DH   print*,'LESSIVAGE =',lessivage
    597851      IF (lessivage) THEN
    598852
    … …  
    675929         ENDDO
    676930      ENDDO
     931
     932
     933
    677934      itra=itra+1
    678935      ndex(1) = 0
     936
     937      i = NINT(zout/zsto)
     938      CALL gr_fi_ecrit(1, klon,iim,jjm+1, paire,zx_tmp_2d)
     939      CALL histwrite(nid_tra,"aire",itra,zx_tmp_2d,
     940     .     iim*(jjm+1),ndex)
     941
     942#ifdef INCA
     943      CALL gr_fi_ecrit(1, klon,iim,jjm+1, ps,zx_tmp_2d)
     944      CALL histwrite(nid_tra,"ps",itra,zx_tmp_2d,
     945     .     iim*(jjm+1),ndex)
     946
     947      CALL gr_fi_ecrit(1, klon,iim,jjm+1, ptrop,zx_tmp_2d)
     948      CALL histwrite(nid_tra,"ptrop",itra,zx_tmp_2d,
     949     .     iim*(jjm+1),ndex)
     950
     951      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,t_seri, zx_tmp_3d)
     952      CALL histwrite(nid_tra,"temp",itra,zx_tmp_3d,
     953     .                                   iim*(jjm+1)*klev,ndex)
     954
     955      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,u, zx_tmp_3d)
     956      CALL histwrite(nid_tra,"u",itra,zx_tmp_3d,
     957     .                                   iim*(jjm+1)*klev,ndex)
     958
     959      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,v, zx_tmp_3d)
     960      CALL histwrite(nid_tra,"v",itra,zx_tmp_3d,
     961     .                                   iim*(jjm+1)*klev,ndex)
     962
     963      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,sh, zx_tmp_3d)
     964      CALL histwrite(nid_tra,"h2o",itra,zx_tmp_3d,
     965     .                                   iim*(jjm+1)*klev,ndex)
     966
     967      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pdel, zx_tmp_3d)
     968      CALL histwrite(nid_tra,"pdel",itra,zx_tmp_3d,
     969     .                                   iim*(jjm+1)*klev,ndex)
     970
     971      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pplay, zx_tmp_3d)
     972      CALL histwrite(nid_tra,"pmid",itra,zx_tmp_3d,
     973     .                                   iim*(jjm+1)*klev,ndex)
     974
     975#ifdef INCA_CH4
     976#ifdef INCAINFO
     977      DO it=1, phtcnt
     978      WRITE(str2,'(i2.2)') it
     979      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,jrates(1,1,it),
     980     .     zx_tmp_3d)
     981      CALL histwrite(nid_tra,"j"//str2,itra,zx_tmp_3d,
     982     .                                   iim*(jjm+1)*klev,ndex)
     983      ENDDO
     984
     985      DO it=1, hetcnt
     986      WRITE(str2,'(i2.2)') it
     987      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,hrates(1,1,it),
     988     .     zx_tmp_3d)
     989      CALL histwrite(nid_tra,"w"//str2,itra,zx_tmp_3d,
     990     .                                   iim*(jjm+1)*klev,ndex)
     991      ENDDO
     992
     993 
     994      DO it=1, extcnt
     995      WRITE(str2,'(i2.2)') it
     996
     997
     998      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,extflx(1,1,it),
     999     .     zx_tmp_3d)
     1000      CALL histwrite(nid_tra,"ext"//str2,itra,zx_tmp_3d,
     1001     .                                   iim*(jjm+1)*klev,ndex)
     1002      ENDDO
     1003
     1004      DO it=1, nfs
     1005      WRITE(str2,'(i2.2)') it
     1006      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,invariants(1,1,it),
     1007     .     zx_tmp_3d)
     1008      CALL histwrite(nid_tra,"INV"//str2,itra,zx_tmp_3d,
     1009     .                                   iim*(jjm+1)*klev,ndex)
     1010      ENDDO
     1011
     1012
     1013#else
     1014      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,jrates(1,1,2),
     1015     .     zx_tmp_3d)
     1016      CALL histwrite(nid_tra,"jO3",itra,zx_tmp_3d,
     1017     .                                   iim*(jjm+1)*klev,ndex)
     1018
     1019      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,jrates(1,1,4),
     1020     .     zx_tmp_3d)
     1021      CALL histwrite(nid_tra,"jNO2",itra,zx_tmp_3d,
     1022     .                                   iim*(jjm+1)*klev,ndex)
     1023
     1024      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,jrates(1,1,13),
     1025     .     zx_tmp_3d)
     1026      CALL histwrite(nid_tra,"jH2O2",itra,zx_tmp_3d,
     1027     .                                   iim*(jjm+1)*klev,ndex)
     1028
     1029      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,hrates(1,1,1),
     1030     .     zx_tmp_3d)
     1031      CALL histwrite(nid_tra,"wHNO3",itra,zx_tmp_3d,
     1032     .                                   iim*(jjm+1)*klev,ndex)
     1033
     1034      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,krates(1,1,1),
     1035     .     zx_tmp_3d)
     1036      CALL histwrite(nid_tra,"kN2O5",itra,zx_tmp_3d,
     1037     .                                   iim*(jjm+1)*klev,ndex)
     1038
     1039      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,extflx(1,1,1),
     1040     .     zx_tmp_3d)
     1041      CALL histwrite(nid_tra,"LghtNO",itra,zx_tmp_3d,
     1042     .                                   iim*(jjm+1)*klev,ndex)
     1043#endif
     1044      DO it=1, grpcnt
     1045      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,nas(1,1,it),zx_tmp_3d)
     1046      zx_tmp_3d = zx_tmp_3d * dry_mass / nadv_mass(it)
     1047      CALL histwrite(nid_tra,grpsym(it),itra,zx_tmp_3d,
     1048     .                                   iim*(jjm+1)*klev,ndex)
     1049      ENDDO
     1050#endif
     1051
     1052#endif
     1053#ifdef INCA_AER
     1054C   3d FIELDS
     1055
     1056      it = id_CIDUSTM
     1057       CALL gr_fi_ecrit(klev,klon,iim,jjm+1,scavcoef_st(1,1,it),
     1058     .                  zx_tmp_3d)
     1059       CALL histwrite(nid_tra,"scavcoef_st",itra,zx_tmp_3d,
     1060     .                  iim*(jjm+1)*klev,ndex)
     1061
     1062       CALL gr_fi_ecrit(klev,klon,iim,jjm+1,scavcoef_cv(1,1,it),
     1063     .                  zx_tmp_3d)
     1064       CALL histwrite(nid_tra,"scavcoef_cv",itra,zx_tmp_3d,
     1065     .                  iim*(jjm+1)*klev,ndex)
     1066
     1067#endif
     1068
    6791069      DO it=1,nqmax
    6801070      IF (it.LE.99) THEN
    6811071       WRITE(str2,'(i2.2)') it
     1072
     1073#ifdef INCA
     1074      IF ( prt_flag_ts(it) == 0 ) CYCLE
     1075C champs 2D
     1076      CALL gr_fi_ecrit(1, klon,iim,jjm+1, eflux(1,it),zx_tmp_2d)
     1077      CALL histwrite(nid_tra,"Emi_"//solsym(it),itra,zx_tmp_2d,
     1078     .     iim*(jjm+1),ndex)
     1079
     1080      CALL gr_fi_ecrit(1, klon,iim,jjm+1, dvel(1,it),zx_tmp_2d)
     1081      CALL histwrite(nid_tra,"Dep_"//solsym(it),itra,zx_tmp_2d,
     1082     .     iim*(jjm+1),ndex)
     1083#ifdef INCA_AER
     1084        IF (solsym(it) .eq. 'CIDUSTM'.or.solsym(it) .eq. 'CIN'
     1085     .  .or.solsym(it) .eq. 'CSSSM'  .or.solsym(it) .eq. 'CSN'
     1086     .  .or.solsym(it) .eq. 'ASSSM'  .or.solsym(it) .eq. 'ASN' ) THEN
     1087      CALL gr_fi_ecrit(1, klon,iim,jjm+1,sflux(1,it),zx_tmp_2d)
     1088      CALL histwrite(nid_tra,"Sed_"//solsym(it),itra,zx_tmp_2d,
     1089     .     iim*(jjm+1),ndex)
     1090        endif
     1091        IF (solsym(it) .eq. 'CIDUSTM') THEN
     1092      CALL gr_fi_ecrit(1, klon,iim,jjm+1,tausum(1),zx_tmp_2d)
     1093      CALL histwrite(nid_tra,"OD_"//solsym(it),itra,zx_tmp_2d,
     1094     .     iim*(jjm+1),ndex)
     1095        endif
     1096
     1097#endif
     1098C champs 3D
     1099       CALL gr_fi_ecrit(klev,klon,iim,jjm+1,tr_seri(1,1,it),zx_tmp_3d)
     1100
     1101       !Prefer vmr to mmr for transported species
     1102       if( adv_mass(it) /= 0. ) then
     1103       zx_tmp_3d = zx_tmp_3d * dry_mass / adv_mass(it)
     1104       else
     1105#ifdef INCA_CH4
     1106       if ( solsym(it) == 'OX' ) then
     1107       zx_tmp_3d = zx_tmp_3d * dry_mass / nadv_mass(id_o3)
     1108       end if
     1109#endif
     1110       end if
     1111
     1112       CALL histwrite(nid_tra,solsym(it),itra,zx_tmp_3d,
     1113     .                                   iim*(jjm+1)*klev,ndex)
     1114#else
     1115
    6821116C champs 3D
    6831117       CALL gr_fi_ecrit(klev,klon,iim,jjm+1,tr_seri(1,1,it),zx_tmp_3d)
    6841118       CALL histwrite(nid_tra,"tr"//str2,itra,zx_tmp_3d,
    6851119     .                                   iim*(jjm+1)*klev,ndex)
    686 c      IF (lessivage) THEN
    687 c      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,flestottr(1,1,it),zx_tmp_3d)
    688 c      CALL histwrite(nid_tra,"fl"//str2,itra,zx_tmp_3d,
    689 c    .                                   iim*(jjm+1)*klev,ndex)
    690 c     ENDIF
     1120       IF (lessivage) THEN
     1121       CALL gr_fi_ecrit(klev,klon,iim,jjm+1,flestottr(1,1,it),zx_tmp_3d)
     1122       CALL histwrite(nid_tra,"fl"//str2,itra,zx_tmp_3d,
     1123     .                                   iim*(jjm+1)*klev,ndex)
     1124       ENDIF
     1125#endif
    6911126      ELSE
    6921127         PRINT*, "Trop de traceurs"
    … …  
    6951130      ENDDO
    6961131
    697         goto 888
    698         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pplay,zx_tmp_3d)
    699         CALL histwrite(nid_tra,"pplay",itra,zx_tmp_3d,
    700      .                  iim*(jjm+1)*klev,ndex)
    701 
    702         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,t_seri,zx_tmp_3d)
    703         CALL histwrite(nid_tra,"t",itra,zx_tmp_3d,
    704      .                  iim*(jjm+1)*klev,ndex)
    705         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pmfu,zx_tmp_3d)
    706         CALL histwrite(nid_tra,"mfu",itra,zx_tmp_3d,
    707      .                  iim*(jjm+1)*klev,ndex)
    708         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pmfd,zx_tmp_3d)
    709         CALL histwrite(nid_tra,"mfd",itra,zx_tmp_3d,
    710      .                  iim*(jjm+1)*klev,ndex)
    711         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pen_u,zx_tmp_3d)
    712         CALL histwrite(nid_tra,"en_u",itra,zx_tmp_3d,
    713      .                  iim*(jjm+1)*klev,ndex)
    714         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pen_d,zx_tmp_3d)
    715         CALL histwrite(nid_tra,"en_d",itra,zx_tmp_3d,
    716      .                  iim*(jjm+1)*klev,ndex)
    717         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pde_d,zx_tmp_3d)
    718         CALL histwrite(nid_tra,"de_d",itra,zx_tmp_3d,
    719      .                  iim*(jjm+1)*klev,ndex)
    720         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,pde_u,zx_tmp_3d)
    721         CALL histwrite(nid_tra,"de_u",itra,zx_tmp_3d,
    722      .                  iim*(jjm+1)*klev,ndex)
    723         CALL gr_fi_ecrit(klev,klon,iim,jjm+1,coefh,zx_tmp_3d)
    724         CALL histwrite(nid_tra,"coefh",itra,zx_tmp_3d,
    725      .                  iim*(jjm+1)*klev,ndex)
    726 
    727 888     continue
    728 
    729 c       print*,'Sortie phytrac'
    730 c      do it=1,nqmax
    731 c         WRITE(itn,'(i1)') it
    732 c        call diagtracphy(tr_seri(:,:,it),paprs,'Fin Phys  '//itn)
    733 c      enddo
     1132#ifdef INCA_CH4
     1133      CALL gr_fi_ecrit(1, klon,iim,jjm+1, o3_tr_col(1), zx_tmp_2d)
     1134      CALL histwrite(nid_tra,"O3_column",itra,zx_tmp_2d,
     1135     .     iim*(jjm+1),ndex)
     1136
     1137      CALL gr_fi_ecrit(1, klon,iim,jjm+1, co_tr_col(1), zx_tmp_2d)
     1138      CALL histwrite(nid_tra,"CO_column",itra,zx_tmp_2d,
     1139     .     iim*(jjm+1),ndex)
     1140
     1141      CALL gr_fi_ecrit(1, klon,iim,jjm+1, ch4_tr_col(1), zx_tmp_2d)
     1142      CALL histwrite(nid_tra,"CH4_column",itra,zx_tmp_2d,
     1143     .     iim*(jjm+1),ndex)
     1144
     1145      CALL gr_fi_ecrit(1, klon,iim,jjm+1, no2_tr_col(1), zx_tmp_2d)
     1146      CALL histwrite(nid_tra,"NO2_column",itra,zx_tmp_2d,
     1147     .     iim*(jjm+1),ndex)
     1148
     1149      CALL gr_fi_ecrit(1, klon,iim,jjm+1, o3_st_flx(1), zx_tmp_2d)
     1150      CALL histwrite(nid_tra,"O3_ste",itra,zx_tmp_2d,
     1151     .     iim*(jjm+1),ndex)
     1152
     1153      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,o3_prod(1,1),
     1154     .     zx_tmp_3d)
     1155      CALL histwrite(nid_tra,"O3_prod",itra,zx_tmp_3d,
     1156     .                                   iim*(jjm+1)*klev,ndex)
     1157
     1158      CALL gr_fi_ecrit(klev,klon,iim,jjm+1,o3_loss(1,1),
     1159     .     zx_tmp_3d)
     1160      CALL histwrite(nid_tra,"O3_loss",itra,zx_tmp_3d,
     1161     .                                   iim*(jjm+1)*klev,ndex)
     1162
     1163!     ... Special section for daytime averaging
     1164!       CALL gr_fi_ecrit(klev,klon,iim,jjm+1,day_cnt(1,1),
     1165!    .       zx_tmp_3d)
     1166!       CALL histwrite(nid_tra,"day_cnt",itra,zx_tmp_3d,
     1167!    .                                  iim*(jjm+1)*klev,ndex)
     1168!       CALL gr_fi_ecrit(klev,klon,iim,jjm+1,no_daytime(1,1),
     1169!    .       zx_tmp_3d)
     1170!       CALL histwrite(nid_tra,"NO_day",itra,zx_tmp_3d,
     1171!    .                                  iim*(jjm+1)*klev,ndex)
     1172
     1173#endif
     1174
     1175      if (ok_sync) call histsync(nid_tra)
    7341176
    7351177      if (lafin) then
    736          print*, 'c est la fin de la physique'
     1178!DH      print*, 'c est la fin de la physique'
    7371179         open (99,file='restarttrac',  form='formatted')
    7381180         do i=1,klon
    … …  
    7421184         close(99)
    7431185      else
    744          print*, 'physique pas fini'
     1186!DH      print*, 'physique pas fini'
    7451187      endif
    7461188
Note: See TracChangeset for help on using the changeset viewer.