Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4GDMLWriteSolids Class Reference

#include <G4GDMLWriteSolids.hh>

+ Inheritance diagram for G4GDMLWriteSolids:

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 DefineWrite (xercesc::DOMElement *)=0
 
virtual void MaterialsWrite (xercesc::DOMElement *)=0
 
virtual void SolidsWrite (xercesc::DOMElement *)=0
 
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 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 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 76 of file G4GDMLWriteSolids.hh.

Constructor & Destructor Documentation

◆ G4GDMLWriteSolids()

G4GDMLWriteSolids::G4GDMLWriteSolids ( )
protected

Definition at line 71 of file G4GDMLWriteSolids.cc.

73{
74}
xercesc::DOMElement * solidsElement

◆ ~G4GDMLWriteSolids()

G4GDMLWriteSolids::~G4GDMLWriteSolids ( )
protectedvirtual

Definition at line 76 of file G4GDMLWriteSolids.cc.

77{
78}

Member Function Documentation

◆ AddSolid()

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

Definition at line 902 of file G4GDMLWriteSolids.cc.

903{
904 for (size_t i=0; i<solidList.size(); i++) // Check if solid is
905 { // already in the list!
906 if (solidList[i] == solidPtr) { return; }
907 }
908
909 solidList.push_back(solidPtr);
910
911 if (const G4BooleanSolid* const booleanPtr
912 = dynamic_cast<const G4BooleanSolid*>(solidPtr))
913 { BooleanWrite(solidsElement,booleanPtr); } else
914 if (const G4Box* const boxPtr
915 = dynamic_cast<const G4Box*>(solidPtr))
916 { BoxWrite(solidsElement,boxPtr); } else
917 if (const G4Cons* const conePtr
918 = dynamic_cast<const G4Cons*>(solidPtr))
919 { ConeWrite(solidsElement,conePtr); } else
920 if (const G4EllipticalCone* const elconePtr
921 = dynamic_cast<const G4EllipticalCone*>(solidPtr))
922 { ElconeWrite(solidsElement,elconePtr); } else
923 if (const G4Ellipsoid* const ellipsoidPtr
924 = dynamic_cast<const G4Ellipsoid*>(solidPtr))
925 { EllipsoidWrite(solidsElement,ellipsoidPtr); } else
926 if (const G4EllipticalTube* const eltubePtr
927 = dynamic_cast<const G4EllipticalTube*>(solidPtr))
928 { EltubeWrite(solidsElement,eltubePtr); } else
929 if (const G4ExtrudedSolid* const xtruPtr
930 = dynamic_cast<const G4ExtrudedSolid*>(solidPtr))
931 { XtruWrite(solidsElement,xtruPtr); } else
932 if (const G4Hype* const hypePtr
933 = dynamic_cast<const G4Hype*>(solidPtr))
934 { HypeWrite(solidsElement,hypePtr); } else
935 if (const G4Orb* const orbPtr
936 = dynamic_cast<const G4Orb*>(solidPtr))
937 { OrbWrite(solidsElement,orbPtr); } else
938 if (const G4Para* const paraPtr
939 = dynamic_cast<const G4Para*>(solidPtr))
940 { ParaWrite(solidsElement,paraPtr); } else
941 if (const G4Paraboloid* const paraboloidPtr
942 = dynamic_cast<const G4Paraboloid*>(solidPtr))
943 { ParaboloidWrite(solidsElement,paraboloidPtr); } else
944 if (const G4Polycone* const polyconePtr
945 = dynamic_cast<const G4Polycone*>(solidPtr))
946 { PolyconeWrite(solidsElement,polyconePtr); } else
947 if (const G4Polyhedra* const polyhedraPtr
948 = dynamic_cast<const G4Polyhedra*>(solidPtr))
949 { PolyhedraWrite(solidsElement,polyhedraPtr); } else
950 if (const G4Sphere* const spherePtr
951 = dynamic_cast<const G4Sphere*>(solidPtr))
952 { SphereWrite(solidsElement,spherePtr); } else
953 if (const G4TessellatedSolid* const tessellatedPtr
954 = dynamic_cast<const G4TessellatedSolid*>(solidPtr))
955 { TessellatedWrite(solidsElement,tessellatedPtr); } else
956 if (const G4Tet* const tetPtr
957 = dynamic_cast<const G4Tet*>(solidPtr))
958 { TetWrite(solidsElement,tetPtr); } else
959 if (const G4Torus* const torusPtr
960 = dynamic_cast<const G4Torus*>(solidPtr))
961 { TorusWrite(solidsElement,torusPtr); } else
962 if (const G4GenericTrap* const gtrapPtr
963 = dynamic_cast<const G4GenericTrap*>(solidPtr))
964 { GenTrapWrite(solidsElement,gtrapPtr); } else
965 if (const G4Trap* const trapPtr
966 = dynamic_cast<const G4Trap*>(solidPtr))
967 { TrapWrite(solidsElement,trapPtr); } else
968 if (const G4Trd* const trdPtr
969 = dynamic_cast<const G4Trd*>(solidPtr))
970 { TrdWrite(solidsElement,trdPtr); } else
971 if (const G4Tubs* const tubePtr
972 = dynamic_cast<const G4Tubs*>(solidPtr))
973 { TubeWrite(solidsElement,tubePtr); } else
974 if (const G4CutTubs* const cuttubePtr
975 = dynamic_cast<const G4CutTubs*>(solidPtr))
976 { CutTubeWrite(solidsElement,cuttubePtr); } else
977 if (const G4TwistedBox* const twistedboxPtr
978 = dynamic_cast<const G4TwistedBox*>(solidPtr))
979 { TwistedboxWrite(solidsElement,twistedboxPtr); } else
980 if (const G4TwistedTrap* const twistedtrapPtr
981 = dynamic_cast<const G4TwistedTrap*>(solidPtr))
982 { TwistedtrapWrite(solidsElement,twistedtrapPtr); } else
983 if (const G4TwistedTrd* const twistedtrdPtr
984 = dynamic_cast<const G4TwistedTrd*>(solidPtr))
985 { TwistedtrdWrite(solidsElement,twistedtrdPtr); } else
986 if (const G4TwistedTubs* const twistedtubsPtr
987 = dynamic_cast<const G4TwistedTubs*>(solidPtr))
988 { TwistedtubsWrite(solidsElement,twistedtubsPtr); }
989 else
990 {
991 G4String error_msg = "Unknown solid: " + solidPtr->GetName()
992 + "; Type: " + solidPtr->GetEntityType();
993 G4Exception("G4GDMLWriteSolids::AddSolid()", "WriteError",
994 FatalException, error_msg);
995 }
996}
@ FatalException
Definition: G4Box.hh:55
Definition: G4Cons.hh:75
void TwistedtrdWrite(xercesc::DOMElement *, const G4TwistedTrd *const)
void TorusWrite(xercesc::DOMElement *, const G4Torus *const)
void TetWrite(xercesc::DOMElement *, const G4Tet *const)
void TrapWrite(xercesc::DOMElement *, const G4Trap *const)
void HypeWrite(xercesc::DOMElement *, const G4Hype *const)
void TwistedtrapWrite(xercesc::DOMElement *, const G4TwistedTrap *const)
void PolyhedraWrite(xercesc::DOMElement *, const G4Polyhedra *const)
void TessellatedWrite(xercesc::DOMElement *, const G4TessellatedSolid *const)
void GenTrapWrite(xercesc::DOMElement *, const G4GenericTrap *const)
void TwistedboxWrite(xercesc::DOMElement *, const G4TwistedBox *const)
void TubeWrite(xercesc::DOMElement *, const G4Tubs *const)
void ParaboloidWrite(xercesc::DOMElement *, const G4Paraboloid *const)
std::vector< const G4VSolid * > solidList
void SphereWrite(xercesc::DOMElement *, const G4Sphere *const)
void BoxWrite(xercesc::DOMElement *, const G4Box *const)
void OrbWrite(xercesc::DOMElement *, const G4Orb *const)
void TwistedtubsWrite(xercesc::DOMElement *, const G4TwistedTubs *const)
void EltubeWrite(xercesc::DOMElement *, const G4EllipticalTube *const)
void BooleanWrite(xercesc::DOMElement *, const G4BooleanSolid *const)
void ElconeWrite(xercesc::DOMElement *, const G4EllipticalCone *const)
void ConeWrite(xercesc::DOMElement *, const G4Cons *const)
void CutTubeWrite(xercesc::DOMElement *, const G4CutTubs *const)
void EllipsoidWrite(xercesc::DOMElement *, const G4Ellipsoid *const)
void PolyconeWrite(xercesc::DOMElement *, const G4Polycone *const)
void ParaWrite(xercesc::DOMElement *, const G4Para *const)
void TrdWrite(xercesc::DOMElement *, const G4Trd *const)
void XtruWrite(xercesc::DOMElement *, const G4ExtrudedSolid *const)
Definition: G4Hype.hh:67
Definition: G4Orb.hh:52
Definition: G4Para.hh:77
Definition: G4Tet.hh:57
Definition: G4Trd.hh:63
Definition: G4Tubs.hh:77
G4String GetName() const
virtual G4GeometryType GetEntityType() const =0
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

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

