Application/src/GFlashEMShowerModel.cc

0001 //
0002 // initial setup : E.Barberio & Joanna Weng
0003 // big changes : Soon Jun & Dongwook Jang
0004 // V.Ivanchenko rename the class, cleanup, and move
0005 //              to SimG4Core/Application - 2012/08/14
0006 //
0007 #include "SimG4Core/Application/interface/SteppingAction.h"
0008 #include "SimG4Core/Application/interface/GFlashEMShowerModel.h"
0009
0010 #include "SimGeneral/GFlash/interface/GflashEMShowerProfile.h"
0011 #include "SimGeneral/GFlash/interface/GflashHit.h"
0012
0013 #include "G4Electron.hh"
0014 #include "G4Positron.hh"
0015 #include "G4VProcess.hh"
0016 #include "G4VPhysicalVolume.hh"
0017 #include "G4LogicalVolume.hh"
0018 #include "G4TransportationManager.hh"
0019 #include "G4EventManager.hh"
0020 #include "G4FastSimulationManager.hh"
0021 #include "G4TouchableHandle.hh"
0022 #include "G4VSensitiveDetector.hh"
0023 #include "G4SystemOfUnits.hh"
0024
0025 GFlashEMShowerModel::GFlashEMShowerModel(const G4String& modelName,
0026                                          G4Envelope* envelope,
0027                                          const edm::ParameterSet& parSet)
0028     : G4VFastSimulationModel(modelName, envelope), theParSet(parSet) {
0029   theWatcherOn = parSet.getParameter<bool>("watcherOn");
0030
0031   theProfile = new GflashEMShowerProfile(parSet);
0032   theRegion = reinterpret_cast<G4Region*>(envelope);
0033
0034   theGflashStep = new G4Step();
0035   theGflashTouchableHandle = new G4TouchableHistory();
0036   theGflashNavigator = new G4Navigator();
0037 }
0038
0039 // --------------------------------------------------------------------------
0040
0041 GFlashEMShowerModel::~GFlashEMShowerModel() {
0042   delete theProfile;
0043   delete theGflashStep;
0044 }
0045
0046 G4bool GFlashEMShowerModel::IsApplicable(const G4ParticleDefinition& particleType) {
0047   return (&particleType == G4Electron::Electron() || &particleType == G4Positron::Positron());
0048 }
0049
0050 // --------------------------------------------------------------------------
0051 G4bool GFlashEMShowerModel::ModelTrigger(const G4FastTrack& fastTrack) {
0052   // Mininum energy cutoff to parameterize
0053   if (fastTrack.GetPrimaryTrack()->GetKineticEnergy() < Gflash::energyCutOff) {
0054     return false;
0055   }
0056   if (excludeDetectorRegion(fastTrack)) {
0057     return false;
0058   }
0059
0060   // This will be changed accordingly when the way
0061   // dealing with CaloRegion changes later.
0062   const G4TouchableHistory* touch = static_cast<const G4TouchableHistory*>(fastTrack.GetPrimaryTrack()->GetTouchable());
0063   const G4VPhysicalVolume* pCurrentVolume = touch->GetVolume();
0064   if (pCurrentVolume == nullptr) {
0065     return false;
0066   }
0067
0068   const G4LogicalVolume* lv = pCurrentVolume->GetLogicalVolume();
0069   if (lv->GetRegion() != theRegion) {
0070     return false;
0071   }
0072   return true;
0073 }
0074
0075 // ---------------------------------------------------------------------------
0076 void GFlashEMShowerModel::DoIt(const G4FastTrack& fastTrack, G4FastStep& fastStep) {
0077   // Kill the parameterised particle:
0078   fastStep.KillPrimaryTrack();
0079   fastStep.ProposePrimaryTrackPathLength(0.0);
0080
0081   // Input variables for GFlashEMShowerProfile with showerType = 1,5
0082   // Shower starts inside crystals
0083   G4double energy = fastTrack.GetPrimaryTrack()->GetKineticEnergy() / GeV;
0084   G4double globalTime = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetGlobalTime();
0085   G4double charge = fastTrack.GetPrimaryTrack()->GetStep()->GetPreStepPoint()->GetCharge();
0086   G4ThreeVector position = fastTrack.GetPrimaryTrack()->GetPosition() / cm;
0087   G4ThreeVector momentum = fastTrack.GetPrimaryTrack()->GetMomentum() / GeV;
0088   G4int showerType = Gflash::findShowerType(position);
0089
0090   // Do actual parameterization
0091   // The result of parameterization is gflashHitList
0092   theProfile->initialize(showerType, energy, globalTime, charge, position, momentum);
0093   theProfile->parameterization();
0094
0095   // Make hits
0096   makeHits(fastTrack);
0097 }
0098
0099 // ---------------------------------------------------------------------------
0100 void GFlashEMShowerModel::makeHits(const G4FastTrack& fastTrack) {
0101   std::vector<GflashHit>& gflashHitList = theProfile->getGflashHitList();
0102
0103   theGflashStep->SetTrack(const_cast<G4Track*>(fastTrack.GetPrimaryTrack()));
0104
0105   theGflashStep->GetPostStepPoint()->SetProcessDefinedStep(
0106       const_cast<G4VProcess*>(fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetProcessDefinedStep()));
0107   theGflashNavigator->SetWorldVolume(
0108       G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking()->GetWorldVolume());
0109
0110   std::vector<GflashHit>::const_iterator spotIter = gflashHitList.begin();
0111   std::vector<GflashHit>::const_iterator spotIterEnd = gflashHitList.end();
0112
0113   for (; spotIter != spotIterEnd; ++spotIter) {
0114     // Put touchable for each hit so that touchable history
0115     //     keeps track of each step.
0116     theGflashNavigator->LocateGlobalPointAndUpdateTouchableHandle(
0117         spotIter->getPosition(), G4ThreeVector(0, 0, 0), theGflashTouchableHandle, false);
0118     updateGflashStep(spotIter->getPosition(), spotIter->getTime());
0119
0120     // If there is a watcher defined in a job and the flag is turned on
0121     if (theWatcherOn) {
0122       SteppingAction* userSteppingAction = (SteppingAction*)G4EventManager::GetEventManager()->GetUserSteppingAction();
0123       userSteppingAction->m_g4StepSignal(theGflashStep);
0124     }
0125
0126     // Send G4Step information to Hit/Digi if the volume is sensitive
0127     // Copied from G4SteppingManager.cc
0128
0129     G4VPhysicalVolume* aCurrentVolume = theGflashStep->GetPreStepPoint()->GetPhysicalVolume();
0130     if (aCurrentVolume == nullptr) {
0131       continue;
0132     }
0133
0134     G4LogicalVolume* lv = aCurrentVolume->GetLogicalVolume();
0135     if (lv->GetRegion() != theRegion) {
0136       continue;
0137     }
0138
0139     theGflashStep->GetPreStepPoint()->SetSensitiveDetector(aCurrentVolume->GetLogicalVolume()->GetSensitiveDetector());
0140     G4VSensitiveDetector* aSensitive = theGflashStep->GetPreStepPoint()->GetSensitiveDetector();
0141
0142     if (aSensitive == nullptr) {
0143       continue;
0144     }
0145
0146     theGflashStep->SetTotalEnergyDeposit(spotIter->getEnergy());
0147     aSensitive->Hit(theGflashStep);
0148   }
0149 }
0150
0151 // ---------------------------------------------------------------------------
0152 void GFlashEMShowerModel::updateGflashStep(const G4ThreeVector& spotPosition, G4double timeGlobal) {
0153   theGflashStep->GetPostStepPoint()->SetGlobalTime(timeGlobal);
0154   theGflashStep->GetPreStepPoint()->SetPosition(spotPosition);
0155   theGflashStep->GetPostStepPoint()->SetPosition(spotPosition);
0156   theGflashStep->GetPreStepPoint()->SetTouchableHandle(theGflashTouchableHandle);
0157 }
0158
0159 // ---------------------------------------------------------------------------
0160 G4bool GFlashEMShowerModel::excludeDetectorRegion(const G4FastTrack& fastTrack) {
0161   //exclude regions where geometry are complicated
0162   //+- one supermodule around the EB/EE boundary: 1.479 +- 0.0174*5
0163   G4double eta = fastTrack.GetPrimaryTrack()->GetPosition().pseudoRapidity();
0164   return (std::fabs(eta) > 1.392 && std::fabs(eta) < 1.566);
0165 }
0166
0167 // ---------------------------------------------------------------------------