Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes
G4GDMLWriteSolids Class Reference

#include <G4GDMLWriteSolids.hh>

Inheritance diagram for G4GDMLWriteSolids:
G4GDMLWriteMaterials G4GDMLWriteDefine G4GDMLWrite G4GDMLWriteSetup G4GDMLWriteParamvol G4GDMLWriteStructure G03ColorWriter

Public Member Functions

virtual void AddSolid (const G4VSolid *const)
 
virtual void SolidsWrite (xercesc::DOMElement *)
 
- Public Member Functions inherited from G4GDMLWriteMaterials
void AddIsotope (const G4Isotope *const)
 
void AddElement (const G4Element *const)
 
void AddMaterial (const G4Material *const)
 
virtual void MaterialsWrite (xercesc::DOMElement *)
 
- Public Member Functions inherited from G4GDMLWriteDefine
G4ThreeVector GetAngles (const G4RotationMatrix &)
 
void ScaleWrite (xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &scl)
 
void RotationWrite (xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &rot)
 
void PositionWrite (xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &pos)
 
void FirstrotationWrite (xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &rot)
 
void FirstpositionWrite (xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &pos)
 
void AddPosition (const G4String &name, const G4ThreeVector &pos)
 
virtual void DefineWrite (xercesc::DOMElement *)
 
- Public Member Functions inherited from G4GDMLWrite
G4Transform3D Write (const G4String &filename, const G4LogicalVolume *const topLog, const G4String &schemaPath, const G4int depth, G4bool storeReferences=true)
 
void AddModule (const G4VPhysicalVolume *const topVol)
 
void AddModule (const G4int depth)
 
virtual void StructureWrite (xercesc::DOMElement *)=0
 
virtual G4Transform3D TraverseVolumeTree (const G4LogicalVolume *const, const G4int)=0
 
virtual void SurfacesWrite ()=0
 
virtual void SetupWrite (xercesc::DOMElement *, const G4LogicalVolume *const)=0
 
virtual void ExtensionWrite (xercesc::DOMElement *)
 
virtual void AddExtension (xercesc::DOMElement *, const G4LogicalVolume *const)
 

Protected Member Functions

 G4GDMLWriteSolids ()
 
virtual ~G4GDMLWriteSolids ()
 
void BooleanWrite (xercesc::DOMElement *, const G4BooleanSolid *const)
 
void BoxWrite (xercesc::DOMElement *, const G4Box *const)
 
void ConeWrite (xercesc::DOMElement *, const G4Cons *const)
 
void ElconeWrite (xercesc::DOMElement *, const G4EllipticalCone *const)
 
void EllipsoidWrite (xercesc::DOMElement *, const G4Ellipsoid *const)
 
void EltubeWrite (xercesc::DOMElement *, const G4EllipticalTube *const)
 
void XtruWrite (xercesc::DOMElement *, const G4ExtrudedSolid *const)
 
void HypeWrite (xercesc::DOMElement *, const G4Hype *const)
 
void OrbWrite (xercesc::DOMElement *, const G4Orb *const)
 
void ParaWrite (xercesc::DOMElement *, const G4Para *const)
 
void ParaboloidWrite (xercesc::DOMElement *, const G4Paraboloid *const)
 
void PolyconeWrite (xercesc::DOMElement *, const G4Polycone *const)
 
void GenericPolyconeWrite (xercesc::DOMElement *, const G4GenericPolycone *const)
 
void PolyhedraWrite (xercesc::DOMElement *, const G4Polyhedra *const)
 
void SphereWrite (xercesc::DOMElement *, const G4Sphere *const)
 
void TessellatedWrite (xercesc::DOMElement *, const G4TessellatedSolid *const)
 
void TetWrite (xercesc::DOMElement *, const G4Tet *const)
 
void TorusWrite (xercesc::DOMElement *, const G4Torus *const)
 
void GenTrapWrite (xercesc::DOMElement *, const G4GenericTrap *const)
 
void TrapWrite (xercesc::DOMElement *, const G4Trap *const)
 
void TrdWrite (xercesc::DOMElement *, const G4Trd *const)
 
void TubeWrite (xercesc::DOMElement *, const G4Tubs *const)
 
void CutTubeWrite (xercesc::DOMElement *, const G4CutTubs *const)
 
void TwistedboxWrite (xercesc::DOMElement *, const G4TwistedBox *const)
 
void TwistedtrapWrite (xercesc::DOMElement *, const G4TwistedTrap *const)
 
void TwistedtrdWrite (xercesc::DOMElement *, const G4TwistedTrd *const)
 
void TwistedtubsWrite (xercesc::DOMElement *, const G4TwistedTubs *const)
 
void ZplaneWrite (xercesc::DOMElement *, const G4double &, const G4double &, const G4double &)
 
void RZPointWrite (xercesc::DOMElement *, const G4double &, const G4double &)
 
void OpticalSurfaceWrite (xercesc::DOMElement *, const G4OpticalSurface *const)
 
- Protected Member Functions inherited from G4GDMLWriteMaterials
 G4GDMLWriteMaterials ()
 
virtual ~G4GDMLWriteMaterials ()
 
void AtomWrite (xercesc::DOMElement *, const G4double &)
 
void DWrite (xercesc::DOMElement *, const G4double &)
 
void PWrite (xercesc::DOMElement *, const G4double &)
 
void TWrite (xercesc::DOMElement *, const G4double &)
 
void MEEWrite (xercesc::DOMElement *, const G4double &)
 
void IsotopeWrite (const G4Isotope *const)
 
void ElementWrite (const G4Element *const)
 
void MaterialWrite (const G4Material *const)
 
void PropertyWrite (xercesc::DOMElement *, const G4Material *const)
 
void PropertyVectorWrite (const G4String &, const G4PhysicsOrderedFreeVector *const)
 
- Protected Member Functions inherited from G4GDMLWriteDefine
 G4GDMLWriteDefine ()
 
virtual ~G4GDMLWriteDefine ()
 
void Scale_vectorWrite (xercesc::DOMElement *, const G4String &, const G4String &, const G4ThreeVector &)
 
void Rotation_vectorWrite (xercesc::DOMElement *, const G4String &, const G4String &, const G4ThreeVector &)
 
void Position_vectorWrite (xercesc::DOMElement *, const G4String &, const G4String &, const G4ThreeVector &)
 
- Protected Member Functions inherited from G4GDMLWrite
 G4GDMLWrite ()
 
virtual ~G4GDMLWrite ()
 
VolumeMapType & VolumeMap ()
 
G4String GenerateName (const G4String &, const void *const)
 
xercesc::DOMAttr * NewAttribute (const G4String &, const G4String &)
 
xercesc::DOMAttr * NewAttribute (const G4String &, const G4double &)
 
xercesc::DOMElement * NewElement (const G4String &)
 
G4String Modularize (const G4VPhysicalVolume *const topvol, const G4int depth)
 
G4bool FileExists (const G4String &) const
 
PhysVolumeMapType & PvolumeMap ()
 
DepthMapType & DepthMap ()
 

Protected Attributes

std::vector< const G4VSolid * > solidList
 
xercesc::DOMElement * solidsElement
 
- Protected Attributes inherited from G4GDMLWriteMaterials
std::vector< const G4Isotope * > isotopeList
 
std::vector< const G4Element * > elementList
 
std::vector< const G4Material * > materialList
 
xercesc::DOMElement * materialsElement
 
- Protected Attributes inherited from G4GDMLWriteDefine
xercesc::DOMElement * defineElement
 
- Protected Attributes inherited from G4GDMLWrite
G4String SchemaLocation
 
xercesc::DOMDocument * doc
 
xercesc::DOMElement * extElement
 
XMLCh tempStr [10000]
 

Static Protected Attributes

static const G4int maxTransforms = 8
 
- Static Protected Attributes inherited from G4GDMLWriteDefine
static const G4double kRelativePrecision = DBL_EPSILON
 
static const G4double kAngularPrecision = DBL_EPSILON
 
static const G4double kLinearPrecision = DBL_EPSILON
 
- Static Protected Attributes inherited from G4GDMLWrite
static G4bool addPointerToName = true
 

Additional Inherited Members

- Static Public Member Functions inherited from G4GDMLWrite
static void SetAddPointerToName (G4bool)
 

Detailed Description

Definition at line 77 of file G4GDMLWriteSolids.hh.

Constructor & Destructor Documentation

G4GDMLWriteSolids::G4GDMLWriteSolids ( )
protected

Definition at line 72 of file G4GDMLWriteSolids.cc.

74 {
75 }
xercesc::DOMElement * solidsElement
G4GDMLWriteSolids::~G4GDMLWriteSolids ( )
protectedvirtual

Definition at line 77 of file G4GDMLWriteSolids.cc.

78 {
79 }

Member Function Documentation

void G4GDMLWriteSolids::AddSolid ( const G4VSolid * const  solidPtr)
virtual

Definition at line 963 of file G4GDMLWriteSolids.cc.

References BooleanWrite(), BoxWrite(), ConeWrite(), CutTubeWrite(), ElconeWrite(), EllipsoidWrite(), EltubeWrite(), FatalException, G4Exception(), GenericPolyconeWrite(), GenTrapWrite(), G4VSolid::GetEntityType(), G4VSolid::GetName(), HypeWrite(), OrbWrite(), ParaboloidWrite(), ParaWrite(), PolyconeWrite(), PolyhedraWrite(), solidList, solidsElement, SphereWrite(), TessellatedWrite(), TetWrite(), TorusWrite(), TrapWrite(), TrdWrite(), TubeWrite(), TwistedboxWrite(), TwistedtrapWrite(), TwistedtrdWrite(), TwistedtubsWrite(), and XtruWrite().

Referenced by BooleanWrite(), and G4GDMLWriteStructure::TraverseVolumeTree().

