#include <GFlashHadronShowerModel.h>

Inheritance diagram for GFlashHadronShowerModel:

Public Member Functions
void	DoIt (const G4FastTrack &, G4FastStep &) override

	GFlashHadronShowerModel (G4String modelName, G4Region *envelope, const edm::ParameterSet &parSet)

G4bool	IsApplicable (const G4ParticleDefinition &) override

G4bool	ModelTrigger (const G4FastTrack &) override

	~GFlashHadronShowerModel () override

Private Member Functions
G4bool	excludeDetectorRegion (const G4FastTrack &fastTrack)

G4bool	isFirstInelasticInteraction (const G4FastTrack &fastTrack)

void	makeHits (const G4FastTrack &fastTrack)

void	updateGflashStep (const G4ThreeVector &position, G4double time)

Private Attributes
GflashAntiProtonShowerProfile *	theAntiProtonProfile

G4Navigator *	theGflashNavigator

G4Step *	theGflashStep

G4TouchableHandle	theGflashTouchableHandle

GflashKaonMinusShowerProfile *	theKaonMinusProfile

GflashKaonPlusShowerProfile *	theKaonPlusProfile

edm::ParameterSet	theParSet

GflashPiKShowerProfile *	thePiKProfile

GflashHadronShowerProfile *	theProfile

GflashProtonShowerProfile *	theProtonProfile

const G4Region *	theRegion

G4bool	theWatcherOn

Detailed Description

Definition at line 20 of file GFlashHadronShowerModel.h.

Constructor & Destructor Documentation

◆ GFlashHadronShowerModel()

GFlashHadronShowerModel::GFlashHadronShowerModel	(	G4String	modelName,
		G4Region *	envelope,
		const edm::ParameterSet &	parSet
	)

Definition at line 35 of file GFlashHadronShowerModel.cc.

     : G4VFastSimulationModel(modelName, envelope), theParSet(parSet) {
   theWatcherOn = parSet.getParameter<bool>("watcherOn");
  
   theProfile = nullptr;
   thePiKProfile = new GflashPiKShowerProfile(parSet);
   theKaonPlusProfile = new GflashKaonPlusShowerProfile(parSet);
   theKaonMinusProfile = new GflashKaonMinusShowerProfile(parSet);
   theProtonProfile = new GflashProtonShowerProfile(parSet);
   theAntiProtonProfile = new GflashAntiProtonShowerProfile(parSet);
   //theHisto = GflashHistogram::instance();
  
   theRegion = const_cast<const G4Region*>(envelope);
  
   theGflashStep = new G4Step();
   theGflashTouchableHandle = new G4TouchableHistory();
   theGflashNavigator = new G4Navigator();
 }

References edm::ParameterSet::getParameter(), theAntiProtonProfile, theGflashNavigator, theGflashStep, theGflashTouchableHandle, theKaonMinusProfile, theKaonPlusProfile, thePiKProfile, theProfile, theProtonProfile, theRegion, and theWatcherOn.

◆ ~GFlashHadronShowerModel()

GFlashHadronShowerModel::~GFlashHadronShowerModel ( )

override

Definition at line 56 of file GFlashHadronShowerModel.cc.

                                                   {
   delete thePiKProfile;
   delete theKaonPlusProfile;
   delete theKaonMinusProfile;
   delete theProtonProfile;
   delete theAntiProtonProfile;
   delete theGflashStep;
 }

References theAntiProtonProfile, theGflashStep, theKaonMinusProfile, theKaonPlusProfile, thePiKProfile, and theProtonProfile.

Member Function Documentation

◆ DoIt()

void GFlashHadronShowerModel::DoIt	(	const G4FastTrack &	fastTrack,
		G4FastStep &	fastStep
	)

override

Definition at line 110 of file GFlashHadronShowerModel.cc.

                                                                                      {
   // kill the particle
   fastStep.KillPrimaryTrack();
   fastStep.ProposePrimaryTrackPathLength(0.0);
  
   // parameterize energy depostion by the particle type
   G4ParticleDefinition* particleType = fastTrack.GetPrimaryTrack()->GetDefinition();
  
   theProfile = thePiKProfile;
   if (particleType == G4KaonMinus::KaonMinusDefinition())
     theProfile = theKaonMinusProfile;
   else if (particleType == G4KaonPlus::KaonPlusDefinition())
     theProfile = theKaonPlusProfile;
   else if (particleType == G4AntiProton::AntiProtonDefinition())
     theProfile = theAntiProtonProfile;
   else if (particleType == G4Proton::ProtonDefinition())
     theProfile = theProtonProfile;
  
   //input variables for GflashHadronShowerProfile
   G4double energy = fastTrack.GetPrimaryTrack()->GetKineticEnergy() / GeV;
   G4double globalTime = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetGlobalTime();
   G4double charge = fastTrack.GetPrimaryTrack()->GetStep()->GetPreStepPoint()->GetCharge();
   G4ThreeVector position = fastTrack.GetPrimaryTrack()->GetPosition() / cm;
   G4ThreeVector momentum = fastTrack.GetPrimaryTrack()->GetMomentum() / GeV;
   G4int showerType = Gflash::findShowerType(position);
  
   theProfile->initialize(showerType, energy, globalTime, charge, position, momentum);
   theProfile->loadParameters();
   theProfile->hadronicParameterization();
  
   // make hits
   makeHits(fastTrack);
 }

