Skip to content

Commit d4468c3

Browse files
SwatiSaha-1997SwatiSaha-1997
committed
Modified event selection and PID cuts
1 parent fb51c55 commit d4468c3

1 file changed

Lines changed: 178 additions & 32 deletions

File tree

‎PWGCF/EbyEFluctuations/Tasks/v0ptHadPiKaProt.cxx‎

Lines changed: 178 additions & 32 deletions
Original file line numberDiff line numberDiff line change
@@ -73,6 +73,8 @@ struct V0ptHadPiKaProt {
7373
Configurable<int64_t> ccdbNoLaterThan{"ccdbNoLaterThan", std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now().time_since_epoch()).count(), "latest acceptable timestamp of creation for the object"};
7474
Configurable<std::string> ccdbUrl{"ccdbUrl", "https://alice-ccdb.cern.ch", "url of the ccdb repository"};
7575
Configurable<std::string> ccdbPath{"ccdbPath", "Users/s/swati/PhiWeight", "CCDB path to ccdb object containing phi weight in a 3D histogram"};
76+
Configurable<int64_t> ccdbNoLaterThanPtEff{"ccdbNoLaterThanPtEff", std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now().time_since_epoch()).count(), "latest acceptable timestamp of creation for the object"};
77+
Configurable<std::string> ccdbPathPtEff{"ccdbPathPtEff", "Users/s/swati/EfficiencyWeight", "CCDB path to ccdb object containing pt-dependent efficiency"};
7678

