Changeset 1160 for trunk


Ignore:
Timestamp:
Jan 16, 2014, 9:49:42 AM (11 years ago)
Author:
slebonnois
Message:

SL: update for tracers management in phytrac, Venus

Location:
trunk/LMDZ.VENUS
Files:
1 added
4 edited
1 moved

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.VENUS/deftank/physiq.def

    r888 r1160  
    6262solaire = 328.
    6363#
     64# parameters for the tracers
     65#
     66# 0: nothing (passive tracers)
     67# 1: pseudo-chemistry relaxation (phytrac_relax)
     68# 2: surface emission (phytrac_emiss)
     69# 3: full chemistry (phytrac_chem)
     70tr_scheme = 0
     71#
    6472# parametres boundary layer
    6573#

  • trunk/LMDZ.VENUS/libf/phyvenus/clesphys.h

    r1120 r1160  
    1010       LOGICAL cycle_diurne,soil_model
    1111       LOGICAL ok_orodr,ok_orolf,ok_gw_nonoro
     12       LOGICAL ok_kzmin
    1213       INTEGER nbapp_rad, nbapp_chim, iflag_con, iflag_ajs
     14       INTEGER lev_histhf, lev_histday, lev_histmth
     15       INTEGER tr_scheme
    1316       REAL    ecriphy
    1417       REAL    solaire
    1518       REAL    z0, lmixmin
    1619       REAL    ksta, inertie
    17        LOGICAL ok_kzmin
    18        INTEGER lev_histhf, lev_histday, lev_histmth
    1920
    20        COMMON/clesphys/cycle_diurne, soil_model,                        &
    21      &     ok_orodr, ok_orolf, ok_gw_nonoro, nbapp_rad, nbapp_chim      &
    22      &     , ecriphy                                                    &
    23      &     , iflag_con, iflag_ajs, solaire, z0, lmixmin, ksta           &
    24      &     , ok_kzmin, lev_histhf, lev_histday, lev_histmth             &
    25      &     , inertie
     21       COMMON/clesphys_l/cycle_diurne, soil_model,                      &
     22     &     ok_orodr, ok_orolf, ok_gw_nonoro, ok_kzmin
    2623
     24       COMMON/clesphys_i/nbapp_rad, nbapp_chim,                         &
     25     &     iflag_con, iflag_ajs,                                        &
     26     &     lev_histhf, lev_histday, lev_histmth, tr_scheme
     27
     28       COMMON/clesphys_r/ecriphy, solaire, z0, lmixmin,                 &
     29     &     ksta, inertie
     30

  • trunk/LMDZ.VENUS/libf/phyvenus/conf_phys.F90

    r1120 r1160  
    190190    call getin('solaire', solaire)
    191191!
     192!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     193! PARAMETER FOR THE TRACERS
     194!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     195!
     196!Config Key  = tr_scheme
     197!Config Desc =
     198!Config Def  = 0
     199!Config Help =
     200!
     201! 0   = Nothing is done (passive tracers)
     202! 1   = pseudo-chemistry with relaxation toward fixed profile
     203!       See Marcq&Lebonnois 2013
     204! 2   = surface emission
     205!       For the moment, inspired from Mars version
     206!       However, the variable 'source' could be used in physiq
     207!       so the call to phytrac_emiss could be to initialise it.
     208! 3   = Full chemistry
     209!       To be added by Aurelien
     210  tr_scheme = 0
     211  call getin('tr_scheme',tr_scheme)
     212
    192213!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    193214! PARAMETER FOR THE PLANETARY BOUNDARY LAYER AND SOIL
    … …  
    323344  write(numout,*)' Equinoxe = ',R_peri
    324345  write(numout,*)' Inclinaison =',R_incl
     346  write(numout,*)' tr_scheme = ', tr_scheme
    325347  write(numout,*)' iflag_pbl = ', iflag_pbl
    326348  write(numout,*)' z0 = ',z0

  • trunk/LMDZ.VENUS/libf/phyvenus/physiq.F

    r1120 r1160  
    669669
    670670      if (iflag_trac.eq.1) then
    671       call phytrac (
    672      I                   itap, gmtime,
    673      I                   debut,lafin,
    674      I                   nqmax,
    675      I                   nlon,nlev,dtime,
    676      I                   u,v,t,paprs,pplay,
    677      I                   rlatd,
    678      I                   rlond,presnivs,pphis,pphi,
    679      I                   falbe,
     671
     672       if (tr_scheme.eq.1) then
     673! Case 1: pseudo-chemistry with relaxation toward fixed profile
     674         call phytrac_relax (debut,lafin,nqmax,
     675     I                   nlon,nlev,dtime,pplay,
    680676     O                   tr_seri)
     677
     678       elseif (tr_scheme.eq.2) then
     679! Case 2: surface emission
     680! For the moment, inspired from Mars version
     681! However, the variable 'source' could be used in physiq
     682! so the call to phytrac_emiss could be to initialise it.
     683         call phytrac_emiss ( (rjourvrai+gmtime)*RDAY,
     684     I                   debut,lafin,nqmax,
     685     I                   nlon,nlev,dtime,paprs,
     686     I                   rlatd,rlond,
     687     O                   tr_seri)
     688       elseif (tr_scheme.eq.3) then
     689! Case 3: Full chemistry
     690!        call phytrac_chem ( ?? )
     691         print*,"Chemistry not yet implemented..."
     692         print*,"See Aurelien Stolzenbach"
     693       endif
    681694      endif
    682695

  • trunk/LMDZ.VENUS/libf/phyvenus/phytrac_relax.F

    r1157 r1160  
    44c
    55c
    6       SUBROUTINE phytrac (nstep,
    7      I                    gmtime,
    8      I                    debutphy,
    9      I                    lafin,
     6      SUBROUTINE phytrac_relax (debutphy,lafin,
    107     I                    nqmax,
    118     I                    nlon,
    129     I                    nlev,
    1310     I                    pdtphys,
    14      I                    u,
    15      I                    v,
    16      I                    t_seri,
    17      I                    paprs,
    1811     I                    pplay,
    19      I                    xlat,
    20      I                    xlon,
    21      I                    presnivs,
    22      I                    pphis,
    23      I                    pphi,
    24      I                    albsol,
    2512     O                    tr_seri)
    2613
    … …  
    3219cAA--------------------
    3320cAA 1/ le call phytrac se fait avec nqmax
     21c
     22c SL: Janvier 2014
     23c Version developed by E. Marcq for pseudo-chemistry relaxation
     24c See Marcq&Lebonnois 2013.
     25c
    3426c======================================================================
    3527      USE ioipsl
    … …  
    4133#include "YOMCST.h"
    4234#include "dimensions.h"
    43 #include "clesphys.h" !///utile?
     35#include "clesphys.h"
    4436#include "temps.h"
    4537#include "paramet.h"
    … …  
    4739
    4840c Arguments:
    49 c
     41
    5042c   EN ENTREE:
    5143c   ==========
    52 c
    53 c   divers:
    54 c   -------
    55 c
     44
     45      logical debutphy       ! le flag de l'initialisation de la physique
     46      logical lafin          ! le flag de la fin de la physique
     47      integer nqmax ! nombre de traceurs auxquels on applique la physique
    5648      integer nlon  ! nombre de points horizontaux
    5749      integer nlev  ! nombre de couches verticales
    58       integer nqmax ! nombre de traceurs auxquels on applique la physique
    59       integer nstep  ! appel physique
    60       integer nseuil ! numero du premier traceur non CV
    61 c      integer julien !jour julien
    62 c      integer itop_con(nlon)
    63 c      integer ibas_con(nlon)
    64       real gmtime
    6550      real pdtphys  ! pas d'integration pour la physique (seconde)
    66       real t_seri(nlon,nlev) ! temperature
     51      real pplay(nlon,nlev)  ! pression pour le mileu de chaque couche (en Pa)
     52
     53c   EN ENTREE/SORTIE:
     54c   =================
     55
    6756      real tr_seri(nlon,nlev,nqmax) ! traceur 
    68       real u(nlon,nlev)
    69       real v(nlon,nlev)
    70       real albsol(nlon)  ! albedo surface
    71       real paprs(nlon,nlev+1)  ! pression pour chaque inter-couche (en Pa)
    72       real ps(nlon)  ! pression surface
    73       real pplay(nlon,nlev)  ! pression pour le mileu de chaque couche (en Pa)
    74       real pphi(nlon,nlev) ! geopotentiel
    75       real pphis(nlon)
    76       REAL xlat(nlon)       ! latitudes pour chaque point
    77       REAL xlon(nlon)       ! longitudes pour chaque point
    78       REAL presnivs(nlev)
    79       logical debutphy       ! le flag de l'initialisation de la physique
    80       logical lafin          ! le flag de la fin de la physique
    81 c      REAL flxmass_w(nlon,nlev)
    8257
    8358cAA ----------------------------
    8459cAA  VARIABLES LOCALES TRACEURS
    8560cAA ----------------------------
    86 cAA
    8761
    8862C les traceurs
    89 C
    90       logical flagCO_OCS
     63
    9164c===================
    9265c it--------indice de traceur
    … …  
    9467c===================
    9568c Variables deja declarees dont on a besoin pour traceurs   
    96 c   k,i,it,tr_seri(klon,klev,nqmax),pplay(nlon,nlev),
     69c   k,i,it,tr_seri(nlon,nlev,nqmax),pplay(nlon,nlev),
    9770      integer nqCO_OCS
    9871c      real pzero,gamma
    … …  
    10073c      parameter (gamma=5000.)
    10174      REAL alpha
    102       real deltatr(klon,klev,nqtot) ! ecart au profil de ref zprof
     75      real deltatr(nlon,nlev,nqtot) ! ecart au profil de ref zprof
    10376      real,save,allocatable :: zprof(:,:)
    10477      real,save,allocatable :: tau(:,:) ! temps de relaxation vers le profil (s)
    10578c======================================================================
    106 c
    107 c Declaration des procedures appelees
    108 c
    109 c--modif convection tiedtke
     79
    11080      INTEGER i, k, it
    111       INTEGER iq, iiq
    112       REAL delp(klon,klev)
    113 c--end modif
    114 c
     81
    11582c Variables liees a l'ecriture de la bande histoire physique
    116 c
     83
    11784c Variables locales pour effectuer les appels en serie
    11885c----------------------------------------------------
    119 c
    120       REAL d_tr(klon,klev), d_trs(klon) ! tendances de traceurs
    121       REAL d_tr_cl(klon,klev,nqmax) ! tendance de traceurs  couche limite
    122       REAL d_tr_cv(klon,klev,nqmax) ! tendance de traceurs  conv pour chq traceur
    123 C
     86
     87      REAL d_tr(nlon,nlev) ! tendances de traceurs
     88
    12489      character*20 modname
    12590      character*80 abort_message
    126 c
    127 c   Controles
    128 c-------------
    129       logical first,couchelimite,convection
    130       save first,couchelimite,convection
    131 c Olivia
    132        data first,couchelimite,convection
    133      s     /.true.,.false.,.false./
    134 
    135       modname = 'phytrac'
    136 c======================================================================
    137 
    138        if(first) then
    139          allocate(zprof(klev,nqtot),tau(klev,nqtot))
    140          first = .false.
    141        endif
    142 
    143          ps(:)=paprs(:,1)
     91
     92c======================================================================
     93
     94      modname = 'phytrac_relax'
    14495c TRACEURS TYPE CO ET OCS
    145          flagCO_OCS = .true.
    146       if (flagCO_OCS) then
    147          nqCO_OCS   = 6
    148       else
    149          nqCO_OCS   = 0
    150       endif  ! flagCO_OCS
     96      nqCO_OCS   = 6
    15197
    15298c---------
    15399c debutphy
    154100c---------
    155          if (debutphy) then
    156                  print*,"DEBUT PHYTRAC"
    157 C         
     101      if (debutphy) then
     102         print*,"DEBUT PHYTRAC"
     103         print*,"PHYTRAC: RELAXATION"
     104         allocate(zprof(nlev,nqtot),tau(nlev,nqtot))
     105
    158106c=============================================================
    159107c=============================================================
    … …  
    163111c=============================================================
    164112c=============================================================
    165 c
    166 c=============================================================
    167 c=============================================================
    168 
    169 C=========================================================================
    170 C=========================================================================
    171       if (flagCO_OCS) then
     113
     114C=========================================================================
     115C=========================================================================
     116
    172117c II) Declaration d'un profil vertical de traceur OK
    173118c
    … …  
    187132
    188133       print*,"INIT TAU"
    189        do k=1,klev
     134       do k=1,nlev
    190135         tau(k,1)=1.e6
    191136         tau(k,2)=1.e7
    … …  
    200145      do it=1,3
    201146       print*,"INIT ZPROF ",tname(it)
    202        do k=1,klev
     147       do k=1,nlev
    203148         zprof(k,it)=0.
    204149c pour l'instant, tau fixe, mais possibilite de le faire varier avec z
    205         if (pplay(klon/2,k) >= 4.8e6) then
     150        if (pplay(nlon/2,k) >= 4.8e6) then
    206151           zprof(k,it)=14.
    207152        endif
    208         if ((pplay(klon/2,k)<=4.8e6).and.(pplay(klon/2,k)>=1.9e6)) then
    209            alpha=(log(pplay(klon/2,k))-log(1.9e6))/
     153        if ((pplay(nlon/2,k)<=4.8e6).and.(pplay(nlon/2,k)>=1.9e6)) then
     154           alpha=(log(pplay(nlon/2,k))-log(1.9e6))/
    210155     .     (log(4.8e6)-log(1.9e6))
    211156           zprof(k,it)=20.*(14./20.)**alpha
    212157        endif
    213         if ((pplay(klon/2,k)<=1.9e6).and.(pplay(klon/2,k)>=1.5e5)) then
    214            alpha=(log(pplay(klon/2,k))-log(1.5e5))/
     158        if ((pplay(nlon/2,k)<=1.9e6).and.(pplay(nlon/2,k)>=1.5e5)) then
     159           alpha=(log(pplay(nlon/2,k))-log(1.5e5))/
    215160     .     (log(1.9e6)-log(1.5e5))
    216161           zprof(k,it)=39.*(20./39.)**alpha
    217162        endif
    218         if ((pplay(klon/2,k)<=1.5e5).and.(pplay(klon/2,k)>=1.1e4)) then
    219            alpha=(log(pplay(klon/2,k))-log(1.1e4))/
     163        if ((pplay(nlon/2,k)<=1.5e5).and.(pplay(nlon/2,k)>=1.1e4)) then
     164           alpha=(log(pplay(nlon/2,k))-log(1.1e4))/
    220165     .     (log(2.73e5)-log(1.1e4))
    221166           zprof(k,it)=50.*(39./50.)**alpha
    222167        endif
    223         if ((pplay(klon/2,k)<=1.1e4).and.(pplay(klon/2,k)>=1.3e3)) then
    224            alpha=(log(pplay(klon/2,k))-log(1.3e3))/
     168        if ((pplay(nlon/2,k)<=1.1e4).and.(pplay(nlon/2,k)>=1.3e3)) then
     169           alpha=(log(pplay(nlon/2,k))-log(1.3e3))/
    225170     .     (log(1.1e4)-log(1.3e3))
    226171           zprof(k,it)=2.*(50./2.)**alpha
    227172        endif
    228         if ((pplay(klon/2,k)<=1.3e3).and.(pplay(klon/2,k)>=2.4)) then
    229            alpha=(log(pplay(klon/2,k))-log(2.4))/
     173        if ((pplay(nlon/2,k)<=1.3e3).and.(pplay(nlon/2,k)>=2.4)) then
     174           alpha=(log(pplay(nlon/2,k))-log(2.4))/
    230175     .     (log(1.3e3)-log(2.4))
    231176           zprof(k,it)=1000.*(2./1000.)**alpha
    232177        endif
    233         if (pplay(klon/2,k) <= 2.4) then
     178        if (pplay(nlon/2,k) <= 2.4) then
    234179           zprof(k,it)=1000.
    235180        endif
    … …  
    239184c OCS
    240185       print*,"INIT ZPROF ",tname(it+3)
    241        do k=1,klev
     186       do k=1,nlev
    242187         zprof(k,it+3)=0.
    243          if (pplay(klon/2,k) >= 4.8e6) then
     188         if (pplay(nlon/2,k) >= 4.8e6) then
    244189           zprof(k,it+3)=30.
    245190         endif
    246          if ((pplay(klon/2,k)<=4.8e6).and.(pplay(klon/2,k)>=9.4e5))
     191         if ((pplay(nlon/2,k)<=4.8e6).and.(pplay(nlon/2,k)>=9.4e5))
    247192     *   then
    248            alpha=(log(pplay(klon/2,k))-log(9.4e5))/
     193           alpha=(log(pplay(nlon/2,k))-log(9.4e5))/
    249194     *     (log(4.8e6)-log(9.4e5))
    250195           zprof(k,it+3)=20.*(30/20.)**alpha
    251196         endif
    252          if ((pplay(klon/2,k)<=9.4e5).and.(pplay(klon/2,k)>=4.724e5))
     197         if ((pplay(nlon/2,k)<=9.4e5).and.(pplay(nlon/2,k)>=4.724e5))
    253198     *   then
    254            alpha=(log(pplay(klon/2,k))-log(4.724e5))/
     199           alpha=(log(pplay(nlon/2,k))-log(4.724e5))/
    255200     *     (log(9.4e5)-log(4.724e5))
    256201           zprof(k,it+3)=0.5*(20/0.5)**alpha
    257202         endif
    258          if ((pplay(klon/2,k)<=4.724e5).and.(pplay(klon/2,k)>=1.1e4))
     203         if ((pplay(nlon/2,k)<=4.724e5).and.(pplay(nlon/2,k)>=1.1e4))
    259204     *   then
    260            alpha=(log(pplay(klon/2,k))-log(1.1e4))/
     205           alpha=(log(pplay(nlon/2,k))-log(1.1e4))/
    261206     *     (log(4.724e5)-log(1.1e4))
    262207           zprof(k,it+3)=0.005*(0.5/0.005)**alpha
    263208         endif
    264          if (pplay(klon/2,k)<=1.1e4) then
     209         if (pplay(nlon/2,k)<=1.1e4) then
    265210           zprof(k,it+3)=0.
    266211         endif
    … …  
    271216c Initialisation du traceur s'il est nul:
    272217       do it=1,nqCO_OCS
    273         if ((tr_seri(klon/2,1,it).eq.0.).and.
    274      .      (tr_seri(klon/2,klev/2,it).eq.0.).and.
    275      .      (tr_seri(klon/2,klev,it).eq.0.)) then
     218        if ((tr_seri(nlon/2,1,it).eq.0.).and.
     219     .      (tr_seri(nlon/2,nlev/2,it).eq.0.).and.
     220     .      (tr_seri(nlon/2,nlev,it).eq.0.)) then
    276221         print*,"INITIALISATION DE ",tname(it)
    277          do k=1,klev
    278            do i=1,klon
     222         do k=1,nlev
     223           do i=1,nlon
    279224             tr_seri(i,k,it) = zprof(k,it)
    280225           enddo
    … …  
    283228       enddo
    284229
    285 C=========================================================================
    286       endif  ! flagCO_OCS
    287230C=========================================================================
    288231C=========================================================================
    … …  
    294237
    295238c======================================================================
    296       if (flagCO_OCS) then
    297239c Rappel vers un profil
    298240c======================================================================
    299241         do it=1,nqCO_OCS
    300            do k=1,klev
    301              do i=1,klon
     242           do k=1,nlev
     243             do i=1,nlon
    302244c VERIF
    303245           if (tr_seri(i,k,it).lt.0) then
    … …  
    322264
    323265c======================================================================
    324       endif  ! flagCO_OCS
    325 c======================================================================
    326 
    327 c======================================================================
    328 c   Calcul de l'effet de la couche limite remis directement dans physiq
    329266c======================================================================
    330267
    … …  
    334271     
    335272     
    336 c=========================================================================
    337 c=========================================================================
    338 c=========================================================================
    339 c ARCHIVES ===============================================================
    340 c=========================================================================
    341 c=========================================================================
    342 c=========================================================================
    343 
    344 c===========
    345 c definition de traceurs idealises
    346 c==========
    347 c
    348 c I) Declaration directe du traceur a altitude fixee
    349 c
    350 c a) traceur en carre OK
    351 c
    352 c         do i=1,klon
    353 c         tr_seri(i,:,1)=0.
    354 c        if ((xlat(i)>=0.).and.(xlat(i)<=-30.)) then
    355 c          if ((xlon(i)>=0.).and.(xlon(i)<=40.)) then
    356 c              tr_seri(i,10,1)=1.
    357 c          endif
    358 c        endif
    359 c      end do
    360 c
    361 c a bis) 2 traceurs en carre lat/alt, uniforme en longitude OK
    362 c
    363 C entre 45-55 km
    364 c
    365 c         do i=1,klon
    366 c         do k=1,klev+1
    367 cc         tr_seri(i,k,1)=0.
    368 c           if ((xlat(i)>=60.).and.(xlat(i)<=80.)) then
    369 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    370 c           if ((pplay(klon/2,k)>=5.e4).and.(pplay(klon/2,k)<=4.e5)) then
    371 c               tr_seri(i,k,1)=1.
    372 c           endif
    373 c           endif
    374 c           endif
    375 c         else
    376 c            tr_seri(i,k,1)=0.
    377 c         end do
    378 c         end do
    379 cc
    380 c         do i=1,klon
    381 c         do k=1,klev+1
    382 cc         tr_seri(i,k,2)=0.
    383 c           if ((xlat(i)>=-60.).and.(xlat(i)<=-80.)) then
    384 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    385 c           if ((pplay(klon/2,k)>=5.e4).and.(pplay(klon/2,k)<=4.e5)) then
    386 c               tr_seri(i,k,2)=1.
    387 c           endif
    388 c           endif
    389 c           endif
    390 c         else
    391 c            tr_seri(i,k,2)=0.
    392 c         end do
    393 c         end do
    394 cc
    395 c         do i=1,klon
    396 c         do k=1,klev+1
    397 cc         tr_seri(i,k,3)=0.
    398 c           if ((xlat(i)>=40.).and.(xlat(i)<=60.)) then
    399 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    400 c           if ((pplay(klon/2,k)>=5.e4).and.(pplay(klon/2,k)<=4.e5)) then
    401 c               tr_seri(i,k,3)=1.
    402 c           endif
    403 c           endif
    404 c           endif
    405 c         else
    406 c            tr_seri(i,k,3)=0.
    407 c         end do
    408 c         end do
    409 cc
    410 c         do i=1,klon
    411 c         do k=1,klev+1
    412 cc         tr_seri(i,k,4)=0.
    413 c           if ((xlat(i)>=-40.).and.(xlat(i)<=-60.)) then
    414 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    415 c           if ((pplay(klon/2,k)>=5.e4).and.(pplay(klon/2,k)<=4.e5)) then
    416 c               tr_seri(i,k,4)=1.
    417 c           endif
    418 c           endif
    419 c           endif
    420 c         else
    421 c            tr_seri(i,k,4)=0.
    422 c         end do
    423 c         end do
    424 cc
    425 c         do i=1,klon
    426 c         do k=1,klev+1
    427 cc         tr_seri(i,k,5)=0.
    428 c           if ((xlat(i)>=-20.).and.(xlat(i)<=20.)) then
    429 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    430 c           if ((pplay(klon/2,k)>=5.e4).and.(pplay(klon/2,k)<=4.e5)) then
    431 c              tr_seri(i,k,5)=1.
    432 c           endif
    433 c           endif
    434 c           endif
    435 c         else
    436 c            tr_seri(i,k,5)=0.
    437 c         end do
    438 c         end do
    439 c
    440 c entre 35-45 km
    441 c
    442 c         do i=1,klon
    443 c         do k=1,klev+1
    444 cc         tr_seri(i,k,6)=0.
    445 c           if ((xlat(i)>=60.).and.(xlat(i)<=80.)) then
    446 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    447 c           if ((pplay(klon/2,k)>=4.e5).and.(pplay(klon/2,k)<=8.e6)) then
    448 c               tr_seri(i,k,6)=1.
    449 c           endif
    450 c           endif
    451 c           endif
    452 c         else
    453 c            tr_seri(i,k,6)=0.
    454 c         end do
    455 c         end do
    456 c
    457 c         do i=1,klon
    458 c         do k=1,klev+1
    459 cc         tr_seri(i,k,7)=0.
    460 c           if ((xlat(i)>=-60.).and.(xlat(i)<=-80.)) then
    461 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    462 c           if ((pplay(klon/2,k)>=4.e5).and.(pplay(klon/2,k)<=8.e6)) then
    463 c               tr_seri(i,k,7)=1.
    464 c           endif
    465 c           endif
    466 c           endif
    467 c         else
    468 c            tr_seri(i,k,7)=0.
    469 c         end do
    470 c         end do
    471 c
    472 C entre 50-60 km
    473 c
    474 c         do i=1,klon
    475 c         do k=1,klev+1
    476 cc         tr_seri(i,k,8)=0.
    477 c           if ((xlat(i)>=60.).and.(xlat(i)<=80.)) then
    478 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    479 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=1.e5)) then
    480 c               tr_seri(i,k,8)=1.
    481 c           endif
    482 c           endif
    483 c          endif
    484 c         else
    485 c            tr_seri(i,k,8)=0.
    486 c         end do
    487 c         end do
    488 c
    489 c         do i=1,klon
    490 c         do k=1,klev+1
    491 cc         tr_seri(i,k,9)=0.
    492 c           if ((xlat(i)>=-80.).and.(xlat(i)<=-60.)) then
    493 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    494 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=1.e5)) then
    495 c               tr_seri(i,k,9)=1.
    496 c           endif
    497 c           endif
    498 c           endif
    499 c         else
    500 c            tr_seri(i,k,9)=0.
    501 c         end do
    502 c         end do
    503 c
    504 c         do i=1,klon
    505 c         do k=1,klev+1
    506 cc         tr_seri(i,k,10)=0.
    507 c           if ((xlat(i)>=40.).and.(xlat(i)<=60.)) then
    508 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    509 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=1.e5)) then
    510 c               tr_seri(i,k,10)=1.
    511 c           endif
    512 c           endif
    513 c           endif
    514 c         else
    515 c            tr_seri(i,k,10)=0.
    516 c         end do
    517 c         end do
    518 c
    519 c         do i=1,klon
    520 c         do k=1,klev+1
    521 cc         tr_seri(i,k,11)=0.
    522 c           if ((xlat(i)>=-60.).and.(xlat(i)<=-40.)) then
    523 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    524 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=1.e5)) then
    525 c               tr_seri(i,k,11)=1.
    526 c           endif
    527 c           endif
    528 c           endif
    529 c         else
    530 c            tr_seri(i,k,11)=0.
    531 c         end do
    532 c         end do
    533 c
    534 c         do i=1,klon
    535 c         do k=1,klev+1
    536 cc         tr_seri(i,k,12)=0.
    537 c           if ((xlat(i)>=-20.).and.(xlat(i)<=20.)) then
    538 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    539 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=1.e5)) then
    540 c               tr_seri(i,k,12)=1.
    541 c           endif
    542 c           endif
    543 c           endif
    544 c         else
    545 c            tr_seri(i,k,12)=0.
    546 c         end do
    547 c         end do
    548 c
    549 c entre 20-30 km
    550 c
    551 c         do i=1,klon
    552 c         do k=1,klev+1
    553 cc         tr_seri(i,k,13)=0.
    554 c           if ((xlat(i)>=60.).and.(xlat(i)<=80.)) then
    555 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    556 c           if ((pplay(klon/2,k)>=1.e6).and.(pplay(klon/2,k)<=2.e6)) then
    557 c               tr_seri(i,k,13)=1.
    558 c           endif
    559 c           endif
    560 c           endif
    561 c         else
    562 c            tr_seri(i,k,13)=0.
    563 c         end do
    564 c         end do
    565 c
    566 c         do i=1,klon
    567 c         do k=1,klev+1
    568 cc         tr_seri(i,k,14)=0.
    569 c           if ((xlat(i)>=-80.).and.(xlat(i)<=-60.)) then
    570 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    571 c           if ((pplay(klon/2,k)>=1.e6).and.(pplay(klon/2,k)<=2.e6)) then
    572 c               tr_seri(i,k,14)=1.
    573 c           endif
    574 c           endif
    575 c           endif
    576 c         else
    577 c            tr_seri(i,k,14)=0.
    578 c         end do
    579 c         end do
    580 c
    581 c         do i=1,klon
    582 c         do k=1,klev+1
    583 cc         tr_seri(i,k,15)=0.
    584 c           if ((xlat(i)>=-20.).and.(xlat(i)<=20.)) then
    585 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    586 c           if ((pplay(klon/2,k)>=1.e6).and.(pplay(klon/2,k)<=2.e6)) then
    587 c               tr_seri(i,k,15)=1.
    588 c           endif
    589 c           endif
    590 c           endif
    591 c         else
    592 c            tr_seri(i,k,15)=0.
    593 c         end do
    594 c         end do
    595 c
    596 c entre 55-65 km
    597 c
    598 c         do i=1,klon
    599 c         do k=1,klev+1
    600 cc         tr_seri(i,k,16)=0.
    601 c           if ((xlat(i)>=60.).and.(xlat(i)<=80.)) then
    602 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    603 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=5.e4)) then
    604 c               tr_seri(i,k,16)=1.
    605 c           endif
    606 c           endif
    607 c           endif
    608 c           endif
    609 c         else
    610 c            tr_seri(i,k,16)=0.
    611 c         end do
    612 c         end do
    613 c
    614 c         do i=1,klon
    615 c         do k=1,klev+1
    616 cc         tr_seri(i,k,17)=0.
    617 c           if ((xlat(i)>=-80.).and.(xlat(i)<=-60.)) then
    618 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    619 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=5.e4)) then
    620 c               tr_seri(i,k,17)=1.
    621 c          endif
    622 c          endif
    623 c           endif
    624 c           endif
    625 c         else
    626 c            tr_seri(i,k,17)=0.
    627 c         end do
    628 c         end do
    629 c
    630 c         do i=1,klon
    631 c         do k=1,klev+1
    632 cc         tr_seri(i,k,18)=0.
    633 c           if ((xlat(i)>=-20.).and.(xlat(i)<=20.)) then
    634 c           if ((xlon(i)>=-180.).and.(xlon(i)<=180.)) then
    635 c           if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=5.e4)) then
    636 c               tr_seri(i,k,18)=1.
    637 c           endif
    638 c           endif
    639 c           endif
    640 c           endif
    641 c         else
    642 c            tr_seri(i,k,18)=0.
    643 c         end do
    644 c         end do
    645 c
    646 c b) traceur a une bande en latitudeOK
    647 c
    648 c a 65km
    649 c
    650 c        do i=1,klon
    651 c         tr_seri(i,:,1)=0.
    652 c        if ((xlat(i)>=60.).and.(xlat(i)<=80.)) then
    653 c              tr_seri(i,20,1)=1.
    654 c        endif
    655 c      end do 
    656 c
    657 c        do i=1,klon
    658 c         tr_seri(i,:,2)=0.
    659 c        if ((xlat(i)>=40.).and.(xlat(i)<=60.)) then
    660 c              tr_seri(i,20,2)=1.
    661 c        endif
    662 c      end do 
    663 c
    664 c        do i=1,klon
    665 c         tr_seri(i,:,3)=0.
    666 c        if ((xlat(i)>=20.).and.(xlat(i)<=40.)) then
    667 c              tr_seri(i,20,3)=1.
    668 c        endif
    669 c      end do 
    670 c
    671 c        do i=1,klon
    672 c         tr_seri(i,:,4)=0.
    673 c        if ((xlat(i)>=0.).and.(xlat(i)<=20.)) then
    674 c              tr_seri(i,20,4)=1.
    675 c        endif
    676 c      end do 
    677 c
    678 c        do i=1,klon
    679 c         tr_seri(i,:,5)=0.
    680 c        if ((xlat(i)>=-20.).and.(xlat(i)<=0.)) then
    681 c              tr_seri(i,20,5)=1.
    682 c        endif
    683 c      end do 
    684 c
    685 c        do i=1,klon
    686 c         tr_seri(i,:,6)=0.
    687 c        if ((xlat(i)>=-40.).and.(xlat(i)<=-20.)) then
    688 c              tr_seri(i,20,6)=1.
    689 c        endif
    690 c      end do 
    691 c
    692 c        do i=1,klon
    693 c         tr_seri(i,:,7)=0.
    694 c        if ((xlat(i)>=-60.).and.(xlat(i)<=-40.)) then
    695 c              tr_seri(i,20,7)=1.
    696 c        endif
    697 c      end do 
    698 c
    699 c        do i=1,klon
    700 c         tr_seri(i,:,8)=0.
    701 c        if ((xlat(i)>=-80.).and.(xlat(i)<=-60.)) then
    702 c              tr_seri(i,20,8)=1.
    703 c        endif
    704 c      end do 
    705 c
    706 c a 50km
    707 c
    708 c        do i=1,klon
    709 c        tr_seri(i,:,1)=0.
    710 c        if ((xlat(i)>=40.).and.(xlat(i)<=60.)) then
    711 c              tr_seri(i,27,1)=1.
    712 c        endif
    713 c      end do 
    714 c
    715 c        do i=1,klon
    716 c         tr_seri(i,:,2)=0.
    717 c        if ((xlat(i)>=60.).and.(xlat(i)<=80.)) then
    718 c              tr_seri(i,27,2)=1.
    719 c        endif
    720 c      end do 
    721 c
    722 c        do i=1,klon
    723 c         tr_seri(i,:,3)=0.
    724 c        if ((xlat(i)>=20.).and.(xlat(i)<=40.)) then
    725 c              tr_seri(i,27,3)=1.
    726 c        endif
    727 c      end do 
    728 c
    729 c        do i=1,klon
    730 c         tr_seri(i,:4)=0.
    731 c        if ((xlat(i)>=0.).and.(xlat(i)<=20.)) then
    732 c              tr_seri(i,27,4)=1.
    733 c       endif
    734 c      end do 
    735 c
    736 c        do i=1,klon
    737 c         tr_seri(i,:,5)=0.
    738 c        if ((xlat(i)>=-20.).and.(xlat(i)<=0.)) then
    739 c              tr_seri(i,27,5)=1.
    740 c        endif
    741 c      end do 
    742 c
    743 c        do i=1,klon
    744 c         tr_seri(i,:,6)=0.
    745 c        if ((xlat(i)>=-40.).and.(xlat(i)<=-20.)) then
    746 c              tr_seri(i,27,6)=1.
    747 c        endif
    748 c      end do 
    749 c
    750 c        do i=1,klon
    751 c         tr_seri(i,:,7)=0.
    752 c        if ((xlat(i)>=-60.).and.(xlat(i)<=-40.)) then
    753 c              tr_seri(i,27,7)=1.
    754 c        endif
    755 c      end do 
    756 c
    757 c        do i=1,klon
    758 c         tr_seri(i,:,8)=0.
    759 c        if ((xlat(i)>=-80.).and.(xlat(i)<=-60.)) then
    760 c              tr_seri(i,27,8)=1.
    761 c        endif
    762 c      end do 
    763 c
    764 c c) traceur a plusieurs bandes en latitude OK
    765 c
    766 c         do i=1,klon
    767 c        tr_seri(i,:,2)=0.
    768 c        if ((xlat(i)>=50.).and.(xlat(i)<=70.)) then
    769 c             tr_seri(i,10,2)=1.
    770 c        endif
    771 c        if ((xlat(i)>=-10.).and.(xlat(i)<=10.)) then
    772 c              tr_seri(i,10,2)=1.
    773 c        endif
    774 c
    775 c        if ((xlat(i)>=-70.).and.(xlat(i)<=-50.)) then
    776 c              tr_seri(i,10,2)=1.
    777 c        endif         
    778 c      end do
    779 c
    780 c d) traceur a une bande en altitude OK
    781 c
    782 c       do k=1,klev+1
    783 c         tr_seri(:,k,1)=0.
    784 c         if ((pplay(klon/2,k)>=1.e5).and.(pplay(klon/2,k)<=1.e6)) then
    785 c             tr_seri(:,k,1)=1.
    786 c         endif
    787 c       end do
    788 c
    789 c dbis) plusieurs traceurs a une bande en altitude OK
    790 c
    791 c bande tres basse tropo
    792 c      do k=1,klev
    793 c        tr_seri(:,k,1)=0.
    794 c        if ((pplay(klon/2,k)>=5.e5).and.(pplay(klon/2,k)<=5.e6)) then
    795 c            tr_seri(:,k,1)=1.
    796 c        endif
    797 c      end do
    798 c bande dans les nuages et un peu en-dessous   
    799 c      do k=1,klev
    800 c         tr_seri(:,k,2)=0.
    801 c         if ((pplay(klon/2,k)>=5.e4).and.(pplay(klon/2,k)<=5.e5)) then
    802 c             tr_seri(:,k,2)=1.
    803 c         endif
    804 c       end do
    805 cune grosse epaisseur: inclue toute la circulation meridienne
    806 c      do k=1,klev
    807 c        tr_seri(:,k,1)=0.
    808 c        if ((pplay(klon/2,k)>=1.e4).and.(pplay(klon/2,k)<=1.e6)) then
    809 c            tr_seri(:,k,1)=1.
    810 c        endif
    811 c      end do
    812 cune grosse epaisseur: inclue la mesosphere
    813 c      do k=1,klev
    814 c         tr_seri(:,k,2)=0.
    815 c         if ((pplay(klon/2,k)>=2.e2).and.(pplay(klon/2,k)<=1.e4)) then
    816 c             tr_seri(:,k,2)=1.
    817 c         endif
    818 c       end do
    819 c
    820 c      do k=1,klev
    821 c        tr_seri(:,k,3)=0.
    822 c        if ((pplay(klon/2,k)>=5.e1).and.(pplay(klon/2,k)<=5.e2)) then
    823 c            tr_seri(:,k,3)=1.
    824 c        endif
    825 c      end do
    826 c
    827 c e) plusieurs couches verticales de traceurs, a plusieurs bandes en latitude???
    828 c       
    829 c au sol
    830 c         do i=1,klon
    831 c        tr_seri(i,:,1)=0.
    832 c        if ((xlat(i)>=50.).and.(xlat(i)<=70.)) then
    833 c             tr_seri(i,5,1)=1.
    834 c        endif
    835 c        if ((xlat(i)>=-10.).and.(xlat(i)<=10.)) then
    836 c              tr_seri(i,5,1)=1.
    837 c        endif
    838 c
    839 c        if ((xlat(i)>=-70.).and.(xlat(i)<=-50.)) then
    840 c              tr_seri(i,5,1)=1.
    841 c        endif         
    842 c      end do
    843 c
    844 c         do i=1,klon
    845 c        tr_seri(i,:,2)=0.
    846 c        if ((xlat(i)>=50.).and.(xlat(i)<=70.)) then
    847 c             tr_seri(i,10,2)=1.
    848 c        endif
    849 c        if ((xlat(i)>=-10.).and.(xlat(i)<=10.)) then
    850 c              tr_seri(i,10,2)=1.
    851 c      endif
    852 c
    853 c        if ((xlat(i)>=-70.).and.(xlat(i)<=-50.)) then
    854 c              tr_seri(i,10,2)=1.
    855 c        endif         
    856 c      end do
    857 c
    858 c         do i=1,klon
    859 c        tr_seri(i,:,3)=0.
    860 c        if ((xlat(i)>=50.).and.(xlat(i)<=70.)) then
    861 c             tr_seri(i,30,3)=1.
    862 c        endif
    863 c        if ((xlat(i)>=-10.).and.(xlat(i)<=10.)) then
    864 c              tr_seri(i,30,3)=1.
    865 c        endif
    866 c
    867 c        if ((xlat(i)>=-70.).and.(xlat(i)<=-50.)) then
    868 c              tr_seri(i,30,3)=1.
    869 c        endif         
    870 c      end do
    871 c       
    872 c        do i=1,klon
    873 c        tr_seri(i,:,4)=0.
    874 c        if ((xlat(i)>=50.).and.(xlat(i)<=70.)) then
    875 c             tr_seri(i,45,4)=1.
    876 c        endif
    877 c        if ((xlat(i)>=-10.).and.(xlat(i)<=10.)) then
    878 c              tr_seri(i,45,4)=1.
    879 c        endif
    880 c
    881 c        if ((xlat(i)>=-70.).and.(xlat(i)<=-50.)) then
    882 c              tr_seri(i,45,4)=1.
    883 c        endif         
    884 c      end do
    885 c     
    886 
Note: See TracChangeset for help on using the changeset viewer.