CMS 3D CMS Logo

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

const G4double fEmax
 
GFlashHitMaker fHitMaker
 
const G4Envelope *const fRegion
 
LowEnergyFastSimParam param
 

Detailed Description

Definition at line 13 of file LowEnergyFastSimModel.h.

Constructor & Destructor Documentation

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

Definition at line 11 of file LowEnergyFastSimModel.cc.

12  : G4VFastSimulationModel(name, region),
13  fEmax(parSet.getParameter<double>("LowEnergyGflashEcalEmax")),
14  fRegion(region) {}
T getParameter(std::string const &) const
const G4Envelope *const fRegion

Member Function Documentation

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

Definition at line 25 of file LowEnergyFastSimModel.cc.

References funct::cos(), Basic3DVector< T >::cross(), cross(), HCALHighEnergyHPDFilter_cfi::energy, fHitMaker, LowEnergyFastSimParam::GetInPointEnergyFraction(), LowEnergyFastSimParam::GetRadius(), LowEnergyFastSimParam::GetZ(), mps_fire::i, createfilelist::int, param, phi, CosmicsPD_Skims::radius, funct::sin(), and z.

25  {
26  fastStep.KillPrimaryTrack();
27  fastStep.SetPrimaryTrackPathLength(0.0);
28  fastStep.SetTotalEnergyDeposited(fastTrack.GetPrimaryTrack()->GetKineticEnergy());
29 
30  const G4double energy = fastTrack.GetPrimaryTrack()->GetKineticEnergy();
31  const G4ThreeVector& pos = fastTrack.GetPrimaryTrack()->GetPosition();
32 
33  G4double inPointEnergy = param.GetInPointEnergyFraction(energy) * energy;
34 
35  const G4ThreeVector& momDir = fastTrack.GetPrimaryTrack()->GetMomentumDirection();
36  const G4ThreeVector& ortho = momDir.orthogonal();
37  const G4ThreeVector& cross = momDir.cross(ortho);
38 
39  // in point energy deposition
40  GFlashEnergySpot spot;
41  spot.SetEnergy(inPointEnergy);
42  spot.SetPosition(pos);
43  fHitMaker.make(&spot, &fastTrack);
44 
45  // tail energy deposition
46  G4double etail = energy - inPointEnergy;
47  const G4int nspots = int(etail) + 1;
48  for (G4int i = 0; i < nspots; ++i) {
49  const G4double radius = param.GetRadius(energy);
50  const G4double z = param.GetZ();
51 
52  const G4double phi = CLHEP::twopi * G4UniformRand();
53  const G4ThreeVector tailPos = pos + z * momDir + radius * std::cos(phi) * ortho + radius * std::sin(phi) * cross;
54 
55  const G4double tailEnergy = etail / nspots;
56 
57  spot.SetEnergy(tailEnergy);
58  spot.SetPosition(tailPos);
59  fHitMaker.make(&spot, &fastTrack);
60  }
61 }
G4double GetRadius(G4double energy) const
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Basic3DVector cross(const Basic3DVector &lh) const
Vector product, or "cross" product, with a vector of same type.
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
G4double GetInPointEnergyFraction(G4double energy) const
LowEnergyFastSimParam param
Basic3DVector cross(const Basic3DVector &v) const
Vector product, or "cross" product, with a vector of same type.
G4bool LowEnergyFastSimModel::IsApplicable ( const G4ParticleDefinition &  particle)
override

Definition at line 16 of file LowEnergyFastSimModel.cc.

16  {
17  return &particle == G4Electron::Definition();
18 }
G4bool LowEnergyFastSimModel::ModelTrigger ( const G4FastTrack &  fastTrack)
override

Definition at line 20 of file LowEnergyFastSimModel.cc.

References HCALHighEnergyHPDFilter_cfi::energy, fEmax, and fRegion.

20  {
21  G4double energy = fastTrack.GetPrimaryTrack()->GetKineticEnergy();
22  return energy < fEmax && fRegion == fastTrack.GetEnvelope();
23 }
const G4Envelope *const fRegion

Member Data Documentation

const G4double LowEnergyFastSimModel::fEmax
private

Definition at line 22 of file LowEnergyFastSimModel.h.

Referenced by ModelTrigger().

GFlashHitMaker LowEnergyFastSimModel::fHitMaker
private

Definition at line 24 of file LowEnergyFastSimModel.h.

Referenced by DoIt().

const G4Envelope* const LowEnergyFastSimModel::fRegion
private

Definition at line 23 of file LowEnergyFastSimModel.h.

Referenced by ModelTrigger().

LowEnergyFastSimParam LowEnergyFastSimModel::param
private

Definition at line 25 of file LowEnergyFastSimModel.h.

Referenced by DoIt().