Belle II Software development
EclData Class Reference

This class contains data for ECLSimHit's and provides several relevant conversion functions for better event display. More...

#include <EclData.h>

Public Types

enum  EclSubsystem {
  ALL ,
  BARR ,
  FORW ,
  BACKW
}
 Subsystems of ECL: ALL all subsystems BARR barrel only FORW forward endcap only BACKW backward endcap only. More...
 

Public Member Functions

 EclData ()
 Default constructor.
 
 EclData (const EclData &data)
 Copy constructor.
 
EclDataoperator= (const EclData &other)
 Assignment operator: utilizes copy constructor.
 
TTree * getTree ()
 Returns data contained in EclDisplay.
 
int * getEventCounts ()
 Returns array of event counts per crystal[getCrystalsCount()].
 
int * getEventCountsPerCrystal ()
 Alias for GetEventCounts()
 
int getEventCountsMax ()
 
float * getEnergySums ()
 Get array of total energy for each channel in the specified time and event range.
 
float * getEnergySumPerCrystal ()
 Alias for GetEnergySums()
 
float getEnergySumsMax ()
 
float getEnergyTotal ()
 
bool isCrystalInSubsystem (int crystal, EclSubsystem subsys)
 
int getTimeRangeMin ()
 Return min time in time range.
 
int getTimeRangeMax ()
 Return max time in time range.
 
void setTimeRange (int time_min, int time_max, bool do_update=true)
 Display only events in the specified time range.
 
int getEventRangeMin ()
 Return min event number to display.
 
int getEventRangeMax ()
 Return max event number to display.
 
void setEventRange (int ev_min, int ev_max, bool do_update=true)
 Sets event range to (ev_min, ev_max)
 
void setEnergyThreshold (int en_min, int en_max, bool do_update=true)
 Sets energy range to (en_min, en_max).
 
int getTimeMin ()
 
int getTimeMax ()
 
int getLastEventId ()
 
int getChannel (int phi_id, int theta_id)
 Converts (phi_id, theta_id) pair to ECL CellId.
 
int getPhiId (int ch)
 ECL CellId -> phi_id.
 
int getThetaId (int ch)
 ECL CellId -> theta_id.
 
void excludeChannel (int ch, bool do_update=false)
 Excludes specific channel from the count of events and energy.
 
void includeChannel (int ch, bool do_update=false)
 Includes specific channel in the count of events and energy.
 
void loadRootFile (const char *path)
 Load root file containing ECLCalDigit data from the specified path.
 
void update (bool reset_event_ranges=false)
 Update time_min, time_max, event_counts and energy_sums.
 
int addEvent (ECLCalDigit *event, int evtn)
 Add ECLDigit event to inner TTree (m_tree).
 
void fillEnergyHistogram (TH1F *hist, int energy_min, int energy_max, EclSubsystem subsys)
 Fill energy per channel histogram for the specified EclSubsystem (Barrel, forward endcap, backward endcap, all of them).
 
void fillEnergySumHistogram (TH1F *hist, int energy_min, int energy_max, EclSubsystem subsys)
 Fill energy per event histogram for the specified EclSubsystem (Barrel, forward endcap, backward endcap, all of them).
 
void fillTimeHistogram (TH1F *hist, int time_min, int time_max, EclSubsystem subsys)
 Fill time histogram for the specified EclSubsystem (Barrel, forward endcap, backward endcap, all of them).
 

Static Public Member Functions

static int getCrystalCount ()
 Get number of crystals in ECL.
 

Private Member Functions

void cloneFrom (const EclData &other)
 Clone attributes from other EclData.
 
void initVariables ()
 Initialization of arrays.
 
void initEventRanges ()
 Set initial values for time and event range.
 

Private Attributes

TTree * m_tree
 Tree with loaded events.
 
int m_branch_ch
 Tree channel field.
 
double m_branch_energy
 Tree energy branch.
 
double m_branch_time
 Tree time branch.
 
int m_branch_evtn
 Tree event number branch.
 
int * m_event_counts
 Number of events for each crystal.
 
int m_event_count_max
 Max value in event_counts array.
 
float * m_energy_sums
 Sum of energies of every event captured by crystal (MeV).
 
float m_energy_sums_max
 Max value in m_energy_sums array.
 
float m_energy_total
 Total energy for last displayed range of events.
 
