1 |
|
|
MODULE tropopause_m |
2 |
|
|
|
3 |
|
|
IMPLICIT NONE |
4 |
|
|
PRIVATE |
5 |
|
|
PUBLIC :: dyn_tropopause |
6 |
|
|
|
7 |
|
|
CONTAINS |
8 |
|
|
|
9 |
|
|
!------------------------------------------------------------------------------- |
10 |
|
|
! |
11 |
|
|
FUNCTION dyn_tropopause(t, ts, paprs, pplay, rot, itrop, thet0, pvor0) |
12 |
|
|
! |
13 |
|
|
!------------------------------------------------------------------------------- |
14 |
|
|
USE assert_m, ONLY: assert |
15 |
|
|
USE assert_eq_m, ONLY: assert_eq |
16 |
|
|
USE dimphy, ONLY: klon, klev |
17 |
|
|
USE geometry_mod, ONLY: latitude_deg, longitude_deg |
18 |
|
|
USE vertical_layers_mod, ONLY: aps, bps, preff |
19 |
|
|
#ifdef REPROBUS |
20 |
|
|
USE chem_rep, ONLY: itroprep |
21 |
|
|
#endif |
22 |
|
|
|
23 |
|
|
!------------------------------------------------------------------------------- |
24 |
|
|
! Arguments: |
25 |
|
|
REAL :: dyn_tropopause(klon) !--- Pressure at tropopause |
26 |
|
|
REAL, INTENT(IN) :: t(:,:) !--- Cells-centers temperature |
27 |
|
|
REAL, INTENT(IN) :: ts(:) !--- Surface temperature |
28 |
|
|
REAL, INTENT(IN) :: paprs(:,:) !--- Cells-edges pressure |
29 |
|
|
REAL, INTENT(IN) :: pplay(:,:) !--- Cells-centers pressure |
30 |
|
|
REAL, INTENT(IN) :: rot(:,:) !--- Cells-centers relative vorticity |
31 |
|
|
INTEGER, INTENT(OUT), OPTIONAL :: itrop(klon) !--- Last tropospheric layer idx |
32 |
|
|
REAL, INTENT(IN), OPTIONAL :: thet0, pvor0 |
33 |
|
|
!------------------------------------------------------------------------------- |
34 |
|
|
! Local variables: |
35 |
|
|
include "YOMCST.h" |
36 |
|
|
REAL, PARAMETER :: DynPTrMin =8.E+3 !--- Thresholds for minimum and maximum |
37 |
|
|
REAL, PARAMETER :: DynPTrMax =4.E+4 ! dynamical tropopause pressure (Pa). |
38 |
|
|
CHARACTER(LEN=80) :: sub |
39 |
|
|
INTEGER :: i, k, kb, kt, kp, ib, ie, nw |
40 |
|
|
REAL :: al, th0, pv0 |
41 |
|
|
REAL, DIMENSION(klon,klev) :: tpot_cen, tpot_edg, pvor_cen |
42 |
|
|
REAL, PARAMETER :: sg0=0.75 !--- Start level for PV=cte search loop |
43 |
|
|
INTEGER, PARAMETER :: nadj=3 !--- Adjacent levs nb for thresholds detection |
44 |
|
|
REAL, PARAMETER :: w(5)=[0.1,0.25,0.3,0.25,0.1] !--- Vertical smoothing |
45 |
|
|
INTEGER, SAVE :: k0 |
46 |
|
|
LOGICAL, SAVE :: first=.TRUE. |
47 |
|
|
!$OMP THREADPRIVATE(k0,first) |
48 |
|
|
!------------------------------------------------------------------------------- |
49 |
|
|
sub='dyn_tropopause' |
50 |
|
|
CALL assert(SIZE(t ,1)==klon, TRIM(sub)//" t klon") |
51 |
|
|
CALL assert(SIZE(t ,2)==klev, TRIM(sub)//" t klev") |
52 |
|
|
CALL assert(SIZE(ts,1)==klon, TRIM(sub)//" ts klon") |
53 |
|
|
CALL assert(SHAPE(paprs)==[klon,klev+1],TRIM(sub)//" paprs shape") |
54 |
|
|
CALL assert(SHAPE(pplay)==[klon,klev ],TRIM(sub)//" pplay shape") |
55 |
|
|
CALL assert(SHAPE(rot) ==[klon,klev ],TRIM(sub)//" rot shape") |
56 |
|
|
|
57 |
|
|
!--- DEFAULT THRESHOLDS |
58 |
|
|
th0=380.; IF(PRESENT(thet0)) th0=thet0 !--- In kelvins |
59 |
|
|
pv0= 2.; IF(PRESENT(pvor0)) pv0=pvor0 !--- In PVU |
60 |
|
|
IF(first) THEN |
61 |
|
|
DO k0=1,klev; IF(aps(k0)/preff+bps(k0)<sg0) EXIT; END DO; first=.FALSE. |
62 |
|
|
END IF |
63 |
|
|
|
64 |
|
|
!--- POTENTIAL TEMPERATURE AT CELLS CENTERS AND INTERFACES |
65 |
|
|
DO i = 1,klon |
66 |
|
|
tpot_cen(i,1) = t(i,1)*(preff/pplay(i,1))**RKAPPA |
67 |
|
|
tpot_edg(i,1) = ts(i) *(preff/paprs(i,1))**RKAPPA |
68 |
|
|
DO k=2,klev |
69 |
|
|
al = LOG(pplay(i,k-1)/paprs(i,k))/LOG(pplay(i,k-1)/pplay(i,k)) |
70 |
|
|
tpot_cen(i,k) = t(i,k) *(preff/pplay(i,k))**RKAPPA |
71 |
|
|
tpot_edg(i,k) = (t(i,k-1)+al*(t(i,k)-t(i,k-1)))*(preff/paprs(i,k))**RKAPPA |
72 |
|
|
!--- FORCE QUANTITIES TO BE GROWING |
73 |
|
|
IF(tpot_edg(i,k)<tpot_edg(i,k-1)) tpot_edg(i,k)=tpot_edg(i,k-1)+1.E-5 |
74 |
|
|
IF(tpot_cen(i,k)<tpot_cen(i,k-1)) tpot_cen(i,k)=tpot_cen(i,k-1)+1.E-5 |
75 |
|
|
END DO |
76 |
|
|
!--- VERTICAL SMOOTHING |
77 |
|
|
tpot_cen(i,:)=smooth(tpot_cen(i,:),w) |
78 |
|
|
tpot_edg(i,:)=smooth(tpot_edg(i,:),w) |
79 |
|
|
END DO |
80 |
|
|
|
81 |
|
|
!--- ERTEL POTENTIAL VORTICITY AT CELLS CENTERS (except in top layer) |
82 |
|
|
DO i = 1, klon |
83 |
|
|
DO k= 1, klev-1 |
84 |
|
|
pvor_cen(i,k)=-1.E6*RG*(rot(i,k)+2.*ROMEGA*SIN(latitude_deg(i)*RPI/180.))& |
85 |
|
|
* (tpot_edg(i,k+1)-tpot_edg(i,k)) / (paprs(i,k+1)-paprs(i,k)) |
86 |
|
|
END DO |
87 |
|
|
!--- VERTICAL SMOOTHING |
88 |
|
|
pvor_cen(i,1:klev-1)=smooth(pvor_cen(i,1:klev-1),w) |
89 |
|
|
END DO |
90 |
|
|
|
91 |
|
|
!--- LOCATE TROPOPAUSE: LOWEST POINT BETWEEN THETA=380K AND PV=2PVU SURFACES. |
92 |
|
|
DO i = 1, klon |
93 |
|
|
!--- UPPER TROPOPAUSE: |PV|=2PVU POINT STARTING FROM TOP |
94 |
|
|
DO kt=klev-1,1,-1; IF(ALL(ABS(pvor_cen(i,kt-nadj:kt))<=pv0)) EXIT; END DO |
95 |
|
|
!--- LOWER TROPOPAUSE: |PV|=2PVU POINT STARTING FROM BOTTOM |
96 |
|
|
DO kb=k0,klev-1; IF(ALL(ABS(pvor_cen(i,kb:kb+nadj))> pv0)) EXIT; END DO; kb=kb-1 |
97 |
|
|
!--- ISO-THETA POINT: THETA=380K STARTING FROM TOP |
98 |
|
|
DO kp=klev-1,1,-1; IF(ALL(ABS(tpot_cen(i,kp-nadj:kp))<=th0)) EXIT; END DO |
99 |
|
|
!--- CHOOSE BETWEEN LOWER AND UPPER TROPOPAUSE |
100 |
|
|
IF(2*COUNT(ABS(pvor_cen(i,kb:kt))>pv0)>kt-kb+1) kt=kb |
101 |
|
|
!--- PV-DEFINED TROPOPAUSE |
102 |
|
|
al = (ABS(pvor_cen(i,kt+1))-pv0)/ABS(pvor_cen(i,kt+1)-pvor_cen(i,kt)) |
103 |
|
|
dyn_tropopause(i) = pplay(i,kt+1)*(pplay(i,kt)/pplay(i,kt+1))**al |
104 |
|
|
!--- THETA=380K IN THE TROPICAL REGION |
105 |
|
|
al = (tpot_cen(i,kp+1)-th0)/(tpot_cen(i,kp+1)-tpot_cen(i,kp)) |
106 |
|
|
dyn_tropopause(i) = MAX( pplay(i,kp+1)*(pplay(i,kp)/pplay(i,kp+1))**al, & |
107 |
|
|
dyn_tropopause(i) ) |
108 |
|
|
!--- UNREALISTIC VALUES DETECTION |
109 |
|
|
IF(dyn_tropopause(i)<DynPTrMin.OR.dyn_tropopause(i)>DynPTrMax) THEN |
110 |
|
|
dyn_tropopause(i)=MIN(MAX(dyn_tropopause(i),DynPTrMax),DynPTrMin) |
111 |
|
|
DO kt=1,klev-1; IF(pplay(i,kt+1)>dyn_tropopause(i)) EXIT; END DO; kp=kt |
112 |
|
|
END IF |
113 |
|
|
#ifdef REPROBUS |
114 |
|
|
itroprep(i)=MAX(kt,kp) |
115 |
|
|
#endif |
116 |
|
|
!--- LAST TROPOSPHERIC LAYER INDEX NEEDED |
117 |
|
|
IF(PRESENT(itrop)) itrop(i)=MAX(kt,kp) |
118 |
|
|
END DO |
119 |
|
|
|
120 |
|
|
END FUNCTION dyn_tropopause |
121 |
|
|
|
122 |
|
|
|
123 |
|
|
!------------------------------------------------------------------------------- |
124 |
|
|
! |
125 |
|
|
FUNCTION smooth(v,w) |
126 |
|
|
! |
127 |
|
|
!------------------------------------------------------------------------------- |
128 |
|
|
! Arguments: |
129 |
|
|
REAL, INTENT(IN) :: v(:), w(:) |
130 |
|
|
REAL, DIMENSION(SIZE(v)) :: smooth |
131 |
|
|
!------------------------------------------------------------------------------- |
132 |
|
|
! Local variables: |
133 |
|
|
INTEGER :: nv, nw, k, kb, ke, lb, le |
134 |
|
|
!------------------------------------------------------------------------------- |
135 |
|
|
nv=SIZE(v); nw=(SIZE(w)-1)/2 |
136 |
|
|
DO k=1,nv |
137 |
|
|
kb=MAX(k-nw,1 ); lb=MAX(2+nw -k,1) |
138 |
|
|
ke=MIN(k+nw,nv); le=MIN(1+nw+nv-k,1+2*nw) |
139 |
|
|
smooth(k)=SUM(v(kb:ke)*w(lb:le))/SUM(w(lb:le)) |
140 |
|
|
END DO |
141 |
|
|
|
142 |
|
|
END FUNCTION smooth |
143 |
|
|
! |
144 |
|
|
!------------------------------------------------------------------------------- |
145 |
|
|
|
146 |
|
|
END MODULE tropopause_m |