Optimize calculation of profile volume difference

Rather than end points

Makefile

@@ -53,7 +53,7 @@ pull-data: ##@data Copies data from data backup directory to ./data/
 # Command for activating our virtual environment and calling the CLI entry point
 PYTHON_CLI = activate ./.venv && python ./src/cli.py
 
-impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv ./data/interim/impacts_forecasted_mean_slope_sto06.csv ./data/interim/impacts_observed.csv ##@products makes obsered and forecasted impacts
+impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv ./data/interim/impacts_forecasted_premean_slope_sto06.csv ./data/interim/impacts_observed.csv ##@products makes obsered and forecasted impacts
 
 
 ### Parses raw matfiles
@@ -67,7 +67,7 @@ impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv ./data/inte
 	--profiles-output-file "./data/interim/profiles.csv" \
 	--sites-output-file "./data/interim/sites.csv"
 
-./data/interim/waves.csv: ./data/interim/sites.csv ./data/raw/processed_shorelines/waves.mat ./data/raw/processed_shorelines/sites_waves.mat
+./data/interim/waves.csv: ./data/interim/sites.csv ./data/raw/processed_shorelines/waves.mat ./data/raw/processed_shorelines/waves.mat
 	$(PYTHON_CLI) create-waves-csv \
 	--waves-mat "./data/raw/processed_shorelines/waves.mat" \
 	--sites-csv "./data/interim/sites.csv" \
@@ -94,52 +94,251 @@ impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv ./data/inte
 
 ### TWLs
 
+twls: ./data/interim/twl_foreshore_slope_sto06.csv
+twls: ./data/interim/twl_premean_slope_sto06.csv
+twls: ./data/interim/twl_premean_slope_hol86.csv
+twls: ./data/interim/twl_premean_slope_nie91.csv
+twls: ./data/interim/twl_premean_slope_pow18.csv
+twls: ./data/interim/twl_postmean_slope_sto06.csv
+twls: ./data/interim/twl_postmean_slope_hol86.csv
+twls: ./data/interim/twl_postmean_slope_nie91.csv
+twls: ./data/interim/twl_postmean_slope_pow18.csv
+twls: ./data/interim/twl_preintertidal_slope_sto06.csv
+twls: ./data/interim/twl_preintertidal_slope_hol86.csv
+twls: ./data/interim/twl_preintertidal_slope_nie91.csv
+twls: ./data/interim/twl_preintertidal_slope_pow18.csv
+twls: ./data/interim/twl_postintertidal_slope_sto06.csv
+twls: ./data/interim/twl_postintertidal_slope_hol86.csv
+twls: ./data/interim/twl_postintertidal_slope_nie91.csv
+twls: ./data/interim/twl_postintertidal_slope_pow18.csv
+
 ./data/interim/twl_foreshore_slope_sto06.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 "sto06" \
 	--slope "foreshore" \
 	--profile-type "prestorm" \
 	--output-file "./data/interim/twl_foreshore_slope_sto06.csv"
 
-./data/interim/twl_mean_slope_sto06.csv: ./data/interim/waves.csv ./data/interim/tides.csv ./data/interim/profiles.csv ./data/interim/sites.csv ./data/interim/profile_features_crest_toes.csv
+./data/interim/twl_premean_slope_sto06.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 "sto06" \
 	--slope "mean" \
 	--profile-type "prestorm" \
-	--output-file "./data/interim/twl_mean_slope_sto06.csv"
+	--output-file "./data/interim/twl_premean_slope_sto06.csv"
 
-./data/interim/twl_mean_slope_hol86.csv: ./data/interim/waves.csv ./data/interim/tides.csv ./data/interim/profiles.csv ./data/interim/sites.csv ./data/interim/profile_features_crest_toes.csv
+./data/interim/twl_premean_slope_hol86.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 "hol86" \
 	--slope "mean" \
 	--profile-type "prestorm" \
-	--output-file "./data/interim/twl_mean_slope_hol86.csv"
+	--output-file "./data/interim/twl_premean_slope_hol86.csv"
+
+./data/interim/twl_premean_slope_nie91.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 "nie91" \
+	--slope "mean" \
+	--profile-type "prestorm" \
+	--output-file "./data/interim/twl_premean_slope_nie91.csv"
+
+./data/interim/twl_premean_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_premean_slope_pow18.csv"
+
+./data/interim/twl_postmean_slope_sto06.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 "sto06" \
+	--slope "mean" \
+	--profile-type "poststorm" \
+	--output-file "./data/interim/twl_postmean_slope_sto06.csv"
+
+./data/interim/twl_postmean_slope_hol86.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 "hol86" \
+	--slope "mean" \
+	--profile-type "poststorm" \
+	--output-file "./data/interim/twl_postmean_slope_hol86.csv"
 
-./data/interim/twl_mean_slope_nie91.csv: ./data/interim/waves.csv ./data/interim/tides.csv ./data/interim/profiles.csv ./data/interim/sites.csv ./data/interim/profile_features_crest_toes.csv
+./data/interim/twl_postmean_slope_nie91.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 "nie91" \
 	--slope "mean" \
+	--profile-type "poststorm" \
+	--output-file "./data/interim/twl_postmean_slope_nie91.csv"
+
+./data/interim/twl_postmean_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 "poststorm" \
+	--output-file "./data/interim/twl_postmean_slope_pow18.csv"
+
+./data/interim/twl_preintertidal_slope_sto06.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 "sto06" \
+	--slope "intertidal" \
 	--profile-type "prestorm" \
-	--output-file "./data/interim/twl_mean_slope_nie91.csv"
+	--output-file "./data/interim/twl_preintertidal_slope_sto06.csv"
+
+./data/interim/twl_preintertidal_slope_hol86.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 "hol86" \
+	--slope "intertidal" \
+	--profile-type "prestorm" \
+	--output-file "./data/interim/twl_preintertidal_slope_hol86.csv"
+
+./data/interim/twl_preintertidal_slope_nie91.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 "nie91" \
+	--slope "intertidal" \
+	--profile-type "prestorm" \
+	--output-file "./data/interim/twl_preintertidal_slope_nie91.csv"
+
+./data/interim/twl_preintertidal_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 "intertidal" \
+	--profile-type "prestorm" \
+	--output-file "./data/interim/twl_preintertidal_slope_pow18.csv"
+
+./data/interim/twl_postintertidal_slope_sto06.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 "sto06" \
+	--slope "intertidal" \
+	--profile-type "poststorm" \
+	--output-file "./data/interim/twl_postintertidal_slope_sto06.csv"
+
+./data/interim/twl_postintertidal_slope_hol86.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 "hol86" \
+	--slope "intertidal" \
+	--profile-type "poststorm" \
+	--output-file "./data/interim/twl_postintertidal_slope_hol86.csv"
+
+./data/interim/twl_postintertidal_slope_nie91.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 "nie91" \
+	--slope "intertidal" \
+	--profile-type "poststorm" \
+	--output-file "./data/interim/twl_postintertidal_slope_nie91.csv"
+
+./data/interim/twl_postintertidal_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 "intertidal" \
+	--profile-type "poststorm" \
+	--output-file "./data/interim/twl_postintertidal_slope_pow18.csv"
+
 
 ### IMPACTS
 
+impacts: ./data/interim/impacts_observed.csv
+impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv
+impacts: ./data/interim/impacts_forecasted_premean_slope_sto06.csv
+impacts: ./data/interim/impacts_forecasted_premean_slope_hol86.csv
+impacts: ./data/interim/impacts_forecasted_premean_slope_nie91.csv
+impacts: ./data/interim/impacts_forecasted_premean_slope_pow18.csv
+impacts: ./data/interim/impacts_forecasted_postmean_slope_sto06.csv
+impacts: ./data/interim/impacts_forecasted_postmean_slope_hol86.csv
+impacts: ./data/interim/impacts_forecasted_postmean_slope_nie91.csv
+impacts: ./data/interim/impacts_forecasted_postmean_slope_pow18.csv
+impacts: ./data/interim/impacts_forecasted_preintertidal_slope_sto06.csv
+impacts: ./data/interim/impacts_forecasted_preintertidal_slope_hol86.csv
+impacts: ./data/interim/impacts_forecasted_preintertidal_slope_nie91.csv
+impacts: ./data/interim/impacts_forecasted_preintertidal_slope_pow18.csv
+impacts: ./data/interim/impacts_forecasted_postintertidal_slope_sto06.csv
+impacts: ./data/interim/impacts_forecasted_postintertidal_slope_hol86.csv
+impacts: ./data/interim/impacts_forecasted_postintertidal_slope_nie91.csv
+impacts: ./data/interim/impacts_forecasted_postintertidal_slope_pow18.csv
+
+
 ./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
 	$(PYTHON_CLI) create-observed-impacts \
 	--profiles-csv "./data/interim/profiles.csv" \
@@ -147,79 +346,368 @@ impacts: ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv ./data/inte
 	--raw-profile-features-csv "./data/raw/profile_features_chris_leaman/profile_features_chris_leaman.csv" \
 	--output-file "./data/interim/impacts_observed.csv"
 
-./data/interim/impacts_forecasted_mean_slope_sto06.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_mean_slope_sto06.csv
-	$(PYTHON_CLI) create-forecasted-impacts \
-	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
-	--forecasted-twl-csv "./data/interim/twl_mean_slope_sto06.csv" \
-	--output-file "./data/interim/impacts_forecasted_mean_slope_sto06.csv"
-
 ./data/interim/impacts_forecasted_foreshore_slope_sto06.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_foreshore_slope_sto06.csv
 	$(PYTHON_CLI) create-forecasted-impacts \
 	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
 	--forecasted-twl-csv "./data/interim/twl_foreshore_slope_sto06.csv" \
 	--output-file "./data/interim/impacts_forecasted_foreshore_slope_sto06.csv"
 
-./data/interim/impacts_forecasted_mean_slope_hol86.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_mean_slope_hol86.csv
+./data/interim/impacts_forecasted_premean_slope_sto06.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_premean_slope_sto06.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_premean_slope_sto06.csv" \
+	--output-file "./data/interim/impacts_forecasted_premean_slope_sto06.csv"
+
+./data/interim/impacts_forecasted_premean_slope_hol86.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_premean_slope_hol86.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_premean_slope_hol86.csv" \
+	--output-file "./data/interim/impacts_forecasted_premean_slope_hol86.csv"
+
+./data/interim/impacts_forecasted_premean_slope_nie91.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_premean_slope_nie91.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_premean_slope_nie91.csv" \
+	--output-file "./data/interim/impacts_forecasted_premean_slope_nie91.csv"
+
+./data/interim/impacts_forecasted_premean_slope_pow18.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_premean_slope_pow18.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_premean_slope_pow18.csv" \
+	--output-file "./data/interim/impacts_forecasted_premean_slope_pow18.csv"
+
+./data/interim/impacts_forecasted_postmean_slope_sto06.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postmean_slope_sto06.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postmean_slope_sto06.csv" \
+	--output-file "./data/interim/impacts_forecasted_postmean_slope_sto06.csv"
+
+./data/interim/impacts_forecasted_postmean_slope_hol86.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postmean_slope_hol86.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postmean_slope_hol86.csv" \
+	--output-file "./data/interim/impacts_forecasted_postmean_slope_hol86.csv"
+
+./data/interim/impacts_forecasted_postmean_slope_nie91.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postmean_slope_nie91.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postmean_slope_nie91.csv" \
+	--output-file "./data/interim/impacts_forecasted_postmean_slope_nie91.csv"
+
+./data/interim/impacts_forecasted_postmean_slope_pow18.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postmean_slope_pow18.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postmean_slope_pow18.csv" \
+	--output-file "./data/interim/impacts_forecasted_postmean_slope_pow18.csv"
+	
+./data/interim/impacts_forecasted_preintertidal_slope_sto06.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_preintertidal_slope_sto06.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_preintertidal_slope_sto06.csv" \
+	--output-file "./data/interim/impacts_forecasted_preintertidal_slope_sto06.csv"
+
+./data/interim/impacts_forecasted_preintertidal_slope_hol86.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_preintertidal_slope_hol86.csv
 	$(PYTHON_CLI) create-forecasted-impacts \
 	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
