Skip to content

Commit

Permalink
Add strange TOF PID QA task (AliceO2Group#4245)
Browse files Browse the repository at this point in the history
* Make pid task independent of orig ao2d

* Add TOF PID QA task for strangeness

* Please consider the following formatting changes (#212)

---------

Co-authored-by: ALICE Builder <[email protected]>
  • Loading branch information
ddobrigk and alibuild authored Jan 7, 2024
1 parent 88a7673 commit 25dae17
Show file tree
Hide file tree
Showing 5 changed files with 154 additions and 53 deletions.
2 changes: 1 addition & 1 deletion PWGLF/DataModel/LFStrangenessTables.h
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@ namespace o2::aod
// this is optional but will ensure full flexibility
// if required (for 2pc, etc)
DECLARE_SOA_TABLE(StraCollisions, "AOD", "STRACOLLISION", //! basic collision properties: position
collision::PosX, collision::PosY, collision::PosZ);
o2::soa::Index<>, collision::PosX, collision::PosY, collision::PosZ);
DECLARE_SOA_TABLE(StraCents, "AOD", "STRACENTS", //! centrality percentiles
cent::CentFT0M, cent::CentFT0A,
cent::CentFT0C, cent::CentFV0A);
Expand Down
76 changes: 26 additions & 50 deletions PWGLF/TableProducer/lambdakzeropid.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -75,14 +75,17 @@ using TaggedV0s = soa::Join<aod::V0s, aod::V0Tags>;
// For MC association in pre-selection
using LabeledTracksExtra = soa::Join<aod::TracksExtra, aod::McTrackLabels>;

// Cores with references and TOF pid
using V0FullCores = soa::Join<aod::V0Cores, aod::V0TOFs, aod::V0CollRefs>;

struct lambdakzeropid {
// TOF pid for strangeness (recalculated with topology)
Produces<aod::V0TOFPIDs> v0tofpid; // table with Nsigmas

Service<o2::ccdb::BasicCCDBManager> ccdb;

// For manual sliceBy
Preslice<aod::V0Datas> perCollision = o2::aod::v0data::collisionId;
Preslice<V0FullCores> perCollision = o2::aod::v0data::straCollisionId;

HistogramRegistry histos{"Histos", {}, OutputObjHandlingPolicy::AnalysisObject};

Expand Down Expand Up @@ -117,7 +120,7 @@ struct lambdakzeropid {
/// \param x2 x of the first point of the detector segment
/// \param y2 y of the first point of the detector segment
/// \param magneticField the magnetic field to use when propagating
float trackLengthToSegment(o2::track::TrackParCov track, float x1, float y1, float x2, float y2, float magneticField)
float trackLengthToSegment(o2::track::TrackPar track, float x1, float y1, float x2, float y2, float magneticField)
{
// don't make use of the track parametrization
float length = -104;
Expand Down Expand Up @@ -216,7 +219,7 @@ struct lambdakzeropid {
/// function to calculate track length of this track up to a certain segmented detector
/// \param track the input track
/// \param magneticField the magnetic field to use when propagating
float findInterceptLength(o2::track::TrackParCov track, float magneticField)
float findInterceptLength(o2::track::TrackPar track, float magneticField)
{
for (int iSeg = 0; iSeg < 18; iSeg++) {
// Detector segmentation loop
Expand Down Expand Up @@ -246,17 +249,13 @@ struct lambdakzeropid {
ccdb->setLocalObjectValidityChecking();
ccdb->setFatalWhenNull(false);

histos.add("h3dMassK0ShortPositive", "h3dMassK0ShortPositive", kTH3F, {axisPtQA, axisDeltaTime, axisK0ShortMass});
histos.add("h3dMassLambdaPositive", "h3dMassLambdaPositive", kTH3F, {axisPtQA, axisDeltaTime, axisLambdaMass});
histos.add("h3dMassAntiLambdaPositive", "h3dMassAntiLambdaPositive", kTH3F, {axisPtQA, axisDeltaTime, axisLambdaMass});
histos.add("h3dMassK0ShortNegative", "h3dMassK0ShortNegative", kTH3F, {axisPtQA, axisDeltaTime, axisK0ShortMass});
histos.add("h3dMassLambdaNegative", "h3dMassLambdaNegative", kTH3F, {axisPtQA, axisDeltaTime, axisLambdaMass});
histos.add("h3dMassAntiLambdaNegative", "h3dMassAntiLambdaNegative", kTH3F, {axisPtQA, axisDeltaTime, axisLambdaMass});
// per event
histos.add("hCandidateCounter", "hCandidateCounter", kTH1F, {{500, -0.5f, 499.5f}});
}

void initCCDB(aod::BCsWithTimestamps::iterator const& bc)
void initCCDB(soa::Join<aod::StraCollisions, aod::StraStamps>::iterator const& collision)
{
if (mRunNumber == bc.runNumber()) {
if (mRunNumber == collision.runNumber()) {
return;
}

Expand All @@ -268,11 +267,11 @@ struct lambdakzeropid {
grpmag.setL3Current(30000.f / (d_bz / 5.0f));
}
o2::base::Propagator::initFieldFromGRP(&grpmag);
mRunNumber = bc.runNumber();
mRunNumber = collision.runNumber();
return;
}

auto run3grp_timestamp = bc.timestamp();
auto run3grp_timestamp = collision.timestamp();
o2::parameters::GRPObject* grpo = ccdb->getForTimeStamp<o2::parameters::GRPObject>(grpPath, run3grp_timestamp);
o2::parameters::GRPMagField* grpmag = 0x0;
if (grpo) {
Expand All @@ -290,7 +289,7 @@ struct lambdakzeropid {
d_bz = std::lround(5.f * grpmag->getL3Current() / 30000.f);
LOG(info) << "Retrieved GRP for timestamp " << run3grp_timestamp << " with magnetic field of " << d_bz << " kZG";
}
mRunNumber = bc.runNumber();
mRunNumber = collision.runNumber();
}

float velocity(float lMomentum, float lMass)
Expand All @@ -301,15 +300,15 @@ struct lambdakzeropid {
return 0.0299792458 * TMath::Sqrt(lA / (1 + lA));
}

void process(aod::Collisions const& collisions, aod::V0Datas const& V0s, FullTracksExtIU const&, aod::BCsWithTimestamps const&, TaggedV0s const& allV0s)
void process(soa::Join<aod::StraCollisions, aod::StraStamps> const& collisions, V0FullCores const& V0s)
{
for (const auto& collision : collisions) {
// Fire up CCDB
auto bc = collision.bc_as<aod::BCsWithTimestamps>();
initCCDB(bc);
// Fire up CCDB - based on StraCollisions for derived analysis
initCCDB(collision);
// Do analysis with collision-grouped V0s, retain full collision information
const uint64_t collIdx = collision.globalIndex();
auto V0Table_thisCollision = V0s.sliceBy(perCollision, collIdx);
histos.fill(HIST("hCandidateCounter"), V0Table_thisCollision.size());
// V0 table sliced
for (auto const& v0 : V0Table_thisCollision) {
// time of V0 segment
Expand All @@ -319,11 +318,9 @@ struct lambdakzeropid {
float timeK0Short = lengthV0 / velocityK0Short; // in picoseconds
float timeLambda = lengthV0 / velocityLambda; // in picoseconds

auto const& posTrackRow = v0.posTrack_as<FullTracksExtIU>();
auto const& negTrackRow = v0.negTrack_as<FullTracksExtIU>();

auto posTrack = getTrackParCov(posTrackRow);
auto negTrack = getTrackParCov(negTrackRow);
// initialize from V0 position and momenta
o2::track::TrackPar posTrack = o2::track::TrackPar({v0.x(), v0.y(), v0.z()}, {v0.pxpos(), v0.pypos(), v0.pzpos()}, +1);
o2::track::TrackPar negTrack = o2::track::TrackPar({v0.x(), v0.y(), v0.z()}, {v0.pxneg(), v0.pyneg(), v0.pzneg()}, -1);

float deltaTimePositiveLambdaPi = -1e+6;
float deltaTimeNegativeLambdaPi = -1e+6;
Expand All @@ -337,46 +334,25 @@ struct lambdakzeropid {
float velocityNegativePr = velocity(negTrack.getP(), o2::constants::physics::MassProton);
float velocityNegativePi = velocity(negTrack.getP(), o2::constants::physics::MassPionCharged);

// propagate to V0 decay vertex
posTrack.propagateTo(v0.posX(), d_bz);
negTrack.propagateTo(v0.negX(), d_bz);

float lengthPositive = findInterceptLength(posTrack, d_bz); // FIXME: tofPosition ok? adjust?
float lengthNegative = findInterceptLength(negTrack, d_bz); // FIXME: tofPosition ok? adjust?
float timePositivePr = lengthPositive / velocityPositivePr;
float timePositivePi = lengthPositive / velocityPositivePi;
float timeNegativePr = lengthNegative / velocityNegativePr;
float timeNegativePi = lengthNegative / velocityNegativePi;

deltaTimePositiveLambdaPr = (posTrackRow.tofSignal() - posTrackRow.tofEvTime()) - (timeLambda + timePositivePr);
deltaTimePositiveLambdaPi = (posTrackRow.tofSignal() - posTrackRow.tofEvTime()) - (timeLambda + timePositivePi);
deltaTimeNegativeLambdaPr = (negTrackRow.tofSignal() - negTrackRow.tofEvTime()) - (timeLambda + timeNegativePr);
deltaTimeNegativeLambdaPi = (negTrackRow.tofSignal() - negTrackRow.tofEvTime()) - (timeLambda + timeNegativePi);
deltaTimePositiveK0ShortPi = (posTrackRow.tofSignal() - posTrackRow.tofEvTime()) - (timeK0Short + timeNegativePi);
deltaTimeNegativeK0ShortPi = (negTrackRow.tofSignal() - negTrackRow.tofEvTime()) - (timeK0Short + timeNegativePi);
deltaTimePositiveLambdaPr = (v0.posTOFSignal() - v0.posTOFEventTime()) - (timeLambda + timePositivePr);
deltaTimePositiveLambdaPi = (v0.posTOFSignal() - v0.posTOFEventTime()) - (timeLambda + timePositivePi);
deltaTimeNegativeLambdaPr = (v0.negTOFSignal() - v0.negTOFEventTime()) - (timeLambda + timeNegativePr);
deltaTimeNegativeLambdaPi = (v0.negTOFSignal() - v0.negTOFEventTime()) - (timeLambda + timeNegativePi);
deltaTimePositiveK0ShortPi = (v0.posTOFSignal() - v0.posTOFEventTime()) - (timeK0Short + timeNegativePi);
deltaTimeNegativeK0ShortPi = (v0.negTOFSignal() - v0.negTOFEventTime()) - (timeK0Short + timeNegativePi);

v0tofpid(lengthPositive, lengthNegative,
deltaTimePositiveLambdaPi, deltaTimePositiveLambdaPr,
deltaTimeNegativeLambdaPi, deltaTimeNegativeLambdaPr,
deltaTimePositiveK0ShortPi, deltaTimeNegativeK0ShortPi,
0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f); // FIXME

auto originalV0 = v0.v0_as<TaggedV0s>(); // this could look confusing, so:
// the first v0 is the v0data row; the getter de-references the v0 (stored indices) row
// the v0 (stored indices) contain the tags of the lambdakzero preselector

if (originalV0.isdEdxK0Short() || !checkTPCCompatibility) {
histos.fill(HIST("h3dMassK0ShortPositive"), v0.pt(), deltaTimePositiveK0ShortPi, v0.mK0Short());
histos.fill(HIST("h3dMassK0ShortNegative"), v0.pt(), deltaTimePositiveK0ShortPi, v0.mK0Short());
}
if (originalV0.isdEdxLambda() || !checkTPCCompatibility) {
histos.fill(HIST("h3dMassLambdaPositive"), v0.pt(), deltaTimePositiveLambdaPr, v0.mLambda());
histos.fill(HIST("h3dMassLambdaNegative"), v0.pt(), deltaTimeNegativeLambdaPi, v0.mLambda());
}
if (originalV0.isdEdxAntiLambda() || !checkTPCCompatibility) {
histos.fill(HIST("h3dMassAntiLambdaPositive"), v0.pt(), deltaTimePositiveK0ShortPi, v0.mAntiLambda());
histos.fill(HIST("h3dMassAntiLambdaNegative"), v0.pt(), deltaTimeNegativeK0ShortPi, v0.mAntiLambda());
}
}
}
}
Expand Down
4 changes: 2 additions & 2 deletions PWGLF/TableProducer/strangederivedbuilder.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -157,7 +157,7 @@ struct strangederivedbuilder {
strangeCents(collision.centFT0M(), collision.centFT0A(),
collision.centFT0C(), collision.centFV0A());
auto bc = collision.bc_as<aod::BCsWithTimestamps>();
strangeStamps(bc.timestamp(), bc.runNumber());
strangeStamps(bc.runNumber(), bc.timestamp());
}
for (int i = 0; i < V0Table_thisColl.size(); i++)
v0collref(strangeColl.lastIndex());
Expand All @@ -182,7 +182,7 @@ struct strangederivedbuilder {
strangeCents(collision.centFT0M(), collision.centFT0A(),
collision.centFT0C(), collision.centFV0A());
auto bc = collision.bc_as<aod::BCsWithTimestamps>();
strangeStamps(bc.timestamp(), bc.runNumber());
strangeStamps(bc.runNumber(), bc.timestamp());
}
for (int i = 0; i < V0Table_thisColl.size(); i++)
v0collref(strangeColl.lastIndex());
Expand Down
5 changes: 5 additions & 0 deletions PWGLF/Tasks/QC/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -68,3 +68,8 @@ o2physics_add_dpl_workflow(strangenessqcpp
SOURCES strangenessQCPP.cxx
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
COMPONENT_NAME Analysis)

o2physics_add_dpl_workflow(strangepidqa
SOURCES strangepidqa.cxx
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
COMPONENT_NAME Analysis)
120 changes: 120 additions & 0 deletions PWGLF/Tasks/QC/strangepidqa.cxx
Original file line number Diff line number Diff line change
@@ -0,0 +1,120 @@
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.
//
// This task is designed to do QA to the TOF PID applied to strangeness
// in the regular framework

#include "Framework/runDataProcessing.h"
#include "Framework/AnalysisTask.h"
#include "Framework/AnalysisDataModel.h"
#include "Framework/ASoAHelpers.h"
#include "ReconstructionDataFormats/Track.h"
#include "Common/Core/RecoDecay.h"
#include "Common/Core/trackUtilities.h"
#include "PWGLF/DataModel/LFStrangenessTables.h"
#include "PWGLF/DataModel/LFStrangenessPIDTables.h"
#include "Common/Core/TrackSelection.h"
#include "Common/DataModel/TrackSelectionTables.h"
#include "Common/DataModel/EventSelection.h"
#include "Common/DataModel/Centrality.h"
#include "Common/DataModel/PIDResponse.h"

#include <TFile.h>
#include <TH2F.h>
#include <TProfile.h>
#include <TLorentzVector.h>
#include <Math/Vector4D.h>
#include <TPDGCode.h>
#include <TDatabasePDG.h>
#include <cmath>
#include <array>
#include <cstdlib>
#include "Framework/ASoAHelpers.h"

using namespace o2;
using namespace o2::framework;
using namespace o2::framework::expressions;
using std::array;
using std::cout;
using std::endl;

struct strangepidqa {
HistogramRegistry histos{"Histos", {}, OutputObjHandlingPolicy::AnalysisObject};
ConfigurableAxis vertexZ{"vertexZ", {30, -15.0f, 15.0f}, ""};

// base properties
ConfigurableAxis axisPt{"axisPt", {VARIABLE_WIDTH, 0.0f, 0.1f, 0.2f, 0.3f, 0.4f, 0.5f, 0.6f, 0.7f, 0.8f, 0.9f, 1.0f, 1.1f, 1.2f, 1.3f, 1.4f, 1.5f, 1.6f, 1.7f, 1.8f, 1.9f, 2.0f, 2.2f, 2.4f, 2.6f, 2.8f, 3.0f, 3.2f, 3.4f, 3.6f, 3.8f, 4.0f, 4.4f, 4.8f, 5.2f, 5.6f, 6.0f, 6.5f, 7.0f, 7.5f, 8.0f, 9.0f, 10.0f, 11.0f, 12.0f, 13.0f, 14.0f, 15.0f, 17.0f, 19.0f, 21.0f, 23.0f, 25.0f, 30.0f, 35.0f, 40.0f, 50.0f}, "p_{T} (GeV/c)"};
ConfigurableAxis axisRadius{"axisRadius", {200, 0.0f, 100.0f}, "V0 radius (cm)"};

// Invariant Mass
ConfigurableAxis axisLambdaMass{"axisLambdaMass", {200, 1.08f, 1.16f}, "M_{#Lambda} (GeV/c^{2})"};

// Delta time
ConfigurableAxis axisDeltaTime{"axisDeltaTime", {200, -1000.0f, +1000.0f}, "#Delta time (ps)"};

// Length axis
ConfigurableAxis axisLength{"axisLength", {600, 0.0f, +600.0f}, "track Length (cm)"};

void init(InitContext const&)
{
// Event counter
histos.add("hEventVertexZ", "hEventVertexZ", kTH1F, {vertexZ});

// V0 Radius
histos.add("h2dLambdaRadiusVsPt", "hLambdaRadiusVsPt", {HistType::kTH2F, {axisPt, axisRadius}});

// Invariant Mass
histos.add("h2dLambdaMassVsPt", "hLambdaMassVsPt", {HistType::kTH2F, {axisPt, axisLambdaMass}});

// radius vs prong length
histos.add("h2dProtonLengthVsRadius", "h2dProtonLengthVsRadius", {HistType::kTH2F, {axisRadius, axisLength}});
histos.add("h2dPionLengthVsRadius", "h2dPionLengthVsRadius", {HistType::kTH2F, {axisRadius, axisLength}});

// recalculated vs topv lengths
histos.add("h2dProtonLengthVsLengthToPV", "h2dProtonLengthVsRadiusToPV", {HistType::kTH2F, {axisRadius, axisRadius}});
histos.add("h2dPionLengthVsLengthToPV", "h2dPionLengthVsLengthToPV", {HistType::kTH2F, {axisRadius, axisRadius}});

// TOF PID testing for prongs
histos.add("h2dProtonDeltaTimeVsPt", "h2dProtonDeltaTimeVsPt", {HistType::kTH2F, {axisPt, axisDeltaTime}});
histos.add("h2dProtonDeltaTimeVsRadius", "h2dProtonDeltaTimeVsRadius", {HistType::kTH2F, {axisRadius, axisDeltaTime}});
histos.add("h2dPionDeltaTimeVsPt", "h2dPionDeltaTimeVsPt", {HistType::kTH2F, {axisPt, axisDeltaTime}});
histos.add("h2dPionDeltaTimeVsRadius", "h2dPionDeltaTimeVsRadius", {HistType::kTH2F, {axisRadius, axisDeltaTime}});
}

void process(aod::StraCollision const& coll, soa::Join<aod::V0Cores, aod::V0CollRefs, aod::V0Extras, aod::V0MCDatas, aod::V0TOFs, aod::V0TOFPIDs> const& v0s)
{
histos.fill(HIST("hEventVertexZ"), coll.posZ());

for (auto& lambda : v0s) { // selecting photons from Sigma0
if (lambda.pdgCode() != 3122)
continue;

histos.fill(HIST("h2dLambdaRadiusVsPt"), lambda.pt(), lambda.v0radius());
histos.fill(HIST("h2dLambdaMassVsPt"), lambda.pt(), lambda.mLambda());

histos.fill(HIST("h2dProtonDeltaTimeVsPt"), lambda.pt(), lambda.posTOFDeltaTLaPr());
histos.fill(HIST("h2dProtonDeltaTimeVsRadius"), lambda.v0radius(), lambda.posTOFDeltaTLaPr());
histos.fill(HIST("h2dPionDeltaTimeVsPt"), lambda.pt(), lambda.negTOFDeltaTLaPi());
histos.fill(HIST("h2dPionDeltaTimeVsRadius"), lambda.v0radius(), lambda.negTOFDeltaTLaPi());

histos.fill(HIST("h2dProtonLengthVsRadius"), lambda.v0radius(), lambda.posTOFLength());
histos.fill(HIST("h2dPionLengthVsRadius"), lambda.v0radius(), lambda.negTOFLength());

histos.fill(HIST("h2dProtonLengthVsLengthToPV"), lambda.posTOFLengthToPV(), lambda.posTOFLength());
histos.fill(HIST("h2dPionLengthVsLengthToPV"), lambda.negTOFLengthToPV(), lambda.negTOFLength());
}
}
};

WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
{
return WorkflowSpec{adaptAnalysisTask<strangepidqa>(cfgc)};
}

0 comments on commit 25dae17

Please sign in to comment.