From 5f59c8f8ee6c0d50e5d5ab2da0af52a9c99f9a23 Mon Sep 17 00:00:00 2001 From: Chris Leaman Date: Fri, 7 Dec 2018 12:01:08 +1100 Subject: [PATCH] Update slope function to accept top x coordinate Needed to ensure correctly mean slope is calculated. If only inputting dune toe z coordinates, sometimes the function can chose the wrong x coordinate if there are multiple crossings of the z coordinate. --- src/analysis/forecast_twl.py | 38 ++++++++++++++++++++++++++++-------- 1 file changed, 30 insertions(+), 8 deletions(-) diff --git a/src/analysis/forecast_twl.py b/src/analysis/forecast_twl.py index 1bb29c1..82efc65 100644 --- a/src/analysis/forecast_twl.py +++ b/src/analysis/forecast_twl.py @@ -23,6 +23,7 @@ def forecast_twl( runup_function, n_processes=MULTIPROCESS_THREADS, slope="foreshore", + profile_type='prestorm' ): # Use df_waves as a base df_twl = df_waves.copy() @@ -45,12 +46,14 @@ def forecast_twl( df_twl["beta"] = pd.concat(results) elif slope == "mean": - df_temp = df_twl.join(df_profile_features.query("profile_type=='prestorm'").reset_index(level='profile_type') + df_temp = df_twl.join(df_profile_features.query("profile_type=='{}'".format(profile_type)).reset_index( + level='profile_type') , how="inner") df_temp["mhw"] = 0.5 with Pool(processes=n_processes) as pool: results = pool.starmap( - mean_slope_for_site_id, [(site_id, df_temp, df_profiles, "dune_toe_z", "mhw") for site_id in site_ids] + mean_slope_for_site_id, [(site_id, df_temp, df_profiles, "dune_toe_z", "dune_toe_x", "mhw") for + site_id in site_ids] ) df_twl["beta"] = pd.concat(results) @@ -71,7 +74,8 @@ def forecast_twl( return df_twl -def mean_slope_for_site_id(site_id, df_twl, df_profiles, top_elevation_col, btm_elevation_col): +def mean_slope_for_site_id(site_id, df_twl, df_profiles, top_elevation_col, top_x_col, btm_elevation_col, + profile_type='prestorm'): """ Calculates the foreshore slope values a given site_id. Returns a series (with same indicies as df_twl) of foreshore slopes. This function is used to parallelize getting foreshore slopes as it is computationally @@ -83,7 +87,7 @@ def mean_slope_for_site_id(site_id, df_twl, df_profiles, top_elevation_col, btm_ """ # Get the prestorm beach profile - profile = df_profiles.query("site_id =='{}' and profile_type == 'prestorm'".format(site_id)) + profile = df_profiles.query("site_id =='{}' and profile_type == '{}'".format(site_id, profile_type)) profile_x = profile.index.get_level_values("x").tolist() profile_z = profile.z.tolist() @@ -96,6 +100,7 @@ def mean_slope_for_site_id(site_id, df_twl, df_profiles, top_elevation_col, btm_ top_elevation=row[top_elevation_col], btm_elevation=row[btm_elevation_col], method="end_points", + top_x= row[top_x_col] ), axis=1, ) @@ -191,7 +196,7 @@ def foreshore_slope_from_profile(profile_x, profile_z, tide, runup_function, **k iteration_count += 1 -def slope_from_profile(profile_x, profile_z, top_elevation, btm_elevation, method="end_points"): +def slope_from_profile(profile_x, profile_z, top_elevation, btm_elevation, method="end_points", top_x=None, btm_x=None): """ Returns a slope (beta) from a bed profile, given the top and bottom elevations of where the slope should be taken. :param x: List of x bed profile coordinates @@ -199,6 +204,9 @@ def slope_from_profile(profile_x, profile_z, top_elevation, btm_elevation, metho :param top_elevation: Top elevation of where to take the slope :param btm_elevation: Bottom elevation of where to take the slope :param method: Method used to calculate slope (end_points or least_squares) + :param top_x: x-coordinate of the top end point. May be needed, as there may be multiple crossings of the + top_elevation. + :param btm_x: x-coordinate of the bottom end point :return: """ @@ -208,7 +216,18 @@ def slope_from_profile(profile_x, profile_z, top_elevation, btm_elevation, metho end_points = {"top": {"z": top_elevation}, "btm": {"z": btm_elevation}} + + for end_type in end_points.keys(): + + # Add x coordinates if they are specified + if top_x and end_type == 'top': + end_points['top']['x'] = top_x + continue + if btm_x and end_type == 'top': + end_points['btm']['x'] = btm_x + continue + elevation = end_points[end_type]["z"] intersection_x = crossings(profile_x, profile_z, elevation) @@ -285,8 +304,10 @@ def crossings(profile_x, profile_z, constant_z): @click.option("--profile-features-csv", required=True, help="") @click.option("--runup-function", required=True, help="", type=click.Choice(["sto06"])) @click.option("--slope", required=True, help="", type=click.Choice(["foreshore", "mean"])) +@click.option("--profile-type", required=True, help="", type=click.Choice(["prestorm", "poststorm"])) @click.option("--output-file", required=True, help="") -def create_twl_forecast(waves_csv, tides_csv, profiles_csv, profile_features_csv, runup_function, slope, output_file): +def create_twl_forecast(waves_csv, tides_csv, profiles_csv, profile_features_csv, runup_function, slope, + profile_type,output_file): logger.info("Creating forecast of total water levels") logger.info("Importing data") df_waves = pd.read_csv(waves_csv, index_col=[0, 1]) @@ -295,15 +316,16 @@ def create_twl_forecast(waves_csv, tides_csv, profiles_csv, profile_features_csv df_profile_features = pd.read_csv(profile_features_csv, index_col=[0,1]) logger.info("Forecasting TWL") - df_twl_foreshore_slope_sto06 = forecast_twl( + df_twl = forecast_twl( df_tides, df_profiles, df_waves, df_profile_features, runup_function=getattr(runup_models, runup_function), slope=slope, + profile_type=profile_type ) - df_twl_foreshore_slope_sto06.to_csv(output_file) + df_twl.to_csv(output_file) logger.info("Saved to %s", output_file) logger.info("Done!")