Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4PSPassageTrackLength Class Reference

#include <G4PSPassageTrackLength.hh>

+ Inheritance diagram for G4PSPassageTrackLength:

Public Member Functions

 G4PSPassageTrackLength (G4String name, G4int depth=0)
 
 G4PSPassageTrackLength (G4String name, const G4String &unit, G4int depth=0)
 
 ~G4PSPassageTrackLength () override=default
 
void Weighted (G4bool flg=true)
 
void Initialize (G4HCofThisEvent *) override
 
void clear () override
 
void PrintAll () override
 
virtual void SetUnit (const G4String &unit)
 
- Public Member Functions inherited from G4VPrimitivePlotter
 ~G4VPrimitivePlotter () override=default
 
void Plot (G4int copyNo, G4int histID)
 
G4int GetNumberOfHist () const
 
 G4VPrimitiveScorer (G4String name, G4int depth=0)
 
- Public Member Functions inherited from G4VPrimitiveScorer
 G4VPrimitiveScorer (G4String name, G4int depth=0)
 
virtual ~G4VPrimitiveScorer ()=default
 
G4int GetCollectionID (G4int)
 
virtual void EndOfEvent (G4HCofThisEvent *)
 
virtual void DrawAll ()
 
void SetUnit (const G4String &unit)
 
const G4StringGetUnit () const
 
G4double GetUnitValue () const
 
void SetMultiFunctionalDetector (G4MultiFunctionalDetector *d)
 
G4MultiFunctionalDetectorGetMultiFunctionalDetector () const
 
G4String GetName () const
 
void SetFilter (G4VSDFilter *f)
 
G4VSDFilterGetFilter () const
 
void SetVerboseLevel (G4int vl)
 
G4int GetVerboseLevel () const
 
void SetNijk (G4int i, G4int j, G4int k)
 

Protected Member Functions

G4bool ProcessHits (G4Step *, G4TouchableHistory *) override
 
G4bool IsPassed (G4Step *)
 
- Protected Member Functions inherited from G4VPrimitiveScorer
G4VSolidComputeSolid (G4Step *aStep, G4int replicaIdx)
 
G4VSolidComputeCurrentSolid (G4Step *aStep)
 
virtual G4int GetIndex (G4Step *)
 
void CheckAndSetUnit (const G4String &unit, const G4String &category)
 

Additional Inherited Members

- Protected Attributes inherited from G4VPrimitivePlotter
std::map< G4int, G4inthitIDMap
 
- Protected Attributes inherited from G4VPrimitiveScorer
G4String primitiveName
 
G4MultiFunctionalDetectordetector {nullptr}
 
G4VSDFilterfilter {nullptr}
 
G4int verboseLevel {0}
 
G4int indexDepth
 
G4String unitName {"NoUnit"}
 
G4double unitValue {1.0}
 
G4int fNi {0}
 
G4int fNj {0}
 
G4int fNk {0}
 

Detailed Description

Definition at line 48 of file G4PSPassageTrackLength.hh.

Constructor & Destructor Documentation

◆ G4PSPassageTrackLength() [1/2]

G4PSPassageTrackLength::G4PSPassageTrackLength ( G4String name,
G4int depth = 0 )

Definition at line 48 of file G4PSPassageTrackLength.cc.

49 : G4PSPassageTrackLength(name, "mm", depth)
50{}
G4PSPassageTrackLength(G4String name, G4int depth=0)

◆ G4PSPassageTrackLength() [2/2]

G4PSPassageTrackLength::G4PSPassageTrackLength ( G4String name,
const G4String & unit,
G4int depth = 0 )

Definition at line 52 of file G4PSPassageTrackLength.cc.

55 : G4VPrimitivePlotter(name, depth)
56 , HCID(-1)
57 , fCurrentTrkID(-1)
58 , fTrackLength(0.)
59 , EvtMap(nullptr)
60 , weighted(false)
61{
62 SetUnit(unit);
63}
virtual void SetUnit(const G4String &unit)

◆ ~G4PSPassageTrackLength()

G4PSPassageTrackLength::~G4PSPassageTrackLength ( )
overridedefault

Member Function Documentation

◆ clear()

void G4PSPassageTrackLength::clear ( )
overridevirtual

Reimplemented from G4VPrimitiveScorer.

Definition at line 142 of file G4PSPassageTrackLength.cc.

142{ EvtMap->clear(); }

◆ Initialize()

void G4PSPassageTrackLength::Initialize ( G4HCofThisEvent * HCE)
overridevirtual

Reimplemented from G4VPrimitiveScorer.

Definition at line 132 of file G4PSPassageTrackLength.cc.

133{
134 fCurrentTrkID = -1;
135
136 EvtMap = new G4THitsMap<G4double>(detector->GetName(), GetName());
137 if(HCID < 0)
138 HCID = GetCollectionID(0);
139 HCE->AddHitsCollection(HCID, EvtMap);
140}
void AddHitsCollection(G4int HCID, G4VHitsCollection *aHC)
G4String GetName() const
G4MultiFunctionalDetector * detector

