Skip to content

Commit f2fb3b5

Browse files
joachimckhjoachimckh
authored
[Common] Delta-pt light ion event cuts (#13980)
1 parent 5c53e2e commit f2fb3b5

File tree

1 file changed

+214
-5
lines changed

1 file changed

+214
-5
lines changed

‎Common/TableProducer/eseTableProducer.cxx‎

Lines changed: 214 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,7 @@
3636
#include <Framework/OutputObjHeader.h>
3737
#include <Framework/runDataProcessing.h>
3838

39+
#include <TF1.h>
3940
#include <TMath.h>
4041

4142
#include <chrono>
@@ -110,15 +111,68 @@ struct EseTableProducer {
110111
bool correctionsLoaded = false;
111112
} cfg;
112113

114+
TF1* fMultPVCutLow = nullptr;
115+
TF1* fMultPVCutHigh = nullptr;
116+
TF1* fMultCutLow = nullptr;
117+
TF1* fMultCutHigh = nullptr;
118+
TF1* fMultPVGlobalCutHigh = nullptr;
119+
TF1* fMultGlobalV0ACutLow = nullptr;
120+
TF1* fMultGlobalV0ACutHigh = nullptr;
121+
TF1* fMultGlobalT0ACutLow = nullptr;
122+
TF1* fMultGlobalT0ACutHigh = nullptr;
123+
124+
TF1* fPtDepDCAxy = nullptr;
125+
126+
float vtxZlow = -10.0, vtxZup = 10.0;
127+
int nchbins = 300;
128+
float nchlow = 0;
129+
float nchup = 3000;
130+
std::vector<double> multGlobalCorrCutPars;
131+
std::vector<double> multPVCorrCutPars;
132+
std::vector<double> multGlobalPVCorrCutPars;
133+
std::vector<double> multGlobalV0ACutPars;
134+
std::vector<double> multGlobalT0ACutPars;
135+
113136
Configurable<float> cfgVtxZ{"cfgVtxZ", 10.0f, "max z vertex position"};
114137
Configurable<float> cfgEta{"cfgEta", 0.8f, "max eta"};
115138
Configurable<float> cfgPtmin{"cfgPtmin", 0.2f, "min pt"};
116-
Configurable<float> cfgPtmax{"cfgPtmax", 5.0f, "max pt"};
139+
Configurable<float> cfgPtmax{"cfgPtmax", 3.0f, "max pt"};
117140
Configurable<float> cfgChi2PrITSCls{"cfgChi2PrITSCls", 4.0f, "max chi2 per ITS cluster"};
118141
Configurable<float> cfgChi2PrTPCCls{"cfgChi2PrTPCCls", 2.5f, "max chi2 per TPC cluster"};
119142
Configurable<float> cfgDCAz{"cfgDCAz", 2.0f, "max DCAz cut"};
143+
Configurable<bool> cfgDoOccupancySel{"cfgDoOccupancySel", true, "enable occupancy selection"};
144+
Configurable<int> cfgOccupancySelection{"cfgOccupancySelection", 2000, "Max occupancy selection, -999 to disable"};
145+
Configurable<bool> cfgUseAdditionalEventCut{"cfgUseAdditionalEventCut", false, "Use additional event cut on mult correlations"};
146+
Configurable<bool> cfgMultCut{"cfgMultCut", true, "Use additional event cut on mult correlations"};
147+
Configurable<bool> cfgTVXinTRD{"cfgTVXinTRD", true, "Use kTVXinTRD (reject TRD triggered events)"};
148+
Configurable<bool> cfgNoSameBunchPileupCut{"cfgNoSameBunchPileupCut", true, "kNoSameBunchPileupCut"};
149+
Configurable<bool> cfgIsGoodZvtxFT0vsPV{"cfgIsGoodZvtxFT0vsPV", true, "kIsGoodZvtxFT0vsPV"};
150+
Configurable<bool> cfgNoCollInTimeRangeStandard{"cfgNoCollInTimeRangeStandard", true, "kNoCollInTimeRangeStandard"};
151+
Configurable<bool> cfgIsVertexITSTPC{"cfgIsVertexITSTPC", true, "Selects collisions with at least one ITS-TPC track"};
152+
Configurable<bool> cfgIsGoodITSLayersAll{"cfgIsGoodITSLayersAll", true, "kIsGoodITSLayersAll"};
120153
Configurable<std::string> cfgEfficiency{"cfgEfficiency", "", "CCDB path to efficiency object"};
121154

155+
// Cut Configurables
156+
Configurable<std::vector<double>> cfgMultGlobalCutPars{"cfgMultGlobalCutPars", std::vector<double>{2272.16, -76.6932, 1.01204, -0.00631545, 1.59868e-05, 136.336, -4.97006, 0.121199, -0.0015921, 7.66197e-06}, "Global vs FT0C multiplicity cut parameter values"};
157+
Configurable<std::vector<double>> cfgMultPVCutPars{"cfgMultPVCutPars", std::vector<double>{3074.43, -106.192, 1.46176, -0.00968364, 2.61923e-05, 182.128, -7.43492, 0.193901, -0.00256715, 1.22594e-05}, "PV vs FT0C multiplicity cut parameter values"};
158+
Configurable<std::vector<double>> cfgMultGlobalPVCutPars{"cfgMultGlobalPVCutPars", std::vector<double>{-0.223013, 0.715849, 0.664242, 0.0829653, -0.000503733, 1.21185e-06}, "Global vs PV multiplicity cut parameter values"};
159+
Configurable<std::string> cfgMultCorrHighCutFunction{"cfgMultCorrHighCutFunction", "[0] + [1]*x + [2]*x*x + [3]*x*x*x + [4]*x*x*x*x + 3.*([5] + [6]*x + [7]*x*x + [8]*x*x*x + [9]*x*x*x*x)", "Functional for multiplicity correlation cut"};
160+
Configurable<std::string> cfgMultCorrLowCutFunction{"cfgMultCorrLowCutFunction", "[0] + [1]*x + [2]*x*x + [3]*x*x*x + [4]*x*x*x*x - 3.*([5] + [6]*x + [7]*x*x + [8]*x*x*x + [9]*x*x*x*x)", "Functional for multiplicity correlation cut"};
161+
Configurable<std::string> cfgMultGlobalPVCorrCutFunction{"cfgMultGlobalPVCorrCutFunction", "[0] + [1]*x + 3*([2] + [3]*x + [4]*x*x + [5]*x*x*x)", "Functional for global vs pv multiplicity correlation cut"};
162+
struct : ConfigurableGroup {
163+
Configurable<std::string> cfgMultGlobalASideCorrCutFunction{"cfgMultGlobalASideCorrCutFunction", "[0]+[1]*x+[2]*x*x+[3]*x*x*x+[4]*x*x*x*x + [10]*([5]+[6]*x+[7]*x*x+[8]*x*x*x+[9]*x*x*x*x)", "Functional for global vs V0A multiplicity low correlation cut"};
164+
Configurable<std::vector<double>> cfgMultGlobalV0ACutPars{"cfgMultGlobalV0ACutPars", std::vector<double>{567.785, 172.715, 0.77888, -0.00693466, 1.40564e-05, 679.853, 66.8068, -0.444332, 0.00115002, -4.92064e-07}, "Global vs FV0A multiplicity cut parameter values"};
165+
Configurable<std::vector<double>> cfgMultGlobalT0ACutPars{"cfgMultGlobalT0ACutPars", std::vector<double>{241.618, 61.8402, 0.348049, -0.00306078, 6.20357e-06, 315.235, 29.1491, -0.188639, 0.00044528, -9.08912e-08}, "Global vs FT0A multiplicity cut parameter values"};
166+
Configurable<float> cfgGlobalV0ALowSigma{"cfgGlobalV0ALowSigma", -3.0f, "Number of sigma deviations below expected value in global vs V0A correlation"};
167+
Configurable<float> cfgGlobalV0AHighSigma{"cfgGlobalV0AHighSigma", 4.0f, "Number of sigma deviations above expected value in global vs V0A correlation"};
168+
Configurable<float> cfgGlobalT0ALowSigma{"cfgGlobalT0ALowSigma", -3.0f, "Number of sigma deviations below expected value in global vs T0A correlation"};
169+
Configurable<float> cfgGlobalT0AHighSigma{"cfgGlobalT0AHighSigma", 4.0f, "Number of sigma deviations above expected value in global vs T0A correlation"};
170+
} cfgGlobalAsideCorrCuts;
171+
Configurable<std::string> cfgDCAxy{"cfgDCAxy", "(0.0026+0.005/(x^1.01))", "Functional form of pt-dependent DCAxy cut"};
172+
Configurable<float> cfgDCAxyNSigma{"cfgDCAxyNSigma", 7.0f, "Cut on number of sigma deviations from expected DCA in the transverse direction"};
173+
174+
//
175+
122176
// o2::framework::expressions::Filter collisionFilter = nabs(aod::collision::posZ) < cfgVtxZ;
123177
o2::framework::expressions::Filter trackFilter = nabs(aod::track::eta) < cfgEta && aod::track::pt > cfgPtmin&& aod::track::pt < cfgPtmax && ((requireGlobalTrackInFilter()) || (aod::track::isGlobalTrackSDD == static_cast<uint8_t>(true))) && (aod::track::itsChi2NCl < cfgChi2PrITSCls) && (aod::track::tpcChi2NCl < cfgChi2PrTPCCls) && nabs(aod::track::dcaZ) < cfgDCAz;
124178

@@ -127,8 +181,23 @@ struct EseTableProducer {
127181

128182
void init(o2::framework::InitContext&)
129183
{
130-
131184
LOGF(info, "ESETable::init()");
185+
EseTableProducer::vtxZlow = -cfgVtxZ;
186+
EseTableProducer::vtxZup = cfgVtxZ;
187+
EseTableProducer::multPVCorrCutPars = cfgMultPVCutPars;
188+
EseTableProducer::multGlobalCorrCutPars = cfgMultGlobalCutPars;
189+
EseTableProducer::multGlobalV0ACutPars = cfgGlobalAsideCorrCuts.cfgMultGlobalV0ACutPars;
190+
EseTableProducer::multGlobalT0ACutPars = cfgGlobalAsideCorrCuts.cfgMultGlobalT0ACutPars;
191+
EseTableProducer::multGlobalPVCorrCutPars = cfgMultGlobalPVCutPars;
192+
193+
AxisSpec t0cAxis = {1000, 0, 10000, "N_{ch} (T0C)"};
194+
AxisSpec t0aAxis = {300, 0, 30000, "N_{ch} (T0A)"};
195+
AxisSpec v0aAxis = {800, 0, 80000, "N_{ch} (V0A)"};
196+
std::vector<double> nchbinning;
197+
int nchskip = (nchup - nchlow) / nchbins;
198+
for (int i = 0; i <= nchbins; ++i) {
199+
nchbinning.push_back(nchskip * i + nchlow + 0.5);
200+
}
132201

133202
registry.add("hEventCounter", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
134203
registry.add("hESEstat", "ese status;ese status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
@@ -151,6 +220,52 @@ struct EseTableProducer {
151220
weightsFFit->setResolution(cfgnResolution);
152221
weightsFFit->setQnType(veccfg);
153222
weightsFFit->init();
223+
224+
fPtDepDCAxy = new TF1("ptDepDCAxy", Form("[0]*%s", cfgDCAxy->c_str()), 0.001, 100);
225+
fPtDepDCAxy->SetParameter(0, cfgDCAxyNSigma);
226+
LOGF(info, "DCAxy pt-dependence function: %s", Form("[0]*%s", cfgDCAxy->c_str()));
227+
if (cfgUseAdditionalEventCut) {
228+
fMultPVCutLow = new TF1("fMultPVCutLow", cfgMultCorrLowCutFunction->c_str(), 0, 100);
229+
fMultPVCutLow->SetParameters(&(EseTableProducer::multPVCorrCutPars[0]));
230+
fMultPVCutHigh = new TF1("fMultPVCutHigh", cfgMultCorrHighCutFunction->c_str(), 0, 100);
231+
fMultPVCutHigh->SetParameters(&(EseTableProducer::multPVCorrCutPars[0]));
232+
fMultCutLow = new TF1("fMultCutLow", cfgMultCorrLowCutFunction->c_str(), 0, 100);
233+
fMultCutLow->SetParameters(&(EseTableProducer::multGlobalCorrCutPars[0]));
234+
fMultCutHigh = new TF1("fMultCutHigh", cfgMultCorrHighCutFunction->c_str(), 0, 100);
235+
fMultCutHigh->SetParameters(&(EseTableProducer::multGlobalCorrCutPars[0]));
236+
fMultPVGlobalCutHigh = new TF1("fMultPVGlobalCutHigh", cfgMultGlobalPVCorrCutFunction->c_str(), 0, nchbinning.back());
237+
fMultPVGlobalCutHigh->SetParameters(&(EseTableProducer::multGlobalPVCorrCutPars[0]));
238+
239+
LOGF(info, "Global V0A function: %s in range 0-%g", cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), v0aAxis.binEdges.back());
240+
fMultGlobalV0ACutLow = new TF1("fMultGlobalV0ACutLow", cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, v0aAxis.binEdges.back());
241+
for (std::size_t i = 0; i < EseTableProducer::multGlobalV0ACutPars.size(); ++i)
242+
fMultGlobalV0ACutLow->SetParameter(i, EseTableProducer::multGlobalV0ACutPars[i]);
243+
fMultGlobalV0ACutLow->SetParameter(EseTableProducer::multGlobalV0ACutPars.size(), cfgGlobalAsideCorrCuts.cfgGlobalV0ALowSigma);
244+
for (int i = 0; i < fMultGlobalV0ACutLow->GetNpar(); ++i)
245+
LOGF(info, "fMultGlobalV0ACutLow par %d = %g", i, fMultGlobalV0ACutLow->GetParameter(i));
246+
247+
fMultGlobalV0ACutHigh = new TF1("fMultGlobalV0ACutHigh", cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, v0aAxis.binEdges.back());
248+
for (std::size_t i = 0; i < EseTableProducer::multGlobalV0ACutPars.size(); ++i)
249+
fMultGlobalV0ACutHigh->SetParameter(i, EseTableProducer::multGlobalV0ACutPars[i]);
250+
fMultGlobalV0ACutHigh->SetParameter(EseTableProducer::multGlobalV0ACutPars.size(), cfgGlobalAsideCorrCuts.cfgGlobalV0AHighSigma);
251+
for (int i = 0; i < fMultGlobalV0ACutHigh->GetNpar(); ++i)
252+
LOGF(info, "fMultGlobalV0ACutHigh par %d = %g", i, fMultGlobalV0ACutHigh->GetParameter(i));
253+
254+
LOGF(info, "Global T0A function: %s", cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str());
255+
fMultGlobalT0ACutLow = new TF1("fMultGlobalT0ACutLow", cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, t0aAxis.binEdges.back());
256+
for (std::size_t i = 0; i < EseTableProducer::multGlobalT0ACutPars.size(); ++i)
257+
fMultGlobalT0ACutLow->SetParameter(i, EseTableProducer::multGlobalT0ACutPars[i]);
258+
fMultGlobalT0ACutLow->SetParameter(EseTableProducer::multGlobalT0ACutPars.size(), cfgGlobalAsideCorrCuts.cfgGlobalT0ALowSigma);
259+
for (int i = 0; i < fMultGlobalT0ACutLow->GetNpar(); ++i)
260+
LOGF(info, "fMultGlobalT0ACutLow par %d = %g", i, fMultGlobalT0ACutLow->GetParameter(i));
261+
262+
fMultGlobalT0ACutHigh = new TF1("fMultGlobalT0ACutHigh", cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, t0aAxis.binEdges.back());
263+
for (std::size_t i = 0; i < EseTableProducer::multGlobalT0ACutPars.size(); ++i)
264+
fMultGlobalT0ACutHigh->SetParameter(i, EseTableProducer::multGlobalT0ACutPars[i]);
265+
fMultGlobalT0ACutHigh->SetParameter(EseTableProducer::multGlobalT0ACutPars.size(), cfgGlobalAsideCorrCuts.cfgGlobalT0AHighSigma);
266+
for (int i = 0; i < fMultGlobalT0ACutHigh->GetNpar(); ++i)
267+
LOGF(info, "fMultGlobalT0ACutHigh par %d = %g", i, fMultGlobalT0ACutHigh->GetParameter(i));
268+
}
154269
}
155270

156271
void initCCDB(aod::BCsWithTimestamps::iterator const& bc)
@@ -329,16 +444,111 @@ struct EseTableProducer {
329444
}
330445
PROCESS_SWITCH(EseTableProducer, processESE, "process q vectors to calculate reduced q-vector", true);
331446

332-
void processMeanPt(soa::Join<aod::Collisions, aod::EvSels, aod::Mults, aod::CentFT0Cs, aod::CentFT0CVariant1s, aod::CentFT0Ms, aod::CentNTPVs, aod::CentNGlobals, aod::CentMFTs>::iterator const& collision, aod::BCsWithTimestamps const&, GFWTracks const& tracks)
447+
struct XAxis {
448+
float centrality;
449+
int64_t multiplicity;
450+
};
451+
452+
template <typename TCollision>
453+
bool eventSelected(TCollision collision, const int& multTrk, const float& centrality)
333454
{
455+
if (cfgTVXinTRD) {
456+
if (collision.alias_bit(kTVXinTRD)) {
457+
return 0;
458+
}
459+
}
460+
if (cfgNoSameBunchPileupCut) {
461+
if (!collision.selection_bit(o2::aod::evsel::kNoSameBunchPileup)) {
462+
return 0;
463+
}
464+
}
465+
if (cfgIsGoodZvtxFT0vsPV) {
466+
if (!collision.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV)) {
467+
return 0;
468+
}
469+
}
470+
if (cfgNoCollInTimeRangeStandard) {
471+
if (!collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)) {
472+
return 0;
473+
}
474+
}
475+
476+
if (cfgIsVertexITSTPC) {
477+
if (!collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC)) {
478+
return 0;
479+
}
480+
}
481+
482+
if (cfgIsGoodITSLayersAll) {
483+
if (!collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) {
484+
return 0;
485+
}
486+
}
487+
488+
float vtxz = -999;
489+
if (collision.numContrib() > 1) {
490+
vtxz = collision.posZ();
491+
float zRes = std::sqrt(collision.covZZ());
492+
float minZRes = 0.25;
493+
int minNContrib = 20;
494+
if (zRes > minZRes && collision.numContrib() < minNContrib)
495+
vtxz = -999;
496+
}
497+
auto multNTracksPV = collision.multNTracksPV();
498+
499+
if (vtxz > EseTableProducer::vtxZup || vtxz < EseTableProducer::vtxZlow)
500+
return 0;
501+
502+
if (cfgMultCut) {
503+
if (multNTracksPV < fMultPVCutLow->Eval(centrality))
504+
return 0;
505+
if (multNTracksPV > fMultPVCutHigh->Eval(centrality))
506+
return 0;
507+
if (multTrk < fMultCutLow->Eval(centrality))
508+
return 0;
509+
if (multTrk > fMultCutHigh->Eval(centrality))
510+
return 0;
511+
if (multTrk > fMultPVGlobalCutHigh->Eval(collision.multNTracksPV()))
512+
return 0;
513+
514+
if (!(cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast<double>(collision.multFV0A()) < fMultGlobalV0ACutLow->Eval(multTrk))
515+
return 0;
516+
if (!(cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast<double>(collision.multFV0A()) > fMultGlobalV0ACutHigh->Eval(multTrk))
517+
return 0;
518+
if (!(cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast<double>(collision.multFT0A()) < fMultGlobalT0ACutLow->Eval(multTrk))
519+
return 0;
520+
if (!(cfgGlobalAsideCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast<double>(collision.multFT0A()) > fMultGlobalT0ACutHigh->Eval(multTrk))
521+
return 0;
522+
}
523+
return 1;
524+
}
334525

335-
std::vector<float> meanPt{-1};
526+
void processMeanPt(soa::Join<aod::Collisions, aod::EvSels, aod::Mults, aod::CentFT0Cs, aod::CentFT0CVariant1s, aod::CentFT0Ms, aod::CentNTPVs, aod::CentNGlobals, aod::CentMFTs>::iterator const& collision, aod::BCsWithTimestamps const&, GFWTracks const& tracks)
527+
{
528+
std::vector<float> meanPt{-999};
336529
std::vector<float> meanPtShape{-1};
337530
if (collision.posZ() < -cfgVtxZ || collision.posZ() > cfgVtxZ) {
338531
meanPts(meanPt);
339532
meanPtShapes(meanPtShape);
340533
return;
341534
}
535+
if (!collision.sel8()) {
536+
meanPts(meanPt);
537+
meanPtShapes(meanPtShape);
538+
return;
539+
}
540+
if (cfgDoOccupancySel) {
541+
int occupancy = collision.trackOccupancyInTimeRange();
542+
if (occupancy < 0 || occupancy > cfgOccupancySelection) {
543+
meanPts(meanPt);
544+
meanPtShapes(meanPtShape);
545+
return;
546+
}
547+
}
548+
549+
const XAxis xaxis{collision.centFT0C(), tracks.size()};
550+
if (cfgUseAdditionalEventCut && !eventSelected(collision, xaxis.multiplicity, xaxis.centrality))
551+
return;
342552

343553
registry.fill(HIST("hMeanPtStat"), 0.5);
344554
const auto centrality = collision.centFT0C();
@@ -363,7 +573,6 @@ struct EseTableProducer {
363573
}
364574
}
365575
}
366-
367576
meanPts(meanPt);
368577
meanPtShapes(meanPtShape);
369578
}

0 commit comments

Comments
 (0)