File indexing completed on 2023-03-17 11:24:56
0001
0002
0003
0004
0005 #include "SimG4Core/GFlash/interface/GflashEMShowerModel.h"
0006
0007 #include "SimGeneral/GFlash/interface/GflashEMShowerProfile.h"
0008 #include "SimGeneral/GFlash/interface/GflashHit.h"
0009
0010 #include "G4Electron.hh"
0011 #include "G4EventManager.hh"
0012 #include "G4FastSimulationManager.hh"
0013 #include "G4LogicalVolume.hh"
0014 #include "G4Positron.hh"
0015 #include "G4TouchableHandle.hh"
0016 #include "G4TransportationManager.hh"
0017 #include "G4VPhysicalVolume.hh"
0018 #include "G4VProcess.hh"
0019 #include "G4VSensitiveDetector.hh"
0020
0021 using namespace CLHEP;
0022
0023 GflashEMShowerModel::GflashEMShowerModel(const G4String &modelName,
0024 G4Envelope *envelope,
0025 const edm::ParameterSet &parSet)
0026 : G4VFastSimulationModel(modelName, envelope), theParSet(parSet) {
0027 theProfile = new GflashEMShowerProfile(parSet);
0028 theRegion = reinterpret_cast<G4Region *>(envelope);
0029
0030 theGflashStep = new G4Step();
0031 theGflashTouchableHandle = new G4TouchableHistory();
0032 theGflashNavigator = new G4Navigator();
0033 }
0034
0035
0036
0037 GflashEMShowerModel::~GflashEMShowerModel() {
0038 delete theProfile;
0039 delete theGflashStep;
0040 }
0041
0042 G4bool GflashEMShowerModel::IsApplicable(const G4ParticleDefinition &particleType) {
0043 return (&particleType == G4Electron::ElectronDefinition() || &particleType == G4Positron::PositronDefinition());
0044 }
0045
0046
0047 G4bool GflashEMShowerModel::ModelTrigger(const G4FastTrack &fastTrack) {
0048
0049 if (fastTrack.GetPrimaryTrack()->GetKineticEnergy() < GeV) {
0050 return false;
0051 }
0052 if (excludeDetectorRegion(fastTrack)) {
0053 return false;
0054 }
0055
0056
0057
0058 const G4VPhysicalVolume *pCurrentVolume = (fastTrack.GetPrimaryTrack()->GetTouchable())->GetVolume();
0059 if (pCurrentVolume == nullptr) {
0060 return false;
0061 }
0062
0063 const G4LogicalVolume *lv = pCurrentVolume->GetLogicalVolume();
0064 if (lv->GetRegion() != theRegion) {
0065 return false;
0066 }
0067 return true;
0068 }
0069
0070
0071 void GflashEMShowerModel::DoIt(const G4FastTrack &fastTrack, G4FastStep &fastStep) {
0072
0073 fastStep.KillPrimaryTrack();
0074 fastStep.ProposePrimaryTrackPathLength(0.0);
0075
0076
0077
0078 G4double energy = fastTrack.GetPrimaryTrack()->GetKineticEnergy() / GeV;
0079 G4double globalTime = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetGlobalTime();
0080 G4double charge = fastTrack.GetPrimaryTrack()->GetStep()->GetPreStepPoint()->GetCharge();
0081 G4ThreeVector position = fastTrack.GetPrimaryTrack()->GetPosition() / cm;
0082 G4ThreeVector momentum = fastTrack.GetPrimaryTrack()->GetMomentum() / GeV;
0083 G4int showerType = Gflash::findShowerType(position);
0084
0085
0086 theProfile->initialize(showerType, energy, globalTime, charge, position, momentum);
0087 theProfile->parameterization();
0088
0089
0090 makeHits(fastTrack);
0091 }
0092
0093 void GflashEMShowerModel::makeHits(const G4FastTrack &fastTrack) {
0094 std::vector<GflashHit> &gflashHitList = theProfile->getGflashHitList();
0095
0096 theGflashStep->SetTrack(const_cast<G4Track *>(fastTrack.GetPrimaryTrack()));
0097
0098 theGflashStep->GetPostStepPoint()->SetProcessDefinedStep(
0099 const_cast<G4VProcess *>(fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetProcessDefinedStep()));
0100 theGflashNavigator->SetWorldVolume(
0101 G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking()->GetWorldVolume());
0102
0103 std::vector<GflashHit>::const_iterator spotIter = gflashHitList.begin();
0104 std::vector<GflashHit>::const_iterator spotIterEnd = gflashHitList.end();
0105
0106 for (; spotIter != spotIterEnd; spotIter++) {
0107
0108
0109 theGflashNavigator->LocateGlobalPointAndUpdateTouchableHandle(
0110 spotIter->getPosition(), G4ThreeVector(0, 0, 0), theGflashTouchableHandle, false);
0111 updateGflashStep(spotIter->getPosition(), spotIter->getTime());
0112
0113
0114
0115
0116 G4VPhysicalVolume *aCurrentVolume = theGflashStep->GetPreStepPoint()->GetPhysicalVolume();
0117 if (aCurrentVolume == nullptr)
0118 continue;
0119
0120 G4LogicalVolume *lv = aCurrentVolume->GetLogicalVolume();
0121 if (lv->GetRegion() != theRegion)
0122 continue;
0123
0124 theGflashStep->GetPreStepPoint()->SetSensitiveDetector(aCurrentVolume->GetLogicalVolume()->GetSensitiveDetector());
0125 G4VSensitiveDetector *aSensitive = theGflashStep->GetPreStepPoint()->GetSensitiveDetector();
0126
0127 if (aSensitive == nullptr)
0128 continue;
0129
0130 theGflashStep->SetTotalEnergyDeposit(spotIter->getEnergy());
0131 aSensitive->Hit(theGflashStep);
0132 }
0133 }
0134
0135 void GflashEMShowerModel::updateGflashStep(const G4ThreeVector &spotPosition, G4double timeGlobal) {
0136 theGflashStep->GetPostStepPoint()->SetGlobalTime(timeGlobal);
0137 theGflashStep->GetPreStepPoint()->SetPosition(spotPosition);
0138 theGflashStep->GetPostStepPoint()->SetPosition(spotPosition);
0139 theGflashStep->GetPreStepPoint()->SetTouchableHandle(theGflashTouchableHandle);
0140 }
0141
0142
0143 G4bool GflashEMShowerModel::excludeDetectorRegion(const G4FastTrack &fastTrack) {
0144
0145
0146 G4double eta = fastTrack.GetPrimaryTrack()->GetPosition().pseudoRapidity();
0147 return (std::fabs(eta) > 1.392 && std::fabs(eta) < 1.566);
0148 }