◆ BooleanWrite()

void G4GDMLWriteSolids::BooleanWrite ( xercesc::DOMElement *  solElement,
const G4BooleanSolid * const  boolean 
)
protected

Definition at line 80 of file G4GDMLWriteSolids.cc.

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

Referenced by AddSolid().

◆ BoxWrite()

void G4GDMLWriteSolids::BoxWrite ( xercesc::DOMElement *  solElement,
const G4Box * const  box 
)
protected

Definition at line 196 of file G4GDMLWriteSolids.cc.

198{
199 const G4String& name = GenerateName(box->GetName(),box);
200
201 xercesc::DOMElement* boxElement = NewElement("box");
202 boxElement->setAttributeNode(NewAttribute("name",name));
203 boxElement->setAttributeNode(NewAttribute("x",2.0*box->GetXHalfLength()/mm));
204 boxElement->setAttributeNode(NewAttribute("y",2.0*box->GetYHalfLength()/mm));
205 boxElement->setAttributeNode(NewAttribute("z",2.0*box->GetZHalfLength()/mm));
206 boxElement->setAttributeNode(NewAttribute("lunit","mm"));
207 solElement->appendChild(boxElement);
208}
G4double GetYHalfLength() const
G4double GetZHalfLength() const
G4double GetXHalfLength() const

Referenced by AddSolid().

◆ ConeWrite()

void G4GDMLWriteSolids::ConeWrite ( xercesc::DOMElement *  solElement,
const G4Cons * const  cone 
)
protected

Definition at line 210 of file G4GDMLWriteSolids.cc.

212{
213 const G4String& name = GenerateName(cone->GetName(),cone);
214
215 xercesc::DOMElement* coneElement = NewElement("cone");
216 coneElement->setAttributeNode(NewAttribute("name",name));
217 coneElement->
218 setAttributeNode(NewAttribute("rmin1",cone->GetInnerRadiusMinusZ()/mm));
219 coneElement->
220 setAttributeNode(NewAttribute("rmax1",cone->GetOuterRadiusMinusZ()/mm));
221 coneElement->
222 setAttributeNode(NewAttribute("rmin2",cone->GetInnerRadiusPlusZ()/mm));
223 coneElement->
224 setAttributeNode(NewAttribute("rmax2",cone->GetOuterRadiusPlusZ()/mm));
225 coneElement->
226 setAttributeNode(NewAttribute("z",2.0*cone->GetZHalfLength()/mm));
227 coneElement->
228 setAttributeNode(NewAttribute("startphi",cone->GetStartPhiAngle()/degree));
229 coneElement->
230 setAttributeNode(NewAttribute("deltaphi",cone->GetDeltaPhiAngle()/degree));
231 coneElement->setAttributeNode(NewAttribute("aunit","deg"));
232 coneElement->setAttributeNode(NewAttribute("lunit","mm"));
233 solElement->appendChild(coneElement);
234}
G4double GetOuterRadiusPlusZ() const
G4double GetStartPhiAngle() const
G4double GetDeltaPhiAngle() const
G4double GetInnerRadiusMinusZ() const
G4double GetInnerRadiusPlusZ() const
G4double GetOuterRadiusMinusZ() const
G4double GetZHalfLength() const

Referenced by AddSolid().

◆ CutTubeWrite()

void G4GDMLWriteSolids::CutTubeWrite ( xercesc::DOMElement *  solElement,
const G4CutTubs * const  cuttube 
)
protected

Definition at line 725 of file G4GDMLWriteSolids.cc.

727{
728 const G4String& name = GenerateName(cuttube->GetName(),cuttube);
729
730 xercesc::DOMElement* cuttubeElement = NewElement("cutTube");
731 cuttubeElement->setAttributeNode(NewAttribute("name",name));
732 cuttubeElement->setAttributeNode(NewAttribute("rmin",
733 cuttube->GetInnerRadius()/mm));
734 cuttubeElement->setAttributeNode(NewAttribute("rmax",
735 cuttube->GetOuterRadius()/mm));
736 cuttubeElement->setAttributeNode(NewAttribute("z",
737 2.0*cuttube->GetZHalfLength()/mm));
738 cuttubeElement->setAttributeNode(NewAttribute("startphi",
739 cuttube->GetStartPhiAngle()/degree));
740 cuttubeElement->setAttributeNode(NewAttribute("deltaphi",
741 cuttube->GetDeltaPhiAngle()/degree));
742 cuttubeElement->setAttributeNode(NewAttribute("lowX",
743 cuttube->GetLowNorm().getX()/mm));
744 cuttubeElement->setAttributeNode(NewAttribute("lowY",
745 cuttube->GetLowNorm().getY()/mm));
746 cuttubeElement->setAttributeNode(NewAttribute("lowZ",
747 cuttube->GetLowNorm().getZ()/mm));
748 cuttubeElement->setAttributeNode(NewAttribute("highX",
749 cuttube->GetHighNorm().getX()/mm));
750 cuttubeElement->setAttributeNode(NewAttribute("highY",
751 cuttube->GetHighNorm().getY()/mm));
752 cuttubeElement->setAttributeNode(NewAttribute("highZ",
753 cuttube->GetHighNorm().getZ()/mm));
754 cuttubeElement->setAttributeNode(NewAttribute("aunit","deg"));
755 cuttubeElement->setAttributeNode(NewAttribute("lunit","mm"));
756 solElement->appendChild(cuttubeElement);
757}
double getZ() const
double getX() const
double getY() const
G4ThreeVector GetHighNorm() const
G4ThreeVector GetLowNorm() const
G4double GetZHalfLength() const
G4double GetInnerRadius() const
G4double GetOuterRadius() const
G4double GetStartPhiAngle() const
G4double GetDeltaPhiAngle() const

Referenced by AddSolid().

◆ ElconeWrite()

void G4GDMLWriteSolids::ElconeWrite ( xercesc::DOMElement *  solElement,
const G4EllipticalCone * const  elcone 
)
protected

Definition at line 236 of file G4GDMLWriteSolids.cc.

239{
240 const G4String& name = GenerateName(elcone->GetName(),elcone);
241
242 xercesc::DOMElement* elconeElement = NewElement("elcone");
243 elconeElement->setAttributeNode(NewAttribute("name",name));
244 elconeElement->setAttributeNode(NewAttribute("dx",elcone->GetSemiAxisX()/mm));
245 elconeElement->setAttributeNode(NewAttribute("dy",elcone->GetSemiAxisY()/mm));
246 elconeElement->setAttributeNode(NewAttribute("zmax",elcone->GetZMax()/mm));
247 elconeElement->setAttributeNode(NewAttribute("zcut",elcone->GetZTopCut()/mm));
248 elconeElement->setAttributeNode(NewAttribute("lunit","mm"));
249 solElement->appendChild(elconeElement);
250}
G4double GetSemiAxisX() const
G4double GetSemiAxisY() const
G4double GetZMax() const
G4double GetZTopCut() const

Referenced by AddSolid().

◆ EllipsoidWrite()

void G4GDMLWriteSolids::EllipsoidWrite ( xercesc::DOMElement *  solElement,
const G4Ellipsoid * const  ellipsoid 
)
protected

Definition at line 252 of file G4GDMLWriteSolids.cc.

255{
256 const G4String& name = GenerateName(ellipsoid->GetName(),ellipsoid);
257
258 xercesc::DOMElement* ellipsoidElement = NewElement("ellipsoid");
259 ellipsoidElement->setAttributeNode(NewAttribute("name",name));
260 ellipsoidElement->
261 setAttributeNode(NewAttribute("ax",ellipsoid->GetSemiAxisMax(0)/mm));
262 ellipsoidElement->
263 setAttributeNode(NewAttribute("by",ellipsoid->GetSemiAxisMax(1)/mm));
264 ellipsoidElement->
265 setAttributeNode(NewAttribute("cz",ellipsoid->GetSemiAxisMax(2)/mm));
266 ellipsoidElement->
267 setAttributeNode(NewAttribute("zcut1",ellipsoid->GetZBottomCut()/mm));
268 ellipsoidElement->
269 setAttributeNode(NewAttribute("zcut2",ellipsoid->GetZTopCut()/mm));
270 ellipsoidElement->
271 setAttributeNode(NewAttribute("lunit","mm"));
272 solElement->appendChild(ellipsoidElement);
273}
G4double GetSemiAxisMax(G4int i) const
G4double GetZTopCut() const
G4double GetZBottomCut() const

Referenced by AddSolid().

◆ EltubeWrite()

void G4GDMLWriteSolids::EltubeWrite ( xercesc::DOMElement *  solElement,
const G4EllipticalTube * const  eltube 
)
protected

Definition at line 275 of file G4GDMLWriteSolids.cc.

278{
279 const G4String& name = GenerateName(eltube->GetName(),eltube);
280
281 xercesc::DOMElement* eltubeElement = NewElement("eltube");
282 eltubeElement->setAttributeNode(NewAttribute("name",name));
283 eltubeElement->setAttributeNode(NewAttribute("dx",eltube->GetDx()/mm));
284 eltubeElement->setAttributeNode(NewAttribute("dy",eltube->GetDy()/mm));
285 eltubeElement->setAttributeNode(NewAttribute("dz",eltube->GetDz()/mm));
286 eltubeElement->setAttributeNode(NewAttribute("lunit","mm"));
287 solElement->appendChild(eltubeElement);
288}
G4double GetDy() const
G4double GetDx() const
G4double GetDz() const

Referenced by AddSolid().

◆ GenTrapWrite()

void G4GDMLWriteSolids::GenTrapWrite ( xercesc::DOMElement *  solElement,
const G4GenericTrap * const  gtrap 
)
protected

Definition at line 613 of file G4GDMLWriteSolids.cc.

616{
617 const G4String& name = GenerateName(gtrap->GetName(),gtrap);
618
619 std::vector<G4TwoVector> vertices = gtrap->GetVertices();
620
621 xercesc::DOMElement* gtrapElement = NewElement("arb8");
622 gtrapElement->setAttributeNode(NewAttribute("name",name));
623 gtrapElement->setAttributeNode(NewAttribute("dz",
624 gtrap->GetZHalfLength()/mm));
625 gtrapElement->setAttributeNode(NewAttribute("v1x", vertices[0].x()));
626 gtrapElement->setAttributeNode(NewAttribute("v1y", vertices[0].y()));
627 gtrapElement->setAttributeNode(NewAttribute("v2x", vertices[1].x()));
628 gtrapElement->setAttributeNode(NewAttribute("v2y", vertices[1].y()));
629 gtrapElement->setAttributeNode(NewAttribute("v3x", vertices[2].x()));
630 gtrapElement->setAttributeNode(NewAttribute("v3y", vertices[2].y()));
631 gtrapElement->setAttributeNode(NewAttribute("v4x", vertices[3].x()));
632 gtrapElement->setAttributeNode(NewAttribute("v4y", vertices[3].y()));
633 gtrapElement->setAttributeNode(NewAttribute("v5x", vertices[4].x()));
634 gtrapElement->setAttributeNode(NewAttribute("v5y", vertices[4].y()));
635 gtrapElement->setAttributeNode(NewAttribute("v6x", vertices[5].x()));
636 gtrapElement->setAttributeNode(NewAttribute("v6y", vertices[5].y()));
637 gtrapElement->setAttributeNode(NewAttribute("v7x", vertices[6].x()));
638 gtrapElement->setAttributeNode(NewAttribute("v7y", vertices[6].y()));
639 gtrapElement->setAttributeNode(NewAttribute("v8x", vertices[7].x()));
640 gtrapElement->setAttributeNode(NewAttribute("v8y", vertices[7].y()));
641 gtrapElement->setAttributeNode(NewAttribute("lunit","mm"));
642 solElement->appendChild(gtrapElement);
643}
G4double GetZHalfLength() const
const std::vector< G4TwoVector > & GetVertices() const

Referenced by AddSolid().

◆ HypeWrite()

void G4GDMLWriteSolids::HypeWrite ( xercesc::DOMElement *  solElement,
const G4Hype * const  hype 
)
protected

Definition at line 334 of file G4GDMLWriteSolids.cc.

336{
337 const G4String& name = GenerateName(hype->GetName(),hype);
338
339 xercesc::DOMElement* hypeElement = NewElement("hype");
340 hypeElement->setAttributeNode(NewAttribute("name",name));
341 hypeElement->setAttributeNode(NewAttribute("rmin",
342 hype->GetInnerRadius()/mm));
343 hypeElement->setAttributeNode(NewAttribute("rmax",
344 hype->GetOuterRadius()/mm));
345 hypeElement->setAttributeNode(NewAttribute("inst",
346 hype->GetInnerStereo()/degree));
347 hypeElement->setAttributeNode(NewAttribute("outst",
348 hype->GetOuterStereo()/degree));
349 hypeElement->setAttributeNode(NewAttribute("z",
350 2.0*hype->GetZHalfLength()/mm));
351 hypeElement->setAttributeNode(NewAttribute("aunit","deg"));
352 hypeElement->setAttributeNode(NewAttribute("lunit","mm"));
353 solElement->appendChild(hypeElement);
354}
G4double GetInnerStereo() const
G4double GetZHalfLength() const
G4double GetOuterStereo() const
G4double GetOuterRadius() const
G4double GetInnerRadius() const

Referenced by AddSolid().

◆ OpticalSurfaceWrite()

void G4GDMLWriteSolids::OpticalSurfaceWrite ( xercesc::DOMElement *  solElement,
const G4OpticalSurface * const  surf 
)
protected

Definition at line 875 of file G4GDMLWriteSolids.cc.

878{
879 xercesc::DOMElement* optElement = NewElement("opticalsurface");
880 G4OpticalSurfaceModel smodel = surf->GetModel();
881 G4double sval = (smodel==glisur) ? surf->GetPolish() : surf->GetSigmaAlpha();
882
883 optElement->setAttributeNode(NewAttribute("name", surf->GetName()));
884 optElement->setAttributeNode(NewAttribute("model", smodel));
885 optElement->setAttributeNode(NewAttribute("finish", surf->GetFinish()));
886 optElement->setAttributeNode(NewAttribute("type", surf->GetType()));
887 optElement->setAttributeNode(NewAttribute("value", sval));
888
889 solElement->appendChild(optElement);
890}
G4OpticalSurfaceModel
@ glisur
double G4double
Definition: G4Types.hh:64
G4OpticalSurfaceModel GetModel() const
G4double GetSigmaAlpha() const
G4OpticalSurfaceFinish GetFinish() const
G4double GetPolish() const
const G4String & GetName() const
const G4SurfaceType & GetType() const

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

◆ OrbWrite()

void G4GDMLWriteSolids::OrbWrite ( xercesc::DOMElement *  solElement,
const G4Orb * const  orb 
)
protected

Definition at line 356 of file G4GDMLWriteSolids.cc.

358{
359 const G4String& name = GenerateName(orb->GetName(),orb);
360
361 xercesc::DOMElement* orbElement = NewElement("orb");
362 orbElement->setAttributeNode(NewAttribute("name",name));
363 orbElement->setAttributeNode(NewAttribute("r",orb->GetRadius()/mm));
364 orbElement->setAttributeNode(NewAttribute("lunit","mm"));
365 solElement->appendChild(orbElement);
366}
G4double GetRadius() const

Referenced by AddSolid().

◆ ParaboloidWrite()

void G4GDMLWriteSolids::ParaboloidWrite ( xercesc::DOMElement *  solElement,
const G4Paraboloid * const  paraboloid 
)
protected

Definition at line 395 of file G4GDMLWriteSolids.cc.

398{
399 const G4String& name = GenerateName(paraboloid->GetName(),paraboloid);
400
401 xercesc::DOMElement* paraboloidElement = NewElement("paraboloid");
402 paraboloidElement->setAttributeNode(NewAttribute("name",name));
403 paraboloidElement->setAttributeNode(NewAttribute("rlo",
404 paraboloid->GetRadiusMinusZ()/mm));
405 paraboloidElement->setAttributeNode(NewAttribute("rhi",
406 paraboloid->GetRadiusPlusZ()/mm));
407 paraboloidElement->setAttributeNode(NewAttribute("dz",
408 paraboloid->GetZHalfLength()/mm));
409 paraboloidElement->setAttributeNode(NewAttribute("lunit","mm"));
410 solElement->appendChild(paraboloidElement);
411}
G4double GetRadiusPlusZ() const
G4double GetRadiusMinusZ() const
G4double GetZHalfLength() const

Referenced by AddSolid().

◆ ParaWrite()

void G4GDMLWriteSolids::ParaWrite ( xercesc::DOMElement *  solElement,
const G4Para * const  para 
)
protected

Definition at line 368 of file G4GDMLWriteSolids.cc.

370{
371 const G4String& name = GenerateName(para->GetName(),para);
372
373 const G4ThreeVector simaxis = para->GetSymAxis();
374 const G4double alpha = std::atan(para->GetTanAlpha());
375 const G4double theta = std::acos(simaxis.z());
376 const G4double phi = (simaxis.z() != 1.0)
377 ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
378
379 xercesc::DOMElement* paraElement = NewElement("para");
380 paraElement->setAttributeNode(NewAttribute("name",name));
381 paraElement->setAttributeNode(NewAttribute("x",
382 2.0*para->GetXHalfLength()/mm));
383 paraElement->setAttributeNode(NewAttribute("y",
384 2.0*para->GetYHalfLength()/mm));
385 paraElement->setAttributeNode(NewAttribute("z",
386 2.0*para->GetZHalfLength()/mm));
387 paraElement->setAttributeNode(NewAttribute("alpha",alpha/degree));
388 paraElement->setAttributeNode(NewAttribute("theta",theta/degree));
389 paraElement->setAttributeNode(NewAttribute("phi",phi/degree));
390 paraElement->setAttributeNode(NewAttribute("aunit","deg"));
391 paraElement->setAttributeNode(NewAttribute("lunit","mm"));
392 solElement->appendChild(paraElement);
393}
G4double GetTanAlpha() const
G4ThreeVector GetSymAxis() const
G4double GetYHalfLength() const
G4double GetZHalfLength() const
G4double GetXHalfLength() const

