Belle II Software  release-05-02-19
TrgEclFAMFit.cc
1 //---------------------------------------------------------------
2 // $Id$
3 //---------------------------------------------------------------
4 // Filename : TrgEclFAM.cc
5 // Section : TRG ECL
6 // Owner : InSoo Lee/Yuuji Unno
7 // Email : islee@hep.hanyang.ac.kr / yunno@post.kek.jp
8 //---------------------------------------------------------------
9 // Description : A class to represent TRG ECL
10 //---------------------------------------------------------------
11 // $Log$
12 // 2017-02-21 : version 1.0
13 //---------------------------------------------------------------
14 #include <framework/datastore/StoreArray.h>
15 
16 #include <framework/logging/Logger.h>
17 #include "trg/ecl/TrgEclFAMFit.h"
18 
19 #include "trg/ecl/dataobjects/TRGECLHit.h"
20 #include "trg/ecl/dataobjects/TRGECLDigi0.h"
21 #include "trg/ecl/dataobjects/TRGECLFAMAna.h"
22 
23 #include <stdlib.h>
24 #include <iostream>
25 #include <math.h>
26 
27 using namespace std;
28 using namespace Belle2;
29 //using namespace TRG;
30 //
31 //
32 //
33 TrgEclFAMFit::TrgEclFAMFit(): _BeamBkgTag(0), _AnaTag(0), EventId(0) //, bin(0)
34 {
35 
36  CoeffSigPDF0.clear();
37  CoeffSigPDF1.clear();
38  CoeffNoise31.clear();
39  CoeffNoise32.clear();
40  CoeffNoise33.clear();
41 
42  _TCMap = new TrgEclMapping();
43  _DataBase = new TrgEclDataBase();
44 
45  Threshold.clear();
46 
47  TCFitEnergy.clear();
48  TCFitTiming.clear();
49  TCRawEnergy.clear();
50  TCRawTiming.clear();
51  BeamBkgTag.clear();
52  TCLatency.clear();
53 
54  Threshold.resize(576, 100.0);
55  TCFitEnergy.resize(576);
56  TCFitTiming.resize(576);
57  TCRawEnergy.resize(576);
58  TCRawTiming.resize(576);
59  BeamBkgTag.resize(576);
60 }
61 //
62 //
63 //
65 {
66 
67  delete _TCMap;
68  delete _DataBase;
69 }
70 //
71 //
72 //
73 
74 void
75 TrgEclFAMFit::setup(int eventId)
76 {
77  //
78  // prepare coefficient for fitting
79  //
80 
81  _DataBase -> getCoeffSigPDF(CoeffSigPDF0 , CoeffSigPDF1);
82  _DataBase -> getCoeffNoise(0, CoeffNoise31 , CoeffNoise32, CoeffNoise33);
83 
84  EventId = eventId;
85  //
86  return;
87 }
88 
89 void
90 TrgEclFAMFit::FAMFit01(std::vector<std::vector<double>> digiEnergy, std::vector<std::vector<double>> digiTiming)
91 {
92  //============================
93  // In this function,
94  // o Energy unit must be [MeV/c2]
95  // o Time unit must be [us]
96  // but
97  // o return energy unit must be [GeV/c2]
98  // o return time unit must be [ns]
99  //============================
100  double* TCFitSample = new double [12]; // MeV/c2
101  double* preped = new double [4];
102 
103  int nbin_pedestal = 4;
104  int fam_sampling_interval = 125;
105  int NSampling = 64;// 8 [us] / 125 [ns]
106 
107  int pedFlag = 0;
108  double CoeffAAA = 0;
109  double CoeffBBB = 0;
110  double CoeffPPP = 0;
111  int dTBin = 0;
112  int ShiftdTBin = 0;
113  int Nsmalldt = 10;
114  int SmallOffset = 1;
115  double IntervaldT = 125 * 0.001 / Nsmalldt;
116  // double EThreshold = _Threshold; //[MeV]
117  int FitSleepCounter = 100; // counter to suspend fit
118  int FitSleepThreshold = 2; // # of clk to suspend fit
119  /* cppcheck-suppress variableScope */
120  double FitE;
121  /* cppcheck-suppress variableScope */
122  double FitT;
123 
124  for (int iTCIdm = 0; iTCIdm < 576; iTCIdm++) {
125 
126  FitE = 0;
127  FitT = 0;
128  for (int iShift = 20; iShift < (NSampling - 12); iShift++) { // In order to avoid BKG shaping effect, iShift start from 20.
129 
130  FitSleepCounter++;
131  if (FitSleepCounter <= FitSleepThreshold) {continue;}
132  for (int iFitSample = 0; iFitSample < 12; iFitSample++) {
133  int iReplace = iFitSample + iShift;
134  TCFitSample[iFitSample] = digiEnergy[iTCIdm][iReplace] * 1000.0;
135  if (0) {
136  if (pedFlag == 1 && iFitSample < 4) {
137  TCFitSample[iFitSample] = preped[iFitSample];
138  pedFlag = 0;
139  }
140  }
141  }
142  //
143  //
144  //
145 
146  dTBin = (int)(ShiftdTBin + Nsmalldt);
147  if (dTBin < 1) {dTBin = 1;}
148  if (dTBin > 20) {dTBin = 20;}
149 
150  CoeffAAA = 0;
151  CoeffBBB = 0;
152  CoeffPPP = 0;
153  for (int iFitSample = 0; iFitSample < 12; iFitSample++) {
154  CoeffAAA += CoeffNoise31[dTBin - 1][iFitSample] * TCFitSample[iFitSample];
155  CoeffBBB += CoeffNoise32[dTBin - 1][iFitSample] * TCFitSample[iFitSample];
156  CoeffPPP += CoeffNoise33[dTBin - 1][iFitSample] * TCFitSample[iFitSample];
157  }
158  double deltaT = CoeffBBB / CoeffAAA; // deltaT [us]
159 
160  ShiftdTBin = int(deltaT / IntervaldT + dTBin);
161 
162  FitE = CoeffAAA;
163 
164  //-------
165  // Require "expected time" is around middle of table = Nsmalldt.
166  //-------
167  double condition_t = -(deltaT + dTBin * IntervaldT - fam_sampling_interval * 0.001);
168 
169  if (fabs(condition_t) < 0.8 * (fam_sampling_interval * 0.001) && FitE > Threshold[iTCIdm]) {
170  FitT = condition_t + (SmallOffset + iShift + nbin_pedestal - 5) * (fam_sampling_interval * 0.001);
171 
172 
173  pedFlag = 1;
174  double rand_sampling_correction = digiTiming[iTCIdm][iShift] + (nbin_pedestal - iShift + 32) * fam_sampling_interval;
175 
176  TCFitEnergy[iTCIdm].push_back(FitE / 1000.0); // [GeV/c2]
177  TCFitTiming[iTCIdm].push_back(FitT * 1000 - 4000 + (_DataBase->GetTCFLatency(iTCIdm + 1)) + rand_sampling_correction);
178  FitSleepCounter = 0;
179  ShiftdTBin = 0;
180 
181  }
182  }
183  }
184  if (_BeamBkgTag == 1 || _AnaTag == 1) {
185  SetBeamBkgTag();
186  }
187  //
188  //
189  //
190  delete [] TCFitSample;
191  delete [] preped;
192  return;
193 }
194 //
195 //
196 //
197 void
198 TrgEclFAMFit::FAMFit02(std::vector<std::vector<double>> TCDigiE, std::vector<std::vector<double>> TCDigiT)
199 {
200 
201  int NSampling = 64;
202 
203  //==================
204  // Peak search
205  //==================
206  //@ T_a and T_b is time at sampling points in which 0.6*E exists.
207  int ta_id[20] = {1000}; //@ id of T_a
208  double ttt_a[20] = {0}; //@ time of T_a
209  double ttt_b[20] = {0}; //@ time of T_b
210  for (int iTCIdm = 0; iTCIdm < 576; iTCIdm++) {
211  int noutput = 0;
212 
213  // double threshold = _Threshold * 0.001; //GeV
214  int maxId[500] = {0};
215  for (int iii = 0 ; iii < 500 ; iii++) {
216 
217  maxId[iii] = 0;
218  }
219  int count_up = 0;
220  int count_down = 0;
221  int flag_up = 3;
222  int flag_down = 3;
223 
224  double max = 0;
225  for (int iSampling = 1; iSampling < NSampling; iSampling++) {
226  //-------------------------------------------------------------------------
227  //Peak finding Method 1
228  //------------------------------------------------------------------------
229 
230  if (TCDigiE[iTCIdm][iSampling] >= max) {
231 
232  max = TCDigiE[iTCIdm][iSampling];
233  maxId[noutput] = iSampling;
234  count_up ++;
235  count_down = 0;
236  } else {
237  count_down++;
238  if (count_down >= flag_down) {
239  if (count_up >= flag_up) {
240  if (Threshold[iTCIdm] * 0.001 < max) {
241  max = 0;
242  count_up = 0;
243  count_down = 0;
244 
245  double NoiseLevel = 0;
246  double NoiseCount = 0;
247  for (int iNoise = 0; iNoise < 5; iNoise++) {
248  int iNoiseReplace = (maxId[noutput] - 10) + iNoise;
249  if (iNoiseReplace >= 0) {
250  NoiseLevel += TCDigiE[iTCIdm][iNoiseReplace];
251  NoiseCount++;
252  }
253  }
254  if (NoiseCount != 0) { NoiseLevel /= NoiseCount; }
255  //** Peak point is the Energy */
256  TCFitEnergy[iTCIdm].push_back(TCDigiE[iTCIdm][maxId[noutput]] - NoiseLevel);
257  if (!(maxId[noutput] - 1)) {
258  for (int jSampling = 1; jSampling < maxId[noutput] + 3; jSampling++) {
259  TCDigiE[iTCIdm][jSampling] -= NoiseLevel;
260  }
261  } else {
262  for (int jSampling = maxId[noutput] - 1; jSampling < maxId[noutput] + 3; jSampling++) {
263  TCDigiE[iTCIdm][jSampling] -= NoiseLevel;
264  }
265  }
266  //@ Search T_a ID
267  for (int iSearch = 0; iSearch < 5; iSearch++) {
268 
269  if (TCDigiE[iTCIdm][maxId[noutput] - iSearch] > 0.6 * TCFitEnergy[iTCIdm][noutput] &&
270  TCDigiE[iTCIdm][maxId[noutput] - iSearch - 1] < 0.6 * TCFitEnergy[iTCIdm][noutput]) {
271  ta_id[noutput] = maxId[noutput] - iSearch - 1;
272  }
273  }
274 
275  //@ Estimate timing of t0
276  if (ta_id[noutput] == 1000) {
277  printf("TrgEclFAMFit::digi02> Cannot find TC Timing (TCId=%5i, E=%8.5f)!!!\n", iTCIdm - 1, TCFitEnergy[iTCIdm][0]);
278  B2ERROR("TrgEclFAMFit::digi02> Cannot find TC Timing");
279  } else {
280  ttt_a[noutput] = TCDigiT[iTCIdm][ta_id[noutput]];
281  ttt_b[noutput] = TCDigiT[iTCIdm][ta_id[noutput] + 1];
282  TCFitTiming[iTCIdm].push_back((ttt_a[noutput] +
283  (0.6 * TCFitEnergy[iTCIdm][noutput] - TCDigiE[iTCIdm][ta_id[noutput]]) * (ttt_b[noutput] -
284  ttt_a[noutput])
285  / (TCDigiE[iTCIdm][ta_id[noutput] + 1] - TCDigiE[iTCIdm][ta_id[noutput]])) - (278.7 + 2) + (_DataBase->GetTCFLatency(iTCIdm + 1)));
286  //@ time between t0 and 0.6*peak_energy
287  //@ Alex's number = 274.4 (how he got this value ?)
288  //@ by my check = 278.7 [ns]
289  //@ here "+2" is a shift due to imperfectness of no-fit method.
290  }
291 
292  }
293  }
294  }
295  }
296  }
297  }
298  if (_BeamBkgTag == 1 || _AnaTag == 1) {
299  SetBeamBkgTag();
300  }
301 
302  //
303  //
304  //
305  return;
306 
307 }
308 //
309 //
310 //
311 void
312 TrgEclFAMFit::FAMFit03(std::vector<std::vector<double>> TCDigiEnergy, std::vector<std::vector<double>> TCDigiTiming)
313 {
314  //===============
315  // (03)Signal digitization (w/ 12ns interval for method-0)
316  //===============
317  // double cut_energy_tot = 0.03; // [GeV]
318  // int nbin_pedestal = 100;
319  // float fam_sampling_interval = 12; // [ns]
320  int NSampling = 666;
321 
322  //==================
323  // (03)Peak search
324  //==================
325  float max_shape_time = 563.48; // [ns], time between peak of PDF and t0.
326  // double threshold = _Threshold * 0.001; //GeV
327  for (int iTCIdm = 0; iTCIdm < 576; iTCIdm++) {
328  int noutput = 0;
329  int maxId[500] = {0};
330  int count_up = 0;
331  int count_down = 0;
332  int flag_up = 30;;
333  int flag_down = 40;
334  double max = 0;
335  for (int iSampling = 1; iSampling < NSampling; iSampling++) {
336 
337  if (TCDigiEnergy[iTCIdm][iSampling] >= max) {
338 
339  max = TCDigiEnergy[iTCIdm][iSampling];
340  maxId[noutput] = iSampling;
341  count_up ++;
342  count_down = 0;
343  } else {
344  count_down++;
345  if (count_down >= flag_down) {
346  if (count_up >= flag_up) {
347  if (Threshold[iTCIdm] * 0.001 < max) {
348  max = 0;
349  count_up = 0;
350  count_down = 0;
351  //@ Remove noise effect
352  float NoiseLevel = 0;
353  float NoiseCount = 0;
354  for (int iNoise = 0; iNoise < 42; iNoise++) {
355  int iNoiseReplace = (maxId[noutput] - 88) + iNoise;
356  if (iNoiseReplace >= 0) {
357  NoiseLevel += TCDigiEnergy[iTCIdm][iNoiseReplace];
358  NoiseCount++;
359  }
360  }
361  if (NoiseCount != 0) { NoiseLevel /= NoiseCount; }
362  TCFitEnergy[iTCIdm].push_back(TCDigiEnergy[iTCIdm][maxId[noutput]] - NoiseLevel);
363  TCFitTiming[iTCIdm].push_back(TCDigiTiming[iTCIdm][maxId[noutput]] - max_shape_time + (_DataBase->GetTCFLatency(iTCIdm)));
364  noutput++;
365  }
366  }
367  }
368  }
369  }
370  }
371 
372  if (_BeamBkgTag == 1 || _AnaTag == 1) {
373  SetBeamBkgTag();
374  }
375 
376 
377  return;
378 }
379 //
380 //
381 //
382 void
384 {
385  std::vector<int> TCId;
386  std::vector<double> RawTCTiming;
387  std::vector<double> RawTCEnergy;
388  std::vector<double> RawBeamBkgTag;
389  RawBeamBkgTag.clear();
390  TCId.clear();
391  RawTCTiming.clear();
392  RawTCEnergy.clear();
393  // BeamBkgTag.resize(576 ,std::vector<int> (size,100));
394  // TCRawEnergy.resize(576,std::vector<double>(size,0.));
395  // TCRawTiming.resize(576,std::vector<double>(size,0.));
396  BeamBkgTag.resize(576);
397  TCRawEnergy.resize(576);
398  TCRawTiming.resize(576);
399 
400 
401 
402  StoreArray<TRGECLDigi0> trgeclDigiArray;
403  for (int ii = 0; ii < trgeclDigiArray.getEntries(); ii++) {
404  TRGECLDigi0* aTRGECLDigi = trgeclDigiArray[ii];
405  int eventid = aTRGECLDigi->getEventId();
406  if (EventId != eventid) {continue;}
407  TCId.push_back(aTRGECLDigi->getTCId());
408  RawTCTiming.push_back(aTRGECLDigi -> getRawTiming());
409  RawTCEnergy.push_back(aTRGECLDigi -> getRawEnergy());
410  RawBeamBkgTag.push_back(aTRGECLDigi -> getBeamBkgTag());
411  }
412 
413  for (int iTCId = 0; iTCId < 576 ; iTCId ++) {
414  const int hitsize = TCFitEnergy[iTCId].size();
415  for (int iHit = 0; iHit < hitsize; iHit++) {
416  const int rawsize = TCId.size();
417  for (int iDigi = 0; iDigi < rawsize; iDigi++) {
418  if (TCId[iDigi] != (iTCId + 1)) {continue;}
419  if (abs(TCFitEnergy[iTCId][iHit] - RawTCEnergy[iDigi]) < 12) {
420  if (abs(TCFitTiming[iTCId][iHit] - RawTCTiming[iDigi]) < 30) {
421  BeamBkgTag[iTCId].push_back(RawBeamBkgTag[iDigi]);
422  TCRawEnergy[iTCId].push_back(RawTCEnergy[iDigi]);
423  TCRawTiming[iTCId].push_back(RawTCTiming[iDigi]);
424  }
425  }
426  }
427  }
428  }
429 
430 
431 
432 }
433 
434 //
435 //
436 //
437 void
438 TrgEclFAMFit::save(int m_nEvent)
439 {
440  //---------------
441  // Root Output
442  //---------------
443  int m_hitNum = 0;
444 
445  for (int iTCIdm = 0; iTCIdm < 576; iTCIdm++) {
446  const int hitsize = TCFitEnergy[iTCIdm].size();
447  for (int iHit = 0; iHit < hitsize; iHit++) {
448  StoreArray<TRGECLHit> TrgEclHitArray;
449  TrgEclHitArray.appendNew();
450  m_hitNum = TrgEclHitArray.getEntries() - 1;
451  TrgEclHitArray[m_hitNum]->setEventId(m_nEvent);
452  TrgEclHitArray[m_hitNum]->setTCId(iTCIdm + 1);
453  TrgEclHitArray[m_hitNum]->setEnergyDep(TCFitEnergy[iTCIdm][iHit]);
454  TrgEclHitArray[m_hitNum] ->setTimeAve(TCFitTiming[iTCIdm][iHit]);
455  if (_BeamBkgTag == 1) {
456  TrgEclHitArray[m_hitNum]->setBeamBkgTag(BeamBkgTag[iTCIdm][iHit]);
457  }
458  }
459  }
460 
461 
462  m_hitNum = 0;
463  if (_AnaTag == 1) {
464  for (int iTCIdm = 0; iTCIdm < 576; iTCIdm++) {
465  if (TCFitEnergy[iTCIdm].size() != TCRawEnergy[iTCIdm].size()) {continue;}
466  const int hitsize = TCFitEnergy[iTCIdm].size();
467  for (int iHit = 0; iHit < hitsize; iHit++) {
468  StoreArray<TRGECLFAMAna> TrgEclAnaArray;
469  TrgEclAnaArray.appendNew();
470  m_hitNum = TrgEclAnaArray.getEntries() - 1;
471  TrgEclAnaArray[m_hitNum]->setEventId(m_nEvent);
472  TrgEclAnaArray[m_hitNum]->setTCId(iTCIdm + 1);
473  TrgEclAnaArray[m_hitNum]->setPhiId(_TCMap->getTCPhiIdFromTCId(iTCIdm + 1));
474  TrgEclAnaArray[m_hitNum]->setThetaId(_TCMap->getTCThetaIdFromTCId(iTCIdm + 1));
475  TrgEclAnaArray[m_hitNum]->setRawEnergy(TCRawEnergy[iTCIdm][iHit]);
476  TrgEclAnaArray[m_hitNum]->setRawTiming(TCRawTiming[iTCIdm][iHit]);
477  TrgEclAnaArray[m_hitNum]->setFitEnergy(TCFitEnergy[iTCIdm][iHit]);
478  TrgEclAnaArray[m_hitNum]->setFitTiming(TCFitTiming[iTCIdm][iHit]);
479  TrgEclAnaArray[m_hitNum]->setBeamBkgTag(BeamBkgTag[iTCIdm][iHit]);
480  }
481  }
482  }
483 
484  return;
485 }
Belle2::TrgEclDataBase::GetTCFLatency
double GetTCFLatency(int)
TC flight time latency
Definition: TrgEclDataBase.cc:234
Belle2::StoreArray::appendNew
T * appendNew()
Construct a new T object at the end of the array.
Definition: StoreArray.h:256
Belle2::TRGECLDigi0::getTCId
int getTCId() const
Get TC id.
Definition: TRGECLDigi0.h:62
Belle2::TrgEclFAMFit::FAMFit03
void FAMFit03(std::vector< std::vector< double >>, std::vector< std::vector< double >>)
function for backup2
Definition: TrgEclFAMFit.cc:312
Belle2::TrgEclFAMFit::CoeffSigPDF1
std::vector< std::vector< double > > CoeffSigPDF1
Coeffisients of signal PDF1.
Definition: TrgEclFAMFit.h:89
Belle2::TrgEclFAMFit::EventId
int EventId
Fill Analysis table.
Definition: TrgEclFAMFit.h:106
Belle2::TrgEclFAMFit::setup
void setup(int)
setup fam module
Definition: TrgEclFAMFit.cc:75
Belle2::TrgEclFAMFit::BeamBkgTag
std::vector< std::vector< int > > BeamBkgTag
fit timing
Definition: TrgEclFAMFit.h:74
Belle2::TrgEclFAMFit::_DataBase
TrgEclDataBase * _DataBase
Object of DataBase.
Definition: TrgEclFAMFit.h:84
Belle2::TrgEclFAMFit::~TrgEclFAMFit
virtual ~TrgEclFAMFit()
Destructor.
Definition: TrgEclFAMFit.cc:64
Belle2::TrgEclMapping::getTCPhiIdFromTCId
int getTCPhiIdFromTCId(int)
get [TC Phi ID] from [TC ID]
Definition: TrgEclMapping.cc:227
Belle2::TrgEclFAMFit::TCFitTiming
std::vector< std::vector< double > > TCFitTiming
fit timing
Definition: TrgEclFAMFit.h:72
Belle2::TRGECLDigi0::getEventId
int getEventId() const
Get event id.
Definition: TRGECLDigi0.h:60
Belle2::TrgEclFAMFit::_TCMap
TrgEclMapping * _TCMap
Object of TC Mapping.
Definition: TrgEclFAMFit.h:82
Belle2::TrgEclFAMFit::CoeffNoise33
std::vector< std::vector< double > > CoeffNoise33
Coeffisient of noise 3
Definition: TrgEclFAMFit.h:95
Belle2::TrgEclDataBase
class TrgEclDataBase;
Definition: TrgEclDataBase.h:24
Belle2::TrgEclMapping
A class of TC Mapping.
Definition: TrgEclMapping.h:31
Belle2::TrgEclFAMFit::Threshold
std::vector< int > Threshold
Threshold (MeV)
Definition: TrgEclFAMFit.h:104
Belle2::TrgEclFAMFit::TCLatency
std::vector< double > TCLatency
TC Latency.
Definition: TrgEclFAMFit.h:97
Belle2::TrgEclFAMFit::_AnaTag
int _AnaTag
Fill Analysis table.
Definition: TrgEclFAMFit.h:102
Belle2::TrgEclFAMFit::FAMFit01
void FAMFit01(std::vector< std::vector< double >>, std::vector< std::vector< double >>)
function for fitting
Definition: TrgEclFAMFit.cc:90
Belle2::TrgEclFAMFit::CoeffNoise31
std::vector< std::vector< double > > CoeffNoise31
Coeffisients of noise 1.
Definition: TrgEclFAMFit.h:91
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrgEclFAMFit::FAMFit02
void FAMFit02(std::vector< std::vector< double >>, std::vector< std::vector< double >>)
function for backup 1
Definition: TrgEclFAMFit.cc:198
Belle2::TRGECLDigi0
Raw TC result nefor digitizing.
Definition: TRGECLDigi0.h:26
Belle2::TrgEclFAMFit::CoeffSigPDF0
std::vector< std::vector< double > > CoeffSigPDF0
Coeffisients of signal PDF0
Definition: TrgEclFAMFit.h:87
Belle2::TrgEclFAMFit::CoeffNoise32
std::vector< std::vector< double > > CoeffNoise32
Coeffisient of noise 2.
Definition: TrgEclFAMFit.h:93
Belle2::TrgEclFAMFit::SetBeamBkgTag
void SetBeamBkgTag()
Set Beam Background Tag.
Definition: TrgEclFAMFit.cc:383
Belle2::TrgEclFAMFit::_BeamBkgTag
int _BeamBkgTag
Add beambkg.
Definition: TrgEclFAMFit.h:100
Belle2::TrgEclFAMFit::save
void save(int)
save fitting result into tables
Definition: TrgEclFAMFit.cc:438
Belle2::TrgEclFAMFit::TCRawEnergy
std::vector< std::vector< double > > TCRawEnergy
Raw energy.
Definition: TrgEclFAMFit.h:76
Belle2::TrgEclMapping::getTCThetaIdFromTCId
int getTCThetaIdFromTCId(int)
get [TC Theta ID] from [TC ID]
Definition: TrgEclMapping.cc:200
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: ECLMatchingPerformanceExpertModule.h:33
Belle2::TrgEclFAMFit::TCRawTiming
std::vector< std::vector< double > > TCRawTiming
Raw timing.
Definition: TrgEclFAMFit.h:78
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226
Belle2::TrgEclFAMFit::TCFitEnergy
std::vector< std::vector< double > > TCFitEnergy
fit energy
Definition: TrgEclFAMFit.h:70