int m_last_event_id
 Id of the event with max recorded event number (evtn).
 
int m_time_max
 Maximum possible time.
 
int m_time_range_min
 Time range (min) for display.
 
int m_time_range_max
 Time range (max) for display.
 
int m_ev_range_min
 Events from ev_range_min will be counted in energy_sums and event_counts.
 
int m_ev_range_max
 Events up to ev_range_max will be counted in energy_sums and event_counts.
 
int m_en_range_min
 Lower boundary of energy threshold.
 
int m_en_range_max
 Upper boundary of energy threshold.
 
std::set< int > m_excluded_ch
 Channels which are filtered out from the count of events and energy.
 
std::vector< int > m_event_entry
 This vector holds the position of each entry which starts an event.
 

Static Private Attributes

static const int ring_start_id [70]
 First crystal id in the beginning of i-th ECL ring.
 

Detailed Description

This class contains data for ECLSimHit's and provides several relevant conversion functions for better event display.

Definition at line 31 of file EclData.h.

Member Enumeration Documentation

◆ EclSubsystem

Subsystems of ECL: ALL all subsystems BARR barrel only FORW forward endcap only BACKW backward endcap only.

Definition at line 98 of file EclData.h.

98{ALL, BARR, FORW, BACKW};

Constructor & Destructor Documentation

◆ EclData() [1/2]

EclData ( )

Default constructor.

Definition at line 41 of file EclData.cc.

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
51}
TTree * m_tree
Tree with loaded events.
Definition: EclData.h:34
int m_branch_evtn
Tree event number branch.
Definition: EclData.h:42
void initVariables()
Initialization of arrays.
Definition: EclData.cc:117
void initEventRanges()
Set initial values for time and event range.
Definition: EclData.cc:136
double m_branch_energy
Tree energy branch.
Definition: EclData.h:38
int m_branch_ch
Tree channel field.
Definition: EclData.h:36
double m_branch_time
Tree time branch.
Definition: EclData.h:40

◆ EclData() [2/2]

EclData ( const EclData data)

Copy constructor.

Resets m_excluded_ch upon copy.

Definition at line 53 of file EclData.cc.

54{
57
58 cloneFrom(data);
59}
void cloneFrom(const EclData &other)
Clone attributes from other EclData.
Definition: EclData.cc:79

◆ ~EclData()

~EclData ( )

Definition at line 72 of file EclData.cc.

73{
74 delete m_event_counts;
75 delete m_energy_sums;
76 delete m_tree;
77}
int * m_event_counts
Number of events for each crystal.
Definition: EclData.h:45
float * m_energy_sums
Sum of energies of every event captured by crystal (MeV).
Definition: EclData.h:49

Member Function Documentation

◆ addEvent()

int addEvent ( ECLCalDigit event,
int  evtn 
)

Add ECLDigit event to inner TTree (m_tree).

Parameters
eventECLDigit event
evtnNumber of event.
Returns
If ECLDigit contains incorrect data, negative values are returned. Otherwise, return value is 0.

Definition at line 428 of file EclData.cc.

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}
std::vector< int > m_event_entry
This vector holds the position of each entry which starts an event.
Definition: EclData.h:84
int m_time_max
Maximum possible time.
Definition: EclData.h:59
int m_last_event_id
Id of the event with max recorded event number (evtn).
Definition: EclData.h:56

◆ cloneFrom()

void cloneFrom ( const EclData other)
private

Clone attributes from other EclData.

Definition at line 79 of file EclData.cc.

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}
int m_event_count_max
Max value in event_counts array.
Definition: EclData.h:47
float m_energy_total
Total energy for last displayed range of events.
Definition: EclData.h:53
int m_time_range_min
Time range (min) for display.
Definition: EclData.h:62
int m_en_range_max
Upper boundary of energy threshold.
Definition: EclData.h:79
int m_ev_range_min
Events from ev_range_min will be counted in energy_sums and event_counts.
Definition: EclData.h:70
static int getCrystalCount()
Get number of crystals in ECL.
Definition: EclData.cc:148
std::set< int > m_excluded_ch
Channels which are filtered out from the count of events and energy.
Definition: EclData.h:82
int m_time_range_max
Time range (max) for display.
Definition: EclData.h:65
int m_en_range_min
Lower boundary of energy threshold.
Definition: EclData.h:77
int m_ev_range_max
Events up to ev_range_max will be counted in energy_sums and event_counts.
Definition: EclData.h:74