7779
enum Particles {
7880
PIONS = 0,
@@ -110,6 +112,10 @@ struct V0ptHadPiKaProt {
110112
Configurable<float> cfgnSigmaOtherParticles{"cfgnSigmaOtherParticles", 3.0f, "PID nSigma cut to remove other particles (default:3)"};
111113
Configurable<float> cfgnSigmaCutTPC{"cfgnSigmaCutTPC", 2.0f, "PID nSigma cut for TPC"};
112114
Configurable<float> cfgnSigmaCutTOF{"cfgnSigmaCutTOF", 2.0f, "PID nSigma cut for TOF"};
115+
Configurable<bool> cfgUseNewSeperationPid{"cfgUseNewSeperationPid", true, "Use seperation based PID cuts (NEW)"};
116+
Configurable<float> cfgnSigmaCutTPCHigherPt{"cfgnSigmaCutTPCHigherPt", 2.0f, "PID nSigma cut for TPC at higher pt"};
117+
Configurable<float> cfgnSigmaCutTOFHigherPt{"cfgnSigmaCutTOFHigherPt", 2.0f, "PID nSigma cut for TOF at higher pt"};
118+
Configurable<float> cfgnSigmaSeperationCut{"cfgnSigmaSeperationCut", 3.5f, "PID nSigma of other species must be greater than the vale"};
113119
Configurable<float> cfgnSigmaCutCombTPCTOF{"cfgnSigmaCutCombTPCTOF", 2.0f, "PID nSigma combined cut for TPC and TOF"};
114120
ConfigurableAxis nchAxis{"nchAxis", {5000, 0.5, 5000.5}, ""};
115121
ConfigurableAxis centAxis{"centAxis", {90, 0., 90.}, "Centrality/Multiplicity percentile bining"};
@@ -121,15 +127,20 @@ struct V0ptHadPiKaProt {
121127
Configurable<float> cfgCutPtUpper{"cfgCutPtUpper", 10.0f, "Higher pT cut for inclusive hadron analysis"};
122128
Configurable<float> cfgCutPtUpperPID{"cfgCutPtUpperPID", 6.0f, "Higher pT cut for identified particle analysis"};
123129
Configurable<float> cfgCutEta{"cfgCutEta", 0.8f, "absolute Eta cut"};
124-
Configurable<float> cfgCutEtaLeft{"cfgCutEtaLeft", 0.8f, "Left end of eta gap"};
125-
Configurable<float> cfgCutEtaRight{"cfgCutEtaRight", 0.8f, "Right end of eta gap"};
130+
Configurable<float> cfgCutEtaLeft{"cfgCutEtaLeft", -0.4f, "Left end of eta gap"};
131+
Configurable<float> cfgCutEtaRight{"cfgCutEtaRight", 0.4f, "Right end of eta gap"};
126132
Configurable<int> cfgNSubsample{"cfgNSubsample", 20, "Number of subsamples"};
127133
Configurable<int> cfgCentralityChoice{"cfgCentralityChoice", 0, "Which centrality estimator? 0-->FT0C, 1-->FT0A, 2-->FT0M, 3-->FV0A"};
128134
Configurable<bool> cfgEvSelkNoSameBunchPileup{"cfgEvSelkNoSameBunchPileup", true, "Pileup removal"};
129135
Configurable<bool> cfgUseGoodITSLayerAllCut{"cfgUseGoodITSLayerAllCut", true, "Remove time interval with dead ITS zone"};
130136
Configurable<bool> cfgEvSelkNoITSROFrameBorder{"cfgEvSelkNoITSROFrameBorder", true, "ITSROFrame border event selection cut"};
131137
Configurable<bool> cfgEvSelkNoTimeFrameBorder{"cfgEvSelkNoTimeFrameBorder", true, "TimeFrame border event selection cut"};
132138
Configurable<bool> cfgEvSelUseGoodZvtxFT0vsPV{"cfgEvSelUseGoodZvtxFT0vsPV", true, "GoodZvertex and FT0 vs PV cut"};
139+
Configurable<bool> cfgEvSelUseOcuppancyTimeCut{"cfgEvSelUseOcuppancyTimeCut", true, "Occupancy Time pattern cut"};
140+
Configurable<bool> cfgEvSelSetOcuppancyRange{"cfgEvSelSetOcuppancyRange", true, "Use cut on occupancy range"};
141+
Configurable<int> cfgMinOccupancy{"cfgMinOccupancy", 0, "min. value of occupancy"};
142+
Configurable<int> cfgMaxOccupancy{"cfgMaxOccupancy", 3000, "max. value of occupancy"};
143+
133144
Configurable<bool> cfgUseItsPID{"cfgUseItsPID", false, "Use ITS PID for particle identification"};
134145
Configurable<float> cfgPtCutTOF{"cfgPtCutTOF", 0.3f, "Minimum pt to use TOF N-sigma"};
135146
Configurable<LabeledArray<float>> nSigmas{"nSigmas", {LongArrayFloat[0], 3, 6, {"TPC", "TOF", "ITS"}, {"pos_pi", "pos_ka", "pos_pr", "neg_pi", "neg_ka", "neg_pr"}}, "Labeled array for n-sigma values for TPC, TOF, ITS for pions, kaons, protons (positive and negative)"};
@@ -138,6 +149,8 @@ struct V0ptHadPiKaProt {
138149
Configurable<float> cfgCutPtMaxForV02{"cfgCutPtMaxForV02", 3.0f, "Max. pT for v02(pT)"};
139150
Configurable<float> cfgCutEtaWindowB{"cfgCutEtaWindowB", 0.4f, "value of x in |eta|<x for window B"};
140151
Configurable<bool> cfgLoadPhiWeights{"cfgLoadPhiWeights", false, "Load phi weights from CCDB to take care of non-uniform acceptance"};
152+
Configurable<bool> cfgLoadPtEffWeights{"cfgLoadPtEffWeights", false, "Load pt-dependent efficiency weights from CCDB to take care of detector inefficiency"};
153+
Configurable<int> cfgMinNoOfParticles{"cfgMinNoOfParticles", 4, "Minimum no. of particles for calculating v02(pT)"};
141154

142155
// pT dep DCAxy and DCAz cuts
143156
Configurable<bool> cfgUsePtDepDCAxy{"cfgUsePtDepDCAxy", true, "Use pt-dependent DCAxy cut"};
@@ -193,7 +206,12 @@ struct V0ptHadPiKaProt {
193206
TRandom3* funRndm = new TRandom3(0);
194207

195208
// Phi weight histograms initialization
196-
TH2F* hWeightPhiFunctionVzEtaPhi = nullptr;
209+
TH3D* hWeightPhiFunctionVzEtaPhi = nullptr;
210+
// Efficiency of diff. particle histograms initialization
211+
TH1D* hEffAllCharged = nullptr;
212+
TH1D* hEffPion = nullptr;
213+
TH1D* hEffKaon = nullptr;
214+
TH1D* hEffProton = nullptr;
197215

198216
// Filter command***********
199217
Filter collisionFilter = nabs(aod::collision::posZ) < cfgCutVertex;
@@ -235,7 +253,7 @@ struct V0ptHadPiKaProt {
235253
// Loading phi weight histograms from CCDB
236254
if (cfgLoadPhiWeights) {
237255

238-
// Accessing eff histograms
256+
// Accessing phi weight histograms
239257
ccdb->setURL(ccdbUrl.value);
240258
// Enabling object caching, otherwise each call goes to the CCDB server
241259
ccdb->setCaching(true);
@@ -245,11 +263,28 @@ struct V0ptHadPiKaProt {
245263
ccdb->setCreatedNotAfter(ccdbNoLaterThan.value);
246264
LOGF(info, "Getting object %s", ccdbPath.value.data());
247265
TList* lst = ccdb->getForTimeStamp<TList>(ccdbPath.value, ccdbNoLaterThan.value);
248-
hWeightPhiFunctionVzEtaPhi = reinterpret_cast<TH2F*>(lst->FindObject("hWeightPhiFunctionVzEtaPhi"));
266+
hWeightPhiFunctionVzEtaPhi = reinterpret_cast<TH3D*>(lst->FindObject("hWeightPhiFunctionVzEtaPhi"));
249267
if (!hWeightPhiFunctionVzEtaPhi)
250268
LOGF(info, "FATAL!! could not get phi weights---------> check");
251269
}
252270

271+
// Loading pT-dependent efficiency histograms from CCDB
272+
if (cfgLoadPtEffWeights) {
273+
274+
ccdb->setURL(ccdbUrl.value);
275+
ccdb->setCaching(true);
276+
ccdb->setLocalObjectValidityChecking();
277+
ccdb->setCreatedNotAfter(ccdbNoLaterThanPtEff.value);
278+
LOGF(info, "Getting object %s", ccdbPathPtEff.value.data());
279+
TList* lst = ccdb->getForTimeStamp<TList>(ccdbPathPtEff.value, ccdbNoLaterThanPtEff.value);
280+
hEffAllCharged = reinterpret_cast<TH1D*>(lst->FindObject("hEffAllCharged"));
281+
hEffPion = reinterpret_cast<TH1D*>(lst->FindObject("hEffPion"));
282+
hEffKaon = reinterpret_cast<TH1D*>(lst->FindObject("hEffKaon"));
283+
hEffProton = reinterpret_cast<TH1D*>(lst->FindObject("hEffProton"));
284+
if (!hEffAllCharged || !hEffPion || !hEffKaon || !hEffProton)
285+
LOGF(info, "FATAL!! could not get efficiency files---------> !!! check !!!");
286+
}
287+
253288
// Define axes
254289
std::vector<double> ptBin = {0.2, 0.4, 0.6, 0.8, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.5, 4.0, 5.0, 6.0, 8.0, 10.0};
255290
AxisSpec ptAxis = {ptBin, "#it{p}_{T} (GeV/#it{c})"};
@@ -458,9 +493,17 @@ struct V0ptHadPiKaProt {
458493
flag2 += 1;
459494
if (combNSigmaKa < cfgnSigmaOtherParticles)
460495
flag2 += 1;
461-
if (!(flag2 > 1) && !(combNSigmaPr > combNSigmaPi) && !(combNSigmaPr > combNSigmaKa)) {
462-
if (combNSigmaPr < cfgnSigmaCutCombTPCTOF) {
463-
flag = 1;
496+
497+
if (cfgUseNewSeperationPid) {
498+
if (std::abs(candidate.tpcNSigmaPr()) < cfgnSigmaCutTPCHigherPt && std::abs(candidate.tofNSigmaPr()) < cfgnSigmaCutTOFHigherPt) {
499+
if (!(flag2 > 1) && std::abs(candidate.tpcNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tpcNSigmaKa()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaKa()) > cfgnSigmaSeperationCut)
500+
flag = 1;
501+
}
502+
} else {
503+
if (!(flag2 > 1) && !(combNSigmaPr > combNSigmaPi) && !(combNSigmaPr > combNSigmaKa)) {
504+
if (combNSigmaPr < cfgnSigmaCutCombTPCTOF) {
505+
flag = 1;
506+
}
464507
}
465508
}
466509
}
@@ -497,9 +540,17 @@ struct V0ptHadPiKaProt {
497540
flag2 += 1;
498541
if (combNSigmaKa < cfgnSigmaOtherParticles)
499542
flag2 += 1;
500-
if (!(flag2 > 1) && !(combNSigmaPi > combNSigmaPr) && !(combNSigmaPi > combNSigmaKa)) {
501-
if (combNSigmaPi < cfgnSigmaCutCombTPCTOF) {
502-
flag = 1;
543+
544+
if (cfgUseNewSeperationPid) {
545+
if (std::abs(candidate.tpcNSigmaPi()) < cfgnSigmaCutTPCHigherPt && std::abs(candidate.tofNSigmaPi()) < cfgnSigmaCutTOFHigherPt) {
546+
if (!(flag2 > 1) && std::abs(candidate.tpcNSigmaKa()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaKa()) > cfgnSigmaSeperationCut && std::abs(candidate.tpcNSigmaPr()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPr()) > cfgnSigmaSeperationCut)
547+
flag = 1;
548+
}
549+
} else {
550+
if (!(flag2 > 1) && !(combNSigmaPi > combNSigmaPr) && !(combNSigmaPi > combNSigmaKa)) {
551+
if (combNSigmaPi < cfgnSigmaCutCombTPCTOF) {
552+
flag = 1;
553+
}
503554
}
504555
}
505556
}
@@ -536,9 +587,17 @@ struct V0ptHadPiKaProt {
536587
flag2 += 1;
537588
if (combNSigmaKa < cfgnSigmaOtherParticles)
538589
flag2 += 1;
539-
if (!(flag2 > 1) && !(combNSigmaKa > combNSigmaPi) && !(combNSigmaKa > combNSigmaPr)) {
540-
if (combNSigmaKa < cfgnSigmaCutCombTPCTOF) {
541-
flag = 1;
590+
591+
if (cfgUseNewSeperationPid) {
592+
if (std::abs(candidate.tpcNSigmaKa()) < cfgnSigmaCutTPCHigherPt && std::abs(candidate.tofNSigmaKa()) < cfgnSigmaCutTOFHigherPt) {
593+
if (!(flag2 > 1) && std::abs(candidate.tpcNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tpcNSigmaPr()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPr()) > cfgnSigmaSeperationCut)
594+
flag = 1;
595+
}
596+
} else {
597+
if (!(flag2 > 1) && !(combNSigmaKa > combNSigmaPi) && !(combNSigmaKa > combNSigmaPr)) {
598+
if (combNSigmaKa < cfgnSigmaCutCombTPCTOF) {
599+
flag = 1;
600+
}
542601
}
543602
}
544603
}
@@ -675,6 +734,18 @@ struct V0ptHadPiKaProt {
675734
}
676735

677736
histos.fill(HIST("hEventStatData"), 6.5);
737+
// events with selection bits based on occupancy time pattern
738+
if (cfgEvSelUseOcuppancyTimeCut && !(coll.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard))) {
739+
return 0;
740+
}
741+
742+
histos.fill(HIST("hEventStatData"), 7.5);
743+
int occupancy = coll.trackOccupancyInTimeRange();
744+
if (cfgEvSelSetOcuppancyRange && (occupancy < cfgMinOccupancy || occupancy > cfgMaxOccupancy)) {
745+
return 0;
746+
}
747+
748+
histos.fill(HIST("hEventStatData"), 8.5);
678749
return 1;
679750
}
680751

@@ -714,6 +785,66 @@ struct V0ptHadPiKaProt {
714785
return weight;
715786
}
716787

788+
template <typename T>
789+
float getEffAllCharged(const T& candidate)
790+
{
791+
if (!cfgLoadPtEffWeights || !hEffAllCharged) {
792+
return 1.0;
793+
}
794+
int bin = hEffAllCharged->FindBin(candidate.pt());
795+
float eff = hEffAllCharged->GetBinContent(bin);
796+
float ptweight = 1.0 / eff;
797+
if (!std::isfinite(ptweight) || ptweight <= 0) {
798+
return 1.0;
799+
}
800+
return ptweight;
801+
}
802+
803+
template <typename T>
804+
float getEffPion(const T& candidate)
805+
{
806+
if (!cfgLoadPtEffWeights || !hEffPion) {
807+
return 1.0;
808+
}
809+
int bin = hEffPion->FindBin(candidate.pt());
810+
float eff = hEffPion->GetBinContent(bin);
811+
float ptweight = 1.0 / eff;
812+
if (!std::isfinite(ptweight) || ptweight <= 0) {
813+
return 1.0;
814+
}
815+
return ptweight;
816+
}
817+
818+
template <typename T>
819+
float getEffKaon(const T& candidate)
820+
{
821+
if (!cfgLoadPtEffWeights || !hEffKaon) {
822+
return 1.0;
823+
}
824+
int bin = hEffKaon->FindBin(candidate.pt());
825+
float eff = hEffKaon->GetBinContent(bin);
826+
float ptweight = 1.0 / eff;
827+
if (!std::isfinite(ptweight) || ptweight <= 0) {
828+
return 1.0;
829+
}
830+
return ptweight;
831+
}
832+
833+
template <typename T>
834+
float getEffProton(const T& candidate)
835+
{
836+
if (!cfgLoadPtEffWeights || !hEffProton) {
837+
return 1.0;
838+
}
839+
int bin = hEffProton->FindBin(candidate.pt());
840+
float eff = hEffProton->GetBinContent(bin);
841+
float ptweight = 1.0 / eff;
842+
if (!std::isfinite(ptweight) || ptweight <= 0) {
843+
return 1.0;
844+
}
845+
return ptweight;
846+
}
847+
717848
// process Data
718849
void process(AodCollisions::iterator const& coll, aod::BCsWithTimestamps const&, AodTracks const& inputTracks)
719850
{
@@ -820,30 +951,37 @@ struct V0ptHadPiKaProt {
820951
}
821952
}
822953

823-
// fill subevent B for f(pT) in v02(pT)
824-
if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) {
825-
if (std::abs(trkEta) < cfgCutEtaWindowB) {
826-
fPtProfileHadInWinB->Fill(trkPt);
827-
nSumInWinB += 1.0;
828-
}
829-
}
830954
double phiweight = 1.0;
831955
if (cfgLoadPhiWeights) {
832-
phiweight = getPhiWeight(track, coll.posZ());
956+
phiweight = getPhiWeight(track, coll.posZ()); // NUA weight
957+
}
958+
double effweight = 1.0;
959+
if (cfgLoadPtEffWeights) {
960+
effweight = 1.0 / getEffAllCharged(track); // NUE weight
833961
}
962+
double weight = phiweight * effweight;
963+
834964
// fill subevent C for v2^2 in v02(pT)
835965
if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) {
836966
histos.fill(HIST("h3DVtxZetaPhi"), coll.posZ(), trkEta, trkPhi);
837-
if (cfgCutEtaWindowB < trkEta && trkEta < 0.8) {
838-
vecQInWinC += phiweight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi));
839-
nSumInWinC += phiweight;
967+
if (cfgCutEtaWindowB < trkEta && trkEta < cfgCutEta) {
968+
vecQInWinC += weight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi));
969+
nSumInWinC += weight;
840970
}
841971
}
842972
// fill subevent A for v2^2 in v02(pT)
843973
if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) {
844-
if (-0.8 < trkEta && trkEta < -1.0 * cfgCutEtaWindowB) {
845-
vecQInWinA += phiweight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi));
846-
nSumInWinA += phiweight;
974+
if (-1.0 * cfgCutEta < trkEta && trkEta < -1.0 * cfgCutEtaWindowB) {
975+
vecQInWinA += weight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi));
976+
nSumInWinA += weight;
977+
}
978+
}
979+
980+
// fill subevent B for f(pT) in v02(pT)
981+
if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) {
982+
if (std::abs(trkEta) < cfgCutEtaWindowB) {
983+
fPtProfileHadInWinB->Fill(trkPt, effweight);
984+
nSumInWinB += 1.0;
847985
}
848986
}
849987

@@ -919,17 +1057,25 @@ struct V0ptHadPiKaProt {
9191057
}
9201058
}
9211059

