Implement G4Trap in Opticks

[ggalgoczi]

No description provided.

[ggalgoczi]

G4Trd does not need to be implemented since it can be translated analytically into a G4Trap.

[ggalgoczi]

Can I assign it to you @plexoos ?

[plexoos]

Yes, I'll see what I can do

[plexoos]

We have both G4Trap and G4Trd (aka trap and trd in GDML) volumes used in the pfRICH geometry. G4Trd is a special case of G4Trap with the z-axis being its central axis.

[ggalgoczi]

That is correct, I provided a way to automatically convert all G4Trds into G4Traps. I can not find it right now but will add a comment here when I do.

Hence only G4Traps need to implemented in this task into Opticks.

[plexoos]

Looks like we can associate up to 6 double parameters with each node https://github.com/BNLNPPS/esi-opticks/blob/esi/sysrap/s_pa.h#L50. But G4Trap requires 11 parameters (9 of which are independent) 🤔

[ggalgoczi]

Yes, I feel like we have to change this sourcecode and add additional variables which will lead to need to change quite a bit of code...

Simon suggested the following steps to do this:

1. study Geant4 implementation of G4Trap (g4-;g4-cls G4Trap)
2. study old Opticks translations with extg4/X4Solid.*
3. study new Opticks translations with u4/U4Solid.*

[plexoos]

We are in luck. All of the trap volumes in pfrich_min have zero angles:

$ grep -r trap geom/pfrich_min.gdml 
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008b62c830" phi="0" theta="0" x1="116.238315670544" x2="62.0667122987972" x3="116.238315670544" x4="62.0667122987972" y1="54.0187869822485" y2="54.0187869822485" z="6.35"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008bf67cc0" phi="0" theta="0" x1="103.069073741954" x2="48.8974703702069" x3="103.069073741954" x4="48.8974703702069" y1="54.0187869822485" y2="54.0187869822485" z="7.35"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008b759510" phi="0" theta="0" x1="116.238315670544" x2="62.0667122987972" x3="116.238315670544" x4="62.0667122987972" y1="54.0187869822485" y2="54.0187869822485" z="431.65"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008b62c3e0" phi="0" theta="0" x1="103.069073741954" x2="48.8974703702069" x3="103.069073741954" x4="48.8974703702069" y1="54.0187869822485" y2="54.0187869822485" z="432.65"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008bf88180" phi="0" theta="0" x1="116.774668179175" x2="62.6030648074284" x3="116.774668179175" x4="62.6030648074284" y1="54.0187869822485" y2="54.0187869822485" z="435.65"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008b62c550" phi="0" theta="0" x1="103.069073741954" x2="48.8974703702069" x3="103.069073741954" x4="48.8974703702069" y1="54.0187869822485" y2="54.0187869822485" z="435.65"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008bf697f0" phi="0" theta="0" x1="114.42186850798" x2="60.2502651362331" x3="114.42186850798" x4="60.2502651362331" y1="54.0187869822485" y2="54.0187869822485" z="6.35"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008bf6a230" phi="0" theta="0" x1="103.069073741954" x2="48.8974703702069" x3="103.069073741954" x4="48.8974703702069" y1="54.0187869822485" y2="54.0187869822485" z="7.35"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008bf68b00" phi="0" theta="0" x1="103.069073741954" x2="48.8974703702069" x3="103.069073741954" x4="48.8974703702069" y1="54.0187869822485" y2="54.0187869822485" z="7.35"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008c098340" phi="0" theta="0" x1="103.069073741954" x2="48.8974703702069" x3="103.069073741954" x4="48.8974703702069" y1="54.0187869822485" y2="54.0187869822485" z="10.7"/>
    <trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="FlangeWedge0x56008b600820" phi="0" theta="0" x1="103.069073741954" x2="48.8974703702069" x3="103.069073741954" x4="48.8974703702069" y1="54.0187869822485" y2="54.0187869822485" z="492"/>

[ggalgoczi]

I saw your commit @plexoos, very interesting! Did you test it? Does it work or does G4Trap need to be converted into Polyhedron and have a function to calculate it's bounding box for GPU use?

By the way do you know what sn and snd stand for?

[plexoos]

AFAICT, sn/snd are the types that hold information about each node in the geometry tree. I implemented the right angular wedge from the trapezoid definitions, and it seems to allow processing a simplified geometry with a union that includes a trapezoid. More tests are needed.

[ggalgoczi]

nconvexpolyhedron* nconvexpolyhedron::CreateTrapezoid(float z, float x1, float y1, float x2, float y2 )
was used to convert G4Trd in previous Opticks version into Polyhedron. Do we want to do the same?

The intersection is implemented for Polyhedron for Optix7!
./CSG/csg_intersect_leaf_convexpolyhedron.h

[ggalgoczi]

Weirdly in the pfRICH geometry when I try to convert trd into a trap I can not get it right.

For

<trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="Dummy0x56008bf728c0"
      phi="0" theta="0"
      x1="123" x2="123" x3="104" x4="104"
      y1="123" y2="104"
      z="30.022"/>

I get:

simg4ox: /esi/opticks/sysrap/sn.h:2820: static sn* sn::Trapezoid(double, double, double, double): Assertion `x > 0 && y > 0 && z > 0 && ltx > 0 && ltx < x' failed.
Aborted (core dumped)

For:

<trap alpha1="0" alpha2="0" aunit="deg" lunit="mm" name="Dummy0x56008bf728c0"
      phi="0" theta="0"
      x1="123" x2="104" x3="123" x4="104"
      y1="123" y2="104"
      z="30.022"/>

I get:

-------- EEEE ------- G4Exception-START -------- EEEE -------
*** G4Exception : GeomSolids0002
      issued by : G4Trap::MakePlanes()
Side face ~-X is not planar for solid: Dummy0x56008bf728c0
Discrepancy: -0.396192 mm
-----------------------------------------------------------
    *** Dump for solid: Dummy0x56008bf728c0 ***
    ===================================================
 Solid type: G4Trap
 Parameters:
    half length Z: 15.011 mm
    half length Y, face -Dz: 61.5 mm
    half length X, face -Dz, side -Dy1: 61.5 mm
    half length X, face -Dz, side +Dy1: 52 mm
    half length Y, face +Dz: 52 mm
    half length X, face +Dz, side -Dy2: 61.5 mm
    half length X, face +Dz, side +Dy2: 52 mm
    theta: 0 degrees
    phi:   0 degrees
    alpha, face -Dz: 0 degrees
    alpha, face +Dz: 0 degrees
-----------------------------------------------------------

*** Fatal Exception *** core dump ***
 **** Track information is not available at this moment
 **** Step information is not available at this moment

-------- EEEE -------- G4Exception-END --------- EEEE -------

@plexoos are you sure assert( x > 0 && y > 0 && z > 0 && ltx > 0 && ltx < x ); is needed?