Belle II Software  release-05-01-25
EventTime.cc
1 //-----------------------------------------------------------------------------
2 // $Id$
3 //-----------------------------------------------------------------------------
4 // Filename : EventTime.cc
5 // Section : TRG CDC
6 // Owner : KyungTae KIM (K.U.)
7 // Email : ktkim@hep.korea.ac.kr
8 //-----------------------------------------------------------------------------
9 // Description : A class to get Event Time information
10 //-----------------------------------------------------------------------------
11 // $Log$
12 //-----------------------------------------------------------------------------
13 
14 #define TRG_SHORT_NAMES
15 #define TRGCDC_SHORT_NAMES
16 
17 #include "trg/trg/Debug.h"
18 #include "trg/cdc/Segment.h"
19 #include "trg/cdc/Track.h"
20 #include "trg/cdc/EventTime.h"
21 #include <iostream>
22 #include <iomanip>
23 #include <fstream>
24 #include <bitset>
25 #include <vector>
26 
27 #include "trg/cdc/TRGCDC.h"
28 #include "trg/cdc/Layer.h"
29 #include "trg/cdc/Wire.h"
30 #include "trg/cdc/WireHit.h"
31 #include "trg/cdc/SegmentHit.h"
32 #include "trg/trg/Utilities.h"
33 
34 #include "TH1D.h"
35 
36 using namespace std;
37 
38 namespace Belle2 {
44  TRGCDCEventTime::TRGCDCEventTime(const TRGCDC& TRGCDC, bool makeRootFile)
45  : _cdc(TRGCDC),
46  cnt{{0}},
47  ft{{{0}}},
48  m_foundT0(0), // 2019/07/31 by ytlai
49  m_minusET(0),
50  m_noET(0),
51  threshold(0),
52  m_eventN(0),
53  m_makeRootFile(makeRootFile),
54  m_ver(0)
55  {
56 
57  if (m_makeRootFile) m_fileEvtTime = new TFile("ETF.root", "RECREATE");
58 
59  m_evtOutputTs = new TTree("m_evtOutsTs", "EvtTime Outputs for TS");
60  m_evtOut = new TTree("m_evtOut", "EvtTime Outputs for 20 TS");
61 
62  m_evtOutputTs->Branch("tsFastestTime", &m_fastestT, "m_fastestT/I");
63  m_evtOutputTs->Branch("tsFoundTime", &m_foundT, "m_foundT/I");
64  m_evtOutputTs->Branch("foundT_fastT_diff", &m_whdiff, "m_whdiff/I");
65 
66  m_evtOut->Branch("histTime", &m_histT, "m_histT/I");
67 
68  h = new TH1D("h", "h", 1000, -500, 499);
69 
70  memset(cnt, 0, sizeof cnt);
71  memset(ft, 0, sizeof ft);
72 
73  // move to constructor
74  //m_minusET = 0;
75  //m_noET = 0;
76  //m_eventN = 0;
77  //m_ver = 0;
78  //m_foundT0 = 0;
79  }
81  {
82  delete m_evtOutputTs;
83  delete m_evtOut;
84  if (m_makeRootFile) delete m_fileEvtTime;
85  delete h;
86  }
87 
89  {
90  m_minusET = 0;
91  m_noET = 0;
92  m_eventN = 0;
93  }
94 
95  void TRGCDCEventTime::doit(int ver, bool print)
96  {
97  TRGDebug::enterStage("Event Time");
98  m_eventN += 1;
99  m_ver = ver;
100  if (ver == 1) {
101  oldVer();
102  getT0();
103  } else if (ver == 0) {
104  hitcount();
105  hist();
106  getT0();
107  } else if (ver == 2) {
108  m_histT = 0;
109  m_foundT0 = 1;
110  getT0();
111  } else {
112  cout << "verETF error" << endl;
113  cout << "choose 0:hist, 1:old, 2:trueT0" << endl;
114  }
115  if (print) printFirm();
116  TRGDebug::leaveStage("Event Time");
117  }
118 
120  {
121  TRGDebug::enterStage("Event Time");
122 
123  memset(cnt, 0, sizeof cnt);
124  memset(ft, 0, sizeof ft);
125 
126  const vector<const TCSHit*> tsh = _cdc.segmentHits();
127 
128  for (int iClk = 0; iClk < 64; iClk++) {
129  for (int iTS = 0; iTS < (int)tsh.size(); iTS++) {
130  const TRGCDCSegment& ts = tsh[iTS]->segment();
131  int fndC = ts.foundTime() / 16;
132  m_fastestT = ts.fastestTime();
133  m_foundT = ts.foundTime();
134  m_whdiff = ts.foundTime() - ts.fastestTime();
135  if (m_whdiff <= 256) {
136  if (iClk == fndC + 31) {
137  for (int iSL = 0; iSL < 9; iSL++) {
138  if ((int)ts.layer().id() == iSL) {
139  cnt[iSL][iClk] += 1;
140  if (cnt[iSL][iClk] <= 10) {
141  h->Fill(m_fastestT);
142  m_evtOutputTs->Fill();
143  ft[iSL][iClk][cnt[iSL][iClk] - 1] = ts.fastestTime();
144  }
145  }
146  }
147  }
148  }
149  }
150  }
151  TRGDebug::leaveStage("Event Time");
152  }//end of TRGCEventTime::hitcount
153 
155  {
156  TRGDebug::enterStage("Event Time");
157 
158  m_foundT0 = 0;
159  threshold = 3;
160  m_histT = 500;
161 
162  for (int i = 450; i < 600; i++) {
163  if (h->GetBinContent(i) > threshold) {
164  m_foundT0 = 1;
165  m_histT = i - 500;
166  m_evtOut->Fill();
167  if (m_histT < 0) {
168  m_minusET += 1;
169  }
170  break;
171  }
172  }//end of loop over hist bins
173  if (m_histT == 500) {
174  m_noET += 1;
175  }//if no bin has contents over threshold no evt time
176  h->Reset();
177 
178  TRGDebug::leaveStage("Event Time");
179  }//end of TRGCDCEventTime::hist
180 
182  {
183  TRGDebug::enterStage("Event Time");
184  m_histT = 65535;
185  for (unsigned i = 0; i < _cdc.nSegmentLayers(); i++) {
186  const Belle2::TRGCDCLayer* l = _cdc.segmentLayer(i);
187  const unsigned nWires = l->nCells();
188  for (unsigned j = 0; j < nWires; j++) {
189  const TCSegment& s = (TCSegment&) * (*l)[j];
190  const vector<const TCWire*>& wires = s.wires();
191  const TRGSignal& timing = s.signal();
192  if (timing.active()) {
193  for (unsigned k = 0; k < wires.size(); k++) {
194  if (wires[k]->hit()) {
195  int dt = wires[k]->signal()[0]->time();
196  if (m_histT > dt) {
197  m_histT = dt;
198  m_foundT0 = 1;
199  m_evtOut->Fill();
200  }
201  }
202  }
203  }
204  }
205  }
206  TRGDebug::leaveStage("Event Time");
207  }//end of TRTGCDCEventTime::oldVer
208 
210  {
211  TRGDebug::enterStage("Event Time");
212 
213  for (int iClk = 0; iClk < 64; iClk++) {
214  for (int iSL = 0; iSL < 9; iSL++) {
215  int hmbw = 160 + 32 * (iSL - 1); //set bit width of hitmap part
216  if (iSL == 0) hmbw = 160;
217  string hs[10];
218  string t[10];
219  char file[30];
220  sprintf(file, "SL%d.coe", iSL);
221  string cnts = string(cnt[iSL][iClk], '1');
222 
223  string clk = bitset<9>(iClk + 64 * (m_eventN - 1)).to_string();
224  string clk2 = bitset<5>(iClk).to_string();
225  for (int i = 0; i < 10; i++) {
226  hs[i] = bitset<4>(ft[iSL][iClk][i]).to_string();
227  if (cnt[iSL][iClk] > i) t[i] = clk2 + hs[i];
228  else t[i] = "00000" + hs[i];
229  }
230  ofstream fout;
231  fout.open(file, ios::app);
232  fout << clk;
233  fout << setfill('0') << setw(hmbw) << cnts;
234  fout << t[9] << t[8] << t[7] << t[6] << t[5] << t[4] << t[3] << t[2] << t[1] << t[0] << endl;
235  fout.close();
236  }
237  }
238  TRGDebug::leaveStage("Event Time");
239  }//end of TRGCDCEventTime::printFirm
240 
241 
242  int TRGCDCEventTime::getT0(void)const
243  {
244  int et = m_histT;
245  if (m_foundT0 == 0) et = 9999;
246  return et;
247  }
248 
250  {
251  delete h;
252  if (m_makeRootFile) {
253  m_fileEvtTime->Write();
254  m_fileEvtTime->Close();
255  }
256  }
258 }
Belle2::TRGCDCEventTime::oldVer
void oldVer(void)
old version of calculation function
Definition: EventTime.cc:181
Belle2::TRGSignal
A class to represent a digitized signal. Unit is nano second.
Definition: Signal.h:28
Belle2::TRGCDCEventTime::m_histT
int m_histT
calculated T0
Definition: EventTime.h:75
Belle2::TRGCDCEventTime::m_foundT
int m_foundT
Found time of TS.
Definition: EventTime.h:83
Belle2::TRGCDC::segmentHits
std::vector< const TRGCDCSegmentHit * > segmentHits(void) const
returns a list of TRGCDCSegmentHit.
Definition: TRGCDC.h:997
Belle2::TRGCDCEventTime::~TRGCDCEventTime
virtual ~TRGCDCEventTime()
destructor of TRGCDCEventTime class
Definition: EventTime.cc:80
Belle2::TRGCDCEventTime::m_whdiff
int m_whdiff
Drift time of TS.
Definition: EventTime.h:85
Belle2::TRGCDCEventTime::m_fileEvtTime
TFile * m_fileEvtTime
TFile pointer.
Definition: EventTime.h:64
Belle2::TRGCDC::nSegmentLayers
unsigned nSegmentLayers(void) const
returns # of track segment layers.
Definition: TRGCDC.h:1071
Belle2::TRGCDCEventTime::hist
void hist(void)
making hostogram
Definition: EventTime.cc:154
Belle2::TRGCDCEventTime::h
TH1 * h
TH1 pointer of the TFile.
Definition: EventTime.h:70
Belle2::TRGCDCLayer
A class to represent a cell layer.
Definition: Layer.h:34
Belle2::TRGCDCEventTime::m_noET
int m_noET
no ET bin is looped
Definition: EventTime.h:89
Belle2::TRGCDCEventTime::m_eventN
int m_eventN
event number
Definition: EventTime.h:96
Belle2::TRGCDCEventTime::m_ver
int m_ver
Version.
Definition: EventTime.h:101
Belle2::TRGCDCEventTime::cnt
int cnt[9][64]
TS counter for each SL and clock.
Definition: EventTime.h:77
Belle2::TRGCDC::segmentLayer
const TRGCDCLayer * segmentLayer(unsigned id) const
returns a pointer to a track segment layer.
Definition: TRGCDC.h:1062
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TRGCDCEventTime::m_makeRootFile
bool m_makeRootFile
make Root file or not
Definition: EventTime.h:99
Belle2::TRGCDCEventTime::threshold
int threshold
Threshold value.
Definition: EventTime.h:93
Belle2::TRGCDCEventTime::m_foundT0
bool m_foundT0
T0 is found or not.
Definition: EventTime.h:81
Belle2::TRGCDCEventTime::_cdc
const TRGCDC & _cdc
TRGCDC class.
Definition: EventTime.h:61
Belle2::TRGCDC
The instance of TRGCDC is a singleton.
Definition: TRGCDC.h:70
Belle2::TRGCDCSegment
A class to represent a wire in CDC.
Definition: Segment.h:40
Belle2::TRGCDCEventTime::ft
int ft[9][64][10]
Fastest time array each SL and clock.
Definition: EventTime.h:79
Belle2::TRGCDCEventTime::printFirm
void printFirm(void)
Print info in firmware level.
Definition: EventTime.cc:209
Belle2::TRGCDCEventTime::hitcount
void hitcount(void)
hit count of TS
Definition: EventTime.cc:119
Belle2::TRGCDCEventTime::terminate
void terminate(void)
terminate function
Definition: EventTime.cc:249
Belle2::TRGDebug::leaveStage
static void leaveStage(const std::string &stageName)
Declare that you leave a stage.
Definition: Debug.cc:39
Belle2::TRGCDCEventTime::m_minusET
int m_minusET
minus ET bin
Definition: EventTime.h:87
Belle2::TRGDebug::enterStage
static void enterStage(const std::string &stageName)
Declare that you enter new stage.
Definition: Debug.cc:29
Belle2::TRGCDCEventTime::getT0
int getT0(void) const
Calculate T0.
Definition: EventTime.cc:242
Belle2::TRGCDCEventTime::m_fastestT
int m_fastestT
The fastest time of TS.
Definition: EventTime.h:73
Belle2::TRGCDCEventTime::doit
void doit(int ver, bool print)
Calculate T0 based on ver.
Definition: EventTime.cc:95
Belle2::TRGCDCEventTime::m_evtOut
TTree * m_evtOut
TTree pointer of the TFile.
Definition: EventTime.h:68
Belle2::TRGCDCEventTime::initialize
void initialize(void)
initialize the class
Definition: EventTime.cc:88
Belle2::TRGCDCEventTime::m_evtOutputTs
TTree * m_evtOutputTs
TTree pointer of the TFile.
Definition: EventTime.h:66