Skip to content

Commit

Permalink
TPC QA and PID information in the AO2D (#12387)
Browse files Browse the repository at this point in the history 
* ATO-647 - adding trackqa table to the AO2D

* ATO-647 - adding trackqa table to the AO2D as in original pull request of David
#12118

* ATO-647 - adding structure for TrackQA  +all dEdx information

* ATO-647 - adding cursor for the tracksQA table

* ATO-647 - formally declaring trackQA table (not filled yet)

* ATO-647 - filling content of the dEdx and ClusterByteMask
* dEdx and Cluster byte mask compression very effective
* trackIndex - not compressed as data are not sorted

* ATO-647 - adding description string for columns, add DCAr and DCAz

* ATO-647 - adding configuration parameter form TrmTrackQAFraction  #12387 (comment)
removing reference to non relevant pull request

* ATO-647 - applying clang-format

* ATO-647 - adding Tsaliss downsamling using configured parameter

* ATO-647 - clang  format

* ATO-647 - get sqrts from CCDB, clang  format

* ATO-647 - making clang hapy

* ATO-647 - disabling initialization part - causing failure
* CCDB can not be used in init
* problem to read fraction

* ATO-647 - fixing setting and  reading of the fraction

* ATO-647 - initializing SqrtS only once as suggested by Ruben comment

* ATO-647 -cleanup the code / remoing copypaste leftover

* ATO-647 -adding brackets  to resolve pull request

* ATO-647 -use simple track without error for propagation

* ATO-647 - speed up dEdx ratio calculation

* ATO-647 - uing only TPC gid (commented by Ruben)

* ATO-647 - adding mTrackQCNTrCut cut to reduce amount of the "short tracks" triggered by sampling

* ATO-647 - make consistent type for the NTrCut+warning removal + changing index

* ATO-647 -use mTableTrID as index
trackQAInfoHolder.trackID = mTableTrID;

* ATO-647 -fixing minor comments of Jan-Fiete - for documentation string

* ATO-647 -adding tpcTime0 for the QA of the V0 and Cascade finder

* ATO-647 -adding tpcTime0 in respect to collision time from AO2D following
#12387 (comment)

* ATO-647 bug fix in diffBCRef

---------

Co-authored-by: miranov25 <marian.ivanov@cern.cg>
  • Loading branch information
@miranov25
miranov25 and miranov25 authored Dec 12, 2023
1 parent d2ffb1f commit fdb237b
Show file tree
Hide file tree
Showing 3 changed files with 168 additions and 8 deletions.
31 changes: 29 additions & 2 deletions Detectors/AOD/include/AODProducerWorkflow/AODProducerWorkflowSpec.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@
#include <set>
#include <string>
#include <vector>

#include <random>
using namespace o2::framework;
using GID = o2::dataformats::GlobalTrackID;
using GIndex = o2::dataformats::VtxTrackIndex;
Expand Down Expand Up @@ -229,6 +229,10 @@ class AODProducerWorkflowDPL : public Task

bool mPropTracks{false};
bool mPropMuons{false};
float mTrackQCFraction{0.00};
int64_t mTrackQCNTrCut{4};
float mSqrtS{13860.};
std::mt19937 mGenerator; ///< random generator for trackQA sampling
o2::base::Propagator::MatCorrType mMatCorr{o2::base::Propagator::MatCorrType::USEMatCorrLUT};
o2::dataformats::MeanVertexObject mVtx;
float mMinPropR{o2::constants::geom::XTPCInnerRef + 0.1f};
Expand Down Expand Up @@ -381,9 +385,26 @@ class AODProducerWorkflowDPL : public Task
float trackPhiEMCAL = -999.f;
float trackTime = -999.f;
float trackTimeRes = -999.f;
int diffBCRef = 0; // offset of time reference BC from the start of the orbit
int bcSlice[2] = {-1, -1};
};

struct TrackQA {
GID trackID;
float tpcTime0;
int16_t tpcdcaR;
int16_t tpcdcaZ;
uint8_t tpcClusterByteMask;
uint8_t tpcdEdxMax0R;
uint8_t tpcdEdxMax1R;
uint8_t tpcdEdxMax2R;
uint8_t tpcdEdxMax3R;
uint8_t tpcdEdxTot0R;
uint8_t tpcdEdxTot1R;
uint8_t tpcdEdxTot2R;
uint8_t tpcdEdxTot3R;
};