964 {
965  for (size_t i=0; i<solidList.size(); i++) // Check if solid is
966  { // already in the list!
967  if (solidList[i] == solidPtr) { return; }
968  }
969 
970  solidList.push_back(solidPtr);
971 
972  if (const G4BooleanSolid* const booleanPtr
973  = dynamic_cast<const G4BooleanSolid*>(solidPtr))
974  { BooleanWrite(solidsElement,booleanPtr); } else
975  if (const G4Box* const boxPtr
976  = dynamic_cast<const G4Box*>(solidPtr))
977  { BoxWrite(solidsElement,boxPtr); } else
978  if (const G4Cons* const conePtr
979  = dynamic_cast<const G4Cons*>(solidPtr))
980  { ConeWrite(solidsElement,conePtr); } else
981  if (const G4EllipticalCone* const elconePtr
982  = dynamic_cast<const G4EllipticalCone*>(solidPtr))
983  { ElconeWrite(solidsElement,elconePtr); } else
984  if (const G4Ellipsoid* const ellipsoidPtr
985  = dynamic_cast<const G4Ellipsoid*>(solidPtr))
986  { EllipsoidWrite(solidsElement,ellipsoidPtr); } else
987  if (const G4EllipticalTube* const eltubePtr
988  = dynamic_cast<const G4EllipticalTube*>(solidPtr))
989  { EltubeWrite(solidsElement,eltubePtr); } else
990  if (const G4ExtrudedSolid* const xtruPtr
991  = dynamic_cast<const G4ExtrudedSolid*>(solidPtr))
992  { XtruWrite(solidsElement,xtruPtr); } else
993  if (const G4Hype* const hypePtr
994  = dynamic_cast<const G4Hype*>(solidPtr))
995  { HypeWrite(solidsElement,hypePtr); } else
996  if (const G4Orb* const orbPtr
997  = dynamic_cast<const G4Orb*>(solidPtr))
998  { OrbWrite(solidsElement,orbPtr); } else
999  if (const G4Para* const paraPtr
1000  = dynamic_cast<const G4Para*>(solidPtr))
1001  { ParaWrite(solidsElement,paraPtr); } else
1002  if (const G4Paraboloid* const paraboloidPtr
1003  = dynamic_cast<const G4Paraboloid*>(solidPtr))
1004  { ParaboloidWrite(solidsElement,paraboloidPtr); } else
1005  if (const G4Polycone* const polyconePtr
1006  = dynamic_cast<const G4Polycone*>(solidPtr))
1007  { PolyconeWrite(solidsElement,polyconePtr); } else
1008  if (const G4GenericPolycone* const genpolyconePtr
1009  = dynamic_cast<const G4GenericPolycone*>(solidPtr))
1010  { GenericPolyconeWrite(solidsElement,genpolyconePtr); } else
1011  if (const G4Polyhedra* const polyhedraPtr
1012  = dynamic_cast<const G4Polyhedra*>(solidPtr))
1013  { PolyhedraWrite(solidsElement,polyhedraPtr); } else
1014  if (const G4Sphere* const spherePtr
1015  = dynamic_cast<const G4Sphere*>(solidPtr))
1016  { SphereWrite(solidsElement,spherePtr); } else
1017  if (const G4TessellatedSolid* const tessellatedPtr
1018  = dynamic_cast<const G4TessellatedSolid*>(solidPtr))
1019  { TessellatedWrite(solidsElement,tessellatedPtr); } else
1020  if (const G4Tet* const tetPtr
1021  = dynamic_cast<const G4Tet*>(solidPtr))
1022  { TetWrite(solidsElement,tetPtr); } else
1023  if (const G4Torus* const torusPtr
1024  = dynamic_cast<const G4Torus*>(solidPtr))
1025  { TorusWrite(solidsElement,torusPtr); } else
1026  if (const G4GenericTrap* const gtrapPtr
1027  = dynamic_cast<const G4GenericTrap*>(solidPtr))
1028  { GenTrapWrite(solidsElement,gtrapPtr); } else
1029  if (const G4Trap* const trapPtr
1030  = dynamic_cast<const G4Trap*>(solidPtr))
1031  { TrapWrite(solidsElement,trapPtr); } else
1032  if (const G4Trd* const trdPtr
1033  = dynamic_cast<const G4Trd*>(solidPtr))
1034  { TrdWrite(solidsElement,trdPtr); } else
1035  if (const G4Tubs* const tubePtr
1036  = dynamic_cast<const G4Tubs*>(solidPtr))
1037  { TubeWrite(solidsElement,tubePtr); } else
1038  if (const G4CutTubs* const cuttubePtr
1039  = dynamic_cast<const G4CutTubs*>(solidPtr))
1040  { CutTubeWrite(solidsElement,cuttubePtr); } else
1041  if (const G4TwistedBox* const twistedboxPtr
1042  = dynamic_cast<const G4TwistedBox*>(solidPtr))
1043  { TwistedboxWrite(solidsElement,twistedboxPtr); } else
1044  if (const G4TwistedTrap* const twistedtrapPtr
1045  = dynamic_cast<const G4TwistedTrap*>(solidPtr))
1046  { TwistedtrapWrite(solidsElement,twistedtrapPtr); } else
1047  if (const G4TwistedTrd* const twistedtrdPtr
1048  = dynamic_cast<const G4TwistedTrd*>(solidPtr))
1049  { TwistedtrdWrite(solidsElement,twistedtrdPtr); } else
1050  if (const G4TwistedTubs* const twistedtubsPtr
1051  = dynamic_cast<const G4TwistedTubs*>(solidPtr))
1052  { TwistedtubsWrite(solidsElement,twistedtubsPtr); }
1053  else
1054  {
1055  G4String error_msg = "Unknown solid: " + solidPtr->GetName()
1056  + "; Type: " + solidPtr->GetEntityType();
1057  G4Exception("G4GDMLWriteSolids::AddSolid()", "WriteError",
1058  FatalException, error_msg);
1059  }
1060 }
G4String GetName() const
Definition: G4Para.hh:76
void BoxWrite(xercesc::DOMElement *, const G4Box *const)
void PolyhedraWrite(xercesc::DOMElement *, const G4Polyhedra *const)
void TessellatedWrite(xercesc::DOMElement *, const G4TessellatedSolid *const)
Definition: G4Box.hh:63
Definition: G4Tubs.hh:84
void TrdWrite(xercesc::DOMElement *, const G4Trd *const)
void XtruWrite(xercesc::DOMElement *, const G4ExtrudedSolid *const)
Definition: G4Trd.hh:71
virtual G4GeometryType GetEntityType() const =0
void TwistedtrdWrite(xercesc::DOMElement *, const G4TwistedTrd *const)
void TwistedboxWrite(xercesc::DOMElement *, const G4TwistedBox *const)
void OrbWrite(xercesc::DOMElement *, const G4Orb *const)
void TwistedtubsWrite(xercesc::DOMElement *, const G4TwistedTubs *const)
void TrapWrite(xercesc::DOMElement *, const G4Trap *const)
void ConeWrite(xercesc::DOMElement *, const G4Cons *const)
Definition: G4Tet.hh:65
void ElconeWrite(xercesc::DOMElement *, const G4EllipticalCone *const)
Definition: G4Hype.hh:66
void ParaboloidWrite(xercesc::DOMElement *, const G4Paraboloid *const)
void PolyconeWrite(xercesc::DOMElement *, const G4Polycone *const)
Definition: G4Cons.hh:82
void GenericPolyconeWrite(xercesc::DOMElement *, const G4GenericPolycone *const)
void CutTubeWrite(xercesc::DOMElement *, const G4CutTubs *const)
Definition: G4Orb.hh:60
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
void HypeWrite(xercesc::DOMElement *, const G4Hype *const)
void EltubeWrite(xercesc::DOMElement *, const G4EllipticalTube *const)
void ParaWrite(xercesc::DOMElement *, const G4Para *const)
void EllipsoidWrite(xercesc::DOMElement *, const G4Ellipsoid *const)
void TetWrite(xercesc::DOMElement *, const G4Tet *const)
void TubeWrite(xercesc::DOMElement *, const G4Tubs *const)
void BooleanWrite(xercesc::DOMElement *, const G4BooleanSolid *const)
void TorusWrite(xercesc::DOMElement *, const G4Torus *const)
void TwistedtrapWrite(xercesc::DOMElement *, const G4TwistedTrap *const)
void GenTrapWrite(xercesc::DOMElement *, const G4GenericTrap *const)
xercesc::DOMElement * solidsElement
void SphereWrite(xercesc::DOMElement *, const G4Sphere *const)
std::vector< const G4VSolid * > solidList
void G4GDMLWriteSolids::BooleanWrite ( xercesc::DOMElement *  solElement,
const G4BooleanSolid * const  boolean 
)
protected

Definition at line 82 of file G4GDMLWriteSolids.cc.

References AddSolid(), FatalException, G4GDMLWriteDefine::FirstpositionWrite(), G4GDMLWriteDefine::FirstrotationWrite(), G4Exception(), G4GDMLWrite::GenerateName(), G4GDMLWriteDefine::GetAngles(), G4VSolid::GetName(), G4GDMLWriteDefine::kAngularPrecision, G4GDMLWriteDefine::kLinearPrecision, maxTransforms, G4GDMLWrite::NewAttribute(), G4GDMLWrite::NewElement(), G4GDMLWriteDefine::PositionWrite(), G4GDMLWriteDefine::RotationWrite(), CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Referenced by AddSolid().

