src/Tools/t_shape_builder.py

   1 # -*- coding: iso-8859-1 -*-
   2
   3 import sys
   4 import salome
   5
   6 import GEOM
   7 from salome.geom import geomBuilder
   8 import math
   9 import SALOMEDS
  10
  11
  12 def demidisk(study, r1, a1, roty=0, solid_thickness=0):
  13   if solid_thickness < 1e-7:
  14     with_solid = False
  15   else:
  16     with_solid = True
  17
  18   geompy = geomBuilder.New(study)
  19
  20   O = geompy.MakeVertex(0, 0, 0)
  21   OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
  22   OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
  23   OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)
  24
  25   v=range(8)
  26   l=range(8)
  27   v0 = geompy.MakeVertex(0, 0, 0)
  28   v[0] = geompy.MakeVertex(0, r1/2.0, 0)
  29   v[1] = geompy.MakeVertex(0, r1, 0)
  30   l[1] = geompy.MakeLineTwoPnt(v[0], v[1])
  31   l[2] = geompy.MakeRotation(l[1], OX, a1*math.pi/180.0)
  32   v[4] = geompy.MakeRotation(v[0], OX, a1*math.pi/180.0)
  33   v[6] = geompy.MakeRotation(v[1], OX, a1*math.pi/180.0)
  34
  35   v[2] = geompy.MakeVertex(0, -r1/2.0, 0)
  36   v[3] = geompy.MakeVertex(0, -r1, 0)
  37   l[3] = geompy.MakeLineTwoPnt(v[2], v[3])
  38   l[4] = geompy.MakeRotation(l[3], OX, -a1*math.pi/180.0)
  39   v[5] = geompy.MakeRotation(v[2], OX, -a1*math.pi/180.0)
  40   v[7] = geompy.MakeRotation(v[3], OX, -a1*math.pi/180.0)
  41
  42   l[5] = geompy.MakeLineTwoPnt(v[4], v[5])
  43   l[6] = geompy.MakeLineTwoPnt(v[0], v[4])
  44   l[7] = geompy.MakeLineTwoPnt(v[2], v[5])
  45
  46   v7 = geompy.MakeVertex(0, 0, r1)
  47   arc1 = geompy.MakeArc(v[1], v7, v[3])
  48   l[0] = geompy.MakeLineTwoPnt(v[1], v[3])
  49   face1 = geompy.MakeFaceWires([arc1, l[0]], 1)
  50
  51   if with_solid:
  52     # Vertices
  53     v0 = geompy.MakeVertex(0, r1 + solid_thickness, 0)
  54     v1 = geompy.MakeRotation(v0, OX, a1*math.pi/180.0)
  55     v2 = geompy.MakeRotation(v0, OX, math.pi - (a1*math.pi/180.0))
  56     v3 = geompy.MakeRotation(v0, OX, math.pi)
  57     v.extend([v0,v1,v3,v2]) # The order is important for use in pointsProjetes
  58     l0 = geompy.MakeLineTwoPnt(v[1], v0)
  59     l2 = geompy.MakeRotation(l0, OX, a1*math.pi/180.0)
  60     l3 = geompy.MakeRotation(l0, OX, math.pi - (a1*math.pi/180.0))
  61     face2 = geompy.MakeRevolution(l0, OX, a1*math.pi/180.0)
  62     face3 = geompy.MakeRevolution(l2, OX, math.pi - 2*a1*math.pi/180.0)
  63     face4 = geompy.MakeRevolution(l3, OX, a1*math.pi/180.0)
  64     part0 = geompy.MakePartition([face1], [l[2], l[4], l[5], l[6], l[7]], [], [], geompy.ShapeType["FACE"], 0, [], 0, True)
  65     compound1 = geompy.MakeCompound([part0, face2, face3, face4])
  66     part1 = geompy.MakeGlueEdges(compound1,1e-7)
  67   else:
  68     part1 = geompy.MakePartition([face1], [l[2], l[4], l[5], l[6], l[7]], [], [], geompy.ShapeType["FACE"], 0, [], 0, True)
  69
  70   if roty != 0:
  71     vrot = [ geompy.MakeRotation(vert, OY, roty*math.pi/180.0) for vert in v ]
  72     lrot = [ geompy.MakeRotation(lin, OY, roty*math.pi/180.0) for lin in l ]
  73     arc = geompy.MakeRotation(arc1, OY, roty*math.pi/180.0)
  74     part = geompy.MakeRotation(part1, OY, roty*math.pi/180.0)
  75     return vrot, lrot, arc, part
  76   else:
  77     return v, l, arc1, part1
  78
  79 def pointsProjetes(study, vref, face):
  80   geompy = geomBuilder.New(study)
  81   vface = geompy.ExtractShapes(face, geompy.ShapeType["VERTEX"], True)
  82   vord = range(len(vref))
  83   plan = geompy.MakePlaneThreePnt(vref[0], vref[1], vref[-1], 10000)
  84   vproj = [ geompy.MakeProjection(vert, plan) for vert in vface ]
  85   for i,v in enumerate(vproj):
  86     dist = [ (geompy.MinDistance(v, vr), j) for j,vr in enumerate(vref) ]
  87     dist.sort()
  88     #print dist
  89     if dist[0][0] < 1.e-3:
  90       vord[dist[0][1]] = vface[i]
  91   return vord
  92
  93 def arcsProjetes(study, vf, face):
  94   geompy = geomBuilder.New(study)
  95   lface = geompy.ExtractShapes(face, geompy.ShapeType["EDGE"], True)
  96   lord = range(3)
  97   ends = [vf[1], vf[6], vf[7], vf[3]]
  98   for i in range(3):
  99     for lf in lface:
 100       pts = geompy.ExtractShapes(lf, geompy.ShapeType["VERTEX"], True)
 101       if (((geompy.MinDistance(pts[0], ends[i]) < 0.001) and (geompy.MinDistance(pts[1], ends[i+1]) < 0.001)) or
 102           ((geompy.MinDistance(pts[1], ends[i]) < 0.001) and (geompy.MinDistance(pts[0], ends[i+1]) < 0.001))):
 103         lord[i] = lf
 104         print "arc_%d OK"%i
 105         break
 106     pass
 107   return lord
 108
 109 def build_shape(study, r1, r2, h1, h2, solid_thickness=0):
 110   if solid_thickness < 1e-7:
 111     with_solid = False
 112   else:
 113     with_solid = True
 114
 115   geompy = geomBuilder.New(study)
 116
 117   O = geompy.MakeVertex(0, 0, 0)
 118   OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
 119   OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
 120   OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)
 121
 122   a1 = 45.0
 123   seuilmax = 0.1
 124   ratio = float(r2)/float(r1)
 125   if ratio > (1.0 -seuilmax):
 126     a1 = 45.0*(1.0 -ratio)/seuilmax
 127
 128   """
 129   res = geompy.MakeCompound([demicyl1,demicyl2])
 130   return res
 131   """
 132
 133   # --- creation des faces de la jonction
 134   [faci, sect45, arc1, l1, lord90, lord45, edges, arcextru] = jonction(study, r1, r2,\
 135                                                                        h1, h2, a1)
 136   if with_solid:
 137     [faci_ext, sect45_ext, arc1_ext, l1_ext, \
 138      lord90_ext, lord45_ext, edges_ext, arcextru_ext] = jonction(study, r1 + solid_thickness, r2 + solid_thickness,\
 139                                                                  h1, h2, a1)
 140     faces_jonction_ext = []
 141     for i,l in enumerate(lord90):
 142       faces_jonction_ext.append(geompy.MakeQuad2Edges(lord90[i],lord90_ext[i]))
 143     for i in [1, 3, 6, 7]:
 144       faces_jonction_ext.append(geompy.MakeQuad2Edges(edges[i],edges_ext[i]))
 145     for i,l in enumerate(lord45):
 146       faces_jonction_ext.append(geompy.MakeQuad2Edges(lord45[i],lord45_ext[i]))
 147
 148     for i,face in enumerate(faces_jonction_ext):
 149       geompy.addToStudy(faces_jonction_ext[i], "faci_ext_%d"%i)
 150
 151   # --- extrusion droite des faces de jonction, pour reconstituer les demi cylindres
 152   # TODO : ajouter les faces nécessaires à sect45 dans le cas avec solide
 153   if with_solid:
 154     sect45 = geompy.MakeCompound([sect45]+faces_jonction_ext[-3:])
 155     sect45 = geompy.MakeGlueEdges(sect45, 1e-7)
 156
 157   #return sect45, faces_jonction_ext[-3:]
 158   extru1 = geompy.MakePrismVecH(sect45, OX, h1+10)
 159
 160   #base2 = geompy.MakeCompound(faci[5:])
 161   #base2 = geompy.MakeGlueEdges(base2, 1e-7)
 162   # RNC : perf
 163   faces_coupe = faci[5:]
 164   if with_solid:
 165     faces_coupe = faci[5:]+faces_jonction_ext[:3]
 166   base2 = geompy.MakePartition(faces_coupe, [], [], [], geompy.ShapeType["FACE"], 0, [], 0, True)
 167   extru2 = geompy.MakePrismVecH(base2, OZ, h2)
 168
 169   # --- partition et coupe
 170
 171   if with_solid:
 172      demiDisque = geompy.MakeFaceWires([arc1_ext, l1_ext[0]], 1)
 173   else:
 174      demiDisque = geompy.MakeFaceWires([arc1, l1[0]], 1)
 175   demiCylindre = geompy.MakePrismVecH(demiDisque, OX, h1)
 176
 177   box = geompy.MakeBox(0, -2*(r1+h1), -2*(r1+h1), 2*(r1+h1), 2*(r1+h1), 2*(r1+h1))
 178   rot = geompy.MakeRotation(box, OY, 45*math.pi/180.0)
 179
 180   garder = geompy.MakeCutList(demiCylindre, [extru2, rot], True)
 181   geompy.addToStudy(garder,"garder")
 182
 183   faces_coupe = faci[:5]
 184   if with_solid:
 185     faces_coupe.extend(faces_jonction_ext[-7:])
 186   raccord = geompy.MakePartition([garder], faces_coupe + [arcextru], [], [], geompy.ShapeType["SOLID"], 0, [], 0, True)
 187   assemblage = geompy.MakeCompound([raccord, extru1, extru2])
 188   assemblage = geompy.MakeGlueFaces(assemblage, 1e-7)
 189   # RNC : perf
 190   #assemblage = geompy.MakePartition([raccord, extru1, extru2], [], [], [], geompy.ShapeType["SOLID"], 0, [], 0, True)
 191
 192   #return extru2, garder, raccord
 193
 194   box = geompy.MakeBox(-1, -(r1+r2+2*solid_thickness), -1, h1, r1+r2+2*solid_thickness, h2)
 195   geompy.addToStudy(box, "box")
 196   final = geompy.MakeCommonList([box, assemblage], True)
 197
 198   # --- Partie inférieure
 199   v3, l3, arc3, part3 = demidisk(study, r1, a1, 180.0, solid_thickness)
 200   geompy.addToStudy(part3,"part3")
 201   extru3 = geompy.MakePrismVecH(part3, OX, h1)
 202   geompy.addToStudy(extru3,"extru3")
 203
 204   # --- Symétrie
 205
 206   compound = geompy.MakeCompound([final, extru3])
 207   plane = geompy.MakePlane(O,OX,2000)
 208   compound_mirrored = geompy.MakeMirrorByPlane(compound, plane)
 209   final = geompy.MakeCompound([compound, compound_mirrored])
 210
 211   return final
 212
 213
 214 def jonction(study, r1, r2, h1, h2, a1):
 215
 216   O = geompy.MakeVertex(0, 0, 0)
 217   OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
 218   OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
 219   OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)
 220
 221   # --- sections droites des deux demi cylindres avec le partionnement
 222   v1, l1, arc1, part1 = demidisk(study, r1, a1, 0.)
 223   v2, l2, arc2, part2 = demidisk(study, r2, a1, 90.0)
 224   #elems_disk1 = [v1, l1, arc1, part1]
 225   #elems_disk2 = [v2, l2, arc2, part2]
 226
 227   # --- extrusion des sections --> demi cylindres de travail, pour en extraire les sections utilisées au niveau du Té
 228   #     et enveloppe cylindrique du cylindre principal
 229
 230   demicyl1 = geompy.MakePrismVecH(part1, OX, h1)
 231   demicyl2 = geompy.MakePrismVecH(part2, OZ, h2)
 232   arcextru = geompy.MakePrismVecH(arc1, OX, h1)
 233
 234   # --- plan de coupe à 45° sur le cylindre principal,
 235   #     section à 45° du cylndre principal,
 236   #     section du cylindre secondaire par l'enveloppe cylindrique du cylindre principal
 237
 238   plan1 = geompy.MakePlane(O, OX, 4*r1)
 239   planr = geompy.MakeRotation(plan1, OY, 45*math.pi/180.0)
 240   geompy.addToStudy(planr, 'planr')
 241
 242   sect45 = geompy.MakeCommonList([demicyl1, planr], True)
 243   geompy.addToStudy(sect45, 'sect45')
 244
 245   sect90 = geompy.MakeCommonList([demicyl2, arcextru], True)
 246   geompy.addToStudy(sect90, 'sect90')
 247
 248   # --- liste ordonnée des points projetés sur les deux sections
 249
 250   vord45 = pointsProjetes(study, v1, sect45)
 251   vord90 = pointsProjetes(study, v2, sect90)
 252   for i,v in enumerate(vord45):
 253     geompy.addToStudyInFather(sect45, v, 'v%d'%i)
 254   for i,v in enumerate(vord90):
 255     geompy.addToStudyInFather(sect90, v, 'v%d'%i)
 256
 257   # --- identification des projections des trois arcs de cercle, sur les deux sections.
 258
 259   lord45 = arcsProjetes(study, vord45, sect45)
 260   lord90 = arcsProjetes(study, vord90, sect90)
 261   for i,l in enumerate(lord45):
 262     geompy.addToStudyInFather(sect45, l, 'l%d'%i)
 263   for i,l in enumerate(lord90):
 264     geompy.addToStudyInFather(sect90, l, 'l%d'%i)
 265
 266   # --- abaissement des quatre points centraux de la section du cylindre secondaire
 267
 268   #if with_solid:
 269     #dz = -(r2 + solid_thickness)/2.0
 270   #else:
 271     #dz = -r2/2.0
 272   dz = -r2/2.0
 273   for i in (0, 2, 4, 5):
 274     vord90[i] = geompy.TranslateDXDYDZ(vord90[i], 0, 0, dz, True)
 275     geompy.addToStudyInFather(sect90, vord90[i], 'vm%d'%i)
 276   #if with_solid:
 277     #for i in (1, 3, 6, 7):
 278       #vord90[i] = geompy.TranslateDXDYDZ(vord90[i], 0, 0, dz*solid_thickness/(r2+solid_thickness), True)
 279
 280   """
 281   res=vord90
 282   return res
 283   """
 284
 285   # --- création des deux arêtes curvilignes sur l'enveloppe cylindrique du cylindre principal, à la jonction
 286
 287   curv = [None for i in range(4)] # liaisons entre les points 1, 3, 6 et 7 des 2 sections
 288
 289   curv[0] = geompy.MakeArcCenter(O, vord90[1] , vord45[1], False)
 290   curv[1] = geompy.MakeArcCenter(O, vord90[3] , vord45[3], False)
 291
 292   lipts = ((6, 6, 4), (7, 7, 5))
 293   for i, ipts in enumerate(lipts):
 294     print i, ipts
 295     p0 = vord90[ipts[0]]
 296     p1 = vord45[ipts[1]]
 297     p2 = vord45[ipts[2]]
 298     plan = geompy.MakePlaneThreePnt(p0, p1, p2, 10000)
 299     #geompy.addToStudy(plan, "plan%d"%i)
 300     section = geompy.MakeSection(plan, arcextru, True)
 301     secpart = geompy.MakePartition([section], [sect45, sect90], [], [], geompy.ShapeType["EDGE"], 0, [], 0, True)
 302     geompy.addToStudy(secpart, "secpart%d"%i)
 303     lsec = geompy.ExtractShapes(secpart, geompy.ShapeType["EDGE"], True)
 304     #print "len(lsec)", len(lsec)
 305
 306     # TODO : revoir ça dans le cas avec solide
 307     for l in lsec:
 308       pts = geompy.ExtractShapes(l, geompy.ShapeType["VERTEX"], True)
 309       if (((geompy.MinDistance(pts[0], p0) < 0.001) and (geompy.MinDistance(pts[1], p1) < 0.001)) or
 310           ((geompy.MinDistance(pts[1], p0) < 0.001) and (geompy.MinDistance(pts[0], p1) < 0.001))):
 311         curv[i+2] =l
 312         print "curv_%d OK"%i
 313         break
 314   # RNC : commente temporairement
 315   #for i,l in enumerate(curv):
 316   #  geompy.addToStudyInFather(arcextru, l, "curv%d"%i)
 317
 318   # --- creation des arêtes droites manquantes, des faces et volumes pour les quatre volumes de la jonction
 319
 320   edges = [None for i in range(8)]
 321   edges[0] = geompy.MakeLineTwoPnt(vord45[0], vord90[0])
 322   edges[1] = curv[0]
 323   edges[2] = geompy.MakeLineTwoPnt(vord45[2], vord90[2])
 324   edges[3] = curv[1]
 325   edges[4] = geompy.MakeLineTwoPnt(vord45[4], vord90[4])
 326   edges[5] = geompy.MakeLineTwoPnt(vord45[5], vord90[5])
 327   edges[6] = curv[2]
 328   edges[7] = curv[3]
 329   for i,l in enumerate(edges):
 330     print i
 331     geompy.addToStudy( l, "edge%d"%i)
 332
 333   ed45 = [None for i in range(8)]
 334   ed45[0] = geompy.MakeLineTwoPnt(vord45[0], vord45[2])
 335   ed45[1] = geompy.MakeLineTwoPnt(vord45[0], vord45[1])
 336   ed45[2] = geompy.MakeLineTwoPnt(vord45[4], vord45[6])
 337   ed45[3] = geompy.MakeLineTwoPnt(vord45[2], vord45[3])
 338   ed45[4] = geompy.MakeLineTwoPnt(vord45[5], vord45[7])
 339   ed45[5] = geompy.MakeLineTwoPnt(vord45[4], vord45[5])
 340   ed45[6] = geompy.MakeLineTwoPnt(vord45[0], vord45[4])
 341   ed45[7] = geompy.MakeLineTwoPnt(vord45[2], vord45[5])
 342   for i,l in enumerate(ed45):
 343     geompy.addToStudyInFather(sect45, l, "ed45_%d"%i)
 344
 345   ed90 = [None for i in range(8)]
 346   ed90[0] = geompy.MakeLineTwoPnt(vord90[0], vord90[2])
 347   ed90[1] = geompy.MakeLineTwoPnt(vord90[0], vord90[1])
 348   ed90[2] = geompy.MakeLineTwoPnt(vord90[4], vord90[6])
 349   ed90[3] = geompy.MakeLineTwoPnt(vord90[2], vord90[3])
 350   ed90[4] = geompy.MakeLineTwoPnt(vord90[5], vord90[7])
 351   ed90[5] = geompy.MakeLineTwoPnt(vord90[4], vord90[5])
 352   ed90[6] = geompy.MakeLineTwoPnt(vord90[0], vord90[4])
 353   ed90[7] = geompy.MakeLineTwoPnt(vord90[2], vord90[5])
 354   for i,l in enumerate(ed90):
 355     geompy.addToStudyInFather(sect90, l, "ed90_%d"%i)
 356
 357   faci = []
 358   faci.append(geompy.MakeFaceWires([ed45[6], edges[0], ed90[6], edges[4]], 0))
 359   faci.append(geompy.MakeFaceWires([ed45[7], edges[2], ed90[7], edges[5]], 0))
 360   faci.append(geompy.MakeFaceWires([ed45[2], edges[4], ed90[2], edges[6]], 0))
 361   faci.append(geompy.MakeFaceWires([ed45[5], edges[4], ed90[5], edges[5]], 0))
 362   faci.append(geompy.MakeFaceWires([ed45[4], edges[5], ed90[4], edges[7]], 0))
 363   faci.append(geompy.MakeFaceWires([ed90[0], ed90[6], ed90[5], ed90[7]], 0))
 364   faci.append(geompy.MakeFaceWires([ed90[1], ed90[6], ed90[2], lord90[0]], 0))
 365   faci.append(geompy.MakeFaceWires([ed90[2], ed90[5], ed90[4], lord90[1]], 0))
 366   faci.append(geompy.MakeFaceWires([ed90[3], ed90[7], ed90[4], lord90[2]], 0))
 367   for i,f in enumerate(faci):
 368     geompy.addToStudy(f, "faci_%d"%i)
 369
 370   return faci, sect45, arc1, l1, lord90, lord45, edges, arcextru
 371
 372 if __name__=="__main__":
 373   """For testing purpose"""
 374   salome.salome_init()
 375   theStudy = salome.myStudy
 376   geompy = geomBuilder.New(theStudy)
 377   res = build_shape(theStudy, 80., 20., 100., 100., 10.)
 378   """
 379   for i,v in enumerate(res):
 380     geompy.addToStudy(v,"v%d"%i)
 381   """
 382
 383   #res = demidisk(theStudy, 80, 45, 0, 10)
 384   #geompy.addToStudy(res[3], "res")
 385   #for i,v in enumerate(res[0]):
 386   #  geompy.addToStudy(v,"v%d"%i)
 387   geompy.addToStudy(res, "res")
 388