// helper struct for addToFwdTracksTable()
struct FwdTrackInfo {
uint8_t trackTypeId = 0;
Expand Down Expand Up @@ -462,6 +483,9 @@ class AODProducerWorkflowDPL : public Task
template <typename TracksExtraCursorType>
void addToTracksExtraTable(TracksExtraCursorType& tracksExtraCursor, TrackExtraInfo& extraInfoHolder);

template <typename TracksQACursorType>
void addToTracksQATable(TracksQACursorType& tracksQACursor, TrackQA& trackQAInfoHolder);

template <typename mftTracksCursorType, typename AmbigMFTTracksCursorType>
void addToMFTTracksTable(mftTracksCursorType& mftTracksCursor, AmbigMFTTracksCursorType& ambigMFTTracksCursor,
GIndex trackID, const o2::globaltracking::RecoContainer& data, int collisionID,
Expand All @@ -472,14 +496,16 @@ class AODProducerWorkflowDPL : public Task
GIndex trackID, const o2::globaltracking::RecoContainer& data, int collisionID, std::uint64_t collisionBC, const std::map<uint64_t, int>& bcsMap);

TrackExtraInfo processBarrelTrack(int collisionID, std::uint64_t collisionBC, GIndex trackIndex, const o2::globaltracking::RecoContainer& data, const std::map<uint64_t, int>& bcsMap);
TrackQA processBarrelTrackQA(int collisionID, std::uint64_t collisionBC, GIndex trackIndex, const o2::globaltracking::RecoContainer& data, const std::map<uint64_t, int>& bcsMap);

bool propagateTrackToPV(o2::track::TrackParametrizationWithError<float>& trackPar, const o2::globaltracking::RecoContainer& data, int colID);
void extrapolateToCalorimeters(TrackExtraInfo& extraInfoHolder, const o2::track::TrackPar& track);
void cacheTriggers(const o2::globaltracking::RecoContainer& recoData);

// helper for track tables
// * fills tables collision by collision
// * interaction time is for TOF information
template <typename TracksCursorType, typename TracksCovCursorType, typename TracksExtraCursorType, typename AmbigTracksCursorType,
template <typename TracksCursorType, typename TracksCovCursorType, typename TracksExtraCursorType, typename TracksQACursorType, typename AmbigTracksCursorType,
typename MFTTracksCursorType, typename AmbigMFTTracksCursorType,
typename FwdTracksCursorType, typename FwdTracksCovCursorType, typename AmbigFwdTracksCursorType>
void fillTrackTablesPerCollision(int collisionID,
Expand All @@ -490,6 +516,7 @@ class AODProducerWorkflowDPL : public Task
TracksCursorType& tracksCursor,
TracksCovCursorType& tracksCovCursor,
TracksExtraCursorType& tracksExtraCursor,
TracksQACursorType& tracksQACursor,
AmbigTracksCursorType& ambigTracksCursor,
MFTTracksCursorType& mftTracksCursor,
AmbigMFTTracksCursorType& ambigMFTTracksCursor,
Expand Down
119 changes: 113 additions & 6 deletions Detectors/AOD/src/AODProducerWorkflowSpec.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -92,6 +92,8 @@
#include <thread>
#include "TLorentzVector.h"
#include "TVector3.h"
#include "MathUtils/Tsallis.h"
#include <random>
#ifdef WITH_OPENMP
#include <omp.h>
#endif
Expand Down Expand Up @@ -340,6 +342,30 @@ void AODProducerWorkflowDPL::addToTracksExtraTable(TracksExtraCursorType& tracks
truncateFloatFraction(extraInfoHolder.trackTimeRes, mTrackTimeError));
}

template <typename TracksQACursorType>
void AODProducerWorkflowDPL::addToTracksQATable(TracksQACursorType& tracksQACursor, TrackQA& trackQAInfoHolder)
{

// trackQA
tracksQACursor(

// truncateFloatFraction(trackQAInfoHolder.tpcdcaR, mTrackChi2),
// truncateFloatFraction(trackQAInfoHolder.tpcdcaZ, mTrackChi2),
trackQAInfoHolder.trackID,
trackQAInfoHolder.tpcTime0,
trackQAInfoHolder.tpcdcaR,
trackQAInfoHolder.tpcdcaZ,
trackQAInfoHolder.tpcClusterByteMask,
trackQAInfoHolder.tpcdEdxMax0R,
trackQAInfoHolder.tpcdEdxMax1R,
trackQAInfoHolder.tpcdEdxMax2R,
trackQAInfoHolder.tpcdEdxMax3R,
trackQAInfoHolder.tpcdEdxTot0R,
trackQAInfoHolder.tpcdEdxTot1R,
trackQAInfoHolder.tpcdEdxTot2R,
trackQAInfoHolder.tpcdEdxTot3R);
}

template <typename mftTracksCursorType, typename AmbigMFTTracksCursorType>
void AODProducerWorkflowDPL::addToMFTTracksTable(mftTracksCursorType& mftTracksCursor, AmbigMFTTracksCursorType& ambigMFTTracksCursor,
GIndex trackID, const o2::globaltracking::RecoContainer& data, int collisionID,
Expand Down Expand Up @@ -381,7 +407,7 @@ void AODProducerWorkflowDPL::addToMFTTracksTable(mftTracksCursorType& mftTracksC
}
}

template <typename TracksCursorType, typename TracksCovCursorType, typename TracksExtraCursorType, typename AmbigTracksCursorType,
template <typename TracksCursorType, typename TracksCovCursorType, typename TracksExtraCursorType, typename TracksQACursorType, typename AmbigTracksCursorType,
typename MFTTracksCursorType, typename AmbigMFTTracksCursorType,
typename FwdTracksCursorType, typename FwdTracksCovCursorType, typename AmbigFwdTracksCursorType>
void AODProducerWorkflowDPL::fillTrackTablesPerCollision(int collisionID,
Expand All @@ -392,6 +418,7 @@ void AODProducerWorkflowDPL::fillTrackTablesPerCollision(int collisionID,
TracksCursorType& tracksCursor,
TracksCovCursorType& tracksCovCursor,
TracksExtraCursorType& tracksExtraCursor,
TracksQACursorType& tracksQACursor,
AmbigTracksCursorType& ambigTracksCursor,
MFTTracksCursorType& mftTracksCursor,
AmbigMFTTracksCursorType& ambigMFTTracksCursor,
Expand Down Expand Up @@ -440,6 +467,21 @@ void AODProducerWorkflowDPL::fillTrackTablesPerCollision(int collisionID,
continue;
}
auto extraInfoHolder = processBarrelTrack(collisionID, collisionBC, trackIndex, data, bcsMap);

float weight = 0;
std::uniform_real_distribution<> distr(0., 1.);
bool writeQAData = o2::math_utils::Tsallis::downsampleTsallisCharged(data.getTrackParam(trackIndex).getPt(), mTrackQCFraction, mSqrtS, weight, distr(mGenerator));
if (writeQAData) {
auto trackQAInfoHolder = processBarrelTrackQA(collisionID, collisionBC, trackIndex, data, bcsMap);
if (std::bitset<8>(trackQAInfoHolder.tpcClusterByteMask).count() >= mTrackQCNTrCut) {
trackQAInfoHolder.trackID = mTableTrID;
// LOGP(info, "orig time0 in bc: {} diffBCRef: {}, ttime: {} -> {}", trackQAInfoHolder.tpcTime0*8, extraInfoHolder.diffBCRef, extraInfoHolder.trackTime, (trackQAInfoHolder.tpcTime0 * 8 - extraInfoHolder.diffBCRef) * o2::constants::lhc::LHCBunchSpacingNS - extraInfoHolder.trackTime);
trackQAInfoHolder.tpcTime0 = (trackQAInfoHolder.tpcTime0 * 8 - extraInfoHolder.diffBCRef) * o2::constants::lhc::LHCBunchSpacingNS - extraInfoHolder.trackTime;
// difference between TPC track time0 and stored track nominal time in ns instead of TF start
addToTracksQATable(tracksQACursor, trackQAInfoHolder);
}
}

if (extraInfoHolder.trackTimeRes < 0.f) { // failed or rejected?
LOG(warning) << "Barrel track " << trackIndex << " has no time set, rejection is not expected : time=" << extraInfoHolder.trackTime
<< " timeErr=" << extraInfoHolder.trackTimeRes << " BCSlice: " << extraInfoHolder.bcSlice[0] << ":" << extraInfoHolder.bcSlice[1];
Expand All @@ -458,7 +500,8 @@ void AODProducerWorkflowDPL::fillTrackTablesPerCollision(int collisionID,
addToTracksTable(tracksCursor, tracksCovCursor, trOrig, collisionID, aod::track::TrackIU);
}
addToTracksExtraTable(tracksExtraCursor, extraInfoHolder);
// collecting table indices of barrel tracks for V0s table
// addToTracksQATable(tracksQACursor, trackQAInfoHolder);
// collecting table indices of barrel tracks for V0s table
if (extraInfoHolder.bcSlice[0] >= 0 && collisionID < 0) {
ambigTracksCursor(mTableTrID, extraInfoHolder.bcSlice);
}
Expand Down Expand Up @@ -1628,6 +1671,9 @@ void AODProducerWorkflowDPL::init(InitContext& ic)
mEMCselectLeading = ic.options().get<bool>("emc-select-leading");
mPropTracks = ic.options().get<bool>("propagate-tracks");
mPropMuons = ic.options().get<bool>("propagate-muons");
mTrackQCFraction = ic.options().get<float>("trackqc-fraction");
mTrackQCNTrCut = ic.options().get<int64_t>("trackqc-NTrCut");
mGenerator = std::mt19937(std::random_device{}());
#ifdef WITH_OPENMP
LOGP(info, "Multi-threaded parts will run with {} OpenMP threads", mNThreads);
#else
Expand Down Expand Up @@ -1761,6 +1807,7 @@ void AODProducerWorkflowDPL::run(ProcessingContext& pc)
auto tracksCursor = createTableCursor<o2::aod::StoredTracksIU>(pc);
auto tracksCovCursor = createTableCursor<o2::aod::StoredTracksCovIU>(pc);
auto tracksExtraCursor = createTableCursor<o2::aod::StoredTracksExtra>(pc);
auto tracksQACursor = createTableCursor<o2::aod::TrackQA>(pc);
auto ambigTracksCursor = createTableCursor<o2::aod::AmbiguousTracks>(pc);
auto ambigMFTTracksCursor = createTableCursor<o2::aod::AmbiguousMFTTracks>(pc);
auto ambigFwdTracksCursor = createTableCursor<o2::aod::AmbiguousFwdTracks>(pc);
Expand Down Expand Up @@ -2065,7 +2112,7 @@ void AODProducerWorkflowDPL::run(ProcessingContext& pc)
// so that all unassigned tracks are stored in the beginning of the table together
auto& trackRef = primVer2TRefs.back(); // references to unassigned tracks are at the end
// fixme: interaction time is undefined for unassigned tracks (?)
fillTrackTablesPerCollision(-1, std::uint64_t(-1), trackRef, primVerGIs, recoData, tracksCursor, tracksCovCursor, tracksExtraCursor,
fillTrackTablesPerCollision(-1, std::uint64_t(-1), trackRef, primVerGIs, recoData, tracksCursor, tracksCovCursor, tracksExtraCursor, tracksQACursor,
ambigTracksCursor, mftTracksCursor, ambigMFTTracksCursor,
fwdTracksCursor, fwdTracksCovCursor, ambigFwdTracksCursor, bcsMap);

Expand Down Expand Up @@ -2107,7 +2154,7 @@ void AODProducerWorkflowDPL::run(ProcessingContext& pc)

auto& trackRef = primVer2TRefs[collisionID];
// passing interaction time in [ps]
fillTrackTablesPerCollision(collisionID, globalBC, trackRef, primVerGIs, recoData, tracksCursor, tracksCovCursor, tracksExtraCursor, ambigTracksCursor,
fillTrackTablesPerCollision(collisionID, globalBC, trackRef, primVerGIs, recoData, tracksCursor, tracksCovCursor, tracksExtraCursor, tracksQACursor, ambigTracksCursor,
mftTracksCursor, ambigMFTTracksCursor,
fwdTracksCursor, fwdTracksCovCursor, ambigFwdTracksCursor, bcsMap);
collisionID++;
Expand Down Expand Up @@ -2341,6 +2388,7 @@ AODProducerWorkflowDPL::TrackExtraInfo AODProducerWorkflowDPL::processBarrelTrac
bcOfTimeRef = fillBCSlice(extraInfoHolder.bcSlice, t - error, t + error, bcsMap);
}
extraInfoHolder.trackTime = float(t - bcOfTimeRef * o2::constants::lhc::LHCBunchSpacingNS);
extraInfoHolder.diffBCRef = int(bcOfTimeRef);
LOGP(debug, "time : {}/{} -> {}/{} -> trunc: {}/{} CollID: {} Amb: {}", t, terr, t - bcOfTimeRef * o2::constants::lhc::LHCBunchSpacingNS, terr,
truncateFloatFraction(extraInfoHolder.trackTime, mTrackTime), truncateFloatFraction(extraInfoHolder.trackTimeRes, mTrackTimeError),
collisionID, trackIndex.isAmbiguous());
Expand Down Expand Up @@ -2426,6 +2474,61 @@ AODProducerWorkflowDPL::TrackExtraInfo AODProducerWorkflowDPL::processBarrelTrac
return extraInfoHolder;
}

AODProducerWorkflowDPL::TrackQA AODProducerWorkflowDPL::processBarrelTrackQA(int collisionID, std::uint64_t collisionBC, GIndex trackIndex,
const o2::globaltracking::RecoContainer& data, const std::map<uint64_t, int>& bcsMap)
{
TrackQA trackQAHolder;
auto contributorsGID = data.getTPCContributorGID(trackIndex);
const auto& trackPar = data.getTrackParam(trackIndex);
// auto src = trackIndex.getSource();
if (contributorsGID.isIndexSet()) {
const auto& tpcOrig = data.getTPCTrack(contributorsGID);
/// getDCA - should be done with the copy of TPC only track
// LOGP(info, "GloIdx: {} TPCIdx: {}, NTPCTracks: {}", trackIndex.asString(), contributorsGID.asString(), data.getTPCTracks().size());
o2::track::TrackParametrization<float> tpcTMP = tpcOrig; /// get backup of the track
o2::base::Propagator::MatCorrType mMatType = o2::base::Propagator::MatCorrType::USEMatCorrLUT; /// should be parameterized
o2::dataformats::VertexBase v = mVtx.getMeanVertex(collisionID < 0 ? 0.f : data.getPrimaryVertex(collisionID).getZ());
o2::gpu::gpustd::array<float, 2> dcaInfo{-999., -999.};
if (o2::base::Propagator::Instance()->propagateToDCABxByBz({v.getX(), v.getY(), v.getZ()}, tpcTMP, 2.f, mMatType, &dcaInfo)) {
trackQAHolder.tpcdcaR = 100. * dcaInfo[0] / sqrt(1. + trackPar.getQ2Pt() * trackPar.getQ2Pt());
trackQAHolder.tpcdcaZ = 100. * dcaInfo[1] / sqrt(1. + trackPar.getQ2Pt() * trackPar.getQ2Pt());
}
/// get tracklet byteMask
uint8_t clusterCounters[8] = {0};
{
uint8_t sectorIndex, rowIndex;
uint32_t clusterIndex;
const auto& tpcClusRefs = data.getTPCTracksClusterRefs();
for (int i = 0; i < tpcOrig.getNClusterReferences(); i++) {
o2::tpc::TrackTPC::getClusterReference(tpcClusRefs, i, sectorIndex, rowIndex, clusterIndex, tpcOrig.getClusterRef());
char indexTracklet = (rowIndex % 152) / 19;
clusterCounters[indexTracklet]++;
}
}
uint8_t byteMask = 0;
for (int i = 0; i < 8; i++) {
if (clusterCounters[i] > 5) {
byteMask |= (1 << i);
}
}
trackQAHolder.tpcTime0 = tpcOrig.getTime0();
trackQAHolder.tpcClusterByteMask = byteMask;
float dEdxNorm = (tpcOrig.getdEdx().dEdxTotTPC > 0) ? 100. / tpcOrig.getdEdx().dEdxTotTPC : 0;
trackQAHolder.tpcdEdxMax0R = uint8_t(tpcOrig.getdEdx().dEdxMaxIROC * dEdxNorm);
trackQAHolder.tpcdEdxMax1R = uint8_t(tpcOrig.getdEdx().dEdxMaxOROC1 * dEdxNorm);
trackQAHolder.tpcdEdxMax2R = uint8_t(tpcOrig.getdEdx().dEdxMaxOROC2 * dEdxNorm);
trackQAHolder.tpcdEdxMax3R = uint8_t(tpcOrig.getdEdx().dEdxMaxOROC3 * dEdxNorm);
//
trackQAHolder.tpcdEdxTot0R = uint8_t(tpcOrig.getdEdx().dEdxTotIROC * dEdxNorm);
trackQAHolder.tpcdEdxTot1R = uint8_t(tpcOrig.getdEdx().dEdxTotOROC1 * dEdxNorm);
trackQAHolder.tpcdEdxTot2R = uint8_t(tpcOrig.getdEdx().dEdxTotOROC2 * dEdxNorm);
trackQAHolder.tpcdEdxTot3R = uint8_t(tpcOrig.getdEdx().dEdxTotOROC3 * dEdxNorm);
///
}

return trackQAHolder;
}

bool AODProducerWorkflowDPL::propagateTrackToPV(o2::track::TrackParametrizationWithError<float>& trackPar,
const o2::globaltracking::RecoContainer& data,
int colID)
Expand Down Expand Up @@ -2551,7 +2654,7 @@ void AODProducerWorkflowDPL::updateTimeDependentParams(ProcessingContext& pc)
if (!initOnceDone) { // this params need to be queried only once
initOnceDone = true;
// Note: DPLAlpideParam for ITS and MFT will be loaded by the RecoContainer

mSqrtS = o2::base::GRPGeomHelper::instance().getGRPLHCIF()->getSqrtS();
// apply settings
auto grpECS = o2::base::GRPGeomHelper::instance().getGRPECS();
o2::BunchFilling bcf = o2::base::GRPGeomHelper::instance().getGRPLHCIF()->getBunchFilling();
Expand Down Expand Up @@ -2829,6 +2932,7 @@ DataProcessorSpec getAODProducerWorkflowSpec(GID::mask_t src, bool enableSV, boo
OutputForTable<StoredTracksIU>::spec(),
OutputForTable<StoredTracksCovIU>::spec(),
OutputForTable<StoredTracksExtra>::spec(),
OutputForTable<TracksQA>::spec(),
OutputForTable<TrackedCascades>::spec(),
OutputForTable<TrackedV0s>::spec(),
OutputForTable<Tracked3Bodys>::spec(),
Expand Down Expand Up @@ -2872,7 +2976,10 @@ DataProcessorSpec getAODProducerWorkflowSpec(GID::mask_t src, bool enableSV, boo
ConfigParamSpec{"ctpreadout-create", VariantType::Int, 0, {"Create CTP digits from detector readout and CTP inputs. !=1 -- off, 1 -- on"}},
ConfigParamSpec{"emc-select-leading", VariantType::Bool, false, {"Flag to select if only the leading contributing particle for an EMCal cell should be stored"}},
ConfigParamSpec{"propagate-tracks", VariantType::Bool, false, {"Propagate tracks (not used for secondary vertices) to IP"}},
ConfigParamSpec{"propagate-muons", VariantType::Bool, false, {"Propagate muons to IP"}}}};
ConfigParamSpec{"propagate-muons", VariantType::Bool, false, {"Propagate muons to IP"}},
ConfigParamSpec{"trackqc-fraction", VariantType::Float, float(0.1), {"Fraction of tracks to QC"}},
ConfigParamSpec{"trackqc-NTrCut", VariantType::Int64, 4L, {"Minimal length of the track - in amount of tracklets"}},
}};
}

} // namespace o2::aodproducer
Loading

0 comments on commit fdb237b

Please sign in to comment.