-	--forecasted-twl-csv "./data/interim/twl_mean_slope_hol86.csv" \
+	--forecasted-twl-csv "./data/interim/twl_preintertidal_slope_hol86.csv" \
-	--output-file "./data/interim/impacts_forecasted_mean_slope_hol86.csv"
+	--output-file "./data/interim/impacts_forecasted_preintertidal_slope_hol86.csv"
 
-./data/interim/impacts_forecasted_mean_slope_nie91.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_mean_slope_nie91.csv
+./data/interim/impacts_forecasted_preintertidal_slope_nie91.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_preintertidal_slope_nie91.csv
 	$(PYTHON_CLI) create-forecasted-impacts \
 	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
-	--forecasted-twl-csv "./data/interim/twl_mean_slope_nie91.csv" \
+	--forecasted-twl-csv "./data/interim/twl_preintertidal_slope_nie91.csv" \
-	--output-file "./data/interim/impacts_forecasted_mean_slope_nie91.csv"
+	--output-file "./data/interim/impacts_forecasted_preintertidal_slope_nie91.csv"
+
+./data/interim/impacts_forecasted_preintertidal_slope_pow18.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_preintertidal_slope_pow18.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_preintertidal_slope_pow18.csv" \
+	--output-file "./data/interim/impacts_forecasted_preintertidal_slope_pow18.csv"
+
+./data/interim/impacts_forecasted_postintertidal_slope_sto06.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postintertidal_slope_sto06.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postintertidal_slope_sto06.csv" \
+	--output-file "./data/interim/impacts_forecasted_postintertidal_slope_sto06.csv"
+
+./data/interim/impacts_forecasted_postintertidal_slope_hol86.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postintertidal_slope_hol86.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postintertidal_slope_hol86.csv" \
+	--output-file "./data/interim/impacts_forecasted_postintertidal_slope_hol86.csv"
+
+./data/interim/impacts_forecasted_postintertidal_slope_nie91.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postintertidal_slope_nie91.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postintertidal_slope_nie91.csv" \
+	--output-file "./data/interim/impacts_forecasted_postintertidal_slope_nie91.csv"
+
+./data/interim/impacts_forecasted_postintertidal_slope_pow18.csv: ./data/interim/profile_features_crest_toes.csv ./data/interim/twl_postintertidal_slope_pow18.csv
+	$(PYTHON_CLI) create-forecasted-impacts \
+	--profile-features-csv "./data/interim/profile_features_crest_toes.csv" \
+	--forecasted-twl-csv "./data/interim/twl_postintertidal_slope_pow18.csv" \
+	--output-file "./data/interim/impacts_forecasted_postintertidal_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
+geojsons: ./data/interim/impacts_forecasted_premean_slope_pow18.geojson ./data/interim/impacts_forecasted_premean_slope_pow18_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_premean_slope_nie91.geojson ./data/interim/impacts_forecasted_premean_slope_nie91_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_premean_slope_hol86.geojson ./data/interim/impacts_forecasted_premean_slope_hol86_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_premean_slope_sto06.geojson ./data/interim/impacts_forecasted_premean_slope_sto06_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postmean_slope_pow18.geojson ./data/interim/impacts_forecasted_postmean_slope_pow18_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postmean_slope_nie91.geojson ./data/interim/impacts_forecasted_postmean_slope_nie91_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postmean_slope_hol86.geojson ./data/interim/impacts_forecasted_postmean_slope_hol86_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postmean_slope_sto06.geojson ./data/interim/impacts_forecasted_postmean_slope_sto06_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_preintertidal_slope_pow18.geojson ./data/interim/impacts_forecasted_preintertidal_slope_pow18_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_preintertidal_slope_nie91.geojson ./data/interim/impacts_forecasted_preintertidal_slope_nie91_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_preintertidal_slope_hol86.geojson ./data/interim/impacts_forecasted_preintertidal_slope_hol86_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_preintertidal_slope_sto06.geojson ./data/interim/impacts_forecasted_preintertidal_slope_sto06_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postintertidal_slope_pow18.geojson ./data/interim/impacts_forecasted_postintertidal_slope_pow18_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postintertidal_slope_nie91.geojson ./data/interim/impacts_forecasted_postintertidal_slope_nie91_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postintertidal_slope_hol86.geojson ./data/interim/impacts_forecasted_postintertidal_slope_hol86_R_high.geojson
+geojsons: ./data/interim/impacts_forecasted_postintertidal_slope_sto06.geojson ./data/interim/impacts_forecasted_postintertidal_slope_sto06_R_high.geojson
+geojsons: ./data/interim/profile_features_crest_toes.geojson ./data/interim/sites.geojson
+
+./data/interim/impacts_forecasted_premean_slope_sto06.geojson: ./data/interim/impacts_forecasted_premean_slope_sto06.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_premean_slope_sto06.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_sto06.geojson"
+
+./data/interim/impacts_forecasted_premean_slope_sto06_R_high.geojson: ./data/interim/impacts_forecasted_premean_slope_sto06.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_premean_slope_sto06.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_sto06_R_high.geojson"
+
+./data/interim/impacts_forecasted_premean_slope_hol86.geojson: ./data/interim/impacts_forecasted_premean_slope_hol86.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_premean_slope_hol86.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_hol86.geojson"
+
+./data/interim/impacts_forecasted_premean_slope_hol86_R_high.geojson: ./data/interim/impacts_forecasted_premean_slope_hol86.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_premean_slope_hol86.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_hol86_R_high.geojson"
+
+./data/interim/impacts_forecasted_premean_slope_nie91.geojson: ./data/interim/impacts_forecasted_premean_slope_nie91.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_premean_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_nie91.geojson"
+
+./data/interim/impacts_forecasted_premean_slope_nie91_R_high.geojson: ./data/interim/impacts_forecasted_premean_slope_nie91.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_premean_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_nie91_R_high.geojson"
+
+./data/interim/impacts_forecasted_premean_slope_pow18.geojson: ./data/interim/impacts_forecasted_premean_slope_pow18.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_premean_slope_pow18.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_pow18.geojson"
+
+./data/interim/impacts_forecasted_premean_slope_pow18_R_high.geojson: ./data/interim/impacts_forecasted_premean_slope_pow18.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_premean_slope_pow18.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_premean_slope_pow18_R_high.geojson"
+
+./data/interim/impacts_forecasted_postmean_slope_sto06.geojson: ./data/interim/impacts_forecasted_postmean_slope_sto06.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postmean_slope_sto06.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_sto06.geojson"
+
+./data/interim/impacts_forecasted_postmean_slope_sto06_R_high.geojson: ./data/interim/impacts_forecasted_postmean_slope_sto06.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postmean_slope_sto06.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_sto06_R_high.geojson"
+
+./data/interim/impacts_forecasted_postmean_slope_hol86.geojson: ./data/interim/impacts_forecasted_postmean_slope_hol86.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postmean_slope_hol86.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_hol86.geojson"
+
+./data/interim/impacts_forecasted_postmean_slope_hol86_R_high.geojson: ./data/interim/impacts_forecasted_postmean_slope_hol86.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postmean_slope_hol86.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_hol86_R_high.geojson"
+
+./data/interim/impacts_forecasted_postmean_slope_nie91.geojson: ./data/interim/impacts_forecasted_postmean_slope_nie91.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postmean_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_nie91.geojson"
 
-./data/interim/impacts_forecasted_mean_slope_sto06.geojson: ./data/interim/impacts_forecasted_mean_slope_sto06.csv ./data/interim/impacts_observed.csv
+./data/interim/impacts_forecasted_postmean_slope_nie91_R_high.geojson: ./data/interim/impacts_forecasted_postmean_slope_nie91.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postmean_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_nie91_R_high.geojson"
+
+./data/interim/impacts_forecasted_postmean_slope_pow18.geojson: ./data/interim/impacts_forecasted_postmean_slope_pow18.csv ./data/interim/impacts_observed.csv
 	$(PYTHON_CLI) impacts-to-geojson \
 	--sites-csv "./data/interim/sites.csv" \
 	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
-	--forecast-impacts-csv "./data/interim/impacts_forecasted_mean_slope_sto06.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postmean_slope_pow18.csv" \
-	--output-geojson "./data/interim/impacts_forecasted_mean_slope_sto06.geojson"
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_pow18.geojson"
 
-./data/interim/impacts_forecasted_mean_slope_sto06_R_high.geojson: ./data/interim/impacts_forecasted_mean_slope_sto06.csv
+./data/interim/impacts_forecasted_postmean_slope_pow18_R_high.geojson: ./data/interim/impacts_forecasted_postmean_slope_pow18.csv
 	$(PYTHON_CLI) r-high-to-geojson \
 	--sites-csv "./data/interim/sites.csv" \
 	--profiles-csv "./data/interim/profiles.csv" \
 	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
-	--impacts-csv "./data/interim/impacts_forecasted_mean_slope_sto06.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postmean_slope_pow18.csv" \
-	--output-geojson "./data/interim/impacts_forecasted_mean_slope_sto06_R_high.geojson"
+	--output-geojson "./data/interim/impacts_forecasted_postmean_slope_pow18_R_high.geojson"
 
-./data/interim/impacts_forecasted_mean_slope_hol86.geojson: ./data/interim/impacts_forecasted_mean_slope_hol86.csv ./data/interim/impacts_observed.csv
+./data/interim/impacts_forecasted_preintertidal_slope_sto06.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_sto06.csv ./data/interim/impacts_observed.csv
 	$(PYTHON_CLI) impacts-to-geojson \
 	--sites-csv "./data/interim/sites.csv" \
 	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
-	--forecast-impacts-csv "./data/interim/impacts_forecasted_mean_slope_hol86.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_sto06.csv" \
-	--output-geojson "./data/interim/impacts_forecasted_mean_slope_hol86.geojson"
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_sto06.geojson"
 
-./data/interim/impacts_forecasted_mean_slope_hol86_R_high.geojson: ./data/interim/impacts_forecasted_mean_slope_hol86.csv
+./data/interim/impacts_forecasted_preintertidal_slope_sto06_R_high.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_sto06.csv
 	$(PYTHON_CLI) r-high-to-geojson \
 	--sites-csv "./data/interim/sites.csv" \
 	--profiles-csv "./data/interim/profiles.csv" \
 	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
-	--impacts-csv "./data/interim/impacts_forecasted_mean_slope_hol86.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_sto06.csv" \
-	--output-geojson "./data/interim/impacts_forecasted_mean_slope_hol86_R_high.geojson"
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_sto06_R_high.geojson"
 
-./data/interim/impacts_forecasted_mean_slope_nie91.geojson: ./data/interim/impacts_forecasted_mean_slope_nie91.csv ./data/interim/impacts_observed.csv
+./data/interim/impacts_forecasted_preintertidal_slope_hol86.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_hol86.csv ./data/interim/impacts_observed.csv
 	$(PYTHON_CLI) impacts-to-geojson \
 	--sites-csv "./data/interim/sites.csv" \
 	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
