Project
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
TPCFourierTransformAggregatorSpec.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 O2_TPCFOURIERTRANSFORMAGGREGATORSPEC_H
18#define O2_TPCFOURIERTRANSFORMAGGREGATORSPEC_H
19
20#include <vector>
21#include <fmt/format.h>
22#include "Framework/Task.h"
24#include "Framework/Logger.h"
26#include "Headers/DataHeader.h"
27#include "CCDB/CcdbApi.h"
31#include "TPCBase/CRU.h"
34
35using namespace o2::framework;
37using namespace o2::tpc;
38
39namespace o2::tpc
40{
41
43 std::array<long, 2> timestamp;
44 std::array<std::vector<float>, 2> idc1;
45};
46
48 std::unordered_map<long, TPCScalerProc> idcs;
49};
50
52{
53 public:
54 // Fourier type
56
57 TPCFourierTransformAggregatorSpec(const unsigned int nFourierCoefficientsStore, const unsigned int rangeIDC, const bool senddebug = false, const bool processSACs = false, const int inputLanes = 1)
58 : mIDCFourierTransform{IDCFType(rangeIDC, nFourierCoefficientsStore), IDCFType(rangeIDC, nFourierCoefficientsStore)}, mSendOutDebug{senddebug}, mProcessSACs{processSACs}, mInputLanes{inputLanes} {};
59
61 {
62 mDumpFFT = ic.options().get<bool>("dump-coefficients-agg");
63 mIntervalsSACs = ic.options().get<int>("intervalsSACs");
64 mLengthIDCScalerSeconds = ic.options().get<float>("tpcScalerLengthS");
65 mDisableScaler = ic.options().get<bool>("disable-scaler");
66 mEnableFFTCCDB = ic.options().get<bool>("enable-fft-CCDB");
67 resizeBuffer(mInputLanes);
68 }
69
71 {
72 const int lane = pc.inputs().get<int>("lane");
73 if (lane >= mInputLanes) {
74 LOGP(error, "Received data from lane {} which is >= than the specified number of expected lanes of {}!", lane, mInputLanes);
75 return;
76 }
77
78 const auto tsTmp = pc.inputs().get<std::vector<long>>("tsccdb");
79 if (tsTmp.front() == 0) {
80 LOGP(warning, "Received dummy data with empty timestamp");
81 return;
82 }
83 mCCDBBuffer[lane] = tsTmp;
84 if (mProcessedTimeStamp > mCCDBBuffer[lane].front()) {
85 LOGP(warning, "Already received data from a later time stamp {} then the currently received time stamp {}! (This might not be an issue)", mProcessedTimeStamp, mCCDBBuffer[lane].front());
86 } else {
87 mProcessedTimeStamp = mCCDBBuffer[lane].front();
88 }
89
90 if (!mProcessSACs) {
91 mIntervalsBuffer[lane] = pc.inputs().get<std::vector<unsigned int>>("intervals");
92 }
93
94 for (auto& ref : InputRecordWalker(pc.inputs(), mFilter[mProcessSACs])) {
95 auto const* dataHeader = o2::framework::DataRefUtils::getHeader<o2::header::DataHeader*>(ref);
96 const int side = dataHeader->subSpecification;
97 mIDCOneBuffer[lane][side].mIDCOne = pc.inputs().get<std::vector<float>>(ref);
98 LOGP(info, "Received {} 1D-IDCs for side {}", mIDCOneBuffer[lane][side].mIDCOne.size(), side);
99
100 if (mProcessSACs && mIntervalsBuffer[lane].empty()) {
101 const auto nValues = mIDCOneBuffer[lane][side].mIDCOne.size();
102 const int nIntervals = nValues / mIntervalsSACs;
103 const int nFirstInterval = nValues % mIntervalsSACs;
104 if (nFirstInterval == 0) {
105 mIntervalsBuffer[lane] = std::vector<unsigned int>(nIntervals, mIntervalsSACs);
106 } else {
107 mIntervalsBuffer[lane] = std::vector<unsigned int>(nIntervals + 1, mIntervalsSACs);
108 mIntervalsBuffer[lane].front() = nFirstInterval;
109 }
110 }
111 }
112
113 // buffer IDCs for TPC scaler
114 if (!mProcessSACs && !mDisableScaler) {
115 const long startTS = mCCDBBuffer[lane].front();
116 TPCScalerProc& scaler = mTPCScalerCont.idcs[startTS];
117 scaler.timestamp[0] = startTS;
118 scaler.timestamp[1] = mCCDBBuffer[lane].back();
119 for (auto& ref : InputRecordWalker(pc.inputs(), mFilterI0)) {
120 auto const* dataHeader = o2::framework::DataRefUtils::getHeader<o2::header::DataHeader*>(ref);
121 const int side = dataHeader->subSpecification;
122 const float idc0mean = pc.inputs().get<float>(ref);
123 LOGP(info, "Received {} IDC0 mean for side {}", idc0mean, side);
124 scaler.idc1[side] = mIDCOneBuffer[lane][side].mIDCOne;
125 auto& vecIDC = scaler.idc1[side];
126
127 // normalize 1D-IDCs to unity as it should be
128 if (vecIDC.size() > 0) {
129 const float mean = std::reduce(vecIDC.begin(), vecIDC.end()) / static_cast<float>(vecIDC.size());
130 LOGP(info, "normalizing by {}", mean);
131 if (std::abs(mean) > 0.001) {
132 std::transform(vecIDC.begin(), vecIDC.end(), vecIDC.begin(), [&mean](auto val) { return val / mean; });
133 }
134 }
135
136 // scale IDC1 with IDC0Mean
137 std::transform(scaler.idc1[side].begin(), scaler.idc1[side].end(), scaler.idc1[side].begin(), [&idc0mean](auto idc) { return idc * idc0mean; });
138 }
139 // check if A- and C-side has the same length!
140 const int lenA = scaler.idc1[0].size();
141 const int lenC = scaler.idc1[1].size();
142 if (lenA != lenC) {
143 // This should never happen
144 LOGP(warning, "Received IDCs have different length! A-side length: {} and C-side length: {}", lenA, lenC);
145 // add dummy to shorter vector
146 const int maxLen = std::max(lenA, lenC);
147 scaler.idc1[0].resize(maxLen);
148 scaler.idc1[1].resize(maxLen);
149 }
150
152 makeTPCScaler(pc.outputs(), false);
153 }
154
155 FourierCoeffSAC coeffSAC;
156 if (lane == mExpectedInputLane) {
157 const int nSides = mIDCOneBuffer[lane][Side::A].mIDCOne.empty() + mIDCOneBuffer[lane][Side::C].mIDCOne.empty();
158 // int iProcessLane = lane;
159 for (int iProcessLaneTmp = 0; iProcessLaneTmp < mInputLanes; ++iProcessLaneTmp) {
160 const int nSidesCurrLane = mIDCOneBuffer[mExpectedInputLane][Side::A].mIDCOne.empty() + mIDCOneBuffer[mExpectedInputLane][Side::C].mIDCOne.empty();
161 if (nSidesCurrLane != nSides) {
162 break;
163 }
164
165 for (int iSide = 0; iSide < SIDES; ++iSide) {
166 const Side side = (iSide == 0) ? A : C;
167 if (mIDCOneBuffer[mExpectedInputLane][side].mIDCOne.empty()) {
168 continue;
169 }
170 LOGP(info, "Processing input lane: {} for Side: {}", mExpectedInputLane, iSide);
171
172 // perform fourier transform of 1D-IDCs
173 mIDCFourierTransform[side].setIDCs(std::move(mIDCOneBuffer[mExpectedInputLane][side]), mIntervalsBuffer[mExpectedInputLane]);
174 mIDCFourierTransform[side].calcFourierCoefficients(mIntervalsBuffer[mExpectedInputLane].size());
175
176 if (!mProcessSACs) {
177 if (mEnableFFTCCDB) {
178 o2::ccdb::CcdbObjectInfo ccdbInfo(CDBTypeMap.at(((side == 0) ? CDBType::CalIDCFourierA : CDBType::CalIDCFourierC)), std::string{}, std::string{}, std::map<std::string, std::string>{}, mCCDBBuffer[mExpectedInputLane].front(), mCCDBBuffer[mExpectedInputLane].back());
179 auto imageFFT = o2::ccdb::CcdbApi::createObjectImage(&mIDCFourierTransform[side].getFourierCoefficients(), &ccdbInfo);
180 LOGP(info, "Sending object {} / {} of size {} bytes, valid for {} : {} ", ccdbInfo.getPath(), ccdbInfo.getFileName(), imageFFT->size(), ccdbInfo.getStartValidityTimestamp(), ccdbInfo.getEndValidityTimestamp());
181 pc.outputs().snapshot(Output{o2::calibration::Utils::gDataOriginCDBPayload, getDataDescriptionCCDBFourier(), 0}, *imageFFT.get());
183 }
184 } else {
185 coeffSAC.mCoeff[side] = mIDCFourierTransform[side].getFourierCoefficients();
186 }
187
188 if (mDumpFFT) {
189 LOGP(info, "dumping FT to file");
190 mIDCFourierTransform[side].dumpToFile(fmt::format("FourierAGG_{:02}_side{}.root", processing_helpers::getCurrentTF(pc), (int)side).data());
191 }
192
193 if (mSendOutDebug) {
194 sendOutput(pc.outputs(), side);
195 }
196 }
197
198 if (mProcessSACs && mEnableFFTCCDB) {
199 o2::ccdb::CcdbObjectInfo ccdbInfo(CDBTypeMap.at(CDBType::CalSACFourier), std::string{}, std::string{}, std::map<std::string, std::string>{}, mCCDBBuffer[mExpectedInputLane].front(), mCCDBBuffer[mExpectedInputLane].back());
200 auto imageFFT = o2::ccdb::CcdbApi::createObjectImage(&coeffSAC, &ccdbInfo);
201 LOGP(info, "Sending object {} / {} of size {} bytes, valid for {} : {} ", ccdbInfo.getPath(), ccdbInfo.getFileName(), imageFFT->size(), ccdbInfo.getStartValidityTimestamp(), ccdbInfo.getEndValidityTimestamp());
202 pc.outputs().snapshot(Output{o2::calibration::Utils::gDataOriginCDBPayload, getDataDescriptionCCDBFourier(), 0}, *imageFFT.get());
204 }
205 mExpectedInputLane = ++mExpectedInputLane % mInputLanes;
206 }
207 }
208 }
209
211 {
212 if (!mDisableScaler && !mProcessSACs) {
213 makeTPCScaler(ec.outputs(), true);
214 }
215 ec.services().get<ControlService>().readyToQuit(QuitRequest::Me);
216 }
217
221
222 private:
223 std::array<IDCFType, SIDES> mIDCFourierTransform{};
224 const bool mSendOutDebug{false};
225 const bool mProcessSACs{false};
226 const int mInputLanes{1};
227 bool mDumpFFT{false};
228 uint64_t mProcessedTimeStamp{0};
229 std::vector<std::vector<long>> mCCDBBuffer{};
230 std::vector<std::vector<unsigned int>> mIntervalsBuffer{};
231 std::vector<std::array<o2::tpc::IDCOne, SIDES>> mIDCOneBuffer{};
232 unsigned int mIntervalsSACs{12};
233 int mExpectedInputLane{0};
234 TPCScalerProcContainer mTPCScalerCont;
235 float mLengthIDCScalerSeconds = 300;
236 long mIDCSCalerEndTSLast = 0;
237 o2::tpc::TPCScaler mScalerLast;
238 bool mDisableScaler{false};
239 bool mEnableFFTCCDB{false};
240 int mRun{};
241 const std::array<std::vector<InputSpec>, 2> mFilter = {std::vector<InputSpec>{{"idcone", ConcreteDataTypeMatcher{o2::header::gDataOriginTPC, TPCFactorizeIDCSpec::getDataDescriptionIDC1()}, Lifetime::Sporadic}},
242 std::vector<InputSpec>{{"sacone", ConcreteDataTypeMatcher{o2::header::gDataOriginTPC, TPCFactorizeSACSpec::getDataDescriptionSAC1()}, Lifetime::Sporadic}}};
243 const std::vector<InputSpec> mFilterI0 = std::vector<InputSpec>{{"idczeromean", ConcreteDataTypeMatcher{o2::header::gDataOriginTPC, TPCFactorizeIDCSpec::getDataDescriptionIDC0Mean()}, Lifetime::Sporadic}};
244
245 void sendOutput(DataAllocator& output, const int side)
246 {
247 output.snapshot(Output{gDataOriginTPC, TPCFourierTransformAggregatorSpec::getDataDescriptionFourier()}, mIDCFourierTransform[side].getFourierCoefficients());
248 }
249
250 void resizeBuffer(const int expectedLanes)
251 {
252 mCCDBBuffer.resize(expectedLanes);
253 mIntervalsBuffer.resize(expectedLanes);
254 mIDCOneBuffer.resize(expectedLanes);
255 }
256
257 void makeTPCScaler(DataAllocator& output, const bool eos)
258 {
259 LOGP(info, "Making TPC scalers");
260 if (mTPCScalerCont.idcs.empty()) {
261 LOGP(warning, "No IDCs received for TPC scaler creation");
262 return;
263 }
264
265 // check if IDC scalers can be created - check length of continous received IDCs
266 std::vector<std::pair<long, long>> times;
267 times.reserve(mTPCScalerCont.idcs.size());
268 for (const auto& idc : mTPCScalerCont.idcs) {
269 times.emplace_back(idc.second.timestamp[0], idc.second.timestamp[1]);
270 }
271
272 // sort received times of the IDCs
273 std::sort(times.begin(), times.end());
274
275 // loop over times and make checks
276 const int checkGapp = 10;
277 // if the diff between end of data[i] and start of data[i+1] is smaller than this, the data is contigous
278 long timesDuration = (times.front().second - times.front().first);
279
280 // make check and store lastValid index in case IDC scalers can be created
281 int lastValidIdx = -1;
282 for (int i = 1; i < times.size(); ++i) {
283 // check time diff between start of current IDCs and end time of last IDCs
284 const auto deltaTime = times[i].first - times[i - 1].second;
285 // check if IDCs are contigous
286 if (deltaTime > (timesDuration / checkGapp)) {
287 // check if the gap is very large - in this case the gapp might be lost, so just write out the TPC scaler until the gap
288 if (deltaTime > (checkGapp * timesDuration)) {
289 lastValidIdx = i - 1;
290 }
291 LOGP(info, "Breaking as big gap between IDCs of {} detected", deltaTime);
292 break;
293 }
294
295 // check if time length is >= than mLengthIDCScalerSeconds
296 if ((times[i].first - times.front().first) / 1000 >= mLengthIDCScalerSeconds) {
297 lastValidIdx = i;
298 }
299 }
300
301 LOGP(info, "Creating IDC scalers with {} IDC objects", lastValidIdx);
302
303 if (eos) {
304 // in case of eos write out everything
305 lastValidIdx = times.empty() ? -1 : times.size() - 1;
306 LOGP(info, "End of stream detected: Creating IDC scalers with {} IDC objects", lastValidIdx);
307 }
308
309 // create IDC scaler in case index is valid
310 if (lastValidIdx >= 0) {
311 o2::tpc::TPCScaler scaler;
312 scaler.setIonDriftTimeMS(170);
313 scaler.setRun(mRun);
314 scaler.setStartTimeStampMS(times.front().first);
315 const auto idcIntegrationTime = 12 /*12 orbits integration interval per IDC*/ * o2::constants::lhc::LHCOrbitMUS / 1000;
316 scaler.setIntegrationTimeMS(idcIntegrationTime);
317
318 std::vector<float> idc1A;
319 std::vector<float> idc1C;
320 long idc1ASize = 0;
321 long idc1CSize = 0;
322
323 // in case already one object is stored add internal overlap
324 if (mIDCSCalerEndTSLast != 0) {
325 const int nOverlap = 500;
326 idc1ASize += nOverlap;
327 idc1CSize += nOverlap;
328 const auto& scalerALast = mScalerLast.getScalers(o2::tpc::Side::A);
329 const auto& scalerCLast = mScalerLast.getScalers(o2::tpc::Side::C);
330 if (scalerALast.size() > nOverlap) {
331 idc1A.insert(idc1A.end(), scalerALast.end() - nOverlap, scalerALast.end());
332 idc1C.insert(idc1C.end(), scalerCLast.end() - nOverlap, scalerCLast.end());
333 // adjust start time
334 scaler.setStartTimeStampMS(scaler.getStartTimeStampMS() - nOverlap * idcIntegrationTime);
335 }
336 } else {
337 // store end timestamp as start time stamp for first object for correct time stamp in CCDB
338 mIDCSCalerEndTSLast = scaler.getStartTimeStampMS();
339 }
340
341 for (int iter = 0; iter < 2; ++iter) {
342 if (iter == 1) {
343 idc1A.reserve(idc1ASize);
344 idc1C.reserve(idc1CSize);
345 }
346 for (int i = 0; i <= lastValidIdx; ++i) {
347 const auto& time = times[i];
348 const auto& idc = mTPCScalerCont.idcs[time.first];
349 if (iter == 0) {
350 idc1ASize += idc.idc1[0].size();
351 idc1CSize += idc.idc1[1].size();
352 } else {
353 idc1A.insert(idc1A.end(), idc.idc1[0].begin(), idc.idc1[0].end());
354 idc1C.insert(idc1C.end(), idc.idc1[1].begin(), idc.idc1[1].end());
355 // in case of eos check if the IDCs are contigous and add dummy values!
356 if (eos && (i < lastValidIdx)) {
357 const float deltaTime = times[i + 1].first - time.second;
358 // if delta time is too large add dummy values
359 if (deltaTime > (timesDuration / checkGapp)) {
360 int nDummyValues = deltaTime / idcIntegrationTime + 0.5;
361 // restrict dummy values
362 const int nMaxDummyValues = checkGapp * timesDuration / idcIntegrationTime;
363 if (nDummyValues > nMaxDummyValues) {
364 nDummyValues = nMaxDummyValues;
365 }
366
367 // add dummy to A
368 if (idc.idc1[0].size() > 0) {
369 float meanA = std::reduce(idc.idc1[0].begin(), idc.idc1[0].end()) / static_cast<float>(idc.idc1[0].size());
370 idc1A.insert(idc1A.end(), nDummyValues, meanA);
371 }
372
373 if (idc.idc1[1].size() > 0) {
374 // add dummy to C
375 float meanC = std::reduce(idc.idc1[1].begin(), idc.idc1[1].end()) / static_cast<float>(idc.idc1[1].size());
376 idc1C.insert(idc1C.end(), nDummyValues, meanC);
377 }
378 }
379 }
380 mTPCScalerCont.idcs.erase(time.first);
381 }
382 }
383 }
384 scaler.setScaler(idc1A, o2::tpc::Side::A);
385 scaler.setScaler(idc1C, o2::tpc::Side::C);
386
387 // store in CCDB
388 TTree tree("ccdb_object", "ccdb_object");
389 tree.Branch("TPCScaler", &scaler);
390 tree.Fill();
391
392 o2::ccdb::CcdbObjectInfo ccdbInfoIDC(CDBTypeMap.at(CDBType::CalScaler), std::string{}, std::string{}, std::map<std::string, std::string>{}, mIDCSCalerEndTSLast, scaler.getEndTimeStampMS(o2::tpc::Side::A));
393 auto imageIDC = o2::ccdb::CcdbApi::createObjectImage(&tree, &ccdbInfoIDC);
394 LOGP(info, "Sending object {} / {} of size {} bytes, valid for {} : {} ", ccdbInfoIDC.getPath(), ccdbInfoIDC.getFileName(), imageIDC->size(), ccdbInfoIDC.getStartValidityTimestamp(), ccdbInfoIDC.getEndValidityTimestamp());
397
398 // store end timestamp
399 mIDCSCalerEndTSLast = scaler.getEndTimeStampMS(o2::tpc::Side::A);
400
401 // for debugging
402 if (mDumpFFT) {
403 static int countwrite = 0;
404 scaler.dumpToFile(fmt::format("TPCScaler_snapshot_{}.root", countwrite++).data(), "ccdb_object");
405 }
406
407 // buffer current scaler object
408 mScalerLast = std::move(scaler);
409 }
410 }
411};
412DataProcessorSpec getTPCFourierTransformAggregatorSpec(const unsigned int rangeIDC, const unsigned int nFourierCoefficientsStore, const bool senddebug, const bool processSACs, const int inputLanes)
413{
414 std::vector<OutputSpec> outputSpecs;
419
420 if (senddebug) {
421 outputSpecs.emplace_back(ConcreteDataTypeMatcher{gDataOriginTPC, TPCFourierTransformAggregatorSpec::getDataDescriptionFourier()}, Lifetime::Sporadic);
422 }
423
424 std::vector<InputSpec> inputSpecs;
425 if (!processSACs) {
426 inputSpecs.emplace_back(InputSpec{"idczeromean", ConcreteDataTypeMatcher{gDataOriginTPC, TPCFactorizeIDCSpec::getDataDescriptionIDC0Mean()}, Lifetime::Sporadic});
427 inputSpecs.emplace_back(InputSpec{"idcone", ConcreteDataTypeMatcher{gDataOriginTPC, TPCFactorizeIDCSpec::getDataDescriptionIDC1()}, Lifetime::Sporadic});
428 inputSpecs.emplace_back(InputSpec{"tsccdb", gDataOriginTPC, TPCFactorizeIDCSpec::getDataDescriptionTimeStamp(), Lifetime::Sporadic});
429 inputSpecs.emplace_back(InputSpec{"intervals", gDataOriginTPC, TPCFactorizeIDCSpec::getDataDescriptionIntervals(), Lifetime::Sporadic});
430 inputSpecs.emplace_back(InputSpec{"lane", gDataOriginTPC, TPCFactorizeIDCSpec::getDataDescriptionLane(), Lifetime::Sporadic});
431 } else {
432 inputSpecs.emplace_back(InputSpec{"sacone", ConcreteDataTypeMatcher{gDataOriginTPC, TPCFactorizeSACSpec::getDataDescriptionSAC1()}, Lifetime::Sporadic});
433 inputSpecs.emplace_back(InputSpec{"tsccdb", gDataOriginTPC, TPCFactorizeSACSpec::getDataDescriptionTimeStamp(), Lifetime::Sporadic});
434 inputSpecs.emplace_back(InputSpec{"lane", gDataOriginTPC, TPCFactorizeSACSpec::getDataDescriptionLane(), Lifetime::Sporadic});
435 }
436
437 std::string processorName = "tpc-aggregator-ft";
438 if (processSACs) {
439 processorName = "tpc-aggregator-ft-sac";
440 }
441
442 return DataProcessorSpec{
443 processorName,
444 inputSpecs,
445 outputSpecs,
446 AlgorithmSpec{adaptFromTask<TPCFourierTransformAggregatorSpec>(nFourierCoefficientsStore, rangeIDC, senddebug, processSACs, inputLanes)},
447 Options{{"intervalsSACs", VariantType::Int, 11, {"Number of integration intervals which will be sampled for the fourier coefficients"}},
448 {"dump-coefficients-agg", VariantType::Bool, false, {"Dump fourier coefficients to file"}},
449 {"tpcScalerLengthS", VariantType::Float, 300.f, {"Length of the TPC scalers in seconds"}},
450 {"disable-scaler", VariantType::Bool, false, {"Disable creation of IDC scaler"}},
451 {"enable-fft-CCDB", VariantType::Bool, false, {"Enable writing of FFT coefficients to CCDB"}}}};
452}
453
454} // namespace o2::tpc
455
456#endif
std::vector< unsigned long > times
int16_t time
Definition RawEventData.h:4
int32_t i
class for calculating the fourier coefficients from 1D-IDCs
void output(const std::map< std::string, ChannelStat > &channels)
Definition rawdump.cxx:197
Definition of the Names Generator class.
uint32_t side
Definition RawData.h:0
TPC factorization of SACs.
Definition A.h:16
static std::unique_ptr< std::vector< char > > createObjectImage(const T *obj, CcdbObjectInfo *info=nullptr)
Definition CcdbApi.h:103
long getEndValidityTimestamp() const
const std::string & getPath() const
long getStartValidityTimestamp() const
const std::string & getFileName() const
A helper class to iteratate over all parts of all input routes.
static constexpr header::DataDescription getDataDescriptionTimeStamp()
static constexpr header::DataDescription getDataDescriptionLane()
static constexpr header::DataDescription getDataDescriptionIDC0Mean()
static constexpr header::DataDescription getDataDescriptionIDC1()
static constexpr header::DataDescription getDataDescriptionIntervals()
static constexpr header::DataDescription getDataDescriptionSAC1()
static constexpr header::DataDescription getDataDescriptionLane()
static constexpr header::DataDescription getDataDescriptionTimeStamp()
void init(o2::framework::InitContext &ic) final
void endOfStream(o2::framework::EndOfStreamContext &ec) final
This is invoked whenever we have an EndOfStream event.
static constexpr header::DataDescription getDataDescriptionCCDBTPCScaler()
static constexpr header::DataDescription getDataDescriptionFourier()
static constexpr header::DataDescription getDataDescriptionCCDBFourier()
TPCFourierTransformAggregatorSpec(const unsigned int nFourierCoefficientsStore, const unsigned int rangeIDC, const bool senddebug=false, const bool processSACs=false, const int inputLanes=1)
void run(o2::framework::ProcessingContext &pc) final
float getScalers(unsigned int idx, o2::tpc::Side side) const
Definition TPCScaler.h:97
void setScaler(const std::vector< float > &values, const o2::tpc::Side side)
Definition TPCScaler.h:55
void setIonDriftTimeMS(float ionDriftTimeMS)
Definition TPCScaler.h:58
double getStartTimeStampMS() const
Definition TPCScaler.h:112
double getEndTimeStampMS(o2::tpc::Side side) const
Definition TPCScaler.h:115
void dumpToFile(const char *file, const char *name)
Definition TPCScaler.cxx:25
void setStartTimeStampMS(double timeStampMS)
Definition TPCScaler.h:67
void setRun(int run)
Definition TPCScaler.h:61
void setIntegrationTimeMS(float integrationTimeMS)
Definition TPCScaler.h:70
GLsizeiptr size
Definition glcorearb.h:659
GLboolean * data
Definition glcorearb.h:298
GLuint GLfloat * val
Definition glcorearb.h:1582
GLint ref
Definition glcorearb.h:291
constexpr o2::header::DataOrigin gDataOriginTPC
Definition DataHeader.h:576
constexpr double LHCOrbitMUS
Defining PrimaryVertex explicitly as messageable.
Definition TFIDInfo.h:20
std::vector< ConfigParamSpec > Options
uint32_t getCurrentTF(o2::framework::ProcessingContext &pc)
uint64_t getRunNumber(o2::framework::ProcessingContext &pc)
Global TPC definitions and constants.
Definition SimTraits.h:167
const std::unordered_map< CDBType, const std::string > CDBTypeMap
Storage name in CCDB for each calibration and parameter type.
Definition CDBTypes.h:95
DataProcessorSpec getTPCFourierTransformAggregatorSpec(const unsigned int rangeIDC, const unsigned int nFourierCoefficientsStore, const bool senddebug, const bool processSACs, const int inputLanes)
constexpr unsigned char SIDES
Definition Defs.h:41
Side
TPC readout sidE.
Definition Defs.h:35
@ A
Definition Defs.h:35
@ C
Definition Defs.h:36
@ CalSACFourier
Fourier coefficients of CalSAC1.
@ CalScaler
Scaler from IDCs or combined estimator.
@ CalIDCFourierC
Fourier coefficients of CalIDC1.
@ CalIDCFourierA
Fourier coefficients of CalIDC1.
void empty(int)
static constexpr o2::header::DataOrigin gDataOriginCDBWrapper
Definition Utils.h:44
static constexpr o2::header::DataOrigin gDataOriginCDBPayload
Definition Utils.h:43
std::array< FourierCoeff, SIDES > mCoeff
std::unordered_map< long, TPCScalerProc > idcs
std::array< long, 2 > timestamp
start -> end timestamp
std::array< std::vector< float >, 2 > idc1
IDC1.
std::unique_ptr< TTree > tree((TTree *) flIn.Get(std::string(o2::base::NameConf::CTFTREENAME).c_str()))