Coverage for src/CSET/operators/radar_filter.py: 10%
60 statements
« prev ^ index » next coverage.py v7.15.1, created at 2026-07-14 13:50 +0000
1# © Crown copyright, Met Office (2022-2026) and CSET 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"""Operators to perform various kind of filtering."""
17import iris
18import iris.cube
19import iris.exceptions
21from CSET._common import iter_maybe
22from CSET.operators.filters import apply_mask, generate_mask
25def mask_list(model_names: list[str]) -> list[str]:
26 """Determine the Nimrod weights files to use.
28 Parameters
29 ----------
30 model_names: list[str]
31 A list of model and Nimrod hourly rainfall accumulation files.
33 Returns
34 -------
35 list[str]
36 A list of the Nimrod weights files to use with each of the input
37 model / observations files.
39 Notes
40 -----
41 At lest one of the entries in the input list must be a Nimrod hourly
42 rainfall accumulation file.
44 If just one Nimrod file is specified, then then the weights file associated
45 with this field is used.
47 If more than one Nimrod file is in the input list, then each of the these
48 Nimrod files is associated with its own weights file e.g. if the input list
49 contains ["Nimrod1km", "Nimrod2km"] then the weights files for these will
50 be ["Nimrod1km_weights", "Nimrod2km_weights"]. Any model fields in the input
51 list will be allocated a weights file according to the order of preference
52 specified in the list nimrod_preference e.g. if the input list is
53 ["UM_model", "Nimrod1km", "Nimrod2km"] then the output weights files list will
54 be ["Nimrod2km_weights", "Nimrod1km_weights", "Nimrod2km_weights"] as the Nimrod
55 weights for 2km data are preferred over those for 1km.
57 Examples
58 --------
59 >>> list_weights = make_list( ["UM_model", "Nimrod1km", "Nimrod2km"] )
60 >>> print(list_weights)
61 ["Nimrod2km_weights", "Nimrod1km_weights", "Nimrod2km_weights"]
63 """
64 # Set the preference order for choosing a Nimrod radar weights source
65 # in order of most to least preferred.
66 nimrod_preference = [
67 "Nimrod2km",
68 "Nimrod_2km",
69 "Nimrodxkm",
70 "Nimrod_xkm",
71 "Nimrod1km",
72 "Nimrod_1km",
73 ]
75 # Define the string that helps form a Nimrod weights file.
76 wei = "_weights"
78 # Determine the preferred Nimrod mask to use.
79 empty_string = ""
80 preferred_nimrod = empty_string
81 for prefer in reversed(nimrod_preference):
82 if any(prefer in model for model in model_names):
83 preferred_nimrod = prefer
85 # Create the list of the required Nimrod masks.
86 mask_names_list = []
87 if preferred_nimrod != empty_string:
88 # Loop over the input model_names.
89 for model in model_names:
90 if any(model in nimrod for nimrod in nimrod_preference):
91 nimrod_mask = model + wei
92 else:
93 nimrod_mask = preferred_nimrod + wei
94 mask_names_list.append(nimrod_mask)
96 return mask_names_list
99def mask_by_weights(
100 cubes: iris.cube.CubeList,
101 model_names: list[str],
102 weights_names: list[str],
103 **kwargs,
104) -> iris.cube.CubeList:
105 """Filter a field using a second field as a mask.
107 Parameters
108 ----------
109 cubes: iris.cube.CubeList
110 Two cubes containing the radar observations and their weights.
112 Returns
113 -------
114 Cube
116 Raises
117 ------
118 ValueError, iris.exceptions.NotYetImplementedError
119 When the cubes are not compatible.
121 Notes
122 -----
123 This is a simple operator designed for combination of diagnostics or
124 creating new diagnostics by using recipes.
126 Examples
127 --------
128 >>> field_filtered = mask_by_weights(cubelist, model_names)
130 """
131 print("model_names are: ", model_names)
132 print("weights_names", weights_names)
134 for cube in cubes:
135 print(" cube.var_name ", cube.var_name)
136 print(" cube.name ", cube.name)
137 print(" cube: ")
138 print(cube)
139 print(" cube.attributes.model_name ", cube.attributes["model_name"])
141 # Check the input unfiltered cubes and the mask cubes are both cubelists
142 # with the same number of cubes. If not, then add extra mask cubes.
143 if len(model_names) != len(weights_names):
144 weights_names = mask_list(model_names)
146 # Create an empty cubelist to hold the filtered fields.
147 filtered_list = iris.cube.CubeList([])
149 # Loop over the fields to filter.
150 var_constraint = iris.NameConstraint(var_name="hourly_rain_accumulation")
151 mask_var_constraint = iris.NameConstraint(var_name="hourly_wts_accumulation")
152 for model, mask in zip(
153 iter_maybe(model_names),
154 iter_maybe(weights_names),
155 strict=True,
156 # iter_maybe(model_names), iter_maybe(weights_names), strict=True
157 ):
158 print(" model, mask ", model, mask)
160 # grab the field to filter
161 model_constraint = iris.AttributeConstraint(model_name=model)
162 unfiltered_field = cubes.extract_cube(var_constraint & model_constraint)
164 # Select the field to use as the mask.
165 # Nice to do - put in support for a static mask.
166 mask_constraint = iris.AttributeConstraint(model_name=mask)
167 mask_field = cubes.extract_cube(mask_var_constraint & mask_constraint)
169 # Create the mask - note that the condition e.g. "ge" can be set by a loader
170 # as can the threshold value.
171 mask_radar_wts = generate_mask(mask_field, "ge", 11)
173 # print(" This is cube radar_obs: ", radar_obs)
174 # print(" This is cube radar_weights: ", radar_wts)
175 # print(" This is cube unfiltered: ", unfiltered_field)
177 # check the coords of the unfiltered field and the mask field.
178 # If these do not match, then regrid the unfiltered field onto
179 # the grid used for the mask field.
180 # For radar weights fields can use the function regrid_onto_xyspacing in regrid.py,
181 # but then might have to extract a subarea to match the mask grid.
182 # Might have to consider serval cases for regridding:
183 # (1) model_field(lat, lon) to radar_weights_field(x, y)
184 # (2) model_field(lat, lon) to other_model_field(lat, lon)
185 # (3) Nimrod_field(x, y) to radar_weights_field(x, y)
186 # (4) Nimrod_field(x, y) to model_field(lat,lon) ?
187 #
189 # Apply the mask.
190 masked_radar_obs = apply_mask(unfiltered_field, mask_radar_wts)
192 # Put the filtered cube into the list of filtered cubes.
193 filtered_list.append(masked_radar_obs)
195 # Preserve returning a cube if only a cube has been supplied to filter.
196 if len(filtered_list) == 1:
197 return filtered_list[0]
198 else:
199 return filtered_list
202def radar_apply_mask(
203 original_field: iris.cube.Cube | iris.cube.CubeList,
204 mask: iris.cube.Cube | iris.cube.CubeList,
205) -> iris.cube.Cube | iris.cube.CubeList:
206 """Apply a mask to given data as a masked array.
208 Parameters
209 ----------
210 original_field: iris.cube.Cube | iris.cube.CubeList
211 The field(s) to be masked.
212 mask: iris.cube.Cube | iris.cube.CubeList
213 The mask(s) being applied to the original field(s).
215 Returns
216 -------
217 masked_field: iris.cube.Cube | iris.cube.CubeList
218 A cube or cubelist of the masked field(s).
220 Notes
221 -----
222 The mask is first converted to 1s and NaNs before multiplication with
223 the original data.
225 As discussed in generate_mask, you can combine multiple masks in a
226 recipe using other functions before applying the mask to the data.
228 Examples
229 --------
230 >>> land_points_only = radar_apply_mask( surface_microphysical_rainfall_rate, Nimrod2km)
231 """
232 # Create an empty cubelist to hold the filtered fields.
233 masked_fields = iris.cube.CubeList([])
235 # Loop over the input mask and field cubes.
236 for M, F in zip(iter_maybe(mask), iter_maybe(original_field), strict=True):
237 masked_field = F.copy()
239 # If the field and mask are on different grids, then regrid the field.
240 if M[0].shape != masked_field[0].shape:
241 masked_field = masked_field.regrid(M, iris.analysis.Linear())
243 # Apply the mask.
244 min_timesteps = min(M.shape[0], masked_field.shape[0])
245 masked_field = apply_mask(masked_field[0:min_timesteps], M[0:min_timesteps])
247 # Attach and attribute to the masked field detailing the mask used.
248 masked_field.attributes["mask"] = f"mask_of_{F.name()}"
250 # Append the masked field to the output list of masked fields.
251 masked_fields.append(masked_field)
253 # Return either a single cube or a cubelist.
254 if len(masked_fields) == 1:
255 return masked_fields[0]
256 else:
257 # return masked_fields
258 return masked_fields.merge()