-	--forecast-impacts-csv "./data/interim/impacts_forecasted_mean_slope_nie91.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_hol86.csv" \
-	--output-geojson "./data/interim/impacts_forecasted_mean_slope_nie91.geojson"
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_hol86.geojson"
 
-./data/interim/impacts_forecasted_mean_slope_nie91_R_high.geojson: ./data/interim/impacts_forecasted_mean_slope_nie91.csv
+./data/interim/impacts_forecasted_preintertidal_slope_hol86_R_high.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_hol86.csv
 	$(PYTHON_CLI) r-high-to-geojson \
 	--sites-csv "./data/interim/sites.csv" \
 	--profiles-csv "./data/interim/profiles.csv" \
 	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
-	--impacts-csv "./data/interim/impacts_forecasted_mean_slope_nie91.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_hol86.csv" \
-	--output-geojson "./data/interim/impacts_forecasted_mean_slope_nie91_R_high.geojson"
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_hol86_R_high.geojson"
+
+./data/interim/impacts_forecasted_preintertidal_slope_nie91.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_nie91.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_nie91.geojson"
 
+./data/interim/impacts_forecasted_preintertidal_slope_nie91_R_high.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_nie91.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_nie91_R_high.geojson"
+
+./data/interim/impacts_forecasted_preintertidal_slope_pow18.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_pow18.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_pow18.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_pow18.geojson"
+
+./data/interim/impacts_forecasted_preintertidal_slope_pow18_R_high.geojson: ./data/interim/impacts_forecasted_preintertidal_slope_pow18.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_preintertidal_slope_pow18.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_preintertidal_slope_pow18_R_high.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_sto06.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_sto06.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_sto06.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_sto06.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_sto06_R_high.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_sto06.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_sto06.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_sto06_R_high.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_hol86.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_hol86.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_hol86.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_hol86.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_hol86_R_high.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_hol86.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_hol86.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_hol86_R_high.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_nie91.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_nie91.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_nie91.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_nie91_R_high.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_nie91.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_nie91.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_nie91_R_high.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_pow18.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_pow18.csv ./data/interim/impacts_observed.csv
+	$(PYTHON_CLI) impacts-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--observed-impacts-csv "./data/interim/impacts_observed.csv" \
+	--forecast-impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_pow18.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_pow18.geojson"
+
+./data/interim/impacts_forecasted_postintertidal_slope_pow18_R_high.geojson: ./data/interim/impacts_forecasted_postintertidal_slope_pow18.csv
+	$(PYTHON_CLI) r-high-to-geojson \
+	--sites-csv "./data/interim/sites.csv" \
+	--profiles-csv "./data/interim/profiles.csv" \
+	--crest-toes-csv "./data/interim/profile_features_crest_toes.csv" \
+	--impacts-csv "./data/interim/impacts_forecasted_postintertidal_slope_pow18.csv" \
+	--output-geojson "./data/interim/impacts_forecasted_postintertidal_slope_pow18_R_high.geojson"
+# Site specific properties
 ./data/interim/profile_features_crest_toes.geojson: ./data/interim/profile_features_crest_toes.csv
 	$(PYTHON_CLI) profile-features-crest-toes-to-geojson \
 	--sites-csv "./data/interim/sites.csv" \

README.md

