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

#include <G4PSPassageCellCurrent.hh>

+ Inheritance diagram for G4PSPassageCellCurrent:

Public Member Functions

 G4PSPassageCellCurrent (G4String name, G4int depth=0)
 
 ~G4PSPassageCellCurrent () 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 (G4String name, G4int depth=0)
 
virtual ~G4VPrimitivePlotter ()
 
void Plot (G4int copyNo, G4int histID)
 
G4int GetNumberOfHist () const
 
- Public Member Functions inherited from G4VPrimitiveScorer
 G4VPrimitiveScorer (G4String name, G4int depth=0)
 
virtual ~G4VPrimitiveScorer ()
 
G4int GetCollectionID (G4int)
 
virtual void Initialize (G4HCofThisEvent *)
 
virtual void EndOfEvent (G4HCofThisEvent *)
 
virtual void clear ()
 
virtual void DrawAll ()
 
virtual void PrintAll ()
 
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
 
virtual G4bool IsPassed (G4Step *)
 
- Protected Member Functions inherited from G4VPrimitiveScorer
virtual G4bool ProcessHits (G4Step *, G4TouchableHistory *)=0
 
virtual G4int GetIndex (G4Step *)
 
void CheckAndSetUnit (const G4String &unit, const G4String &category)
 
G4VSolidComputeSolid (G4Step *aStep, G4int replicaIdx)
 
G4VSolidComputeCurrentSolid (G4Step *aStep)
 

Additional Inherited Members

- Protected Attributes inherited from G4VPrimitivePlotter
std::map< G4int, G4inthitIDMap
 
- Protected Attributes inherited from G4VPrimitiveScorer
G4String primitiveName
 
G4MultiFunctionalDetectordetector
 
G4VSDFilterfilter
 
G4int verboseLevel
 
G4int indexDepth
 
G4String unitName
 
G4double unitValue
 
G4int fNi
 
G4int fNj
 
G4int fNk
 

Detailed Description

Definition at line 49 of file G4PSPassageCellCurrent.hh.

Constructor & Destructor Documentation

◆ G4PSPassageCellCurrent()

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

Definition at line 50 of file G4PSPassageCellCurrent.cc.

51 : G4VPrimitivePlotter(name, depth)
52{
53 SetUnit("");
54}
virtual void SetUnit(const G4String &unit)

◆ ~G4PSPassageCellCurrent()

G4PSPassageCellCurrent::~G4PSPassageCellCurrent ( )
overridedefault

Member Function Documentation

◆ clear()

void G4PSPassageCellCurrent::clear ( )
overridevirtual

Reimplemented from G4VPrimitiveScorer.

Definition at line 129 of file G4PSPassageCellCurrent.cc.

129{ EvtMap->clear(); }
void clear()
Definition: G4THitsMap.hh:524

◆ Initialize()

void G4PSPassageCellCurrent::Initialize ( G4HCofThisEvent HCE)
overridevirtual

Reimplemented from G4VPrimitiveScorer.

Definition at line 119 of file G4PSPassageCellCurrent.cc.

120{
121 fCurrentTrkID = -1;
122
123 EvtMap = new G4THitsMap<G4double>(detector->GetName(), GetName());
124 if(HCID < 0)
125 HCID = GetCollectionID(0);
126 HCE->AddHitsCollection(HCID, EvtMap);
127}
void AddHitsCollection(G4int HCID, G4VHitsCollection *aHC)
G4String GetName() const
G4MultiFunctionalDetector * detector
G4int GetCollectionID(G4int)

◆ IsPassed()

G4bool G4PSPassageCellCurrent::IsPassed ( G4Step aStep)
protectedvirtual

Definition at line 86 of file G4PSPassageCellCurrent.cc.

