Project
Loading...
Searching...
No Matches
IntegratedClusterCalibrator.h
Go to the documentation of this file.
1// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3// All rights not expressly granted are reserved.
4//
5// This software is distributed under the terms of the GNU General Public
6// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7//
8// In applying this license CERN does not waive the privileges and immunities
9// granted to it by virtue of its status as an Intergovernmental Organization
10// or submit itself to any jurisdiction.
11
16
17#ifndef INTEGRATEDCLUSTERCALIBRATOR_H_
18#define INTEGRATEDCLUSTERCALIBRATOR_H_
19
22
23class TTree;
24
25namespace o2
26{
27
28// see https://stackoverflow.com/questions/56483053/if-condition-checking-multiple-vector-sizes-are-equal
29// comparing if all vector have the same size
30template <typename T0, typename... Ts>
31bool sameSize(T0 const& first, Ts const&... rest)
32{
33 return ((first.size() == rest.size()) && ...);
34}
35
36namespace tof
37{
38
40struct ITOFC {
41 std::vector<float> mITOFCNCl;
42 std::vector<float> mITOFCQ;
43 long mTimeMS{};
44 bool areSameSize() const { return sameSize(mITOFCNCl, mITOFCQ); }
45 bool isEmpty() const { return mITOFCNCl.empty(); }
46 size_t getEntries() const { return mITOFCNCl.size(); }
47 void setStartTime(long timeMS) { mTimeMS = timeMS; }
48
52 void fill(const unsigned int posIndex, const ITOFC& data)
53 {
54 std::copy(data.mITOFCNCl.begin(), data.mITOFCNCl.end(), mITOFCNCl.begin() + posIndex);
55 std::copy(data.mITOFCQ.begin(), data.mITOFCQ.end(), mITOFCQ.begin() + posIndex);
56 }
57
59 void insert(const unsigned int nDummyValues)
60 {
61 std::vector<float> vecTmp(nDummyValues, 0);
62 mITOFCNCl.insert(mITOFCNCl.begin(), vecTmp.begin(), vecTmp.end());
63 mITOFCQ.insert(mITOFCQ.begin(), vecTmp.begin(), vecTmp.end());
64 }
65
67 void resize(const unsigned int nTotal)
68 {
69 mITOFCNCl.resize(nTotal);
70 mITOFCQ.resize(nTotal);
71 }
72
74};
75} // end namespace tof
76
77namespace tpc
78{
79
81struct ITPCC {
82 std::vector<float> mIQMaxA;
83 std::vector<float> mIQMaxC;
84 std::vector<float> mIQTotA;
85 std::vector<float> mIQTotC;
86 std::vector<float> mINClA;
87 std::vector<float> mINClC;
88 long mTimeMS{};
89
90 float compression(float value, const int nBits) const
91 {
92 const int shiftN = std::pow(2, nBits);
93 int exp2;
94 const auto mantissa = std::frexp(value, &exp2);
95 const auto mantissaRounded = std::round(mantissa * shiftN) / shiftN;
96 return std::ldexp(mantissaRounded, exp2);
97 }
98
99 void compress(const int nBits)
100 {
101 std::transform(mIQMaxA.begin(), mIQMaxA.end(), mIQMaxA.begin(), [this, nBits](float val) { return compression(val, nBits); });
102 std::transform(mIQMaxC.begin(), mIQMaxC.end(), mIQMaxC.begin(), [this, nBits](float val) { return compression(val, nBits); });
103 std::transform(mIQTotA.begin(), mIQTotA.end(), mIQTotA.begin(), [this, nBits](float val) { return compression(val, nBits); });
104 std::transform(mIQTotC.begin(), mIQTotC.end(), mIQTotC.begin(), [this, nBits](float val) { return compression(val, nBits); });
105 std::transform(mINClA.begin(), mINClA.end(), mINClA.begin(), [this, nBits](float val) { return compression(val, nBits); });
106 std::transform(mINClC.begin(), mINClC.end(), mINClC.begin(), [this, nBits](float val) { return compression(val, nBits); });
107 }
108
110 bool isEmpty() const { return mIQMaxA.empty(); }
111 size_t getEntries() const { return mIQMaxA.size(); }
112 void setStartTime(long timeMS) { mTimeMS = timeMS; }
113
117 void fill(const unsigned int posIndex, const ITPCC& data)
118 {
119 std::copy(data.mIQMaxA.begin(), data.mIQMaxA.end(), mIQMaxA.begin() + posIndex);
120 std::copy(data.mIQMaxC.begin(), data.mIQMaxC.end(), mIQMaxC.begin() + posIndex);
121 std::copy(data.mIQTotA.begin(), data.mIQTotA.end(), mIQTotA.begin() + posIndex);
122 std::copy(data.mIQTotC.begin(), data.mIQTotC.end(), mIQTotC.begin() + posIndex);
123 std::copy(data.mINClA.begin(), data.mINClA.end(), mINClA.begin() + posIndex);
124 std::copy(data.mINClC.begin(), data.mINClC.end(), mINClC.begin() + posIndex);
125 }
126
128 void insert(const unsigned int nDummyValues)
129 {
130 std::vector<float> vecTmp(nDummyValues, 0);
131 mIQMaxA.insert(mIQMaxA.begin(), vecTmp.begin(), vecTmp.end());
132 mIQMaxC.insert(mIQMaxC.begin(), vecTmp.begin(), vecTmp.end());
133 mIQTotA.insert(mIQTotA.begin(), vecTmp.begin(), vecTmp.end());
134 mIQTotC.insert(mIQTotC.begin(), vecTmp.begin(), vecTmp.end());
135 mINClA.insert(mINClA.begin(), vecTmp.begin(), vecTmp.end());
136 mINClC.insert(mINClC.begin(), vecTmp.begin(), vecTmp.end());
137 }
138
140 void resize(const unsigned int nTotal)
141 {
142 mIQMaxA.resize(nTotal);
143 mIQMaxC.resize(nTotal);
144 mIQTotA.resize(nTotal);
145 mIQTotC.resize(nTotal);
146 mINClA.resize(nTotal);
147 mINClC.resize(nTotal);
148 }
149
151 void reset()
152 {
153 std::fill(mIQMaxA.begin(), mIQMaxA.end(), 0);
154 std::fill(mIQMaxC.begin(), mIQMaxC.end(), 0);
155 std::fill(mIQTotA.begin(), mIQTotA.end(), 0);
156 std::fill(mIQTotC.begin(), mIQTotC.end(), 0);
157 std::fill(mINClA.begin(), mINClA.end(), 0);
158 std::fill(mINClC.begin(), mINClC.end(), 0);
159 }
160
162 void normalize(const float factor)
163 {
164 std::transform(mIQMaxA.begin(), mIQMaxA.end(), mIQMaxA.begin(), [factor](const float val) { return val * factor; });
165 std::transform(mIQMaxC.begin(), mIQMaxC.end(), mIQMaxC.begin(), [factor](const float val) { return val * factor; });
166 std::transform(mIQTotA.begin(), mIQTotA.end(), mIQTotA.begin(), [factor](const float val) { return val * factor; });
167 std::transform(mIQTotC.begin(), mIQTotC.end(), mIQTotC.begin(), [factor](const float val) { return val * factor; });
168 std::transform(mINClA.begin(), mINClA.end(), mINClA.begin(), [factor](const float val) { return val * factor; });
169 std::transform(mINClC.begin(), mINClC.end(), mINClC.begin(), [factor](const float val) { return val * factor; });
170 }
171
173};
174
177 std::vector<float> mDCAr_A_Median;
178 std::vector<float> mDCAr_C_Median;
179 std::vector<float> mDCAr_A_WeightedMean;
180 std::vector<float> mDCAr_C_WeightedMean;
181 std::vector<float> mDCAr_A_RMS;
182 std::vector<float> mDCAr_C_RMS;
183 std::vector<float> mDCAr_A_NTracks;
184 std::vector<float> mDCAr_C_NTracks;
185 std::vector<float> mDCAz_A_Median;
186 std::vector<float> mDCAz_C_Median;
187 std::vector<float> mDCAz_A_WeightedMean;
188 std::vector<float> mDCAz_C_WeightedMean;
189 std::vector<float> mDCAz_A_RMS;
190 std::vector<float> mDCAz_C_RMS;
191 std::vector<float> mDCAz_A_NTracks;
192 std::vector<float> mDCAz_C_NTracks;
193 std::vector<float> mMIPdEdxRatioQMaxA;
194 std::vector<float> mMIPdEdxRatioQMaxC;
195 std::vector<float> mMIPdEdxRatioQTotA;
196 std::vector<float> mMIPdEdxRatioQTotC;
197 std::vector<float> mTPCChi2A;
198 std::vector<float> mTPCChi2C;
199 std::vector<float> mTPCNClA;
200 std::vector<float> mTPCNClC;
201
202 unsigned char mNBinsPhi{};
203 unsigned char mNBinsTgl{};
204 float mTglMax{};
205 unsigned char mNBinsqPt{};
206 float mQPtMax{};
207 unsigned char mMultBins{};
208 float mMultMax{};
209 long mTimeMS{};
210
214 void setStartTime(long timeMS) { mTimeMS = timeMS; }
215
217 int getNBins() const { return mNBinsTgl + mNBinsPhi + mNBinsqPt + mMultBins + 1; }
218
220 int getIndexPhi(const int iPhi, int slice = 0) const { return iPhi + slice * getNBins(); }
221
223 int getIndexTgl(const int iTgl, int slice = 0) const { return mNBinsPhi + iTgl + slice * getNBins(); }
224
226 int getIndexqPt(const int iqPt, int slice = 0) const { return mNBinsPhi + mNBinsTgl + iqPt + slice * getNBins(); }
227
229 int getIndexMult(const int iMult, int slice = 0) const { return mNBinsPhi + mNBinsTgl + mNBinsqPt + iMult + slice * getNBins(); }
230
232 int getIndexInt(int slice = 0) const { return getNBins() - 1 + slice * getNBins(); }
233
235 void resize(const unsigned int nTotal)
236 {
237 mDCAr_A_Median.resize(nTotal);
238 mDCAr_C_Median.resize(nTotal);
239 mDCAr_A_RMS.resize(nTotal);
240 mDCAr_C_RMS.resize(nTotal);
241 mDCAz_A_Median.resize(nTotal);
242 mDCAz_C_Median.resize(nTotal);
243 mDCAz_A_RMS.resize(nTotal);
244 mDCAz_C_RMS.resize(nTotal);
245 mDCAr_C_NTracks.resize(nTotal);
246 mDCAr_A_NTracks.resize(nTotal);
247 mDCAz_A_NTracks.resize(nTotal);
248 mDCAz_C_NTracks.resize(nTotal);
249 mDCAr_A_WeightedMean.resize(nTotal);
250 mDCAr_C_WeightedMean.resize(nTotal);
251 mDCAz_A_WeightedMean.resize(nTotal);
252 mDCAz_C_WeightedMean.resize(nTotal);
253 mMIPdEdxRatioQMaxA.resize(nTotal);
254 mMIPdEdxRatioQMaxC.resize(nTotal);
255 mMIPdEdxRatioQTotA.resize(nTotal);
256 mMIPdEdxRatioQTotC.resize(nTotal);
257 mTPCChi2A.resize(nTotal);
258 mTPCChi2C.resize(nTotal);
259 mTPCNClA.resize(nTotal);
260 mTPCNClC.resize(nTotal);
261 }
262
264};
265
267 std::vector<float> mITSTPC_A_MatchEff;
268 std::vector<float> mITSTPC_C_MatchEff;
269 std::vector<float> mITSTPC_A_Chi2Match;
270 std::vector<float> mITSTPC_C_Chi2Match;
271
273 void resize(const unsigned int nTotal)
274 {
275 mITSTPC_A_MatchEff.resize(nTotal);
276 mITSTPC_C_MatchEff.resize(nTotal);
277 mITSTPC_A_Chi2Match.resize(nTotal);
278 mITSTPC_C_Chi2Match.resize(nTotal);
279 }
280
282};
283
285 std::vector<float> mLogdEdx_A_Median;
286 std::vector<float> mLogdEdx_A_RMS;
287 std::vector<float> mLogdEdx_A_IROC_Median;
288 std::vector<float> mLogdEdx_A_IROC_RMS;
289 std::vector<float> mLogdEdx_A_OROC1_Median;
290 std::vector<float> mLogdEdx_A_OROC1_RMS;
291 std::vector<float> mLogdEdx_A_OROC2_Median;
292 std::vector<float> mLogdEdx_A_OROC2_RMS;
293 std::vector<float> mLogdEdx_A_OROC3_Median;
294 std::vector<float> mLogdEdx_A_OROC3_RMS;
295 std::vector<float> mLogdEdx_C_Median;
296 std::vector<float> mLogdEdx_C_RMS;
297 std::vector<float> mLogdEdx_C_IROC_Median;
298 std::vector<float> mLogdEdx_C_IROC_RMS;
299 std::vector<float> mLogdEdx_C_OROC1_Median;
300 std::vector<float> mLogdEdx_C_OROC1_RMS;
301 std::vector<float> mLogdEdx_C_OROC2_Median;
302 std::vector<float> mLogdEdx_C_OROC2_RMS;
303 std::vector<float> mLogdEdx_C_OROC3_Median;
304 std::vector<float> mLogdEdx_C_OROC3_RMS;
305
306 void resize(const unsigned int nTotal)
307 {
308 mLogdEdx_A_Median.resize(nTotal);
309 mLogdEdx_A_RMS.resize(nTotal);
310 mLogdEdx_A_IROC_Median.resize(nTotal);
311 mLogdEdx_A_IROC_RMS.resize(nTotal);
312 mLogdEdx_A_OROC1_Median.resize(nTotal);
313 mLogdEdx_A_OROC1_RMS.resize(nTotal);
314 mLogdEdx_A_OROC2_Median.resize(nTotal);
315 mLogdEdx_A_OROC2_RMS.resize(nTotal);
316 mLogdEdx_A_OROC3_Median.resize(nTotal);
317 mLogdEdx_A_OROC3_RMS.resize(nTotal);
318 mLogdEdx_C_Median.resize(nTotal);
319 mLogdEdx_C_RMS.resize(nTotal);
320 mLogdEdx_C_IROC_Median.resize(nTotal);
321 mLogdEdx_C_IROC_RMS.resize(nTotal);
322 mLogdEdx_C_OROC1_Median.resize(nTotal);
323 mLogdEdx_C_OROC1_RMS.resize(nTotal);
324 mLogdEdx_C_OROC2_Median.resize(nTotal);
325 mLogdEdx_C_OROC2_RMS.resize(nTotal);
326 mLogdEdx_C_OROC3_Median.resize(nTotal);
327 mLogdEdx_C_OROC3_RMS.resize(nTotal);
328 }
329
331};
332
341 std::vector<unsigned int> mOccupancyMapTPC;
342
343 std::vector<float> nPrimVertices;
344 std::vector<float> nPrimVertices_ITS;
346 std::vector<float> nVertexContributors_ITS_RMS;
347 std::vector<float> vertexX_ITS_Median;
348 std::vector<float> vertexY_ITS_Median;
349 std::vector<float> vertexZ_ITS_Median;
350 std::vector<float> vertexX_ITS_RMS;
351 std::vector<float> vertexY_ITS_RMS;
352 std::vector<float> vertexZ_ITS_RMS;
353
354 std::vector<float> nPrimVertices_ITSTPC;
357 std::vector<float> vertexX_ITSTPC_Median;
358 std::vector<float> vertexY_ITSTPC_Median;
359 std::vector<float> vertexZ_ITSTPC_Median;
360 std::vector<float> vertexX_ITSTPC_RMS;
361 std::vector<float> vertexY_ITSTPC_RMS;
362 std::vector<float> vertexZ_ITSTPC_RMS;
363
364 int quantileValues = 23;
365 std::vector<float> nVertexContributors_Quantiles;
366
367 std::vector<float> mDCAr_comb_A_Median;
368 std::vector<float> mDCAz_comb_A_Median;
369 std::vector<float> mDCAr_comb_A_RMS;
370 std::vector<float> mDCAz_comb_A_RMS;
371 std::vector<float> mDCAr_comb_C_Median;
372 std::vector<float> mDCAz_comb_C_Median;
373 std::vector<float> mDCAr_comb_C_RMS;
374 std::vector<float> mDCAz_comb_C_RMS;
375 std::vector<float> mITS_A_NCl_Median;
376 std::vector<float> mITS_A_NCl_RMS;
377 std::vector<float> mITS_C_NCl_Median;
378 std::vector<float> mITS_C_NCl_RMS;
379 std::vector<float> mSqrtITSChi2_Ncl_A_Median;
380 std::vector<float> mSqrtITSChi2_Ncl_C_Median;
381 std::vector<float> mSqrtITSChi2_Ncl_A_RMS;
382 std::vector<float> mSqrtITSChi2_Ncl_C_RMS;
383
384 std::vector<float> mITSTPCDeltaP2_A_Median;
385 std::vector<float> mITSTPCDeltaP3_A_Median;
386 std::vector<float> mITSTPCDeltaP4_A_Median;
387 std::vector<float> mITSTPCDeltaP2_C_Median;
388 std::vector<float> mITSTPCDeltaP3_C_Median;
389 std::vector<float> mITSTPCDeltaP4_C_Median;
390 std::vector<float> mITSTPCDeltaP2_A_RMS;
391 std::vector<float> mITSTPCDeltaP3_A_RMS;
392 std::vector<float> mITSTPCDeltaP4_A_RMS;
393 std::vector<float> mITSTPCDeltaP2_C_RMS;
394 std::vector<float> mITSTPCDeltaP3_C_RMS;
395 std::vector<float> mITSTPCDeltaP4_C_RMS;
396
397 std::vector<float> mTPCSigmaY2A_Median;
398 std::vector<float> mTPCSigmaZ2A_Median;
399 std::vector<float> mTPCSigmaY2C_Median;
400 std::vector<float> mTPCSigmaZ2C_Median;
401 std::vector<float> mTPCSigmaY2A_RMS;
402 std::vector<float> mTPCSigmaZ2A_RMS;
403 std::vector<float> mTPCSigmaY2C_RMS;
404 std::vector<float> mTPCSigmaZ2C_RMS;
405
406 void setStartTime(long timeMS)
407 {
408 mTSTPC.setStartTime(timeMS);
409 mTSITSTPC.setStartTime(timeMS);
410 }
411
412 void setBinning(const int nBinsPhi, const int nBinsTgl, const int qPtBins, const int nBinsMult, float tglMax, float qPtMax, float multMax)
413 {
414 mTSTPC.mNBinsPhi = nBinsPhi;
415 mTSTPC.mNBinsTgl = nBinsTgl;
416 mTSTPC.mNBinsqPt = qPtBins;
417 mTSTPC.mMultBins = nBinsMult;
418 mTSTPC.mTglMax = tglMax;
419 mTSTPC.mQPtMax = qPtMax;
420 mTSTPC.mMultMax = multMax;
421 mTSITSTPC.mNBinsPhi = nBinsPhi;
422 mTSITSTPC.mNBinsTgl = nBinsTgl;
423 mTSITSTPC.mNBinsqPt = qPtBins;
424 mTSITSTPC.mMultBins = nBinsMult;
425 mTSITSTPC.mTglMax = tglMax;
426 mTSITSTPC.mQPtMax = qPtMax;
427 mTSITSTPC.mMultMax = multMax;
428 }
429
431 void resize(const unsigned int nTotal)
432 {
433 mTSTPC.resize(nTotal);
434 mTSITSTPC.resize(nTotal);
435 mITSTPCAll.resize(nTotal);
438 mdEdxQTot.resize(nTotal);
439 mdEdxQMax.resize(nTotal);
440 mDCAr_comb_A_Median.resize(nTotal);
441 mDCAz_comb_A_Median.resize(nTotal);
442 mDCAr_comb_A_RMS.resize(nTotal);
443 mDCAz_comb_A_RMS.resize(nTotal);
444 mDCAr_comb_C_Median.resize(nTotal);
445 mDCAz_comb_C_Median.resize(nTotal);
446 mDCAr_comb_C_RMS.resize(nTotal);
447 mDCAz_comb_C_RMS.resize(nTotal);
448 mITS_A_NCl_Median.resize(nTotal);
449 mITS_A_NCl_RMS.resize(nTotal);
450 mITS_C_NCl_Median.resize(nTotal);
451 mITS_C_NCl_RMS.resize(nTotal);
452 mSqrtITSChi2_Ncl_A_Median.resize(nTotal);
453 mSqrtITSChi2_Ncl_C_Median.resize(nTotal);
454 mSqrtITSChi2_Ncl_A_RMS.resize(nTotal);
455 mSqrtITSChi2_Ncl_C_RMS.resize(nTotal);
456 mITSTPCDeltaP2_A_Median.resize(nTotal);
457 mITSTPCDeltaP3_A_Median.resize(nTotal);
458 mITSTPCDeltaP4_A_Median.resize(nTotal);
459 mITSTPCDeltaP2_C_Median.resize(nTotal);
460 mITSTPCDeltaP3_C_Median.resize(nTotal);
461 mITSTPCDeltaP4_C_Median.resize(nTotal);
462 mITSTPCDeltaP2_A_RMS.resize(nTotal);
463 mITSTPCDeltaP3_A_RMS.resize(nTotal);
464 mITSTPCDeltaP4_A_RMS.resize(nTotal);
465 mITSTPCDeltaP2_C_RMS.resize(nTotal);
466 mITSTPCDeltaP3_C_RMS.resize(nTotal);
467 mITSTPCDeltaP4_C_RMS.resize(nTotal);
468 mTPCSigmaY2A_Median.resize(nTotal);
469 mTPCSigmaZ2A_Median.resize(nTotal);
470 mTPCSigmaY2C_Median.resize(nTotal);
471 mTPCSigmaZ2C_Median.resize(nTotal);
472 mTPCSigmaY2A_RMS.resize(nTotal);
473 mTPCSigmaZ2A_RMS.resize(nTotal);
474 mTPCSigmaY2C_RMS.resize(nTotal);
475 mTPCSigmaZ2C_RMS.resize(nTotal);
476
477 const int nTotalVtx = nTotal / mTSTPC.getNBins();
478 nPrimVertices.resize(nTotalVtx);
479 nPrimVertices_ITS.resize(nTotalVtx);
480 nVertexContributors_ITS_Median.resize(nTotalVtx);
481 nVertexContributors_ITS_RMS.resize(nTotalVtx);
482 vertexX_ITS_Median.resize(nTotalVtx);
483 vertexY_ITS_Median.resize(nTotalVtx);
484 vertexZ_ITS_Median.resize(nTotalVtx);
485 vertexX_ITS_RMS.resize(nTotalVtx);
486 vertexY_ITS_RMS.resize(nTotalVtx);
487 vertexZ_ITS_RMS.resize(nTotalVtx);
488 nPrimVertices_ITSTPC.resize(nTotalVtx);
489 nVertexContributors_ITSTPC_Median.resize(nTotalVtx);
490 nVertexContributors_ITSTPC_RMS.resize(nTotalVtx);
491 vertexX_ITSTPC_Median.resize(nTotalVtx);
492 vertexY_ITSTPC_Median.resize(nTotalVtx);
493 vertexZ_ITSTPC_Median.resize(nTotalVtx);
494 vertexX_ITSTPC_RMS.resize(nTotalVtx);
495 vertexY_ITSTPC_RMS.resize(nTotalVtx);
496 vertexZ_ITSTPC_RMS.resize(nTotalVtx);
497
498 const int nTotalQ = quantileValues * nTotal / mTSTPC.getNBins();
499 nVertexContributors_Quantiles.resize(nTotalQ);
500 }
501
503};
504
505} // end namespace tpc
506
507namespace fit
508{
509
511struct IFT0C {
512 std::vector<float> mINChanA;
513 std::vector<float> mINChanC;
514 std::vector<float> mIAmplA;
515 std::vector<float> mIAmplC;
516 long mTimeMS{};
517 bool areSameSize() const { return sameSize(mINChanA, mINChanC, mIAmplA, mIAmplC); }
518 bool isEmpty() const { return mINChanA.empty(); }
519 size_t getEntries() const { return mINChanA.size(); }
520 void setStartTime(long timeMS) { mTimeMS = timeMS; }
521
525 void fill(const unsigned int posIndex, const IFT0C& data)
526 {
527 std::copy(data.mINChanA.begin(), data.mINChanA.end(), mINChanA.begin() + posIndex);
528 std::copy(data.mINChanC.begin(), data.mINChanC.end(), mINChanC.begin() + posIndex);
529 std::copy(data.mIAmplA.begin(), data.mIAmplA.end(), mIAmplA.begin() + posIndex);
530 std::copy(data.mIAmplC.begin(), data.mIAmplC.end(), mIAmplC.begin() + posIndex);
531 }
532
534 void insert(const unsigned int nDummyValues)
535 {
536 std::vector<float> vecTmp(nDummyValues, 0);
537 mINChanA.insert(mINChanA.begin(), vecTmp.begin(), vecTmp.end());
538 mINChanC.insert(mINChanC.begin(), vecTmp.begin(), vecTmp.end());
539 mIAmplA.insert(mIAmplA.begin(), vecTmp.begin(), vecTmp.end());
540 mIAmplC.insert(mIAmplC.begin(), vecTmp.begin(), vecTmp.end());
541 }
542
544 void resize(const unsigned int nTotal)
545 {
546 mINChanA.resize(nTotal);
547 mINChanC.resize(nTotal);
548 mIAmplA.resize(nTotal);
549 mIAmplC.resize(nTotal);
550 }
551
553 void reset()
554 {
555 std::fill(mINChanA.begin(), mINChanA.end(), 0);
556 std::fill(mINChanC.begin(), mINChanC.end(), 0);
557 std::fill(mIAmplA.begin(), mIAmplA.end(), 0);
558 std::fill(mIAmplC.begin(), mIAmplC.end(), 0);
559 }
560
562 void normalize(const float factor)
563 {
564 std::transform(mINChanA.begin(), mINChanA.end(), mINChanA.begin(), [factor](const float val) { return val * factor; });
565 std::transform(mINChanC.begin(), mINChanC.end(), mINChanC.begin(), [factor](const float val) { return val * factor; });
566 std::transform(mIAmplA.begin(), mIAmplA.end(), mIAmplA.begin(), [factor](const float val) { return val * factor; });
567 std::transform(mIAmplC.begin(), mIAmplC.end(), mIAmplC.begin(), [factor](const float val) { return val * factor; });
568 }
569
571};
572
574struct IFV0C {
575 std::vector<float> mINChanA;
576 std::vector<float> mIAmplA;
577 long mTimeMS{};
578 bool areSameSize() const { return sameSize(mINChanA, mIAmplA); }
579 bool isEmpty() const { return mINChanA.empty(); }
580 size_t getEntries() const { return mINChanA.size(); }
581 void setStartTime(long timeMS) { mTimeMS = timeMS; }
582
586 void fill(const unsigned int posIndex, const IFV0C& data)
587 {
588 std::copy(data.mINChanA.begin(), data.mINChanA.end(), mINChanA.begin() + posIndex);
589 std::copy(data.mIAmplA.begin(), data.mIAmplA.end(), mIAmplA.begin() + posIndex);
590 }
591
593 void insert(const unsigned int nDummyValues)
594 {
595 std::vector<float> vecTmp(nDummyValues, 0);
596 mINChanA.insert(mINChanA.begin(), vecTmp.begin(), vecTmp.end());
597 mIAmplA.insert(mIAmplA.begin(), vecTmp.begin(), vecTmp.end());
598 }
599
601 void resize(const unsigned int nTotal)
602 {
603 mINChanA.resize(nTotal);
604 mIAmplA.resize(nTotal);
605 }
606
608 void reset()
609 {
610 std::fill(mINChanA.begin(), mINChanA.end(), 0);
611 std::fill(mIAmplA.begin(), mIAmplA.end(), 0);
612 }
613
615 void normalize(const float factor)
616 {
617 std::transform(mINChanA.begin(), mINChanA.end(), mINChanA.begin(), [factor](const float val) { return val * factor; });
618 std::transform(mIAmplA.begin(), mIAmplA.end(), mIAmplA.begin(), [factor](const float val) { return val * factor; });
619 }
620
622};
623
625struct IFDDC {
626 std::vector<float> mINChanA;
627 std::vector<float> mINChanC;
628 std::vector<float> mIAmplA;
629 std::vector<float> mIAmplC;
630 long mTimeMS{};
631 bool areSameSize() const { return sameSize(mINChanA, mINChanC, mIAmplA, mIAmplC); }
632 bool isEmpty() const { return mINChanA.empty(); }
633 size_t getEntries() const { return mINChanA.size(); }
634 void setStartTime(long timeMS) { mTimeMS = timeMS; }
635
639 void fill(const unsigned int posIndex, const IFDDC& data)
640 {
641 std::copy(data.mINChanA.begin(), data.mINChanA.end(), mINChanA.begin() + posIndex);
642 std::copy(data.mINChanC.begin(), data.mINChanC.end(), mINChanC.begin() + posIndex);
643 std::copy(data.mIAmplA.begin(), data.mIAmplA.end(), mIAmplA.begin() + posIndex);
644 std::copy(data.mIAmplC.begin(), data.mIAmplC.end(), mIAmplC.begin() + posIndex);
645 }
646
648 void insert(const unsigned int nDummyValues)
649 {
650 std::vector<float> vecTmp(nDummyValues, 0);
651 mINChanA.insert(mINChanA.begin(), vecTmp.begin(), vecTmp.end());
652 mINChanC.insert(mINChanC.begin(), vecTmp.begin(), vecTmp.end());
653 mIAmplA.insert(mIAmplA.begin(), vecTmp.begin(), vecTmp.end());
654 mIAmplC.insert(mIAmplC.begin(), vecTmp.begin(), vecTmp.end());
655 }
656
658 void resize(const unsigned int nTotal)
659 {
660 mINChanA.resize(nTotal);
661 mINChanC.resize(nTotal);
662 mIAmplA.resize(nTotal);
663 mIAmplC.resize(nTotal);
664 }
665
667 void reset()
668 {
669 std::fill(mINChanA.begin(), mINChanA.end(), 0);
670 std::fill(mINChanC.begin(), mINChanC.end(), 0);
671 std::fill(mIAmplA.begin(), mIAmplA.end(), 0);
672 std::fill(mIAmplC.begin(), mIAmplC.end(), 0);
673 }
674
676 void normalize(const float factor)
677 {
678 std::transform(mINChanA.begin(), mINChanA.end(), mINChanA.begin(), [factor](const float val) { return val * factor; });
679 std::transform(mINChanC.begin(), mINChanC.end(), mINChanC.begin(), [factor](const float val) { return val * factor; });
680 std::transform(mIAmplA.begin(), mIAmplA.end(), mIAmplA.begin(), [factor](const float val) { return val * factor; });
681 std::transform(mIAmplC.begin(), mIAmplC.end(), mIAmplC.begin(), [factor](const float val) { return val * factor; });
682 }
683
685};
686
687} // end namespace fit
688} // end namespace o2
689
690namespace o2
691{
692namespace calibration
693{
694
696template <typename DataT>
698{
699 public:
703 IntegratedClusters(o2::calibration::TFType tFirst, o2::calibration::TFType tLast) : mTFFirst{tFirst}, mTFLast{tLast} {};
704
707
709 void print() const { LOGP(info, "TF Range from {} to {} with {} of remaining data", mTFFirst, mTFLast, mRemainingData); }
710
714 void fill(const o2::calibration::TFType tfID, const DataT& currentsContainer);
715
717 void merge(const IntegratedClusters* prev);
718
720 bool hasEnoughData() const { return (mRemainingData != -1); }
721
723 const auto& getCurrents() const& { return mCurrents; }
724
726 auto getCurrents() && { return std::move(mCurrents); }
727
729 void setCurrents(const DataT& currents) { mCurrents = currents; }
730
734 void dumpToFile(const char* outFileName = "IntegratedClusters.root", const char* outName = "IC") const;
735
738 void dumpToTree(const char* outFileName = "ICTree.root");
739
741 void setStartTime(long timeMS) { mCurrents.setStartTime(timeMS); }
742
743 private:
744 DataT mCurrents;
745 o2::calibration::TFType mTFFirst{};
746 o2::calibration::TFType mTFLast{};
747 o2::calibration::TFType mRemainingData = -1;
748 unsigned int mNValuesPerTF{};
749 bool mInitialize{true};
750
753 void initData(const unsigned int valuesPerTF);
754
755 ClassDefNV(IntegratedClusters, 1);
756};
757
758template <typename DataT>
760{
761 using TFType = o2::calibration::TFType;
763 using CalibVector = std::vector<DataT>;
764 using TFinterval = std::vector<std::pair<TFType, TFType>>;
765 using TimeInterval = std::vector<std::pair<long, long>>;
766
767 public:
770
773
775 bool hasEnoughData(const Slot& slot) const final { return slot.getContainer()->hasEnoughData(); }
776
778 void initOutput() final;
779
781 void finalizeSlot(Slot& slot) final;
782
784 Slot& emplaceNewSlot(bool front, TFType tstart, TFType tend) final;
785
787 const TFinterval& getTFinterval() const { return mIntervals; }
788
790 const TimeInterval& getTimeIntervals() const { return mTimeIntervals; }
791
793 auto getCalibs() && { return std::move(mCalibs); }
794
796 bool hasCalibrationData() const { return mCalibs.size() > 0; }
797
799 void setDebug(const bool debug) { mDebug = debug; }
800
801 private:
802 TFinterval mIntervals;
803 TimeInterval mTimeIntervals;
804 CalibVector mCalibs;
805 bool mDebug{false};
806
807 ClassDefOverride(IntegratedClusterCalibrator, 1);
808};
809
810} // end namespace calibration
811} // end namespace o2
812
813#endif
std::ostringstream debug
void initOutput() final
clearing all calibration objects in the output buffer
~IntegratedClusterCalibrator() final=default
default destructor
bool hasCalibrationData() const
check if calibration data is available
Slot & emplaceNewSlot(bool front, TFType tstart, TFType tend) final
Creates new time slot.
IntegratedClusterCalibrator()=default
default constructor
bool hasEnoughData(const Slot &slot) const final
check if given slot has already enough data
void setDebug(const bool debug)
set if debug objects will be created
void finalizeSlot(Slot &slot) final
storing the integrated currents for given slot
class for accumulating integrated currents
IntegratedClusters()=default
\default constructor for ROOT I/O
void dumpToTree(const char *outFileName="ICTree.root")
void merge(const IntegratedClusters *prev)
merging TOF currents with previous interval
IntegratedClusters(o2::calibration::TFType tFirst, o2::calibration::TFType tLast)
void print() const
print summary informations
void setStartTime(long timeMS)
setting the start time
void fill(const o2::calibration::TFType tfID, const DataT &currentsContainer)
void dumpToFile(const char *outFileName="IntegratedClusters.root", const char *outName="IC") const
GLsizei const GLfloat * value
Definition glcorearb.h:819
GLboolean * data
Definition glcorearb.h:298
GLuint GLfloat * val
Definition glcorearb.h:1582
uint32_t TFType
Definition TimeSlot.h:29
a couple of static helper functions to create timestamp values for CCDB queries or override obsolete ...
bool sameSize(T0 const &first, Ts const &... rest)
struct containing the integrated FDD currents
ClassDefNV(IFDDC, 2)
std::vector< float > mINChanA
integrated 1D FIT currents for NChan A
void normalize(const float factor)
normalize currents
long mTimeMS
start time in ms
std::vector< float > mIAmplC
integrated 1D FIT currents for Ampl C
bool isEmpty() const
check if values are empty
void insert(const unsigned int nDummyValues)
void fill(const unsigned int posIndex, const IFDDC &data)
std::vector< float > mINChanC
integrated 1D FIT currents for NChan C
void reset()
reset buffered currents
std::vector< float > mIAmplA
integrated 1D FIT currents for Ampl A
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
bool areSameSize() const
check if stored currents have same number of entries
void setStartTime(long timeMS)
struct containing the integrated FT0 currents
bool isEmpty() const
check if values are empty
ClassDefNV(IFT0C, 2)
void setStartTime(long timeMS)
std::vector< float > mIAmplA
integrated 1D FIT currents for Ampl A
void fill(const unsigned int posIndex, const IFT0C &data)
std::vector< float > mINChanA
integrated 1D FIT currents for NChan A
void normalize(const float factor)
normalize currents
bool areSameSize() const
check if stored currents have same number of entries
long mTimeMS
start time in ms
std::vector< float > mIAmplC
integrated 1D FIT currents for Ampl C
std::vector< float > mINChanC
integrated 1D FIT currents for NChan C
void reset()
reset buffered currents
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
void insert(const unsigned int nDummyValues)
struct containing the integrated FV0 currents
bool isEmpty() const
check if values are empty
bool areSameSize() const
check if stored currents have same number of entries
ClassDefNV(IFV0C, 2)
void fill(const unsigned int posIndex, const IFV0C &data)
long mTimeMS
start time in ms
void reset()
reset buffered currents
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
std::vector< float > mIAmplA
integrated 1D FIT currents for Ampl A
void normalize(const float factor)
normalize currents
void setStartTime(long timeMS)
std::vector< float > mINChanA
integrated 1D FIT currents for NChan A
void insert(const unsigned int nDummyValues)
struct containing the integrated TOF currents
bool isEmpty() const
check if values are empty
std::vector< float > mITOFCNCl
integrated 1D TOF cluster currents
void fill(const unsigned int posIndex, const ITOFC &data)
void setStartTime(long timeMS)
ClassDefNV(ITOFC, 2)
void insert(const unsigned int nDummyValues)
long mTimeMS
start time in ms
std::vector< float > mITOFCQ
integrated 1D TOF qTot currents
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
bool areSameSize() const
check if stored currents have same number of entries
struct containing the integrated TPC currents
ClassDefNV(ITPCC, 2)
std::vector< float > mINClA
integrated 1D-currents for NCl A-side
void fill(const unsigned int posIndex, const ITPCC &data)
std::vector< float > mIQMaxC
integrated 1D-currents for QMax C-side
void insert(const unsigned int nDummyValues)
bool areSameSize() const
check if stored currents have same number of entries
void setStartTime(long timeMS)
std::vector< float > mIQTotA
integrated 1D-currents for QTot A-side
float compression(float value, const int nBits) const
void normalize(const float factor)
normalize currents
std::vector< float > mIQMaxA
integrated 1D-currents for QMax A-side
void compress(const int nBits)
std::vector< float > mINClC
integrated 1D-currents for NCl A-side
std::vector< float > mIQTotC
integrated 1D-currents for QTot A-side
bool isEmpty() const
check if values are empty
void reset()
reset buffered currents
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
long mTimeMS
start time in ms
std::vector< float > mITSTPC_A_MatchEff
matching efficiency of ITS-TPC tracks A-side
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
std::vector< float > mITSTPC_C_MatchEff
matching efficiency of ITS-TPC tracks C-side
ClassDefNV(ITSTPC_Matching, 2)
std::vector< float > mITSTPC_A_Chi2Match
ITS-TPC chi2 A-side.
std::vector< float > mITSTPC_C_Chi2Match
ITS-TPC chi2 C-side.
std::vector< float > mDCAr_comb_A_RMS
DCAr RMS for ITS-TPC track - A-side.
std::vector< float > mTPCSigmaZ2A_RMS
sigmaZ2 RMS at vertex
std::vector< float > mSqrtITSChi2_Ncl_C_Median
sqrt(ITC chi2 / ncl)
ITSTPC_Matching mITSTPCAll
ITS-TPC matching efficiency for ITS standalone + afterburner.
std::vector< float > nPrimVertices
number of primary vertices
std::vector< float > mITSTPCDeltaP4_C_RMS
RMS of track param TPC - track param ITS-TPC for param 4 - A-side.
std::vector< float > mTPCSigmaZ2C_RMS
sigmaZ2 RMS at vertex
TimeSeries mTSTPC
TPC standalone DCAs.
std::vector< float > mITSTPCDeltaP4_A_Median
track param TPC - track param ITS-TPC for param 4 - A-side
ITSTPC_Matching mITSTPCStandalone
ITS-TPC matching efficiency for ITS standalone.
std::vector< float > vertexY_ITS_RMS
vertex y RMS selected with ITS cut 0.2<nContributorsITS/nContributors<0.8
std::vector< float > mITSTPCDeltaP4_C_Median
track param TPC - track param ITS-TPC for param 4 - A-side
std::vector< float > nPrimVertices_ITS
number of primary vertices selected with ITS cut 0.2<nContributorsITS/nContributors<0....
std::vector< float > vertexZ_ITSTPC_Median
vertex z position with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0.95
std::vector< float > mITSTPCDeltaP3_A_RMS
RMS of track param TPC - track param ITS-TPC for param 3 - A-side.
std::vector< float > mITS_A_NCl_Median
its number of clusters
std::vector< float > mDCAr_comb_C_Median
DCAr for ITS-TPC track - C-side.
std::vector< float > nVertexContributors_ITS_Median
number of primary vertices selected with ITS cut 0.2<nContributorsITS/nContributors<0....
std::vector< float > mSqrtITSChi2_Ncl_A_Median
sqrt(ITC chi2 / ncl)
std::vector< float > vertexZ_ITS_Median
vertex z position selected with ITS cut 0.2<nContributorsITS/nContributors<0.8
std::vector< float > mSqrtITSChi2_Ncl_C_RMS
sqrt(ITC chi2 / ncl)
std::vector< float > mDCAz_comb_A_RMS
DCAz RMS for ITS-TPC track - A-side.
std::vector< float > vertexX_ITSTPC_RMS
vertex x RMS with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0.95
std::vector< float > mITSTPCDeltaP2_C_RMS
RMS of track param TPC - track param ITS-TPC for param 2 - A-side.
std::vector< float > mITSTPCDeltaP3_A_Median
track param TPC - track param ITS-TPC for param 3 - A-side
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
std::vector< float > mITSTPCDeltaP2_A_Median
track param TPC - track param ITS-TPC for param 2 - A-side
std::vector< float > vertexZ_ITS_RMS
vertex z RMS selected with ITS cut 0.2<nContributorsITS/nContributors<0.8
std::vector< float > mSqrtITSChi2_Ncl_A_RMS
sqrt(ITC chi2 / ncl)
ClassDefNV(TimeSeriesITSTPC, 5)
std::vector< float > vertexX_ITS_Median
vertex x position selected with ITS cut 0.2<nContributorsITS/nContributors<0.8
std::vector< float > mITSTPCDeltaP2_C_Median
track param TPC - track param ITS-TPC for param 2 - A-side
std::vector< unsigned int > mOccupancyMapTPC
cluster occupancy map
std::vector< float > mITSTPCDeltaP4_A_RMS
RMS of track param TPC - track param ITS-TPC for param 4 - A-side.
std::vector< float > mTPCSigmaY2A_Median
sigmaY2 at vertex
std::vector< float > vertexZ_ITSTPC_RMS
vertex z RMS with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0.95
std::vector< float > nVertexContributors_ITSTPC_Median
number of primary vertices with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0....
std::vector< float > vertexX_ITS_RMS
vertex x RMS selected with ITS cut 0.2<nContributorsITS/nContributors<0.8
std::vector< float > nVertexContributors_ITS_RMS
number of primary vertices selected with ITS cut 0.2<nContributorsITS/nContributors<0....
std::vector< float > mITSTPCDeltaP3_C_RMS
RMS of track param TPC - track param ITS-TPC for param 3 - A-side.
std::vector< float > nVertexContributors_Quantiles
number of primary vertices for quantiles 0.1, 0.2, ... 0.9 and truncated mean values 0....
TimeSeries mTSITSTPC
ITS-TPC standalone DCAs.
std::vector< float > mTPCSigmaY2C_Median
sigmaY2 at vertex
std::vector< float > vertexY_ITSTPC_Median
vertex y position with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0.95
std::vector< float > nPrimVertices_ITSTPC
number of primary vertices with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0....
std::vector< float > mTPCSigmaZ2C_Median
sigmaZ2 at vertex
std::vector< float > mITS_C_NCl_RMS
its number of clusters
std::vector< float > vertexY_ITS_Median
vertex y position selected with ITS cut 0.2<nContributorsITS/nContributors<0.8
std::vector< float > nVertexContributors_ITSTPC_RMS
number of primary vertices with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0....
std::vector< float > mDCAz_comb_A_Median
DCAz for ITS-TPC track - A-side.
std::vector< float > mDCAz_comb_C_RMS
DCAz RMS for ITS-TPC track - C-side.
ITSTPC_Matching mITSTPCAfterburner
ITS-TPC matchin efficiency fir ITS afterburner.
std::vector< float > mTPCSigmaY2C_RMS
sigmaY2 RMS at vertex
int quantileValues
! number of values in quantiles + truncated mean (hardcoded for the moment)
std::vector< float > mITSTPCDeltaP3_C_Median
track param TPC - track param ITS-TPC for param 3 - A-side
std::vector< float > mDCAz_comb_C_Median
DCAz for ITS-TPC track - C-side.
std::vector< float > mTPCSigmaZ2A_Median
sigmaZ2 at vertex
std::vector< float > mDCAr_comb_C_RMS
DCAr RMS for ITS-TPC track - C-side.
std::vector< float > vertexY_ITSTPC_RMS
vertex y RMS with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0.95
std::vector< float > mITS_C_NCl_Median
its number of clusters
TimeSeriesdEdx mdEdxQTot
time series for dE/dx qTot monitoring
std::vector< float > mITS_A_NCl_RMS
its number of clusters
std::vector< float > mITSTPCDeltaP2_A_RMS
RMS of track param TPC - track param ITS-TPC for param 2 - A-side.
std::vector< float > mTPCSigmaY2A_RMS
sigmaY2 RMS at vertex
std::vector< float > mDCAr_comb_A_Median
DCAr for ITS-TPC track - A-side.
std::vector< float > vertexX_ITSTPC_Median
vertex x position with ITS-TPC cut (nContributorsITS + nContributorsITSTPC)<0.95
void setBinning(const int nBinsPhi, const int nBinsTgl, const int qPtBins, const int nBinsMult, float tglMax, float qPtMax, float multMax)
TimeSeriesdEdx mdEdxQMax
time series for dE/dx qMax monitoring
struct containing time series values
unsigned char mNBinsTgl
number of phi bins
std::vector< float > mDCAr_A_NTracks
number of tracks used to calculate the DCAs
std::vector< float > mMIPdEdxRatioQTotA
ratio of MIP/dEdx - qTot -
unsigned char mMultBins
multiplicity bins
std::vector< float > mMIPdEdxRatioQTotC
ratio of MIP/dEdx - qTot -
int getIndexInt(int slice=0) const
unsigned char mNBinsPhi
number of tgl bins
std::vector< float > mDCAr_C_NTracks
number of tracks used to calculate the DCAs
std::vector< float > mDCAr_C_Median
integrated 1D DCAr for C-side weighted mean in phi/tgl slices
std::vector< float > mTPCNClC
number of TPC cluster
ClassDefNV(TimeSeries, 1)
std::vector< float > mMIPdEdxRatioQMaxC
ratio of MIP/dEdx - qMax -
std::vector< float > mDCAr_A_WeightedMean
integrated 1D DCAr for A-side weighted mean in phi/tgl slices
std::vector< float > mDCAz_A_RMS
integrated 1D DCAz for A-side RMS in phi/tgl slices
float mMultMax
max local multiplicity
std::vector< float > mTPCNClA
number of TPC cluster
std::vector< float > mDCAz_C_Median
integrated 1D DCAz for C-side median in phi/tgl slices
int getIndexPhi(const int iPhi, int slice=0) const
std::vector< float > mTPCChi2C
Chi2 of TPC tracks.
std::vector< float > mDCAz_A_Median
integrated 1D DCAz for A-side median in phi/tgl slices
std::vector< float > mDCAr_C_RMS
integrated 1D DCAr for C-side RMS in phi/tgl slices
std::vector< float > mDCAz_A_NTracks
number of tracks used to calculate the DCAs
std::vector< float > mDCAr_A_Median
integrated 1D DCAr for A-side median in phi/tgl slices
std::vector< float > mDCAz_C_NTracks
number of tracks used to calculate the DCAs
int getIndexTgl(const int iTgl, int slice=0) const
std::vector< float > mMIPdEdxRatioQMaxA
ratio of MIP/dEdx - qMax -
void resize(const unsigned int nTotal)
resize buffer for accumulated currents
std::vector< float > mDCAz_C_WeightedMean
integrated 1D DCAz for C-side weighted mean in phi/tgl slices
std::vector< float > mTPCChi2A
Chi2 of TPC tracks.
std::vector< float > mDCAz_C_RMS
integrated 1D DCAz for C-side RMS in phi/tgl slices
int getIndexqPt(const int iqPt, int slice=0) const
std::vector< float > mDCAr_A_RMS
integrated 1D DCAr for A-side RMS in phi/tgl slices
unsigned char mNBinsqPt
number of qPt bins
int getIndexMult(const int iMult, int slice=0) const
std::vector< float > mDCAz_A_WeightedMean
integrated 1D DCAz for A-side weighted mean in phi/tgl slices
std::vector< float > mDCAr_C_WeightedMean
integrated 1D DCAr for C-side median in phi/tgl slices
std::vector< float > mLogdEdx_C_OROC2_Median
log(dedxOROC2 / dEdx) - C-side
void resize(const unsigned int nTotal)
std::vector< float > mLogdEdx_A_IROC_RMS
log(dedxIROC / dEdx) - A-side
std::vector< float > mLogdEdx_C_IROC_Median
log(dedxIROC / dEdx) - C-side
std::vector< float > mLogdEdx_A_OROC3_Median
log(dedxOROC3 / dEdx) - A-side
std::vector< float > mLogdEdx_A_OROC2_RMS
log(dedxOROC2 / dEdx) - A-side
std::vector< float > mLogdEdx_A_OROC3_RMS
log(dedxOROC3 / dEdx) - A-side
ClassDefNV(TimeSeriesdEdx, 1)
std::vector< float > mLogdEdx_A_Median
log(dEdx_exp(pion)/dEdx) - A-side
std::vector< float > mLogdEdx_C_OROC1_Median
log(dedxOROC1 / dEdx) - C-side
std::vector< float > mLogdEdx_C_Median
log(dEdx_exp(pion)/dEdx) - C-side
std::vector< float > mLogdEdx_C_RMS
log(dEdx_exp(pion)/dEdx) - C-side
std::vector< float > mLogdEdx_A_OROC1_RMS
log(dedxOROC1 / dEdx) - A-side
std::vector< float > mLogdEdx_C_OROC3_Median
log(dedxOROC3 / dEdx) - C-side
std::vector< float > mLogdEdx_C_OROC1_RMS
log(dedxOROC1 / dEdx) - C-side
std::vector< float > mLogdEdx_C_IROC_RMS
log(dedxIROC / dEdx) - C-side
std::vector< float > mLogdEdx_A_IROC_Median
log(dedxIROC / dEdx) - A-side
std::vector< float > mLogdEdx_A_OROC2_Median
log(dedxOROC2 / dEdx) - A-side
std::vector< float > mLogdEdx_C_OROC2_RMS
log(dedxOROC2 / dEdx) - C-side
std::vector< float > mLogdEdx_A_OROC1_Median
log(dedxOROC1 / dEdx) - A-side
std::vector< float > mLogdEdx_A_RMS
log(dEdx_exp(pion)/dEdx) - A-side
std::vector< float > mLogdEdx_C_OROC3_RMS
log(dedxOROC3 / dEdx) - C-side