Belle II Software development
KLMTimeAlgorithm.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9/* Own header. */
10#include <klm/calibration/KLMTimeAlgorithm.h>
11
12/* KLM headers. */
13#include <klm/dataobjects/bklm/BKLMElementNumbers.h>
14
15/* Basf2 headers. */
16#include <framework/database/Database.h>
17#include <framework/database/DBObjPtr.h>
18#include <framework/database/DBStore.h>
19#include <framework/dataobjects/EventMetaData.h>
20#include <framework/gearbox/Const.h>
21#include <framework/logging/Logger.h>
22
23/* ROOT headers. */
24#include <Math/MinimizerOptions.h>
25#include <Math/Vector3D.h>
26#include <TFile.h>
27#include <TFitResult.h>
28#include <TMinuit.h>
29#include <TROOT.h>
30#include <TString.h>
31#include <TTree.h>
32
33using namespace Belle2;
34using namespace ROOT::Math;
35
37const int c_NBinsTime = 100;
38
40const int c_NBinsDistance = 100;
41
43static double s_BinnedData[c_NBinsTime][c_NBinsDistance];
44
46static double s_LowerTimeBoundary = 0;
47
49static double s_UpperTimeBoundary = 0;
50
52static double s_StripLength = 0;
53
54static bool compareEventNumber(const std::pair<KLMChannelNumber, unsigned int>& pair1,
55 const std::pair<KLMChannelNumber, unsigned int>& pair2)
56{
57 return pair1.second < pair2.second;
58}
59
60static double timeDensity(const double x[2], const double* par)
61{
62 double polynomial, t0, gauss;
63 polynomial = par[0];
64 t0 = par[2] + par[4] * x[1];
65 gauss = par[1] / (sqrt(2.0 * M_PI) * par[3]) *
66 exp(-0.5 * pow((x[0] - t0) / par[3], 2));
67 return fabs(polynomial + gauss);
68}
69
70static void fcn(int& npar, double* grad, double& fval, double* par, int iflag)
71{
72 (void)npar;
73 (void)grad;
74 (void)iflag;
75 double x[2];
76 fval = 0;
77 for (int i = 0; i < c_NBinsTime; ++i) {
78 x[0] = s_LowerTimeBoundary +
79 (s_UpperTimeBoundary - s_LowerTimeBoundary) *
80 (double(i) + 0.5) / c_NBinsTime;
81 for (int j = 0; j < c_NBinsDistance; ++j) {
82 x[1] = s_StripLength * (double(j) + 0.5) / c_NBinsDistance;
83 double f = timeDensity(x, par);
84 if (s_BinnedData[i][j] == 0)
85 fval = fval + 2.0 * f;
86 else
87 fval = fval + 2.0 * (f - s_BinnedData[i][j] *
88 (1.0 - log(s_BinnedData[i][j] / f)));
89 }
90 }
91}
92
93KLMTimeAlgorithm::KLMTimeAlgorithm() :
94 CalibrationAlgorithm("KLMTimeCollector")
95{
96 m_ElementNumbers = &(KLMElementNumbers::Instance());
97 m_minimizerOptions = ROOT::Math::MinimizerOptions();
98}
99
100KLMTimeAlgorithm::~KLMTimeAlgorithm()
101{
102}
103
104void KLMTimeAlgorithm::setupDatabase()
105{
106 const std::vector<Calibration::ExpRun>& runs = getRunList();
107 int firstExperiment = runs[0].first;
108 int lastExperiment = runs[runs.size() - 1].first;
109 if (firstExperiment != lastExperiment) {
110 B2FATAL("Runs from different experiments are used "
111 "for KLM time calibration (single algorithm run).");
112 }
113 /* DataStore. */
114 DataStore::Instance().setInitializeActive(true);
115 StoreObjPtr<EventMetaData> eventMetaData;
116 eventMetaData.registerInDataStore();
117 DataStore::Instance().setInitializeActive(false);
118 /* Database. */
119 if (eventMetaData.isValid()) {
120 if (eventMetaData->getExperiment() != firstExperiment) {
121 B2FATAL("Runs from different experiments are used "
122 "for KLM time calibration (consecutive algorithm runs).");
123 }
124 eventMetaData->setExperiment(firstExperiment);
125 eventMetaData->setRun(runs[0].second);
126 } else {
127 eventMetaData.construct(1, runs[0].second, firstExperiment);
128 }
129 DBStore& dbStore = DBStore::Instance();
130 dbStore.update();
131 dbStore.updateEvent();
132 /*
133 * For calibration algorithms, the database is not initialized on class
134 * creation. Do not move the database object to class members.
135 */
136 DBObjPtr<BKLMGeometryPar> bklmGeometry;
137 m_BKLMGeometry = bklm::GeometryPar::instance(*bklmGeometry);
138 m_EKLMGeometry = &(EKLM::GeometryData::Instance());
139}
140
141CalibrationAlgorithm::EResult KLMTimeAlgorithm::readCalibrationData()
142{
143 B2INFO("Read tree entries and separate events by module id.");
144 Event event;
145 std::shared_ptr<TTree> timeCalibrationData;
146 timeCalibrationData = getObjectPtr<TTree>("time_calibration_data");
147 timeCalibrationData->SetBranchAddress("t0", &event.t0);
148 timeCalibrationData->SetBranchAddress("flyTime", &event.flyTime);
149 timeCalibrationData->SetBranchAddress("recTime", &event.recTime);
150 timeCalibrationData->SetBranchAddress("dist", &event.dist);
151 timeCalibrationData->SetBranchAddress("diffDistX", &event.diffDistX);
152 timeCalibrationData->SetBranchAddress("diffDistY", &event.diffDistY);
153 timeCalibrationData->SetBranchAddress("diffDistZ", &event.diffDistZ);
154 timeCalibrationData->SetBranchAddress("eDep", &event.eDep);
155 timeCalibrationData->SetBranchAddress("nPE", &event.nPE);
156 timeCalibrationData->SetBranchAddress("channelId", &event.channelId);
157 timeCalibrationData->SetBranchAddress("inRPC", &event.inRPC);
158 timeCalibrationData->SetBranchAddress("isFlipped", &event.isFlipped);
159
160 B2INFO(LogVar("Total number of digits:", timeCalibrationData->GetEntries()));
161 m_evts.clear();
162
163 int n = timeCalibrationData->GetEntries();
164 if (n < m_MinimalDigitNumber)
165 return CalibrationAlgorithm::c_NotEnoughData;
166 for (int i = 0; i < n; ++i) {
167 timeCalibrationData->GetEntry(i);
168 m_evts[event.channelId].push_back(event);
169 }
170 B2INFO("Events packing finish.");
171 return CalibrationAlgorithm::c_OK;
172}
173
174void KLMTimeAlgorithm::createHistograms()
175{
176 if (m_mc) {
177 m_LowerTimeBoundaryRPC = -10.0;
178 m_UpperTimeBoundaryRPC = 10.0;
179 m_LowerTimeBoundaryScintillatorsBKLM = 20.0;
180 m_UpperTimeBoundaryScintillatorsBKLM = 70.0;
181 m_LowerTimeBoundaryScintillatorsEKLM = 20.0;
182 m_UpperTimeBoundaryScintillatorsEKLM = 70.0;
183 } else {
184 m_LowerTimeBoundaryRPC = -800.0;
185 m_UpperTimeBoundaryRPC = -600.0;
186 m_LowerTimeBoundaryScintillatorsBKLM = -4800.0;
187 m_UpperTimeBoundaryScintillatorsBKLM = -4400.0;
188 m_LowerTimeBoundaryScintillatorsEKLM = -4950.0;
189 m_UpperTimeBoundaryScintillatorsEKLM = -4650.0;
190
191 }
192 int nBin = 200;
193 int nBin_scint = 400;
194
195 TString iFstring[2] = {"Backward", "Forward"};
196 TString iPstring[2] = {"ZReadout", "PhiReadout"};
197 TString hn, ht;
198
199 h_diff = new TH1F("h_diff", "Position difference between bklmHit2d and extHit;position difference", 100, 0, 10);
200 h_calibrated = new TH1I("h_calibrated_summary", "h_calibrated_summary;calibrated or not", 3, 0, 3);
201 hc_calibrated = new TH1I("hc_calibrated_summary", "hc_calibrated_summary;calibrated or not", 3, 0, 3);
202
203 gre_time_channel_scint = new TGraphErrors();
204 gre_time_channel_rpc = new TGraphErrors();
205 gre_time_channel_scint_end = new TGraphErrors();
206
207 gr_timeShift_channel_scint = new TGraph();
208 gr_timeShift_channel_rpc = new TGraph();
209 gr_timeShift_channel_scint_end = new TGraph();
210
211 gre_ctime_channel_scint = new TGraphErrors();
212 gre_ctime_channel_rpc = new TGraphErrors();
213 gre_ctime_channel_scint_end = new TGraphErrors();
214
215 gr_timeRes_channel_scint = new TGraph();
216 gr_timeRes_channel_rpc = new TGraph();
217 gr_timeRes_channel_scint_end = new TGraph();
218
219 double maximalPhiStripLengthBKLM =
220 m_BKLMGeometry->getMaximalPhiStripLength();
221 double maximalZStripLengthBKLM =
222 m_BKLMGeometry->getMaximalZStripLength();
223 double maximalStripLengthEKLM =
224 m_EKLMGeometry->getMaximalStripLength() / CLHEP::cm * Unit::cm;
225
226 m_ProfileRpcPhi = new TProfile("hprf_rpc_phi_effC",
227 "Time over propagation length for RPCs (Phi_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]", 400, 0.0,
228 400.0);
229 m_ProfileRpcZ = new TProfile("hprf_rpc_z_effC",
230 "Time over propagation length for RPCs (Z_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]", 400, 0.0,
231 400.0);
232 m_ProfileBKLMScintillatorPhi = new TProfile("hprf_scint_phi_effC",
233 "Time over propagation length for scintillators (Phi_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
234 200, 0.0, maximalPhiStripLengthBKLM);
235 m_ProfileBKLMScintillatorZ = new TProfile("hprf_scint_z_effC",
236 "Time over propagation length for scintillators (Z_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
237 200, 0.0, maximalZStripLengthBKLM);
238 m_ProfileEKLMScintillatorPlane1 = new TProfile("hprf_scint_plane1_effC_end",
239 "Time over propagation length for scintillators (plane1, Endcap); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
240 200, 0.0, maximalStripLengthEKLM);
241 m_ProfileEKLMScintillatorPlane2 = new TProfile("hprf_scint_plane2_effC_end",
242 "Time over propagation length for scintillators (plane2, Endcap); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
243 200, 0.0, maximalStripLengthEKLM);
244
245 m_Profile2RpcPhi = new TProfile("hprf2_rpc_phi_effC",
246 "Time over propagation length for RPCs (Phi_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]", 400, 0.0,
247 400.0);
248 m_Profile2RpcZ = new TProfile("hprf2_rpc_z_effC",
249 "Time over propagation length for RPCs (Z_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]", 400, 0.0,
250 400.0);
251 m_Profile2BKLMScintillatorPhi = new TProfile("hprf2_scint_phi_effC",
252 "Time over propagation length for scintillators (Phi_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
253 200, 0.0, maximalPhiStripLengthBKLM);
254 m_Profile2BKLMScintillatorZ = new TProfile("hprf2_scint_z_effC",
255 "Time over propagation length for scintillators (Z_Readout); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
256 200, 0.0, maximalZStripLengthBKLM);
257 m_Profile2EKLMScintillatorPlane1 = new TProfile("hprf2_scint_plane1_effC_end",
258 "Time over propagation length for scintillators (plane1, Endcap); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
259 200, 0.0, maximalStripLengthEKLM);
260 m_Profile2EKLMScintillatorPlane2 = new TProfile("hprf2_scint_plane2_effC_end",
261 "Time over propagation length for scintillators (plane2, Endcap); propagation distance[cm]; T_rec-T_0-T_fly-'T_calibration'[ns]",
262 200, 0.0, maximalStripLengthEKLM);
263
264 h_time_rpc_tc = new TH1F("h_time_rpc_tc", "time distribution for RPC", nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
265 h_time_scint_tc = new TH1F("h_time_scint_tc", "time distribution for Scintillator", nBin_scint,
266 m_LowerTimeBoundaryScintillatorsBKLM, m_UpperTimeBoundaryScintillatorsBKLM);
267 h_time_scint_tc_end = new TH1F("h_time_scint_tc_end", "time distribution for Scintillator (Endcap)", nBin_scint,
268 m_LowerTimeBoundaryScintillatorsEKLM,
269 m_UpperTimeBoundaryScintillatorsEKLM);
270
272 h_time_rpc = new TH1F("h_time_rpc", "time distribution for RPC; T_rec-T_0-T_fly-T_propagation[ns]", nBin, m_LowerTimeBoundaryRPC,
273 m_UpperTimeBoundaryRPC);
274 h_time_scint = new TH1F("h_time_scint", "time distribution for Scintillator; T_rec-T_0-T_fly-T_propagation[ns]", nBin_scint,
275 m_LowerTimeBoundaryScintillatorsBKLM, m_UpperTimeBoundaryScintillatorsBKLM);
276 h_time_scint_end = new TH1F("h_time_scint_end", "time distribution for Scintillator (Endcap); T_rec-T_0-T_fly-T_propagation[ns]",
277 nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM, m_UpperTimeBoundaryScintillatorsEKLM);
278
279 hc_time_rpc = new TH1F("hc_time_rpc", "Calibrated time distribution for RPC; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
280 nBin, m_LowerTimeBoundaryCalibratedRPC, m_UpperTimeBoundaryCalibratedRPC);
281 hc_time_scint = new TH1F("hc_time_scint",
282 "Calibrated time distribution for Scintillator; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]", nBin_scint,
283 m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
284 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
285 hc_time_scint_end = new TH1F("hc_time_scint_end",
286 "Calibrated time distribution for Scintillator (Endcap); T_rec-T_0-T_fly-T_propagation-T_calibration[ns]", nBin_scint,
287 m_LowerTimeBoundaryCalibratedScintillatorsEKLM, m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
288
289 for (int iF = 0; iF < 2; ++iF) {
290 hn = Form("h_timeF%d_rpc", iF);
291 ht = Form("Time distribution for RPC of %s; T_rec-T_0-T_fly-T_propagation[ns]", iFstring[iF].Data());
292 h_timeF_rpc[iF] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
293 hn = Form("h_timeF%d_scint", iF);
294 ht = Form("Time distribution for Scintillator of %s; T_rec-T_0-T_fly-T_propagation[ns]", iFstring[iF].Data());
295 h_timeF_scint[iF] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
296 m_UpperTimeBoundaryScintillatorsBKLM);
297 hn = Form("h_timeF%d_scint_end", iF);
298 ht = Form("Time distribution for Scintillator of %s (Endcap); T_rec-T_0-T_fly-T_propagation[ns]", iFstring[iF].Data());
299 h_timeF_scint_end[iF] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
300 m_UpperTimeBoundaryScintillatorsEKLM);
301
302 hn = Form("h2_timeF%d_rpc", iF);
303 ht = Form("Time distribution for RPC of %s; Sector Index; T_rec-T_0-T_fly-T_propagation[ns]", iFstring[iF].Data());
304 h2_timeF_rpc[iF] = new TH2F(hn.Data(), ht.Data(), 8, 0, 8, nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
305 hn = Form("h2_timeF%d_scint", iF);
306 ht = Form("Time distribution for Scintillator of %s; Sector Index; T_rec-T_0-T_fly-T_propagation[ns]", iFstring[iF].Data());
307 h2_timeF_scint[iF] = new TH2F(hn.Data(), ht.Data(), 8, 0, 8, nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
308 m_UpperTimeBoundaryScintillatorsBKLM);
309 hn = Form("h2_timeF%d_scint_end", iF);
310 ht = Form("Time distribution for Scintillator of %s (Endcap); Sector Index; T_rec-T_0-T_fly-T_propagation[ns]",
311 iFstring[iF].Data());
312 h2_timeF_scint_end[iF] = new TH2F(hn.Data(), ht.Data(), 4, 0, 4, nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
313 m_UpperTimeBoundaryScintillatorsEKLM);
314
315 hn = Form("hc_timeF%d_rpc", iF);
316 ht = Form("Calibrated time distribution for RPC of %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]", iFstring[iF].Data());
317 hc_timeF_rpc[iF] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryCalibratedRPC, m_UpperTimeBoundaryCalibratedRPC);
318 hn = Form("hc_timeF%d_scint", iF);
319 ht = Form("Calibrated time distribution for Scintillator of %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
320 iFstring[iF].Data());
321 hc_timeF_scint[iF] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
322 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
323 hn = Form("hc_timeF%d_scint_end", iF);
324 ht = Form("Calibrated time distribution for Scintillator of %s (Endcap); T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
325 iFstring[iF].Data());
326 hc_timeF_scint_end[iF] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
327 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
328
329 hn = Form("h2c_timeF%d_rpc", iF);
330 ht = Form("Calibrated time distribution for RPC of %s; Sector Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
331 iFstring[iF].Data());
332 h2c_timeF_rpc[iF] = new TH2F(hn.Data(), ht.Data(), 8, 0, 8, nBin, m_LowerTimeBoundaryCalibratedRPC,
333 m_UpperTimeBoundaryCalibratedRPC);
334 hn = Form("h2c_timeF%d_scint", iF);
335 ht = Form("Calibrated time distribution for Scintillator of %s; Sector Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
336 iFstring[iF].Data());
337 h2c_timeF_scint[iF] = new TH2F(hn.Data(), ht.Data(), 8, 0, 8, nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
338 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
339 hn = Form("h2c_timeF%d_scint_end", iF);
340 ht = Form("Calibrated time distribution for Scintillator of %s (Endcap) ; Sector Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
341 iFstring[iF].Data());
342 h2c_timeF_scint_end[iF] = new TH2F(hn.Data(), ht.Data(), 4, 0, 4, nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
343 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
344
345 for (int iS = 0; iS < 8; ++iS) {
346 // Barrel parts
347 hn = Form("h_timeF%d_S%d_scint", iF, iS);
348 ht = Form("Time distribution for Scintillator of Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iS, iFstring[iF].Data());
349 h_timeFS_scint[iF][iS] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
350 m_UpperTimeBoundaryScintillatorsBKLM);
351 hn = Form("h_timeF%d_S%d_rpc", iF, iS);
352 ht = Form("Time distribution for RPC of Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iS, iFstring[iF].Data());
353 h_timeFS_rpc[iF][iS] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
354 hn = Form("h2_timeF%d_S%d", iF, iS);
355 ht = Form("Time distribution of Sector%d, %s; Layer Index; T_rec-T_0-T_fly-T_propagation[ns]", iS, iFstring[iF].Data());
356 h2_timeFS[iF][iS] = new TH2F(hn.Data(), ht.Data(), 15, 0, 15, nBin_scint, m_LowerTimeBoundaryRPC,
357 m_UpperTimeBoundaryScintillatorsBKLM);
358
359 hn = Form("hc_timeF%d_S%d_scint", iF, iS);
360 ht = Form("Calibrated time distribution for Scintillator of Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]", iS,
361 iFstring[iF].Data());
362 hc_timeFS_scint[iF][iS] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
363 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
364 hn = Form("hc_timeF%d_S%d_rpc", iF, iS);
365 ht = Form("Calibrated time distribution for RPC of Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]", iS,
366 iFstring[iF].Data());
367 hc_timeFS_rpc[iF][iS] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedRPC,
368 m_UpperTimeBoundaryCalibratedRPC);
369 hn = Form("h2c_timeF%d_S%d", iF, iS);
370 ht = Form("Calibrated time distribution of Sector%d, %s; Layer Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]", iS,
371 iFstring[iF].Data());
372 h2c_timeFS[iF][iS] = new TH2F(hn.Data(), ht.Data(), 15, 0, 15, nBin_scint, m_LowerTimeBoundaryCalibratedRPC,
373 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
374
375 // Inner 2 layers --> Scintillators
376 for (int iL = 0; iL < 2; ++iL) {
377 hn = Form("h_timeF%d_S%d_L%d", iF, iS, iL);
378 ht = Form("Time distribution for Scintillator of Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iL, iS,
379 iFstring[iF].Data());
380 h_timeFSL[iF][iS][iL] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
381 m_UpperTimeBoundaryScintillatorsBKLM);
382 hn = Form("hc_timeF%d_S%d_L%d", iF, iS, iL);
383 ht = Form("Calibrated time distribution for Scintillator of Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
384 iL, iS, iFstring[iF].Data());
385 hc_timeFSL[iF][iS][iL] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
386 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
387
388 for (int iP = 0; iP < 2; ++iP) {
389 hn = Form("h_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
390 ht = Form("Time distribution for Scintillator of %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]",
391 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
392 h_timeFSLP[iF][iS][iL][iP] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
393 m_UpperTimeBoundaryScintillatorsBKLM);
394 hn = Form("h2_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
395 ht = Form("Time distribution for Scintillator of %s, Layer%d, Sector%d, %s; Channel Index; T_rec-T_0-T_fly-T_propagation[ns]",
396 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
397 h2_timeFSLP[iF][iS][iL][iP] = new TH2F(hn.Data(), ht.Data(), 54, 0, 54, nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
398 m_UpperTimeBoundaryScintillatorsBKLM);
399
400 hn = Form("hc_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
401 ht = Form("Calibrated time distribution for Scintillator of %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
402 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
403 hc_timeFSLP[iF][iS][iL][iP] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
404 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
405 hn = Form("h2c_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
406 ht = Form("Calibrated time distribution for Scintillator of %s, Layer%d, Sector%d, %s; Channel Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
407 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
408 h2c_timeFSLP[iF][iS][iL][iP] = new TH2F(hn.Data(), ht.Data(), 54, 0, 54, nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
409 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
410
411 int nchannel_max = BKLMElementNumbers::getNStrips(iF, iS + 1, iL + 1, iP);
412 for (int iC = 0; iC < nchannel_max; ++iC) {
413 hn = Form("h_timeF%d_S%d_L%d_P%d_C%d_tc", iF, iS, iL, iP, iC);
414 ht = Form("time distribution for Scintillator of Channel%d, %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iC,
415 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
416 h_timeFSLPC_tc[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
417 m_UpperTimeBoundaryScintillatorsBKLM);
418
419 hn = Form("h_timeF%d_S%d_L%d_P%d_C%d", iF, iS, iL, iP, iC);
420 ht = Form("time distribution for Scintillator of Channel%d, %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iC,
421 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
422 h_timeFSLPC[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsBKLM,
423 m_UpperTimeBoundaryScintillatorsBKLM);
424
425 hn = Form("hc_timeF%d_S%d_L%d_P%d_C%d", iF, iS, iL, iP, iC);
426 ht = Form("Calibrated time distribution for Scintillator of Channel%d, %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
427 iC, iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
428 hc_timeFSLPC[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
429 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
430 hn = Form("time_length_bklm_F%d_S%d_L%d_P%d_C%d", iF, iS, iL, iP, iC);
431 double stripLength = 200;
432 m_HistTimeLengthBKLM[iF][iS][iL][iP][iC] =
433 new TH2F(hn.Data(),
434 "Time versus propagation length; "
435 "propagation distance[cm]; "
436 "T_rec-T_0-T_fly-'T_calibration'[ns]",
437 200, 0.0, stripLength,
438 400, m_LowerTimeBoundaryCalibratedScintillatorsBKLM,
439 m_UpperTimeBoundaryCalibratedScintillatorsBKLM);
440 }
441 }
442 }
443
444 for (int iL = 2; iL < 15; ++iL) {
445 hn = Form("h_timeF%d_S%d_L%d", iF, iS, iL);
446 ht = Form("time distribution for RPC of Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iL, iS, iFstring[iF].Data());
447 h_timeFSL[iF][iS][iL] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
448
449 hn = Form("hc_timeF%d_S%d_L%d", iF, iS, iL);
450 ht = Form("Calibrated time distribution for RPC of Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]", iL, iS,
451 iFstring[iF].Data());
452 hc_timeFSL[iF][iS][iL] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryCalibratedRPC, m_UpperTimeBoundaryCalibratedRPC);
453
454 for (int iP = 0; iP < 2; ++iP) {
455 hn = Form("h_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
456 ht = Form("time distribution for RPC of %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iPstring[iP].Data(), iL, iS,
457 iFstring[iF].Data());
458 h_timeFSLP[iF][iS][iL][iP] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
459
460 hn = Form("h2_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
461 ht = Form("time distribution for RPC of %s, Layer%d, Sector%d, %s; Channel Index; T_rec-T_0-T_fly-T_propagation[ns]",
462 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
463 h2_timeFSLP[iF][iS][iL][iP] = new TH2F(hn.Data(), ht.Data(), 48, 0, 48, nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
464
465 hn = Form("hc_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
466 ht = Form("Calibrated time distribution for RPC of %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
467 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
468 hc_timeFSLP[iF][iS][iL][iP] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryCalibratedRPC,
469 m_UpperTimeBoundaryCalibratedRPC);
470
471 hn = Form("h2c_timeF%d_S%d_L%d_P%d", iF, iS, iL, iP);
472 ht = Form("Calibrated time distribution for RPC of %s, Layer%d, Sector%d, %s; Channel Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
473 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
474 h2c_timeFSLP[iF][iS][iL][iP] = new TH2F(hn.Data(), ht.Data(), 48, 0, 48, nBin, m_LowerTimeBoundaryCalibratedRPC,
475 m_UpperTimeBoundaryCalibratedRPC);
476
477 int nchannel_max = BKLMElementNumbers::getNStrips(iF, iS + 1, iL + 1, iP);
478 for (int iC = 0; iC < nchannel_max; ++iC) {
479 hn = Form("h_timeF%d_S%d_L%d_P%d_C%d_tc", iF, iS, iL, iP, iC);
480 ht = Form("Time distribution for RPC of Channel%d, %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iC,
481 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
482 h_timeFSLPC_tc[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
483
484 hn = Form("h_timeF%d_S%d_L%d_P%d_C%d", iF, iS, iL, iP, iC);
485 ht = Form("Time distribution for RPC of Channel%d, %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation[ns]", iC,
486 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
487 h_timeFSLPC[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryRPC, m_UpperTimeBoundaryRPC);
488
489 hn = Form("hc_timeF%d_S%d_L%d_P%d_C%d", iF, iS, iL, iP, iC);
490 ht = Form("Calibrated time distribution for RPC of Channel%d, %s, Layer%d, Sector%d, %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
491 iC, iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
492 hc_timeFSLPC[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin, m_LowerTimeBoundaryCalibratedRPC,
493 m_UpperTimeBoundaryCalibratedRPC);
494 }
495 }
496 }
497 }
498 // Endcap part
499 int maxLay = 12 + 2 * iF;
500 for (int iS = 0; iS < 4; ++iS) {
501 hn = Form("h_timeF%d_S%d_scint_end", iF, iS);
502 ht = Form("Time distribution for Scintillator of Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation[ns]", iS,
503 iFstring[iF].Data());
504 h_timeFS_scint_end[iF][iS] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
505 m_UpperTimeBoundaryScintillatorsEKLM);
506 hn = Form("h2_timeF%d_S%d_end", iF, iS);
507 ht = Form("Time distribution of Sector%d, %s (Endcap); Layer Index; T_rec-T_0-T_fly-T_propagation[ns]", iS, iFstring[iF].Data());
508 h2_timeFS_end[iF][iS] = new TH2F(hn.Data(), ht.Data(), maxLay, 0, maxLay, nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
509 m_UpperTimeBoundaryScintillatorsEKLM);
510 hn = Form("hc_timeF%d_S%d_scint_end", iF, iS);
511 ht = Form("Calibrated time distribution for Scintillator of Sector%d (Endcap), %s; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
512 iS, iFstring[iF].Data());
513 hc_timeFS_scint_end[iF][iS] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
514 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
515 hn = Form("h2c_timeF%d_S%d_end", iF, iS);
516 ht = Form("Calibrated time distribution of Sector%d, %s (Endcap); Layer Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
517 iS, iFstring[iF].Data());
518 h2c_timeFS_end[iF][iS] = new TH2F(hn.Data(), ht.Data(), maxLay, 0, maxLay, nBin_scint,
519 m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
520 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
521
522 for (int iL = 0; iL < maxLay; ++iL) {
523 hn = Form("h_timeF%d_S%d_L%d_end", iF, iS, iL);
524 ht = Form("Time distribution for Scintillator of Layer%d, Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation[ns]", iL, iS,
525 iFstring[iF].Data());
526 h_timeFSL_end[iF][iS][iL] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
527 m_UpperTimeBoundaryScintillatorsEKLM);
528 hn = Form("hc_timeF%d_S%d_L%d_end", iF, iS, iL);
529 ht = Form("Calibrated time distribution for Scintillator of Layer%d, Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
530 iL, iS, iFstring[iF].Data());
531 hc_timeFSL_end[iF][iS][iL] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
532 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
533
534 for (int iP = 0; iP < 2; ++iP) {
535 hn = Form("h_timeF%d_S%d_L%d_P%d_end", iF, iS, iL, iP);
536 ht = Form("Time distribution for Scintillator of %s, Layer%d, Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation[ns]",
537 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
538 h_timeFSLP_end[iF][iS][iL][iP] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
539 m_UpperTimeBoundaryScintillatorsEKLM);
540
541 hn = Form("h2_timeF%d_S%d_L%d_P%d_end", iF, iS, iL, iP);
542 ht = Form("Time distribution for Scintillator of %s, Layer%d, Sector%d, %s (Endcap); Channel Index; T_rec-T_0-T_fly-T_propagation[ns]",
543 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
544 h2_timeFSLP_end[iF][iS][iL][iP] = new TH2F(hn.Data(), ht.Data(), 75, 0, 75, nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
545 m_UpperTimeBoundaryScintillatorsEKLM);
546
547 hn = Form("hc_timeF%d_S%d_L%d_P%d_end", iF, iS, iL, iP);
548 ht = Form("Calibrated time distribution for Scintillator of %s, Layer%d, Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
549 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
550 hc_timeFSLP_end[iF][iS][iL][iP] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
551 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
552
553 hn = Form("h2c_timeF%d_S%d_L%d_P%d_end", iF, iS, iL, iP);
554 ht = Form("Calibrated time distribution for Scintillator of %s, Layer%d, Sector%d, %s (Endcap); Channel Index; T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
555 iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
556 h2c_timeFSLP_end[iF][iS][iL][iP] = new TH2F(hn.Data(), ht.Data(), 75, 0, 75, nBin_scint,
557 m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
558 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
559
560 for (int iC = 0; iC < 75; ++iC) {
561 hn = Form("h_timeF%d_S%d_L%d_P%d_C%d_tc_end", iF, iS, iL, iP, iC);
562 ht = Form("Time distribution for Scintillator of Channel%d, %s, Layer%d, Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation[ns]",
563 iC, iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
564 h_timeFSLPC_tc_end[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
565 m_UpperTimeBoundaryScintillatorsEKLM);
566
567 hn = Form("h_timeF%d_S%d_L%d_P%d_C%d_end", iF, iS, iL, iP, iC);
568 ht = Form("Time distribution for Scintillator of Channel%d, %s, Layer%d, Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation[ns]",
569 iC, iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
570 h_timeFSLPC_end[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryScintillatorsEKLM,
571 m_UpperTimeBoundaryScintillatorsEKLM);
572 hn = Form("hc_timeF%d_S%d_L%d_P%d_C%d_end", iF, iS, iL, iP, iC);
573 ht = Form("Calibrated time distribution for Scintillator of Channel%d, %s, Layer%d, Sector%d, %s (Endcap); T_rec-T_0-T_fly-T_propagation-T_calibration[ns]",
574 iC, iPstring[iP].Data(), iL, iS, iFstring[iF].Data());
575 hc_timeFSLPC_end[iF][iS][iL][iP][iC] = new TH1F(hn.Data(), ht.Data(), nBin_scint, m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
576 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
577 hn = Form("time_length_eklm_F%d_S%d_L%d_P%d_C%d", iF, iS, iL, iP, iC);
578 double stripLength = m_EKLMGeometry->getStripLength(iC + 1) /
579 CLHEP::cm * Unit::cm;
580 m_HistTimeLengthEKLM[iF][iS][iL][iP][iC] =
581 new TH2F(hn.Data(),
582 "Time versus propagation length; "
583 "propagation distance[cm]; "
584 "T_rec-T_0-T_fly-'T_calibration'[ns]",
585 200, 0.0, stripLength,
586 400, m_LowerTimeBoundaryCalibratedScintillatorsEKLM,
587 m_UpperTimeBoundaryCalibratedScintillatorsEKLM);
588 }
589 }
590 }
591 }
592 }
593}
594
595void KLMTimeAlgorithm::fillTimeDistanceProfiles(
596 TProfile* profileRpcPhi, TProfile* profileRpcZ,
597 TProfile* profileBKLMScintillatorPhi, TProfile* profileBKLMScintillatorZ,
598 TProfile* profileEKLMScintillatorPlane1,
599 TProfile* profileEKLMScintillatorPlane2, bool fill2dHistograms)
600{
601 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
602 KLMChannelNumber channel = klmChannel.getKLMChannelNumber();
603 if (m_cFlag[channel] == ChannelCalibrationStatus::c_NotEnoughData)
604 continue;
605
606 std::vector<struct Event> eventsChannel;
607 eventsChannel = m_evts[channel];
608 int iSub = klmChannel.getSubdetector();
609
610 for (const Event& event : eventsChannel) {
611 double timeHit = event.time() - m_timeShift[channel];
612 if (m_useEventT0)
613 timeHit = timeHit - event.t0;
614 double distHit = event.dist;
615
616 if (iSub == KLMElementNumbers::c_BKLM) {
617 int iF = klmChannel.getSection();
618 int iS = klmChannel.getSector() - 1;
619 int iL = klmChannel.getLayer() - 1;
620 int iP = klmChannel.getPlane();
621 int iC = klmChannel.getStrip() - 1;
622 if (iL > 1) {
623 if (iP) {
624 profileRpcPhi->Fill(distHit, timeHit);
625 } else {
626 profileRpcZ->Fill(distHit, timeHit);
627 }
628 } else {
629 if (fill2dHistograms)
630 m_HistTimeLengthBKLM[iF][iS][iL][iP][iC]->Fill(distHit, timeHit);
631 if (iP) {
632 profileBKLMScintillatorPhi->Fill(distHit, timeHit);
633 } else {
634 profileBKLMScintillatorZ->Fill(distHit, timeHit);
635 }
636 }
637 } else {
638 int iF = klmChannel.getSection() - 1;
639 int iS = klmChannel.getSector() - 1;
640 int iL = klmChannel.getLayer() - 1;
641 int iP = klmChannel.getPlane() - 1;
642 int iC = klmChannel.getStrip() - 1;
643 if (fill2dHistograms)
644 m_HistTimeLengthEKLM[iF][iS][iL][iP][iC]->Fill(distHit, timeHit);
645 if (iP) {
646 profileEKLMScintillatorPlane1->Fill(distHit, timeHit);
647 } else {
648 profileEKLMScintillatorPlane2->Fill(distHit, timeHit);
649 }
650 }
651 }
652 }
653}
654
655void KLMTimeAlgorithm::timeDistance2dFit(
656 const std::vector< std::pair<KLMChannelNumber, unsigned int> >& channels,
657 double& delay, double& delayError)
658{
659 std::vector<struct Event> eventsChannel;
660 int nFits = 1000;
661 int nConvergedFits = 0;
662 delay = 0;
663 delayError = 0;
664 if (nFits > (int)channels.size())
665 nFits = channels.size();
666 for (int i = 0; i < nFits; ++i) {
667 int subdetector, section, sector, layer, plane, strip;
668 m_ElementNumbers->channelNumberToElementNumbers(
669 channels[i].first, &subdetector, &section, &sector, &layer, &plane,
670 &strip);
671 if (subdetector == KLMElementNumbers::c_BKLM) {
672 s_LowerTimeBoundary = m_LowerTimeBoundaryScintillatorsBKLM;
673 s_UpperTimeBoundary = m_UpperTimeBoundaryScintillatorsBKLM;
674 const bklm::Module* module =
675 m_BKLMGeometry->findModule(section, sector, layer);
676 s_StripLength = module->getStripLength(plane, strip);
677 } else {
678 s_LowerTimeBoundary = m_LowerTimeBoundaryScintillatorsEKLM;
679 s_UpperTimeBoundary = m_UpperTimeBoundaryScintillatorsEKLM;
680 s_StripLength = m_EKLMGeometry->getStripLength(strip) /
681 CLHEP::cm * Unit::cm;
682 }
683 for (int j = 0; j < c_NBinsTime; ++j) {
684 for (int k = 0; k < c_NBinsDistance; ++k)
685 s_BinnedData[j][k] = 0;
686 }
687 eventsChannel = m_evts[channels[i].first];
688 double averageTime = 0;
689 for (const Event& event : eventsChannel) {
690 double timeHit = event.time();
691 if (m_useEventT0)
692 timeHit = timeHit - event.t0;
693 averageTime = averageTime + timeHit;
694 int timeBin = std::floor((timeHit - s_LowerTimeBoundary) * c_NBinsTime /
695 (s_UpperTimeBoundary - s_LowerTimeBoundary));
696 if (timeBin < 0 || timeBin >= c_NBinsTime)
697 continue;
698 int distanceBin = std::floor(event.dist * c_NBinsDistance / s_StripLength);
699 if (distanceBin < 0 || distanceBin >= c_NBinsDistance) {
700 B2ERROR("The distance to SiPM is greater than the strip length.");
701 continue;
702 }
703 s_BinnedData[timeBin][distanceBin] += 1;
704 }
705 averageTime = averageTime / eventsChannel.size();
706 TMinuit minuit(5);
707 minuit.SetPrintLevel(-1);
708 int minuitResult;
709 minuit.mnparm(0, "P0", 1, 0.001, 0, 0, minuitResult);
710 minuit.mnparm(1, "N", 10, 0.001, 0, 0, minuitResult);
711 minuit.mnparm(2, "T0", averageTime, 0.001, 0, 0, minuitResult);
712 minuit.mnparm(3, "SIGMA", 10, 0.001, 0, 0, minuitResult);
713 minuit.mnparm(4, "DELAY", 0.0, 0.001, 0, 0, minuitResult);
714 minuit.SetFCN(fcn);
715 minuit.mncomd("FIX 2 3 4 5", minuitResult);
716 minuit.mncomd("MIGRAD 10000", minuitResult);
717 minuit.mncomd("RELEASE 2", minuitResult);
718 minuit.mncomd("MIGRAD 10000", minuitResult);
719 minuit.mncomd("RELEASE 3", minuitResult);
720 minuit.mncomd("MIGRAD 10000", minuitResult);
721 minuit.mncomd("RELEASE 4", minuitResult);
722 minuit.mncomd("MIGRAD 10000", minuitResult);
723 minuit.mncomd("RELEASE 5", minuitResult);
724 minuit.mncomd("MIGRAD 10000", minuitResult);
725 /* Require converged fit with accurate error matrix. */
726 if (minuit.fISW[1] != 3)
727 continue;
728 nConvergedFits++;
729 double channelDelay, channelDelayError;
730 minuit.GetParameter(4, channelDelay, channelDelayError);
731 delay = delay + channelDelay;
732 delayError = delayError + channelDelayError * channelDelayError;
733 }
734 delay = delay / nConvergedFits;
735 delayError = sqrt(delayError) / (nConvergedFits - 1);
736}
737
738CalibrationAlgorithm::EResult KLMTimeAlgorithm::calibrate()
739{
740 int channelId;
741 gROOT->SetBatch(kTRUE);
742 setupDatabase();
743 m_timeCableDelay = new KLMTimeCableDelay();
744 m_timeConstants = new KLMTimeConstants();
745 m_timeResolution = new KLMTimeResolution();
746
747 fcn_gaus = new TF1("fcn_gaus", "gaus");
748 fcn_land = new TF1("fcn_land", "landau");
749 fcn_pol1 = new TF1("fcn_pol1", "pol1");
750 fcn_const = new TF1("fcn_const", "pol0");
751
752 CalibrationAlgorithm::EResult result = readCalibrationData();
753 if (result != CalibrationAlgorithm::c_OK)
754 return result;
755
756 /* Choose non-existing file name. */
757 std::string name = "time_calibration.root";
758 int i = 1;
759 while (1) {
760 struct stat buffer;
761 if (stat(name.c_str(), &buffer) != 0)
762 break;
763 name = "time_calibration_" + std::to_string(i) + ".root";
764 i = i + 1;
765 /* Overflow. */
766 if (i < 0)
767 break;
768 }
769 m_outFile = new TFile(name.c_str(), "recreate");
770 createHistograms();
771
772 std::vector<struct Event> eventsChannel;
773
774 eventsChannel.clear();
775 m_cFlag.clear();
776 m_minimizerOptions.SetDefaultStrategy(2);
777
778 /* Sort channels by number of events. */
779 std::vector< std::pair<KLMChannelNumber, unsigned int> > channelsBKLM;
780 std::vector< std::pair<KLMChannelNumber, unsigned int> > channelsEKLM;
781 KLMChannelIndex klmChannels;
782 for (KLMChannelIndex& klmChannel : klmChannels) {
783 KLMChannelNumber channel = klmChannel.getKLMChannelNumber();
784 m_cFlag[channel] = ChannelCalibrationStatus::c_NotEnoughData;
785 if (m_evts.find(channel) == m_evts.end())
786 continue;
787 int nEvents = m_evts[channel].size();
788 if (nEvents < m_lower_limit_counts) {
789 B2WARNING("Not enough calibration data collected."
790 << LogVar("channel", channel)
791 << LogVar("number of digit", nEvents));
792 continue;
793 }
794 m_cFlag[channel] = ChannelCalibrationStatus::c_FailedFit;
795 if (klmChannel.getSubdetector() == KLMElementNumbers::c_BKLM &&
796 klmChannel.getLayer() < BKLMElementNumbers::c_FirstRPCLayer) {
797 channelsBKLM.push_back(
798 std::pair<KLMChannelNumber, unsigned int>(channel, nEvents));
799 }
800 if (klmChannel.getSubdetector() == KLMElementNumbers::c_EKLM) {
801 channelsEKLM.push_back(
802 std::pair<KLMChannelNumber, unsigned int>(channel, nEvents));
803 }
804 }
805 std::sort(channelsBKLM.begin(), channelsBKLM.end(), compareEventNumber);
806 std::sort(channelsEKLM.begin(), channelsEKLM.end(), compareEventNumber);
807
808 /* Two-dimensional fit for the channel with the maximal number of events. */
809 double delayBKLM, delayBKLMError;
810 double delayEKLM, delayEKLMError;
811 timeDistance2dFit(channelsBKLM, delayBKLM, delayBKLMError);
812 timeDistance2dFit(channelsEKLM, delayEKLM, delayEKLMError);
813
814 /**********************************
815 * First loop
816 * Estimation of effective light speed for Scintillators and RPCs, separately.
817 **********************************/
818 B2INFO("Effective light speed Estimation.");
819 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
820 channelId = klmChannel.getKLMChannelNumber();
821 if (m_cFlag[channelId] == ChannelCalibrationStatus::c_NotEnoughData)
822 continue;
823
824 eventsChannel = m_evts[channelId];
825 int iSub = klmChannel.getSubdetector();
826
827 for (const Event& event : eventsChannel) {
828 XYZVector diffD = XYZVector(event.diffDistX, event.diffDistY, event.diffDistZ);
829 h_diff->Fill(diffD.R());
830 double timeHit = event.time();
831 if (m_useEventT0)
832 timeHit = timeHit - event.t0;
833 if (iSub == KLMElementNumbers::c_BKLM) {
834 int iF = klmChannel.getSection();
835 int iS = klmChannel.getSector() - 1;
836 int iL = klmChannel.getLayer() - 1;
837 int iP = klmChannel.getPlane();
838 int iC = klmChannel.getStrip() - 1;
839 h_timeFSLPC_tc[iF][iS][iL][iP][iC]->Fill(timeHit);
840 if (iL > 1) {
841 h_time_rpc_tc->Fill(timeHit);
842 } else {
843 h_time_scint_tc->Fill(timeHit);
844 }
845 } else {
846 int iF = klmChannel.getSection() - 1;
847 int iS = klmChannel.getSector() - 1;
848 int iL = klmChannel.getLayer() - 1;
849 int iP = klmChannel.getPlane() - 1;
850 int iC = klmChannel.getStrip() - 1;
851 h_timeFSLPC_tc_end[iF][iS][iL][iP][iC]->Fill(timeHit);
852 h_time_scint_tc_end->Fill(timeHit);
853 }
854 }
855 }
856 B2INFO("Effective light speed Estimation! Hists and Graph filling done.");
857
858 m_timeShift.clear();
859
860 double tmpMean_rpc_global = h_time_rpc_tc->GetMean();
861 double tmpMean_scint_global = h_time_scint_tc->GetMean();
862 double tmpMean_scint_global_end = h_time_scint_tc_end->GetMean();
863
864 B2INFO("Global Mean for Raw." << LogVar("RPC", tmpMean_rpc_global) << LogVar("Scint BKLM",
865 tmpMean_scint_global) << LogVar("Scint EKLM", tmpMean_scint_global_end));
866
867 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
868 channelId = klmChannel.getKLMChannelNumber();
869 if (m_cFlag[channelId] == ChannelCalibrationStatus::c_NotEnoughData)
870 continue;
871
872 int iSub = klmChannel.getSubdetector();
873 if (iSub == KLMElementNumbers::c_BKLM) {
874 int iF = klmChannel.getSection();
875 int iS = klmChannel.getSector() - 1;
876 int iL = klmChannel.getLayer() - 1;
877 int iP = klmChannel.getPlane();
878 int iC = klmChannel.getStrip() - 1;
879 h_timeFSLPC_tc[iF][iS][iL][iP][iC]->Fit(fcn_gaus, "LESQ");
880 double tmpMean_channel = fcn_gaus->GetParameter(1);
881 if (iL > 1) {
882 m_timeShift[channelId] = tmpMean_channel - tmpMean_rpc_global;
883 } else {
884 m_timeShift[channelId] = tmpMean_channel - tmpMean_scint_global;
885 }
886 } else {
887 int iF = klmChannel.getSection() - 1;
888 int iS = klmChannel.getSector() - 1;
889 int iL = klmChannel.getLayer() - 1;
890 int iP = klmChannel.getPlane() - 1;
891 int iC = klmChannel.getStrip() - 1;
892 h_timeFSLPC_tc_end[iF][iS][iL][iP][iC]->Fit(fcn_gaus, "LESQ");
893 double tmpMean_channel = fcn_gaus->GetParameter(1);
894 m_timeShift[channelId] = tmpMean_channel - tmpMean_scint_global_end;
895 }
896 }
897
898 delete h_time_scint_tc;
899 delete h_time_scint_tc_end;
900 delete h_time_rpc_tc;
901 B2INFO("Effective Light m_timeShift obtained. done.");
902
903 fillTimeDistanceProfiles(
904 m_ProfileRpcPhi, m_ProfileRpcZ,
905 m_ProfileBKLMScintillatorPhi, m_ProfileBKLMScintillatorZ,
906 m_ProfileEKLMScintillatorPlane1, m_ProfileEKLMScintillatorPlane2, false);
907
908 B2INFO("Effective light speed fitting.");
909 m_ProfileRpcPhi->Fit("fcn_pol1", "EMQ");
910 double delayRPCPhi = fcn_pol1->GetParameter(1);
911 double e_slope_rpc_phi = fcn_pol1->GetParError(1);
912
913 m_ProfileRpcZ->Fit("fcn_pol1", "EMQ");
914 double delayRPCZ = fcn_pol1->GetParameter(1);
915 double e_slope_rpc_z = fcn_pol1->GetParError(1);
916
917 m_ProfileBKLMScintillatorPhi->Fit("fcn_pol1", "EMQ");
918 double slope_scint_phi = fcn_pol1->GetParameter(1);
919 double e_slope_scint_phi = fcn_pol1->GetParError(1);
920
921 m_ProfileBKLMScintillatorZ->Fit("fcn_pol1", "EMQ");
922 double slope_scint_z = fcn_pol1->GetParameter(1);
923 double e_slope_scint_z = fcn_pol1->GetParError(1);
924
925 m_ProfileEKLMScintillatorPlane1->Fit("fcn_pol1", "EMQ");
926 double slope_scint_plane1_end = fcn_pol1->GetParameter(1);
927 double e_slope_scint_plane1_end = fcn_pol1->GetParError(1);
928
929 m_ProfileEKLMScintillatorPlane2->Fit("fcn_pol1", "EMQ");
930 double slope_scint_plane2_end = fcn_pol1->GetParameter(1);
931 double e_slope_scint_plane2_end = fcn_pol1->GetParError(1);
932
933 TString logStr_phi, logStr_z;
934 logStr_phi = Form("%.4f +/- %.4f ns/cm", delayRPCPhi, e_slope_rpc_phi);
935 logStr_z = Form("%.4f +/- %.4f ns/cm", delayRPCZ, e_slope_rpc_z);
936 B2INFO("Delay in RPCs:"
937 << LogVar("Fitted Value (phi readout) ", logStr_phi.Data())
938 << LogVar("Fitted Value (z readout) ", logStr_z.Data()));
939 logStr_phi = Form("%.4f +/- %.4f ns/cm", slope_scint_phi, e_slope_scint_phi);
940 logStr_z = Form("%.4f +/- %.4f ns/cm", slope_scint_z, e_slope_scint_z);
941 B2INFO("Delay in BKLM scintillators:"
942 << LogVar("Fitted Value (phi readout) ", logStr_phi.Data())
943 << LogVar("Fitted Value (z readout) ", logStr_z.Data()));
944 logStr_phi = Form("%.4f +/- %.4f ns/cm", slope_scint_plane1_end,
945 e_slope_scint_plane1_end);
946 logStr_z = Form("%.4f +/- %.4f ns/cm", slope_scint_plane2_end,
947 e_slope_scint_plane2_end);
948 B2INFO("Delay in EKLM scintillators:"
949 << LogVar("Fitted Value (plane1 readout) ", logStr_phi.Data())
950 << LogVar("Fitted Value (plane2 readout) ", logStr_z.Data()));
951
952 logStr_z = Form("%.4f +/- %.4f ns/cm", delayBKLM, delayBKLMError);
953 B2INFO("Delay in BKLM scintillators:"
954 << LogVar("Fitted Value (2d fit) ", logStr_z.Data()));
955 logStr_z = Form("%.4f +/- %.4f ns/cm", delayEKLM, delayEKLMError);
956 B2INFO("Delay in EKLM scintillators:"
957 << LogVar("Fitted Value (2d fit) ", logStr_z.Data()));
958
959 // Default Effective light speed in current Database
960 //delayEKLM = 0.5 * (slope_scint_plane1_end + slope_scint_plane2_end);
961 //delayBKLM = 0.5 * (slope_scint_phi + slope_scint_z);
962
963 m_timeConstants->setDelay(delayEKLM, KLMTimeConstants::c_EKLM);
964 m_timeConstants->setDelay(delayBKLM, KLMTimeConstants::c_BKLM);
965 m_timeConstants->setDelay(delayRPCPhi, KLMTimeConstants::c_RPCPhi);
966 m_timeConstants->setDelay(delayRPCZ, KLMTimeConstants::c_RPCZ);
967
969 B2INFO("Time distribution filling begins.");
970 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
971 channelId = klmChannel.getKLMChannelNumber();
972 int iSub = klmChannel.getSubdetector();
973
974 if (m_cFlag[channelId] == ChannelCalibrationStatus::c_NotEnoughData)
975 continue;
976 eventsChannel = m_evts[channelId];
977
978 for (const Event& event : eventsChannel) {
979 double timeHit = event.time();
980 if (m_useEventT0)
981 timeHit = timeHit - event.t0;
982 if (iSub == KLMElementNumbers::c_BKLM) {
983 int iF = klmChannel.getSection();
984 int iS = klmChannel.getSector() - 1;
985 int iL = klmChannel.getLayer() - 1;
986 int iP = klmChannel.getPlane();
987 int iC = klmChannel.getStrip() - 1;
988 if (iL > 1) {
989 double propgationT;
990 if (iP == BKLMElementNumbers::c_ZPlane)
991 propgationT = event.dist * delayRPCZ;
992 else
993 propgationT = event.dist * delayRPCPhi;
994 double time = timeHit - propgationT;
995 h_time_rpc->Fill(time);
996 h_timeF_rpc[iF]->Fill(time);
997 h_timeFS_rpc[iF][iS]->Fill(time);
998 h_timeFSL[iF][iS][iL]->Fill(time);
999 h_timeFSLP[iF][iS][iL][iP]->Fill(time);
1000 h_timeFSLPC[iF][iS][iL][iP][iC]->Fill(time);
1001 h2_timeF_rpc[iF]->Fill(iS, time);
1002 h2_timeFS[iF][iS]->Fill(iL, time);
1003 h2_timeFSLP[iF][iS][iL][iP]->Fill(iC, time);
1004 } else {
1005 double propgationT = event.dist * delayBKLM;
1006 double time = timeHit - propgationT;
1007 h_time_scint->Fill(time);
1008 h_timeF_scint[iF]->Fill(time);
1009 h_timeFS_scint[iF][iS]->Fill(time);
1010 h_timeFSL[iF][iS][iL]->Fill(time);
1011 h_timeFSLP[iF][iS][iL][iP]->Fill(time);
1012 h_timeFSLPC[iF][iS][iL][iP][iC]->Fill(time);
1013 h2_timeF_scint[iF]->Fill(iS, time);
1014 h2_timeFS[iF][iS]->Fill(iL, time);
1015 h2_timeFSLP[iF][iS][iL][iP]->Fill(iC, time);
1016 }
1017 } else {
1018 int iF = klmChannel.getSection() - 1;
1019 int iS = klmChannel.getSector() - 1;
1020 int iL = klmChannel.getLayer() - 1;
1021 int iP = klmChannel.getPlane() - 1;
1022 int iC = klmChannel.getStrip() - 1;
1023 double propgationT = event.dist * delayEKLM;
1024 double time = timeHit - propgationT;
1025 h_time_scint_end->Fill(time);
1026 h_timeF_scint_end[iF]->Fill(time);
1027 h_timeFS_scint_end[iF][iS]->Fill(time);
1028 h_timeFSL_end[iF][iS][iL]->Fill(time);
1029 h_timeFSLP_end[iF][iS][iL][iP]->Fill(time);
1030 h_timeFSLPC_end[iF][iS][iL][iP][iC]->Fill(time);
1031 h2_timeF_scint_end[iF]->Fill(iS, time);
1032 h2_timeFS_end[iF][iS]->Fill(iL, time);
1033 h2_timeFSLP_end[iF][iS][iL][iP]->Fill(iC, time);
1034 }
1035 }
1036 }
1037
1038 B2INFO("Original filling done.");
1039
1040 int iChannel_rpc = 0;
1041 int iChannel = 0;
1042 int iChannel_end = 0;
1043 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1044 channelId = klmChannel.getKLMChannelNumber();
1045 if (m_cFlag[channelId] == ChannelCalibrationStatus::c_NotEnoughData)
1046 continue;
1047 int iSub = klmChannel.getSubdetector();
1048
1049 if (iSub == KLMElementNumbers::c_BKLM) {
1050 int iF = klmChannel.getSection();
1051 int iS = klmChannel.getSector() - 1;
1052 int iL = klmChannel.getLayer() - 1;
1053 int iP = klmChannel.getPlane();
1054 int iC = klmChannel.getStrip() - 1;
1055
1056 TFitResultPtr r = h_timeFSLPC[iF][iS][iL][iP][iC]->Fit(fcn_gaus, "LESQ");
1057 if (int(r) != 0)
1058 continue;
1059 if (int(r) == 0)
1060 m_cFlag[channelId] = ChannelCalibrationStatus::c_SuccessfulCalibration;
1061 m_time_channel[channelId] = fcn_gaus->GetParameter(1);
1062 m_etime_channel[channelId] = fcn_gaus->GetParError(1);
1063 if (iL > 1) {
1064 gre_time_channel_rpc->SetPoint(iChannel_rpc, channelId, m_time_channel[channelId]);
1065 gre_time_channel_rpc->SetPointError(iChannel_rpc, 0., m_etime_channel[channelId]);
1066 iChannel_rpc++;
1067 } else {
1068 gre_time_channel_scint->SetPoint(iChannel, channelId, m_time_channel[channelId]);
1069 gre_time_channel_scint->SetPointError(iChannel, 0., m_etime_channel[channelId]);
1070 iChannel++;
1071 }
1072 } else {
1073 int iF = klmChannel.getSection() - 1;
1074 int iS = klmChannel.getSector() - 1;
1075 int iL = klmChannel.getLayer() - 1;
1076 int iP = klmChannel.getPlane() - 1;
1077 int iC = klmChannel.getStrip() - 1;
1078
1079 TFitResultPtr r = h_timeFSLPC_end[iF][iS][iL][iP][iC]->Fit(fcn_gaus, "LESQ");
1080 if (int(r) != 0)
1081 continue;
1082 if (int(r) == 0)
1083 m_cFlag[channelId] = ChannelCalibrationStatus::c_SuccessfulCalibration;
1084 m_time_channel[channelId] = fcn_gaus->GetParameter(1);
1085 m_etime_channel[channelId] = fcn_gaus->GetParError(1);
1086 gre_time_channel_scint_end->SetPoint(iChannel_end, channelId, m_time_channel[channelId]);
1087 gre_time_channel_scint_end->SetPointError(iChannel_end, 0., m_etime_channel[channelId]);
1088 iChannel_end++;
1089 }
1090 }
1091
1092 gre_time_channel_scint->Fit("fcn_const", "EMQ");
1093 m_time_channelAvg_scint = fcn_const->GetParameter(0);
1094 m_etime_channelAvg_scint = fcn_const->GetParError(0);
1095
1096 gre_time_channel_scint_end->Fit("fcn_const", "EMQ");
1097 m_time_channelAvg_scint_end = fcn_const->GetParameter(0);
1098 m_etime_channelAvg_scint_end = fcn_const->GetParError(0);
1099
1100 gre_time_channel_rpc->Fit("fcn_const", "EMQ");
1101 m_time_channelAvg_rpc = fcn_const->GetParameter(0);
1102 m_etime_channelAvg_rpc = fcn_const->GetParError(0);
1103
1104 B2INFO("Channel's time distribution fitting done.");
1105 B2DEBUG(20, LogVar("Average time (RPC)", m_time_channelAvg_rpc)
1106 << LogVar("Average time (BKLM scintillators)", m_time_channelAvg_scint)
1107 << LogVar("Average time (EKLM scintillators)", m_time_channelAvg_scint_end));
1108
1109 B2INFO("Calibrated channel's time distribution filling begins.");
1110
1111 m_timeShift.clear();
1112 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1113 channelId = klmChannel.getKLMChannelNumber();
1114 h_calibrated->Fill(m_cFlag[channelId]);
1115 if (m_time_channel.find(channelId) == m_time_channel.end())
1116 continue;
1117 double timeShift = m_time_channel[channelId];
1118 int iSub = klmChannel.getSubdetector();
1119 if (iSub == KLMElementNumbers::c_BKLM) {
1120 int iL = klmChannel.getLayer() - 1;
1121 if (iL > 1)
1122 m_timeShift[channelId] = timeShift;
1123 else
1124 m_timeShift[channelId] = timeShift;
1125 } else {
1126 m_timeShift[channelId] = timeShift;
1127 }
1128 m_timeCableDelay->setTimeDelay(channelId, m_timeShift[channelId]);
1129 }
1130
1131 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1132 channelId = klmChannel.getKLMChannelNumber();
1133 if (m_timeShift.find(channelId) != m_timeShift.end())
1134 continue;
1135 m_timeShift[channelId] = esti_timeShift(klmChannel);
1136 m_timeCableDelay->setTimeDelay(channelId, m_timeShift[channelId]);
1137 B2DEBUG(20, "Uncalibrated Estimation " << LogVar("Channel", channelId) << LogVar("Estimated value", m_timeShift[channelId]));
1138 }
1139
1140 iChannel_rpc = 0;
1141 iChannel = 0;
1142 iChannel_end = 0;
1143 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1144 channelId = klmChannel.getKLMChannelNumber();
1145 if (m_timeShift.find(channelId) == m_timeShift.end()) {
1146 B2ERROR("!!! Not All Channels Calibration Constant Set. Error Happened on " << LogVar("Channel", channelId));
1147 continue;
1148 }
1149 int iSub = klmChannel.getSubdetector();
1150 if (iSub == KLMElementNumbers::c_BKLM) {
1151 int iL = klmChannel.getLayer();
1152 if (iL > 2) {
1153 gr_timeShift_channel_rpc->SetPoint(iChannel_rpc, channelId, m_timeShift[channelId]);
1154 iChannel_rpc++;
1155 } else {
1156 gr_timeShift_channel_scint->SetPoint(iChannel, channelId, m_timeShift[channelId]);
1157 iChannel++;
1158 }
1159 } else {
1160 gr_timeShift_channel_scint_end->SetPoint(iChannel_end, channelId, m_timeShift[channelId]);
1161 iChannel_end++;
1162 }
1163 }
1164
1165 fillTimeDistanceProfiles(
1166 m_Profile2RpcPhi, m_Profile2RpcZ,
1167 m_Profile2BKLMScintillatorPhi, m_Profile2BKLMScintillatorZ,
1168 m_Profile2EKLMScintillatorPlane1, m_Profile2EKLMScintillatorPlane2, true);
1169
1170 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1171 channelId = klmChannel.getKLMChannelNumber();
1172 int iSub = klmChannel.getSubdetector();
1173 eventsChannel = m_evts[channelId];
1174 for (const Event& event : eventsChannel) {
1175 double timeHit = event.time();
1176 if (m_useEventT0)
1177 timeHit = timeHit - event.t0;
1178 if (iSub == KLMElementNumbers::c_BKLM) {
1179 int iF = klmChannel.getSection();
1180 int iS = klmChannel.getSector() - 1;
1181 int iL = klmChannel.getLayer() - 1;
1182 int iP = klmChannel.getPlane();
1183 int iC = klmChannel.getStrip() - 1;
1184 if (iL > 1) {
1185 double propgationT;
1186 if (iP == BKLMElementNumbers::c_ZPlane)
1187 propgationT = event.dist * delayRPCZ;
1188 else
1189 propgationT = event.dist * delayRPCPhi;
1190 double time = timeHit - propgationT - m_timeShift[channelId];
1191 hc_time_rpc->Fill(time);
1192 hc_timeF_rpc[iF]->Fill(time);
1193 hc_timeFS_rpc[iF][iS]->Fill(time);
1194 hc_timeFSL[iF][iS][iL]->Fill(time);
1195 hc_timeFSLP[iF][iS][iL][iP]->Fill(time);
1196 hc_timeFSLPC[iF][iS][iL][iP][iC]->Fill(time);
1197 h2c_timeF_rpc[iF]->Fill(iS, time);
1198 h2c_timeFS[iF][iS]->Fill(iL, time);
1199 h2c_timeFSLP[iF][iS][iL][iP]->Fill(iC, time);
1200 } else {
1201 double propgationT = event.dist * delayBKLM;
1202 double time = timeHit - propgationT - m_timeShift[channelId];
1203 hc_time_scint->Fill(time);
1204 hc_timeF_scint[iF]->Fill(time);
1205 hc_timeFS_scint[iF][iS]->Fill(time);
1206 hc_timeFSL[iF][iS][iL]->Fill(time);
1207 hc_timeFSLP[iF][iS][iL][iP]->Fill(time);
1208 hc_timeFSLPC[iF][iS][iL][iP][iC]->Fill(time);
1209 h2c_timeF_scint[iF]->Fill(iS, time);
1210 h2c_timeFS[iF][iS]->Fill(iL, time);
1211 h2c_timeFSLP[iF][iS][iL][iP]->Fill(iC, time);
1212 }
1213 } else {
1214 int iF = klmChannel.getSection() - 1;
1215 int iS = klmChannel.getSector() - 1;
1216 int iL = klmChannel.getLayer() - 1;
1217 int iP = klmChannel.getPlane() - 1;
1218 int iC = klmChannel.getStrip() - 1;
1219 double propgationT = event.dist * delayEKLM;
1220 double time = timeHit - propgationT - m_timeShift[channelId];
1221 hc_time_scint_end->Fill(time);
1222 hc_timeF_scint_end[iF]->Fill(time);
1223 hc_timeFS_scint_end[iF][iS]->Fill(time);
1224 hc_timeFSL_end[iF][iS][iL]->Fill(time);
1225 hc_timeFSLP_end[iF][iS][iL][iP]->Fill(time);
1226 hc_timeFSLPC_end[iF][iS][iL][iP][iC]->Fill(time);
1227 h2c_timeF_scint_end[iF]->Fill(iS, time);
1228 h2c_timeFS_end[iF][iS]->Fill(iL, time);
1229 h2c_timeFSLP_end[iF][iS][iL][iP]->Fill(iC, time);
1230 }
1231 }
1232 }
1233
1234 int icChannel_rpc = 0;
1235 int icChannel = 0;
1236 int icChannel_end = 0;
1237 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1238 channelId = klmChannel.getKLMChannelNumber();
1239 if (m_cFlag[channelId] == ChannelCalibrationStatus::c_NotEnoughData)
1240 continue;
1241 int iSub = klmChannel.getSubdetector();
1242
1243 if (iSub == KLMElementNumbers::c_BKLM) {
1244 int iF = klmChannel.getSection();
1245 int iS = klmChannel.getSector() - 1;
1246 int iL = klmChannel.getLayer() - 1;
1247 int iP = klmChannel.getPlane();
1248 int iC = klmChannel.getStrip() - 1;
1249
1250 TFitResultPtr rc = hc_timeFSLPC[iF][iS][iL][iP][iC]->Fit(fcn_gaus, "LESQ");
1251 if (int(rc) != 0)
1252 continue;
1253 if (int(rc) == 0)
1254 m_cFlag[channelId] = ChannelCalibrationStatus::c_SuccessfulCalibration;
1255 m_ctime_channel[channelId] = fcn_gaus->GetParameter(1);
1256 mc_etime_channel[channelId] = fcn_gaus->GetParError(1);
1257 if (iL > 1) {
1258 gre_ctime_channel_rpc->SetPoint(icChannel_rpc, channelId, m_ctime_channel[channelId]);
1259 gre_ctime_channel_rpc->SetPointError(icChannel_rpc, 0., mc_etime_channel[channelId]);
1260 icChannel_rpc++;
1261 } else {
1262 gre_ctime_channel_scint->SetPoint(icChannel, channelId, m_ctime_channel[channelId]);
1263 gre_ctime_channel_scint->SetPointError(icChannel, 0., mc_etime_channel[channelId]);
1264 icChannel++;
1265 }
1266 } else {
1267 int iF = klmChannel.getSection() - 1;
1268 int iS = klmChannel.getSector() - 1;
1269 int iL = klmChannel.getLayer() - 1;
1270 int iP = klmChannel.getPlane() - 1;
1271 int iC = klmChannel.getStrip() - 1;
1272
1273 TFitResultPtr rc = hc_timeFSLPC_end[iF][iS][iL][iP][iC]->Fit(fcn_gaus, "LESQ");
1274 if (int(rc) != 0)
1275 continue;
1276 if (int(rc) == 0)
1277 m_cFlag[channelId] = ChannelCalibrationStatus::c_SuccessfulCalibration;
1278 m_ctime_channel[channelId] = fcn_gaus->GetParameter(1);
1279 mc_etime_channel[channelId] = fcn_gaus->GetParError(1);
1280 gre_ctime_channel_scint_end->SetPoint(icChannel_end, channelId, m_ctime_channel[channelId]);
1281 gre_ctime_channel_scint_end->SetPointError(icChannel_end, 0., mc_etime_channel[channelId]);
1282 icChannel_end++;
1283 }
1284 }
1285
1286 gre_ctime_channel_scint->Fit("fcn_const", "EMQ");
1287 m_ctime_channelAvg_scint = fcn_const->GetParameter(0);
1288 mc_etime_channelAvg_scint = fcn_const->GetParError(0);
1289
1290 gre_ctime_channel_scint_end->Fit("fcn_const", "EMQ");
1291 m_ctime_channelAvg_scint_end = fcn_const->GetParameter(0);
1292 mc_etime_channelAvg_scint_end = fcn_const->GetParError(0);
1293
1294 gre_ctime_channel_rpc->Fit("fcn_const", "EMQ");
1295 m_ctime_channelAvg_rpc = fcn_const->GetParameter(0);
1296 mc_etime_channelAvg_rpc = fcn_const->GetParError(0);
1297
1298 B2INFO("Channel's time distribution fitting done.");
1299 B2DEBUG(20, LogVar("Average calibrated time (RPC)", m_ctime_channelAvg_rpc)
1300 << LogVar("Average calibrated time (BKLM scintillators)", m_ctime_channelAvg_scint)
1301 << LogVar("Average calibrated time (EKLM scintillators)", m_ctime_channelAvg_scint_end));
1302
1303 B2INFO("Calibrated channel's time distribution filling begins.");
1304
1305 m_timeRes.clear();
1306 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1307 channelId = klmChannel.getKLMChannelNumber();
1308 hc_calibrated->Fill(m_cFlag[channelId]);
1309 if (m_ctime_channel.find(channelId) == m_ctime_channel.end())
1310 continue;
1311 double timeRes = m_ctime_channel[channelId];
1312 int iSub = klmChannel.getSubdetector();
1313 if (iSub == KLMElementNumbers::c_BKLM) {
1314 int iL = klmChannel.getLayer() - 1;
1315 if (iL > 1)
1316 m_timeRes[channelId] = timeRes;
1317 else
1318 m_timeRes[channelId] = timeRes;
1319 } else {
1320 m_timeRes[channelId] = timeRes;
1321 }
1322 m_timeResolution->setTimeResolution(channelId, m_timeRes[channelId]);
1323 }
1324
1325 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1326 channelId = klmChannel.getKLMChannelNumber();
1327 if (m_timeRes.find(channelId) != m_timeRes.end())
1328 continue;
1329 m_timeRes[channelId] = esti_timeRes(klmChannel);
1330 m_timeResolution->setTimeResolution(channelId, m_timeRes[channelId]);
1331 B2DEBUG(20, "Calibrated Estimation " << LogVar("Channel", channelId) << LogVar("Estimated value", m_timeRes[channelId]));
1332 }
1333
1334 icChannel_rpc = 0;
1335 icChannel = 0;
1336 icChannel_end = 0;
1337 for (KLMChannelIndex klmChannel = m_klmChannels.begin(); klmChannel != m_klmChannels.end(); ++klmChannel) {
1338 channelId = klmChannel.getKLMChannelNumber();
1339 if (m_timeRes.find(channelId) == m_timeRes.end()) {
1340 B2ERROR("!!! Not All Channels Calibration Constant Set. Error Happened on " << LogVar("Channel", channelId));
1341 continue;
1342 }
1343 int iSub = klmChannel.getSubdetector();
1344 if (iSub == KLMElementNumbers::c_BKLM) {
1345 int iL = klmChannel.getLayer();
1346 if (iL > 2) {
1347 gr_timeRes_channel_rpc->SetPoint(icChannel_rpc, channelId, m_timeRes[channelId]);
1348 icChannel_rpc++;
1349 } else {
1350 gr_timeRes_channel_scint->SetPoint(icChannel, channelId, m_timeRes[channelId]);
1351 icChannel++;
1352 }
1353 } else {
1354 gr_timeRes_channel_scint_end->SetPoint(icChannel_end, channelId, m_timeRes[channelId]);
1355 icChannel_end++;
1356 }
1357 }
1358 saveHist();
1359
1360 delete fcn_const;
1361 m_evts.clear();
1362 m_timeShift.clear();
1363 m_timeRes.clear();
1364 m_cFlag.clear();
1365
1366 saveCalibration(m_timeCableDelay, "KLMTimeCableDelay");
1367 saveCalibration(m_timeConstants, "KLMTimeConstants");
1368 saveCalibration(m_timeResolution, "KLMTimeResolution");
1369 return CalibrationAlgorithm::c_OK;
1370}
1371
1372
1373void KLMTimeAlgorithm::saveHist()
1374{
1375 m_outFile->cd();
1376 B2INFO("Save Histograms into Files.");
1377 TDirectory* dir_monitor = m_outFile->mkdir("monitor_Hists", "", true);
1378 dir_monitor->cd();
1379 h_calibrated->SetDirectory(dir_monitor);
1380 hc_calibrated->SetDirectory(dir_monitor);
1381 h_diff->SetDirectory(dir_monitor);
1382
1383 m_outFile->cd();
1384 TDirectory* dir_effC = m_outFile->mkdir("effC_Hists", "", true);
1385 dir_effC->cd();
1386 m_ProfileRpcPhi->SetDirectory(dir_effC);
1387 m_ProfileRpcZ->SetDirectory(dir_effC);
1388 m_ProfileBKLMScintillatorPhi->SetDirectory(dir_effC);
1389 m_ProfileBKLMScintillatorZ->SetDirectory(dir_effC);
1390 m_ProfileEKLMScintillatorPlane1->SetDirectory(dir_effC);
1391 m_ProfileEKLMScintillatorPlane2->SetDirectory(dir_effC);
1392 m_Profile2RpcPhi->SetDirectory(dir_effC);
1393 m_Profile2RpcZ->SetDirectory(dir_effC);
1394 m_Profile2BKLMScintillatorPhi->SetDirectory(dir_effC);
1395 m_Profile2BKLMScintillatorZ->SetDirectory(dir_effC);
1396 m_Profile2EKLMScintillatorPlane1->SetDirectory(dir_effC);
1397 m_Profile2EKLMScintillatorPlane2->SetDirectory(dir_effC);
1398
1399 m_outFile->cd();
1400 TDirectory* dir_time = m_outFile->mkdir("time", "", true);
1401 dir_time->cd();
1402
1403 h_time_scint->SetDirectory(dir_time);
1404 hc_time_scint->SetDirectory(dir_time);
1405
1406 h_time_scint_end->SetDirectory(dir_time);
1407 hc_time_scint_end->SetDirectory(dir_time);
1408
1409 h_time_rpc->SetDirectory(dir_time);
1410 hc_time_rpc->SetDirectory(dir_time);
1411
1412 gre_time_channel_rpc->Write("gre_time_channel_rpc");
1413 gre_time_channel_scint->Write("gre_time_channel_scint");
1414 gre_time_channel_scint_end->Write("gre_time_channel_scint_end");
1415 gr_timeShift_channel_rpc->Write("gr_timeShift_channel_rpc");
1416 gr_timeShift_channel_scint->Write("gr_timeShift_channel_scint");
1417 gr_timeShift_channel_scint_end->Write("gr_timeShift_channel_scint_end");
1418 gre_ctime_channel_rpc->Write("gre_ctime_channel_rpc");
1419 gre_ctime_channel_scint->Write("gre_ctime_channel_scint");
1420 gre_ctime_channel_scint_end->Write("gre_ctime_channel_scint_end");
1421 gr_timeRes_channel_rpc->Write("gr_timeRes_channel_rpc");
1422 gr_timeRes_channel_scint->Write("gr_timeRes_channel_scint");
1423 gr_timeRes_channel_scint_end->Write("gr_timeRes_channel_scint_end");
1424
1425 B2INFO("Top file setup Done.");
1426
1427 TDirectory* dir_time_F[2];
1428 TDirectory* dir_time_FS[2][8];
1429 TDirectory* dir_time_FSL[2][8][15];
1430 TDirectory* dir_time_FSLP[2][8][15][2];
1431 TDirectory* dir_time_F_end[2];
1432 TDirectory* dir_time_FS_end[2][4];
1433 TDirectory* dir_time_FSL_end[2][4][14];
1434 TDirectory* dir_time_FSLP_end[2][4][14][2];
1435 char dirname[50];
1436 B2INFO("Sub files declare Done.");
1437 for (int iF = 0; iF < 2; ++iF) {
1438 h_timeF_rpc[iF]->SetDirectory(dir_time);
1439 hc_timeF_rpc[iF]->SetDirectory(dir_time);
1440
1441 h2_timeF_rpc[iF]->SetDirectory(dir_time);
1442 h2c_timeF_rpc[iF]->SetDirectory(dir_time);
1443
1444 h_timeF_scint[iF]->SetDirectory(dir_time);
1445 hc_timeF_scint[iF]->SetDirectory(dir_time);
1446
1447 h2_timeF_scint[iF]->SetDirectory(dir_time);
1448 h2c_timeF_scint[iF]->SetDirectory(dir_time);
1449
1450 h_timeF_scint_end[iF]->SetDirectory(dir_time);
1451 hc_timeF_scint_end[iF]->SetDirectory(dir_time);
1452
1453 h2_timeF_scint_end[iF]->SetDirectory(dir_time);
1454 h2c_timeF_scint_end[iF]->SetDirectory(dir_time);
1455
1456 sprintf(dirname, "isForward_%d", iF);
1457 dir_time_F[iF] = dir_time->mkdir(dirname, "", true);
1458 dir_time_F[iF]->cd();
1459
1460 for (int iS = 0; iS < 8; ++iS) {
1461 h_timeFS_rpc[iF][iS]->SetDirectory(dir_time_F[iF]);
1462 hc_timeFS_rpc[iF][iS]->SetDirectory(dir_time_F[iF]);
1463
1464 h_timeFS_scint[iF][iS]->SetDirectory(dir_time_F[iF]);
1465 hc_timeFS_scint[iF][iS]->SetDirectory(dir_time_F[iF]);
1466
1467 h2_timeFS[iF][iS]->SetDirectory(dir_time_F[iF]);
1468 h2c_timeFS[iF][iS]->SetDirectory(dir_time_F[iF]);
1469
1470 sprintf(dirname, "Sector_%d", iS + 1);
1471 dir_time_FS[iF][iS] = dir_time_F[iF]->mkdir(dirname, "", true);
1472 dir_time_FS[iF][iS]->cd();
1473
1474 for (int iL = 0; iL < 15; ++iL) {
1475 h_timeFSL[iF][iS][iL]->SetDirectory(dir_time_FS[iF][iS]);
1476 hc_timeFSL[iF][iS][iL]->SetDirectory(dir_time_FS[iF][iS]);
1477
1478 sprintf(dirname, "Layer_%d", iL + 1);
1479 dir_time_FSL[iF][iS][iL] = dir_time_FS[iF][iS]->mkdir(dirname, "", true);
1480 dir_time_FSL[iF][iS][iL]->cd();
1481 for (int iP = 0; iP < 2; ++iP) {
1482 h_timeFSLP[iF][iS][iL][iP]->SetDirectory(dir_time_FSL[iF][iS][iL]);
1483 hc_timeFSLP[iF][iS][iL][iP]->SetDirectory(dir_time_FSL[iF][iS][iL]);
1484 h2_timeFSLP[iF][iS][iL][iP]->SetDirectory(dir_time_FSL[iF][iS][iL]);
1485 h2c_timeFSLP[iF][iS][iL][iP]->SetDirectory(dir_time_FSL[iF][iS][iL]);
1486
1487 sprintf(dirname, "Plane_%d", iP);
1488 dir_time_FSLP[iF][iS][iL][iP] = dir_time_FSL[iF][iS][iL]->mkdir(dirname, "", true);
1489 dir_time_FSLP[iF][iS][iL][iP]->cd();
1490
1491 int nchannel_max = BKLMElementNumbers::getNStrips(iF, iS + 1, iL + 1, iP);
1492 for (int iC = 0; iC < nchannel_max; ++iC) {
1493 if (iL < 2)
1494 m_HistTimeLengthBKLM[iF][iS][iL][iP][iC]->SetDirectory(dir_time_FSLP[iF][iS][iL][iP]);
1495 h_timeFSLPC[iF][iS][iL][iP][iC]->SetDirectory(dir_time_FSLP[iF][iS][iL][iP]);
1496 hc_timeFSLPC[iF][iS][iL][iP][iC]->SetDirectory(dir_time_FSLP[iF][iS][iL][iP]);
1497 delete h_timeFSLPC_tc[iF][iS][iL][iP][iC];
1498 }
1499 }
1500 }
1501 }
1502
1503 sprintf(dirname, "isForward_%d_end", iF + 1);
1504 dir_time_F_end[iF] = dir_time->mkdir(dirname, "", true);
1505 dir_time_F_end[iF]->cd();
1506 int maxLayer = 12 + 2 * iF;
1507 for (int iS = 0; iS < 4; ++iS) {
1508 h_timeFS_scint_end[iF][iS]->SetDirectory(dir_time_F_end[iF]);
1509 hc_timeFS_scint_end[iF][iS]->SetDirectory(dir_time_F_end[iF]);
1510
1511 h2_timeFS_end[iF][iS]->SetDirectory(dir_time_F_end[iF]);
1512 h2c_timeFS_end[iF][iS]->SetDirectory(dir_time_F_end[iF]);
1513
1514 sprintf(dirname, "Sector_%d_end", iS + 1);
1515 dir_time_FS_end[iF][iS] = dir_time_F_end[iF]->mkdir(dirname, "", true);
1516 dir_time_FS_end[iF][iS]->cd();
1517 for (int iL = 0; iL < maxLayer; ++iL) {
1518 h_timeFSL_end[iF][iS][iL]->SetDirectory(dir_time_FS_end[iF][iS]);
1519 hc_timeFSL_end[iF][iS][iL]->SetDirectory(dir_time_FS_end[iF][iS]);
1520
1521 sprintf(dirname, "Layer_%d_end", iL + 1);
1522 dir_time_FSL_end[iF][iS][iL] = dir_time_FS_end[iF][iS]->mkdir(dirname, "", true);
1523 dir_time_FSL_end[iF][iS][iL]->cd();
1524 for (int iP = 0; iP < 2; ++iP) {
1525 h_timeFSLP_end[iF][iS][iL][iP]->SetDirectory(dir_time_FSL_end[iF][iS][iL]);
1526 hc_timeFSLP_end[iF][iS][iL][iP]->SetDirectory(dir_time_FSL_end[iF][iS][iL]);
1527 h2_timeFSLP_end[iF][iS][iL][iP]->SetDirectory(dir_time_FSL_end[iF][iS][iL]);
1528 h2c_timeFSLP_end[iF][iS][iL][iP]->SetDirectory(dir_time_FSL_end[iF][iS][iL]);
1529
1530 sprintf(dirname, "plane_%d_end", iP);
1531 dir_time_FSLP_end[iF][iS][iL][iP] = dir_time_FSL_end[iF][iS][iL]->mkdir(dirname, "", true);
1532 dir_time_FSLP_end[iF][iS][iL][iP]->cd();
1533
1534 for (int iC = 0; iC < 75; ++iC) {
1535 h_timeFSLPC_end[iF][iS][iL][iP][iC]->SetDirectory(dir_time_FSLP_end[iF][iS][iL][iP]);
1536 hc_timeFSLPC_end[iF][iS][iL][iP][iC]->SetDirectory(dir_time_FSLP_end[iF][iS][iL][iP]);
1537 m_HistTimeLengthEKLM[iF][iS][iL][iP][iC]->SetDirectory(dir_time_FSLP_end[iF][iS][iL][iP]);
1538 delete h_timeFSLPC_tc_end[iF][iS][iL][iP][iC];
1539 }
1540 }
1541 }
1542 }
1543 }
1544 m_outFile->cd();
1545 m_outFile->Write();
1546 m_outFile->Close();
1547 B2INFO("File Write and Close. Done.");
1548}
1549
1550double KLMTimeAlgorithm::esti_timeShift(const KLMChannelIndex& klmChannel)
1551{
1552 double tS = 0.0;
1553 int iSub = klmChannel.getSubdetector();
1554 int iF = klmChannel.getSection();
1555 int iS = klmChannel.getSector();
1556 int iL = klmChannel.getLayer();
1557 int iP = klmChannel.getPlane();
1558 int iC = klmChannel.getStrip();
1559 int totNStrips = EKLMElementNumbers::getMaximalStripNumber();
1560 if (iSub == KLMElementNumbers::c_BKLM)
1561 totNStrips = BKLMElementNumbers::getNStrips(iF, iS, iL, iP);
1562 if (iC == 1) {
1563 KLMChannelIndex kCIndex_upper(iSub, iF, iS, iL, iP, iC + 1);
1564 tS = tS_upperStrip(kCIndex_upper).second;
1565 } else if (iC == totNStrips) {
1566 KLMChannelIndex kCIndex_lower(iSub, iF, iS, iL, iP, iC - 1);
1567 tS = tS_lowerStrip(kCIndex_lower).second;
1568 } else {
1569 KLMChannelIndex kCIndex_upper(iSub, iF, iS, iL, iP, iC + 1);
1570 KLMChannelIndex kCIndex_lower(iSub, iF, iS, iL, iP, iC - 1);
1571 std::pair<int, double> tS_upper = tS_upperStrip(kCIndex_upper);
1572 std::pair<int, double> tS_lower = tS_lowerStrip(kCIndex_lower);
1573 unsigned int td_upper = tS_upper.first - iC;
1574 unsigned int td_lower = iC - tS_lower.first;
1575 unsigned int td = tS_upper.first - tS_lower.first;
1576 tS = (double(td_upper) * tS_lower.second + double(td_lower) * tS_upper.second) / double(td);
1577 }
1578 return tS;
1579}
1580
1581std::pair<int, double> KLMTimeAlgorithm::tS_upperStrip(const KLMChannelIndex& klmChannel)
1582{
1583 std::pair<int, double> tS;
1584 int cId = klmChannel.getKLMChannelNumber();
1585 int iSub = klmChannel.getSubdetector();
1586 int iF = klmChannel.getSection();
1587 int iS = klmChannel.getSector();
1588 int iL = klmChannel.getLayer();
1589 int iP = klmChannel.getPlane();
1590 int iC = klmChannel.getStrip();
1591 int totNStrips = EKLMElementNumbers::getMaximalStripNumber();
1592 if (iSub == KLMElementNumbers::c_BKLM)
1593 totNStrips = BKLMElementNumbers::getNStrips(iF, iS, iL, iP);
1594 if (m_timeShift.find(cId) != m_timeShift.end()) {
1595 tS.first = iC;
1596 tS.second = m_timeShift[cId];
1597 } else if (iC == totNStrips) {
1598 tS.first = iC;
1599 tS.second = 0.0;
1600 } else {
1601 KLMChannelIndex kCIndex(iSub, iF, iS, iL, iP, iC + 1);
1602 tS = tS_upperStrip(kCIndex);
1603 }
1604 return tS;
1605}
1606
1607std::pair<int, double> KLMTimeAlgorithm::tS_lowerStrip(const KLMChannelIndex& klmChannel)
1608{
1609 std::pair<int, double> tS;
1610 int cId = klmChannel.getKLMChannelNumber();
1611 int iSub = klmChannel.getSubdetector();
1612 int iF = klmChannel.getSection();
1613 int iS = klmChannel.getSector();
1614 int iL = klmChannel.getLayer();
1615 int iP = klmChannel.getPlane();
1616 int iC = klmChannel.getStrip();
1617 if (m_timeShift.find(cId) != m_timeShift.end()) {
1618 tS.first = iC;
1619 tS.second = m_timeShift[cId];
1620 } else if (iC == 1) {
1621 tS.first = iC;
1622 tS.second = 0.0;
1623 } else {
1624 KLMChannelIndex kCIndex(iSub, iF, iS, iL, iP, iC - 1);
1625 tS = tS_lowerStrip(kCIndex);
1626 }
1627 return tS;
1628}
1629
1630double KLMTimeAlgorithm::esti_timeRes(const KLMChannelIndex& klmChannel)
1631{
1632 double tR = 0.0;
1633 int iSub = klmChannel.getSubdetector();
1634 int iF = klmChannel.getSection();
1635 int iS = klmChannel.getSector();
1636 int iL = klmChannel.getLayer();
1637 int iP = klmChannel.getPlane();
1638 int iC = klmChannel.getStrip();
1639 int totNStrips = EKLMElementNumbers::getMaximalStripNumber();
1640 if (iSub == KLMElementNumbers::c_BKLM)
1641 totNStrips = BKLMElementNumbers::getNStrips(iF, iS, iL, iP);
1642 if (iC == 1) {
1643 KLMChannelIndex kCIndex_upper(iSub, iF, iS, iL, iP, iC + 1);
1644 tR = tR_upperStrip(kCIndex_upper).second;
1645 } else if (iC == totNStrips) {
1646 KLMChannelIndex kCIndex_lower(iSub, iF, iS, iL, iP, iC - 1);
1647 tR = tR_lowerStrip(kCIndex_lower).second;
1648 } else {
1649 KLMChannelIndex kCIndex_upper(iSub, iF, iS, iL, iP, iC + 1);
1650 KLMChannelIndex kCIndex_lower(iSub, iF, iS, iL, iP, iC - 1);
1651 std::pair<int, double> tR_upper = tR_upperStrip(kCIndex_upper);
1652 std::pair<int, double> tR_lower = tR_lowerStrip(kCIndex_lower);
1653 unsigned int tr_upper = tR_upper.first - iC;
1654 unsigned int tr_lower = iC - tR_lower.first;
1655 unsigned int tr = tR_upper.first - tR_lower.first;
1656 tR = (double(tr_upper) * tR_lower.second + double(tr_lower) * tR_upper.second) / double(tr);
1657 }
1658 return tR;
1659}
1660
1661std::pair<int, double> KLMTimeAlgorithm::tR_upperStrip(const KLMChannelIndex& klmChannel)
1662{
1663 std::pair<int, double> tR;
1664 int cId = klmChannel.getKLMChannelNumber();
1665 int iSub = klmChannel.getSubdetector();
1666 int iF = klmChannel.getSection();
1667 int iS = klmChannel.getSector();
1668 int iL = klmChannel.getLayer();
1669 int iP = klmChannel.getPlane();
1670 int iC = klmChannel.getStrip();
1671 int totNStrips = EKLMElementNumbers::getMaximalStripNumber();
1672 if (iSub == KLMElementNumbers::c_BKLM)
1673 totNStrips = BKLMElementNumbers::getNStrips(iF, iS, iL, iP);
1674 if (m_timeRes.find(cId) != m_timeRes.end()) {
1675 tR.first = iC;
1676 tR.second = m_timeRes[cId];
1677 } else if (iC == totNStrips) {
1678 tR.first = iC;
1679 tR.second = 0.0;
1680 } else {
1681 KLMChannelIndex kCIndex(iSub, iF, iS, iL, iP, iC + 1);
1682 tR = tR_upperStrip(kCIndex);
1683 }
1684 return tR;
1685}
1686
1687std::pair<int, double> KLMTimeAlgorithm::tR_lowerStrip(const KLMChannelIndex& klmChannel)
1688{
1689 std::pair<int, double> tR;
1690 int cId = klmChannel.getKLMChannelNumber();
1691 int iSub = klmChannel.getSubdetector();
1692 int iF = klmChannel.getSection();
1693 int iS = klmChannel.getSector();
1694 int iL = klmChannel.getLayer();
1695 int iP = klmChannel.getPlane();
1696 int iC = klmChannel.getStrip();
1697 if (m_timeRes.find(cId) != m_timeRes.end()) {
1698 tR.first = iC;
1699 tR.second = m_timeRes[cId];
1700 } else if (iC == 1) {
1701 tR.first = iC;
1702 tR.second = 0.0;
1703 } else {
1704 KLMChannelIndex kCIndex(iSub, iF, iS, iL, iP, iC - 1);
1705 tR = tR_lowerStrip(kCIndex);
1706 }
1707 return tR;
1708}
1709
Belle2::BKLMElementNumbers::c_FirstRPCLayer
@ c_FirstRPCLayer
First RPC layer.
Definition: BKLMElementNumbers.h:64
Belle2::BKLMElementNumbers::getNStrips
static int getNStrips(int section, int sector, int layer, int plane)
Get number of strips.
Definition: BKLMElementNumbers.cc:107
Belle2::CalibrationAlgorithm
Base class for calibration algorithms.
Definition: CalibrationAlgorithm.h:37
Belle2::CalibrationAlgorithm::saveCalibration
void saveCalibration(TClonesArray *data, const std::string &name)
Store DBArray payload with given name with default IOV.
Definition: CalibrationAlgorithm.cc:297
Belle2::CalibrationAlgorithm::getRunList
const std::vector< Calibration::ExpRun > & getRunList() const
Get the list of runs for which calibration is called.
Definition: CalibrationAlgorithm.h:266
Belle2::CalibrationAlgorithm::EResult
EResult
The result of calibration.
Definition: CalibrationAlgorithm.h:40
Belle2::CalibrationAlgorithm::c_OK
@ c_OK
Finished successfully =0 in Python.
Definition: CalibrationAlgorithm.h:41
Belle2::CalibrationAlgorithm::c_NotEnoughData
@ c_NotEnoughData
Needs more data =2 in Python.
Definition: CalibrationAlgorithm.h:43
Belle2::DBObjPtr
Class for accessing objects in the database.
Definition: DBObjPtr.h:21
Belle2::DBStore
Singleton class to cache database objects.
Definition: DBStore.h:31
Belle2::DataStore::Instance
static DataStore & Instance()
Instance of singleton Store.
Definition: DataStore.cc:54
Belle2::DataStore::setInitializeActive
void setInitializeActive(bool active)
Setter for m_initializeActive.
Definition: DataStore.cc:94
Belle2::EKLMElementNumbers::getMaximalStripNumber
static constexpr int getMaximalStripNumber()
Get maximal strip number.
Definition: EKLMElementNumbers.h:319
Belle2::EKLM::GeometryData::Instance
static const GeometryData & Instance(enum DataSource dataSource=c_Database, const GearDir *gearDir=nullptr)
Instantiation.
Definition: GeometryData.cc:33
Belle2::EKLM::GeometryData::getStripLength
double getStripLength(int strip) const
Get strip length.
Definition: GeometryData.h:71
Belle2::EKLM::GeometryData::getMaximalStripLength
double getMaximalStripLength() const
Get maximal strip length.
Definition: GeometryData.h:79
Belle2::KLMChannelIndex
KLM channel index.
Definition: KLMChannelIndex.h:23
Belle2::KLMChannelIndex::getSubdetector
int getSubdetector() const
Get subdetector.
Definition: KLMChannelIndex.h:106
Belle2::KLMChannelIndex::getLayer
int getLayer() const
Get layer.
Definition: KLMChannelIndex.h:130
Belle2::KLMChannelIndex::begin
KLMChannelIndex begin()
First channel.
Definition: KLMChannelIndex.h:199
Belle2::KLMChannelIndex::getSection
int getSection() const
Get section.
Definition: KLMChannelIndex.h:114
Belle2::KLMChannelIndex::getPlane
int getPlane() const
Get plane.
Definition: KLMChannelIndex.h:138
Belle2::KLMChannelIndex::getStrip
int getStrip() const
Get strip.
Definition: KLMChannelIndex.h:146
Belle2::KLMChannelIndex::getSector
int getSector() const
Get sector.
Definition: KLMChannelIndex.h:122
Belle2::KLMChannelIndex::getKLMChannelNumber
KLMChannelNumber getKLMChannelNumber() const
Get KLM channel number.
Definition: KLMChannelIndex.cc:143
Belle2::KLMChannelIndex::end
KLMChannelIndex & end()
Last channel.
Definition: KLMChannelIndex.h:207
Belle2::KLMElementNumbers::Instance
static const KLMElementNumbers & Instance()
Instantiation.
Definition: KLMElementNumbers.cc:29
Belle2::KLMElementNumbers::channelNumberToElementNumbers
void channelNumberToElementNumbers(KLMChannelNumber channel, int *subdetector, int *section, int *sector, int *layer, int *plane, int *strip) const
Get element numbers by channel number.
Definition: KLMElementNumbers.cc:104
Belle2::KLMTimeAlgorithm::m_Profile2EKLMScintillatorPlane2
TProfile * m_Profile2EKLMScintillatorPlane2
For EKLM scintillator plane2.
Definition: KLMTimeAlgorithm.h:518
Belle2::KLMTimeAlgorithm::mc_etime_channelAvg_rpc
double mc_etime_channelAvg_rpc
Calibrated central value error of the global time distribution (BKLM RPC part).
Definition: KLMTimeAlgorithm.h:367
Belle2::KLMTimeAlgorithm::m_HistTimeLengthEKLM
TH2F * m_HistTimeLengthEKLM[2][4][14][2][75]
Two-dimensional distributions of time versus propagation length.
Definition: KLMTimeAlgorithm.h:728
Belle2::KLMTimeAlgorithm::h_timeFSLPC_tc
TH1F * h_timeFSLPC_tc[2][8][15][2][54]
BKLM part, used for effective light speed estimation.
Definition: KLMTimeAlgorithm.h:713
Belle2::KLMTimeAlgorithm::h2c_timeF_scint_end
TH2F * h2c_timeF_scint_end[2]
EKLM part.
Definition: KLMTimeAlgorithm.h:610
Belle2::KLMTimeAlgorithm::h_timeFSLPC_tc_end
TH1F * h_timeFSLPC_tc_end[2][4][14][2][75]
EKLM part, used for effective light speed estimation.
Definition: KLMTimeAlgorithm.h:722
Belle2::KLMTimeAlgorithm::m_timeResolution
KLMTimeResolution * m_timeResolution
DBObject of time resolution.
Definition: KLMTimeAlgorithm.h:416
Belle2::KLMTimeAlgorithm::h_time_scint_tc_end
TH1F * h_time_scint_tc_end
EKLM part.
Definition: KLMTimeAlgorithm.h:532
Belle2::KLMTimeAlgorithm::createHistograms
void createHistograms()
Create histograms.
Definition: KLMTimeAlgorithm.cc:174
Belle2::KLMTimeAlgorithm::gre_time_channel_scint
TGraphErrors * gre_time_channel_scint
BKLM Scintillator.
Definition: KLMTimeAlgorithm.h:442
Belle2::KLMTimeAlgorithm::h_timeFSL
TH1F * h_timeFSL[2][8][15]
BKLM part.
Definition: KLMTimeAlgorithm.h:659
Belle2::KLMTimeAlgorithm::hc_timeFSLPC_end
TH1F * hc_timeFSLPC_end[2][4][14][2][75]
EKLM part.
Definition: KLMTimeAlgorithm.h:736
Belle2::KLMTimeAlgorithm::hc_timeFSL_end
TH1F * hc_timeFSL_end[2][4][14]
EKLM part.
Definition: KLMTimeAlgorithm.h:670
Belle2::KLMTimeAlgorithm::h_timeFSLP_end
TH1F * h_timeFSLP_end[2][4][14][2]
EKLM part.
Definition: KLMTimeAlgorithm.h:678
Belle2::KLMTimeAlgorithm::gr_timeRes_channel_rpc
TGraph * gr_timeRes_channel_rpc
BKLM RPC.
Definition: KLMTimeAlgorithm.h:472
Belle2::KLMTimeAlgorithm::hc_timeFSLP_end
TH1F * hc_timeFSLP_end[2][4][14][2]
EKLM part.
Definition: KLMTimeAlgorithm.h:686
Belle2::KLMTimeAlgorithm::h_timeFSLP
TH1F * h_timeFSLP[2][8][15][2]
BKLM part.
Definition: KLMTimeAlgorithm.h:675
Belle2::KLMTimeAlgorithm::hc_timeF_scint_end
TH1F * hc_timeF_scint_end[2]
EKLM part.
Definition: KLMTimeAlgorithm.h:582
Belle2::KLMTimeAlgorithm::m_timeShift
std::map< KLMChannelNumber, double > m_timeShift
Shift values of each channel.
Definition: KLMTimeAlgorithm.h:288
Belle2::KLMTimeAlgorithm::h_time_scint
TH1F * h_time_scint
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:540
Belle2::KLMTimeAlgorithm::m_time_channelAvg_scint
double m_time_channelAvg_scint
Central value of the global time distribution (BKLM scintillator part).
Definition: KLMTimeAlgorithm.h:349
Belle2::KLMTimeAlgorithm::hc_timeFS_scint_end
TH1F * hc_timeFS_scint_end[2][4]
EKLM part.
Definition: KLMTimeAlgorithm.h:632
Belle2::KLMTimeAlgorithm::esti_timeRes
double esti_timeRes(const KLMChannelIndex &klmChannel)
Estimate value of calibration constant for calibrated channels.
Definition: KLMTimeAlgorithm.cc:1630
Belle2::KLMTimeAlgorithm::m_UpperTimeBoundaryCalibratedRPC
double m_UpperTimeBoundaryCalibratedRPC
Upper time boundary for RPC (calibrated data).
Definition: KLMTimeAlgorithm.h:327
Belle2::KLMTimeAlgorithm::m_ctime_channelAvg_rpc
double m_ctime_channelAvg_rpc
Calibrated central value of the global time distribution (BKLM RPC part).
Definition: KLMTimeAlgorithm.h:364
Belle2::KLMTimeAlgorithm::m_timeConstants
KLMTimeConstants * m_timeConstants
DBObject of time cost on some parts of the detector.
Definition: KLMTimeAlgorithm.h:407
Belle2::KLMTimeAlgorithm::setupDatabase
void setupDatabase()
Setup the database.
Definition: KLMTimeAlgorithm.cc:104
Belle2::KLMTimeAlgorithm::hc_timeFS_scint
TH1F * hc_timeFS_scint[2][8]
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:629
Belle2::KLMTimeAlgorithm::m_time_channel
std::map< KLMChannelNumber, double > m_time_channel
Time distribution central value of each channel.
Definition: KLMTimeAlgorithm.h:294
Belle2::KLMTimeAlgorithm::m_ctime_channelAvg_scint_end
double m_ctime_channelAvg_scint_end
Calibrated central value of the global time distribution (EKLM scintillator part).
Definition: KLMTimeAlgorithm.h:379
Belle2::KLMTimeAlgorithm::readCalibrationData
CalibrationAlgorithm::EResult readCalibrationData()
Read calibration data.
Definition: KLMTimeAlgorithm.cc:141
Belle2::KLMTimeAlgorithm::gr_timeShift_channel_scint_end
TGraph * gr_timeShift_channel_scint_end
EKLM.
Definition: KLMTimeAlgorithm.h:467
Belle2::KLMTimeAlgorithm::gr_timeRes_channel_scint
TGraph * gr_timeRes_channel_scint
BKLM scintillator.
Definition: KLMTimeAlgorithm.h:475
Belle2::KLMTimeAlgorithm::hc_timeF_scint
TH1F * hc_timeF_scint[2]
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:579
Belle2::KLMTimeAlgorithm::h_timeFS_scint
TH1F * h_timeFS_scint[2][8]
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:618
Belle2::KLMTimeAlgorithm::m_UpperTimeBoundaryScintillatorsBKLM
double m_UpperTimeBoundaryScintillatorsBKLM
Upper time boundary for BKLM scintillators.
Definition: KLMTimeAlgorithm.h:315
Belle2::KLMTimeAlgorithm::m_ElementNumbers
const KLMElementNumbers * m_ElementNumbers
Element numbers.
Definition: KLMTimeAlgorithm.h:392
Belle2::KLMTimeAlgorithm::tR_upperStrip
std::pair< int, double > tR_upperStrip(const KLMChannelIndex &klmChannel)
Tracing available channels with increasing strip number.
Definition: KLMTimeAlgorithm.cc:1661
Belle2::KLMTimeAlgorithm::hc_timeF_rpc
TH1F * hc_timeF_rpc[2]
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:576
Belle2::KLMTimeAlgorithm::h2c_timeFS_end
TH2F * h2c_timeFS_end[2][4]
EKLM part.
Definition: KLMTimeAlgorithm.h:654
Belle2::KLMTimeAlgorithm::h_timeFSLPC_end
TH1F * h_timeFSLPC_end[2][4][14][2][75]
EKLM part.
Definition: KLMTimeAlgorithm.h:725
Belle2::KLMTimeAlgorithm::m_EKLMGeometry
const EKLM::GeometryData * m_EKLMGeometry
EKLM geometry data.
Definition: KLMTimeAlgorithm.h:398
Belle2::KLMTimeAlgorithm::gre_ctime_channel_scint_end
TGraphErrors * gre_ctime_channel_scint_end
EKLM.
Definition: KLMTimeAlgorithm.h:456
Belle2::KLMTimeAlgorithm::m_Profile2BKLMScintillatorPhi
TProfile * m_Profile2BKLMScintillatorPhi
For BKLM scintillator phi plane.
Definition: KLMTimeAlgorithm.h:509
Belle2::KLMTimeAlgorithm::hc_time_scint_end
TH1F * hc_time_scint_end
EKLM part.
Definition: KLMTimeAlgorithm.h:554
Belle2::KLMTimeAlgorithm::gre_time_channel_scint_end
TGraphErrors * gre_time_channel_scint_end
EKLM.
Definition: KLMTimeAlgorithm.h:445
Belle2::KLMTimeAlgorithm::h2_timeFSLP
TH2F * h2_timeFSLP[2][8][15][2]
BKLM part.
Definition: KLMTimeAlgorithm.h:694
Belle2::KLMTimeAlgorithm::m_UpperTimeBoundaryCalibratedScintillatorsEKLM
double m_UpperTimeBoundaryCalibratedScintillatorsEKLM
Upper time boundary for BKLM scintillators (calibrated data).
Definition: KLMTimeAlgorithm.h:339
Belle2::KLMTimeAlgorithm::gr_timeShift_channel_scint
TGraph * gr_timeShift_channel_scint
BKLM scintillator.
Definition: KLMTimeAlgorithm.h:464
Belle2::KLMTimeAlgorithm::m_time_channelAvg_scint_end
double m_time_channelAvg_scint_end
Central value of the global time distribution (EKLM scintillator part).
Definition: KLMTimeAlgorithm.h:357
Belle2::KLMTimeAlgorithm::m_Profile2EKLMScintillatorPlane1
TProfile * m_Profile2EKLMScintillatorPlane1
For EKLM scintillator plane1.
Definition: KLMTimeAlgorithm.h:515
Belle2::KLMTimeAlgorithm::hc_timeFS_rpc
TH1F * hc_timeFS_rpc[2][8]
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:626
Belle2::KLMTimeAlgorithm::m_UpperTimeBoundaryScintillatorsEKLM
double m_UpperTimeBoundaryScintillatorsEKLM
Upper time boundary for BKLM scintillators.
Definition: KLMTimeAlgorithm.h:321
Belle2::KLMTimeAlgorithm::fillTimeDistanceProfiles
void fillTimeDistanceProfiles(TProfile *profileRpcPhi, TProfile *profileRpcZ, TProfile *profileBKLMScintillatorPhi, TProfile *profileBKLMScintillatorZ, TProfile *profileEKLMScintillatorPlane1, TProfile *profileEKLMScintillatorPlane2, bool fill2dHistograms)
Fill profiles of time versus distance.
Definition: KLMTimeAlgorithm.cc:595
Belle2::KLMTimeAlgorithm::m_outFile
TFile * m_outFile
Output file.
Definition: KLMTimeAlgorithm.h:753
Belle2::KLMTimeAlgorithm::esti_timeShift
double esti_timeShift(const KLMChannelIndex &klmChannel)
Estimate value of calibration constant for uncalibrated channels.
Definition: KLMTimeAlgorithm.cc:1550
Belle2::KLMTimeAlgorithm::tS_upperStrip
std::pair< int, double > tS_upperStrip(const KLMChannelIndex &klmChannel)
Tracing available channels with increasing strip number.
Definition: KLMTimeAlgorithm.cc:1581
Belle2::KLMTimeAlgorithm::m_LowerTimeBoundaryCalibratedScintillatorsEKLM
double m_LowerTimeBoundaryCalibratedScintillatorsEKLM
Lower time boundary for EKLM scintillators (calibrated data).
Definition: KLMTimeAlgorithm.h:336
Belle2::KLMTimeAlgorithm::hc_timeFSLP
TH1F * hc_timeFSLP[2][8][15][2]
BKLM part.
Definition: KLMTimeAlgorithm.h:683
Belle2::KLMTimeAlgorithm::gre_ctime_channel_rpc
TGraphErrors * gre_ctime_channel_rpc
BKLM RPC.
Definition: KLMTimeAlgorithm.h:450
Belle2::KLMTimeAlgorithm::saveHist
void saveHist()
Save histograms to file.
Definition: KLMTimeAlgorithm.cc:1373
Belle2::KLMTimeAlgorithm::m_BKLMGeometry
const bklm::GeometryPar * m_BKLMGeometry
BKLM geometry data.
Definition: KLMTimeAlgorithm.h:395
Belle2::KLMTimeAlgorithm::h2c_timeFSLP
TH2F * h2c_timeFSLP[2][8][15][2]
BKLM part.
Definition: KLMTimeAlgorithm.h:705
Belle2::KLMTimeAlgorithm::m_ctime_channelAvg_scint
double m_ctime_channelAvg_scint
Calibrated central value of the global time distribution (BKLM scintillator part).
Definition: KLMTimeAlgorithm.h:371
Belle2::KLMTimeAlgorithm::fcn_const
TF1 * fcn_const
Const function.
Definition: KLMTimeAlgorithm.h:744
Belle2::KLMTimeAlgorithm::m_UpperTimeBoundaryCalibratedScintillatorsBKLM
double m_UpperTimeBoundaryCalibratedScintillatorsBKLM
Upper time boundary for BKLM scintillators (calibrated data).
Definition: KLMTimeAlgorithm.h:333
Belle2::KLMTimeAlgorithm::m_Profile2RpcZ
TProfile * m_Profile2RpcZ
For BKLM RPC z plane.
Definition: KLMTimeAlgorithm.h:506
Belle2::KLMTimeAlgorithm::h_timeFSLPC
TH1F * h_timeFSLPC[2][8][15][2][54]
BKLM part.
Definition: KLMTimeAlgorithm.h:716
Belle2::KLMTimeAlgorithm::m_HistTimeLengthBKLM
TH2F * m_HistTimeLengthBKLM[2][8][15][2][54]
Two-dimensional distributions of time versus propagation length.
Definition: KLMTimeAlgorithm.h:719
Belle2::KLMTimeAlgorithm::h2_timeFSLP_end
TH2F * h2_timeFSLP_end[2][4][14][2]
EKLM part.
Definition: KLMTimeAlgorithm.h:697
Belle2::KLMTimeAlgorithm::hc_calibrated
TH1I * hc_calibrated
Calibration statistics for each channel.
Definition: KLMTimeAlgorithm.h:431
Belle2::KLMTimeAlgorithm::timeDistance2dFit
void timeDistance2dFit(const std::vector< std::pair< KLMChannelNumber, unsigned int > > &channels, double &delay, double &delayError)
Two-dimensional fit for individual channels.
Definition: KLMTimeAlgorithm.cc:655
Belle2::KLMTimeAlgorithm::gr_timeRes_channel_scint_end
TGraph * gr_timeRes_channel_scint_end
EKLM.
Definition: KLMTimeAlgorithm.h:478
Belle2::KLMTimeAlgorithm::h_calibrated
TH1I * h_calibrated
Calibration statistics for each channel.
Definition: KLMTimeAlgorithm.h:428
Belle2::KLMTimeAlgorithm::m_ProfileBKLMScintillatorZ
TProfile * m_ProfileBKLMScintillatorZ
For BKLM scintillator z plane.
Definition: KLMTimeAlgorithm.h:492
Belle2::KLMTimeAlgorithm::m_LowerTimeBoundaryScintillatorsBKLM
double m_LowerTimeBoundaryScintillatorsBKLM
Lower time boundary for BKLM scintillators.
Definition: KLMTimeAlgorithm.h:312
Belle2::KLMTimeAlgorithm::h_time_rpc_tc
TH1F * h_time_rpc_tc
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:526
Belle2::KLMTimeAlgorithm::h_time_scint_end
TH1F * h_time_scint_end
EKLM part.
Definition: KLMTimeAlgorithm.h:543
Belle2::KLMTimeAlgorithm::h2c_timeF_scint
TH2F * h2c_timeF_scint[2]
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:607
Belle2::KLMTimeAlgorithm::fcn_pol1
TF1 * fcn_pol1
Pol1 function.
Definition: KLMTimeAlgorithm.h:741
Belle2::KLMTimeAlgorithm::m_etime_channelAvg_scint_end
double m_etime_channelAvg_scint_end
Central value error of the global time distribution (EKLM scintillator part).
Definition: KLMTimeAlgorithm.h:361
Belle2::KLMTimeAlgorithm::hc_time_rpc
TH1F * hc_time_rpc
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:548
Belle2::KLMTimeAlgorithm::mc_etime_channelAvg_scint
double mc_etime_channelAvg_scint
Calibrated central value error of the global time distribution (BKLM scintillator part).
Definition: KLMTimeAlgorithm.h:375
Belle2::KLMTimeAlgorithm::h2c_timeFSLP_end
TH2F * h2c_timeFSLP_end[2][4][14][2]
EKLM part.
Definition: KLMTimeAlgorithm.h:708
Belle2::KLMTimeAlgorithm::mc_etime_channelAvg_scint_end
double mc_etime_channelAvg_scint_end
Calibrated central value error of the global time distribution (EKLM scintillator part).
Definition: KLMTimeAlgorithm.h:383
Belle2::KLMTimeAlgorithm::h2_timeF_scint_end
TH2F * h2_timeF_scint_end[2]
EKLM part.
Definition: KLMTimeAlgorithm.h:596
Belle2::KLMTimeAlgorithm::fcn_gaus
TF1 * fcn_gaus
Gaussian function.
Definition: KLMTimeAlgorithm.h:747
Belle2::KLMTimeAlgorithm::m_LowerTimeBoundaryCalibratedScintillatorsBKLM
double m_LowerTimeBoundaryCalibratedScintillatorsBKLM
Lower time boundary for BKLM scintillators (calibrated data).
Definition: KLMTimeAlgorithm.h:330
Belle2::KLMTimeAlgorithm::h_timeF_rpc
TH1F * h_timeF_rpc[2]
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:562
Belle2::KLMTimeAlgorithm::m_ProfileBKLMScintillatorPhi
TProfile * m_ProfileBKLMScintillatorPhi
For BKLM scintillator phi plane.
Definition: KLMTimeAlgorithm.h:489
Belle2::KLMTimeAlgorithm::hc_time_scint
TH1F * hc_time_scint
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:551
Belle2::KLMTimeAlgorithm::h2_timeFS
TH2F * h2_timeFS[2][8]
BKLM part.
Definition: KLMTimeAlgorithm.h:640
Belle2::KLMTimeAlgorithm::m_etime_channelAvg_scint
double m_etime_channelAvg_scint
Central value error of the global time distribution (BKLM scintillator part).
Definition: KLMTimeAlgorithm.h:353
Belle2::KLMTimeAlgorithm::h_timeF_scint_end
TH1F * h_timeF_scint_end[2]
EKLM part.
Definition: KLMTimeAlgorithm.h:568
Belle2::KLMTimeAlgorithm::m_ProfileRpcPhi
TProfile * m_ProfileRpcPhi
For BKLM RPC phi plane.
Definition: KLMTimeAlgorithm.h:483
Belle2::KLMTimeAlgorithm::gre_ctime_channel_scint
TGraphErrors * gre_ctime_channel_scint
BKLM Scintillator.
Definition: KLMTimeAlgorithm.h:453
Belle2::KLMTimeAlgorithm::m_ProfileEKLMScintillatorPlane2
TProfile * m_ProfileEKLMScintillatorPlane2
For EKLM scintillator plane2.
Definition: KLMTimeAlgorithm.h:498
Belle2::KLMTimeAlgorithm::h_time_rpc
TH1F * h_time_rpc
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:537
Belle2::KLMTimeAlgorithm::h_timeFSL_end
TH1F * h_timeFSL_end[2][4][14]
EKLM part.
Definition: KLMTimeAlgorithm.h:662
Belle2::KLMTimeAlgorithm::m_Profile2RpcPhi
TProfile * m_Profile2RpcPhi
For BKLM RPC phi plane.
Definition: KLMTimeAlgorithm.h:503
Belle2::KLMTimeAlgorithm::h_timeF_scint
TH1F * h_timeF_scint[2]
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:565
Belle2::KLMTimeAlgorithm::hc_timeFSL
TH1F * hc_timeFSL[2][8][15]
BKLM part.
Definition: KLMTimeAlgorithm.h:667
Belle2::KLMTimeAlgorithm::m_Profile2BKLMScintillatorZ
TProfile * m_Profile2BKLMScintillatorZ
For BKLM scintillator z plane.
Definition: KLMTimeAlgorithm.h:512
Belle2::KLMTimeAlgorithm::h_timeFS_rpc
TH1F * h_timeFS_rpc[2][8]
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:615
Belle2::KLMTimeAlgorithm::m_klmChannels
KLMChannelIndex m_klmChannels
KLM ChannelIndex object.
Definition: KLMTimeAlgorithm.h:401
Belle2::KLMTimeAlgorithm::gr_timeShift_channel_rpc
TGraph * gr_timeShift_channel_rpc
BKLM RPC.
Definition: KLMTimeAlgorithm.h:461
Belle2::KLMTimeAlgorithm::m_timeRes
std::map< KLMChannelNumber, double > m_timeRes
Resolution values of each channel.
Definition: KLMTimeAlgorithm.h:291
Belle2::KLMTimeAlgorithm::h_diff
TH1F * h_diff
Distance between global and local position.
Definition: KLMTimeAlgorithm.h:434
Belle2::KLMTimeAlgorithm::h2_timeF_scint
TH2F * h2_timeF_scint[2]
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:593
Belle2::KLMTimeAlgorithm::h_time_scint_tc
TH1F * h_time_scint_tc
BKLM scintillator part.
Definition: KLMTimeAlgorithm.h:529
Belle2::KLMTimeAlgorithm::m_LowerTimeBoundaryRPC
double m_LowerTimeBoundaryRPC
Lower time boundary for RPC.
Definition: KLMTimeAlgorithm.h:306
Belle2::KLMTimeAlgorithm::calibrate
virtual EResult calibrate() override
Run algorithm on data.
Definition: KLMTimeAlgorithm.cc:738
Belle2::KLMTimeAlgorithm::m_ctime_channel
std::map< KLMChannelNumber, double > m_ctime_channel
Calibrated time distribution central value of each channel.
Definition: KLMTimeAlgorithm.h:300
Belle2::KLMTimeAlgorithm::m_LowerTimeBoundaryCalibratedRPC
double m_LowerTimeBoundaryCalibratedRPC
Lower time boundary for RPC (calibrated data).
Definition: KLMTimeAlgorithm.h:324
Belle2::KLMTimeAlgorithm::m_useEventT0
bool m_useEventT0
Whether to use event T0 from CDC.
Definition: KLMTimeAlgorithm.h:425
Belle2::KLMTimeAlgorithm::m_MinimalDigitNumber
int m_MinimalDigitNumber
Minimal digit number (total).
Definition: KLMTimeAlgorithm.h:386
Belle2::KLMTimeAlgorithm::m_ProfileEKLMScintillatorPlane1
TProfile * m_ProfileEKLMScintillatorPlane1
For EKLM scintillator plane1.
Definition: KLMTimeAlgorithm.h:495
Belle2::KLMTimeAlgorithm::m_UpperTimeBoundaryRPC
double m_UpperTimeBoundaryRPC
Upper time boundary for RPC.
Definition: KLMTimeAlgorithm.h:309
Belle2::KLMTimeAlgorithm::h2c_timeF_rpc
TH2F * h2c_timeF_rpc[2]
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:604
Belle2::KLMTimeAlgorithm::fcn_land
TF1 * fcn_land
Landau function.
Definition: KLMTimeAlgorithm.h:750
Belle2::KLMTimeAlgorithm::m_timeCableDelay
KLMTimeCableDelay * m_timeCableDelay
DBObject of the calibration constant of each channel due to cable decay.
Definition: KLMTimeAlgorithm.h:413
Belle2::KLMTimeAlgorithm::tS_lowerStrip
std::pair< int, double > tS_lowerStrip(const KLMChannelIndex &klmChannel)
Tracing available channels with decreasing strip number.
Definition: KLMTimeAlgorithm.cc:1607
Belle2::KLMTimeAlgorithm::m_ProfileRpcZ
TProfile * m_ProfileRpcZ
For BKLM RPC z plane.
Definition: KLMTimeAlgorithm.h:486
Belle2::KLMTimeAlgorithm::mc_etime_channel
std::map< KLMChannelNumber, double > mc_etime_channel
Calibrated time distribution central value Error of each channel.
Definition: KLMTimeAlgorithm.h:303
Belle2::KLMTimeAlgorithm::tR_lowerStrip
std::pair< int, double > tR_lowerStrip(const KLMChannelIndex &klmChannel)
Tracing available channels with decreasing strip number.
Definition: KLMTimeAlgorithm.cc:1687
Belle2::KLMTimeAlgorithm::hc_timeFSLPC
TH1F * hc_timeFSLPC[2][8][15][2][54]
BKLM part.
Definition: KLMTimeAlgorithm.h:733
Belle2::KLMTimeAlgorithm::m_minimizerOptions
ROOT::Math::MinimizerOptions m_minimizerOptions
Minimization options.
Definition: KLMTimeAlgorithm.h:404
Belle2::KLMTimeAlgorithm::m_cFlag
std::map< KLMChannelNumber, int > m_cFlag
Calibration flag if the channel has enough hits collected and fitted OK.
Definition: KLMTimeAlgorithm.h:285
Belle2::KLMTimeAlgorithm::gre_time_channel_rpc
TGraphErrors * gre_time_channel_rpc
BKLM RPC.
Definition: KLMTimeAlgorithm.h:439
Belle2::KLMTimeAlgorithm::h2_timeFS_end
TH2F * h2_timeFS_end[2][4]
EKLM part.
Definition: KLMTimeAlgorithm.h:643
Belle2::KLMTimeAlgorithm::h2_timeF_rpc
TH2F * h2_timeF_rpc[2]
BKLM RPC part.
Definition: KLMTimeAlgorithm.h:590
Belle2::KLMTimeAlgorithm::m_evts
std::map< KLMChannelNumber, std::vector< struct Event > > m_evts
Container of hit information.
Definition: KLMTimeAlgorithm.h:279
Belle2::KLMTimeAlgorithm::h_timeFS_scint_end
TH1F * h_timeFS_scint_end[2][4]
EKLM part.
Definition: KLMTimeAlgorithm.h:621
Belle2::KLMTimeAlgorithm::m_time_channelAvg_rpc
double m_time_channelAvg_rpc
Central value of the global time distribution (BKLM RPC part).
Definition: KLMTimeAlgorithm.h:342
Belle2::KLMTimeAlgorithm::h2c_timeFS
TH2F * h2c_timeFS[2][8]
BKLM part.
Definition: KLMTimeAlgorithm.h:651
Belle2::KLMTimeAlgorithm::m_etime_channelAvg_rpc
double m_etime_channelAvg_rpc
Central value error of the global time distribution (BKLM RPC part).
Definition: KLMTimeAlgorithm.h:345
Belle2::KLMTimeAlgorithm::m_LowerTimeBoundaryScintillatorsEKLM
double m_LowerTimeBoundaryScintillatorsEKLM
Lower time boundary for EKLM scintillators.
Definition: KLMTimeAlgorithm.h:318
Belle2::KLMTimeAlgorithm::m_etime_channel
std::map< KLMChannelNumber, double > m_etime_channel
Time distribution central value Error of each channel.
Definition: KLMTimeAlgorithm.h:297
Belle2::KLMTimeAlgorithm::m_lower_limit_counts
int m_lower_limit_counts
Lower limit of hits collected for on single channel.
Definition: KLMTimeAlgorithm.h:389
Belle2::KLMTimeCableDelay
Class to store BKLM delay time coused by cable in the database.
Definition: KLMTimeCableDelay.h:29
Belle2::KLMTimeCableDelay::setTimeDelay
void setTimeDelay(KLMChannelNumber channel, float timeDelay)
Set the time calibration constant.
Definition: KLMTimeCableDelay.cc:25
Belle2::KLMTimeConstants
Class to store KLM constants related to time.
Definition: KLMTimeConstants.h:23
Belle2::KLMTimeConstants::c_BKLM
@ c_BKLM
BKLM scintillator.
Definition: KLMTimeConstants.h:36
Belle2::KLMTimeConstants::c_EKLM
@ c_EKLM
EKLM scintillator.
Definition: KLMTimeConstants.h:33
Belle2::KLMTimeConstants::setDelay
void setDelay(float delay, int cType)
Set effective light speed of scintillators.
Definition: KLMTimeConstants.cc:32
Belle2::KLMTimeResolution
Class to store KLM time resolution in the database.
Definition: KLMTimeResolution.h:30
Belle2::KLMTimeResolution::setTimeResolution
void setTimeResolution(KLMChannelNumber channel, float resolution)
Set time resolution.
Definition: KLMTimeResolution.cc:26
Belle2::StoreAccessorBase::registerInDataStore
bool registerInDataStore(DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut)
Register the object/array in the DataStore.
Definition: StoreAccessorBase.h:52
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: StoreObjPtr.h:96
Belle2::StoreObjPtr::isValid
bool isValid() const
Check whether the object was created.
Definition: StoreObjPtr.h:111
Belle2::StoreObjPtr::construct
bool construct(Args &&... params)
Construct an object of type T in this StoreObjPtr, using the provided constructor arguments.
Definition: StoreObjPtr.h:119
Belle2::Unit::cm
static const double cm
Standard units with the value = 1.
Definition: Unit.h:47
Belle2::bklm::GeometryPar::findModule
const Module * findModule(int section, int sector, int layer) const
Get the pointer to the definition of a module.
Definition: GeometryPar.cc:721
Belle2::bklm::GeometryPar::getMaximalZStripLength
double getMaximalZStripLength() const
Get maximal Z strip length (for scintillators).
Definition: GeometryPar.h:289
Belle2::bklm::GeometryPar::instance
static GeometryPar * instance(void)
Static method to get a reference to the singleton GeometryPar instance.
Definition: GeometryPar.cc:27
Belle2::bklm::GeometryPar::getMaximalPhiStripLength
double getMaximalPhiStripLength() const
Get maximal phi strip length (for scintillators).
Definition: GeometryPar.h:283
Belle2::bklm::Module
Define the geometry of a BKLM module Each sector [octant] contains Modules.
Definition: Module.h:76
LogVar
Class to store variables with their name which were sent to the logging service.
Definition: LogVariableStream.h:24
Belle2::DBStore::Instance
static DBStore & Instance()
Instance of a singleton DBStore.
Definition: DBStore.cc:28
Belle2::DBStore::updateEvent
void updateEvent()
Updates all intra-run dependent objects.
Definition: DBStore.cc:142
Belle2::DBStore::update
void update()
Updates all objects that are outside their interval of validity.
Definition: DBStore.cc:79
Belle2::sqrt
double sqrt(double a)
sqrt for double
Definition: beamHelpers.h:28
Belle2::KLMChannelNumber
uint16_t KLMChannelNumber
Channel number.
Definition: KLMElementNumberDefinitions.h:36
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:17