Line data    Source code

       1             : 
       2       58824 : subroutine cpslec (ncol, pmid, phis, ps, t, psl, gravit, rair)
       3             : 
       4             : !----------------------------------------------------------------------- 
       5             : ! 
       6             : ! Purpose: 
       7             : ! Hybrid coord version:  Compute sea level pressure for a latitude line
       8             : ! 
       9             : ! Method: 
      10             : ! CCM2 hybrid coord version using ECMWF formulation
      11             : ! Algorithm: See section 3.1.b in NCAR NT-396 "Vertical 
      12             : ! Interpolation and Truncation of Model-Coordinate Data
      13             : !
      14             : ! Author: Stolen from the Processor by Erik Kluzek
      15             : ! 
      16             : !-----------------------------------------------------------------------
      17             : !
      18             : ! $Id$
      19             : ! $Author$
      20             : !
      21             : !-----------------------------------------------------------------------
      22             : 
      23             :   use shr_kind_mod, only: r8 => shr_kind_r8
      24             :   use ppgrid, only: pcols, pver
      25             : 
      26             :   implicit none
      27             : 
      28             : !-----------------------------Arguments---------------------------------
      29             :   integer , intent(in) :: ncol             ! longitude dimension
      30             : 
      31             :   real(r8), intent(in) :: pmid(pcols,pver) ! Atmospheric pressure (pascals)
      32             :   real(r8), intent(in) :: phis(pcols)      ! Surface geopotential (m**2/sec**2)
      33             :   real(r8), intent(in) :: ps(pcols)        ! Surface pressure (pascals)
      34             :   real(r8), intent(in) :: T(pcols,pver)    ! Vertical slice of temperature (top to bot)
      35             :   real(r8), intent(in) :: gravit           ! Gravitational acceleration
      36             :   real(r8), intent(in) :: rair             ! gas constant for dry air
      37             : 
      38             :   real(r8), intent(out):: psl(pcols)       ! Sea level pressures (pascals)
      39             : !-----------------------------------------------------------------------
      40             : 
      41             : !-----------------------------Parameters--------------------------------
      42             :   real(r8), parameter :: xlapse = 6.5e-3_r8   ! Temperature lapse rate (K/m)
      43             : !-----------------------------------------------------------------------
      44             : 
      45             : !-----------------------------Local Variables---------------------------
      46             :   integer i              ! Loop index
      47             :   real(r8) alpha         ! Temperature lapse rate in terms of pressure ratio (unitless)
      48             :   real(r8) Tstar         ! Computed surface temperature
      49             :   real(r8) TT0           ! Computed temperature at sea-level
      50             :   real(r8) alph          ! Power to raise P/Ps to get rate of increase of T with pressure
      51             :   real(r8) beta          ! alpha*phis/(R*T) term used in approximation of PSL
      52             : !-----------------------------------------------------------------------
      53             : !
      54       58824 :   alpha = rair*xlapse/gravit
      55      982224 :   do i=1,ncol
      56      982224 :      if ( abs(phis(i)/gravit) < 1.e-4_r8 )then
      57      555427 :         psl(i)=ps(i)
      58             :      else
      59      367973 :         Tstar=T(i,pver)*(1._r8+alpha*(ps(i)/pmid(i,pver)-1._r8)) ! pg 7 eq 5
      60             : 
      61      367973 :         TT0=Tstar + xlapse*phis(i)/gravit                  ! pg 8 eq 13
      62             : 
      63      367973 :         if ( Tstar<=290.5_r8 .and. TT0>290.5_r8 ) then           ! pg 8 eq 14.1
      64       17550 :            alph=rair/phis(i)*(290.5_r8-Tstar)  
      65      350423 :         else if (Tstar>290.5_r8  .and. TT0>290.5_r8) then        ! pg 8 eq 14.2
      66      139463 :            alph=0._r8
      67      139463 :            Tstar= 0.5_r8 * (290.5_r8 + Tstar)  
      68             :         else  
      69      210960 :            alph=alpha  
      70      210960 :            if (Tstar<255._r8) then  
      71       56447 :               Tstar= 0.5_r8 * (255._r8 + Tstar)                  ! pg 8 eq 14.3
      72             :            endif
      73             :         endif
      74             : 
      75      367973 :         beta = phis(i)/(rair*Tstar)
      76      367973 :         psl(i)=ps(i)*exp( beta*(1._r8-alph*beta/2._r8+((alph*beta)**2)/3._r8))
      77             :      end if
      78             :   enddo
      79             : 
      80       58824 :   return
      81             : end subroutine cpslec