add benthic oxygen consumption

ecoval/__init__.py

@@ -587,6 +587,8 @@ def validate(title="Automated model evaluation", author=None, variables = "all",
                     Variable = "pCO2"
                 if vv == "benbio":
                     Variable = "macrobenthos biomass"
+                if vv == "oxycons":
+                    Variable = "oxygen consumption"
                 if vv == "susfrac":
                     Variable = "suspension feeding fraction"
                 if vv == "pft":

ecoval/data/point_template.ipynb

@@ -48,7 +48,7 @@
    "source": [
     "variable = \"point_variable\".lower()\n",
     "vv_name = variable\n",
-    "if vv_name in [\"benbio\", \"carbon\", \"susfrac\"]:\n",
+    "if vv_name in [\"benbio\", \"carbon\", \"susfrac\", \"oxycons\"]:\n",
     "    compact = True\n",
     "if vv_name.lower() == \"ph\":\n",
     "    vv_name = \"pH\"\n",
@@ -58,6 +58,8 @@
     "    vv_name = \"biomass of macrobenthos\"\n",
     "if vv_name == \"susfrac\":\n",
     "    vv_name = \"Suspension feeding fraction\"\n",
+    "if vv_name == \"oxycons\":\n",
+    "    vv_name = \"Benthic oxygen consumption\"\n",
     "layer = \"point_layer\"\n",
     "# get the units. File inspection could be randomized in case people have put loose files in there...\n",
     "import glob\n",
@@ -132,7 +134,7 @@
     "df[\"lat\"] = df[\"lat\"].apply(lambda x: bin_value(x, 0.5))\n",
     "if variable == \"benbio\":\n",
     "    df = df.assign(observation = lambda x: 1000 * 0.45 * x.observation) \n",
-    "if variable not in  [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in  [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    df = df.groupby([\"lon\", \"lat\", \"year\", \"month\"]).mean().reset_index()\n",
     "else:\n",
     "    df = df.groupby([\"lon\", \"lat\"]).mean().reset_index()"
@@ -382,7 +384,7 @@
     "  )+\n",
     "    theme(\n",
     "\n",
-    "    legend.position = \"bottom\", legend.direction = \"horizontal\", legend.box = \"horizontal\", legend.key.width = unit(6.0, \"cm\"),\n",
+    "    legend.position = \"bottom\", legend.direction = \"horizontal\", legend.box = \"horizontal\", legend.key.width = unit(3.0, \"cm\"),\n",
     "    legend.key.height = unit(1.0, \"cm\"))\n",
     "    # use ggtext to ensure things are superscripted\n",
     "    #theme(legend.title = element_markdown()) \n",
@@ -783,7 +785,7 @@
    },
    "outputs": [],
    "source": [
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    if compact is False:\n",
     "        if layer == \"surface\":\n",
     "            md(f\"**Figure {chapter}{i_figure}**: Simulated versus observed {vv_name} in the top 5 m of the water column. The blue curve is a generalized additive model (GAM) fit to the data, and the black line represents 1-1 relationship between the simulation and observations. The data has been binned to 0.5 degree resolution.\") \n",
@@ -890,7 +892,7 @@
    },
    "outputs": [],
    "source": [
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    df_bias = (\n",
     "        df_raw\n",
     "        .assign(bias = lambda x: x.model - x.observation)\n",
@@ -910,7 +912,7 @@
     "else:\n",
     "    # only want annual\n",
     "    df_bias = pd.DataFrame({\"month\": [\"All\"], \"bias\": [df_raw.model.mean() - df_raw.observation.mean()]})\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    # now create an rmse dataframe\n",
     "    df_rmse = (\n",
     "        df_raw\n",
@@ -933,7 +935,7 @@
     "df_table = copy.deepcopy(df_bias).merge(df_rmse)\n",
     "df_table = df_table.round(2)\n",
     "# create df_corr\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    df_corr = (\n",
     "        df_raw\n",
     "        .groupby(\"month\")\n",
@@ -960,7 +962,7 @@
     "# change Month to Period\n",
     "df_table = df_table.rename(columns={\"Month\": \"Time period\"})\n",
     "\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    # add commas to bias and rmse\n",
     "    df_number = df_raw.groupby(\"month\").count().reset_index().loc[:,[\"month\", \"observation\"]]\n",
     "# convert month number to name\n",
@@ -971,7 +973,7 @@
     "\n",
     "# add total number of observations\n",
     "annual_number = len(df_raw)\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    df_number = pd.concat([df_number, pd.DataFrame({\"Time period\": [\"All\"], \"Number of observations\": [annual_number]})])\n",
     "# df_number = df_number.append({\"Time period\": \"All\", \"Number of observations\": annual_number}, ignore_index=True)\n",
     "df_table = df_table.merge(df_number)\n",
@@ -1008,7 +1010,7 @@
    },
    "outputs": [],
    "source": [
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    df_bias = (\n",
     "        df\n",
     "        .assign(bias = lambda x: x.model - x.observation)\n",
@@ -1028,7 +1030,7 @@
     "else:\n",
     "    # only want annual\n",
     "    df_bias = pd.DataFrame({\"month\": [\"All\"], \"bias\": [df.model.mean() - df.observation.mean()]})\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    # now create an rmse dataframe\n",
     "    df_rmse = (\n",
     "        df\n",
@@ -1051,7 +1053,7 @@
     "df_table = copy.deepcopy(df_bias).merge(df_rmse)\n",
     "df_table = df_table.round(2)\n",
     "# create df_corr\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    df_corr = (\n",
     "        df\n",
     "        .groupby(\"month\")\n",
@@ -1078,7 +1080,7 @@
     "# change Month to Period\n",
     "df_table = df_table.rename(columns={\"Month\": \"Time period\"})\n",
     "\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    # add commas to bias and rmse\n",
     "    df_number = df.groupby(\"month\").count().reset_index().loc[:,[\"month\", \"observation\"]]\n",
     "# convert month number to name\n",
@@ -1089,7 +1091,7 @@
     "\n",
     "# add total number of observations\n",
     "annual_number = len(df)\n",
-    "if variable not in [\"carbon\", \"benbio\", \"susfrac\"]:\n",
+    "if variable not in [\"carbon\", \"benbio\", \"susfrac\", \"oxycons\"]:\n",
     "    df_number = pd.concat([df_number, pd.DataFrame({\"Time period\": [\"All\"], \"Number of observations\": [annual_number]})])\n",
     "# df_number = df_number.append({\"Time period\": \"All\", \"Number of observations\": annual_number}, ignore_index=True)\n",
     "df_table = df_table.merge(df_number)\n",
@@ -1319,7 +1321,10 @@
    "outputs": [],
    "source": [
     "if variable in [\"benbio\", \"susfrac\"]:\n",
-    "    md(\"URL: https://www.vliz.be/vmdcdata/nsbs/about.php\")"
+    "    md(\"URL: https://www.vliz.be/vmdcdata/nsbs/about.php\")\n",
+    "if variable == \"oxycons\":\n",
+    "    md(\"Stratmann, Tanja; Soetaert, Karline; Wei, Chih-Lin et al. (2022). Data from: The SCOC database – a large, open and global database with sediment community oxygen consumption rates [Dataset]. Dryad. https://doi.org/10.5061/dryad.25nd083\")\n",
+    "\n"
    ]
   },
   {

ecoval/matchall.py

@@ -38,7 +38,8 @@
     "carbon",
     "benbio",
     "benthic_carbon_flux",
-    "mesozoo"
+    "mesozoo",
+    "oxycons"
 
 ]
 
@@ -152,7 +153,11 @@ def mm_match(
 
                         ds.set_units({"total":"kg/m3"})
                     else:
-                        ds.sum_all()
+                        if variable == "oxycons":
+                            #Y2_fYG3c/12.011 + Y3_fYG3c/12.011 + Y4_fYG3c/12.011 + H1_fHG3c/12.011 + H2_fHG3c/12.011 + 2.0 * ben_nit_nrate
+                            ds.assign(total = lambda x: x.Ymacro_fYG3c_result/12.011 + x.Y4_fYG3c/12.011 + x.H1_fHG3c/12.011 + x.H2_fHG3c/12.011 + 2.0 * x.ben_nit_nrate, drop = True) 
+                        else:
+                            ds.sum_all()
 
             if len(df_locs) > 0:
                 if top_layer:
@@ -285,7 +290,7 @@ def extract_variable_mapping(folder, exclude=[], n_check=None, fvcom = False):
                         pass
         options = glob.glob(new_directory + "/**.nc")
         if True:
-            options = [x for x in options if "part" not in os.path.basename(x)]
+            #options = [x for x in options if "part" not in os.path.basename(x)]
             options = [x for x in options if "restart" not in os.path.basename(x)]
 
         if len([x for x in options if ".nc" in x]) > 0:
@@ -1288,7 +1293,7 @@ def write_report(x):
                 return None
 
     # check if benbio is in all_df variable column
-    if "benbio" in list(all_df.variable):
+    if "benbio" in point_benthic:
         # add another column using benbio row
         df_benbio = all_df.query("variable == 'benbio'").reset_index(drop=True)
         df_benbio["variable"] = ["susfrac"]
@@ -1617,6 +1622,8 @@ def point_match(variable, layer="all", ds_depths=None, df_times = None):
                                     "variable == @point_variable"
                                 ).model_variable
                             )[0]
+
+
                             if session_info["user_dir"]:
                                 paths = glob.glob(
                                     f"{obs_dir}/point/user/**/{variable}/**{variable}**.feather"
@@ -1629,6 +1636,7 @@ def point_match(variable, layer="all", ds_depths=None, df_times = None):
                                 paths = glob.glob(
                                     f"{obs_dir}/point/nws/**/{variable}/**{variable}**.feather"
                                 )
+
                             if variable == "pft":
                                 point_variable = "pft"
                             source = os.path.basename(paths[0]).split("_")[0]
@@ -1707,7 +1715,7 @@ def point_match(variable, layer="all", ds_depths=None, df_times = None):
                                 df = df.merge(df_include).reset_index(drop=True)
                             sel_these = point_time_res
                             sel_these = [x for x in df.columns if x in sel_these]
-                            if variable not in ["carbon", "benbio", "susfrac"]:
+                            if variable not in ["carbon", "benbio", "susfrac", "oxycons"]:
                                 paths = list(
                                     set(
                                         df.loc[:, sel_these]

ecoval/parsers.py

@@ -56,7 +56,8 @@ def generate_mapping(ds, fvcom = False):
         "nano",
         "pico",
         "benthic_carbon_flux",
-        "mesozoo"
+        "mesozoo",
+        "oxycons"
     ]
     if fvcom is False:
         ds1 = nc.open_data(ds[0], checks=False)
@@ -221,6 +222,14 @@ def generate_mapping(ds, fvcom = False):
                 for x in ds_contents.long_name
                 if "chloroph" in x.lower() and "pico" in x
             ]
+
+        if vv == "oxycons":
+            oxy_con_vars = list(set(["Ymacro_fYG3c_result", "Y4_fYG3c", "H1_fHG3c", "H2_fHG3c", "ben_nit_nrate"]))
+            the_vars = [x for x in ds_contents_top.variable if x in oxy_con_vars]
+            print(the_vars)
+            if len(the_vars) != len(oxy_con_vars):
+                the_vars = []
+
 
         if vv in ["carbon", "benbio", "benthic_carbon_flux"]:
             model_vars = ds_contents_top.query("long_name in @the_vars").variable
@@ -241,13 +250,16 @@ def generate_mapping(ds, fvcom = False):
                     ).variable
 
         add = True
+        if vv == "oxycons":
+            model_vars = the_vars
 
         if len(model_vars) > 1 and vv not in [
             "doc",
             "chlorophyll",
             "carbon",
             "benbio",
             "micro",
+            "oxycons"
         ]:
             add = False