-
Notifications
You must be signed in to change notification settings - Fork 1
/
precip_class.py
91 lines (80 loc) · 2.77 KB
/
precip_class.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
#
# Column-based precipitation classification algorithm designed for application on
# numerical model output.
#
# It has been designed using WRF model output using the Thompson and Eidhammer
# (2014, JAS) microphysics scheme, which has 2 liquid and 3 frozen categories as
# listed and expected below.
#
# Input:
#
# Q_INT: n-D array of vertically integrated hydrometeors as f(q, X), where
# q(5) is the hydrometeor dimension, arranged as
# ['QCLOUD', 'QRAIN', 'QICE', 'QSNOW', 'QGRAUP'] and X includes the
# remaining (time and) spatial dimensions.
# Returns:
#
# C_TYPE: (n-2)-D array as f(X) with classification results:
# 0: non-cloud
# Convective:
# 1: deep convective
# 2: congestus
# 3: shallow
# Layered:
# 4: stratiform
# 5: anvil (weaker rainfall)
#
# Emily Luschen - [email protected]
# James Ruppert - [email protected]
# 5/19/23
import numpy as np
def precip_class(q_int):
shape = q_int.shape
ndims=len(shape)
shape_out = shape[1:ndims]
# Integrated water variables
LWP = q_int[0] + q_int[1] # Liquid water path = cloud + rain
IWP = q_int[2] + q_int[3] + q_int[4] # Ice water path = ice + snow + graupel
TWP = LWP + IWP # Total water path [mm]
IWP = np.ma.masked_where((LWP == 0), IWP, copy=False)
# Threshold p]arameters
twp_thresh = 1e-1
cr_thresh = 2
# ice_thresh = 1e-8
graup_thresh = 1e-4
rain_thresh_conv = 1e-1
rain_thresh_strat = 1e-2
# Initialize output array
if np.ma.is_masked(q_int):
c_type = np.ma.zeros(shape_out)
domask=True
else:
c_type = np.zeros(shape_out)
domask=False
cr = IWP/LWP
# Deep convection
c_type[( ((LWP != 0) & (TWP > twp_thresh)) &
(cr <= cr_thresh) &
(q_int[1] >= rain_thresh_conv) &
(q_int[4] >= graup_thresh) ).nonzero() ] = 1
# Congestus
c_type[( ((LWP != 0) & (TWP > twp_thresh)) &
(cr <= cr_thresh) &
(q_int[1] >= rain_thresh_conv) &
(q_int[4] < graup_thresh) ).nonzero() ] = 2
# Shallow
c_type[( ((LWP != 0) & (TWP > twp_thresh)) &
(cr <= cr_thresh) &
(q_int[1] < rain_thresh_conv) ).nonzero() ] = 3
# Stratiform
c_type[( ((LWP != 0) & (TWP > twp_thresh)) &
(cr > cr_thresh) &
(q_int[1] >= rain_thresh_strat) ).nonzero() ] = 4
# Anvil
c_type[( ((LWP != 0) & (TWP > twp_thresh)) &
(cr > cr_thresh) &
(q_int[1] < rain_thresh_strat) ).nonzero() ] = 5
# Fill in mask of original array if it exists
if domask:
c_type.mask = q_int.mask[0,...]
return c_type