◆ excludeChannel()

void excludeChannel ( int  ch,
bool  do_update = false 
)

Excludes specific channel from the count of events and energy.

Definition at line 298 of file EclData.cc.

299{
300 m_excluded_ch.insert(_ch);
301 if (do_update)
302 update();
303}
void update(bool reset_event_ranges=false)
Update time_min, time_max, event_counts and energy_sums.
Definition: EclData.cc:368

◆ 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, backward endcap, all of them).

Definition at line 453 of file EclData.cc.

454{
455 int start, end;
456
457 if (m_ev_range_min < 0)
458 return;
460 return;
461
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}
bool isCrystalInSubsystem(int crystal, EclSubsystem subsys)
Definition: EclData.cc:193

◆ fillEnergySumHistogram()

void fillEnergySumHistogram ( TH1F *  hist,
int  energy_min,
int  energy_max,
EclData::EclSubsystem  subsys 
)

Fill energy per event histogram for the specified EclSubsystem (Barrel, forward endcap, backward endcap, all of them).

Definition at line 477 of file EclData.cc.

478{
479 int start, end;
480
481 if (m_ev_range_min < 0)
482 return;
484 return;
485
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}

◆ fillTimeHistogram()

void fillTimeHistogram ( TH1F *  hist,
int  time_min,
int  time_max,
EclData::EclSubsystem  subsys 
)

Fill time histogram for the specified EclSubsystem (Barrel, forward endcap, backward endcap, all of them).

Definition at line 517 of file EclData.cc.

518{
519 int start, end;
520
521 if (m_ev_range_min < 0)
522 return;
524 return;
525
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}

◆ getChannel()

int getChannel ( int  phi_id,
int  theta_id 
)

Converts (phi_id, theta_id) pair to ECL CellId.

Definition at line 273 of file EclData.cc.

274{
275 return ring_start_id[theta_id] + phi_id;
276}
static const int ring_start_id[70]
First crystal id in the beginning of i-th ECL ring.
Definition: EclData.h:89

◆ getCrystalCount()

int getCrystalCount ( )
static

Get number of crystals in ECL.

Definition at line 148 of file EclData.cc.

149{
151}
const int c_NCrystals
Number of crystals.

◆ getEnergySumPerCrystal()

float * getEnergySumPerCrystal ( )

Alias for GetEnergySums()

Definition at line 178 of file EclData.cc.

179{
180 return getEnergySums();
181}
float * getEnergySums()
Get array of total energy for each channel in the specified time and event range.
Definition: EclData.cc:173

◆ getEnergySums()

float * getEnergySums ( )

Get array of total energy for each channel in the specified time and event range.

Definition at line 173 of file EclData.cc.

174{
175 return m_energy_sums;
176}

◆ getEnergySumsMax()

float getEnergySumsMax ( )
Returns
Maximum value in array of energy sums.

Definition at line 183 of file EclData.cc.

184{
185 return m_energy_sums_max;
186}
float m_energy_sums_max
Max value in m_energy_sums array.
Definition: EclData.h:51

◆ getEnergyTotal()

float getEnergyTotal ( )
Returns
Get total energy in currently set range of time and/or events.

Definition at line 188 of file EclData.cc.

189{
190 return m_energy_total;
191}

◆ getEventCounts()

int * getEventCounts ( )

Returns array of event counts per crystal[getCrystalsCount()].

Definition at line 158 of file EclData.cc.

159{
160 return m_event_counts;
161}

◆ getEventCountsMax()

int getEventCountsMax ( )
Returns
Maximum value in array of event counts.

Definition at line 168 of file EclData.cc.

169{
170 return m_event_count_max;
171}

◆ getEventCountsPerCrystal()

int * getEventCountsPerCrystal ( )

Alias for GetEventCounts()

Definition at line 163 of file EclData.cc.

164{
165 return getEventCounts();
166}
int * getEventCounts()
Returns array of event counts per crystal[getCrystalsCount()].
Definition: EclData.cc:158

◆ getEventRangeMax()

int getEventRangeMax ( )

Return max event number to display.

Definition at line 233 of file EclData.cc.

234{
235 return m_ev_range_max;
236}

◆ getEventRangeMin()

int getEventRangeMin ( )

Return min event number to display.

