DDG4/src/Geant4VolumeManager.cpp

0001 //==========================================================================
0002 //  AIDA Detector description implementation
0003 //--------------------------------------------------------------------------
0004 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
0005 // All rights reserved.
0006 //
0007 // For the licensing terms see $DD4hepINSTALL/LICENSE.
0008 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
0009 //
0010 // Author     : M.Frank
0011 //
0012 //==========================================================================
0013
0014 // Framework include files
0015 #include <DD4hep/Printout.h>
0016 #include <DD4hep/Volumes.h>
0017 #include <DD4hep/DetElement.h>
0018 #include <DD4hep/DetectorTools.h>
0019 #include <DDG4/Geant4VolumeManager.h>
0020 #include <DDG4/Geant4TouchableHandler.h>
0021 #include <DDG4/Geant4Mapping.h>
0022
0023 // Geant4 include files
0024 #include <G4VTouchable.hh>
0025 #include <G4LogicalVolume.hh>
0026 #include <G4VPhysicalVolume.hh>
0027
0028 // C/C++ include files
0029 #include <sstream>
0030
0031 using namespace dd4hep::sim;
0032 using namespace dd4hep;
0033
0034 #include <DDG4/Geant4AssemblyVolume.h>
0035 using VolIDDescriptor = std::pair<VolumeID,std::vector<std::pair<const BitFieldElement*, VolumeID> > >;
0036
0037 namespace  {
0038
0039   /// Helper class to populate the Geant4 volume manager
0040   struct Populator {
0041
0042     typedef std::vector<const TGeoNode*> Chain;
0043     typedef std::map<VolumeID,Geant4GeometryInfo::Geant4PlacementPath> Registries;
0044
0045     /// Reference to the Detector instance
0046     const Detector&     m_detDesc;
0047     /// Set of already added entries
0048     Registries          m_entries;
0049     /// Reference to Geant4 translation information
0050     Geant4GeometryInfo& m_geo;
0051
0052     /// Default constructor
0053     Populator(const Detector& description, Geant4GeometryInfo& g)
0054       : m_detDesc(description), m_geo(g) {
0055     }
0056
0057     /// Populate the Volume manager
0058     void populate(DetElement e) {
0059       const DetElement::Children& c = e.children();
0060       for (const auto& i : c)  {
0061         DetElement de = i.second;
0062         PlacedVolume pv = de.placement();
0063         if (pv.isValid()) {
0064           Chain chain;
0065           SensitiveDetector sd;
0066           PlacedVolume::VolIDs ids;
0067           m_entries.clear();
0068           chain.emplace_back(m_detDesc.world().placement().ptr());
0069           scanPhysicalVolume(pv.ptr(), std::move(ids), sd, chain);
0070           continue;
0071         }
0072         printout(WARNING, "Geant4VolumeManager",
0073                  "++ Detector element %s of type %s has no placement.", de.name(), de.type().c_str());
0074       }
0075       /// Needed to compute the cellID of parameterized volumes
0076       for( const auto& pv : m_geo.g4Placements )  {
0077         if ( pv.second->IsParameterised() )
0078           m_geo.g4Parameterised[pv.second] = pv.first;
0079         if ( pv.second->IsReplicated() )
0080           m_geo.g4Replicated[pv.second] = pv.first;
0081       }
0082     }
0083
0084     /// Scan a single physical volume and look for sensitive elements below
0085     void scanPhysicalVolume(const TGeoNode* node, PlacedVolume::VolIDs ids, SensitiveDetector& sd, Chain& chain) {
0086       PlacedVolume pv = node;
0087       Volume vol = pv.volume();
0088       PlacedVolume::VolIDs pv_ids = pv.volIDs();
0089
0090       chain.emplace_back(node);
0091       ids.PlacedVolume::VolIDs::Base::insert(ids.end(), pv_ids.begin(), pv_ids.end());
0092       if (vol.isSensitive()) {
0093         sd = vol.sensitiveDetector();
0094         if (sd.readout().isValid()) {
0095           add_entry(sd, node, ids, chain);
0096         }
0097         else {
0098           printout(WARNING, "Geant4VolumeManager",
0099                    "populate: Strange constellation volume %s is sensitive, but has no readout! sd:%p", pv.volume().name(),
0100                    sd.ptr());
0101         }
0102       }
0103       for (Int_t idau = 0, ndau = node->GetNdaughters(); idau < ndau; ++idau) {
0104         TGeoNode* daughter = node->GetDaughter(idau);
0105         PlacedVolume placement(daughter);
0106         if ( placement.data() ) {
0107           scanPhysicalVolume(daughter, ids, sd, chain);
0108         }
0109       }
0110       chain.pop_back();
0111     }
0112
0113     void add_entry(SensitiveDetector sd, const TGeoNode* n, const PlacedVolume::VolIDs& ids, const Chain& nodes) {
0114       Chain control;
0115       Volume       vol;
0116       const TGeoNode* node = nullptr;
0117       Readout      ro = sd.readout();
0118       IDDescriptor iddesc = ro.idSpec();
0119       VolumeID     code = iddesc.encode(ids);
0120       Geant4GeometryInfo::Geant4PlacementPath path;
0121       Registries::const_iterator i = m_entries.find(code);
0122       PrintLevel print_level  = m_geo.printLevel;
0123       PrintLevel print_action = print_level;
0124       PrintLevel print_chain  = print_level;
0125       PrintLevel print_res    = print_level;
0126
0127       printout(print_action,"Geant4VolumeManager","+++ Add path:%s vid:%016X",
0128                detail::tools::placementPath(nodes,false).c_str(),code);
0129
0130       if( i == m_entries.end() ) {
0131         path.reserve(nodes.size());
0132         for( Chain::const_reverse_iterator k = nodes.rbegin(), kend=nodes.rend(); k != kend; ++k ) {
0133           node = *(k);
0134           auto g4pit = m_geo.g4Placements.find(node);
0135           if( g4pit != m_geo.g4Placements.end() )  {
0136             G4VPhysicalVolume* phys = g4pit->second;
0137             if( phys->IsParameterised() )  {
0138               PlacedVolume pv(n);
0139               PlacedVolumeExtension* ext = pv.data();
0140               if( nullptr == ext->params->field )  {
0141                 ext->params->field = iddesc.field(ext->volIDs.at(0).first);
0142               }
0143             }
0144             path.emplace_back(phys);
0145             printout(print_chain, "Geant4VolumeManager", "+++     Chain: Node OK: %s [%s]",
0146                      node->GetName(), phys->GetName().c_str());
0147             continue;
0148           }
0149           control.insert(control.begin(),node);
0150           vol = Volume(node->GetVolume());
0151           auto iVolImp = m_geo.g4VolumeImprints.find(vol);
0152           if ( iVolImp != m_geo.g4VolumeImprints.end() )  {
0153             for(const auto& imp : iVolImp->second )  {
0154               const auto& c = imp.first;
0155               if ( c.size() <= control.size() && control == c )  {
0156                 path.emplace_back(imp.second);
0157                 printout(print_chain, "Geant4VolumeManager", "+++     Chain: Node OK: %s %s -> %s",
0158                          node->GetName(), detail::tools::placementPath(c,false).c_str(),
0159                          imp.second->GetName().c_str());
0160                 control.clear();
0161                 break;
0162               }
0163             }
0164           }
0165         }
0166         if ( control.empty() )  {
0167           printout(print_res, "Geant4VolumeManager", "+++     Volume  IDs:%s",
0168                    detail::tools::toString(ro.idSpec(),ids,code).c_str());
0169           path.erase(path.begin()+path.size()-1);
0170           printout(print_res, "Geant4VolumeManager", "+++     Map %016X to Geant4 Path:%s",
0171                    (void*)code, Geant4GeometryInfo::placementPath(path).c_str());
0172           if ( m_geo.g4Paths.find(path) == m_geo.g4Paths.end() ) {
0173             Geant4GeometryInfo::PlacementFlags opt;
0174             for(const auto* phys : path)  {
0175               opt.flags.path_has_parametrised = phys->IsParameterised() ? 1 : 0;
0176               opt.flags.path_has_replicated   = phys->IsReplicated()    ? 1 : 0;
0177             }
0178             opt.flags.parametrised = path.front()->IsParameterised() ? 1 : 0;
0179             opt.flags.replicated   = path.front()->IsReplicated()    ? 1 : 0;
0180             m_geo.g4Paths[path] = { code, opt.value };
0181             m_entries.emplace(code,path);
0182             return;
0183           }
0184           /// This is a normal case for parametrized volumes and no error
0185           if ( !path.empty() && (path.front()->IsParameterised() || path.front()->IsReplicated()) )  {
0186             return;
0187           }
0188           printout(ERROR, "Geant4VolumeManager", "populate: Severe error: Duplicated Geant4 path!!!! %s %s",
0189                    " [THIS SHOULD NEVER HAPPEN]",Geant4GeometryInfo::placementPath(path).c_str());
0190           goto Err;
0191         }
0192         printout(INFO, "Geant4VolumeManager", "Control block has still %d entries:%s",
0193                  int(control.size()),detail::tools::placementPath(control,true).c_str());
0194         goto Err;
0195       }
0196       else  {
0197         /// This is a normal case for parametrized volumes and no error
0198         if ( !path.empty() && (path.front()->IsParameterised() || path.front()->IsReplicated()) )  {
0199           return;
0200         }
0201       }
0202       printout(ERROR, "Geant4VolumeManager", "populate: Severe error: Duplicated Volume entry: 0x%X"
0203                " [THIS SHOULD NEVER HAPPEN]", code);
0204
0205     Err:
0206       if ( i != m_entries.end() )
0207         printout(ERROR,"Geant4VolumeManager"," Known G4 path: %s",Geant4GeometryInfo::placementPath((*i).second).c_str());
0208       if ( !path.empty() )
0209         printout(ERROR,"Geant4VolumeManager"," New   G4 path: %s",Geant4GeometryInfo::placementPath(path).c_str());
0210       if ( !nodes.empty() )
0211         printout(ERROR,"Geant4VolumeManager","     TGeo path: %s",detail::tools::placementPath(nodes,false).c_str());
0212       printout(ERROR,"Geant4VolumeManager",  " Offend.VolIDs: %s",detail::tools::toString(ro.idSpec(),ids,code).c_str());
0213       throw std::runtime_error("Failed to populate Geant4 volume manager!");
0214     }
0215   };
0216 }
0217
0218 /// Initializing constructor. The tree will automatically be built if possible
0219 Geant4VolumeManager::Geant4VolumeManager(const Detector& description, Geant4GeometryInfo* info)
0220   : Handle<Geant4GeometryInfo>(info)  {
0221   if (info && info->valid && info->g4Paths.empty()) {
0222     Populator p(description, *info);
0223     p.populate(description.world());
0224     return;
0225   }
0226   except("Geant4VolumeManager", "Attempt populate from invalid Geant4 geometry info [Invalid-Info]");
0227 }
0228
0229 /// Helper: Generate placement path from touchable object
0230 std::vector<const G4VPhysicalVolume*>
0231 Geant4VolumeManager::placementPath(const G4VTouchable* touchable, bool exception) const {
0232   Geant4TouchableHandler handler(touchable);
0233   return handler.placementPath(exception);
0234 }
0235
0236 /// Check the validity of the information before accessing it.
0237 bool Geant4VolumeManager::checkValidity() const {
0238   if (!isValid()) {
0239     except("Geant4VolumeManager", "Attempt to use invalid Geant4 volume manager [Invalid-Handle]");
0240   }
0241   else if (!ptr()->valid) {
0242     except("Geant4VolumeManager", "Attempt to use invalid Geant4 geometry info [Invalid-Info]");
0243   }
0244   return true;
0245 }
0246
0247 #if 0
0248 /// Access CELLID by placement path
0249 VolumeID Geant4VolumeManager::volumeID(const std::vector<const G4VPhysicalVolume*>& path) const {
0250   if (!path.empty() && checkValidity()) {
0251     const auto& mapping = ptr()->g4Paths;
0252     auto i = mapping.find(path);
0253     if ( i != mapping.end() )  {
0254       return (*i).second.first;
0255     }
0256     if (!path[0])
0257       return InvalidPath;
0258     else if (!path[0]->GetLogicalVolume()->GetSensitiveDetector())
0259       return Insensitive;
0260   }
0261   printout(INFO, "Geant4VolumeManager","+++   Bad volume Geant4 Path: %s",
0262            Geant4GeometryInfo::placementPath(path).c_str());
0263   return NonExisting;
0264 }
0265 #endif
0266
0267 /// Access CELLID by Geant4 touchable object
0268 VolumeID Geant4VolumeManager::volumeID(const G4VTouchable* touchable) const {
0269   Geant4TouchableHandler handler(touchable);
0270   std::vector<const G4VPhysicalVolume*> path = handler.placementPath();
0271   if( !path.empty() && checkValidity() )  {
0272     const auto& mapping = ptr()->g4Paths;
0273     auto i = mapping.find(path);
0274     if( i != mapping.end() )  {
0275       const auto& e = (*i).second;
0276       /// No parametrization or replication.
0277       if( e.flags == 0 )  {
0278         return e.volumeID;
0279       }
0280       VolumeID volid = e.volumeID;
0281       const auto& paramterised = ptr()->g4Parameterised;
0282       const auto& replicated   = ptr()->g4Replicated;
0283       /// This is incredibly slow .... but what can I do ? Need a better idea.
0284       for( std::size_t j=0; j < path.size(); ++j )  {
0285         const auto* phys = path[j];
0286         if( phys->IsParameterised() )  {
0287           int copy_no = touchable->GetCopyNumber(j);
0288           const auto it = paramterised.find(phys);
0289           if( it != paramterised.end() )  {
0290             //printout(INFO,"Geant4VolumeManager",
0291             //         "Copy number:   %ld  <--> %ld", copy_no, long(phys->GetCopyNo()));
0292             const auto* field = (*it).second.data()->params->field;
0293             volid |= IDDescriptor::encode(field, copy_no);
0294             continue;
0295           }
0296           except("Geant4VolumeManager","Error  Geant4VolumeManager::volumeID(const G4VTouchable* touchable)");
0297         }
0298         else if( phys->IsReplicated() )   {
0299           int copy_no = touchable->GetCopyNumber(j);
0300           const auto it = replicated.find(phys);
0301           if( it != replicated.end() )  {
0302             const auto* field = (*it).second.data()->params->field;
0303             volid |= IDDescriptor::encode(field, copy_no);
0304             continue;
0305           }
0306           except("Geant4VolumeManager","Error  Geant4VolumeManager::volumeID(const G4VTouchable* touchable)");
0307         }
0308       }
0309       return volid;
0310     }
0311     if( !path[0] )
0312       return InvalidPath;
0313     else if( !path[0]->GetLogicalVolume()->GetSensitiveDetector() )
0314       return Insensitive;
0315   }
0316   printout(INFO, "Geant4VolumeManager","+++   Bad volume Geant4 Path: %s",
0317            Geant4GeometryInfo::placementPath(path).c_str());
0318   return NonExisting;
0319 }
0320
0321 /// Accessfully decoded volume fields  by placement path
0322 void Geant4VolumeManager::volumeDescriptor(const std::vector<const G4VPhysicalVolume*>& path,
0323                                            VolIDDescriptor& vol_desc) const
0324 {
0325   vol_desc.second.clear();
0326   vol_desc.first = NonExisting;
0327   if( !path.empty() && checkValidity() )  {
0328     const auto& mapping = ptr()->g4Paths;
0329     auto i = mapping.find(path);
0330     if( i != mapping.end() ) {
0331       VolumeID vid = (*i).second.volumeID;
0332       G4LogicalVolume* lvol = path[0]->GetLogicalVolume();
0333       if( lvol->GetSensitiveDetector() ) {
0334         const auto* node = path[0];
0335         const auto& pm = ptr()->g4Placements;
0336         for( const auto& ipm : pm )  {
0337           if ( ipm.second == node )  {
0338             PlacedVolume      pv  = ipm.first;
0339             SensitiveDetector sd  = pv.volume().sensitiveDetector();
0340             IDDescriptor      dsc = sd.readout().idSpec();
0341             vol_desc.first = vid;
0342             dsc.decodeFields(vid, vol_desc.second);
0343             return;
0344           }
0345         }
0346       }
0347       vol_desc.first = Insensitive;
0348       return;
0349     }
0350     if( !path[0] )
0351       vol_desc.first = InvalidPath;
0352     else if( !path[0]->GetLogicalVolume()->GetSensitiveDetector() )
0353       vol_desc.first = Insensitive;
0354     else
0355       vol_desc.first = NonExisting;
0356   }
0357 }
0358
0359 /// Access fully decoded volume fields by Geant4 touchable object
0360 void Geant4VolumeManager::volumeDescriptor(const G4VTouchable* touchable,
0361                                            VolIDDescriptor&    vol_desc)  const  {
0362   volumeDescriptor(placementPath(touchable), vol_desc);
0363 }