Introduction of Icing Multivariate Regression #2187
Conversation
A previous merge created extraneous lines for integrate-time-bounds test file which are not even part of the icing multivariate regression. Put the correct lines in place to match the master branch,
There was a problem hiding this comment.
Thanks @alexproschitsky
I've made some suggestions that I think will improve what is already good code.
| IF RH is greater than 70% and T is between 250.0K and 273.15K THEN: | ||
| AISI=100*TANH(0.06*RH-4.0)[TANH(0.1(T-247.0)] | ||
|
|
||
| Categorical icing severity levels are defined as | ||
| AISI < 58 : "No Icing" | ||
| 58 <= AISI < 85 : "Light Icing" | ||
| 85 <= AISI < 92 : "Moderate Icing" | ||
| 92 <= AISI : "Severe Icing" |
There was a problem hiding this comment.
Is there a published paper that we can refer to here?
| output_cubes.append(self._extract_input(cubes, input_name)) | ||
| if not output_cubes[-1].units in units: | ||
| expected_unit_string = " or ".join(map(str, units)) | ||
| received_unit_string = str(output_cubes[-1].units) | ||
| raise ValueError( | ||
| f"The {output_cubes[-1].name()} units are incorrect, expected " | ||
| f"units as {expected_unit_string} but received {received_unit_string})." | ||
| ) |
There was a problem hiding this comment.
If you attempt to convert the cube to the units you want, it will either work, and you can use whatever you've got (Fahrenheit becomes Kelvin) or an error message is raised (Metres cannot become Kelvin).
Also, do this before appending the result to get nicer (well, I think it's nicer) code.
| output_cubes.append(self._extract_input(cubes, input_name)) | |
| if not output_cubes[-1].units in units: | |
| expected_unit_string = " or ".join(map(str, units)) | |
| received_unit_string = str(output_cubes[-1].units) | |
| raise ValueError( | |
| f"The {output_cubes[-1].name()} units are incorrect, expected " | |
| f"units as {expected_unit_string} but received {received_unit_string})." | |
| ) | |
| cube = self._extract_input(cubes, input_name) | |
| cube.convert_units(units) | |
| output_cubes.append(cube) |
| f"units as {expected_unit_string} but received {received_unit_string})." | ||
| ) | ||
|
|
||
| t,rh = output_cubes |
There was a problem hiding this comment.
This isn't safe, as dictionaries do not guarantee ordering. As you only have two inputs, probably easiest and clearest to unroll this loop. Also t and rh are very brief and open to misinterpretation. I recommend spelling out temperature and relative_humidity for clarity.
| if t.coord("forecast_reference_time") != rh.coord("forecast_reference_time"): | ||
| raise ValueError( | ||
| f"{t.name()} and {rh.name()} do not have the same forecast reference time" | ||
| ) |
There was a problem hiding this comment.
Is this test required?
What I'm thinking is that "forecast_reference_time" is not the only way of defining this. Some IMPROVER cubes have "blend_time" instead. While I admit that I would not recommend applying this function to blended data, I think it would be more important to check that both have the same scalar coords, rather than specifying forecast reference time and time specifically. How about height or pressure? Should these be present and should they match? Would it affect the algorithm if they did not?
| """ | ||
| try: | ||
| cube = cubes.extract_cube(iris.Constraint(cube_name)) | ||
| except iris.exceptions.ConstraintMismatchError: |
There was a problem hiding this comment.
Is it necessary to trap this exception? Is the error message from iris not sufficient? If the iris message is adequate, this whole method can be removed and CubeList.extract_cube can be called directly.
| # Regression equation if RH is greater than 70% and T is between 250.0K and 273.15K | ||
| aisi = 100*np.tanh(0.06*rh-4.0)*(np.tanh(0.1*(t-247.0))) | ||
| aisi[np.where(rh<70)] = 0 | ||
| aisi[np.where(t<250.0)] = 0 | ||
| aisi[np.where(t>273.15)] = 0 |
There was a problem hiding this comment.
I recommend putting the actual calculations into a separate class-method. Then the numbers can be tested in isolation. Also, I think the equation isn't Ruff-standard. Try installing and running the pre-commit hooks described in README.md.
| ) | ||
| from improver.utilities.cube_checker import spatial_coords_match | ||
|
|
||
| class IcingSeverityMultivariateRegression_USAF2024(PostProcessingPlugin): |
There was a problem hiding this comment.
I recommend against including the source and year in the class name. Github tracks the authors already.
| class IcingSeverityMultivariateRegression_USAF2024(PostProcessingPlugin): | |
| class IcingSeverityMultivariateRegression(PostProcessingPlugin): |
| @@ -0,0 +1,71 @@ | |||
| #!/usr/bin/env python | |||
There was a problem hiding this comment.
I recommend removing usaf2024 from this file name too.
| assert result.xml().splitlines(keepends=True) == expected_cube.xml().splitlines( | ||
| keepends=True | ||
| ) | ||
| assert np.allclose(result.data, expected_cube.data) |
There was a problem hiding this comment.
What tolerance is sensible here? We usually specify something to avoid false-positive errors when underlying packages change the behaviour slightly.
| assert np.allclose(result.data, expected_cube.data) | |
| assert np.allclose(result.data, expected_cube.data, rtol=1e-4) |
|
@MoseleyS Thank you for taking your time to review. I will address these suggestions. Alex. |
2 participants
Addresses #GitHubissuenum
Adding MIT Lincoln Laboratory code to produce Icing Severity Multivariate Regression.
Data for acceptance tests have been provided via PR metoppv/improver_test_data#106
Please advise of any necessary changes.
Alex Proschitsky
Testing:
CLA