ATLAS Offline Software
Public Member Functions | Protected Member Functions | Protected Attributes | Private Attributes | List of all members
G4UA::RadiationMapsMakerTool Class Reference

Tool which manages the RadiationMapsMaker action. More...

#include <RadiationMapsMakerTool.h>

Inheritance diagram for G4UA::RadiationMapsMakerTool:
Collaboration diagram for G4UA::RadiationMapsMakerTool:

Public Member Functions

 RadiationMapsMakerTool (const std::string &type, const std::string &name, const IInterface *parent)
 standard tool ctor More...
 
virtual StatusCode initialize () override final
 Initialize configurable properties. More...
 
virtual StatusCode finalize () override final
 Finalize and merge results from all threads. More...
 
virtual StatusCode fillUserAction (G4AtlasUserActions &actionLists) override final
 Fill the user action lists. More...
 

Protected Member Functions

virtual std::unique_ptr< RadiationMapsMakermakeAndFillAction (G4AtlasUserActions &) override final
 create action for this thread More...
 

Protected Attributes

ThreadSpecificUserAction< RadiationMapsMakerm_actions
 Thread-specific storage of the user action. More...
 

Private Attributes

std::string m_radMapsFileName
 Output Filename for the Radiation Maps. More...
 
RadiationMapsMaker::Config m_config
 Radiation Map ranges and granularities. More...
 

Detailed Description

Tool which manages the RadiationMapsMaker action.

Create the RadiationMapsMaker for each worker thread

