Skip to content

Commit bed5dbd

Browse files
author
jokonig
committed
[PWGEM] Add QA histograms to Dalitz Skimmer
- In order to monitor the various cuts performed to the electron tracks, the basic track QA histograms can now be filled also before any track selection. In addition, in MC the histograms can also be filled for true primary electrons only
1 parent e82d63e commit bed5dbd

1 file changed

Lines changed: 158 additions & 58 deletions

File tree

‎PWGEM/PhotonMeson/TableProducer/skimmerPrimaryElectronFromDalitzEE.h‎

Lines changed: 158 additions & 58 deletions
Original file line numberDiff line numberDiff line change
@@ -121,6 +121,33 @@ struct skimmerPrimaryElectronFromDalitzEE {
121121
EpEmPairsAndPhoton = 3
122122
};
123123

124+
enum RejectReason {
125+
kAllTracks = 0,
126+
kPassed,
127+
kNoMcParticle,
128+
kNotTrueElectron,
129+
kTOFRequired,
130+
kNoITS,
131+
kITSChi2,
132+
kITSNCls,
133+
kITSNClsInnerBarrel,
134+
kNoTPC,
135+
kTPCChi2,
136+
kTPCNClsFound,
137+
kTPCCrossedRows,
138+
kTPCCrossedRowsOverFindable,
139+
kTPCSharedFraction,
140+
kPropagation,
141+
kDCAxyz,
142+
kDCA3D,
143+
kPtEta,
144+
kITSClusterSize,
145+
kITSsaPt,
146+
kNotElectron,
147+
kNRejectReasons
148+
};
149+
static constexpr std::array<const char*, kNRejectReasons> RejectReasonLabels = {"all tracks", "accepted tracks", "no MC particle", "not true electron", "TOF required", "no ITS", "ITS chi2", "ITS Ncls", "ITS Ncls IB", "no TPC", "TPC chi2", "TPC Ncls found", "TPC crossed rows", "TPC Ncr/Nf", "TPC shared frac", "propagation", "DCAxy/z", "DCA 3D", "pt/eta", "ITS cluster size", "ITSsa pt", "not electron"};
150+
124151
// Configurables
125152
Configurable<std::string> ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
126153
Configurable<std::string> grpPath{"grpPath", "GLO/GRP/GRP", "Path of the grp file"};
@@ -172,9 +199,11 @@ struct skimmerPrimaryElectronFromDalitzEE {
172199
Configurable<float> intercept{"intercept", -0.0380, "intercept for m vs. phiv"};
173200
Configurable<bool> useTOFNSigmaDeltaBC{"useTOFNSigmaDeltaBC", false, "Flag to shift delta BC for TOF n sigma (only with TTCA)"};
174201
Configurable<bool> storeOnlyTrueElectronMC{"storeOnlyTrueElectronMC", false, "Flag to store only true electron in MC"};
202+
Configurable<bool> fillHistogramsBeforeSelection{"fillHistogramsBeforeSelection", false, "Flag to fill track histograms before track selections (only for fillingMode 1)"};
175203

176204
o2::framework::HistogramRegistry fRegistry{"output", {}, o2::framework::OutputObjHandlingPolicy::AnalysisObject, false, false};
177205
static constexpr std::array<std::string_view, 3> DileptonSigns = {"uls/", "lspp/", "lsmm/"};
206+
static constexpr std::array<std::string_view, 4> TrackDirs = {"Track/", "TrackTrueElec/", "TrackNoSel/", "TrackNoSelTrueElec/"};
178207

179208
int mRunNumber = 0;
180209
float dBz = 0.;
@@ -237,6 +266,24 @@ struct skimmerPrimaryElectronFromDalitzEE {
237266
fRegistry.add("Track/hTOFNsigmaEl", "TOF n sigma el;p_{pv} (GeV/c);n #sigma_{e}^{TOF}", kTH2F, {{1000, 0, 10}, {100, -5, +5}}, false);
238267
fRegistry.add("Track/hTOFNsigmaPi", "TOF n sigma pi;p_{pv} (GeV/c);n #sigma_{#pi}^{TOF}", kTH2F, {{1000, 0, 10}, {100, -5, +5}}, false);
239268

269+
fRegistry.add("Track/hTrackRejectReason", "track rejection reason vs. p_{T};;p_{T} (GeV/c)", kTH2F, {{static_cast<int>(kNRejectReasons), -0.5, static_cast<double>(kNRejectReasons) - 0.5}, {200, 0, 20}}, false);
270+
auto hRejectReason = fRegistry.get<TH2>(HIST("Track/hTrackRejectReason"));
271+
for (int i = 0; i < kNRejectReasons; i++) {
272+
hRejectReason->GetXaxis()->SetBinLabel(i + 1, RejectReasonLabels[i]);
273+
}
274+
275+
if (doprocessMC) {
276+
LOG(info) << "adding TrackTrueElec folder";
277+
fRegistry.addClone("Track/", "TrackTrueElec/");
278+
}
279+
280+
if (fillHistogramsBeforeSelection) {
281+
fRegistry.addClone("Track/", "TrackNoSel/");
282+
if (doprocessMC) {
283+
fRegistry.addClone("Track/", "TrackNoSelTrueElec/");
284+
}
285+
}
286+
240287
// pair
241288
fRegistry.add("Pair/uls/hTrackMvsPt", "m_{ee} vs. p_{T,ee};m_{ee} (GeV/c^{2});p_{T,ee} (GeV/c)", kTH2F, {{100, 0, 0.1}, {200, 0, 2}}, false);
242289
fRegistry.add("Pair/uls/hCheckEMvsPt", "m_{ee} vs. p_{T,ee};m_{ee} (GeV/c^{2});p_{T,ee} (GeV/c)", kTH2F, {{100, 0, 0.1}, {200, 0, 2}}, false);
@@ -354,63 +401,72 @@ struct skimmerPrimaryElectronFromDalitzEE {
354401
}
355402

356403
template <bool isMC, o2::soa::is_iterator TCollision, o2::soa::is_iterator TTrack>
357-
bool checkTrack(TCollision const& collision, TTrack const& track)
404+
bool checkTrack(TCollision const& collision, TTrack const& track, int* rejectReason = nullptr)
358405
{
406+
// Lambda function: Sets the reject reason to the correct value and returns false.
407+
// Used for every return statement in the checkTrack function where false is returned
408+
auto reject = [rejectReason](RejectReason reason) {
409+
if (rejectReason != nullptr) {
410+
*rejectReason = reason;
411+
}
412+
return false;
413+
};
414+
359415
if constexpr (isMC) {
360416
if (!track.has_mcParticle()) {
361-
return false;
417+
return reject(kNoMcParticle);
362418
}
363419
if (storeOnlyTrueElectronMC) {
364420
const auto& mcParticle = track.template mcParticle_as<o2::aod::McParticles>();
365421
if (std::abs(mcParticle.pdgCode()) != 11) {
366-
return false;
422+
return reject(kNotTrueElectron);
367423
}
368424
}
369425
}
370426

371427
float tofNSigmaEl = mapTOFNsigmaReassociated[std::make_pair(collision.globalIndex(), track.globalIndex())];
372428
if (requireTOF && !(track.hasTOF() && std::fabs(tofNSigmaEl) < maxTOFNsigmaEl)) {
373-
return false;
429+
return reject(kTOFRequired);
374430
}
375431

376432
if (!track.hasITS()) {
377-
return false;
433+
return reject(kNoITS);
378434
}
379435

380436
if (track.itsChi2NCl() < minchi2its || maxchi2its < track.itsChi2NCl()) {
381-
return false;
437+
return reject(kITSChi2);
382438
}
383439

384440
if (track.itsNCls() < min_ncluster_its) {
385-
return false;
441+
return reject(kITSNCls);
386442
}
387443
if (track.itsNClsInnerBarrel() < min_ncluster_itsib) {
388-
return false;
444+
return reject(kITSNClsInnerBarrel);
389445
}
390446

391447
if (!includeITSsa && !track.hasTPC()) {
392-
return false;
448+
return reject(kNoTPC);
393449
}
394450

395451
if (track.hasTPC()) {
396452
if (track.tpcChi2NCl() < minchi2tpc || maxchi2tpc < track.tpcChi2NCl()) {
397-
return false;
453+
return reject(kTPCChi2);
398454
}
399455

400456
if (track.tpcNClsFound() < min_ncluster_tpc) {
401-
return false;
457+
return reject(kTPCNClsFound);
402458
}
403459

404460
if (track.tpcNClsCrossedRows() < mincrossedrows) {
405-
return false;
461+
return reject(kTPCCrossedRows);
406462
}
407463

408464
if (track.tpcCrossedRowsOverFindableCls() < min_tpc_cr_findable_ratio) {
409-
return false;
465+
return reject(kTPCCrossedRowsOverFindable);
410466
}
411467

412468
if (track.tpcFractionSharedCls() > max_frac_shared_clusters_tpc) {
413-
return false;
469+
return reject(kTPCSharedFraction);
414470
}
415471
}
416472

@@ -422,13 +478,13 @@ struct skimmerPrimaryElectronFromDalitzEE {
422478
mVtx.setCov(collision.covXX(), collision.covXY(), collision.covYY(), collision.covXZ(), collision.covYZ(), collision.covZZ());
423479
bool isPropOK = o2::base::Propagator::Instance()->propagateToDCABxByBz(mVtx, trackParCov, 2.f, matCorr, &mDcaInfoCov);
424480
if (!isPropOK) {
425-
return false;
481+
return reject(kPropagation);
426482
}
427483
float dcaXY = mDcaInfoCov.getY();
428484
float dcaZ = mDcaInfoCov.getZ();
429485

430486
if (std::fabs(dcaXY) > dca_xy_max || std::fabs(dcaZ) > dca_z_max) {
431-
return false;
487+
return reject(kDCAxyz);
432488
}
433489

434490
float dca3D = 999.f;
@@ -440,11 +496,11 @@ struct skimmerPrimaryElectronFromDalitzEE {
440496
dca3D = std::sqrt(std::fabs(chi2) / 2.);
441497
}
442498
if (dca3D > dca_3d_sigma_max) {
443-
return false;
499+
return reject(kDCA3D);
444500
}
445501

446502
if (trackParCov.getPt() < minpt || std::fabs(trackParCov.getEta()) > maxeta) {
447-
return false;
503+
return reject(kPtEta);
448504
}
449505

450506
int total_cluster_size = 0, nl = 0;
@@ -457,11 +513,11 @@ struct skimmerPrimaryElectronFromDalitzEE {
457513
}
458514

459515
if (maxMeanITSClusterSize < static_cast<float>(total_cluster_size) / static_cast<float>(nl) * std::cos(std::atan(trackParCov.getTgl()))) {
460-
return false;
516+
return reject(kITSClusterSize);
461517
}
462518

463519
if ((track.hasITS() && !track.hasTPC() && !track.hasTRD() && !track.hasTOF()) && maxpt_itssa < track.pt()) {
464-
return false;
520+
return reject(kITSsaPt);
465521
}
466522

467523
// if ((track.hasITS() && !track.hasTPC() && !track.hasTRD() && !track.hasTOF()) && maxpt_itssa < track.pt()) {
@@ -528,11 +584,24 @@ struct skimmerPrimaryElectronFromDalitzEE {
528584
{
529585

530586
for (const auto& track : tracks) {
531-
if (!checkTrack<isMC>(collision, track)) {
532-
continue;
587+
588+
if (fillHistogramsBeforeSelection) {
589+
fillTrackHistograms<isMC, true>(track);
533590
}
534591

535-
if (!isElectron(collision, track)) {
592+
fRegistry.fill(HIST("Track/hTrackRejectReason"), kAllTracks, track.pt());
593+
int rejectReason = kPassed;
594+
if (checkTrack<isMC>(collision, track, &rejectReason) && !isElectron(collision, track)) {
595+
rejectReason = kNotElectron;
596+
}
597+
fRegistry.fill(HIST("Track/hTrackRejectReason"), rejectReason, track.pt());
598+
if constexpr (isMC) {
599+
if (isTruePrimaryElectron(track)) {
600+
fRegistry.fill(HIST("TrackTrueElec/hTrackRejectReason"), kAllTracks, track.pt());
601+
fRegistry.fill(HIST("TrackTrueElec/hTrackRejectReason"), rejectReason, track.pt());
602+
}
603+
}
604+
if (rejectReason != kPassed) {
536605
continue;
537606
}
538607

@@ -544,7 +613,16 @@ struct skimmerPrimaryElectronFromDalitzEE {
544613
template <bool isMC, int pairtype, o2::soa::is_iterator TCollision, o2::soa::is_table TTracks1, o2::soa::is_table TTracks2>
545614
void fillPairInfo(TCollision const& collision, TTracks1 const& tracks1, TTracks2 const& tracks2)
546615
{
547-
if constexpr (pairtype == 0) { // ULS
616+
if constexpr (pairtype == 0) { // ULS
617+
if (fillHistogramsBeforeSelection) { // we cannot fill the track histograms in the pair loop as they would get filled multiple times. Hence we loop over all tracks here
618+
for (const auto& track : tracks1) {
619+
fillTrackHistograms<isMC, true>(track);
620+
}
621+
for (const auto& track : tracks2) {
622+
fillTrackHistograms<isMC, true>(track);
623+
}
624+
}
625+
548626
for (const auto& [t1, t2] : combinations(o2::soa::CombinationsFullIndexPolicy(tracks1, tracks2))) {
549627

550628
if (!checkTrack<isMC>(collision, t1) || !checkTrack<isMC>(collision, t2)) {
@@ -621,44 +699,66 @@ struct skimmerPrimaryElectronFromDalitzEE {
621699
}
622700
}
623701

624-
template <bool isMC, o2::soa::is_iterator TTrack>
702+
template <o2::soa::is_iterator TTrack>
703+
bool isTruePrimaryElectron(TTrack const& track)
704+
{
705+
if (!track.has_mcParticle()) {
706+
return false;
707+
}
708+
const auto mcParticle = track.template mcParticle_as<o2::aod::McParticles>();
709+
return std::abs(mcParticle.pdgCode()) == 11 && mcParticle.isPhysicalPrimary(); // o2-linter: disable=pdg/explicit-code (electron)
710+
}
711+
712+
template <bool isMC, bool beforeSelection = false, o2::soa::is_iterator TTrack>
625713
void fillTrackHistograms(TTrack const& track)
714+
{
715+
constexpr int kDirAllTracks = beforeSelection ? 2 : 0;
716+
fillTrackHistogramsImpl<isMC, kDirAllTracks>(track);
717+
if constexpr (isMC) {
718+
if (isTruePrimaryElectron(track)) {
719+
fillTrackHistogramsImpl<isMC, kDirAllTracks + 1>(track);
720+
}
721+
}
722+
}
723+
724+
template <bool isMC, int trackType, o2::soa::is_iterator TTrack>
725+
void fillTrackHistogramsImpl(TTrack const& track)
626726
{
627727
float mcTunedTPCSignal = 0.f;
628728
if constexpr (isMC) {
629729
mcTunedTPCSignal = track.mcTunedTPCSignal();
630730
}
631731

632-
fRegistry.fill(HIST("Track/hPt"), track.pt());
633-
fRegistry.fill(HIST("Track/hEtaPhi"), track.phi(), track.eta());
634-
fRegistry.fill(HIST("Track/hQoverPt"), track.sign() / track.pt());
635-
fRegistry.fill(HIST("Track/hRelDeltaPt"), track.pt(), track.sigma1Pt() * track.pt());
636-
fRegistry.fill(HIST("Track/hDCAxyz"), track.dcaXY(), track.dcaZ());
637-
fRegistry.fill(HIST("Track/hDCAxy_Pt"), track.pt(), track.dcaXY());
638-
fRegistry.fill(HIST("Track/hDCAz_Pt"), track.pt(), track.dcaZ());
639-
fRegistry.fill(HIST("Track/hDCAxyzSigma"), track.dcaXY() / std::sqrt(track.cYY()), track.dcaZ() / std::sqrt(track.cZZ()));
640-
fRegistry.fill(HIST("Track/hDCAxyRes_Pt"), track.pt(), std::sqrt(track.cYY()) * 1e+4); // convert cm to um
641-
fRegistry.fill(HIST("Track/hDCAzRes_Pt"), track.pt(), std::sqrt(track.cZZ()) * 1e+4); // convert cm to um
642-
643-
fRegistry.fill(HIST("Track/hNclsTPC"), track.tpcNClsFound());
644-
fRegistry.fill(HIST("Track/hNcrTPC"), track.tpcNClsCrossedRows());
645-
fRegistry.fill(HIST("Track/hChi2TPC"), track.tpcChi2NCl());
646-
fRegistry.fill(HIST("Track/hTPCNcr2Nf"), track.tpcCrossedRowsOverFindableCls());
647-
fRegistry.fill(HIST("Track/hTPCNcls2Nf"), track.tpcFoundOverFindableCls());
648-
fRegistry.fill(HIST("Track/hTPCNclsShared"), track.pt(), track.tpcFractionSharedCls());
649-
fRegistry.fill(HIST("Track/hTPCdEdx"), track.tpcInnerParam(), track.tpcSignal());
650-
fRegistry.fill(HIST("Track/hTPCdEdxMC"), track.tpcInnerParam(), mcTunedTPCSignal);
651-
fRegistry.fill(HIST("Track/hTPCNsigmaEl"), track.tpcInnerParam(), track.tpcNSigmaEl());
652-
fRegistry.fill(HIST("Track/hTPCNsigmaPi"), track.tpcInnerParam(), track.tpcNSigmaPi());
653-
654-
fRegistry.fill(HIST("Track/hChi2TOF"), track.tofChi2());
655-
fRegistry.fill(HIST("Track/hTOFbeta"), track.p(), track.beta());
656-
fRegistry.fill(HIST("Track/hTOFNsigmaEl"), track.p(), track.tofNSigmaEl());
657-
fRegistry.fill(HIST("Track/hTOFNsigmaPi"), track.p(), track.tofNSigmaPi());
658-
659-
fRegistry.fill(HIST("Track/hNclsITS"), track.itsNCls());
660-
fRegistry.fill(HIST("Track/hChi2ITS"), track.itsChi2NCl());
661-
fRegistry.fill(HIST("Track/hITSClusterMap"), track.itsClusterMap());
732+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hPt"), track.pt());
733+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hEtaPhi"), track.phi(), track.eta());
734+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hQoverPt"), track.sign() / track.pt());
735+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hRelDeltaPt"), track.pt(), track.sigma1Pt() * track.pt());
736+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hDCAxyz"), track.dcaXY(), track.dcaZ());
737+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hDCAxy_Pt"), track.pt(), track.dcaXY());
738+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hDCAz_Pt"), track.pt(), track.dcaZ());
739+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hDCAxyzSigma"), track.dcaXY() / std::sqrt(track.cYY()), track.dcaZ() / std::sqrt(track.cZZ()));
740+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hDCAxyRes_Pt"), track.pt(), std::sqrt(track.cYY()) * 1e+4); // convert cm to um
741+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hDCAzRes_Pt"), track.pt(), std::sqrt(track.cZZ()) * 1e+4); // convert cm to um
742+
743+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hNclsTPC"), track.tpcNClsFound());
744+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hNcrTPC"), track.tpcNClsCrossedRows());
745+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hChi2TPC"), track.tpcChi2NCl());
746+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTPCNcr2Nf"), track.tpcCrossedRowsOverFindableCls());
747+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTPCNcls2Nf"), track.tpcFoundOverFindableCls());
748+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTPCNclsShared"), track.pt(), track.tpcFractionSharedCls());
749+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTPCdEdx"), track.tpcInnerParam(), track.tpcSignal());
750+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTPCdEdxMC"), track.tpcInnerParam(), mcTunedTPCSignal);
751+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTPCNsigmaEl"), track.tpcInnerParam(), track.tpcNSigmaEl());
752+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTPCNsigmaPi"), track.tpcInnerParam(), track.tpcNSigmaPi());
753+
754+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hChi2TOF"), track.tofChi2());
755+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTOFbeta"), track.p(), track.beta());
756+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTOFNsigmaEl"), track.p(), track.tofNSigmaEl());
757+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hTOFNsigmaPi"), track.p(), track.tofNSigmaPi());
758+
759+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hNclsITS"), track.itsNCls());
760+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hChi2ITS"), track.itsChi2NCl());
761+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hITSClusterMap"), track.itsClusterMap());
662762

663763
int totalClusterSize = 0, nl = 0;
664764
for (unsigned int layer = 0; layer < 7; layer++) { // o2-linter: disable=magic-number (number of ITS layers)
@@ -686,9 +786,9 @@ struct skimmerPrimaryElectronFromDalitzEE {
686786
}
687787
totalClusterSizeOb += clusterSizePerLayer;
688788
}
689-
fRegistry.fill(HIST("Track/hMeanClusterSizeITS"), track.p(), static_cast<float>(totalClusterSize) / static_cast<float>(nl) * std::cos(std::atan(track.tgl())));
690-
fRegistry.fill(HIST("Track/hMeanClusterSizeITSib"), track.p(), static_cast<float>(totalClusterSizeIb) / static_cast<float>(nl_ib) * std::cos(std::atan(track.tgl())));
691-
fRegistry.fill(HIST("Track/hMeanClusterSizeITSob"), track.p(), static_cast<float>(totalClusterSizeOb) / static_cast<float>(nlOb) * std::cos(std::atan(track.tgl())));
789+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hMeanClusterSizeITS"), track.p(), static_cast<float>(totalClusterSize) / static_cast<float>(nl) * std::cos(std::atan(track.tgl())));
790+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hMeanClusterSizeITSib"), track.p(), static_cast<float>(totalClusterSizeIb) / static_cast<float>(nl_ib) * std::cos(std::atan(track.tgl())));
791+
fRegistry.fill(HIST(TrackDirs[trackType]) + HIST("hMeanClusterSizeITSob"), track.p(), static_cast<float>(totalClusterSizeOb) / static_cast<float>(nlOb) * std::cos(std::atan(track.tgl())));
692792
}
693793

694794
template <bool isMC, o2::soa::is_iterator TCollision, o2::soa::is_iterator TTrack>

0 commit comments

Comments
 (0)