Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-12-17 09:21:55

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 //
0027 // $Id: PHG4GDMLWriteSolids.cc 81843 2014-06-06 09:11:11Z gcosmo $
0028 //
0029 // class PHG4GDMLWriteSolids Implementation
0030 //
0031 // Original author: Zoltan Torzsok, November 2007
0032 //
0033 // --------------------------------------------------------------------
0034 
0035 #include "PHG4GDMLWriteSolids.hh"
0036 
0037 #include <Geant4/G4BooleanSolid.hh>
0038 #include <Geant4/G4Box.hh>
0039 #include <Geant4/G4Cons.hh>
0040 #include <Geant4/G4CutTubs.hh>
0041 #include <Geant4/G4Ellipsoid.hh>
0042 #include <Geant4/G4EllipticalCone.hh>
0043 #include <Geant4/G4EllipticalTube.hh>
0044 #include <Geant4/G4ExtrudedSolid.hh>
0045 #include <Geant4/G4GenericPolycone.hh>
0046 #include <Geant4/G4GenericTrap.hh>
0047 #include <Geant4/G4Hype.hh>
0048 #include <Geant4/G4IntersectionSolid.hh>
0049 #include <Geant4/G4OpticalSurface.hh>
0050 #include <Geant4/G4Orb.hh>
0051 #include <Geant4/G4Para.hh>
0052 #include <Geant4/G4Paraboloid.hh>
0053 #include <Geant4/G4Polycone.hh>
0054 #include <Geant4/G4Polyhedra.hh>
0055 #include <Geant4/G4ReflectedSolid.hh>
0056 #include <Geant4/G4Sphere.hh>
0057 #include <Geant4/G4SubtractionSolid.hh>
0058 #include <Geant4/G4SurfaceProperty.hh>
0059 #include <Geant4/G4SystemOfUnits.hh>
0060 #include <Geant4/G4TessellatedSolid.hh>
0061 #include <Geant4/G4Tet.hh>
0062 #include <Geant4/G4Torus.hh>
0063 #include <Geant4/G4Trap.hh>
0064 #include <Geant4/G4Trd.hh>
0065 #include <Geant4/G4Tubs.hh>
0066 #include <Geant4/G4TwistedBox.hh>
0067 #include <Geant4/G4TwistedTrap.hh>
0068 #include <Geant4/G4TwistedTrd.hh>
0069 #include <Geant4/G4TwistedTubs.hh>
0070 #include <Geant4/G4UnionSolid.hh>
0071 
0072 PHG4GDMLWriteSolids::PHG4GDMLWriteSolids()
0073   : solidsElement(nullptr)
0074 {
0075 }
0076 
0077 PHG4GDMLWriteSolids::~PHG4GDMLWriteSolids()
0078 {
0079 }
0080 
0081 #if !defined(G4GEOM_USE_USOLIDS)
0082 void PHG4GDMLWriteSolids::
0083     MultiUnionWrite(xercesc::DOMElement* /*unused*/,
0084                     const G4MultiUnion* const /*unused*/)
0085 {
0086   G4Exception("PHG4GDMLWriteSolids::MultiUnionWrite()",
0087               "InvalidSetup", FatalException,
0088               "Installation with USolids primitives required!");
0089   return;
0090 }
0091 #else
0092 void PHG4GDMLWriteSolids::
0093     MultiUnionWrite(xercesc::DOMElement* solElement,
0094                     const G4MultiUnion* const munionSolid)
0095 {
0096   G4int numSolids = munionSolid->GetNumberOfSolids();
0097   G4String tag("multiUnion");
0098 
0099   const G4String& name = GenerateName(munionSolid->GetName(), munionSolid);
0100   xercesc::DOMElement* multiUnionElement = NewElement(tag);
0101   multiUnionElement->setAttributeNode(NewAttribute("name", name));
0102 
0103   for (G4int i = 0; i < numSolids; ++i)
0104   {
0105     G4VSolid* solid = munionSolid->GetSolid(i);
0106     G4Transform3D* transform = munionSolid->GetTransformation(i);
0107 
0108     HepGeom::Rotate3D rot3d;
0109     HepGeom::Translate3D transl;
0110     HepGeom::Scale3D scale;
0111     transform->getDecomposition(scale, rot3d, transl);
0112 
0113     G4ThreeVector pos = transl.getTranslation();
0114     G4RotationMatrix
0115         rotm(CLHEP::HepRep3x3(rot3d.xx(), rot3d.xy(), rot3d.xz(),
0116                               rot3d.yx(), rot3d.yy(), rot3d.yz(),
0117                               rot3d.zx(), rot3d.zy(), rot3d.zz()));
0118     G4ThreeVector rot = GetAngles(rotm);
0119 
0120     AddSolid(solid);
0121     const G4String& solidref = GenerateName(solid->GetName(), solid);
0122     std::ostringstream os;
0123     os << i + 1;
0124     const G4String& nodeName = "Node-" + G4String(os.str());
0125     xercesc::DOMElement* solidElement = NewElement("solid");
0126     solidElement->setAttributeNode(NewAttribute("ref", solidref));
0127     xercesc::DOMElement* multiUnionNodeElement = NewElement("multiUnionNode");
0128     multiUnionNodeElement->setAttributeNode(NewAttribute("name", nodeName));
0129     multiUnionNodeElement->appendChild(solidElement);  // Append solid to node
0130     if ((std::fabs(pos.x()) > kLinearPrecision) || (std::fabs(pos.y()) > kLinearPrecision) || (std::fabs(pos.z()) > kLinearPrecision))
0131     {
0132       PositionWrite(multiUnionNodeElement, name + "_pos", pos);
0133     }
0134     if ((std::fabs(rot.x()) > kAngularPrecision) || (std::fabs(rot.y()) > kAngularPrecision) || (std::fabs(rot.z()) > kAngularPrecision))
0135     {
0136       RotationWrite(multiUnionNodeElement, name + "_rot", rot);
0137     }
0138     multiUnionElement->appendChild(multiUnionNodeElement);  // Append node
0139   }
0140 
0141   solElement->appendChild(multiUnionElement);
0142   // Add the multiUnion solid AFTER the constituent nodes!
0143 }
0144 #endif
0145 
0146 void PHG4GDMLWriteSolids::
0147     BooleanWrite(xercesc::DOMElement* solElement,
0148                  const G4BooleanSolid* const boolean)
0149 {
0150   G4int displaced = 0;
0151 
0152   G4String tag("undefined");
0153   if (dynamic_cast<const G4IntersectionSolid*>(boolean))
0154   {
0155     tag = "intersection";
0156   }
0157   else if (dynamic_cast<const G4SubtractionSolid*>(boolean))
0158   {
0159     tag = "subtraction";
0160   }
0161   else if (dynamic_cast<const G4UnionSolid*>(boolean))
0162   {
0163     tag = "union";
0164   }
0165 
0166   G4VSolid* firstPtr = const_cast<G4VSolid*>(boolean->GetConstituentSolid(0));
0167   G4VSolid* secondPtr = const_cast<G4VSolid*>(boolean->GetConstituentSolid(1));
0168 
0169   G4ThreeVector firstpos;
0170   G4ThreeVector firstrot;
0171   G4ThreeVector pos;
0172   G4ThreeVector rot;
0173 
0174   // Solve possible displacement of referenced solids!
0175   //
0176   while (true)
0177   {
0178     if (displaced > 8)
0179     {
0180       G4String ErrorMessage = "The referenced solid '" + firstPtr->GetName() +
0181                               +"in the Boolean shape '" +
0182                               +boolean->GetName() +
0183                               +"' was displaced too many times!";
0184       G4Exception("PHG4GDMLWriteSolids::BooleanWrite()",
0185                   "InvalidSetup", FatalException, ErrorMessage);
0186     }
0187 
0188     if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(firstPtr))
0189     {
0190       firstpos += disp->GetObjectTranslation();
0191       firstrot += GetAngles(disp->GetObjectRotation());
0192       firstPtr = disp->GetConstituentMovedSolid();
0193       displaced++;
0194       continue;
0195     }
0196     break;
0197   }
0198   displaced = 0;
0199   while (true)
0200   {
0201     if (displaced > maxTransforms)
0202     {
0203       G4String ErrorMessage = "The referenced solid '" + secondPtr->GetName() +
0204                               +"in the Boolean shape '" +
0205                               +boolean->GetName() +
0206                               +"' was displaced too many times!";
0207       G4Exception("PHG4GDMLWriteSolids::BooleanWrite()",
0208                   "InvalidSetup", FatalException, ErrorMessage);
0209     }
0210 
0211     if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(secondPtr))
0212     {
0213       pos += disp->GetObjectTranslation();
0214       rot += GetAngles(disp->GetObjectRotation());
0215       secondPtr = disp->GetConstituentMovedSolid();
0216       displaced++;
0217       continue;
0218     }
0219     break;
0220   }
0221 
0222   AddSolid(firstPtr);  // At first add the constituent solids!
0223   AddSolid(secondPtr);
0224 
0225   const G4String& name = GenerateName(boolean->GetName(), boolean);
0226   const G4String& firstref = GenerateName(firstPtr->GetName(), firstPtr);
0227   const G4String& secondref = GenerateName(secondPtr->GetName(), secondPtr);
0228 
0229   xercesc::DOMElement* booleanElement = NewElement(tag);
0230   booleanElement->setAttributeNode(NewAttribute("name", name));
0231   xercesc::DOMElement* firstElement = NewElement("first");
0232   firstElement->setAttributeNode(NewAttribute("ref", firstref));
0233   booleanElement->appendChild(firstElement);
0234   xercesc::DOMElement* secondElement = NewElement("second");
0235   secondElement->setAttributeNode(NewAttribute("ref", secondref));
0236   booleanElement->appendChild(secondElement);
0237   solElement->appendChild(booleanElement);
0238   // Add the boolean solid AFTER the constituent solids!
0239 
0240   if ((std::fabs(pos.x()) > kLinearPrecision) || (std::fabs(pos.y()) > kLinearPrecision) || (std::fabs(pos.z()) > kLinearPrecision))
0241   {
0242     PositionWrite(booleanElement, name + "_pos", pos);
0243   }
0244 
0245   if ((std::fabs(rot.x()) > kAngularPrecision) || (std::fabs(rot.y()) > kAngularPrecision) || (std::fabs(rot.z()) > kAngularPrecision))
0246   {
0247     RotationWrite(booleanElement, name + "_rot", rot);
0248   }
0249 
0250   if ((std::fabs(firstpos.x()) > kLinearPrecision) || (std::fabs(firstpos.y()) > kLinearPrecision) || (std::fabs(firstpos.z()) > kLinearPrecision))
0251   {
0252     FirstpositionWrite(booleanElement, name + "_fpos", firstpos);
0253   }
0254 
0255   if ((std::fabs(firstrot.x()) > kAngularPrecision) || (std::fabs(firstrot.y()) > kAngularPrecision) || (std::fabs(firstrot.z()) > kAngularPrecision))
0256   {
0257     FirstrotationWrite(booleanElement, name + "_frot", firstrot);
0258   }
0259 }
0260 
0261 void PHG4GDMLWriteSolids::
0262     BoxWrite(xercesc::DOMElement* solElement, const G4Box* const box)
0263 {
0264   const G4String& name = GenerateName(box->GetName(), box);
0265 
0266   xercesc::DOMElement* boxElement = NewElement("box");
0267   boxElement->setAttributeNode(NewAttribute("name", name));
0268   boxElement->setAttributeNode(NewAttribute("x", 2.0 * box->GetXHalfLength() / mm));
0269   boxElement->setAttributeNode(NewAttribute("y", 2.0 * box->GetYHalfLength() / mm));
0270   boxElement->setAttributeNode(NewAttribute("z", 2.0 * box->GetZHalfLength() / mm));
0271   boxElement->setAttributeNode(NewAttribute("lunit", "mm"));
0272   solElement->appendChild(boxElement);
0273 }
0274 
0275 void PHG4GDMLWriteSolids::
0276     ConeWrite(xercesc::DOMElement* solElement, const G4Cons* const cone)
0277 {
0278   const G4String& name = GenerateName(cone->GetName(), cone);
0279 
0280   xercesc::DOMElement* coneElement = NewElement("cone");
0281   coneElement->setAttributeNode(NewAttribute("name", name));
0282   coneElement->setAttributeNode(NewAttribute("rmin1", cone->GetInnerRadiusMinusZ() / mm));
0283   coneElement->setAttributeNode(NewAttribute("rmax1", cone->GetOuterRadiusMinusZ() / mm));
0284   coneElement->setAttributeNode(NewAttribute("rmin2", cone->GetInnerRadiusPlusZ() / mm));
0285   coneElement->setAttributeNode(NewAttribute("rmax2", cone->GetOuterRadiusPlusZ() / mm));
0286   coneElement->setAttributeNode(NewAttribute("z", 2.0 * cone->GetZHalfLength() / mm));
0287   coneElement->setAttributeNode(NewAttribute("startphi", cone->GetStartPhiAngle() / degree));
0288   coneElement->setAttributeNode(NewAttribute("deltaphi", cone->GetDeltaPhiAngle() / degree));
0289   coneElement->setAttributeNode(NewAttribute("aunit", "deg"));
0290   coneElement->setAttributeNode(NewAttribute("lunit", "mm"));
0291   solElement->appendChild(coneElement);
0292 }
0293 
0294 void PHG4GDMLWriteSolids::
0295     ElconeWrite(xercesc::DOMElement* solElement,
0296                 const G4EllipticalCone* const elcone)
0297 {
0298   const G4String& name = GenerateName(elcone->GetName(), elcone);
0299 
0300   xercesc::DOMElement* elconeElement = NewElement("elcone");
0301   elconeElement->setAttributeNode(NewAttribute("name", name));
0302   elconeElement->setAttributeNode(NewAttribute("dx", elcone->GetSemiAxisX() / mm));
0303   elconeElement->setAttributeNode(NewAttribute("dy", elcone->GetSemiAxisY() / mm));
0304   elconeElement->setAttributeNode(NewAttribute("zmax", elcone->GetZMax() / mm));
0305   elconeElement->setAttributeNode(NewAttribute("zcut", elcone->GetZTopCut() / mm));
0306   elconeElement->setAttributeNode(NewAttribute("lunit", "mm"));
0307   solElement->appendChild(elconeElement);
0308 }
0309 
0310 void PHG4GDMLWriteSolids::
0311     EllipsoidWrite(xercesc::DOMElement* solElement,
0312                    const G4Ellipsoid* const ellipsoid)
0313 {
0314   const G4String& name = GenerateName(ellipsoid->GetName(), ellipsoid);
0315 
0316   xercesc::DOMElement* ellipsoidElement = NewElement("ellipsoid");
0317   ellipsoidElement->setAttributeNode(NewAttribute("name", name));
0318   ellipsoidElement->setAttributeNode(NewAttribute("ax", ellipsoid->GetSemiAxisMax(0) / mm));
0319   ellipsoidElement->setAttributeNode(NewAttribute("by", ellipsoid->GetSemiAxisMax(1) / mm));
0320   ellipsoidElement->setAttributeNode(NewAttribute("cz", ellipsoid->GetSemiAxisMax(2) / mm));
0321   ellipsoidElement->setAttributeNode(NewAttribute("zcut1", ellipsoid->GetZBottomCut() / mm));
0322   ellipsoidElement->setAttributeNode(NewAttribute("zcut2", ellipsoid->GetZTopCut() / mm));
0323   ellipsoidElement->setAttributeNode(NewAttribute("lunit", "mm"));
0324   solElement->appendChild(ellipsoidElement);
0325 }
0326 
0327 void PHG4GDMLWriteSolids::
0328     EltubeWrite(xercesc::DOMElement* solElement,
0329                 const G4EllipticalTube* const eltube)
0330 {
0331   const G4String& name = GenerateName(eltube->GetName(), eltube);
0332 
0333   xercesc::DOMElement* eltubeElement = NewElement("eltube");
0334   eltubeElement->setAttributeNode(NewAttribute("name", name));
0335   eltubeElement->setAttributeNode(NewAttribute("dx", eltube->GetDx() / mm));
0336   eltubeElement->setAttributeNode(NewAttribute("dy", eltube->GetDy() / mm));
0337   eltubeElement->setAttributeNode(NewAttribute("dz", eltube->GetDz() / mm));
0338   eltubeElement->setAttributeNode(NewAttribute("lunit", "mm"));
0339   solElement->appendChild(eltubeElement);
0340 }
0341 
0342 void PHG4GDMLWriteSolids::
0343     XtruWrite(xercesc::DOMElement* solElement,
0344               const G4ExtrudedSolid* const xtru)
0345 {
0346   const G4String& name = GenerateName(xtru->GetName(), xtru);
0347 
0348   xercesc::DOMElement* xtruElement = NewElement("xtru");
0349   xtruElement->setAttributeNode(NewAttribute("name", name));
0350   xtruElement->setAttributeNode(NewAttribute("lunit", "mm"));
0351   solElement->appendChild(xtruElement);
0352 
0353   const G4int NumVertex = xtru->GetNofVertices();
0354 
0355   for (G4int i = 0; i < NumVertex; i++)
0356   {
0357     xercesc::DOMElement* twoDimVertexElement = NewElement("twoDimVertex");
0358     xtruElement->appendChild(twoDimVertexElement);
0359 
0360     const G4TwoVector& vertex = xtru->GetVertex(i);
0361 
0362     twoDimVertexElement->setAttributeNode(NewAttribute("x", vertex.x() / mm));
0363     twoDimVertexElement->setAttributeNode(NewAttribute("y", vertex.y() / mm));
0364   }
0365 
0366   const G4int NumSection = xtru->GetNofZSections();
0367 
0368   for (G4int i = 0; i < NumSection; i++)
0369   {
0370     xercesc::DOMElement* sectionElement = NewElement("section");
0371     xtruElement->appendChild(sectionElement);
0372 
0373     const G4ExtrudedSolid::ZSection section = xtru->GetZSection(i);
0374 
0375     sectionElement->setAttributeNode(NewAttribute("zOrder", i));
0376     sectionElement->setAttributeNode(NewAttribute("zPosition", section.fZ / mm));
0377     sectionElement->setAttributeNode(NewAttribute("xOffset", section.fOffset.x() / mm));
0378     sectionElement->setAttributeNode(NewAttribute("yOffset", section.fOffset.y() / mm));
0379     sectionElement->setAttributeNode(NewAttribute("scalingFactor", section.fScale));
0380   }
0381 }
0382 
0383 void PHG4GDMLWriteSolids::
0384     HypeWrite(xercesc::DOMElement* solElement, const G4Hype* const hype)
0385 {
0386   const G4String& name = GenerateName(hype->GetName(), hype);
0387 
0388   xercesc::DOMElement* hypeElement = NewElement("hype");
0389   hypeElement->setAttributeNode(NewAttribute("name", name));
0390   hypeElement->setAttributeNode(NewAttribute("rmin",
0391                                              hype->GetInnerRadius() / mm));
0392   hypeElement->setAttributeNode(NewAttribute("rmax",
0393                                              hype->GetOuterRadius() / mm));
0394   hypeElement->setAttributeNode(NewAttribute("inst",
0395                                              hype->GetInnerStereo() / degree));
0396   hypeElement->setAttributeNode(NewAttribute("outst",
0397                                              hype->GetOuterStereo() / degree));
0398   hypeElement->setAttributeNode(NewAttribute("z",
0399                                              2.0 * hype->GetZHalfLength() / mm));
0400   hypeElement->setAttributeNode(NewAttribute("aunit", "deg"));
0401   hypeElement->setAttributeNode(NewAttribute("lunit", "mm"));
0402   solElement->appendChild(hypeElement);
0403 }
0404 
0405 void PHG4GDMLWriteSolids::
0406     OrbWrite(xercesc::DOMElement* solElement, const G4Orb* const orb)
0407 {
0408   const G4String& name = GenerateName(orb->GetName(), orb);
0409 
0410   xercesc::DOMElement* orbElement = NewElement("orb");
0411   orbElement->setAttributeNode(NewAttribute("name", name));
0412   orbElement->setAttributeNode(NewAttribute("r", orb->GetRadius() / mm));
0413   orbElement->setAttributeNode(NewAttribute("lunit", "mm"));
0414   solElement->appendChild(orbElement);
0415 }
0416 
0417 void PHG4GDMLWriteSolids::
0418     ParaWrite(xercesc::DOMElement* solElement, const G4Para* const para)
0419 {
0420   const G4String& name = GenerateName(para->GetName(), para);
0421 
0422   const G4ThreeVector simaxis = para->GetSymAxis();
0423   const G4double alpha = std::atan(para->GetTanAlpha());
0424   const G4double phi = simaxis.phi();
0425   const G4double theta = simaxis.theta();
0426 
0427   xercesc::DOMElement* paraElement = NewElement("para");
0428   paraElement->setAttributeNode(NewAttribute("name", name));
0429   paraElement->setAttributeNode(NewAttribute("x",
0430                                              2.0 * para->GetXHalfLength() / mm));
0431   paraElement->setAttributeNode(NewAttribute("y",
0432                                              2.0 * para->GetYHalfLength() / mm));
0433   paraElement->setAttributeNode(NewAttribute("z",
0434                                              2.0 * para->GetZHalfLength() / mm));
0435   paraElement->setAttributeNode(NewAttribute("alpha", alpha / degree));
0436   paraElement->setAttributeNode(NewAttribute("theta", theta / degree));
0437   paraElement->setAttributeNode(NewAttribute("phi", phi / degree));
0438   paraElement->setAttributeNode(NewAttribute("aunit", "deg"));
0439   paraElement->setAttributeNode(NewAttribute("lunit", "mm"));
0440   solElement->appendChild(paraElement);
0441 }
0442 
0443 void PHG4GDMLWriteSolids::
0444     ParaboloidWrite(xercesc::DOMElement* solElement,
0445                     const G4Paraboloid* const paraboloid)
0446 {
0447   const G4String& name = GenerateName(paraboloid->GetName(), paraboloid);
0448 
0449   xercesc::DOMElement* paraboloidElement = NewElement("paraboloid");
0450   paraboloidElement->setAttributeNode(NewAttribute("name", name));
0451   paraboloidElement->setAttributeNode(NewAttribute("rlo",
0452                                                    paraboloid->GetRadiusMinusZ() / mm));
0453   paraboloidElement->setAttributeNode(NewAttribute("rhi",
0454                                                    paraboloid->GetRadiusPlusZ() / mm));
0455   paraboloidElement->setAttributeNode(NewAttribute("dz",
0456                                                    paraboloid->GetZHalfLength() / mm));
0457   paraboloidElement->setAttributeNode(NewAttribute("lunit", "mm"));
0458   solElement->appendChild(paraboloidElement);
0459 }
0460 void PHG4GDMLWriteSolids::
0461     PolyconeWrite(xercesc::DOMElement* solElement,
0462                   const G4Polycone* const polycone)
0463 {
0464   const G4String& name = GenerateName(polycone->GetName(), polycone);
0465 
0466   xercesc::DOMElement* polyconeElement = NewElement("polycone");
0467   polyconeElement->setAttributeNode(NewAttribute("name", name));
0468   polyconeElement->setAttributeNode(NewAttribute("startphi",
0469                                                  polycone->GetOriginalParameters()->Start_angle / degree));
0470   polyconeElement->setAttributeNode(NewAttribute("deltaphi",
0471                                                  polycone->GetOriginalParameters()->Opening_angle / degree));
0472   polyconeElement->setAttributeNode(NewAttribute("aunit", "deg"));
0473   polyconeElement->setAttributeNode(NewAttribute("lunit", "mm"));
0474   solElement->appendChild(polyconeElement);
0475 
0476   const size_t num_zplanes = polycone->GetOriginalParameters()->Num_z_planes;
0477   const G4double* z_array = polycone->GetOriginalParameters()->Z_values;
0478   const G4double* rmin_array = polycone->GetOriginalParameters()->Rmin;
0479   const G4double* rmax_array = polycone->GetOriginalParameters()->Rmax;
0480 
0481   for (size_t i = 0; i < num_zplanes; i++)
0482   {
0483     ZplaneWrite(polyconeElement, z_array[i], rmin_array[i], rmax_array[i]);
0484   }
0485 }
0486 
0487 void PHG4GDMLWriteSolids::
0488     GenericPolyconeWrite(xercesc::DOMElement* solElement,
0489                          const G4GenericPolycone* const polycone)
0490 {
0491   const G4String& name = GenerateName(polycone->GetName(), polycone);
0492   xercesc::DOMElement* polyconeElement = NewElement("genericPolycone");
0493   const G4double startPhi = polycone->GetStartPhi();
0494   polyconeElement->setAttributeNode(NewAttribute("name", name));
0495   polyconeElement->setAttributeNode(NewAttribute("startphi",
0496                                                  startPhi / degree));
0497   polyconeElement->setAttributeNode(NewAttribute("deltaphi",
0498                                                  (polycone->GetEndPhi() - startPhi) / degree));
0499   polyconeElement->setAttributeNode(NewAttribute("aunit", "deg"));
0500   polyconeElement->setAttributeNode(NewAttribute("lunit", "mm"));
0501   solElement->appendChild(polyconeElement);
0502 
0503   const size_t num_rzpoints = polycone->GetNumRZCorner();
0504   for (size_t i = 0; i < num_rzpoints; i++)
0505   {
0506     const G4double r_point = polycone->GetCorner(i).r;
0507     const G4double z_point = polycone->GetCorner(i).z;
0508     RZPointWrite(polyconeElement, r_point, z_point);
0509   }
0510 }
0511 
0512 void PHG4GDMLWriteSolids::
0513     PolyhedraWrite(xercesc::DOMElement* solElement,
0514                    const G4Polyhedra* const polyhedra)
0515 {
0516   const G4String& name = GenerateName(polyhedra->GetName(), polyhedra);
0517   if (polyhedra->IsGeneric() == false)
0518   {
0519     xercesc::DOMElement* polyhedraElement = NewElement("polyhedra");
0520     polyhedraElement->setAttributeNode(NewAttribute("name", name));
0521     polyhedraElement->setAttributeNode(NewAttribute("startphi",
0522                                                     polyhedra->GetOriginalParameters()->Start_angle / degree));
0523     polyhedraElement->setAttributeNode(NewAttribute("deltaphi",
0524                                                     polyhedra->GetOriginalParameters()->Opening_angle / degree));
0525     polyhedraElement->setAttributeNode(NewAttribute("numsides",
0526                                                     polyhedra->GetOriginalParameters()->numSide));
0527     polyhedraElement->setAttributeNode(NewAttribute("aunit", "deg"));
0528     polyhedraElement->setAttributeNode(NewAttribute("lunit", "mm"));
0529     solElement->appendChild(polyhedraElement);
0530 
0531     const size_t num_zplanes = polyhedra->GetOriginalParameters()->Num_z_planes;
0532     const G4double* z_array = polyhedra->GetOriginalParameters()->Z_values;
0533     const G4double* rmin_array = polyhedra->GetOriginalParameters()->Rmin;
0534     const G4double* rmax_array = polyhedra->GetOriginalParameters()->Rmax;
0535 
0536     const G4double convertRad =
0537         std::cos(0.5 * polyhedra->GetOriginalParameters()->Opening_angle / polyhedra->GetOriginalParameters()->numSide);
0538 
0539     for (size_t i = 0; i < num_zplanes; i++)
0540     {
0541       ZplaneWrite(polyhedraElement, z_array[i],
0542                   rmin_array[i] * convertRad, rmax_array[i] * convertRad);
0543     }
0544   }
0545   else
0546   {  // generic polyhedra
0547 
0548     xercesc::DOMElement* polyhedraElement = NewElement("genericPolyhedra");
0549     polyhedraElement->setAttributeNode(NewAttribute("name", name));
0550     polyhedraElement->setAttributeNode(NewAttribute("startphi",
0551                                                     polyhedra->GetOriginalParameters()->Start_angle / degree));
0552     polyhedraElement->setAttributeNode(NewAttribute("deltaphi",
0553                                                     polyhedra->GetOriginalParameters()->Opening_angle / degree));
0554     polyhedraElement->setAttributeNode(NewAttribute("numsides",
0555                                                     polyhedra->GetOriginalParameters()->numSide));
0556     polyhedraElement->setAttributeNode(NewAttribute("aunit", "deg"));
0557     polyhedraElement->setAttributeNode(NewAttribute("lunit", "mm"));
0558     solElement->appendChild(polyhedraElement);
0559 
0560     const size_t num_rzpoints = polyhedra->GetNumRZCorner();
0561 
0562     for (size_t i = 0; i < num_rzpoints; i++)
0563     {
0564       const G4double r_point = polyhedra->GetCorner(i).r;
0565       const G4double z_point = polyhedra->GetCorner(i).z;
0566       RZPointWrite(polyhedraElement, r_point, z_point);
0567     }
0568   }
0569 }
0570 
0571 void PHG4GDMLWriteSolids::
0572     SphereWrite(xercesc::DOMElement* solElement, const G4Sphere* const sphere)
0573 {
0574   const G4String& name = GenerateName(sphere->GetName(), sphere);
0575 
0576   xercesc::DOMElement* sphereElement = NewElement("sphere");
0577   sphereElement->setAttributeNode(NewAttribute("name", name));
0578   sphereElement->setAttributeNode(NewAttribute("rmin",
0579                                                sphere->GetInnerRadius() / mm));
0580   sphereElement->setAttributeNode(NewAttribute("rmax",
0581                                                sphere->GetOuterRadius() / mm));
0582   sphereElement->setAttributeNode(NewAttribute("startphi",
0583                                                sphere->GetStartPhiAngle() / degree));
0584   sphereElement->setAttributeNode(NewAttribute("deltaphi",
0585                                                sphere->GetDeltaPhiAngle() / degree));
0586   sphereElement->setAttributeNode(NewAttribute("starttheta",
0587                                                sphere->GetStartThetaAngle() / degree));
0588   sphereElement->setAttributeNode(NewAttribute("deltatheta",
0589                                                sphere->GetDeltaThetaAngle() / degree));
0590   sphereElement->setAttributeNode(NewAttribute("aunit", "deg"));
0591   sphereElement->setAttributeNode(NewAttribute("lunit", "mm"));
0592   solElement->appendChild(sphereElement);
0593 }
0594 
0595 void PHG4GDMLWriteSolids::
0596     TessellatedWrite(xercesc::DOMElement* solElement,
0597                      const G4TessellatedSolid* const tessellated)
0598 {
0599   const G4String& solid_name = tessellated->GetName();
0600   const G4String& name = GenerateName(solid_name, tessellated);
0601 
0602   xercesc::DOMElement* tessellatedElement = NewElement("tessellated");
0603   tessellatedElement->setAttributeNode(NewAttribute("name", name));
0604   tessellatedElement->setAttributeNode(NewAttribute("aunit", "deg"));
0605   tessellatedElement->setAttributeNode(NewAttribute("lunit", "mm"));
0606   solElement->appendChild(tessellatedElement);
0607 
0608   std::map<G4ThreeVector, G4String, G4ThreeVectorCompare> vertexMap;
0609 
0610   const size_t NumFacets = tessellated->GetNumberOfFacets();
0611   size_t NumVertex = 0;
0612 
0613   for (size_t i = 0; i < NumFacets; i++)
0614   {
0615     const G4VFacet* facet = tessellated->GetFacet(i);
0616     const size_t NumVertexPerFacet = facet->GetNumberOfVertices();
0617 
0618     G4String FacetTag;
0619 
0620     if (NumVertexPerFacet == 3)
0621     {
0622       FacetTag = "triangular";
0623     }
0624     else
0625     {
0626       if (NumVertexPerFacet == 4)
0627       {
0628         FacetTag = "quadrangular";
0629       }
0630       else
0631       {
0632         G4Exception("PHG4GDMLWriteSolids::TessellatedWrite()", "InvalidSetup",
0633                     FatalException, "Facet should contain 3 or 4 vertices!");
0634       }
0635     }
0636     xercesc::DOMElement* facetElement = NewElement(FacetTag);
0637     tessellatedElement->appendChild(facetElement);
0638 
0639     for (size_t j = 0; j < NumVertexPerFacet; j++)
0640     {
0641       std::stringstream name_stream;
0642       std::stringstream ref_stream;
0643 
0644       name_stream << "vertex" << (j + 1);
0645       ref_stream << solid_name << "_v" << NumVertex;
0646 
0647       const G4String& fname = name_stream.str();  // facet's tag variable
0648       G4String ref = ref_stream.str();            // vertex tag to be associated
0649 
0650       // Now search for the existance of the current vertex in the
0651       // map of cached vertices. If existing, do NOT store it as
0652       // position in the GDML file, so avoiding duplication; otherwise
0653       // cache it in the local map and add it as position in the
0654       // "define" section of the GDML file.
0655 
0656       const G4ThreeVector& vertex = facet->GetVertex(j);
0657 
0658       if (vertexMap.contains(vertex))  // Vertex is cached
0659       {
0660         ref = vertexMap[vertex];  // Set the proper tag for it
0661       }
0662       else  // Vertex not found
0663       {
0664         vertexMap.insert(std::make_pair(vertex, ref));  // Cache vertex and ...
0665         AddPosition(ref, vertex);                       // ... add it to define section!
0666         NumVertex++;
0667       }
0668 
0669       // Now create association of the vertex with its facet
0670       //
0671       facetElement->setAttributeNode(NewAttribute(fname, ref));
0672     }
0673   }
0674 }
0675 
0676 void PHG4GDMLWriteSolids::
0677     TetWrite(xercesc::DOMElement* solElement, const G4Tet* const tet)
0678 {
0679   const G4String& solid_name = tet->GetName();
0680   const G4String& name = GenerateName(solid_name, tet);
0681 
0682   std::vector<G4ThreeVector> vertexList = tet->GetVertices();
0683 
0684   xercesc::DOMElement* tetElement = NewElement("tet");
0685   tetElement->setAttributeNode(NewAttribute("name", name));
0686   tetElement->setAttributeNode(NewAttribute("vertex1", solid_name + "_v1"));
0687   tetElement->setAttributeNode(NewAttribute("vertex2", solid_name + "_v2"));
0688   tetElement->setAttributeNode(NewAttribute("vertex3", solid_name + "_v3"));
0689   tetElement->setAttributeNode(NewAttribute("vertex4", solid_name + "_v4"));
0690   tetElement->setAttributeNode(NewAttribute("lunit", "mm"));
0691   solElement->appendChild(tetElement);
0692 
0693   AddPosition(solid_name + "_v1", vertexList[0]);
0694   AddPosition(solid_name + "_v2", vertexList[1]);
0695   AddPosition(solid_name + "_v3", vertexList[2]);
0696   AddPosition(solid_name + "_v4", vertexList[3]);
0697 }
0698 
0699 void PHG4GDMLWriteSolids::
0700     TorusWrite(xercesc::DOMElement* solElement, const G4Torus* const torus)
0701 {
0702   const G4String& name = GenerateName(torus->GetName(), torus);
0703 
0704   xercesc::DOMElement* torusElement = NewElement("torus");
0705   torusElement->setAttributeNode(NewAttribute("name", name));
0706   torusElement->setAttributeNode(NewAttribute("rmin", torus->GetRmin() / mm));
0707   torusElement->setAttributeNode(NewAttribute("rmax", torus->GetRmax() / mm));
0708   torusElement->setAttributeNode(NewAttribute("rtor", torus->GetRtor() / mm));
0709   torusElement->setAttributeNode(NewAttribute("startphi", torus->GetSPhi() / degree));
0710   torusElement->setAttributeNode(NewAttribute("deltaphi", torus->GetDPhi() / degree));
0711   torusElement->setAttributeNode(NewAttribute("aunit", "deg"));
0712   torusElement->setAttributeNode(NewAttribute("lunit", "mm"));
0713   solElement->appendChild(torusElement);
0714 }
0715 
0716 void PHG4GDMLWriteSolids::
0717     GenTrapWrite(xercesc::DOMElement* solElement,
0718                  const G4GenericTrap* const gtrap)
0719 {
0720   const G4String& name = GenerateName(gtrap->GetName(), gtrap);
0721 
0722   const std::vector<G4TwoVector>& vertices = gtrap->GetVertices();
0723 
0724   xercesc::DOMElement* gtrapElement = NewElement("arb8");
0725   gtrapElement->setAttributeNode(NewAttribute("name", name));
0726   gtrapElement->setAttributeNode(NewAttribute("dz",
0727                                               gtrap->GetZHalfLength() / mm));
0728   gtrapElement->setAttributeNode(NewAttribute("v1x", vertices[0].x()));
0729   gtrapElement->setAttributeNode(NewAttribute("v1y", vertices[0].y()));
0730   gtrapElement->setAttributeNode(NewAttribute("v2x", vertices[1].x()));
0731   gtrapElement->setAttributeNode(NewAttribute("v2y", vertices[1].y()));
0732   gtrapElement->setAttributeNode(NewAttribute("v3x", vertices[2].x()));
0733   gtrapElement->setAttributeNode(NewAttribute("v3y", vertices[2].y()));
0734   gtrapElement->setAttributeNode(NewAttribute("v4x", vertices[3].x()));
0735   gtrapElement->setAttributeNode(NewAttribute("v4y", vertices[3].y()));
0736   gtrapElement->setAttributeNode(NewAttribute("v5x", vertices[4].x()));
0737   gtrapElement->setAttributeNode(NewAttribute("v5y", vertices[4].y()));
0738   gtrapElement->setAttributeNode(NewAttribute("v6x", vertices[5].x()));
0739   gtrapElement->setAttributeNode(NewAttribute("v6y", vertices[5].y()));
0740   gtrapElement->setAttributeNode(NewAttribute("v7x", vertices[6].x()));
0741   gtrapElement->setAttributeNode(NewAttribute("v7y", vertices[6].y()));
0742   gtrapElement->setAttributeNode(NewAttribute("v8x", vertices[7].x()));
0743   gtrapElement->setAttributeNode(NewAttribute("v8y", vertices[7].y()));
0744   gtrapElement->setAttributeNode(NewAttribute("lunit", "mm"));
0745   solElement->appendChild(gtrapElement);
0746 }
0747 
0748 void PHG4GDMLWriteSolids::
0749     TrapWrite(xercesc::DOMElement* solElement, const G4Trap* const trap)
0750 {
0751   const G4String& name = GenerateName(trap->GetName(), trap);
0752 
0753   const G4ThreeVector& simaxis = trap->GetSymAxis();
0754   const G4double phi = simaxis.phi();
0755   const G4double theta = simaxis.theta();
0756   const G4double alpha1 = std::atan(trap->GetTanAlpha1());
0757   const G4double alpha2 = std::atan(trap->GetTanAlpha2());
0758 
0759   xercesc::DOMElement* trapElement = NewElement("trap");
0760   trapElement->setAttributeNode(NewAttribute("name", name));
0761   trapElement->setAttributeNode(NewAttribute("z",
0762                                              2.0 * trap->GetZHalfLength() / mm));
0763   trapElement->setAttributeNode(NewAttribute("theta", theta / degree));
0764   trapElement->setAttributeNode(NewAttribute("phi", phi / degree));
0765   trapElement->setAttributeNode(NewAttribute("y1",
0766                                              2.0 * trap->GetYHalfLength1() / mm));
0767   trapElement->setAttributeNode(NewAttribute("x1",
0768                                              2.0 * trap->GetXHalfLength1() / mm));
0769   trapElement->setAttributeNode(NewAttribute("x2",
0770                                              2.0 * trap->GetXHalfLength2() / mm));
0771   trapElement->setAttributeNode(NewAttribute("alpha1", alpha1 / degree));
0772   trapElement->setAttributeNode(NewAttribute("y2",
0773                                              2.0 * trap->GetYHalfLength2() / mm));
0774   trapElement->setAttributeNode(NewAttribute("x3",
0775                                              2.0 * trap->GetXHalfLength3() / mm));
0776   trapElement->setAttributeNode(NewAttribute("x4",
0777                                              2.0 * trap->GetXHalfLength4() / mm));
0778   trapElement->setAttributeNode(NewAttribute("alpha2", alpha2 / degree));
0779   trapElement->setAttributeNode(NewAttribute("aunit", "deg"));
0780   trapElement->setAttributeNode(NewAttribute("lunit", "mm"));
0781   solElement->appendChild(trapElement);
0782 }
0783 
0784 void PHG4GDMLWriteSolids::
0785     TrdWrite(xercesc::DOMElement* solElement, const G4Trd* const trd)
0786 {
0787   const G4String& name = GenerateName(trd->GetName(), trd);
0788 
0789   xercesc::DOMElement* trdElement = NewElement("trd");
0790   trdElement->setAttributeNode(NewAttribute("name", name));
0791   trdElement->setAttributeNode(NewAttribute("x1",
0792                                             2.0 * trd->GetXHalfLength1() / mm));
0793   trdElement->setAttributeNode(NewAttribute("x2",
0794                                             2.0 * trd->GetXHalfLength2() / mm));
0795   trdElement->setAttributeNode(NewAttribute("y1",
0796                                             2.0 * trd->GetYHalfLength1() / mm));
0797   trdElement->setAttributeNode(NewAttribute("y2",
0798                                             2.0 * trd->GetYHalfLength2() / mm));
0799   trdElement->setAttributeNode(NewAttribute("z",
0800                                             2.0 * trd->GetZHalfLength() / mm));
0801   trdElement->setAttributeNode(NewAttribute("lunit", "mm"));
0802   solElement->appendChild(trdElement);
0803 }
0804 
0805 void PHG4GDMLWriteSolids::
0806     TubeWrite(xercesc::DOMElement* solElement, const G4Tubs* const tube)
0807 {
0808   const G4String& name = GenerateName(tube->GetName(), tube);
0809 
0810   xercesc::DOMElement* tubeElement = NewElement("tube");
0811   tubeElement->setAttributeNode(NewAttribute("name", name));
0812   tubeElement->setAttributeNode(NewAttribute("rmin",
0813                                              tube->GetInnerRadius() / mm));
0814   tubeElement->setAttributeNode(NewAttribute("rmax",
0815                                              tube->GetOuterRadius() / mm));
0816   tubeElement->setAttributeNode(NewAttribute("z",
0817                                              2.0 * tube->GetZHalfLength() / mm));
0818   tubeElement->setAttributeNode(NewAttribute("startphi",
0819                                              tube->GetStartPhiAngle() / degree));
0820   tubeElement->setAttributeNode(NewAttribute("deltaphi",
0821                                              tube->GetDeltaPhiAngle() / degree));
0822   tubeElement->setAttributeNode(NewAttribute("aunit", "deg"));
0823   tubeElement->setAttributeNode(NewAttribute("lunit", "mm"));
0824   solElement->appendChild(tubeElement);
0825 }
0826 
0827 void PHG4GDMLWriteSolids::
0828     CutTubeWrite(xercesc::DOMElement* solElement, const G4CutTubs* const cuttube)
0829 {
0830   const G4String& name = GenerateName(cuttube->GetName(), cuttube);
0831 
0832   xercesc::DOMElement* cuttubeElement = NewElement("cutTube");
0833   cuttubeElement->setAttributeNode(NewAttribute("name", name));
0834   cuttubeElement->setAttributeNode(NewAttribute("rmin",
0835                                                 cuttube->GetInnerRadius() / mm));
0836   cuttubeElement->setAttributeNode(NewAttribute("rmax",
0837                                                 cuttube->GetOuterRadius() / mm));
0838   cuttubeElement->setAttributeNode(NewAttribute("z",
0839                                                 2.0 * cuttube->GetZHalfLength() / mm));
0840   cuttubeElement->setAttributeNode(NewAttribute("startphi",
0841                                                 cuttube->GetStartPhiAngle() / degree));
0842   cuttubeElement->setAttributeNode(NewAttribute("deltaphi",
0843                                                 cuttube->GetDeltaPhiAngle() / degree));
0844   cuttubeElement->setAttributeNode(NewAttribute("lowX",
0845                                                 cuttube->GetLowNorm().getX() / mm));
0846   cuttubeElement->setAttributeNode(NewAttribute("lowY",
0847                                                 cuttube->GetLowNorm().getY() / mm));
0848   cuttubeElement->setAttributeNode(NewAttribute("lowZ",
0849                                                 cuttube->GetLowNorm().getZ() / mm));
0850   cuttubeElement->setAttributeNode(NewAttribute("highX",
0851                                                 cuttube->GetHighNorm().getX() / mm));
0852   cuttubeElement->setAttributeNode(NewAttribute("highY",
0853                                                 cuttube->GetHighNorm().getY() / mm));
0854   cuttubeElement->setAttributeNode(NewAttribute("highZ",
0855                                                 cuttube->GetHighNorm().getZ() / mm));
0856   cuttubeElement->setAttributeNode(NewAttribute("aunit", "deg"));
0857   cuttubeElement->setAttributeNode(NewAttribute("lunit", "mm"));
0858   solElement->appendChild(cuttubeElement);
0859 }
0860 
0861 void PHG4GDMLWriteSolids::
0862     TwistedboxWrite(xercesc::DOMElement* solElement,
0863                     const G4TwistedBox* const twistedbox)
0864 {
0865   const G4String& name = GenerateName(twistedbox->GetName(), twistedbox);
0866 
0867   xercesc::DOMElement* twistedboxElement = NewElement("twistedbox");
0868   twistedboxElement->setAttributeNode(NewAttribute("name", name));
0869   twistedboxElement->setAttributeNode(NewAttribute("x",
0870                                                    2.0 * twistedbox->GetXHalfLength() / mm));
0871   twistedboxElement->setAttributeNode(NewAttribute("y",
0872                                                    2.0 * twistedbox->GetYHalfLength() / mm));
0873   twistedboxElement->setAttributeNode(NewAttribute("z",
0874                                                    2.0 * twistedbox->GetZHalfLength() / mm));
0875   twistedboxElement->setAttributeNode(NewAttribute("PhiTwist",
0876                                                    twistedbox->GetPhiTwist() / degree));
0877   twistedboxElement->setAttributeNode(NewAttribute("aunit", "deg"));
0878   twistedboxElement->setAttributeNode(NewAttribute("lunit", "mm"));
0879   solElement->appendChild(twistedboxElement);
0880 }
0881 
0882 void PHG4GDMLWriteSolids::
0883     TwistedtrapWrite(xercesc::DOMElement* solElement,
0884                      const G4TwistedTrap* const twistedtrap)
0885 {
0886   const G4String& name = GenerateName(twistedtrap->GetName(), twistedtrap);
0887 
0888   xercesc::DOMElement* twistedtrapElement = NewElement("twistedtrap");
0889   twistedtrapElement->setAttributeNode(NewAttribute("name", name));
0890   twistedtrapElement->setAttributeNode(NewAttribute("y1",
0891                                                     2.0 * twistedtrap->GetY1HalfLength() / mm));
0892   twistedtrapElement->setAttributeNode(NewAttribute("x1",
0893                                                     2.0 * twistedtrap->GetX1HalfLength() / mm));
0894   twistedtrapElement->setAttributeNode(NewAttribute("x2",
0895                                                     2.0 * twistedtrap->GetX2HalfLength() / mm));
0896   twistedtrapElement->setAttributeNode(NewAttribute("y2",
0897                                                     2.0 * twistedtrap->GetY2HalfLength() / mm));
0898   twistedtrapElement->setAttributeNode(NewAttribute("x3",
0899                                                     2.0 * twistedtrap->GetX3HalfLength() / mm));
0900   twistedtrapElement->setAttributeNode(NewAttribute("x4",
0901                                                     2.0 * twistedtrap->GetX4HalfLength() / mm));
0902   twistedtrapElement->setAttributeNode(NewAttribute("z",
0903                                                     2.0 * twistedtrap->GetZHalfLength() / mm));
0904   twistedtrapElement->setAttributeNode(NewAttribute("Alph",
0905                                                     twistedtrap->GetTiltAngleAlpha() / degree));
0906   twistedtrapElement->setAttributeNode(NewAttribute("Theta",
0907                                                     twistedtrap->GetPolarAngleTheta() / degree));
0908   twistedtrapElement->setAttributeNode(NewAttribute("Phi",
0909                                                     twistedtrap->GetAzimuthalAnglePhi() / degree));
0910   twistedtrapElement->setAttributeNode(NewAttribute("PhiTwist",
0911                                                     twistedtrap->GetPhiTwist() / degree));
0912   twistedtrapElement->setAttributeNode(NewAttribute("aunit", "deg"));
0913   twistedtrapElement->setAttributeNode(NewAttribute("lunit", "mm"));
0914 
0915   solElement->appendChild(twistedtrapElement);
0916 }
0917 
0918 void PHG4GDMLWriteSolids::
0919     TwistedtrdWrite(xercesc::DOMElement* solElement,
0920                     const G4TwistedTrd* const twistedtrd)
0921 {
0922   const G4String& name = GenerateName(twistedtrd->GetName(), twistedtrd);
0923 
0924   xercesc::DOMElement* twistedtrdElement = NewElement("twistedtrd");
0925   twistedtrdElement->setAttributeNode(NewAttribute("name", name));
0926   twistedtrdElement->setAttributeNode(NewAttribute("x1",
0927                                                    2.0 * twistedtrd->GetX1HalfLength() / mm));
0928   twistedtrdElement->setAttributeNode(NewAttribute("x2",
0929                                                    2.0 * twistedtrd->GetX2HalfLength() / mm));
0930   twistedtrdElement->setAttributeNode(NewAttribute("y1",
0931                                                    2.0 * twistedtrd->GetY1HalfLength() / mm));
0932   twistedtrdElement->setAttributeNode(NewAttribute("y2",
0933                                                    2.0 * twistedtrd->GetY2HalfLength() / mm));
0934   twistedtrdElement->setAttributeNode(NewAttribute("z",
0935                                                    2.0 * twistedtrd->GetZHalfLength() / mm));
0936   twistedtrdElement->setAttributeNode(NewAttribute("PhiTwist",
0937                                                    twistedtrd->GetPhiTwist() / degree));
0938   twistedtrdElement->setAttributeNode(NewAttribute("aunit", "deg"));
0939   twistedtrdElement->setAttributeNode(NewAttribute("lunit", "mm"));
0940   solElement->appendChild(twistedtrdElement);
0941 }
0942 
0943 void PHG4GDMLWriteSolids::
0944     TwistedtubsWrite(xercesc::DOMElement* solElement,
0945                      const G4TwistedTubs* const twistedtubs)
0946 {
0947   const G4String& name = GenerateName(twistedtubs->GetName(), twistedtubs);
0948 
0949   xercesc::DOMElement* twistedtubsElement = NewElement("twistedtubs");
0950   twistedtubsElement->setAttributeNode(NewAttribute("name", name));
0951   twistedtubsElement->setAttributeNode(NewAttribute("twistedangle",
0952                                                     twistedtubs->GetPhiTwist() / degree));
0953   twistedtubsElement->setAttributeNode(NewAttribute("endinnerrad",
0954                                                     twistedtubs->GetInnerRadius() / mm));
0955   twistedtubsElement->setAttributeNode(NewAttribute("endouterrad",
0956                                                     twistedtubs->GetOuterRadius() / mm));
0957   twistedtubsElement->setAttributeNode(NewAttribute("zlen",
0958                                                     2.0 * twistedtubs->GetZHalfLength() / mm));
0959   twistedtubsElement->setAttributeNode(NewAttribute("phi",
0960                                                     twistedtubs->GetDPhi() / degree));
0961   twistedtubsElement->setAttributeNode(NewAttribute("aunit", "deg"));
0962   twistedtubsElement->setAttributeNode(NewAttribute("lunit", "mm"));
0963   solElement->appendChild(twistedtubsElement);
0964 }
0965 
0966 void PHG4GDMLWriteSolids::
0967     ZplaneWrite(xercesc::DOMElement* element, const G4double& z,
0968                 const G4double& rmin, const G4double& rmax)
0969 {
0970   xercesc::DOMElement* zplaneElement = NewElement("zplane");
0971   zplaneElement->setAttributeNode(NewAttribute("z", z / mm));
0972   zplaneElement->setAttributeNode(NewAttribute("rmin", rmin / mm));
0973   zplaneElement->setAttributeNode(NewAttribute("rmax", rmax / mm));
0974   element->appendChild(zplaneElement);
0975 }
0976 
0977 void PHG4GDMLWriteSolids::
0978     RZPointWrite(xercesc::DOMElement* element, const G4double& r,
0979                  const G4double& z)
0980 {
0981   xercesc::DOMElement* rzpointElement = NewElement("rzpoint");
0982   rzpointElement->setAttributeNode(NewAttribute("r", r / mm));
0983   rzpointElement->setAttributeNode(NewAttribute("z", z / mm));
0984   element->appendChild(rzpointElement);
0985 }
0986 
0987 void PHG4GDMLWriteSolids::
0988     OpticalSurfaceWrite(xercesc::DOMElement* solElement,
0989                         const G4OpticalSurface* const surf)
0990 {
0991   xercesc::DOMElement* optElement = NewElement("opticalsurface");
0992   G4OpticalSurfaceModel smodel = surf->GetModel();
0993   G4double sval = (smodel == glisur) ? surf->GetPolish() : surf->GetSigmaAlpha();
0994 
0995   optElement->setAttributeNode(NewAttribute("name", surf->GetName()));
0996   optElement->setAttributeNode(NewAttribute("model", smodel));
0997   optElement->setAttributeNode(NewAttribute("finish", surf->GetFinish()));
0998   optElement->setAttributeNode(NewAttribute("type", surf->GetType()));
0999   optElement->setAttributeNode(NewAttribute("value", sval));
1000 
1001   solElement->appendChild(optElement);
1002 }
1003 
1004 void PHG4GDMLWriteSolids::SolidsWrite(xercesc::DOMElement* gdmlElement)
1005 {
1006   std::cout << "PHG4GDML: Writing solids..." << std::endl;
1007 
1008   solidsElement = NewElement("solids");
1009   gdmlElement->appendChild(solidsElement);
1010 
1011   solidList.clear();
1012 }
1013 
1014 void PHG4GDMLWriteSolids::AddSolid(const G4VSolid* const solidPtr)
1015 {
1016   for (auto& i : solidList)  // Check if solid is
1017   {                          // already in the list!
1018     if (i == solidPtr)
1019     {
1020       return;
1021     }
1022   }
1023 
1024   solidList.push_back(solidPtr);
1025 
1026   if (const G4BooleanSolid* const booleanPtr = dynamic_cast<const G4BooleanSolid*>(solidPtr))
1027   {
1028     BooleanWrite(solidsElement, booleanPtr);
1029   }
1030   else if (solidPtr->GetEntityType() == "G4MultiUnion")
1031   {
1032     const G4MultiUnion* const munionPtr = static_cast<const G4MultiUnion*>(solidPtr);
1033     MultiUnionWrite(solidsElement, munionPtr);
1034   }
1035   else if (solidPtr->GetEntityType() == "G4Box")
1036   {
1037     const G4Box* const boxPtr = static_cast<const G4Box*>(solidPtr);
1038     BoxWrite(solidsElement, boxPtr);
1039   }
1040   else if (solidPtr->GetEntityType() == "G4Cons")
1041   {
1042     const G4Cons* const conePtr = static_cast<const G4Cons*>(solidPtr);
1043     ConeWrite(solidsElement, conePtr);
1044   }
1045   else if (solidPtr->GetEntityType() == "G4EllipticalCone")
1046   {
1047     const G4EllipticalCone* const elconePtr = static_cast<const G4EllipticalCone*>(solidPtr);
1048     ElconeWrite(solidsElement, elconePtr);
1049   }
1050   else if (solidPtr->GetEntityType() == "G4Ellipsoid")
1051   {
1052     const G4Ellipsoid* const ellipsoidPtr = static_cast<const G4Ellipsoid*>(solidPtr);
1053     EllipsoidWrite(solidsElement, ellipsoidPtr);
1054   }
1055   else if (solidPtr->GetEntityType() == "G4EllipticalTube")
1056   {
1057     const G4EllipticalTube* const eltubePtr = static_cast<const G4EllipticalTube*>(solidPtr);
1058     EltubeWrite(solidsElement, eltubePtr);
1059   }
1060   else if (solidPtr->GetEntityType() == "G4ExtrudedSolid")
1061   {
1062     const G4ExtrudedSolid* const xtruPtr = static_cast<const G4ExtrudedSolid*>(solidPtr);
1063     XtruWrite(solidsElement, xtruPtr);
1064   }
1065   else if (solidPtr->GetEntityType() == "G4Hype")
1066   {
1067     const G4Hype* const hypePtr = static_cast<const G4Hype*>(solidPtr);
1068     HypeWrite(solidsElement, hypePtr);
1069   }
1070   else if (solidPtr->GetEntityType() == "G4Orb")
1071   {
1072     const G4Orb* const orbPtr = static_cast<const G4Orb*>(solidPtr);
1073     OrbWrite(solidsElement, orbPtr);
1074   }
1075   else if (solidPtr->GetEntityType() == "G4Para")
1076   {
1077     const G4Para* const paraPtr = static_cast<const G4Para*>(solidPtr);
1078     ParaWrite(solidsElement, paraPtr);
1079   }
1080   else if (solidPtr->GetEntityType() == "G4Paraboloid")
1081   {
1082     const G4Paraboloid* const paraboloidPtr = static_cast<const G4Paraboloid*>(solidPtr);
1083     ParaboloidWrite(solidsElement, paraboloidPtr);
1084   }
1085   else if (solidPtr->GetEntityType() == "G4Polycone")
1086   {
1087     const G4Polycone* const polyconePtr = static_cast<const G4Polycone*>(solidPtr);
1088     PolyconeWrite(solidsElement, polyconePtr);
1089   }
1090   else if (solidPtr->GetEntityType() == "G4GenericPolycone")
1091   {
1092     const G4GenericPolycone* const genpolyconePtr = static_cast<const G4GenericPolycone*>(solidPtr);
1093     GenericPolyconeWrite(solidsElement, genpolyconePtr);
1094   }
1095   else if (solidPtr->GetEntityType() == "G4Polyhedra")
1096   {
1097     const G4Polyhedra* const polyhedraPtr = static_cast<const G4Polyhedra*>(solidPtr);
1098     PolyhedraWrite(solidsElement, polyhedraPtr);
1099   }
1100   else if (solidPtr->GetEntityType() == "G4Sphere")
1101   {
1102     const G4Sphere* const spherePtr = static_cast<const G4Sphere*>(solidPtr);
1103     SphereWrite(solidsElement, spherePtr);
1104   }
1105   else if (solidPtr->GetEntityType() == "G4TessellatedSolid")
1106   {
1107     const G4TessellatedSolid* const tessellatedPtr = static_cast<const G4TessellatedSolid*>(solidPtr);
1108     TessellatedWrite(solidsElement, tessellatedPtr);
1109   }
1110   else if (solidPtr->GetEntityType() == "G4Tet")
1111   {
1112     const G4Tet* const tetPtr = static_cast<const G4Tet*>(solidPtr);
1113     TetWrite(solidsElement, tetPtr);
1114   }
1115   else if (solidPtr->GetEntityType() == "G4Torus")
1116   {
1117     const G4Torus* const torusPtr = static_cast<const G4Torus*>(solidPtr);
1118     TorusWrite(solidsElement, torusPtr);
1119   }
1120   else if (solidPtr->GetEntityType() == "G4GenericTrap")
1121   {
1122     const G4GenericTrap* const gtrapPtr = static_cast<const G4GenericTrap*>(solidPtr);
1123     GenTrapWrite(solidsElement, gtrapPtr);
1124   }
1125   else if (solidPtr->GetEntityType() == "G4Trap")
1126   {
1127     const G4Trap* const trapPtr = static_cast<const G4Trap*>(solidPtr);
1128     TrapWrite(solidsElement, trapPtr);
1129   }
1130   else if (solidPtr->GetEntityType() == "G4Trd")
1131   {
1132     const G4Trd* const trdPtr = static_cast<const G4Trd*>(solidPtr);
1133     TrdWrite(solidsElement, trdPtr);
1134   }
1135   else if (solidPtr->GetEntityType() == "G4Tubs")
1136   {
1137     const G4Tubs* const tubePtr = static_cast<const G4Tubs*>(solidPtr);
1138     TubeWrite(solidsElement, tubePtr);
1139   }
1140   else if (solidPtr->GetEntityType() == "G4CutTubs")
1141   {
1142     const G4CutTubs* const cuttubePtr = static_cast<const G4CutTubs*>(solidPtr);
1143     CutTubeWrite(solidsElement, cuttubePtr);
1144   }
1145   else if (solidPtr->GetEntityType() == "G4TwistedBox")
1146   {
1147     const G4TwistedBox* const twistedboxPtr = static_cast<const G4TwistedBox*>(solidPtr);
1148     TwistedboxWrite(solidsElement, twistedboxPtr);
1149   }
1150   else if (solidPtr->GetEntityType() == "G4TwistedTrap")
1151   {
1152     const G4TwistedTrap* const twistedtrapPtr = static_cast<const G4TwistedTrap*>(solidPtr);
1153     TwistedtrapWrite(solidsElement, twistedtrapPtr);
1154   }
1155   else if (solidPtr->GetEntityType() == "G4TwistedTrd")
1156   {
1157     const G4TwistedTrd* const twistedtrdPtr = static_cast<const G4TwistedTrd*>(solidPtr);
1158     TwistedtrdWrite(solidsElement, twistedtrdPtr);
1159   }
1160   else if (solidPtr->GetEntityType() == "G4TwistedTubs")
1161   {
1162     const G4TwistedTubs* const twistedtubsPtr = static_cast<const G4TwistedTubs*>(solidPtr);
1163     TwistedtubsWrite(solidsElement, twistedtubsPtr);
1164   }
1165   else
1166   {
1167     G4String error_msg = "Unknown solid: " + solidPtr->GetName() + "; Type: " + solidPtr->GetEntityType();
1168     G4Exception("PHG4GDMLWriteSolids::AddSolid()", "WriteError",
1169                 FatalException, error_msg);
1170   }
1171 }