CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
GFlashHadronShowerModel Class Reference

#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

GflashAntiProtonShowerProfiletheAntiProtonProfile
 
G4Navigator * theGflashNavigator
 
G4Step * theGflashStep
 
G4TouchableHandle theGflashTouchableHandle
 
GflashKaonMinusShowerProfiletheKaonMinusProfile
 
GflashKaonPlusShowerProfiletheKaonPlusProfile
 
edm::ParameterSet theParSet
 
GflashPiKShowerProfilethePiKProfile
 
GflashHadronShowerProfiletheProfile
 
GflashProtonShowerProfiletheProtonProfile
 
const G4Region * theRegion
 
G4bool theWatcherOn
 

Detailed Description

Definition at line 20 of file GFlashHadronShowerModel.h.

Constructor & Destructor Documentation

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

Definition at line 35 of file GFlashHadronShowerModel.cc.

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

38  : G4VFastSimulationModel(modelName, envelope), theParSet(parSet) {
39  theWatcherOn = parSet.getParameter<bool>("watcherOn");
40 
41  theProfile = nullptr;
47  //theHisto = GflashHistogram::instance();
48 
49  theRegion = const_cast<const G4Region*>(envelope);
50 
51  theGflashStep = new G4Step();
52  theGflashTouchableHandle = new G4TouchableHistory();
53  theGflashNavigator = new G4Navigator();
54 }
T getParameter(std::string const &) const
GflashPiKShowerProfile * thePiKProfile
GflashKaonMinusShowerProfile * theKaonMinusProfile
G4TouchableHandle theGflashTouchableHandle
GflashKaonPlusShowerProfile * theKaonPlusProfile
GflashAntiProtonShowerProfile * theAntiProtonProfile
GflashHadronShowerProfile * theProfile
GflashProtonShowerProfile * theProtonProfile
GFlashHadronShowerModel::~GFlashHadronShowerModel ( )
override

Definition at line 56 of file GFlashHadronShowerModel.cc.

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

56  {
57  delete thePiKProfile;
58  delete theKaonPlusProfile;
59  delete theKaonMinusProfile;
60  delete theProtonProfile;
61  delete theAntiProtonProfile;
62  delete theGflashStep;
63 }
GflashPiKShowerProfile * thePiKProfile
GflashKaonMinusShowerProfile * theKaonMinusProfile
GflashKaonPlusShowerProfile * theKaonPlusProfile
GflashAntiProtonShowerProfile * theAntiProtonProfile
GflashProtonShowerProfile * theProtonProfile

Member Function Documentation

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

Definition at line 110 of file GFlashHadronShowerModel.cc.

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.

110  {
111  // kill the particle
112  fastStep.KillPrimaryTrack();
113  fastStep.ProposePrimaryTrackPathLength(0.0);
114 
115  // parameterize energy depostion by the particle type
116  G4ParticleDefinition* particleType = fastTrack.GetPrimaryTrack()->GetDefinition();
117 
119  if (particleType == G4KaonMinus::KaonMinusDefinition())
121  else if (particleType == G4KaonPlus::KaonPlusDefinition())
123  else if (particleType == G4AntiProton::AntiProtonDefinition())
125  else if (particleType == G4Proton::ProtonDefinition())
127 
128  //input variables for GflashHadronShowerProfile
129  G4double energy = fastTrack.GetPrimaryTrack()->GetKineticEnergy() / GeV;
130  G4double globalTime = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetGlobalTime();
131  G4double charge = fastTrack.GetPrimaryTrack()->GetStep()->GetPreStepPoint()->GetCharge();
132  G4ThreeVector position = fastTrack.GetPrimaryTrack()->GetPosition() / cm;
133  G4ThreeVector momentum = fastTrack.GetPrimaryTrack()->GetMomentum() / GeV;
134  G4int showerType = Gflash::findShowerType(position);
135 
136  theProfile->initialize(showerType, energy, globalTime, charge, position, momentum);
139 
140  // make hits
141  makeHits(fastTrack);
142 }
const double GeV
Definition: MathUtil.h:16
GflashPiKShowerProfile * thePiKProfile
GflashKaonMinusShowerProfile * theKaonMinusProfile
GflashKaonPlusShowerProfile * theKaonPlusProfile
GflashAntiProtonShowerProfile * theAntiProtonProfile
GflashHadronShowerProfile * theProfile
void makeHits(const G4FastTrack &fastTrack)
GflashProtonShowerProfile * theProtonProfile
int findShowerType(const Gflash3Vector &position)
static int position[264][3]
Definition: ReadPGInfo.cc:289
void initialize(int showerType, double energy, double globalTime, double charge, Gflash3Vector &position, Gflash3Vector &momentum)
G4bool GFlashHadronShowerModel::excludeDetectorRegion ( const G4FastTrack &  fastTrack)
private

