Belle II Software  release-05-02-19
EclData.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2015 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Milkail Remnev, Dmitry Matvienko *
7  * *
8  * This software is provided "as is" without any warranty. *
9  ***************************************************************************/
10 
11 #pragma once
12 
13 //STL
14 #include <vector>
15 #include <set>
16 
17 //Root
18 #include <TTree.h>
19 #include <TH1F.h>
20 
21 namespace Belle2 {
27  class ECLCalDigit;
28 
33  class EclData {
34  private:
36  TTree* m_tree;
38  int m_branch_ch;
40  double m_branch_energy;
42  double m_branch_time;
45 
47  int* m_event_counts;
51  float* m_energy_sums;
53  float m_energy_sums_max;
56 
58  int m_last_event_id;
59 
64  int m_time_range_min;
67  int m_time_range_max;
68 
76  int m_ev_range_max;
77 
79  int m_en_range_min;
81  int m_en_range_max;
82 
84  std::set<int> m_excluded_ch;
86  std::vector<int> m_event_entry;
87 
91  static const int ring_start_id[70];
92 
93  public:
100  enum EclSubsystem {ALL, BARR, FORW, BACKW};
101 
105  EclData();
106 
110  EclData(const EclData& data);
114  EclData& operator=(const EclData& other);
115 
116 
117  ~EclData();
118 
119  private:
123  void cloneFrom(const EclData& other);
124 
128  void initVariables();
132  void initEventRanges();
133 
134  public:
138  static int getCrystalCount();
139 
143  TTree* getTree();
144 
148  int* getEventCounts();
156  int getEventCountsMax();
157 
162  float* getEnergySums();
166  float* getEnergySumPerCrystal();
170  float getEnergySumsMax();
174  float getEnergyTotal();
175 
180  bool isCrystalInSubsystem(int crystal, EclSubsystem subsys);
181 
185  int getTimeRangeMin();
189  int getTimeRangeMax();
194  void setTimeRange(int time_min, int time_max, bool do_update = true);
195 
199  int getEventRangeMin();
203  int getEventRangeMax();
209  void setEventRange(int ev_min, int ev_max, bool do_update = true);
210 
216  void setEnergyThreshold(int en_min, int en_max, bool do_update = true);
217 
221  int getTimeMin();
225  int getTimeMax();
226 
230  int getLastEventId();
231 
235  int getChannel(int phi_id, int theta_id);
239  int getPhiId(int ch);
243  int getThetaId(int ch);
244 
248  void excludeChannel(int ch, bool do_update = false);
252  void includeChannel(int ch, bool do_update = false);
253 
257  void loadRootFile(const char* path);
258 
262  void update(bool reset_event_ranges = false);
263 
270  int addEvent(ECLCalDigit* event, int evtn);
275  void fillEnergyHistogram(TH1F* hist, int energy_min, int energy_max, EclSubsystem subsys);
280  void fillEnergySumHistogram(TH1F* hist, int energy_min, int energy_max, EclSubsystem subsys);
285  void fillTimeHistogram(TH1F* hist, int time_min, int time_max, EclSubsystem subsys);
286  };
288 }
Belle2::EclData::loadRootFile
void loadRootFile(const char *path)
Load root file containing ECLCalDigit data from the specified path.
Definition: EclData.cc:312
Belle2::EclData::getChannel
int getChannel(int phi_id, int theta_id)
Converts (phi_id, theta_id) pair to ECL CellId.
Definition: EclData.cc:273
Belle2::EclData::getLastEventId
int getLastEventId()
Definition: EclData.cc:268
Belle2::EclData::getEventRangeMin
int getEventRangeMin()
Return min event number to display.
Definition: EclData.cc:229
Belle2::EclData::initVariables
void initVariables()
Initialization of arrays.
Definition: EclData.cc:117
Belle2::EclData::includeChannel
void includeChannel(int ch, bool do_update=false)
Includes specific channel in the count of events and energy.
Definition: EclData.cc:305
Belle2::ECLCalDigit
Class to store calibrated ECLDigits: ECLCalDigits.
Definition: ECLCalDigit.h:38
Belle2::EclData::fillEnergySumHistogram
void fillEnergySumHistogram(TH1F *hist, int energy_min, int energy_max, EclSubsystem subsys)
Fill energy per event histogram for the specified EclSubsystem (Barrel, forward endcap,...
Definition: EclData.cc:477
Belle2::EclData::getEnergySumPerCrystal
float * getEnergySumPerCrystal()
Alias for GetEnergySums()
Definition: EclData.cc:178
Belle2::EclData::m_event_entry
std::vector< int > m_event_entry
This vector holds the position of each entry which starts an event.
Definition: EclData.h:94
Belle2::EclData::update
void update(bool reset_event_ranges=false)
Update time_min, time_max, event_counts and energy_sums.
Definition: EclData.cc:368
Belle2::EclData::initEventRanges
void initEventRanges()
Set initial values for time and event range.
Definition: EclData.cc:136
Belle2::EclData::getEventCounts
int * getEventCounts()
Returns array of event counts per crystal[getCrystalsCount()].
Definition: EclData.cc:158
Belle2::EclData::m_event_count_max
int m_event_count_max
Max value in event_counts array.
Definition: EclData.h:57
Belle2::EclData::m_time_range_max
int m_time_range_max
Time range (max) for display.
Definition: EclData.h:75
Belle2::EclData::addEvent
int addEvent(ECLCalDigit *event, int evtn)
Add ECLDigit event to inner TTree (m_tree).
Definition: EclData.cc:428
Belle2::EclData::m_ev_range_min
int m_ev_range_min
Events from ev_range_min will be counted in energy_sums and event_counts.
Definition: EclData.h:80
Belle2::EclData::excludeChannel
void excludeChannel(int ch, bool do_update=false)
Excludes specific channel from the count of events and energy.
Definition: EclData.cc:298
Belle2::EclData::m_tree
TTree * m_tree
Tree with loaded events.
Definition: EclData.h:44
Belle2::EclData::operator=
EclData & operator=(const EclData &other)
Assignment operator: utilizes copy constructor.
Definition: EclData.cc:61
Belle2::EclData::getPhiId
int getPhiId(int ch)
ECL CellId -> phi_id.
Definition: EclData.cc:278
Belle2::EclData::m_ev_range_max
int m_ev_range_max
Events up to ev_range_max will be counted in energy_sums and event_counts.
Definition: EclData.h:84
Belle2::EclData::setEnergyThreshold
void setEnergyThreshold(int en_min, int en_max, bool do_update=true)
Sets energy range to (en_min, en_max).
Definition: EclData.cc:248
Belle2::EclData::getTimeRangeMin
int getTimeRangeMin()
Return min time in time range.
Definition: EclData.cc:209
Belle2::EclData::m_energy_sums
float * m_energy_sums
Sum of energies of every event captured by crystal (MeV).
Definition: EclData.h:59
Belle2::EclData::EclSubsystem
EclSubsystem
Subsystems of ECL: ALL all subsystems BARR barrel only FORW forward endcap only BACKW backward endcap...
Definition: EclData.h:108
Belle2::EclData::m_energy_total
float m_energy_total
Total energy for last displayed range of events.
Definition: EclData.h:63
Belle2::EclData::m_event_counts
int * m_event_counts
Number of events for each crystal.
Definition: EclData.h:55
Belle2::EclData::m_branch_evtn
int m_branch_evtn
Tree event number branch.
Definition: EclData.h:52
Belle2::EclData::cloneFrom
void cloneFrom(const EclData &other)
Clone attributes from other EclData.
Definition: EclData.cc:79
Belle2::EclData::fillEnergyHistogram
void fillEnergyHistogram(TH1F *hist, int energy_min, int energy_max, EclSubsystem subsys)
Fill energy per channel histogram for the specified EclSubsystem (Barrel, forward endcap,...
Definition: EclData.cc:453
Belle2::EclData::getTimeRangeMax
int getTimeRangeMax()
Return max time in time range.
Definition: EclData.cc:213
Belle2::EclData::EclData
EclData()
Default constructor.
Definition: EclData.cc:41
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::EclData::m_en_range_max
int m_en_range_max
Upper boundary of energy threshold.
Definition: EclData.h:89
Belle2::EclData::m_energy_sums_max
float m_energy_sums_max
Max value in m_energy_sums array.
Definition: EclData.h:61
Belle2::EclData::setTimeRange
void setTimeRange(int time_min, int time_max, bool do_update=true)
Display only events in the specified time range.
Definition: EclData.cc:218
Belle2::EclData::isCrystalInSubsystem
bool isCrystalInSubsystem(int crystal, EclSubsystem subsys)
Definition: EclData.cc:193
Belle2::EclData::getCrystalCount
static int getCrystalCount()
Get number of crystals in ECL.
Definition: EclData.cc:148
Belle2::EclData::m_branch_ch
int m_branch_ch
Tree channel field.
Definition: EclData.h:46
Belle2::EclData::m_branch_time
double m_branch_time
Tree time branch.
Definition: EclData.h:50
Belle2::EclData::setEventRange
void setEventRange(int ev_min, int ev_max, bool do_update=true)
Sets event range to (ev_min, ev_max)
Definition: EclData.cc:237
Belle2::EclData::getEventCountsMax
int getEventCountsMax()
Definition: EclData.cc:168
Belle2::EclData::getThetaId
int getThetaId(int ch)
ECL CellId -> theta_id.
Definition: EclData.cc:288
Belle2::EclData::fillTimeHistogram
void fillTimeHistogram(TH1F *hist, int time_min, int time_max, EclSubsystem subsys)
Fill time histogram for the specified EclSubsystem (Barrel, forward endcap, backward endcap,...
Definition: EclData.cc:517
Belle2::EclData::m_en_range_min
int m_en_range_min
Lower boundary of energy threshold.
Definition: EclData.h:87
Belle2::EclData
This class contains data for ECLSimHit's and provides several relevant conversion functions for bette...
Definition: EclData.h:41
Belle2::EclData::getTimeMin
int getTimeMin()
Definition: EclData.cc:259
Belle2::EclData::getTimeMax
int getTimeMax()
Definition: EclData.cc:263
Belle2::EclData::getEventCountsPerCrystal
int * getEventCountsPerCrystal()
Alias for GetEventCounts()
Definition: EclData.cc:163
Belle2::EclData::m_last_event_id
int m_last_event_id
Id of the event with max recorded event number (evtn).
Definition: EclData.h:66
Belle2::EclData::m_time_max
int m_time_max
Maximum possible time.
Definition: EclData.h:69
Belle2::EclData::ring_start_id
static const int ring_start_id[70]
First crystal id in the beginning of i-th ECL ring.
Definition: EclData.h:99
Belle2::EclData::m_excluded_ch
std::set< int > m_excluded_ch
Channels which are filtered out from the count of events and energy.
Definition: EclData.h:92
Belle2::EclData::getTree
TTree * getTree()
Returns data contained in EclDisplay.
Definition: EclData.cc:153
Belle2::EclData::getEnergyTotal
float getEnergyTotal()
Definition: EclData.cc:188
Belle2::EclData::getEventRangeMax
int getEventRangeMax()
Return max event number to display.
Definition: EclData.cc:233
Belle2::EclData::getEnergySumsMax
float getEnergySumsMax()
Definition: EclData.cc:183
Belle2::EclData::getEnergySums
float * getEnergySums()
Get array of total energy for each channel in the specified time and event range.
Definition: EclData.cc:173
Belle2::EclData::m_time_range_min
int m_time_range_min
Time range (min) for display.
Definition: EclData.h:72
Belle2::EclData::m_branch_energy
double m_branch_energy
Tree energy branch.
Definition: EclData.h:48