Initial implementation of bulk shear

lib/xmet/sounding.py

@@ -1,3 +1,4 @@
+import math
 import datetime
 import shapely
 
@@ -11,7 +12,11 @@ from xmet.thermo import follow_dry_adiabat, \
                         follow_saturated_mixing_ratio, \
                         pressure_height, \
                         virtual_temp, \
-                        kelvin
+                        kelvin, \
+                        wind_u, \
+                        wind_v
+
+from xmet.units import deg
 
 LAPSE_RATE_DRY   = 9.8 # degrees C per 1000m
 LAPSE_RATE_MOIST = 4.0
@@ -204,6 +209,43 @@ class Sounding(DatabaseTable):
 
         return series
 
+    def shear(self, level: int):
+        s1 = self.samples[level]
+        s2 = self.samples[level+1]
+
+        u1 = wind_u(s1.wind_speed, s1.wind_dir)
+        u2 = wind_u(s2.wind_speed, s2.wind_dir)
+
+        v1 = wind_v(s1.wind_speed, s1.wind_dir)
+        v2 = wind_v(s2.wind_speed, s2.wind_dir)
+
+        zd = s2.height - s1.height
+
+        ud = (u2 - u1) / zd
+        vd = (v2 - v1) / zd
+
+        return (
+            ((ud**2 + vd**2)**0.5) / zd,
+            90 - deg(math.atan(vd / ud))
+        )
+
+    def bulk_shear(self):
+        levels = len(self.samples) - 1
+
+        shear_speed = 0.0
+        shear_dir   = 0.0
+
+        for level in range(0, levels):
+            shear = self.shear(level)
+
+            shear_speed += shear[0]
+            shear_dir   += shear[1]
+
+        return (
+            shear_speed / levels,
+            shear_dir   / levels
+        )
+
 def between(n, a, b):
     return n > a and n < b
 

lib/xmet/thermo.py

@@ -64,8 +64,9 @@ def mixing_ratio_temp(ratio: float, pressure: float) -> float:
     w = ratio * 1000.0
     x = math.log10(w * pressure / (622.0 + w))
 
-    return celsius(math.pow(10.0, ((c1 * x) + c2)) - c3 + \
-             (c4 * math.pow((math.pow(10, (c5 * x)) - c6), 2)))
+    temp = (10.0 ** ((c1 * x) + c2)) - c3 + (c4 * ((10 ** (c5 * x)) - c6) ** 2)
+
+    return celsius(temp)
 
 def lcl(temp: float, dewpoint: float) -> float:
     """