Changeset 3318 for trunk


Ignore:
Timestamp:
Apr 26, 2024, 4:27:26 PM (14 months ago)
Author:
slebonnois
Message:

Titan PCM update : optics + microphysics

Location:
trunk/LMDZ.TITAN/libf
Files:
20 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.TITAN/libf/muphytitan/fsystem.F90

    r3090 r3318  
    428428    CHARACTER(len=:), ALLOCATABLE :: opath
    429429      !! A Fortran allocated string with the parent directory path or an empty string if method fails
     430    CHARACTER(len=:), ALLOCATABLE :: cpath
    430431    TYPE(C_PTR) :: zpath
    431432    IF (LEN_TRIM(path) == 0) THEN
     
    433434      RETURN
    434435    ENDIF
    435     zpath = dirname_c(cstring(ADJUSTL(path)))
     436    cpath = cstring(ADJUSTL(path))
     437    zpath = dirname_c(cpath)
    436438    IF (.NOT.C_ASSOCIATED(zpath)) THEN
    437439      opath = ""
     
    448450    CHARACTER(len=:), ALLOCATABLE :: opath
    449451      !! The basename of the path or an empty string if method fails
     452    CHARACTER(len=:), ALLOCATABLE :: cpath
    450453    TYPE(C_PTR) :: zpath
    451454    IF (LEN_TRIM(path) == 0) THEN
     
    453456      RETURN
    454457    ENDIF
    455     zpath = basename_c(cstring(ADJUSTL(path)))
     458    cpath = cstring(ADJUSTL(path))
     459    zpath = basename_c(cpath)
    456460    IF (.NOT.C_ASSOCIATED(zpath)) THEN
    457461      opath = ""
     
    472476    CHARACTER(len=:), ALLOCATABLE :: opath
    473477      !! The absolute of the path or an empty string if method fails
     478    CHARACTER(len=:), ALLOCATABLE :: cpath
    474479    TYPE(C_PTR) :: zpath
    475     zpath = realpath_c(cstring(ADJUSTL(path)))
     480    cpath = cstring(ADJUSTL(path))
     481    zpath = realpath_c(cpath)
    476482    IF (.NOT.C_ASSOCIATED(zpath)) THEN
    477483      opath = ""
     
    490496                                    reldir  !! A directory path from which output should be relative to
    491497    CHARACTER(len=:), ALLOCATABLE :: res    !! An allocated string with the resulting path
     498    CHARACTER(len=:), ALLOCATABLE :: cpath1,cpath2
    492499    TYPE(C_PTR) :: zpath
    493     zpath = relpath_c(cstring(ADJUSTL(path)),cstring(ADJUSTL(reldir)))
     500    cpath1 = cstring(ADJUSTL(path))
     501    cpath2 = cstring(ADJUSTL(reldir))
     502    zpath = relpath_c(cpath1,cpath2)
    494503    IF (.NOT.C_ASSOCIATED(zpath)) THEN
    495504      res = TRIM(ADJUSTL(path))
     
    520529    CHARACTER(len=*), INTENT(in)  :: output !! Output file path destination.
    521530    LOGICAL :: ret                          !! True on success, false otherwise.
     531    CHARACTER(len=:), ALLOCATABLE :: cpath1,cpath2
     532
    522533    IF (LEN_TRIM(input) == 0 .OR. LEN_TRIM(output) == 0 .OR. input == output) THEN
    523534      ret = .false.
    524535    ELSE
    525       ret = INT(copy_c(cstring(ADJUSTL(output)),cstring(ADJUSTL(input)))) == 0
     536      cpath1 = cstring(ADJUSTL(output))
     537      cpath2 = cstring(ADJUSTL(input))
     538      ret = INT(copy_c(cpath1,cpath2)) == 0
    526539    ENDIF
    527540    RETURN
     
    532545    CHARACTER(len=*), INTENT(in)  :: path !! A string with the (valid) file path to delete
    533546    LOGICAL :: ret                        !! True on success, false otherwise.
     547    CHARACTER(len=:), ALLOCATABLE :: cpath
    534548    IF (LEN_TRIM(path) == 0) THEN
    535549      ret = .false.
    536550    ELSE
    537       ret = INT(remove_c(cstring(ADJUSTL(path)))) == 0
     551      cpath = cstring(ADJUSTL(path))
     552      ret = INT(remove_c(cpath)) == 0
    538553    ENDIF
    539554    RETURN
     
    545560                                    new    !! A string with the new name of the path
    546561    LOGICAL :: ret                         !! True on success, false otherwise.
     562    CHARACTER(len=:), ALLOCATABLE :: cpath1,cpath2
    547563    IF (LEN_TRIM(old) == 0.OR.LEN_TRIM(new) == 0) THEN
    548564      ret = .false.
    549565    ELSE
    550       ret = INT(rename_c(cstring(ADJUSTL(old)),cstring(ADJUSTL(new)))) == 0
     566      cpath1 = cstring(ADJUSTL(old))
     567      cpath2 = cstring(ADJUSTL(new))
     568      ret = INT(rename_c(cpath1,cpath2)) == 0
    551569    ENDIF
    552570    RETURN
     
    559577    LOGICAL  :: ret                      !! True on success, false otherwise.
    560578    INTEGER(kind=C_INT) :: zmode
     579    CHARACTER(len=:), ALLOCATABLE :: cpath
    561580    IF (LEN_TRIM(path) == 0) THEN
    562581      ret = .false.
    563582    ELSE
    564583      zmode = INT(oct_2_dec(mode),kind=C_INT)
    565       ret = INT(chmod_c(cstring(ADJUSTL(path)), zmode)) == 0
     584      cpath = cstring(ADJUSTL(path))
     585      ret = INT(chmod_c(cpath, zmode)) == 0
    566586    ENDIF
    567587    RETURN
     
    572592    CHARACTER(len=*), INTENT(in) :: path !! Path of the new working directory
    573593    LOGICAL :: ret                       !! True on success, false otherwise.
     594    CHARACTER(len=:), ALLOCATABLE :: cpath
    574595    IF (LEN_TRIM(path) == 0) THEN
    575596      ret = .false.
    576597    ELSE
    577       ret = INT(chdir_c(cstring(ADJUSTL(path)))) == 0
     598      cpath = cstring(ADJUSTL(path))
     599      ret = INT(chdir_c(cpath)) == 0
    578600    ENDIF
    579601    RETURN
     
    595617    INTEGER :: zmode
    596618    LOGICAL :: zperm
     619    CHARACTER(len=:), ALLOCATABLE :: cpath
     620
    597621    IF (LEN_TRIM(path) == 0) THEN
    598622      ret = .false.
     
    605629        zmode = oct_2_dec(mode)
    606630      ENDIF
     631      cpath = cstring(ADJUSTL(path))
    607632      zperm = .false. ; IF (PRESENT(permissive)) zperm = permissive
    608633      IF (zperm) THEN
    609         ret = INT(mkdirp_c(cstring(ADJUSTL(path)),INT(zmode,kind=C_INT))) == 0
     634        ret = INT(mkdirp_c(cpath,INT(zmode,kind=C_INT))) == 0
    610635      ELSE
    611         ret = INT(mkdir_c(cstring(ADJUSTL(path)),INT(zmode,kind=C_INT))) == 0
     636        ret = INT(mkdir_c(cpath,INT(zmode,kind=C_INT))) == 0
    612637      ENDIF
    613638    ENDIF
     
    627652      !! True on success, false otherwise.
    628653    LOGICAL :: zforce
     654    CHARACTER(len=:), ALLOCATABLE :: cpath
    629655    IF (LEN_TRIM(path) == 0) THEN
    630656      ret = .false.
    631657    ELSE
    632658      zforce = .false. ; IF (PRESENT(forced)) zforce = forced
     659      cpath = cstring(ADJUSTL(path))
    633660      IF (.NOT.zforce) THEN
    634         ret = INT(rmdir_c(cstring(ADJUSTL(path)))) == 0
     661        ret = INT(rmdir_c(cpath)) == 0
    635662      ELSE
    636         ret = INT(rmdirf_c(cstring(ADJUSTL(path)))) == 0
     663        ret = INT(rmdirf_c(cpath)) == 0
    637664      ENDIF
    638665    ENDIF
     
    668695    INTEGER(kind=c_long)          :: f
    669696    CHARACTER(len=20,kind=C_CHAR) :: ta,tm,tc
     697    CHARACTER(len=:), ALLOCATABLE :: cpath
    670698    IF (LEN_TRIM(path) == 0) THEN
    671699      ret = .false.; RETURN
     
    677705      ! set default values
    678706      pe=-1 ; ty=-1 ; ud=-1 ; gd=-1 ; fs=-1 ; at="" ; mt="" ; ct=""
    679       ret = INT(fstat_c(cstring(ADJUSTL(path)),p,l,t,u,g,f,ta,tm,tc)) == 0
     707      cpath = cstring(ADJUSTL(path))
     708      ret = INT(fstat_c(cpath,p,l,t,u,g,f,ta,tm,tc)) == 0
    680709      IF (ret) THEN
    681710        pe=INT(p) ; ln=INT(l) ; ty=INT(t) ; ud=INT(u) ; gd=INT(g)
     
    752781    LOGICAL :: ret                              !! True on success, false otherwise.
    753782    INTEGER(kind=C_INT) :: zp
     783    CHARACTER(len=:), ALLOCATABLE :: cpath
    754784    IF (LEN_TRIM(path) == 0) THEN
    755785      ret = .false.
     
    757787      zp = 0 ; IF (PRESENT(permission)) zp = INT(permission,kind=C_INT)
    758788      ! Defaults are set in the C function.
    759       ret = INT(access_c(cstring(ADJUSTL(path)),zp)) == 0
     789      cpath = cstring(ADJUSTL(path))
     790      ret = INT(access_c(cpath,zp)) == 0
    760791    ENDIF
    761792    RETURN
     
    822853    INTEGER                       :: zmd,zt,zp
    823854    CHARACTER(len=:), ALLOCATABLE :: b,e
     855    CHARACTER(len=:), ALLOCATABLE :: cpath
    824856    ret = .false.
    825857    ! Checking for existence
     
    856888    ENDIF
    857889    zp = 0 ; IF(PRESENT(permissive)) THEN ; IF(permissive) zp=1 ; ENDIF
    858     ret = INT(create_c(cstring(ADJUSTL(path)),INT(zmd,kind=C_INT),INT(zt,kind=C_INT),INT(zp,kind=C_INT))) == 0
     890
     891    cpath = cstring(ADJUSTL(path))
     892    ret = INT(create_c(cpath,INT(zmd,kind=C_INT),INT(zt,kind=C_INT),INT(zp,kind=C_INT))) == 0
    859893    RETURN
    860894  END FUNCTION fs_create
  • trunk/LMDZ.TITAN/libf/muphytitan/mm_clouds.f90

    r3090 r3318  
    112112      mm_ccn_vsed(:) = wsettle(mm_play,mm_temp,mm_zlay,mm_drho,mm_drad)
    113113
    114       ! Computes flux [kg.m-2.s-1] and precipitation [m.iphysiq] of ccn
     114      ! Computes flux [kg.m-2.s-1] and precipitation [kg.m-2.iphysiq] of ccn
    115115      mm_ccn_flux(:) = get_mass_flux(mm_rhoaer,mm_m3ccn(:))
    116       mm_ccn_prec = SUM(zdm3n*mm_dzlev)
    117 
    118       ! Computes flux [kg.m-2.s-1] and precipitation [m.iphysiq] of ices
     116      mm_ccn_prec = SUM(zdm3n*mm_dzlev*mm_rhoaer)
     117
     118      ! Computes flux [kg.m-2.s-1] and precipitation [kg.m-2.iphysiq] of ices
    119119      DO i = 1, mm_nesp
    120120        mm_ice_fluxes(:,i) = get_mass_flux(mm_xESPS(i)%rho,(3._mm_wp*mm_m3ice(:,i))/(4._mm_wp*mm_pi))
    121         mm_ice_prec(i) = SUM(zdm3i(:,i)*mm_dzlev)
     121        mm_ice_prec(i) = SUM(zdm3i(:,i)*mm_dzlev*mm_xESPS(i)%rho)
    122122      ENDDO
    123123
     
    258258    ! Saturation ratio
    259259    Xsat = zvapX / qsat
     260
    260261   
    261262    ! Gets nucleation rate (ccn radius is the monomer !)
     
    740741    Us = (2._mm_wp * rad**2 * rho * mm_effg(z)) / (9._mm_wp * mm_eta_g(t))
    741742   
    742     ! Computes settling velocity (correction factor : x2.0)
    743     w = Us * Fc * 2._mm_wp
     743    ! Computes settling velocity (correction factor : x3.0)
     744    w = Us * Fc * 3._mm_wp
    744745  END FUNCTION wsettle
    745746
  • trunk/LMDZ.TITAN/libf/muphytitan/mm_globals.f90

    r3090 r3318  
    227227  REAL(kind=mm_wp), PARAMETER :: mm_rgas = mm_kboltz * mm_navo
    228228  !> Desorption energy (\(J\)) (nucleation).
    229   REAL(kind=mm_wp), PARAMETER :: mm_fdes = 0.288e-19_mm_wp
     229  REAL(kind=mm_wp), PARAMETER :: mm_fdes = 1.519e-20_mm_wp
    230230  !> Surface diffusion energy (\(J\)) (nucleation).
    231   REAL(kind=mm_wp), PARAMETER :: mm_fdif = 0.288e-20_mm_wp
     231  REAL(kind=mm_wp), PARAMETER :: mm_fdif = 1.519e-21_mm_wp
    232232  !> Jump frequency (\(s^{-1}\)) (nucleation).
    233233  REAL(kind=mm_wp), PARAMETER :: mm_nus = 1.e+13_mm_wp
     
    429429  REAL(kind=mm_wp), DIMENSION(:), ALLOCATABLE, SAVE :: mm_drho
    430430
    431   !> Aerosols precipitations (m).
     431  !> Aerosols precipitations (kg.m-2.s-1).
    432432  !!
    433433  !! Aerosols precipitations take into account both spherical and fractal modes.
     
    435435  REAL(kind=mm_wp), SAVE :: mm_aer_prec = 0._mm_wp
    436436
    437   !> CCN precipitations (m).
     437  !> CCN precipitations (kg.m-2.s-1).
    438438  !! It is updated in [[mm_clouds(module):mm_cloud_microphysics(subroutine)]].
    439439  REAL(kind=mm_wp), SAVE :: mm_ccn_prec = 0._mm_wp
     
    505505  REAL(kind=mm_wp), DIMENSION(:), ALLOCATABLE, SAVE :: mm_ccn_flux
    506506
    507   !> Ice components precipitations (m).
     507  !> Ice components precipitations (kg.m-2.s-1).
    508508  !!
    509509  !! It is a vector of [[mm_globals(module):mm_nesp(variable)]] values which share the same indexing
     
    14311431    ! Initialization :
    14321432    Ntot = m0ccn
    1433     Vtot = pifac*m3ccn + SUM(m3ice)
    1434     Wtot = pifac*m3ccn*mm_rhoaer + SUM(m3ice*mm_xESPS(:)%rho)
     1433    Vtot = pifac*m3ccn + pifac*SUM(m3ice)
     1434    Wtot = pifac*m3ccn*mm_rhoaer + pifac*SUM(m3ice*mm_xESPS(:)%rho)
    14351435
    14361436    IF (Ntot <= mm_m0n_min .OR. Vtot <= mm_m3cld_min) THEN
  • trunk/LMDZ.TITAN/libf/muphytitan/mm_haze.f90

    r3090 r3318  
    117117
    118118      ! Computes precipitations
    119       mm_aer_prec = SUM(zdm3as*mm_dzlev) + SUM(zdm3af*mm_dzlev)
     119      mm_aer_prec = SUM(zdm3as*mm_dzlev*mm_rhoaer) + SUM(zdm3af*mm_dzlev*mm_rhoaer)
    120120
    121121      ! Updates tendencies
  • trunk/LMDZ.TITAN/libf/muphytitan/mm_microphysic.f90

    r3090 r3318  
    163163  END FUNCTION muphys_nocld
    164164
    165   SUBROUTINE mm_diagnostics(aer_prec,aer_s_w,aer_f_w,aer_s_flux,aer_f_flux,ccn_prec,ccn_w,ccn_flux,ice_prec,ice_fluxes,gazs_sat)
     165  SUBROUTINE mm_diagnostics(dt,aer_prec,aer_s_w,aer_f_w,aer_s_flux,aer_f_flux,ccn_prec,ccn_w,ccn_flux,ice_prec,ice_fluxes,gazs_sat)
    166166    !! Get various diagnostic fields of the microphysics.
    167167    !!
     
    185185    !! __ccnprec__, __iceprec__, __icefluxes__ and __gazsat__ are always set to 0 if clouds
    186186    !! microphysics is disabled (see [[mm_globals(module):mm_w_clouds(variable)]] documentation).
    187     REAL(kind=mm_wp), INTENT(out), OPTIONAL                 :: aer_prec   !! Aerosols precipitations (both modes) (m).
    188     REAL(kind=mm_wp), INTENT(out), OPTIONAL                 :: ccn_prec   !! CCN precipitations (m).
     187    REAL(kind=8), INTENT(IN)                                :: dt         !! Physics timestep (s).
     188    REAL(kind=mm_wp), INTENT(out), OPTIONAL                 :: aer_prec   !! Aerosols precipitations (both modes) (kg.m-2.s-1).
     189    REAL(kind=mm_wp), INTENT(out), OPTIONAL                 :: ccn_prec   !! CCN precipitations (kg.m-2.s-1).
    189190    REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: aer_s_w    !! Spherical aerosol settling velocity (\(m.s^{-1}\)).
    190191    REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: aer_f_w    !! Fractal aerosol settling velocity (\(m.s^{-1}\)).
     
    195196    REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:,:) :: ice_fluxes !! Ice sedimentation fluxes (\(kg.m^{-2}.s^{-1}\)).
    196197    REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:,:) :: gazs_sat   !! Condensible gaz saturation ratios (--).
    197     REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: ice_prec   !! Ice precipitations (m).
    198 
    199     IF (PRESENT(aer_prec))   aer_prec   = ABS(mm_aer_prec)
     198    REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: ice_prec   !! Ice precipitations (kg.m-2.s-1).
     199
     200    IF (PRESENT(aer_prec))   aer_prec   = ABS(mm_aer_prec) / dt
    200201    IF (PRESENT(aer_s_w))    aer_s_w    = -mm_m3as_vsed(mm_nla:1:-1)
    201202    IF (PRESENT(aer_f_w))    aer_f_w    = -mm_m3af_vsed(mm_nla:1:-1)
     
    204205
    205206    IF (mm_w_clouds) THEN
    206       IF (PRESENT(ccn_prec))   ccn_prec   = ABS(mm_ccn_prec)
    207       IF (PRESENT(ice_prec))   ice_prec   = ABS(mm_ice_prec)
     207      IF (PRESENT(ccn_prec))   ccn_prec   = ABS(mm_ccn_prec) / dt
     208      IF (PRESENT(ice_prec))   ice_prec   = ABS(mm_ice_prec) / dt
    208209      IF (PRESENT(ccn_w))      ccn_w      = mm_ccn_vsed(mm_nla:1:-1)
    209210      IF (PRESENT(ccn_flux))   ccn_flux   = mm_ccn_flux(mm_nla:1:-1)
  • trunk/LMDZ.TITAN/libf/phytitan/calc_ysat.F90

    r3090 r3318  
    8282      else if(trim(cnames(ic)).eq."C2H4") then
    8383         ! not far from Fray and Schmidt (2009)
    84          where (temp(:,:).lt.89.0)
    85             ysat(:,:,ic) = 10.0**(1.5477e0 + (1.0e0/temp(:,:) - 0.011e0)                     &
    86                   * (16537.0e0*(1.0e0/temp(:,:) - 0.011e0) - 1038.1e0))                      &
    87                   / press(:,:) * 1.01325e0 / 760.0
    88          elsewhere (temp(:,:).lt.104.0)
    89             ysat(:,:,ic) = 10.0**(8.724e0 - 901.6e0/(temp(:,:) - 2.555e0) )                  &
    90                   / press(:,:) * 1013.25e0 / 760.0
    91          elsewhere (temp(:,:).lt.120.0)
    92             ysat(:,:,ic) = 10.0**(50.79e0 - 1703.0e0/temp(:,:) - 17.141e0                    &
    93                   * alog10(temp(:,:)) ) / press(:,:) * 1013.25e0 / 760.0
    94          elsewhere (temp(:,:).lt.155.0)
    95             ysat(:,:,ic) = 10.0**(6.74756e0 - 585.0e0/(temp(:,:) - 18.16e0) )                &
    96                   / press(:,:) * 1013.25e0 / 760.0
    97          endwhere
     84         !where (temp(:,:).lt.89.0)
     85         !   ysat(:,:,ic) = 10.0**(1.5477e0 + (1.0e0/temp(:,:) - 0.011e0)                     &
     86         !         * (16537.0e0*(1.0e0/temp(:,:) - 0.011e0) - 1038.1e0))                      &
     87         !         / press(:,:) * 1.01325e0 / 760.0
     88         !elsewhere (temp(:,:).lt.104.0)
     89         !   ysat(:,:,ic) = 10.0**(8.724e0 - 901.6e0/(temp(:,:) - 2.555e0) )                  &
     90         !         / press(:,:) * 1013.25e0 / 760.0
     91         !elsewhere (temp(:,:).lt.120.0)
     92         !   ysat(:,:,ic) = 10.0**(50.79e0 - 1703.0e0/temp(:,:) - 17.141e0                    &
     93         !         * alog10(temp(:,:)) ) / press(:,:) * 1013.25e0 / 760.0
     94         !elsewhere (temp(:,:).lt.155.0)
     95         !   ysat(:,:,ic) = 10.0**(6.74756e0 - 585.0e0/(temp(:,:) - 18.16e0) )                &
     96         !         / press(:,:) * 1013.25e0 / 760.0
     97         !endwhere
     98
     99         ! Fray and Schmidt (2009)
     100         ysat(:,:,ic) = (1.e3 / press(:,:)) * exp(1.540e1 - 2.206e3/temp(:,:)                &
     101                         - 1.216e4/temp(:,:)**2 + 2.843e5/temp(:,:)**3 - 2.203e6/temp(:,:)**4)
    98102     
    99103      ! C2H6 :
  • trunk/LMDZ.TITAN/libf/phytitan/calchim.F90

    r2368 r3318  
    6262
    6363  USE, INTRINSIC :: iso_c_binding
     64  use tracer_h
    6465  USE comchem_h
    6566  USE dimphy
     
    334335       ENDDO
    335336     endif
     337
     338     ! [BdBdT : On force la chimie... kedd nul dans le PCM]
     339     kedd(:klev) = 1.e3 ! Default value =/= zero
    336340 
    337341     firstcall = .FALSE.
     
    521525        DO ic=1,nkim
    522526           DO l=1,klev
    523               dqyc(ig,l,ic) = (cqy(l,ic) - cqy0(l,ic))/dtchim ! (mol/mol/s)
     527            dqyc(ig,l,ic) = (cqy(l,ic) - cqy0(l,ic))/dtchim ! (mol/mol/s)
    524528           ENDDO
    525529        ENDDO
  • trunk/LMDZ.TITAN/libf/phytitan/callcorrk.F90

    r3090 r3318  
    66          fluxabs_sw,fluxtop_dn,                               &
    77          OLR_nu,OSR_nu,                                       &
    8           int_dtaui,int_dtauv,popthi,popthv,poptci,poptcv,     &
     8          int_dtaui,int_dtauv,popthi,popthv,poptti,popttv,     &
    99          lastcall)
    1010
     
    8888      REAL,INTENT(OUT) :: int_dtauv(ngrid,nlayer,L_NSPECTV) ! VI optical thickness of layers within narrowbands for diags ().
    8989      ! Diagnostics
    90       REAL,INTENT(OUT) :: popthi(ngrid,nlayer,L_NSPECTI,5)  ! IR optical properties of haze within narrowbands (dtau,tau,k,w,g).
    91       REAL,INTENT(OUT) :: popthv(ngrid,nlayer,L_NSPECTV,5)  ! VI optical properties of haze within narrowbands (dtau,tau,k,w,g).
    92       REAL,INTENT(OUT) :: poptci(ngrid,nlayer,L_NSPECTI,5)  ! IR optical properties of haze and clouds within narrowbands (dtau,tau,k,w,g).
    93       REAL,INTENT(OUT) :: poptcv(ngrid,nlayer,L_NSPECTV,5)  ! VI optical properties of haze and clouds within narrowbands (dtau,tau,k,w,g).
     90      REAL,INTENT(OUT) :: popthi(ngrid,nlayer,L_NSPECTI,8)  ! IR optical properties of haze within narrowbands (dtau,tau,k,w,g,drayaer,taugaz,dcont).
     91      REAL,INTENT(OUT) :: popthv(ngrid,nlayer,L_NSPECTV,8)  ! VI optical properties of haze within narrowbands (dtau,tau,k,w,g,drayaer,taugaz,dcont).
     92      REAL,INTENT(OUT) :: poptti(ngrid,nlayer,L_NSPECTI,8)  ! IR optical properties of haze and clouds within narrowbands (dtau,tau,k,w,g,drayaer,taugaz,dcont).
     93      REAL,INTENT(OUT) :: popttv(ngrid,nlayer,L_NSPECTV,8)  ! VI optical properties of haze and clouds within narrowbands (dtau,tau,k,w,g,drayaer,taugaz,dcont).
    9494     
    9595     
     
    139139 
    140140      ! Optical diagnostics
    141       REAL*8 popt_hazei(L_LEVELS,L_NSPECTI,3)
    142       REAL*8 popt_hazev(L_LEVELS,L_NSPECTV,3)
    143       REAL*8 popt_cloudsi(L_LEVELS,L_NSPECTI,3)
    144       REAL*8 popt_cloudsv(L_LEVELS,L_NSPECTV,3)
     141      ! Haze
     142      REAL*8 diag_opthi(L_LEVELS,L_NSPECTI,6)
     143      REAL*8 diag_opthv(L_LEVELS,L_NSPECTV,6)
     144      ! Clouds
     145      REAL*8 diag_optti(L_LEVELS,L_NSPECTI,6)
     146      REAL*8 diag_opttv(L_LEVELS,L_NSPECTV,6)
    145147      ! Temporary optical diagnostics (clear column)
    146       REAL*8 popt_hazei_cc(L_LEVELS,L_NSPECTI,3)
    147       REAL*8 popt_hazev_cc(L_LEVELS,L_NSPECTV,3)
    148       REAL*8 popt_cloudsi_cc(L_LEVELS,L_NSPECTI,3)
    149       REAL*8 popt_cloudsv_cc(L_LEVELS,L_NSPECTV,3)
     148      REAL*8 diag_optti_cc(L_LEVELS,L_NSPECTI,6)
     149      REAL*8 diag_opttv_cc(L_LEVELS,L_NSPECTV,6)
    150150      ! Temporary optical diagnostics (dark column)
    151       REAL*8 popt_hazei_dc(L_LEVELS,L_NSPECTI,3)
    152       REAL*8 popt_hazev_dc(L_LEVELS,L_NSPECTV,3)
    153       REAL*8 popt_cloudsi_dc(L_LEVELS,L_NSPECTI,3)
    154       REAL*8 popt_cloudsv_dc(L_LEVELS,L_NSPECTV,3)
     151      REAL*8 diag_optti_dc(L_LEVELS,L_NSPECTI,6)
     152      REAL*8 diag_opttv_dc(L_LEVELS,L_NSPECTV,6)
    155153
    156154
     
    508506         call optcv(pqmo(ig,:,:),nlayer,zzlev(ig,:),plevrad,tmid,pmid,    &
    509507              dtauv_cc,tauv_cc,taucumv_cc,wbarv_cc,cosbv_cc,tauray,taugsurf,seashazefact,&
    510               popt_hazev_cc,popt_cloudsv_cc,cdcolumn)
     508              diag_opthv,diag_opttv_cc,cdcolumn)
    511509         ! Dark column :
    512510         cdcolumn = 1
    513511         call optcv(pqmo(ig,:,:),nlayer,zzlev(ig,:),plevrad,tmid,pmid,    &
    514512              dtauv_dc,tauv_dc,taucumv_dc,wbarv_dc,cosbv_dc,tauray,taugsurf,seashazefact,&
    515               popt_hazev_dc,popt_cloudsv_dc,cdcolumn)
     513              diag_opthv,diag_opttv_dc,cdcolumn)
    516514         
    517515         ! Mean opacity, ssa and asf :
     
    537535            end do
    538536         end do
    539          wbarv = (dtauv_cc*wbarv_cc + dtauv_dc*wbarv_dc) / (dtauv_cc + dtauv_dc)
    540          cosbv = (dtauv_cc*wbarv_cc*cosbv_cc + dtauv_dc*wbarv_dc*cosbv_dc) / (dtauv_cc*wbarv_cc + dtauv_dc*wbarv_dc)
     537         wbarv = (1-Fcloudy) * wbarv_cc + Fcloudy * wbarv_dc
     538         cosbv = (1-Fcloudy) * cosbv_cc + Fcloudy * cosbv_dc
    541539         
    542          ! Diagnostics for haze :
    543          where (popt_hazev_cc(:,:,1) .lt. 1.d-30)
    544             popt_hazev_cc(:,:,1) = 1.d-30
    545          endwhere
    546          where (popt_hazev_dc(:,:,1) .lt. 1.d-30)
    547             popt_hazev_dc(:,:,1) = 1.d-30
    548          endwhere
    549          popt_hazev(:,:,1) = -log( (1-Fcloudy)*exp(-popt_hazev_cc(:,:,1)) + Fcloudy*exp(-popt_hazev_dc(:,:,1)) )
    550          popt_hazev(:,:,2) = ((popt_hazev_cc(:,:,1)*popt_hazev_cc(:,:,2) + popt_hazev_dc(:,:,1)*popt_hazev_dc(:,:,2)) &
    551                              / (popt_hazev_cc(:,:,1) + popt_hazev_dc(:,:,1)))
    552          where ((popt_hazev_cc(:,:,2).gt.0.0D0) .or. (popt_hazev_dc(:,:,2).gt.0.0D0))
    553             popt_hazev(:,:,3) = (popt_hazev_cc(:,:,1)*popt_hazev_cc(:,:,2)*popt_hazev_cc(:,:,3) + popt_hazev_dc(:,:,1)*popt_hazev_dc(:,:,2)*popt_hazev_dc(:,:,3)) &
    554                              / (popt_hazev_cc(:,:,1)*popt_hazev_cc(:,:,2) + popt_hazev_dc(:,:,1)*popt_hazev_dc(:,:,2))
    555          elsewhere
    556             popt_hazev(:,:,3) = 0.0D0
    557          endwhere
    558540         ! Diagnostics for clouds :
    559541         if (callclouds) then
    560             where (popt_cloudsv_cc(:,:,1) .lt. 1.d-30)
    561                popt_cloudsv_cc(:,:,1) = 1.d-30
     542            where (diag_opttv_cc(:,:,1) .lt. 1.d-30)
     543               diag_opttv_cc(:,:,1) = 1.d-30
    562544            endwhere
    563             where (popt_cloudsv_dc(:,:,1) .lt. 1.d-30)
    564                popt_cloudsv_dc(:,:,1) = 1.d-30
     545            where (diag_opttv_dc(:,:,1) .lt. 1.d-30)
     546               diag_opttv_dc(:,:,1) = 1.d-30
    565547            endwhere
    566             popt_cloudsv(:,:,1) = -log( (1-Fcloudy)*exp(-popt_cloudsv_cc(:,:,1)) + Fcloudy*exp(-popt_cloudsv_dc(:,:,1)) )
    567             popt_cloudsv(:,:,2) = ((popt_cloudsv_cc(:,:,1)*popt_cloudsv_cc(:,:,2) + popt_cloudsv_dc(:,:,1)*popt_cloudsv_dc(:,:,2)) &
    568                               / (popt_cloudsv_cc(:,:,1) + popt_cloudsv_dc(:,:,1)))
    569             where ((popt_cloudsv_cc(:,:,2).gt.0.0D0) .or. (popt_cloudsv_dc(:,:,2).gt.0.0D0))
    570                popt_cloudsv(:,:,3) = (popt_cloudsv_cc(:,:,1)*popt_cloudsv_cc(:,:,2)*popt_cloudsv_cc(:,:,3) + popt_cloudsv_dc(:,:,1)*popt_cloudsv_dc(:,:,2)*popt_cloudsv_dc(:,:,3)) &
    571                               / (popt_cloudsv_cc(:,:,1)*popt_cloudsv_cc(:,:,2) + popt_cloudsv_dc(:,:,1)*popt_cloudsv_dc(:,:,2))
    572             elsewhere
    573                popt_cloudsv(:,:,3) = 0.0D0
    574             endwhere
     548            diag_opttv(:,:,1) = -log( (1-Fcloudy)*exp(-diag_opttv_cc(:,:,1)) + Fcloudy*exp(-diag_opttv_dc(:,:,1)) )
     549            diag_opttv(:,:,2) = (1-Fcloudy) * diag_opttv_cc(:,:,2) + Fcloudy * diag_opttv_dc(:,:,2)
     550            diag_opttv(:,:,3) = (1-Fcloudy) * diag_opttv_cc(:,:,3) + Fcloudy * diag_opttv_dc(:,:,3)
     551            diag_opttv(:,:,4) = -log( (1-Fcloudy)*exp(-diag_opttv_cc(:,:,4)) + Fcloudy*exp(-diag_opttv_dc(:,:,4)) )
     552            diag_opttv(:,:,5) = -log( (1-Fcloudy)*exp(-diag_opttv_cc(:,:,5)) + Fcloudy*exp(-diag_opttv_dc(:,:,5)) )
     553            diag_opttv(:,:,6) = -log( (1-Fcloudy)*exp(-diag_opttv_cc(:,:,6)) + Fcloudy*exp(-diag_opttv_dc(:,:,6)) )
    575554         endif
    576555         
     
    629608         call optci(pqmo(ig,:,:),nlayer,zzlev(ig,:),plevrad,tlevrad,tmid,pmid,&
    630609              dtaui_cc,taucumi_cc,cosbi_cc,wbari_cc,taugsurfi,seashazefact,   &
    631               popt_hazei_cc,popt_cloudsi_cc,cdcolumn)
     610              diag_opthi,diag_optti_cc,cdcolumn)
    632611         ! Dark column :
    633612         cdcolumn = 1
    634613         call optci(pqmo(ig,:,:),nlayer,zzlev(ig,:),plevrad,tlevrad,tmid,pmid,&
    635614              dtaui_dc,taucumi_dc,cosbi_dc,wbari_dc,taugsurfi,seashazefact,   &
    636               popt_hazei_dc,popt_cloudsi_dc,cdcolumn)
     615              diag_opthi,diag_optti_dc,cdcolumn)
    637616
    638617         ! Mean opacity, ssa and asf :
     
    654633            end do
    655634         end do
    656          where ((wbari_cc.gt.0.0D0) .or. (wbari_dc.gt.0.0D0))
    657             wbari = (dtaui_cc*wbari_cc + dtaui_dc*wbari_dc) / (dtaui_cc + dtaui_dc)
    658             cosbi = (dtaui_cc*wbari_cc*cosbi_cc + dtaui_dc*wbari_dc*cosbi_dc) / (dtaui_cc*wbari_cc + dtaui_dc*wbari_dc)
    659          elsewhere
    660             wbari = 0.0D0
    661             cosbi = 0.0D0
    662          endwhere
    663 
    664          ! Diagnostics for haze :
    665          where (popt_hazei_cc(:,:,1) .lt. 1.d-30)
    666             popt_hazei_cc(:,:,1) = 1.d-30
    667          endwhere
    668          where (popt_hazei_dc(:,:,1) .lt. 1.d-30)
    669             popt_hazei_dc(:,:,1) = 1.d-30
    670          endwhere
    671          popt_hazei(:,:,1) = -log( (1-Fcloudy)*exp(-popt_hazei_cc(:,:,1)) + Fcloudy*exp(-popt_hazei_dc(:,:,1)) )
    672          popt_hazei(:,:,2) = (popt_hazei_cc(:,:,1)*popt_hazei_cc(:,:,2) + popt_hazei_dc(:,:,1)*popt_hazei_dc(:,:,2)) &
    673                              / (popt_hazei_cc(:,:,1) + popt_hazei_dc(:,:,1))
    674          where ((popt_hazei_cc(:,:,2).gt.0.0D0) .or. (popt_hazei_dc(:,:,2).gt.0.0D0))
    675             popt_hazei(:,:,3) = (popt_hazei_cc(:,:,1)*popt_hazei_cc(:,:,2)*popt_hazei_cc(:,:,3) + popt_hazei_dc(:,:,1)*popt_hazei_dc(:,:,2)*popt_hazei_dc(:,:,3)) &
    676                              / (popt_hazei_cc(:,:,1)*popt_hazei_cc(:,:,2) + popt_hazei_dc(:,:,1)*popt_hazei_dc(:,:,2))
    677          elsewhere
    678             popt_hazei(:,:,3) = 0.0D0
    679          endwhere
     635         wbari = (1-Fcloudy) * wbari_cc + Fcloudy * wbari_dc
     636         cosbi = (1-Fcloudy) * cosbi_cc + Fcloudy * cosbi_dc
     637
    680638         ! Diagnostics for clouds :
    681639         if (callclouds) then
    682             where (popt_cloudsi_cc(:,:,1) .lt. 1.d-30)
    683                popt_cloudsi_cc(:,:,1) = 1.d-30
     640            where (diag_optti_cc(:,:,1) .lt. 1.d-30)
     641               diag_optti_cc(:,:,1) = 1.d-30
    684642            endwhere
    685             where (popt_cloudsi_dc(:,:,1) .lt. 1.d-30)
    686                popt_cloudsi_dc(:,:,1) = 1.d-30
     643            where (diag_optti_dc(:,:,1) .lt. 1.d-30)
     644               diag_optti_dc(:,:,1) = 1.d-30
    687645            endwhere
    688             popt_cloudsi(:,:,1) = -log( (1-Fcloudy)*exp(-popt_cloudsi_cc(:,:,1)) + Fcloudy*exp(-popt_cloudsi_dc(:,:,1)) )
    689             popt_cloudsi(:,:,2) = ((popt_cloudsi_cc(:,:,1)*popt_cloudsi_cc(:,:,2) + popt_cloudsi_dc(:,:,1)*popt_cloudsi_dc(:,:,2)) &
    690                               / (popt_cloudsi_cc(:,:,1) + popt_cloudsi_dc(:,:,1)))
    691             where ((popt_cloudsi_cc(:,:,2).gt.0.0D0) .or. (popt_cloudsi_dc(:,:,2).gt.0.0D0))
    692                popt_cloudsi(:,:,3) = (popt_cloudsi_cc(:,:,1)*popt_cloudsi_cc(:,:,2)*popt_cloudsi_cc(:,:,3) + popt_cloudsi_dc(:,:,1)*popt_cloudsi_dc(:,:,2)*popt_cloudsi_dc(:,:,3)) &
    693                               / (popt_cloudsi_cc(:,:,1)*popt_cloudsi_cc(:,:,2) + popt_cloudsi_dc(:,:,1)*popt_cloudsi_dc(:,:,2))
    694             elsewhere
    695                popt_cloudsi(:,:,3) = 0.0D0
    696             endwhere
     646            diag_optti(:,:,1) = -log( (1-Fcloudy)*exp(-diag_optti_cc(:,:,1)) + Fcloudy*exp(-diag_optti_dc(:,:,1)) )
     647            diag_optti(:,:,2) = (1-Fcloudy) * diag_optti_cc(:,:,2) + Fcloudy * diag_optti_dc(:,:,2)
     648            diag_optti(:,:,3) = (1-Fcloudy) * diag_optti_cc(:,:,3) + Fcloudy * diag_optti_dc(:,:,3)
     649            diag_optti(:,:,4) = -log( (1-Fcloudy)*exp(-diag_optti_cc(:,:,4)) + Fcloudy*exp(-diag_optti_dc(:,:,4)) )
     650            diag_optti(:,:,5) = -log( (1-Fcloudy)*exp(-diag_optti_cc(:,:,5)) + Fcloudy*exp(-diag_optti_dc(:,:,5)) )
     651            diag_optti(:,:,6) = -log( (1-Fcloudy)*exp(-diag_optti_cc(:,:,6)) + Fcloudy*exp(-diag_optti_dc(:,:,6)) )
    697652         endif
    698653
     
    779734         do nw = 1, L_NSPECTV
    780735            popthv(ig,lev2lay,nw,:) = 0.0d0
    781             poptcv(ig,lev2lay,nw,:) = 0.0d0
     736            popttv(ig,lev2lay,nw,:) = 0.0d0
    782737            ! Optical thickness (dtau) :
    783             popthv(ig,lev2lay,nw,1) = (popt_hazev(l,nw,1) + popt_hazev(l+1,nw,1)) / 2.0
     738            popthv(ig,lev2lay,nw,1) = (diag_opthv(l,nw,1) + diag_opthv(l+1,nw,1)) / 2.0
    784739            if (callclouds) then
    785                poptcv(ig,lev2lay,nw,1) = (popt_cloudsv(l,nw,1) + popt_cloudsv(l+1,nw,1)) / 2.0
     740               popttv(ig,lev2lay,nw,1) = (diag_opttv(l,nw,1) + diag_opttv(l+1,nw,1)) / 2.0
    786741            endif
    787742            ! Opacity (tau) :
     
    789744               popthv(ig,lev2lay,nw,2) = popthv(ig,lev2lay,nw,2) + popthv(ig,k,nw,1)
    790745               if (callclouds) then
    791                   poptcv(ig,lev2lay,nw,2) = poptcv(ig,lev2lay,nw,2) + poptcv(ig,k,nw,1)
     746                  popttv(ig,lev2lay,nw,2) = popttv(ig,lev2lay,nw,2) + popttv(ig,k,nw,1)
    792747               endif
    793748            enddo
     
    795750            popthv(ig,lev2lay,nw,3) = popthv(ig,lev2lay,nw,1) / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
    796751            if (callclouds) then
    797                poptcv(ig,lev2lay,nw,3) = poptcv(ig,lev2lay,nw,1) / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     752               popttv(ig,lev2lay,nw,3) = popttv(ig,lev2lay,nw,1) / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
    798753            endif
    799754            ! Simple Scattering Albedo (w) and Asymmetry Parameter (g) :
    800             popthv(ig,lev2lay,nw,4) = (popt_hazev(l,nw,2) + popt_hazev(l+1,nw,2)) / 2.0
    801             popthv(ig,lev2lay,nw,5) = (popt_hazev(l,nw,3) + popt_hazev(l+1,nw,3)) / 2.0
     755            popthv(ig,lev2lay,nw,4) = (diag_opthv(l,nw,2) + diag_opthv(l+1,nw,2)) / 2.0
     756            popthv(ig,lev2lay,nw,5) = (diag_opthv(l,nw,3) + diag_opthv(l+1,nw,3)) / 2.0
    802757            if (callclouds) then
    803                poptcv(ig,lev2lay,nw,4) = (popt_cloudsv(l,nw,2) + popt_cloudsv(l+1,nw,2)) / 2.0
    804                poptcv(ig,lev2lay,nw,5) = (popt_cloudsv(l,nw,3) + popt_cloudsv(l+1,nw,3)) / 2.0
     758               popttv(ig,lev2lay,nw,4) = (diag_opttv(l,nw,2) + diag_opttv(l+1,nw,2)) / 2.0
     759               popttv(ig,lev2lay,nw,5) = (diag_opttv(l,nw,3) + diag_opttv(l+1,nw,3)) / 2.0
     760            endif
     761            ! DRAYAER, TAUGAS, DCONT :
     762            popthv(ig,lev2lay,nw,6) = (diag_opthv(l,nw,4) + diag_opthv(l+1,nw,4)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     763            popthv(ig,lev2lay,nw,7) = (diag_opthv(l,nw,5) + diag_opthv(l+1,nw,5)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     764            popthv(ig,lev2lay,nw,8) = (diag_opthv(l,nw,6) + diag_opthv(l+1,nw,6)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     765            if (callclouds) then
     766               popttv(ig,lev2lay,nw,6) = (diag_opttv(l,nw,4) + diag_opttv(l+1,nw,4)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     767               popttv(ig,lev2lay,nw,7) = (diag_opttv(l,nw,5) + diag_opttv(l+1,nw,5)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     768               popttv(ig,lev2lay,nw,8) = (diag_opttv(l,nw,6) + diag_opttv(l+1,nw,6)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
    805769            endif
    806770         enddo
     
    808772         do nw=1,L_NSPECTI
    809773            popthi(ig,lev2lay,nw,:) = 0.0d0
    810             poptci(ig,lev2lay,nw,:) = 0.0d0
     774            poptti(ig,lev2lay,nw,:) = 0.0d0
    811775            ! Optical thickness (dtau) :
    812             popthi(ig,lev2lay,nw,1) = (popt_hazei(l,nw,1) + popt_hazei(l+1,nw,1)) / 2.0
     776            popthi(ig,lev2lay,nw,1) = (diag_opthi(l,nw,1) + diag_opthi(l+1,nw,1)) / 2.0
    813777            if (callclouds) then
    814                poptci(ig,lev2lay,nw,1) = (popt_cloudsi(l,nw,1) + popt_cloudsi(l+1,nw,1)) / 2.0
     778               poptti(ig,lev2lay,nw,1) = (diag_optti(l,nw,1) + diag_optti(l+1,nw,1)) / 2.0
    815779            endif
    816780            ! Opacity (tau) :
     
    818782               popthi(ig,lev2lay,nw,2) = popthi(ig,lev2lay,nw,2) + popthi(ig,k,nw,1)
    819783               if (callclouds) then
    820                   poptci(ig,lev2lay,nw,2) = poptci(ig,lev2lay,nw,2) + poptci(ig,k,nw,1)
     784                  poptti(ig,lev2lay,nw,2) = poptti(ig,lev2lay,nw,2) + poptti(ig,k,nw,1)
    821785               endif
    822786            enddo
     
    824788            popthi(ig,lev2lay,nw,3) = popthi(ig,lev2lay,nw,1) / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
    825789            if (callclouds) then
    826                poptci(ig,lev2lay,nw,3) = poptci(ig,lev2lay,nw,1) / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     790               poptti(ig,lev2lay,nw,3) = poptti(ig,lev2lay,nw,1) / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
    827791            endif
    828792            ! Simple Scattering Albedo (w) and Asymmetry Parameter (g) :
    829             popthi(ig,lev2lay,nw,4) = (popt_hazei(l,nw,2) + popt_hazei(l+1,nw,2)) / 2.0
    830             popthi(ig,lev2lay,nw,5) = (popt_hazei(l,nw,3) + popt_hazei(l+1,nw,3)) / 2.0
     793            popthi(ig,lev2lay,nw,4) = (diag_opthi(l,nw,2) + diag_opthi(l+1,nw,2)) / 2.0
     794            popthi(ig,lev2lay,nw,5) = (diag_opthi(l,nw,3) + diag_opthi(l+1,nw,3)) / 2.0
    831795            if (callclouds) then
    832                poptci(ig,lev2lay,nw,4) = (popt_cloudsi(l,nw,2) + popt_cloudsi(l+1,nw,2)) / 2.0
    833                poptci(ig,lev2lay,nw,5) = (popt_cloudsi(l,nw,3) + popt_cloudsi(l+1,nw,3)) / 2.0
     796               poptti(ig,lev2lay,nw,4) = (diag_optti(l,nw,2) + diag_optti(l+1,nw,2)) / 2.0
     797               poptti(ig,lev2lay,nw,5) = (diag_optti(l,nw,3) + diag_optti(l+1,nw,3)) / 2.0
     798            endif
     799            ! DRAYAER, TAUGAS, DCONT :
     800            popthi(ig,lev2lay,nw,6) = (diag_opthi(l,nw,4) + diag_opthi(l+1,nw,4)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     801            popthi(ig,lev2lay,nw,7) = (diag_opthi(l,nw,5) + diag_opthi(l+1,nw,5)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     802            popthi(ig,lev2lay,nw,8) = (diag_opthi(l,nw,6) + diag_opthi(l+1,nw,6)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     803            if (callclouds) then
     804               poptti(ig,lev2lay,nw,6) = (diag_optti(l,nw,4) + diag_optti(l+1,nw,4)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     805               poptti(ig,lev2lay,nw,7) = (diag_optti(l,nw,5) + diag_optti(l+1,nw,5)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
     806               poptti(ig,lev2lay,nw,8) = (diag_optti(l,nw,6) + diag_optti(l+1,nw,6)) / 2.0 / (zzlev(ig,lev2lay+1) - zzlev(ig,lev2lay))
    834807            endif
    835808         enddo
  • trunk/LMDZ.TITAN/libf/phytitan/callkeys_mod.F90

    r3090 r3318  
    5454!$OMP THREADPRIVATE(opt4clouds)
    5555      real,save    :: Fcloudy
    56       real,save    :: Fssa
     56      real,save    :: FCVIS
     57      real,save    :: FCIR
    5758      real,save    :: FHVIS
    5859      real,save    :: FHIR
    59 !$OMP THREADPRIVATE(Fcloudy,Fssa,FHVIS,FHIR)
     60!$OMP THREADPRIVATE(Fcloudy,FCVIS,FCIR,FHVIS,FHIR)
    6061     
    6162      logical,save :: resCH4_inf
  • trunk/LMDZ.TITAN/libf/phytitan/calmufi.F90

    r3090 r3318  
    118118      m3as(:) = zq(ilon,:,2) * int2ext(ilon,:)                  ! then mixed, even though both are above the threshold, their ratio can be nonsense
    119119    ELSEWHERE
    120       m0as(:)=0.D0
    121       m3as(:)=0.D0
     120      m0as(:) = 0.D0
     121      m3as(:) = 0.D0
    122122    ENDWHERE
    123123    WHERE (pq(ilon,:,3) > 2.D-200 .AND. pq(ilon,:,4) > 2.D-200)
     
    125125      m3af(:) = zq(ilon,:,4) * int2ext(ilon,:)
    126126    ELSEWHERE
    127       m0af(:)=0.D0
    128       m3af(:)=0.D0
     127      m0af(:) = 0.D0
     128      m3af(:) = 0.D0
    129129    ENDWHERE
    130130   
     
    168168    ENDIF
    169169    ! save diags (if no clouds, relevant arrays will be set to 0 !)
    170     call mm_diagnostics(mmd_aer_prec(ilon),mmd_aer_s_w(ilon,:),mmd_aer_f_w(ilon,:),mmd_aer_s_flux(ilon,:),mmd_aer_f_flux(ilon,:),  &
     170    call mm_diagnostics(dt,mmd_aer_prec(ilon),mmd_aer_s_w(ilon,:),mmd_aer_f_w(ilon,:),mmd_aer_s_flux(ilon,:),mmd_aer_f_flux(ilon,:),  &
    171171                        mmd_ccn_prec(ilon),mmd_ccn_w(ilon,:),mmd_ccn_flux(ilon,:),mmd_ice_prec(ilon,:),   &
    172172                        mmd_ice_fluxes(ilon,:,:),mmd_gazs_sat(ilon,:,:))
  • trunk/LMDZ.TITAN/libf/phytitan/datafile_mod.F90

    r3083 r3318  
    3333      ! Set in physiq_mod
    3434      character(LEN=100),save :: nudging_file ='datagcm/SuperNudging.dat'
     35      ! character(LEN=100),save :: nudging_file ='datagcm/NewNudging.dat'
    3536!$OMP THREADPRIVATE(nudging_file)
    3637
  • trunk/LMDZ.TITAN/libf/phytitan/dyn1d/rcm1d.F

    r2366 r3318  
    728728! Initialize methane surface tank
    729729! -------------------------------
    730       tankCH4=2. ! default value for tankCH4
     730      tankCH4=200. ! default value for tankCH4
    731731     
    732732! Initialize tracers
  • trunk/LMDZ.TITAN/libf/phytitan/evapCH4.F90

    r3090 r3318  
    7676! Parameters :
    7777!-------------
    78 real, parameter :: karman = 0.4     ! Karman constant [-]
    79 real, parameter :: humCH4 = 0.5     ! Imposed surface humidity for CH4 [-]
    80 
    81 real, parameter :: Flnp = 0.05      ! Fraction occupied by lakes (North Pole)
    82 real, parameter :: Flsp = 0.005     ! Fraction occupied by lakes (South Pole)
     78real, parameter :: karman = 0.4 ! Karman constant [-]
     79real, parameter :: humCH4 = 0.4 ! Imposed surface humidity for CH4 [-]
     80
     81real, parameter :: Flnp = 0.08  ! Fraction occupied by lakes (North Pole)
     82real, parameter :: Flsp = 0.02  ! Fraction occupied by lakes (South Pole)
     83real, parameter :: Flml = 0.02  ! Fraction not infiltrated into the ground (Mid latitudes)
    8384
    8485real, parameter :: mmolair = 28.e-3 ! Molar mass of air [kg.mol-1]
     
    9495
    9596! Initialisation :
     97real*8  :: Tlake      ! Lake temperature [K]
    9698real*8  :: rhoair     ! Density of air [kg.m-3]
    9799real*8  :: ws         ! Horizontal wind at the surface [m.s-1]
     
    127129
    128130  ! Saturation profile of CH4 [mol.mol-1] (Fray and Schmidt 2009) :
    129   qsatCH4 = (1.0e5 / pplev(ig,1)) * exp(1.051e1 - 1.110e3/tsurf(ig) - 4.341e3/tsurf(ig)**2 + 1.035e5/tsurf(ig)**3 - 7.910e5/tsurf(ig)**4)
    130   qsatCH4 = humCH4 * qsatCH4
     131  Tlake = tsurf(ig) - 7 ! Lakes are 2-7K less than surface.
     132  qsatCH4 = (1.0e5 / pplev(ig,1)) * exp(1.051e1 - 1.110e3/Tlake - 4.341e3/Tlake**2 + 1.035e5/Tlake**3 - 7.910e5/Tlake**4)
    131133  ! CH4 : 0.80 * qsat because of dissolution in N2
    132134  qsatCH4 = 0.80 * qsatCH4
     
    134136  ! Flux at the surface [kg.m-2.s-1] :
    135137  flux = rhoair * Cd * ws
     138 
     139  ! Surface humidity :
     140  qsatCH4 = humCH4 * qsatCH4
    136141
    137142  ! <North Pole>
    138   if (REAL(latitude_deg(ig)) .ge. 70.) then
     143  if (REAL(latitude_deg(ig)) .ge. 70.0) then
    139144    flux = Flnp * flux
    140145  ! <South Pole>
    141   else if (REAL(latitude_deg(ig)) .le. -70.) then
     146  else if (REAL(latitude_deg(ig)) .le. -70.0) then
    142147    flux = Flsp * flux
    143148  ! <Mid Latitudes>
    144149  else
    145     flux = 0.0
     150    flux = Flml * flux
    146151  endif
    147152
     
    151156    tankCH4(ig) = 1.e-30
    152157  endif
    153  
     158
     159!------------------------------------
     160! 2. IMPLICIT SCHEME
     161!------------------------------------
     162
    154163  ! Flux of CH4 at the surface [kg.m-2.s-1.mol.mol-1] :
    155164  fluxCH4 = flux * (qsatCH4 - pqCH4(ig,1))
    156165
    157 
    158 !------------------------------------
    159 ! 2. IMPLICIT SCHEME
    160 !------------------------------------
    161 
    162166  ! Flux at the surface [kg.m-2.s-1] --> [s-1] :
    163167  flux = flux / rhoair / (zzlev(ig,2) - zzlev(ig,1))
     
    170174
    171175  ! New tank depth of CH4 [m]
    172   if ((tankCH4(ig) + dtankCH4) .ge. 0.) then
     176  if ((tankCH4(ig) + dtankCH4) .ge. 1.e-30) then
    173177    tankCH4(ig) = tankCH4(ig) + dtankCH4
    174178  else
     179    !write(*,*) 'Evaporation CH4 : Empty lakes...'
    175180    newpqCH4 = tankCH4(ig) / (rhoair * (zzlev(ig,2) - zzlev(ig,1)) * mmolCH4 / mmolair / rhoiCH4) + pqCH4(ig,1)
    176181    tankCH4(ig) = 1.e-30
     
    185190!------------------------------------
    186191
    187   ftm = (1. - tsurf(ig) / TcCH4)
     192  ftm = (1. - Tlake / TcCH4)
    188193  if(ftm .le. 1.e-3) then
    189194    ftm = 1.e-3
     
    191196 
    192197  ! Latente heat of CH4 vaporisation [J.kg-1.mol.mol-1]
    193   LvCH4 = 8.314 * 190.4 * (7.08 * ftm**0.354 + 10.95 * 1.1e-2 * ftm**0.456) / mmolCH4
     198  LvCH4 = 8.314 * TcCH4 * (7.08 * ftm**0.354 + 10.95 * 1.1e-2 * ftm**0.456) / mmolCH4
    194199
    195200  ! Evaporation heating rate [K.s-1]
  • trunk/LMDZ.TITAN/libf/phytitan/inifis_mod.F90

    r3090 r3318  
    218218     write(*,*) "Fcloudy = ", Fcloudy
    219219
    220      write(*,*) "Adjustment of single scattering albedo for cloudy dropplet"
    221      Fssa = 1.0
    222      call getin_p("Fssa", Fssa)
    223      write(*,*) "Fssa = ", Fssa
     220     write(*,*) "Adjustment of droplet properties in the VIS"
     221     FCVIS = 1.0
     222     call getin_p("FCVIS", FCVIS)
     223     write(*,*) "FCVIS = ", FCVIS
     224
     225     write(*,*) "Adjustment of droplet properties in the IR"
     226     FCIR = 1.0
     227     call getin_p("FCIR", FCIR)
     228     write(*,*) "FCIR = ", FCIR
    224229     
    225230     write(*,*) "Adjustment of aerosol properties in the VIS"
  • trunk/LMDZ.TITAN/libf/phytitan/optci.F90

    r3090 r3318  
    11subroutine optci(PQMO,NLAY,ZLEV,PLEV,TLEV,TMID,PMID, &
    22     DTAUI,TAUCUMI,COSBI,WBARI,TAUGSURF,SEASHAZEFACT,&
    3      POPT_HAZE,POPT_CLOUDS,CDCOLUMN)
     3     DIAG_OPTH,DIAG_OPTT,CDCOLUMN)
    44
    55  use radinc_h
     
    1111  use callkeys_mod, only: continuum,graybody,corrk_recombin,              &
    1212                          callclouds,callmufi,seashaze,uncoupl_optic_haze,&
    13                           opt4clouds,FHIR
     13                          opt4clouds,FHIR,FCIR
    1414  use tracer_h, only: nmicro,nice,ices_indx
    1515
     
    6161  REAL*8, INTENT(OUT) :: WBARI(L_NLAYRAD,L_NSPECTI,L_NGAUSS)
    6262  REAL*8, INTENT(OUT) :: TAUGSURF(L_NSPECTI,L_NGAUSS-1)
    63   REAL*8, INTENT(OUT) :: POPT_HAZE(L_LEVELS,L_NSPECTI,3)
    64   REAL*8, INTENT(OUT) :: POPT_CLOUDS(L_LEVELS,L_NSPECTI,3)
     63  REAL*8, INTENT(OUT) :: DIAG_OPTH(L_LEVELS,L_NSPECTI,6)
     64  REAL*8, INTENT(OUT) :: DIAG_OPTT(L_LEVELS,L_NSPECTI,6)
    6565  ! ==========================================================
    6666 
     
    230230            ENDIF
    231231
     232            ! Tuning of optical properties for haze :
     233            !dtauaer_s = dtauaer_s * (FHIR * (1-ssa_s(nw)) + ssa_s(nw))
     234            !ssa_s(nw) = ssa_s(nw) / (FHIR * (1-ssa_s(nw)) + ssa_s(nw))
     235            dtauaer_f = dtauaer_f * (FHIR * (1-ssa_f(nw)) + ssa_f(nw))
     236            ssa_f(nw) = ssa_f(nw) / (FHIR * (1-ssa_f(nw)) + ssa_f(nw))
     237
    232238            ! Total of Haze opacity (dtau), SSA (w) and ASF (COS) :
    233239            DHAZE_T(k,nw) = dtauaer_s + dtauaer_f
     
    241247               ASF_T(k,nw)   = 1.0
    242248            ENDIF
     249           
    243250            ! Diagnostics for the haze :
    244             POPT_HAZE(k,nw,1) = DHAZE_T(k,nw) ! dtau
    245             POPT_HAZE(k,nw,2) = SSA_T(k,nw)   ! wbar
    246             POPT_HAZE(k,nw,3) = ASF_T(k,nw)   ! gbar
     251            DIAG_OPTH(k,nw,1) = DHAZE_T(k,nw) ! dtau
     252            DIAG_OPTH(k,nw,2) = SSA_T(k,nw)   ! wbar
     253            DIAG_OPTH(k,nw,3) = ASF_T(k,nw)   ! gbar
    247254
    248255            !---------------------
     
    253260               m0ccn = pqmo(ilay,5) / 2.0
    254261               m3ccn = pqmo(ilay,6) / 2.0
    255                m3cld  = 0.0d0
     262               m3cld  = pqmo(ilay,6) / 2.0
    256263               
    257264               ! Clear / Dark column method :
     
    285292               
    286293               ! For small dropplets, opacity of nucleus dominates...
     294               dtau_cld = (dtau_cld*m3ccn + dtau_cld*m3cld) / (m3ccn + m3cld)
    287295               ssa_cld(nw) = (ssa_ccn(nw)*m3ccn + ssa_cld(nw)*m3cld) / (m3ccn + m3cld)
    288 
     296               
     297               ! Tuning of optical properties for clouds :
     298               dtau_cld = dtau_cld * (FCIR * (1-ssa_cld(nw)) + ssa_cld(nw))
     299               ssa_cld(nw) = ssa_cld(nw) / (FCIR * (1-ssa_cld(nw)) + ssa_cld(nw))
     300               
    289301               ! Total of Haze + Clouds opacity (dtau), SSA (w) and ASF (COS) :
    290302               DHAZE_T(k,nw) = dtauaer_s + dtauaer_f + dtau_cld
     
    299311               ENDIF
    300312
    301                ! Tuning of optical properties (now useless...) :
    302                DHAZE_T(k,nw) = DHAZE_T(k,nw) * (FHIR * (1-SSA_T(k,nw)) + SSA_T(k,nw))
    303                SSA_T(k,nw)   = SSA_T(k,nw) / (FHIR * (1-SSA_T(k,nw)) + SSA_T(k,nw))
    304 
    305313               ! Diagnostics for clouds :
    306                POPT_CLOUDS(k,nw,1) = DHAZE_T(k,nw) ! dtau
    307                POPT_CLOUDS(k,nw,2) = SSA_T(k,nw)   ! wbar
    308                POPT_CLOUDS(k,nw,3) = ASF_T(k,nw)   ! gbar
     314               DIAG_OPTT(k,nw,1) = DHAZE_T(k,nw) ! dtau
     315               DIAG_OPTT(k,nw,2) = SSA_T(k,nw)   ! wbar
     316               DIAG_OPTT(k,nw,3) = ASF_T(k,nw)   ! gbar
    309317
    310318            ELSE
    311319               ! Diagnostics for clouds :
    312                POPT_CLOUDS(k,nw,1) = 0.D0 ! dtau
    313                POPT_CLOUDS(k,nw,2) = 1.0  ! wbar
    314                POPT_CLOUDS(k,nw,3) = 1.0  ! gbar
     320               DIAG_OPTT(k,nw,1) = 0.D0 ! dtau
     321               DIAG_OPTT(k,nw,2) = 1.0  ! wbar
     322               DIAG_OPTT(k,nw,3) = 1.0  ! gbar
    315323            ENDIF
    316324           
     
    321329            if (seashaze) DHAZE_T(k,nw) = DHAZE_T(k,nw)*seashazefact(k)
    322330            ! Diagnostics for the haze :
    323             POPT_HAZE(k,nw,1) = DHAZE_T(k,nw) ! dtau
    324             POPT_HAZE(k,nw,2) = SSA_T(k,nw)   ! wbar
    325             POPT_HAZE(k,nw,3) = ASF_T(k,nw)   ! gbar
     331            DIAG_OPTH(k,nw,1) = DHAZE_T(k,nw) ! dtau
     332            DIAG_OPTH(k,nw,2) = SSA_T(k,nw)   ! wbar
     333            DIAG_OPTH(k,nw,3) = ASF_T(k,nw)   ! gbar
    326334            ! Diagnostics for clouds :
    327             POPT_CLOUDS(k,nw,1) = 0.D0 ! dtau
    328             POPT_CLOUDS(k,nw,2) = 1.0  ! wbar
    329             POPT_CLOUDS(k,nw,3) = 1.0  ! gbar
     335            DIAG_OPTT(k,nw,1) = 0.D0 ! dtau
     336            DIAG_OPTT(k,nw,2) = 1.0  ! wbar
     337            DIAG_OPTT(k,nw,3) = 1.0  ! gbar
    330338         ENDIF ! ENDIF callmufi
    331339         
     
    436444                           + DCONT   & ! For parameterized continuum absorption
    437445                        + DHAZE_T(K,NW)     ! For Titan Haze
     446
     447         DIAG_OPTH(K,nw,4) = 0.d0
     448         DIAG_OPTH(K,nw,5) = TAUGAS
     449         DIAG_OPTH(K,nw,6) = DCONT
     450         DIAG_OPTT(K,nw,4) = 0.d0
     451         DIAG_OPTT(K,nw,5) = TAUGAS
     452         DIAG_OPTT(K,nw,6) = DCONT
    438453
    439454      end do ! k = L_LEVELS
  • trunk/LMDZ.TITAN/libf/phytitan/optcv.F90

    r3090 r3318  
    11SUBROUTINE OPTCV(PQMO,NLAY,ZLEV,PLEV,TMID,PMID,  &
    22     DTAUV,TAUV,TAUCUMV,WBARV,COSBV,TAURAY,TAUGSURF,SEASHAZEFACT,&
    3      POPT_HAZE,POPT_CLOUDS,CDCOLUMN)
     3     DIAG_OPTH,DIAG_OPTT,CDCOLUMN)
    44
    55  use radinc_h
     
    1111  use callkeys_mod, only: continuum,graybody,callgasvis,corrk_recombin,   &
    1212                          callclouds,callmufi,seashaze,uncoupl_optic_haze,&
    13                           opt4clouds,FHVIS, Fssa
     13                          opt4clouds,FHVIS,FCVIS
    1414  use tracer_h, only: nmicro,nice,ices_indx
    1515
     
    6969  REAL*8, INTENT(OUT) :: WBARV(L_NLAYRAD,L_NSPECTV,L_NGAUSS)
    7070  REAL*8, INTENT(OUT) :: TAUGSURF(L_NSPECTV,L_NGAUSS-1)
    71   REAL*8, INTENT(OUT) :: POPT_HAZE(L_LEVELS,L_NSPECTI,3)
    72   REAL*8, INTENT(OUT) :: POPT_CLOUDS(L_LEVELS,L_NSPECTI,3)
     71  REAL*8, INTENT(OUT) :: DIAG_OPTH(L_LEVELS,L_NSPECTV,6)
     72  REAL*8, INTENT(OUT) :: DIAG_OPTT(L_LEVELS,L_NSPECTV,6)
    7373  ! ==========================================================
    7474 
     
    124124  real*8,save :: rhoaer_s(L_NSPECTV),ssa_s(L_NSPECTV),asf_s(L_NSPECTV)
    125125  real*8,save :: rhoaer_f(L_NSPECTV),ssa_f(L_NSPECTV),asf_f(L_NSPECTV)
    126   real*8,save :: ssa_ccn(L_NSPECTI),asf_ccn(L_NSPECTI)
     126  real*8,save :: ssa_ccn(L_NSPECTV),asf_ccn(L_NSPECTV)
    127127  real*8,save :: ssa_cld(L_NSPECTV),asf_cld(L_NSPECTV)
    128128!$OMP THREADPRIVATE(rhoaer_s,rhoaer_f,ssa_s,ssa_f,ssa_cld,asf_s,asf_f,asf_cld)
     
    191191   end do
    192192 
     193   DIAG_OPTH(:,:,:) = 0.D0
     194   DIAG_OPTT(:,:,:) = 0.D0
     195
    193196   do NW=1,L_NSPECTV
    194197      ! We ignore K=1...
     
    254257            ENDIF
    255258
     259            ! Tuning of optical properties for haze :
     260            !dtauaer_s = dtauaer_s * (FHVIS * (1-ssa_s(nw)) + ssa_s(nw))
     261            !ssa_s(nw) = ssa_s(nw) / (FHVIS * (1-ssa_s(nw)) + ssa_s(nw))
     262            dtauaer_f = dtauaer_f * (FHVIS * (1-ssa_f(nw)) + ssa_f(nw))
     263            ssa_f(nw) = ssa_f(nw) / (FHVIS * (1-ssa_f(nw)) + ssa_f(nw))
     264
    256265            ! Total of Haze opacity (dtau), SSA (w) and ASF (COS) :
    257266            DHAZE_T(k,nw) = dtauaer_s + dtauaer_f
     
    265274               ASF_T(k,nw)   = 1.0
    266275            ENDIF
     276
    267277            ! Diagnostics for the haze :
    268             POPT_HAZE(k,nw,1) = DHAZE_T(k,nw) ! dtau
    269             POPT_HAZE(k,nw,2) = SSA_T(k,nw)   ! wbar
    270             POPT_HAZE(k,nw,3) = ASF_T(k,nw)   ! gbar
     278            DIAG_OPTH(k,nw,1) = DHAZE_T(k,nw) ! dtau
     279            DIAG_OPTH(k,nw,2) = SSA_T(k,nw)   ! wbar
     280            DIAG_OPTH(k,nw,3) = ASF_T(k,nw)   ! gbar
    271281
    272282            !---------------------
     
    277287               m0ccn = pqmo(ilay,5) / 2.0
    278288               m3ccn = pqmo(ilay,6) / 2.0
    279                m3cld = 0.0d0
     289               m3cld  = pqmo(ilay,6) / 2.0
    280290               
    281291               ! Clear / Dark column method :
     
    309319
    310320               ! For small dropplets, opacity of nucleus dominates
     321               dtau_cld = (dtau_cld*m3ccn + dtau_cld*m3cld) / (m3ccn + m3cld)
    311322               ssa_cld(nw) = (ssa_ccn(nw)*m3ccn + ssa_cld(nw)*m3cld) / (m3ccn + m3cld)
    312                ssa_cld(nw) = ssa_cld(nw) * Fssa
    313 
     323
     324               ! Tuning of optical properties for clouds :
     325               dtau_cld = dtau_cld * (FCVIS * (1-ssa_cld(nw)) + ssa_cld(nw))
     326               ssa_cld(nw) = ssa_cld(nw) / (FCVIS * (1-ssa_cld(nw)) + ssa_cld(nw))
     327               
    314328               ! Total of Haze + Clouds opacity (dtau), SSA (w) and ASF (COS) :
    315329               DHAZE_T(k,nw) = dtauaer_s + dtauaer_f + dtau_cld
     
    324338               ENDIF
    325339               
    326                ! Tuning of optical properties (now useless...) :
    327                DHAZE_T(k,nw) = DHAZE_T(k,nw) * (FHVIS * (1-SSA_T(k,nw)) + SSA_T(k,nw))
    328                SSA_T(k,nw)   = SSA_T(k,nw) / (FHVIS * (1-SSA_T(k,nw)) + SSA_T(k,nw))
    329 
    330340               ! Diagnostics for clouds :
    331                POPT_CLOUDS(k,nw,1) = DHAZE_T(k,nw) ! dtau
    332                POPT_CLOUDS(k,nw,2) = SSA_T(k,nw)   ! wbar
    333                POPT_CLOUDS(k,nw,3) = ASF_T(k,nw)   ! gbar
     341               DIAG_OPTT(k,nw,1) = DHAZE_T(k,nw) ! dtau
     342               DIAG_OPTT(k,nw,2) = SSA_T(k,nw)   ! wbar
     343               DIAG_OPTT(k,nw,3) = ASF_T(k,nw)   ! gbar
    334344           
    335345            ELSE
    336346               ! Diagnostics for clouds :
    337                POPT_CLOUDS(k,nw,1) = 0.D0 ! dtau
    338                POPT_CLOUDS(k,nw,2) = 1.0  ! wbar
    339                POPT_CLOUDS(k,nw,3) = 1.0  ! gbar
     347               DIAG_OPTT(k,nw,1) = 0.D0 ! dtau
     348               DIAG_OPTT(k,nw,2) = 1.0  ! wbar
     349               DIAG_OPTT(k,nw,3) = 1.0  ! gbar
    340350            ENDIF
    341351           
     
    346356            if (seashaze) DHAZE_T(k,nw) = DHAZE_T(k,nw)*seashazefact(k)
    347357            ! Diagnostics for the haze :
    348             POPT_HAZE(k,nw,1) = DHAZE_T(k,nw) ! dtau
    349             POPT_HAZE(k,nw,2) = SSA_T(k,nw)   ! wbar
    350             POPT_HAZE(k,nw,3) = ASF_T(k,nw)   ! gbar
     358            DIAG_OPTH(k,nw,1) = DHAZE_T(k,nw) ! dtau
     359            DIAG_OPTH(k,nw,2) = SSA_T(k,nw)   ! wbar
     360            DIAG_OPTH(k,nw,3) = ASF_T(k,nw)   ! gbar
    351361            ! Diagnostics for clouds :
    352             POPT_CLOUDS(k,nw,1) = 0.D0 ! dtau
    353             POPT_CLOUDS(k,nw,2) = 1.0  ! wbar
    354             POPT_CLOUDS(k,nw,3) = 1.0  ! gbar
     362            DIAG_OPTT(k,nw,1) = 0.D0 ! dtau
     363            DIAG_OPTT(k,nw,2) = 1.0  ! wbar
     364            DIAG_OPTT(k,nw,3) = 1.0  ! gbar
    355365         ENDIF ! ENDIF callmufi
    356366         
     
    472482         DTAUKV(K,nw,ng) = DRAYAER + DCONT ! Scattering + parameterized continuum absorption, including Titan's haze
    473483
    474       end do ! k = L_LEVELS
    475    end do ! nw = L_NSPECTI
     484         DIAG_OPTH(K,nw,4) = DRAYAER
     485         DIAG_OPTH(K,nw,5) = TAUGAS
     486         DIAG_OPTH(K,nw,6) = DCONT
     487         DIAG_OPTT(K,nw,4) = DRAYAER
     488         DIAG_OPTT(K,nw,5) = TAUGAS
     489         DIAG_OPTT(K,nw,6) = DCONT
     490
     491      end do ! K = L_LEVELS
     492   end do ! nw = L_NSPECTV
    476493
    477494  !=======================================================================
  • trunk/LMDZ.TITAN/libf/phytitan/phyetat0_mod.F90

    r1943 r3318  
    2020                     inquire_dimension, inquire_dimension_length
    2121  use callkeys_mod, only: surfalbedo,surfemis
     22  use geometry_mod, only: latitude_deg
    2223
    2324  implicit none
     
    219220             minval(emis), maxval(emis)
    220221
    221 ! Depth of methane tank (added by JVO 2017)
     222! Depth of methane tank (added by BdBdT 2023)
    222223 if (startphy_file) then
    223224   call get_field("tankCH4",tankCH4,found,indextime)
     
    225226     write(*,*) "phyetat0: Failed loading <tankCH4>"
    226227     !  call abort
    227      tankCH4(:)=2.
     228     DO ig = 1, ngrid
     229      if (REAL(latitude_deg(ig)) .ge. 70.0) then
     230        tankCH4(ig) = 200.0 ! [m]
     231      else if (REAL(latitude_deg(ig)) .le. -70.0) then
     232        tankCH4(ig) = 50.0  ! [m]
     233      else
     234        tankCH4(ig) = 0.0   ! [m]
     235      endif
     236     ENDDO
    228237   endif
    229238 else
    230    tankCH4(:)=2.
     239  DO ig = 1, ngrid
     240    if (REAL(latitude_deg(ig)) .ge. 70.0) then
     241      tankCH4(ig) = 200.0 ! [m]
     242    else if (REAL(latitude_deg(ig)) .le. -70.0) then
     243      tankCH4(ig) = 50.0  ! [m]
     244    else
     245      tankCH4(ig) = 0.0   ! [m]
     246    endif
     247   ENDDO
    231248 endif ! of if (startphy_file)
    232249 write(*,*) "phyetat0: Depth of methane tank <tankCH4> range:", &
  • trunk/LMDZ.TITAN/libf/phytitan/phys_state_var_mod.F90

    r3083 r3318  
    7272!$OMP THREADPRIVATE(int_dtaui,int_dtauv)
    7373
    74       real,dimension(:,:,:,:),allocatable,save :: zpopthi ! IR optical properties of haze within narrowbands for diags (dtau,tau,k,wbar,gbar).
    75       real,dimension(:,:,:,:),allocatable,save :: zpopthv ! VI optical properties of haze within narrowbands for diags (dtau,tau,k,wbar,gbar).
    76       real,dimension(:,:,:,:),allocatable,save :: zpoptci ! IR optical properties of clouds within narrowbands for diags (dtau,tau,k,wbar,gbar).
    77       real,dimension(:,:,:,:),allocatable,save :: zpoptcv ! VI optical properties of clouds within narrowbands for diags (dtau,tau,k,wbar,gbar).
    78 !$OMP THREADPRIVATE(zpopthi,zpopthv,zpoptci,zpoptcv)
     74      real,dimension(:,:,:,:),allocatable,save :: zpopthi ! IR optical properties [haze] within narrowbands for diags (dtau,tau,k,wbar,gbar,drayaer,taugaz,dcont).
     75      real,dimension(:,:,:,:),allocatable,save :: zpopthv ! VI optical properties [haze] within narrowbands for diags (dtau,tau,k,wbar,gbar,drayaer,taugaz,dcont).
     76      real,dimension(:,:,:,:),allocatable,save :: zpoptti ! IR optical properties [haze+clouds] within narrowbands for diags (dtau,tau,k,wbar,gbar,drayaer,taugaz,dcont).
     77      real,dimension(:,:,:,:),allocatable,save :: zpopttv ! VI optical properties [haze+clouds] within narrowbands for diags (dtau,tau,k,wbar,gbar,drayaer,taugaz,dcont).
     78!$OMP THREADPRIVATE(zpopthi,zpopthv,zpoptti,zpopttv)
    7979
    8080      real,dimension(:,:),allocatable,save :: u_ref ! Reference profile for the zonal wind nudging (m.s-1).
     
    139139        ALLOCATE(int_dtaui(klon,klev,L_NSPECTI))
    140140        ALLOCATE(int_dtauv(klon,klev,L_NSPECTV))
    141         ALLOCATE(zpopthi(klon,klev,L_NSPECTI,5))
    142         ALLOCATE(zpopthv(klon,klev,L_NSPECTV,5))
    143         ALLOCATE(zpoptci(klon,klev,L_NSPECTI,5))
    144         ALLOCATE(zpoptcv(klon,klev,L_NSPECTV,5))
     141        ALLOCATE(zpopthi(klon,klev,L_NSPECTI,8))
     142        ALLOCATE(zpopthv(klon,klev,L_NSPECTV,8))
     143        ALLOCATE(zpoptti(klon,klev,L_NSPECTI,8))
     144        ALLOCATE(zpopttv(klon,klev,L_NSPECTV,8))
    145145        ALLOCATE(u_ref(49,24))
    146146        allocate(dycchi(klon,klev,nkim))
     
    207207        DEALLOCATE(zpopthi)
    208208        DEALLOCATE(zpopthv)
    209         DEALLOCATE(zpoptci)
    210         DEALLOCATE(zpoptcv)
     209        DEALLOCATE(zpoptti)
     210        DEALLOCATE(zpopttv)
    211211        DEALLOCATE(u_ref)
    212212        DEALLOCATE(dycchi)
  • trunk/LMDZ.TITAN/libf/phytitan/physiq_mod.F90

    r3090 r3318  
    1515 
    1616      use radinc_h, only : L_NSPECTI,L_NSPECTV
    17       use radcommon_h, only: sigma, gzlat, grav, BWNV, WNOI, WNOV
     17      use radcommon_h, only: sigma, gzlat, grav, BWNV, WAVEI, WAVEV
    1818      use comchem_h, only: nkim, cnames, nlaykim_up, ykim_up, ykim_tot, botCH4
    1919      use comdiurn_h, only: coslat, sinlat, coslon, sinlon
     
    161161!                + clean of all too-generic (ocean, water, co2 ...) routines
    162162!                + Titan's chemistry
    163 !           Microphysical moment model - J.Burgalat / B. de Batz de Trenquelléon (2022-2023)
    164 !           Modified : B. de Batz de Trenquelléon (2023)
    165 !                + Optical improvements (haze and cloud)
    166 !                + Nudging of troposphere !!
     163!           Microphysical moment model: J.Burgalat / B. de Batz de Trenquelléon (2022-2023)
     164!           Optics for haze and clouds: B. de Batz de Trenquelléon (2023)
    167165!============================================================================================
    168166
     
    279277      real zdpsrfmr(ngrid)        ! Pressure tendency for mass_redistribution routine (Pa/s).
    280278
     279      ! Temporal tracers :
     280      real ppq(ngrid,nlayer,nq)   ! Tracers (kg/kg_of_air).
    281281     
    282282   
     
    422422         zpopthi(:,:,:,:) = 0.D0
    423423         zpopthv(:,:,:,:) = 0.D0
    424          zpoptci(:,:,:,:) = 0.D0
    425          zpoptcv(:,:,:,:) = 0.D0
     424         zpoptti(:,:,:,:) = 0.D0
     425         zpopttv(:,:,:,:) = 0.D0
    426426
    427427!        Initialize setup for correlated-k radiative transfer
     
    589589         endif
    590590
    591          ! Read SuperNudging.dat and initialize the nudging for pu
    592          ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     591         ! Read NewNudging.dat and initialize the nudging for pu
     592         ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    593593         if (nudging_u) then
    594             nudging_file=trim(datadir)//'/SuperNudging.dat'
     594            nudging_file=trim(datadir)//'/NewNudging.dat'
    595595            inquire(file=trim(nudging_file),exist=file_ok)
    596596            if(.not.file_ok) then
     
    620620         write(*,*) "physiq: call initialize_xios_output"
    621621         call initialize_xios_output(pday,ptime,ptimestep,daysec, &
    622                                      presnivs,pseudoalt,wnoi,wnov)
     622                                     presnivs,pseudoalt,wavei,wavev)
    623623#endif
    624624         write(*,*) "physiq: end of firstcall"
     
    840840
    841841               ! standard callcorrk
     842               !call callcorrk(ngrid,nlayer,ppq,nq,qsurf,zday,                      &
    842843               call callcorrk(ngrid,nlayer,pq,nq,qsurf,zday,                      &
    843844                              albedo,albedo_equivalent,emis,mu0,pplev,pplay,zzlev,&
     
    846847                              fluxsurfabs_sw,fluxtop_lw,                          &
    847848                              fluxabs_sw,fluxtop_dn,OLR_nu,OSR_nu,                &
    848                               int_dtaui,int_dtauv,zpopthi,zpopthv,zpoptci,zpoptcv,&
     849                              int_dtaui,int_dtauv,zpopthi,zpopthv,zpoptti,zpopttv,&
    849850                              lastcall)
    850851
     
    932933!  III. Vertical diffusion (turbulent mixing) :
    933934!  --------------------------------------------
    934 
    935935      if (calldifv) then
    936936     
     
    12331233            ! iv. Surface pseudo-evaporation
    12341234            ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    1235             ! Infinite tank of CH4 (global)
     1235            ! Infinite tank of CH4
    12361236            if (resCH4_inf) then
    12371237               do ig=1,ngrid
     
    12431243               enddo
    12441244
    1245             ! Tank at the pole (for the moment infinit...)
    12461245            else
    1247                tankCH4(:) = 2.
     1246               ! Fill lakes with precipitation :
     1247               tankCH4(:) = tankCH4(:) + (mmd_ice_prec(:,1) / 425. * ptimestep) ! [m]
     1248               
     1249               ! Evaporation of lakes :
    12481250               if (moyzon_ch) then
    12491251                 tpq_CH4(:,:) = ychimbar(:,:,7) + dycchi(:,:,7)*ptimestep ! + dycchi because ychim not yet updated [mol/mol]
     
    12561258            endif
    12571259
    1258             pdq(:,1,7+nmicro) = pdq(:,1,7+nmicro) + dycevapCH4(:)*rat_mmol(7+nmicro)
     1260
     1261            pdq(:,1,7+nmicro) = pdq(:,1,7+nmicro) + dycevapCH4(:)*rat_mmol(7+nmicro) ! convert tendencies to mass mixing ratio
    12591262           
    12601263            ! v. Updates and positivity check
     
    13911394      enddo
    13921395
    1393       ! Tracers : for now because problem with re-evaporation in microphysics ...
    1394       if (callclouds) then
    1395          pdq(:,1:13,gazs_indx(2)) = (2e-8 - pq(:,1:13,gazs_indx(2))) / (ptimestep * 100)  ! C2H6
    1396          pdq(:,1:12,gazs_indx(3)) = (1e-9 - pq(:,1:12,gazs_indx(3))) / (ptimestep * 100)  ! C2H2
    1397          pdq(:,1:11,gazs_indx(4)) = (1e-20 - pq(:,1:11,gazs_indx(4))) / (ptimestep * 100) ! HCN
    1398       endif
     1396      ! [BBT TEMPO] Firstcall : Initialize Tracers
     1397      !if (firstcall) then
     1398      !   do ig = 1, ngrid
     1399      !      ! Aerosols :
     1400      !      !-----------
     1401      !      WHERE (pq(ig,:24,3) < 3.0e7) &
     1402      !               pdq(ig,:24,3) = (3.0e7 - pq(ig,:24,3)) / ptimestep   ! M0 aerf
     1403      !      WHERE (pq(ig,:24,4) < 1.5e-12) &
     1404      !               pdq(ig,:24,4) = (1.5e-12 - pq(ig,:24,4)) / ptimestep ! M3 aerf
     1405      !      ! Lakes :
     1406      !      !--------
     1407      !      if (REAL(latitude_deg(ig)) .ge. 70.0) then
     1408      !         tankCH4(ig) = 100.0 ! [m]
     1409      !      else if (REAL(latitude_deg(ig)) .le. -70.0) then
     1410      !         tankCH4(ig) = 20.0  ! [m]
     1411      !      else
     1412      !         tankCH4(ig) = 0.0   ! [m]
     1413      !      endif
     1414      !      ! Species :
     1415      !      !----------
     1416      !      do iq = 1, size(ices_indx)
     1417      !         ! CH4 :
     1418      !         !------
     1419      !         if(trim(nameOfTracer(gazs_indx(iq))) .eq. "CH4") then
     1420      !            WHERE (pq(ig,:,gazs_indx(iq))/rat_mmol(gazs_indx(iq)) > 0.05) &
     1421      !               pdq(ig,:,gazs_indx(iq)) = (0.05 * rat_mmol(gazs_indx(iq)) - pq(ig,:,gazs_indx(iq))) / ptimestep
     1422      !            WHERE (pq(ig,:,gazs_indx(iq))/rat_mmol(gazs_indx(iq)) < 0.014) &
     1423      !               pdq(ig,:,gazs_indx(iq)) = (0.014 * rat_mmol(gazs_indx(iq)) - pq(ig,:,gazs_indx(iq))) / ptimestep
     1424      !         endif
     1425      !         ! C2H2 :
     1426      !         !-------
     1427      !         if(trim(nameOfTracer(gazs_indx(iq))) .eq. "C2H2") then
     1428      !            WHERE (pq(ig,26:,gazs_indx(iq))/rat_mmol(gazs_indx(iq)) < 3.0e-6) &
     1429      !               pdq(ig,26:,gazs_indx(iq)) = (5.0e-6 * rat_mmol(gazs_indx(iq)) - pq(ig,26:,gazs_indx(iq))) / ptimestep
     1430      !         endif
     1431      !         ! C2H6 :
     1432      !         !-------
     1433      !         if(trim(nameOfTracer(gazs_indx(iq))) .eq. "C2H6") then
     1434      !            WHERE (pq(ig,26:,gazs_indx(iq))/rat_mmol(gazs_indx(iq)) < 2.0e-5) &
     1435      !               pdq(ig,26:,gazs_indx(iq)) = (2.0e-5 * rat_mmol(gazs_indx(iq)) - pq(ig,26:,gazs_indx(iq))) / ptimestep
     1436      !         endif
     1437      !         ! HCN :
     1438      !         !------
     1439      !         if(trim(nameOfTracer(gazs_indx(iq))) .eq. "HCN") then
     1440      !            WHERE (pq(ig,26:,gazs_indx(iq))/rat_mmol(gazs_indx(iq)) < 8.0e-6) &
     1441      !               pdq(ig,26:,gazs_indx(iq)) = (8.0e-6 * rat_mmol(gazs_indx(iq)) - pq(ig,26:,gazs_indx(iq))) / ptimestep
     1442      !         endif
     1443      !      enddo
     1444      !   enddo
     1445      !endif
     1446      ! [BdBdT : Forcage de la photochimie]
     1447      do ig = 1, ngrid
     1448         do iq = 1, size(ices_indx)
     1449            ! C2H2 :
     1450            !-------
     1451            if(trim(nameOfTracer(gazs_indx(iq))) .eq. "C2H2") then
     1452               pdq(ig,nlayer-3:,gazs_indx(iq)) = (4.0e-5 * rat_mmol(gazs_indx(iq)) - pq(ig,nlayer-3:,gazs_indx(iq))) / ptimestep
     1453            endif
     1454            ! C2H6 :
     1455            !-------
     1456            if(trim(nameOfTracer(gazs_indx(iq))) .eq. "C2H6") then
     1457               pdq(ig,nlayer-3:,gazs_indx(iq)) = (1.0e-4 * rat_mmol(gazs_indx(iq)) - pq(ig,nlayer-3:,gazs_indx(iq))) / ptimestep
     1458            endif
     1459            ! HCN :
     1460            !------
     1461            if(trim(nameOfTracer(gazs_indx(iq))) .eq. "HCN") then
     1462               pdq(ig,nlayer-3:,gazs_indx(iq)) = (4.0e-5 * rat_mmol(gazs_indx(iq)) - pq(ig,nlayer-3:,gazs_indx(iq))) / ptimestep
     1463            endif
     1464         enddo
     1465      enddo
    13991466      zq(:,:,:) = pq(:,:,:) + pdq(:,:,:)*ptimestep
    14001467
     
    16841751      ENDIF
    16851752
    1686       ! Subsurface temperatures (3D)
    1687       !CALL send_xios_field("tempsoil",tsoil)
    1688       ! Total energy balance diagnostics (2D)
    1689       !CALL send_xios_field("ALB",albedo_equivalent)
    1690 
    16911753      !--------------------------------------------------------
    16921754      ! Winds trends :
     
    17471809      ! Optical diagnostics :
    17481810      !--------------------------------------------------------
    1749       ! Diagnostics for haze :
    1750       ! IR
     1811      ! Haze opacity :
     1812      CALL send_xios_field('ttauhv_14',zpopthv(:,:,14,2)) ! 14 --> 1.076 um
     1813      CALL send_xios_field('ttauhv_20',zpopthv(:,:,20,2)) ! 20 --> 0.671 um
     1814      CALL send_xios_field('ttauhv_23',zpopthv(:,:,23,2)) ! 23 --> 0.346 um
     1815      CALL send_xios_field('ttauhi_02',zpopthi(:,:,2,2))  ! 02 --> 175.3 um
     1816      CALL send_xios_field('ttauhi_17',zpopthi(:,:,17,2)) ! 17 --> 11.00 um
     1817      CALL send_xios_field('ttauhi_23',zpopthi(:,:,23,2)) ! 23 --> 4.849 um
     1818      ! Haze extinction :
     1819      CALL send_xios_field('kkhv_14',zpopthv(:,:,14,3))
     1820      CALL send_xios_field('kkhv_20',zpopthv(:,:,20,3))
     1821      CALL send_xios_field('kkhv_23',zpopthv(:,:,23,3))
     1822      CALL send_xios_field('kkhi_02',zpopthi(:,:,2,3))
     1823      CALL send_xios_field('kkhi_17',zpopthi(:,:,17,3))
     1824      CALL send_xios_field('kkhi_23',zpopthi(:,:,23,3))
     1825      ! Haze single scattering albedo :
     1826      CALL send_xios_field('wwhv_14',zpopthv(:,:,14,4))
     1827      CALL send_xios_field('wwhv_20',zpopthv(:,:,20,4))
     1828      CALL send_xios_field('wwhv_23',zpopthv(:,:,23,4))
     1829      CALL send_xios_field('wwhi_02',zpopthi(:,:,2,4))
     1830      CALL send_xios_field('wwhi_17',zpopthi(:,:,17,4))
     1831      CALL send_xios_field('wwhi_23',zpopthi(:,:,23,4))
     1832      ! Haze asymmetry parameter :
     1833      CALL send_xios_field('gghv_14',zpopthv(:,:,14,5))
     1834      CALL send_xios_field('gghv_20',zpopthv(:,:,20,5))
     1835      CALL send_xios_field('gghv_23',zpopthv(:,:,23,5))
     1836      CALL send_xios_field('gghi_02',zpopthi(:,:,2,5))
     1837      CALL send_xios_field('gghi_17',zpopthi(:,:,17,5))
     1838      CALL send_xios_field('gghi_23',zpopthi(:,:,23,5))
     1839
     1840      ! Diagnostics for haze and clouds :
     1841      IF (callclouds) THEN
     1842         ! Opacity :
     1843         CALL send_xios_field('ttauv_14',zpopttv(:,:,14,2)) ! 14 --> 1.076 um
     1844         CALL send_xios_field('ttauv_20',zpopttv(:,:,20,2)) ! 20 --> 0.671 um
     1845         CALL send_xios_field('ttauv_23',zpopttv(:,:,23,2)) ! 23 --> 0.346 um
     1846         CALL send_xios_field('ttaui_02',zpoptti(:,:,2,2))  ! 02 --> 175.3 um
     1847         CALL send_xios_field('ttaui_17',zpoptti(:,:,17,2)) ! 17 --> 11.00 um
     1848         CALL send_xios_field('ttaui_23',zpoptti(:,:,23,2)) ! 23 --> 4.849 um
     1849         ! Extinction :
     1850         CALL send_xios_field('kkv_14',zpopttv(:,:,14,3))
     1851         CALL send_xios_field('kkv_20',zpopttv(:,:,20,3))
     1852         CALL send_xios_field('kkv_23',zpopttv(:,:,23,3))
     1853         CALL send_xios_field('kki_02',zpoptti(:,:,2,3))
     1854         CALL send_xios_field('kki_17',zpoptti(:,:,17,3))
     1855         CALL send_xios_field('kki_23',zpoptti(:,:,23,3))
     1856         ! Single scattering albedo :
     1857         CALL send_xios_field('wwv_14',zpopttv(:,:,14,4))
     1858         CALL send_xios_field('wwv_20',zpopttv(:,:,20,4))
     1859         CALL send_xios_field('wwv_23',zpopttv(:,:,23,4))
     1860         CALL send_xios_field('wwi_02',zpoptti(:,:,2,4))
     1861         CALL send_xios_field('wwi_17',zpoptti(:,:,17,4))
     1862         CALL send_xios_field('wwi_23',zpoptti(:,:,23,4))
     1863         ! Asymmetry parameter :
     1864         CALL send_xios_field('ggv_14',zpopttv(:,:,14,5))
     1865         CALL send_xios_field('ggv_20',zpopttv(:,:,20,5))
     1866         CALL send_xios_field('ggv_23',zpopttv(:,:,23,5))
     1867         CALL send_xios_field('ggi_02',zpoptti(:,:,2,5))
     1868         CALL send_xios_field('ggi_17',zpoptti(:,:,17,5))
     1869         CALL send_xios_field('ggi_23',zpoptti(:,:,23,5))
     1870         ! DRAYAER, TAUGAS, DCONT :
     1871         CALL send_xios_field('drayaerv_20',zpopttv(:,:,20,6)) ! 20 --> 0.671um
     1872         CALL send_xios_field('taugasv_20',zpopttv(:,:,20,7))
     1873         CALL send_xios_field('dcontv_20',zpopttv(:,:,20,8))
     1874         CALL send_xios_field('drayaeri_17',zpoptti(:,:,17,6)) ! 17 --> 11.00um
     1875         CALL send_xios_field('taugasi_17',zpoptti(:,:,17,7))
     1876         CALL send_xios_field('dconti_17',zpoptti(:,:,17,8))
     1877      ENDIF     
     1878
     1879      ! Diagnostics for haze and clouds (4D) :
    17511880      CALL send_xios_field('dtauhi',zpopthi(:,:,:,1))
    17521881      CALL send_xios_field('tauhi',zpopthi(:,:,:,2))
     
    17541883      CALL send_xios_field('whi',zpopthi(:,:,:,4))
    17551884      CALL send_xios_field('ghi',zpopthi(:,:,:,5))
    1756       ! VI
    17571885      CALL send_xios_field('dtauhv',zpopthv(:,:,:,1))
    17581886      CALL send_xios_field('tauhv',zpopthv(:,:,:,2))
     
    17601888      CALL send_xios_field('whv',zpopthv(:,:,:,4))
    17611889      CALL send_xios_field('ghv',zpopthv(:,:,:,5))
    1762 
    1763       ! Diagnostics for haze and clouds :
    17641890      IF (callclouds) THEN
    1765          ! IR
    1766          CALL send_xios_field('dtaui',zpoptci(:,:,:,1))
    1767          CALL send_xios_field('taui',zpoptci(:,:,:,2))
    1768          CALL send_xios_field('ki',zpoptci(:,:,:,3))
    1769          CALL send_xios_field('wi',zpoptci(:,:,:,4))
    1770          CALL send_xios_field('gi',zpoptci(:,:,:,5))
    1771          ! VI
    1772          CALL send_xios_field('dtauv',zpoptcv(:,:,:,1))
    1773          CALL send_xios_field('tauv',zpoptcv(:,:,:,2))
    1774          CALL send_xios_field('kv',zpoptcv(:,:,:,3))
    1775          CALL send_xios_field('wv',zpoptcv(:,:,:,4))
    1776          CALL send_xios_field('gv',zpoptcv(:,:,:,5))
    1777       ENDIF     
     1891         CALL send_xios_field('dtaui',zpoptti(:,:,:,1))
     1892         CALL send_xios_field('taui',zpoptti(:,:,:,2))
     1893         CALL send_xios_field('ki',zpoptti(:,:,:,3))
     1894         CALL send_xios_field('wi',zpoptti(:,:,:,4))
     1895         CALL send_xios_field('gi',zpoptti(:,:,:,5))
     1896         CALL send_xios_field('dtauv',zpopttv(:,:,:,1))
     1897         CALL send_xios_field('tauv',zpopttv(:,:,:,2))
     1898         CALL send_xios_field('kv',zpopttv(:,:,:,3))
     1899         CALL send_xios_field('wv',zpopttv(:,:,:,4))
     1900         CALL send_xios_field('gv',zpopttv(:,:,:,5))
     1901      ENDIF
    17781902
    17791903      !--------------------------------------------------------
     
    18261950      !--------------------------------------------------------
    18271951      IF (callchim) THEN
     1952         ! Surface (2D) :
     1953         CALL send_xios_field("evapCH4",dycevapCH4(:)) ! Pseudo-evaporation flux (mol/mol/s)
     1954         CALL send_xios_field("tankCH4",tankCH4(:))    ! CH4 tank at the surface (m)
     1955
    18281956         ! Atmosphere (3D) :
    18291957         ! Chemical species :
     
    18772005            CALL send_xios_field(trim(cnames(iq))//"_tot",ykim_tot(iq,:,:)) ! mol/mol
    18782006         ENDDO
    1879 
    1880          ! Surface (2D) :
    1881          CALL send_xios_field("evapCH4",dycevapCH4(:)) ! Pseudo-evaporation flux (mol/mol/s)
    1882          CALL send_xios_field("tankCH4",tankCH4(:))    ! CH4 tank at the surface (m)
    18832007      ENDIF ! of 'if callchim'
    18842008
  • trunk/LMDZ.TITAN/libf/phytitan/xios_output_mod.F90

    r3090 r3318  
    1818
    1919  SUBROUTINE initialize_xios_output(day,timeofday,dtphys,daysec,&
    20                                     presnivs,pseudoalt,wnoi,wnov)
     20                                    presnivs,pseudoalt,wavei,wavev)
    2121!  USE mod_phys_lmdz_para, only: gather, bcast, &
    2222!                                jj_nb, jj_begin, jj_end, ii_begin, ii_end, &
     
    4040  IMPLICIT NONE
    4141 
    42   REAL,INTENT(IN) :: day ! Number of elapsed sols since reference Ls=0.
    43   REAL,INTENT(IN) :: timeofday ! "Universal time", given as fraction of sol (e.g.: 0.5 for noon).
    44   REAL,INTENT(IN) :: dtphys ! physics time step (s)
    45   REAL,INTENT(IN) :: daysec ! lengthof a standard day (s)
    46   REAL,INTENT(IN) :: presnivs(:) ! vertical grid approximate pressure (Pa)
     42  REAL,INTENT(IN) :: day          ! Number of elapsed sols since reference Ls=0.
     43  REAL,INTENT(IN) :: timeofday    ! "Universal time", given as fraction of sol (e.g.: 0.5 for noon).
     44  REAL,INTENT(IN) :: dtphys       ! physics time step (s)
     45  REAL,INTENT(IN) :: daysec       ! lengthof a standard day (s)
     46  REAL,INTENT(IN) :: presnivs(:)  ! vertical grid approximate pressure (Pa)
    4747  REAL,INTENT(IN) :: pseudoalt(:) ! vertical grid approximate altitude (km)
    48   REAL,INTENT(IN) :: wnoi(:) ! Array of wavelength channels for the infrared.
    49   REAL,INTENT(IN) :: wnov(:) ! Array of wavelength channels for the visible.
     48  REAL,INTENT(IN) :: wavei(:)     ! Array of the wavelenght (in microns) at the center of each IR spectral interval
     49  REAL,INTENT(IN) :: wavev(:)     ! Array of the wavelenght (in microns) at the center of each VI spectral interval
    5050 
    5151  INTEGER :: i
     
    7474                            unit="Pa",positive="down")
    7575   
    76     IF (prt_level>=10) WRITE(lunout,*) "initialize_xios_output: call xios_set_axis_attr for IR_Wavelength"
    77     CALL xios_set_axis_attr("IR_Wavelength", n_glo=size(wnoi), value=wnoi,unit="m",positive="up")
    78 
    79     IF (prt_level>=10) WRITE(lunout,*) "initialize_xios_output: call xios_set_axis_attr for VI_Wavelength"
    80     CALL xios_set_axis_attr("VI_Wavelength", n_glo=size(wnov), value=wnov,unit="m",positive="up")
     76    IF (prt_level>=10) WRITE(lunout,*) "initialize_xios_output: call xios_set_axis_attr for wavelength_ir"
     77    CALL xios_set_axis_attr("wavelength_ir", n_glo=size(wavei), value=wavei,unit="micrometer",positive="down")
     78
     79    IF (prt_level>=10) WRITE(lunout,*) "initialize_xios_output: call xios_set_axis_attr for wavelength_vi"
     80    CALL xios_set_axis_attr("wavelength_vi", n_glo=size(wavev), value=wavev,unit="micrometer",positive="down")
    8181                 
    8282    ! Calculation of pseudo-altitudes is still to be done
Note: See TracChangeset for help on using the changeset viewer.