diff --git a/src/Model/Geometry/Profile.py b/src/Model/Geometry/Profile.py index d28da7d652040ac802f4848b3459071947b0e484..ad1b82252ab790872796e8f402ac428aea0ef438 100644 --- a/src/Model/Geometry/Profile.py +++ b/src/Model/Geometry/Profile.py @@ -345,3 +345,9 @@ class Profile(object): # Abstract method, must be implemented for in non abstract class def get_station(self): raise NotImplementedMethodeError(self, self.get_station) + + # Computation method + + # Abstract method for width approximation + def width_approximation(self): + raise NotImplementedMethodeError(self, self.width_approximation) diff --git a/src/Model/Geometry/ProfileXYZ.py b/src/Model/Geometry/ProfileXYZ.py index 4375d18cba7a600bf1e3ad1262cd69cad3a2403b..dfe187024d3c8cee425c54adfa07e09985a47152 100644 --- a/src/Model/Geometry/ProfileXYZ.py +++ b/src/Model/Geometry/ProfileXYZ.py @@ -287,6 +287,27 @@ class ProfileXYZ(Profile, SQLSubModel): except IndexError: raise IndexError(f"Invalid point index: {index}") + def get_point_by_name(self, name: str) -> PointXYZ: + """Get point by name. + + Args: + name: Point name. + + Returns: + The point. + """ + try: + n_name = name.lower().strip() + return next( + filter( + lambda p: p.name.lower().strip() == n_name, + self.points + ) + ) + except Exception as e: + logger.debug(f"{e}") + raise IndexError(f"Invalid point name: {name}") + def add(self): """Add a new PointXYZ to profile. @@ -415,3 +436,9 @@ class ProfileXYZ(Profile, SQLSubModel): constant = station[index_profile_z_min[0]] return list(map(lambda s: s - constant, station)) + + def width_approximation(self): + rg = self.get_point_by_name("rg") + rd = self.get_point_by_name("rd") + + return abs(rg.dist(rd)) diff --git a/src/Model/Geometry/Reach.py b/src/Model/Geometry/Reach.py index c3ee413cd6430e3649f12f3a372b1252e0571412..5351a9cbbabeb1dd7cf9a62012ccb8542b779168 100644 --- a/src/Model/Geometry/Reach.py +++ b/src/Model/Geometry/Reach.py @@ -578,3 +578,103 @@ class Reach(SQLSubModel): file_st.write(" 999.9990 999.9990 999.9990") file_st.write("\n") + + def get_incline(self): + first = self.profile(0) + last = self.profile(len(self) - 1) + + z_first = first.z_min() + kp_first = first.kp + + z_last = last.z_min() + kp_last = last.kp + + return ( + (z_last - z_first) + / + (kp_last - kp_first) + ) + + def get_incline_mean(self): + profiles = self.profiles + previous = profiles[0] + + incline_acc = 0 + + for profile in profiles[1:]: + z_first = previous.z_min() + kp_first = previous.kp + + z_last = profile.z_min() + kp_last = profile.kp + + incline_acc += ( + (z_last - z_first) + / + (kp_last - kp_first) + ) + + previous = profile + + return (incline_acc / (len(profiles) - 1)) + + def get_incline_median(self): + profiles = self.profiles + previous = profiles[0] + + incline_acc = [] + + for profile in profiles[1:]: + z_first = previous.z_min() + kp_first = previous.kp + + z_last = profile.z_min() + kp_last = profile.kp + + incline_acc += [ + (z_last - z_first) + / + (kp_last - kp_first) + ] + + previous = profile + + incline_acc.sort() + return incline_acc[round(len(profiles)/2)] + + def get_incline_median_mean(self): + profiles = self.profiles + previous = profiles[0] + + incline_acc = [] + + for profile in profiles[1:]: + z_first = previous.z_min() + kp_first = previous.kp + + z_last = profile.z_min() + kp_last = profile.kp + + incline_acc += [ + (z_last - z_first) + / + (kp_last - kp_first) + ] + + previous = profile + + incline_acc.sort() + + marge = round(len(incline_acc) * 0.1) + incline_set = incline_acc[marge:-marge] + + logger.debug(f"+{incline_acc}") + logger.debug(f"-{incline_set}") + + return ( + reduce( + lambda acc, x: acc + x, + incline_set, + 0.0 + ) / (len(incline_set)) + ) diff --git a/src/Model/InitialConditions/InitialConditions.py b/src/Model/InitialConditions/InitialConditions.py index 63d9b40a5fd3c769d0bed9595fd7d4e7900eb91e..c81f940c119b0c5ad1958513f9c56e642d99a1ec 100644 --- a/src/Model/InitialConditions/InitialConditions.py +++ b/src/Model/InitialConditions/InitialConditions.py @@ -16,12 +16,16 @@ # -*- coding: utf-8 -*- +import logging + from copy import copy, deepcopy from tools import trace, timer from functools import reduce from Model.Tools.PamhyrDB import SQLSubModel +logger = logging.getLogger() + class Data(SQLSubModel): def __init__(self, name: str = "", @@ -358,37 +362,79 @@ class InitialConditions(SQLSubModel): profiles = self._reach.reach.profiles.copy() self._sort_by_z_and_kp(profiles) - prev = None + incline = self._reach.reach.get_incline_median_mean() + logger.debug(f"incline = {incline}") + + previous_elevation = -99999.99 for profile in profiles: + width = profile.width_approximation() + strickler = 25 + discharge = ( + ((width * 0.8) + * strickler + * (height ** (5/3)) + * (abs(incline) ** (0.5))) + ) + + elevation= max( + profile.z_min() + height, + previous_elevation + ) + + logger.debug(f"({profile.kp}):") + logger.debug(f" width = {width}") + logger.debug(f" strickler = {strickler}") + logger.debug(f" discharge = {discharge}") + new = Data(reach=self._reach, status=self._status) new["kp"] = profile.kp + new["discharge"] = discharge - if prev is None: - new["elevation"] = profile.z_min() + height - else: - new["elevation"] = max( - profile.z_min() + height, - prev["elevation"] - ) + new["elevation"] = elevation self._data.append(new) - prev = new + previous_elevation = elevation + - is_reverse = False - if profiles[0].kp > profiles[-1].kp: - is_reverse = True + def generate_discharge(self, discharge: float): + self._data = [] - self._data.sort( - reverse=not is_reverse, - key=lambda d: d['kp'] + self._generate_height_estimation_from_discharge( + discharge ) - def generate_discharge(self, discharge: float): - self._new = [] + def _generate_height_estimation_from_discharge(self, discharge: float): + profiles = self._reach.reach.profiles.copy() + self._sort_by_z_and_kp(profiles) - for d in self._data: - n = d.copy() - n['discharge'] = discharge - self._new.append(n) + previous_elevation = -99999.99 + + incline = self._reach.reach.get_incline_median_mean() + logger.debug(f"incline = {incline}") + + for profile in profiles: + width = profile.width_approximation() + strickler = 25 + height = ( + discharge + / + ((width * 0.8) * strickler * (abs(incline) ** (0.5))) + ) ** (0.6) + + elevation= max( + profile.z_min() + height, + previous_elevation + ) + + logger.debug(f"({profile.kp}):") + logger.debug(f" width = {width}") + logger.debug(f" strickler = {strickler}") + logger.debug(f" hieght = {height}") + + new = Data(reach=self._reach, status=self._status) + new["kp"] = profile.kp + new["discharge"] = discharge + new["elevation"] = elevation - self._data = self._new + previous_elevation = elevation + self._data.append(new)