Fix Open-MSS#2930 (Open-MSS#2931)

Implement a standard atmosphere up to 100km altitude, raising the
ceiling from the current 71km.

Author: joernu76

mslib/utils/thermolib.py

@@ -126,14 +126,17 @@ def omega_to_w(omega, p, t):
 
 # Values according to the 1976 U.S. Standard atmosphere [NOAA1976]_.
 # List of tuples (height, temperature, pressure, temperature gradient)
+# Added a value for 100km roughly based on CIRA
 _STANDARD_ATMOSPHERE = [
     (0 * units.km, 288.15 * units.K, 101325 * units.Pa, 0.0065 * units.K / units.m),
     (11 * units.km, 216.65 * units.K, 22632.1 * units.Pa, 0 * units.K / units.m),
     (20 * units.km, 216.65 * units.K, 5474.89 * units.Pa, -0.001 * units.K / units.m),
     (32 * units.km, 228.65 * units.K, 868.019 * units.Pa, -0.0028 * units.K / units.m),
     (47 * units.km, 270.65 * units.K, 110.906 * units.Pa, 0 * units.K / units.m),
     (51 * units.km, 270.65 * units.K, 66.9389 * units.Pa, 0.0028 * units.K / units.m),
-    (71 * units.km, 214.65 * units.K, 3.95642 * units.Pa, np.nan * units.K / units.m)
+    (71 * units.km, 214.65 * units.K, 3.95642 * units.Pa, 0.002 * units.K / units.m),
+    (84.852 * units.km, 186.95 * units.K, 0.3734 * units.Pa, 0 * units.K / units.m),
+    (100 * units.km, 186.95 * units.K, 0.025641 * units.Pa, np.nan * units.K / units.m),
 ]
 _HEIGHT, _TEMPERATURE, _PRESSURE, _TEMPERATURE_GRADIENT = 0, 1, 2, 3
 
@@ -235,8 +238,8 @@ def pressure2flightlevel(pressure):
             z[indices] = z0 - (Rd * t0) / g * np.log(pressure[indices] / p0)
 
     if np.isnan(z).any():
-        raise ValueError("flight level to pressure conversion not "
-                         "implemented for z > 71km")
+        raise ValueError("pressure to flight level conversion not "
+                         "implemented for p < 0.025641 Pa")
 
     return z if is_array else z[0]
 
@@ -262,7 +265,7 @@ def isa_temperature(height):
             return t0 - gamma * (height - z0)
 
     raise ValueError("ISA temperature from flight level not "
-                     "implemented for z > 71km")
+                     "implemented for z > 100 km")
 
 
 def convert_pressure_to_vertical_axis_measure(vertical_axis, pressure):

tests/_test_utils/test_thermolib.py

@@ -33,14 +33,14 @@
 
 
 def test_flightlevel2pressure2flightlevel():
-    fs = (np.arange(1, 71000, 1000.) * units.m).to(units.hft)
+    fs = (np.arange(1, 95000, 1000.) * units.m).to(units.hft)
     ps = thermolib.flightlevel2pressure(fs)
     fs_p = thermolib.pressure2flightlevel(ps).magnitude
     assert fs.magnitude == pytest.approx(fs_p)
 
 
 def test_pressure2flightlevel2pressure():
-    ps = np.arange(5, 105000, 1.)[::-1] * units.Pa
+    ps = np.arange(1, 105000, 1.)[::-1] * units.Pa
     fs = thermolib.pressure2flightlevel(ps)
     ps_p = thermolib.flightlevel2pressure(fs).magnitude
     assert ps.magnitude == pytest.approx(ps_p)
@@ -54,7 +54,7 @@ def test_flightlevel2pressure():
     assert thermolib.flightlevel2pressure(1626.8966 * units.hft).magnitude == pytest.approx(80)
     assert thermolib.flightlevel2pressure(1804.2727 * units.hft).magnitude == pytest.approx(40)
     with pytest.raises(ValueError):
-        thermolib.flightlevel2pressure(72000 * units.m)
+        thermolib.flightlevel2pressure(102000 * units.m)
 
 
 def test_pressure2flightlevel():
@@ -67,7 +67,7 @@ def test_pressure2flightlevel():
     assert thermolib.pressure2flightlevel(80 * units.Pa).magnitude == pytest.approx(1626.8966)
     assert thermolib.pressure2flightlevel(40 * units.Pa).magnitude == pytest.approx(1804.2727)
     with pytest.raises(ValueError):
-        thermolib.pressure2flightlevel(3.9 * units.Pa)
+        thermolib.pressure2flightlevel(0.02 * units.Pa)
 
 
 def test_isa_temperature():
@@ -81,12 +81,14 @@ def test_isa_temperature():
     assert thermolib.isa_temperature(800 * units.hft).magnitude == pytest.approx(221.034)
     assert thermolib.isa_temperature(1000 * units.hft).magnitude == pytest.approx(227.13)
     with pytest.raises(ValueError):
-        thermolib.isa_temperature(71001 * units.m)
+        thermolib.isa_temperature(100001 * units.m)
     assert thermolib.isa_temperature(11000 * units.m).magnitude == pytest.approx(216.65)
     assert thermolib.isa_temperature(20000 * units.m).magnitude == pytest.approx(216.65)
     assert thermolib.isa_temperature(32000 * units.m).magnitude == pytest.approx(228.65)
     assert thermolib.isa_temperature(47000 * units.m).magnitude == pytest.approx(270.65)
     assert thermolib.isa_temperature(51000 * units.m).magnitude == pytest.approx(270.65)
+    assert thermolib.isa_temperature(71000 * units.m).magnitude == pytest.approx(214.65)
+    assert thermolib.isa_temperature(84852 * units.m).magnitude == pytest.approx(186.95)
 
 
 class TestConverter: