38 #include <Geant4/G4Material.hh>
39 #include <Geant4/G4ReflectedSolid.hh>
40 #include <Geant4/G4DisplacedSolid.hh>
41 #include <Geant4/G4LogicalVolumeStore.hh>
42 #include <Geant4/G4PhysicalVolumeStore.hh>
43 #include <Geant4/G4ReflectionFactory.hh>
44 #include <Geant4/G4PVDivision.hh>
45 #include <Geant4/G4PVReplica.hh>
46 #include <Geant4/G4Region.hh>
47 #include <Geant4/G4OpticalSurface.hh>
48 #include <Geant4/G4LogicalSkinSurface.hh>
49 #include <Geant4/G4LogicalBorderSurface.hh>
51 #include <Geant4/G4ProductionCuts.hh>
52 #include <Geant4/G4ProductionCutsTable.hh>
53 #include <Geant4/G4Gamma.hh>
54 #include <Geant4/G4Electron.hh>
55 #include <Geant4/G4Positron.hh>
56 #include <Geant4/G4Proton.hh>
76 const G4PVDivision*
const divisionvol)
78 EAxis axis = kUndefined;
82 G4bool consuming =
false;
84 divisionvol->GetReplicationData(axis,number,width,offset,consuming);
85 axis = divisionvol->GetDivisionAxis();
87 G4String unitString(
"mm");
88 G4String axisString(
"kUndefined");
89 if (axis==kXAxis) { axisString =
"kXAxis"; }
90 else if (axis==kYAxis) { axisString =
"kYAxis"; }
91 else if (axis==kZAxis) { axisString =
"kZAxis"; }
92 else if (axis==kRho) { axisString =
"kRho"; }
93 else if (axis==
kPhi) { axisString =
"kPhi"; unitString =
"rad"; }
97 const G4String volumeref
98 =
GenerateName(divisionvol->GetLogicalVolume()->GetName(),
99 divisionvol->GetLogicalVolume());
101 xercesc::DOMElement* divisionvolElement =
NewElement(
"divisionvol");
102 divisionvolElement->setAttributeNode(
NewAttribute(
"axis",axisString));
103 divisionvolElement->setAttributeNode(
NewAttribute(
"number",number));
104 divisionvolElement->setAttributeNode(
NewAttribute(
"width",width));
105 divisionvolElement->setAttributeNode(
NewAttribute(
"offset",offset));
106 divisionvolElement->setAttributeNode(
NewAttribute(
"unit",unitString));
107 xercesc::DOMElement* volumerefElement =
NewElement(
"volumeref");
108 volumerefElement->setAttributeNode(
NewAttribute(
"ref",volumeref));
109 divisionvolElement->appendChild(volumerefElement);
110 volumeElement->appendChild(divisionvolElement);
114 const G4VPhysicalVolume*
const physvol,
115 const G4Transform3D&
T,
116 const G4String& ModuleName)
118 HepGeom::Scale3D scale;
119 HepGeom::Rotate3D rotate;
120 HepGeom::Translate3D translate;
122 T.getDecomposition(scale,rotate,translate);
124 const G4ThreeVector scl(scale(0,0),scale(1,1),scale(2,2));
125 const G4ThreeVector rot =
GetAngles(rotate.getRotation());
126 const G4ThreeVector
pos = T.getTranslation();
129 const G4int copynumber = physvol->GetCopyNo();
131 xercesc::DOMElement* physvolElement =
NewElement(
"physvol");
132 physvolElement->setAttributeNode(
NewAttribute(
"name",name));
133 if (copynumber) physvolElement->setAttributeNode(
NewAttribute(
"copynumber", copynumber));
135 volumeElement->appendChild(physvolElement);
139 if (
reflFactory->IsReflected(physvol->GetLogicalVolume()))
141 lv =
reflFactory->GetConstituentLV(physvol->GetLogicalVolume());
145 lv = physvol->GetLogicalVolume();
148 const G4String volumeref =
GenerateName(lv->GetName(), lv);
150 if (ModuleName.empty())
152 xercesc::DOMElement* volumerefElement =
NewElement(
"volumeref");
153 volumerefElement->setAttributeNode(
NewAttribute(
"ref",volumeref));
154 physvolElement->appendChild(volumerefElement);
158 xercesc::DOMElement* fileElement =
NewElement(
"file");
159 fileElement->setAttributeNode(
NewAttribute(
"name",ModuleName));
160 fileElement->setAttributeNode(
NewAttribute(
"volname",volumeref));
161 physvolElement->appendChild(fileElement);
185 const G4VPhysicalVolume*
const replicavol)
187 EAxis axis = kUndefined;
189 G4double
width = 0.0;
191 G4bool consuming =
false;
192 G4String unitString(
"mm");
194 replicavol->GetReplicationData(axis,number,width,offset,consuming);
196 const G4String volumeref
197 =
GenerateName(replicavol->GetLogicalVolume()->GetName(),
198 replicavol->GetLogicalVolume());
200 xercesc::DOMElement* replicavolElement =
NewElement(
"replicavol");
201 replicavolElement->setAttributeNode(
NewAttribute(
"number",number));
202 xercesc::DOMElement* volumerefElement =
NewElement(
"volumeref");
203 volumerefElement->setAttributeNode(
NewAttribute(
"ref",volumeref));
204 replicavolElement->appendChild(volumerefElement);
205 xercesc::DOMElement* replicateElement =
NewElement(
"replicate_along_axis");
206 replicavolElement->appendChild(replicateElement);
208 xercesc::DOMElement* dirElement =
NewElement(
"direction");
210 { dirElement->setAttributeNode(
NewAttribute(
"x",
"1")); }
211 else if(axis==kYAxis)
212 { dirElement->setAttributeNode(
NewAttribute(
"y",
"1")); }
213 else if(axis==kZAxis)
214 { dirElement->setAttributeNode(
NewAttribute(
"z",
"1")); }
216 { dirElement->setAttributeNode(
NewAttribute(
"rho",
"1")); }
218 { dirElement->setAttributeNode(
NewAttribute(
"phi",
"1"));
220 replicateElement->appendChild(dirElement);
222 xercesc::DOMElement* widthElement =
NewElement(
"width");
223 widthElement->setAttributeNode(
NewAttribute(
"value",width));
224 widthElement->setAttributeNode(
NewAttribute(
"unit",unitString));
225 replicateElement->appendChild(widthElement);
227 xercesc::DOMElement* offsetElement =
NewElement(
"offset");
228 offsetElement->setAttributeNode(
NewAttribute(
"value",offset));
229 offsetElement->setAttributeNode(
NewAttribute(
"unit",unitString));
230 replicateElement->appendChild(offsetElement);
232 volumeElement->appendChild(replicavolElement);
238 if (!bsurf) {
return; }
240 const G4SurfaceProperty* psurf = bsurf->GetSurfaceProperty();
244 xercesc::DOMElement* borderElement =
NewElement(
"bordersurface");
245 borderElement->setAttributeNode(
NewAttribute(
"name", bsurf->GetName()));
246 borderElement->setAttributeNode(
NewAttribute(
"surfaceproperty",
249 const G4String volumeref1 =
GenerateName(bsurf->GetVolume1()->GetName(),
250 bsurf->GetVolume1());
251 const G4String volumeref2 =
GenerateName(bsurf->GetVolume2()->GetName(),
252 bsurf->GetVolume2());
253 xercesc::DOMElement* volumerefElement1 =
NewElement(
"physvolref");
254 xercesc::DOMElement* volumerefElement2 =
NewElement(
"physvolref");
255 volumerefElement1->setAttributeNode(
NewAttribute(
"ref",volumeref1));
256 volumerefElement2->setAttributeNode(
NewAttribute(
"ref",volumeref2));
257 borderElement->appendChild(volumerefElement1);
258 borderElement->appendChild(volumerefElement2);
262 const G4OpticalSurface* opsurf =
263 dynamic_cast<const G4OpticalSurface*
>(psurf);
266 G4Exception(
"PHG4GDMLWriteStructure::BorderSurfaceCache()",
267 "InvalidSetup", FatalException,
"No optical surface found!");
279 if (!ssurf) {
return; }
281 const G4SurfaceProperty* psurf = ssurf->GetSurfaceProperty();
285 xercesc::DOMElement* skinElement =
NewElement(
"skinsurface");
286 skinElement->setAttributeNode(
NewAttribute(
"name", ssurf->GetName()));
287 skinElement->setAttributeNode(
NewAttribute(
"surfaceproperty",
290 const G4String volumeref =
GenerateName(ssurf->GetLogicalVolume()->GetName(),
291 ssurf->GetLogicalVolume());
292 xercesc::DOMElement* volumerefElement =
NewElement(
"volumeref");
293 volumerefElement->setAttributeNode(
NewAttribute(
"ref",volumeref));
294 skinElement->appendChild(volumerefElement);
298 const G4OpticalSurface* opsurf =
299 dynamic_cast<const G4OpticalSurface*
>(psurf);
302 G4Exception(
"PHG4GDMLWriteStructure::SkinSurfaceCache()",
303 "InvalidSetup", FatalException,
"No optical surface found!");
314 const G4OpticalSurface* osurf =
dynamic_cast<const G4OpticalSurface*
>(psurf);
315 std::vector<const G4OpticalSurface*>::const_iterator
pos;
317 if (pos !=
opt_vec.end()) {
return false; }
323 const G4LogicalSkinSurface*
326 G4LogicalSkinSurface* surf = 0;
327 G4int nsurf = G4LogicalSkinSurface::GetNumberOfSkinSurfaces();
330 const G4LogicalSkinSurfaceTable* stable =
331 G4LogicalSkinSurface::GetSurfaceTable();
332 std::vector<G4LogicalSkinSurface*>::const_iterator
pos;
333 for (pos = stable->begin(); pos != stable->end(); ++
pos)
335 if (lvol == (*pos)->GetLogicalVolume())
344 const G4LogicalBorderSurface*
347 G4LogicalBorderSurface* surf = 0;
348 G4int nsurf = G4LogicalBorderSurface::GetNumberOfBorderSurfaces();
351 const G4LogicalBorderSurfaceTable* btable =
352 G4LogicalBorderSurface::GetSurfaceTable();
353 std::vector<G4LogicalBorderSurface*>::const_iterator
pos;
354 for (pos = btable->begin(); pos != btable->end(); ++
pos)
356 if (pvol == (*pos)->GetVolume1())
367 std::cout <<
"PHG4GDML: Writing surfaces..." << std::endl;
369 std::vector<xercesc::DOMElement*>::const_iterator
pos;
382 std::cout <<
"PHG4GDML: Writing structure..." << std::endl;
400 return G4Transform3D::Identity;
404 G4VSolid* solidPtr = volumePtr->GetSolid();
405 G4Transform3D
R,invR;
408 std::map<const G4LogicalVolume*, PHG4GDMLAuxListType>::iterator auxiter;
414 G4String ErrorMessage =
"Referenced solid in volume '"
415 + volumePtr->GetName()
416 +
"' was displaced/reflected too many times!";
417 G4Exception(
"PHG4GDMLWriteStructure::TraverseVolumeTree()",
418 "InvalidSetup", FatalException, ErrorMessage);
421 if (G4ReflectedSolid* refl = dynamic_cast<G4ReflectedSolid*>(solidPtr))
423 R = R*refl->GetTransform3D();
424 solidPtr = refl->GetConstituentMovedSolid();
429 if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(solidPtr))
431 R = R*G4Transform3D(disp->GetObjectRotation(),
432 disp->GetObjectTranslation());
433 solidPtr = disp->GetConstituentMovedSolid();
442 G4LogicalVolume* tmplv =
const_cast<G4LogicalVolume*
>(volumePtr);
455 if (trans>0) { invR = R.inverse(); }
459 const G4String materialref
461 volumePtr->GetMaterial());
462 const G4String solidref
465 xercesc::DOMElement* volumeElement =
NewElement(
"volume");
466 volumeElement->setAttributeNode(
NewAttribute(
"name",name));
467 xercesc::DOMElement* materialrefElement =
NewElement(
"materialref");
468 materialrefElement->setAttributeNode(
NewAttribute(
"ref",materialref));
469 volumeElement->appendChild(materialrefElement);
470 xercesc::DOMElement* solidrefElement =
NewElement(
"solidref");
471 solidrefElement->setAttributeNode(
NewAttribute(
"ref",solidref));
472 volumeElement->appendChild(solidrefElement);
474 const G4int daughterCount = volumePtr->GetNoDaughters();
476 for (G4int i=0;i<daughterCount;i++)
478 const G4VPhysicalVolume*
const physvol = volumePtr->GetDaughter(i);
487 const G4String ModuleName =
Modularize(physvol,depth);
489 G4Transform3D daughterR;
491 if (ModuleName.empty())
498 daughterR = writer.
Write(ModuleName,physvol->GetLogicalVolume(),
502 if (
const G4PVDivision*
const divisionvol
503 = dynamic_cast<const G4PVDivision*>(physvol))
507 G4String ErrorMessage =
"Division volume in '" + name
508 +
"' can not be related to reflected solid!";
509 G4Exception(
"PHG4GDMLWriteStructure::TraverseVolumeTree()",
510 "InvalidSetup", FatalException, ErrorMessage);
514 if (physvol->IsParameterised())
518 G4String ErrorMessage =
"Parameterised volume in '" + name
519 +
"' can not be related to reflected solid!";
520 G4Exception(
"PHG4GDMLWriteStructure::TraverseVolumeTree()",
521 "InvalidSetup", FatalException, ErrorMessage);
525 if (physvol->IsReplicated())
529 G4String ErrorMessage =
"Replica volume in '" + name
530 +
"' can not be related to reflected solid!";
531 G4Exception(
"PHG4GDMLWriteStructure::TraverseVolumeTree()",
532 "InvalidSetup", FatalException, ErrorMessage);
538 G4RotationMatrix rot;
539 if (physvol->GetFrameRotation() != 0)
541 rot = *(physvol->GetFrameRotation());
543 G4Transform3D P(rot,physvol->GetObjectTranslation());
545 PhysvolWrite(volumeElement,physvol,invR*P*daughterR,ModuleName);
555 auxiter =
auxmap.find(volumePtr);
556 if (auxiter !=
auxmap.end())
558 AddAuxInfo(&(auxiter->second), volumeElement);
582 const G4LogicalVolume*
const lvol)
584 std::map<
const G4LogicalVolume*,
589 auxmap[lvol].push_back(myaux);
602 G4ProductionCuts* pcuts = lvol->GetRegion()->GetProductionCuts();
603 G4ProductionCutsTable* ctab = G4ProductionCutsTable::GetProductionCutsTable();
604 G4Gamma*
gamma = G4Gamma::Gamma();
605 G4Electron* eminus = G4Electron::Electron();
606 G4Positron* eplus = G4Positron::Positron();
607 G4Proton* proton = G4Proton::Proton();
609 G4double gamma_cut = ctab->ConvertRangeToEnergy(gamma, lvol->GetMaterial(),
610 pcuts->GetProductionCut(
"gamma"));
611 G4double eminus_cut = ctab->ConvertRangeToEnergy(eminus, lvol->GetMaterial(),
612 pcuts->GetProductionCut(
"e-"));
613 G4double eplus_cut = ctab->ConvertRangeToEnergy(eplus, lvol->GetMaterial(),
614 pcuts->GetProductionCut(
"e+"));
615 G4double proton_cut = ctab->ConvertRangeToEnergy(proton, lvol->GetMaterial(),
616 pcuts->GetProductionCut(
"proton"));
635 std::ostringstream os;
637 G4String vl = os.str();