References ALCARECOTkAlJpsiMuMu_cff::charge, HCALHighEnergyHPDFilter_cfi::energy, Gflash::findShowerType(), GeV, GflashHadronShowerProfile::hadronicParameterization(), GflashHadronShowerProfile::initialize(), GflashHadronShowerProfile::loadParameters(), makeHits(), PbPb_ZMuSkimMuonDPG_cff::particleType, position, theAntiProtonProfile, theKaonMinusProfile, theKaonPlusProfile, thePiKProfile, theProfile, and theProtonProfile.

◆ excludeDetectorRegion()

G4bool GFlashHadronShowerModel::excludeDetectorRegion ( const G4FastTrack & fastTrack )

private

Definition at line 260 of file GFlashHadronShowerModel.cc.

                                                                                   {
   G4bool isExcluded = false;
   int verbosity = theParSet.getUntrackedParameter<int>("Verbosity");
  
   //exclude regions where geometry are complicated
   //+- one supermodule around the EB/EE boundary: 1.479 +- 0.0174*5
   G4double eta = fastTrack.GetPrimaryTrack()->GetPosition().pseudoRapidity();
   if (std::fabs(eta) > 1.392 && std::fabs(eta) < 1.566) {
     if (verbosity > 0) {
       edm::LogInfo("SimG4CoreApplication") << "GFlashHadronShowerModel: excluding region of eta = " << eta;
     }
     return true;
   } else {
     G4StepPoint* postStep = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint();
  
     Gflash::CalorimeterNumber kCalor = Gflash::getCalorimeterNumber(postStep->GetPosition() / cm);
     G4double distOut = 9999.0;
  
     //exclude the region where the shower starting point is inside the preshower
     if (std::fabs(eta) > Gflash::EtaMin[Gflash::kENCA] &&
         std::fabs((postStep->GetPosition()).getZ() / cm) < Gflash::Zmin[Gflash::kENCA]) {
       return true;
     }
  
     //<---the shower starting point is always inside envelopes
     //@@@exclude the region where the shower starting point is too close to the end of
     //the hadronic envelopes (may need to be optimized further!)
     //@@@if we extend parameterization including Magnet/HO, we need to change this strategy
     if (kCalor == Gflash::kHB) {
       distOut = Gflash::Rmax[Gflash::kHB] - postStep->GetPosition().getRho() / cm;
       if (distOut < Gflash::MinDistanceToOut)
         isExcluded = true;
     } else if (kCalor == Gflash::kHE) {
       distOut = Gflash::Zmax[Gflash::kHE] - std::fabs(postStep->GetPosition().getZ() / cm);
       if (distOut < Gflash::MinDistanceToOut)
         isExcluded = true;
     }
  
     //@@@remove this print statement later
     if (isExcluded && verbosity > 0) {
       edm::LogInfo("SimG4CoreApplication")
           << "GFlashHadronShowerModel: skipping kCalor = " << kCalor << " DistanceToOut " << distOut << " from ("
           << (postStep->GetPosition()).getRho() / cm << ":" << (postStep->GetPosition()).getZ() / cm
           << ") of KE = " << fastTrack.GetPrimaryTrack()->GetKineticEnergy() / GeV;
     }
   }
  
   return isExcluded;
 }

References PVValHelper::eta, Gflash::EtaMin, Gflash::getCalorimeterNumber(), edm::ParameterSet::getUntrackedParameter(), GeV, Gflash::kENCA, Gflash::kHB, Gflash::kHE, Gflash::MinDistanceToOut, Gflash::Rmax, theParSet, HIPAlignmentAlgorithm_cfi::verbosity, Gflash::Zmax, and Gflash::Zmin.

Referenced by ModelTrigger().

◆ IsApplicable()

G4bool GFlashHadronShowerModel::IsApplicable ( const G4ParticleDefinition & particleType )

override

Definition at line 65 of file GFlashHadronShowerModel.cc.

                                                                                      {
   return &particleType == G4PionMinus::PionMinusDefinition() || &particleType == G4PionPlus::PionPlusDefinition() ||
          &particleType == G4KaonMinus::KaonMinusDefinition() || &particleType == G4KaonPlus::KaonPlusDefinition() ||
          &particleType == G4AntiProton::AntiProtonDefinition() || &particleType == G4Proton::ProtonDefinition();
 }

