Skip to content

Commit c3f1440

Browse files
wang-yuanzhewang-yuanzhe
authored
[PWGLF] Fix bug in getting ITSNSigma and update datamodel for hyphe4s (#11706)
1 parent 20d1e4e commit c3f1440

File tree

2 files changed

+41
-21
lines changed

2 files changed

+41
-21
lines changed

‎PWGLF/DataModel/LFHyperhelium4sigmaTables.h‎

Lines changed: 4 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -13,8 +13,8 @@
1313
/// \brief Slim hyperhelium4sigma tables
1414
/// \author Yuanzhe Wang <yuanzhe.wang@cern.ch>
1515

16-
#include "Framework/AnalysisDataModel.h"
1716
#include "Framework/ASoAHelpers.h"
17+
#include "Framework/AnalysisDataModel.h"
1818

1919
#ifndef PWGLF_DATAMODEL_LFHYPERHELIUM4SIGMATABLES_H_
2020
#define PWGLF_DATAMODEL_LFHYPERHELIUM4SIGMATABLES_H_
@@ -44,6 +44,7 @@ DECLARE_SOA_COLUMN(ItsChi2Moth, itsChi2Moth, float); // ITS
4444
DECLARE_SOA_COLUMN(ItsClusterSizesMoth, itsClusterSizesMoth, uint32_t); // ITS cluster size of the mother track
4545
DECLARE_SOA_COLUMN(ItsClusterSizesAlpha, itsClusterSizesAlpha, uint32_t); // ITS cluster size of the daughter alpha track
4646
DECLARE_SOA_COLUMN(NSigmaTPCAlpha, nSigmaTPCAlpha, float); // Number of tpc sigmas of the daughter alpha track
47+
DECLARE_SOA_COLUMN(NSigmaITSAlpha, nSigmaITSAlpha, float); // Number of ITS sigmas of the daughter alpha track
4748

4849
DECLARE_SOA_COLUMN(IsSignal, isSignal, bool); // bool: true for hyperhelium4signal
4950
DECLARE_SOA_COLUMN(IsSignalReco, isSignalReco, bool); // bool: true if the signal is reconstructed
@@ -71,7 +72,7 @@ DECLARE_SOA_TABLE(He4S2BCands, "AOD", "HE4S2BCANDS",
7172
he4scand::PxAlpha, he4scand::PyAlpha, he4scand::PzAlpha,
7273
he4scand::DcaMothPv, he4scand::DcaAlphaPv, he4scand::DcaKinkTopo,
7374
he4scand::ItsChi2Moth, he4scand::ItsClusterSizesMoth, he4scand::ItsClusterSizesAlpha,
74-
he4scand::NSigmaTPCAlpha);
75+
he4scand::NSigmaTPCAlpha, he4scand::NSigmaITSAlpha);
7576

