g4simulation/g4gdml/PHG4GDMLWriteStructure.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 //
0027 // $Id: PHG4GDMLWriteStructure.cc 68053 2013-03-13 14:39:51Z gcosmo $
0028 //
0029 // class PHG4GDMLWriteStructure Implementation
0030 //
0031 // Original author: Zoltan Torzsok, November 2007
0032 //
0033 // --------------------------------------------------------------------
0034
0035 #include "PHG4GDMLWriteStructure.hh"
0036 #include "PHG4GDMLConfig.hh"
0037
0038 #include <Geant4/G4DisplacedSolid.hh>
0039 #include <Geant4/G4LogicalBorderSurface.hh>
0040 #include <Geant4/G4LogicalSkinSurface.hh>
0041 #include <Geant4/G4LogicalVolumeStore.hh>
0042 #include <Geant4/G4Material.hh>
0043 #include <Geant4/G4OpticalSurface.hh>
0044 #include <Geant4/G4PVDivision.hh>
0045 #include <Geant4/G4PVReplica.hh>
0046 #include <Geant4/G4PhysicalVolumeStore.hh>
0047 #include <Geant4/G4ReflectedSolid.hh>
0048 #include <Geant4/G4ReflectionFactory.hh>
0049 #include <Geant4/G4Region.hh>
0050
0051 #include <Geant4/G4Electron.hh>
0052 #include <Geant4/G4Gamma.hh>
0053 #include <Geant4/G4Positron.hh>
0054 #include <Geant4/G4ProductionCuts.hh>
0055 #include <Geant4/G4ProductionCutsTable.hh>
0056 #include <Geant4/G4Proton.hh>
0057 #include <Geant4/G4Version.hh>
0058
0059 #include <cassert>
0060
0061 PHG4GDMLWriteStructure::PHG4GDMLWriteStructure(const PHG4GDMLConfig* config_input)
0062   : structureElement(nullptr)
0063   , cexport(false)
0064   , config(config_input)
0065 {
0066   assert(config);
0067   reflFactory = G4ReflectionFactory::Instance();
0068 }
0069
0070 void PHG4GDMLWriteStructure::DivisionvolWrite(xercesc::DOMElement* volumeElement,
0071                                               const G4PVDivision* const divisionvol)
0072 {
0073   EAxis axis = kUndefined;
0074   G4int number = 0;
0075   G4double width = 0.0;
0076   G4double offset = 0.0;
0077   G4bool consuming = false;
0078
0079   divisionvol->GetReplicationData(axis, number, width, offset, consuming);
0080   axis = divisionvol->GetDivisionAxis();
0081
0082   G4String unitString("mm");
0083   G4String axisString("kUndefined");
0084   if (axis == kXAxis)
0085   {
0086     axisString = "kXAxis";
0087   }
0088   else if (axis == kYAxis)
0089   {
0090     axisString = "kYAxis";
0091   }
0092   else if (axis == kZAxis)
0093   {
0094     axisString = "kZAxis";
0095   }
0096   else if (axis == kRho)
0097   {
0098     axisString = "kRho";
0099   }
0100   else if (axis == kPhi)
0101   {
0102     axisString = "kPhi";
0103     unitString = "rad";
0104   }
0105
0106   const G4String name = GenerateName(divisionvol->GetName(), divisionvol);
0107   const G4String volumeref = GenerateName(divisionvol->GetLogicalVolume()->GetName(),
0108                                           divisionvol->GetLogicalVolume());
0109
0110   xercesc::DOMElement* divisionvolElement = NewElement("divisionvol");
0111   divisionvolElement->setAttributeNode(NewAttribute("axis", axisString));
0112   divisionvolElement->setAttributeNode(NewAttribute("number", number));
0113   divisionvolElement->setAttributeNode(NewAttribute("width", width));
0114   divisionvolElement->setAttributeNode(NewAttribute("offset", offset));
0115   divisionvolElement->setAttributeNode(NewAttribute("unit", unitString));
0116   xercesc::DOMElement* volumerefElement = NewElement("volumeref");
0117   volumerefElement->setAttributeNode(NewAttribute("ref", volumeref));
0118   divisionvolElement->appendChild(volumerefElement);
0119   volumeElement->appendChild(divisionvolElement);
0120 }
0121
0122 void PHG4GDMLWriteStructure::PhysvolWrite(xercesc::DOMElement* volumeElement,
0123                                           const G4VPhysicalVolume* const physvol,
0124                                           const G4Transform3D& T,
0125                                           const G4String& ModuleName)
0126 {
0127   HepGeom::Scale3D scale;
0128   HepGeom::Rotate3D rotate;
0129   HepGeom::Translate3D translate;
0130
0131   T.getDecomposition(scale, rotate, translate);
0132
0133   const G4ThreeVector scl(scale(0, 0), scale(1, 1), scale(2, 2));
0134   const G4ThreeVector rot = GetAngles(rotate.getRotation());
0135   const G4ThreeVector pos = T.getTranslation();
0136
0137   const G4String name = GenerateName(physvol->GetName(), physvol);
0138   const G4int copynumber = physvol->GetCopyNo();
0139
0140   xercesc::DOMElement* physvolElement = NewElement("physvol");
0141   physvolElement->setAttributeNode(NewAttribute("name", name));
0142   if (copynumber) { physvolElement->setAttributeNode(NewAttribute("copynumber", copynumber));
0143 }
0144
0145   volumeElement->appendChild(physvolElement);
0146
0147   G4LogicalVolume* lv;
0148   // Is it reflected?
0149   if (reflFactory->IsReflected(physvol->GetLogicalVolume()))
0150   {
0151     lv = reflFactory->GetConstituentLV(physvol->GetLogicalVolume());
0152   }
0153   else
0154   {
0155     lv = physvol->GetLogicalVolume();
0156   }
0157
0158   const G4String volumeref = GenerateName(lv->GetName(), lv);
0159
0160   if (ModuleName.empty())
0161   {
0162     xercesc::DOMElement* volumerefElement = NewElement("volumeref");
0163     volumerefElement->setAttributeNode(NewAttribute("ref", volumeref));
0164     physvolElement->appendChild(volumerefElement);
0165   }
0166   else
0167   {
0168     xercesc::DOMElement* fileElement = NewElement("file");
0169     fileElement->setAttributeNode(NewAttribute("name", ModuleName));
0170     fileElement->setAttributeNode(NewAttribute("volname", volumeref));
0171     physvolElement->appendChild(fileElement);
0172   }
0173
0174   if (std::fabs(pos.x()) > kLinearPrecision || std::fabs(pos.y()) > kLinearPrecision || std::fabs(pos.z()) > kLinearPrecision)
0175   {
0176     PositionWrite(physvolElement, name + "_pos", pos);
0177   }
0178   if (std::fabs(rot.x()) > kAngularPrecision || std::fabs(rot.y()) > kAngularPrecision || std::fabs(rot.z()) > kAngularPrecision)
0179   {
0180     RotationWrite(physvolElement, name + "_rot", rot);
0181   }
0182   if (std::fabs(scl.x() - 1.0) > kRelativePrecision || std::fabs(scl.y() - 1.0) > kRelativePrecision || std::fabs(scl.z() - 1.0) > kRelativePrecision)
0183   {
0184     ScaleWrite(physvolElement, name + "_scl", scl);
0185   }
0186 }
0187
0188 void PHG4GDMLWriteStructure::ReplicavolWrite(xercesc::DOMElement* volumeElement,
0189                                              const G4VPhysicalVolume* const replicavol)
0190 {
0191   EAxis axis = kUndefined;
0192   G4int number = 0;
0193   G4double width = 0.0;
0194   G4double offset = 0.0;
0195   G4bool consuming = false;
0196   G4String unitString("mm");
0197
0198   replicavol->GetReplicationData(axis, number, width, offset, consuming);
0199
0200   const G4String volumeref = GenerateName(replicavol->GetLogicalVolume()->GetName(),
0201                                           replicavol->GetLogicalVolume());
0202
0203   xercesc::DOMElement* replicavolElement = NewElement("replicavol");
0204   replicavolElement->setAttributeNode(NewAttribute("number", number));
0205   xercesc::DOMElement* volumerefElement = NewElement("volumeref");
0206   volumerefElement->setAttributeNode(NewAttribute("ref", volumeref));
0207   replicavolElement->appendChild(volumerefElement);
0208   xercesc::DOMElement* replicateElement = NewElement("replicate_along_axis");
0209   replicavolElement->appendChild(replicateElement);
0210
0211   xercesc::DOMElement* dirElement = NewElement("direction");
0212   if (axis == kXAxis)
0213   {
0214     dirElement->setAttributeNode(NewAttribute("x", "1"));
0215   }
0216   else if (axis == kYAxis)
0217   {
0218     dirElement->setAttributeNode(NewAttribute("y", "1"));
0219   }
0220   else if (axis == kZAxis)
0221   {
0222     dirElement->setAttributeNode(NewAttribute("z", "1"));
0223   }
0224   else if (axis == kRho)
0225   {
0226     dirElement->setAttributeNode(NewAttribute("rho", "1"));
0227   }
0228   else if (axis == kPhi)
0229   {
0230     dirElement->setAttributeNode(NewAttribute("phi", "1"));
0231     unitString = "rad";
0232   }
0233   replicateElement->appendChild(dirElement);
0234
0235   xercesc::DOMElement* widthElement = NewElement("width");
0236   widthElement->setAttributeNode(NewAttribute("value", width));
0237   widthElement->setAttributeNode(NewAttribute("unit", unitString));
0238   replicateElement->appendChild(widthElement);
0239
0240   xercesc::DOMElement* offsetElement = NewElement("offset");
0241   offsetElement->setAttributeNode(NewAttribute("value", offset));
0242   offsetElement->setAttributeNode(NewAttribute("unit", unitString));
0243   replicateElement->appendChild(offsetElement);
0244
0245   volumeElement->appendChild(replicavolElement);
0246 }
0247
0248 void PHG4GDMLWriteStructure::
0249     BorderSurfaceCache(const G4LogicalBorderSurface* const bsurf)
0250 {
0251   if (!bsurf)
0252   {
0253     return;
0254   }
0255
0256   const G4SurfaceProperty* psurf = bsurf->GetSurfaceProperty();
0257
0258   // Generate the new element for border-surface
0259   //
0260   xercesc::DOMElement* borderElement = NewElement("bordersurface");
0261   borderElement->setAttributeNode(NewAttribute("name", bsurf->GetName()));
0262   borderElement->setAttributeNode(NewAttribute("surfaceproperty",
0263                                                psurf->GetName()));
0264
0265   const G4String volumeref1 = GenerateName(bsurf->GetVolume1()->GetName(),
0266                                            bsurf->GetVolume1());
0267   const G4String volumeref2 = GenerateName(bsurf->GetVolume2()->GetName(),
0268                                            bsurf->GetVolume2());
0269   xercesc::DOMElement* volumerefElement1 = NewElement("physvolref");
0270   xercesc::DOMElement* volumerefElement2 = NewElement("physvolref");
0271   volumerefElement1->setAttributeNode(NewAttribute("ref", volumeref1));
0272   volumerefElement2->setAttributeNode(NewAttribute("ref", volumeref2));
0273   borderElement->appendChild(volumerefElement1);
0274   borderElement->appendChild(volumerefElement2);
0275
0276   if (FindOpticalSurface(psurf))
0277   {
0278     const G4OpticalSurface* opsurf =
0279         dynamic_cast<const G4OpticalSurface*>(psurf);
0280     if (!opsurf)
0281     {
0282       G4Exception("PHG4GDMLWriteStructure::BorderSurfaceCache()",
0283                   "InvalidSetup", FatalException, "No optical surface found!");
0284       return;
0285     }
0286     OpticalSurfaceWrite(solidsElement, opsurf);
0287   }
0288
0289   borderElementVec.push_back(borderElement);
0290 }
0291
0292 void PHG4GDMLWriteStructure::
0293     SkinSurfaceCache(const G4LogicalSkinSurface* const ssurf)
0294 {
0295   if (!ssurf)
0296   {
0297     return;
0298   }
0299
0300   const G4SurfaceProperty* psurf = ssurf->GetSurfaceProperty();
0301
0302   // Generate the new element for border-surface
0303   //
0304   xercesc::DOMElement* skinElement = NewElement("skinsurface");
0305   skinElement->setAttributeNode(NewAttribute("name", ssurf->GetName()));
0306   skinElement->setAttributeNode(NewAttribute("surfaceproperty",
0307                                              psurf->GetName()));
0308
0309   const G4String volumeref = GenerateName(ssurf->GetLogicalVolume()->GetName(),
0310                                           ssurf->GetLogicalVolume());
0311   xercesc::DOMElement* volumerefElement = NewElement("volumeref");
0312   volumerefElement->setAttributeNode(NewAttribute("ref", volumeref));
0313   skinElement->appendChild(volumerefElement);
0314
0315   if (FindOpticalSurface(psurf))
0316   {
0317     const G4OpticalSurface* opsurf =
0318         dynamic_cast<const G4OpticalSurface*>(psurf);
0319     if (!opsurf)
0320     {
0321       G4Exception("PHG4GDMLWriteStructure::SkinSurfaceCache()",
0322                   "InvalidSetup", FatalException, "No optical surface found!");
0323       return;
0324     }
0325     OpticalSurfaceWrite(solidsElement, opsurf);
0326   }
0327
0328   skinElementVec.push_back(skinElement);
0329 }
0330
0331 G4bool PHG4GDMLWriteStructure::FindOpticalSurface(const G4SurfaceProperty* psurf)
0332 {
0333   const G4OpticalSurface* osurf = dynamic_cast<const G4OpticalSurface*>(psurf);
0334   std::vector<const G4OpticalSurface*>::const_iterator pos;
0335   pos = std::find(opt_vec.begin(), opt_vec.end(), osurf);
0336   if (pos != opt_vec.end())
0337   {
0338     return false;
0339   }  // item already created!
0340
0341   opt_vec.push_back(osurf);  // cache it for future reference
0342   return true;
0343 }
0344
0345 const G4LogicalSkinSurface*
0346 PHG4GDMLWriteStructure::GetSkinSurface(const G4LogicalVolume* const lvol)
0347 {
0348   G4LogicalSkinSurface* surf = nullptr;
0349   G4int nsurf = G4LogicalSkinSurface::GetNumberOfSkinSurfaces();
0350   if (nsurf)
0351   {
0352     const G4LogicalSkinSurfaceTable* stable =
0353         G4LogicalSkinSurface::GetSurfaceTable();
0354     std::vector<G4LogicalSkinSurface*>::const_iterator pos;
0355     for (pos = stable->begin(); pos != stable->end(); ++pos)
0356     {
0357       if (lvol == (*pos)->GetLogicalVolume())
0358       {
0359         surf = *pos;
0360         break;
0361       }
0362     }
0363   }
0364   return surf;
0365 }
0366
0367 #if G4VERSION_NUMBER >= 1007
0368
0369 const G4LogicalBorderSurface* PHG4GDMLWriteStructure::GetBorderSurface(
0370     const G4VPhysicalVolume* const pvol)
0371 {
0372   G4LogicalBorderSurface* surf = nullptr;
0373   G4int nsurf = G4LogicalBorderSurface::GetNumberOfBorderSurfaces();
0374   if (nsurf)
0375   {
0376     const G4LogicalBorderSurfaceTable* btable =
0377         G4LogicalBorderSurface::GetSurfaceTable();
0378     for (const auto & pos : *btable)
0379     {
0380       if (pvol == pos.first.first)  // just the first in the couple
0381       {                              // could be enough?
0382         surf = pos.second;          // break;
0383         BorderSurfaceCache(surf);
0384       }
0385     }
0386   }
0387   return surf;
0388 }
0389
0390 #else
0391
0392 const G4LogicalBorderSurface*
0393 PHG4GDMLWriteStructure::GetBorderSurface(const G4VPhysicalVolume* const pvol)
0394 {
0395   G4LogicalBorderSurface* surf = 0;
0396   G4int nsurf = G4LogicalBorderSurface::GetNumberOfBorderSurfaces();
0397   if (nsurf)
0398   {
0399     const G4LogicalBorderSurfaceTable* btable =
0400         G4LogicalBorderSurface::GetSurfaceTable();
0401     std::vector<G4LogicalBorderSurface*>::const_iterator pos;
0402     for (pos = btable->begin(); pos != btable->end(); ++pos)
0403     {
0404       if (pvol == (*pos)->GetVolume1())  // just the first in the couple
0405       {                                  // is enough
0406         surf = *pos;
0407         break;
0408       }
0409     }
0410   }
0411   return surf;
0412 }
0413
0414 #endif
0415
0416 void PHG4GDMLWriteStructure::SurfacesWrite()
0417 {
0418   std::cout << "PHG4GDML: Writing surfaces..." << std::endl;
0419
0420   std::vector<xercesc::DOMElement*>::const_iterator pos;
0421   for (pos = skinElementVec.begin(); pos != skinElementVec.end(); ++pos)
0422   {
0423     structureElement->appendChild(*pos);
0424   }
0425   for (pos = borderElementVec.begin(); pos != borderElementVec.end(); ++pos)
0426   {
0427     structureElement->appendChild(*pos);
0428   }
0429 }
0430
0431 void PHG4GDMLWriteStructure::StructureWrite(xercesc::DOMElement* gdmlElement)
0432 {
0433   std::cout << "PHG4GDML: Writing structure..." << std::endl;
0434
0435   structureElement = NewElement("structure");
0436   gdmlElement->appendChild(structureElement);
0437 }
0438
0439 G4Transform3D PHG4GDMLWriteStructure::
0440     TraverseVolumeTree(const G4LogicalVolume* const volumePtr, const G4int depth)
0441 {
0442   if (VolumeMap().contains(volumePtr))
0443   {
0444     return VolumeMap()[volumePtr];  // Volume is already processed
0445   }
0446
0447   //jump over the exclusions
0448   assert(config);
0449   if (config->get_excluded_logical_vol().contains(volumePtr))
0450   {
0451     return G4Transform3D::Identity;
0452   }
0453
0454   G4VSolid* solidPtr = volumePtr->GetSolid();
0455   G4Transform3D R;
0456   G4Transform3D invR;
0457   G4int trans = 0;
0458
0459   std::map<const G4LogicalVolume*, PHG4GDMLAuxListType>::iterator auxiter;
0460
0461   while (true)  // Solve possible displacement/reflection
0462   {             // of the referenced solid!
0463     if (trans > maxTransforms)
0464     {
0465       G4String ErrorMessage = "Referenced solid in volume '" + volumePtr->GetName() + "' was displaced/reflected too many times!";
0466       G4Exception("PHG4GDMLWriteStructure::TraverseVolumeTree()",
0467                   "InvalidSetup", FatalException, ErrorMessage);
0468     }
0469
0470     if (G4ReflectedSolid* refl = dynamic_cast<G4ReflectedSolid*>(solidPtr))
0471     {
0472       R = R * refl->GetTransform3D();
0473       solidPtr = refl->GetConstituentMovedSolid();
0474       trans++;
0475       continue;
0476     }
0477
0478     if (G4DisplacedSolid* disp = dynamic_cast<G4DisplacedSolid*>(solidPtr))
0479     {
0480       R = R * G4Transform3D(disp->GetObjectRotation(),
0481                             disp->GetObjectTranslation());
0482       solidPtr = disp->GetConstituentMovedSolid();
0483       trans++;
0484       continue;
0485     }
0486
0487     break;
0488   }
0489
0490   //check if it is a reflected volume
0491   G4LogicalVolume* tmplv = const_cast<G4LogicalVolume*>(volumePtr);
0492
0493   if (reflFactory->IsReflected(tmplv))
0494   {
0495     tmplv = reflFactory->GetConstituentLV(tmplv);
0496     if (VolumeMap().contains(tmplv))
0497     {
0498       return R;  // Volume is already processed
0499     }
0500   }
0501
0502   // Only compute the inverse when necessary!
0503   //
0504   if (trans > 0)
0505   {
0506     invR = R.inverse();
0507   }
0508
0509   const G4String name = GenerateName(tmplv->GetName(), tmplv);
0510   const G4String materialref = GenerateName(volumePtr->GetMaterial()->GetName(),
0511                                             volumePtr->GetMaterial());
0512   const G4String solidref = GenerateName(solidPtr->GetName(), solidPtr);
0513
0514   xercesc::DOMElement* volumeElement = NewElement("volume");
0515   volumeElement->setAttributeNode(NewAttribute("name", name));
0516   xercesc::DOMElement* materialrefElement = NewElement("materialref");
0517   materialrefElement->setAttributeNode(NewAttribute("ref", materialref));
0518   volumeElement->appendChild(materialrefElement);
0519   xercesc::DOMElement* solidrefElement = NewElement("solidref");
0520   solidrefElement->setAttributeNode(NewAttribute("ref", solidref));
0521   volumeElement->appendChild(solidrefElement);
0522
0523   const G4int daughterCount = volumePtr->GetNoDaughters();
0524
0525   for (G4int i = 0; i < daughterCount; i++)  // Traverse all the children!
0526   {
0527     const G4VPhysicalVolume* const physvol = volumePtr->GetDaughter(i);
0528
0529     //jump over the exclusions
0530     assert(config);
0531     if (config->get_excluded_physical_vol().contains(physvol)) {
0532       continue;
0533 }
0534
0535     const G4String ModuleName = Modularize(physvol, depth);
0536
0537     G4Transform3D daughterR;
0538
0539     if (ModuleName.empty())  // Check if subtree requested to be
0540     {                        // a separate module!
0541       daughterR = TraverseVolumeTree(physvol->GetLogicalVolume(), depth + 1);
0542     }
0543     else
0544     {
0545       PHG4GDMLWriteStructure writer(config);
0546       daughterR = writer.Write(ModuleName, physvol->GetLogicalVolume(),
0547                                SchemaLocation, depth + 1);
0548     }
0549
0550     if (const G4PVDivision* const divisionvol = dynamic_cast<const G4PVDivision*>(physvol))  // Is it division?
0551     {
0552       if (!G4Transform3D::Identity.isNear(invR * daughterR, kRelativePrecision))
0553       {
0554         G4String ErrorMessage = "Division volume in '" + name + "' can not be related to reflected solid!";
0555         G4Exception("PHG4GDMLWriteStructure::TraverseVolumeTree()",
0556                     "InvalidSetup", FatalException, ErrorMessage);
0557       }
0558       DivisionvolWrite(volumeElement, divisionvol);
0559     }
0560     else if (physvol->IsParameterised())  // Is it a paramvol?
0561     {
0562       if (!G4Transform3D::Identity.isNear(invR * daughterR, kRelativePrecision))
0563       {
0564         G4String ErrorMessage = "Parameterised volume in '" + name + "' can not be related to reflected solid!";
0565         G4Exception("PHG4GDMLWriteStructure::TraverseVolumeTree()",
0566                     "InvalidSetup", FatalException, ErrorMessage);
0567       }
0568       ParamvolWrite(volumeElement, physvol);
0569     }
0570     else if (physvol->IsReplicated())  // Is it a replicavol?
0571     {
0572       if (!G4Transform3D::Identity.isNear(invR * daughterR, kRelativePrecision))
0573       {
0574         G4String ErrorMessage = "Replica volume in '" + name + "' can not be related to reflected solid!";
0575         G4Exception("PHG4GDMLWriteStructure::TraverseVolumeTree()",
0576                     "InvalidSetup", FatalException, ErrorMessage);
0577       }
0578       ReplicavolWrite(volumeElement, physvol);
0579     }
0580     else  // Is it a physvol?
0581     {
0582       G4RotationMatrix rot;
0583       if (physvol->GetFrameRotation() != nullptr)
0584       {
0585         rot = *(physvol->GetFrameRotation());
0586       }
0587       G4Transform3D P(rot, physvol->GetObjectTranslation());
0588
0589       PhysvolWrite(volumeElement, physvol, invR * P * daughterR, ModuleName);
0590     }
0591     BorderSurfaceCache(GetBorderSurface(physvol));
0592   }
0593
0594   if (cexport)
0595   {
0596     ExportEnergyCuts(volumePtr);
0597   }
0598   // Add optional energy cuts
0599
0600   // Here write the auxiliary info
0601   //
0602   auxiter = auxmap.find(volumePtr);
0603   if (auxiter != auxmap.end())
0604   {
0605     AddAuxInfo(&(auxiter->second), volumeElement);
0606   }
0607
0608   structureElement->appendChild(volumeElement);
0609   // Append the volume AFTER traversing the children so that
0610   // the order of volumes will be correct!
0611
0612   VolumeMap()[tmplv] = R;
0613
0614   AddExtension(volumeElement, volumePtr);
0615   // Add any possible user defined extension attached to a volume
0616
0617   AddMaterial(volumePtr->GetMaterial());
0618   // Add the involved materials and solids!
0619
0620   AddSolid(solidPtr);
0621
0622   SkinSurfaceCache(GetSkinSurface(volumePtr));
0623
0624   return R;
0625 }
0626
0627 void PHG4GDMLWriteStructure::AddVolumeAuxiliary(const PHG4GDMLAuxStructType& myaux,
0628                                                 const G4LogicalVolume* const lvol)
0629 {
0630   std::map<const G4LogicalVolume*,
0631            PHG4GDMLAuxListType>::iterator pos = auxmap.find(lvol);
0632
0633   if (pos == auxmap.end())
0634   {
0635     auxmap[lvol] = PHG4GDMLAuxListType();
0636   }
0637
0638   auxmap[lvol].push_back(myaux);
0639 }
0640
0641 void PHG4GDMLWriteStructure::SetEnergyCutsExport(G4bool fcuts)
0642 {
0643   cexport = fcuts;
0644 }
0645
0646 void PHG4GDMLWriteStructure::ExportEnergyCuts(const G4LogicalVolume* const lvol)
0647 {
0648   //  PHG4GDMLEvaluator eval;
0649   G4ProductionCuts* pcuts = lvol->GetRegion()->GetProductionCuts();
0650   G4ProductionCutsTable* ctab = G4ProductionCutsTable::GetProductionCutsTable();
0651   G4Gamma* gamma = G4Gamma::Gamma();
0652   G4Electron* eminus = G4Electron::Electron();
0653   G4Positron* eplus = G4Positron::Positron();
0654   G4Proton* proton = G4Proton::Proton();
0655
0656   G4double gamma_cut = ctab->ConvertRangeToEnergy(gamma, lvol->GetMaterial(),
0657                                                   pcuts->GetProductionCut("gamma"));
0658   G4double eminus_cut = ctab->ConvertRangeToEnergy(eminus, lvol->GetMaterial(),
0659                                                    pcuts->GetProductionCut("e-"));
0660   G4double eplus_cut = ctab->ConvertRangeToEnergy(eplus, lvol->GetMaterial(),
0661                                                   pcuts->GetProductionCut("e+"));
0662   G4double proton_cut = ctab->ConvertRangeToEnergy(proton, lvol->GetMaterial(),
0663                                                    pcuts->GetProductionCut("proton"));
0664
0665   PHG4GDMLAuxStructType gammainfo = {"gammaECut",
0666                                      ConvertToString(gamma_cut), "MeV", nullptr};
0667   PHG4GDMLAuxStructType eminusinfo = {"electronECut",
0668                                       ConvertToString(eminus_cut), "MeV", nullptr};
0669   PHG4GDMLAuxStructType eplusinfo = {"positronECut",
0670                                      ConvertToString(eplus_cut), "MeV", nullptr};
0671   PHG4GDMLAuxStructType protinfo = {"protonECut",
0672                                     ConvertToString(proton_cut), "MeV", nullptr};
0673
0674   AddVolumeAuxiliary(gammainfo, lvol);
0675   AddVolumeAuxiliary(eminusinfo, lvol);
0676   AddVolumeAuxiliary(eplusinfo, lvol);
0677   AddVolumeAuxiliary(protinfo, lvol);
0678 }
0679
0680 G4String PHG4GDMLWriteStructure::ConvertToString(G4double dval)
0681 {
0682   std::ostringstream os;
0683   os << dval;
0684   G4String vl = os.str();
0685   return vl;
0686 }