@@ -78,3 +78,6 @@ been corrected for systematic errors, so actual elevations should be taken from
 - [ ] Investigate using [modin](https://github.com/modin-project/modin) to help speed up analysis.
 - [X] Need to think about how relative imports are handled, see [here](https://chrisyeh96.github.io/2017/08/08/definitive-guide-python-imports.html). Maybe the click CLI interface should be moved to the `./src/` folder and it can import all the other packages?
 - [ ] Simplify runup_models in Stockdon06 - we should really only have one function for each runup model. Need to make it work with individual values or entire dataframe. Use [np.maskedarray](https://docs.scipy.org/doc/numpy-1.15.0/reference/maskedarray.generic.html)
+
+## Misc
+If Qt5Agg error comes up when trying to plot in PyCharm: https://stackoverflow.com/a/42231526

notebooks/01_exploration.ipynb

@@ -90,8 +90,7 @@
     "# Note that the forecasted data sets should be in the same order for impacts and twls\n",
     "impacts = {\n",
     "    'forecasted': {\n",
-    "        'foreshore_slope_sto06': df_from_csv('impacts_forecasted_foreshore_slope_sto06.csv', index_col=[0]),\n",
+    "        'postintertidal_slope_sto06': df_from_csv('impacts_forecasted_postintertidal_slope_sto06.csv', index_col=[0]),\n",
-    "        'mean_slope_sto06': df_from_csv('impacts_forecasted_mean_slope_sto06.csv', index_col=[0]),\n",
     "        },\n",
     "    'observed': df_from_csv('impacts_observed.csv', index_col=[0])\n",
     "    }\n",
@@ -99,8 +98,7 @@
     "\n",
     "twls = {\n",
     "    'forecasted': {\n",
-    "        'foreshore_slope_sto06': df_from_csv('twl_foreshore_slope_sto06.csv', index_col=[0, 1]),\n",
+    "        'postintertidal_slope_sto06':df_from_csv('twl_postintertidal_slope_sto06.csv', index_col=[0, 1]),\n",
-    "        'mean_slope_sto06':df_from_csv('twl_mean_slope_sto06.csv', index_col=[0, 1]),\n",
     "    }\n",
     "}\n",
     "print('Done!')"
@@ -260,6 +258,7 @@
     "    title='Bed Profiles',\n",
     "    height=300,\n",
     "    legend=dict(font={'size': 10}),\n",
+    "    showlegend=False,\n",
     "    margin=dict(t=50, b=50, l=50, r=20),\n",
     "    xaxis=dict(\n",
     "        title='x (m)',\n",

notebooks/07_evaluate_model_performance.ipynb

@@ -0,0 +1,592 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Evaluate prediction metrics \n",
+    "- This notebook is used to check out each of our storm impact prediction models performed in comparison to our observed storm impacts."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Setup notebook\n",
+    "Import our required packages and set default plotting options."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Enable autoreloading of our modules. \n",
+    "# Most of the code will be located in the /src/ folder, \n",
+    "# and then called from the notebook.\n",
+    "%matplotlib inline\n",
+    "%reload_ext autoreload\n",
+    "%autoreload"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from IPython.core.debugger import set_trace\n",
+    "\n",
+    "import pandas as pd\n",
+    "import numpy as np\n",
+    "import os\n",
+    "import decimal\n",
+    "import plotly\n",
+    "import plotly.graph_objs as go\n",
+    "import plotly.plotly as py\n",
+    "import plotly.tools as tls\n",
+    "import plotly.figure_factory as ff\n",
+    "from plotly import tools\n",
+    "import plotly.io as pio\n",
+    "from scipy import stats\n",
+    "import math\n",
+    "import matplotlib\n",
+    "from matplotlib import cm\n",
+    "import colorlover as cl\n",
+    "from tqdm import tqdm_notebook\n",
+    "from ipywidgets import widgets, Output\n",
+    "from IPython.display import display, clear_output, Image, HTML\n",
+    "from scipy import stats\n",
+    "from sklearn.metrics import confusion_matrix\n",
+    "import matplotlib.pyplot as plt\n",
+    "from matplotlib.ticker import MultipleLocator\n",
+    "from matplotlib.lines import Line2D\n",
+    "from cycler import cycler\n",
+    "from scipy.interpolate import interp1d\n",
+    "from pandas.api.types import CategoricalDtype"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Matplot lib default settings\n",
+    "plt.rcParams[\"figure.figsize\"] = (10,6)\n",
+    "plt.rcParams['axes.grid']=True\n",
+    "plt.rcParams['grid.alpha'] = 0.5\n",
+    "plt.rcParams['grid.color'] = \"grey\"\n",
+    "plt.rcParams['grid.linestyle'] = \"--\"\n",
+    "plt.rcParams['axes.grid']=True\n",
+    "\n",
+    "# https://stackoverflow.com/a/20709149\n",
+    "matplotlib.rcParams['text.usetex'] = True\n",
+    "\n",
+    "matplotlib.rcParams['text.latex.preamble'] = [\n",
+    "       r'\\usepackage{siunitx}',   # i need upright \\micro symbols, but you need...\n",
+    "       r'\\sisetup{detect-all}',   # ...this to force siunitx to actually use your fonts\n",
+    "       r'\\usepackage{helvet}',    # set the normal font here\n",
+    "       r'\\usepackage{amsmath}',\n",
+    "       r'\\usepackage{sansmath}',  # load up the sansmath so that math -> helvet\n",
+    "       r'\\sansmath',              # <- tricky! -- gotta actually tell tex to use!\n",
+    "]  "
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Import data\n",
+    "Import our data from the `./data/interim/` folder and load it into pandas dataframes. "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def df_from_csv(csv, index_col, data_folder='../data/interim'):\n",
+    "    print('Importing {}'.format(csv))\n",
+    "    return pd.read_csv(os.path.join(data_folder,csv), index_col=index_col)\n",
+    "\n",
+    "df_waves = df_from_csv('waves.csv', index_col=[0, 1])\n",
+    "df_tides = df_from_csv('tides.csv', index_col=[0, 1])\n",
+    "df_profiles = df_from_csv('profiles.csv', index_col=[0, 1, 2])\n",
+    "df_sites = df_from_csv('sites.csv', index_col=[0])\n",
+    "df_profile_features_crest_toes = df_from_csv('profile_features_crest_toes.csv', index_col=[0,1])\n",
+    "\n",
+    "# Note that the forecasted data sets should be in the same order for impacts and twls\n",
+    "impacts = {\n",
+    "    'forecasted': {\n",
+    "    'postintertidal_slope_hol86': df_from_csv('impacts_forecasted_postintertidal_slope_hol86.csv', index_col=[0]),\n",
+    "    'postintertidal_slope_nie91': df_from_csv('impacts_forecasted_postintertidal_slope_nie91.csv', index_col=[0]),\n",
+    "    'postintertidal_slope_pow18': df_from_csv('impacts_forecasted_postintertidal_slope_pow18.csv', index_col=[0]),\n",
+    "    'postintertidal_slope_sto06': df_from_csv('impacts_forecasted_postintertidal_slope_sto06.csv', index_col=[0]),\n",
+    "    'postmean_slope_hol86': df_from_csv('impacts_forecasted_postmean_slope_hol86.csv', index_col=[0]),\n",
+    "    'postmean_slope_nie91': df_from_csv('impacts_forecasted_postmean_slope_nie91.csv', index_col=[0]),\n",
+    "    'postmean_slope_pow18': df_from_csv('impacts_forecasted_postmean_slope_pow18.csv', index_col=[0]),\n",
+    "    'postmean_slope_sto06': df_from_csv('impacts_forecasted_postmean_slope_sto06.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_hol86': df_from_csv('impacts_forecasted_preintertidal_slope_hol86.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_nie91': df_from_csv('impacts_forecasted_preintertidal_slope_nie91.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_pow18': df_from_csv('impacts_forecasted_preintertidal_slope_pow18.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_sto06': df_from_csv('impacts_forecasted_preintertidal_slope_sto06.csv', index_col=[0]),\n",
+    "    'premean_slope_hol86': df_from_csv('impacts_forecasted_premean_slope_hol86.csv', index_col=[0]),\n",
+    "    'premean_slope_nie91': df_from_csv('impacts_forecasted_premean_slope_nie91.csv', index_col=[0]),\n",
+    "    'premean_slope_pow18': df_from_csv('impacts_forecasted_premean_slope_pow18.csv', index_col=[0]),\n",
+    "    'premean_slope_sto06': df_from_csv('impacts_forecasted_premean_slope_sto06.csv', index_col=[0]),\n",
+    "        },\n",
+    "    'observed': df_from_csv('impacts_observed.csv', index_col=[0])\n",
+    "    }\n",
+    "\n",
+    "\n",
+    "twls = {\n",
+    "    'forecasted': {\n",
+    "    'postintertidal_slope_hol86': df_from_csv('twl_postintertidal_slope_hol86.csv', index_col=[0,1]),\n",
+    "    'postintertidal_slope_nie91': df_from_csv('twl_postintertidal_slope_nie91.csv', index_col=[0,1]),\n",
+    "    'postintertidal_slope_pow18': df_from_csv('twl_postintertidal_slope_pow18.csv', index_col=[0,1]),\n",
+    "    'postintertidal_slope_sto06': df_from_csv('twl_postintertidal_slope_sto06.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_hol86': df_from_csv('twl_postmean_slope_hol86.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_nie91': df_from_csv('twl_postmean_slope_nie91.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_pow18': df_from_csv('twl_postmean_slope_pow18.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_sto06': df_from_csv('twl_postmean_slope_sto06.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_hol86': df_from_csv('twl_preintertidal_slope_hol86.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_nie91': df_from_csv('twl_preintertidal_slope_nie91.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_pow18': df_from_csv('twl_preintertidal_slope_pow18.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_sto06': df_from_csv('twl_preintertidal_slope_sto06.csv', index_col=[0,1]),\n",
+    "    'premean_slope_hol86': df_from_csv('twl_premean_slope_hol86.csv', index_col=[0,1]),\n",
+    "    'premean_slope_nie91': df_from_csv('twl_premean_slope_nie91.csv', index_col=[0,1]),\n",
+    "    'premean_slope_pow18': df_from_csv('twl_premean_slope_pow18.csv', index_col=[0,1]),\n",
+    "    'premean_slope_sto06': df_from_csv('twl_premean_slope_sto06.csv', index_col=[0,1]),\n",
+    "    }\n",
+    "}\n",
+    "print('Done!')"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Generate longshore plots for each beach"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "code_folding": []
+   },
+   "outputs": [],
+   "source": [
+    "beaches = list(\n",
+    "    set([\n",
+    "        x[:-4] for x in df_profiles.index.get_level_values('site_id').unique()\n",
+    "    ]))\n",
+    "\n",
+    "for beach in beaches:\n",
+    "   \n",
+    "    df_obs_impacts = impacts['observed'].loc[impacts['observed'].index.str.\n",
+    "                                             contains(beach)]\n",
+    "\n",
+    "    # Get index for each site on the beach\n",
+    "    n = [x for x in range(len(df_obs_impacts))][::-1]\n",
+    "    n_sites = [x for x in df_obs_impacts.index][::-1]\n",
+    "\n",
+    "    # Convert storm regimes to categorical datatype\n",
+    "    cat_type = CategoricalDtype(\n",
+    "        categories=['swash', 'collision', 'overwash', 'inundation'],\n",
+    "        ordered=True)\n",
+    "    df_obs_impacts.storm_regime = df_obs_impacts.storm_regime.astype(cat_type)\n",
+    "\n",
+    "    # Create figure\n",
+    "    \n",
+    "    # Determine the height of the figure, based on the number of sites.\n",
+    "    fig_height = max(6, 0.18 * len(n_sites))\n",
+    "    f, (ax1, ax2, ax3, ax4, ax5, ax6, ax7, ax8) = plt.subplots(\n",
+    "        1,\n",
+    "        8,\n",
+    "        sharey=True,\n",
+    "        figsize=(18, fig_height),\n",
+    "        gridspec_kw={'width_ratios': [4, 4, 2, 2, 2, 2, 2, 2]})\n",
+    "\n",
+    "    # ax1: Impacts\n",
+    "\n",
+    "    # Define colors for storm regime\n",
+    "    cmap = {'swash': '#1a9850', 'collision': '#fee08b', 'overwash': '#d73027'}\n",
+    "\n",
+    "    # Common marker style\n",
+    "    marker_style = {\n",
+    "        's': 60,\n",
+    "        'linewidths': 0.7,\n",
+    "        'alpha': 1,\n",
+    "        'edgecolors': 'k',\n",
+    "        'marker': 'o',\n",
+    "    }\n",
+    "\n",
+    "    # Plot observed impacts\n",
+    "    colors = [cmap.get(x) for x in df_obs_impacts.storm_regime]\n",
+    "    colors = ['#aaaaaa' if c is None else c for c in colors]\n",
+    "    ax1.scatter([0 for x in n], n, color=colors, **marker_style)\n",
+    "\n",
+    "    # Plot model impacts\n",
+    "    for i, model in enumerate(impacts['forecasted']):\n",
+    "\n",
+    "        # Only get model results for this beach\n",
+    "        df_model = impacts['forecasted'][model].loc[\n",
+    "            impacts['forecasted'][model].index.str.contains(beach)]\n",
+    "\n",
+    "        # Recast storm regimes as categorical data\n",
+    "        df_model.storm_regime = df_model.storm_regime.astype(cat_type)\n",
+    "\n",
+    "        # Assign colors\n",
+    "        colors = [cmap.get(x) for x in df_model.storm_regime]\n",
+    "        colors = ['#aaaaaa' if c is None else c for c in colors]\n",
+    "\n",
+    "        # Only plot markers which are different to the observed storm regime. \n",
+    "        # This makes it easier to find where model predictions differ\n",
+    "        y_coords = []\n",
+    "        for obs_impact, for_impact in zip(df_model.storm_regime,\n",
+    "                                          df_obs_impacts.storm_regime):\n",
+    "            if obs_impact == for_impact:\n",
+    "                y_coords.append(None)\n",
+    "            else:\n",
+    "                y_coords.append(i + 1)\n",
+    "\n",
+    "        ax1.scatter(y_coords, n, color=colors, **marker_style)\n",
+    "\n",
+    "    # Add model names to each impact on x axis\n",
+    "    ax1.set_xticks(range(len(impacts['forecasted']) + 1))\n",
+    "    ax1.set_xticklabels(['observed'] +\n",
+    "                        [x.replace('_', '\\_') for x in impacts['forecasted']])\n",
+    "    ax1.xaxis.set_tick_params(rotation=90)\n",
+    "\n",
+    "    # Add title\n",
+    "    ax1.set_title('Storm regime')\n",
+    "\n",
+    "    # Create custom legend\n",
+    "    legend_elements = [\n",
+    "        Line2D([0], [0],\n",
+    "               marker='o',\n",
+    "               color='w',\n",
+    "               label='Swash',\n",
+    "               markerfacecolor='#1a9850',\n",
+    "               markersize=8,\n",
+    "               markeredgewidth=1.0,\n",
+    "               markeredgecolor='k'),\n",
+    "        Line2D([0], [0],\n",
+    "               marker='o',\n",
+    "               color='w',\n",
+    "               label='Collision',\n",
+    "               markerfacecolor='#fee08b',\n",
+    "               markersize=8,\n",
+    "               markeredgewidth=1.0,\n",
+    "               markeredgecolor='k'),\n",
+    "        Line2D([0], [0],\n",
+    "               marker='o',\n",
+    "               color='w',\n",
+    "               label='Overwash',\n",
+    "               markerfacecolor='#d73027',\n",
+    "               markersize=8,\n",
+    "               markeredgewidth=1.0,\n",
+    "               markeredgecolor='k'),\n",
+    "    ]\n",
+    "    ax1.legend(\n",
+    "        handles=legend_elements, loc='lower center', bbox_to_anchor=(0.5, 1.1))\n",
+    "\n",
+    "    # Replace axis ticks with names of site ids\n",
+    "    ytick_labels = ax1.get_yticks().tolist()\n",
+    "    yticks = [\n",
+    "        n_sites[int(y)] if all([y >= 0, y < len(n_sites)]) else ''\n",
+    "        for y in ytick_labels\n",
+    "    ]\n",
+    "    yticks = [x.replace('_', '\\_') for x in yticks]\n",
+    "    ax1.set_yticklabels(yticks)\n",
+    "\n",
+    "    # ax2: elevations\n",
+    "\n",
+    "    # Dune elevations\n",
+    "    df_feats = df_profile_features_crest_toes.xs(['prestorm'],\n",
+    "                                                 level=['profile_type'])\n",
+    "    df_feats = df_feats.loc[df_feats.index.str.contains(beach)]\n",
+    "\n",
+    "    ax2.plot(df_feats.dune_crest_z, n, color='#fdae61')\n",
+    "    ax2.plot(df_feats.dune_toe_z, n, color='#fdae61')\n",
+    "    ax2.fill_betweenx(\n",
+    "        n,\n",
+    "        df_feats.dune_toe_z,\n",
+    "        df_feats.dune_crest_z,\n",
+    "        alpha=0.2,\n",
+    "        color='#fdae61',\n",
+    "        label='$D_{low}$ to $D_{high}$')\n",
+    "\n",
+    "    model_colors = [\n",
+    "        '#1f78b4',\n",
+    "        '#33a02c',\n",
+    "        '#e31a1c',\n",
+    "        '#6a3d9a',\n",
+    "        '#a6cee3',\n",
+    "        '#b2df8a',\n",
+    "        '#fb9a99',\n",
+    "        '#cab2d6',\n",
+    "        '#ffff99',\n",
+    "    ]\n",
+    "\n",
+    "    # Define colors to cycle through for our R_high\n",
+    "    ax2.set_prop_cycle(cycler('color', model_colors))\n",
+    "\n",
+    "    # For TWL elevations, Rhigh-Dlow and R2 axis, only plot a few models\n",
+    "    models_to_plot = [\n",
+    "        'premean_slope_hol86',\n",
+    "        'premean_slope_sto06',\n",
+    "        'preintertidal_slope_hol86',\n",
+    "        'preintertidal_slope_sto06',\n",
+    "    ]\n",
+    "    models_linewidth = 0.8\n",
+    "\n",
+    "    # Plot R_high values\n",
+    "    for model in models_to_plot:\n",
+    "\n",
+    "        # Only get model results for this beach\n",
+    "        df_model = impacts['forecasted'][model].loc[\n",
+    "            impacts['forecasted'][model].index.str.contains(beach)]\n",
+    "\n",
+    "        # Recast storm regimes as categorical data\n",
+    "        ax2.plot(\n",
+    "            df_model.R_high,\n",
+    "            n,\n",
+    "            label=model.replace('_', '\\_'),\n",
+    "            linewidth=models_linewidth)\n",
+    "\n",
+    "    # Set title, legend and labels\n",
+    "    ax2.set_title('TWL \\& dune\\nelevations')\n",
+    "    ax2.legend(loc='lower center', bbox_to_anchor=(0.5, 1.1))\n",
+    "    ax2.set_xlabel('Elevation (m AHD)')\n",
+    "#     ax2.set_xlim([0, max(df_feats.dune_crest_z)])\n",
+    "\n",
+    "    # ax3: Plot R_high - D_low\n",
+    "\n",
+    "    # Define colors to cycle through for our R_high\n",
+    "    ax3.set_prop_cycle(cycler('color', model_colors))\n",
+    "\n",
+    "    # Plot R_high values\n",
+    "    for model in models_to_plot:\n",
+    "\n",
+    "        df_model = impacts['forecasted'][model].loc[\n",
+    "            impacts['forecasted'][model].index.str.contains(beach)]\n",
+    "        # R_high - D_low\n",
+    "        ax3.plot(\n",
+    "            df_model.R_high - df_feats.dune_toe_z,\n",
+    "            n,\n",
+    "            label=model.replace('_', '\\_'),\n",
+    "            linewidth=models_linewidth)\n",
+    "\n",
+    "    ax3.axvline(x=0, color='black', linestyle=':')\n",
+    "    ax3.set_title('$R_{high}$ - $D_{low}$')\n",
+    "    ax3.set_xlabel('Height (m)')\n",
+    "#     ax3.set_xlim([-2, 2])\n",
+    "\n",
+    "    # Define colors to cycle through for our R2\n",
+    "    ax4.set_prop_cycle(cycler('color', model_colors))\n",
+    "\n",
+    "    # R_high - D_low\n",
+    "    for model in models_to_plot:\n",
+    "        df_R2 = impacts['forecasted'][model].merge(\n",
+    "            twls['forecasted'][model], on=['site_id', 'datetime'], how='left')\n",
+    "        df_R2 = df_R2.loc[df_R2.index.str.contains(beach)]\n",
+    "        ax4.plot(\n",
+    "            df_R2.R2,\n",
+    "            n,\n",
+    "            label=model.replace('_', '\\_'),\n",
+    "            linewidth=models_linewidth)\n",
+    "\n",
+    "    ax4.set_title(r'$R_{2\\%}$')\n",
+    "    ax4.set_xlabel('Height (m)')\n",
+    "#     ax4.set_xlim([0, 10])\n",
+    "\n",
+    "    # Beach slope\n",
+    "    slope_colors = [\n",
+    "        '#bebada',\n",
+    "        '#bc80bd',\n",
+    "        '#ffed6f',\n",
+    "        '#fdb462',\n",
+    "    ]\n",
+    "    ax5.set_prop_cycle(cycler('color', slope_colors))\n",
+    "    slope_models = {\n",
+    "        'prestorm mean': 'premean_slope_sto06',\n",
+    "        'poststorm mean': 'postmean_slope_sto06',\n",
+    "        'prestorm intertidal': 'preintertidal_slope_sto06',\n",
+    "        'poststorm intertidal': 'postintertidal_slope_sto06',\n",
+    "    }\n",
+    "\n",
+    "    for label in slope_models:\n",
+    "        model = slope_models[label]\n",
+    "        df_beta = impacts['forecasted'][model].merge(\n",
+    "            twls['forecasted'][model], on=['site_id', 'datetime'], how='left')\n",
+    "        df_beta = df_beta.loc[df_beta.index.str.contains(beach)]\n",
+    "        ax5.plot(df_beta.beta, n, label=label, linewidth=models_linewidth)\n",
+    "\n",
+    "    ax5.set_title(r'$\\beta$')\n",
+    "    ax5.set_xlabel('Beach slope')\n",
+    "    ax5.legend(loc='lower center', bbox_to_anchor=(0.5, 1.1))\n",
+    "    # ax5.set_xlim([0, 0.15])\n",
+    "\n",
+    "    # Need to chose a model to extract environmental parameters at maximum R_high time\n",
+    "    model = 'premean_slope_sto06'\n",
+    "    df_beach = impacts['forecasted'][model].merge(\n",
+    "        twls['forecasted'][model], on=['site_id', 'datetime'], how='left')\n",
+    "    df_beach = df_beach.loc[df_beach.index.str.contains(beach)]\n",
+    "\n",
+    "    # Wave height, wave period\n",
+    "    ax6.plot(df_beach.Hs0, n, color='#999999')\n",
+    "    ax6.set_title('$H_{s0}$')\n",
+    "    ax6.set_xlabel('Sig. wave height (m)')\n",
+    "    ax6.set_xlim([3, 5])\n",
+    "\n",
+    "    ax7.plot(df_beach.Tp, n, color='#999999')\n",
+    "    ax7.set_title('$T_{p}$')\n",
+    "    ax7.set_xlabel('Peak wave period (s)')\n",
+    "    ax7.set_xlim([8, 14])\n",
+    "\n",
+    "    ax8.plot(df_beach.tide, n, color='#999999')\n",
+    "    ax8.set_title('Tide \\& surge')\n",
+    "    ax8.set_xlabel('Elevation (m AHD)')\n",
+    "    ax8.set_xlim([0, 2])\n",
+    "\n",
+    "    plt.tight_layout()\n",
+    "    f.subplots_adjust(top=0.88)\n",
+    "    f.suptitle(beach.replace('_', '\\_'))\n",
+    "\n",
+    "    # Set minor axis ticks on each plot\n",
+    "    ax1.yaxis.set_minor_locator(MultipleLocator(1))\n",
+    "    ax1.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "    ax2.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "    ax3.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "    ax4.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "    ax5.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "    ax6.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "    ax7.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "    ax8.yaxis.grid(True, which='minor', linestyle='--', alpha=0.1)\n",
+    "\n",
+    "    # # Print to figure\n",
+    "    plt.savefig('07_{}.png'.format(beach), dpi=600, bbox_inches='tight')\n",
+    "\n",
+    "#     plt.show()\n",
+    "    plt.close()\n",
+    "    print('Done: {}'.format(beach))\n",
+    "print('Done!')"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Generate classification reports for each model\n",
+    "Use sklearn metrics to generate classification reports for each forecasting model."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import sklearn.metrics\n",
+    "\n",
+    "# Get observed impacts\n",
+    "df_obs = impacts['observed']\n",
+    "\n",
+    "# Convert storm regimes to categorical datatype\n",
+    "cat_type = CategoricalDtype(\n",
+    "    categories=['swash', 'collision', 'overwash', 'inundation'], ordered=True)\n",
+    "df_obs.storm_regime = df_obs.storm_regime.astype(cat_type)\n",
+    "\n",
+    "for model in impacts['forecasted']:\n",
+    "    df_for = impacts['forecasted'][model]\n",
+    "    df_for.storm_regime = df_for.storm_regime.astype(cat_type)\n",
+    "\n",
+    "    m = sklearn.metrics.classification_report(\n",
+    "            df_obs.storm_regime.astype(cat_type).cat.codes.values,\n",
+    "            df_for.storm_regime.astype(cat_type).cat.codes.values,\n",
+    "            labels=[0, 1, 2, 3],\n",
+    "            target_names=['swash', 'collision', 'overwash', 'inundation'])\n",
+    "    print(model)\n",
+    "    print(m)\n",
+    "    print()\n",
+    "    "
+   ]
+  }
+ ],
+ "metadata": {
+  "hide_input": false,
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.6.6"
+  },
+  "toc": {
+   "base_numbering": 1,
+   "nav_menu": {},
+   "number_sections": true,
+   "sideBar": true,
+   "skip_h1_title": false,
+   "title_cell": "Table of Contents",
+   "title_sidebar": "Contents",
+   "toc_cell": false,
+   "toc_position": {},
+   "toc_section_display": true,
+   "toc_window_display": false
+  },
+  "varInspector": {
+   "cols": {
+    "lenName": 16,
+    "lenType": 16,
+    "lenVar": 40
+   },
+   "kernels_config": {
+    "python": {
+     "delete_cmd_postfix": "",
+     "delete_cmd_prefix": "del ",
+     "library": "var_list.py",
+     "varRefreshCmd": "print(var_dic_list())"
+    },
+    "r": {
+     "delete_cmd_postfix": ") ",
+     "delete_cmd_prefix": "rm(",
+     "library": "var_list.r",
+     "varRefreshCmd": "cat(var_dic_list()) "
+    }
+   },
+   "types_to_exclude": [
+    "module",
+    "function",
+    "builtin_function_or_method",
+    "instance",
+    "_Feature"
+   ],
+   "window_display": false
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}