Definition at line 229 of file EclData.cc.

230{
231 return m_ev_range_min;
232}

◆ getLastEventId()

int getLastEventId ( )
Returns
Id of the last event.

Definition at line 268 of file EclData.cc.

269{
270 return m_last_event_id;
271}

◆ getPhiId()

int getPhiId ( int  ch)

ECL CellId -> phi_id.

Definition at line 278 of file EclData.cc.

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}

◆ getThetaId()

int getThetaId ( int  ch)

ECL CellId -> theta_id.

Definition at line 288 of file EclData.cc.

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}

◆ getTimeMax()

int getTimeMax ( )
Returns
Time of the last event in TTrees.

Definition at line 263 of file EclData.cc.

264{
265 return m_time_max;
266}

◆ getTimeMin()

int getTimeMin ( )
Returns
Time of the first event in TTree.

Definition at line 259 of file EclData.cc.

260{
261 return m_time_range_min;
262}

◆ getTimeRangeMax()

int getTimeRangeMax ( )

Return max time in time range.

Definition at line 213 of file EclData.cc.

214{
215 return m_time_range_max;
216}

◆ getTimeRangeMin()

int getTimeRangeMin ( )

Return min time in time range.

Definition at line 209 of file EclData.cc.

210{
211 return m_time_range_min;
212}

◆ getTree()

TTree * getTree ( )

Returns data contained in EclDisplay.

Definition at line 153 of file EclData.cc.

154{
155 return m_tree;
156}

◆ includeChannel()

void includeChannel ( int  ch,
bool  do_update = false 
)

Includes specific channel in the count of events and energy.

Definition at line 305 of file EclData.cc.

306{
307 m_excluded_ch.erase(_ch);
308 if (do_update)
309 update();
310}

◆ initEventRanges()

void initEventRanges ( )
private

Set initial values for time and event range.

Definition at line 136 of file EclData.cc.

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}

◆ initVariables()

void initVariables ( )
private

Initialization of arrays.

Definition at line 117 of file EclData.cc.

118{
119 // Setting
121
122 m_event_counts = new int[getCrystalCount() + 1];
123 m_energy_sums = new float[getCrystalCount() + 1];
124
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}

◆ isCrystalInSubsystem()

bool isCrystalInSubsystem ( int  crystal,
EclData::EclSubsystem  subsys 
)
Parameters
crystalcrystal number
subsysID of the ECL subsystem.
Returns
Whether crystal is in barrel or forward/backward endcap.

Definition at line 193 of file EclData.cc.

194{
195 switch (subsys) {
196 case ALL:
197 return true;
198 case BARR:
199 return ECLElementNumbers::isBarrel(crystal);
200 case FORW:
201 return ECLElementNumbers::isForward(crystal);
202 case BACKW:
203 return ECLElementNumbers::isBackward(crystal);
204 default:
205 return false;
206 }
207}
bool isForward(int cellId)
Check whether the crystal is in forward ECL.
bool isBarrel(int cellId)
Check whether the crystal is in barrel ECL.
bool isBackward(int cellId)
Check whether the crystal is in backward ECL.

◆ loadRootFile()

void loadRootFile ( const char *  path)

Load root file containing ECLCalDigit data from the specified path.

Definition at line 312 of file EclData.cc.

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) {
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}

◆ operator=()

EclData & operator= ( const EclData other)

Assignment operator: utilizes copy constructor.

Definition at line 61 of file EclData.cc.

62{
63 if (this != &other) { // self-assignment check
66
67 cloneFrom(other);
68 }
69 return *this;
70}

◆ setEnergyThreshold()

void setEnergyThreshold ( int  en_min,
int  en_max,
bool  do_update = true 
)

Sets energy range to (en_min, en_max).

Parameters
en_minminimal energy
en_maxmaximal energy
do_updateImmediately update sums of energy and event count for each channel.

Definition at line 248 of file EclData.cc.

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}

◆ setEventRange()

void setEventRange ( int  ev_min,
int  ev_max,
bool  do_update = true 
)

Sets event range to (ev_min, ev_max)

Parameters
ev_minminimal event number
ev_maxmaximal event number
do_updateImmediately update sums of energy and event count for each channel.

Definition at line 237 of file EclData.cc.

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}

◆ setTimeRange()

void setTimeRange ( int  time_min,
int  time_max,
bool  do_update = true 
)

