Belle II Software  release-05-01-25
TrgEclTimingCalibration.cc
1 //---------------------------------------------------------------
2 // $Id$
3 //---------------------------------------------------------------
4 // Filename : TRGECLTimingCalibration.cc
5 // Section : TRG ECL
6 // Owner : YoungJun Kim
7 // Email : rladudwns118@korea.ac.kr
8 //---------------------------------------------------------------
9 // Description : TRGECL TC time offset calibration Module for TRG ECL
10 //
11 //---------------------------------------------------------------
12 // 1.00 : 2019/06/26 : First version
13 //---------------------------------------------------------------
14 
15 #include "trg/ecl/modules/trgeclTimingCal/TrgEclTimingCalibration.h"
16 
17 #include "TDecompLU.h"
18 #include "TFile.h"
19 #include "TTree.h"
20 
21 using namespace Belle2;
22 using namespace std;
23 
24 REG_MODULE(TRGECLTimingCal); // Register module
25 
27 {
28 
29  setDescription("example module for TRGECL timing calibration");
31 
32  addParam("TRGECLCalSim", fSimulation,
33  "0 : Raw data 1 : Simulation data", 0);
34 
35  addParam("TRGECLCalType", fCalType,
36  "0 : Beam data 1 : Cosmic ray data", 0);
37 
38  addParam("TRGECLCal3DBhabhaVeto", f3Dbhabha_veto,
39  "0 : Not using 3D bhabha veto bit 1 : Using 3D bhabha veto bit ", 1);
40 
41  addParam("TRGECLCalTCRef", TC_ref,
42  "Reference TCID. Default setting : 184", 184);
43 
44  addParam("TCEnergyCalibrationConstant", TCEnergyCalibrationConstant,
45  "A TC energy calibration factor (GeV/ADC), default 0.00525", 0.00525);
46 
47  addParam("TRGECLCalnTC", cut_ntc,
48  "Number of TC cut : Events will be used < ntc", 999);
49 
50  addParam("TRGECLCalHighEnergyCut", cut_high_energy,
51  "High energy cut", 9999.0);
52 
53  addParam("TRGECLCalLowEnergyCut", cut_low_energy,
54  "Low energy cut", 0.0);
55 
56  addParam("TRGECLCalFWD", fInclude_FWD,
57  "Flag including FWD end cap", 1);
58 
59  addParam("TRGECLCalBR", fInclude_BR,
60  "Flag including Barrel", 1);
61 
62  addParam("TRGECLCalBWD", fInclude_BWD,
63  "Flag including BWD end cap", 1);
64 
65  addParam("TRGECLCalofname", str_ofilename,
66  "Output file name", str_default_name);
67 
68  addParam("TRGECLCalIteration", fIteration,
69  "iteration", 0);
70 
71  addParam("TRGECLCalChisqCut", cut_chisq,
72  "iteration chisq cut", 1000.0);
73 
74  addParam("TRGECLCalOffsetFname", str_timeoffset_fname,
75  "Input offset file name", str_timeoffset_fname);
76 
77 }
78 
80 {
81 }
82 
84 {
86  StoreArray<TRGECLUnpackerStore> TCHitUnpacker;
89  //StoreArray<TRGECLCluster> TCCluster;
90 
91  if (fSimulation == 1) {
92  TCHit.isRequired();
93  if (f3Dbhabha_veto == 1) {
94  TRGECLTrg.isRequired();
95  }
96  //StoreArray<TRGECLHit>::registerPersistent();
97  } else if (fSimulation == 0) {
98  TCHitUnpacker.isRequired();
99  if (f3Dbhabha_veto == 1) {
100  TRGECLEvtStore.isRequired();
101  }
102  //StoreArray<TRGECLUnpackerStore>::registerPersistent();
103  }
104 
105  InitParams();
106 }
107 
109 {
110 }
111 
113 {
114 }
115 
117 {
118 
120  // -1 : No event for reference TC
121  // 0 : Singular matrix
122  // 1 : Solvable matrix
123  Calculate_Chisq();
124  Save_Result();
125 
126 }
127 
129 {
130  nevt_read++;
131 
132  if (f3Dbhabha_veto == 1) {
133  if (fSimulation == 0) {
134  StoreArray<TRGECLUnpackerEvtStore> TRGECLEvtStore;
135  int nEvtStore = TRGECLEvtStore.getEntries();
136  if (nEvtStore <= 0) return;
137  if (TRGECLEvtStore[0] -> getCheckSum() != 0) return;
138  if (TRGECLEvtStore[0] -> get3DBhabhaV() != 1) return;
139  }
140 
141  if (fSimulation == 1) {
143  int nEvtStore = TRGECLTrg.getEntries();
144  if (nEvtStore <= 0) return;
145  if (TRGECLTrg[0] -> get3DBhabha() != 1) return;
146  }
147  }
149 
150  if (ReadData() == 0) return;
151  if (FillMatrix() == 1)
152  nevt_selected++;
153 
154  return;
155 }
156 
157 // 0: Skip this event 1: Use this event
159 {
160 
161  int ntc;
162  if (fSimulation == 1) { // Simulation data : TRGECLHit table
163  StoreArray<TRGECLHit> TCHit;
164  ntc = TCHit.getEntries();
165  if (ntc < 2 || ntc > cut_ntc) return 0; // ntc cut
166  for (int itc = 0; itc < ntc; itc++) {
167  TCId.push_back(TCHit[itc]->getTCId());
168  TCEnergy.push_back(TCHit[itc]->getEnergyDep());
169  TCTiming.push_back((double) TCHit[itc]->getTimeAve());
170  }
171  } else { // Raw data : TRGECLUnpackerStore table
173  ntc = TCHit.getEntries();
174  if (ntc < 2 || ntc > cut_ntc) return 0; // ntc cut
175  for (int itc = 0; itc < ntc; itc++) {
176  if (TCHit[itc] -> getChecksum() != 0) return 0;
177  if (!(TCHit[itc] -> getHitWin() == 3 || TCHit[itc] -> getHitWin() == 4)) continue;
178  TCId.push_back(TCHit[itc]->getTCId());
179  TCEnergy.push_back((double) TCHit[itc]->getTCEnergy() * TCEnergyCalibrationConstant);
180  TCTiming.push_back((double) TCHit[itc]->getTCTime());
181  }
182  }
183  return 1;
184 }
185 
186 // 0: Skip this event 1: Use this event
188 {
189  int flag_event_used = 0;
190  int ntc = TCId.size();
191 
192  for (int itc1 = 0; itc1 < ntc; itc1++) {
193  if (TCEnergy[itc1] > cut_high_energy || TCEnergy[itc1] < cut_low_energy) continue;
194 
195  int m_TCID1 = TCId[itc1] - 1;
196  if (b_Exclude_TC[m_TCID1] == true) continue;
197 
198  for (int itc2 = itc1 + 1; itc2 < ntc; itc2++) {
199  int m_TCID2 = TCId[itc2] - 1;
200  if (b_Exclude_TC[m_TCID2] == true) continue;
201  // energy cut
202  if (TCEnergy[itc2] > cut_high_energy || TCEnergy[itc2] < cut_low_energy) continue;
203  //double sigma2 = 1.0/(TCEnergy[itc1]*TCEnergy[itc1]) + 1.0/(TCEnergy[itc2]*TCEnergy[itc2]); // sigma i,j or j,i
204 
205  double sigma2 = 5.0 * 5.0; // temporal resolution
206  double tmp_chisq = pow((TCTiming[itc1] - TimeOffset[m_TCID1]) - (TCTiming[itc2] - TimeOffset[m_TCID2]), 2) / sigma2;
207  if (fIteration == 1 && tmp_chisq > cut_chisq) continue;
208 
209  // Fill symmetric matrix
210  _Matrix[m_TCID1][m_TCID2] += (-1.0 / sigma2);
211  _Matrix[m_TCID2][m_TCID1] += (-1.0 / sigma2);
212  // Fill diagonal component
213  _Matrix[m_TCID1][m_TCID1] += (1.0 / sigma2);
214  _Matrix[m_TCID2][m_TCID2] += (1.0 / sigma2);
215 
216  // Fill Chisq matrix1
217  Chisq1[m_TCID1][m_TCID2] += (1.0 / sigma2);
218  Chisq1[m_TCID2][m_TCID1] += (1.0 / sigma2);
219 
220 
221  // TOF correction & Fill vector component
222  double Vector_component = 0;
223  if (fCalType == 0) { //no TOF correction for beam calibration
224  Vector_component = (TCTiming[itc1] - TCTiming[itc2]);
225  }
226 
227  else if (fCalType == 1) { //cosmic ray data
228  double TOF_TC2TC = (TCPosition[m_TCID1] - TCPosition[m_TCID2]).Mag() / 29.9792458;
229  if (TCPosition[m_TCID1].y() > TCPosition[m_TCID2].y()) {
230  Vector_component = (TCTiming[itc1] - (TCTiming[itc2] - TOF_TC2TC));
231  } else Vector_component = ((TCTiming[itc1] - TOF_TC2TC) - TCTiming[itc2]);
232  }
233  Vector_component /= sigma2;
234  _Vector[m_TCID1] += Vector_component;
235  _Vector[m_TCID2] += -Vector_component;
236  // Fill Chisq matrix2,3
237  Chisq2[m_TCID1][m_TCID2] += -2.0 * Vector_component;
238  Chisq2[m_TCID2][m_TCID1] += 2.0 * Vector_component;
239 
240  Chisq3[m_TCID1][m_TCID2] += Vector_component * Vector_component * sigma2;
241  Chisq3[m_TCID2][m_TCID1] += Vector_component * Vector_component * sigma2;
242 
243  Nevent_TC[m_TCID1]++;
244  Nevent_TC[m_TCID2]++;
245  flag_event_used = 1;
246  if (b_Inc_TC[m_TCID2] == false) {b_Inc_TC[m_TCID2] = true;}
247  }// for itc2
248  if (b_Inc_TC[m_TCID1] == false) {b_Inc_TC[m_TCID1] = true;}
249  }// for itc1
250  return flag_event_used;
251 }
252 
254 {
255  tcal_result.clear();
256  tcal_result_err.clear();
257  tcal_result.resize(576);
258  tcal_result_err.resize(576);
259  const int nTC_tot = 576;
260  TMatrixDSym mat_A(575);
261  mat_A.ResizeTo(575, 575, 0.0);
262  TVectorD vec_b(575);
263  vec_b.ResizeTo(575);
264 
265  // copy matrix excluding reference TC
266  for (int i = 0; i < 576; i++) {
267  if (i == TC_ref - 1) continue;
268  int indexi = i;
269  if (i > TC_ref - 1) indexi = i - 1;
270  vec_b[indexi] = _Vector[i];
271  for (int j = 0; j < 576; j++) {
272  if (j == TC_ref - 1) continue;
273  int indexj = j;
274  if (j > TC_ref - 1) indexj = j - 1;
275  mat_A[indexi][indexj] = _Matrix[i][j];
276  }
277  }
278 
279  // error calculation
280  //**** Not supported yet ***********//
281 
282  // (A^T * A)^-1_ii : error^2 of t_i
283  //TMatrixDSym mat_offset_err(mat_A);
284  //TMatrixDSym mat_offset_errT = mat_offset_err.Transpose();
285  //mat_offset_err.ResizeTo(576,576,0.0);
286  //mat_offset_errT.ResizeTo(576,576,0.0);
287 
288  //mat_offset_errT.TMult(mat_offset_err);
289  //mat_offset_errT.Invert();
290  //double *element_mat_offset_err = mat_offset_err.GetMatrixArray();
291  for (int i = 0; i < 576; i++) {
292  //tcal_result_err[i] = mat_offset_errT[i][i];
293  tcal_result_err[i] = 0.0;
294  }
295 
296  TDecompLU lu(mat_A);
297  bool b_solve; // Is this matrix solvable? 0 : no, singular 1 : yes, solvable
298 
299  // Solve Ax = b
300  TVectorD x = lu.Solve(vec_b, b_solve);
301  /* cppcheck-suppress variableScope */
302  int ix = 0;
303  if (b_solve) {
304  for (int i = 0; i < nTC_tot; i++) {
305  if (i != TC_ref - 1) {
306  tcal_result[i] = x[ix];
307  ix++;
308  } else if (i == TC_ref - 1) {
309  tcal_result[i] = 0.0;
310  }
311  }
312  return 1;
313  } else {
314  for (int i = 0; i < nTC_tot; i++) {
315  tcal_result[i] = -999.0;
316  tcal_result_err[i] = -999.0;
317  }
318  return 0;
319  }
320 }
321 
323 {
324  for (int iTCID = 0; iTCID < 576; iTCID++) {
325  TrgEclMapping* trgeclMap = new TrgEclMapping();
326  TCPosition.push_back(trgeclMap->getTCPosition(iTCID + 1));
327  delete trgeclMap;
328  }
329 }
330 
332 {
333  nevt_selected = 0;
334  nevt_read = 0;
335 
336  if (fCalType == 1) {
337  Set_TCposition();
338  }
339 
340  // initialize matrix and vector
341  for (int i = 0; i < 576; i++) {
342  b_Inc_TC.push_back(false);
343  }
344  TimeOffset.assign(576, 0.0);
345  _Matrix = TMatrixDSym(576);
346  _Vector = TVectorD(576);
347  _Matrix.ResizeTo(576, 576, 0.0);
348  _Vector.ResizeTo(576);
349  for (int i = 0; i < 576; i++) {
350  _Vector[i] = 0.0;
351  for (int j = 0; j < 576; j++) {
352  _Matrix[i][j] = 0.0;
353  }
354  }
355  // initialize chisq matrix
356  Chisq1.resize(576);
357  Chisq2.resize(576);
358  Chisq3.resize(576);
359  Nevent_TC.assign(576, 0);
360  for (int i = 0; i < 576; i++) {
361  Chisq1[i].assign(576, 0.0);
362  Chisq2[i].assign(576, 0.0);
363  Chisq3[i].assign(576, 0.0);
364  }
365 
366  // Exclude TC setup
367  b_Exclude_TC.resize(576);
368  for (int i = 0; i < 576; i++) {
369  b_Exclude_TC[i] = false;
370  }
371 
372  if (fInclude_FWD == 0) {
373  for (int i = 0; i < 80; i++) {
374  b_Exclude_TC[i] = true;
375  }
376  }
377  if (fInclude_BR == 0) {
378  for (int i = 80; i < 512; i++) {
379  b_Exclude_TC[i] = true;
380  }
381  }
382  if (fInclude_BWD == 0) {
383  for (int i = 512; i < 576; i++) {
384  b_Exclude_TC[i] = true;
385  }
386  }
387 
388  GetTimeOffset();
389 
390 }
391 
393 {
394  TString fname = str_ofilename;
395  TFile* tf = new TFile(fname, "RECREATE");
396  TTree* tr = new TTree("tree", "TC time calibration result");
397  vector<double> time_offset;
398  vector<double> time_offset_err;
399  vector<int> TCID;
400  vector<double> chisq;
401  vector<int> Nevent;
402  time_offset.resize(576);
403  time_offset_err.resize(576);
404  TCID.resize(576);
405  chisq.resize(576);
406  Nevent.resize(576);
407  double timeoffsetave = 0.0;
408  tr -> Branch("TCID", &TCID);
409  tr -> Branch("TimeOffset", &time_offset);
410  tr -> Branch("TimeOffAve", &timeoffsetave);
411  tr -> Branch("TimeOffsetErr", &time_offset_err);
412  tr -> Branch("Chisq", &chisq);
413  tr -> Branch("Nevent", &Nevent);
414  tr -> Branch("FlagSolved", &FlagMatrixSolved);
415  tr -> Branch("Matrix", &_Matrix);
416  tr -> Branch("Vector", &_Vector);
417  tr -> Branch("ChisqComponent1", &Chisq1);
418  tr -> Branch("ChisqComponent2", &Chisq2);
419  tr -> Branch("ChisqComponent3", &Chisq3);
420  tr -> Branch("NeventRead", &nevt_read);
421  tr -> Branch("NeventUsed", &nevt_selected);
422  for (int i = 0; i < 576; i++) {
423  TCID[i] = i + 1;
424  time_offset[i] = tcal_result[i];
425  timeoffsetave += tcal_result[i];
426  time_offset_err[i] = tcal_result_err[i];
427  chisq[i] = chisq_result[i];
428  Nevent[i] = Nevent_TC[i];
429  }
430  timeoffsetave /= 576.0;
431 
432  tr -> Fill();
433  tf -> cd();
434  tr -> Write();
435  tf -> Close();
436 
437 }
438 
440 {
441  TCId.clear();
442  TCEnergy.clear();
443  TCTiming.clear();
444 
445 }
446 
448 {
449  chisq_result.assign(576, 0.0);
450  if (FlagMatrixSolved == 1) {
451  for (int i = 0; i < 576; i++) {
452  if (b_Inc_TC[i] == false && i != TC_ref - 1) continue;
453  std::vector<double> offsetdiff;
454  offsetdiff.assign(576, 0.0);
455  for (int j = 0; j < 576; j++) {
456  if (i == j || (b_Inc_TC[j] == false && j != TC_ref - 1)) continue;
457  offsetdiff[j] = tcal_result[i] - tcal_result[j];
458  chisq_result[i] += Chisq1[i][j] * offsetdiff[j] * offsetdiff[j] + Chisq2[i][j] * offsetdiff[j] + Chisq3[i][j];
459  }
460  }
461  } else {
462  for (int i = 0; i < 576; i++) {
463  chisq_result[i] = -1.0;
464  }
465  }
466 }
467 
469 {
470  if (fIteration != 0) {
471  TString str_fname_tmp = str_timeoffset_fname;
472  TFile* tf = new TFile(str_fname_tmp, "READ");
473  auto tr = (TTree*) tf -> Get("tree");
474  vector<double>* TimeOffset_tmp = 0;
475  tr -> SetBranchAddress("TimeOffset", &TimeOffset_tmp);
476  TimeOffset = *TimeOffset_tmp;
477  tr -> GetEntry(0);
478  tf -> Close();
479  }
480 }
Belle2::TRGECLTimingCalModule::Chisq2
std::vector< std::vector< double > > Chisq2
chisq component2
Definition: TrgEclTimingCalibration.h:137
Belle2::TRGECLTimingCalModule::fInclude_BR
int fInclude_BR
Flag include barrel 0 : exclude BR 1 : include BR.
Definition: TrgEclTimingCalibration.h:90
Belle2::TRGECLTimingCalModule::fCalType
int fCalType
Calibration type 0 : beam 1 : cosmic.
Definition: TrgEclTimingCalibration.h:84
Belle2::TRGECLTimingCalModule::beginRun
virtual void beginRun() override
Begin Run function.
Definition: TrgEclTimingCalibration.cc:108
Belle2::TRGECLTimingCalModule::chisq_result
std::vector< double > chisq_result
Chisq values based on time offset obtained from "Solve_Matrix" function.
Definition: TrgEclTimingCalibration.h:144
Belle2::TRGECLTimingCalModule::initialize
virtual void initialize() override
initialize function
Definition: TrgEclTimingCalibration.cc:83
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
Belle2::TRGECLTimingCalModule::cut_chisq
double cut_chisq
TC chisq cut (iteration mode)
Definition: TrgEclTimingCalibration.h:112
Belle2::TRGECLTimingCalModule::Solve_Matrix
int Solve_Matrix()
Solve matrix.
Definition: TrgEclTimingCalibration.cc:253
Belle2::TRGECLTimingCalModule::cut_ntc
int cut_ntc
The number of TC cut.
Definition: TrgEclTimingCalibration.h:101
Belle2::TRGECLTimingCalModule::Save_Result
void Save_Result()
Save time offset and chisq as a root file.
Definition: TrgEclTimingCalibration.cc:392
Belle2::TRGECLTimingCalModule::InitParams
void InitParams()
Parameters initialization.
Definition: TrgEclTimingCalibration.cc:331
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::Module::c_ParallelProcessingCertified
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition: Module.h:82
Belle2::TRGECLTimingCalModule::b_Inc_TC
std::vector< bool > b_Inc_TC
Flag TC included.
Definition: TrgEclTimingCalibration.h:148
Belle2::TRGECLTimingCalModule::Nevent_TC
std::vector< int > Nevent_TC
The number of uesed TC hit (TC by TC)
Definition: TrgEclTimingCalibration.h:142
Belle2::TRGECLTimingCalModule::fInclude_FWD
int fInclude_FWD
Flag include forward endcap 0 : exclude FWD 1 : include FWD.
Definition: TrgEclTimingCalibration.h:88
Belle2::TRGECLTimingCalModule::TimeOffset
std::vector< double > TimeOffset
Time offset Input time offset (iteration mode)
Definition: TrgEclTimingCalibration.h:161
Belle2::TRGECLTimingCalModule::tcal_result
std::vector< double > tcal_result
Time offset calibration result.
Definition: TrgEclTimingCalibration.h:152
Belle2::TRGECLTimingCalModule::str_default_name
std::string str_default_name
Default output root file name.
Definition: TrgEclTimingCalibration.h:115
Belle2::TRGECLTimingCalModule::cut_low_energy
double cut_low_energy
TC energy cut lower limit.
Definition: TrgEclTimingCalibration.h:99
Belle2::TRGECLTimingCalModule::fInclude_BWD
int fInclude_BWD
Flag include backward endcap 0 : exclude BWD 1 : include BWD.
Definition: TrgEclTimingCalibration.h:92
Belle2::TRGECLTimingCalModule::TC_ref
int TC_ref
Reference TC whose time offset set to be 0.
Definition: TrgEclTimingCalibration.h:103
Belle2::TrgEclMapping
A class of TC Mapping.
Definition: TrgEclMapping.h:31
Belle2::TRGECLTimingCalModule::cut_high_energy
double cut_high_energy
TC energy cut upper limit.
Definition: TrgEclTimingCalibration.h:97
Belle2::TRGECLTimingCalModule::event
virtual void event() override
Event function.
Definition: TrgEclTimingCalibration.cc:128
Belle2::TRGECLTimingCalModule::f3Dbhabha_veto
int f3Dbhabha_veto
Flag 3D bhabha veto 0 : not use 1 : use.
Definition: TrgEclTimingCalibration.h:86
Belle2::TRGECLTimingCalModule::Chisq1
std::vector< std::vector< double > > Chisq1
chisq component1
Definition: TrgEclTimingCalibration.h:135
Belle2::TRGECLTimingCalModule::_Vector
TVectorD _Vector
Vector component.
Definition: TrgEclTimingCalibration.h:146
Belle2::TRGECLTimingCalModule::TCTiming
std::vector< double > TCTiming
TC time.
Definition: TrgEclTimingCalibration.h:127
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::Module::setPropertyFlags
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition: Module.cc:210
Belle2::TRGECLTimingCalModule::Chisq3
std::vector< std::vector< double > > Chisq3
chisq component3
Definition: TrgEclTimingCalibration.h:139
Belle2::TRGECLTimingCalModule::TRGECLTimingCalClear
void TRGECLTimingCalClear()
Clear vector every event.
Definition: TrgEclTimingCalibration.cc:439
Belle2::TRGECLTimingCalModule::fSimulation
int fSimulation
Flag simulation 0 : real data 1 : simulation data.
Definition: TrgEclTimingCalibration.h:82
Belle2::TRGECLTimingCalModule::GetTimeOffset
void GetTimeOffset()
Get previous time offsets (iteration mode)
Definition: TrgEclTimingCalibration.cc:468
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TRGECLTimingCalModule::TCPosition
std::vector< TVector3 > TCPosition
TC position from TRGECLMap (cosmic calibration)
Definition: TrgEclTimingCalibration.h:129
Belle2::TRGECLTimingCalModule::TCEnergyCalibrationConstant
double TCEnergyCalibrationConstant
ADC to GeV energy conversion factor.
Definition: TrgEclTimingCalibration.h:110
Belle2::TrgEclMapping::getTCPosition
TVector3 getTCPosition(int)
TC position (cm)
Definition: TrgEclMapping.cc:251
Belle2::TRGECLTimingCalModule::FlagMatrixSolved
int FlagMatrixSolved
Flag matrix solved (0 : not solved 1 : solved)
Definition: TrgEclTimingCalibration.h:157
Belle2::TRGECLTimingCalModule::endRun
virtual void endRun() override
End Run function.
Definition: TrgEclTimingCalibration.cc:112
Belle2::TRGECLTimingCalModule::fIteration
int fIteration
Flag iteration 0 : first calbiration 1 : iteration.
Definition: TrgEclTimingCalibration.h:94
Belle2::TRGECLTimingCalModule::_Matrix
TMatrixDSym _Matrix
Matrix Marix component.
Definition: TrgEclTimingCalibration.h:133
Belle2::TRGECLTimingCalModule::TRGECLTimingCalModule
TRGECLTimingCalModule()
Constructor.
Definition: TrgEclTimingCalibration.cc:26
Belle2::TRGECLTimingCalModule::str_timeoffset_fname
std::string str_timeoffset_fname
Input time offset file name (iteration mode)
Definition: TrgEclTimingCalibration.h:119
Belle2::TRGECLTimingCalModule::terminate
virtual void terminate() override
Terminate function.
Definition: TrgEclTimingCalibration.cc:116
Belle2::Module::addParam
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition: Module.h:562
Belle2::TRGECLTimingCalModule::Set_TCposition
void Set_TCposition()
Set TC position from TRGECLMap.
Definition: TrgEclTimingCalibration.cc:322
Belle2::TRGECLTimingCalModule::nevt_selected
int nevt_selected
The number of selected events.
Definition: TrgEclTimingCalibration.h:106
Belle2::TRGECLTrg
Example Detector.
Definition: TRGECLTrg.h:27
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: ECLMatchingPerformanceExpertModule.h:33
Belle2::TRGECLTimingCalModule::FillMatrix
int FillMatrix()
Fill matrix and vector components.
Definition: TrgEclTimingCalibration.cc:187
Belle2::TRGECLTimingCalModule::TCId
std::vector< int > TCId
Trigger Cell TCID.
Definition: TrgEclTimingCalibration.h:123
Belle2::TRGECLTimingCalModule::str_ofilename
std::string str_ofilename
Output root file name.
Definition: TrgEclTimingCalibration.h:117
Belle2::TRGECLTimingCalModule::Calculate_Chisq
void Calculate_Chisq()
Calculate chisq based on time offset obtained from "Solve_Matrix" function.
Definition: TrgEclTimingCalibration.cc:447
Belle2::TRGECLTimingCalModule::b_Exclude_TC
std::vector< bool > b_Exclude_TC
Flag TC excluded.
Definition: TrgEclTimingCalibration.h:150
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226
Belle2::TRGECLTimingCalModule::nevt_read
int nevt_read
The number of read events.
Definition: TrgEclTimingCalibration.h:108
Belle2::TRGECLTimingCalModule::ReadData
int ReadData()
Fill data from.
Definition: TrgEclTimingCalibration.cc:158
Belle2::TRGECLTimingCalModule::tcal_result_err
std::vector< double > tcal_result_err
Time offset error.
Definition: TrgEclTimingCalibration.h:154
Belle2::TRGECLTimingCalModule::TCEnergy
std::vector< double > TCEnergy
TC energy.
Definition: TrgEclTimingCalibration.h:125
Belle2::TRGECLTimingCalModule::~TRGECLTimingCalModule
virtual ~TRGECLTimingCalModule()
Destructor.
Definition: TrgEclTimingCalibration.cc:79