notebooks/07_longshore_plotting.ipynb

@@ -1,348 +0,0 @@
-{
- "cells": [
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "# Longshore plots of each beach\n",
-    "- Need to create a longshore plot of each beach to see how the variables change alongshore."
-   ]
-  },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## Setup notebook"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# Enable autoreloading of our modules. \n",
-    "# Most of the code will be located in the /src/ folder, \n",
-    "# and then called from the notebook.\n",
-    "%matplotlib inline\n",
-    "%reload_ext autoreload\n",
-    "%autoreload"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "from IPython.core.debugger import set_trace\n",
-    "\n",
-    "import pandas as pd\n",
-    "import numpy as np\n",
-    "import os\n",
-    "import decimal\n",
-    "import plotly\n",
-    "import plotly.graph_objs as go\n",
-    "import plotly.plotly as py\n",
-    "import plotly.tools as tls\n",
-    "import plotly.figure_factory as ff\n",
-    "from plotly import tools\n",
-    "import plotly.io as pio\n",
-    "from scipy import stats\n",
-    "import math\n",
-    "import matplotlib\n",
-    "from matplotlib import cm\n",
-    "import colorlover as cl\n",
-    "from tqdm import tqdm_notebook\n",
-    "from ipywidgets import widgets, Output\n",
-    "from IPython.display import display, clear_output, Image, HTML\n",
-    "from scipy import stats\n",
-    "from sklearn.metrics import confusion_matrix\n",
-    "import matplotlib.pyplot as plt\n",
-    "from scipy.interpolate import interp1d\n",
-    "from pandas.api.types import CategoricalDtype"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# Matplot lib default settings\n",
-    "plt.rcParams[\"figure.figsize\"] = (10,6)\n",
-    "plt.rcParams['axes.grid']=True\n",
-    "plt.rcParams['grid.alpha'] = 0.5\n",
-    "plt.rcParams['grid.color'] = \"grey\"\n",
-    "plt.rcParams['grid.linestyle'] = \"--\"\n",
-    "plt.rcParams['axes.grid']=True\n",
-    "\n",
-    "# https://stackoverflow.com/a/20709149\n",
-    "matplotlib.rcParams['text.usetex'] = True\n",
-    "\n",
-    "matplotlib.rcParams['text.latex.preamble'] = [\n",
-    "       r'\\usepackage{siunitx}',   # i need upright \\micro symbols, but you need...\n",
-    "       r'\\sisetup{detect-all}',   # ...this to force siunitx to actually use your fonts\n",
-    "       r'\\usepackage{helvet}',    # set the normal font here\n",
-    "       r'\\usepackage{amsmath}',\n",
-    "       r'\\usepackage{sansmath}',  # load up the sansmath so that math -> helvet\n",
-    "       r'\\sansmath',              # <- tricky! -- gotta actually tell tex to use!\n",
-    "]  "
-   ]
-  },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## Import data"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "def df_from_csv(csv, index_col, data_folder='../data/interim'):\n",
-    "    print('Importing {}'.format(csv))\n",
-    "    return pd.read_csv(os.path.join(data_folder,csv), index_col=index_col)\n",
-    "\n",
-    "df_waves = df_from_csv('waves.csv', index_col=[0, 1])\n",
-    "df_tides = df_from_csv('tides.csv', index_col=[0, 1])\n",
-    "df_profiles = df_from_csv('profiles.csv', index_col=[0, 1, 2])\n",
-    "df_sites = df_from_csv('sites.csv', index_col=[0])\n",
-    "df_profile_features_crest_toes = df_from_csv('profile_features_crest_toes.csv', index_col=[0,1])\n",
-    "\n",
-    "# Note that the forecasted data sets should be in the same order for impacts and twls\n",
-    "impacts = {\n",
-    "    'forecasted': {\n",
-    "        'foreshore_slope_sto06': df_from_csv('impacts_forecasted_foreshore_slope_sto06.csv', index_col=[0]),\n",
-    "        'mean_slope_sto06': df_from_csv('impacts_forecasted_mean_slope_sto06.csv', index_col=[0]),\n",
-    "        },\n",
-    "    'observed': df_from_csv('impacts_observed.csv', index_col=[0])\n",
-    "    }\n",
-    "\n",
-    "\n",
-    "twls = {\n",
-    "    'forecasted': {\n",
-    "        'foreshore_slope_sto06': df_from_csv('twl_foreshore_slope_sto06.csv', index_col=[0, 1]),\n",
-    "        'mean_slope_sto06':df_from_csv('twl_mean_slope_sto06.csv', index_col=[0, 1]),\n",
-    "    }\n",
-    "}\n",
-    "print('Done!')"
-   ]
-  },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## Generate plot for each beach"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "beach = 'NARRA'\n",
-    "\n",
-    "# Get the dataframe\n",
-    "df = impacts['forecasted']['mean_slope_sto06']\n",
-    "df = df.rename(columns={'storm_regime': 'forecasted_regime'})\n",
-    "\n",
-    "df_beach = df.loc[df.index.str.contains(beach)]\n",
-    "\n",
-    "# Add information about hydrodynamics at max(R_high) time\n",
-    "df_beach = df_beach.merge(\n",
-    "    twls['forecasted']['mean_slope_sto06'].drop(columns=['R_high', 'R_low']),\n",
-    "    left_on=['site_id', 'datetime'],\n",
-    "    right_on=['site_id', 'datetime'])\n",
-    "\n",
-    "# Add information about observed impacts\n",
-    "obs_impacts = impacts['observed'].rename(columns={\n",
-    "    'storm_regime': 'observed_regime'\n",
-    "}).observed_regime.to_frame()\n",
-    "df_beach = df_beach.merge(obs_impacts, left_on='site_id', right_on='site_id')\n",
-    "\n",
-    "# Convert storm regimes to categorical datatype\n",
-    "cat_type = CategoricalDtype(\n",
-    "    categories=['swash', 'collision', 'overwash', 'inundation'], ordered=True)\n",
-    "df_beach.forecasted_regime = df_beach.forecasted_regime.astype(cat_type)\n",
-    "df_beach.observed_regime = df_beach.observed_regime.astype(cat_type)\n",
-    "\n",
-    "# Get index\n",
-    "n = [x for x in range(len(df_beach))][::-1]\n",
-    "n_sites = [x for x in df_beach.index][::-1]\n",
-    "\n",
-    "f, (ax1, ax2, ax3, ax4, ax5, ax6, ax7, ax8) = plt.subplots(\n",
-    "    1,\n",
-    "    8,\n",
-    "    sharey=True,\n",
-    "    figsize=(14, 8),\n",
-    "    gridspec_kw={'width_ratios': [4, 4, 2, 2, 2, 2,2,2]})\n",
-    "\n",
-    "# Specify colors for storm regimes\n",
-    "cmap = {\n",
-    "    'swash': '#1a9850',\n",
-    "    'collision': '#fee08b',\n",
-    "    'overwash': '#d73027'\n",
-    "}\n",
-    "\n",
-    "colors = [cmap.get(x) for x in df_beach.observed_regime]\n",
-    "colors = ['#d73027' if c is None else c for c in colors]\n",
-    "\n",
-    "# Plot forecasted and observed storm regime\n",
-    "ax1.scatter(\n",
-    "    df_beach.observed_regime.cat.codes.replace(-1,np.NaN),\n",
-    "    n,\n",
-    "    color=colors,\n",
-    "    marker='o',\n",
-    "    label='Observed regime')\n",
-    "\n",
-    "ax1.scatter(\n",
-    "    df_beach.forecasted_regime.cat.codes.replace(-1,np.NaN),\n",
-    "    n,\n",
-    "    color='b',\n",
-    "    marker='o',\n",
-    "    edgecolors='black',\n",
-    "    facecolors='none',\n",
-    "    label='Forecasted regime')\n",
-    "\n",
-    "ax1.set_title('Storm\\nregime')\n",
-    "ax1.set_xticks([0,1,2,3])\n",
-    "ax1.set_xticklabels(['swash','collision','overwash','inundation'])\n",
-    "ax1.tick_params(axis='x', rotation=45)\n",
-    "ax1.legend(loc='center', bbox_to_anchor=(0.5, -0.15))\n",
-    "\n",
-    "# Replace yticks with site_ids\n",
-    "yticks = ax1.get_yticks().tolist()\n",
-    "yticks = [n_sites[int(y)] if 0 <= y <= len(n_sites) else y for y in yticks ]\n",
-    "ax1.set_yticklabels(yticks)\n",
-    "\n",
-    "# Water levels\n",
-    "ax2.plot(df_beach.R_high, n, color='#2c7bb6')\n",
-    "ax2.plot(df_beach.R_low, n, color='#2c7bb6')\n",
-    "ax2.fill_betweenx(\n",
-    "    n, df_beach.R_low, df_beach.R_high, alpha=0.2, color='#2c7bb6', label='$R_{low}$ to $R_{high}$')\n",
-    "\n",
-    "# Dune elevations\n",
-    "ax2.plot(df_beach.dune_crest_z, n, color='#fdae61')\n",
-    "ax2.plot(df_beach.dune_toe_z, n, color='#fdae61')\n",
-    "ax2.fill_betweenx(\n",
-    "    n, df_beach.dune_toe_z, df_beach.dune_crest_z, alpha=0.2, color='#fdae61', label='$D_{low}$ to $D_{high}$')\n",
-    "\n",
-    "ax2.set_title('TWL \\& Dune\\nElevations')\n",
-    "ax2.legend(loc='center',bbox_to_anchor=(0.5,-0.15))\n",
-    "ax2.set_xlabel('Elevation (m AHD)')\n",
-    "\n",
-    "# Plot R_high - D_low\n",
-    "ax3.plot(df_beach.R_high - df_beach.dune_toe_z,n,color='#999999')\n",
-    "ax3.axvline(x=0,color='black',linestyle=':')\n",
-    "ax3.set_title('$R_{high}$ - $D_{low}$')\n",
-    "ax3.set_xlabel('Height (m)')\n",
-    "ax3.set_xlim([-2,2])\n",
-    "\n",
-    "# Wave height, wave period, beach slope\n",
-    "ax4.plot(df_beach.Hs0, n,color='#377eb8')\n",
-    "ax4.set_title('$H_{s0}$')\n",
-    "ax4.set_xlabel('Sig. wave height (m)')\n",
-    "ax4.set_xlim([3,5])\n",
-    "\n",
-    "ax5.plot(df_beach.Tp, n,color='#e41a1c')\n",
-    "ax5.set_title('$T_{p}$')\n",
-    "ax5.set_xlabel('Peak wave period (s)')\n",
-    "ax5.set_xlim([8,14])\n",
-    "\n",
-    "ax6.plot(df_beach.tide, n,color='#a6cee3')\n",
-    "ax6.set_title('Tide')\n",
-    "ax6.set_xlabel('Elevation (m AHD)')\n",
-    "ax6.set_xlim([0,2])\n",
-    "\n",
-    "ax7.plot(df_beach.beta, n,color='#4daf4a')\n",
-    "ax7.set_title(r'$\\beta$')\n",
-    "ax7.set_xlabel('Mean prestorm\\nbeach slope')\n",
-    "ax7.set_xlim([0,0.15])\n",
-    "\n",
-    "ax8.plot(df_beach.R2, n,color='#6a3d9a')\n",
-    "ax8.set_title(r'$R_{2\\%}$')\n",
-    "ax8.set_xlabel('Height (m)')\n",
-    "\n",
-    "plt.tight_layout()\n",
-    "f.subplots_adjust(top=0.88)\n",
-    "f.suptitle(beach)\n",
-    "\n",
-    "\n",
-    "# Print to figure\n",
-    "plt.savefig('07-{}.png'.format(beach), dpi=600, bbox_inches='tight') \n",
-    "\n",
-    "plt.show()\n",
-    "plt.close()"
-   ]
-  }
- ],
- "metadata": {
-  "hide_input": false,
-  "kernelspec": {
-   "display_name": "Python 3",
-   "language": "python",
-   "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.6"
-  },
-  "toc": {
-   "base_numbering": 1,
-   "nav_menu": {},
-   "number_sections": true,
-   "sideBar": true,
-   "skip_h1_title": false,
-   "title_cell": "Table of Contents",
-   "title_sidebar": "Contents",
-   "toc_cell": false,
-   "toc_position": {},
-   "toc_section_display": true,
-   "toc_window_display": false
-  },
-  "varInspector": {
-   "cols": {
-    "lenName": 16,
-    "lenType": 16,
-    "lenVar": 40
-   },
-   "kernels_config": {
-    "python": {
-     "delete_cmd_postfix": "",
-     "delete_cmd_prefix": "del ",
-     "library": "var_list.py",
-     "varRefreshCmd": "print(var_dic_list())"
-    },
-    "r": {
-     "delete_cmd_postfix": ") ",
-     "delete_cmd_prefix": "rm(",
-     "library": "var_list.r",
-     "varRefreshCmd": "cat(var_dic_list()) "
-    }
-   },
-   "types_to_exclude": [
-    "module",
-    "function",
-    "builtin_function_or_method",
-    "instance",
-    "_Feature"
-   ],
-   "window_display": false
-  }
- },
- "nbformat": 4,
- "nbformat_minor": 2
-}