Definition at line 260 of file GFlashHadronShowerModel.cc.

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().

260  {
261  G4bool isExcluded = false;
262  int verbosity = theParSet.getUntrackedParameter<int>("Verbosity");
263 
264  //exclude regions where geometry are complicated
265  //+- one supermodule around the EB/EE boundary: 1.479 +- 0.0174*5
266  G4double eta = fastTrack.GetPrimaryTrack()->GetPosition().pseudoRapidity();
267  if (std::fabs(eta) > 1.392 && std::fabs(eta) < 1.566) {
268  if (verbosity > 0) {
269  edm::LogInfo("SimG4CoreApplication") << "GFlashHadronShowerModel: excluding region of eta = " << eta;
270  }
271  return true;
272  } else {
273  G4StepPoint* postStep = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint();
274 
275  Gflash::CalorimeterNumber kCalor = Gflash::getCalorimeterNumber(postStep->GetPosition() / cm);
276  G4double distOut = 9999.0;
277 
278  //exclude the region where the shower starting point is inside the preshower
279  if (std::fabs(eta) > Gflash::EtaMin[Gflash::kENCA] &&
280  std::fabs((postStep->GetPosition()).getZ() / cm) < Gflash::Zmin[Gflash::kENCA]) {
281  return true;
282  }
283 
284  //<---the shower starting point is always inside envelopes
285  //@@@exclude the region where the shower starting point is too close to the end of
286  //the hadronic envelopes (may need to be optimized further!)
287  //@@@if we extend parameterization including Magnet/HO, we need to change this strategy
288  if (kCalor == Gflash::kHB) {
289  distOut = Gflash::Rmax[Gflash::kHB] - postStep->GetPosition().getRho() / cm;
290  if (distOut < Gflash::MinDistanceToOut)
291  isExcluded = true;
292  } else if (kCalor == Gflash::kHE) {
293  distOut = Gflash::Zmax[Gflash::kHE] - std::fabs(postStep->GetPosition().getZ() / cm);
294  if (distOut < Gflash::MinDistanceToOut)
295  isExcluded = true;
296  }
297 
298  //@@@remove this print statement later
299  if (isExcluded && verbosity > 0) {
300  edm::LogInfo("SimG4CoreApplication")
301  << "GFlashHadronShowerModel: skipping kCalor = " << kCalor << " DistanceToOut " << distOut << " from ("
302  << (postStep->GetPosition()).getRho() / cm << ":" << (postStep->GetPosition()).getZ() / cm
303  << ") of KE = " << fastTrack.GetPrimaryTrack()->GetKineticEnergy() / GeV;
304  }
305  }
306 
307  return isExcluded;
308 }
T getUntrackedParameter(std::string const &, T const &) const
const double MinDistanceToOut
const double GeV
Definition: MathUtil.h:16
const double Zmax[kNumberCalorimeter]
CalorimeterNumber getCalorimeterNumber(const Gflash3Vector &position)
const double EtaMin[kNumberCalorimeter]
const double Zmin[kNumberCalorimeter]
const double Rmax[kNumberCalorimeter]
G4bool GFlashHadronShowerModel::IsApplicable ( const G4ParticleDefinition &  particleType)
override

Definition at line 65 of file GFlashHadronShowerModel.cc.

