ISFG4GeoHelper.cxx

Go to the documentation of this file.

/*
  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
*/
 
#include "ISF_Geant4Event/ISFG4GeoHelper.h"
 
// Athena includes
 
//Geant4
#include "G4LogicalVolumeStore.hh"
#include "G4Step.hh"
#include "globals.hh"
#include "ISF_Interfaces/IGeoIDSvc.h"
#include "G4TransportationManager.hh"
#include "SimHelpers/StepHelper.h"
 
 
//________________________________________________________________________
AtlasDetDescr::AtlasRegion
iGeant4::ISFG4GeoHelper::nextGeoId(const G4Step* aStep, int truthVolLevel, ISF::IGeoIDSvc *geoIDSvc)
{
  struct VolumeHolder {
    G4LogicalVolume *BP, *ID, *CALO, *MU, *TTR;
  };
 
  // one time initialization of volumes
  static const VolumeHolder vh = [&]() {
    VolumeHolder h{};
    G4LogicalVolumeStore * lvs = G4LogicalVolumeStore::GetInstance();
    for (size_t i=0;i<lvs->size();++i) {
 
      if ( !(*lvs)[i] ) { continue; }
      if ( (*lvs)[i]->GetName() == "BeamPipe::BeamPipe" ) { h.BP = (*lvs)[i]; }
      else if ( (*lvs)[i]->GetName() == "IDET::IDET" || (*lvs)[i]->GetName() == "ITK::ITK" ) { h.ID = (*lvs)[i]; }
      else if ( (*lvs)[i]->GetName() == "CALO::CALO" ) { h.CALO = (*lvs)[i]; }
      else if ( (*lvs)[i]->GetName() == "MUONQ02::MUONQ02" ) { h.MU = (*lvs)[i]; }
      else if ( (*lvs)[i]->GetName() == "TTR_BARREL::TTR_BARREL" ) { h.TTR = (*lvs)[i]; }
    }
 
    const auto& worldVolume = G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking()->GetWorldVolume()->GetLogicalVolume();
    ISFG4GeoHelper::checkVolumeDepth(worldVolume, truthVolLevel);
    return h;
  }();
 
  AtlasDetDescr::AtlasRegion nextGeoID = truthVolLevel > 1 ? AtlasDetDescr::fAtlasCavern
                                                           : AtlasDetDescr::fUndefinedAtlasRegion;
 
  const G4StepPoint *postStep = aStep->GetPostStepPoint();
  bool leavingG4World       = postStep->GetStepStatus()==fWorldBoundary;
 
  if ( leavingG4World ) {
    nextGeoID = AtlasDetDescr::fAtlasCavern;
    return nextGeoID;
  }
 
  // If in mother volume, use the ISF_GeoIDSvc to resolve the geoID
  if (G4StepHelper::postStepBranchDepth(aStep)<truthVolLevel) {
    nextGeoID = getNextGeoIDFromSvc(*postStep, *geoIDSvc);
    return nextGeoID;
  }
 
  // Ordering inside out (most truth in the ID anyway...)
  if ( vh.ID==G4StepHelper::getPostStepLogicalVolume(aStep,truthVolLevel) ) {
    nextGeoID = AtlasDetDescr::fAtlasID;
  }
  else if ( vh.CALO==G4StepHelper::getPostStepLogicalVolume(aStep,truthVolLevel) ) {
    nextGeoID = AtlasDetDescr::fAtlasCalo;
  }
  else if ( vh.MU==G4StepHelper::getPostStepLogicalVolume(aStep,truthVolLevel) ) {
    nextGeoID = AtlasDetDescr::fAtlasMS;
  }
  else if ( vh.BP==G4StepHelper::getPostStepLogicalVolume(aStep,truthVolLevel) ) {
    nextGeoID = (G4StepHelper::postStepBranchDepth(aStep)>truthVolLevel && G4StepHelper::getPostStepLogicalVolumeName(aStep,truthVolLevel+1)=="BeamPipe::BeamPipeCentral")?AtlasDetDescr::fAtlasID:AtlasDetDescr::fAtlasForward;
  }
  else if ( vh.TTR==G4StepHelper::getPostStepLogicalVolume(aStep,truthVolLevel) ) {
    nextGeoID = AtlasDetDescr::fAtlasCavern;
  }
  else if (truthVolLevel>0 && G4StepHelper::getPostStepLogicalVolumeName(aStep,truthVolLevel-1).find("CavernInfra")!=std::string::npos) {
    nextGeoID = AtlasDetDescr::fAtlasCavern;
  }
  else {
    nextGeoID = getNextGeoIDFromSvc(*postStep, *geoIDSvc);
  }
 
  return nextGeoID;
}
 