Referenced by AddSolid().

◆ PolyconeWrite()

void G4GDMLWriteSolids::PolyconeWrite ( xercesc::DOMElement *  solElement,
const G4Polycone * const  polycone 
)
protected

Definition at line 413 of file G4GDMLWriteSolids.cc.

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",
422 polycone->GetOriginalParameters()->Start_angle/degree));
423 polyconeElement->setAttributeNode(NewAttribute("deltaphi",
424 polycone->GetOriginalParameters()->Opening_angle/degree));
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}
void ZplaneWrite(xercesc::DOMElement *, const G4double &, const G4double &, const G4double &)
G4double Opening_angle
Definition: G4Polycone.hh:78
G4double * Z_values
Definition: G4Polycone.hh:80
G4PolyconeHistorical * GetOriginalParameters() const

Referenced by AddSolid().

◆ PolyhedraWrite()

void G4GDMLWriteSolids::PolyhedraWrite ( xercesc::DOMElement *  solElement,
const G4Polyhedra * const  polyhedra 
)
protected

Definition at line 440 of file G4GDMLWriteSolids.cc.

443{
444 const G4String& name = GenerateName(polyhedra->GetName(),polyhedra);
445
446 xercesc::DOMElement* polyhedraElement = NewElement("polyhedra");
447 polyhedraElement->setAttributeNode(NewAttribute("name",name));
448 polyhedraElement->setAttributeNode(NewAttribute("startphi",
449 polyhedra->GetOriginalParameters()->Start_angle/degree));
450 polyhedraElement->setAttributeNode(NewAttribute("deltaphi",
451 polyhedra->GetOriginalParameters()->Opening_angle/degree));
452 polyhedraElement->setAttributeNode(NewAttribute("numsides",
453 polyhedra->GetOriginalParameters()->numSide));
454 polyhedraElement->setAttributeNode(NewAttribute("aunit","deg"));
455 polyhedraElement->setAttributeNode(NewAttribute("lunit","mm"));
456 solElement->appendChild(polyhedraElement);
457
458 const size_t num_zplanes = polyhedra->GetOriginalParameters()->Num_z_planes;
459 const G4double* z_array = polyhedra->GetOriginalParameters()->Z_values;
460 const G4double* rmin_array = polyhedra->GetOriginalParameters()->Rmin;
461 const G4double* rmax_array = polyhedra->GetOriginalParameters()->Rmax;
462
463 const G4double convertRad =
464 std::cos(0.5*polyhedra->GetOriginalParameters()->Opening_angle
465 / polyhedra->GetOriginalParameters()->numSide);
466
467 for (size_t i=0;i<num_zplanes;i++)
468 {
469 ZplaneWrite(polyhedraElement,z_array[i],
470 rmin_array[i]*convertRad, rmax_array[i]*convertRad);
471 }
472}
G4PolyhedraHistorical * GetOriginalParameters() const

Referenced by AddSolid().

◆ SolidsWrite()

void G4GDMLWriteSolids::SolidsWrite ( xercesc::DOMElement *  gdmlElement)
virtual

Implements G4GDMLWrite.

Definition at line 892 of file G4GDMLWriteSolids.cc.

893{
894 G4cout << "G4GDML: Writing solids..." << G4endl;
895
896 solidsElement = NewElement("solids");
897 gdmlElement->appendChild(solidsElement);
898
899 solidList.clear();
900}
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout

◆ SphereWrite()

void G4GDMLWriteSolids::SphereWrite ( xercesc::DOMElement *  solElement,
const G4Sphere * const  sphere 
)
protected

Definition at line 474 of file G4GDMLWriteSolids.cc.

476{
477 const G4String& name = GenerateName(sphere->GetName(),sphere);
478
479 xercesc::DOMElement* sphereElement = NewElement("sphere");
480 sphereElement->setAttributeNode(NewAttribute("name",name));
481 sphereElement->setAttributeNode(NewAttribute("rmin",
482 sphere->GetInsideRadius()/mm));
483 sphereElement->setAttributeNode(NewAttribute("rmax",
484 sphere->GetOuterRadius()/mm));
485 sphereElement->setAttributeNode(NewAttribute("startphi",
486 sphere->GetStartPhiAngle()/degree));
487 sphereElement->setAttributeNode(NewAttribute("deltaphi",
488 sphere->GetDeltaPhiAngle()/degree));
489 sphereElement->setAttributeNode(NewAttribute("starttheta",
490 sphere->GetStartThetaAngle()/degree));
491 sphereElement->setAttributeNode(NewAttribute("deltatheta",
492 sphere->GetDeltaThetaAngle()/degree));
493 sphereElement->setAttributeNode(NewAttribute("aunit","deg"));
494 sphereElement->setAttributeNode(NewAttribute("lunit","mm"));
495 solElement->appendChild(sphereElement);
496}
G4double GetStartPhiAngle() const
G4double GetDeltaPhiAngle() const
G4double GetOuterRadius() const
G4double GetDeltaThetaAngle() const
G4double GetStartThetaAngle() const
G4double GetInsideRadius() const

Referenced by AddSolid().

◆ TessellatedWrite()

void G4GDMLWriteSolids::TessellatedWrite ( xercesc::DOMElement *  solElement,
const G4TessellatedSolid * const  tessellated 
)
protected

Definition at line 498 of file G4GDMLWriteSolids.cc.

