Coverage for src/CSET/operators/convection.py: 78%
42 statements
« prev ^ index » next coverage.py v7.5.4, created at 2024-07-01 10:31 +0000
1# Copyright 2022-2023 Met Office and contributors.
2#
3# Licensed under the Apache License, Version 2.0 (the "License");
4# you may not use this file except in compliance with the License.
5# You may obtain a copy of the License at
6#
7# http://www.apache.org/licenses/LICENSE-2.0
8#
9# Unless required by applicable law or agreed to in writing, software
10# distributed under the License is distributed on an "AS IS" BASIS,
11# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12# See the License for the specific language governing permissions and
13# limitations under the License.
15"""A module containing different diagnostics for convection.
17The diagnostics are calculated from output from the Unified Model, although
18precalculated values in the required input form may also be used.
20"""
22import copy
23import logging
24import warnings
26import iris
27import iris.cube
28import numpy as np
31def cape_ratio(SBCAPE, MUCAPE, MUCIN, MUCIN_thresh=-75.0):
32 r"""Ratio of two fields, one filtered to allow physical values to be output.
34 Parameters
35 ----------
36 SBCAPE: Cube
37 Surface-based convective available potential energy as calculated by the
38 model. If using the UM please use STASH ``m01s20i114``
39 MUCAPE: Cube
40 Most-unstable convective available potential energy as calculated by the
41 model. If using the UM please use STASH ``m01s20i112``
42 MUCIN: Cube
43 Most-unstable convective inhibition associated with the most-unstable
44 ascent as calculated by the model. If using the UM please use STASH
45 ``m01s20i113``
46 MUCIN_thresh: float, optional, default is -75. J/kg.
47 Threshold to filter the MUCAPE by values are realistically realisable.
49 Returns
50 -------
51 Cube
53 Notes
54 -----
55 This diagnostic is based on Clark et al. (2012) [Clarketal2012]_. It is
56 based around the idea that for elevated convection the convective
57 instability is not based at the surface. This utilises two flavours of CAPE:
58 the surface-based CAPE (SBCAPE) and the most-unstable CAPE (MUCAPE). The
59 MUCAPE is filtered by the MUCIN associated with that parcel's ascent to
60 ensure that any CAPE can at least theoretically be released. The default
61 value is set at -75 J/kg but it can be changes depending on location and
62 users requirements.
64 .. math:: 1 - (\frac{SBCAPE}{MUCAPE})
66 The ratio is defined in this way such that if SBCAPE=MUCAPE the ratio will
67 equal 1. If the ratio was reversed when MUCAPE exists and SBCAPE is zero the
68 ratio would be undefined.
70 The diagnostic varies smoothly between zero and unity. A value of 0 implies
71 an environment is suitable for surface-based convection. A value of 1
72 implies an environment is suitable for elevated convection. Values between
73 imply transition convection with values closer to one imply elevated
74 convection is more likely and values closer to zero implying that
75 surface-based convection is more likely.
77 Further details about this diagnostic for elevated convection identification
78 can be found in Flack et al. (2023) [FlackCAPE2023]_.
80 Expected applicability ranges: Convective-scale models will be noisier than
81 parametrized models as they are more responsive to the convection, and thus
82 it may be more sensible to view as a larger spatial average rather than on
83 the native resolution.
85 Interpretation notes: UM stash for CAPE and CIN are calculated at the end of
86 the timestep. Therefore this diagnostic is applicable after precipitation
87 has occurred, not before as is the usual interpretation of CAPE related
88 diagnostics.
90 References
91 ----------
92 .. [Clarketal2012] Clark, A. J., Kain J. S., Marsh P. T., Correia J., Xue
93 M., and Kong F., (2012) "Forecasting tornado pathlengths using a
94 three-dimensional object identification algorithm applied to
95 convection-allowing forecasts." Weather and Forecasting, vol. 27,
96 1090–1113, doi: 10.1175/WAF-D-11-00147.1
97 .. [FlackCAPE2023] Flack, D.L.A., Lehnert, M., Lean, H.W., and Willington,
98 S. (2023) "Characteristics of Diagnostics for Identifying Elevated
99 Convection over the British Isles in a Convection-Allowing Model."
100 Weather and Forecasting, vol. 30, 1079-1094, doi: 10.1175/WAF-D-22-0219.1
102 Examples
103 --------
104 >>> CAPE_ratios=convection.cape_ratio(
105 SBCAPE,MUCAPE,MUCIN)
106 >>> iplt.pcolormesh(CAPE_ratios[0,:,:],cmap=mpl.cm.RdBu)
107 >>> plt.gca().coastlines('10m')
108 >>> plt.colorbar()
109 >>> plt.clim(0,1)
110 >>> plt.show()
112 >>> CAPE_ratios=convection.cape_ratio(
113 SBCAPE,MUCAPE,MUCIN,MUCIN_thresh=-1.5)
114 >>> iplt.pcolormesh(CAPE_ratios[0,:,:],cmap=mpl.cm.RdBu)
115 >>> plt.gca().coastlines('10m')
116 >>> plt.clim(0,1)
117 >>> plt.colorbar()
118 >>> plt.show()
119 """
120 # Load in the data into the new arrays.
121 MUCAPE_data = copy.deepcopy(MUCAPE.data)
122 if isinstance(MUCAPE_data, np.ma.MaskedArray): 122 ↛ 125line 122 didn't jump to line 125 because the condition on line 122 was always true
123 MUCAPE_data = MUCAPE_data.filled(np.nan)
124 # Remove all MUCAPE below MUCIN threshold.
125 MUCAPE_data[MUCIN.data <= MUCIN_thresh] = np.nan
126 with warnings.catch_warnings():
127 # Ignore possible divide by zero warnings, as they are replaced by NaNs.
128 warnings.filterwarnings("ignore", category=RuntimeWarning)
129 # Now calculate the main diagnostic.
130 EC_Flagb = 1 - (SBCAPE.data / MUCAPE_data)
131 if isinstance(EC_Flagb, np.ma.MaskedArray): 131 ↛ 126line 131 didn't jump to line 126
132 EC_Flagb = EC_Flagb.filled(np.nan)
133 # Filter to reduce NaN values and -inf values for plotting ease.
134 # There are multiple types of NaN values so need to convert them all to same type.
135 EC_Flagb[np.isnan(EC_Flagb)] = np.nan
136 EC_Flagb[np.isinf(EC_Flagb)] = np.nan
137 # Take the coordinates from an existing cube and replace the data.
138 cape_ratio_cube = SBCAPE.copy()
139 cape_ratio_cube.data = EC_Flagb
140 # Rename and remove STASH code.
141 cape_ratio_cube.var_name = "cape_ratio"
142 cape_ratio_cube.attributes.pop("STASH", None)
143 return cape_ratio_cube
146def inflow_layer_properties(EIB, BLheight, Orography):
147 r"""Filter to create a binary mask identifying elevated convection.
149 Parameters
150 ----------
151 EIB: Cube
152 Effective inflow layer base (precalculated or as identified by the
153 model). If using the UM please use STASH ``m01s20i119``.
154 BLheight: Cube
155 Boundary layer height (precalculated or as identified by the model). If
156 using the UM please use STASH ``m01s00i025``.
157 Orography: Cube
158 Model or actual orography, expected to be 2 dimensional. If 3 or 4
159 dimensional cube given converts to 2 dimensions assuming static
160 orography field in ensemble realization and time. If using the UM please
161 use STASH ``m01s00i033``.
163 Returns
164 -------
165 Cube
167 Notes
168 -----
169 This diagnostic is based on the concept of an effective inflow layer. This
170 concept was first introduced by Thompson et al. (2007) [Thompsonetal2007]_.
171 The inflow layer defined the region of air that is most likely to be
172 ingested into the convective event. It is defined by thresholding the CAPE
173 and CIN values: CAPE > 100 J/kg and \|CIN\| < 250 J/kg.
175 To turn this into a diagnostic for elevated convection the inflow layer base
176 is filtered against the boundary layer height. The model orography is added
177 to the boundary layer height to ensure reference height consistency as the
178 BL height is defined above ground level and the inflow layer base is defined
179 above sea level in the model output.
181 .. math:: EIB > BLheight + Orography
183 This is a binary diagnostic. It has a value of 0 to imply the environment is
184 suitable for surface-based convection. It has a value of 1 to indicate the
185 environment is suitable to produce elevated convection.
187 Further details about this diagnostic for elevated convection identification
188 can be found in Flack et al. (2023) [Flackinf2023]_.
190 Expected applicability ranges: Convective-scale models will be noisier than
191 parametrized models as they are more responsive to the convection, and thus
192 it may be more sensible to view as a larger spatial average rather than at
193 native resolution.
195 Interpretation notes: The effective inflow layer base diagnostic from UM
196 STASH is dependent upon the UM CAPE and CIN diagnostics. These diagnostics
197 are calculated at the end of the timestep. Therefore this diagnostic is
198 applicable after precipitation has occurred, not before as is the usual
199 interpretation of CAPE related diagnostics.
201 You might encounter warnings with the following text ``Orography assumed not
202 to vary with time or ensemble member.`` or ``Orography assumed not to vary
203 with time and ensemble member.`` these warnings are expected when the
204 orography files are not 2-dimensional, and do not cause any problems unless
205 ordering is not as expected.
207 References
208 ----------
209 .. [Thompsonetal2007] Thompson, R. L. Mead, C. M., and Edwards, R., (2007)
210 "Effective Storm-Relative Helicity and Bulk Shear in Supercell
211 Thunderstorm Environments." Weather and Forecasting, vol. 22, 102-115,
212 doi: 10.1175/WAF969.1
213 .. [Flackinf2023] Flack, D.L.A., Lehnert, M., Lean, H.W., and Willington, S.
214 (2023) "Characteristics of Diagnostics for Identifying Elevated
215 Convection over the British Isles in a Convection-Allowing Model."
216 Weather and Forecasting, vol. 30, 1079-1094, doi: 10.1175/WAF-D-22-0219.1
218 Examples
219 --------
220 >>> Inflow_properties=convection.inflow_layer_properties(EIB,BLheight,Orography)
221 >>> iplt.pcolormesh(Inflow_properties[0,:,:],cmap=mpl.cm.Purples)
222 >>> plt.gca().coastlines('10m')
223 >>> plt.colorbar()
224 >>> plt.clim(0,1)
225 >>> plt.show()
227 """
228 # Setup new array for output of the diagnostic.
229 EC_Flagd = np.zeros(EIB.shape)
230 # Check dimensions for Orography cube and replace with 2D array if not 2D.
231 if Orography.ndim == 3: 231 ↛ 232line 231 didn't jump to line 232 because the condition on line 231 was never true
232 try:
233 Orography = Orography.slices_over("realization").next()
234 except iris.exceptions.CoordinateNotFoundError:
235 Orography = Orography.slices_over("time").next()
236 logging.warning("Orography assumed not to vary with time or ensemble member")
237 elif Orography.ndim == 4: 237 ↛ 238line 237 didn't jump to line 238 because the condition on line 237 was never true
238 Orography = Orography.slices_over(("time", "realization")).next()
239 logging.warning("Orography assumed not to vary with time or ensemble member. ")
240 # Masked arrays are not respected, so convert masked values into NaNs.
241 if isinstance(EIB.data, np.ma.MaskedArray): 241 ↛ 245line 241 didn't jump to line 245 because the condition on line 241 was always true
242 EIB.data = EIB.data.filled(np.nan)
243 # Change points where Effective inflow layer base is larger than boundary
244 # layer height to 1 implying elevated convection.
245 EC_Flagd[EIB.data > (BLheight.data + Orography.data)] = 1.0
246 # Take the coordinates from an existing cube and replace the data.
247 inflow_properties_cube = EIB.copy()
248 inflow_properties_cube.data = EC_Flagd
249 # Rename and remove STASH code.
250 inflow_properties_cube.var_name = "inflow_layer_properties"
251 inflow_properties_cube.attributes.pop("STASH", None)
252 return inflow_properties_cube