ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4GDMLWriteStructure.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4GDMLWriteStructure.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 //
28 // class G4GDMLWriteStructure Implementation
29 //
30 // Original author: Zoltan Torzsok, November 2007
31 //
32 // --------------------------------------------------------------------
33 
34 #include "G4GDMLWriteStructure.hh"
35 #include "G4GDMLEvaluator.hh"
36 
37 #include "G4Material.hh"
38 #include "G4ReflectedSolid.hh"
39 #include "G4DisplacedSolid.hh"
40 #include "G4LogicalVolumeStore.hh"
41 #include "G4PhysicalVolumeStore.hh"
42 #include "G4ReflectionFactory.hh"
43 #include "G4PVDivision.hh"
44 #include "G4PVReplica.hh"
45 #include "G4Region.hh"
46 #include "G4OpticalSurface.hh"
47 #include "G4LogicalSkinSurface.hh"
48 #include "G4LogicalBorderSurface.hh"
49 
50 #include "G4ProductionCuts.hh"
51 #include "G4ProductionCutsTable.hh"
52 #include "G4Gamma.hh"
53 #include "G4Electron.hh"
54 #include "G4Positron.hh"
55 #include "G4Proton.hh"
56 
57 #include "G4VSensitiveDetector.hh"
58 #include "G4AssemblyStore.hh"
59 #include "G4AssemblyVolume.hh"
60 
61 G4int G4GDMLWriteStructure::levelNo = 0; // Counter for level being exported
62 
63 G4GDMLWriteStructure::G4GDMLWriteStructure()
64  : G4GDMLWriteParamvol(), cexport(false), maxLevel(INT_MAX)
65 {
66  reflFactory = G4ReflectionFactory::Instance();
67 }
68 
69 G4GDMLWriteStructure::~G4GDMLWriteStructure()
70 {
71 }
72 
73 void
74 G4GDMLWriteStructure::DivisionvolWrite(xercesc::DOMElement* volumeElement,
75  const G4PVDivision* const divisionvol)
76 {
77  EAxis axis = kUndefined;
78  G4int number = 0;
79  G4double width = 0.0;
80  G4double offset = 0.0;
81  G4bool consuming = false;
82 
83  divisionvol->GetReplicationData(axis,number,width,offset,consuming);
84  axis = divisionvol->GetDivisionAxis();
85 
86  G4String unitString("mm");
87  G4String axisString("kUndefined");
88  if (axis==kXAxis) { axisString = "kXAxis"; }
89  else if (axis==kYAxis) { axisString = "kYAxis"; }
90  else if (axis==kZAxis) { axisString = "kZAxis"; }
91  else if (axis==kRho) { axisString = "kRho"; }
92  else if (axis==kPhi) { axisString = "kPhi"; unitString = "rad"; }
93 
94  const G4String name
95  = GenerateName(divisionvol->GetName(),divisionvol);
96  const G4String volumeref
97  = GenerateName(divisionvol->GetLogicalVolume()->GetName(),
98  divisionvol->GetLogicalVolume());
99 
100  xercesc::DOMElement* divisionvolElement = NewElement("divisionvol");
101  divisionvolElement->setAttributeNode(NewAttribute("axis",axisString));
102  divisionvolElement->setAttributeNode(NewAttribute("number",number));
103  divisionvolElement->setAttributeNode(NewAttribute("width",width));
104  divisionvolElement->setAttributeNode(NewAttribute("offset",offset));
105  divisionvolElement->setAttributeNode(NewAttribute("unit",unitString));
106  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
107  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
108  divisionvolElement->appendChild(volumerefElement);
109  volumeElement->appendChild(divisionvolElement);
110 }
111 
112 void G4GDMLWriteStructure::PhysvolWrite(xercesc::DOMElement* volumeElement,
113  const G4VPhysicalVolume* const physvol,
114  const G4Transform3D& T,
115  const G4String& ModuleName)
116 {
117  HepGeom::Scale3D scale;
118  HepGeom::Rotate3D rotate;
119  HepGeom::Translate3D translate;
120 
121  T.getDecomposition(scale,rotate,translate);
122 
123  const G4ThreeVector scl(scale(0,0),scale(1,1),scale(2,2));
124  const G4ThreeVector rot = GetAngles(rotate.getRotation());
125  const G4ThreeVector pos = T.getTranslation();
126 
127  const G4String name = GenerateName(physvol->GetName(),physvol);
128  const G4int copynumber = physvol->GetCopyNo();
129 
130  xercesc::DOMElement* physvolElement = NewElement("physvol");
131  physvolElement->setAttributeNode(NewAttribute("name",name));
132  if (copynumber) physvolElement->setAttributeNode(NewAttribute("copynumber", copynumber));
133 
134  volumeElement->appendChild(physvolElement);
135 
136  G4LogicalVolume* lv;
137  // Is it reflected?
138  if (reflFactory->IsReflected(physvol->GetLogicalVolume()))
139  {
140  lv = reflFactory->GetConstituentLV(physvol->GetLogicalVolume());
141  }
142  else
143  {
144  lv = physvol->GetLogicalVolume();
145  }
146 
147  const G4String volumeref = GenerateName(lv->GetName(), lv);
148 
149  if (ModuleName.empty())
150  {
151  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
152  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
153  physvolElement->appendChild(volumerefElement);
154  }
155  else
156  {
157  xercesc::DOMElement* fileElement = NewElement("file");
158  fileElement->setAttributeNode(NewAttribute("name",ModuleName));
159  fileElement->setAttributeNode(NewAttribute("volname",volumeref));
160  physvolElement->appendChild(fileElement);
161  }
162 
163  if (std::fabs(pos.x()) > kLinearPrecision
164  || std::fabs(pos.y()) > kLinearPrecision
165  || std::fabs(pos.z()) > kLinearPrecision)
166  {
167  PositionWrite(physvolElement,name+"_pos",pos);
168  }
169  if (std::fabs(rot.x()) > kAngularPrecision
170  || std::fabs(rot.y()) > kAngularPrecision
171  || std::fabs(rot.z()) > kAngularPrecision)
172  {
173  RotationWrite(physvolElement,name+"_rot",rot);
174  }
175  if (std::fabs(scl.x()-1.0) > kRelativePrecision
176  || std::fabs(scl.y()-1.0) > kRelativePrecision
177  || std::fabs(scl.z()-1.0) > kRelativePrecision)
178  {
179  ScaleWrite(physvolElement,name+"_scl",scl);
180  }
181 }
182 
183 void G4GDMLWriteStructure::ReplicavolWrite(xercesc::DOMElement* volumeElement,
184  const G4VPhysicalVolume* const replicavol)
185 {
186  EAxis axis = kUndefined;
187  G4int number = 0;
188  G4double width = 0.0;
189  G4double offset = 0.0;
190  G4bool consuming = false;
191  G4String unitString("mm");
192 
193  replicavol->GetReplicationData(axis,number,width,offset,consuming);
194 
195  const G4String volumeref
196  = GenerateName(replicavol->GetLogicalVolume()->GetName(),
197  replicavol->GetLogicalVolume());
198 
199  xercesc::DOMElement* replicavolElement = NewElement("replicavol");
200  replicavolElement->setAttributeNode(NewAttribute("number",number));
201  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
202  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
203  replicavolElement->appendChild(volumerefElement);
204  xercesc::DOMElement* replicateElement = NewElement("replicate_along_axis");
205  replicavolElement->appendChild(replicateElement);
206 
207  xercesc::DOMElement* dirElement = NewElement("direction");
208  if(axis==kXAxis)
209  { dirElement->setAttributeNode(NewAttribute("x","1")); }
210  else if(axis==kYAxis)
211  { dirElement->setAttributeNode(NewAttribute("y","1")); }
212  else if(axis==kZAxis)
213  { dirElement->setAttributeNode(NewAttribute("z","1")); }
214  else if(axis==kRho)
215  { dirElement->setAttributeNode(NewAttribute("rho","1")); }
216  else if(axis==kPhi)
217  { dirElement->setAttributeNode(NewAttribute("phi","1"));
218  unitString="rad"; }
219  replicateElement->appendChild(dirElement);
220 
221  xercesc::DOMElement* widthElement = NewElement("width");
222  widthElement->setAttributeNode(NewAttribute("value",width));
223  widthElement->setAttributeNode(NewAttribute("unit",unitString));
224  replicateElement->appendChild(widthElement);
225 
226  xercesc::DOMElement* offsetElement = NewElement("offset");
227  offsetElement->setAttributeNode(NewAttribute("value",offset));
228  offsetElement->setAttributeNode(NewAttribute("unit",unitString));
229  replicateElement->appendChild(offsetElement);
230 
231  volumeElement->appendChild(replicavolElement);
232 }
233 
234 
235 void G4GDMLWriteStructure::AssemblyWrite(xercesc::DOMElement* volumeElement,
236  const int assemblyID)
237 {
238 
239  G4AssemblyStore* assemblies = G4AssemblyStore::GetInstance();
240  G4AssemblyVolume* myassembly = assemblies->GetAssembly(assemblyID);
241 
242  xercesc::DOMElement* assemblyElement = NewElement("assembly");
243  G4String name = "Assembly_" + std::to_string(assemblyID);
244 
245  assemblyElement->setAttributeNode(NewAttribute("name",name));
246 
247  std::vector<G4AssemblyTriplet>::iterator vit = myassembly->GetTripletsIterator();
248 
249  int depth = 0;
250  const G4String ModuleName;
251 
252  for (size_t i5=0; i5<myassembly->TotalTriplets(); i5++)
253  {
254  TraverseVolumeTree((*vit).GetVolume(),depth+1);
255 
256  const G4ThreeVector rot = GetAngles((*vit).GetRotation()->inverse());
257  const G4ThreeVector pos = (*vit).GetTranslation();
258 
259  const G4String pname = GenerateName((*vit).GetVolume()->GetName()+"_pv", &(*vit));
260 
261  xercesc::DOMElement* physvolElement = NewElement("physvol");
262  physvolElement->setAttributeNode(NewAttribute("name",pname));
263 
264  assemblyElement->appendChild(physvolElement);
265 
266  const G4String volumeref = GenerateName((*vit).GetVolume()->GetName(), (*vit).GetVolume());
267 
268  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
269  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
270  physvolElement->appendChild(volumerefElement);
271 
272  if (std::fabs(pos.x()) > kLinearPrecision
273  || std::fabs(pos.y()) > kLinearPrecision
274  || std::fabs(pos.z()) > kLinearPrecision)
275  {
276  PositionWrite(physvolElement,"Position_" + std::to_string(i5), pos);
277  }
278  if (std::fabs(rot.x()) > kAngularPrecision
279  || std::fabs(rot.y()) > kAngularPrecision
280  || std::fabs(rot.z()) > kAngularPrecision)
281  {
282  RotationWrite(physvolElement,"Rotation_" + std::to_string(i5), rot);
283  }
284  vit++;
285  }
286 
287  volumeElement->appendChild(assemblyElement);
288 }
289 
290 
291 void G4GDMLWriteStructure::
292 BorderSurfaceCache(const G4LogicalBorderSurface* const bsurf)
293 {
294  if (!bsurf) { return; }
295 
296  const G4SurfaceProperty* psurf = bsurf->GetSurfaceProperty();
297 
298  // Generate the new element for border-surface
299  //
300  const G4String& bsname = GenerateName(bsurf->GetName(), bsurf);
301  const G4String& psname = GenerateName(psurf->GetName(), psurf);
302  xercesc::DOMElement* borderElement = NewElement("bordersurface");
303  borderElement->setAttributeNode(NewAttribute("name", bsname));
304  borderElement->setAttributeNode(NewAttribute("surfaceproperty", psname));
305 
306  const G4String volumeref1 = GenerateName(bsurf->GetVolume1()->GetName(),
307  bsurf->GetVolume1());
308  const G4String volumeref2 = GenerateName(bsurf->GetVolume2()->GetName(),
309  bsurf->GetVolume2());
310  xercesc::DOMElement* volumerefElement1 = NewElement("physvolref");
311  xercesc::DOMElement* volumerefElement2 = NewElement("physvolref");
312  volumerefElement1->setAttributeNode(NewAttribute("ref",volumeref1));
313  volumerefElement2->setAttributeNode(NewAttribute("ref",volumeref2));
314  borderElement->appendChild(volumerefElement1);
315  borderElement->appendChild(volumerefElement2);
316 
317  if (FindOpticalSurface(psurf))
318  {
319  const G4OpticalSurface* opsurf =
320  dynamic_cast<const G4OpticalSurface*>(psurf);
321  if (!opsurf)
322  {
323  G4Exception("G4GDMLWriteStructure::BorderSurfaceCache()",
324  "InvalidSetup", FatalException, "No optical surface found!");
325  return;
326  }
327  OpticalSurfaceWrite(solidsElement, opsurf);
328  }
329 
330  borderElementVec.push_back(borderElement);
331 }
332 
333 void G4GDMLWriteStructure::
334 SkinSurfaceCache(const G4LogicalSkinSurface* const ssurf)
335 {
336  if (!ssurf) { return; }
337 
338  const G4SurfaceProperty* psurf = ssurf->GetSurfaceProperty();
339 
340  // Generate the new element for border-surface
341  //
342  const G4String& ssname = GenerateName(ssurf->GetName(), ssurf);
343  const G4String& psname = GenerateName(psurf->GetName(), psurf);
344  xercesc::DOMElement* skinElement = NewElement("skinsurface");
345  skinElement->setAttributeNode(NewAttribute("name", ssname));
346  skinElement->setAttributeNode(NewAttribute("surfaceproperty", psname));
347 
348  const G4String volumeref = GenerateName(ssurf->GetLogicalVolume()->GetName(),
349  ssurf->GetLogicalVolume());
350  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
351  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
352  skinElement->appendChild(volumerefElement);
353 
354  if (FindOpticalSurface(psurf))
355  {
356  const G4OpticalSurface* opsurf =
357  dynamic_cast<const G4OpticalSurface*>(psurf);
358  if (!opsurf)
359  {
360  G4Exception("G4GDMLWriteStructure::SkinSurfaceCache()",
361  "InvalidSetup", FatalException, "No optical surface found!");
362  return;
363  }
364  OpticalSurfaceWrite(solidsElement, opsurf);
365  }
366 
367  skinElementVec.push_back(skinElement);
368 }
369 
370 G4bool G4GDMLWriteStructure::FindOpticalSurface(const G4SurfaceProperty* psurf)
371 {
372  const G4OpticalSurface* osurf = dynamic_cast<const G4OpticalSurface*>(psurf);
373  std::vector<const G4OpticalSurface*>::const_iterator pos;
374  pos = std::find(opt_vec.begin(), opt_vec.end(), osurf);
375  if (pos != opt_vec.end()) { return false; } // item already created!
376 
377  opt_vec.push_back(osurf); // cache it for future reference
378  return true;
379 }
380 
381 const G4LogicalSkinSurface*
382 G4GDMLWriteStructure::GetSkinSurface(const G4LogicalVolume* const lvol)
383 {
384  G4LogicalSkinSurface* surf = 0;
385  G4int nsurf = G4LogicalSkinSurface::GetNumberOfSkinSurfaces();
386  if (nsurf)
387  {
388  const G4LogicalSkinSurfaceTable* stable =
389  G4LogicalSkinSurface::GetSurfaceTable();
390  std::vector<G4LogicalSkinSurface*>::const_iterator pos;
391  for (pos = stable->begin(); pos != stable->end(); pos++)
392  {
393  if (lvol == (*pos)->GetLogicalVolume())
394  {
395  surf = *pos; break;
396  }
397  }
398  }
399  return surf;
400 }
401 
402 const G4LogicalBorderSurface*
403 G4GDMLWriteStructure::GetBorderSurface(const G4VPhysicalVolume* const pvol)
404 {
405  G4LogicalBorderSurface* surf = 0;
406  G4int nsurf = G4LogicalBorderSurface::GetNumberOfBorderSurfaces();
407  if (nsurf)
408  {
409  const G4LogicalBorderSurfaceTable* btable =
410  G4LogicalBorderSurface::GetSurfaceTable();
411  std::vector<G4LogicalBorderSurface*>::const_iterator pos;
412  for (pos = btable->begin(); pos != btable->end(); pos++)
413  {
414  if (pvol == (*pos)->GetVolume1()) // just the first in the couple
415  { // could be enough?
416  surf = *pos; // break;
417  BorderSurfaceCache(surf);
418  }
419  }
420  }
421  return surf;
422 }
423 
424 void G4GDMLWriteStructure::SurfacesWrite()
425 {
426 #ifdef G4VERBOSE
427  G4cout << "G4GDML: Writing surfaces..." << G4endl;
428 #endif
429  std::vector<xercesc::DOMElement*>::const_iterator pos;
430  for (pos = skinElementVec.begin(); pos != skinElementVec.end(); pos++)
431  {
432  structureElement->appendChild(*pos);
433  }
434  for (pos = borderElementVec.begin(); pos != borderElementVec.end(); pos++)
435  {
436  structureElement->appendChild(*pos);
437  }
438 }
439 
440 void G4GDMLWriteStructure::StructureWrite(xercesc::DOMElement* gdmlElement)
441 {
442 #ifdef G4VERBOSE
443  G4cout << "G4GDML: Writing structure..." << G4endl;
444 #endif
445 
446  // filling the list of phys volumes that are parts of assemblies
447 
448  G4AssemblyStore* assemblies = G4AssemblyStore::GetInstance();
449 
450  for(G4AssemblyStore::iterator it=assemblies->begin(); it!=assemblies->end(); it++)
451  {
452 
453  std::vector<G4VPhysicalVolume*>::iterator vit = (*it)->GetVolumesIterator();
454 
455  for (size_t i5=0; i5<(*it)->TotalImprintedVolumes(); i5++)
456  {
457  G4String pvname = (*vit)->GetName();
458  std::size_t pos = pvname.find("_impr_") + 6;
459  G4String impID = pvname.substr(pos);
460 
461  pos = impID.find("_");
462  impID = impID.substr(0, pos);
463 
464  assemblyVolMap[*vit] = (*it)->GetAssemblyID();
465  imprintsMap[*vit] = std::atoi(impID.c_str());
466  vit++;
467  }
468  }
469 
470  structureElement = NewElement("structure");
471  gdmlElement->appendChild(structureElement);
472 }
473 
474 G4Transform3D G4GDMLWriteStructure::
475 TraverseVolumeTree(const G4LogicalVolume* const volumePtr, const G4int depth)
476 {
477  if (VolumeMap().find(volumePtr) != VolumeMap().end())
478  {
479  return VolumeMap()[volumePtr]; // Volume is already processed
480  }
481 
482  G4VSolid* solidPtr = volumePtr->GetSolid();
483  G4Transform3D R,invR;
484  G4int trans=0;
485 
486  std::map<const G4LogicalVolume*, G4GDMLAuxListType>::iterator auxiter;
487 
488  levelNo++;
489 
490  while (true) // Solve possible displacement/reflection
491  { // of the referenced solid!
492  if (trans>maxTransforms)
493  {
494  G4String ErrorMessage = "Referenced solid in volume '"
495  + volumePtr->GetName()
496  + "' was displaced/reflected too many times!";
497  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
498  "InvalidSetup", FatalException, ErrorMessage);
499  }
500 
501  if (G4ReflectedSolid* refl = dynamic_cast<G4ReflectedSolid*>(solidPtr))
502  {
503  R = R*refl->GetTransform3D();
504  solidPtr = refl->GetConstituentMovedSolid();
505  trans++;
506  continue;
507  }
508 
509  if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(solidPtr))
510  {
511  R = R*G4Transform3D(disp->GetObjectRotation(),
512  disp->GetObjectTranslation());
513  solidPtr = disp->GetConstituentMovedSolid();
514  trans++;
515  continue;
516  }
517 
518  break;
519  }
520 
521  //check if it is a reflected volume
522  G4LogicalVolume* tmplv = const_cast<G4LogicalVolume*>(volumePtr);
523 
524  if (reflFactory->IsReflected(tmplv))
525  {
526  tmplv = reflFactory->GetConstituentLV(tmplv);
527  if (VolumeMap().find(tmplv) != VolumeMap().end())
528  {
529  return R; // Volume is already processed
530  }
531  }
532 
533  // Only compute the inverse when necessary!
534  //
535  if (trans>0) { invR = R.inverse(); }
536 
537  const G4String name
538  = GenerateName(tmplv->GetName(), tmplv);
539 
540  G4String materialref = "NULL";
541 
542  if(volumePtr->GetMaterial())
543  {
544  materialref = GenerateName(volumePtr->GetMaterial()->GetName(),
545  volumePtr->GetMaterial());
546  }
547 
548  const G4String solidref
549  = GenerateName(solidPtr->GetName(),solidPtr);
550 
551  xercesc::DOMElement* volumeElement = NewElement("volume");
552  volumeElement->setAttributeNode(NewAttribute("name",name));
553  xercesc::DOMElement* materialrefElement = NewElement("materialref");
554  materialrefElement->setAttributeNode(NewAttribute("ref",materialref));
555  volumeElement->appendChild(materialrefElement);
556  xercesc::DOMElement* solidrefElement = NewElement("solidref");
557  solidrefElement->setAttributeNode(NewAttribute("ref",solidref));
558  volumeElement->appendChild(solidrefElement);
559 
560  G4int daughterCount = volumePtr->GetNoDaughters();
561 
562  if (levelNo == maxLevel) // Stop exporting if reached levels limit
563  {
564  daughterCount = 0;
565  }
566 
567  std::vector<int> addedImprints;
568 
569  for (G4int i=0;i<daughterCount;i++) // Traverse all the children!
570  {
571  const G4VPhysicalVolume* const physvol = volumePtr->GetDaughter(i);
572  const G4String ModuleName = Modularize(physvol,depth);
573 
574  G4Transform3D daughterR;
575 
576  if (ModuleName.empty()) // Check if subtree requested to be
577  { // a separate module!
578  daughterR = TraverseVolumeTree(physvol->GetLogicalVolume(),depth+1);
579  }
580  else
581  {
582  G4GDMLWriteStructure writer;
583  daughterR = writer.Write(ModuleName,physvol->GetLogicalVolume(),
584  SchemaLocation,depth+1);
585  }
586 
587  if (const G4PVDivision* const divisionvol
588  = dynamic_cast<const G4PVDivision*>(physvol)) // Is it division?
589  {
590  if (!G4Transform3D::Identity.isNear(invR*daughterR,kRelativePrecision))
591  {
592  G4String ErrorMessage = "Division volume in '" + name
593  + "' can not be related to reflected solid!";
594  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
595  "InvalidSetup", FatalException, ErrorMessage);
596  }
597  DivisionvolWrite(volumeElement,divisionvol);
598  } else
599  if (physvol->IsParameterised()) // Is it a paramvol?
600  {
601  if (!G4Transform3D::Identity.isNear(invR*daughterR,kRelativePrecision))
602  {
603  G4String ErrorMessage = "Parameterised volume in '" + name
604  + "' can not be related to reflected solid!";
605  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
606  "InvalidSetup", FatalException, ErrorMessage);
607  }
608  ParamvolWrite(volumeElement,physvol);
609  } else
610  if (physvol->IsReplicated()) // Is it a replicavol?
611  {
612  if (!G4Transform3D::Identity.isNear(invR*daughterR,kRelativePrecision))
613  {
614  G4String ErrorMessage = "Replica volume in '" + name
615  + "' can not be related to reflected solid!";
616  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
617  "InvalidSetup", FatalException, ErrorMessage);
618  }
619  ReplicavolWrite(volumeElement,physvol);
620  }
621  else // Is it a physvol or an assembly?
622  {
623  if(assemblyVolMap.find(physvol) != assemblyVolMap.end())
624  {
625  int assemblyID = assemblyVolMap[physvol];
626 
627  G4String assemblyref = "Assembly_" + std::to_string(assemblyID);
628 
629  // here I need to retrieve the imprint ID
630 
631  G4int imprintID = imprintsMap[physvol];
632 
633  // there are 2 steps:
634  //
635  // 1) add assembly to the structure if that has not yet been done
636  // (but after the constituents volumes have been added)
637  //
638 
639  if(std::find(addedAssemblies.begin(), addedAssemblies.end(), assemblyID) == addedAssemblies.end())
640  {
641  AssemblyWrite(structureElement, assemblyID);
642  addedAssemblies.push_back(assemblyID);
643  }
644 
645  // 2) add the assembly (as physical volume) to the mother volume (but only once),
646  // using it's original position and rotation.
647  //
648 
649  // here I need a check if assembly has been already added to the mother volume
650  if(std::find(addedImprints.begin(), addedImprints.end(), imprintID) == addedImprints.end())
651  {
652  G4String imprintname = "Imprint_" + std::to_string(imprintID);
653  imprintname = GenerateName(imprintname, physvol);
654 
655  // I need to get those two from the container of imprints from the assembly
656  // I have the imprint ID, I need to get pos and rot
657  //
658  G4Transform3D& transf =
659  G4AssemblyStore::GetInstance()->GetAssembly(assemblyID)->GetImprintTransformation(imprintID);
660 
661  HepGeom::Scale3D scale;
662  HepGeom::Rotate3D rotate;
663  HepGeom::Translate3D translate;
664 
665  transf.getDecomposition(scale,rotate,translate);
666 
667  const G4ThreeVector scl(scale(0,0),scale(1,1),scale(2,2));
668  const G4ThreeVector rot = GetAngles(rotate.getRotation().inverse());
669  const G4ThreeVector pos = transf.getTranslation();
670 
671  // here I need a normal physvol referencing to my assemblyref
672 
673  xercesc::DOMElement* physvolElement = NewElement("physvol");
674  physvolElement->setAttributeNode(NewAttribute("name",imprintname));
675 
676  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
677  volumerefElement->setAttributeNode(NewAttribute("ref",assemblyref));
678  physvolElement->appendChild(volumerefElement);
679 
680  if (std::fabs(pos.x()) > kLinearPrecision
681  || std::fabs(pos.y()) > kLinearPrecision
682  || std::fabs(pos.z()) > kLinearPrecision)
683  {
684  PositionWrite(physvolElement,imprintname+"_pos",pos);
685  }
686  if (std::fabs(rot.x()) > kAngularPrecision
687  || std::fabs(rot.y()) > kAngularPrecision
688  || std::fabs(rot.z()) > kAngularPrecision)
689  {
690  RotationWrite(physvolElement,imprintname+"_rot",rot);
691  }
692  if (std::fabs(scl.x()-1.0) > kRelativePrecision
693  || std::fabs(scl.y()-1.0) > kRelativePrecision
694  || std::fabs(scl.z()-1.0) > kRelativePrecision)
695  {
696  ScaleWrite(physvolElement,name+"_scl",scl);
697  }
698 
699  volumeElement->appendChild(physvolElement);
700  //
701  addedImprints.push_back(imprintID);
702  }
703  }
704  else // not part of assembly, so a normal physical volume
705  {
706  G4RotationMatrix rot;
707  if (physvol->GetFrameRotation() != 0)
708  {
709  rot = *(physvol->GetFrameRotation());
710  }
711  G4Transform3D P(rot,physvol->GetObjectTranslation());
712 
713  PhysvolWrite(volumeElement,physvol,invR*P*daughterR,ModuleName);
714  }
715  }
716  // BorderSurfaceCache(GetBorderSurface(physvol));
717  GetBorderSurface(physvol);
718  }
719 
720  if (cexport) { ExportEnergyCuts(volumePtr); }
721  // Add optional energy cuts
722 
723  if (sdexport) { ExportSD(volumePtr); }
724  // Add optional SDs
725 
726  // Here write the auxiliary info
727  //
728  auxiter = auxmap.find(volumePtr);
729  if (auxiter != auxmap.end())
730  {
731  AddAuxInfo(&(auxiter->second), volumeElement);
732  }
733 
734  structureElement->appendChild(volumeElement);
735  // Append the volume AFTER traversing the children so that
736  // the order of volumes will be correct!
737 
738  VolumeMap()[tmplv] = R;
739 
740  AddExtension(volumeElement, volumePtr);
741  // Add any possible user defined extension attached to a volume
742 
743  AddMaterial(volumePtr->GetMaterial());
744  // Add the involved materials and solids!
745 
746  AddSolid(solidPtr);
747 
748  SkinSurfaceCache(GetSkinSurface(volumePtr));
749 
750  return R;
751 }
752 
753 void
754 G4GDMLWriteStructure::AddVolumeAuxiliary(G4GDMLAuxStructType myaux,
755  const G4LogicalVolume* const lvol)
756 {
757  std::map<const G4LogicalVolume*,
758  G4GDMLAuxListType>::iterator pos = auxmap.find(lvol);
759 
760  if (pos == auxmap.end()) { auxmap[lvol] = G4GDMLAuxListType(); }
761 
762  auxmap[lvol].push_back(myaux);
763 }
764 
765 void
766 G4GDMLWriteStructure::SetEnergyCutsExport(G4bool fcuts)
767 {
768  cexport = fcuts;
769 }
770 
771 void
772 G4GDMLWriteStructure::ExportEnergyCuts(const G4LogicalVolume* const lvol)
773 {
774  G4GDMLEvaluator eval;
775  G4ProductionCuts* pcuts = lvol->GetRegion()->GetProductionCuts();
776  G4ProductionCutsTable* ctab = G4ProductionCutsTable::GetProductionCutsTable();
777  G4Gamma* gamma = G4Gamma::Gamma();
778  G4Electron* eminus = G4Electron::Electron();
779  G4Positron* eplus = G4Positron::Positron();
780  G4Proton* proton = G4Proton::Proton();
781 
782  G4double gamma_cut = ctab->ConvertRangeToEnergy(gamma, lvol->GetMaterial(),
783  pcuts->GetProductionCut("gamma"));
784  G4double eminus_cut = ctab->ConvertRangeToEnergy(eminus, lvol->GetMaterial(),
785  pcuts->GetProductionCut("e-"));
786  G4double eplus_cut = ctab->ConvertRangeToEnergy(eplus, lvol->GetMaterial(),
787  pcuts->GetProductionCut("e+"));
788  G4double proton_cut = ctab->ConvertRangeToEnergy(proton, lvol->GetMaterial(),
789  pcuts->GetProductionCut("proton"));
790 
791  G4GDMLAuxStructType gammainfo = {"gammaECut",
792  eval.ConvertToString(gamma_cut), "MeV", 0};
793  G4GDMLAuxStructType eminusinfo = {"electronECut",
794  eval.ConvertToString(eminus_cut), "MeV", 0};
795  G4GDMLAuxStructType eplusinfo = {"positronECut",
796  eval.ConvertToString(eplus_cut), "MeV", 0};
797  G4GDMLAuxStructType protinfo = {"protonECut",
798  eval.ConvertToString(proton_cut), "MeV", 0};
799 
800  AddVolumeAuxiliary(gammainfo, lvol);
801  AddVolumeAuxiliary(eminusinfo, lvol);
802  AddVolumeAuxiliary(eplusinfo, lvol);
803  AddVolumeAuxiliary(protinfo, lvol);
804 }
805 
806 void
807 G4GDMLWriteStructure::SetSDExport(G4bool fsd)
808 {
809  sdexport = fsd;
810 }
811 
812 
813 void
814 G4GDMLWriteStructure::ExportSD(const G4LogicalVolume* const lvol)
815 {
816  G4VSensitiveDetector* sd = lvol->GetSensitiveDetector();
817 
818  if(sd)
819  {
820  G4String SDname = sd->GetName();
821 
822  G4GDMLAuxStructType SDinfo = {"SensDet", SDname, "", 0};
823  AddVolumeAuxiliary(SDinfo, lvol);
824  }
825 }
826 
827 G4int
828 G4GDMLWriteStructure::GetMaxExportLevel() const
829 {
830  return maxLevel;
831 }
832 
833 void
834 G4GDMLWriteStructure::SetMaxExportLevel(G4int level)
835 {
836  if (level <= 0)
837  {
838  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
839  "InvalidSetup", FatalException,
840  "Levels to export must be greater than zero!");
841  return;
842  }
843  maxLevel = level;
844  levelNo = 0;
845 }