Skip to content

Commit 2c0f6f2

Browse files
dajones2dajones2
authored
[PWGJE] Changing outlier rejection on response to the event level (#11481)
1 parent a093163 commit 2c0f6f2

File tree

1 file changed

+60
-58
lines changed

1 file changed

+60
-58
lines changed

‎PWGJE/Tasks/jetHadronRecoil.cxx‎

Lines changed: 60 additions & 58 deletions
Original file line numberDiff line numberDiff line change
@@ -108,7 +108,6 @@ struct JetHadronRecoil {
108108
0.57, 0.60};
109109

110110
AxisSpec dRAxis = {dRBinning, "#Delta R"};
111-
112111
AxisSpec ptAxisDet = {ptBinningDet, "#it{p}_{T,det} (GeV/c)"};
113112
AxisSpec ptAxisPart = {ptBinningPart, "#it{p}_{T,part} (GeV/c)"};
114113
AxisSpec phiAxisDet = {100, 0.0, o2::constants::math::TwoPI, "#phi_{det}"};
@@ -145,6 +144,8 @@ struct JetHadronRecoil {
145144
{"hPtPart", "Particle p_{T};p_{T};entries", {HistType::kTH1F, {{200, 0, 200}}}},
146145
{"hEtaPart", "Particle #eta;#eta;entries", {HistType::kTH1F, {{100, -1.0, 1.0}}}},
147146
{"hPhiPart", "Particle #phi;#phi;entries", {HistType::kTH1F, {{100, 0.0, o2::constants::math::TwoPI}}}},
147+
{"hPart3D", "3D tracks histogram;p_{T};#eta;#phi", {HistType::kTH3F, {{200, 0, 200}, {100, -1.0, 1.0}, {100, 0.0, o2::constants::math::TwoPI}}}},
148+
{"hPtPartPtHard", "Track p_{T} vs #hat{p};p_{T};#frac{p_{T}}{#hat{p}}", {HistType::kTH2F, {{200, 0, 200}, {20, 0, 5}}}},
148149
{"hDeltaR", "#DeltaR;#DeltaR;#frac{dN_{jets}}{d#DeltaR}", {HistType::kTH1F, {dRAxis}}},
149150
{"hDeltaRPart", "Particle #DeltaR;#DeltaR;#frac{1}{N_{jets}}#frac{dN_{jets}}{d#DeltaR}", {HistType::kTH1F, {dRAxis}}},
150151
{"hDeltaRpT", "jet p_{T} vs #DeltaR;p_{T,jet};#DeltaR", {HistType::kTH2F, {{500, -100, 400}, dRAxis}}},
@@ -330,7 +331,9 @@ struct JetHadronRecoil {
330331
{
331332
bool isSigCol;
332333
std::vector<double> phiTTAr;
334+
std::vector<double> ptTTAr;
333335
double phiTT = 0;
336+
double ptTT = 0;
334337
int trigNumber = 0;
335338
int nTT = 0;
336339
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
@@ -354,27 +357,34 @@ struct JetHadronRecoil {
354357
}
355358
if (isSigCol && particle.pt() < ptTTsigMax && particle.pt() > ptTTsigMin) {
356359
phiTTAr.push_back(particle.phi());
360+
ptTTAr.push_back(particle.pt());
357361
nTT++;
358362
}
359363
if (!isSigCol && particle.pt() < ptTTrefMax && particle.pt() > ptTTrefMin) {
360364
phiTTAr.push_back(particle.phi());
365+
ptTTAr.push_back(particle.pt());
361366
nTT++;
362367
}
363368
registry.fill(HIST("hPtPart"), particle.pt(), weight);
364369
registry.fill(HIST("hEtaPart"), particle.eta(), weight);
365370
registry.fill(HIST("hPhiPart"), particle.phi(), weight);
371+
registry.fill(HIST("hPart3D"), particle.pt(), particle.eta(), particle.phi(), weight);
372+
registry.fill(HIST("hPtPartPtHard"), particle.pt(), particle.pt() / pTHat, weight);
366373
}
367374

368375
if (nTT > 0) {
369376
trigNumber = rand->Integer(nTT);
370377
phiTT = phiTTAr[trigNumber];
378+
ptTT = ptTTAr[trigNumber];
371379
if (isSigCol) {
372380
registry.fill(HIST("hNtrig"), 1.5, weight);
373381
registry.fill(HIST("hSigEventTriggers"), nTT, weight);
382+
registry.fill(HIST("hSignalTriggersPtHard"), ptTT / pTHat, weight);
374383
}
375384
if (!isSigCol) {
376385
registry.fill(HIST("hNtrig"), 0.5, weight);
377386
registry.fill(HIST("hRefEventTriggers"), nTT, weight);
387+
registry.fill(HIST("hReferenceTriggersPtHard"), ptTT / pTHat, weight);
378388
}
379389
}
380390

@@ -405,6 +415,7 @@ struct JetHadronRecoil {
405415
registry.fill(HIST("hSignalPtDPhi"), dphi, jet.pt(), weight);
406416
if (std::abs(dphi - o2::constants::math::PI) < 0.6) {
407417
registry.fill(HIST("hSignalPt"), jet.pt(), weight);
418+
registry.fill(HIST("hSignalPtHard"), jet.pt(), ptTT / pTHat, weight);
408419
}
409420
}
410421
if (!isSigCol) {
@@ -416,53 +427,56 @@ struct JetHadronRecoil {
416427
registry.fill(HIST("hReferencePtDPhi"), dphi, jet.pt(), weight);
417428
if (std::abs(dphi - o2::constants::math::PI) < 0.6) {
418429
registry.fill(HIST("hReferencePt"), jet.pt(), weight);
430+
registry.fill(HIST("hReferencePtHard"), jet.pt(), ptTT / pTHat, weight);
419431
}
420432
}
421433
}
422434
}
423435
}
424436

425437
template <typename T, typename V, typename W, typename U, typename X, typename Y>
426-
void fillMatchedHistograms(T const& jetBase, V const& mcdjetsWTA, W const& mcpjetsWTA, U const&, X const& tracks, Y const& particles, float weight = 1.0, float rho = 0.0)
438+
void fillMatchedHistograms(T const& jetsBase, V const& mcdjetsWTA, W const& mcpjetsWTA, U const&, X const& tracks, Y const& particles, float weight = 1.0, float rho = 0.0)
427439
{
428-
double dR = 0;
429-
double dRp = 0;
430-
431-
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
432-
if (jetBase.pt() > pTHatMaxMCD * pTHat) {
433-
return;
434-
}
440+
for (const auto& jetBase : jetsBase) {
441+
double dR = 0;
442+
double dRp = 0;
435443

436-
for (const auto& mcdjetWTA : mcdjetsWTA) {
437-
double djet = RecoDecay::sqrtSumOfSquares(RecoDecay::constrainAngle(jetBase.phi() - mcdjetWTA.phi(), -o2::constants::math::PI), jetBase.eta() - mcdjetWTA.eta());
438-
if (mcdjetWTA.pt() > pTHatMaxMCD * pTHat) {
439-
continue;
440-
}
441-
if (djet < 0.6 * jetR) {
442-
dR = djet;
443-
break;
444+
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
445+
if (jetBase.pt() > pTHatMaxMCD * pTHat) {
446+
return;
444447
}
445-
}
446-
447-
dR = getWTAaxisDifference(jetBase, mcdjetsWTA, tracks, true);
448448

449-
if (jetBase.has_matchedJetGeo()) {
450-
for (const auto& jetTag : jetBase.template matchedJetGeo_as<std::decay_t<U>>()) {
451-
if (jetTag.pt() > pTHatMaxMCP * pTHat) {
449+
for (const auto& mcdjetWTA : mcdjetsWTA) {
450+
double djet = RecoDecay::sqrtSumOfSquares(RecoDecay::constrainAngle(jetBase.phi() - mcdjetWTA.phi(), -o2::constants::math::PI), jetBase.eta() - mcdjetWTA.eta());
451+
if (mcdjetWTA.pt() > pTHatMaxMCD * pTHat) {
452452
continue;
453453
}
454+
if (djet < 0.6 * jetR) {
455+
dR = djet;
456+
break;
457+
}
458+
}
459+
460+
dR = getWTAaxisDifference(jetBase, mcdjetsWTA, tracks, true);
461+
462+
if (jetBase.has_matchedJetGeo()) {
463+
for (const auto& jetTag : jetBase.template matchedJetGeo_as<std::decay_t<U>>()) {
464+
if (jetTag.pt() > pTHatMaxMCP * pTHat) {
465+
return;
466+
}
454467

455-
dRp = getWTAaxisDifference(jetTag, mcpjetsWTA, particles, true);
456-
457-
registry.fill(HIST("hPtMatched"), jetBase.pt() - (rho * jetBase.area()), jetTag.pt(), weight);
458-
registry.fill(HIST("hPhiMatched"), jetBase.phi(), jetTag.phi(), weight);
459-
registry.fill(HIST("hPtResolution"), jetTag.pt(), (jetTag.pt() - (jetBase.pt() - (rho * jetBase.area()))) / jetTag.pt(), weight);
460-
registry.fill(HIST("hPhiResolution"), jetTag.pt(), jetTag.phi() - jetBase.phi(), weight);
461-
registry.fill(HIST("hDeltaRMatched"), dR, dRp, weight);
462-
registry.fill(HIST("hDeltaRResolution"), jetTag.pt(), dRp - dR, weight);
463-
registry.fill(HIST("hFullMatching"), jetBase.pt() - (rho * jetBase.area()), jetTag.pt(), jetBase.phi(), jetTag.phi(), dR, dRp, weight);
464-
registry.fill(HIST("hPtMatched1d"), jetTag.pt(), weight);
465-
registry.fill(HIST("hDeltaRMatched1d"), dRp, weight);
468+
dRp = getWTAaxisDifference(jetTag, mcpjetsWTA, particles, true);
469+
470+
registry.fill(HIST("hPtMatched"), jetBase.pt() - (rho * jetBase.area()), jetTag.pt(), weight);
471+
registry.fill(HIST("hPhiMatched"), jetBase.phi(), jetTag.phi(), weight);
472+
registry.fill(HIST("hPtResolution"), jetTag.pt(), (jetTag.pt() - (jetBase.pt() - (rho * jetBase.area()))) / jetTag.pt(), weight);
473+
registry.fill(HIST("hPhiResolution"), jetTag.pt(), jetTag.phi() - jetBase.phi(), weight);
474+
registry.fill(HIST("hDeltaRMatched"), dR, dRp, weight);
475+
registry.fill(HIST("hDeltaRResolution"), jetTag.pt(), dRp - dR, weight);
476+
registry.fill(HIST("hFullMatching"), jetBase.pt() - (rho * jetBase.area()), jetTag.pt(), jetBase.phi(), jetTag.phi(), dR, dRp, weight);
477+
registry.fill(HIST("hPtMatched1d"), jetTag.pt(), weight);
478+
registry.fill(HIST("hDeltaRMatched1d"), dRp, weight);
479+
}
466480
}
467481
}
468482
}
@@ -626,9 +640,7 @@ struct JetHadronRecoil {
626640
}
627641
registry.fill(HIST("hZvtxSelected"), collision.posZ());
628642
const auto& mcpjetsWTACut = mcpjetsWTA.sliceBy(partJetsPerCollision, collision.mcCollisionId());
629-
for (const auto& mcdjet : mcdjets) {
630-
fillMatchedHistograms(mcdjet, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles);
631-
}
643+
fillMatchedHistograms(mcdjets, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles);
632644
}
633645
PROCESS_SWITCH(JetHadronRecoil, processJetsMCPMCDMatched, "process MC matched (inc jets)", false);
634646

@@ -649,20 +661,18 @@ struct JetHadronRecoil {
649661
}
650662
registry.fill(HIST("hZvtxSelected"), collision.posZ());
651663
const auto& mcpjetsWTACut = mcpjetsWTA.sliceBy(partJetsPerCollision, collision.mcCollisionId());
652-
for (const auto& mcdjet : mcdjets) {
653-
fillMatchedHistograms(mcdjet, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles, 1.0, collision.rho());
654-
}
664+
fillMatchedHistograms(mcdjets, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles);
655665
}
656666
PROCESS_SWITCH(JetHadronRecoil, processJetsMCPMCDMatchedWithRhoSubtraction, "process MC matched (inc jets) with rho subtraction", false);
657667

658668
void processJetsMCPMCDMatchedWeighted(soa::Filtered<aod::JetCollisionsMCD>::iterator const& collision,
659-
soa::Filtered<soa::Join<aod::ChargedMCDetectorLevelJets, aod::ChargedMCDetectorLevelJetConstituents, aod::ChargedMCDetectorLevelJetsMatchedToChargedMCParticleLevelJets, aod::ChargedMCDetectorLevelJetEventWeights>> const& mcdjets,
669+
soa::Filtered<soa::Join<aod::ChargedMCDetectorLevelJets, aod::ChargedMCDetectorLevelJetConstituents, aod::ChargedMCDetectorLevelJetsMatchedToChargedMCParticleLevelJets>> const& mcdjets,
660670
soa::Filtered<soa::Join<aod::Charged1MCDetectorLevelJets, aod::Charged1MCDetectorLevelJetConstituents>> const& mcdjetsWTA,
661671
soa::Filtered<soa::Join<aod::Charged1MCParticleLevelJets, aod::Charged1MCParticleLevelJetConstituents>> const& mcpjetsWTA,
662672
aod::JetTracks const& tracks,
663673
aod::JetParticles const& particles,
664674
aod::JetMcCollisions const&,
665-
soa::Filtered<soa::Join<aod::ChargedMCParticleLevelJets, aod::ChargedMCParticleLevelJetConstituents, aod::ChargedMCParticleLevelJetsMatchedToChargedMCDetectorLevelJets, aod::ChargedMCParticleLevelJetEventWeights>> const& mcpjets)
675+
soa::Filtered<soa::Join<aod::ChargedMCParticleLevelJets, aod::ChargedMCParticleLevelJetConstituents, aod::ChargedMCParticleLevelJetsMatchedToChargedMCDetectorLevelJets>> const& mcpjets)
666676
{
667677
if (!jetderiveddatautilities::selectCollision(collision, eventSelectionBits)) {
668678
return;
@@ -672,20 +682,18 @@ struct JetHadronRecoil {
672682
}
673683
registry.fill(HIST("hZvtxSelected"), collision.posZ());
674684
const auto& mcpjetsWTACut = mcpjetsWTA.sliceBy(partJetsPerCollision, collision.mcCollisionId());
675-
for (const auto& mcdjet : mcdjets) {
676-
fillMatchedHistograms(mcdjet, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles, mcdjet.eventWeight());
677-
}
685+
fillMatchedHistograms(mcdjets, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles, collision.mcCollision().weight());
678686
}
679687
PROCESS_SWITCH(JetHadronRecoil, processJetsMCPMCDMatchedWeighted, "process MC matched with event weights (inc jets)", false);
680688

681689
void processJetsMCPMCDMatchedWeightedWithRhoSubtraction(soa::Filtered<soa::Join<aod::JetCollisionsMCD, aod::BkgChargedRhos>>::iterator const& collision,
682-
soa::Filtered<soa::Join<aod::ChargedMCDetectorLevelJets, aod::ChargedMCDetectorLevelJetConstituents, aod::ChargedMCDetectorLevelJetsMatchedToChargedMCParticleLevelJets, aod::ChargedMCDetectorLevelJetEventWeights>> const& mcdjets,
690+
soa::Filtered<soa::Join<aod::ChargedMCDetectorLevelJets, aod::ChargedMCDetectorLevelJetConstituents, aod::ChargedMCDetectorLevelJetsMatchedToChargedMCParticleLevelJets>> const& mcdjets,
683691
soa::Filtered<soa::Join<aod::Charged1MCDetectorLevelJets, aod::Charged1MCDetectorLevelJetConstituents>> const& mcdjetsWTA,
684692
soa::Filtered<soa::Join<aod::Charged1MCParticleLevelJets, aod::Charged1MCParticleLevelJetConstituents>> const& mcpjetsWTA,
685693
aod::JetTracks const& tracks,
686694
aod::JetParticles const& particles,
687695
aod::JetMcCollisions const&,
688-
soa::Filtered<soa::Join<aod::ChargedMCParticleLevelJets, aod::ChargedMCParticleLevelJetConstituents, aod::ChargedMCParticleLevelJetsMatchedToChargedMCDetectorLevelJets, aod::ChargedMCParticleLevelJetEventWeights>> const& mcpjets)
696+
soa::Filtered<soa::Join<aod::ChargedMCParticleLevelJets, aod::ChargedMCParticleLevelJetConstituents, aod::ChargedMCParticleLevelJetsMatchedToChargedMCDetectorLevelJets>> const& mcpjets)
689697
{
690698
if (!jetderiveddatautilities::selectCollision(collision, eventSelectionBits)) {
691699
return;
@@ -695,9 +703,7 @@ struct JetHadronRecoil {
695703
}
696704
registry.fill(HIST("hZvtxSelected"), collision.posZ());
697705
const auto& mcpjetsWTACut = mcpjetsWTA.sliceBy(partJetsPerCollision, collision.mcCollisionId());
698-
for (const auto& mcdjet : mcdjets) {
699-
fillMatchedHistograms(mcdjet, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles, mcdjet.eventWeight(), collision.rho());
700-
}
706+
fillMatchedHistograms(mcdjets, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles, collision.mcCollision().weight(), collision.rho());
701707
}
702708
PROCESS_SWITCH(JetHadronRecoil, processJetsMCPMCDMatchedWeightedWithRhoSubtraction, "process MC matched with event weights (inc jets) and rho subtraction", false);
703709

@@ -726,21 +732,19 @@ struct JetHadronRecoil {
726732
}
727733
}
728734
if (ishJetEvent) {
729-
for (const auto& mcdjet : mcdjets) {
730-
fillMatchedHistograms(mcdjet, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles);
731-
}
735+
fillMatchedHistograms(mcdjets, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles);
732736
}
733737
}
734738
PROCESS_SWITCH(JetHadronRecoil, processRecoilJetsMCPMCDMatched, "process MC matched (recoil jets)", false);
735739

736740
void processRecoilJetsMCPMCDMatchedWeighted(soa::Filtered<aod::JetCollisionsMCD>::iterator const& collision,
737-
soa::Filtered<soa::Join<aod::ChargedMCDetectorLevelJets, aod::ChargedMCDetectorLevelJetConstituents, aod::ChargedMCDetectorLevelJetsMatchedToChargedMCParticleLevelJets, aod::ChargedMCDetectorLevelJetEventWeights>> const& mcdjets,
741+
soa::Filtered<soa::Join<aod::ChargedMCDetectorLevelJets, aod::ChargedMCDetectorLevelJetConstituents, aod::ChargedMCDetectorLevelJetsMatchedToChargedMCParticleLevelJets>> const& mcdjets,
738742
soa::Filtered<soa::Join<aod::Charged1MCDetectorLevelJets, aod::Charged1MCDetectorLevelJetConstituents>> const& mcdjetsWTA,
739743
soa::Filtered<soa::Join<aod::Charged1MCParticleLevelJets, aod::Charged1MCParticleLevelJetConstituents>> const& mcpjetsWTA,
740744
soa::Filtered<aod::JetTracks> const& tracks,
741745
soa::Filtered<aod::JetParticles> const& particles,
742746
aod::JetMcCollisions const&,
743-
soa::Filtered<soa::Join<aod::ChargedMCParticleLevelJets, aod::ChargedMCParticleLevelJetConstituents, aod::ChargedMCParticleLevelJetsMatchedToChargedMCDetectorLevelJets, aod::ChargedMCParticleLevelJetEventWeights>> const& mcpjets)
747+
soa::Filtered<soa::Join<aod::ChargedMCParticleLevelJets, aod::ChargedMCParticleLevelJetConstituents, aod::ChargedMCParticleLevelJetsMatchedToChargedMCDetectorLevelJets>> const& mcpjets)
744748
{
745749
if (!jetderiveddatautilities::selectCollision(collision, eventSelectionBits)) {
746750
return;
@@ -758,9 +762,7 @@ struct JetHadronRecoil {
758762
}
759763
}
760764
if (ishJetEvent) {
761-
for (const auto& mcdjet : mcdjets) {
762-
fillMatchedHistograms(mcdjet, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles, mcdjet.eventWeight());
763-
}
765+
fillMatchedHistograms(mcdjets, mcdjetsWTA, mcpjetsWTACut, mcpjets, tracks, particles, collision.mcCollision().weight());
764766
}
765767
}
766768
PROCESS_SWITCH(JetHadronRecoil, processRecoilJetsMCPMCDMatchedWeighted, "process MC matched with event weights (recoil jets)", false);

0 commit comments

Comments
 (0)