Use virtual temperature to fix CAPE calculation

Implement Series.nearest() to allow finding value by nearest key in a
Series object

lib/xmet/series.py

@@ -50,3 +50,8 @@ class Series(dict):
 
             if sign is intersection:
                 return v1, pair[0]
+
+    def nearest(self, key):
+        k = sorted(self.keys(), key=lambda k: abs(k-key))[0]
+
+        return self[k]

lib/xmet/skew_t.py

@@ -260,7 +260,7 @@ class SkewTGraph():
 
         if self.opts.draw_cape_mixing_ratio:
             cr.set_source_rgb(0.8, 0.8, 1.0)
-            self.draw_series(cr, x, y, params.saturated_mr_line)
+            self.draw_series(cr, x, y, params.mixing_ratio_line)
 
     def draw(self,
              cr: cairo.Context,

lib/xmet/sounding.py

@@ -10,6 +10,7 @@ from xmet.thermo import follow_dry_adiabat, \
                         follow_moist_adiabat, \
                         follow_saturated_mixing_ratio, \
                         pressure_height, \
+                        virtual_temp, \
                         kelvin
 
 LAPSE_RATE_DRY   = 9.8 # degrees C per 1000m
@@ -233,19 +234,31 @@ class SoundingParams():
         return 9.8076 * cape
 
     def load_sounding(self, sounding: Sounding):
+        surface = sounding.samples[0]
+
+        tv = virtual_temp(surface.temp,
+                          surface.dewpoint,
+                          surface.pressure)
+
         temp_line     = sounding.follow_temp()
-        dry_adiabat = follow_dry_adiabat(sounding.samples[0].temp,
-                                         sounding.samples[0].pressure)
+        dewpoint_line = sounding.follow_dewpoint()
+        dry_adiabat   = follow_dry_adiabat(tv, surface.pressure)
 
-        saturated_mr_line = follow_saturated_mixing_ratio(sounding.samples[0].dewpoint,
-                                                          sounding.samples[0].pressure)
+        mixing_ratio_line = follow_saturated_mixing_ratio(surface.dewpoint,
+                                                          surface.pressure)
 
-        lcl = dry_adiabat.intersect(saturated_mr_line, SeriesIntersection.LESSER)
+        lcl = dry_adiabat.intersect(mixing_ratio_line,
+                                    SeriesIntersection.LESSER)
 
-        moist_adiabat = follow_moist_adiabat(*lcl)
+        tv2 = virtual_temp(temp_line.nearest(lcl[1]),
+                           dewpoint_line.nearest(lcl[1]),
+                           lcl[1])
+
+        moist_adiabat = follow_moist_adiabat(tv2, lcl[1])
 
         lfc = temp_line.intersect(moist_adiabat, SeriesIntersection.LESSER)
-        el  = temp_line.intersect(moist_adiabat, SeriesIntersection.GREATER, lfc[1]-100)
+        el  = temp_line.intersect(moist_adiabat, SeriesIntersection.GREATER,
+                                  lfc[1])
 
         self.lcl  = lcl
         self.lfc  = lfc
@@ -254,7 +267,7 @@ class SoundingParams():
 
         self.temp_line         = temp_line
         self.dry_adiabat       = dry_adiabat
-        self.saturated_mr_line = saturated_mr_line
+        self.mixing_ratio_line = mixing_ratio_line
         self.moist_adiabat     = moist_adiabat
 
     @staticmethod