notebooks/09_superseded_run_comparison.ipynb

@@ -118,10 +118,22 @@
     "# Note that the forecasted data sets should be in the same order for impacts and twls\n",
     "impacts = {\n",
     "    'forecasted': {\n",
-    "        'foreshore_slope_sto06': df_from_csv('impacts_forecasted_foreshore_slope_sto06.csv', index_col=[0]),\n",
+    "    'postintertidal_slope_hol86': df_from_csv('impacts_forecasted_postintertidal_slope_hol86.csv', index_col=[0]),\n",
-    "        'mean_slope_sto06': df_from_csv('impacts_forecasted_mean_slope_sto06.csv', index_col=[0]),\n",
+    "    'postintertidal_slope_nie91': df_from_csv('impacts_forecasted_postintertidal_slope_nie91.csv', index_col=[0]),\n",
-    "        'mean_slope_nie91': df_from_csv('impacts_forecasted_mean_slope_nie91.csv', index_col=[0]),\n",
+    "    'postintertidal_slope_pow18': df_from_csv('impacts_forecasted_postintertidal_slope_pow18.csv', index_col=[0]),\n",
-    "        'mean_slope_hol86': df_from_csv('impacts_forecasted_mean_slope_hol86.csv', index_col=[0]),\n",
+    "    'postintertidal_slope_sto06': df_from_csv('impacts_forecasted_postintertidal_slope_sto06.csv', index_col=[0]),\n",
+    "    'postmean_slope_hol86': df_from_csv('impacts_forecasted_postmean_slope_hol86.csv', index_col=[0]),\n",
+    "    'postmean_slope_nie91': df_from_csv('impacts_forecasted_postmean_slope_nie91.csv', index_col=[0]),\n",
+    "    'postmean_slope_pow18': df_from_csv('impacts_forecasted_postmean_slope_pow18.csv', index_col=[0]),\n",
+    "    'postmean_slope_sto06': df_from_csv('impacts_forecasted_postmean_slope_sto06.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_hol86': df_from_csv('impacts_forecasted_preintertidal_slope_hol86.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_nie91': df_from_csv('impacts_forecasted_preintertidal_slope_nie91.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_pow18': df_from_csv('impacts_forecasted_preintertidal_slope_pow18.csv', index_col=[0]),\n",
+    "    'preintertidal_slope_sto06': df_from_csv('impacts_forecasted_preintertidal_slope_sto06.csv', index_col=[0]),\n",
+    "    'premean_slope_hol86': df_from_csv('impacts_forecasted_premean_slope_hol86.csv', index_col=[0]),\n",
+    "    'premean_slope_nie91': df_from_csv('impacts_forecasted_premean_slope_nie91.csv', index_col=[0]),\n",
+    "    'premean_slope_pow18': df_from_csv('impacts_forecasted_premean_slope_pow18.csv', index_col=[0]),\n",
+    "    'premean_slope_sto06': df_from_csv('impacts_forecasted_premean_slope_sto06.csv', index_col=[0]),\n",
     "        },\n",
     "    'observed': df_from_csv('impacts_observed.csv', index_col=[0])\n",
     "    }\n",
