source: trunk/LMDZ.PLUTO/libf/phypluto/cooling_hcn_c2h2.F90 @ 3586

      subroutine cooling_hcn_c2h2(ngrid,nlayer,pplay,pt,dtlw)

      use callkeys_mod, only: alpha_top, deltap, pref 
      
          
      implicit none

!==================================================================
!     Purpose
!     -------
!     Calculation of cooling rate for C2H2-HCN
!     = f(pplay) * B(lambda,T)
!
!     Inputs
!     ------ 
!     ngrid                 Number of vertical columns
!     nlayer                Number of layers
!     pt
!
!     Outputs
!     -------
!
!     dtlw       ! cooling rate
!
!     Authors
!     ------- 
!     Tanguy Bertrand (2016)
!     FF (2016)
!==================================================================
!     Arguments
      INTEGER ngrid, nlayer
      REAL pplay(ngrid,nlayer)            ! pres. level in GCM mid of layer
      REAL pt(ngrid,nlayer)
      REAL dtlw(ngrid,nlayer)
!     Local variables
      INTEGER l,ig
      REAL lonw
      REAL alpha !, alpha_top
      REAL transition
      REAL BB 
!-----------------------------------------------------------------------

      lonw = 14.e-6  ! 14um
      ! alpha_top=5.e-11  ! 1.e-13  ! cooling constant at top of atmosphere
      ! pref = 0.02      ! pressure at mid transition fo alpha_top (Pa)
      ! deltap = 0.1    ! width of transition to alpha_top (Pa)  
      !lonw = 11.e-6  ! 14um
      
!     transition = 0 if p>pref+deltap/2 and 1 if p< pref-deltap/2
      DO l = 1, nlayer
         DO ig = 1, ngrid
          transition  =  0.5*(1-tanh((pplay(ig,l)-pref)/deltap))
          dtlw(ig,l)=-transition*alpha_top*BB(lonw,pt(ig,l))
         ENDDO
      ENDDO   
      end 

!******************************************************
!   FUNCTION Blackbody  (Planck)
!**********************************************************
      function BB (lw, T)
! Variable declaration
! --------------------
! wavelenght (m), Temperature (K)
      real lw,T
! constant
      real c1,c2
      parameter ( c1=1.19103E-16 )
      parameter (c2=1.43887E-2 )
! function
!---------
      BB= (c1/lw**5)/(-1.+exp(c2/(lw*T)))
      return
      end