da_tpq_to_slp.inc

1 subroutine da_tpq_to_slp ( T, Q, P, TERR, PSFC, SLP, xp )
2 
3    !-----------------------------------------------------------------------
4    ! Purpose:  computes sea level pressure from the rule                
5    !              t1/t2=(p1/p2)**(gamma*r/g).                              
6    !                                                                       
7    !     input       t        temperature
8    !                 q        mixing ratio
9    !                 p        pressure
10    !                 terr     terrain
11    !                 psfc     surface pressure
12    !                                                                       
13    !     output      slp      sea level pressure    
14    !-----------------------------------------------------------------------        
15 
16    implicit none
17 
18    type (xpose_type), intent(in)                 :: xp
19    real,                intent(in)               :: TERR, PSFC
20    real, dimension(xp%kms:xp%kme), intent(in)    :: T, Q, P
21    real,                intent(inout)            :: SLP
22 
23    integer              :: K, KLO, KHI
24    real                 :: PL, T0, TS, XTERM,    &
25                            TLO, THI, TL
26                                                
27    real, parameter      :: GAMMA = 6.5E-3,  &
28                            TC=t_kelvin+17.5,  &
29                            PCONST=10000. ,  &
30                            EPS   = 0.622
31                                                                        
32    ! SEA LEVEL PRESSURE                                            
33                                                                          
34    XTERM=GAMMA* gas_constant / gravity                                                   
35                                                                        
36    ! COMPUTE PRESSURE AT PCONST MB ABOVE SURFACE (PL)              
37                                                                         
38    if (terr <= 0.) then
39       slp = psfc
40       return
41    end if
42 
43    PL  = psfc - PCONST                                        
44    klo = 0
45 
46    ! FinD 2 LEVELS ON SIGMA SURFACES SURROUNDinG PL AT EACH I,J    
47 
48    do k=xp%kts, xp%kte-1
49       if ((p(k) >= pl) .and. (p(k+1) < pl)) then
50          khi = k+1
51          klo = k
52          exit
53       end if
54    end do
55 
56    if (klo < 1) then                                      
57       write(unit=message(1),fmt='(A,F11.3,A)') &
58          'ERROR FinDinG PRESSURE LEVEL ',PCONST,' MB ABOVE THE SURFACE'
59       write(unit=message(2),fmt='(A,F11.3,2X,A,F11.3)') 'PL=',PL,'  PSFC=',psfc
60       call da_error(__FILE__,__LINE__,message(1:2))                                
61    end if                                                         
62 
63    ! GET TEMPERATURE AT PL (TL), EXTRAPOLATE T AT SURFACE (TS)     
64    ! AND T AT SEA LEVEL (T0) WITH 6.5 K/KM LAPSE RATE              
65 
66    TLO=t(KLO) * (EPS+q(KLO))/(EPS*(1.+q(KLO)))
67    THI=t(KHI) * (EPS+q(KHI))/(EPS*(1.+q(KHI)))
68    TL=THI-(THI-TLO)*LOG(PL/p(KHI)) &
69                       /LOG(p(KLO)/p(KHI))               
70    TS=TL*(psfc/PL)**XTERM                           
71    T0=TS +GAMMA*terr
72 
73    ! CORRECT SEA LEVEL TEMPERATURE if TOO HOT                      
74 
75    if ( t0 >= tc ) then
76       if ( ts <= tc ) then
77         t0 = tc
78       else
79         t0 = tc-0.005*(ts-tc)**2
80       end if
81    end if
82 
83    ! COMPUTE SEA LEVEL PRESSURE                                    
84 
85    slp=psfc*EXP(2.*gravity*terr/(gas_constant*(TS+T0)))
86 
87 end subroutine da_tpq_to_slp
88 
89