65  {
66  return &particleType == G4PionMinus::PionMinusDefinition() || &particleType == G4PionPlus::PionPlusDefinition() ||
67  &particleType == G4KaonMinus::KaonMinusDefinition() || &particleType == G4KaonPlus::KaonPlusDefinition() ||
68  &particleType == G4AntiProton::AntiProtonDefinition() || &particleType == G4Proton::ProtonDefinition();
69 }
G4bool GFlashHadronShowerModel::isFirstInelasticInteraction ( const G4FastTrack &  fastTrack)
private

Definition at line 205 of file GFlashHadronShowerModel.cc.

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

205  {
206  G4bool isFirst = false;
207 
208  G4StepPoint* preStep = fastTrack.GetPrimaryTrack()->GetStep()->GetPreStepPoint();
209  G4StepPoint* postStep = fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint();
210 
211  G4String procName = postStep->GetProcessDefinedStep()->GetProcessName();
212  G4ParticleDefinition* particleType = fastTrack.GetPrimaryTrack()->GetDefinition();
213 
214  //@@@ this part is still temporary and the cut for the variable ratio should be optimized later
215 
216  if ((particleType == G4PionPlus::PionPlusDefinition() && procName == "WrappedPionPlusInelastic") ||
217  (particleType == G4PionMinus::PionMinusDefinition() && procName == "WrappedPionMinusInelastic") ||
218  (particleType == G4KaonPlus::KaonPlusDefinition() && procName == "WrappedKaonPlusInelastic") ||
219  (particleType == G4KaonMinus::KaonMinusDefinition() && procName == "WrappedKaonMinusInelastic") ||
220  (particleType == G4AntiProton::AntiProtonDefinition() && procName == "WrappedAntiProtonInelastic") ||
221  (particleType == G4Proton::ProtonDefinition() && procName == "WrappedProtonInelastic")) {
222  //skip to the second interaction if the first inelastic is a quasi-elastic like interaction
223  //@@@ the cut may be optimized later
224 
225  const G4TrackVector* fSecondaryVector = fastTrack.GetPrimaryTrack()->GetStep()->GetSecondary();
226  G4double leadingEnergy = 0.0;
227 
228  //loop over 'all' secondaries including those produced by continuous processes.
229  //@@@may require an additional condition only for hadron interaction with the process name,
230  //but it will not change the result anyway
231 
232  for (unsigned int isec = 0; isec < fSecondaryVector->size(); isec++) {
233  G4Track* fSecondaryTrack = (*fSecondaryVector)[isec];
234  G4double secondaryEnergy = fSecondaryTrack->GetKineticEnergy();
235 
236  if (secondaryEnergy > leadingEnergy) {
237  leadingEnergy = secondaryEnergy;
238  }
239  }
240 
241  if ((preStep->GetTotalEnergy() != 0) && (leadingEnergy / preStep->GetTotalEnergy() < Gflash::QuasiElasticLike))
242  isFirst = true;
243 
244  //Fill debugging histograms and check information on secondaries -
245  //remove after final implimentation
246  /*
247  if(theHisto->getStoreFlag()) {
248  theHisto->preStepPosition->Fill(preStep->GetPosition().getRho()/cm);
249  theHisto->postStepPosition->Fill(postStep->GetPosition().getRho()/cm);
250  theHisto->deltaStep->Fill((postStep->GetPosition() - preStep->GetPosition()).getRho()/cm);
251  theHisto->kineticEnergy->Fill(fastTrack.GetPrimaryTrack()->GetKineticEnergy()/GeV);
252  theHisto->energyLoss->Fill(fabs(fastTrack.GetPrimaryTrack()->GetStep()->GetDeltaEnergy()/GeV));
253  theHisto->energyRatio->Fill(leadingEnergy/preStep->GetTotalEnergy());
254  }
255  */
256  }
257  return isFirst;
258 }
isFirst
Definition: cuy.py:419
const double QuasiElasticLike
void GFlashHadronShowerModel::makeHits ( const G4FastTrack &  fastTrack)
private

