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