1060+
double effweightPion = 1.0;
1061+
double effweightKaon = 1.0;
1062+
double effweightProton = 1.0;
1063+
if (cfgLoadPtEffWeights) {
1064+
effweightPion = 1.0 / getEffPion(track); // NUE weight for pion
1065+
effweightKaon = 1.0 / getEffKaon(track); // NUE weight for kaon
1066+
effweightProton = 1.0 / getEffProton(track); // NUE weight for proton
1067+
}
9221068
// fill subevent B for ***identified particles'*** f(pT) in v02(pT)
9231069
if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) {
9241070
if (std::abs(trkEta) < cfgCutEtaWindowB) {
9251071
if (isPion) {
926-
fPtProfilePiInWinB->Fill(trkPt);
1072+
fPtProfilePiInWinB->Fill(trkPt, effweightPion);
9271073
}
9281074
if (isKaon) {
929-
fPtProfileKaInWinB->Fill(trkPt);
1075+
fPtProfileKaInWinB->Fill(trkPt, effweightKaon);
9301076
}
9311077
if (isProton && trkPt > cfgCutPtLowerProt) {
932-
fPtProfileProtInWinB->Fill(trkPt);
1078+
fPtProfileProtInWinB->Fill(trkPt, effweightProton);
9331079
}
9341080
}
9351081
}
@@ -1003,7 +1149,7 @@ struct V0ptHadPiKaProt {
10031149
}
10041150
}
10051151

1006-
if (nSumInWinA > 4 && nSumInWinB > 4 && nSumInWinC > 4) {
1152+
if (nSumInWinA > cfgMinNoOfParticles && nSumInWinB > cfgMinNoOfParticles && nSumInWinC > cfgMinNoOfParticles) {
10071153
double twoParCorr = (vecQInWinA * TComplex::Conjugate(vecQInWinC)).Re();
10081154
twoParCorr *= 1.0 / (nSumInWinA * nSumInWinC);
10091155
histos.get<TProfile2D>(HIST("Prof_XY"))->Fill(cent, 0.5, twoParCorr);

0 commit comments

Comments
 (0)