CMS 3D CMS Logo

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

#include <LowEnergyFastSimModel.h>

Inheritance diagram for LowEnergyFastSimModel:

Public Member Functions

void DoIt (const G4FastTrack &fastTrack, G4FastStep &fastStep) override
 
G4bool IsApplicable (const G4ParticleDefinition &particle) override
 
 LowEnergyFastSimModel (const G4String &name, G4Region *region, const edm::ParameterSet &parSet)
 
G4bool ModelTrigger (const G4FastTrack &fastTrack) override
 

Private Attributes

G4bool fCheck
 
G4double fEmax
 
GFlashHitMaker fHitMaker
 
const G4Material * fMaterial
 
LowEnergyFastSimParam fParam
 
const G4ParticleDefinition * fPositron
 
const G4Envelope * fRegion
 
G4ThreeVector fTailPos
 
const TrackingActionfTrackingAction
 

Detailed Description

Definition at line 18 of file LowEnergyFastSimModel.h.

Constructor & Destructor Documentation

◆ LowEnergyFastSimModel()

LowEnergyFastSimModel::LowEnergyFastSimModel ( const G4String &  name,
G4Region *  region,
const edm::ParameterSet parSet 
)

Definition at line 22 of file LowEnergyFastSimModel.cc.

References fEmax, fMaterial, fPositron, edm::ParameterSet::getParameter(), dqmiodumpmetadata::n, DD4hep2DDDName::noNameSpace(), and TableParser::table.

23  : G4VFastSimulationModel(name, region),
24  fRegion(region),
25  fTrackingAction(nullptr),
26  fCheck(true),
27  fTailPos(0., 0., 0.) {
28  fEmax = parSet.getParameter<double>("LowEnergyGflashEcalEmax") * CLHEP::GeV;
29  fPositron = G4Positron::Positron();
30  fMaterial = nullptr;
31  auto table = G4Material::GetMaterialTable();
32  for (auto const& mat : *table) {
33  G4String nam = (G4String)(DD4hep2DDDName::noNameSpace(mat->GetName()));
34  size_t n = nam.size();
35  if (n > 4) {
36  G4String sn = nam.substr(n - 5, 5);
37  if (sn == "PbWO4") {
38  fMaterial = mat;
39  break;
40  }
41  }
42  }
43  G4String nm = (nullptr == fMaterial) ? "not found!" : (G4String)(DD4hep2DDDName::noNameSpace(fMaterial->GetName()));
44  edm::LogVerbatim("LowEnergyFastSimModel") << "LowEGFlash material: <" << nm << ">";
45 }
Log< level::Info, true > LogVerbatim
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
const G4ParticleDefinition * fPositron
const G4Envelope * fRegion
const G4Material * fMaterial
const TrackingAction * fTrackingAction
std::string noNameSpace(const std::string &name)

Member Function Documentation

◆ DoIt()

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

Definition at line 74 of file LowEnergyFastSimModel.cc.

References funct::cos(), hcalRecHitTable_cff::energy, fHitMaker, fParam, fPositron, fTailPos, LowEnergyFastSimParam::GetInPointEnergyFraction(), LowEnergyFastSimParam::GetRadius(), LowEnergyFastSimParam::GetZ(), mps_fire::i, phi, alignCSCRings::r, scaleFactor, funct::sin(), HLT_2023v12_cff::track, twomass, and z.

74  {
75  fastStep.KillPrimaryTrack();
76  fastStep.ProposePrimaryTrackPathLength(0.0);
77  auto track = fastTrack.GetPrimaryTrack();
78  G4double energy = track->GetKineticEnergy() * scaleFactor;
79 
80  const G4ThreeVector& pos = track->GetPosition();
81 
82  G4double inPointEnergy = fParam.GetInPointEnergyFraction(energy) * energy;
83 
84  // take into account positron annihilation (not included in in-point)
85  if (fPositron == track->GetDefinition())
86  energy += twomass;
87 
88  const G4ThreeVector& momDir = track->GetMomentumDirection();
89 
90  // in point energy deposition
91  GFlashEnergySpot spot;
92  spot.SetEnergy(inPointEnergy);
93  spot.SetPosition(pos);
94  fHitMaker.make(&spot, &fastTrack);
95 
96  // Russian roulette
97  double wt2 = track->GetWeight();
98  if (wt2 <= 0.0) {
99  wt2 = 1.0;
100  }
101 
102  // tail energy deposition
103  const G4double etail = energy - inPointEnergy;
104  const G4int nspots = etail;
105  const G4double tailEnergy = etail * wt2 / (nspots + 1);
106  /*
107  edm::LogVerbatim("LowEnergyFastSimModel") << track->GetDefinition()->GetParticleName()
108  << " Ekin(MeV)=" << energy << " material: <"
109  << track->GetMaterial()->GetName()
110  << "> Elocal=" << inPointEnergy
111  << " Etail=" << tailEnergy
112  << " Nspots=" << nspots+1;
113  */
114  for (G4int i = 0; i <= nspots; ++i) {
115  const G4double r = fParam.GetRadius(energy);
116  const G4double z = fParam.GetZ();
117 
118  const G4double phi = CLHEP::twopi * G4UniformRand();
119  fTailPos.set(r * std::cos(phi), r * std::sin(phi), z);
120  fTailPos.rotateUz(momDir);
121  fTailPos += pos;
122 
123  spot.SetEnergy(tailEnergy);
124  spot.SetPosition(fTailPos);
125  fHitMaker.make(&spot, &fastTrack);
126  }
127 }
const G4ParticleDefinition * fPositron
constexpr G4double twomass
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
LowEnergyFastSimParam fParam
G4double GetInPointEnergyFraction(G4double energy) const
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
constexpr G4double scaleFactor
G4double GetRadius(G4double energy) const

