Belle II Software  release-05-02-19
MaterialScan.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010-2018 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Martin Ritter *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #pragma once
12 
13 #include <framework/core/Module.h>
14 #include <framework/gearbox/Unit.h>
15 
16 #include <string>
17 #include <vector>
18 #include <set>
19 
20 #include <G4UserSteppingAction.hh>
21 #include <G4ThreeVector.hh>
22 
23 #include <TFile.h>
24 class TH1D;
25 class TH2D;
26 
27 #include <boost/optional.hpp>
28 
29 namespace Belle2 {
36  class MaterialScanBase: public G4UserSteppingAction {
37  public:
40  MaterialScanBase(TFile* rootFile, const std::string& name, const std::string& axisLabel):
41  m_rootFile(rootFile), m_name(name), m_axisLabel(axisLabel)
42  {
43  m_rootFile->mkdir(name.c_str());
44  }
46  std::string getName() const { return m_name; }
48  virtual int getNRays() const = 0;
49 
57  virtual bool createNext(G4ThreeVector& origin, G4ThreeVector& direction) = 0;
58 
59  protected:
61  bool checkStep(const G4Step* step);
63  TFile* m_rootFile;
65  std::string m_name;
67  std::string m_axisLabel;
68  private:
70  static constexpr double c_zeroTolerance = 1e-6;
72  static constexpr int c_maxZeroStepsNudge = 10;
74  static constexpr int c_maxZeroStepsKill = 20;
76  int m_zeroSteps {0};
77  };
78 
84  public:
86  struct ScanParams {
88  ScanParams(): nU(0), nV(0), minU(0), maxU(0), minV(0), maxV(0), maxDepth(-1), splitByMaterials(false) {}
90  int nU;
92  int nV;
94  double minU;
96  double maxU;
98  double minV;
100  double maxV;
102  double maxDepth;
104  std::vector<std::string> ignoredMaterials;
107  };
108 
115  MaterialScan2D(TFile* rootFile, const std::string& name, const std::string& axisLabel, const ScanParams& params);
116 
122  bool createNext(G4ThreeVector& origin, G4ThreeVector& direction) override;
123 
125  int getNRays() const override { return m_params.nU * m_params.nV; }
126 
128  void UserSteppingAction(const G4Step* step) override;
129  protected:
135  virtual void getRay(G4ThreeVector& origin, G4ThreeVector& direction) = 0;
136 
140  TH2D* getHistogram(const std::string& name);
141 
146  void fillValue(const std::string& name, double value);
147 
151  double m_curU;
153  double m_stepU;
155  double m_curV;
157  double m_stepV;
159  double m_curDepth;
161  std::map<std::string, std::unique_ptr<TH2D>> m_regions;
162  };
163 
168  public:
175  MaterialScanSpherical(TFile* rootFile, const G4ThreeVector& origin, const ScanParams& params, bool doCosTheta):
176  MaterialScan2D(rootFile, "Spherical", doCosTheta ? "cos(#theta);#phi [deg]" : "#theta [deg];#phi [deg]", params), m_origin(origin),
177  m_doCosTheta(doCosTheta)
178  {
179  if (doCosTheta) {
182  if (m_params.minU > m_params.maxU) std::swap(m_params.minU, m_params.maxU);
184  m_curU = m_params.minU - m_stepU / 2.;
185  }
186  }
187  protected:
189  void getRay(G4ThreeVector& origin, G4ThreeVector& direction) override;
190 
192  G4ThreeVector m_origin;
193 
196  };
197 
206  public:
214  MaterialScanPlanar(TFile* rootFile, const G4ThreeVector& origin, const G4ThreeVector& dirU, const G4ThreeVector& dirV,
215  const ScanParams& params):
216  MaterialScan2D(rootFile, "Planar", "u [cm];v [cm]", params), m_origin(origin), m_dirU(dirU.unit()), m_dirV(dirV.unit()),
217  m_dirW(m_dirU.cross(m_dirV))
218  {
219  }
220  protected:
222  void getRay(G4ThreeVector& origin, G4ThreeVector& direction) override;
223  protected:
225  G4ThreeVector m_origin;
227  G4ThreeVector m_dirU;
229  G4ThreeVector m_dirV;
231  G4ThreeVector m_dirW;
232  };
233 
243  public:
245  MaterialScanRay(TFile* rootFile, const G4ThreeVector& origin, const G4ThreeVector& dir, double opening,
246  int count, double sampleDepth, double maxDepth, bool splitByMaterials,
247  const std::vector<std::string>& ignoredMaterials):
248  MaterialScanBase(rootFile, "Ray", "ray length [cm]; material budget [X_0/cm]"),
249  m_ignoredMaterials(ignoredMaterials.begin(), ignoredMaterials.end()),
250  m_origin(origin), m_dir(dir), m_opening(opening), m_sampleDepth(sampleDepth),
251  m_maxDepth(maxDepth), m_count(std::max(1, count)), m_splitByMaterials(splitByMaterials)
252  {
253  if (m_opening <= 0) m_count = 1;
254  }
255 
257  int getNRays() const override { return m_count + 1; }
259  void UserSteppingAction(const G4Step* step) override;
261  bool createNext(G4ThreeVector& origin, G4ThreeVector& direction) override;
262  private:
266  TH1D* getHistogram(const std::string& name);
272  void fillValue(const std::string& name, double value, double steplength);
273 
275  std::set<std::string> m_ignoredMaterials;
277  std::map<std::string, std::unique_ptr<TH1D>> m_regions;
279  G4ThreeVector m_origin;
281  G4ThreeVector m_dir;
283  double m_opening;
287  double m_maxDepth{0};
289  double m_curDepth{0};
295  double m_scanDepth{0};
297  int m_count;
299  int m_curRay{ -1};
302  };
303 
314  class MaterialScanModule : public Module {
315 
316  public:
317 
318  /* The constructor of the module.
319  * Sets the description and the parameters of the module.
320  */
322 
324  void initialize() override;
325 
327  void beginRun() override;
328  private:
331  G4ThreeVector getAxis(char name);
332 
334  TFile* m_rootFile;
336  std::string m_filename;
338  std::string m_planeName;
344  std::vector<double> m_customPlane;
346  std::vector<double> m_sphericalOrigin;
348  std::vector<double> m_rayOrigin{0, 0, 0};
350  std::vector<std::string> m_rayIgnoredMaterials;
352  boost::optional<std::vector<double>> m_rayDirection;
354  double m_rayTheta{0};
356  double m_rayPhi{0};
358  double m_rayOpening{0};
360  double m_raySampleDepth{0.1};
362  double m_rayMaxDepth{1000};
364  unsigned int m_rayCount{0};
373  bool m_doRay;
376  };
378 }
Belle2::MaterialScan2D::ScanParams::ScanParams
ScanParams()
Default Constructor.
Definition: MaterialScan.h:88
Belle2::MaterialScanModule::m_doRay
bool m_doRay
Perform a material ray scan: scan along a certain direction with a given opening angle and plot mater...
Definition: MaterialScan.h:373
Belle2::MaterialScan2D::getHistogram
TH2D * getHistogram(const std::string &name)
get histogram for a given name, create if needed.
Definition: MaterialScan.cc:127
Belle2::MaterialScanRay::m_opening
double m_opening
Opening angle in radian.
Definition: MaterialScan.h:283
Belle2::MaterialScan2D::fillValue
void fillValue(const std::string &name, double value)
Fill the recorded material budget into the corresponding histogram.
Definition: MaterialScan.cc:140
Belle2::MaterialScanModule::m_raySplitByMaterials
bool m_raySplitByMaterials
If true Split by materials instead of regions when doing ray scan.
Definition: MaterialScan.h:375
Belle2::MaterialScan2D::getRay
virtual void getRay(G4ThreeVector &origin, G4ThreeVector &direction)=0
Get the origin and direction for the next scan particle.
Belle2::MaterialScanRay::m_regions
std::map< std::string, std::unique_ptr< TH1D > > m_regions
Map holding pointers to all created histograms.
Definition: MaterialScan.h:277
Belle2::MaterialScanModule::m_rayMaxDepth
double m_rayMaxDepth
Max depth of the ray before stopping.
Definition: MaterialScan.h:362
Belle2::MaterialScan2D::ScanParams::nU
int nU
Number of rays along u coordinate.
Definition: MaterialScan.h:90
Belle2::MaterialScan2D::m_regions
std::map< std::string, std::unique_ptr< TH2D > > m_regions
Map holding pointers to all created histograms.
Definition: MaterialScan.h:161
Belle2::MaterialScanModule::m_rayOrigin
std::vector< double > m_rayOrigin
Origin of the ray if a ray scan is performed
Definition: MaterialScan.h:348
Belle2::MaterialScanBase::m_zeroSteps
int m_zeroSteps
Count the number of steps with (almost) zero length.
Definition: MaterialScan.h:76
Belle2::MaterialScanModule::m_rayDirection
boost::optional< std::vector< double > > m_rayDirection
Direction of the ray, cannot be set at the same time as Theta and Phi.
Definition: MaterialScan.h:352
Belle2::MaterialScanModule::m_rayCount
unsigned int m_rayCount
Number of rays to shoot (if opening angle is >0)
Definition: MaterialScan.h:364
Belle2::MaterialScanPlanar::MaterialScanPlanar
MaterialScanPlanar(TFile *rootFile, const G4ThreeVector &origin, const G4ThreeVector &dirU, const G4ThreeVector &dirV, const ScanParams &params)
Create a Planar Scan object with the given parameters.
Definition: MaterialScan.h:214
Belle2::MaterialScan2D::ScanParams::maxV
double maxV
Maximum v value to scan.
Definition: MaterialScan.h:100
Belle2::MaterialScanRay::getHistogram
TH1D * getHistogram(const std::string &name)
get histogram for a given name, create if needed.
Definition: MaterialScan.cc:256
Belle2::MaterialScan2D
Base class to create a Material Scan of the detector geometry.
Definition: MaterialScan.h:83
Belle2::MaterialScanModule::beginRun
void beginRun() override
Do the actual scanning.
Definition: MaterialScan.cc:480
Belle2::MaterialScan2D::ScanParams::nV
int nV
Number of rays along v coordinate.
Definition: MaterialScan.h:92
Belle2::MaterialScanModule::m_spherical
MaterialScan2D::ScanParams m_spherical
Scan parameters for the spherical scan.
Definition: MaterialScan.h:340
Belle2::MaterialScanBase::c_maxZeroStepsKill
static constexpr int c_maxZeroStepsKill
maximum number of consecutive zero steps before killing the track
Definition: MaterialScan.h:74
Belle2::MaterialScanRay::m_scanDepth
double m_scanDepth
The first ray does not record any material but just checks for the maximum useful depth to not get a ...
Definition: MaterialScan.h:295
Belle2::MaterialScanRay
MaterialScan implementation to shoot one ray along a defined direction and record the Material as a f...
Definition: MaterialScan.h:242
Belle2::MaterialScan2D::ScanParams::splitByMaterials
bool splitByMaterials
If true, split output by Materials (otherwise by region)
Definition: MaterialScan.h:106
Belle2::MaterialScanPlanar::m_origin
G4ThreeVector m_origin
Origin of the scan plane.
Definition: MaterialScan.h:225
Belle2::MaterialScanSpherical::m_doCosTheta
bool m_doCosTheta
Flag to indicate if polar-angular sampling is uniform in cos(theta) rather than theta.
Definition: MaterialScan.h:195
Belle2::MaterialScanSpherical
Specific implementation of MaterialScan to do Spherical scanning.
Definition: MaterialScan.h:167
Belle2::MaterialScan2D::m_stepV
double m_stepV
Stepsize for the parameter v.
Definition: MaterialScan.h:157
Belle2::MaterialScanModule::m_raySampleDepth
double m_raySampleDepth
Sample depth of the ray: create one bin every x cm.
Definition: MaterialScan.h:360
Belle2::MaterialScanSpherical::getRay
void getRay(G4ThreeVector &origin, G4ThreeVector &direction) override
Create a ray with the current parameter values according to a spherical distribution.
Definition: MaterialScan.cc:196
Belle2::MaterialScanBase::getName
std::string getName() const
Return the name of the scan.
Definition: MaterialScan.h:46
Belle2::MaterialScanBase::m_axisLabel
std::string m_axisLabel
Labels for the coordinate axes.
Definition: MaterialScan.h:67
Belle2::MaterialScanRay::m_curRay
int m_curRay
Current Ray number: 0 = scan for maximum depth, 1..N = record materials.
Definition: MaterialScan.h:299
Belle2::MaterialScanBase::MaterialScanBase
MaterialScanBase(TFile *rootFile, const std::string &name, const std::string &axisLabel)
This is indeed the constructor and it creates a TDirectory in the output root file and sets all varia...
Definition: MaterialScan.h:40
Belle2::MaterialScanRay::m_curDepth
double m_curDepth
Current depth of the current ray.
Definition: MaterialScan.h:289
Belle2::MaterialScan2D::m_stepU
double m_stepU
Stepsize for the parameter u.
Definition: MaterialScan.h:153
Belle2::MaterialScanModule::m_sphericalOrigin
std::vector< double > m_sphericalOrigin
Origin for spherical scan.
Definition: MaterialScan.h:346
Belle2::MaterialScanRay::m_maxDepth
double m_maxDepth
Maximum depth for each ray after which it will be stopped.
Definition: MaterialScan.h:287
Belle2::MaterialScanPlanar::m_dirV
G4ThreeVector m_dirV
v direction of the scan plane
Definition: MaterialScan.h:229
Belle2::MaterialScanRay::m_count
int m_count
Amount of rays to shoot.
Definition: MaterialScan.h:297
Belle2::MaterialScan2D::getNRays
int getNRays() const override
Return the number of rays in this scan.
Definition: MaterialScan.h:125
Belle2::MaterialScanModule::m_rayPhi
double m_rayPhi
Phi direction of the ray if custom direction is not set.
Definition: MaterialScan.h:356
Belle2::MaterialScanBase
Base class for Material Scans.
Definition: MaterialScan.h:36
Belle2::MaterialScanSpherical::m_origin
G4ThreeVector m_origin
Origin for the spherical scan.
Definition: MaterialScan.h:192
Belle2::MaterialScanModule::m_customPlane
std::vector< double > m_customPlane
Custom plane definition if m_planName is "custom".
Definition: MaterialScan.h:344
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::MaterialScan2D::createNext
bool createNext(G4ThreeVector &origin, G4ThreeVector &direction) override
Get the origin and direction for the next scan particle.
Definition: MaterialScan.cc:99
Belle2::MaterialScanModule::m_rayOpening
double m_rayOpening
Opening angle of the ray.
Definition: MaterialScan.h:358
Belle2::MaterialScan2D::UserSteppingAction
void UserSteppingAction(const G4Step *step) override
Record the material budget for each step of the particles.
Definition: MaterialScan.cc:146
Belle2::MaterialScanModule::m_rootFile
TFile * m_rootFile
Pointer to the ROOT file for the histograms.
Definition: MaterialScan.h:334
Belle2::MaterialScan2D::ScanParams::minU
double minU
Minimum u value to scan.
Definition: MaterialScan.h:94
Belle2::MaterialScanModule::m_planeName
std::string m_planeName
Name of the plane to use for scanning.
Definition: MaterialScan.h:338
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::MaterialScanPlanar::m_dirU
G4ThreeVector m_dirU
u direction of the scan plane
Definition: MaterialScan.h:227
Belle2::MaterialScanModule::m_rayIgnoredMaterials
std::vector< std::string > m_rayIgnoredMaterials
Materials ignored when ray scanning.
Definition: MaterialScan.h:350
Belle2::MaterialScanModule::m_doSpherical
bool m_doSpherical
Wether or not to do a spherical scan.
Definition: MaterialScan.h:366
Belle2::MaterialScan2D::MaterialScan2D
MaterialScan2D(TFile *rootFile, const std::string &name, const std::string &axisLabel, const ScanParams &params)
Constructor.
Definition: MaterialScan.cc:80
Belle2::MaterialScanPlanar::getRay
void getRay(G4ThreeVector &origin, G4ThreeVector &direction) override
Create a ray with the current parameter values according to a planar distribution.
Definition: MaterialScan.cc:210
Belle2::MaterialScanBase::m_name
std::string m_name
Name of the scan, will be prefixed to all histogram names.
Definition: MaterialScan.h:65
Belle2::MaterialScanSpherical::MaterialScanSpherical
MaterialScanSpherical(TFile *rootFile, const G4ThreeVector &origin, const ScanParams &params, bool doCosTheta)
Create a Spherical Scan object with the given parameters.
Definition: MaterialScan.h:175
Belle2::MaterialScanPlanar::m_dirW
G4ThreeVector m_dirW
direction perpendicluar to u and v
Definition: MaterialScan.h:231
Belle2::MaterialScan2D::m_curU
double m_curU
Current value of the parametetr u.
Definition: MaterialScan.h:151
Belle2::MaterialScan2D::ScanParams::maxDepth
double maxDepth
Maximum depth of the scan.
Definition: MaterialScan.h:102
Belle2::MaterialScanRay::m_ignoredMaterials
std::set< std::string > m_ignoredMaterials
Materials ignored when scanning.
Definition: MaterialScan.h:275
Belle2::MaterialScanBase::checkStep
bool checkStep(const G4Step *step)
check for stuck tracks by looking at the step length
Definition: MaterialScan.cc:46
Belle2::MaterialScanModule::m_doCosTheta
bool m_doCosTheta
Perform the spherical scan uniform in cos(theta) instead of theta.
Definition: MaterialScan.h:370
Belle2::MaterialScanModule
The MaterialScan module.
Definition: MaterialScan.h:314
Belle2::Unit::deg
static const double deg
degree to radians
Definition: Unit.h:119
Belle2::MaterialScan2D::m_curV
double m_curV
Current value of the parametetr v.
Definition: MaterialScan.h:155
Belle2::MaterialScanRay::m_sampleDepth
double m_sampleDepth
The ray length after which to sample.
Definition: MaterialScan.h:285
Belle2::MaterialScan2D::m_curDepth
double m_curDepth
Tracklength of the current Ray.
Definition: MaterialScan.h:159
Belle2::MaterialScan2D::ScanParams::minV
double minV
Minimum v value to scan.
Definition: MaterialScan.h:98
Belle2::MaterialScanBase::c_maxZeroStepsNudge
static constexpr int c_maxZeroStepsNudge
maximum number of consecutive zero steps before nudging the track along
Definition: MaterialScan.h:72
Belle2::MaterialScanPlanar
Specific implementaion of MaterialScan to scan parallel to a given plane.
Definition: MaterialScan.h:205
Belle2::MaterialScanRay::m_origin
G4ThreeVector m_origin
Origin of the scan.
Definition: MaterialScan.h:279
Belle2::MaterialScanModule::initialize
void initialize() override
Initialize the output file and check the parameters.
Definition: MaterialScan.cc:441
Belle2::MaterialScan2D::ScanParams
Helper struct to Store Parameters of a Scan.
Definition: MaterialScan.h:86
Belle2::MaterialScanRay::UserSteppingAction
void UserSteppingAction(const G4Step *step) override
Record the material budget for each step of the particles.
Definition: MaterialScan.cc:284
Belle2::MaterialScanRay::fillValue
void fillValue(const std::string &name, double value, double steplength)
Fill the recorded material budget into the corresponding histogram.
Definition: MaterialScan.cc:268
Belle2::MaterialScanRay::createNext
bool createNext(G4ThreeVector &origin, G4ThreeVector &direction) override
Implement shooting along the ray.
Definition: MaterialScan.cc:217
Belle2::MaterialScanModule::m_doPlanar
bool m_doPlanar
Wether or not to do a planar scan.
Definition: MaterialScan.h:368
Belle2::MaterialScanBase::getNRays
virtual int getNRays() const =0
Return the number of rays necessary to perform the scan.
Belle2::MaterialScanModule::getAxis
G4ThreeVector getAxis(char name)
Return a vector along the axis with the given name.
Definition: MaterialScan.cc:579
Belle2::MaterialScanModule::m_rayTheta
double m_rayTheta
Theta direction of the ray if custom direction is not set.
Definition: MaterialScan.h:354
Belle2::MaterialScanRay::m_splitByMaterials
bool m_splitByMaterials
If true Split by materials instead of regions.
Definition: MaterialScan.h:301
Belle2::MaterialScanModule::m_planar
MaterialScan2D::ScanParams m_planar
Scan parameters for the planar scan.
Definition: MaterialScan.h:342
Belle2::MaterialScanRay::getNRays
int getNRays() const override
Return the number of rays.
Definition: MaterialScan.h:257
Belle2::MaterialScanBase::m_rootFile
TFile * m_rootFile
Pointer to the root file for the histograms.
Definition: MaterialScan.h:63
Belle2::MaterialScanBase::c_zeroTolerance
static constexpr double c_zeroTolerance
maximum Step length to be considered zero
Definition: MaterialScan.h:70
Belle2::MaterialScanRay::m_dir
G4ThreeVector m_dir
Direction of the ray.
Definition: MaterialScan.h:281
Belle2::MaterialScan2D::ScanParams::ignoredMaterials
std::vector< std::string > ignoredMaterials
Names of ignored Materials.
Definition: MaterialScan.h:104
Belle2::MaterialScan2D::ScanParams::maxU
double maxU
Maximum u value to scan.
Definition: MaterialScan.h:96
Belle2::MaterialScan2D::m_params
ScanParams m_params
Parameters for the scan.
Definition: MaterialScan.h:149
Belle2::MaterialScanModule::m_filename
std::string m_filename
Filename for the ROOT output.
Definition: MaterialScan.h:336
Belle2::MaterialScanBase::createNext
virtual bool createNext(G4ThreeVector &origin, G4ThreeVector &direction)=0
Belle2::MaterialScanBase::createNext() implemention Get the origin and direction for the next scan pa...
Belle2::MaterialScanRay::MaterialScanRay
MaterialScanRay(TFile *rootFile, const G4ThreeVector &origin, const G4ThreeVector &dir, double opening, int count, double sampleDepth, double maxDepth, bool splitByMaterials, const std::vector< std::string > &ignoredMaterials)
Construct a new instance and set all parameters.
Definition: MaterialScan.h:245