Skip to content

Commit f2e941b

Browse files
alcalivaalcaliva
authored
[PWGLF] Remove jet area cut and use uniform area for UE estimate (#16285)
1 parent 234717a commit f2e941b

1 file changed

Lines changed: 18 additions & 59 deletions

File tree

‎PWGLF/Tasks/Nuspex/antinucleiInJets.cxx‎

Lines changed: 18 additions & 59 deletions
Original file line numberDiff line numberDiff line change
@@ -132,7 +132,7 @@ struct JetMatching {
132132

133133
struct AntinucleiInJets {
134134

135-
// Random engine
135+
// Random engine (Mersenne Twister)
136136
std::mt19937 rng;
137137
std::uniform_int_distribution<int> generateRandomNr{0, 1};
138138

@@ -171,8 +171,7 @@ struct AntinucleiInJets {
171171
Configurable<double> ptLeadingMin{"ptLeadingMin", 5.0, "pt Leading Min"};
172172
Configurable<double> rJet{"rJet", 0.4, "Jet resolution parameter R"};
173173
Configurable<double> zVtx{"zVtx", 10.0, "Maximum zVertex"};
174-
Configurable<bool> applyAreaCut{"applyAreaCut", true, "apply area cut"};
175-
Configurable<double> maxNormalizedJetArea{"maxNormalizedJetArea", 1.0, "area cut"};
174+
Configurable<bool> applyAreaCut{"applyAreaCut", false, "apply area cut A > f * pi R^2"};
176175
Configurable<double> deltaEtaEdge{"deltaEtaEdge", 0.05, "eta gap from the edge"};
177176
Configurable<int> nSyst{"nSyst", 50, "number of systematic variations"};
178177
Configurable<int> nSubsamples{"nSubsamples", 50, "number of subsamples"};
@@ -1296,26 +1295,22 @@ struct AntinucleiInJets {
12961295
continue;
12971296

12981297
// Apply area cut if required
1299-
double normalizedJetArea = jet.area() / (PI * rJet * rJet);
1300-
if (applyAreaCut && (!isppRefAnalysis) && normalizedJetArea > maxNormalizedJetArea)
1301-
continue;
1302-
if (isppRefAnalysis && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
1298+
if (applyAreaCut && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
13031299
continue;
13041300
isAtLeastOneJetSelected = true;
13051301

1302+
// Fill histograms with jet effective area / piR^2 for normalization
1303+
registryData.fill(HIST("jetEffectiveAreaOverPiR2"), jet.area() / (PI * rJet * rJet));
1304+
registryData.fill(HIST("jetArea"), jet.area());
1305+
13061306
// Perpendicular cones
1307-
double coneRadius = std::sqrt(jet.area() / PI);
13081307
TVector3 jetAxis(jet.px(), jet.py(), jet.pz());
13091308
TVector3 ueAxis1(0, 0, 0), ueAxis2(0, 0, 0);
13101309
getPerpendicularDirections(jetAxis, ueAxis1, ueAxis2);
13111310
if (ueAxis1.Mag() == 0 || ueAxis2.Mag() == 0) {
13121311
continue;
13131312
}
13141313

1315-
// Fill histogram with jet effective area / piR^2
1316-
registryData.fill(HIST("jetEffectiveAreaOverPiR2"), jet.area() / (PI * rJet * rJet));
1317-
registryData.fill(HIST("jetArea"), jet.area());
1318-
13191314
// Get jet constituents
13201315
std::vector<fastjet::PseudoJet> jetConstituents = jet.constituents();
13211316

@@ -1427,14 +1422,8 @@ struct AntinucleiInJets {
14271422
double deltaPhiUe2 = getDeltaPhi(track.phi(), ueAxis2.Phi());
14281423
double deltaRUe2 = std::sqrt(deltaEtaUe2 * deltaEtaUe2 + deltaPhiUe2 * deltaPhiUe2);
14291424

1430-
// Determine the maximum allowed distance from UE axes for particle selection
1431-
double maxConeRadius = coneRadius;
1432-
if (applyAreaCut) {
1433-
maxConeRadius = std::sqrt(maxNormalizedJetArea) * rJet;
1434-
}
1435-
14361425
// Reject tracks that lie outside the maxConeRadius from both UE axes
1437-
if (deltaRUe1 > maxConeRadius && deltaRUe2 > maxConeRadius)
1426+
if (deltaRUe1 > rJet && deltaRUe2 > rJet)
14381427
continue;
14391428

14401429
// Define variables
@@ -1590,8 +1579,7 @@ struct AntinucleiInJets {
15901579
continue;
15911580

15921581
// Apply area cut if required
1593-
double normalizedJetArea = jet.area() / (PI * rJet * rJet);
1594-
if (applyAreaCut && normalizedJetArea > maxNormalizedJetArea)
1582+
if (applyAreaCut && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
15951583
continue;
15961584
isAtLeastOneJetSelected = true;
15971585
}
@@ -1672,7 +1660,6 @@ struct AntinucleiInJets {
16721660
// Jet properties and perpendicular cone
16731661
std::vector<fastjet::PseudoJet> jetConstituents = jet.constituents();
16741662
TVector3 jetAxis(jet.px(), jet.py(), jet.pz());
1675-
double coneRadius = std::sqrt(jet.area() / PI);
16761663
TVector3 ueAxis1(0, 0, 0), ueAxis2(0, 0, 0);
16771664
getPerpendicularDirections(jetAxis, ueAxis1, ueAxis2);
16781665
if (ueAxis1.Mag() == 0 || ueAxis2.Mag() == 0) {
@@ -1704,7 +1691,7 @@ struct AntinucleiInJets {
17041691
double deltaEtaUe2 = track.eta() - ueAxis2.Eta();
17051692
double deltaPhiUe2 = getDeltaPhi(track.phi(), ueAxis2.Phi());
17061693
double deltaRUe2 = std::sqrt(deltaEtaUe2 * deltaEtaUe2 + deltaPhiUe2 * deltaPhiUe2);
1707-
if (deltaRUe1 > coneRadius && deltaRUe2 > coneRadius)
1694+
if (deltaRUe1 > rJet && deltaRUe2 > rJet)
17081695
continue;
17091696

17101697
ptPerp = ptPerp + track.pt();
@@ -2136,10 +2123,7 @@ struct AntinucleiInJets {
21362123
continue;
21372124

21382125
// Apply area cut if required
2139-
double normalizedJetArea = jet.area() / (PI * rJet * rJet);
2140-
if (applyAreaCut && (!isppRefAnalysis) && normalizedJetArea > maxNormalizedJetArea)
2141-
continue;
2142-
if (isppRefAnalysis && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
2126+
if (applyAreaCut && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
21432127
continue;
21442128
isAtLeastOneJetSelected = true;
21452129

@@ -2181,7 +2165,6 @@ struct AntinucleiInJets {
21812165

21822166
// Set up two perpendicular cone axes for underlying event estimation
21832167
TVector3 jetAxis(jet.px(), jet.py(), jet.pz());
2184-
double coneRadius = std::sqrt(jet.area() / PI);
21852168
TVector3 ueAxis1(0, 0, 0), ueAxis2(0, 0, 0);
21862169
getPerpendicularDirections(jetAxis, ueAxis1, ueAxis2);
21872170
if (ueAxis1.Mag() == 0 || ueAxis2.Mag() == 0) {
@@ -2199,14 +2182,8 @@ struct AntinucleiInJets {
21992182
double deltaPhiUe2 = getDeltaPhi(protonVec.Phi(), ueAxis2.Phi());
22002183
double deltaRUe2 = std::sqrt(deltaEtaUe2 * deltaEtaUe2 + deltaPhiUe2 * deltaPhiUe2);
22012184

2202-
// Determine the maximum allowed distance from UE axes for particle selection
2203-
double maxConeRadius = coneRadius;
2204-
if (applyAreaCut) {
2205-
maxConeRadius = std::sqrt(maxNormalizedJetArea) * rJet;
2206-
}
2207-
22082185
// Reject tracks that lie outside the maxConeRadius from both UE axes
2209-
if (deltaRUe1 > maxConeRadius && deltaRUe2 > maxConeRadius)
2186+
if (deltaRUe1 > rJet && deltaRUe2 > rJet)
22102187
continue;
22112188

22122189
// Fill normalization histogram
@@ -2393,18 +2370,14 @@ struct AntinucleiInJets {
23932370
continue;
23942371

23952372
// Apply area cut if required
2396-
double normalizedJetArea = jet.area() / (PI * rJet * rJet);
2397-
if (applyAreaCut && (!isppRefAnalysis) && normalizedJetArea > maxNormalizedJetArea)
2398-
continue;
2399-
if (isppRefAnalysis && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
2373+
if (applyAreaCut && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
24002374
continue;
24012375
isAtLeastOneJetSelected = true;
24022376

24032377
// Reconstructed jets
24042378
registryMC.fill(HIST("recJets"), 0.5);
24052379

24062380
// Set up two perpendicular cone axes for underlying event estimation
2407-
double coneRadius = std::sqrt(jet.area() / PI);
24082381
TVector3 jetAxis(jet.px(), jet.py(), jet.pz());
24092382
TVector3 ueAxis1(0, 0, 0), ueAxis2(0, 0, 0);
24102383
getPerpendicularDirections(jetAxis, ueAxis1, ueAxis2);
@@ -2550,14 +2523,8 @@ struct AntinucleiInJets {
25502523
double deltaPhiUe2 = getDeltaPhi(track.phi(), ueAxis2.Phi());
25512524
double deltaRUe2 = std::sqrt(deltaEtaUe2 * deltaEtaUe2 + deltaPhiUe2 * deltaPhiUe2);
25522525

2553-
// Determine the maximum allowed distance from UE axes for particle selection
2554-
double maxConeRadius = coneRadius;
2555-
if (applyAreaCut) {
2556-
maxConeRadius = std::sqrt(maxNormalizedJetArea) * rJet;
2557-
}
2558-
25592526
// Reject tracks that lie outside the maxConeRadius from both UE axes
2560-
if (deltaRUe1 > maxConeRadius && deltaRUe2 > maxConeRadius)
2527+
if (deltaRUe1 > rJet && deltaRUe2 > rJet)
25612528
continue;
25622529

25632530
// Particle identification using the ITS cluster size
@@ -3753,8 +3720,7 @@ struct AntinucleiInJets {
37533720
continue;
37543721

37553722
// Apply area cut if required
3756-
double normalizedJetArea = jet.area() / (PI * rJet * rJet);
3757-
if (applyAreaCut && normalizedJetArea > maxNormalizedJetArea)
3723+
if (applyAreaCut && (jet.area() < cfgAreaFrac * PI * rJet * rJet))
37583724
continue;
37593725
isAtLeastOneJetSelected = true;
37603726

@@ -3769,7 +3735,6 @@ struct AntinucleiInJets {
37693735

37703736
// Set up two perpendicular cone axes for underlying event estimation
37713737
TVector3 jetAxis(jet.px(), jet.py(), jet.pz());
3772-
double coneRadius = std::sqrt(jet.area() / PI);
37733738
TVector3 ueAxis1(0, 0, 0), ueAxis2(0, 0, 0);
37743739
getPerpendicularDirections(jetAxis, ueAxis1, ueAxis2);
37753740
if (ueAxis1.Mag() == 0 || ueAxis2.Mag() == 0) {
@@ -3791,14 +3756,8 @@ struct AntinucleiInJets {
37913756
double deltaPhiUe2 = getDeltaPhi(chParticle.phi(), ueAxis2.Phi());
37923757
double deltaRUe2 = std::sqrt(deltaEtaUe2 * deltaEtaUe2 + deltaPhiUe2 * deltaPhiUe2);
37933758

3794-
// Determine the maximum allowed distance from UE axes for particle selection
3795-
double maxConeRadius = coneRadius;
3796-
if (applyAreaCut) {
3797-
maxConeRadius = std::sqrt(maxNormalizedJetArea) * rJet;
3798-
}
3799-
38003759
// Reject tracks that lie outside the maxConeRadius from both UE axes
3801-
if (deltaRUe1 > maxConeRadius && deltaRUe2 > maxConeRadius)
3760+
if (deltaRUe1 > rJet && deltaRUe2 > rJet)
38023761
continue;
38033762

38043763
// Fill histograms for UE
@@ -4121,7 +4080,7 @@ struct AntinucleiInJets {
41214080
continue;
41224081

41234082
// Apply area cut if required
4124-
if (applyAreaCut && (jetRec.area() / (PI * rJet * rJet)) > maxNormalizedJetArea)
4083+
if (applyAreaCut && (jetRec.area() < cfgAreaFrac * PI * rJet * rJet))
41254084
continue;
41264085

41274086
// Clear jet-pair container
@@ -4134,7 +4093,7 @@ struct AntinucleiInJets {
41344093
continue;
41354094

41364095
// Apply area cut if required
4137-
if (applyAreaCut && (jetGen.area() / (PI * rJet * rJet)) > maxNormalizedJetArea)
4096+
if (applyAreaCut && (jetGen.area() < cfgAreaFrac * PI * rJet * rJet))
41384097
continue;
41394098

41404099
double deltaEta = jetGen.eta() - jetRec.eta();

0 commit comments

Comments
 (0)