7677
DECLARE_SOA_TABLE(MCHe4S2BCands, "AOD", "MCHE4S2BCANDS",
7778
o2::soa::Index<>,
@@ -81,7 +82,7 @@ DECLARE_SOA_TABLE(MCHe4S2BCands, "AOD", "MCHE4S2BCANDS",
8182
he4scand::PxAlpha, he4scand::PyAlpha, he4scand::PzAlpha,
8283
he4scand::DcaMothPv, he4scand::DcaAlphaPv, he4scand::DcaKinkTopo,
8384
he4scand::ItsChi2Moth, he4scand::ItsClusterSizesMoth, he4scand::ItsClusterSizesAlpha,
84-
he4scand::NSigmaTPCAlpha,
85+
he4scand::NSigmaTPCAlpha, he4scand::NSigmaITSAlpha,
8586
he4scand::IsSignal, he4scand::IsSignalReco, he4scand::IsCollReco, he4scand::IsSurvEvSelection,
8687
he4scand::TrueXDecVtx, he4scand::TrueYDecVtx, he4scand::TrueZDecVtx,
8788
he4scand::GenPxMoth, he4scand::GenPyMoth, he4scand::GenPzMoth,

‎PWGLF/TableProducer/Nuspex/hyperhelium4sigmaRecoTask.cxx‎

Lines changed: 37 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -242,6 +242,7 @@ struct Hyphe4sCandidate {
242242
uint32_t itsClusterSizeMoth = 0u;
243243
uint32_t itsClusterSizeDau = 0u;
244244
float nSigmaTPCDau = -999.f;
245+
float nSigmaITSDau = -999.f;
245246

246247
// mc information
247248
bool isSignal = false;
@@ -290,13 +291,17 @@ struct Hyperhelium4sigmaRecoTask {
290291
float mBz;
291292
o2::base::MatLayerCylSet* lut = nullptr;
292293

294+
o2::aod::ITSResponse itsResponse;
295+
293296
void init(InitContext const&)
294297
{
295298
// Axes
296299
const AxisSpec vertexZAxis{100, -15., 15., "vtx_{Z} [cm]"};
297300
const AxisSpec ptAxis{50, -10, 10, "#it{p}_{T} (GeV/#it{c})"};
298301
const AxisSpec nSigmaAxis{120, -6.f, 6.f, "n#sigma_{#alpha}"};
299302
const AxisSpec massAxis{100, 3.85, 4.25, "m (GeV/#it{c}^{2})"};
303+
const AxisSpec diffPxAxis{200, -10.f, 10.f, "#Delta #it{p}_{x} (GeV/#it{c})"};
304+
const AxisSpec diffPyAxis{200, -10.f, 10.f, "#Delta #it{p}_{y} (GeV/#it{c})"};
300305
const AxisSpec diffPtAxis{200, -10.f, 10.f, "#Delta #it{p}_{T} (GeV/#it{c})"};
301306
const AxisSpec diffPzAxis{200, -10.f, 10.f, "#Delta #it{p}_{z} (GeV/#it{c})"};
302307
const AxisSpec radiusAxis{40, 0.f, 40.f, "R (cm)"};
@@ -306,11 +311,15 @@ struct Hyperhelium4sigmaRecoTask {
306311
registry.add<TH1>("hCandidateCounter", "hCandidateCounter", HistType::kTH1F, {{3, 0, 3}});
307312

308313
if (doprocessMC == true) {
314+
itsResponse.setMCDefaultParameters();
315+
309316
registry.add<TH1>("hTrueCandidateCounter", "hTrueCandidateCounter", HistType::kTH1F, {{3, 0, 3}});
310317
registry.add<TH1>("hDiffSVx", ";#Delta x (cm);", HistType::kTH1F, {{200, -10, 10}});
311318
registry.add<TH1>("hDiffSVy", ";#Delta y (cm);", HistType::kTH1F, {{200, -10, 10}});
312319
registry.add<TH1>("hDiffSVz", ";#Delta z (cm);", HistType::kTH1F, {{200, -10, 10}});
313320
registry.add<TH2>("h2RecSVRVsTrueSVR", ";Reconstruced SV R (cm);True SV R (cm);", HistType::kTH2F, {radiusAxis, radiusAxis});
321+
registry.add<TH2>("h2TrueMotherDiffPxVsRecSVR", ";Reconstruced SV R (cm);#Delta #it{p}_{T} (GeV/#it{c});", HistType::kTH2F, {radiusAxis, diffPxAxis});
322+
registry.add<TH2>("h2TrueMotherDiffPyVsRecSVR", ";Reconstruced SV R (cm);#Delta #it{p}_{T} (GeV/#it{c});", HistType::kTH2F, {radiusAxis, diffPyAxis});
314323
registry.add<TH2>("h2TrueMotherDiffPtVsRecSVR", ";Reconstruced SV R (cm);#Delta #it{p}_{T} (GeV/#it{c});", HistType::kTH2F, {radiusAxis, diffPtAxis});
315324
registry.add<TH2>("h2TrueMotherDiffPzVsRecSVR", ";Reconstruced SV R (cm);#Delta #it{p}_{z} (GeV/#it{c});", HistType::kTH2F, {radiusAxis, diffPzAxis});
316325
registry.add<TH2>("h2TrueMotherDiffTglVsRecSVR", ";Reconstruced SV R (cm);#Delta tan#lambda;", HistType::kTH2F, {radiusAxis, {200, -1.f, 1.f}});
@@ -384,6 +393,7 @@ struct Hyperhelium4sigmaRecoTask {
384393

385394
hyphe4sCand.itsClusterSizeDau = trackDau.itsClusterSizes();
386395
hyphe4sCand.nSigmaTPCDau = trackDau.tpcNSigmaAl();
396+
hyphe4sCand.nSigmaITSDau = itsResponse.nSigmaITS<o2::track::PID::Alpha>(trackDau);
387397
}
388398

389399
template <typename TTrack>
@@ -449,7 +459,7 @@ struct Hyperhelium4sigmaRecoTask {
449459
hyphe4sCand.momDaug[0], hyphe4sCand.momDaug[1], hyphe4sCand.momDaug[2],
450460
hyphe4sCand.dcaXYMothPv, hyphe4sCand.dcaXYDauPv, hyphe4sCand.dcaKinkTopo,
451461
hyphe4sCand.chi2ITSMoth, hyphe4sCand.itsClusterSizeMoth, hyphe4sCand.itsClusterSizeDau,
452-
hyphe4sCand.nSigmaTPCDau);
462+
hyphe4sCand.nSigmaTPCDau, hyphe4sCand.nSigmaITSDau);
453463
}
454464
}
455465
PROCESS_SWITCH(Hyperhelium4sigmaRecoTask, processData, "process data", true);
@@ -550,22 +560,26 @@ struct Hyperhelium4sigmaRecoTask {
550560
o2::base::Propagator::Instance()->propagateToDCABxByBz({kinkCand.xDecVtx(), kinkCand.yDecVtx(), kinkCand.zDecVtx()}, mcMotherTrackPar, 2.f, matCorr, &dcaInfo);
551561
registry.fill(HIST("hDCAXYMothToRecSV"), dcaInfo[0]);
552562
registry.fill(HIST("hDCAZMothToRecSV"), dcaInfo[1]);
563+
std::array<float, 3> pMotherAtSV = {-999.f, -999.f, -999.f};
564+
mcMotherTrackPar.getPxPyPzGlo(pMotherAtSV);
565+
registry.fill(HIST("h2TrueMotherDiffPxVsRecSVR"), recSVR, pMotherAtSV[0] - kinkCand.pxMoth());
566+
registry.fill(HIST("h2TrueMotherDiffPyVsRecSVR"), recSVR, pMotherAtSV[1] - kinkCand.pyMoth());
553567

554-
float spKink = kinkCand.pxMoth() * kinkCand.pxDaug() + kinkCand.pyMoth() * kinkCand.pyDaug() + kinkCand.pzMoth() * kinkCand.pzDaug();
555-
float sptKink = kinkCand.pxMoth() * kinkCand.pxDaug() + kinkCand.pyMoth() * kinkCand.pyDaug();
568+
float spKinkXY = kinkCand.pxMoth() * kinkCand.pxDaug() + kinkCand.pyMoth() * kinkCand.pyDaug();
569+
float spKink = spKinkXY + kinkCand.pzMoth() * kinkCand.pzDaug();
556570
float pMoth = std::hypot(kinkCand.pxMoth(), kinkCand.pyMoth(), kinkCand.pzMoth());
557571
float pDaug = std::hypot(kinkCand.pxDaug(), kinkCand.pyDaug(), kinkCand.pzDaug());
558572

559-
float pMothDirVector[3] = {kinkCand.xDecVtx() - collision.posX(), kinkCand.yDecVtx() - collision.posY(), kinkCand.zDecVtx() - collision.posZ()};
560-
float epMothDirVector = std::hypot(pMothDirVector[0], pMothDirVector[1], pMothDirVector[2]);
561-
float eptMothDirVector = std::hypot(pMothDirVector[0], pMothDirVector[1]);
562-
float spKinkSV = pMothDirVector[0] * kinkCand.pxDaug() + pMothDirVector[1] * kinkCand.pyDaug() + pMothDirVector[2] * kinkCand.pzDaug();
563-
float sptKinkSV = pMothDirVector[0] * kinkCand.pxDaug() + pMothDirVector[1] * kinkCand.pyDaug();
573+
float mothPDir[3] = {kinkCand.xDecVtx() - collision.posX(), kinkCand.yDecVtx() - collision.posY(), kinkCand.zDecVtx() - collision.posZ()};
574+
float magMothPDirXY = std::hypot(mothPDir[0], mothPDir[1]);
575+
float magMothPDir = std::hypot(mothPDir[0], mothPDir[1], mothPDir[2]);
576+
float spKinkSV = mothPDir[0] * kinkCand.pxDaug() + mothPDir[1] * kinkCand.pyDaug() + mothPDir[2] * kinkCand.pzDaug();
577+
float sptKinkSV = mothPDir[0] * kinkCand.pxDaug() + mothPDir[1] * kinkCand.pyDaug();
564578

565579
float kinkAngle = spKink / (pMoth * pDaug);
566-
float kinkAngleSV = spKinkSV / (epMothDirVector * pDaug);
567-
float kinkAngleXY = sptKink / (kinkCand.ptMoth() * kinkCand.ptDaug());
568-
float kinkAngleXYSV = sptKinkSV / (eptMothDirVector * kinkCand.ptDaug());
580+
float kinkAngleSV = spKinkSV / (magMothPDir * pDaug);
581+
float kinkAngleXY = spKinkXY / (kinkCand.ptMoth() * kinkCand.ptDaug());
582+
float kinkAngleXYSV = sptKinkSV / (magMothPDirXY * kinkCand.ptDaug());
569583

570584
registry.fill(HIST("h2TrueMotherDiffPtVsKinkAngle"), kinkAngle, (mcMotherTrack.pt() - kinkCand.ptMoth()) / kinkCand.ptMoth());
571585
registry.fill(HIST("h2TrueMotherDiffPtVsKinkAngleSV"), kinkAngleSV, (mcMotherTrack.pt() - kinkCand.ptMoth()) / kinkCand.ptMoth());
@@ -580,7 +594,7 @@ struct Hyperhelium4sigmaRecoTask {
580594
hyphe4sCand.momDaug[0], hyphe4sCand.momDaug[1], hyphe4sCand.momDaug[2],
581595
hyphe4sCand.dcaXYMothPv, hyphe4sCand.dcaXYDauPv, hyphe4sCand.dcaKinkTopo,
582596
hyphe4sCand.chi2ITSMoth, hyphe4sCand.itsClusterSizeMoth, hyphe4sCand.itsClusterSizeDau,
583-
hyphe4sCand.nSigmaTPCDau,
597+
hyphe4sCand.nSigmaTPCDau, hyphe4sCand.nSigmaITSDau,
584598
hyphe4sCand.isSignal, hyphe4sCand.isSignalReco, hyphe4sCand.isCollReco, hyphe4sCand.isSurvEvSelection,
585599
hyphe4sCand.trueDecVtx[0], hyphe4sCand.trueDecVtx[1], hyphe4sCand.trueDecVtx[2],
586600
hyphe4sCand.gMomMoth[0], hyphe4sCand.gMomMoth[1], hyphe4sCand.gMomMoth[2],
@@ -614,7 +628,7 @@ struct Hyperhelium4sigmaRecoTask {
614628
-1, -1, -1,
615629
-1, -1, -1,
616630
-1, -1, -1,
617-
-1,
631+
-1, -1,
618632
true, false, isReconstructedMCCollisions[mcparticle.mcCollisionId()], isSelectedMCCollisions[mcparticle.mcCollisionId()],
619633
hyphe4sCand.trueDecVtx[0], hyphe4sCand.trueDecVtx[1], hyphe4sCand.trueDecVtx[2],
620634
hyphe4sCand.gMomMoth[0], hyphe4sCand.gMomMoth[1], hyphe4sCand.gMomMoth[2],
@@ -651,6 +665,7 @@ struct Hyperhelium4sigmaQa {
651665
const AxisSpec ctAxis{ctBins, "c#it{t} (cm)"};
652666
const AxisSpec rigidityAxis{rigidityBins, "p/z (GeV/#it{c})"};
653667
const AxisSpec nsigmaAxis{nsigmaBins, "TPC n#sigma"};
668+
const AxisSpec itsnsigmaAxis{nsigmaBins, "ITS n#sigma"};
654669
const AxisSpec invMassAxis{invMassBins, "Inv Mass (GeV/#it{c}^{2})"};
655670
const AxisSpec diffPtAxis{200, -10.f, 10.f, "#Delta p_{T} (GeV/#it{c})"};
656671
const AxisSpec diffPzAxis{200, -10.f, 10.f, "#Delta p_{z} (GeV/#it{c})"};
@@ -731,9 +746,11 @@ struct Hyperhelium4sigmaQa {
731746
}
732747

733748
recoQAHist.add<TH1>("hMotherIsPVContributer", "", HistType::kTH1F, {{2, 0.f, 2.f}});
749+
recoQAHist.add<TH1>("hMotherITSCls", "", HistType::kTH1F, {{8, 0.f, 8.f}});
734750
recoQAHist.add<TH1>("hDauAlphaIsPVContributer", "", HistType::kTH1F, {{2, 0.f, 2.f}});
735-
recoQAHist.add<TH2>("hDauAlphaPVsITSNSigma", "", HistType::kTH2F, {pAxis, nsigmaAxis});
736-
recoQAHist.add<TH2>("h2BCandDauAlphaPVsITSNSigma", "", HistType::kTH2F, {pAxis, nsigmaAxis});
751+
recoQAHist.add<TH1>("hDauAlphaITSCls", "", HistType::kTH1F, {{8, 0.f, 8.f}});
752+
recoQAHist.add<TH2>("hDauAlphaPVsITSNSigma", "", HistType::kTH2F, {rigidityAxis, itsnsigmaAxis});
753+
recoQAHist.add<TH2>("h2BCandDauAlphaPVsITSNSigma", "", HistType::kTH2F, {rigidityAxis, nsigmaAxis});
737754
recoQAHist.add<TH1>("hReco2BCandidateCount", "", HistType::kTH1F, {{4, 0.f, 4.f}});
738755
}
739756
}
@@ -975,12 +992,14 @@ struct Hyperhelium4sigmaQa {
975992

976993
recoQAHist.fill(HIST("hReco2BCandidateCount"), 0.5);
977994
recoQAHist.fill(HIST("hRecoMotherCounter"), 0.5);
995+
recoQAHist.fill(HIST("hMotherITSCls"), motherTrack.itsNCls());
978996
recoQAHist.fill(HIST("hRecoDauAlphaCounter"), 0.5);
979997
recoQAHist.fill(HIST("hMotherIsPVContributer"), motherTrack.isPVContributor() ? 1.5 : 0.5);
980998
recoQAHist.fill(HIST("hDauAlphaIsPVContributer"), daughterTrack.isPVContributor() ? 1.5 : 0.5);
981999

982-
float itsNSigma = itsResponse.nSigmaITS<o2::track::PID::Alpha>(daughterTrack.itsClusterSizes(), 2 * daughterTrack.p(), daughterTrack.eta());
983-
recoQAHist.fill(HIST("hDauAlphaPVsITSNSigma"), 2 * daughterTrack.p(), itsNSigma);
1000+
float itsNSigma = itsResponse.nSigmaITS<o2::track::PID::Alpha>(daughterTrack);
1001+
recoQAHist.fill(HIST("hDauAlphaPVsITSNSigma"), daughterTrack.sign() * daughterTrack.p(), itsNSigma);
1002+
recoQAHist.fill(HIST("hDauAlphaITSCls"), daughterTrack.itsNCls());
9841003

9851004
if (motherTrack.has_collision() && daughterTrack.has_collision()) {
9861005
recoQAHist.fill(HIST("hReco2BCandidateCount"), 1.5);
@@ -991,7 +1010,7 @@ struct Hyperhelium4sigmaQa {
9911010

9921011
if (isMoth && isDaug) {
9931012
recoQAHist.fill(HIST("hReco2BCandidateCount"), 3.5);
994-
recoQAHist.fill(HIST("h2BCandDauAlphaPVsITSNSigma"), 2 * daughterTrack.p(), itsNSigma);
1013+
recoQAHist.fill(HIST("h2BCandDauAlphaPVsITSNSigma"), daughterTrack.sign() * daughterTrack.p(), itsNSigma);
9951014
}
9961015
}
9971016
}

0 commit comments

Comments
 (0)