Add Power et al. (2018) runup model

Makefile

@@ -142,6 +142,18 @@ impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv ./data/inte
 	--profile-type "prestorm" \
 	--output-file "./data/interim/twl_mean_slope_nie91.csv"
 
+./data/interim/twl_mean_slope_pow18.csv: ./data/interim/waves.csv ./data/interim/tides.csv ./data/interim/profiles.csv ./data/interim/sites.csv ./data/interim/profile_features_crest_toes.csv
+	$(PYTHON_CLI) create-twl-forecast \
+	--waves-csv "./data/interim/waves.csv" \
+	--tides-csv "./data/interim/tides.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--grain-size-csv "./data/interim/sites_grain_size.csv" \
+	--runup-function "pow18" \
+	--slope "mean" \
+	--profile-type "prestorm" \
+	--output-file "./data/interim/twl_mean_slope_pow18.csv"
+
 ### IMPACTS
 
 ./data/interim/impacts_observed.csv: ./data/interim/profiles.csv ./data/interim/profile_features_crest_toes.csv ./data/raw/profile_features_chris_leaman/profile_features_chris_leaman.csv
@@ -175,6 +187,12 @@ impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv ./data/inte
 	--forecasted-twl-csv "./data/interim/twl_mean_slope_nie91.csv" \
 	--output-file "./data/interim/impacts_forecasted_mean_slope_nie91.csv"
 
+./data/interim/impacts_forecasted_mean_slope_pow18.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_mean_slope_pow18.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_mean_slope_pow18.csv" \
+	--output-file "./data/interim/impacts_forecasted_mean_slope_pow18.csv"
+
 ### GEOJSONs
 
 geojsons: ./data/interim/impacts_forecasted_mean_slope_nie91.geojson ./data/interim/impacts_forecasted_mean_slope_nie91_R_high.geojson ./data/interim/impacts_forecasted_mean_slope_hol86.geojson ./data/interim/impacts_forecasted_mean_slope_hol86_R_high.geojson ./data/interim/impacts_forecasted_mean_slope_sto06.geojson ./data/interim/impacts_forecasted_mean_slope_sto06_R_high.geojson ./data/interim/profile_features_crest_toes.geojson ./data/interim/sites.geojson

src/analysis/forecast_twl.py

@@ -330,7 +330,7 @@ def slope_from_profile(
     "--runup-function",
     required=True,
     help="",
-    type=click.Choice(["sto06", "hol86", "nie91"]),
+    type=click.Choice(["sto06", "hol86", "nie91", "pow18"]),
 )
 @click.option(
     "--slope", required=True, help="", type=click.Choice(["foreshore", "mean"])

src/analysis/runup_models.py

@@ -1,6 +1,8 @@
 import numpy as np
 import pandas as pd
+from logs import setup_logging
 
+logger = setup_logging()
 
 
 def sto06(Hs0, Tp, beta, **kwargs):
@@ -104,8 +106,98 @@ def atk18(Hs0, Tp, beta):
     pass
 
 
-def pow18(Hs0, Tp, beta):
-    pass
+def pow18(Hs0, Tp, beta, r, **kwargs):
+    logger.info("Calculating runup using Power et al. (2018)")
+    df = pd.DataFrame(
+        {"Hs0": Hs0, "Tp": Tp, "beta": beta, "r": r},
+        index=[x for x in range(0, np.size(Hs0))],
+    )
+
+    df["Lp"] = 9.8 * df["Tp"] ** 2 / 2 / np.pi
+    df["x1"] = df["Hs0"] / df["Lp"]
+    df["x2"] = df["beta"]
+    df["x3"] = df["r"] / df["Hs0"]
+
+    df["R2"] = df.Hs0 * (
+        (df.x2 + (((df.x3 * 3) / np.exp(-5)) * ((3 * df.x3) * df.x3)))
+        + ((((df.x1 + df.x3) - 2) - (df.x3 - df.x2)) + ((df.x2 - df.x1) - df.x3))
+        + (((df.x3 ** df.x1) - (df.x3 ** (1 / 3))) - (np.exp(df.x2) ** (df.x1 * 3)))
+        + np.sqrt((((df.x3 + df.x1) - df.x2) - (df.x2 + np.log10(df.x3))))
+        + ((((df.x2 ** 2) / (df.x1 ** (1 / 3))) ** (df.x1 ** (1 / 3))) - np.sqrt(df.x3))
+        + (
+            (df.x2 + ((df.x3 / df.x1) ** (1 / 3)))
+            + (np.log(2) - (1 / (1 + np.exp(-(df.x2 + df.x3)))))
+        )
+        + ((np.sqrt(df.x3) - (((3 ** 2) + 3) * (df.x2 ** 2))) ** 2)
+        + ((((df.x3 * -5) ** 2) ** 2) + (((df.x3 + df.x3) * df.x1) / (df.x2 ** 2)))
+        + np.log(
+            (np.sqrt(((df.x2 ** 2) + (df.x3 ** (1 / 3)))) + ((df.x2 + 3) ** (1 / 3)))
+        )
+        + (
+            (((df.x1 / df.x3) * (-5 ** 2)) * (df.x3 ** 2))
+            - np.log10((1 / (1 + np.exp(-(df.x2 + df.x3)))))
+        )
+        + (df.x1 ** df.x3)
+        + np.exp(-((((df.x3 / df.x1) ** np.exp(4)) + (np.exp(df.x3) ** 3)) ** 2))
+        + np.exp((np.log((df.x2 - df.x3)) - np.log(np.exp(-((-1 + df.x1) ** 2)))))
+        + ((np.sqrt(4) * (((df.x3 / df.x2) - df.x2) - (0 - df.x1))) ** 2)
+        + (2 * ((((-5 * df.x3) + df.x1) * (2 - df.x3)) - 2))
+        + ((np.sqrt(4) * (((df.x3 / df.x2) - df.x2) - (0 - df.x1))) ** 2)
+        + ((((-5 + df.x1) - df.x2) * (df.x2 - df.x3)) * ((df.x1 - df.x2) - (-4 ** -5)))
+        + (np.exp(-((df.x2 + (-5 - df.x1)) ** 2)) + ((df.x2 + 5) * (df.x3 ** 2)))
+        + np.sqrt(
+            1
+            / (
+                1
+                + np.exp(
+                    -(
+                        (np.exp(df.x1) - np.exp(-((df.x3 + df.x3) ** 2)))
+                        + ((df.x1 ** df.x3) - (df.x3 * 4))
+                    )
+                )
+            )
+        )
+        + (
+            (
+                np.exp(
+                    -(
+                        (
+                            (
+                                (
+                                    np.exp(
+                                        -(
+                                            (
+                                                (np.sqrt(df.x3) * 4)
+                                                + (1 / (1 + np.exp(-(df.x2 + 2))))
+                                            )
+                                            ** 2
+                                        )
+                                    )
+                                )
+                                ** 2
+                            )
+                            + df.x1
+                        )
+                        ** 2
+                    )
+                )
+            )
+            ** 3
+        )
+    )
+
+    df["setup"] = np.nan
+    df["S_ig"] = np.nan
+    df["S_inc"] = np.nan
+    df["S_total"] = np.nan
+
+    return (
+        float_or_list(df["R2"].tolist()),
+        float_or_list(df["setup"].tolist()),
+        float_or_list(df["S_total"].tolist()),
+        float_or_list(df["S_inc"].tolist()),
+        float_or_list(df["S_ig"].tolist()),
+    )
 
 
 def beu(Hs0, Tp, beta):