GDML polyhedra

[lambdam94]

Hello Keith

With the polyhedra: `

  <zplane rmin="70" rmax="100" z="80" />
</polyhedra>`

We don't see 11 sides for phi between 1 rad and 5 rad (4+1).

With the polyhedra: `

  <zplane rmin="70" rmax="100" z="80" />
</polyhedra>`

We see the 11 sides for phi between 0 and 360 deg.

The polyhedra should have N (number of numsides) "same" elements for phi between startphi and thetaphi. Do you think so?

Regards Dam

[KeithSloan]

Another good error spot. I ran gdmlview and it shows there should be 11 sides even when the polyhedra is not a 360 revolution. Should not be too hard to fix, Need to adjust the number of sides before the cut.

[KeithSloan]

Have made a change `numsides = int(numsides * 360 / getAngleDeg(fp.aunit,fp.deltaphi)

Not sure about rounding though.`

[lambdam94]

Hello Keith,

I tried to code a simple createGeometry for GDMLpolyhedra. It works well on Freecad 0.18 Could you take a look? Thanks!

Regards Dam

` def createGeometry(self,fp):

GDMLShared.setTrace(True)

   #GDMLShared.trace("Execute Polyhedra")

   parms = fp.OutList
   numsides = fp.numsides
   if numsides<1: numsides=1
   mul = GDMLShared.getMult(fp)
   startphi_rad = getAngleRad(fp.aunit, fp.startphi)
   deltaphi_rad = getAngleRad(fp.aunit, fp.deltaphi)
   if(deltaphi_rad<0.0): deltaphi_rad=0.0
   elif(deltaphi_rad>2.0*math.pi): deltaphi_rad=2.0*math.pi

   deltaPhiSide_rad = deltaphi_rad/float(numsides)

   faceList = []

   if deltaphi_rad<=0.0: return

   for indiceParm in range(len(parms)-1):
       zmin    = parms[indiceParm].z * mul
       zmax    = parms[indiceParm+1].z * mul
       rmin_i = 0.0
       rmin_ip1 = 0.0
       if hasattr(parms[indiceParm], 'rmin'):
           rmin_i = max( parms[indiceParm].rmin * mul, 0.0)
       if hasattr(parms[indiceParm+1], 'rmin'):
           rmin_ip1 = max( parms[indiceParm+1].rmin * mul, 0.0)
       rmax_i = max(parms[indiceParm].rmax * mul, 0.0)
       rmax_ip1 = max(parms[indiceParm+1].rmax * mul, 0.0)

       for indiceSide in range(numsides):
           cosmin = math.cos(float(indiceSide)*deltaPhiSide_rad+startphi_rad)
           sinmin = math.sin(float(indiceSide)*deltaPhiSide_rad+startphi_rad)
           cosmax = math.cos(float(indiceSide+1)*deltaPhiSide_rad+startphi_rad)
           sinmax = math.sin(float(indiceSide+1)*deltaPhiSide_rad+startphi_rad)

           # def of the 8 corners
           rmin_zmin_phimin = FreeCAD.Vector(rmin_i*cosmin, rmin_i*sinmin, zmin)
           rmax_zmin_phimin = FreeCAD.Vector(rmax_i*cosmin, rmax_i*sinmin, zmin)
           rminp1_zmax_phimin = FreeCAD.Vector(rmin_ip1*cosmin, rmin_ip1*sinmin, zmax)
           rmaxp1_zmax_phimin = FreeCAD.Vector(rmax_ip1*cosmin, rmax_ip1*sinmin, zmax)

           rmin_zmin_phimax = FreeCAD.Vector(rmin_i*cosmax, rmin_i*sinmax, zmin)
           rmax_zmin_phimax = FreeCAD.Vector(rmax_i*cosmax, rmax_i*sinmax, zmin)
           rminp1_zmax_phimax = FreeCAD.Vector(rmin_ip1*cosmax, rmin_ip1*sinmax, zmax)
           rmaxp1_zmax_phimax = FreeCAD.Vector(rmax_ip1*cosmax, rmax_ip1*sinmax, zmax)

           # face Z max
           if indiceParm ==(len(parms)-2):
               faceList.append(Part.Face(Part.makePolygon([rminp1_zmax_phimin,
                                                  rmaxp1_zmax_phimin,
                                                  rmaxp1_zmax_phimax,
                                                  rminp1_zmax_phimax,
                                                  rminp1_zmax_phimin])))
           #face Z min
           if indiceParm ==0:
               faceList.append(Part.Face(Part.makePolygon([rmin_zmin_phimin,
                                                  rmax_zmin_phimin,
                                                  rmax_zmin_phimax,
                                                  rmin_zmin_phimax,
                                                  rmin_zmin_phimin])))

           # R out
           faceList.append(Part.Face(Part.makePolygon([rmax_zmin_phimin,
                                                  rmaxp1_zmax_phimin,
                                                  rmaxp1_zmax_phimax,
                                                  rmax_zmin_phimax,
                                                  rmax_zmin_phimin])))

           # R min 
           if( rmin_i > 0.0 and rmin_ip1 > 0.0):
               faceList.append(Part.Face(Part.makePolygon([rmin_zmin_phimin,
                                                  rminp1_zmax_phimin,
                                                  rminp1_zmax_phimax,
                                                  rmin_zmin_phimax,
                                                  rmin_zmin_phimin])))

           if ( deltaphi_rad < 2*math.pi ):
               # border at startphi
               if indiceSide == 0 :
                   faceList.append(Part.Face(Part.makePolygon([rmin_zmin_phimin,
                                                  rmax_zmin_phimin,
                                                  rmaxp1_zmax_phimin,
                                                  rminp1_zmax_phimin,
                                                  rmin_zmin_phimin])))

               # border at startphi+deltaphi
               if indiceSide == numsides-1 :
                   faceList.append(Part.Face(Part.makePolygon([rmin_zmin_phimax,
                                                  rmax_zmin_phimax,
                                                  rmaxp1_zmax_phimax,
                                                  rminp1_zmax_phimax,
                                                  rmin_zmin_phimax])))

           #faceList.append(faceUp)
           #faceList.append(faceDown)

   fp.Shape = Part.makeShell(faceList)
   del faceList[:]
   del faceList

`

[KeithSloan]

Seems a lot more complicated and use a number of trig functions.

Are you saying the current version does not work?

[lambdam94]

I have done some tests with pb on Freecad 0.18 with previous version. Now, with current version the classic cases works!

I saw a small pb with:

where each subparts are not equal  The suggested function is very badly commented. It is a loop on each sub-volume: the 8 cornes are calculated then the outer faces are added to a list. On my old pc, the suggested function seems faster. But Seems is not a measurement! Perhaps, some comparison tests need to be done (with loops or others ideas). With or without del function etc. Have you some ideas to do them? Regards Dam

[KeithSloan]

You could always create a branch and remove all the GDMLShared.trace calls. I am sure that would speed things up a lot.

i.e. git branch notrace git checkout notrace

edit files and remove all GDMLShared.trace calls.

Or try in Utils - gdml2fc.py

I cannot remember if it works or not, so see if it works with a small file and if yes just leave it running overnight on the big CERN files.