◆ IsPassed()

G4bool G4PSPassageTrackLength::IsPassed ( G4Step * aStep)
protected

Definition at line 91 of file G4PSPassageTrackLength.cc.

92{
93 G4bool Passed = false;
94
95 G4bool IsEnter = aStep->GetPreStepPoint()->GetStepStatus() == fGeomBoundary;
96 G4bool IsExit = aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary;
97
98 G4int trkid = aStep->GetTrack()->GetTrackID();
99 G4double trklength = aStep->GetStepLength();
100 if(weighted)
101 trklength *= aStep->GetPreStepPoint()->GetWeight();
102
103 if(IsEnter && IsExit)
104 { // Passed at one step
105 fTrackLength = trklength; // Track length is absolutely given.
106 Passed = true;
107 }
108 else if(IsEnter)
109 { // Enter a new geometry
110 fCurrentTrkID = trkid; // Resetting the current track.
111 fTrackLength = trklength;
112 }
113 else if(IsExit)
114 { // Exit a current geometry
115 if(fCurrentTrkID == trkid)
116 {
117 fTrackLength += trklength; // Adding the track length to current one,
118 Passed = true; // if the track is same as entered.
119 }
120 }
121 else
122 { // Inside geometry
123 if(fCurrentTrkID == trkid)
124 { // Adding the track length to current one ,
125 fTrackLength += trklength; // if the track is same as entered.
126 }
127 }
128
129 return Passed;
130}
@ fGeomBoundary
double G4double
Definition G4Types.hh:83
bool G4bool
Definition G4Types.hh:86
int G4int
Definition G4Types.hh:85
G4StepStatus GetStepStatus() const
G4double GetWeight() const
G4Track * GetTrack() const
G4StepPoint * GetPreStepPoint() const
G4double GetStepLength() const
G4StepPoint * GetPostStepPoint() const
G4int GetTrackID() const

Referenced by ProcessHits().

◆ PrintAll()

void G4PSPassageTrackLength::PrintAll ( )
overridevirtual

Reimplemented from G4VPrimitiveScorer.

Definition at line 144 of file G4PSPassageTrackLength.cc.

145{
146 G4cout << " MultiFunctionalDet " << detector->GetName() << G4endl;
147 G4cout << " PrimitiveSenstivity " << GetName() << G4endl;
148 G4cout << " Number of entries " << EvtMap->entries() << G4endl;
149 for(const auto& [copy, length] : *(EvtMap->GetMap()))
150 {
151 G4cout << " copy no.: " << copy
152 << " track length : " << *(length) / GetUnitValue() << " ["
153 << GetUnit() << "]" << G4endl;
154 }
155}
#define G4endl
Definition G4ios.hh:67
G4GLOB_DLL std::ostream G4cout
const G4String & GetUnit() const
G4double GetUnitValue() const
Map_t * GetMap() const
size_t entries() const
void copy(G4double dst[], const G4double src[], std::size_t size=G4FieldTrack::ncompSVEC)

◆ ProcessHits()

G4bool G4PSPassageTrackLength::ProcessHits ( G4Step * aStep,
G4TouchableHistory *  )
overrideprotectedvirtual

Implements G4VPrimitiveScorer.

Definition at line 65 of file G4PSPassageTrackLength.cc.

66{
67 if(IsPassed(aStep))
68 {
69 G4int index = GetIndex(aStep);
70 EvtMap->add(index, fTrackLength);
71
72 if(!hitIDMap.empty() && hitIDMap.find(index) != hitIDMap.cend())
73 {
74 auto filler = G4VScoreHistFiller::Instance();
75 if(filler == nullptr)
76 {
78 "G4PSPassageTrackLength::ProcessHits", "SCORER0123", JustWarning,
79 "G4TScoreHistFiller is not instantiated!! Histogram is not filled.");
80 }
81 else
82 {
83 filler->FillH1(hitIDMap[index], fTrackLength, 1.);
84 }
85 }
86 }
87
88 return true;
89}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
std::map< G4int, G4int > hitIDMap
virtual G4int GetIndex(G4Step *)
static G4VScoreHistFiller * Instance()
size_t add(const G4int &key, U *&aHit) const

◆ SetUnit()

void G4PSPassageTrackLength::SetUnit ( const G4String & unit)
virtual

Definition at line 157 of file G4PSPassageTrackLength.cc.

158{
159 CheckAndSetUnit(unit, "Length");
160}
void CheckAndSetUnit(const G4String &unit, const G4String &category)

Referenced by G4PSPassageTrackLength(), and G4PSPassageTrackLength3D::G4PSPassageTrackLength3D().

◆ Weighted()

void G4PSPassageTrackLength::Weighted ( G4bool flg = true)
inline

Definition at line 55 of file G4PSPassageTrackLength.hh.

55{ weighted = flg; }

Referenced by G4ScoreQuantityMessenger::SetNewValue().


The documentation for this class was generated from the following files: