Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
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 // $Id: G4GDMLWriteStructure.cc 108895 2018-03-15 10:27:25Z gcosmo $
28 //
29 // class G4GDMLWriteStructure Implementation
30 //
31 // Original author: Zoltan Torzsok, November 2007
32 //
33 // --------------------------------------------------------------------
34 
35 #include "G4GDMLWriteStructure.hh"
36 #include "G4GDMLEvaluator.hh"
37 
38 #include "G4Material.hh"
39 #include "G4ReflectedSolid.hh"
40 #include "G4DisplacedSolid.hh"
41 #include "G4LogicalVolumeStore.hh"
42 #include "G4PhysicalVolumeStore.hh"
43 #include "G4ReflectionFactory.hh"
44 #include "G4PVDivision.hh"
45 #include "G4PVReplica.hh"
46 #include "G4Region.hh"
47 #include "G4OpticalSurface.hh"
48 #include "G4LogicalSkinSurface.hh"
50 
51 #include "G4ProductionCuts.hh"
52 #include "G4ProductionCutsTable.hh"
53 #include "G4Gamma.hh"
54 #include "G4Electron.hh"
55 #include "G4Positron.hh"
56 #include "G4Proton.hh"
57 
58 #include "G4VSensitiveDetector.hh"
59 
60 G4int G4GDMLWriteStructure::levelNo = 0; // Counter for level being exported
61 
63  : G4GDMLWriteParamvol(), cexport(false), maxLevel(INT_MAX)
64 {
66 }
67 
69 {
70 }
71 
72 void
73 G4GDMLWriteStructure::DivisionvolWrite(xercesc::DOMElement* volumeElement,
74  const G4PVDivision* const divisionvol)
75 {
76  EAxis axis = kUndefined;
77  G4int number = 0;
78  G4double width = 0.0;
79  G4double offset = 0.0;
80  G4bool consuming = false;
81 
82  divisionvol->GetReplicationData(axis,number,width,offset,consuming);
83  axis = divisionvol->GetDivisionAxis();
84 
85  G4String unitString("mm");
86  G4String axisString("kUndefined");
87  if (axis==kXAxis) { axisString = "kXAxis"; }
88  else if (axis==kYAxis) { axisString = "kYAxis"; }
89  else if (axis==kZAxis) { axisString = "kZAxis"; }
90  else if (axis==kRho) { axisString = "kRho"; }
91  else if (axis==kPhi) { axisString = "kPhi"; unitString = "rad"; }
92 
93  const G4String name
94  = GenerateName(divisionvol->GetName(),divisionvol);
95  const G4String volumeref
96  = GenerateName(divisionvol->GetLogicalVolume()->GetName(),
97  divisionvol->GetLogicalVolume());
98 
99  xercesc::DOMElement* divisionvolElement = NewElement("divisionvol");
100  divisionvolElement->setAttributeNode(NewAttribute("axis",axisString));
101  divisionvolElement->setAttributeNode(NewAttribute("number",number));
102  divisionvolElement->setAttributeNode(NewAttribute("width",width));
103  divisionvolElement->setAttributeNode(NewAttribute("offset",offset));
104  divisionvolElement->setAttributeNode(NewAttribute("unit",unitString));
105  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
106  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
107  divisionvolElement->appendChild(volumerefElement);
108  volumeElement->appendChild(divisionvolElement);
109 }
110 
111 void G4GDMLWriteStructure::PhysvolWrite(xercesc::DOMElement* volumeElement,
112  const G4VPhysicalVolume* const physvol,
113  const G4Transform3D& T,
114  const G4String& ModuleName)
115 {
117  HepGeom::Rotate3D rotate;
118  HepGeom::Translate3D translate;
119 
120  T.getDecomposition(scale,rotate,translate);
121 
122  const G4ThreeVector scl(scale(0,0),scale(1,1),scale(2,2));
123  const G4ThreeVector rot = GetAngles(rotate.getRotation());
124  const G4ThreeVector pos = T.getTranslation();
125 
126  const G4String name = GenerateName(physvol->GetName(),physvol);
127  const G4int copynumber = physvol->GetCopyNo();
128 
129  xercesc::DOMElement* physvolElement = NewElement("physvol");
130  physvolElement->setAttributeNode(NewAttribute("name",name));
131  if (copynumber) physvolElement->setAttributeNode(NewAttribute("copynumber", copynumber));
132 
133  volumeElement->appendChild(physvolElement);
134 
135  G4LogicalVolume* lv;
136  // Is it reflected?
137  if (reflFactory->IsReflected(physvol->GetLogicalVolume()))
138  {
140  }
141  else
142  {
143  lv = physvol->GetLogicalVolume();
144  }
145 
146  const G4String volumeref = GenerateName(lv->GetName(), lv);
147 
148  if (ModuleName.empty())
149  {
150  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
151  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
152  physvolElement->appendChild(volumerefElement);
153  }
154  else
155  {
156  xercesc::DOMElement* fileElement = NewElement("file");
157  fileElement->setAttributeNode(NewAttribute("name",ModuleName));
158  fileElement->setAttributeNode(NewAttribute("volname",volumeref));
159  physvolElement->appendChild(fileElement);
160  }
161 
162  if (std::fabs(pos.x()) > kLinearPrecision
163  || std::fabs(pos.y()) > kLinearPrecision
164  || std::fabs(pos.z()) > kLinearPrecision)
165  {
166  PositionWrite(physvolElement,name+"_pos",pos);
167  }
168  if (std::fabs(rot.x()) > kAngularPrecision
169  || std::fabs(rot.y()) > kAngularPrecision
170  || std::fabs(rot.z()) > kAngularPrecision)
171  {
172  RotationWrite(physvolElement,name+"_rot",rot);
173  }
174  if (std::fabs(scl.x()-1.0) > kRelativePrecision
175  || std::fabs(scl.y()-1.0) > kRelativePrecision
176  || std::fabs(scl.z()-1.0) > kRelativePrecision)
177  {
178  ScaleWrite(physvolElement,name+"_scl",scl);
179  }
180 }
181 
182 void G4GDMLWriteStructure::ReplicavolWrite(xercesc::DOMElement* volumeElement,
183  const G4VPhysicalVolume* const replicavol)
184 {
185  EAxis axis = kUndefined;
186  G4int number = 0;
187  G4double width = 0.0;
188  G4double offset = 0.0;
189  G4bool consuming = false;
190  G4String unitString("mm");
191 
192  replicavol->GetReplicationData(axis,number,width,offset,consuming);
193 
194  const G4String volumeref
195  = GenerateName(replicavol->GetLogicalVolume()->GetName(),
196  replicavol->GetLogicalVolume());
197 
198  xercesc::DOMElement* replicavolElement = NewElement("replicavol");
199  replicavolElement->setAttributeNode(NewAttribute("number",number));
200  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
201  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
202  replicavolElement->appendChild(volumerefElement);
203  xercesc::DOMElement* replicateElement = NewElement("replicate_along_axis");
204  replicavolElement->appendChild(replicateElement);
205 
206  xercesc::DOMElement* dirElement = NewElement("direction");
207  if(axis==kXAxis)
208  { dirElement->setAttributeNode(NewAttribute("x","1")); }
209  else if(axis==kYAxis)
210  { dirElement->setAttributeNode(NewAttribute("y","1")); }
211  else if(axis==kZAxis)
212  { dirElement->setAttributeNode(NewAttribute("z","1")); }
213  else if(axis==kRho)
214  { dirElement->setAttributeNode(NewAttribute("rho","1")); }
215  else if(axis==kPhi)
216  { dirElement->setAttributeNode(NewAttribute("phi","1"));
217  unitString="rad"; }
218  replicateElement->appendChild(dirElement);
219 
220  xercesc::DOMElement* widthElement = NewElement("width");
221  widthElement->setAttributeNode(NewAttribute("value",width));
222  widthElement->setAttributeNode(NewAttribute("unit",unitString));
223  replicateElement->appendChild(widthElement);
224 
225  xercesc::DOMElement* offsetElement = NewElement("offset");
226  offsetElement->setAttributeNode(NewAttribute("value",offset));
227  offsetElement->setAttributeNode(NewAttribute("unit",unitString));
228  replicateElement->appendChild(offsetElement);
229 
230  volumeElement->appendChild(replicavolElement);
231 }
232 
235 {
236  if (!bsurf) { return; }
237 
238  const G4SurfaceProperty* psurf = bsurf->GetSurfaceProperty();
239 
240  // Generate the new element for border-surface
241  //
242  xercesc::DOMElement* borderElement = NewElement("bordersurface");
243  borderElement->setAttributeNode(NewAttribute("name", bsurf->GetName()));
244  borderElement->setAttributeNode(NewAttribute("surfaceproperty",
245  psurf->GetName()));
246 
247  const G4String volumeref1 = GenerateName(bsurf->GetVolume1()->GetName(),
248  bsurf->GetVolume1());
249  const G4String volumeref2 = GenerateName(bsurf->GetVolume2()->GetName(),
250  bsurf->GetVolume2());
251  xercesc::DOMElement* volumerefElement1 = NewElement("physvolref");
252  xercesc::DOMElement* volumerefElement2 = NewElement("physvolref");
253  volumerefElement1->setAttributeNode(NewAttribute("ref",volumeref1));
254  volumerefElement2->setAttributeNode(NewAttribute("ref",volumeref2));
255  borderElement->appendChild(volumerefElement1);
256  borderElement->appendChild(volumerefElement2);
257 
258  if (FindOpticalSurface(psurf))
259  {
260  const G4OpticalSurface* opsurf =
261  dynamic_cast<const G4OpticalSurface*>(psurf);
262  if (!opsurf)
263  {
264  G4Exception("G4GDMLWriteStructure::BorderSurfaceCache()",
265  "InvalidSetup", FatalException, "No optical surface found!");
266  return;
267  }
269  }
270 
271  borderElementVec.push_back(borderElement);
272 }
273 
276 {
277  if (!ssurf) { return; }
278 
279  const G4SurfaceProperty* psurf = ssurf->GetSurfaceProperty();
280 
281  // Generate the new element for border-surface
282  //
283  xercesc::DOMElement* skinElement = NewElement("skinsurface");
284  skinElement->setAttributeNode(NewAttribute("name", ssurf->GetName()));
285  skinElement->setAttributeNode(NewAttribute("surfaceproperty",
286  psurf->GetName()));
287 
288  const G4String volumeref = GenerateName(ssurf->GetLogicalVolume()->GetName(),
289  ssurf->GetLogicalVolume());
290  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
291  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
292  skinElement->appendChild(volumerefElement);
293 
294  if (FindOpticalSurface(psurf))
295  {
296  const G4OpticalSurface* opsurf =
297  dynamic_cast<const G4OpticalSurface*>(psurf);
298  if (!opsurf)
299  {
300  G4Exception("G4GDMLWriteStructure::SkinSurfaceCache()",
301  "InvalidSetup", FatalException, "No optical surface found!");
302  return;
303  }
305  }
306 
307  skinElementVec.push_back(skinElement);
308 }
309 
311 {
312  const G4OpticalSurface* osurf = dynamic_cast<const G4OpticalSurface*>(psurf);
313  std::vector<const G4OpticalSurface*>::const_iterator pos;
314  pos = std::find(opt_vec.begin(), opt_vec.end(), osurf);
315  if (pos != opt_vec.end()) { return false; } // item already created!
316 
317  opt_vec.push_back(osurf); // cache it for future reference
318  return true;
319 }
320 
323 {
324  G4LogicalSkinSurface* surf = 0;
326  if (nsurf)
327  {
328  const G4LogicalSkinSurfaceTable* stable =
330  std::vector<G4LogicalSkinSurface*>::const_iterator pos;
331  for (pos = stable->begin(); pos != stable->end(); pos++)
332  {
333  if (lvol == (*pos)->GetLogicalVolume())
334  {
335  surf = *pos; break;
336  }
337  }
338  }
339  return surf;
340 }
341 
344 {
345  G4LogicalBorderSurface* surf = 0;
347  if (nsurf)
348  {
349  const G4LogicalBorderSurfaceTable* btable =
351  std::vector<G4LogicalBorderSurface*>::const_iterator pos;
352  for (pos = btable->begin(); pos != btable->end(); pos++)
353  {
354  if (pvol == (*pos)->GetVolume1()) // just the first in the couple
355  { // is enough
356  surf = *pos; break;
357  }
358  }
359  }
360  return surf;
361 }
362 
364 {
365 #ifdef G4VERBOSE
366  G4cout << "G4GDML: Writing surfaces..." << G4endl;
367 #endif
368  std::vector<xercesc::DOMElement*>::const_iterator pos;
369  for (pos = skinElementVec.begin(); pos != skinElementVec.end(); pos++)
370  {
371  structureElement->appendChild(*pos);
372  }
373  for (pos = borderElementVec.begin(); pos != borderElementVec.end(); pos++)
374  {
375  structureElement->appendChild(*pos);
376  }
377 }
378 
379 void G4GDMLWriteStructure::StructureWrite(xercesc::DOMElement* gdmlElement)
380 {
381 #ifdef G4VERBOSE
382  G4cout << "G4GDML: Writing structure..." << G4endl;
383 #endif
384  structureElement = NewElement("structure");
385  gdmlElement->appendChild(structureElement);
386 }
387 
389 TraverseVolumeTree(const G4LogicalVolume* const volumePtr, const G4int depth)
390 {
391  if (VolumeMap().find(volumePtr) != VolumeMap().end())
392  {
393  return VolumeMap()[volumePtr]; // Volume is already processed
394  }
395 
396  G4VSolid* solidPtr = volumePtr->GetSolid();
397  G4Transform3D R,invR;
398  G4int trans=0;
399 
400  std::map<const G4LogicalVolume*, G4GDMLAuxListType>::iterator auxiter;
401 
402  levelNo++;
403 
404  while (true) // Solve possible displacement/reflection
405  { // of the referenced solid!
406  if (trans>maxTransforms)
407  {
408  G4String ErrorMessage = "Referenced solid in volume '"
409  + volumePtr->GetName()
410  + "' was displaced/reflected too many times!";
411  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
412  "InvalidSetup", FatalException, ErrorMessage);
413  }
414 
415  if (G4ReflectedSolid* refl = dynamic_cast<G4ReflectedSolid*>(solidPtr))
416  {
417  R = R*refl->GetTransform3D();
418  solidPtr = refl->GetConstituentMovedSolid();
419  trans++;
420  continue;
421  }
422 
423  if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(solidPtr))
424  {
425  R = R*G4Transform3D(disp->GetObjectRotation(),
426  disp->GetObjectTranslation());
427  solidPtr = disp->GetConstituentMovedSolid();
428  trans++;
429  continue;
430  }
431 
432  break;
433  }
434 
435  //check if it is a reflected volume
436  G4LogicalVolume* tmplv = const_cast<G4LogicalVolume*>(volumePtr);
437 
438  if (reflFactory->IsReflected(tmplv))
439  {
440  tmplv = reflFactory->GetConstituentLV(tmplv);
441  if (VolumeMap().find(tmplv) != VolumeMap().end())
442  {
443  return R; // Volume is already processed
444  }
445  }
446 
447  // Only compute the inverse when necessary!
448  //
449  if (trans>0) { invR = R.inverse(); }
450 
451  const G4String name
452  = GenerateName(tmplv->GetName(), tmplv);
453 
454  G4String materialref = "NULL";
455 
456  if(volumePtr->GetMaterial())
457  {
458  materialref = GenerateName(volumePtr->GetMaterial()->GetName(),
459  volumePtr->GetMaterial());
460  }
461 
462  const G4String solidref
463  = GenerateName(solidPtr->GetName(),solidPtr);
464 
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);
473 
474  G4int daughterCount = volumePtr->GetNoDaughters();
475 
476  if (levelNo == maxLevel) // Stop exporting if reached levels limit
477  {
478  daughterCount = 0;
479  }
480 
481  for (G4int i=0;i<daughterCount;i++) // Traverse all the children!
482  {
483  const G4VPhysicalVolume* const physvol = volumePtr->GetDaughter(i);
484  const G4String ModuleName = Modularize(physvol,depth);
485 
486  G4Transform3D daughterR;
487 
488  if (ModuleName.empty()) // Check if subtree requested to be
489  { // a separate module!
490  daughterR = TraverseVolumeTree(physvol->GetLogicalVolume(),depth+1);
491  }
492  else
493  {
494  G4GDMLWriteStructure writer;
495  daughterR = writer.Write(ModuleName,physvol->GetLogicalVolume(),
496  SchemaLocation,depth+1);
497  }
498 
499  if (const G4PVDivision* const divisionvol
500  = dynamic_cast<const G4PVDivision*>(physvol)) // Is it division?
501  {
502  if (!G4Transform3D::Identity.isNear(invR*daughterR,kRelativePrecision))
503  {
504  G4String ErrorMessage = "Division volume in '" + name
505  + "' can not be related to reflected solid!";
506  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
507  "InvalidSetup", FatalException, ErrorMessage);
508  }
509  DivisionvolWrite(volumeElement,divisionvol);
510  } else
511  if (physvol->IsParameterised()) // Is it a paramvol?
512  {
513  if (!G4Transform3D::Identity.isNear(invR*daughterR,kRelativePrecision))
514  {
515  G4String ErrorMessage = "Parameterised volume in '" + name
516  + "' can not be related to reflected solid!";
517  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
518  "InvalidSetup", FatalException, ErrorMessage);
519  }
520  ParamvolWrite(volumeElement,physvol);
521  } else
522  if (physvol->IsReplicated()) // Is it a replicavol?
523  {
524  if (!G4Transform3D::Identity.isNear(invR*daughterR,kRelativePrecision))
525  {
526  G4String ErrorMessage = "Replica volume in '" + name
527  + "' can not be related to reflected solid!";
528  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
529  "InvalidSetup", FatalException, ErrorMessage);
530  }
531  ReplicavolWrite(volumeElement,physvol);
532  }
533  else // Is it a physvol?
534  {
535  G4RotationMatrix rot;
536  if (physvol->GetFrameRotation() != 0)
537  {
538  rot = *(physvol->GetFrameRotation());
539  }
540  G4Transform3D P(rot,physvol->GetObjectTranslation());
541 
542  PhysvolWrite(volumeElement,physvol,invR*P*daughterR,ModuleName);
543  }
545  }
546 
547  if (cexport) { ExportEnergyCuts(volumePtr); }
548  // Add optional energy cuts
549 
550  if (sdexport) { ExportSD(volumePtr); }
551  // Add optional SDs
552 
553  // Here write the auxiliary info
554  //
555  auxiter = auxmap.find(volumePtr);
556  if (auxiter != auxmap.end())
557  {
558  AddAuxInfo(&(auxiter->second), volumeElement);
559  }
560 
561  structureElement->appendChild(volumeElement);
562  // Append the volume AFTER traversing the children so that
563  // the order of volumes will be correct!
564 
565  VolumeMap()[tmplv] = R;
566 
567  AddExtension(volumeElement, volumePtr);
568  // Add any possible user defined extension attached to a volume
569 
570  AddMaterial(volumePtr->GetMaterial());
571  // Add the involved materials and solids!
572 
573  AddSolid(solidPtr);
574 
575  SkinSurfaceCache(GetSkinSurface(volumePtr));
576 
577  return R;
578 }
579 
580 void
582  const G4LogicalVolume* const lvol)
583 {
584  std::map<const G4LogicalVolume*,
585  G4GDMLAuxListType>::iterator pos = auxmap.find(lvol);
586 
587  if (pos == auxmap.end()) { auxmap[lvol] = G4GDMLAuxListType(); }
588 
589  auxmap[lvol].push_back(myaux);
590 }
591 
592 void
594 {
595  cexport = fcuts;
596 }
597 
598 void
600 {
601  G4GDMLEvaluator eval;
602  G4ProductionCuts* pcuts = lvol->GetRegion()->GetProductionCuts();
604  G4Gamma* gamma = G4Gamma::Gamma();
605  G4Electron* eminus = G4Electron::Electron();
608 
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"));
617 
618  G4GDMLAuxStructType gammainfo = {"gammaECut",
619  eval.ConvertToString(gamma_cut), "MeV", 0};
620  G4GDMLAuxStructType eminusinfo = {"electronECut",
621  eval.ConvertToString(eminus_cut), "MeV", 0};
622  G4GDMLAuxStructType eplusinfo = {"positronECut",
623  eval.ConvertToString(eplus_cut), "MeV", 0};
624  G4GDMLAuxStructType protinfo = {"protonECut",
625  eval.ConvertToString(proton_cut), "MeV", 0};
626 
627  AddVolumeAuxiliary(gammainfo, lvol);
628  AddVolumeAuxiliary(eminusinfo, lvol);
629  AddVolumeAuxiliary(eplusinfo, lvol);
630  AddVolumeAuxiliary(protinfo, lvol);
631 }
632 
633 void
635 {
636  sdexport = fsd;
637 }
638 
639 
640 void
642 {
644 
645  if(sd)
646  {
647  G4String SDname = sd->GetName();
648 
649  G4GDMLAuxStructType SDinfo = {"SensDet", SDname, "", 0};
650  AddVolumeAuxiliary(SDinfo, lvol);
651  }
652 }
653 
654 G4int
656 {
657  return maxLevel;
658 }
659 
660 void
662 {
663  if (level <= 0)
664  {
665  G4Exception("G4GDMLWriteStructure::TraverseVolumeTree()",
666  "InvalidSetup", FatalException,
667  "Levels to export must be greater than zero!");
668  return;
669  }
670  maxLevel = level;
671  levelNo = 0;
672 }
Definition: geomdefs.hh:54
const G4RotationMatrix * GetFrameRotation() const
void ExportSD(const G4LogicalVolume *const)
CLHEP::HepRotation getRotation() const
const XML_Char * name
Definition: expat.h:151
G4String Modularize(const G4VPhysicalVolume *const topvol, const G4int depth)
Definition: G4GDMLWrite.cc:332
CLHEP::Hep3Vector getTranslation() const
static const G4double pos
G4LogicalVolume * GetLogicalVolume() const
virtual G4bool IsParameterised() const =0
G4double ConvertRangeToEnergy(const G4ParticleDefinition *particle, const G4Material *material, G4double range)
G4VSensitiveDetector * GetSensitiveDetector() const
static const G4LogicalSkinSurfaceTable * GetSurfaceTable()
const G4VPhysicalVolume * GetVolume1() const
HepGeom::Transform3D G4Transform3D
const G4LogicalVolume * GetLogicalVolume() const
void ReplicavolWrite(xercesc::DOMElement *, const G4VPhysicalVolume *const)
G4String GenerateName(const G4String &, const void *const)
Definition: G4GDMLWrite.cc:124
void AddVolumeAuxiliary(G4GDMLAuxStructType myaux, const G4LogicalVolume *const)
#define G4endl
Definition: G4ios.hh:61
const G4LogicalBorderSurface * GetBorderSurface(const G4VPhysicalVolume *const)
Transform3D inverse() const
Definition: Transform3D.cc:142
static size_t GetNumberOfSkinSurfaces()
static const G4LogicalBorderSurfaceTable * GetSurfaceTable()
const G4VPhysicalVolume * GetVolume2() const
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:161
static const G4double kLinearPrecision
G4Transform3D TraverseVolumeTree(const G4LogicalVolume *const topVol, const G4int depth)
G4int GetNoDaughters() const
double z() const
std::vector< const G4OpticalSurface * > opt_vec
G4String ConvertToString(G4int ival)
static G4Proton * Proton()
Definition: G4Proton.cc:93
G4Material * GetMaterial() const
G4double GetProductionCut(G4int index) const
const G4String & GetName() const
virtual void GetReplicationData(EAxis &axis, G4int &nReplicas, G4double &width, G4double &offset, G4bool &consuming) const
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
void DivisionvolWrite(xercesc::DOMElement *, const G4PVDivision *const)
std::vector< xercesc::DOMElement * > borderElementVec
const G4String & GetName() const
Definition: G4Material.hh:179
G4ThreeVector GetAngles(const G4RotationMatrix &)
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
Double_t scale
const G4String & GetName() const
static const G4double kAngularPrecision
#define width
static const G4int maxTransforms
void AddAuxInfo(G4GDMLAuxListType *auxInfoList, xercesc::DOMElement *element)
Definition: G4GDMLWrite.cc:90
void PhysvolWrite(xercesc::DOMElement *, const G4VPhysicalVolume *const topVol, const G4Transform3D &transform, const G4String &moduleName)
virtual void AddSolid(const G4VSolid *const)
EAxis GetDivisionAxis() const
void AddMaterial(const G4Material *const)
Double_t R
static G4Positron * Positron()
Definition: G4Positron.cc:94
G4SurfaceProperty * GetSurfaceProperty() const
virtual G4int GetCopyNo() const =0
void BorderSurfaceCache(const G4LogicalBorderSurface *const)
virtual void ParamvolWrite(xercesc::DOMElement *, const G4VPhysicalVolume *const)
G4bool FindOpticalSurface(const G4SurfaceProperty *)
static G4Electron * Electron()
Definition: G4Electron.cc:94
std::vector< G4LogicalSkinSurface * > G4LogicalSkinSurfaceTable
std::vector< G4LogicalBorderSurface * > G4LogicalBorderSurfaceTable
virtual void StructureWrite(xercesc::DOMElement *)
#define INT_MAX
Definition: templates.hh:111
static DLL_API const Transform3D Identity
Definition: Transform3D.h:197
G4ProductionCuts * GetProductionCuts() const
G4int GetMaxExportLevel() const
static double P[]
static G4ProductionCutsTable * GetProductionCutsTable()
static constexpr double eplus
Definition: G4SIunits.hh:199
const G4LogicalSkinSurface * GetSkinSurface(const G4LogicalVolume *const)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
int G4int
Definition: G4Types.hh:78
VolumeMapType & VolumeMap()
Definition: G4GDMLWrite.cc:60
EAxis
Definition: geomdefs.hh:54
G4String GetName() const
void OpticalSurfaceWrite(xercesc::DOMElement *, const G4OpticalSurface *const)
std::map< const G4LogicalVolume *, G4GDMLAuxListType > auxmap
virtual void GetReplicationData(EAxis &axis, G4int &nReplicas, G4double &width, G4double &offset, G4bool &consuming) const =0
std::vector< G4GDMLAuxStructType > G4GDMLAuxListType
void SkinSurfaceCache(const G4LogicalSkinSurface *const)
static const G4double kRelativePrecision
virtual void AddExtension(xercesc::DOMElement *, const G4LogicalVolume *const)
Definition: G4GDMLWrite.cc:78
G4GLOB_DLL std::ostream G4cout
double x() const
xercesc::DOMElement * structureElement
virtual G4bool IsReplicated() const =0
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:137
G4VSolid * GetSolid() const
G4bool IsReflected(G4LogicalVolume *lv) const
G4String SchemaLocation
Definition: G4GDMLWrite.hh:131
void RotationWrite(xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &rot)
void getDecomposition(Scale3D &scale, Rotate3D &rotation, Translate3D &translation) const
Definition: Transform3D.cc:174
Definition: geomdefs.hh:54
double y() const
void ExportEnergyCuts(const G4LogicalVolume *const)
G4ReflectionFactory * reflFactory
G4ThreeVector GetObjectTranslation() const
std::vector< xercesc::DOMElement * > skinElementVec
G4Transform3D Write(const G4String &filename, const G4LogicalVolume *const topLog, const G4String &schemaPath, const G4int depth, G4bool storeReferences=true)
Definition: G4GDMLWrite.cc:167
xercesc::DOMElement * solidsElement
static size_t GetNumberOfBorderSurfaces()
G4Region * GetRegion() const
void ScaleWrite(xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &scl)
void PositionWrite(xercesc::DOMElement *element, const G4String &name, const G4ThreeVector &pos)
const G4String & GetName() const
G4VPhysicalVolume * GetDaughter(const G4int i) const
G4LogicalVolume * GetConstituentLV(G4LogicalVolume *reflLV) const
const G4String & GetName() const
static G4ReflectionFactory * Instance()