Definition at line 144 of file GFlashHadronShowerModel.cc.

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

Referenced by DoIt().

144  {
145  std::vector<GflashHit>& gflashHitList = theProfile->getGflashHitList();
146 
147  theGflashStep->SetTrack(const_cast<G4Track*>(fastTrack.GetPrimaryTrack()));
148  theGflashStep->GetPostStepPoint()->SetProcessDefinedStep(
149  const_cast<G4VProcess*>(fastTrack.GetPrimaryTrack()->GetStep()->GetPostStepPoint()->GetProcessDefinedStep()));
150  theGflashNavigator->SetWorldVolume(
151  G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking()->GetWorldVolume());
152 
153  std::vector<GflashHit>::const_iterator spotIter = gflashHitList.begin();
154  std::vector<GflashHit>::const_iterator spotIterEnd = gflashHitList.end();
155 
156  for (; spotIter != spotIterEnd; spotIter++) {
157  theGflashNavigator->LocateGlobalPointAndUpdateTouchableHandle(
158  spotIter->getPosition(), G4ThreeVector(0, 0, 0), theGflashTouchableHandle, false);
159  updateGflashStep(spotIter->getPosition(), spotIter->getTime());
160 
161  // if there is a watcher defined in a job and the flag is turned on
162  if (theWatcherOn) {
163  theGflashStep->SetTotalEnergyDeposit(spotIter->getEnergy());
164  SteppingAction* userSteppingAction = (SteppingAction*)G4EventManager::GetEventManager()->GetUserSteppingAction();
165  userSteppingAction->m_g4StepSignal(theGflashStep);
166  }
167 
168  G4VPhysicalVolume* aCurrentVolume = theGflashStep->GetPreStepPoint()->GetPhysicalVolume();
169  if (aCurrentVolume == nullptr) {
170  continue;
171  }
172 
173  G4LogicalVolume* lv = aCurrentVolume->GetLogicalVolume();
174  if (lv->GetRegion() != theRegion) {
175  continue;
176  }
177 
178  theGflashStep->GetPreStepPoint()->SetSensitiveDetector(aCurrentVolume->GetLogicalVolume()->GetSensitiveDetector());
179  G4VSensitiveDetector* aSensitive = theGflashStep->GetPreStepPoint()->GetSensitiveDetector();
180 
181  if (aSensitive == nullptr)
182  continue;
183 
184  G4String nameCalor = aCurrentVolume->GetName();
185  nameCalor.assign(nameCalor, 0, 2);
186  G4double samplingWeight = 1.0;
187  if (nameCalor == "HB") {
188  samplingWeight = Gflash::scaleSensitiveHB;
189  } else if (nameCalor == "HE" || nameCalor == "HT") {
190  samplingWeight = Gflash::scaleSensitiveHE;
191  }
192  theGflashStep->SetTotalEnergyDeposit(spotIter->getEnergy() * samplingWeight);
193 
194  aSensitive->Hit(theGflashStep);
195  }
196 }
SimActivityRegistry::G4StepSignal m_g4StepSignal
G4TouchableHandle theGflashTouchableHandle
GflashHadronShowerProfile * theProfile
void updateGflashStep(const G4ThreeVector &position, G4double time)
const double scaleSensitiveHE
std::vector< GflashHit > & getGflashHitList()
const double scaleSensitiveHB
G4bool GFlashHadronShowerModel::ModelTrigger ( const G4FastTrack &  fastTrack)
override

Definition at line 71 of file GFlashHadronShowerModel.cc.

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

