diff --git a/src/Meshing/Internal.py b/src/Meshing/Internal.py new file mode 100644 index 00000000..3cdacf0a --- /dev/null +++ b/src/Meshing/Internal.py @@ -0,0 +1,326 @@ +# Internal.py -- Pamhyr +# Copyright (C) 2023-2025 INRAE +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . + +# -*- coding: utf-8 -*- + +import logging +from numpy import sign +from copy import deepcopy + +from tools import logger_color_red, logger_color_reset, logger_exception +from Meshing.AMeshingTool import AMeshingTool + +from Model.Geometry.ProfileXYZ import ProfileXYZ + +logger = logging.getLogger() + + +class InternalMeshing(AMeshingTool): + def __init__(self): + super(InternalMeshing, self).__init__() + + ########### + # Meshing # + ########### + + def meshing(self, reach, + step: float = 50, + limites=[-1, -1], + linear: bool = True): + if reach is None or len(reach.profiles) == 0: + return reach + + if limites[0] > limites[1]: + l = limites[1] + limites[1] = limites[0] + limites[0] = l + elif limites[0] == limites[1]: + return reach + + self.st_to_m(reach, limites) + self.interpolate_transversal_step(reach, + limites, + step) + return reach + + def st_to_m(self, reach, limites): + + gl = reach.compute_guidelines() + profiles = reach.profiles[limites[0]:limites[1]+1] + + print(profiles) + + # we make sure that the lines are in the left-to-right order + guide_list = [ + x.name + for x in profiles[0].named_points() + if x.name in gl[0] + ] + + gl1 = [""] + guide_list + gl2 = guide_list + [""] + max_values = [0] * (len(guide_list) + 1) + max_values_index = [0] * (len(guide_list) + 1) + # we work between each couple of guidelines + for j, p in enumerate(profiles): + index1 = 0 + index2 = p.named_point_index(guide_list[0]) + if index2 - index1 > max_values[0]: + max_values[0] = index2 - index1 + max_values_index[0] = j + for i in range(len(guide_list) - 1): + index1 = p.named_point_index(guide_list[i]) + index2 = p.named_point_index(guide_list[i+1]) + if index2 - index1 > max_values[i + 1]: + max_values[i + 1] = index2 - index1 + max_values_index[i + 1] = j + index1 = p.named_point_index(guide_list[-1]) + index2 = len(p) - 1 + if index2 - index1 > max_values[-1]: + max_values[-1] = index2 - index1 + max_values_index[-1] = j + + for i in range(len(max_values)): + for isect in range(max_values_index[i], len(profiles)-1): + self.compl_sect(profiles[isect], + profiles[isect+1], + gl1[i], gl2[i]) + for isect in reversed(range(0, max_values_index[i])): + self.compl_sect(profiles[isect+1], + profiles[isect], + gl1[i], gl2[i]) + + def interpolate_transversal_step(self, + reach, + limites, + step): + # calcul number of intermediate profiles + np = [] + for i in range(limites[0], limites[1]): + np.append(int((reach.profiles[i+1].rk - reach.profiles[i].rk)/step) -1) + if np[-1] < 0: + np[-1] = 0 + + + d = [] # ratios + p = [] # new profiles + ptr = int(limites[0]) + for i in range(limites[1] - limites[0]): + d = 1.0/float(np[i]+1) + ptr0 = ptr + for j in range(np[i]): + #p = reach.profiles[ptr0].copy() + p = reach.profiles[ptr0].copy() + # RATIO entre les deux sections initiales + dj=float(j+1)*d + ptr += 1 # next profile, original + for k in range(len(reach.profiles[ptr0].points)): + p.points[k].x = reach.profiles[ptr0].points[k].x + \ + dj*(reach.profiles[ptr].points[k].x - \ + reach.profiles[ptr0].points[k].x) + p.points[k].y = reach.profiles[ptr0].points[k].y + \ + dj*(reach.profiles[ptr].points[k].y - \ + reach.profiles[ptr0].points[k].y) + p.points[k].z = reach.profiles[ptr0].points[k].z + \ + dj*(reach.profiles[ptr].points[k].z - \ + reach.profiles[ptr0].points[k].z) + p.rk = reach.profiles[ptr0].rk + \ + dj*(reach.profiles[ptr].rk - \ + reach.profiles[ptr0].rk) + p.name = f'interpol{p.rk}' + reach.insert_profile(ptr, p) + ptr += 1 # next profile, original + + def compl_sect(self, sect1, sect2, tag1, tag2): + # sect1: reference section + # sect2: section to complete + # start1 and 2: indices of the first point of the bed to complete + # end = start + len + # len1 and 2: number of intervals + # len1: target len + + if tag1 == '': # left end point + start1 = 0 + start2 = 0 + else: + start1 = sect1.named_point_index(tag1) + start2 = sect2.named_point_index(tag1) + + if tag2 == '': # right end point + end1 = sect1.nb_points-1 + end2 = sect2.nb_points-1 + else: + end1 = sect1.named_point_index(tag2) + end2 = sect2.named_point_index(tag2) + + len1 = end1 - start1 + len2 = end2 - start2 + + if len1 < len2: + self.compl_sect(sect2, sect1, start2, start1, + len2, len1, tag1, tag2) + return + elif len1 == len2: + return + + ltot1 = 0.0 + ltot2 = 0.0 + alpha = [] + beta = [] + # remove names + sect2.points[start2].name = '' + sect2.points[end2].name = '' + for i in range(len1): + ltot1 += sect1.point(start1+i).dist(sect1.point(start1+i+1)) + alpha.append(ltot1) + alpha = list(map(lambda x: x/ltot1, alpha)) # target ratios + for i in range(len2): + ltot2 += sect2.point(start2+i).dist(sect2.point(start2+i+1)) + beta.append(ltot2) + beta = list(map(lambda x: x/ltot2, beta)) # current ratios + for i in range(len1 - len2): + beta.append(1.0) + beta.append(0.0) # beta[-1] + if len2 < 1 or ltot2 < 0.0001: + # unique point (copy len1 times) + # we are at an edge of the profile + for i in range(len1-len2): + p = sect2.point(start2).copy() + sect2.insert_point(start2, p) + if tag1 == '': # left end point + sect2.point(start2).name = '' + if tag2 == '': # left end point + sect2.point(start2+1).name = '' + elif ltot1 < 0.0001: + sect2.add_npoints(len1-len2) + else: # regular case + # loop on the points to insert + for k in range(len1-len2): + trouve = False + for j0 in range(len2): + if trouve: + break + undeplus = True + for j1 in range(j0, len2): + if undeplus: + if beta[j1] <= alpha[j1]: + undeplus = False + break + if undeplus: + if j0 == len2: + trouve = True + else: + j1 = j0+1 + if beta[j0]-alpha[j0] > alpha[j1]-beta[j0]: + trouve = True + + if trouve: + len2 += 1 + for i in reversed(range(j0, len2)): + beta[i+1] = beta[i] + if (beta[j0] == beta[j0-1]): + # rien + pass + else: + ratio = (alpha[j0] - beta[j0-1]) \ + / (beta[j0] - beta[j0-1]) + if ratio < 0.0: + print(f"ratio négatif {ratio}") + # on double le point a gauche + p = sect2.point(start2+j0-1).copy() + sect2.insert_point(start2+j0-1, p) + sect2.point(start2+j0-1).name = '' + beta[j0] = beta[j0-1] + else: + p = sect2.point(start2+j0).copy() + sect2.insert_point(start2+j0, p) + sect2.point(start2+j0+1).name = '' + sect2.point(start2+j0+1).x = \ + (1.0-ratio) * sect2.point(start2+j0).x + \ + ratio * sect2.point(start2+j0+2).x + sect2.point(start2+j0+1).y = \ + (1.0-ratio) * sect2.point(start2+j0).y + \ + ratio * sect2.point(start2+j0+2).y + sect2.point(start2+j0+1).z = \ + (1.0-ratio) * sect2.point(start2+j0).z + \ + ratio * sect2.point(start2+j0+2).z + beta[j0] = (1.-ratio)*beta[j0-1] \ + + ratio*beta[j0+1] + + if len1 > len2: + for i in range(len1 - len2): + p = sect2.point(start2+len2).copy() + sect2.insert_point(start2+len2, p) + sect2.point(start2+len2).name = '' + sect2.point(start2).name = tag1 + sect2.point(start2+len2).name = tag2 + sect2.modified() + + def update_rk(self, reach, + step: float = 50, + limites=[-1, -1], + origin=0, + directrices=['un', 'np'], + lplan: bool = False, + lm: int = 3, + linear: bool = False, + origin_value=0.0, + orientation=0): + if reach is None or len(reach.profiles) == 0: + return reach + + nprof = reach.number_profiles + if origin >= nprof or origin < 0: + logger.warning( + f"Origin outside reach" + ) + return reach + + # orientation is the orintation of the bief: + # 1: up -> downstream + # 2: down -> upstream + # else: keep current orientation + if orientation == 1: + sgn = 1.0 + elif orientation == 2: + sgn = -1.0 + else: + sgn = sign(reach.profiles[-1].rk - reach.profiles[0].rk) + + reach.profiles[origin].rk = origin_value + + if origin < nprof - 1: + for i in range(origin+1, nprof): + # 2D + d1 = reach.profiles[i-1].named_point( + directrices[0] + ).dist_2d(reach.profiles[i].named_point(directrices[0])) + d2 = reach.profiles[i-1].named_point( + directrices[1] + ).dist_2d(reach.profiles[i].named_point(directrices[1])) + reach.profiles[i].rk = reach.profiles[i-1].rk \ + + (sgn * (d1 + d2) / 2) + if origin > 0: + for i in reversed(range(0, origin)): + # 2D + d1 = reach.profiles[i+1].named_point( + directrices[0] + ).dist_2d(reach.profiles[i].named_point(directrices[0])) + d2 = reach.profiles[i+1].named_point( + directrices[1] + ).dist_2d(reach.profiles[i].named_point(directrices[1])) + reach.profiles[i].rk = reach.profiles[i+1].rk \ + - (sgn * (d1 + d2) / 2) diff --git a/src/Meshing/Mage.py b/src/Meshing/Mage.py index 2d28fb04..0dd271ff 100644 --- a/src/Meshing/Mage.py +++ b/src/Meshing/Mage.py @@ -19,6 +19,7 @@ import os import logging import tempfile +from numpy import sign from ctypes import ( cdll, diff --git a/src/Model/Geometry/PointXYZ.py b/src/Model/Geometry/PointXYZ.py index 9ba79e66..4d16f08f 100644 --- a/src/Model/Geometry/PointXYZ.py +++ b/src/Model/Geometry/PointXYZ.py @@ -228,3 +228,13 @@ class PointXYZ(Point, SQLSubModel): c = PointXYZ.distance(p3, p1) s = (a + b + c) / 2 return (s*(s-a) * (s-b)*(s-c)) ** 0.5 + + def copy(self): + p = PointXYZ(self.x, + self.y, + self.z, + self.name, + self._profile, + self._status) + p.sl = self.sl + return p diff --git a/src/Model/Geometry/Profile.py b/src/Model/Geometry/Profile.py index 6a4b45ef..7f9d5373 100644 --- a/src/Model/Geometry/Profile.py +++ b/src/Model/Geometry/Profile.py @@ -88,6 +88,14 @@ class Profile(object): def point(self, index): return self.points[index] + def named_point(self, name): + return next((p for p in self.points if p.name == name), None) + + def named_point_index(self, name): + return next( + (p for p in enumerate(self.points) if p[1].name == name), None + )[0] + @property def reach(self): return self._reach @@ -206,7 +214,7 @@ class Profile(object): """Insert point at index. Args: - index: The index of new profile. + index: The index of new point. point: The point. Returns: diff --git a/src/Model/Geometry/ProfileXYZ.py b/src/Model/Geometry/ProfileXYZ.py index 6a4f196b..d0f44ea9 100644 --- a/src/Model/Geometry/ProfileXYZ.py +++ b/src/Model/Geometry/ProfileXYZ.py @@ -882,3 +882,33 @@ class ProfileXYZ(Profile, SQLSubModel): def modified(self): self.tab_up_to_date = False self.station_up_to_date = False + + def add_npoints(self, npoints): + # add npoints in a profile + for k in range(npoints): + disti = 0.0 + i = 1 + for j in range(len(self.points)-1): + distj = self.point(j).dist(self.point(j+1)) + if distj > disti: + disti = distj + i = j + self.insert_point(i, self.point(i).copy()) + self.point(i+1).name = '' + self.point(i+1).x = 0.5 * self.point(i).x + 0.5 * self.point(i+2).x + self.point(i+1).y = 0.5 * self.point(i).y + 0.5 * self.point(i+2).y + self.point(i+1).z = 0.5 * self.point(i).z + 0.5 * self.point(i+2).z + + def copy(self): + p=ProfileXYZ(self.id, + self.name, + self.rk, + self.reach, + self.num, + 0, + 0, 0, + self._status) + for i, k in enumerate(self.points): + p.insert_point(i, k.copy()) + + return p diff --git a/src/View/Geometry/MeshingDialog.py b/src/View/Geometry/MeshingDialog.py index d33f66b6..e5928472 100644 --- a/src/View/Geometry/MeshingDialog.py +++ b/src/View/Geometry/MeshingDialog.py @@ -23,12 +23,12 @@ from PyQt5.QtGui import ( ) from PyQt5.QtCore import ( - Qt, QVariant, QAbstractTableModel, + Qt, QVariant, QAbstractTableModel, QItemSelectionModel, ) from PyQt5.QtWidgets import ( QDialogButtonBox, QComboBox, QUndoStack, QShortcut, - QDoubleSpinBox, + QDoubleSpinBox, QRadioButton, ) @@ -71,7 +71,7 @@ class MeshingDialog(PamhyrDialog): self._origin = self._reach.profile(0) self._init_default_values_profiles() - self._init_default_values_guidelines() + # self._init_default_values_guidelines() self.set_double_spin_box( "doubleSpinBox_space_step", @@ -83,19 +83,22 @@ class MeshingDialog(PamhyrDialog): else: self.set_radio_button("radioButton_spline", True) + self.find(QRadioButton, "radioButton_spline").setEnabled(False) + self.find(QRadioButton, "radioButton_linear").setEnabled(False) + def _init_default_values_profiles(self): profiles = self.profiles self.combobox_add_items("comboBox_begin_rk", profiles) self.combobox_add_items("comboBox_end_rk", profiles) - if len(self.parent.tableView.selectedIndexes()) <= 1: + r = self.parent.tableView\ + .selectionModel()\ + .selectedRows() + if len(r) <= 1: self.set_combobox_text("comboBox_begin_rk", profiles[0]) self.set_combobox_text("comboBox_end_rk", profiles[-1]) else: - r = self.parent.tableView\ - .selectionModel()\ - .selectedRows() self.set_combobox_text("comboBox_begin_rk", profiles[r[0].row()]) self.set_combobox_text("comboBox_end_rk", profiles[r[-1].row()]) @@ -174,8 +177,10 @@ class MeshingDialog(PamhyrDialog): self._begin_cs = self.profiles.index(p1) self._end_cs = self.profiles.index(p2) - self._begin_dir = self.get_combobox_text("comboBox_begin_gl") - self._end_dir = self.get_combobox_text("comboBox_end_gl") + # self._begin_dir = self.get_combobox_text("comboBox_begin_gl") + # self._end_dir = self.get_combobox_text("comboBox_end_gl") + + self.parent.tableView.selectionModel().clearSelection() super().accept() diff --git a/src/View/Geometry/UndoCommand.py b/src/View/Geometry/UndoCommand.py index f031c008..e60ad848 100644 --- a/src/View/Geometry/UndoCommand.py +++ b/src/View/Geometry/UndoCommand.py @@ -28,8 +28,6 @@ from PyQt5.QtWidgets import ( from Model.Geometry import Reach from Model.Except import exception_message_box -from Meshing.Mage import MeshingWithMage - logger = logging.getLogger() @@ -240,7 +238,7 @@ class MeshingCommand(QUndoCommand): self._mesher = mesher self._command = command - self._profiles = reach.profiles.copy() + self._profiles = [p.copy() for p in reach.profiles] self._profiles.reverse() self._new_profiles = None @@ -264,7 +262,7 @@ class MeshingCommand(QUndoCommand): **self._data ) - self._new_profiles = self._reach.profiles.copy() + self._new_profiles = [p.copy() for p in self._reach.profiles] self._new_profiles.reverse() else: self._reach.purge() diff --git a/src/View/Geometry/Window.py b/src/View/Geometry/Window.py index 4356755b..237f70ed 100644 --- a/src/View/Geometry/Window.py +++ b/src/View/Geometry/Window.py @@ -49,9 +49,7 @@ from View.Tools.Plot.PamhyrCanvas import MplCanvas from View.SelectReach.Window import SelectReachWindow -from Meshing.Mage import ( - MeshingWithMage, MeshingWithMageMailleurTT -) +from Meshing.Internal import InternalMeshing from View.Geometry.Table import GeometryReachTableModel from View.Geometry.PlotXY import PlotXY @@ -298,8 +296,6 @@ class GeometryWindow(PamhyrWindow): data = { "step": dlg.space_step, "limites": [dlg.begin_cs, dlg.end_cs], - "directrices": [dlg.begin_dir, dlg.end_dir], - "lplan": dlg.lplan, "linear": dlg.linear, } self._edit_meshing(data) @@ -309,16 +305,10 @@ class GeometryWindow(PamhyrWindow): def _edit_meshing(self, data): try: - mesher = MeshingWithMageMailleurTT() + mesher = InternalMeshing() self._table.meshing(mesher, data) except Exception as e: logger_exception(e) - raise ExternFileMissingError( - module="mage", - filename="MailleurTT", - path=MeshingWithMageMailleurTT._path(), - src_except=e - ) def update_rk(self): try: diff --git a/src/View/ui/MeshingOptions.ui b/src/View/ui/MeshingOptions.ui index de48342a..b7eab409 100644 --- a/src/View/ui/MeshingOptions.ui +++ b/src/View/ui/MeshingOptions.ui @@ -7,24 +7,14 @@ 0 0 520 - 341 + 245 Dialog - - - - Qt::Horizontal - - - QDialogButtonBox::Cancel|QDialogButtonBox::Ok - - - - + Parameters @@ -79,7 +69,7 @@ Spline - true + false @@ -88,6 +78,9 @@ Linear + + true + @@ -108,47 +101,14 @@ - - - - Guideline used for distance computation + + + + Qt::Horizontal + + + QDialogButtonBox::Cancel|QDialogButtonBox::Ok - - - - - true - - - Second guideline - - - - - - - true - - - First guideline - - - - - - - true - - - - - - - true - - - -