◆ IsApplicable()

G4bool LowEnergyFastSimModel::IsApplicable ( const G4ParticleDefinition &  particle)
override

Definition at line 47 of file LowEnergyFastSimModel.cc.

References funct::abs().

47  {
48  return (11 == std::abs(particle.GetPDGEncoding()));
49 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22

◆ ModelTrigger()

G4bool LowEnergyFastSimModel::ModelTrigger ( const G4FastTrack &  fastTrack)
override

Definition at line 51 of file LowEnergyFastSimModel.cc.

References funct::abs(), hcalRecHitTable_cff::energy, fCheck, fEmax, fMaterial, fTrackingAction, TrackInformation::genParticlePID(), and HLT_2023v12_cff::track.

51  {
52  const G4Track* track = fastTrack.GetPrimaryTrack();
53  G4double energy = track->GetKineticEnergy();
54  if (fMaterial != track->GetMaterial() || energy >= fEmax)
55  return false;
56 
57  /*
58  edm::LogVerbatim("LowEnergyFastSimModel") << track->GetDefinition()->GetParticleName()
59  << " Ekin(MeV)=" << energy << " material: <"
60  << track->GetMaterial()->GetName() << ">";
61  */
62  if (fCheck) {
63  if (nullptr == fTrackingAction) {
64  fTrackingAction = static_cast<const TrackingAction*>(G4EventManager::GetEventManager()->GetUserTrackingAction());
65  }
66  const TrackInformation* ptrti = static_cast<TrackInformation*>(track->GetUserInformation());
67  int pdg = ptrti->genParticlePID();
68  if (std::abs(pdg) == 11 || pdg == 22)
69  return false;
70  }
71  return true;
72 }
const G4Material * fMaterial
const TrackingAction * fTrackingAction
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int genParticlePID() const

Member Data Documentation

◆ fCheck

G4bool LowEnergyFastSimModel::fCheck
private

Definition at line 32 of file LowEnergyFastSimModel.h.

Referenced by ModelTrigger().

◆ fEmax

G4double LowEnergyFastSimModel::fEmax
private

Definition at line 27 of file LowEnergyFastSimModel.h.

Referenced by LowEnergyFastSimModel(), and ModelTrigger().

◆ fHitMaker

GFlashHitMaker LowEnergyFastSimModel::fHitMaker
private

Definition at line 34 of file LowEnergyFastSimModel.h.

Referenced by DoIt().

◆ fMaterial

const G4Material* LowEnergyFastSimModel::fMaterial
private

Definition at line 31 of file LowEnergyFastSimModel.h.

Referenced by LowEnergyFastSimModel(), and ModelTrigger().

◆ fParam

LowEnergyFastSimParam LowEnergyFastSimModel::fParam
private

Definition at line 35 of file LowEnergyFastSimModel.h.

Referenced by DoIt().

◆ fPositron

const G4ParticleDefinition* LowEnergyFastSimModel::fPositron
private

Definition at line 30 of file LowEnergyFastSimModel.h.

Referenced by DoIt(), and LowEnergyFastSimModel().

◆ fRegion

const G4Envelope* LowEnergyFastSimModel::fRegion
private

Definition at line 28 of file LowEnergyFastSimModel.h.

◆ fTailPos

G4ThreeVector LowEnergyFastSimModel::fTailPos
private

Definition at line 33 of file LowEnergyFastSimModel.h.

Referenced by DoIt().

◆ fTrackingAction

const TrackingAction* LowEnergyFastSimModel::fTrackingAction
private

Definition at line 29 of file LowEnergyFastSimModel.h.

Referenced by ModelTrigger().