@@ -129,10 +141,22 @@
     "\n",
     "twls = {\n",
     "    'forecasted': {\n",
-    "        'foreshore_slope_sto06': df_from_csv('twl_foreshore_slope_sto06.csv', index_col=[0, 1]),\n",
+    "    'postintertidal_slope_hol86.csv': df_from_csv('twl_postintertidal_slope_hol86.csv', index_col=[0,1]),\n",
-    "        'mean_slope_sto06':df_from_csv('twl_mean_slope_sto06.csv', index_col=[0, 1]),\n",
+    "    'postintertidal_slope_nie91.csv': df_from_csv('twl_postintertidal_slope_nie91.csv', index_col=[0,1]),\n",
-    "        'mean_slope_nie91':df_from_csv('twl_mean_slope_nie91.csv', index_col=[0, 1]),\n",
+    "    'postintertidal_slope_pow18.csv': df_from_csv('twl_postintertidal_slope_pow18.csv', index_col=[0,1]),\n",
-    "        'mean_slope_hol86':df_from_csv('twl_mean_slope_hol86.csv', index_col=[0, 1]),\n",
+    "    'postintertidal_slope_sto06.csv': df_from_csv('twl_postintertidal_slope_sto06.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_hol86.csv': df_from_csv('twl_postmean_slope_hol86.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_nie91.csv': df_from_csv('twl_postmean_slope_nie91.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_pow18.csv': df_from_csv('twl_postmean_slope_pow18.csv', index_col=[0,1]),\n",
+    "    'postmean_slope_sto06.csv': df_from_csv('twl_postmean_slope_sto06.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_hol86.csv': df_from_csv('twl_preintertidal_slope_hol86.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_nie91.csv': df_from_csv('twl_preintertidal_slope_nie91.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_pow18.csv': df_from_csv('twl_preintertidal_slope_pow18.csv', index_col=[0,1]),\n",
+    "    'preintertidal_slope_sto06.csv': df_from_csv('twl_preintertidal_slope_sto06.csv', index_col=[0,1]),\n",
+    "    'premean_slope_hol86.csv': df_from_csv('twl_premean_slope_hol86.csv', index_col=[0,1]),\n",
+    "    'premean_slope_nie91.csv': df_from_csv('twl_premean_slope_nie91.csv', index_col=[0,1]),\n",
+    "    'premean_slope_pow18.csv': df_from_csv('twl_premean_slope_pow18.csv', index_col=[0,1]),\n",
+    "    'premean_slope_sto06.csv': df_from_csv('twl_premean_slope_sto06.csv', index_col=[0,1]),\n",
     "    }\n",
     "}\n",
     "print('Done!')"

notebooks/10_profile_clustering.ipynb

@@ -0,0 +1,186 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import pandas as pd\n",
+    "import matplotlib.pyplot as plt\n",
+    "import os\n",
+    "\n",
+    "from dtaidistance import dtw\n",
+    "from dtaidistance import dtw_visualisation as dtwvis"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Import data\n",
+    "def df_from_csv(csv, index_col, data_folder='../data/interim'):\n",
+    "    print('Importing {}'.format(csv))\n",
+    "    return pd.read_csv(os.path.join(data_folder,csv), index_col=index_col)\n",
+    "\n",
+    "df_profiles = df_from_csv('profiles.csv', index_col=[0, 1, 2])\n",
+    "\n",
+    "print('Done!')"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Use dtaidistance package"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "p1 = df_profiles.dropna(subset=['z']).xs(['AVOCAn0003','prestorm'],level=['site_id','profile_type']).z.values\n",
+    "p2 = df_profiles.dropna(subset=['z']).xs(['AVOCAn0004','prestorm'],level=['site_id','profile_type']).z.values\n",
+    "path = dtw.warping_path(p1,p2)\n",
+    "dtwvis.plot_warping(p1,p2,path)\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Use kshape package"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "profiles = df_profiles.dropna(subset=['z'])\\\n",
+    "            .xs(['prestorm'],level=['profile_type'])\\\n",
+    "            .groupby('site_id').z\\\n",
+    "            .apply(list).tolist()\n",
+    "\n",
+    "\n",
+    "# profiles = [x[-50:] for x in profiles]\n",
+    "# print(min(len(x) for x in profiles))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from kshape.core import kshape, zscore\n",
+    "\n",
+    "time_series = [[1,2,3,4], [0,1,2,3], [0,1,2,3], [1,2,2,3]]\n",
+    "cluster_num = 4\n",
+    "clusters = kshape(zscore(profiles, axis=1), cluster_num)\n",
+    "# print(clusters)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "cluster_no = 0\n",
+    "\n",
+    "# Plot shape of all clusters\n",
+    "plt.figure(0)\n",
+    "for n,cluster in enumerate(clusters):\n",
+    "    plt.plot(cluster[0],label=n)\n",
+    "plt.legend()\n",
+    "\n",
+    "plt.figure(1)\n",
+    "# Plot all profiles in partiuclar cluster\n",
+    "for profile_no in clusters[cluster_no][1]:\n",
+    "    plt.plot(profiles[profile_no])\n",
+    "\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "a = [1,2,3,4,5,6]\n",
+    "a[-1:]"
+   ]
+  }
+ ],
+ "metadata": {
+  "hide_input": false,
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.6.6"
+  },
+  "toc": {
+   "base_numbering": 1,
+   "nav_menu": {},
+   "number_sections": true,
+   "sideBar": true,
+   "skip_h1_title": false,
+   "title_cell": "Table of Contents",
+   "title_sidebar": "Contents",
+   "toc_cell": false,
+   "toc_position": {},
+   "toc_section_display": true,
+   "toc_window_display": false
+  },
+  "varInspector": {
+   "cols": {
+    "lenName": 16,
+    "lenType": 16,
+    "lenVar": 40
+   },
+   "kernels_config": {
+    "python": {
+     "delete_cmd_postfix": "",
+     "delete_cmd_prefix": "del ",
+     "library": "var_list.py",
+     "varRefreshCmd": "print(var_dic_list())"
+    },
+    "r": {
+     "delete_cmd_postfix": ") ",
+     "delete_cmd_prefix": "rm(",
+     "library": "var_list.r",
+     "varRefreshCmd": "cat(var_dic_list()) "
+    }
+   },
+   "types_to_exclude": [
+    "module",
+    "function",
+    "builtin_function_or_method",
+    "instance",
+    "_Feature"
+   ],
+   "window_display": false
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}

src/analysis/forecast_twl.py