References PbPb_ZMuSkimMuonDPG_cff::particleType.

◆ isFirstInelasticInteraction()

G4bool GFlashHadronShowerModel::isFirstInelasticInteraction ( const G4FastTrack & fastTrack )

private

Definition at line 205 of file GFlashHadronShowerModel.cc.

                                                                                         {
   G4bool isFirst = false;
  
   G4StepPoint* preStep = fastTrack.GetPrimaryTrack()->GetStep()->GetPreStepPoint();
   G4StepPoint* postStep = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint();
  
   G4String procName = postStep->GetProcessDefinedStep()->GetProcessName();
   G4ParticleDefinition* particleType = fastTrack.GetPrimaryTrack()->GetDefinition();
  
   //@@@ this part is still temporary and the cut for the variable ratio should be optimized later
  
   if ((particleType == G4PionPlus::PionPlusDefinition() && procName == "WrappedPionPlusInelastic") ||
       (particleType == G4PionMinus::PionMinusDefinition() && procName == "WrappedPionMinusInelastic") ||
       (particleType == G4KaonPlus::KaonPlusDefinition() && procName == "WrappedKaonPlusInelastic") ||
       (particleType == G4KaonMinus::KaonMinusDefinition() && procName == "WrappedKaonMinusInelastic") ||
       (particleType == G4AntiProton::AntiProtonDefinition() && procName == "WrappedAntiProtonInelastic") ||
       (particleType == G4Proton::ProtonDefinition() && procName == "WrappedProtonInelastic")) {
     //skip to the second interaction if the first inelastic is a quasi-elastic like interaction
     //@@@ the cut may be optimized later
  
     const G4TrackVector* fSecondaryVector = fastTrack.GetPrimaryTrack()->GetStep()->GetSecondary();
     G4double leadingEnergy = 0.0;
  
     //loop over 'all' secondaries including those produced by continuous processes.
     //@@@may require an additional condition only for hadron interaction with the process name,
     //but it will not change the result anyway
  
     for (unsigned int isec = 0; isec < fSecondaryVector->size(); isec++) {
       G4Track* fSecondaryTrack = (*fSecondaryVector)[isec];
       G4double secondaryEnergy = fSecondaryTrack->GetKineticEnergy();
  
       if (secondaryEnergy > leadingEnergy) {
         leadingEnergy = secondaryEnergy;
       }
     }
  
     if ((preStep->GetTotalEnergy() != 0) && (leadingEnergy / preStep->GetTotalEnergy() < Gflash::QuasiElasticLike))
       isFirst = true;
  
     //Fill debugging histograms and check information on secondaries -
     //remove after final implimentation
     /*
     if(theHisto->getStoreFlag()) {
       theHisto->preStepPosition->Fill(preStep->GetPosition().getRho()/cm);
       theHisto->postStepPosition->Fill(postStep->GetPosition().getRho()/cm);
       theHisto->deltaStep->Fill((postStep->GetPosition() - preStep->GetPosition()).getRho()/cm);
       theHisto->kineticEnergy->Fill(fastTrack.GetPrimaryTrack()->GetKineticEnergy()/GeV);
       theHisto->energyLoss->Fill(fabs(fastTrack.GetPrimaryTrack()->GetStep()->GetDeltaEnergy()/GeV));
       theHisto->energyRatio->Fill(leadingEnergy/preStep->GetTotalEnergy());
     }
     */
   }
   return isFirst;
 }

References cuy::isFirst, PbPb_ZMuSkimMuonDPG_cff::particleType, and Gflash::QuasiElasticLike.

◆ makeHits()

void GFlashHadronShowerModel::makeHits ( const G4FastTrack & fastTrack )

private

