Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-05 08:18:03 

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

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!