Skip to content

Commit e41b4fa

Browse files
committed
Added a filter to prevent track mixing effects based on the radial distance method
1 parent b2a1496 commit e41b4fa

File tree

2 files changed

+127
-27
lines changed

2 files changed

+127
-27
lines changed

PWGCF/MultiparticleCorrelations/Tasks/CMakeLists.txt

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -19,7 +19,7 @@ o2physics_add_dpl_workflow(multiparticle-correlations-ar
1919
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::PWGCFCore
2020
COMPONENT_NAME Analysis)
2121

22-
o2physics_add_dpl_workflow(threeparticle-correlations
22+
o2physics_add_dpl_workflow(three-particle-correlations
2323
SOURCES threeParticleCorrelations.cxx
2424
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::PWGCFCore
2525
COMPONENT_NAME Analysis)

PWGCF/MultiparticleCorrelations/Tasks/threeParticleCorrelations.cxx

Lines changed: 126 additions & 26 deletions
Original file line numberDiff line numberDiff line change
@@ -20,10 +20,11 @@
2020
#include "CCDB/BasicCCDBManager.h"
2121
#include "Common/DataModel/Centrality.h"
2222
#include "Common/DataModel/PIDResponse.h"
23+
#include "DataFormatsParameters/GRPMagField.h"
2324
#include "PWGLF/DataModel/LFStrangenessTables.h"
2425

25-
#include "TPDGCode.h"
2626
#include "RecoDecay.h"
27+
#include "TPDGCode.h"
2728

2829
using namespace o2;
2930
using namespace o2::framework;
@@ -36,6 +37,7 @@ struct ThreeParticleCorrelations {
3637
HistogramRegistry rMECorrRegistry{"MECorrRegistry", {}, OutputObjHandlingPolicy::AnalysisObject, false, true};
3738
HistogramRegistry rSECorrRegistry{"SECorrRegistry", {}, OutputObjHandlingPolicy::AnalysisObject, false, true};
3839
HistogramRegistry rMCRegistry{"MCRegistry", {}, OutputObjHandlingPolicy::AnalysisObject, false, true};
40+
HistogramRegistry rPhiStarRegistry{"PhiStarRegistry", {}, OutputObjHandlingPolicy::AnalysisObject, false, true};
3941
HistogramRegistry rQARegistry{"QARegistry", {}, OutputObjHandlingPolicy::AnalysisObject, false, true};
4042

4143
// Collision & Event filters
@@ -93,7 +95,7 @@ struct ThreeParticleCorrelations {
9395
TH1D** hEffPions = new TH1D*[2];
9496
TH1D** hEffKaons = new TH1D*[2];
9597
TH1D** hEffProtons = new TH1D*[2];
96-
98+
9799
// Correlation variables
98100
int triggSign;
99101
double candMass;
@@ -134,13 +136,26 @@ struct ThreeParticleCorrelations {
134136
rQARegistry.add("hNSigmaKaon", "hNSigmaKaon", {HistType::kTH2D, {{201, -5.025, 5.025}, {201, -5.025, 5.025}}});
135137
rQARegistry.add("hNSigmaProton", "hNSigmaProton", {HistType::kTH2D, {{201, -5.025, 5.025}, {201, -5.025, 5.025}}});
136138

139+
rQARegistry.add("hTPCPion", "hTPCPion", {HistType::kTH2D, {{trackPtAxis}, {1000, -50, 50}}});
140+
rQARegistry.add("hTPCKaon", "hTPCKaon", {HistType::kTH2D, {{trackPtAxis}, {1000, -50, 50}}});
141+
rQARegistry.add("hTPCProton", "hTPCProton", {HistType::kTH2D, {{trackPtAxis}, {1000, -50, 50}}});
137142
rQARegistry.add("hTOFPion", "hTOFPion", {HistType::kTH2D, {{trackPtAxis}, {1000, -50, 50}}});
138143
rQARegistry.add("hTOFKaon", "hTOFKaon", {HistType::kTH2D, {{trackPtAxis}, {1000, -50, 50}}});
139144
rQARegistry.add("hTOFProton", "hTOFProton", {HistType::kTH2D, {{trackPtAxis}, {1000, -50, 50}}});
140145

141146
rQARegistry.add("hInvMassLambda", "hInvMassLambda", {HistType::kTH3D, {{lambdaInvMassAxis}, {v0PtAxis}, {centralityAxis}}});
142147
rQARegistry.add("hInvMassAntiLambda", "hInvMassAntiLambda", {HistType::kTH3D, {{lambdaInvMassAxis}, {v0PtAxis}, {centralityAxis}}});
143148

149+
// PhiStar
150+
rPhiStarRegistry.add("hSEProtonPreCut", "hSEProtonPreCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
151+
rPhiStarRegistry.add("hSEProtonPostCut", "hSEProtonPostCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
152+
rPhiStarRegistry.add("hSEPionPreCut", "hSEPionPreCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
153+
rPhiStarRegistry.add("hSEPionPostCut", "hSEPionPostCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
154+
rPhiStarRegistry.add("hMEProtonPreCut", "hMEProtonPreCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
155+
rPhiStarRegistry.add("hMEProtonPostCut", "hMEProtonPostCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
156+
rPhiStarRegistry.add("hMEPionPreCut", "hMEPionPreCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
157+
rPhiStarRegistry.add("hMEPionPostCut", "hMEPionPostCut", {HistType::kTH2D, {{80, -0.2, 0.2}, {40, -0.1, 0.1}}});
158+
144159
// Efficiency
145160
rMCRegistry.add("hGenerated", "hGenerated", {HistType::kTH1D, {trackPtAxis}});
146161
rMCRegistry.add("hGenPionP", "hGenPionP", {HistType::kTH1D, {trackPtAxis}});
@@ -199,14 +214,19 @@ struct ThreeParticleCorrelations {
199214

200215
//================================================================================================================================================================================================================
201216

202-
void processSame(MyFilteredCollision const& collision, MyFilteredV0s const& v0s, MyFilteredTracks const& tracks)
217+
void processSame(MyFilteredCollision const& collision, MyFilteredV0s const& v0s, MyFilteredTracks const& tracks, aod::BCsWithTimestamps const&)
203218
{
204219

220+
auto bc = collision.bc_as<aod::BCsWithTimestamps>();
221+
auto bField = getMagneticField(bc.timestamp());
205222
rQARegistry.fill(HIST("hEventCentrality"), collision.centFT0C());
206223
rQARegistry.fill(HIST("hEventZvtx"), collision.posZ());
207224

208225
// Start of the Track QA
209226
for (const auto& track : tracks) {
227+
rQARegistry.fill(HIST("hTPCPion"), track.pt(), track.tpcNSigmaPi());
228+
rQARegistry.fill(HIST("hTPCKaon"), track.pt(), track.tpcNSigmaKa());
229+
rQARegistry.fill(HIST("hTPCProton"), track.pt(), track.tpcNSigmaPr());
210230
if (track.hasTOF()) {
211231
rQARegistry.fill(HIST("hTOFPion"), track.pt(), track.tofNSigmaPi());
212232
rQARegistry.fill(HIST("hTOFKaon"), track.pt(), track.tofNSigmaKa());
@@ -252,7 +272,7 @@ struct ThreeParticleCorrelations {
252272

253273
for (const auto& associate : tracks) {
254274
if (trackFilters(associate)) {
255-
if (correlationFilters(trigger, associate) && fakeV0Filter(trigger, associate)) {
275+
if (correlationFilters(trigger, associate) && radialDistanceFilter(trigger, associate, bField, false) && fakeV0Filter(trigger, associate)) {
256276

257277
assocPID = trackPID(associate);
258278
deltaPhi = RecoDecay::constrainAngle(trigger.phi() - associate.phi(), -constants::math::PIHalf);
@@ -283,14 +303,17 @@ struct ThreeParticleCorrelations {
283303
// End of the V0-Track Correlations
284304
}
285305

286-
void processMixed(MyFilteredCollisions const&, MyFilteredV0s const&, MyFilteredTracks const&)
306+
void processMixed(MyFilteredCollisions const&, MyFilteredV0s const&, MyFilteredTracks const&, aod::BCsWithTimestamps const&)
287307
{
288-
308+
289309
// Start of the Mixed-events Correlations
290310
for (const auto& [coll_1, v0_1, coll_2, track_2] : pair) {
311+
312+
auto bc = coll_1.bc_as<aod::BCsWithTimestamps>();
313+
auto bField = getMagneticField(bc.timestamp());
291314
for (const auto& [trigger, associate] : soa::combinations(soa::CombinationsFullIndexPolicy(v0_1, track_2))) {
292315
if (v0Filters(trigger) && trackFilters(associate)) {
293-
if (correlationFilters(trigger, associate) && fakeV0Filter(trigger, associate)) {
316+
if (radialDistanceFilter(trigger, associate, bField, true) && fakeV0Filter(trigger, associate)) {
294317

295318
triggSign = v0Sign(trigger);
296319
if (triggSign == 1) {
@@ -438,6 +461,22 @@ struct ThreeParticleCorrelations {
438461

439462
//================================================================================================================================================================================================================
440463

464+
double getMagneticField(uint64_t timestamp)
465+
{
466+
static parameters::GRPMagField* grpo = nullptr;
467+
if (grpo == nullptr) {
468+
grpo = ccdb->getForTimeStamp<parameters::GRPMagField>("GLO/Config/GRPMagField", timestamp);
469+
if (grpo == nullptr) {
470+
LOGF(fatal, "GRP object not found for timestamp %llu", timestamp);
471+
return 0;
472+
}
473+
LOGF(info, "Retrieved GRP for timestamp %llu with magnetic field of %d kG", timestamp, grpo->getNominalL3Field());
474+
}
475+
476+
return 0.1 * (grpo->getNominalL3Field()); // 1 T = 10 kG
477+
}
478+
479+
441480
double trackEff(TH1D** efficiencies, int sign, double pT)
442481
{
443482

@@ -470,25 +509,30 @@ struct ThreeParticleCorrelations {
470509

471510
static double pid[2]; // {PID, NSigma}
472511

473-
double nSigma[3];
474-
double nSigmaTOF[3];
475-
nSigmaTOF[0] = track.tofNSigmaPi();
476-
nSigmaTOF[1] = track.tofNSigmaKa();
477-
nSigmaTOF[2] = track.tofNSigmaPr();
512+
double nSigmaTPC[3], nSigma[3];
513+
double nSigmaTOF[3] = {0.0, 0.0, 0.0};
514+
nSigmaTPC[0] = track.tpcNSigmaPi();
515+
nSigmaTPC[1] = track.tpcNSigmaKa();
516+
nSigmaTPC[2] = track.tpcNSigmaPr();
517+
if (track.hasTOF()) {
518+
nSigmaTOF[0] = track.tofNSigmaPi();
519+
nSigmaTOF[1] = track.tofNSigmaKa();
520+
nSigmaTOF[2] = track.tofNSigmaPr();
521+
}
478522

479-
nSigma[0] = std::abs(nSigmaTOF[0]);
480-
nSigma[1] = std::abs(nSigmaTOF[1]);
481-
nSigma[2] = std::abs(nSigmaTOF[2]);
523+
nSigma[0] = std::sqrt(std::pow(nSigmaTPC[0], 2) + std::pow(nSigmaTOF[0], 2));
524+
nSigma[1] = std::sqrt(std::pow(nSigmaTPC[1], 2) + std::pow(nSigmaTOF[1], 2));
525+
nSigma[2] = std::sqrt(std::pow(nSigmaTPC[2], 2) + std::pow(nSigmaTOF[2], 2));
482526

483527
if (nSigma[0] <= std::min(nSigma[1], nSigma[2])) { // Pions
484528
pid[0] = 0.0;
485-
pid[1] = nSigmaTOF[0];
529+
pid[1] = nSigma[0];
486530
} else if (nSigma[1] <= std::min(nSigma[0], nSigma[2])) { // Kaons
487531
pid[0] = 1.0;
488-
pid[1] = nSigmaTOF[1];
532+
pid[1] = nSigma[1];
489533
} else if (nSigma[2] < std::min(nSigma[0], nSigma[1])) { // Protons
490534
pid[0] = 2.0;
491-
pid[1] = nSigmaTOF[2];
535+
pid[1] = nSigma[2];
492536
}
493537

494538
return pid;
@@ -530,11 +574,11 @@ struct ThreeParticleCorrelations {
530574
if (track.pt() < 0.3) {
531575
return kFALSE;
532576
} else if (track.pt() > 0.3 && track.pt() < 1.5) {
533-
if (std::abs(trackPID(track)[1]) > 4.0) {
577+
if (std::abs(track.tofNSigmaPi()) > 4.0) {
534578
return kFALSE;
535579
}
536580
} else if (track.pt() > 1.5 && track.pt() < 2.3) {
537-
if (trackPID(track)[1] < -4.0 || trackPID(track)[1] > 0.0) {
581+
if (track.tofNSigmaPi() < -4.0 ||track.tofNSigmaPi() > 0.0) {
538582
return kFALSE;
539583
}
540584
} else if (track.pt() > 2.3) {
@@ -548,15 +592,15 @@ struct ThreeParticleCorrelations {
548592
if (track.pt() < 0.5) {
549593
return kFALSE;
550594
} else if (track.pt() > 0.5 && track.pt() < 1.5) {
551-
if (std::abs(trackPID(track)[1]) > 4.0) {
595+
if (std::abs(track.tofNSigmaKa()) > 4.0) {
552596
return kFALSE;
553597
}
554598
} else if (track.pt() > 1.5 && track.pt() < 2.0) {
555-
if (trackPID(track)[1] < -2.0 || trackPID(track)[1] > 4.0) {
599+
if (track.tofNSigmaKa() < -2.0 || track.tofNSigmaKa() > 4.0) {
556600
return kFALSE;
557601
}
558602
} else if (track.pt() > 2.0 && track.pt() < 2.5) {
559-
if (trackPID(track)[1] < 0.0 || trackPID(track)[1] > 4.0) {
603+
if (track.tofNSigmaKa() < 0.0 || track.tofNSigmaKa() > 4.0) {
560604
return kFALSE;
561605
}
562606
} else if (track.pt() > 2.5) {
@@ -570,15 +614,15 @@ struct ThreeParticleCorrelations {
570614
if (track.pt() < 0.5) {
571615
return kFALSE;
572616
} else if (track.pt() > 0.5 && track.pt() < 0.7) {
573-
if (trackPID(track)[1] < -2.0 || trackPID(track)[1] > 4.0) {
617+
if (track.tofNSigmaPr() < -2.0 || track.tofNSigmaPr() > 4.0) {
574618
return kFALSE;
575619
}
576620
} else if (track.pt() > 0.7 && track.pt() < 2.5) {
577-
if (std::abs(trackPID(track)[1]) > 4.0) {
621+
if (std::abs(track.tofNSigmaPr()) > 4.0) {
578622
return kFALSE;
579623
}
580624
} else if (track.pt() > 2.5) {
581-
if (trackPID(track)[1] < -2.0 || trackPID(track)[1] > 4.0) {
625+
if (track.tofNSigmaPr() < -2.0 || track.tofNSigmaPr() > 4.0) {
582626
return kFALSE;
583627
}
584628
}
@@ -641,6 +685,62 @@ struct ThreeParticleCorrelations {
641685

642686
return kTRUE;
643687
}
688+
689+
template <class V0Cand, class TrackCand>
690+
bool radialDistanceFilter(const V0Cand& v0, const TrackCand& track, double B, bool Mix)
691+
{
692+
693+
const auto& dProton = v0.template posTrack_as<MyFilteredTracks>();
694+
const auto& dPion = v0.template negTrack_as<MyFilteredTracks>();
695+
if (v0Sign(v0) == -1) {
696+
const auto& dProton = v0.template negTrack_as<MyFilteredTracks>();
697+
const auto& dPion = v0.template posTrack_as<MyFilteredTracks>();
698+
}
699+
700+
double dEtaProton = dProton.eta() - track.eta();
701+
double dEtaPion = dPion.eta() - track.eta();
702+
if (std::abs(dEtaProton) > 0.02 && std::abs(dEtaPion) > 0.02) {
703+
return kTRUE;
704+
}
705+
706+
double phiStarProton, phiStarPion;
707+
double dPhiProton = dProton.phi() - track.phi();
708+
double dPhiPion = dPion.phi() - track.phi();
709+
double phaseProton = (-0.3*B*dProton.sign()) / (2*dProton.pt());
710+
double phasePion = (-0.3*B*dPion.sign()) / (2*dPion.pt());
711+
double phaseTrack = (-0.3*B*track.sign()) / (2*track.pt());
712+
713+
for (double r = 0.8; r <= 2.5; r += 0.01) {
714+
phiStarProton = RecoDecay::constrainAngle(dPhiProton + std::asin(phaseProton*r) - std::asin(phaseTrack*r), -constants::math::PIHalf);
715+
phiStarPion = RecoDecay::constrainAngle(dPhiPion + std::asin(phasePion*r) - std::asin(phaseTrack*r), -constants::math::PIHalf);
716+
717+
if (r == 0.8) {
718+
if (!Mix) {
719+
rPhiStarRegistry.fill(HIST("hSEProtonPreCut"), phiStarProton, dEtaProton);
720+
rPhiStarRegistry.fill(HIST("hSEPionPreCut"), phiStarPion, dEtaPion);
721+
} else {
722+
rPhiStarRegistry.fill(HIST("hMEProtonPreCut"), phiStarProton, dEtaProton);
723+
rPhiStarRegistry.fill(HIST("hMEPionPreCut"), phiStarPion, dEtaPion);
724+
}
725+
}
726+
727+
if ((std::abs(dEtaProton) <= 0.02 && std::abs(phiStarProton) <= 0.045) || (std::abs(dEtaPion) <= 0.02 && std::abs(phiStarPion) <= 0.045)) {
728+
return kFALSE;
729+
}
730+
731+
if (r == 0.8) {
732+
if (!Mix) {
733+
rPhiStarRegistry.fill(HIST("hSEProtonPostCut"), phiStarProton, dEtaProton);
734+
rPhiStarRegistry.fill(HIST("hSEPionPostCut"), phiStarPion, dEtaPion);
735+
} else {
736+
rPhiStarRegistry.fill(HIST("hMEProtonPostCut"), phiStarProton, dEtaProton);
737+
rPhiStarRegistry.fill(HIST("hMEPionPostCut"), phiStarPion, dEtaPion);
738+
}
739+
}
740+
}
741+
742+
return kTRUE;
743+
}
644744
};
645745

646746
//==================================================================================================================================================================================================================

0 commit comments

Comments
 (0)