CMS 3D CMS Logo
CMSSteppingVerbose.cc
Go to the documentation of this file.
1 
2 #include "SimG4Core/Notification/interface/CMSSteppingVerbose.h"
3 #include "G4Event.hh"
4 #include "G4EventManager.hh"
5 #include "G4Track.hh"
6 #include "G4TrackStatus.hh"
7 #include "G4Step.hh"
8 #include "G4SteppingManager.hh"
9 #include "G4SteppingVerbose.hh"
10 #include "G4ParticleDefinition.hh"
11 #include "G4VProcess.hh"
12 #include "G4SystemOfUnits.hh"
13 
14 CMSSteppingVerbose::CMSSteppingVerbose(
15  G4int verb, G4double ekin, std::vector<G4int>& evtNum, std::vector<G4int>& primV, std::vector<G4int>& trNum)
16  : m_PrintEvent(false),
17  m_PrintTrack(false),
18  m_smInitialized(false),
19  m_verbose(verb),
20  m_EventNumbers(evtNum),
21  m_PrimaryVertex(primV),
22  m_TrackNumbers(trNum),
23  m_EkinThreshold(ekin) {
24  m_nEvents = m_EventNumbers.size();
25  m_nVertex = m_PrimaryVertex.size();
26  m_nTracks = m_TrackNumbers.size();
27  m_g4SteppingVerbose = new G4SteppingVerbose();
28  G4VSteppingVerbose::SetInstance(m_g4SteppingVerbose);
29  m_g4SteppingVerbose->SetSilent(1);
30 }
31 
32 CMSSteppingVerbose::~CMSSteppingVerbose() {}
33 
34 void CMSSteppingVerbose::BeginOfEvent(const G4Event* evt) {
35  m_PrintEvent = false;
36  if (0 >= m_verbose) {
37  return;
38  }
39  if (m_nEvents == 0) {
40  m_PrintEvent = true;
41  } else {
42  for (G4int i = 0; i < m_nEvents; ++i) {
43  // check event number
44  if (evt->GetEventID() == m_EventNumbers[i]) {
45  // check number of vertex
46  if (m_nVertex == m_nEvents && evt->GetNumberOfPrimaryVertex() != m_PrimaryVertex[i]) {
47  continue;
48  }
49  m_PrintEvent = true;
50  break;
51  }
52  }
53  }
54  if (!m_PrintEvent) {
55  return;
56  }
57  G4cout << "========== Event #" << evt->GetEventID() << " " << evt->GetNumberOfPrimaryVertex()
58  << " primary vertexes ======" << G4endl;
59  G4cout << G4endl;
60 }
61 
62 void CMSSteppingVerbose::TrackStarted(const G4Track* track, bool isKilled) {
63  m_PrintTrack = false;
64  if (!m_PrintEvent) {
65  return;
66  }
67 
68  if (!m_smInitialized) {
69  G4SteppingManager* stepman = G4EventManager::GetEventManager()->GetTrackingManager()->GetSteppingManager();
70  m_g4SteppingVerbose->SetManager(stepman);
71  stepman->SetVerboseLevel(m_verbose);
72  m_smInitialized = true;
73  }
74 
75  if (m_nTracks == 0) {
76  if (track->GetKineticEnergy() >= m_EkinThreshold) {
77  m_PrintTrack = true;
78  }
79 
80  } else {
81  for (G4int i = 0; i < m_nTracks; ++i) {
82  if (track->GetTrackID() == m_TrackNumbers[i]) {
83  m_PrintTrack = true;
84  break;
85  }
86  }
87  }
88  if (!m_PrintTrack) {
89  return;
90  }
91 
92  G4cout << "*********************************************************************************************************"
93  << G4endl;
94  const G4ParticleDefinition* pd = track->GetDefinition();
95  G4cout << "* G4Track Information: Particle = ";
96  if (pd) {
97  G4cout << pd->GetParticleName();
98  }
99  G4cout << ", Track ID = " << track->GetTrackID() << ", Parent ID = " << track->GetParentID() << G4endl;
100  G4cout << "*********************************************************************************************************"
101  << G4endl;
102 
103  G4cout << std::setw(5) << "Step#"
104  << " " << std::setw(8) << "X(cm)"
105  << " " << std::setw(8) << "Y(cm)"
106  << " " << std::setw(8) << "Z(cm)"
107  << " " << std::setw(9) << "KinE(GeV)"
108  << " " << std::setw(8) << "dE(MeV)"
109  << " " << std::setw(8) << "Step(mm)"
110  << " " << std::setw(9) << "TrackL(cm)"
111  << " " << std::setw(30) << "PhysVolume"
112  << " " << std::setw(8) << "ProcName" << G4endl;
113 
114  G4int prec = G4cout.precision(4);
115 
116  G4cout << std::setw(5) << track->GetCurrentStepNumber() << " " << std::setw(8) << track->GetPosition().x() / CLHEP::cm
117  << " " << std::setw(8) << track->GetPosition().y() / CLHEP::cm << " " << std::setw(8)
118  << track->GetPosition().z() / CLHEP::cm << " " << std::setw(9) << track->GetKineticEnergy() / CLHEP::GeV << " "
119  << std::setw(8) << " " << std::setw(8) << " " << std::setw(9) << " ";
120  if (track->GetVolume() != nullptr) {
121  G4cout << std::setw(30) << track->GetVolume()->GetName() << " ";
122  }
123  if (isKilled) {
124  G4cout << "isKilled";
125  }
126  G4cout << G4endl;
127  G4cout.precision(prec);
128 }
129 
130 void CMSSteppingVerbose::TrackEnded(const G4Track* track) const {
131  if (!m_PrintTrack || 1 == m_verbose) {
132  return;
133  }
134 }
135 
136 void CMSSteppingVerbose::StackFilled(const G4Track* track, bool isKilled) const {
137  if (2 >= m_verbose || !m_PrintTrack || track->GetKineticEnergy() < m_EkinThreshold) {
138  return;
139  }
140  G4int prec = G4cout.precision(4);
141 
142  G4cout << std::setw(10) << track->GetTrackID() << " " << std::setw(8) << track->GetPosition().x() / CLHEP::cm << " "
143  << std::setw(8) << track->GetPosition().y() / CLHEP::cm << " " << std::setw(8)
144  << track->GetPosition().z() / CLHEP::cm << " " << std::setw(9) << track->GetKineticEnergy() / CLHEP::GeV
145  << " ";
146  if (track->GetVolume() != nullptr) {
147  G4cout << std::setw(24) << track->GetVolume()->GetName() << " ";
148  }
149  if (isKilled) {
150  G4cout << "isKilled";
151  }
152  G4cout << G4endl;
153  G4cout.precision(prec);
154 }
155 
156 void CMSSteppingVerbose::NextStep(const G4Step* step, const G4SteppingManager* sManager, bool isKilled) const {
157  if (!m_PrintTrack) {
158  return;
159  }
160 
161  G4int prec;
162  const G4Track* track = step->GetTrack();
163  const G4StepPoint* preStep = step->GetPreStepPoint();
164  const G4StepPoint* postStep = step->GetPostStepPoint();
165 
166  if (3 <= m_verbose) {
167  m_g4SteppingVerbose->SetSilent(0);
168  m_g4SteppingVerbose->DPSLStarted();
169  m_g4SteppingVerbose->DPSLAlongStep();
170  m_g4SteppingVerbose->DPSLPostStep();
171  if (4 <= m_verbose) {
172  m_g4SteppingVerbose->AlongStepDoItAllDone();
173  m_g4SteppingVerbose->PostStepDoItAllDone();
174  }
175  m_g4SteppingVerbose->SetSilent(1);
176 
177  prec = G4cout.precision(16);
178 
179  G4cout << G4endl;
180  G4cout << " ++G4Step Information " << G4endl;
181  G4cout << " Address of G4Track : " << track << G4endl;
182  G4cout << " Step Length (mm) : " << track->GetStepLength() << G4endl;
183  G4cout << " Energy Deposit (MeV) : " << step->GetTotalEnergyDeposit() << G4endl;
184 
185  G4cout << " -------------------------------------------------------"
186  << "-------------------------------" << G4endl;
187  G4cout << " StepPoint Information " << std::setw(30) << "PreStep" << std::setw(30) << "PostStep" << G4endl;
188  G4cout << " -------------------------------------------------------"
189  << "-------------------------------" << G4endl;
190  G4cout << " Position - x (cm) : " << std::setw(30) << preStep->GetPosition().x() / CLHEP::cm << std::setw(30)
191  << postStep->GetPosition().x() / CLHEP::cm << G4endl;
192  G4cout << " Position - y (cm) : " << std::setw(30) << preStep->GetPosition().y() / CLHEP::cm << std::setw(30)
193  << postStep->GetPosition().y() / CLHEP::cm << G4endl;
194  G4cout << " Position - z (cm) : " << std::setw(30) << preStep->GetPosition().z() / CLHEP::cm << std::setw(30)
195  << postStep->GetPosition().z() / CLHEP::cm << G4endl;
196  G4cout << " Global Time (ns) : " << std::setw(30) << preStep->GetGlobalTime() / CLHEP::ns << std::setw(30)
197  << postStep->GetGlobalTime() / CLHEP::ns << G4endl;
198  G4cout << " Local Time (ns) : " << std::setw(30) << preStep->GetLocalTime() / CLHEP::ns << std::setw(30)
199  << postStep->GetLocalTime() / CLHEP::ns << G4endl;
200  G4cout << " Proper Time (ns) : " << std::setw(30) << preStep->GetProperTime() / CLHEP::ns << std::setw(30)
201  << postStep->GetProperTime() / CLHEP::ns << G4endl;
202  G4cout << " Momentum Direct - x : " << std::setw(30) << preStep->GetMomentumDirection().x() << std::setw(30)
203  << postStep->GetMomentumDirection().x() << G4endl;
204  G4cout << " Momentum Direct - y : " << std::setw(30) << preStep->GetMomentumDirection().y() << std::setw(30)
205  << postStep->GetMomentumDirection().y() << G4endl;
206  G4cout << " Momentum Direct - z : " << std::setw(30) << preStep->GetMomentumDirection().z() << std::setw(30)
207  << postStep->GetMomentumDirection().z() << G4endl;
208  G4cout << " Momentum - x (GeV/c): " << std::setw(30) << preStep->GetMomentum().x() / CLHEP::GeV
209  << std::setw(30) << postStep->GetMomentum().x() / CLHEP::GeV << G4endl;
210  G4cout << " Momentum - y (GeV/c): " << std::setw(30) << preStep->GetMomentum().y() / CLHEP::GeV
211  << std::setw(30) << postStep->GetMomentum().y() / CLHEP::GeV << G4endl;
212  G4cout << " Momentum - z (GeV/c): " << std::setw(30) << preStep->GetMomentum().z() / CLHEP::GeV
213  << std::setw(30) << postStep->GetMomentum().z() / CLHEP::GeV << G4endl;
214  G4cout << " Total Energy (GeV) : " << std::setw(30) << preStep->GetTotalEnergy() / CLHEP::GeV << std::setw(30)
215  << postStep->GetTotalEnergy() / CLHEP::GeV << G4endl;
216  G4cout << " Kinetic Energy (GeV): " << std::setw(30) << preStep->GetKineticEnergy() / CLHEP::GeV
217  << std::setw(30) << postStep->GetKineticEnergy() / CLHEP::GeV << G4endl;
218  G4cout << " Velocity (mm/ns) : " << std::setw(30) << preStep->GetVelocity() << std::setw(30)
219  << postStep->GetVelocity() << G4endl;
220  G4cout << " Volume Name : " << std::setw(30) << preStep->GetPhysicalVolume()->GetName();
221  G4String volName = (postStep->GetPhysicalVolume()) ? postStep->GetPhysicalVolume()->GetName() : "OutOfWorld";
222 
223  G4cout << std::setw(30) << volName << G4endl;
224  G4cout << " Safety (mm) : " << std::setw(30) << preStep->GetSafety() << std::setw(30)
225  << postStep->GetSafety() << G4endl;
226  G4cout << " Polarization - x : " << std::setw(30) << preStep->GetPolarization().x() << std::setw(30)
227  << postStep->GetPolarization().x() << G4endl;
228  G4cout << " Polarization - y : " << std::setw(30) << preStep->GetPolarization().y() << std::setw(30)
229  << postStep->GetPolarization().y() << G4endl;
230  G4cout << " Polarization - Z : " << std::setw(30) << preStep->GetPolarization().z() << std::setw(30)
231  << postStep->GetPolarization().z() << G4endl;
232  G4cout << " Weight : " << std::setw(30) << preStep->GetWeight() << std::setw(30)
233  << postStep->GetWeight() << G4endl;
234  G4cout << " Step Status : ";
235  G4StepStatus tStepStatus = preStep->GetStepStatus();
236  if (tStepStatus == fGeomBoundary) {
237  G4cout << std::setw(30) << "Geom Limit";
238  } else if (tStepStatus == fAlongStepDoItProc) {
239  G4cout << std::setw(30) << "AlongStep Proc.";
240  } else if (tStepStatus == fPostStepDoItProc) {
241  G4cout << std::setw(30) << "PostStep Proc";
242  } else if (tStepStatus == fAtRestDoItProc) {
243  G4cout << std::setw(30) << "AtRest Proc";
244  } else if (tStepStatus == fUndefined) {
245  G4cout << std::setw(30) << "Undefined";
246  }
247 
248  tStepStatus = postStep->GetStepStatus();
249  if (tStepStatus == fGeomBoundary) {
250  G4cout << std::setw(30) << "Geom Limit";
251  } else if (tStepStatus == fAlongStepDoItProc) {
252  G4cout << std::setw(30) << "AlongStep Proc.";
253  } else if (tStepStatus == fPostStepDoItProc) {
254  G4cout << std::setw(30) << "PostStep Proc";
255  } else if (tStepStatus == fAtRestDoItProc) {
256  G4cout << std::setw(30) << "AtRest Proc";
257  } else if (tStepStatus == fUndefined) {
258  G4cout << std::setw(30) << "Undefined";
259  }
260 
261  G4cout << G4endl;
262  G4cout << " Process defined Step: ";
263  if (preStep->GetProcessDefinedStep() == nullptr) {
264  G4cout << std::setw(30) << "Undefined";
265  } else {
266  G4cout << std::setw(30) << preStep->GetProcessDefinedStep()->GetProcessName();
267  }
268  if (postStep->GetProcessDefinedStep() == nullptr) {
269  G4cout << std::setw(30) << "Undefined";
270  } else {
271  G4cout << std::setw(30) << postStep->GetProcessDefinedStep()->GetProcessName();
272  }
273  G4cout.precision(prec);
274 
275  G4cout << G4endl;
276  G4cout << " -------------------------------------------------------"
277  << "-------------------------------" << G4endl;
278  }
279 
280  // geometry information
281  if (4 <= m_verbose) {
282  const G4VTouchable* tch1 = preStep->GetTouchable();
283  const G4VTouchable* tch2 = postStep->GetTouchable();
284 
285  G4double x = postStep->GetPosition().x();
286  G4double y = postStep->GetPosition().y();
287  G4cout << "Touchable depth pre= " << tch1->GetHistoryDepth() << " post= " << tch2->GetHistoryDepth()
288  << " trans1= " << tch1->GetTranslation(tch1->GetHistoryDepth())
289  << " trans2= " << tch2->GetTranslation(tch2->GetHistoryDepth()) << " r= " << std::sqrt(x * x + y * y)
290  << G4endl;
291  const G4VPhysicalVolume* pv1 = preStep->GetPhysicalVolume();
292  const G4VPhysicalVolume* pv2 = postStep->GetPhysicalVolume();
293  const G4RotationMatrix* rotm = pv1->GetFrameRotation();
294  G4cout << "PreStepVolume: " << pv1->GetName() << G4endl;
295  G4cout << " Translation: " << pv1->GetObjectTranslation() << G4endl;
296  if (nullptr != rotm)
297  G4cout << " Rotation: " << *rotm << G4endl;
298  const G4VSolid* sv1 = pv1->GetLogicalVolume()->GetSolid();
299  sv1->StreamInfo(G4cout);
300  G4cout << G4endl;
301  if (pv2 && pv2 != pv1) {
302  G4cout << "PostStepVolume: " << pv2->GetName() << G4endl;
303  G4cout << " Translation: " << pv2->GetObjectTranslation() << G4endl;
304  rotm = pv2->GetFrameRotation();
305  if (nullptr != rotm)
306  G4cout << " Rotation: " << *rotm << G4endl;
307  const G4VSolid* sv2 = pv2->GetLogicalVolume()->GetSolid();
308  sv2->StreamInfo(G4cout);
309  }
310  G4cout << G4endl;
311 
312  if (5 <= m_verbose) {
313  G4cout << "##### Geometry Depth Analysis" << G4endl;
314  for (G4int k = 1; k < tch1->GetHistoryDepth(); ++k) {
315  const G4VPhysicalVolume* pv = tch1->GetVolume(k);
316  if (pv) {
317  const G4VSolid* sol = pv->GetLogicalVolume()->GetSolid();
318  G4cout << " Depth # " << k << " PhysVolume " << pv->GetName() << G4endl;
319  G4cout << " Translation: " << pv->GetObjectTranslation() << G4endl;
320  const G4RotationMatrix* rotm = pv->GetFrameRotation();
321  if (nullptr != rotm)
322  G4cout << " Rotation: " << *rotm << G4endl;
323  sol->StreamInfo(G4cout);
324  }
325  }
326  }
327  G4cout << G4endl;
328  }
329 
330  // verbose == 1
331  prec = G4cout.precision(4);
332  G4cout << std::setw(5) << track->GetCurrentStepNumber() << " " << std::setw(8) << track->GetPosition().x() / CLHEP::cm
333  << " " << std::setw(8) << track->GetPosition().y() / CLHEP::cm << " " << std::setw(8)
334  << track->GetPosition().z() / CLHEP::cm << " " << std::setw(9) << track->GetKineticEnergy() / CLHEP::GeV
335  << " ";
336  G4cout.precision(prec);
337  prec = G4cout.precision(3);
338  G4cout << std::setw(8) << step->GetTotalEnergyDeposit() / CLHEP::MeV << " " << std::setw(8)
339  << step->GetStepLength() / CLHEP::mm << " " << std::setw(9) << track->GetTrackLength() / CLHEP::cm << " ";
340 
341  G4bool endTracking = false;
342  if (track->GetNextVolume() != nullptr) {
343  G4cout << std::setw(30) << track->GetVolume()->GetName() << " ";
344  } else {
345  G4cout << std::setw(30) << "OutOfWorld"
346  << " ";
347  endTracking = true;
348  }
349  if (isKilled) {
350  G4cout << "isKilled";
351  endTracking = true;
352  } else if (postStep->GetProcessDefinedStep() != nullptr) {
353  G4cout << postStep->GetProcessDefinedStep()->GetProcessName();
354  }
355  G4cout << G4endl;
356  G4cout.precision(prec);
357 
358  if (1 >= m_verbose) {
359  return;
360  }
361  // verbose > 1
362  if (!endTracking && fStopAndKill != track->GetTrackStatus()) {
363  return;
364  }
365 
366  prec = G4cout.precision(4);
367  const G4TrackVector* tv = step->GetSecondary();
368  G4int nt = tv->size();
369  if (nt > 0) {
370  G4cout << " ++List of " << nt << " secondaries generated; "
371  << " Ekin > " << m_EkinThreshold << " MeV are shown:" << G4endl;
372  }
373  for (G4int i = 0; i < nt; ++i) {
374  if ((*tv)[i]->GetKineticEnergy() < m_EkinThreshold) {
375  continue;
376  }
377  G4cout << " (" << std::setw(9) << (*tv)[i]->GetPosition().x() / CLHEP::cm << " " << std::setw(9)
378  << (*tv)[i]->GetPosition().y() / CLHEP::cm << " " << std::setw(9) << (*tv)[i]->GetPosition().z() / CLHEP::cm
379  << ") cm, " << std::setw(9) << (*tv)[i]->GetKineticEnergy() / CLHEP::GeV << " GeV, " << std::setw(9)
380  << (*tv)[i]->GetGlobalTime() / CLHEP::ns << " ns, " << std::setw(18)
381  << (*tv)[i]->GetDefinition()->GetParticleName() << G4endl;
382  }
383  G4cout.precision(prec);
384 }
mps_fire.i
i
Definition: mps_fire.py:429
CMSSteppingVerbose::m_EventNumbers
std::vector< G4int > m_EventNumbers
Definition: CMSSteppingVerbose.h:45
ecalTB2006H4_GenSimDigiReco_cfg.G4cout
G4cout
Definition: ecalTB2006H4_GenSimDigiReco_cfg.py:285
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
CMSSteppingVerbose::m_TrackNumbers
std::vector< G4int > m_TrackNumbers
Definition: CMSSteppingVerbose.h:47
CMSSteppingVerbose::CMSSteppingVerbose
CMSSteppingVerbose(G4int verb, G4double ekin, std::vector< G4int > &evtNum, std::vector< G4int > &primV, std::vector< G4int > &trNum)
Definition: CMSSteppingVerbose.cc:14
nt
int nt
Definition: AMPTWrapper.h:42
CMSSteppingVerbose::m_g4SteppingVerbose
G4SteppingVerbose * m_g4SteppingVerbose
Definition: CMSSteppingVerbose.h:49
CMSSteppingVerbose::BeginOfEvent
void BeginOfEvent(const G4Event *)
Definition: CMSSteppingVerbose.cc:34
mkfit::Const::sol
constexpr float sol
Definition: Config.h:48
CMSSteppingVerbose::m_PrimaryVertex
std::vector< G4int > m_PrimaryVertex
Definition: CMSSteppingVerbose.h:46
HLT_2022v15_cff.track
track
Definition: HLT_2022v15_cff.py:9644
CMSSteppingVerbose::StackFilled
void StackFilled(const G4Track *, bool isKilled) const
Definition: CMSSteppingVerbose.cc:136
CMSSteppingVerbose::TrackEnded
void TrackEnded(const G4Track *) const
Definition: CMSSteppingVerbose.cc:130
CMSSteppingVerbose::TrackStarted
void TrackStarted(const G4Track *, bool isKilled)
Definition: CMSSteppingVerbose.cc:62
step
step
Definition: StallMonitor.cc:98
CMSSteppingVerbose::NextStep
void NextStep(const G4Step *, const G4SteppingManager *ptr, bool isKilled) const
Definition: CMSSteppingVerbose.cc:156
Generated for CMSSW Reference Manual by doxygen 1.8.14