501{
502 const G4String& solid_name = tessellated->GetName();
503 const G4String& name = GenerateName(solid_name, tessellated);
504
505 xercesc::DOMElement* tessellatedElement = NewElement("tessellated");
506 tessellatedElement->setAttributeNode(NewAttribute("name",name));
507 tessellatedElement->setAttributeNode(NewAttribute("aunit","deg"));
508 tessellatedElement->setAttributeNode(NewAttribute("lunit","mm"));
509 solElement->appendChild(tessellatedElement);
510
511 std::map<G4ThreeVector, G4String> vertexMap;
512
513 const size_t NumFacets = tessellated->GetNumberOfFacets();
514 size_t NumVertex = 0;
515
516 for (size_t i=0;i<NumFacets;i++)
517 {
518 const G4VFacet* facet = tessellated->GetFacet(i);
519 const size_t NumVertexPerFacet = facet->GetNumberOfVertices();
520
521 G4String FacetTag;
522
523 if (NumVertexPerFacet==3) { FacetTag="triangular"; } else
524 if (NumVertexPerFacet==4) { FacetTag="quadrangular"; }
525 else
526 {
527 G4Exception("G4GDMLWriteSolids::TessellatedWrite()", "InvalidSetup",
528 FatalException, "Facet should contain 3 or 4 vertices!");
529 }
530
531 xercesc::DOMElement* facetElement = NewElement(FacetTag);
532 tessellatedElement->appendChild(facetElement);
533
534 for (size_t j=0; j<NumVertexPerFacet; j++)
535 {
536 std::stringstream name_stream;
537 std::stringstream ref_stream;
538
539 name_stream << "vertex" << (j+1);
540 ref_stream << solid_name << "_v" << NumVertex;
541
542 const G4String& fname = name_stream.str(); // facet's tag variable
543 G4String ref = ref_stream.str(); // vertex tag to be associated
544
545 // Now search for the existance of the current vertex in the
546 // map of cached vertices. If existing, do NOT store it as
547 // position in the GDML file, so avoiding duplication; otherwise
548 // cache it in the local map and add it as position in the
549 // "define" section of the GDML file.
550
551 const G4ThreeVector& vertex = facet->GetVertex(j);
552
553 if(vertexMap.find(vertex) != vertexMap.end()) // Vertex is cached
554 {
555 ref = vertexMap[vertex]; // Set the proper tag for it
556 }
557 else // Vertex not found
558 {
559 vertexMap.insert(std::make_pair(vertex,ref)); // Cache vertex and ...
560 AddPosition(ref, vertex); // ... add it to define section!
561 NumVertex++;
562 }
563
564 // Now create association of the vertex with its facet
565 //
566 facetElement->setAttributeNode(NewAttribute(fname,ref));
567 }
568 }
569}
void AddPosition(const G4String &name, const G4ThreeVector &pos)
G4int GetNumberOfFacets() const
G4VFacet * GetFacet(G4int i) const
virtual G4ThreeVector GetVertex(G4int i) const =0
virtual G4int GetNumberOfVertices() const =0

Referenced by AddSolid().

◆ TetWrite()

void G4GDMLWriteSolids::TetWrite ( xercesc::DOMElement *  solElement,
const G4Tet * const  tet 
)
protected

Definition at line 571 of file G4GDMLWriteSolids.cc.

573{
574 const G4String& solid_name = tet->GetName();
575 const G4String& name = GenerateName(solid_name, tet);
576
577 std::vector<G4ThreeVector> vertexList = tet->GetVertices();
578
579 xercesc::DOMElement* tetElement = NewElement("tet");
580 tetElement->setAttributeNode(NewAttribute("name",name));
581 tetElement->setAttributeNode(NewAttribute("vertex1",solid_name+"_v1"));
582 tetElement->setAttributeNode(NewAttribute("vertex2",solid_name+"_v2"));
583 tetElement->setAttributeNode(NewAttribute("vertex3",solid_name+"_v3"));
584 tetElement->setAttributeNode(NewAttribute("vertex4",solid_name+"_v4"));
585 tetElement->setAttributeNode(NewAttribute("lunit","mm"));
586 solElement->appendChild(tetElement);
587
588 AddPosition(solid_name+"_v1",vertexList[0]);
589 AddPosition(solid_name+"_v2",vertexList[1]);
590 AddPosition(solid_name+"_v3",vertexList[2]);
591 AddPosition(solid_name+"_v4",vertexList[3]);
592}
std::vector< G4ThreeVector > GetVertices() const
Definition: G4Tet.cc:739

Referenced by AddSolid().

◆ TorusWrite()

void G4GDMLWriteSolids::TorusWrite ( xercesc::DOMElement *  solElement,
const G4Torus * const  torus 
)
protected

Definition at line 594 of file G4GDMLWriteSolids.cc.

596{
597 const G4String& name = GenerateName(torus->GetName(),torus);
598
599 xercesc::DOMElement* torusElement = NewElement("torus");
600 torusElement->setAttributeNode(NewAttribute("name",name));
601 torusElement->setAttributeNode(NewAttribute("rmin",torus->GetRmin()/mm));
602 torusElement->setAttributeNode(NewAttribute("rmax",torus->GetRmax()/mm));
603 torusElement->setAttributeNode(NewAttribute("rtor",torus->GetRtor()/mm));
604 torusElement->
605 setAttributeNode(NewAttribute("startphi",torus->GetSPhi()/degree));
606 torusElement->
607 setAttributeNode(NewAttribute("deltaphi",torus->GetDPhi()/degree));
608 torusElement->setAttributeNode(NewAttribute("aunit","deg"));
609 torusElement->setAttributeNode(NewAttribute("lunit","mm"));
610 solElement->appendChild(torusElement);
611}
G4double GetDPhi() const
G4double GetRmin() const
G4double GetRtor() const
G4double GetRmax() const
G4double GetSPhi() const

Referenced by AddSolid().

◆ TrapWrite()

void G4GDMLWriteSolids::TrapWrite ( xercesc::DOMElement *  solElement,
const G4Trap * const  trap 
)
protected

Definition at line 645 of file G4GDMLWriteSolids.cc.

647{
648 const G4String& name = GenerateName(trap->GetName(),trap);
649
650 const G4ThreeVector& simaxis = trap->GetSymAxis();
651 const G4double phi = (simaxis.z() != 1.0)
652 ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
653 const G4double theta = std::acos(simaxis.z());
654 const G4double alpha1 = std::atan(trap->GetTanAlpha1());
655 const G4double alpha2 = std::atan(trap->GetTanAlpha2());
656
657 xercesc::DOMElement* trapElement = NewElement("trap");
658 trapElement->setAttributeNode(NewAttribute("name",name));
659 trapElement->setAttributeNode(NewAttribute("z",
660 2.0*trap->GetZHalfLength()/mm));
661 trapElement->setAttributeNode(NewAttribute("theta",theta/degree));
662 trapElement->setAttributeNode(NewAttribute("phi",phi/degree));
663 trapElement->setAttributeNode(NewAttribute("y1",
664 2.0*trap->GetYHalfLength1()/mm));
665 trapElement->setAttributeNode(NewAttribute("x1",
666 2.0*trap->GetXHalfLength1()/mm));
667 trapElement->setAttributeNode(NewAttribute("x2",
668 2.0*trap->GetXHalfLength2()/mm));
669 trapElement->setAttributeNode(NewAttribute("alpha1",alpha1/degree));
670 trapElement->setAttributeNode(NewAttribute("y2",
671 2.0*trap->GetYHalfLength2()/mm));
672 trapElement->setAttributeNode(NewAttribute("x3",
673 2.0*trap->GetXHalfLength3()/mm));
674 trapElement->setAttributeNode(NewAttribute("x4",
675 2.0*trap->GetXHalfLength4()/mm));
676 trapElement->setAttributeNode(NewAttribute("alpha2",alpha2/degree));
677 trapElement->setAttributeNode(NewAttribute("aunit","deg"));
678 trapElement->setAttributeNode(NewAttribute("lunit","mm"));
679 solElement->appendChild(trapElement);
680}
G4double GetYHalfLength1() const
G4double GetTanAlpha2() const
G4double GetXHalfLength2() const
G4ThreeVector GetSymAxis() const
G4double GetXHalfLength4() const
G4double GetZHalfLength() const
G4double GetYHalfLength2() const
G4double GetTanAlpha1() const
G4double GetXHalfLength3() const
G4double GetXHalfLength1() const

Referenced by AddSolid().

◆ TrdWrite()

void G4GDMLWriteSolids::TrdWrite ( xercesc::DOMElement *  solElement,
const G4Trd * const  trd 
)
protected

Definition at line 682 of file G4GDMLWriteSolids.cc.

684{
685 const G4String& name = GenerateName(trd->GetName(),trd);
686
687 xercesc::DOMElement* trdElement = NewElement("trd");
688 trdElement->setAttributeNode(NewAttribute("name",name));
689 trdElement->setAttributeNode(NewAttribute("x1",
690 2.0*trd->GetXHalfLength1()/mm));
691 trdElement->setAttributeNode(NewAttribute("x2",
692 2.0*trd->GetXHalfLength2()/mm));
693 trdElement->setAttributeNode(NewAttribute("y1",
694 2.0*trd->GetYHalfLength1()/mm));
695 trdElement->setAttributeNode(NewAttribute("y2",
696 2.0*trd->GetYHalfLength2()/mm));
697 trdElement->setAttributeNode(NewAttribute("z",
698 2.0*trd->GetZHalfLength()/mm));
699 trdElement->setAttributeNode(NewAttribute("lunit","mm"));
700 solElement->appendChild(trdElement);
701}
G4double GetXHalfLength2() const
G4double GetYHalfLength2() const
G4double GetXHalfLength1() const
G4double GetYHalfLength1() const
G4double GetZHalfLength() const

Referenced by AddSolid().

◆ TubeWrite()

void G4GDMLWriteSolids::TubeWrite ( xercesc::DOMElement *  solElement,
const G4Tubs * const  tube 
)
protected

Definition at line 703 of file G4GDMLWriteSolids.cc.

705{
706 const G4String& name = GenerateName(tube->GetName(),tube);
707
708 xercesc::DOMElement* tubeElement = NewElement("tube");
709 tubeElement->setAttributeNode(NewAttribute("name",name));
710 tubeElement->setAttributeNode(NewAttribute("rmin",
711 tube->GetInnerRadius()/mm));
712 tubeElement->setAttributeNode(NewAttribute("rmax",
713 tube->GetOuterRadius()/mm));
714 tubeElement->setAttributeNode(NewAttribute("z",
715 2.0*tube->GetZHalfLength()/mm));
716 tubeElement->setAttributeNode(NewAttribute("startphi",
717 tube->GetStartPhiAngle()/degree));
718 tubeElement->setAttributeNode(NewAttribute("deltaphi",
719 tube->GetDeltaPhiAngle()/degree));
720 tubeElement->setAttributeNode(NewAttribute("aunit","deg"));
721 tubeElement->setAttributeNode(NewAttribute("lunit","mm"));
722 solElement->appendChild(tubeElement);
723}

Referenced by AddSolid().

◆ TwistedboxWrite()

void G4GDMLWriteSolids::TwistedboxWrite ( xercesc::DOMElement *  solElement,
const G4TwistedBox * const  twistedbox 
)
protected

Definition at line 759 of file G4GDMLWriteSolids.cc.

762{
763 const G4String& name = GenerateName(twistedbox->GetName(),twistedbox);
764
765 xercesc::DOMElement* twistedboxElement = NewElement("twistedbox");
766 twistedboxElement->setAttributeNode(NewAttribute("name",name));
767 twistedboxElement->setAttributeNode(NewAttribute("x",
768 2.0*twistedbox->GetXHalfLength()/mm));
769 twistedboxElement->setAttributeNode(NewAttribute("y",
770 2.0*twistedbox->GetYHalfLength()/mm));
771 twistedboxElement->setAttributeNode(NewAttribute("z",
772 2.0*twistedbox->GetZHalfLength()/mm));
773 twistedboxElement->setAttributeNode(NewAttribute("PhiTwist",
774 twistedbox->GetPhiTwist()/degree));
775 twistedboxElement->setAttributeNode(NewAttribute("aunit","deg"));
776 twistedboxElement->setAttributeNode(NewAttribute("lunit","mm"));
777 solElement->appendChild(twistedboxElement);
778}
G4double GetPhiTwist() const
Definition: G4TwistedBox.hh:76
G4double GetXHalfLength() const
Definition: G4TwistedBox.hh:73
G4double GetZHalfLength() const
Definition: G4TwistedBox.hh:75
G4double GetYHalfLength() const
Definition: G4TwistedBox.hh:74

Referenced by AddSolid().

◆ TwistedtrapWrite()

void G4GDMLWriteSolids::TwistedtrapWrite ( xercesc::DOMElement *  solElement,
const G4TwistedTrap * const  twistedtrap 
)
protected

Definition at line 780 of file G4GDMLWriteSolids.cc.

