Created a new function to calculate the SOMO metrics.

toarstats/ozone_metrics.py

@@ -42,8 +42,8 @@ import pandas as pd
 from toarstats.stats_utils import (aot40stat, dma8_processor,
                                    get_elevation_angle, kth_highest,
                                    prepare_data, resample_with_date,
-                                   stat_processor_1, threshold_processor,
-                                   w126stat)
+                                   somo_processor, stat_processor_1,
+                                   threshold_processor, w126stat)
 
 
 def max1h_values(ser, ref, mtype, metadata, seasons, data_capture):
@@ -301,17 +301,8 @@ def somo10(ser, ref, mtype, metadata, seasons, data_capture):
 
     :return: A list of dictionaries containing the processed data
     """
-    tmp = dma8eu(ser, ref, "daily", metadata, seasons, data_capture)
-    if mtype == "daily":
-        res = tmp
-        res[0]["name"] = "somo10"
-    else:
-        new_ser = tmp[0]["ser"] - 10.
-        new_ser[new_ser < 0.] = 0.
-        new_ref = pd.Series(1., tmp[0]["ref"].resample("D").asfreq().index)
-        res = stat_processor_1("somo10", new_ser, new_ref, mtype, seasons,
-                               how="sum", minfrac=data_capture)
-    return res
+    return somo_processor(ser, ref, mtype, metadata, seasons, data_capture,
+                          10.)
 
 
 def somo35(ser, ref, mtype, metadata, seasons, data_capture):
@@ -333,17 +324,8 @@ def somo35(ser, ref, mtype, metadata, seasons, data_capture):
 
     :return: A list of dictionaries containing the processed data
     """
-    tmp = dma8eu(ser, ref, "daily", metadata, seasons, data_capture)
-    if mtype == "daily":
-        res = tmp
-        res[0]["name"] = "somo35"
-    else:
-        new_ser = tmp[0]["ser"] - 35.
-        new_ser[new_ser < 0.] = 0.
-        new_ref = pd.Series(1., tmp[0]["ref"].resample("D").asfreq().index)
-        res = stat_processor_1("somo35", new_ser, new_ref, mtype, seasons,
-                               how="sum", minfrac=data_capture)
-    return res
+    return somo_processor(ser, ref, mtype, metadata, seasons, data_capture,
+                          35.)
 
 
 def somo10_strict(ser, ref, mtype, metadata, seasons, data_capture):
@@ -365,17 +347,8 @@ def somo10_strict(ser, ref, mtype, metadata, seasons, data_capture):
 
     :return: A list of dictionaries containing the processed data
     """
-    tmp = dma8eu_strict(ser, ref, "daily", metadata, seasons, data_capture)
-    if mtype == "daily":
-        res = tmp
-        res[0]["name"] = "somo10_strict"
-    else:
-        new_ser = tmp[0]["ser"] - 10.
-        new_ser[new_ser < 0.] = 0.
-        new_ref = pd.Series(1., tmp[0]["ref"].resample("D").asfreq().index)
-        res = stat_processor_1("somo10_strict", new_ser, new_ref, mtype,
-                               seasons, how="sum", minfrac=data_capture)
-    return res
+    return somo_processor(ser, ref, mtype, metadata, seasons, data_capture,
+                          10., strict=True)
 
 
 def somo35_strict(ser, ref, mtype, metadata, seasons, data_capture):
@@ -397,17 +370,8 @@ def somo35_strict(ser, ref, mtype, metadata, seasons, data_capture):
 
     :return: A list of dictionaries containing the processed data
     """
-    tmp = dma8eu_strict(ser, ref, "daily", metadata, seasons, data_capture)
-    if mtype == "daily":
-        res = tmp
-        res[0]["name"] = "somo35_strict"
-    else:
-        new_ser = tmp[0]["ser"] - 35.
-        new_ser[new_ser < 0.] = 0.
-        new_ref = pd.Series(1., tmp[0]["ref"].resample("D").asfreq().index)
-        res = stat_processor_1("somo35_strict", new_ser, new_ref, mtype,
-                               seasons, how="sum", minfrac=data_capture)
-    return res
+    return somo_processor(ser, ref, mtype, metadata, seasons, data_capture,
+                          35., strict=True)
 
 
 def w90(ser, ref, mtype, metadata, seasons, data_capture):

toarstats/stats_utils.py

@@ -14,6 +14,7 @@ perc - calculate percentiles
 prepare_data - prepare the result list
 resample - resample the given series
 resample_with_date - resample the given series
+somo_processor - calculate the sum of excess of dma8
 stat_processor_1 - calculate statistics
 stat_processor_2 - calculate statistics
 threshold_processor - evaluate counts above the given threshold
@@ -409,6 +410,37 @@ def resample_with_date(ser, ref, sampling, how, mincount=0, minfrac=None,
     return dfres
 
 
+def somo_processor(ser, ref, mtype, metadata, seasons, data_capture, threshold,
+                   strict=False):
+    """Calculate the sum of excess of daily maximum 8-hour means.
+
+    :param ser: a series containing the values
+    :param ref: reference series
+    :param mtype: aggregation type
+    :param metadata: metadata needed for some statistics
+    :param seasons: season names
+    :param data_capture: a fractional coverage value
+    :param threshold: the cut-off value
+    :param strict: a boolean denoting whether to use the strict rules
+                   when to save a running mean
+
+    :return: A list of dictionaries containing the processed data
+    """
+    label = f"somo{int(threshold)}{'_strict' if strict else ''}"
+    tmp = dma8_processor(ser, ref, "daily", metadata, seasons, data_capture,
+                         dname="EU", strict=strict)
+    if mtype == "daily":
+        res = tmp
+        res[0]["name"] = label
+    else:
+        new_ser = tmp[0]["ser"] - threshold
+        new_ser[new_ser < 0.] = 0.
+        new_ref = pd.Series(1., tmp[0]["ref"].resample("D").asfreq().index)
+        res = stat_processor_1(label, new_ser, new_ref, mtype, seasons,
+                               how="sum", minfrac=data_capture)
+    return res
+
+
 def stat_processor_1(label, ser, ref, mtype, seasons, func=resample,
                      how="mean", mincount=0, minfrac=None, **kwargs):
     """Calculate statistics.