@@ -20,6 +20,7 @@ def forecast_twl(
     df_profiles,
     df_waves,
     df_profile_features,
+    df_grain_size,
     runup_function,
     n_processes=MULTIPROCESS_THREADS,
     slope="foreshore",
@@ -34,10 +35,11 @@ def forecast_twl(
 
     # Estimate foreshore slope. Do the analysis per site_id. This is so we only have to query the x and z
     # cross-section profiles once per site.
-    logger.info("Calculating beach slopes")
+
     site_ids = df_twl.index.get_level_values("site_id").unique()
 
     if slope == "foreshore":
+        logger.info("Calculating foreshore slopes")
         # Process each site_id with a different process and combine results at the end
         with Pool(processes=n_processes) as pool:
             results = pool.starmap(
@@ -47,41 +49,79 @@ def forecast_twl(
         df_twl["beta"] = pd.concat(results)
 
     elif slope == "mean":
-        df_temp = df_twl.join(
+        logger.info("Calculating mean (dune toe to MHW) slopes")
-            df_profile_features.query(
+        btm_z = 0.5  # m AHD
-                "profile_type=='{}'".format(profile_type)
-            ).reset_index(level="profile_type"),
-            how="inner",
-        )
-        df_temp["mhw"] = 0.5
 
         # When calculating mean slope, we go from the dune toe to mhw. However, in some profiles, the dune toe is not
-        # defined. In these cases, we should go to the dune crest
+        # defined. In these cases, we should go to the dune crest. Let's make a temporary dataframe which has this
-        df_temp["top_elevation"] = df_temp["dune_toe_z"]
+        # already calculated.
-        df_temp.loc[df_temp.dune_toe_z.isnull(), "top_elevation"] = df_temp.loc[
+        df_top_ele = df_profile_features.xs(profile_type, level="profile_type").copy()
-            df_temp.dune_toe_z.isnull(), "dune_crest_z"
+        df_top_ele.loc[:, "top_ele"] = df_top_ele.dune_toe_z
-        ]
+        df_top_ele.loc[
-        df_temp["top_x"] = df_temp["dune_toe_x"]
+            df_top_ele.top_ele.isnull().values, "top_ele"
-        df_temp.loc[df_temp.dune_toe_x.isnull(), "top_x"] = df_temp.loc[
+        ] = df_top_ele.dune_crest_z
-            df_temp.dune_toe_x.isnull(), "dune_crest_x"
+
-        ]
+        n_no_top_ele = len(df_top_ele[df_top_ele.top_ele.isnull()].index)
+        if n_no_top_ele != 0:
+            logger.warning(
+                "{} sites do not have dune toes/crests to calculate mean slope".format(
+                    n_no_top_ele
+                )
+            )
 
-        with Pool(processes=n_processes) as pool:
+        df_slopes = (
-            results = pool.starmap(
+            df_profiles.xs(profile_type, level="profile_type")
-                mean_slope_for_site_id,
+            .dropna(subset=["z"])
-                [
+            .groupby("site_id")
-                    (site_id, df_temp, df_profiles, "top_elevation", "top_x", "mhw")
+            .apply(
-                    for site_id in site_ids
+                lambda x: slope_from_profile(
-                ],
+                    profile_x=x.index.get_level_values("x").tolist(),
+                    profile_z=x.z.tolist(),
+                    top_elevation=df_top_ele.loc[x.index[0][0], :].top_ele,
+                    btm_elevation=btm_z,
+                    method="least_squares",
+                )
+            )
+            .rename("beta")
+            .to_frame()
         )
 
-        df_twl["beta"] = pd.concat(results)
+        # Merge calculated slopes onto each twl timestep
+        df_twl = df_twl.merge(df_slopes, left_index=True, right_index=True)
+
+    elif slope == "intertidal":
+
+        logger.info("Calculating intertidal slopes")
+        top_z = 1.15  # m AHD = HAT from MHL annual ocean tides summary report
+        btm_z = -0.9  # m AHD = HAT from MHL annual ocean tides summary report
+
+        # Calculate slopes for each profile
+        df_slopes = (
+            df_profiles.xs(profile_type, level="profile_type")
+            .dropna(subset=["z"])
+            .groupby("site_id")
+            .apply(
+                lambda x: slope_from_profile(
+                    profile_x=x.index.get_level_values("x").tolist(),
+                    profile_z=x.z.tolist(),
+                    top_elevation=top_z,
+                    btm_elevation=max(min(x.z), btm_z),
+                    method="least_squares",
+                )
+            )
+            .rename("beta")
+            .to_frame()
+        )
+
+        # Merge calculated slopes onto each twl timestep
+        df_twl = df_twl.merge(df_slopes, left_index=True, right_index=True)
 
     # Estimate runup
     R2, setup, S_total, S_inc, S_ig = runup_function(
         Hs0=df_twl["Hs0"].tolist(),
         Tp=df_twl["Tp"].tolist(),
         beta=df_twl["beta"].tolist(),
+        r=df_twl.merge(df_grain_size, on="site_id").r.tolist(),
     )
 
     df_twl["R2"] = R2
@@ -133,7 +173,7 @@ def mean_slope_for_site_id(
             profile_z=profile_z,
             top_elevation=row[top_elevation_col],
             btm_elevation=row[btm_elevation_col],
-            method="end_points",
+            method="least_squares",
             top_x=row[top_x_col],
         ),
         axis=1,
@@ -209,7 +249,7 @@ def foreshore_slope_from_profile(profile_x, profile_z, tide, runup_function, **k
         beta_new = slope_from_profile(
             profile_x=profile_x,
             profile_z=profile_z,
-            method="end_points",
+            method="least_squares",
             top_elevation=tide + setup + S_total / 2,
             btm_elevation=tide + setup - S_total / 2,
         )
@@ -297,10 +337,16 @@ def slope_from_profile(
                 ]
                 if len(end_point_btm) == 0:
                     # If there doesn't seem to be an intersection seaward of the top elevation, return none.
+                    logger.warning("No intersections found seaward of top elevation")
                     return None
                 else:
                     end_points[end_type]["x"] = end_point_btm[0]
 
+    # Ensure that top point is landward of bottom point
+    if end_points["top"]["x"] > end_points["btm"]["x"]:
+        logger.warning("Top point is not landward of bottom point")
+        return None
+
     if method == "end_points":
         x_top = end_points["top"]["x"]
         x_btm = end_points["btm"]["x"]
@@ -309,6 +355,7 @@ def slope_from_profile(
         return -(z_top - z_btm) / (x_top - x_btm)
 
     elif method == "least_squares":
+
         profile_mask = [
             True if end_points["top"]["x"] < pts < end_points["btm"]["x"] else False
             for pts in profile_x
@@ -328,10 +375,13 @@ 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"])
+    "--slope",
+    required=True,
+    help="",
+    type=click.Choice(["foreshore", "mean", "intertidal"]),
 )
 @click.option(
     "--profile-type",
@@ -340,6 +390,7 @@ def slope_from_profile(
     type=click.Choice(["prestorm", "poststorm"]),
 )
 @click.option("--output-file", required=True, help="")
+@click.option("--grain-size-csv", required=False, help="")
 def create_twl_forecast(
     waves_csv,
     tides_csv,
@@ -349,6 +400,7 @@ def create_twl_forecast(
     slope,
     profile_type,
     output_file,
+    grain_size_csv,
 ):
     logger.info("Creating forecast of total water levels")
     logger.info("Importing data")
@@ -356,6 +408,7 @@ def create_twl_forecast(
     df_tides = pd.read_csv(tides_csv, index_col=[0, 1])
     df_profiles = pd.read_csv(profiles_csv, index_col=[0, 1, 2])
     df_profile_features = pd.read_csv(profile_features_csv, index_col=[0, 1])
+    df_grain_size = pd.read_csv(grain_size_csv, index_col=[0])
 
     logger.info("Forecasting TWL")
     df_twl = forecast_twl(
@@ -363,6 +416,7 @@ def create_twl_forecast(
         df_profiles,
         df_waves,
         df_profile_features,
+        df_grain_size,
         runup_function=getattr(runup_models, runup_function),
         slope=slope,
         profile_type=profile_type,

src/analysis/observed_storm_impacts.py

@@ -4,6 +4,7 @@ import pandas as pd
 from scipy.integrate import simps
 
 from logs import setup_logging
+from utils import crossings
 
 logger = setup_logging()
 
@@ -41,14 +42,9 @@ def volume_change(df_profiles, df_profile_features, zone):
             "Calculating change in beach volume at {} in {} zone".format(site_id, zone)
         )
 
-        # TODO Change this query to an index
+        prestorm_row = df_profile_features.loc[(site_id, "prestorm")]
-        query = "site_id=='{}'&profile_type=='prestorm'".format(site_id)
+        prestorm_dune_toe_x = prestorm_row.dune_toe_x
-        prestorm_dune_toe_x = df_profile_features.query(query).dune_toe_x.tolist()
+        prestorm_dune_crest_x = prestorm_row.dune_crest_x
-        prestorm_dune_crest_x = df_profile_features.query(query).dune_crest_x.tolist()
-
-        # We may not have a dune toe or crest defined, or there may be multiple defined.
-        prestorm_dune_crest_x = return_first_or_nan(prestorm_dune_crest_x)
-        prestorm_dune_toe_x = return_first_or_nan(prestorm_dune_toe_x)
 
         # If no dune toe has been defined, Dlow = Dhigh. Refer to Sallenger (2000).
         if np.isnan(prestorm_dune_toe_x):
@@ -56,8 +52,7 @@ def volume_change(df_profiles, df_profile_features, zone):
 
         # If no prestorm and poststorm profiles, skip site and continue
         profile_lengths = [
-            len(df_site.query("profile_type == '{}'".format(x)))
+            len(df_site.xs(x, level="profile_type")) for x in ["prestorm", "poststorm"]
-            for x in ["prestorm", "poststorm"]
         ]
         if any([length == 0 for length in profile_lengths]):
             continue
@@ -113,20 +108,52 @@ def volume_change(df_profiles, df_profile_features, zone):
             x_max=x_max,
         )
 
-        # Volume change needs to be calculated including a tolerance for LIDAR accuracy. If difference between
+        # Identify the x location where our pre and post storm profiles start to differ. This is so changes no due to
-        # profiles is less than 20 cm, consider them as zero difference.
+        # the storm are not included when calculating volume.
-        prestorm_z = (
+        df_prestorm = (
-            df_zone.query("profile_type=='prestorm'").reset_index("profile_type").z
+            df_site.xs("prestorm", level="profile_type").z.rename("z_pre").to_frame()
         )
-        poststorm_z = (
+        df_poststorm = (
-            df_zone.query("profile_type=='poststorm'").reset_index("profile_type").z
+            df_site.xs("poststorm", level="profile_type").z.rename("z_post").to_frame()
         )
-        diff_z = prestorm_z - poststorm_z
+        df_diff = df_prestorm.merge(df_poststorm, on=["site_id", "x"])
-        diff_z[abs(diff_z) < 0.2] = 0
+        df_diff["z_diff"] = df_diff.z_pre - df_diff.z_post
+
+        # Find all locations where the difference in pre and post storm is zero. Take the most seaward location as the
+        # x location where our profiles are the same.
+        x_crossings = crossings(df_diff.index.get_level_values("x"), df_diff.z_diff, 0)
+
+        if len(x_crossings) != 0:
+            x_change_point = x_crossings[-1]
+        else:
+            x_change_point = np.nan
+
+        # # For debugging
+        # import matplotlib.pyplot as plt
+        # f,(ax1,ax2) = plt.subplots(2,1,sharex=True)
+        # ax1.plot(df_prestorm.index.get_level_values('x'), df_prestorm.z_pre,label='prestorm')
+        # ax1.plot(df_poststorm.index.get_level_values('x'), df_poststorm.z_post,label='poststorm')
+        # ax1.axvline(x_crossings[-1], color='red', linestyle='--', linewidth=0.5, label='Change point')
+        # ax1.legend()
+        # ax1.set_title(site_id)
+        # ax1.set_ylabel('elevation (m AHD)')
+        # ax2.plot(df_diff.index.get_level_values('x'), df_diff.z_diff)
+        # ax2.set_xlabel('x coordinate (m)')
+        # ax2.set_ylabel('elevation diff (m)')
+        # ax2.axvline(x_crossings[-1],color='red',linestyle='--',linewidth=0.5)
+        # plt.show()
+
         diff_vol = beach_volume(
-            x=diff_z.index.get_level_values("x"), z=diff_z, x_min=x_min, x_max=x_max
+            x=df_diff.index.get_level_values("x"),
+            z=df_diff.z_diff,
+            x_min=np.nanmax([x_min, x_change_point]),
+            x_max=np.nanmax([x_max, x_change_point]),
         )
 
+        # Here, if cannot calculate the difference volume, assume no volume change
+        if np.isnan(diff_vol):
+            diff_vol = 0
+
         df_vol_changes.loc[site_id, "prestorm_{}_vol".format(zone)] = prestorm_vol
         df_vol_changes.loc[site_id, "poststorm_{}_vol".format(zone)] = poststorm_vol
         df_vol_changes.loc[site_id, "{}_vol_change".format(zone)] = diff_vol

src/analysis/runup_models.py

@@ -1,8 +1,11 @@
 import numpy as np
 import pandas as pd
+from logs import setup_logging
 
+logger = setup_logging()
 
-def sto06(Hs0, Tp, beta):
+
+def sto06(Hs0, Tp, beta, **kwargs):
     """
     :param Hs0: List or float of offshore significant wave height values
     :param Tp: List or float of peak wave period
@@ -43,8 +46,7 @@ def sto06(Hs0, Tp, beta):
     )
 
 
-def hol86(Hs0, Tp, beta):
+def hol86(Hs0, Tp, beta, **kwargs):
-
     df = pd.DataFrame(
         {"Hs0": Hs0, "Tp": Tp, "beta": beta}, index=[x for x in range(0, np.size(Hs0))]
     )
@@ -67,8 +69,7 @@ def hol86(Hs0, Tp, beta):
     )
 
 
-def nie91(Hs0, Tp, beta):
+def nie91(Hs0, Tp, beta, **kwargs):
-
     df = pd.DataFrame(
         {"Hs0": Hs0, "Tp": Tp, "beta": beta}, index=[x for x in range(0, np.size(Hs0))]
     )
@@ -105,8 +106,98 @@ def atk18(Hs0, Tp, beta):
     pass
 
 
-def pow18(Hs0, Tp, beta):
+def pow18(Hs0, Tp, beta, r, **kwargs):
-    pass
+    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):

src/utils.py

@@ -23,6 +23,10 @@ def crossings(profile_x, profile_z, constant_z):
     profile_z = np.array(profile_z)[valid]
     profile_x = np.array(profile_x)[valid]
 
+    # Return empty list if mask removes all values
+    if profile_x.size == 0:
+        return []
+
     # Normalize the 'signal' to zero.
     # Use np.subtract rather than a list comprehension for performance reasons
     z = np.subtract(profile_z, constant_z)
@@ -32,11 +36,20 @@ def crossings(profile_x, profile_z, constant_z):
     indicies = np.where(z[:-1] * z[1:] < 0)[0]
 
     # Use linear interpolation to find intersample crossings.
-    return [
+    x_crossings = [
         profile_x[i] - (profile_x[i] - profile_x[i + 1]) / (z[i] - z[i + 1]) * (z[i])
         for i in indicies
     ]
 
+    # Also need to check the end points as the above will not include them if they are close
+    x_crossings += [
+        profile_x[i]
+        for i in [0, len(profile_z) - 1]
+        if np.isclose(constant_z, profile_z[i])
+    ]
+
+    return sorted(x_crossings)
+
 
 # TODO Think of a better way to do this than having to manually specify the coordinate systems
 def convert_coord_systems(