AtlasDetDescr::AtlasRegion iGeant4::ISFG4GeoHelper::getNextGeoIDFromSvc(const G4StepPoint& postStep,
                                                                        const ISF::IGeoIDSvc &geoIDSvc) {
  const G4ThreeVector &postPos = postStep.GetPosition();
  AtlasDetDescr::AtlasRegion nextGeoID = geoIDSvc.identifyGeoID(postPos.x(),
                                                                postPos.y(),
                                                                postPos.z());
 
  // if we ever run into problems with the current approach, the following
  // takes the particle's traveling direction into account for finding the
  // *next* volume it enters
  //const G4ThreeVector     &postMom  = postStep->GetMomentum();
  //nextGeoID = m_geoIDSvcQuick->identifyNextGeoID( postPos.x(),
  //                                                postPos.y(),
  //                                                postPos.z(),
  //                                                postMom.x(),
  //                                                postMom.y(),
  //                                                postMom.z() );
 
  return nextGeoID;
}
 
bool iGeant4::ISFG4GeoHelper::checkVolumeDepth(G4LogicalVolume* lv, int volLevel, int depth) {
 
  //FIXME - can replace all this code with similar methods to those in MCTruthBase/src/RecordingEnvelope.cxx
  if (lv==nullptr) { return false; }
 
  bool Cavern = false;
  // Check the volumes rather explicitly
 
  if ( lv->GetName() == "BeamPipe::BeamPipe" ||
       lv->GetName() == "IDET::IDET" ||
       lv->GetName() == "ITK::ITK" ||
       lv->GetName() == "CALO::CALO" ||
       lv->GetName() == "MUONQ02::MUONQ02" ||
       lv->GetName() == "TTR_BARREL::TTR_BARREL" ) {
    if(depth!=volLevel) {
      G4ExceptionDescription description;
      description << "Volume " << lv->GetName() << " at depth " << depth << " instead of depth " << volLevel;
      G4Exception("iGeant4::ISFG4GeoHelper", "checkVolumeDepth", FatalException, description);
      return false; //The G4Exception call above should abort the job, but Coverity does not seem to pick this up.
    }
    //ATH_MSG_DEBUG("Volume " << lv->GetName() << " is correctly registered at depth " << depth);
  }
  else if ( lv->GetName()=="BeamPipe::BeamPipeCentral" ) { // Things that are supposed to be one deeper
    if (depth!=volLevel+1) {
      G4ExceptionDescription description;
      description << "Volume " << lv->GetName() << " at depth " << depth << " instead of depth " << volLevel+1;
      G4Exception("iGeant4::ISFG4GeoHelper", "checkVolumeDepth", FatalException, description);
      return false; //The G4Exception call above should abort the job, but Coverity does not seem to pick this up.
    }
    //ATH_MSG_DEBUG("Volume " << lv->GetName() << " is correctly registered at depth " << depth);
  }
  else if ( lv->GetName().find("CavernInfra")!=std::string::npos ) { // Things that are supposed to be one shallower
    if (depth==volLevel-1) {
      Cavern=true;
      //ATH_MSG_DEBUG("Volume " << lv->GetName() << " is correctly registered at depth " << depth);
      // Note: a number of volumes exist with "CavernInfra" in the name at the wrong depth, so we just need to
      //   check that there's at least one at the right depth
    } // Check of volume level
  } // Check of volume name
 
  // Going through the volume depth
  for (unsigned int i=0; i<lv->GetNoDaughters(); ++i) {
    Cavern = Cavern || checkVolumeDepth( lv->GetDaughter(i)->GetLogicalVolume() , volLevel , depth+1 );
  }
  if (depth==0 && !Cavern && volLevel>1) {
    G4ExceptionDescription description;
    description << "No CavernInfra volume registered at depth " << volLevel-1;
    G4Exception("iGeant4::ISFG4GeoHelper", "checkVolumeDepth", FatalException, description);
    return false; //The G4Exception call above should abort the job, but Coverity does not seem to pick this up.
  }
 
  return Cavern;
}