87{
88 G4bool Passed = false;
89
90 G4bool IsEnter = aStep->GetPreStepPoint()->GetStepStatus() == fGeomBoundary;
91 G4bool IsExit = aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary;
92
93 G4int trkid = aStep->GetTrack()->GetTrackID();
94
95 if(IsEnter && IsExit)
96 { // Passed at one step
97 Passed = true;
98 }
99 else if(IsEnter)
100 { // Enter a new geometry
101 fCurrentTrkID = trkid; // Resetting the current track.
102 }
103 else if(IsExit)
104 { // Exit a current geometry
105 if(fCurrentTrkID == trkid)
106 {
107 Passed = true; // if the track is same as entered.
108 }
109 }
110 else
111 { // Inside geometry
112 if(fCurrentTrkID == trkid)
113 { // Adding the track length to current one ,
114 }
115 }
116 return Passed;
117}
@ fGeomBoundary
Definition: G4StepStatus.hh:43
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
G4StepStatus GetStepStatus() const
G4Track * GetTrack() const
G4StepPoint * GetPreStepPoint() const
G4StepPoint * GetPostStepPoint() const
G4int GetTrackID() const

Referenced by ProcessHits().

◆ PrintAll()

void G4PSPassageCellCurrent::PrintAll ( )
overridevirtual

Reimplemented from G4VPrimitiveScorer.

Definition at line 131 of file G4PSPassageCellCurrent.cc.

132{
133 G4cout << " MultiFunctionalDet " << detector->GetName() << G4endl;
134 G4cout << " PrimitiveScorer " << GetName() << G4endl;
135 G4cout << " Number of entries " << EvtMap->entries() << G4endl;
136 for(const auto& [copy, current] : *(EvtMap->GetMap()))
137 {
138 G4cout << " copy no.: " << copy
139 << " cell current : " << *(current) << " [tracks] " << G4endl;
140 }
141}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
Map_t * GetMap() const
Definition: G4THitsMap.hh:155
size_t entries() const
Definition: G4THitsMap.hh:158
void copy(G4double dst[], const G4double src[], std::size_t size=G4FieldTrack::ncompSVEC)
Definition: G4FieldUtils.cc:98

◆ ProcessHits()

G4bool G4PSPassageCellCurrent::ProcessHits ( G4Step aStep,
G4TouchableHistory  
)
overrideprotectedvirtual

Implements G4VPrimitiveScorer.

Definition at line 56 of file G4PSPassageCellCurrent.cc.

57{
58 if(IsPassed(aStep))
59 {
60 fCurrent = 1.;
61 if(weighted)
62 fCurrent = aStep->GetPreStepPoint()->GetWeight();
63 G4int index = GetIndex(aStep);
64 EvtMap->add(index, fCurrent);
65
66 if(!hitIDMap.empty() && hitIDMap.find(index) != hitIDMap.cend())
67 {
68 auto filler = G4VScoreHistFiller::Instance();
69 if(filler == nullptr)
70 {
72 "G4PSVolumeFlux::ProcessHits", "SCORER0123", JustWarning,
73 "G4TScoreHistFiller is not instantiated!! Histogram is not filled.");
74 }
75 else
76 {
77 filler->FillH1(hitIDMap[index],
78 aStep->GetPreStepPoint()->GetKineticEnergy(), fCurrent);
79 }
80 }
81 }
82
83 return true;
84}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
virtual G4bool IsPassed(G4Step *)
G4double GetKineticEnergy() const
G4double GetWeight() const
std::map< G4int, G4int > hitIDMap
virtual G4int GetIndex(G4Step *)
static G4VScoreHistFiller * Instance()
size_t add(const G4int &key, U *&aHit) const
Definition: G4THitsMap.hh:177

◆ SetUnit()

void G4PSPassageCellCurrent::SetUnit ( const G4String unit)
virtual

Definition at line 143 of file G4PSPassageCellCurrent.cc.

144{
145 if(unit.empty())
146 {
147 unitName = unit;
148 unitValue = 1.0;
149 }
150 else
151 {
152 G4String msg = "Invalid unit [" + unit + "] (Current unit is [" +
153 GetUnit() + "] ) for " + GetName();
154 G4Exception("G4PSPassageCellCurrent::SetUnit", "DetPS0012", JustWarning,
155 msg);
156 }
157}
const G4String & GetUnit() const

Referenced by G4PSPassageCellCurrent().

◆ Weighted()

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

Definition at line 55 of file G4PSPassageCellCurrent.hh.

55{ weighted = flg; }

Referenced by G4ScoreQuantityMessenger::SetNewValue().


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