ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
TSRunAction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file TSRunAction.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
28 //
29 //
30 //
31 //
32 //
33 //
34 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
35 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
36 
37 
38 #include "TSRunAction.hh"
40 
41 #include "G4RunManager.hh"
42 #include "G4Run.hh"
43 #include "G4Timer.hh"
44 
45 #include "TSRun.hh"
47 #include "G4StatAnalysis.hh"
48 
49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
50 
52 : fDetector(TSDetectorConstruction::Instance()),
53  fName(fDetector->GetMFDName())
54 { }
55 
56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
57 
59 { }
60 
61 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
62 
64 {
65  return new TSRun(fName);
66 }
67 
68 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
69 
71 {
72  G4int evts_to_process = aRun->GetNumberOfEventToBeProcessed();
73  G4RunManager::GetRunManager()->SetPrintProgress((evts_to_process > 100)
74  ? evts_to_process/100
75  : 1);
76  if(IsMaster())
77  G4PrintEnv();
78 }
79 
80 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
81 
83 {
84 
85  if(IsMaster())
86  {
87  G4cout << " ###### EndOfTSRunAction ###### " << G4endl;
88 
89  aRun->GetNumberOfEvent();
90  std::ofstream fileout;
91  G4String fname = "";
92  std::stringstream separator;
93 
94  separator
95  << "============================================================";
96 
97  typedef std::set<G4int> IDSet_t;
98  IDSet_t IDs;
99 
100  //- TSRun object.
101  const TSRun* tsRun = static_cast<const TSRun*>(aRun);
102  //--- Dump all scored quantities involved in TSRun.
103 
104  //---------------------------------------------
105  // Dump accumulated quantities for this RUN.
106  //---------------------------------------------
107  std::vector<G4String> primScorerNames { "EnergyDeposit",
108  "NumberOfSteps" };
109  std::vector<G4String> fnames { "mfd_tl", "mfd_tg" };
110  std::vector<G4double> units { CLHEP::eV, CLHEP::keV, 1, 1 };
111  std::vector<G4String> unitstr { "keV", "steps" };
112 
113  //----------------------------------------------------------------------//
114  // lambda to print double value
115  auto print = [] (std::ostream& fout,
116  G4int first, G4double second,
117  G4double unit1, G4double unit2, G4String unit2str)
118  {
119  if(fout)
120  fout << first
121  << " " << second/unit1
122  << G4endl;
123 
124  G4cout
125  << " " << std::setw(10) << first
126  << " " << std::setw(15) << std::setprecision(6)
127  << std::fixed << second/unit2 << " " << unit2str
128  << G4endl;
129  G4cout.unsetf(std::ios::fixed);
130  };
131  //----------------------------------------------------------------------//
132  // lambda to print statistics
133  auto stat_print = [] (std::ostream& fout,
134  G4int first, G4StatAnalysis* stat, G4ConvergenceTester* conv,
135  G4double unit1, G4double unit2, G4String unit2str)
136  {
137  if(!stat || !conv)
138  return;
139  auto fsecond = (*stat);
140  auto psecond = (*stat);
141  fsecond /= unit1;
142  psecond /= unit2;
143  if(fout)
144  {
145  fout << first << " " << fsecond << G4endl;
146  conv->ShowResult(fout);
147  }
148  std::stringstream ss;
149  ss << " " << std::setw(10) << first
150  << " " << std::setw(15) << std::setprecision(6)
151  << std::fixed << psecond << " " << unit2str;
152  // skip print of ConvergenceTester to stdout
153  G4cout << ss.str() << G4endl;
154  };
155  //----------------------------------------------------------------------//
156 
157  for(unsigned i = 0; i < primScorerNames.size(); ++i)
158  {
159  for(unsigned j = 0; j < fnames.size(); ++j)
160  {
161  fname = fnames.at(j) + "_" + primScorerNames.at(i) + ".out";
162  fileout.open(fname);
163  G4cout << separator.str() << G4endl;
164  G4cout << " opened file " << fname << " for output" << G4endl;
165  G4cout << separator.str() << G4endl;
166 
167  G4bool valid = true;
168  if(j == 0)
169  {
170  G4THitsMap<G4double>* hitmap
171  = tsRun->GetHitsMap(fName + "/" + primScorerNames.at(i));
173  = tsRun->GetStatMap(fName + "/" + primScorerNames.at(i));
175  = tsRun->GetConvMap(fName + "/" + primScorerNames.at(i));
176 
177  if(hitmap && hitmap->size() != 0)
178  {
179  for(auto itr = hitmap->begin(); itr != hitmap->end(); itr++)
180  {
181  if(!hitmap->GetObject(itr))
182  continue;
183  IDs.insert(itr->first);
184  std::get<0>(fTypeCompare[primScorerNames.at(i)][itr->first])
185  = *itr->second/units.at(i);
186  print(fileout, itr->first, *itr->second,
187  units.at(i), units.at(i+1), unitstr.at(i));
188  }
189  }
190  else
191  {
192  valid = false;
193  }
194 
195  if(statmap && statmap->size() != 0 &&
196  convmap && convmap->size() != 0)
197  {
198  auto stat_fname = "stat_" + fname;
199  std::ofstream statout;
200  statout.open(stat_fname);
201  for(auto itr = statmap->begin(); itr != statmap->end(); itr++)
202  {
203  G4int _f = statmap->GetIndex(itr);
204  G4StatAnalysis* _s = statmap->GetObject(itr);
205  G4ConvergenceTester* _c = convmap->GetObject(_f);
206  stat_print(statout, _f, _s, _c,
207  units.at(i), units.at(i+1), unitstr.at(i));
208  }
209  statout.close();
210  }
211  else
212  {
213  std::stringstream ss;
214  ss << " StatMap/ConvMap is either not "
215  << "created or the StatMap/ConvMap was empty";
216  if(statmap)
217  ss << " (StatMap size == " << statmap->size() << ")";
218  if(convmap)
219  ss << " (ConvMap size == " << convmap->size() << ")";
220 
221  G4Exception("TSRunAction", "002", JustWarning,
222  G4String(primScorerNames.at(i) +
223  ss.str()).c_str());
224  }
225 
226  if(!valid)
227  {
228  G4Exception("TSRunAction", "000", JustWarning,
229  G4String(primScorerNames.at(i) +
230  " HitsMap is either not "
231  "created or the HitsMap was empty").c_str());
232  }
233  }
234  else
235  {
237  = tsRun->GetAtomicHitsMap(fName + "/" +
238  primScorerNames.at(i));
239  if(hitmap && hitmap->size() != 0)
240  {
241  for(auto itr = hitmap->begin(); itr != hitmap->end(); itr++)
242  {
243  IDs.insert(itr->first);
244  std::get<1>(fTypeCompare[primScorerNames.at(i)][itr->first])
245  = *itr->second/units.at(i);
246  print(fileout, itr->first, *itr->second,
247  units.at(i), units.at(i+1), unitstr.at(i));
248  }
249  }
250  else
251  {
252  valid = false;
253  }
254 
255  if(!valid)
256  {
257  G4Exception("TSRunAction", "001", JustWarning,
258  G4String(primScorerNames.at(i) +
259  " HitsMap is either not "
260  "created or the HitsMap was empty").c_str());
261  }
262  }
263 
264  fileout.close();
265  G4cout << separator.str() << G4endl;
266  G4cout << " closed file " << fname << " for output" << G4endl;
267  }
268  // add the mutex data
269  TSRun::MutexHitsMap_t* hitmap
270  = tsRun->GetMutexHitsMap(fName + "/" +
271  primScorerNames.at(i));
272  if(hitmap && hitmap->size() != 0)
273  {
274  for(auto itr = hitmap->begin();
275  itr != hitmap->end(); itr++)
276  {
277  IDs.insert(itr->first);
278  std::get<2>(fTypeCompare[primScorerNames.at(i)][itr->first])
279  = itr->second/units.at(i);
280  }
281  }
282 
283  }
284 
285  //--------------------------------------------------------------------//
286  // Check that the values are equivalent and there are no
287  // IDs in one container that aren't in another
288  //--------------------------------------------------------------------//
289 
290  fname = "mfd_diff.out";
291  fileout.open(fname);
292 
293  G4cout << separator.str() << G4endl;
294  G4cout << " opened file " << fname << " for difference output" << G4endl;
295  G4cout << separator.str() << G4endl;
296 
297  fileout << " " << std::setw(10) << "ID"
298  << " "
299  << std::setw(30) << std::setprecision(12) << std::fixed
300  << "MFD value"
301  << " "
302  << std::setw(30) << std::setprecision(12) << std::fixed
303  << "Atomic Hits Map value"
304  << " "
305  << std::setw(30) << std::setprecision(8) << std::scientific
306  << "Difference"
307  << " "
308  << std::setw(30) << std::setprecision(8) << std::scientific
309  << "Diff (MFD - MUTEXED)"
310  << " "
311  << std::setw(30) << std::setprecision(8) << std::scientific
312  << "Diff (ATOM_HIT_MAP - MUTEXED)"
313  << G4endl << G4endl;
314 
315  for(auto itr1 = fTypeCompare.begin();
316  itr1 != fTypeCompare.end(); ++itr1)
317  {
318  fileout << "\n\nType = " << itr1->first << "\n" << G4endl;
319  for(auto itr2 = itr1->second.begin();
320  itr2 != itr1->second.end(); ++itr2)
321  {
322  G4double d01
323  = std::fabs(std::get<0>(itr2->second) -
324  std::get<1>(itr2->second));
325  G4double d02
326  = std::fabs(std::get<0>(itr2->second) -
327  std::get<2>(itr2->second));
328  G4double d03
329  = std::fabs(std::get<1>(itr2->second) -
330  std::get<2>(itr2->second));
331 
332 
333  auto _print_diff = [&] (const G4double& _dval)
334  {
335  if(_dval > 0.0)
336  fileout << std::setprecision(8) << std::scientific
337  << std::setw(30) << _dval << " ";
338  else
339  fileout << std::setprecision(1) << std::fixed
340  << std::setw(30) << _dval << " ";
341  };
342 
343  fileout
344  << " " << std::setw(10) << itr2->first
345  << " "
346  << std::setw(30) << std::setprecision(12) << std::fixed
347  << std::get<0>(itr2->second)
348  << " "
349  << std::setw(30) << std::setprecision(12) << std::fixed
350  << std::get<1>(itr2->second)
351  << " ";
352 
353  _print_diff(d01);
354  _print_diff(d02);
355  _print_diff(d03);
356 
357  fileout << G4endl;
358 
359  }
360  }
361 
362  fileout.close();
363  G4cout << " closed file " << fname << " for difference output" << G4endl;
364  G4cout << separator.str() << G4endl;
365  }
366 
367 
368 }
369 
370 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
371 
372 
373 
374