71  {
72  // ModelTrigger returns false for Gflash Hadronic Shower Model if it is not
73  // tested from the corresponding wrapper process, GflashHadronWrapperProcess.
74  // Temporarily the track status is set to fPostponeToNextEvent at the wrapper
75  // process before ModelTrigger is really tested for the second time through
76  // PostStepGPIL of enviaged parameterization processes. The better
77  // implmentation may be using via G4VUserTrackInformation of each track, which
78  // requires to modify a geant source code of stepping (G4SteppingManager2)
79 
80  // mininum energy cutoff to parameterize
81  if (fastTrack.GetPrimaryTrack()->GetKineticEnergy() < GeV * Gflash::energyCutOff) {
82  return false;
83  }
84 
85  // This will be changed accordingly when the way
86  // dealing with CaloRegion changes later.
87  G4TouchableHistory* touch = (G4TouchableHistory*)(fastTrack.GetPrimaryTrack()->GetTouchable());
88  G4VPhysicalVolume* pCurrentVolume = touch->GetVolume();
89  if (pCurrentVolume == nullptr) {
90  return false;
91  }
92 
93  G4LogicalVolume* lv = pCurrentVolume->GetLogicalVolume();
94  if (lv->GetRegion() != theRegion) {
95  return false;
96  }
97 
98  // check whether this is called from the normal GPIL or the wrapper process GPIL
99  //if(fastTrack.GetPrimaryTrack()->GetTrackStatus() == fPostponeToNextEvent ) {
100 
101  // Shower pameterization start at the first inelastic interaction point
102  //if(isFirstInelasticInteraction(fastTrack) && excludeDetectorRegion(fastTrack))
103  if (excludeDetectorRegion(fastTrack)) {
104  return false;
105  }
106 
107  return true;
108 }
const double GeV
Definition: MathUtil.h:16
G4bool excludeDetectorRegion(const G4FastTrack &fastTrack)
const double energyCutOff
void GFlashHadronShowerModel::updateGflashStep ( const G4ThreeVector &  position,
G4double  time 
)
private

Definition at line 198 of file GFlashHadronShowerModel.cc.

References theGflashStep, and theGflashTouchableHandle.

Referenced by makeHits().

198  {
199  theGflashStep->GetPostStepPoint()->SetGlobalTime(timeGlobal);
200  theGflashStep->GetPreStepPoint()->SetPosition(spotPosition);
201  theGflashStep->GetPostStepPoint()->SetPosition(spotPosition);
202  theGflashStep->GetPreStepPoint()->SetTouchableHandle(theGflashTouchableHandle);
203 }
G4TouchableHandle theGflashTouchableHandle

Member Data Documentation

GflashAntiProtonShowerProfile* GFlashHadronShowerModel::theAntiProtonProfile
private
G4Navigator* GFlashHadronShowerModel::theGflashNavigator
private

Definition at line 55 of file GFlashHadronShowerModel.h.

Referenced by GFlashHadronShowerModel(), and makeHits().

G4Step* GFlashHadronShowerModel::theGflashStep
private
G4TouchableHandle GFlashHadronShowerModel::theGflashTouchableHandle
private

Definition at line 56 of file GFlashHadronShowerModel.h.

Referenced by GFlashHadronShowerModel(), makeHits(), and updateGflashStep().

GflashKaonMinusShowerProfile* GFlashHadronShowerModel::theKaonMinusProfile
private
GflashKaonPlusShowerProfile* GFlashHadronShowerModel::theKaonPlusProfile
private
edm::ParameterSet GFlashHadronShowerModel::theParSet
private

Definition at line 44 of file GFlashHadronShowerModel.h.

Referenced by excludeDetectorRegion().

GflashPiKShowerProfile* GFlashHadronShowerModel::thePiKProfile
private
GflashHadronShowerProfile* GFlashHadronShowerModel::theProfile
private

Definition at line 45 of file GFlashHadronShowerModel.h.

Referenced by DoIt(), GFlashHadronShowerModel(), and makeHits().

GflashProtonShowerProfile* GFlashHadronShowerModel::theProtonProfile
private
const G4Region* GFlashHadronShowerModel::theRegion
private

Definition at line 52 of file GFlashHadronShowerModel.h.

Referenced by GFlashHadronShowerModel(), makeHits(), and ModelTrigger().

G4bool GFlashHadronShowerModel::theWatcherOn
private

Definition at line 43 of file GFlashHadronShowerModel.h.

Referenced by GFlashHadronShowerModel(), and makeHits().