28 subroutine maximum_hourly_diagnostics_run(im, levs, reset, lradar, imp_physics, &
29 imp_physics_gfdl, imp_physics_thompson, &
30 imp_physics_fer_hires, imp_physics_nssl, &
32 gt0, refl_10cm, refdmax, refdmax263k, u10m, v10m, &
33 u10max, v10max, spd10max, pgr, t2m, q2m, t02max, &
34 t02min, rh02max, rh02min, dtp, rain, pratemax, &
35 lightning_threat, ltg1_max,ltg2_max,ltg3_max, &
36 wgrs, prsi, qgraupel, qsnowwat, qicewat, tgrs, con_rd,&
37 prsl, kdt, errmsg, errflg)
40 integer,
intent(in) :: im, levs, kdt
41 logical,
intent(in) :: reset, lradar, lightning_threat
42 integer,
intent(in) :: imp_physics, imp_physics_gfdl, imp_physics_thompson, imp_physics_fer_hires, &
44 real(kind_phys),
intent(in ) :: con_g
45 real(kind_phys),
intent(in ) :: con_rd
46 real(kind_phys),
intent(in ) :: phil(:,:)
47 real(kind_phys),
intent(in ) :: gt0(:,:)
48 real(kind_phys),
intent(in ) :: refl_10cm(:,:)
49 real(kind_phys),
intent(inout) :: refdmax(:)
50 real(kind_phys),
intent(inout) :: refdmax263k(:)
51 real(kind_phys),
intent(in ) :: u10m(:)
52 real(kind_phys),
intent(in ) :: v10m(:)
53 real(kind_phys),
intent(inout) :: u10max(:)
54 real(kind_phys),
intent(inout) :: v10max(:)
55 real(kind_phys),
intent(inout) :: spd10max(:)
56 real(kind_phys),
intent(in ) :: pgr(:)
57 real(kind_phys),
intent(in ) :: t2m(:)
58 real(kind_phys),
intent(in ) :: q2m(:)
59 real(kind_phys),
intent(inout) :: t02max(:)
60 real(kind_phys),
intent(inout) :: t02min(:)
61 real(kind_phys),
intent(inout) :: rh02max(:)
62 real(kind_phys),
intent(inout) :: rh02min(:)
63 real(kind_phys),
intent(in ) :: dtp
64 real(kind_phys),
intent(in ) :: rain(:)
65 real(kind_phys),
intent(in ) :: tgrs(:,:)
66 real(kind_phys),
intent(in ) :: prsl(:,:)
67 real(kind_phys),
intent(inout) :: pratemax(:)
69 real(kind_phys),
intent(in),
dimension(:,:) :: prsi, qgraupel, qsnowwat, qicewat
70 real(kind_phys),
intent(in),
dimension(:,:),
optional :: wgrs
71 real(kind_phys),
intent(inout),
dimension(:),
optional :: ltg1_max, ltg2_max, ltg3_max
72 character(len=*),
intent(out) :: errmsg
73 integer,
intent(out) :: errflg
76 real(kind_phys),
dimension(:),
allocatable :: refd, refd263k
77 real(kind_phys) :: tem, pshltr, qcq, rh02, dp, q
85 if (lightning_threat)
then
86 call lightning_threat_indices
90 if (lradar .and. (imp_physics == imp_physics_gfdl .or. &
91 imp_physics == imp_physics_thompson .or. &
92 imp_physics == imp_physics_fer_hires .or. &
93 imp_physics == imp_physics_nssl ))
then
95 allocate(refd263k(im))
96 call max_fields(phil,refl_10cm,con_g,im,levs,refd,gt0,refd263k)
98 IF ( imp_physics == imp_physics_nssl )
THEN
106 refdmax263k(i) = -35.
111 refdmax(i) = max(refdmax(i),refd(i))
112 refdmax263k(i) = max(refdmax263k(i),refd263k(i))
115 deallocate (refd263k)
132 tem = sqrt(u10m(i)*u10m(i) + v10m(i)*v10m(i))
133 if (tem > spd10max(i))
then
138 pshltr=pgr(i)*exp(-0.068283/gt0(i,1))
139 qcq=pq0/pshltr*exp(a2a*(t2m(i)-a3)/(t2m(i)-a4))
144 IF (rh02 < rhmin)
THEN
147 rh02max(i) = max(rh02max(i),rh02)
148 rh02min(i) = min(rh02min(i),rh02)
149 t02max(i) = max(t02max(i),t2m(i))
150 t02min(i) = min(t02min(i),t2m(i))
151 pratemax(i) = max(pratemax(i),(3.6e6/dtp)*rain(i))
156 subroutine lightning_threat_indices
158 REAL(kind_phys),
PARAMETER :: clim1=1.50
159 REAL(kind_phys),
PARAMETER :: clim2=0.40*1.22
160 REAL(kind_phys),
PARAMETER :: clim3=0.02*1.22
167 REAL(kind_phys),
PARAMETER :: coef1=0.042*1000.*1.22
168 REAL(kind_phys),
PARAMETER :: coef2=0.20*1.22
170 REAL(kind_phys) :: totice_colint(im), ltg1, ltg2, high_ltg1, high_wgrs, high_graupel, rho
171 LOGICAL :: ltg1_calc(im)
172 integer :: k, i,
count
183 dp = prsi(i,k) - prsi(i,k+1)
184 q = qgraupel(i,k) + qsnowwat(i,k) + qicewat(i,k)
185 rho = prsl(i,k) / (con_rd * tgrs(i,k))
186 totice_colint(i) = totice_colint(i) + q * rho * dp / con_g
188 IF ( .not.ltg1_calc(i) )
THEN
189 IF ( 0.5*(tgrs(i,k+1) + tgrs(i,k)) < 258.15 )
THEN
191 ltg1_calc(i) = .true.
193 ltg1 = coef1*wgrs(i,k)* &
194 (( qgraupel(i,k+1) + qgraupel(i,k) )*0.5 )
195 if(ltg1 > high_ltg1)
then
197 high_graupel = qgraupel(i,k)
198 high_wgrs = wgrs(i,k)
201 IF ( ltg1 .LT. clim1 ) ltg1 = 0.
204 ltg1 = ltg1 * scaling_factor
206 IF ( ltg1 .GT. ltg1_max(i) )
THEN
215 ltg2 = coef2 * totice_colint(i)
217 IF ( ltg2 .LT. clim2 ) ltg2 = 0.
220 ltg2 = ltg2 * scaling_factor
222 IF ( ltg2 .GT. ltg2_max(i) )
THEN
227 ltg3_max(i) = 0.95 * ltg1_max(i) + 0.05 * ltg2_max(i)
230 IF ( ltg3_max(i) .LT. clim3 * scaling_factor ) ltg3_max(i) = 0.
233 end subroutine lightning_threat_indices
237 subroutine max_fields(phil,ref3D,grav,im,levs,refd,tk,refd263k)
238 integer,
intent(in) :: im,levs
239 real (kind=kind_phys),
intent(in) :: grav
240 real (kind=kind_phys),
intent(in),
dimension(:,:) :: phil,ref3d,tk
241 integer :: i,k,ll,ipt,kpt
242 real :: dbz1avg,zmidp1,zmidloc,refl,fact
243 real,
dimension(im,levs) :: z
244 real,
dimension(im) :: zintsfc
245 real,
dimension(:),
intent(inout) :: refd,refd263k
246 REAL :: dbz1(2),dbzk,dbzk1
250 z(i,k)=phil(i,k)/grav
256 if ( z(i,k+1) >= 1000. .and. z(i,k) <= 1000.)
then
274 if (dbz1(ll)>-35.) refl=10.**(0.1*dbz1(ll))
279 fact=(1000.-zmidloc)/(zmidloc-zmidp1)
280 dbz1avg=dbz1(2)+(dbz1(2)-dbz1(1))*fact
282 if (dbz1avg>0.01)
then
283 dbz1avg=10.*alog10(dbz1avg)
287 refd(i)=max(refd(i),dbz1avg)
294 vloopm10:
do k=1,levs-1
295 if (tk(i,k+1) .le. 263.15 .and. tk(i,k) .ge. 263.15)
then
304 if (dbz1(ll)>-35.) refl=10.**(0.1*dbz1(ll))
311 if (dbz1avg>0.01)
then
312 dbz1avg=10.*alog10(dbz1avg)