GRSISort "v4.0.0.5"
An extension of the ROOT analysis Framework
Loading...
Searching...
No Matches
TGriffinAngles.h
Go to the documentation of this file.
1#ifndef TGRIFFINANGLES_H
2#define TGRIFFINANGLES_H
3
4/** \addtogroup Analysis
5 * @{
6 */
7
8////////////////////////////////////////////////////////////////////////////////
9///
10/// \class TGriffinAngles
11///
12/// This class provides the number of angles, their values, and the index a
13/// specific angle is at, as well as how many combinations contribute to this
14/// angle.
15/// It is meant to be used for the angular correlation analysis.
16///
17////////////////////////////////////////////////////////////////////////////////
18
19#include <set>
20#include <map>
21#include <functional>
22
23#include "TCollection.h"
24#include "TNamed.h"
25#include "TGraphErrors.h"
26
27#include "Globals.h"
28
29class TGriffinAngles : public TNamed {
30public:
31 explicit TGriffinAngles(double distance = 145., bool folding = false, bool grouping = false, bool addback = true);
32 TGriffinAngles(const TGriffinAngles&) = default;
33 TGriffinAngles(TGriffinAngles&&) noexcept = default;
34 TGriffinAngles& operator=(const TGriffinAngles&) = default;
35 TGriffinAngles& operator=(TGriffinAngles&&) noexcept = default;
36 ~TGriffinAngles() = default;
37
38 double Distance() const { return fDistance; }
39 bool Folding() const { return fFolding; }
40 bool Grouping() const { return fGrouping; }
41 bool Addback() const { return fAddback; }
42
43 int Index(double angle);
44 int NumberOfAngles() const { return fAngles.size(); }
45 double Angle(int index) const
46 {
47 auto it = fAngles.begin();
48 std::advance(it, index);
49 return *it;
50 }
51 double AverageAngle(int index) const;
52 int Count(double angle)
53 {
54 /// If the angle is in our map, report how often it exists, otherwise return zero.
55 if(fAngleCount.count(static_cast<int>(std::round(angle / fRounding))) == 1) { return fAngleCount.at(static_cast<int>(std::round(angle / fRounding))); }
56 return 0;
57 }
58
59 void FoldOrGroup(TGraphErrors* z0, TGraphErrors* z2, TGraphErrors* z4, bool verbose = false) const;
60 std::set<double>::iterator begin() const { return fAngles.begin(); }
61 std::set<double>::iterator end() const { return fAngles.end(); }
62
63 bool ExcludeDetector(int detector) const;
64 bool ExcludeCrystal(int detector, int crystal) const;
65
66 void Print(Option_t* = "") const override;
67
68 Long64_t Merge(TCollection* list)
69 {
70 for(auto* griffinAngles : *list) { Add(static_cast<TGriffinAngles*>(griffinAngles)); }
71 return 0;
72 }
73
74 static double Rounding() { return fRounding; }
75 static void Rounding(const double& val) { fRounding = val; }
76
77 static EVerbosity Verbosity() { return fVerbosity; }
78 static void Verbosity(const EVerbosity& val) { fVerbosity = val; }
79
80private:
81 void Add(TGriffinAngles* griffinAngles);
82
83 static EVerbosity fVerbosity; ///< verbosity level
84 double fDistance{145.}; ///< distance of detector from center of array in mmm
85 bool fFolding{false}; ///< flag indicating whether we fold our distribution around 90 degree
86 bool fGrouping{false}; ///< flag indicating whether we group close angles together
87 bool fAddback{true}; ///< flag indicating whether we use addback
88 static double fRounding; ///< we consider any angles whose difference is less than this to be equal
89 std::vector<int> fExcludedDetectors; ///< list of detectors that are excluded in calculating the angles
90 std::vector<int> fExcludedCrystals; ///< list of crystals that are excluded in calculating the angles, the crystals are numbered as 4*(det-1)+cry, so start at 0 and go up to 63
91 std::vector<int> fCustomGrouping; ///< list of custom groups
92 std::set<double> fAngles; ///< set of unique angles, when grouping is used, the largest angle of the group is used!
93 std::map<double, int> fAngleMap; ///< Maps angles to indices. This is fairly straight forward without grouping, but if grouping is used multiple angles can be mapped to the same index.
94 std::map<int, int> fAngleCount; ///< Maps angles (divided by rounding and cast to integers) to number of combinations contributing to it.
95
96 /// \cond CLASSIMP
97 ClassDefOverride(TGriffinAngles, 5) // NOLINT(readability-else-after-return)
98 /// \endcond
99};
100/*! @} */
101#endif
EVerbosity
Definition Globals.h:143
void Add(TGriffinAngles *griffinAngles)
std::set< double > fAngles
set of unique angles, when grouping is used, the largest angle of the group is used!
static void Rounding(const double &val)
TGriffinAngles(TGriffinAngles &&) noexcept=default
Long64_t Merge(TCollection *list)
bool Addback() const
bool Grouping() const
static double fRounding
we consider any angles whose difference is less than this to be equal
void FoldOrGroup(TGraphErrors *z0, TGraphErrors *z2, TGraphErrors *z4, bool verbose=false) const
bool ExcludeDetector(int detector) const
bool fGrouping
flag indicating whether we group close angles together
double Distance() const
std::set< double >::iterator end() const
std::vector< int > fCustomGrouping
list of custom groups
bool Folding() const
void Print(Option_t *="") const override
double fDistance
distance of detector from center of array in mmm
bool fFolding
flag indicating whether we fold our distribution around 90 degree
bool fAddback
flag indicating whether we use addback
int NumberOfAngles() const
TGriffinAngles(double distance=145., bool folding=false, bool grouping=false, bool addback=true)
TGriffinAngles(const TGriffinAngles &)=default
double Angle(int index) const
static void Verbosity(const EVerbosity &val)
int Index(double angle)
std::map< int, int > fAngleCount
Maps angles (divided by rounding and cast to integers) to number of combinations contributing to it.
static EVerbosity Verbosity()
std::vector< int > fExcludedDetectors
list of detectors that are excluded in calculating the angles
std::map< double, int > fAngleMap
Maps angles to indices. This is fairly straight forward without grouping, but if grouping is used mul...
static EVerbosity fVerbosity
verbosity level
std::vector< int > fExcludedCrystals
list of crystals that are excluded in calculating the angles, the crystals are numbered as 4*(det-1)+...
static double Rounding()
std::set< double >::iterator begin() const
bool ExcludeCrystal(int detector, int crystal) const
int Count(double angle)
double AverageAngle(int index) const