84 {
85  G4int displaced=0;
86 
87  G4String tag("undefined");
88  if (dynamic_cast<const G4IntersectionSolid*>(boolean))
89  { tag = "intersection"; } else
90  if (dynamic_cast<const G4SubtractionSolid*>(boolean))
91  { tag = "subtraction"; } else
92  if (dynamic_cast<const G4UnionSolid*>(boolean))
93  { tag = "union"; }
94 
95  G4VSolid* firstPtr = const_cast<G4VSolid*>(boolean->GetConstituentSolid(0));
96  G4VSolid* secondPtr = const_cast<G4VSolid*>(boolean->GetConstituentSolid(1));
97 
98  G4ThreeVector firstpos,firstrot,pos,rot;
99 
100  // Solve possible displacement of referenced solids!
101  //
102  while (true)
103  {
104  if ( displaced>8 )
105  {
106  G4String ErrorMessage = "The referenced solid '"
107  + firstPtr->GetName() +
108  + "in the Boolean shape '" +
109  + boolean->GetName() +
110  + "' was displaced too many times!";
111  G4Exception("G4GDMLWriteSolids::BooleanWrite()",
112  "InvalidSetup", FatalException, ErrorMessage);
113  }
114 
115  if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(firstPtr))
116  {
117  firstpos += disp->GetObjectTranslation();
118  firstrot += firstrot + GetAngles(disp->GetObjectRotation());
119  firstPtr = disp->GetConstituentMovedSolid();
120  displaced++;
121  continue;
122  }
123  break;
124  }
125  displaced = 0;
126  while (true)
127  {
128  if ( displaced>maxTransforms )
129  {
130  G4String ErrorMessage = "The referenced solid '"
131  + secondPtr->GetName() +
132  + "in the Boolean shape '" +
133  + boolean->GetName() +
134  + "' was displaced too many times!";
135  G4Exception("G4GDMLWriteSolids::BooleanWrite()",
136  "InvalidSetup", FatalException, ErrorMessage);
137  }
138 
139  if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(secondPtr))
140  {
141  pos += disp->GetObjectTranslation();
142  rot += GetAngles(disp->GetObjectRotation());
143  secondPtr = disp->GetConstituentMovedSolid();
144  displaced++;
145  continue;
146  }
147  break;
148  }
149 
150  AddSolid(firstPtr); // At first add the constituent solids!
151  AddSolid(secondPtr);
152 
153  const G4String& name = GenerateName(boolean->GetName(),boolean);
154  const G4String& firstref = GenerateName(firstPtr->GetName(),firstPtr);
155  const G4String& secondref = GenerateName(secondPtr->GetName(),secondPtr);
156 
157  xercesc::DOMElement* booleanElement = NewElement(tag);
158  booleanElement->setAttributeNode(NewAttribute("name",name));
159  xercesc::DOMElement* firstElement = NewElement("first");
160  firstElement->setAttributeNode(NewAttribute("ref",firstref));
161  booleanElement->appendChild(firstElement);
162  xercesc::DOMElement* secondElement = NewElement("second");
163  secondElement->setAttributeNode(NewAttribute("ref",secondref));
164  booleanElement->appendChild(secondElement);
165  solElement->appendChild(booleanElement);
166  // Add the boolean solid AFTER the constituent solids!
167 
168  if ( (std::fabs(pos.x()) > kLinearPrecision)
169  || (std::fabs(pos.y()) > kLinearPrecision)
170  || (std::fabs(pos.z()) > kLinearPrecision) )
171  {
172  PositionWrite(booleanElement,name+"_pos",pos);
173  }
174 
175  if ( (std::fabs(rot.x()) > kAngularPrecision)
176  || (std::fabs(rot.y()) > kAngularPrecision)
177  || (std::fabs(rot.z()) > kAngularPrecision) )
178  {
179  RotationWrite(booleanElement,name+"_rot",rot);
180  }
181 
182  if ( (std::fabs(firstpos.x()) > kLinearPrecision)
183  || (std::fabs(firstpos.y()) > kLinearPrecision)
184  || (std::fabs(firstpos.z()) > kLinearPrecision) )
185  {
186  FirstpositionWrite(booleanElement,name+"_fpos",firstpos);
187  }
188 
189  if ( (std::fabs(firstrot.x()) > kAngularPrecision)
190  || (std::fabs(firstrot.y()) > kAngularPrecision)
191  || (std::fabs(firstrot.z()) > kAngularPrecision) )
192  {
193  FirstrotationWrite(booleanElement,name+"_frot",firstrot);
194  }
195 }
static const G4double kAngularPrecision
G4String GetName() const
double x() const
Definition: xmlparse.cc:179
const XML_Char * name
void FirstpositionWrite(xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &pos)
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
int G4int
Definition: G4Types.hh:78
double z() const
static const G4int maxTransforms
void FirstrotationWrite(xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &rot)
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
void RotationWrite(xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &rot)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
static const G4double kLinearPrecision
G4ThreeVector GetAngles(const G4RotationMatrix &)
double y() const
void PositionWrite(xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &pos)
virtual void AddSolid(const G4VSolid *const)
void G4GDMLWriteSolids::BoxWrite ( xercesc::DOMElement *  solElement,
const G4Box * const  box 
)
protected

Definition at line 198 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Box::GetXHalfLength(), G4Box::GetYHalfLength(), G4Box::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

199 {
200  const G4String& name = GenerateName(box->GetName(),box);
201 
202  xercesc::DOMElement* boxElement = NewElement("box");
203  boxElement->setAttributeNode(NewAttribute("name",name));
204  boxElement->setAttributeNode(NewAttribute("x",2.0*box->GetXHalfLength()/mm));
205  boxElement->setAttributeNode(NewAttribute("y",2.0*box->GetYHalfLength()/mm));
206  boxElement->setAttributeNode(NewAttribute("z",2.0*box->GetZHalfLength()/mm));
207  boxElement->setAttributeNode(NewAttribute("lunit","mm"));
208  solElement->appendChild(boxElement);
209 }
G4String GetName() const
G4double GetXHalfLength() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetZHalfLength() const
G4double GetYHalfLength() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::ConeWrite ( xercesc::DOMElement *  solElement,
const G4Cons * const  cone 
)
protected

Definition at line 212 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4Cons::GetDeltaPhiAngle(), G4Cons::GetInnerRadiusMinusZ(), G4Cons::GetInnerRadiusPlusZ(), G4VSolid::GetName(), G4Cons::GetOuterRadiusMinusZ(), G4Cons::GetOuterRadiusPlusZ(), G4Cons::GetStartPhiAngle(), G4Cons::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

213 {
214  const G4String& name = GenerateName(cone->GetName(),cone);
215 
216  xercesc::DOMElement* coneElement = NewElement("cone");
217  coneElement->setAttributeNode(NewAttribute("name",name));
218  coneElement->
219  setAttributeNode(NewAttribute("rmin1",cone->GetInnerRadiusMinusZ()/mm));
220  coneElement->
221  setAttributeNode(NewAttribute("rmax1",cone->GetOuterRadiusMinusZ()/mm));
222  coneElement->
223  setAttributeNode(NewAttribute("rmin2",cone->GetInnerRadiusPlusZ()/mm));
224  coneElement->
225  setAttributeNode(NewAttribute("rmax2",cone->GetOuterRadiusPlusZ()/mm));
226  coneElement->
227  setAttributeNode(NewAttribute("z",2.0*cone->GetZHalfLength()/mm));
228  coneElement->
229  setAttributeNode(NewAttribute("startphi",cone->GetStartPhiAngle()/degree));
230  coneElement->
231  setAttributeNode(NewAttribute("deltaphi",cone->GetDeltaPhiAngle()/degree));
232  coneElement->setAttributeNode(NewAttribute("aunit","deg"));
233  coneElement->setAttributeNode(NewAttribute("lunit","mm"));
234  solElement->appendChild(coneElement);
235 }
G4String GetName() const
const XML_Char * name
G4double GetOuterRadiusMinusZ() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
tuple degree
Definition: hepunit.py:69
G4double GetStartPhiAngle() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetInnerRadiusPlusZ() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetInnerRadiusMinusZ() const
G4double GetOuterRadiusPlusZ() const
G4double GetZHalfLength() const
G4double GetDeltaPhiAngle() const
void G4GDMLWriteSolids::CutTubeWrite ( xercesc::DOMElement *  solElement,
const G4CutTubs * const  cuttube 
)
protected

Definition at line 777 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4OTubs::GetDeltaPhiAngle(), G4CutTubs::GetHighNorm(), G4OTubs::GetInnerRadius(), G4CutTubs::GetLowNorm(), G4VSolid::GetName(), G4OTubs::GetOuterRadius(), G4OTubs::GetStartPhiAngle(), CLHEP::Hep3Vector::getX(), CLHEP::Hep3Vector::getY(), CLHEP::Hep3Vector::getZ(), G4OTubs::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

778 {
779  const G4String& name = GenerateName(cuttube->GetName(),cuttube);
780 
781  xercesc::DOMElement* cuttubeElement = NewElement("cutTube");
782  cuttubeElement->setAttributeNode(NewAttribute("name",name));
783  cuttubeElement->setAttributeNode(NewAttribute("rmin",
784  cuttube->GetInnerRadius()/mm));
785  cuttubeElement->setAttributeNode(NewAttribute("rmax",
786  cuttube->GetOuterRadius()/mm));
787  cuttubeElement->setAttributeNode(NewAttribute("z",
788  2.0*cuttube->GetZHalfLength()/mm));
789  cuttubeElement->setAttributeNode(NewAttribute("startphi",
790  cuttube->GetStartPhiAngle()/degree));
791  cuttubeElement->setAttributeNode(NewAttribute("deltaphi",
792  cuttube->GetDeltaPhiAngle()/degree));
793  cuttubeElement->setAttributeNode(NewAttribute("lowX",
794  cuttube->GetLowNorm().getX()/mm));
795  cuttubeElement->setAttributeNode(NewAttribute("lowY",
796  cuttube->GetLowNorm().getY()/mm));
797  cuttubeElement->setAttributeNode(NewAttribute("lowZ",
798  cuttube->GetLowNorm().getZ()/mm));
799  cuttubeElement->setAttributeNode(NewAttribute("highX",
800  cuttube->GetHighNorm().getX()/mm));
801  cuttubeElement->setAttributeNode(NewAttribute("highY",
802  cuttube->GetHighNorm().getY()/mm));
803  cuttubeElement->setAttributeNode(NewAttribute("highZ",
804  cuttube->GetHighNorm().getZ()/mm));
805  cuttubeElement->setAttributeNode(NewAttribute("aunit","deg"));
806  cuttubeElement->setAttributeNode(NewAttribute("lunit","mm"));
807  solElement->appendChild(cuttubeElement);
808 }
G4String GetName() const
G4ThreeVector GetLowNorm() const
const XML_Char * name
double getY() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetStartPhiAngle() const
double getX() const
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetZHalfLength() const
G4double GetOuterRadius() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
double getZ() const
G4ThreeVector GetHighNorm() const
G4double GetInnerRadius() const
G4double GetDeltaPhiAngle() const
void G4GDMLWriteSolids::ElconeWrite ( xercesc::DOMElement *  solElement,
const G4EllipticalCone * const  elcone 
)
protected