Author
Sven Menke (based on Andrea Di Simone's FluxRecorderTool)

Definition at line 26 of file RadiationMapsMakerTool.h.

Constructor & Destructor Documentation

◆ RadiationMapsMakerTool()

G4UA::RadiationMapsMakerTool::RadiationMapsMakerTool ( const std::string &  type,
const std::string &  name,
const IInterface *  parent 
)

standard tool ctor

Output Filename for the Radiation Maps

Name of the materials to make radiation maps for (take all if empty)

Name of the file to store activated unstable isotopes (don't save if empty)

If true consider hits with y>0 only – useful for shafts

map granularities number of bins in r and z for all 2D maps

number of bins in logE for energy spectra of neutrons in 2D grids

number of bins in logE for energy spectra of other particles in 2D grids

number of bins in dphi for dphi x theta dependent energy spectra

number of bins in theta for dphi x theta dependent energy spectra

number of bins in r, z and phi for all 3D maps

number of bins in logTimeCut for time dependent TID and H_T 2D maps

map ranges for Zoomed area in 2D and 3D

for Full detector in 2D

for Zoomed area in 3D

for logE of neutrons in 2D spectra

for logE of other particles in 2D spectra

for Theta in 2D spectra

for logT in time-dependent TID and H_T 2D maps

for elements mass fracion 2D maps

Definition at line 11 of file RadiationMapsMakerTool.cxx.

14  : UserActionToolBase<RadiationMapsMaker>(type, name, parent),
15  m_radMapsFileName("RadMaps.root")
16  {
18  declareProperty("RadMapsFileName" , m_radMapsFileName);
20  declareProperty("Materials" , m_config.materials);
22  declareProperty("ActivationFileName", m_config.activationFileName);
24  declareProperty("PositiveYOnly" , m_config.posYOnly);
27  declareProperty("NBinsR" , m_config.nBinsr);
28  declareProperty("NBinsZ" , m_config.nBinsz);
30  declareProperty("NBinsLogEn" , m_config.nBinslogEn);
32  declareProperty("NBinsLogEo" , m_config.nBinslogEo);
34  declareProperty("NBinsDPhi" , m_config.nBinsdphi);
36  declareProperty("NBinsTheta" , m_config.nBinstheta);
38  declareProperty("NBinsR3D" , m_config.nBinsr3d);
39  declareProperty("NBinsZ3D" , m_config.nBinsz3d);
40  declareProperty("NBinsPhi3D" , m_config.nBinsphi3d);
42  declareProperty("NBinsLogTimeCut" , m_config.nBinslogT);
45  declareProperty("RMinZoom" , m_config.rMinZoom);
46  declareProperty("RMaxZoom" , m_config.rMaxZoom);
47  declareProperty("ZMinZoom" , m_config.zMinZoom);
48  declareProperty("ZMaxZoom" , m_config.zMaxZoom);
50  declareProperty("RMinFull" , m_config.rMinFull);
51  declareProperty("RMaxFull" , m_config.rMaxFull);
52  declareProperty("ZMinFull" , m_config.zMinFull);
53  declareProperty("ZMaxFull" , m_config.zMaxFull);
55  declareProperty("PhiMinZoom" , m_config.phiMinZoom);
56  declareProperty("PhiMaxZoom" , m_config.phiMaxZoom);
58  declareProperty("LogEMinn" , m_config.logEMinn);
59  declareProperty("LogEMaxn" , m_config.logEMaxn);
61  declareProperty("LogEMino" , m_config.logEMino);
62  declareProperty("LogEMaxo" , m_config.logEMaxo);
64  declareProperty("ThetaMin" , m_config.thetaMin);
65  declareProperty("ThetaMax" , m_config.thetaMax);
67  declareProperty("LogTMin" , m_config.logTMin);
68  declareProperty("LogTMax" , m_config.logTMax);
70  declareProperty("ElemZMin" , m_config.elemZMin);
71  declareProperty("ElemZMax" , m_config.elemZMax);
72  }

Member Function Documentation

◆ fillUserAction()

virtual StatusCode G4UA::UserActionToolBase< RadiationMapsMaker >::fillUserAction ( G4AtlasUserActions actionLists)
inlinefinaloverridevirtualinherited

Fill the user action lists.

Definition at line 45 of file UserActionToolBase.h.

46  {
47  auto myAction = makeAndFillAction(actionLists);
48  if(myAction == nullptr) {
49  ATH_MSG_ERROR( "Failed to construct user action in " << name() );
50  return StatusCode::FAILURE;
51  }
52  m_actions.set( std::move(myAction) );
53  return StatusCode::SUCCESS;
54  }

◆ finalize()

StatusCode G4UA::RadiationMapsMakerTool::finalize ( )
finaloverridevirtual

Finalize and merge results from all threads.

Definition at line 134 of file RadiationMapsMakerTool.cxx.

135  {
136  ATH_MSG_DEBUG( "Finalizing " << name() );
137 
138  RadiationMapsMaker::Report maps;
139 
140  // vector of pointers to vectors of double to save space
141  std::vector<std::vector<double> *> pVal = {
142  &(maps.m_rz_tid),&(maps.m_rz_eion),&(maps.m_rz_niel),&(maps.m_rz_h20),&(maps.m_rz_neut),&(maps.m_rz_chad),
143  &(maps.m_full_rz_tid),&(maps.m_full_rz_eion),&(maps.m_full_rz_niel),&(maps.m_full_rz_h20),&(maps.m_full_rz_neut),&(maps.m_full_rz_chad)
144  };
145 
146  std::vector<std::vector<double> *> pVal3d = {
147  &(maps.m_3d_tid),&(maps.m_3d_eion),&(maps.m_3d_niel),&(maps.m_3d_h20),&(maps.m_3d_neut),&(maps.m_3d_chad)
148  };
149 
150  std::vector<std::vector<double> *> pValSpec = {
151  &(maps.m_rz_neut_spec),&(maps.m_rz_gamm_spec),&(maps.m_rz_elec_spec),&(maps.m_rz_muon_spec),&(maps.m_rz_pion_spec),&(maps.m_rz_prot_spec),&(maps.m_rz_rest_spec),
152  &(maps.m_full_rz_neut_spec),&(maps.m_full_rz_gamm_spec),&(maps.m_full_rz_elec_spec),&(maps.m_full_rz_muon_spec),&(maps.m_full_rz_pion_spec),&(maps.m_full_rz_prot_spec),&(maps.m_full_rz_rest_spec)
153  };
154 
155  std::vector<int> nBinslogE = {
158  };
159 
160  std::vector<std::vector<double> *> pValThetaSpec = {
161  &(maps.m_theta_full_rz_neut_spec),&(maps.m_theta_full_rz_gamm_spec),&(maps.m_theta_full_rz_elec_spec),&(maps.m_theta_full_rz_muon_spec),&(maps.m_theta_full_rz_pion_spec),&(maps.m_theta_full_rz_prot_spec),&(maps.m_theta_full_rz_rchgd_spec),&(maps.m_theta_full_rz_rneut_spec)
162  };
163 
164  std::vector<int> nBinsThetalogE = {
166  };
167 
168  for(unsigned int hi=0;hi<pVal.size();hi++) {
169 
170  // first make sure the vectors are empty
171 
172  (pVal[hi])->resize(0);
173 
174  // then resize to proper size and initialize with 0's
175  // all maps are needed for the merge - so have to do resize
176  // for all before merging any ...
177 
178  (pVal[hi])->resize(m_config.nBinsz*m_config.nBinsr,0.0);
179  }
180 
181  // same for 3d
182 
183  for(unsigned int hi=0;hi<pVal3d.size();hi++) {
184  (pVal3d[hi])->resize(0);
185  (pVal3d[hi])->resize(m_config.nBinsz3d*m_config.nBinsr3d*m_config.nBinsphi3d,0.0);
186  }
187 
188  // same for spectra
189 
190  for(unsigned int hi=0;hi<pValSpec.size();hi++) {
191  (pValSpec[hi])->resize(0);
192  (pValSpec[hi])->resize(m_config.nBinsz*m_config.nBinsr*nBinslogE[hi],0.0);
193  }
194 
195  // same for theta x dphi spectra
196 
197  for(unsigned int hi=0;hi<pValThetaSpec.size();hi++) {
198  (pValThetaSpec[hi])->resize(0);
199  (pValThetaSpec[hi])->resize(m_config.nBinsdphi*m_config.nBinstheta*m_config.nBinsz*m_config.nBinsr*nBinsThetalogE[hi],0.0);
200  }
201 
202  // same for misc individual maps
203 
204  maps.m_rz_tid_time .resize(0);
205  maps.m_rz_ht_time .resize(0);
206  maps.m_full_rz_tid_time .resize(0);
207  maps.m_full_rz_ht_time .resize(0);
208 
209  maps.m_rz_element .resize(0);
210  maps.m_full_rz_element .resize(0);
211 
212  if (!m_config.materials.empty()) {
213  // need volume fraction only if particular materials are selected
214  // 2d zoom
215  maps.m_rz_vol .resize(0);
216  maps.m_rz_norm.resize(0);
217  // 2d full
218  maps.m_full_rz_vol.resize(0);
219  maps.m_full_rz_norm.resize(0);
220  // 3d
221  maps.m_3d_vol .resize(0);
222  maps.m_3d_norm.resize(0);
223  }
224 
225  // then resize to proper size and initialize with 0's
226  // all maps are needed for the merge - so have to do resize
227  // for all first ...
228 
229  maps.m_rz_tid_time .resize(m_config.nBinsz*m_config.nBinsr*m_config.nBinslogT,0.0);
230  maps.m_rz_ht_time .resize(m_config.nBinsz*m_config.nBinsr*m_config.nBinslogT,0.0);
231  maps.m_full_rz_tid_time .resize(m_config.nBinsz*m_config.nBinsr*m_config.nBinslogT,0.0);
232  maps.m_full_rz_ht_time .resize(m_config.nBinsz*m_config.nBinsr*m_config.nBinslogT,0.0);
233 
234  maps.m_rz_element .resize(m_config.nBinsz*m_config.nBinsr*(m_config.elemZMax-m_config.elemZMin+1),0.0);
235  maps.m_full_rz_element .resize(m_config.nBinsz*m_config.nBinsr*(m_config.elemZMax-m_config.elemZMin+1),0.0);
236 
237  if (!m_config.materials.empty()) {
238  // need volume fraction only if particular materials are selected
239  // 2d zoom
240  maps.m_rz_vol .resize(m_config.nBinsz*m_config.nBinsr,0.0);
241  maps.m_rz_norm.resize(m_config.nBinsz*m_config.nBinsr,0.0);
242  // 2d full
243  maps.m_full_rz_vol .resize(m_config.nBinsz*m_config.nBinsr,0.0);
244  maps.m_full_rz_norm.resize(m_config.nBinsz*m_config.nBinsr,0.0);
245  // 3d
246  maps.m_3d_vol .resize(m_config.nBinsz3d*m_config.nBinsr3d*m_config.nBinsphi3d,0.0);
247  maps.m_3d_norm.resize(m_config.nBinsz3d*m_config.nBinsr3d*m_config.nBinsphi3d,0.0);
248  }
249 
250  // merge radiation map vectors from threads
251  // Accumulate the results across threads
254 
255  TFile * f = new TFile(m_radMapsFileName.c_str(),"RECREATE");
256 
257  // histograms can be booked, filled, written and deleted in a loop
258  // in order to save memory. Need vectors to loop over them ...
259 
260  std::vector<const char *> h2dNames = {
261  "rz_tid","rz_eion","rz_niel","rz_h20","rz_neut","rz_chad",
262  "full_rz_tid","full_rz_eion","full_rz_niel","full_rz_h20","full_rz_neut","full_rz_chad"
263  };
264 
265  std::vector<const char *> h2dZTitles = {
266  "TID [Gy]","E_{ion}/V [MeV/cm^{3}]","NIEL [n_{eq}/cm^{2}]","SEE [h_{>20 MeV}/cm^{2}]","FLUX [n_{>100 keV}/cm^{2}]","FLUX [h_{charged}/cm^{2}]",
267  "TID [Gy]","E_{ion}/V [MeV/cm^{3}]","NIEL [n_{eq}/cm^{2}]","SEE [h_{>20 MeV}/cm^{2}]","FLUX [n_{>100 keV}/cm^{2}]","FLUX [h_{charged}/cm^{2}]"
268  };
269 
270  std::vector<double> zMin = {
273  };
274 
275  std::vector<double> zMax = {
278  };
279 
280  std::vector<double> rMin = {
283  };
284 
285  std::vector<double> rMax = {
288  };
289 
290  const char * xtitle = (m_config.zMinFull<0?"z [cm]":"|z| [cm]");
291 
292  TH2D * h_rz_vol = 0;
293  TH2D * h_rz_norm = 0;
294  TH2D * h_full_rz_vol = 0;
295  TH2D * h_full_rz_norm = 0;
296 
297  for(unsigned int hi=0;hi<h2dNames.size();hi++) {
298  TH2D *h2d = new TH2D(h2dNames[hi],h2dNames[hi],m_config.nBinsz,zMin[hi],zMax[hi],m_config.nBinsr,rMin[hi],rMax[hi]);
299  h2d->SetXTitle(xtitle);
300  h2d->SetYTitle("r [cm]");
301  h2d->SetZTitle(h2dZTitles[hi]);
302  if (hi==0 && !m_config.materials.empty()) {
303  // need volume fraction only if particular materials are selected
304  //
305  // the maps for TID, NIEL and SEE need to be divided by the ratio of (vol/norm) in order to get
306  // the proper estimate per volume bin for the selected material.
307  // This is *not* done in the tool directly and left to the user after having summed the histograms
308  // from many individual jobs.
309  //
312  h_full_rz_vol = new TH2D("full_rz_vol" ,"full_rz_vol" ,m_config.nBinsz,m_config.zMinFull,m_config.zMaxFull,m_config.nBinsr,m_config.rMinFull,m_config.rMaxFull);
313  h_full_rz_norm = new TH2D("full_rz_norm","full_rz_norm",m_config.nBinsz,m_config.zMinFull,m_config.zMaxFull,m_config.nBinsr,m_config.rMinFull,m_config.rMaxFull);
314 
315  h_rz_vol ->SetXTitle(xtitle);
316  h_rz_norm ->SetXTitle(xtitle);
317  h_rz_vol ->SetYTitle("r [cm]");
318  h_rz_norm ->SetYTitle("r [cm]");
319  std::string hname("Volume fraction of");
320  char c = ' ';
321  for(const std::string& matName : m_config.materials) {
322  hname += c;
323  hname += matName;
324  c = ',';
325  }
326  h_rz_vol ->SetZTitle(hname.data());
327  h_rz_norm ->SetZTitle("Volume norm");
328 
329  h_full_rz_vol ->SetXTitle(xtitle);
330  h_full_rz_norm ->SetXTitle(xtitle);
331  h_full_rz_vol ->SetYTitle("r [cm]");
332  h_full_rz_norm ->SetYTitle("r [cm]");
333  h_full_rz_vol ->SetZTitle(hname.data());
334  h_full_rz_norm ->SetZTitle("Volume norm");
335  }
336 
337  // normalize to volume element per bin
338  for(int i=0;i<h2d->GetNbinsX();i++) {
339  for(int j=0;j<h2d->GetNbinsY();j++) {
340  int iBin = h2d->GetBin(i+1,j+1);
341  int vBin = m_config.nBinsr*i+j;
342  double r0=h2d->GetYaxis()->GetBinLowEdge(j+1);
343  double r1=h2d->GetYaxis()->GetBinUpEdge(j+1);
344  double z0=h2d->GetXaxis()->GetBinLowEdge(i+1);
345  double z1=h2d->GetXaxis()->GetBinUpEdge(i+1);
346  double vol=(z1-z0)*M_PI*(r1*r1-r0*r0);
347  // if |z| instead of z double the volume
348  if ( m_config.zMinFull >= 0 ) vol *= 2;
349  // if positive y hemisphere is used only -- half the volume
350  if ( m_config.posYOnly ) vol *= 0.5;
351  double val = (*(pVal[hi]))[vBin];
352  h2d->SetBinContent(iBin,val/vol);
353  if (hi ==0 && !m_config.materials.empty()) {
354  // need volume fraction only if particular materials are selected
355  // VOL
356  val =maps.m_rz_vol[vBin];
357  h_rz_vol->SetBinContent(iBin,val/vol);
358  // NORM
359  val =maps.m_rz_norm[vBin];
360  h_rz_norm->SetBinContent(iBin,val/vol);
361  }
362  if (hi ==h2dNames.size()/2 && !m_config.materials.empty()) {
363  // need volume fraction only if particular materials are selected
364  // VOL
365  val =maps.m_full_rz_vol[vBin];
366  h_full_rz_vol->SetBinContent(iBin,val/vol);
367  // NORM
368  val =maps.m_full_rz_norm[vBin];
369  h_full_rz_norm->SetBinContent(iBin,val/vol);
370  }
371  }
372  }
373  h2d->Write();
374  h2d->Delete();
375  if (hi ==h2dNames.size()/2 && !m_config.materials.empty()) {
376  h_rz_vol->Write();
377  h_rz_vol->Delete();
378  h_rz_norm->Write();
379  h_rz_norm->Delete();
380  h_full_rz_vol->Write();
381  h_full_rz_vol->Delete();
382  h_full_rz_norm->Write();
383  h_full_rz_norm->Delete();
384  }
385  }
386 
387  std::vector<const char *> h3dNames = {
388  "h3d_tid","h3d_eion","h3d_niel","h3d_h20","h3d_neut","h3d_chad"
389  };
390 
391  std::vector<const char *> h3dTitles = {
392  "TID [Gy]","E_{ion}/V [MeV/cm^{3}]","NIEL [n_{eq}/cm^{2}]","SEE [h_{>20 MeV}/cm^{2}]","FLUX [n_{>100 keV}/cm^{2}]","FLUX [h_{charged}/cm^{2}]"
393  };
394 
395  TH3D * h_3d_vol = 0;
396  TH3D * h_3d_norm = 0;
397 
398  for(unsigned int hi=0;hi<h3dNames.size();hi++) {
400  h3d->SetXTitle(xtitle);
401  h3d->SetYTitle("r [cm]");
402  h3d->SetZTitle("#phi [#circ]");
403  h3d->SetTitle(h3dTitles[hi]);
404  if (hi == 0 && !m_config.materials.empty()) {
407  h_3d_vol ->SetXTitle(xtitle);
408  h_3d_norm ->SetXTitle(xtitle);
409  h_3d_vol ->SetYTitle("r [cm]");
410  h_3d_norm ->SetYTitle("r [cm]");
411  h_3d_vol ->SetZTitle("#phi [#circ]");
412  h_3d_norm ->SetZTitle("#phi [#circ]");
413  std::string hname("Volume fraction of");
414  char c = ' ';
415  for(const std::string& matName : m_config.materials) {
416  hname += c;
417  hname += matName;
418  c = ',';
419  }
420  h_3d_vol ->SetTitle(hname.data());
421  h_3d_norm ->SetTitle("Volume norm");
422  }
423  // normalize to volume element per bin
424  for(int i=0;i<h3d->GetNbinsX();i++) {
425  for(int j=0;j<h3d->GetNbinsY();j++) {
426  for(int k=0;k<h3d->GetNbinsZ();k++) {
428  int iBin = h3d->GetBin(i+1,j+1,k+1);
429  double phi0=h3d->GetZaxis()->GetBinLowEdge(k+1);
430  double phi1=h3d->GetZaxis()->GetBinUpEdge(k+1);
431  double r0=h3d->GetYaxis()->GetBinLowEdge(j+1);
432  double r1=h3d->GetYaxis()->GetBinUpEdge(j+1);
433  double z0=h3d->GetXaxis()->GetBinLowEdge(i+1);
434  double z1=h3d->GetXaxis()->GetBinUpEdge(i+1);
435  double vol=(z1-z0)*M_PI*(r1*r1-r0*r0)*(phi1-phi0)/360.;
436  // if |z| instead of z double the volume
437  if ( m_config.zMinFull >= 0 ) vol *= 2;
438  // assume that phi-range corresponds to full 360 degrees in case
439  // lower phi boundary is 0 - i.e. all phi-segments mapped to first
440  if ( m_config.phiMinZoom == 0 ) {
441  vol *= 360./m_config.phiMaxZoom;
442  // if positive y hemisphere is used only -- half the volume
443  if ( m_config.posYOnly )
444  vol *= 0.5;
445  }
446  double val = (*(pVal3d[hi]))[vBin];
447  h3d->SetBinContent(iBin,val/vol);
448  if (hi == 0 && !m_config.materials.empty()) {
449  // VOL
450  val =maps.m_3d_vol[vBin];
451  h_3d_vol->SetBinContent(iBin,val/vol);
452  // NORM
453  val =maps.m_3d_norm[vBin];
454  h_3d_norm->SetBinContent(iBin,val/vol);
455  }
456  }
457  }
458  }
459  h3d->Write();
460  h3d->Delete();
461  if (hi == 0 && !m_config.materials.empty()) {
462  h_3d_vol->Write();
463  h_3d_vol->Delete();
464  h_3d_norm->Write();
465  h_3d_norm->Delete();
466  }
467  }
468 
469  // spectra
470 
471  std::vector<const char *> hSpecNames = {
472  "rz_neut_spec","rz_gamm_spec","rz_elec_spec","rz_muon_spec","rz_pion_spec","rz_prot_spec","rz_rest_spec",
473  "full_rz_neut_spec","full_rz_gamm_spec","full_rz_elec_spec","full_rz_muon_spec","full_rz_pion_spec","full_rz_prot_spec","full_rz_rest_spec"
474  };
475 
476  std::vector<const char *> hSpecTitles = {
477  "FLUX [n(log_{10}(E/MeV))/cm^{2}]","FLUX [#gamma(log_{10}(E/MeV))/cm^{2}]","FLUX [e^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [#mu^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [#pi^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [p(log_{10}(E/MeV))/cm^{2}]","FLUX [rest(log_{10}(E/MeV))/cm^{2}]",
478  "FLUX [n(log_{10}(E/MeV))/cm^{2}]","FLUX [#gamma(log_{10}(E/MeV))/cm^{2}]","FLUX [e^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [#mu^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [#pi^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [p(log_{10}(E/MeV))/cm^{2}]","FLUX [rest(log_{10}(E/MeV))/cm^{2}]"
479  };
480 
481  std::vector<double> zMinSpec = {
484  };
485 
486  std::vector<double> zMaxSpec = {
489  };
490 
491  std::vector<double> rMinSpec = {
494  };
495 
496  std::vector<double> rMaxSpec = {
499  };
500 
501  std::vector<double> logEMin = {
504  };
505 
506  std::vector<double> logEMax = {
509  };
510 
511  for(unsigned int hi=0;hi<hSpecNames.size();hi++) {
512  TH3D *hSpec = new TH3D(hSpecNames[hi],hSpecNames[hi],m_config.nBinsz,zMinSpec[hi],zMaxSpec[hi],m_config.nBinsr,rMinSpec[hi],rMaxSpec[hi],nBinslogE[hi],logEMin[hi],logEMax[hi]);
513  hSpec->SetXTitle(xtitle);
514  hSpec->SetYTitle("r [cm]");
515  hSpec->SetZTitle("log_{10}(E/MeV)");
516  hSpec->SetTitle(hSpecTitles[hi]);
517  // normalize to volume element per bin
518  for(int i=0;i<hSpec->GetNbinsX();i++) {
519  for(int j=0;j<hSpec->GetNbinsY();j++) {
520  double r0=hSpec->GetYaxis()->GetBinLowEdge(j+1);
521  double r1=hSpec->GetYaxis()->GetBinUpEdge(j+1);
522  double z0=hSpec->GetXaxis()->GetBinLowEdge(i+1);
523  double z1=hSpec->GetXaxis()->GetBinUpEdge(i+1);
524  double vol=(z1-z0)*M_PI*(r1*r1-r0*r0);
525  // if |z| instead of z double the volume
526  if ( m_config.zMinFull >= 0 ) vol *= 2;
527  // if positive y hemisphere is used only -- half the volume
528  if ( m_config.posYOnly ) vol *= 0.5;
529  double val;
530  // Energy Spectra
531  for(int k=0;k<hSpec->GetNbinsZ();k++) {
532  int kBin = hSpec->GetBin(i+1,j+1,k+1);
533  int vBinE = m_config.nBinsr*nBinslogE[hi]*i+j*nBinslogE[hi]+k;
534  val = (*(pValSpec[hi]))[vBinE];
535  hSpec->SetBinContent(kBin,val/vol);
536  }
537  }
538  }
539  hSpec->Write();
540  hSpec->Delete();
541  }
542 
543  std::vector<const char *> hThetaSpecNames = {
544  "theta_dphi_full_rz_neut_spec","theta_dphi_full_rz_gamm_spec","theta_dphi_full_rz_elec_spec","theta_dphi_full_rz_muon_spec","theta_dphi_full_rz_pion_spec","theta_dphi_full_rz_prot_spec","theta_dphi_full_rz_rchgd_spec","theta_dphi_full_rz_rneut_spec"
545  };
546 
547  std::vector<const char *> hThetaSpecTitles = {
548  "FLUX [n(log_{10}(E/MeV))/cm^{2}]","FLUX [#gamma(log_{10}(E/MeV))/cm^{2}]","FLUX [e^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [#mu^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [#pi^{#pm}(log_{10}(E/MeV))/cm^{2}]","FLUX [p(log_{10}(E/MeV))/cm^{2}]","FLUX [rest^{chgd}(log_{10}(E/MeV))/cm^{2}]","FLUX [rest^{neut}(log_{10}(E/MeV))/cm^{2}]"
549  };
550 
551  std::vector<double> thetalogEMin = {
553  };
554 
555  std::vector<double> thetalogEMax = {
557  };
558 
559  for(unsigned int hi=0;hi<hThetaSpecNames.size();hi++) {
560  for(int i=0;i<m_config.nBinstheta;i++) {
561  for(int j=0;j<m_config.nBinsdphi;j++) {
562  TString hname(hThetaSpecNames[hi]);
563  hname += "_";
565  hname += "_";
567  hname += "_";
568  hname += (int)(j*360./m_config.nBinsdphi);
569  hname += "_";
570  hname += (int)((j+1)*360./m_config.nBinsdphi);
571  TH3D * hThetaSpec = new TH3D(hname.Data(),hname.Data(),m_config.nBinsz,m_config.zMinFull,m_config.zMaxFull,m_config.nBinsr,m_config.rMinFull,m_config.rMaxFull,nBinsThetalogE[hi],thetalogEMin[hi],thetalogEMax[hi]);
572 
573  hThetaSpec->SetXTitle(xtitle);
574  hThetaSpec->SetYTitle("r [cm]");
575  hThetaSpec->SetZTitle("log_{10}(E/MeV)");
576  hname = TString(hThetaSpecTitles[hi]);
577  hname += " Theta: ";
579  hname += " - ";
581  hname += " DPhi: ";
582  hname += (int)(j*360./m_config.nBinsdphi);
583  hname += " - ";
584  hname += (int)((j+1)*360./m_config.nBinsdphi);
585  hThetaSpec->SetTitle(hname.Data());
586  // normalize to volume element per bin
587  for(int k=0;k<hThetaSpec->GetNbinsX();k++) {
588  for(int l=0;l<hThetaSpec->GetNbinsY();l++) {
589  double r0=hThetaSpec->GetYaxis()->GetBinLowEdge(l+1);
590  double r1=hThetaSpec->GetYaxis()->GetBinUpEdge(l+1);
591  double z0=hThetaSpec->GetXaxis()->GetBinLowEdge(k+1);
592  double z1=hThetaSpec->GetXaxis()->GetBinUpEdge(k+1);
593  double vol=(z1-z0)*M_PI*(r1*r1-r0*r0);
594  // if |z| instead of z double the volume
595  if ( m_config.zMinFull >= 0 ) vol *= 2;
596  // if positive y hemisphere is used only -- half the volume
597  if ( m_config.posYOnly ) vol *= 0.5;
598  double val;
599  // Energy Spectra
600  for(int m=0;m<hThetaSpec->GetNbinsZ();m++) {
601  int mBin = hThetaSpec->GetBin(k+1,l+1,m+1);
603  val = (*(pValThetaSpec[hi]))[vBinETh];
604  hThetaSpec->SetBinContent(mBin,val/vol);
605  }
606  }
607  }
608  hThetaSpec->Write();
609  hThetaSpec->Delete();
610  }
611  }
612  }
613 
614  // time dependent TID and H_T maps zoom
616  h_rz_tid_time ->SetXTitle(xtitle);
617  h_rz_tid_time ->SetYTitle("r [cm]");
618  h_rz_tid_time ->SetZTitle("log_{10}(t_{cut}/s)");
619 
621  h_rz_ht_time ->SetXTitle(xtitle);
622  h_rz_ht_time ->SetYTitle("r [cm]");
623  h_rz_ht_time ->SetZTitle("log_{10}(t_{cut}/s)");
624 
625  // element mass fraction maps zoom
627  h_rz_element ->SetXTitle(xtitle);
628  h_rz_element ->SetYTitle("r [cm]");
629  h_rz_element ->SetZTitle("Z");
630 
631  // normalize to volume element per bin
632  for(int i=0;i<h_rz_tid_time->GetNbinsX();i++) {
633  for(int j=0;j<h_rz_tid_time->GetNbinsY();j++) {
634  double r0=h_rz_tid_time->GetYaxis()->GetBinLowEdge(j+1);
635  double r1=h_rz_tid_time->GetYaxis()->GetBinUpEdge(j+1);
636  double z0=h_rz_tid_time->GetXaxis()->GetBinLowEdge(i+1);
637  double z1=h_rz_tid_time->GetXaxis()->GetBinUpEdge(i+1);
638  double vol=(z1-z0)*M_PI*(r1*r1-r0*r0);
639  // if |z| instead of z double the volume
640  if ( m_config.zMinFull >= 0 ) vol *= 2;
641  // if positive y hemisphere is used only -- half the volume
642  if ( m_config.posYOnly ) vol *= 0.5;
643  double val;
644  // Time dependent TID and H_T maps
645  for(int k=0;k<h_rz_tid_time->GetNbinsZ();k++) {
646  int kBin = h_rz_tid_time->GetBin(i+1,j+1,k+1);
648  // TID
649  val =maps.m_rz_tid_time[vBinT];
650  h_rz_tid_time->SetBinContent(kBin,val/vol);
651  // H_T
652  val =maps.m_rz_ht_time[vBinT];
653  h_rz_ht_time->SetBinContent(kBin,val/vol);
654  }
655  // Element mass fraction maps
656  for(int k=0;k<h_rz_element->GetNbinsZ();k++) {
657  int kBin = h_rz_element->GetBin(i+1,j+1,k+1);
659  val =maps.m_rz_element[vBinElem];
660  h_rz_element->SetBinContent(kBin,val/vol);
661  }
662  }
663  }
664 
665  h_rz_tid_time->Write();
666  h_rz_tid_time->Delete();
667  h_rz_ht_time->Write();
668  h_rz_ht_time->Delete();
669  h_rz_element->Write();
670  h_rz_element->Delete();
671 
672  // time dependent TID and H_T maps full
674  h_full_rz_tid_time->SetXTitle(xtitle);
675  h_full_rz_tid_time->SetYTitle("r [cm]");
676  h_full_rz_tid_time->SetZTitle("log_{10}(t_{cut}/s)");
677 
679  h_full_rz_ht_time->SetXTitle(xtitle);
680  h_full_rz_ht_time->SetYTitle("r [cm]");
681  h_full_rz_ht_time->SetZTitle("log_{10}(t_{cut}/s)");
682 
683  // element mass fraction maps full
685  h_full_rz_element->SetXTitle(xtitle);
686  h_full_rz_element->SetYTitle("r [cm]");
687  h_full_rz_element->SetZTitle("Z");
688 
689  // normalize to volume element per bin
690  for(int i=0;i<h_full_rz_tid_time->GetNbinsX();i++) {
691  for(int j=0;j<h_full_rz_tid_time->GetNbinsY();j++) {
692  double r0=h_full_rz_tid_time->GetYaxis()->GetBinLowEdge(j+1);
693  double r1=h_full_rz_tid_time->GetYaxis()->GetBinUpEdge(j+1);
694  double z0=h_full_rz_tid_time->GetXaxis()->GetBinLowEdge(i+1);
695  double z1=h_full_rz_tid_time->GetXaxis()->GetBinUpEdge(i+1);
696  double vol=(z1-z0)*M_PI*(r1*r1-r0*r0);
697  // if |z| instead of z double the volume
698  if ( m_config.zMinFull >= 0 ) vol *= 2;
699  // if positive y hemisphere is used only -- half the volume
700  if ( m_config.posYOnly ) vol *= 0.5;
701  double val;
702  // Time dependent TID and H_T maps
703  for(int k=0;k<h_full_rz_tid_time->GetNbinsZ();k++) {
704  int kBin = h_full_rz_tid_time->GetBin(i+1,j+1,k+1);
706  // TID
707  val =maps.m_full_rz_tid_time[vBinT];
708  h_full_rz_tid_time->SetBinContent(kBin,val/vol);
709  // H_T
710  val =maps.m_full_rz_ht_time[vBinT];
711  h_full_rz_ht_time->SetBinContent(kBin,val/vol);
712  }
713  // Element mass fraction maps
714  for(int k=0;k<h_full_rz_element->GetNbinsZ();k++) {
715  int kBin = h_full_rz_element->GetBin(i+1,j+1,k+1);
717  val =maps.m_full_rz_element[vBinElem];
718  h_full_rz_element->SetBinContent(kBin,val/vol);
719  }
720  }
721  }
722  h_full_rz_tid_time->Write();
723  h_full_rz_tid_time->Delete();
724  h_full_rz_ht_time->Write();
725  h_full_rz_ht_time->Delete();
726  h_full_rz_element->Write();
727  h_full_rz_element->Delete();
728 
729  f->Close();
730 
731  return StatusCode::SUCCESS;
732  }

◆ initialize()

StatusCode G4UA::RadiationMapsMakerTool::initialize ( )
finaloverridevirtual

Initialize configurable properties.

Definition at line 77 of file RadiationMapsMakerTool.cxx.

78  {
79  msg(MSG::INFO)
80  << "Initializing " << name() << "\n"
81  << "OutputFile: " << m_radMapsFileName << "\n"
82  << "ActivationFile: " << m_config.activationFileName << "\n"
83  << "Materials: ";
84  char c=' ';
85  for(const std::string& matName : m_config.materials) {
86  msg() << c << matName;
87  c=',';
88  }
89  msg()
90  << "\n"
91  << "PositiveYOnly: " << m_config.posYOnly << "\n"
92  << "2D Maps: " << m_config.nBinsz << (m_config.zMinFull<0?" z-bins, ":" |z|-bins, ")
93  << m_config.nBinsr << " r-bins" << "\n"
94  << "Zoom: " << m_config.zMinZoom << (m_config.zMinFull<0?" < z/cm < ":" < |z|/cm < ")
95  << m_config.zMaxZoom << ", "
96  << m_config.rMinZoom << " < r/cm < "
97  << m_config.rMaxZoom << "\n"
98  << "Full: " << m_config.zMinFull << (m_config.zMinFull<0?" < z/cm < ":" < |z|/cm < ")
99  << m_config.zMaxFull << ", "
100  << m_config.rMinFull << " < r/cm < "
101  << m_config.rMaxFull << "\n"
102  << "Neutron Spectra: " << m_config.nBinslogEn << " log10E-bins, "
103  << m_config.logEMinn << " < log10(E/MeV) < " << m_config.logEMaxn << "\n"
104  << "Other Spectra: " << m_config.nBinslogEo << " log10E-bins, "
105  << m_config.logEMino << " < log10(E/MeV) < " << m_config.logEMaxo << "\n"
106  << "DPhi-Bins: " << m_config.nBinsdphi << " dphi-bins, 0 < dphi < 360 \n"
107  << "Theta-Bins: " << m_config.nBinstheta << " theta-bins, "
108  << m_config.thetaMin << " < theta < " << m_config.thetaMax << "\n"
109  << "3D Maps: " << m_config.nBinsz3d << (m_config.zMinFull<0?" z-bins, ":" |z|-bins, ")
110  << m_config.nBinsr3d << " r-bins, " << m_config.nBinsphi3d << " phi-bins" << "\n"
111  << "Zoom: " << m_config.zMinZoom << (m_config.zMinFull<0?" < z/cm < ":" < |z|/cm < ")
112  << m_config.zMaxZoom << ", "
113  << m_config.rMinZoom << " < r/cm < " << m_config.rMaxZoom << ", "
114  << m_config.phiMinZoom << " < phi/degrees < " << m_config.phiMaxZoom << "\n"
115  << "Time TID/H Maps: " << m_config.nBinslogT << " Time-cut bins, "
116  << m_config.logTMin << " < log10(t_cut/s) < "<< m_config.logTMax << "\n"
117  << "Mass frac. Maps: " << m_config.elemZMax-m_config.elemZMin+1 << " Element bins, "
118  << m_config.elemZMin << " <= Z <= < " << m_config.elemZMax
119  << endmsg;
120 
121  // clear the activation file if requested
122  if ( !m_config.activationFileName.empty() ) {
123  std::ofstream actf;
124  actf.open(m_config.activationFileName);
125  actf.close();
126  }
127 
128  return StatusCode::SUCCESS;
129  }

◆ makeAndFillAction()

std::unique_ptr< RadiationMapsMaker > G4UA::RadiationMapsMakerTool::makeAndFillAction ( G4AtlasUserActions actionList)
finaloverrideprotectedvirtual

create action for this thread

Implements G4UA::UserActionToolBase< RadiationMapsMaker >.

Definition at line 734 of file RadiationMapsMakerTool.cxx.

734  {
735  ATH_MSG_INFO("Making a RadiationMapsMaker action");
736  auto action = std::make_unique<RadiationMapsMaker>(m_config);
737  actionList.runActions.push_back( action.get() );
738  actionList.steppingActions.push_back( action.get() );
739  return action;
740  }

Member Data Documentation

◆ m_actions

Thread-specific storage of the user action.

Definition at line 63 of file UserActionToolBase.h.

◆ m_config

RadiationMapsMaker::Config G4UA::RadiationMapsMakerTool::m_config
private

Radiation Map ranges and granularities.

Definition at line 51 of file RadiationMapsMakerTool.h.

◆ m_radMapsFileName

std::string G4UA::RadiationMapsMakerTool::m_radMapsFileName
private

Output Filename for the Radiation Maps.

Definition at line 49 of file RadiationMapsMakerTool.h.


The documentation for this class was generated from the following files:
G4UA::RadiationMapsMaker::Config::nBinslogT
int nBinslogT
Definition: RadiationMapsMaker.h:79
python.CaloRecoConfig.f
f
Definition: CaloRecoConfig.py:127
python.SystemOfUnits.m
int m
Definition: SystemOfUnits.py:91
G4UA::RadiationMapsMaker::Config::logEMino
double logEMino
Definition: RadiationMapsMaker.h:75
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
TH3D
Definition: rootspy.cxx:505
CaloCellPos2Ntuple.int
int
Definition: CaloCellPos2Ntuple.py:24
G4UA::RadiationMapsMaker::Config::nBinslogEn
int nBinslogEn
Definition: RadiationMapsMaker.h:69
dqt_zlumi_pandas.hname
string hname
Definition: dqt_zlumi_pandas.py:272
G4UA::ThreadSpecificUserAction::set
void set(std::unique_ptr< ActionType > action)
Assign the object of the current thread.
Definition: ThreadSpecificUserAction.h:61
G4UA::RadiationMapsMaker::Config::logEMinn
double logEMinn
Definition: RadiationMapsMaker.h:70
G4UA::RadiationMapsMaker::Config::activationFileName
std::string activationFileName
Definition: RadiationMapsMaker.h:34
InDetAccessor::phi0
@ phi0
Definition: InDetAccessor.h:33
G4UA::RadiationMapsMaker::Config::rMinZoom
double rMinZoom
Definition: RadiationMapsMaker.h:41
TH2D::SetBinContent
void SetBinContent(int, double)
Definition: rootspy.cxx:436
TH3D::SetBinContent
void SetBinContent(int, double)
Definition: rootspy.cxx:511
G4UA::RadiationMapsMaker::Config::logTMax
double logTMax
Definition: RadiationMapsMaker.h:81
M_PI
#define M_PI
Definition: ActiveFraction.h:11
UploadAMITag.l
list l
Definition: UploadAMITag.larcaf.py:158
G4UA::ThreadSpecificUserAction::accumulate
void accumulate(ResultType &result, Mapper mapOp, Reducer reduceOp)
Accumulate results across user actions with specified operations.
Definition: ThreadSpecificUserAction.h:88
G4UA::RadiationMapsMaker::Config::nBinstheta
int nBinstheta
Definition: RadiationMapsMaker.h:63
G4UA::RadiationMapsMaker::Config::nBinsz3d
int nBinsz3d
Definition: RadiationMapsMaker.h:54
G4UA::RadiationMapsMaker::Config::zMinZoom
double zMinZoom
Definition: RadiationMapsMaker.h:47
G4UA::RadiationMapsMaker::Config::materials
std::vector< std::string > materials
bin sizes and ranges match the requirements for the Radiation Estimate Web tool for the default value...
Definition: RadiationMapsMaker.h:32
G4UA::RadiationMapsMaker::Config::thetaMin
double thetaMin
Definition: RadiationMapsMaker.h:65
G4UA::RadiationMapsMaker::Config::rMaxFull
double rMaxFull
Definition: RadiationMapsMaker.h:45
G4UA::RadiationMapsMaker::Config::zMaxFull
double zMaxFull
Definition: RadiationMapsMaker.h:51
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
lumiFormat.i
int i
Definition: lumiFormat.py:92
G4UA::RadiationMapsMakerTool::m_config
RadiationMapsMaker::Config m_config
Radiation Map ranges and granularities.
Definition: RadiationMapsMakerTool.h:51
endmsg
#define endmsg
Definition: AnalysisConfig_Ntuple.cxx:63
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
G4UA::UserActionToolBase< RadiationMapsMaker >::m_actions
ThreadSpecificUserAction< RadiationMapsMaker > m_actions
Thread-specific storage of the user action.
Definition: UserActionToolBase.h:63
G4UA::RadiationMapsMaker::Config::nBinsz
int nBinsz
Definition: RadiationMapsMaker.h:39
G4UA::RadiationMapsMaker::Config::zMinFull
double zMinFull
Definition: RadiationMapsMaker.h:48
G4UA::RadiationMapsMaker::Config::phiMaxZoom
double phiMaxZoom
Definition: RadiationMapsMaker.h:58
G4UA::RadiationMapsMaker::Config::nBinsr
int nBinsr
Definition: RadiationMapsMaker.h:38
TRT_PAI_physicsConstants::r0
const double r0
electron radius{cm}
Definition: TRT_PAI_physicsConstants.h:20
test_pyathena.parent
parent
Definition: test_pyathena.py:15
G4UA::RadiationMapsMaker::Report::merge
void merge(const Report &maps)
Definition: RadiationMapsMaker.cxx:155
TH2D
Definition: rootspy.cxx:430
TRT::Track::z0
@ z0
Definition: InnerDetector/InDetCalibEvent/TRT_CalibData/TRT_CalibData/TrackInfo.h:63
G4UA::RadiationMapsMaker::Config::elemZMax
int elemZMax
Definition: RadiationMapsMaker.h:85
G4UA::RadiationMapsMakerTool::m_radMapsFileName
std::string m_radMapsFileName
Output Filename for the Radiation Maps.
Definition: RadiationMapsMakerTool.h:49
G4UA::RadiationMapsMaker::Config::rMinFull
double rMinFull
Definition: RadiationMapsMaker.h:42
G4UA::RadiationMapsMaker::Config::logEMaxo
double logEMaxo
Definition: RadiationMapsMaker.h:76
G4UA::RadiationMapsMaker::Config::rMaxZoom
double rMaxZoom
Definition: RadiationMapsMaker.h:44
G4UA::RadiationMapsMaker::Config::logEMaxn
double logEMaxn
Definition: RadiationMapsMaker.h:71
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:195
G4UA::RadiationMapsMaker::Config::thetaMax
double thetaMax
Definition: RadiationMapsMaker.h:66
G4UA::UserActionToolBase< RadiationMapsMaker >::makeAndFillAction
virtual std::unique_ptr< RadiationMapsMaker > makeAndFillAction(G4AtlasUserActions &actionLists)=0
Make the action and push onto the lists.
G4UA::RadiationMapsMaker::Config::logTMin
double logTMin
Definition: RadiationMapsMaker.h:80
G4UA::RadiationMapsMaker::Config::nBinslogEo
int nBinslogEo
Definition: RadiationMapsMaker.h:74
G4UA::RadiationMapsMaker::getReport
const Report & getReport() const
Retrieve my maps.
Definition: RadiationMapsMaker.h:252
python.CaloScaleNoiseConfig.action
action
Definition: CaloScaleNoiseConfig.py:77
Pythia8_RapidityOrderMPI.val
val
Definition: Pythia8_RapidityOrderMPI.py:14
python.CaloScaleNoiseConfig.type
type
Definition: CaloScaleNoiseConfig.py:78
plotting.yearwise_efficiency.xtitle
string xtitle
Definition: yearwise_efficiency.py:38
G4UA::RadiationMapsMaker::Config::posYOnly
bool posYOnly
Definition: RadiationMapsMaker.h:36
if
if(febId1==febId2)
Definition: LArRodBlockPhysicsV0.cxx:569
G4UA::RadiationMapsMaker::Config::nBinsr3d
int nBinsr3d
Definition: RadiationMapsMaker.h:53
declareProperty
#define declareProperty(n, p, h)
Definition: BaseFakeBkgTool.cxx:15
G4UA::RadiationMapsMaker::Config::phiMinZoom
double phiMinZoom
Definition: RadiationMapsMaker.h:57
G4UA::RadiationMapsMaker::Config::nBinsphi3d
int nBinsphi3d
Definition: RadiationMapsMaker.h:55
MCP::ScaleSmearParam::r1
@ r1
G4UA::RadiationMapsMaker::Config::zMaxZoom
double zMaxZoom
Definition: RadiationMapsMaker.h:50
python.compressB64.c
def c
Definition: compressB64.py:93
G4UA::RadiationMapsMaker::Config::nBinsdphi
int nBinsdphi
Definition: RadiationMapsMaker.h:62
G4UA::RadiationMapsMaker::Config::elemZMin
int elemZMin
Definition: RadiationMapsMaker.h:84
python.AutoConfigFlags.msg
msg
Definition: AutoConfigFlags.py:7
fitman.k
k
Definition: fitman.py:528