diff --git a/src/Model/Geometry/Profile.py b/src/Model/Geometry/Profile.py index 1a36164d3535760a79839164f1540d91651336c2..33741eaa63699730e10247b285223e55f6b6f1da 100644 --- a/src/Model/Geometry/Profile.py +++ b/src/Model/Geometry/Profile.py @@ -19,6 +19,7 @@ import logging from tools import timer +from shapely import geometry from Model.Geometry.Point import Point from Model.Except import NotImplementedMethodeError @@ -341,3 +342,12 @@ class Profile(object): # Abstract method for width approximation def get_water_limits(self, z): raise NotImplementedMethodeError(self, self.get_water_limits) + + def wet_points(self, z): + raise NotImplementedMethodeError(self, self.wet_point) + + def wet_perimeter(self, z): + raise NotImplementedMethodeError(self, self.wet_perimeter) + + def wet_area(self, z): + raise NotImplementedMethodeError(self, self.wet_area) diff --git a/src/Model/Geometry/ProfileXYZ.py b/src/Model/Geometry/ProfileXYZ.py index f3ebe88fb220b018332ac02ed7a1033b1629be8f..a5e24b99702701de4eec13a5b45bc39a8e0d2c1d 100644 --- a/src/Model/Geometry/ProfileXYZ.py +++ b/src/Model/Geometry/ProfileXYZ.py @@ -22,6 +22,7 @@ from typing import List from functools import reduce from tools import timer +from shapely import geometry from Model.Tools.PamhyrDB import SQLSubModel from Model.Except import ClipboardFormatError @@ -373,109 +374,60 @@ class ProfileXYZ(Profile, SQLSubModel): """ return [x for x in lst if not np.isnan(x)] - def _first_point_not_nan(self): - first_point = self.points[0] + def speed(self, q, z): + area = self.wet_area(z) - for point in self.points: - if not point.is_nan(): - first_point = point - break + if area == 0: + return 0 - return first_point + return q / area - def _last_point_not_nan(self): - last_point = self.points[-1] + def width_approximation(self): + if self.has_standard_named_points(): + rg = self.get_point_by_name("rg") + rd = self.get_point_by_name("rd") + else: + rg = self.points[0] + rd = self.points[-1] - for point in self.points[::-1]: - if not point.is_nan(): - last_point = point - break + return abs(rg.dist(rd)) - return last_point + def wet_perimeter(self, z): + poly = self.wet_polygon(z) - @timer - def get_station(self) -> np.ndarray: - """Projection of the points of the profile on a plane. + if poly is None: + return 0 - Args: - profile: The profile + return poly.length - Returns: - Projection of the points of the profile on a plane. - """ - if self.nb_points < 3: + def wet_area(self, z): + poly = self.wet_polygon(z) + + if poly is None: + return 0 + + return poly.area + + def wet_polygon(self, z): + points = self.wet_points(z) + if len(points) < 3: return None - else: - first_named_point = None - index_first_named_point = None - last_named_point = None - - first_point_not_nan = self._first_point_not_nan() - last_point_not_nan = self._last_point_not_nan() - - for index, point in enumerate(self.points): - if point.point_is_named(): - index_first_named_point = index - first_named_point = point - break - - for point in reversed(self.points): - if point.point_is_named(): - last_named_point = point - break - - station = [] - constant = 0.0 - - if (first_named_point is not None and - last_named_point is not None): - if (first_named_point != last_named_point and - first_named_point.x != last_named_point.x): - vector = Vector1d(first_named_point, last_named_point) - norm_dir_vec = vector.normalized_direction_vector() - else: - vector = Vector1d(first_point_not_nan, last_point_not_nan) - norm_dir_vec = vector.normalized_direction_vector() - - for point in self.points: - xi = point.x - first_named_point.x - yi = point.y - first_named_point.y - station_i = (norm_dir_vec[0] * xi + - norm_dir_vec[1] * yi) - station.append(station_i) - - constant = station[index_first_named_point] - elif first_named_point is None: - vector = Vector1d(first_point_not_nan, last_point_not_nan) - norm_dir_vec = vector.normalized_direction_vector() - for point in self.points: - xi = point.x - first_point_not_nan.x - yi = point.y - first_point_not_nan.y - station_i = (norm_dir_vec[0] * xi + - norm_dir_vec[1] * yi) - station.append(station_i) + z = map(lambda p: p.z, points) + station = self._get_station(points) - z_min = self.z_min() - index_profile_z_min = list( - filter( - lambda z: z[1] == z_min, - enumerate(self.z()) - ) - )[0] - constant = station[index_profile_z_min[0]] + poly = geometry.Polygon(list(zip(station, z))) - return list(map(lambda s: s - constant, station)) + return poly - def width_approximation(self): - if self.has_standard_named_points(): - rg = self.get_point_by_name("rg") - rd = self.get_point_by_name("rd") - else: - rg = self.points[0] - rd = self.points[-1] + def wet_points(self, z): + left, right = self.get_water_limits(z) - return abs(rg.dist(rd)) + points = list(filter(lambda p: p.z < z, self._points)) + points = [left] + points + [right] + points = sorted(points, key=lambda p: p.x) + + return points def get_water_limits(self, z): """ @@ -509,7 +461,7 @@ class ProfileXYZ(Profile, SQLSubModel): [self.point(i_left).z, self.point(i_left - 1).z], [self.point(i_left).y, self.point(i_left - 1).y] ) - pt_left = PointXYZ(x, y, z) + pt_left = PointXYZ(x, y, z, name="wl_left") else: pt_left = self.point(0) @@ -525,8 +477,105 @@ class ProfileXYZ(Profile, SQLSubModel): [self.point(i_right).z, self.point(i_right + 1).z], [self.point(i_right).y, self.point(i_right + 1).y] ) - pt_right = PointXYZ(x, y, z) + pt_right = PointXYZ(x, y, z, name="wl_right") else: pt_right = self.point(self.number_points - 1) return pt_left, pt_right + + def get_station(self): + """Projection of the points of the profile on a plane. + + Args: + self: The profile + + Returns: + Projection of the points of the profile on a plane. + """ + if self.nb_points < 3: + return None + else: + return self._get_station(self.points) + + @timer + def _get_station(self, points): + first_named_point = None + index_first_named_point = None + last_named_point = None + + first_point_not_nan = self._first_point_not_nan(points) + last_point_not_nan = self._last_point_not_nan(points) + + for index, point in enumerate(points): + if point.point_is_named(): + index_first_named_point = index + first_named_point = point + break + + for point in reversed(points): + if point.point_is_named(): + last_named_point = point + break + + station = [] + constant = 0.0 + + if (first_named_point is not None and + last_named_point is not None): + if (first_named_point != last_named_point and + first_named_point.x != last_named_point.x): + vector = Vector1d(first_named_point, last_named_point) + norm_dir_vec = vector.normalized_direction_vector() + else: + vector = Vector1d(first_point_not_nan, last_point_not_nan) + norm_dir_vec = vector.normalized_direction_vector() + + for point in points: + xi = point.x - first_named_point.x + yi = point.y - first_named_point.y + station_i = (norm_dir_vec[0] * xi + + norm_dir_vec[1] * yi) + station.append(station_i) + + constant = station[index_first_named_point] + elif first_named_point is None: + vector = Vector1d(first_point_not_nan, last_point_not_nan) + norm_dir_vec = vector.normalized_direction_vector() + + for point in points: + xi = point.x - first_point_not_nan.x + yi = point.y - first_point_not_nan.y + station_i = (norm_dir_vec[0] * xi + + norm_dir_vec[1] * yi) + station.append(station_i) + + z_min = self.z_min() + index_profile_z_min = list( + filter( + lambda z: z[1] == z_min, + enumerate(self.z()) + ) + )[0] + constant = station[index_profile_z_min[0]] + + return list(map(lambda s: s - constant, station)) + + def _first_point_not_nan(self, points): + first_point = None + + for point in points: + if not point.is_nan(): + first_point = point + break + + return first_point + + def _last_point_not_nan(self, points): + last_point = None + + for point in reversed(points): + if not point.is_nan(): + last_point = point + break + + return last_point