Definition at line 144 of file GFlashHadronShowerModel.cc.

                                                                    {
   std::vector<GflashHit>& gflashHitList = theProfile->getGflashHitList();
  
   theGflashStep->SetTrack(const_cast<G4Track*>(fastTrack.GetPrimaryTrack()));
   theGflashStep->GetPostStepPoint()->SetProcessDefinedStep(
       const_cast<G4VProcess*>(fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetProcessDefinedStep()));
   theGflashNavigator->SetWorldVolume(
       G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking()->GetWorldVolume());
  
   std::vector<GflashHit>::const_iterator spotIter = gflashHitList.begin();
   std::vector<GflashHit>::const_iterator spotIterEnd = gflashHitList.end();
  
   for (; spotIter != spotIterEnd; spotIter++) {
     theGflashNavigator->LocateGlobalPointAndUpdateTouchableHandle(
         spotIter->getPosition(), G4ThreeVector(0, 0, 0), theGflashTouchableHandle, false);
     updateGflashStep(spotIter->getPosition(), spotIter->getTime());
  
     // if there is a watcher defined in a job and the flag is turned on
     if (theWatcherOn) {
       theGflashStep->SetTotalEnergyDeposit(spotIter->getEnergy());
       SteppingAction* userSteppingAction = (SteppingAction*)G4EventManager::GetEventManager()->GetUserSteppingAction();
       userSteppingAction->m_g4StepSignal(theGflashStep);
     }
  
     G4VPhysicalVolume* aCurrentVolume = theGflashStep->GetPreStepPoint()->GetPhysicalVolume();
     if (aCurrentVolume == nullptr) {
       continue;
     }
  
     G4LogicalVolume* lv = aCurrentVolume->GetLogicalVolume();
     if (lv->GetRegion() != theRegion) {
       continue;
     }
  
     theGflashStep->GetPreStepPoint()->SetSensitiveDetector(aCurrentVolume->GetLogicalVolume()->GetSensitiveDetector());
     G4VSensitiveDetector* aSensitive = theGflashStep->GetPreStepPoint()->GetSensitiveDetector();
  
     if (aSensitive == nullptr)
       continue;
  
     G4String nameCalor = aCurrentVolume->GetName();
     nameCalor.assign(nameCalor, 0, 2);
     G4double samplingWeight = 1.0;
     if (nameCalor == "HB") {
       samplingWeight = Gflash::scaleSensitiveHB;
     } else if (nameCalor == "HE" || nameCalor == "HT") {
       samplingWeight = Gflash::scaleSensitiveHE;
     }
     theGflashStep->SetTotalEnergyDeposit(spotIter->getEnergy() * samplingWeight);
  
     aSensitive->Hit(theGflashStep);
   }
 }

References GflashHadronShowerProfile::getGflashHitList(), SteppingAction::m_g4StepSignal, Gflash::scaleSensitiveHB, Gflash::scaleSensitiveHE, theGflashNavigator, theGflashStep, theGflashTouchableHandle, theProfile, theRegion, theWatcherOn, and updateGflashStep().

Referenced by DoIt().

◆ ModelTrigger()

G4bool GFlashHadronShowerModel::ModelTrigger ( const G4FastTrack & fastTrack )

override

Definition at line 71 of file GFlashHadronShowerModel.cc.

                                                                          {
   // ModelTrigger returns false for Gflash Hadronic Shower Model if it is not
   // tested from the corresponding wrapper process, GflashHadronWrapperProcess.
   // Temporarily the track status is set to fPostponeToNextEvent at the wrapper
   // process before ModelTrigger is really tested for the second time through
   // PostStepGPIL of enviaged parameterization processes.  The better
   // implmentation may be using via G4VUserTrackInformation of each track, which
   // requires to modify a geant source code of stepping (G4SteppingManager2)
  
   // mininum energy cutoff to parameterize
   if (fastTrack.GetPrimaryTrack()->GetKineticEnergy() < GeV * Gflash::energyCutOff) {
     return false;
   }
  
   // This will be changed accordingly when the way
   // dealing with CaloRegion changes later.
   G4TouchableHistory* touch = (G4TouchableHistory*)(fastTrack.GetPrimaryTrack()->GetTouchable());
   G4VPhysicalVolume* pCurrentVolume = touch->GetVolume();
   if (pCurrentVolume == nullptr) {
     return false;
   }
  
   G4LogicalVolume* lv = pCurrentVolume->GetLogicalVolume();
   if (lv->GetRegion() != theRegion) {
     return false;
   }
  
   // check whether this is called from the normal GPIL or the wrapper process GPIL
   //if(fastTrack.GetPrimaryTrack()->GetTrackStatus() == fPostponeToNextEvent ) {
  
   // Shower pameterization start at the first inelastic interaction point
   //if(isFirstInelasticInteraction(fastTrack) && excludeDetectorRegion(fastTrack))
   if (excludeDetectorRegion(fastTrack)) {
     return false;
   }
  
   return true;
 }

References Gflash::energyCutOff, excludeDetectorRegion(), GeV, and theRegion.

◆ updateGflashStep()

void GFlashHadronShowerModel::updateGflashStep	(	const G4ThreeVector &	position,
		G4double	time
	)

private

Definition at line 198 of file GFlashHadronShowerModel.cc.

                                                                                                      {
   theGflashStep->GetPostStepPoint()->SetGlobalTime(timeGlobal);
   theGflashStep->GetPreStepPoint()->SetPosition(spotPosition);
   theGflashStep->GetPostStepPoint()->SetPosition(spotPosition);
   theGflashStep->GetPreStepPoint()->SetTouchableHandle(theGflashTouchableHandle);
 }