783{
784 const G4String& name = GenerateName(twistedtrap->GetName(),twistedtrap);
785
786 xercesc::DOMElement* twistedtrapElement = NewElement("twistedtrap");
787 twistedtrapElement->setAttributeNode(NewAttribute("name",name));
788 twistedtrapElement->setAttributeNode(NewAttribute("y1",
789 2.0*twistedtrap->GetY1HalfLength()/mm));
790 twistedtrapElement->setAttributeNode(NewAttribute("x1",
791 2.0*twistedtrap->GetX1HalfLength()/mm));
792 twistedtrapElement->setAttributeNode(NewAttribute("x2",
793 2.0*twistedtrap->GetX2HalfLength()/mm));
794 twistedtrapElement->setAttributeNode(NewAttribute("y2",
795 2.0*twistedtrap->GetY2HalfLength()/mm));
796 twistedtrapElement->setAttributeNode(NewAttribute("x3",
797 2.0*twistedtrap->GetX3HalfLength()/mm));
798 twistedtrapElement->setAttributeNode(NewAttribute("x4",
799 2.0*twistedtrap->GetX4HalfLength()/mm));
800 twistedtrapElement->setAttributeNode(NewAttribute("z",
801 2.0*twistedtrap->GetZHalfLength()/mm));
802 twistedtrapElement->setAttributeNode(NewAttribute("Alph",
803 twistedtrap->GetTiltAngleAlpha()/degree));
804 twistedtrapElement->setAttributeNode(NewAttribute("Theta",
805 twistedtrap->GetPolarAngleTheta()/degree));
806 twistedtrapElement->setAttributeNode(NewAttribute("Phi",
807 twistedtrap->GetAzimuthalAnglePhi()/degree));
808 twistedtrapElement->setAttributeNode(NewAttribute("PhiTwist",
809 twistedtrap->GetPhiTwist()/degree));
810 twistedtrapElement->setAttributeNode(NewAttribute("aunit","deg"));
811 twistedtrapElement->setAttributeNode(NewAttribute("lunit","mm"));
812
813 solElement->appendChild(twistedtrapElement);
814}
G4double GetPolarAngleTheta() const
G4double GetAzimuthalAnglePhi() const
G4double GetTiltAngleAlpha() const
G4double GetZHalfLength() const
G4double GetX1HalfLength() const
G4double GetX2HalfLength() const
G4double GetX3HalfLength() const
G4double GetX4HalfLength() const
G4double GetY2HalfLength() const
G4double GetPhiTwist() const
G4double GetY1HalfLength() const

Referenced by AddSolid().

◆ TwistedtrdWrite()

void G4GDMLWriteSolids::TwistedtrdWrite ( xercesc::DOMElement *  solElement,
const G4TwistedTrd * const  twistedtrd 
)
protected

Definition at line 816 of file G4GDMLWriteSolids.cc.

819{
820 const G4String& name = GenerateName(twistedtrd->GetName(),twistedtrd);
821
822 xercesc::DOMElement* twistedtrdElement = NewElement("twistedtrd");
823 twistedtrdElement->setAttributeNode(NewAttribute("name",name));
824 twistedtrdElement->setAttributeNode(NewAttribute("x1",
825 2.0*twistedtrd->GetX1HalfLength()/mm));
826 twistedtrdElement->setAttributeNode(NewAttribute("x2",
827 2.0*twistedtrd->GetX2HalfLength()/mm));
828 twistedtrdElement->setAttributeNode(NewAttribute("y1",
829 2.0*twistedtrd->GetY1HalfLength()/mm));
830 twistedtrdElement->setAttributeNode(NewAttribute("y2",
831 2.0*twistedtrd->GetY2HalfLength()/mm));
832 twistedtrdElement->setAttributeNode(NewAttribute("z",
833 2.0*twistedtrd->GetZHalfLength()/mm));
834 twistedtrdElement->setAttributeNode(NewAttribute("PhiTwist",
835 twistedtrd->GetPhiTwist()/degree));
836 twistedtrdElement->setAttributeNode(NewAttribute("aunit","deg"));
837 twistedtrdElement->setAttributeNode(NewAttribute("lunit","mm"));
838 solElement->appendChild(twistedtrdElement);
839}
G4double GetX2HalfLength() const
Definition: G4TwistedTrd.hh:78
G4double GetY2HalfLength() const
Definition: G4TwistedTrd.hh:80
G4double GetY1HalfLength() const
Definition: G4TwistedTrd.hh:79
G4double GetZHalfLength() const
Definition: G4TwistedTrd.hh:81
G4double GetPhiTwist() const
Definition: G4TwistedTrd.hh:82
G4double GetX1HalfLength() const
Definition: G4TwistedTrd.hh:77

Referenced by AddSolid().

◆ TwistedtubsWrite()

void G4GDMLWriteSolids::TwistedtubsWrite ( xercesc::DOMElement *  solElement,
const G4TwistedTubs * const  twistedtubs 
)
protected

Definition at line 841 of file G4GDMLWriteSolids.cc.

844{
845 const G4String& name = GenerateName(twistedtubs->GetName(),twistedtubs);
846
847 xercesc::DOMElement* twistedtubsElement = NewElement("twistedtubs");
848 twistedtubsElement->setAttributeNode(NewAttribute("name",name));
849 twistedtubsElement->setAttributeNode(NewAttribute("twistedangle",
850 twistedtubs->GetPhiTwist()/degree));
851 twistedtubsElement->setAttributeNode(NewAttribute("endinnerrad",
852 twistedtubs->GetInnerRadius()/mm));
853 twistedtubsElement->setAttributeNode(NewAttribute("endouterrad",
854 twistedtubs->GetOuterRadius()/mm));
855 twistedtubsElement->setAttributeNode(NewAttribute("zlen",
856 2.0*twistedtubs->GetZHalfLength()/mm));
857 twistedtubsElement->setAttributeNode(NewAttribute("phi",
858 twistedtubs->GetDPhi()/degree));
859 twistedtubsElement->setAttributeNode(NewAttribute("aunit","deg"));
860 twistedtubsElement->setAttributeNode(NewAttribute("lunit","mm"));
861 solElement->appendChild(twistedtubsElement);
862}
G4double GetOuterRadius() const
G4double GetZHalfLength() const
G4double GetPhiTwist() const
G4double GetInnerRadius() const
G4double GetDPhi() const

Referenced by AddSolid().

◆ XtruWrite()

void G4GDMLWriteSolids::XtruWrite ( xercesc::DOMElement *  solElement,
const G4ExtrudedSolid * const  xtru 
)
protected

Definition at line 290 of file G4GDMLWriteSolids.cc.

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

Referenced by AddSolid().

◆ ZplaneWrite()

void G4GDMLWriteSolids::ZplaneWrite ( xercesc::DOMElement *  element,
const G4double z,
const G4double rmin,
const G4double rmax 
)
protected

Definition at line 864 of file G4GDMLWriteSolids.cc.

867{
868 xercesc::DOMElement* zplaneElement = NewElement("zplane");
869 zplaneElement->setAttributeNode(NewAttribute("z",z/mm));
870 zplaneElement->setAttributeNode(NewAttribute("rmin",rmin/mm));
871 zplaneElement->setAttributeNode(NewAttribute("rmax",rmax/mm));
872 element->appendChild(zplaneElement);
873}

Referenced by PolyconeWrite(), and PolyhedraWrite().

Member Data Documentation

◆ maxTransforms

const G4int G4GDMLWriteSolids::maxTransforms = 8
staticprotected

Definition at line 123 of file G4GDMLWriteSolids.hh.

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

◆ solidList

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

Definition at line 121 of file G4GDMLWriteSolids.hh.

Referenced by AddSolid(), and SolidsWrite().

◆ solidsElement

xercesc::DOMElement* G4GDMLWriteSolids::solidsElement
protected

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