src/Tools/MacMesh/MacMesh/CutnGroup.py

   1 # Copyright (C) 2007-2014  CEA/DEN, EDF R&D, OPEN CASCADE
   2
   3 # Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 # CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5
   6 # This library is free software; you can redistribute it and/or
   7 # modify it under the terms of the GNU Lesser General Public
   8 # License as published by the Free Software Foundation; either
   9 # version 2.1 of the License, or (at your option) any later version.
  10
  11 # This library is distributed in the hope that it will be useful,
  12 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 # Lesser General Public License for more details.
  15
  16 # You should have received a copy of the GNU Lesser General Public
  17 # License along with this library; if not, write to the Free Software
  18 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19
  20 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21
  22
  23
  24
  25 # This module allows cutting and grouping geometries by defining plane sections, with level of cutting as well as customizable Prefixes.
  26
  27 import math, Config
  28
  29 from salome.geom import geomBuilder
  30 geompy = geomBuilder.New( Config.theStudy )
  31
  32 def Go(GeoObj, CutPlnLst, OutLvlLst, PrefixLst, Publish):
  33
  34         """
  35         This function cuts any geometry (with infinite trim !) into several subgeometries that are cleanly saved inside the navigation tree. (Check GoTrim for the same functionality with custom trim size)
  36           - GeoObj is the geometrical object to be cut and grouped
  37           - CutPlnLst is a list of plane definitions. Each plane is a 6-tuple (contains 6 values). The first three represent the coordinates of the origin point and the second three represent the coordinates of the normal vector to the plane
  38                 Example 1: [(0,0,0,1,0,0)]: cut along a plane passing through the origin and normal to the x-axis
  39                 Example 2: [(0,0,0,1,0,0),(50,0,0,0,1,0)]: in addition to the first plane cut, cut through a plane passing by (50,0,0) and normal to the y axis.
  40                 Note that the plane size us determined automatically from the size of the geometry in question (using a very big trim size = 100 x length of geometry!)
  41           - OutLvlLst is a list containing integers that represent the inner sectioning level with respect to the original geometry type
  42                 A value of 1 means that the section will provide elements of one level lower than the original type. For exemple a solid sectioned at level 1 will produce faces. A Face sectioned at level 1 will produce edges.
  43                 A value of 2 means that a deeper sectioning will be applied. A solid sectioned with level 2 will give faces and edges. A face will give edges and vertices. An edge will give only vertices
  44                 The number of elements in this list should be (this is verified in the code) equal to the number of elements in the plane cut list. This is logical.
  45                 Example 1: [1]
  46                 Example 2: [1, 2], This means that the cut over the second plane will produce two types of elements unlike the first cut which will only output the first level objects.
  47           - PrefixLst is a list of strings that contains the naming Prefixes that are used by the script to generate the subshape names. This is very useful for relating the results to the sectioning requested.
  48                 Example 1: ['Entry']
  49                 Example 2: ['Entry','Exit']     The resulting groups from the sectioning with plane no.1 will then be saved as "Entry_FACE" and/or "Entry_EDGE" according to the original geometry object type and the cutting level
  50
  51         Imagine that we have a solid called ExampleSolid, an example command will be:
  52         CutnGroup.Go(ExampleSolid,[(0,0,0,1,0,0),(50,0,0,0,1,0)],[1, 2],['Entry','Exit'])
  53         """
  54
  55         NumCuts = CheckInput(CutPlnLst, OutLvlLst, PrefixLst, 1)
  56         OrigType = FindStandType(GeoObj,0)
  57         InvDictionary = dict((v,k) for k, v in geompy.ShapeType.iteritems())    # Give geometry type name as a function of standard type numbering, ex: 4=FACE, 6=EDGE, 7=VERTEX
  58         TrimSize = geompy.BasicProperties(GeoObj)[0]*100
  59         CutPlane = [] ; Sections = [] ; Parts = []
  60
  61         if NumCuts:
  62                 for i in range(0, NumCuts):             # Loop over the cutting planes to create them one by one
  63                         CutPlane.append(CreatePlane(CutPlnLst[i],TrimSize))
  64                 OutFather = geompy.MakePartition([GeoObj],CutPlane, [], [],FindStandType(GeoObj,1), 0, [], 0)   #Creating the partition object
  65                 if Publish: geompy.addToStudy(OutFather,'SectionedObject')
  66                 for i in range(0, NumCuts):
  67                         for j in range(OrigType+1+2, OrigType+1+2*(OutLvlLst[i]+1),2):
  68                                 if j == 8 : j = 7;      # Exception for the vertex case (=7)
  69                                 PossSubShapesID = geompy.SubShapeAllIDs(OutFather,j)    # List of subshape IDs than correspond to the required cutting level (section type : face/wire/vertex)
  70                                 PossSubShapes = geompy.ExtractShapes(OutFather,j)               # and the corresponding objects
  71                                 Accepted = []
  72                                 for k in range(0,len(PossSubShapesID)):         # Loop over all the subshapes checking if they belong to the cutting plane! if yes add them to current list
  73                                         if  IsOnPlane(PossSubShapes[k], CutPlnLst[i], 1e-7):
  74                                                 Accepted.append(PossSubShapesID[k])
  75                                 if Accepted :                                   # If some element is found, save it as a group with the prescribed Prefix
  76                                         dummyObj = geompy.CreateGroup(OutFather, j)
  77                                         geompy.UnionIDs(dummyObj, Accepted)
  78                                         Sections.append(dummyObj)
  79                                         if Publish:geompy.addToStudyInFather(OutFather, dummyObj, PrefixLst[i]+"_"+InvDictionary[j][0:2])
  80                                 else :
  81                                         print "Warning: For the section no.", i, ", No intersection of type " + InvDictionary[j] + " was found. Hence, no corresponding groups were created"
  82
  83                 SubShapesID = geompy.SubShapeAllIDs(OutFather,OrigType+1)               # Saving also the groups corresponding to the sectioned item of the same type: ex. A solid into n sub-solids due to the sections
  84                 for i in range(0,len(SubShapesID)):
  85                         dummyObj = geompy.CreateGroup(OutFather, OrigType+1)
  86                         geompy.UnionIDs(dummyObj, [SubShapesID[i]])
  87                         Parts.append(dummyObj)
  88                         if Publish: geompy.addToStudyInFather(OutFather, dummyObj, "SB"+"_"+InvDictionary[OrigType+1][0:3]+"_"+str(i+1))
  89
  90                 return OutFather, Sections, Parts
  91         else:
  92                 print("Fatal error, the routine cannot continue any further, check your input variables")
  93                 return -1
  94
  95 def GoTrim(GeoObj, CutPlnLst, OutLvlLst, PrefixLst, Publish):
  96
  97         """
  98         This function cuts any geometry into several subgeometries that are cleanly saved inside the navigation tree with a fully customizable trim size.
  99           - GeoObj is the geometrical object to be cut and grouped
 100           - CutPlnLst is a list of plane definitions. Each plane is a 7-tuple (contains 7 values). The first three represent the coordinates of the origin point and the second three represent the coordinates of the normal vector to the plane, the last value corresponds to the trim size of the planes
 101                 Example 1: [(0,0,0,1,0,0,5)]: cut along a plane passing through the origin and normal to the x-axis with a trim size of 5
 102                 Example 2: [(0,0,0,1,0,0,5),(50,0,0,0,1,0,10)]: in addition to the first plane cut, cut through a plane passing by (50,0,0) and normal to the y axis with a trim size of 10
 103           - OutLvlLst is a list containing integers that represent the inner sectioning level with respect to the original geometry type
 104                 A value of 1 means that the section will provide elements of one level lower than the original type. For exemple a solid sectioned at level 1 will produce faces. A Face sectioned at level 1 will produce edges.
 105                 A value of 2 means that a deeper sectioning will be applied. A solid sectioned with level 2 will give faces and edges. A face will give edges and vertices. An edge will give only vertices
 106                 The number of elements in this list should be (this is verified in the code) equal to the number of elements in the plane cut list. This is logical.
 107                 Example 1: [1]
 108                 Example 2: [1, 2], This means that the cut over the second plane will produce two types of elements unlike the first cut which will only output the first level objects.
 109           - PrefixLst is a list of strings that contains the naming Prefixes that are used by the script to generate the subshape names. This is very useful for relating the results to the sectioning requested.
 110                 Example 1: ['Entry']
 111                 Example 2: ['Entry','Exit']     The resulting groups from the sectioning with plane no.1 will then be saved as "Entry_FACE" and/or "Entry_EDGE" according to the original geometry object type and the cutting level
 112
 113         Imagine that we have a solid called ExampleSolid, an example command will be:
 114         CutnGroup.Go(ExampleSolid,[(0,0,0,1,0,0,5),(50,0,0,0,1,0,10)],[1, 2],['Entry','Exit'])
 115         """
 116
 117         NumCuts = CheckInput(CutPlnLst, OutLvlLst, PrefixLst, 0)
 118         OrigType = FindStandType(GeoObj,0)
 119         InvDictionary = dict((v,k) for k, v in geompy.ShapeType.iteritems())    # Give geometry type name as a function of standard type numbering, ex: 4=FACE, 6=EDGE, 7=VERTEX
 120         CutPlane = [] ; Sections = [] ; Parts = []
 121         if NumCuts:
 122                 for i in range(0, NumCuts):             # Loop over the cutting planes to create them one by one
 123                         CutPlane.append(CreatePlane(CutPlnLst[i][0:6],CutPlnLst[i][6]))
 124                 OutFather = geompy.MakePartition([GeoObj],CutPlane, [], [],FindStandType(GeoObj,1), 0, [], 0)   #Creating the partition object
 125                 if Publish: geompy.addToStudy(OutFather,'SectionedObject')
 126                 for i in range(0, NumCuts):
 127                         for j in range(OrigType+1+2, OrigType+1+2*(OutLvlLst[i]+1),2):
 128                                 if j == 8 : j = 7;      # Exception for the vertex case (=7)
 129                                 PossSubShapesID = geompy.SubShapeAllIDs(OutFather,j)    # List of subshape IDs than correspond to the required cutting level (section type : face/wire/vertex)
 130                                 PossSubShapes = geompy.ExtractShapes(OutFather,j)               # and the corresponding objects
 131                                 Accepted = []
 132                                 for k in range(0,len(PossSubShapesID)):         # Loop over all the subshapes checking if they belong to the cutting plane WITH THE TRIM SIZE CONDITION! if yes add them to current list
 133                                         if  IsOnPlane(PossSubShapes[k], CutPlnLst[i], 1e-7) and Distance2Pt(geompy.PointCoordinates(geompy.MakeCDG(PossSubShapes[k])),CutPlnLst[i][0:3])<=CutPlnLst[i][-1]:
 134                                                 Accepted.append(PossSubShapesID[k])
 135                                 if Accepted :                                   # If some element is found, save it as a group with the prescribed Prefix
 136                                         dummyObj = geompy.CreateGroup(OutFather, j)
 137                                         geompy.UnionIDs(dummyObj, Accepted)
 138                                         Sections.append(dummyObj)
 139                                         if Publish: geompy.addToStudyInFather(OutFather, dummyObj, PrefixLst[i]+"_"+InvDictionary[j][0:2])
 140                                 else :
 141                                         print "Warning: For the section no.", i, ", No intersection of type " + InvDictionary[j] + " was found. Hence, no corresponding groups were created"
 142
 143                 SubShapesID = geompy.SubShapeAllIDs(OutFather,OrigType+1)               # Saving also the groups corresponding to the sectioned item of the same type: ex. A solid into n sub-solids due to the sections
 144                 for i in range(0,len(SubShapesID)):
 145                         dummyObj = geompy.CreateGroup(OutFather, OrigType+1)
 146                         geompy.UnionIDs(dummyObj, [SubShapesID[i]])
 147                         Parts.append(dummyObj)
 148                         if Publish: geompy.addToStudyInFather(OutFather, dummyObj, "SB"+"_"+InvDictionary[OrigType+1][0:3]+"_"+str(i+1))
 149
 150                 return OutFather, Sections, Parts
 151         else:
 152                 print("Fatal error, the routine cannot continue any further, check your input variables")
 153                 return -1
 154 def FindStandType(GeoObj, method):
 155         """
 156         Find the standard index for the Geometrical object/compound type input, see dictionary in geompy.ShapeType
 157         """
 158         TopType = GeoObj.GetMaxShapeType().__str__()
 159         UnModType = geompy.ShapeType[TopType]
 160         if method == 0 :
 161                 StandType = UnModType-int(not(UnModType%2))             # So that wires and edges and considered the same, faces and shells, and so on
 162         else :
 163                 StandType = UnModType
 164
 165         return(StandType)
 166
 167 def CreatePlane(CutPlnVar,Trim):
 168         """
 169         Creates a temporary point and vector in salome in order to build the sectioning planes needed
 170         """
 171         Temp_Vtx = geompy.MakeVertex(CutPlnVar[0], CutPlnVar[1], CutPlnVar[2])
 172         Temp_Vec = geompy.MakeVectorDXDYDZ(CutPlnVar[3], CutPlnVar[4], CutPlnVar[5])
 173         CutPlane = geompy.MakePlane(Temp_Vtx, Temp_Vec, Trim)
 174         return(CutPlane)
 175
 176 def CheckInput(CutPlnLst, OutLvlLst, PrefixLst, AutoTrim):
 177         """
 178         Checks the user input specifically if all needed parameters are provided
 179         """
 180         if not ((len(CutPlnLst) == len(OutLvlLst)) and (len(CutPlnLst) == len(PrefixLst))):
 181                 print("Missing information about one or more of the cut planes")
 182                 return 0
 183         elif not ((len(CutPlnLst[0]) == 6+int(not AutoTrim))):
 184                 print("For each cutting plane you need to specify 6 parameters = 2 x 3 coordinates")
 185                 return 0
 186         else:
 187                 return len(CutPlnLst)
 188
 189 def IsOnPlane(GeoSubObj, CutPlnVar, tolerance):
 190         """
 191         Checks whether a geometry (vertex, segment, or face) belongs *completely* to the plane defined as a point and a normal vector
 192         """
 193         # lambda function that represents the plane equation, function = 0 <=> Pt defined with Coor belongs to plane
 194         PlaneEq = lambda Coor: CutPlnVar[3]*(Coor[0]-CutPlnVar[0])+CutPlnVar[4]*(Coor[1]-CutPlnVar[1])+CutPlnVar[5]*(Coor[2]-CutPlnVar[2])
 195         OrigType = FindStandType(GeoSubObj,0)
 196         if (OrigType >= 7):             # Vertex
 197                 NonTrimDecision = abs(PlaneEq(geompy.PointCoordinates(GeoSubObj))) < tolerance
 198                 if len(CutPlnVar) == 6 : return NonTrimDecision # No trim condition used
 199                 else : return (NonTrimDecision and Distance2Pt(CutPlnVar[0:3],geompy.PointCoordinates(GeoSubObj))<=CutPlnVar[6]/2)
 200         elif (OrigType >= 5):           # Line, decompose into two points then call recursively IsOnPlane function!
 201                 Verdict = True
 202                 for i in range(0,2):
 203                         Verdict = Verdict and IsOnPlane(geompy.GetVertexByIndex(GeoSubObj,i), CutPlnVar, tolerance)
 204                 return Verdict
 205         elif (OrigType >= 3):           # Face, decompose into three points then call recursively IsOnPlane function!
 206                 if IsOnPlane(geompy.MakeCDG(GeoSubObj),CutPlnVar, tolerance) : # Center of gravity belongs to plane, check if normal is parallel to plane
 207                         NormalP1Coor = geompy.PointCoordinates(geompy.GetVertexByIndex(geompy.GetNormal(GeoSubObj),0))
 208                         NormalP2Coor = geompy.PointCoordinates(geompy.GetVertexByIndex(geompy.GetNormal(GeoSubObj),1))
 209                         Normal = [NormalP1Coor[0]-NormalP2Coor[0], NormalP1Coor[1]-NormalP2Coor[1], NormalP1Coor[2]-NormalP2Coor[2]]
 210                         CrossP = CrossProd(CutPlnVar[3:6],Normal)               # Checks whether normals (of section plane and of face) are parallel or not
 211                         if (abs(CrossP[0])<tolerance and abs(CrossP[1])<tolerance and abs(CrossP[2])<tolerance):        # meaning zero cross product => parallel
 212                                 return True
 213                         else :
 214                                 return False
 215                 else :
 216                         return False
 217
 218
 219 def CrossProd(V1,V2):
 220         """
 221         Determines the cross product of two 3D vectors
 222         """
 223         return ([V1[1]*V2[2]-V1[2]*V2[1], V1[2]*V2[0]-V1[0]*V2[2], V1[0]*V2[1]-V1[1]*V2[0]])
 224
 225 def Distance2Pt(P1,P2):
 226         """
 227         Returns the distance between two points
 228         """
 229         return (math.sqrt((P1[0]-P2[0])**2+(P1[1]-P2[1])**2+(P1[2]-P2[2])**2))