ATLAS Offline Software
Public Member Functions | Protected Attributes | List of all members
StoppedParticleFastSim Class Reference

#include <StoppedParticleFastSim.h>

Inheritance diagram for StoppedParticleFastSim:
Collaboration diagram for StoppedParticleFastSim:

Public Member Functions

 StoppedParticleFastSim (const std::string &name, const std::string &fsSDname)
 
 ~StoppedParticleFastSim ()
 
G4bool IsApplicable (const G4ParticleDefinition &) override final
 
virtual G4bool ModelTrigger (const G4FastTrack &) override final
 
void DoIt (const G4FastTrack &, G4FastStep &) override final
 

Protected Attributes

TrackFastSimSDm_fsSD {}
 
bool m_init {false}
 
std::string m_fsSDname {""}
 

Detailed Description

Definition at line 16 of file StoppedParticleFastSim.h.

Constructor & Destructor Documentation

◆ StoppedParticleFastSim()

StoppedParticleFastSim::StoppedParticleFastSim ( const std::string &  name,
const std::string &  fsSDname 
)

Definition at line 23 of file StoppedParticleFastSim.cxx.

24  : G4VFastSimulationModel(name)
25  , m_fsSDname(fsSDname)
26 {
27 }

◆ ~StoppedParticleFastSim()

StoppedParticleFastSim::~StoppedParticleFastSim ( )
inline

Definition at line 21 of file StoppedParticleFastSim.h.

21 {}

Member Function Documentation

◆ DoIt()

void StoppedParticleFastSim::DoIt ( const G4FastTrack &  fastTrack,
G4FastStep &  fastStep 
)
finaloverride

Definition at line 56 of file StoppedParticleFastSim.cxx.

57 {
58  if (!m_init){
59  m_init = true;
60 
61  G4SDManager *sdm = G4SDManager::GetSDMpointer();
62  G4VSensitiveDetector * vsd = sdm->FindSensitiveDetector( m_fsSDname );
63  if (vsd) {
64  m_fsSD = dynamic_cast<TrackFastSimSD*>(vsd);
65  if (!m_fsSD) {
66  G4ExceptionDescription description;
67  description << "DoIt: Could not cast the SD into an instance of TrackFasSimSD.";
68  G4Exception("StoppedParticleFastSim", "MissingTrackFastSimSD", FatalException, description);
69  abort();
70  }
71  }
72  else {
73  G4cout << "StoppedParticleFastSim::DoIt INFO Could not get TrackFastSimSD sensitive detector. If you are not writing track records this is expected." << G4endl;
74  }
75  // found the SD
76  } // End of lazy init
77  const int id = fastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetDefinition()->GetPDGEncoding();
78  if (m_fsSD &&
79  (MC::isSquarkLH(id) ||
80  id == 1000021 || // gluino
81  MC::isRHadron(id))) {
82  m_fsSD->WriteTrack( fastTrack.GetPrimaryTrack() , false , true );
83  }
84  fastStep.KillPrimaryTrack();
85 }

◆ IsApplicable()

G4bool StoppedParticleFastSim::IsApplicable ( const G4ParticleDefinition &  )
finaloverride

Definition at line 29 of file StoppedParticleFastSim.cxx.

30 {
31  return true;
32 }

◆ ModelTrigger()

G4bool StoppedParticleFastSim::ModelTrigger ( const G4FastTrack &  fastTrack)
finaloverridevirtual

Definition at line 34 of file StoppedParticleFastSim.cxx.

35 {
36  // Trigger if the energy is below our threshold or if the time is over 150 ns
37  int id = fastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetDefinition()->GetPDGEncoding();
38  if (id<1000000 || id>1100000) return true; // skip SM particles and super-partners of RH-fermions
39  if (MC::isSquarkLH(id) ||
40  id == 1000021 || // gluino
41  MC::isRHadron(id)) {
42  G4Material * mat = fastTrack.GetPrimaryTrack()->GetMaterial();
43  double minA=1500000.;
44  for (unsigned int i=0;i<mat->GetNumberOfElements();++i){
45  if (mat->GetElement(i) &&
46  minA>mat->GetElement(i)->GetN()){
47  minA=mat->GetElement(i)->GetN();
48  }
49  }
50  if (fastTrack.GetPrimaryTrack()->GetVelocity()<0.15*std::pow(minA,-2./3.)*CLHEP::c_light) return true;
51  return false;
52  }
53  return true;
54 }

Member Data Documentation

◆ m_fsSD

TrackFastSimSD* StoppedParticleFastSim::m_fsSD {}
protected

Definition at line 30 of file StoppedParticleFastSim.h.

◆ m_fsSDname

std::string StoppedParticleFastSim::m_fsSDname {""}
protected

Definition at line 32 of file StoppedParticleFastSim.h.

◆ m_init

bool StoppedParticleFastSim::m_init {false}
protected

Definition at line 31 of file StoppedParticleFastSim.h.


The documentation for this class was generated from the following files:
StoppedParticleFastSim::m_fsSD
TrackFastSimSD * m_fsSD
Definition: StoppedParticleFastSim.h:30
TrackFastSimSD
Definition: TrackFastSimSD.h:24
mat
GeoMaterial * mat
Definition: LArDetectorConstructionTBEC.cxx:55
StoppedParticleFastSim::m_fsSDname
std::string m_fsSDname
Definition: StoppedParticleFastSim.h:32
StoppedParticleFastSim::m_init
bool m_init
Definition: StoppedParticleFastSim.h:31
lumiFormat.i
int i
Definition: lumiFormat.py:85
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:228
python.PhysicalConstants.c_light
float c_light
Definition: PhysicalConstants.py:63
TrackFastSimSD::WriteTrack
void WriteTrack(const G4Track *, const bool, const bool)
Definition: TrackFastSimSD.cxx:97
isSquarkLH
bool isSquarkLH(const T &p)
Definition: AtlasPID.h:402
isRHadron
bool isRHadron(const T &p)
Definition: AtlasPID.h:871
pow
constexpr int pow(int base, int exp) noexcept
Definition: ap_fixedTest.cxx:15
description
std::string description
glabal timer - how long have I taken so far?
Definition: hcg.cxx:88