Definition at line 238 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4EllipticalCone::GetSemiAxisX(), G4EllipticalCone::GetSemiAxisY(), G4EllipticalCone::GetZMax(), G4EllipticalCone::GetZTopCut(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

240 {
241  const G4String& name = GenerateName(elcone->GetName(),elcone);
242 
243  xercesc::DOMElement* elconeElement = NewElement("elcone");
244  elconeElement->setAttributeNode(NewAttribute("name",name));
245  elconeElement->setAttributeNode(NewAttribute("dx",elcone->GetSemiAxisX()/mm));
246  elconeElement->setAttributeNode(NewAttribute("dy",elcone->GetSemiAxisY()/mm));
247  elconeElement->setAttributeNode(NewAttribute("zmax",elcone->GetZMax()/mm));
248  elconeElement->setAttributeNode(NewAttribute("zcut",elcone->GetZTopCut()/mm));
249  elconeElement->setAttributeNode(NewAttribute("lunit","mm"));
250  solElement->appendChild(elconeElement);
251 }
G4String GetName() const
G4double GetSemiAxisX() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetZMax() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetSemiAxisY() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetZTopCut() const
void G4GDMLWriteSolids::EllipsoidWrite ( xercesc::DOMElement *  solElement,
const G4Ellipsoid * const  ellipsoid 
)
protected

Definition at line 254 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Ellipsoid::GetSemiAxisMax(), G4Ellipsoid::GetZBottomCut(), G4Ellipsoid::GetZTopCut(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

256 {
257  const G4String& name = GenerateName(ellipsoid->GetName(),ellipsoid);
258 
259  xercesc::DOMElement* ellipsoidElement = NewElement("ellipsoid");
260  ellipsoidElement->setAttributeNode(NewAttribute("name",name));
261  ellipsoidElement->
262  setAttributeNode(NewAttribute("ax",ellipsoid->GetSemiAxisMax(0)/mm));
263  ellipsoidElement->
264  setAttributeNode(NewAttribute("by",ellipsoid->GetSemiAxisMax(1)/mm));
265  ellipsoidElement->
266  setAttributeNode(NewAttribute("cz",ellipsoid->GetSemiAxisMax(2)/mm));
267  ellipsoidElement->
268  setAttributeNode(NewAttribute("zcut1",ellipsoid->GetZBottomCut()/mm));
269  ellipsoidElement->
270  setAttributeNode(NewAttribute("zcut2",ellipsoid->GetZTopCut()/mm));
271  ellipsoidElement->
272  setAttributeNode(NewAttribute("lunit","mm"));
273  solElement->appendChild(ellipsoidElement);
274 }
G4String GetName() const
G4double GetZTopCut() const
G4double GetZBottomCut() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetSemiAxisMax(G4int i) const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::EltubeWrite ( xercesc::DOMElement *  solElement,
const G4EllipticalTube * const  eltube 
)
protected

Definition at line 277 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4EllipticalTube::GetDx(), G4EllipticalTube::GetDy(), G4EllipticalTube::GetDz(), G4VSolid::GetName(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

279 {
280  const G4String& name = GenerateName(eltube->GetName(),eltube);
281 
282  xercesc::DOMElement* eltubeElement = NewElement("eltube");
283  eltubeElement->setAttributeNode(NewAttribute("name",name));
284  eltubeElement->setAttributeNode(NewAttribute("dx",eltube->GetDx()/mm));
285  eltubeElement->setAttributeNode(NewAttribute("dy",eltube->GetDy()/mm));
286  eltubeElement->setAttributeNode(NewAttribute("dz",eltube->GetDz()/mm));
287  eltubeElement->setAttributeNode(NewAttribute("lunit","mm"));
288  solElement->appendChild(eltubeElement);
289 }
G4String GetName() const
G4double GetDy() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetDz() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetDx() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::GenericPolyconeWrite ( xercesc::DOMElement *  solElement,
const G4GenericPolycone * const  polycone 
)
protected

Definition at line 443 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4GenericPolycone::GetCorner(), G4GenericPolycone::GetEndPhi(), G4VSolid::GetName(), G4GenericPolycone::GetNumRZCorner(), G4GenericPolycone::GetStartPhi(), G4GDMLWrite::NewAttribute(), G4GDMLWrite::NewElement(), G4PolyconeSideRZ::r, RZPointWrite(), and G4PolyconeSideRZ::z.

Referenced by AddSolid().

445 {
446  const G4String& name = GenerateName(polycone->GetName(),polycone);
447  xercesc::DOMElement* polyconeElement = NewElement("genericPolycone");
448  const G4double startPhi=polycone->GetStartPhi();
449  polyconeElement->setAttributeNode(NewAttribute("name",name));
450  polyconeElement->setAttributeNode(NewAttribute("startphi",
451  startPhi/degree));
452  polyconeElement->setAttributeNode(NewAttribute("deltaphi",
453  (polycone->GetEndPhi()-startPhi)/degree));
454  polyconeElement->setAttributeNode(NewAttribute("aunit","deg"));
455  polyconeElement->setAttributeNode(NewAttribute("lunit","mm"));
456  solElement->appendChild(polyconeElement);
457 
458  const size_t num_rzpoints = polycone->GetNumRZCorner();
459  for (size_t i=0; i<num_rzpoints; i++)
460  {
461  const G4double r_point=polycone->GetCorner(i).r;
462  const G4double z_point=polycone->GetCorner(i).z;
463  RZPointWrite(polyconeElement,r_point,z_point);
464  }
465 
466 }
G4PolyconeSideRZ GetCorner(G4int index) const
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetEndPhi() const
void RZPointWrite(xercesc::DOMElement *, const G4double &, const G4double &)
tuple degree
Definition: hepunit.py:69
G4double GetStartPhi() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
double G4double
Definition: G4Types.hh:76
G4int GetNumRZCorner() const
void G4GDMLWriteSolids::GenTrapWrite ( xercesc::DOMElement *  solElement,
const G4GenericTrap * const  gtrap 
)
protected

Definition at line 665 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4GenericTrap::GetVertices(), G4GenericTrap::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), G4GDMLWrite::NewElement(), and test::x.

Referenced by AddSolid().

667 {
668  const G4String& name = GenerateName(gtrap->GetName(),gtrap);
669 
670  std::vector<G4TwoVector> vertices = gtrap->GetVertices();
671 
672  xercesc::DOMElement* gtrapElement = NewElement("arb8");
673  gtrapElement->setAttributeNode(NewAttribute("name",name));
674  gtrapElement->setAttributeNode(NewAttribute("dz",
675  gtrap->GetZHalfLength()/mm));
676  gtrapElement->setAttributeNode(NewAttribute("v1x", vertices[0].x()));
677  gtrapElement->setAttributeNode(NewAttribute("v1y", vertices[0].y()));
678  gtrapElement->setAttributeNode(NewAttribute("v2x", vertices[1].x()));
679  gtrapElement->setAttributeNode(NewAttribute("v2y", vertices[1].y()));
680  gtrapElement->setAttributeNode(NewAttribute("v3x", vertices[2].x()));
681  gtrapElement->setAttributeNode(NewAttribute("v3y", vertices[2].y()));
682  gtrapElement->setAttributeNode(NewAttribute("v4x", vertices[3].x()));
683  gtrapElement->setAttributeNode(NewAttribute("v4y", vertices[3].y()));
684  gtrapElement->setAttributeNode(NewAttribute("v5x", vertices[4].x()));
685  gtrapElement->setAttributeNode(NewAttribute("v5y", vertices[4].y()));
686  gtrapElement->setAttributeNode(NewAttribute("v6x", vertices[5].x()));
687  gtrapElement->setAttributeNode(NewAttribute("v6y", vertices[5].y()));
688  gtrapElement->setAttributeNode(NewAttribute("v7x", vertices[6].x()));
689  gtrapElement->setAttributeNode(NewAttribute("v7y", vertices[6].y()));
690  gtrapElement->setAttributeNode(NewAttribute("v8x", vertices[7].x()));
691  gtrapElement->setAttributeNode(NewAttribute("v8y", vertices[7].y()));
692  gtrapElement->setAttributeNode(NewAttribute("lunit","mm"));
693  solElement->appendChild(gtrapElement);
694 }
G4String GetName() const
const XML_Char * name
G4double GetZHalfLength() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
const std::vector< G4TwoVector > & GetVertices() const
void G4GDMLWriteSolids::HypeWrite ( xercesc::DOMElement *  solElement,
const G4Hype * const  hype 
)
protected

Definition at line 336 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4Hype::GetInnerRadius(), G4Hype::GetInnerStereo(), G4VSolid::GetName(), G4Hype::GetOuterRadius(), G4Hype::GetOuterStereo(), G4Hype::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

337 {
338  const G4String& name = GenerateName(hype->GetName(),hype);
339 
340  xercesc::DOMElement* hypeElement = NewElement("hype");
341  hypeElement->setAttributeNode(NewAttribute("name",name));
342  hypeElement->setAttributeNode(NewAttribute("rmin",
343  hype->GetInnerRadius()/mm));
344  hypeElement->setAttributeNode(NewAttribute("rmax",
345  hype->GetOuterRadius()/mm));
346  hypeElement->setAttributeNode(NewAttribute("inst",
347  hype->GetInnerStereo()/degree));
348  hypeElement->setAttributeNode(NewAttribute("outst",
349  hype->GetOuterStereo()/degree));
350  hypeElement->setAttributeNode(NewAttribute("z",
351  2.0*hype->GetZHalfLength()/mm));
352  hypeElement->setAttributeNode(NewAttribute("aunit","deg"));
353  hypeElement->setAttributeNode(NewAttribute("lunit","mm"));
354  solElement->appendChild(hypeElement);
355 }
G4String GetName() const
G4double GetOuterStereo() const
G4double GetInnerStereo() const
G4double GetInnerRadius() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
tuple degree
Definition: hepunit.py:69
G4double GetOuterRadius() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetZHalfLength() const
void G4GDMLWriteSolids::OpticalSurfaceWrite ( xercesc::DOMElement *  solElement,
const G4OpticalSurface * const  surf 
)
protected

Definition at line 937 of file G4GDMLWriteSolids.cc.

References G4OpticalSurface::GetFinish(), G4OpticalSurface::GetModel(), G4SurfaceProperty::GetName(), G4OpticalSurface::GetPolish(), G4OpticalSurface::GetSigmaAlpha(), G4SurfaceProperty::GetType(), glisur, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by G4GDMLWriteStructure::BorderSurfaceCache(), and G4GDMLWriteStructure::SkinSurfaceCache().

939 {
940  xercesc::DOMElement* optElement = NewElement("opticalsurface");
941  G4OpticalSurfaceModel smodel = surf->GetModel();
942  G4double sval = (smodel==glisur) ? surf->GetPolish() : surf->GetSigmaAlpha();
943 
944  optElement->setAttributeNode(NewAttribute("name", surf->GetName()));
945  optElement->setAttributeNode(NewAttribute("model", smodel));
946  optElement->setAttributeNode(NewAttribute("finish", surf->GetFinish()));
947  optElement->setAttributeNode(NewAttribute("type", surf->GetType()));
948  optElement->setAttributeNode(NewAttribute("value", sval));
949 
950  solElement->appendChild(optElement);
951 }
G4OpticalSurfaceModel GetModel() const
const G4String & GetName() const
const G4SurfaceType & GetType() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetPolish() const
G4OpticalSurfaceFinish GetFinish() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4OpticalSurfaceModel
double G4double
Definition: G4Types.hh:76
G4double GetSigmaAlpha() const
void G4GDMLWriteSolids::OrbWrite ( xercesc::DOMElement *  solElement,
const G4Orb * const  orb 
)
protected

Definition at line 358 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Orb::GetRadius(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

359 {
360  const G4String& name = GenerateName(orb->GetName(),orb);
361 
362  xercesc::DOMElement* orbElement = NewElement("orb");
363  orbElement->setAttributeNode(NewAttribute("name",name));
364  orbElement->setAttributeNode(NewAttribute("r",orb->GetRadius()/mm));
365  orbElement->setAttributeNode(NewAttribute("lunit","mm"));
366  solElement->appendChild(orbElement);
367 }
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetRadius() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::ParaboloidWrite ( xercesc::DOMElement *  solElement,
const G4Paraboloid * const  paraboloid 
)
protected

Definition at line 397 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Paraboloid::GetRadiusMinusZ(), G4Paraboloid::GetRadiusPlusZ(), G4Paraboloid::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

399 {
400  const G4String& name = GenerateName(paraboloid->GetName(),paraboloid);
401 
402  xercesc::DOMElement* paraboloidElement = NewElement("paraboloid");
403  paraboloidElement->setAttributeNode(NewAttribute("name",name));
404  paraboloidElement->setAttributeNode(NewAttribute("rlo",
405  paraboloid->GetRadiusMinusZ()/mm));
406  paraboloidElement->setAttributeNode(NewAttribute("rhi",
407  paraboloid->GetRadiusPlusZ()/mm));
408  paraboloidElement->setAttributeNode(NewAttribute("dz",
409  paraboloid->GetZHalfLength()/mm));
410  paraboloidElement->setAttributeNode(NewAttribute("lunit","mm"));
411  solElement->appendChild(paraboloidElement);
412 }
G4String GetName() const
const XML_Char * name
G4double GetRadiusMinusZ() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetZHalfLength() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetRadiusPlusZ() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::ParaWrite ( xercesc::DOMElement *  solElement,
const G4Para * const  para 
)
protected

Definition at line 370 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Para::GetSymAxis(), G4Para::GetTanAlpha(), G4Para::GetXHalfLength(), G4Para::GetYHalfLength(), G4Para::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

371 {
372  const G4String& name = GenerateName(para->GetName(),para);
373 
374  const G4ThreeVector simaxis = para->GetSymAxis();
375  const G4double alpha = std::atan(para->GetTanAlpha());
376  const G4double theta = std::acos(simaxis.z());
377  const G4double phi = (simaxis.z() != 1.0)
378  ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
379 
380  xercesc::DOMElement* paraElement = NewElement("para");
381  paraElement->setAttributeNode(NewAttribute("name",name));
382  paraElement->setAttributeNode(NewAttribute("x",
383  2.0*para->GetXHalfLength()/mm));
384  paraElement->setAttributeNode(NewAttribute("y",
385  2.0*para->GetYHalfLength()/mm));
386  paraElement->setAttributeNode(NewAttribute("z",
387  2.0*para->GetZHalfLength()/mm));
388  paraElement->setAttributeNode(NewAttribute("alpha",alpha/degree));
389  paraElement->setAttributeNode(NewAttribute("theta",theta/degree));
390  paraElement->setAttributeNode(NewAttribute("phi",phi/degree));
391  paraElement->setAttributeNode(NewAttribute("aunit","deg"));
392  paraElement->setAttributeNode(NewAttribute("lunit","mm"));
393  solElement->appendChild(paraElement);
394 }
G4ThreeVector GetSymAxis() const
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetXHalfLength() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetTanAlpha() const
G4double GetZHalfLength() const
double G4double
Definition: G4Types.hh:76
G4double GetYHalfLength() const
void G4GDMLWriteSolids::PolyconeWrite ( xercesc::DOMElement *  solElement,
const G4Polycone * const  polycone 
)
protected

Definition at line 414 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Polycone::GetOriginalParameters(), G4GDMLWrite::NewAttribute(), G4GDMLWrite::NewElement(), G4PolyconeHistorical::Num_z_planes, G4PolyconeHistorical::Opening_angle, G4PolyconeHistorical::Rmax, G4PolyconeHistorical::Rmin, G4PolyconeHistorical::Start_angle, G4PolyconeHistorical::Z_values, and ZplaneWrite().

Referenced by AddSolid().

416 {
417  const G4String& name = GenerateName(polycone->GetName(),polycone);
418 
419  xercesc::DOMElement* polyconeElement = NewElement("polycone");
420  polyconeElement->setAttributeNode(NewAttribute("name",name));
421  polyconeElement->setAttributeNode(NewAttribute("startphi",
423  polyconeElement->setAttributeNode(NewAttribute("deltaphi",
425  polyconeElement->setAttributeNode(NewAttribute("aunit","deg"));
426  polyconeElement->setAttributeNode(NewAttribute("lunit","mm"));
427  solElement->appendChild(polyconeElement);
428 
429  const size_t num_zplanes = polycone->GetOriginalParameters()->Num_z_planes;
430  const G4double* z_array = polycone->GetOriginalParameters()->Z_values;
431  const G4double* rmin_array = polycone->GetOriginalParameters()->Rmin;
432  const G4double* rmax_array = polycone->GetOriginalParameters()->Rmax;
433 
434  for (size_t i=0; i<num_zplanes; i++)
435  {
436  ZplaneWrite(polyconeElement,z_array[i],rmin_array[i],rmax_array[i]);
437  }
438 
439 
440 }
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4PolyconeHistorical * GetOriginalParameters() const
void ZplaneWrite(xercesc::DOMElement *, const G4double &, const G4double &, const G4double &)
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
double G4double
Definition: G4Types.hh:76
void G4GDMLWriteSolids::PolyhedraWrite ( xercesc::DOMElement *  solElement,
const G4Polyhedra * const  polyhedra 
)
protected

Definition at line 469 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4Polyhedra::GetCorner(), G4VSolid::GetName(), G4Polyhedra::GetNumRZCorner(), G4Polyhedra::GetOriginalParameters(), G4Polyhedra::IsGeneric(), G4GDMLWrite::NewAttribute(), G4GDMLWrite::NewElement(), G4PolyhedraHistorical::Num_z_planes, G4PolyhedraHistorical::numSide, G4PolyhedraHistorical::Opening_angle, G4PolyhedraSideRZ::r, G4PolyhedraHistorical::Rmax, G4PolyhedraHistorical::Rmin, RZPointWrite(), G4PolyhedraHistorical::Start_angle, G4PolyhedraSideRZ::z, G4PolyhedraHistorical::Z_values, and ZplaneWrite().

Referenced by AddSolid().

471 {
472  const G4String& name = GenerateName(polyhedra->GetName(),polyhedra);
473  if(polyhedra->IsGeneric() == false){
474  xercesc::DOMElement* polyhedraElement = NewElement("polyhedra");
475  polyhedraElement->setAttributeNode(NewAttribute("name",name));
476  polyhedraElement->setAttributeNode(NewAttribute("startphi",
477  polyhedra->GetOriginalParameters()->Start_angle/degree));
478  polyhedraElement->setAttributeNode(NewAttribute("deltaphi",
480  polyhedraElement->setAttributeNode(NewAttribute("numsides",
481  polyhedra->GetOriginalParameters()->numSide));
482  polyhedraElement->setAttributeNode(NewAttribute("aunit","deg"));
483  polyhedraElement->setAttributeNode(NewAttribute("lunit","mm"));
484  solElement->appendChild(polyhedraElement);
485 
486  const size_t num_zplanes = polyhedra->GetOriginalParameters()->Num_z_planes;
487  const G4double* z_array = polyhedra->GetOriginalParameters()->Z_values;
488  const G4double* rmin_array = polyhedra->GetOriginalParameters()->Rmin;
489  const G4double* rmax_array = polyhedra->GetOriginalParameters()->Rmax;
490 
491  const G4double convertRad =
492  std::cos(0.5*polyhedra->GetOriginalParameters()->Opening_angle
493  / polyhedra->GetOriginalParameters()->numSide);
494 
495  for (size_t i=0;i<num_zplanes;i++)
496  {
497  ZplaneWrite(polyhedraElement,z_array[i],
498  rmin_array[i]*convertRad, rmax_array[i]*convertRad);
499  }
500  }else{//generic polyhedra
501 
502  xercesc::DOMElement* polyhedraElement = NewElement("genericPolyhedra");
503  polyhedraElement->setAttributeNode(NewAttribute("name",name));
504  polyhedraElement->setAttributeNode(NewAttribute("startphi",
505  polyhedra->GetOriginalParameters()->Start_angle/degree));
506  polyhedraElement->setAttributeNode(NewAttribute("deltaphi",
508  polyhedraElement->setAttributeNode(NewAttribute("numsides",
509  polyhedra->GetOriginalParameters()->numSide));
510  polyhedraElement->setAttributeNode(NewAttribute("aunit","deg"));
511  polyhedraElement->setAttributeNode(NewAttribute("lunit","mm"));
512  solElement->appendChild(polyhedraElement);
513 
514  const size_t num_rzpoints = polyhedra->GetNumRZCorner();
515 
516  for (size_t i=0;i<num_rzpoints;i++)
517  {
518  const G4double r_point = polyhedra->GetCorner(i).r;
519  const G4double z_point = polyhedra->GetCorner(i).z;
520  RZPointWrite(polyhedraElement,r_point,z_point);
521  }
522  }
523 }
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4int GetNumRZCorner() const
void RZPointWrite(xercesc::DOMElement *, const G4double &, const G4double &)
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
void ZplaneWrite(xercesc::DOMElement *, const G4double &, const G4double &, const G4double &)
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4bool IsGeneric() const
G4PolyhedraHistorical * GetOriginalParameters() const
G4PolyhedraSideRZ GetCorner(const G4int index) const
double G4double
Definition: G4Types.hh:76
void G4GDMLWriteSolids::RZPointWrite ( xercesc::DOMElement *  element,
const G4double r,
const G4double z 
)
protected

Definition at line 927 of file G4GDMLWriteSolids.cc.

References python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by GenericPolyconeWrite(), and PolyhedraWrite().

929 {
930  xercesc::DOMElement* rzpointElement = NewElement("rzpoint");
931  rzpointElement->setAttributeNode(NewAttribute("r",r/mm));
932  rzpointElement->setAttributeNode(NewAttribute("z",z/mm));
933  element->appendChild(rzpointElement);
934 }
G4double z
Definition: TRTMaterials.hh:39
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::SolidsWrite ( xercesc::DOMElement *  gdmlElement)
virtual

Implements G4GDMLWrite.

Definition at line 953 of file G4GDMLWriteSolids.cc.

References G4cout, G4endl, G4GDMLWrite::NewElement(), solidList, and solidsElement.

954 {
955  G4cout << "G4GDML: Writing solids..." << G4endl;
956 
957  solidsElement = NewElement("solids");
958  gdmlElement->appendChild(solidsElement);
959 
960  solidList.clear();
961 }
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
xercesc::DOMElement * solidsElement
std::vector< const G4VSolid * > solidList
void G4GDMLWriteSolids::SphereWrite ( xercesc::DOMElement *  solElement,
const G4Sphere * const  sphere 
)
protected

Definition at line 526 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4Sphere::GetDeltaPhiAngle(), G4Sphere::GetDeltaThetaAngle(), G4Sphere::GetInnerRadius(), G4VSolid::GetName(), G4Sphere::GetOuterRadius(), G4Sphere::GetStartPhiAngle(), G4Sphere::GetStartThetaAngle(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

527 {
528  const G4String& name = GenerateName(sphere->GetName(),sphere);
529 
530  xercesc::DOMElement* sphereElement = NewElement("sphere");
531  sphereElement->setAttributeNode(NewAttribute("name",name));
532  sphereElement->setAttributeNode(NewAttribute("rmin",
533  sphere->GetInnerRadius()/mm));
534  sphereElement->setAttributeNode(NewAttribute("rmax",
535  sphere->GetOuterRadius()/mm));
536  sphereElement->setAttributeNode(NewAttribute("startphi",
537  sphere->GetStartPhiAngle()/degree));
538  sphereElement->setAttributeNode(NewAttribute("deltaphi",
539  sphere->GetDeltaPhiAngle()/degree));
540  sphereElement->setAttributeNode(NewAttribute("starttheta",
541  sphere->GetStartThetaAngle()/degree));
542  sphereElement->setAttributeNode(NewAttribute("deltatheta",
543  sphere->GetDeltaThetaAngle()/degree));
544  sphereElement->setAttributeNode(NewAttribute("aunit","deg"));
545  sphereElement->setAttributeNode(NewAttribute("lunit","mm"));
546  solElement->appendChild(sphereElement);
547 }
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetDeltaPhiAngle() const
G4double GetStartThetaAngle() const
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetInnerRadius() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetOuterRadius() const
G4double GetStartPhiAngle() const
G4double GetDeltaThetaAngle() const
void G4GDMLWriteSolids::TessellatedWrite ( xercesc::DOMElement *  solElement,
const G4TessellatedSolid * const  tessellated 
)
protected

Definition at line 550 of file G4GDMLWriteSolids.cc.

References G4GDMLWriteDefine::AddPosition(), FatalException, test::fname, G4Exception(), G4GDMLWrite::GenerateName(), G4TessellatedSolid::GetFacet(), G4VSolid::GetName(), G4TessellatedSolid::GetNumberOfFacets(), G4VFacet::GetNumberOfVertices(), G4VFacet::GetVertex(), G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

552 {
553  const G4String& solid_name = tessellated->GetName();
554  const G4String& name = GenerateName(solid_name, tessellated);
555 
556  xercesc::DOMElement* tessellatedElement = NewElement("tessellated");
557  tessellatedElement->setAttributeNode(NewAttribute("name",name));
558  tessellatedElement->setAttributeNode(NewAttribute("aunit","deg"));
559  tessellatedElement->setAttributeNode(NewAttribute("lunit","mm"));
560  solElement->appendChild(tessellatedElement);
561 
562  std::map<G4ThreeVector, G4String> vertexMap;
563 
564  const size_t NumFacets = tessellated->GetNumberOfFacets();
565  size_t NumVertex = 0;
566 
567  for (size_t i=0;i<NumFacets;i++)
568  {
569  const G4VFacet* facet = tessellated->GetFacet(i);
570  const size_t NumVertexPerFacet = facet->GetNumberOfVertices();
571 
572  G4String FacetTag;
573 
574  if (NumVertexPerFacet==3) { FacetTag="triangular"; } else
575  if (NumVertexPerFacet==4) { FacetTag="quadrangular"; }
576  else
577  {
578  G4Exception("G4GDMLWriteSolids::TessellatedWrite()", "InvalidSetup",
579  FatalException, "Facet should contain 3 or 4 vertices!");
580  }
581 
582  xercesc::DOMElement* facetElement = NewElement(FacetTag);
583  tessellatedElement->appendChild(facetElement);
584 
585  for (size_t j=0; j<NumVertexPerFacet; j++)
586  {
587  std::stringstream name_stream;
588  std::stringstream ref_stream;
589 
590  name_stream << "vertex" << (j+1);
591  ref_stream << solid_name << "_v" << NumVertex;
592 
593  const G4String& fname = name_stream.str(); // facet's tag variable
594  G4String ref = ref_stream.str(); // vertex tag to be associated
595 
596  // Now search for the existance of the current vertex in the
597  // map of cached vertices. If existing, do NOT store it as
598  // position in the GDML file, so avoiding duplication; otherwise
599  // cache it in the local map and add it as position in the
600  // "define" section of the GDML file.
601 
602  const G4ThreeVector& vertex = facet->GetVertex(j);
603 
604  if(vertexMap.find(vertex) != vertexMap.end()) // Vertex is cached
605  {
606  ref = vertexMap[vertex]; // Set the proper tag for it
607  }
608  else // Vertex not found
609  {
610  vertexMap.insert(std::make_pair(vertex,ref)); // Cache vertex and ...
611  AddPosition(ref, vertex); // ... add it to define section!
612  NumVertex++;
613  }
614 
615  // Now create association of the vertex with its facet
616  //
617  facetElement->setAttributeNode(NewAttribute(fname,ref));
618  }
619  }
620 }
G4String GetName() const
virtual G4int GetNumberOfVertices() const =0
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4VFacet * GetFacet(G4int i) const
void AddPosition(const G4String &name, const G4ThreeVector &pos)
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4int GetNumberOfFacets() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
virtual G4ThreeVector GetVertex(G4int i) const =0
void G4GDMLWriteSolids::TetWrite ( xercesc::DOMElement *  solElement,
const G4Tet * const  tet 
)
protected

Definition at line 623 of file G4GDMLWriteSolids.cc.

References G4GDMLWriteDefine::AddPosition(), G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Tet::GetVertices(), G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

624 {
625  const G4String& solid_name = tet->GetName();
626  const G4String& name = GenerateName(solid_name, tet);
627 
628  std::vector<G4ThreeVector> vertexList = tet->GetVertices();
629 
630  xercesc::DOMElement* tetElement = NewElement("tet");
631  tetElement->setAttributeNode(NewAttribute("name",name));
632  tetElement->setAttributeNode(NewAttribute("vertex1",solid_name+"_v1"));
633  tetElement->setAttributeNode(NewAttribute("vertex2",solid_name+"_v2"));
634  tetElement->setAttributeNode(NewAttribute("vertex3",solid_name+"_v3"));
635  tetElement->setAttributeNode(NewAttribute("vertex4",solid_name+"_v4"));
636  tetElement->setAttributeNode(NewAttribute("lunit","mm"));
637  solElement->appendChild(tetElement);
638 
639  AddPosition(solid_name+"_v1",vertexList[0]);
640  AddPosition(solid_name+"_v2",vertexList[1]);
641  AddPosition(solid_name+"_v3",vertexList[2]);
642  AddPosition(solid_name+"_v4",vertexList[3]);
643 }
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
std::vector< G4ThreeVector > GetVertices() const
Definition: G4Tet.cc:780
void AddPosition(const G4String &name, const G4ThreeVector &pos)
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::TorusWrite ( xercesc::DOMElement *  solElement,
const G4Torus * const  torus 
)
protected

Definition at line 646 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4Torus::GetDPhi(), G4VSolid::GetName(), G4Torus::GetRmax(), G4Torus::GetRmin(), G4Torus::GetRtor(), G4Torus::GetSPhi(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

647 {
648  const G4String& name = GenerateName(torus->GetName(),torus);
649 
650  xercesc::DOMElement* torusElement = NewElement("torus");
651  torusElement->setAttributeNode(NewAttribute("name",name));
652  torusElement->setAttributeNode(NewAttribute("rmin",torus->GetRmin()/mm));
653  torusElement->setAttributeNode(NewAttribute("rmax",torus->GetRmax()/mm));
654  torusElement->setAttributeNode(NewAttribute("rtor",torus->GetRtor()/mm));
655  torusElement->
656  setAttributeNode(NewAttribute("startphi",torus->GetSPhi()/degree));
657  torusElement->
658  setAttributeNode(NewAttribute("deltaphi",torus->GetDPhi()/degree));
659  torusElement->setAttributeNode(NewAttribute("aunit","deg"));
660  torusElement->setAttributeNode(NewAttribute("lunit","mm"));
661  solElement->appendChild(torusElement);
662 }
G4double GetSPhi() const
G4String GetName() const
G4double GetRmax() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetRtor() const
G4double GetRmin() const
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetDPhi() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::TrapWrite ( xercesc::DOMElement *  solElement,
const G4Trap * const  trap 
)
protected

Definition at line 697 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Trap::GetSymAxis(), G4Trap::GetTanAlpha1(), G4Trap::GetTanAlpha2(), G4Trap::GetXHalfLength1(), G4Trap::GetXHalfLength2(), G4Trap::GetXHalfLength3(), G4Trap::GetXHalfLength4(), G4Trap::GetYHalfLength1(), G4Trap::GetYHalfLength2(), G4Trap::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

698 {
699  const G4String& name = GenerateName(trap->GetName(),trap);
700 
701  const G4ThreeVector& simaxis = trap->GetSymAxis();
702  const G4double phi = (simaxis.z() != 1.0)
703  ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
704  const G4double theta = std::acos(simaxis.z());
705  const G4double alpha1 = std::atan(trap->GetTanAlpha1());
706  const G4double alpha2 = std::atan(trap->GetTanAlpha2());
707 
708  xercesc::DOMElement* trapElement = NewElement("trap");
709  trapElement->setAttributeNode(NewAttribute("name",name));
710  trapElement->setAttributeNode(NewAttribute("z",
711  2.0*trap->GetZHalfLength()/mm));
712  trapElement->setAttributeNode(NewAttribute("theta",theta/degree));
713  trapElement->setAttributeNode(NewAttribute("phi",phi/degree));
714  trapElement->setAttributeNode(NewAttribute("y1",
715  2.0*trap->GetYHalfLength1()/mm));
716  trapElement->setAttributeNode(NewAttribute("x1",
717  2.0*trap->GetXHalfLength1()/mm));
718  trapElement->setAttributeNode(NewAttribute("x2",
719  2.0*trap->GetXHalfLength2()/mm));
720  trapElement->setAttributeNode(NewAttribute("alpha1",alpha1/degree));
721  trapElement->setAttributeNode(NewAttribute("y2",
722  2.0*trap->GetYHalfLength2()/mm));
723  trapElement->setAttributeNode(NewAttribute("x3",
724  2.0*trap->GetXHalfLength3()/mm));
725  trapElement->setAttributeNode(NewAttribute("x4",
726  2.0*trap->GetXHalfLength4()/mm));
727  trapElement->setAttributeNode(NewAttribute("alpha2",alpha2/degree));
728  trapElement->setAttributeNode(NewAttribute("aunit","deg"));
729  trapElement->setAttributeNode(NewAttribute("lunit","mm"));
730  solElement->appendChild(trapElement);
731 }
G4String GetName() const
G4double GetXHalfLength4() const
G4double GetYHalfLength2() const
G4double GetZHalfLength() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetXHalfLength2() const
G4double GetTanAlpha2() const
G4double GetXHalfLength1() const
G4double GetXHalfLength3() const
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4ThreeVector GetSymAxis() const
G4double GetYHalfLength1() const
double G4double
Definition: G4Types.hh:76
G4double GetTanAlpha1() const
void G4GDMLWriteSolids::TrdWrite ( xercesc::DOMElement *  solElement,
const G4Trd * const  trd 
)
protected

Definition at line 734 of file G4GDMLWriteSolids.cc.

References G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4Trd::GetXHalfLength1(), G4Trd::GetXHalfLength2(), G4Trd::GetYHalfLength1(), G4Trd::GetYHalfLength2(), G4Trd::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

735 {
736  const G4String& name = GenerateName(trd->GetName(),trd);
737 
738  xercesc::DOMElement* trdElement = NewElement("trd");
739  trdElement->setAttributeNode(NewAttribute("name",name));
740  trdElement->setAttributeNode(NewAttribute("x1",
741  2.0*trd->GetXHalfLength1()/mm));
742  trdElement->setAttributeNode(NewAttribute("x2",
743  2.0*trd->GetXHalfLength2()/mm));
744  trdElement->setAttributeNode(NewAttribute("y1",
745  2.0*trd->GetYHalfLength1()/mm));
746  trdElement->setAttributeNode(NewAttribute("y2",
747  2.0*trd->GetYHalfLength2()/mm));
748  trdElement->setAttributeNode(NewAttribute("z",
749  2.0*trd->GetZHalfLength()/mm));
750  trdElement->setAttributeNode(NewAttribute("lunit","mm"));
751  solElement->appendChild(trdElement);
752 }
G4String GetName() const
G4double GetYHalfLength1() const
const XML_Char * name
G4double GetZHalfLength() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetXHalfLength2() const
G4double GetYHalfLength2() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetXHalfLength1() const
void G4GDMLWriteSolids::TubeWrite ( xercesc::DOMElement *  solElement,
const G4Tubs * const  tube 
)
protected

Definition at line 755 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4Tubs::GetDeltaPhiAngle(), G4Tubs::GetInnerRadius(), G4VSolid::GetName(), G4Tubs::GetOuterRadius(), G4Tubs::GetStartPhiAngle(), G4Tubs::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

756 {
757  const G4String& name = GenerateName(tube->GetName(),tube);
758 
759  xercesc::DOMElement* tubeElement = NewElement("tube");
760  tubeElement->setAttributeNode(NewAttribute("name",name));
761  tubeElement->setAttributeNode(NewAttribute("rmin",
762  tube->GetInnerRadius()/mm));
763  tubeElement->setAttributeNode(NewAttribute("rmax",
764  tube->GetOuterRadius()/mm));
765  tubeElement->setAttributeNode(NewAttribute("z",
766  2.0*tube->GetZHalfLength()/mm));
767  tubeElement->setAttributeNode(NewAttribute("startphi",
768  tube->GetStartPhiAngle()/degree));
769  tubeElement->setAttributeNode(NewAttribute("deltaphi",
770  tube->GetDeltaPhiAngle()/degree));
771  tubeElement->setAttributeNode(NewAttribute("aunit","deg"));
772  tubeElement->setAttributeNode(NewAttribute("lunit","mm"));
773  solElement->appendChild(tubeElement);
774 }
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetDeltaPhiAngle() const
tuple degree
Definition: hepunit.py:69
G4double GetStartPhiAngle() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetInnerRadius() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetZHalfLength() const
G4double GetOuterRadius() const
void G4GDMLWriteSolids::TwistedboxWrite ( xercesc::DOMElement *  solElement,
const G4TwistedBox * const  twistedbox 
)
protected

Definition at line 811 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4TwistedBox::GetPhiTwist(), G4TwistedBox::GetXHalfLength(), G4TwistedBox::GetYHalfLength(), G4TwistedBox::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

813 {
814  const G4String& name = GenerateName(twistedbox->GetName(),twistedbox);
815 
816  xercesc::DOMElement* twistedboxElement = NewElement("twistedbox");
817  twistedboxElement->setAttributeNode(NewAttribute("name",name));
818  twistedboxElement->setAttributeNode(NewAttribute("x",
819  2.0*twistedbox->GetXHalfLength()/mm));
820  twistedboxElement->setAttributeNode(NewAttribute("y",
821  2.0*twistedbox->GetYHalfLength()/mm));
822  twistedboxElement->setAttributeNode(NewAttribute("z",
823  2.0*twistedbox->GetZHalfLength()/mm));
824  twistedboxElement->setAttributeNode(NewAttribute("PhiTwist",
825  twistedbox->GetPhiTwist()/degree));
826  twistedboxElement->setAttributeNode(NewAttribute("aunit","deg"));
827  twistedboxElement->setAttributeNode(NewAttribute("lunit","mm"));
828  solElement->appendChild(twistedboxElement);
829 }
G4String GetName() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetZHalfLength() const
Definition: G4TwistedBox.hh:75
tuple degree
Definition: hepunit.py:69
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetYHalfLength() const
Definition: G4TwistedBox.hh:74
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetPhiTwist() const
Definition: G4TwistedBox.hh:76
G4double GetXHalfLength() const
Definition: G4TwistedBox.hh:73
void G4GDMLWriteSolids::TwistedtrapWrite ( xercesc::DOMElement *  solElement,
const G4TwistedTrap * const  twistedtrap 
)
protected

Definition at line 832 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4TwistedTrap::GetAzimuthalAnglePhi(), G4VSolid::GetName(), G4TwistedTrap::GetPhiTwist(), G4TwistedTrap::GetPolarAngleTheta(), G4TwistedTrap::GetTiltAngleAlpha(), G4TwistedTrap::GetX1HalfLength(), G4TwistedTrap::GetX2HalfLength(), G4TwistedTrap::GetX3HalfLength(), G4TwistedTrap::GetX4HalfLength(), G4TwistedTrap::GetY1HalfLength(), G4TwistedTrap::GetY2HalfLength(), G4TwistedTrap::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

834 {
835  const G4String& name = GenerateName(twistedtrap->GetName(),twistedtrap);
836 
837  xercesc::DOMElement* twistedtrapElement = NewElement("twistedtrap");
838  twistedtrapElement->setAttributeNode(NewAttribute("name",name));
839  twistedtrapElement->setAttributeNode(NewAttribute("y1",
840  2.0*twistedtrap->GetY1HalfLength()/mm));
841  twistedtrapElement->setAttributeNode(NewAttribute("x1",
842  2.0*twistedtrap->GetX1HalfLength()/mm));
843  twistedtrapElement->setAttributeNode(NewAttribute("x2",
844  2.0*twistedtrap->GetX2HalfLength()/mm));
845  twistedtrapElement->setAttributeNode(NewAttribute("y2",
846  2.0*twistedtrap->GetY2HalfLength()/mm));
847  twistedtrapElement->setAttributeNode(NewAttribute("x3",
848  2.0*twistedtrap->GetX3HalfLength()/mm));
849  twistedtrapElement->setAttributeNode(NewAttribute("x4",
850  2.0*twistedtrap->GetX4HalfLength()/mm));
851  twistedtrapElement->setAttributeNode(NewAttribute("z",
852  2.0*twistedtrap->GetZHalfLength()/mm));
853  twistedtrapElement->setAttributeNode(NewAttribute("Alph",
854  twistedtrap->GetTiltAngleAlpha()/degree));
855  twistedtrapElement->setAttributeNode(NewAttribute("Theta",
856  twistedtrap->GetPolarAngleTheta()/degree));
857  twistedtrapElement->setAttributeNode(NewAttribute("Phi",
858  twistedtrap->GetAzimuthalAnglePhi()/degree));
859  twistedtrapElement->setAttributeNode(NewAttribute("PhiTwist",
860  twistedtrap->GetPhiTwist()/degree));
861  twistedtrapElement->setAttributeNode(NewAttribute("aunit","deg"));
862  twistedtrapElement->setAttributeNode(NewAttribute("lunit","mm"));
863 
864  solElement->appendChild(twistedtrapElement);
865 }
G4String GetName() const
G4double GetPolarAngleTheta() const
G4double GetZHalfLength() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetY1HalfLength() const
tuple degree
Definition: hepunit.py:69
G4double GetX3HalfLength() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetX4HalfLength() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetAzimuthalAnglePhi() const
G4double GetX1HalfLength() const
G4double GetTiltAngleAlpha() const
G4double GetX2HalfLength() const
G4double GetPhiTwist() const
G4double GetY2HalfLength() const
void G4GDMLWriteSolids::TwistedtrdWrite ( xercesc::DOMElement *  solElement,
const G4TwistedTrd * const  twistedtrd 
)
protected

Definition at line 868 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4TwistedTrd::GetPhiTwist(), G4TwistedTrd::GetX1HalfLength(), G4TwistedTrd::GetX2HalfLength(), G4TwistedTrd::GetY1HalfLength(), G4TwistedTrd::GetY2HalfLength(), G4TwistedTrd::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

870 {
871  const G4String& name = GenerateName(twistedtrd->GetName(),twistedtrd);
872 
873  xercesc::DOMElement* twistedtrdElement = NewElement("twistedtrd");
874  twistedtrdElement->setAttributeNode(NewAttribute("name",name));
875  twistedtrdElement->setAttributeNode(NewAttribute("x1",
876  2.0*twistedtrd->GetX1HalfLength()/mm));
877  twistedtrdElement->setAttributeNode(NewAttribute("x2",
878  2.0*twistedtrd->GetX2HalfLength()/mm));
879  twistedtrdElement->setAttributeNode(NewAttribute("y1",
880  2.0*twistedtrd->GetY1HalfLength()/mm));
881  twistedtrdElement->setAttributeNode(NewAttribute("y2",
882  2.0*twistedtrd->GetY2HalfLength()/mm));
883  twistedtrdElement->setAttributeNode(NewAttribute("z",
884  2.0*twistedtrd->GetZHalfLength()/mm));
885  twistedtrdElement->setAttributeNode(NewAttribute("PhiTwist",
886  twistedtrd->GetPhiTwist()/degree));
887  twistedtrdElement->setAttributeNode(NewAttribute("aunit","deg"));
888  twistedtrdElement->setAttributeNode(NewAttribute("lunit","mm"));
889  solElement->appendChild(twistedtrdElement);
890 }
G4String GetName() const
G4double GetX1HalfLength() const
Definition: G4TwistedTrd.hh:77
G4double GetY1HalfLength() const
Definition: G4TwistedTrd.hh:79
const XML_Char * name
G4double GetX2HalfLength() const
Definition: G4TwistedTrd.hh:78
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetZHalfLength() const
Definition: G4TwistedTrd.hh:81
tuple degree
Definition: hepunit.py:69
G4double GetY2HalfLength() const
Definition: G4TwistedTrd.hh:80
G4double GetPhiTwist() const
Definition: G4TwistedTrd.hh:82
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
void G4GDMLWriteSolids::TwistedtubsWrite ( xercesc::DOMElement *  solElement,
const G4TwistedTubs * const  twistedtubs 
)
protected

Definition at line 893 of file G4GDMLWriteSolids.cc.

References python.hepunit::degree, G4GDMLWrite::GenerateName(), G4TwistedTubs::GetDPhi(), G4TwistedTubs::GetInnerRadius(), G4VSolid::GetName(), G4TwistedTubs::GetOuterRadius(), G4TwistedTubs::GetPhiTwist(), G4TwistedTubs::GetZHalfLength(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by AddSolid().

895 {
896  const G4String& name = GenerateName(twistedtubs->GetName(),twistedtubs);
897 
898  xercesc::DOMElement* twistedtubsElement = NewElement("twistedtubs");
899  twistedtubsElement->setAttributeNode(NewAttribute("name",name));
900  twistedtubsElement->setAttributeNode(NewAttribute("twistedangle",
901  twistedtubs->GetPhiTwist()/degree));
902  twistedtubsElement->setAttributeNode(NewAttribute("endinnerrad",
903  twistedtubs->GetInnerRadius()/mm));
904  twistedtubsElement->setAttributeNode(NewAttribute("endouterrad",
905  twistedtubs->GetOuterRadius()/mm));
906  twistedtubsElement->setAttributeNode(NewAttribute("zlen",
907  2.0*twistedtubs->GetZHalfLength()/mm));
908  twistedtubsElement->setAttributeNode(NewAttribute("phi",
909  twistedtubs->GetDPhi()/degree));
910  twistedtubsElement->setAttributeNode(NewAttribute("aunit","deg"));
911  twistedtubsElement->setAttributeNode(NewAttribute("lunit","mm"));
912  solElement->appendChild(twistedtubsElement);
913 }
G4String GetName() const
const XML_Char * name
G4double GetInnerRadius() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
G4double GetPhiTwist() const
tuple degree
Definition: hepunit.py:69
G4double GetDPhi() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4double GetZHalfLength() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
G4double GetOuterRadius() const
void G4GDMLWriteSolids::XtruWrite ( xercesc::DOMElement *  solElement,
const G4ExtrudedSolid * const  xtru 
)
protected

Definition at line 292 of file G4GDMLWriteSolids.cc.

References G4ExtrudedSolid::ZSection::fOffset, G4ExtrudedSolid::ZSection::fScale, G4ExtrudedSolid::ZSection::fZ, G4GDMLWrite::GenerateName(), G4VSolid::GetName(), G4ExtrudedSolid::GetNofVertices(), G4ExtrudedSolid::GetNofZSections(), G4ExtrudedSolid::GetVertex(), G4ExtrudedSolid::GetZSection(), python.hepunit::mm, G4GDMLWrite::NewAttribute(), G4GDMLWrite::NewElement(), CLHEP::Hep2Vector::x(), and CLHEP::Hep2Vector::y().

Referenced by AddSolid().

294 {
295  const G4String& name = GenerateName(xtru->GetName(),xtru);
296 
297  xercesc::DOMElement* xtruElement = NewElement("xtru");
298  xtruElement->setAttributeNode(NewAttribute("name",name));
299  xtruElement->setAttributeNode(NewAttribute("lunit","mm"));
300  solElement->appendChild(xtruElement);
301 
302  const G4int NumVertex = xtru->GetNofVertices();
303 
304  for (G4int i=0;i<NumVertex;i++)
305  {
306  xercesc::DOMElement* twoDimVertexElement = NewElement("twoDimVertex");
307  xtruElement->appendChild(twoDimVertexElement);
308 
309  const G4TwoVector& vertex = xtru->GetVertex(i);
310 
311  twoDimVertexElement->setAttributeNode(NewAttribute("x",vertex.x()/mm));
312  twoDimVertexElement->setAttributeNode(NewAttribute("y",vertex.y()/mm));
313  }
314 
315  const G4int NumSection = xtru->GetNofZSections();
316 
317  for (G4int i=0;i<NumSection;i++)
318  {
319  xercesc::DOMElement* sectionElement = NewElement("section");
320  xtruElement->appendChild(sectionElement);
321 
322  const G4ExtrudedSolid::ZSection section = xtru->GetZSection(i);
323 
324  sectionElement->setAttributeNode(NewAttribute("zOrder",i));
325  sectionElement->setAttributeNode(NewAttribute("zPosition",section.fZ/mm));
326  sectionElement->
327  setAttributeNode(NewAttribute("xOffset",section.fOffset.x()/mm));
328  sectionElement->
329  setAttributeNode(NewAttribute("yOffset",section.fOffset.y()/mm));
330  sectionElement->
331  setAttributeNode(NewAttribute("scalingFactor",section.fScale));
332  }
333 }
G4String GetName() const
double y() const
double x() const
const XML_Char * name
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
int G4int
Definition: G4Types.hh:78
G4TwoVector GetVertex(G4int index) const
G4int GetNofVertices() const
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:90
G4int GetNofZSections() const
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
ZSection GetZSection(G4int index) const
void G4GDMLWriteSolids::ZplaneWrite ( xercesc::DOMElement *  element,
const G4double z,
const G4double rmin,
const G4double rmax 
)
protected

Definition at line 916 of file G4GDMLWriteSolids.cc.

References python.hepunit::mm, G4GDMLWrite::NewAttribute(), and G4GDMLWrite::NewElement().

Referenced by G4GDMLWriteParamvol::Polycone_dimensionsWrite(), PolyconeWrite(), G4GDMLWriteParamvol::Polyhedra_dimensionsWrite(), and PolyhedraWrite().

918 {
919  xercesc::DOMElement* zplaneElement = NewElement("zplane");
920  zplaneElement->setAttributeNode(NewAttribute("z",z/mm));
921  zplaneElement->setAttributeNode(NewAttribute("rmin",rmin/mm));
922  zplaneElement->setAttributeNode(NewAttribute("rmax",rmax/mm));
923  element->appendChild(zplaneElement);
924 }
G4double z
Definition: TRTMaterials.hh:39
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103

Field Documentation

const G4int G4GDMLWriteSolids::maxTransforms = 8
staticprotected

Definition at line 126 of file G4GDMLWriteSolids.hh.

Referenced by BooleanWrite(), and G4GDMLWriteStructure::TraverseVolumeTree().

std::vector<const G4VSolid*> G4GDMLWriteSolids::solidList
protected

Definition at line 124 of file G4GDMLWriteSolids.hh.

Referenced by AddSolid(), and SolidsWrite().

xercesc::DOMElement* G4GDMLWriteSolids::solidsElement
protected

The documentation for this class was generated from the following files: