6 from salome.smesh import smeshBuilder
10 import medcoupling as mc
17 def freeBordersGroup(meshFile):
18 print(" === freeBordersGroup", meshFile)
19 smesh = smeshBuilder.New()
20 smesh.SetEnablePublish( False ) # Set to False to avoid publish in study if not needed
21 ([MESH], status) = smesh.CreateMeshesFromMED(meshFile)
22 groups = MESH.GetGroups()
24 aCriterion = smesh.GetCriterion(SMESH.EDGE,SMESH.FT_FreeBorders,SMESH.FT_Undefined,0)
25 aCriteria.append(aCriterion)
26 aFilter = smesh.GetFilterFromCriteria(aCriteria)
27 aFilter.SetMesh(MESH.GetMesh())
28 FreeBorders = MESH.GroupOnFilter( SMESH.EDGE, 'FreeBorders', aFilter )
29 smesh.SetName(MESH, 'MESH')
30 newMeshName = '/tmp/freeBorders.med'
32 MESH.ExportMED(newMeshName,auto_groups=0,minor=41,overwrite=1,meshPart=None,autoDimension=1)
35 print('ExportMED() failed. Invalid file name?')
39 def explodeGroup(grp, grpName):
40 print(" === explodeGroup", grpName)
41 nbCells=grp.getNumberOfCells()
43 dicReverse = {} # id noeud --> id edges
44 for i in range(nbCells):
45 nodcell = grp.getNodeIdsOfCell(i)
46 for j in range(len(nodcell)):
47 if nodcell[j] in dicReverse:
48 dicReverse[nodcell[j]].append(i)
50 dicReverse[nodcell[j]] = [i]
53 usedCells = [False] * nbCells
54 while False in usedCells:
55 icell = usedCells.index(False)
56 usedCells[icell] = True
57 nodcell = grp.getNodeIdsOfCell(icell)
59 chain = [nodcell[0], nodcell[1]]
62 while nextnode in dicReverse:
63 nextcells = dicReverse[nextnode]
64 if len(nextcells) != 2: # end of chain(1) or "edges connector"(>2): stop
69 for i in range(len(nextcells)):
71 if not usedCells[ncell]: # the chain of nodes grows
72 usedCells[ncell] = True
74 nodcell = grp.getNodeIdsOfCell(ncell)
75 if nodcell[0] == nextnode:
76 nextnode = nodcell[1] # forward edge
78 nextnode = nodcell[0] # reversed edge ?
79 chain.append(nextnode)
80 if not newcell: # end of chain, closed
83 while prevnode in dicReverse:
84 prevcells = dicReverse[prevnode]
85 if len(prevcells) != 2: # end of chain(1) or "edges connector"(>2): stop
90 for i in range(len(prevcells)):
92 if not usedCells[ncell]: # the chain of nodes grows
93 usedCells[ncell] = True
95 nodcell = grp.getNodeIdsOfCell(ncell)
96 if nodcell[1] == prevnode:
97 prevnode = nodcell[0] # forward edge
99 prevnode = nodcell[1] # reversed edge ?
100 chain.insert(0, prevnode)
101 if not newcell: # end of chain, closed
105 chainDesc = (chain, grpName +"_%s" % len(nodeChains), closed)
106 nodeChains.append(chainDesc)
111 def writeShapeLines(mcMesh, grpName, nodeChains, offsetX=0., offsetY=0.):
112 print(" === writeShapeLines", grpName)
113 coords = mcMesh.getCoords()
114 w = shapefile.Writer(grpName)
117 for (chain, chainName, closed) in nodeChains:
118 print(" --- ", chainName)
121 coord = coords[node].getValues()
122 coordLb93=[coord[0] + offsetX, coord[1] + offsetY]
123 #print(" ", coordLb93)
124 chaincoords.append(coordLb93)
125 w.line([chaincoords])
130 def writeShapePoints(mcMesh, grpName, nodeChains, offsetX=0., offsetY=0.):
131 print(" === writeShapePoints", grpName)
132 coords = mcMesh.getCoords()
133 w = shapefile.Writer(grpName + '_pts')
136 for (chain, chainName, closed) in nodeChains:
137 print(" --- ", chainName)
140 coord = coords[node].getValues()
141 coordLb93=[coord[0] + offsetX, coord[1] + offsetY]
142 #print(" ", coordLb93)
143 chaincoords.append(coordLb93)
144 w.multipoint(chaincoords)
149 def exploreEdgeGroups(meshFile, offsetX=0., offsetY=0.):
150 print(" === exploreEdgeGroups", meshFile)
151 mcMesh = ml.MEDFileMesh.New(meshFile)
152 dim = mcMesh.getSpaceDimension()
153 d1=-1 # when dimension 2, edges are dim -1
154 if dim == 3: # when dimension 3, edges are dim -2
157 grp_names = mcMesh.getGroupsOnSpecifiedLev(d1) #names of edges groups
158 groups = [mcMesh.getGroup(d1, name) for name in grp_names]
159 for (grp, grpName) in zip(groups, grp_names):
160 nodeChains = explodeGroup(grp, grpName)
161 writeShapeLines(mcMesh, grpName, nodeChains, offsetX, offsetY)
162 writeShapePoints(mcMesh, grpName, nodeChains, offsetX, offsetY)
165 def fitShapePointsToMesh(freeBorderShapefile, shapefileToAdjust):
166 print(" === fitShapePointsToMesh", freeBorderShapefile, shapefileToAdjust)
167 fb = shapefile.Reader(freeBorderShapefile)
168 fbShapes = fb.shapes()
169 maxbbox=[1.e30, 1.e30, -1.e30, -1.e30]
171 for i,s in enumerate(fbShapes):
173 if (bb[0] < maxbbox[0] and bb[1] < maxbbox[1] and bb[2] > maxbbox[2] and bb[3] > maxbbox[3]):
176 fbs = fbShapes[outerBboxIndex]
178 sf = shapefile.Reader(shapefileToAdjust)
182 for i,p in enumerate(sfta.points):
185 for j,pfb in enumerate(fbs.points):
186 d = math.sqrt((pfb[0]-p[0])*(pfb[0]-p[0]) + (pfb[1]-p[1])*(pfb[1]-p[1]))
190 pdist.append((dmin, jmin, i)) # distance, index in freeBorder, index in shapeToAdjust
191 pdist.sort() # the 2 closest points must be set on the mesh freeBorder
193 i1 = min(pdist[0][2], pdist[1][2]) # index of first adjusted point in shapeToAdjust
194 i2 = max(pdist[0][2], pdist[1][2]) # index of second adjusted point in shapeToAdjust
195 if i1 == pdist[0][2]:
196 ifb1 = pdist[0][1] # index of first adjusted point in freeBorder
197 ifb2 = pdist[1][1] # index of second adjusted point in freeBorder
202 a = os.path.splitext(shapefileToAdjust)
203 shapefileAdjusted = a[0] + '_adj' + a[1]
204 chainName = os.path.basename(a[0])
206 w = shapefile.Writer(shapefileAdjusted)
211 chaincoords.append(fbs.points[ifb1])
212 for i in range(i1+1, i2):
213 chaincoords.append(sfta.points[i])
214 chaincoords.append(fbs.points[ifb2])
215 w.line([chaincoords])
216 w.record(chainName + '_0')
219 chaincoords.append(fbs.points[ifb2])
220 for i in range(i2+1, len(sfta.points)):
221 chaincoords.append(sfta.points[i])
223 chaincoords.append(sfta.points[i])
224 chaincoords.append(fbs.points[ifb1])
225 w.line([chaincoords])
226 w.record(chainName + '_1')