subroutine da_tpq_to_slp ( t, q, p, terr, psfc, slp) !----------------------------------------------------------------------- ! purpose: computes sea level pressure from the rule ! t1/t2=(p1/p2)**(gamma*r/g). ! ! input t temperature ! q mixing ratio ! p pressure ! terr terrain ! psfc surface pressure ! ! output slp sea level pressure !----------------------------------------------------------------------- implicit none real, intent(in) :: terr, psfc real, intent(in) :: t(kms:kme) real, intent(in) :: q(kms:kme) real, intent(in) :: p(kms:kme) real, intent(inout) :: slp integer :: k, klo, khi real :: pl, t0, ts, xterm,tlo, thi, tl real, parameter :: gamma = 6.5e-3 real, parameter :: tc = t_kelvin+17.5 real, parameter :: pconst = 10000.0 real, parameter :: eps = 0.622 if (trace_use) call da_trace_entry("da_tpq_to_slp") ! sea level pressure xterm=gamma* gas_constant / gravity ! compute pressure at pconst mb above surface (pl) if (terr <= 0.0) then slp = psfc return end if pl = psfc - pconst klo = 0 ! find 2 levels on sigma surfaces surrounding pl at each i,j do k=kts, kte-1 if ((p(k) >= pl) .and. (p(k+1) < pl)) then khi = k+1 klo = k exit end if end do if (klo < 1) then write(unit=message(1),fmt='(a,f11.3,a)') & 'error finding pressure level ',pconst,' mb above the surface' write(unit=message(2),fmt='(a,f11.3,2x,a,f11.3)') 'pl=',pl,' psfc=',psfc call da_error(__FILE__,__LINE__,message(1:2)) end if ! get temperature at pl (tl), extrapolate t at surface (ts) ! and t at sea level (t0) with 6.5 k/km lapse rate tlo=t(klo) * (eps+q(klo))/(eps*(1.0+q(klo))) thi=t(khi) * (eps+q(khi))/(eps*(1.0+q(khi))) tl=thi-(thi-tlo)*log(pl/p(khi)) & /log(p(klo)/p(khi)) ts=tl*(psfc/pl)**xterm t0=ts +gamma*terr ! correct sea level temperature if too hot if ( t0 >= tc ) then if ( ts <= tc ) then t0 = tc else t0 = tc-0.005*(ts-tc)**2 end if end if ! compute sea level pressure slp=psfc*exp(2.0*gravity*terr/(gas_constant*(ts+t0))) if (trace_use) call da_trace_exit("da_tpq_to_slp") end subroutine da_tpq_to_slp