Project
Loading...
Searching...
No Matches
TPCTrackStudy.cxx
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
12#include <vector>
13#include <TStopwatch.h>
30#include "GPUO2InterfaceRefit.h"
34
35namespace o2::trackstudy
36{
37
38using namespace o2::framework;
41
47
49
50class TPCTrackStudySpec : public Task
51{
52 public:
53 TPCTrackStudySpec(std::shared_ptr<DataRequest> dr, std::shared_ptr<o2::base::GRPGeomRequest> gr, const o2::tpc::CorrectionMapsLoaderGloOpts& sclOpts, GTrackID::mask_t src, bool useMC)
54 : mDataRequest(dr), mGGCCDBRequest(gr), mTracksSrc(src), mUseMC(useMC)
55 {
56 mTPCCorrMapsLoader.setLumiScaleType(sclOpts.lumiType);
57 mTPCCorrMapsLoader.setLumiScaleMode(sclOpts.lumiMode);
58 mTPCCorrMapsLoader.setCheckCTPIDCConsistency(sclOpts.checkCTPIDCconsistency);
59 }
60 ~TPCTrackStudySpec() final = default;
61 void init(InitContext& ic) final;
62 void run(ProcessingContext& pc) final;
63 void endOfStream(EndOfStreamContext& ec) final;
64 void finaliseCCDB(ConcreteDataMatcher& matcher, void* obj) final;
65 void process(o2::globaltracking::RecoContainer& recoData);
66
67 private:
68 void updateTimeDependentParams(ProcessingContext& pc);
69 std::shared_ptr<DataRequest> mDataRequest;
70 std::shared_ptr<o2::base::GRPGeomRequest> mGGCCDBRequest;
71 o2::tpc::VDriftHelper mTPCVDriftHelper{};
72 o2::tpc::CorrectionMapsLoader mTPCCorrMapsLoader{};
73 bool mUseMC{false};
74 bool mUseGPUModel{false};
75 float mXRef = 0.;
76 int mNMoves = 6;
77 int mTFStart = 0;
78 int mTFEnd = 999999999;
79 int mTFCount = -1;
80 bool mUseR = false;
81 bool mRepRef = false;
82 std::unique_ptr<o2::utils::TreeStreamRedirector> mDBGOut;
83 std::unique_ptr<o2::utils::TreeStreamRedirector> mDBGOutCl;
84 float mITSROFrameLengthMUS = 0.;
85 GTrackID::mask_t mTracksSrc{};
86 o2::steer::MCKinematicsReader mcReader; // reader of MC information
87 //
88 // TPC data
89 gsl::span<const o2::tpc::TPCClRefElem> mTPCTrackClusIdx;
90 gsl::span<const o2::tpc::TrackTPC> mTPCTracksArray;
91 gsl::span<const unsigned char> mTPCRefitterShMap;
92 gsl::span<const unsigned int> mTPCRefitterOccMap;
93 const o2::tpc::ClusterNativeAccess* mTPCClusterIdxStruct = nullptr;
94 gsl::span<const o2::MCCompLabel> mTPCTrkLabels;
95 std::unique_ptr<o2::gpu::GPUO2InterfaceRefit> mTPCRefitter;
96};
97
99{
101 mXRef = ic.options().get<float>("target-x");
102 mNMoves = std::max(2, ic.options().get<int>("n-moves"));
103 mUseR = ic.options().get<bool>("use-r-as-x");
104 mRepRef = ic.options().get<bool>("repeat-ini-ref");
105 mUseGPUModel = ic.options().get<bool>("use-gpu-fitter");
106 mTFStart = ic.options().get<int>("tf-start");
107 mTFEnd = ic.options().get<int>("tf-end");
108 if (mXRef < 0.) {
109 mXRef = 0.;
110 }
111 mTPCCorrMapsLoader.init(ic);
112 mDBGOut = std::make_unique<o2::utils::TreeStreamRedirector>("tpc-trackStudy.root", "recreate");
113 if (ic.options().get<bool>("dump-clusters")) {
114 mDBGOutCl = std::make_unique<o2::utils::TreeStreamRedirector>("tpc-trackStudy-cl.root", "recreate");
115 }
116}
117
119{
120 mTFCount++;
121 if (mTFCount < mTFStart || mTFCount > mTFEnd) {
122 LOGP(info, "Skipping TF {}", mTFCount);
123 return;
124 }
125
127 recoData.collectData(pc, *mDataRequest.get()); // select tracks of needed type, with minimal cuts, the real selected will be done in the vertexer
128 updateTimeDependentParams(pc); // Make sure this is called after recoData.collectData, which may load some conditions
129 process(recoData);
130
131 if (mTFCount > mTFEnd) {
132 LOGP(info, "Stopping processing after TF {}", mTFCount);
134 return;
135 }
136}
137
138void TPCTrackStudySpec::updateTimeDependentParams(ProcessingContext& pc)
139{
141 mTPCVDriftHelper.extractCCDBInputs(pc);
142 mTPCCorrMapsLoader.extractCCDBInputs(pc);
143 static bool initOnceDone = false;
144 if (!initOnceDone) { // this params need to be queried only once
145 initOnceDone = true;
146 // none at the moment
147 }
148 // we may have other params which need to be queried regularly
149 bool updateMaps = false;
150 if (mTPCCorrMapsLoader.isUpdated()) {
151 mTPCCorrMapsLoader.acknowledgeUpdate();
152 updateMaps = true;
153 }
154 if (mTPCVDriftHelper.isUpdated()) {
155 LOGP(info, "Updating TPC fast transform map with new VDrift factor of {} wrt reference {} and DriftTimeOffset correction {} wrt {} from source {}",
156 mTPCVDriftHelper.getVDriftObject().corrFact, mTPCVDriftHelper.getVDriftObject().refVDrift,
157 mTPCVDriftHelper.getVDriftObject().timeOffsetCorr, mTPCVDriftHelper.getVDriftObject().refTimeOffset,
158 mTPCVDriftHelper.getSourceName());
159 mTPCVDriftHelper.acknowledgeUpdate();
160 updateMaps = true;
161 }
162 if (updateMaps) {
163 mTPCCorrMapsLoader.updateVDrift(mTPCVDriftHelper.getVDriftObject().corrFact, mTPCVDriftHelper.getVDriftObject().refVDrift, mTPCVDriftHelper.getVDriftObject().getTimeOffset());
164 }
165}
166
168{
169 static long counter = -1;
170 auto prop = o2::base::Propagator::Instance();
171
172 mTPCTracksArray = recoData.getTPCTracks();
173 mTPCTrackClusIdx = recoData.getTPCTracksClusterRefs();
174 mTPCClusterIdxStruct = &recoData.inputsTPCclusters->clusterIndex;
175 mTPCRefitterShMap = recoData.clusterShMapTPC;
176 mTPCRefitterOccMap = recoData.occupancyMapTPC;
177
178 std::vector<o2::InteractionTimeRecord> intRecs;
179 if (mUseMC) { // extract MC tracks
180 const o2::steer::DigitizationContext* digCont = nullptr;
181 if (!mcReader.initFromDigitContext("collisioncontext.root")) {
182 throw std::invalid_argument("initialization of MCKinematicsReader failed");
183 }
184 digCont = mcReader.getDigitizationContext();
185 intRecs = digCont->getEventRecords();
186 mTPCTrkLabels = recoData.getTPCTracksMCLabels();
187 }
188 if (mTPCTracksArray.size()) {
189 LOGP(info, "Found {} TPC tracks", mTPCTracksArray.size());
190 mTPCRefitter = std::make_unique<o2::gpu::GPUO2InterfaceRefit>(mTPCClusterIdxStruct, &mTPCCorrMapsLoader, prop->getNominalBz(), mTPCTrackClusIdx.data(), 0, mTPCRefitterShMap.data(), mTPCRefitterOccMap.data(), mTPCRefitterOccMap.size(), nullptr, o2::base::Propagator::Instance());
191 mTPCRefitter->setTrackReferenceX(900); // disable propagation after refit by setting reference to value > 500
192 }
193 float vdriftTB = mTPCVDriftHelper.getVDriftObject().getVDrift() * o2::tpc::ParameterElectronics::Instance().ZbinWidth; // VDrift expressed in cm/TimeBin
194 float tpcTBBias = mTPCVDriftHelper.getVDriftObject().getTimeOffset() / (8 * o2::constants::lhc::LHCBunchSpacingMUS);
195 std::vector<short> clSector, clRow;
196 std::vector<float> clX, clY, clZ;
197
198 auto dumpClusters = [this] {
199 static int tf = 0;
200 const auto* corrMap = this->mTPCCorrMapsLoader.getCorrMap();
201 for (int sector = 0; sector < 36; sector++) {
202 float alp = ((sector % 18) * 20 + 10) * TMath::DegToRad();
203 float sn = TMath::Sin(alp), cs = TMath::Cos(alp);
204 for (int row = 0; row < 152; row++) {
205 for (int ic = 0; ic < this->mTPCClusterIdxStruct->nClusters[sector][row]; ic++) {
206 const auto cl = this->mTPCClusterIdxStruct->clusters[sector][row][ic];
207 float x, y, z, xG, yG;
208 corrMap->TransformIdeal(sector, row, cl.getPad(), cl.getTime(), x, y, z, 0);
209 o2::math_utils::detail::rotateZ(x, y, xG, yG, sn, cs);
210 LOGP(debug, "tf:{} s:{} r:{} p:{} t:{} qm:{} qt:{} f:{} x:{} y:{} z:{}", tf, sector, row, cl.getPad(), cl.getTime(), cl.getQmax(), cl.getQtot(), cl.getFlags(), x, y, z);
211 (*mDBGOutCl) << "tpccl"
212 << "tf=" << tf << "sect=" << sector << "row=" << row << "pad=" << cl.getPad() << "time=" << cl.getTime() << "qmax=" << cl.getQmax() << "qtot=" << cl.getQtot()
213 << "sigT=" << cl.getSigmaTime() << "sigP=" << cl.getSigmaPad()
214 << "flags=" << cl.getFlags()
215 << "x=" << x << "y=" << y << "z=" << z << "xg=" << xG << "yg=" << yG
216 << "\n";
217 }
218 }
219 }
220 tf++;
221 };
222
223 if (mDBGOutCl) {
224 dumpClusters();
225 }
226
227 for (size_t itr = 0; itr < mTPCTracksArray.size(); itr++) {
228 auto tr = mTPCTracksArray[itr]; // create track copy
229 int side = 0;
230 if (tr.hasBothSidesClusters()) {
231 continue;
232 }
233 side = tr.hasASideClustersOnly() ? 1 : -1;
234 //=========================================================================
235 // create refitted copy
236 auto trackRefit = [itr, this](o2::track::TrackParCov& trc, float t) -> bool {
237 float chi2Out = 0;
238 int retVal = mUseGPUModel ? this->mTPCRefitter->RefitTrackAsGPU(trc, this->mTPCTracksArray[itr].getClusterRef(), t, &chi2Out, false, true) : this->mTPCRefitter->RefitTrackAsTrackParCov(trc, this->mTPCTracksArray[itr].getClusterRef(), t, &chi2Out, false, true);
239 if (retVal < 0) {
240 LOGP(warn, "Refit failed ({}) with time={}: track#{}[{}]", retVal, t, counter, trc.asString());
241 return false;
242 }
243 return true;
244 };
245
246 auto trackProp = [&tr, itr, prop, this](o2::track::TrackParCov& trc) -> bool {
247 if (!trc.rotate(tr.getAlpha())) {
248 LOGP(warn, "Rotation to original track alpha {} failed, track#{}[{}]", tr.getAlpha(), counter, trc.asString());
249 return false;
250 }
251 float xtgt = this->mXRef;
252 if (mUseR && !trc.getXatLabR(this->mXRef, xtgt, prop->getNominalBz(), o2::track::DirInward)) {
253 xtgt = 0;
254 return false;
255 }
256 if (!prop->PropagateToXBxByBz(trc, xtgt)) {
257 LOGP(warn, "Propagation to X={} failed, track#{}[{}]", xtgt, counter, trc.asString());
258 return false;
259 }
260 return true;
261 };
262
263 auto prepClus = [this, &tr, &clSector, &clRow, &clX, &clY, &clZ](float t) { // extract cluster info
264 clSector.clear();
265 clRow.clear();
266 clX.clear();
267 clY.clear();
268 clZ.clear();
269 int count = tr.getNClusters();
270 const o2::tpc::ClusterNative* cl = nullptr;
271 for (int ic = count; ic--;) {
272 uint8_t sector, row;
273 cl = &tr.getCluster(this->mTPCTrackClusIdx, ic, *this->mTPCClusterIdxStruct, sector, row);
274 clSector.push_back(sector);
275 clRow.push_back(row);
276 float x, y, z;
277 mTPCCorrMapsLoader.Transform(sector, row, cl->getPad(), cl->getTime(), x, y, z, t); // nominal time of the track
278 clX.push_back(x);
279 clY.push_back(y);
280 clZ.push_back(z);
281 }
282 };
283
284 //=========================================================================
285
286 auto trf = tr.getOuterParam(); // we refit inward original track
287 if (!trackRefit(trf, tr.getTime0()) || !trackProp(trf)) {
288 continue;
289 }
290
291 // propagate original track
292 if (!trackProp(tr)) {
293 continue;
294 }
295
296 prepClus(tr.getTime0()); // original clusters
297 counter++;
298 // store results
299 (*mDBGOut) << "tpcIni"
300 << "counter=" << counter
301 << "iniTrack=" << tr
302 << "iniTrackRef=" << trf
303 << "side=" << side
304 << "time=" << tr.getTime0()
305 << "clSector=" << clSector
306 << "clRow=" << clRow
307 << "clX=" << clX
308 << "clY=" << clY
309 << "clZ=" << clZ
310 << "\n";
311
312 float dz = 0;
313
314 while (mUseMC) { // impose MC time in TPC timebin and refit inward after resetted covariance
315 // extract MC truth
316 const o2::MCTrack* mcTrack = nullptr;
317 auto lbl = mTPCTrkLabels[itr];
318 if (!lbl.isValid() || !(mcTrack = mcReader.getTrack(lbl))) {
319 break;
320 }
321 long bc = intRecs[lbl.getEventID()].toLong(); // bunch crossing of the interaction
322 float bcTB = bc / 8. + tpcTBBias; // the same in TPC timebins, accounting for the TPC time bias
323 // create MC truth track in O2 format
324 std::array<float, 3> xyz{(float)mcTrack->GetStartVertexCoordinatesX(), (float)mcTrack->GetStartVertexCoordinatesY(), (float)mcTrack->GetStartVertexCoordinatesZ()},
325 pxyz{(float)mcTrack->GetStartVertexMomentumX(), (float)mcTrack->GetStartVertexMomentumY(), (float)mcTrack->GetStartVertexMomentumZ()};
326 TParticlePDG* pPDG = TDatabasePDG::Instance()->GetParticle(mcTrack->GetPdgCode());
327 if (!pPDG) {
328 break;
329 }
330 o2::track::TrackPar mctrO2(xyz, pxyz, TMath::Nint(pPDG->Charge() / 3), false);
331 //
332 // propagate it to the alpha/X of the reconstructed track
333 if (!mctrO2.rotate(tr.getAlpha()) || !prop->PropagateToXBxByBz(mctrO2, tr.getX())) {
334 break;
335 }
336 // now create a properly refitted track with correct time and distortions correction
337 {
338 auto trfm = tr.getOuterParam(); // we refit inward
339 // impose MC time in TPC timebin and refit inward after resetted covariance
340 if (!trackRefit(trfm, bcTB) || !trfm.rotate(tr.getAlpha()) || !prop->PropagateToXBxByBz(trfm, tr.getX())) {
341 LOGP(warn, "Failed to propagate MC-time refitted track#{} [{}] to X/alpha of original track [{}]", counter, trfm.asString(), tr.asString());
342 break;
343 }
344 // estimate Z shift in case of no-distortions
345 dz = (tr.getTime0() - bcTB) * vdriftTB;
346 if (tr.hasCSideClustersOnly()) {
347 dz = -dz;
348 }
349 //
350 prepClus(bcTB); // clusters for MC time
351 (*mDBGOut) << "tpcMC"
352 << "counter=" << counter
353 << "movTrackRef=" << trfm
354 << "mcTrack=" << mctrO2
355 << "side=" << side
356 << "imposedTB=" << bcTB
357 << "dz=" << dz
358 << "clX=" << clX
359 << "clY=" << clY
360 << "clZ=" << clZ
361 << "\n";
362 }
363 break;
364 }
365 // refit and store the same track for a few compatible times
366 float tmin = tr.getTime0() - tr.getDeltaTBwd();
367 float tmax = tr.getTime0() + tr.getDeltaTFwd();
368 for (int it = 0; it < mNMoves; it++) {
369 float tb = tmin + it * (tmax - tmin) / (mNMoves - 1);
370 auto trfm = tr.getOuterParam(); // we refit inward
371 // impose time in TPC timebin and refit inward after resetted covariance
372 if (!trackRefit(trfm, tb) || !trfm.rotate(tr.getAlpha()) || !prop->PropagateToXBxByBz(trfm, tr.getX())) {
373 LOGP(warn, "Failed to propagate time={} refitted track#{} [{}] to X/alpha of original track [{}]", tb, counter, trfm.asString(), tr.asString());
374 continue;
375 }
376 // estimate Z shift in case of no-distortions
377 dz = (tr.getTime0() - tb) * vdriftTB;
378 if (tr.hasCSideClustersOnly()) {
379 dz = -dz;
380 }
381 //
382 int mnm = mNMoves - 1;
383 prepClus(tb); // clusters for MC time
384 (*mDBGOut) << "tpcMov"
385 << "counter=" << counter
386 << "copy=" << it
387 << "maxCopy=" << mnm
388 << "movTrackRef=" << trfm;
389 if (mRepRef) {
390 (*mDBGOut) << "tpcMov"
391 << "iniTrackRef=" << trf << "time=" << tr.getTime0();
392 }
393 (*mDBGOut) << "tpcMov"
394 << "side=" << side
395 << "imposedTB=" << tb
396 << "dz=" << dz
397 << "clX=" << clX
398 << "clY=" << clY
399 << "clZ=" << clZ
400 << "\n";
401 }
402 }
403}
404
406{
407 mDBGOut.reset();
408 mDBGOutCl.reset();
409}
410
412{
414 return;
415 }
416 if (mTPCVDriftHelper.accountCCDBInputs(matcher, obj)) {
417 return;
418 }
419 if (mTPCCorrMapsLoader.accountCCDBInputs(matcher, obj)) {
420 return;
421 }
422}
423
425{
426 std::vector<OutputSpec> outputs;
427 Options opts{
428 {"target-x", VariantType::Float, 70.f, {"Try to propagate to this radius"}},
429 {"n-moves", VariantType::Int, 6, {"Number of moves in allow range"}},
430 {"dump-clusters", VariantType::Bool, false, {"dump clusters"}},
431 {"tf-start", VariantType::Int, 0, {"1st TF to process"}},
432 {"tf-end", VariantType::Int, 999999999, {"last TF to process"}},
433 {"use-gpu-fitter", VariantType::Bool, false, {"use GPU track model for refit instead of TrackParCov"}},
434 {"repeat-ini-ref", VariantType::Bool, false, {"store ini-refit param with every moved track"}},
435 {"use-r-as-x", VariantType::Bool, false, {"Use radius instead of target sector X"}}};
436 auto dataRequest = std::make_shared<DataRequest>();
437
438 dataRequest->requestTracks(srcTracks, useMC);
439 dataRequest->requestClusters(srcClusters, useMC);
440 auto ggRequest = std::make_shared<o2::base::GRPGeomRequest>(false, // orbitResetTime
441 true, // GRPECS=true
442 false, // GRPLHCIF
443 true, // GRPMagField
444 true, // askMatLUT
446 dataRequest->inputs,
447 true);
448 o2::tpc::VDriftHelper::requestCCDBInputs(dataRequest->inputs);
449 o2::tpc::CorrectionMapsLoader::requestCCDBInputs(dataRequest->inputs, opts, sclOpts);
450
451 return DataProcessorSpec{
452 "tpc-track-study",
453 dataRequest->inputs,
454 outputs,
455 AlgorithmSpec{adaptFromTask<TPCTrackStudySpec>(dataRequest, ggRequest, sclOpts, srcTracks, useMC)},
456 opts};
457}
458
459} // namespace o2::trackstudy
Helper class to access load maps from CCDB.
Wrapper container for different reconstructed object types.
Definition of the GeometryManager class.
uint64_t bc
Definition RawEventData.h:5
int32_t retVal
Helper for geometry and GRP related CCDB requests.
Global index for barrel track: provides provenance (detectors combination), index in respective array...
Utility functions for MC particles.
Definition of the Names Generator class.
Definition of the parameter class for the detector electronics.
uint32_t side
Definition RawData.h:0
Wrapper container for different reconstructed object types.
Helper class to extract VDrift from different sources.
std::ostringstream debug
Double_t GetStartVertexMomentumZ() const
Definition MCTrack.h:81
Double_t GetStartVertexMomentumX() const
Definition MCTrack.h:79
Double_t GetStartVertexCoordinatesY() const
Definition MCTrack.h:83
Double_t GetStartVertexCoordinatesZ() const
Definition MCTrack.h:84
Double_t GetStartVertexMomentumY() const
Definition MCTrack.h:80
Double_t GetStartVertexCoordinatesX() const
Definition MCTrack.h:82
Int_t GetPdgCode() const
Accessors.
Definition MCTrack.h:72
void checkUpdates(o2::framework::ProcessingContext &pc)
static GRPGeomHelper & instance()
void setRequest(std::shared_ptr< GRPGeomRequest > req)
GPUd() value_type estimateLTFast(o2 static GPUd() float estimateLTIncrement(const o2 PropagatorImpl * Instance(bool uninitialized=false)
Definition Propagator.h:147
Static class with identifiers, bitmasks and names for ALICE detectors.
Definition DetID.h:58
ConfigParamRegistry const & options()
Definition InitContext.h:33
ServiceRegistryRef services()
The services registry associated with this processing context.
std::vector< o2::InteractionTimeRecord > & getEventRecords(bool withQED=false)
bool initFromDigitContext(std::string_view filename)
DigitizationContext const * getDigitizationContext() const
MCTrack const * getTrack(o2::MCCompLabel const &) const
void extractCCDBInputs(o2::framework::ProcessingContext &pc)
void updateVDrift(float vdriftCorr, float vdrifRef, float driftTimeOffset=0)
bool accountCCDBInputs(const o2::framework::ConcreteDataMatcher &matcher, void *obj)
static void requestCCDBInputs(std::vector< o2::framework::InputSpec > &inputs, std::vector< o2::framework::ConfigParamSpec > &options, const CorrectionMapsLoaderGloOpts &gloOpts)
void init(o2::framework::InitContext &ic)
static void requestCCDBInputs(std::vector< o2::framework::InputSpec > &inputs, bool laser=true, bool itstpcTgl=true)
void extractCCDBInputs(o2::framework::ProcessingContext &pc, bool laser=true, bool itstpcTgl=true)
const VDriftCorrFact & getVDriftObject() const
bool accountCCDBInputs(const o2::framework::ConcreteDataMatcher &matcher, void *obj)
static std::string_view getSourceName(Source s)
bool isUpdated() const
void endOfStream(EndOfStreamContext &ec) final
This is invoked whenever we have an EndOfStream event.
TPCTrackStudySpec(std::shared_ptr< DataRequest > dr, std::shared_ptr< o2::base::GRPGeomRequest > gr, const o2::tpc::CorrectionMapsLoaderGloOpts &sclOpts, GTrackID::mask_t src, bool useMC)
void init(InitContext &ic) final
void run(ProcessingContext &pc) final
void process(o2::globaltracking::RecoContainer &recoData)
void finaliseCCDB(ConcreteDataMatcher &matcher, void *obj) final
GLint GLenum GLint x
Definition glcorearb.h:403
GLenum src
Definition glcorearb.h:1767
GLint GLsizei count
Definition glcorearb.h:399
GLint y
Definition glcorearb.h:270
GLdouble GLdouble GLdouble z
Definition glcorearb.h:843
GLuint counter
Definition glcorearb.h:3987
constexpr double LHCBunchSpacingMUS
Defining PrimaryVertex explicitly as messageable.
Definition TFIDInfo.h:20
std::vector< ConfigParamSpec > Options
std::tuple< T, T > rotateZ(T xL, T yL, T snAlp, T csAlp)
detail::Bracket< float > Bracketf_t
Definition Primitive2D.h:40
TrackParCovF TrackParCov
Definition Track.h:33
o2::math_utils::Bracketf_t TBracket
o2::framework::DataProcessorSpec getTPCTrackStudySpec(o2::dataformats::GlobalTrackID::mask_t srcTracks, o2::dataformats::GlobalTrackID::mask_t srcClus, bool useMC, const o2::tpc::CorrectionMapsLoaderGloOpts &sclOpts)
create a processor spec
o2::dataformats::VtxTrackRef V2TRef
o2::dataformats::VtxTrackIndex VTIndex
a couple of static helper functions to create timestamp values for CCDB queries or override obsolete ...
Defining DataPointCompositeObject explicitly as copiable.
std::unique_ptr< GPUReconstructionTimeframe > tf
gsl::span< const unsigned char > clusterShMapTPC
externally set TPC clusters sharing map
void collectData(o2::framework::ProcessingContext &pc, const DataRequest &request)
std::unique_ptr< o2::tpc::internal::getWorkflowTPCInput_ret > inputsTPCclusters
gsl::span< const unsigned int > occupancyMapTPC
externally set TPC clusters occupancy map
unsigned int nClusters[constants::MAXSECTOR][constants::MAXGLOBALPADROW]
const ClusterNative * clusters[constants::MAXSECTOR][constants::MAXGLOBALPADROW]
int lumiType
what estimator to used for corrections scaling: 0: no scaling, 1: CTP, 2: IDC
int lumiMode
what corrections method to use: 0: classical scaling, 1: Using of the derivative map,...
float refTimeOffset
additive time offset reference (\mus)
float refVDrift
reference vdrift for which factor was extracted
float getTimeOffset() const
float timeOffsetCorr
additive time offset correction (\mus)
float corrFact
drift velocity correction factor (multiplicative)
std::vector< int > row