Display only events in the specified time range.

This method currently can't be called from the GUI.

Definition at line 218 of file EclData.cc.

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}

◆ update()

void update ( bool  reset_event_ranges = false)

Update time_min, time_max, event_counts and energy_sums.

Definition at line 368 of file EclData.cc.

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
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 }
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}

Member Data Documentation

◆ m_branch_ch

int m_branch_ch
private

Tree channel field.

Definition at line 36 of file EclData.h.

◆ m_branch_energy

double m_branch_energy
private

Tree energy branch.

Definition at line 38 of file EclData.h.

◆ m_branch_evtn

int m_branch_evtn
private

Tree event number branch.

Definition at line 42 of file EclData.h.

◆ m_branch_time

double m_branch_time
private

Tree time branch.

Definition at line 40 of file EclData.h.

◆ m_en_range_max

int m_en_range_max
private

Upper boundary of energy threshold.

Definition at line 79 of file EclData.h.

◆ m_en_range_min

int m_en_range_min
private

Lower boundary of energy threshold.

Definition at line 77 of file EclData.h.

◆ m_energy_sums

float* m_energy_sums
private

Sum of energies of every event captured by crystal (MeV).

Definition at line 49 of file EclData.h.

◆ m_energy_sums_max

float m_energy_sums_max
private

Max value in m_energy_sums array.

Definition at line 51 of file EclData.h.

◆ m_energy_total

float m_energy_total
private

Total energy for last displayed range of events.

Definition at line 53 of file EclData.h.

◆ m_ev_range_max

int m_ev_range_max
private

Events up to ev_range_max will be counted in energy_sums and event_counts.

Setting time ranges means EclData will filter out excluded events in most of the data acquisition methods.

Definition at line 74 of file EclData.h.

◆ m_ev_range_min

int m_ev_range_min
private

Events from ev_range_min will be counted in energy_sums and event_counts.

Setting time ranges means EclData will filter out excluded events in most of the data acquisition methods.

Definition at line 70 of file EclData.h.

◆ m_event_count_max

int m_event_count_max
private

Max value in event_counts array.

Definition at line 47 of file EclData.h.

◆ m_event_counts

int* m_event_counts
private

Number of events for each crystal.

Definition at line 45 of file EclData.h.

◆ m_event_entry

std::vector<int> m_event_entry
private

This vector holds the position of each entry which starts an event.

Definition at line 84 of file EclData.h.

◆ m_excluded_ch

std::set<int> m_excluded_ch
private

Channels which are filtered out from the count of events and energy.

Definition at line 82 of file EclData.h.

◆ m_last_event_id

int m_last_event_id
private

Id of the event with max recorded event number (evtn).

Definition at line 56 of file EclData.h.

◆ m_time_max

int m_time_max
private

Maximum possible time.

Definition at line 59 of file EclData.h.

◆ m_time_range_max

int m_time_range_max
private

Time range (max) for display.

Setting time ranges means EclData will filter out excluded events in most of the data acquisition methods.

Definition at line 65 of file EclData.h.

◆ m_time_range_min

int m_time_range_min
private

Time range (min) for display.

Setting time ranges means EclData will filter out excluded events in most of the data acquisition methods.

Definition at line 62 of file EclData.h.

◆ m_tree

TTree* m_tree
private

Tree with loaded events.

Definition at line 34 of file EclData.h.

◆ ring_start_id

const int ring_start_id
staticprivate
Initial value:
= {
1, 49, 97, 161, 225, 289, 385, 481, 577, 673,
769, 865, 1009,
1153, 1297, 1441, 1585, 1729, 1873, 2017, 2161, 2305, 2449,
2593, 2737, 2881, 3025, 3169, 3313, 3457, 3601, 3745, 3889,
4033, 4177, 4321, 4465, 4609, 4753, 4897, 5041, 5185, 5329,
5473, 5617, 5761, 5905, 6049, 6193, 6337, 6481, 6625, 6769,
6913, 7057, 7201, 7345, 7489, 7633,
7777, 7921, 8065, 8161, 8257, 8353, 8449, 8545, 8609, 8673,
}

First crystal id in the beginning of i-th ECL ring.

Taken from basf2/ecl/data/ecl_channels_map (two last columns)

Definition at line 89 of file EclData.h.


The documentation for this class was generated from the following files: