Belle II Software  release-05-02-19
EclData.cc
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 //This module
11 #include <ecl/modules/eclDisplay/EclData.h>
12 
13 //ROOT
14 #include <TLeaf.h>
15 #include <TFile.h>
16 
17 //Framework
18 #include <framework/logging/Logger.h>
19 
20 //ECL
21 #include <ecl/dataobjects/ECLCalDigit.h>
22 
23 using namespace Belle2;
24 
25 const int EclData::ring_start_id[70] = {
26  // forward (0-12)
27  1, 49, 97, 161, 225, 289, 385, 481, 577, 673,
28  769, 865, 1009,
29  // barrel (13-58)
30  1153, 1297, 1441, 1585, 1729, 1873, 2017, 2161, 2305, 2449,
31  2593, 2737, 2881, 3025, 3169, 3313, 3457, 3601, 3745, 3889,
32  4033, 4177, 4321, 4465, 4609, 4753, 4897, 5041, 5185, 5329,
33  5473, 5617, 5761, 5905, 6049, 6193, 6337, 6481, 6625, 6769,
34  6913, 7057, 7201, 7345, 7489, 7633,
35  // forward (59-68)
36  7777, 7921, 8065, 8161, 8257, 8353, 8449, 8545, 8609, 8673,
37  // last_crystal+1
39 };
40 
42 {
43  m_tree = new TTree("tree", "tree");
44  m_tree->Branch("ch", &m_branch_ch, "ch/I");
45  m_tree->Branch("energy", &m_branch_energy, "energy/D");
46  m_tree->Branch("time", &m_branch_time, "time/D");
47  m_tree->Branch("evtn", &m_branch_evtn, "evtn/I");
48 
49  initVariables();
51 }
52 
54 {
55  initVariables();
57 
58  cloneFrom(data);
59 }
60 
62 {
63  if (this != &other) { // self-assignment check
64  initVariables();
66 
67  cloneFrom(other);
68  }
69  return *this;
70 }
71 
72 EclData::~EclData()
73 {
74  delete m_event_counts;
75  delete m_energy_sums;
76  delete m_tree;
77 }
78 
79 void EclData::cloneFrom(const EclData& other)
80 {
81  m_tree = other.m_tree->CloneTree();
82 
83  m_tree = new TTree("tree", "tree");
84  m_tree->Branch("ch", &m_branch_ch, "ch/I");
85  m_tree->Branch("energy", &m_branch_energy, "energy/I");
86  m_tree->Branch("time", &m_branch_time, "time/I");
87  m_tree->Branch("evtn", &m_branch_evtn, "evtn/I");
88 
90  m_time_max = other.m_time_max;
91  // m_excluded_ch is not copied.
92  m_excluded_ch.clear();
93 
96 
97  for (int i = 0; i < getCrystalCount() + 1; i++) {
98  m_event_counts[i] = other.m_event_counts[i];
99  m_energy_sums[i] = other.m_energy_sums[i];
100  }
101 
102  // TODO: Vectors surely can be copied in more simple way.
103  for (unsigned int i = 0; i < other.m_event_entry.size(); i++) {
104  m_event_entry.push_back(other.m_event_entry[i]);
105  }
106 
109 
112 
115 }
116 
118 {
119  // Setting
121 
122  m_event_counts = new int[getCrystalCount() + 1];
123  m_energy_sums = new float[getCrystalCount() + 1];
124 
125  m_event_count_max = 0;
126  m_energy_sums_max = 0;
127  m_energy_total = 0;
128 
129  m_last_event_id = -1;
130 
131  m_time_max = 0;
132 
133  m_excluded_ch.clear();
134 }
135 
137 {
138  m_en_range_min = 0;
139  m_en_range_max = -1;
140 
141  m_time_range_min = -2048;
143 
144  m_ev_range_min = 0;
146 }
147 
149 {
150  return 8736;
151 }
152 
154 {
155  return m_tree;
156 }
157 
159 {
160  return m_event_counts;
161 }
162 // Alias for getEventCounts()
164 {
165  return getEventCounts();
166 }
167 
169 {
170  return m_event_count_max;
171 }
172 
174 {
175  return m_energy_sums;
176 }
177 // Alias for getEnergySums()
179 {
180  return getEnergySums();
181 }
182 
184 {
185  return m_energy_sums_max;
186 }
187 
189 {
190  return m_energy_total;
191 }
192 
194 {
195  switch (subsys) {
196  case ALL:
197  return true;
198  case BARR:
199  return crystal >= 1153 && crystal <= 7776;
200  case FORW:
201  return crystal < 1153;
202  case BACKW:
203  return crystal > 7776;
204  default:
205  return false;
206  }
207 }
208 
210 {
211  return m_time_range_min;
212 }
214 {
215  return m_time_range_max;
216 }
217 
218 void EclData::setTimeRange(int time_min, int time_max, bool do_update)
219 {
220  if (m_time_range_min == time_min && m_time_range_max == time_max)
221  return;
222 
223  m_time_range_min = time_min;
224  m_time_range_max = time_max;
225  if (do_update)
226  update();
227 }
228 
230 {
231  return m_ev_range_min;
232 }
234 {
235  return m_ev_range_max;
236 }
237 void EclData::setEventRange(int ev_min, int ev_max, bool do_update)
238 {
239  if (ev_min == m_ev_range_min && ev_max == m_ev_range_max)
240  return;
241 
242  m_ev_range_min = ev_min;
243  m_ev_range_max = ev_max;
244  if (do_update)
245  update();
246 }
247 
248 void EclData::setEnergyThreshold(int en_min, int en_max, bool do_update)
249 {
250  if (en_min == m_en_range_min && en_max == m_en_range_max)
251  return;
252 
253  m_en_range_min = en_min;
254  m_en_range_max = en_max;
255  if (do_update)
256  update();
257 }
258 
260 {
261  return m_time_range_min;
262 }
264 {
265  return m_time_max;
266 }
267 
269 {
270  return m_last_event_id;
271 }
272 
273 int EclData::getChannel(int phi_id, int theta_id)
274 {
275  return ring_start_id[theta_id] + phi_id;
276 }
277 
278 int EclData::getPhiId(int _ch)
279 {
280  for (int i = 0; i < 69; i++) {
281  if (_ch < ring_start_id[i + 1])
282  return _ch - ring_start_id[i];
283  }
284 
285  return -1;
286 }
287 
289 {
290  for (int i = 0; i < 69; i++) {
291  if (_ch < ring_start_id[i + 1])
292  return i;
293  }
294 
295  return -1;
296 }
297 
298 void EclData::excludeChannel(int _ch, bool do_update)
299 {
300  m_excluded_ch.insert(_ch);
301  if (do_update)
302  update();
303 }
304 
305 void EclData::includeChannel(int _ch, bool do_update)
306 {
307  m_excluded_ch.erase(_ch);
308  if (do_update)
309  update();
310 }
311 
312 void EclData::loadRootFile(const char* path)
313 {
314  // TODO: Move this process to some new reset() method.
315  m_excluded_ch.clear();
316  m_event_entry.clear();
317  m_tree->Reset();
318  m_last_event_id = -1;
319  m_en_range_max = -1;
320  m_time_max = 0;
321 
322  TFile* file = new TFile(path, "READ");
323  TTree* tree = 0;
324  file->GetObject("tree", tree);
325  long nentries = tree->GetEntries();
326 
327  tree->SetBranchAddress("ECLCalDigits.m_CellId", &m_branch_ch);
328  tree->SetBranchAddress("ECLCalDigits.m_Energy", &m_branch_energy);
329  tree->SetBranchAddress("ECLCalDigits.m_Time" , &m_branch_time);
330 
331  TLeaf* leafCellId;
332  TLeaf* leafEnergy;
333  TLeaf* leafTimeFit;
334 
335  leafCellId = tree->GetLeaf("ECLCalDigits.m_CellId");
336  leafEnergy = tree->GetLeaf("ECLCalDigits.m_Energy");
337  leafTimeFit = tree->GetLeaf("ECLCalDigits.m_Time");
338 
339  for (long i = 0; i < nentries; i++) {
340  tree->GetEntry(i);
341 
342  const int len = leafCellId->GetLen();
343 
344  if (len > 0) {
345  m_last_event_id++;
347 
348  m_event_entry.push_back(m_tree->GetEntries());
349  }
350 
351  for (int j = 0; j < len; j++) {
352  m_branch_ch = leafCellId ->GetTypedValue<int>(j);
353  m_branch_energy = leafEnergy ->GetTypedValue<double>(j) * 1e3;
354  m_branch_time = leafTimeFit->GetTypedValue<double>(j);
355 
357 
358  B2DEBUG(10, "Added event ECLCalDigits for crystal " << m_branch_ch);
359 
360  // TODO: Add energy cut.
361  m_tree->Fill();
362  }
363  }
364 
365  update(true);
366 }
367 
368 void EclData::update(bool reset_event_ranges)
369 {
370  int start, end;
371 
372  if (reset_event_ranges) initEventRanges();
373 
374  if (m_ev_range_min < 0)
375  return;
377  return;
378 
379  start = m_event_entry[m_ev_range_min];
381  end = m_event_entry[m_ev_range_max + 1];
382  else
383  end = m_tree->GetEntries();
384 
385  for (int i = 1; i <= getCrystalCount(); i++) {
386  m_event_counts[i] = 0;
387  m_energy_sums[i] = 0;
388  }
389  m_event_count_max = 0;
390  m_energy_sums_max = 0;
391  m_energy_total = 0;
392 
393  for (int i = start; i < end; i++) {
394  m_tree->GetEntry(i);
395 
396  if (m_excluded_ch.count(m_branch_ch) > 0) continue;
397 
398  if (m_branch_ch >= 1 && m_branch_ch <= getCrystalCount()) {
399  // Check if current time belongs to time_range, set by user.
400  if (m_time_range_max >= 0)
401  if (m_branch_time < m_time_range_min || m_branch_time > m_time_range_max)
402  continue;
403 
404  if (m_en_range_max >= 0)
405  if (m_branch_energy < m_en_range_min || m_branch_energy > m_en_range_max)
406  continue;
407 
411 
414 
417  }
418  }
419 
420  if (m_time_range_max < 0)
422  if (m_ev_range_max < 0)
424 
425  B2DEBUG(250, end - start << " events handled.");
426 }
427 
428 int EclData::addEvent(ECLCalDigit* event, int _evtn)
429 {
430  if (event->getEnergy() <= 0 || event->getCellId() <= 0) {
431  return -1;
432  }
433 
434  m_branch_ch = event->getCellId();
435  m_branch_energy = event->getEnergy() * 1e3; // To MeV
436  m_branch_time = event->getTime();
437  m_branch_evtn = _evtn;
438 
440  if (m_last_event_id < _evtn) {
441  m_last_event_id = _evtn;
442  m_event_entry.push_back(m_tree->GetEntries());
443  }
444 
445  B2DEBUG(200, "Added event #" << _evtn << ", ch:" << m_branch_ch
446  << ", energy:" << m_branch_energy << ", time:" << m_branch_time);
447 
448  m_tree->Fill();
449 
450  return 0;
451 }
452 
453 void EclData::fillEnergyHistogram(TH1F* hist, int energy_min, int energy_max, EclSubsystem subsys)
454 {
455  int start, end;
456 
457  if (m_ev_range_min < 0)
458  return;
460  return;
461 
462  start = m_event_entry[m_ev_range_min];
464  end = m_event_entry[m_ev_range_max + 1];
465  else
466  end = m_tree->GetEntries();
467 
468  for (int i = start; i < end; i++) {
469  m_tree->GetEntry(i);
470  if (isCrystalInSubsystem(m_branch_ch, subsys)) {
471  if (m_branch_energy >= energy_min && m_branch_energy <= energy_max)
472  hist->Fill(m_branch_energy);
473  }
474  }
475 }
476 
477 void EclData::fillEnergySumHistogram(TH1F* hist, int energy_min, int energy_max, EclData::EclSubsystem subsys)
478 {
479  int start, end;
480 
481  if (m_ev_range_min < 0)
482  return;
484  return;
485 
486  start = m_event_entry[m_ev_range_min];
488  end = m_event_entry[m_ev_range_max + 1];
489  else
490  end = m_tree->GetEntries();
491 
492  m_tree->GetEntry(start);
493  int cur_evtn = m_branch_evtn;
494  int energy_sum = 0;
495 
496  for (int i = start; i < end; i++) {
497  m_tree->GetEntry(i);
498 
499  // After reading all of the data from event cur_evtn,
500  // save the sum in histogram.
501  if (m_branch_evtn > cur_evtn) {
502  if (energy_sum > 0) hist->Fill(energy_sum);
503  energy_sum = 0;
504  cur_evtn = m_branch_evtn;
505  }
506 
507  if (isCrystalInSubsystem(m_branch_ch, subsys)) {
508  if (m_branch_energy >= energy_min && m_branch_energy <= energy_max)
509  energy_sum += m_branch_energy;
510  }
511  }
512 
513  // Add last selected event.
514  if (energy_sum > 0) hist->Fill(energy_sum);
515 }
516 
517 void EclData::fillTimeHistogram(TH1F* hist, int time_min, int time_max, EclData::EclSubsystem subsys)
518 {
519  int start, end;
520 
521  if (m_ev_range_min < 0)
522  return;
524  return;
525 
526  start = m_event_entry[m_ev_range_min];
528  end = m_event_entry[m_ev_range_max + 1];
529  else
530  end = m_tree->GetEntries();
531 
532  for (int i = start; i < end; i++) {
533  m_tree->GetEntry(i);
534 
535  if (m_en_range_max >= 0)
536  if (m_branch_energy < m_en_range_min || m_branch_energy > m_en_range_max)
537  continue;
538 
539  if (isCrystalInSubsystem(m_branch_ch, subsys)) {
540  if (m_branch_time >= time_min && m_branch_time <= time_max)
541  hist->Fill(m_branch_time);
542  }
543  }
544 }
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