diff --git a/PWGJE/Tasks/jetFragmentation.cxx b/PWGJE/Tasks/jetFragmentation.cxx index 5a6104ad872..997a475400f 100644 --- a/PWGJE/Tasks/jetFragmentation.cxx +++ b/PWGJE/Tasks/jetFragmentation.cxx @@ -48,34 +48,14 @@ using namespace o2; using namespace o2::framework; using namespace o2::framework::expressions; -// Charged jets -using MyTracks = soa::Join; -using ChargedJetsWithConstituents = soa::Join; - -using MCDJets = aod::ChargedMCDetectorLevelJets; -using MCDJetsWithConstituents = soa::Join; -using MatchedMCDJets = soa::Join; -using MatchedMCDJetsWithConstituents = soa::Join; - -using MCPJets = aod::ChargedMCParticleLevelJets; -using MatchedMCPJets = soa::Join; -using MCPJetsWithConstituents = soa::Join; -using MatchedMCPJetsWithConstituents = soa::Join; - -// V0 jets using DataV0JetsWithConstituents = soa::Join; -using CandidatesV0DataWithFlags = aod::CandidatesV0Data; -using CandidatesV0MCDWithLabelsAndFlags = soa::Join; -using MCDV0Jets = aod::V0ChargedMCDetectorLevelJets; -using MCDV0JetsWithConstituents = soa::Join; -using MatchedMCDV0Jets = soa::Join; -using MatchedMCDV0JetsWithConstituents = soa::Join; +using CandidatesV0MCDWithLabels = soa::Join; +using MatchedMCDV0Jets = soa::Join; +using MatchedMCDV0JetsWithConstituents = soa::Join; -using MCPV0Jets = aod::V0ChargedMCParticleLevelJets; -using MCPV0JetsWithConstituents = soa::Join; -using MatchedMCPV0Jets = soa::Join; -using MatchedMCPV0JetsWithConstituents = soa::Join; +using MatchedMCPV0Jets = soa::Join; +using MatchedMCPV0JetsWithConstituents = soa::Join; struct JetFragmentation { HistogramRegistry registry{"registry"}; // CallSumw2 = false? @@ -240,6 +220,7 @@ struct JetFragmentation { AxisSpec antiLambdaMassCutAxis = {binAntiLambdaMassCut, "Inv. mass (GeV/#it{c}^{2})"}; if (doprocessDataV0) { + registry.add("data/hEvents", "hEvents", {HistType::kTH1D, {{2, 0.0f, 2.0f}}}); registry.add("data/V0/nV0sEvent", "nV0sEvent", HistType::kTH1D, {v0Count}); registry.add("data/V0/nV0sEventAcc", "nV0s per event (accepted)", HistType::kTH1D, {v0Count}); registry.add("data/V0/nV0sEventAccWeighted", "nV0s per event (accepted, weighted)", HistType::kTH1D, {v0Weight}, true); @@ -342,11 +323,10 @@ struct JetFragmentation { // Cut variation registry.add("data/V0/V0CutVariation", "V0CutVariation", HistType::kTHnSparseD, {ptCutAxis, k0SMassCutAxis, lambdaMassCutAxis, antiLambdaMassCutAxis, rCutAxis, ctauCutAxis, cosPACutAxis, dcapCutAxis, dcanCutAxis, dcadCutAxis}); - } // doprocessDataV0 - if (doprocessDataV0JetsFrag) { + // Jets registry.add("data/jets/jetPtEtaPhi", "Jet #it{p}_{T}, #eta, #phi", HistType::kTH3D, {jetPtAxis, etaAxis, phiAxis}); - + registry.add("data/jets/inclJetPtEtaPhi", "incljetPtEtaPhi", HistType::kTH3D, {jetPtAxis, etaAxis, phiAxis}); registry.add("data/jets/V0/jetPtV0TrackProj", "jetPtV0TrackProj", HistType::kTH2D, {jetPtAxis, zAxis}); registry.add("data/jets/V0/jetPtnV0nK0SnLambdanAntiLambda", "jetPtnV0nK0SnLambdanAntiLambda", HistType::kTHnSparseD, {jetPtAxis, v0Count, v0Count, v0Count, v0Count}); @@ -547,9 +527,51 @@ struct JetFragmentation { registry.add("data/jets/weighted/V0/jetPtBkgTrackProjDCAposAntiLambdaMass", "jetPtBkgTrackProjDCAposAntiLambdaMass", HistType::kTHnSparseD, {jetPtAxis, zAxis, v0DCApAxis, lambdaMassAxis}, true); registry.add("data/jets/weighted/V0/jetPtBkgTrackProjDCAnegAntiLambdaMass", "jetPtBkgTrackProjDCAnegAntiLambdaMass", HistType::kTHnSparseD, {jetPtAxis, zAxis, v0DCAnAxis, lambdaMassAxis}, true); registry.add("data/jets/weighted/V0/jetPtBkgTrackProjDCAdAntiLambdaMass", "jetPtBkgTrackProjDCAdAntiLambdaMass", HistType::kTHnSparseD, {jetPtAxis, zAxis, v0DCAdAxis, lambdaMassAxis}, true); - } - if (doprocessMcMatchedV0) { + // Perp Cone + registry.add("data/PC/JetPtEtaV0Pt", "JetPtEtaV0Pt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); + registry.add("data/PC/V0PtEtaPhi", "V0 #it{p}_{T}, #eta, #phi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); + registry.add("data/PC/V0PtCtau", "V0PtCtau", HistType::kTHnSparseD, {v0PtAxis, v0CtauAxis, v0CtauAxis, v0CtauAxis}); + registry.add("data/PC/V0PtMass", "V0PtMass", HistType::kTHnSparseD, {v0PtAxis, k0SMassAxis, lambdaMassAxis, lambdaMassAxis}); + registry.add("data/PC/V0PtMassWide", "V0PtMassWide", HistType::kTHnSparseD, {v0PtAxis, k0SWideAxis, lambdaMassAxis, lambdaMassAxis}); + registry.add("data/PC/V0PtRadiusCosPA", "V0PtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); + registry.add("data/PC/V0PtDCAposneg", "V0PtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); + registry.add("data/PC/V0PtDCAd", "V0PtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); + + // K0S + registry.add("data/PC/JetPtK0SPtMass", "JetPtK0SPtMass", HistType::kTH3D, {jetPtAxis, v0PtAxis, k0SMassAxis}); + registry.add("data/PC/JetPtEtaK0SPt", "JetPtEtaK0SPt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); + registry.add("data/PC/K0SPtEtaPhi", "K0SPtEtaPhi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); + registry.add("data/PC/K0SPtCtauMass", "K0SPtCtauMass", HistType::kTH3D, {v0PtAxis, v0CtauAxis, k0SMassAxis}); + registry.add("data/PC/K0SPtRadiusCosPA", "K0SPtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); + registry.add("data/PC/K0SPtDCAposneg", "K0SPtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); + registry.add("data/PC/K0SPtDCAd", "K0SPtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); + + // Lambda + registry.add("data/PC/JetPtLambdaPtMass", "JetPtLambdaPtMass", HistType::kTH3D, {jetPtAxis, v0PtAxis, lambdaMassAxis}); + registry.add("data/PC/JetPtEtaLambdaPt", "JetPtEtaLambdaPt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); + registry.add("data/PC/LambdaPtEtaPhi", "LambdaPtEtaPhi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); + registry.add("data/PC/LambdaPtCtauMass", "LambdaPtCtauMass", HistType::kTH3D, {v0PtAxis, v0CtauAxis, lambdaMassAxis}); + registry.add("data/PC/LambdaPtRadiusCosPA", "LambdaPtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); + registry.add("data/PC/LambdaPtDCAposneg", "LambdaPtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); + registry.add("data/PC/LambdaPtDCAd", "LambdaPtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); + + // AntiLambda + registry.add("data/PC/JetPtAntiLambdaPtMass", "JetPtAntiLambdaPtMass", HistType::kTH3D, {jetPtAxis, v0PtAxis, lambdaMassAxis}); + registry.add("data/PC/JetPtEtaAntiLambdaPt", "JetPtEtaAntiLambdaPt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); + registry.add("data/PC/AntiLambdaPtEtaPhi", "AntiLambdaPtEtaPhi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); + registry.add("data/PC/AntiLambdaPtCtauMass", "AntiLambdaPtCtauMass", HistType::kTH3D, {v0PtAxis, v0CtauAxis, lambdaMassAxis}); + registry.add("data/PC/AntiLambdaPtRadiusCosPA", "AntiLambdaPtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); + registry.add("data/PC/AntiLambdaPtDCAposneg", "AntiLambdaPtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); + registry.add("data/PC/AntiLambdaPtDCAd", "AntiLambdaPtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); + + registry.add("data/PC/nV0sConePtEta", "nV0sConePtEta", HistType::kTH3D, {v0Count, jetPtAxis, etaAxis}); + registry.add("data/PC/ConePtEtaPhi", "ConePtEtaPhi", HistType::kTH3D, {jetPtAxis, etaAxis, phiAxis}); + registry.add("data/PC/JetPtEtaConePt", "JetPtEtaConePt", HistType::kTH3D, {jetPtAxis, etaAxis, jetPtAxis}); + } // doprocessDataV0 + + if (doprocessMcV0) { + registry.add("matching/hEvents", "hEvents", {HistType::kTH1D, {{3, 0.0f, 3.0f}}}); // MCP registry.add("mcp/V0/nV0sEventAcc", "nV0s per event (accepted)", HistType::kTH1D, {v0Count}, true); registry.add("mcp/V0/nV0sEventAccWeighted", "nV0s per event (accepted, weighted)", HistType::kTH1D, {v0Weight}, true); @@ -716,13 +738,18 @@ struct JetFragmentation { registry.add("matching/V0/K0SPosNegPtMass", "K0SPosNegPtMass", HistType::kTHnSparseD, {v0partPtAxis, trackPtAxis, trackPtAxis, k0SMassAxis}, true); registry.add("matching/V0/LambdaPosNegPtMass", "LambdaPosNegPtMass", HistType::kTHnSparseD, {v0partPtAxis, trackPtAxis, trackPtAxis, lambdaMassAxis}, true); registry.add("matching/V0/AntiLambdaPosNegPtMass", "AntiLambdaPosNegPtMass", HistType::kTHnSparseD, {v0partPtAxis, trackPtAxis, trackPtAxis, lambdaMassAxis}, true); - } // doprocessMcMatchedV0 - if (doprocessMcMatchedV0JetsFrag) { // Matching - Jets + registry.add("mcp/jets/inclPartJetPtEtaPhi", "inclPartJetPtEtaPhi", HistType::kTH3D, {partJetPtAxis, partEtaAxis, partPhiAxis}, true); registry.add("mcp/jets/partJetPtEtaPhi", "Particle level jet #it{p}_{T}, #eta, #phi", HistType::kTH3D, {partJetPtAxis, partEtaAxis, partPhiAxis}, true); + registry.add("mcd/jets/inclDetJetPtEtaPhi", "inclDetJetPtEtaPhi", HistType::kTH3D, {detJetPtAxis, detEtaAxis, detPhiAxis}, true); registry.add("mcd/jets/detJetPtEtaPhi", "Detector level jet #it{p}_{T}, #eta, #phi", HistType::kTH3D, {detJetPtAxis, detEtaAxis, detPhiAxis}, true); + registry.add("matching/jets/fakeJetPtEtaPhi", "fakeJetPtEtaPhi", HistType::kTH3D, {detJetPtAxis, detEtaAxis, detPhiAxis}, true); + registry.add("matching/jets/inclFakeJetPtEtaPhi", "inclFakeJetPtEtaPhi", HistType::kTH3D, {detJetPtAxis, detEtaAxis, detPhiAxis}, true); + registry.add("matching/jets/inclMissJetPtEtaPhi", "inclMissJetPtEtaPhi", HistType::kTH3D, {partJetPtAxis, partEtaAxis, partPhiAxis}, true); + registry.add("matching/jets/missJetPtEtaPhi", "missJetPtEtaPhi", HistType::kTH3D, {partJetPtAxis, partEtaAxis, partPhiAxis}, true); + registry.add("matching/jets/matchDetJetPtEtaPhi", "Matched detector level jet #it{p}_{T}, #eta, #phi", HistType::kTH3D, {detJetPtAxis, detEtaAxis, detPhiAxis}, true); registry.add("matching/jets/matchPartJetPtEtaPhi", "Matched particle level jet #it{p}_{T}, #eta, #phi", HistType::kTH3D, {partJetPtAxis, partEtaAxis, partPhiAxis}, true); registry.add("matching/jets/matchPartJetPtEtaPhiMatchDist", "matchJetMatchDist", HistType::kTHnSparseD, {partJetPtAxis, partEtaAxis, partPhiAxis, matchDistAxis}, true); @@ -735,6 +762,18 @@ struct JetFragmentation { registry.add("matching/jets/matchPartJetPtResolutionPhi", "#phi^{jet, det} - #phi^{jet, part}", HistType::kTH3D, {partJetPtAxis, partPhiAxis, phiDiffAxis}, true); registry.add("matching/jets/matchPartJetPtRelDiffPt", "#it{p}_{T}^{jet, det} - #it{p}_{T}^{jet, part}", HistType::kTH2D, {partJetPtAxis, ptJetRelDiffAxis}, true); + registry.add("matching/jets/inclMatchDetJetPtEtaPhi", "inclMatchDetJetPtEtaPhi", HistType::kTH3D, {detJetPtAxis, detEtaAxis, detPhiAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtEtaPhi", "inclMatchPartJetPtEtaPhi", HistType::kTH3D, {partJetPtAxis, partEtaAxis, partPhiAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtEtaPhiMatchDist", "matchJetMatchDist", HistType::kTHnSparseD, {partJetPtAxis, partEtaAxis, partPhiAxis, matchDistAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtEnergyScale", "jetEnergyScale", HistType::kTH2D, {partJetPtAxis, ptRatioAxis}, true); + registry.add("matching/jets/inclMatchDetJetPtPartJetPt", "matchDetJetPtPartJetPt", HistType::kTH2D, {detJetPtAxis, partJetPtAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtDetJetEtaPartJetEta", "matchPartJetPtDetJetEtaPartJetEta", HistType::kTH3D, {partJetPtAxis, detEtaAxis, partEtaAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtDetJetPhiPartJetPhi", "matchPartJetPtDetJetPhiPartJetPhi", HistType::kTH3D, {partJetPtAxis, detPhiAxis, partPhiAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtResolutionPt", "inclMatchPartJetPtResolutionPt", HistType::kTH2D, {partJetPtAxis, ptDiffAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtResolutionEta", "inclMatchPartJetPtResolutionEta", HistType::kTH3D, {partJetPtAxis, partEtaAxis, etaDiffAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtResolutionPhi", "inclMatchPartJetPtResolutionPhi", HistType::kTH3D, {partJetPtAxis, partPhiAxis, phiDiffAxis}, true); + registry.add("matching/jets/inclMatchPartJetPtRelDiffPt", "inclMatchPartJetPtRelDiffPt", HistType::kTH2D, {partJetPtAxis, ptJetRelDiffAxis}, true); + registry.add("matching/jets/V0/jetPtnV0MatchednK0SnLambdanAntiLambda", "jet Pt, nV0 matched, nK0S nLambdan AntiLambda", HistType::kTHnSparseD, {detJetPtAxis, v0Count, v0Count, v0Count, v0Count}, true); registry.add("matching/jets/V0/partJetPtV0PtDetPt", "V0PartPtDetPt", HistType::kTH3D, {partJetPtAxis, v0partPtAxis, v0detPtAxis}, true); registry.add("matching/jets/V0/partJetPtDetJetPtPartV0PtRatioPtRelDiffPt", "V0PartPtRatioRelDiffPt", HistType::kTHnSparseD, {partJetPtAxis, detJetPtAxis, v0partPtAxis, v0PtRatioAxis, v0PtRelDiffAxis}, true); @@ -980,51 +1019,7 @@ struct JetFragmentation { registry.add("matching/V0/decayedOtherPtV0PtMass", "decayedOtherPtV0PtMass", HistType::kTHnSparseD, {v0partPtAxis, v0detPtAxis, k0SMassAxis, lambdaMassAxis, lambdaMassAxis}, true); registry.add("matching/jets/V0/decayedOtherPtV0PtMass", "decayedOtherPtV0PtMass", HistType::kTHnSparseD, {partJetPtAxis, detJetPtAxis, v0partPtAxis, v0detPtAxis, k0SMassAxis, lambdaMassAxis, lambdaMassAxis}, true); registry.add("matching/jets/V0/decayedOtherPtV0TrackProjMass", "decayedOtherPtV0TrackProjMass", HistType::kTHnSparseD, {partJetPtAxis, detJetPtAxis, partZAxis, detZAxis, k0SMassAxis, lambdaMassAxis, lambdaMassAxis}, true); - } // doprocessMcMatchedV0JetsFrag - - if (doprocessDataV0PerpCone) { - registry.add("data/PC/JetPtEtaV0Pt", "JetPtEtaV0Pt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); - registry.add("data/PC/V0PtEtaPhi", "V0 #it{p}_{T}, #eta, #phi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); - registry.add("data/PC/V0PtCtau", "V0PtCtau", HistType::kTHnSparseD, {v0PtAxis, v0CtauAxis, v0CtauAxis, v0CtauAxis}); - registry.add("data/PC/V0PtMass", "V0PtMass", HistType::kTHnSparseD, {v0PtAxis, k0SMassAxis, lambdaMassAxis, lambdaMassAxis}); - registry.add("data/PC/V0PtMassWide", "V0PtMassWide", HistType::kTHnSparseD, {v0PtAxis, k0SWideAxis, lambdaMassAxis, lambdaMassAxis}); - registry.add("data/PC/V0PtRadiusCosPA", "V0PtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); - registry.add("data/PC/V0PtDCAposneg", "V0PtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); - registry.add("data/PC/V0PtDCAd", "V0PtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); - // K0S - registry.add("data/PC/JetPtK0SPtMass", "JetPtK0SPtMass", HistType::kTH3D, {jetPtAxis, v0PtAxis, k0SMassAxis}); - registry.add("data/PC/JetPtEtaK0SPt", "JetPtEtaK0SPt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); - registry.add("data/PC/K0SPtEtaPhi", "K0SPtEtaPhi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); - registry.add("data/PC/K0SPtCtauMass", "K0SPtCtauMass", HistType::kTH3D, {v0PtAxis, v0CtauAxis, k0SMassAxis}); - registry.add("data/PC/K0SPtRadiusCosPA", "K0SPtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); - registry.add("data/PC/K0SPtDCAposneg", "K0SPtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); - registry.add("data/PC/K0SPtDCAd", "K0SPtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); - - // Lambda - registry.add("data/PC/JetPtLambdaPtMass", "JetPtLambdaPtMass", HistType::kTH3D, {jetPtAxis, v0PtAxis, lambdaMassAxis}); - registry.add("data/PC/JetPtEtaLambdaPt", "JetPtEtaLambdaPt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); - registry.add("data/PC/LambdaPtEtaPhi", "LambdaPtEtaPhi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); - registry.add("data/PC/LambdaPtCtauMass", "LambdaPtCtauMass", HistType::kTH3D, {v0PtAxis, v0CtauAxis, lambdaMassAxis}); - registry.add("data/PC/LambdaPtRadiusCosPA", "LambdaPtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); - registry.add("data/PC/LambdaPtDCAposneg", "LambdaPtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); - registry.add("data/PC/LambdaPtDCAd", "LambdaPtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); - - // AntiLambda - registry.add("data/PC/JetPtAntiLambdaPtMass", "JetPtAntiLambdaPtMass", HistType::kTH3D, {jetPtAxis, v0PtAxis, lambdaMassAxis}); - registry.add("data/PC/JetPtEtaAntiLambdaPt", "JetPtEtaAntiLambdaPt", HistType::kTH3D, {jetPtAxis, etaAxis, v0PtAxis}); - registry.add("data/PC/AntiLambdaPtEtaPhi", "AntiLambdaPtEtaPhi", HistType::kTH3D, {v0PtAxis, v0EtaAxis, v0PhiAxis}); - registry.add("data/PC/AntiLambdaPtCtauMass", "AntiLambdaPtCtauMass", HistType::kTH3D, {v0PtAxis, v0CtauAxis, lambdaMassAxis}); - registry.add("data/PC/AntiLambdaPtRadiusCosPA", "AntiLambdaPtRadiusCosPA", HistType::kTH3D, {v0PtAxis, v0RadiusAxis, v0CosPAAxis}); - registry.add("data/PC/AntiLambdaPtDCAposneg", "AntiLambdaPtDCAposneg", HistType::kTH3D, {v0PtAxis, v0DCApAxis, v0DCAnAxis}); - registry.add("data/PC/AntiLambdaPtDCAd", "AntiLambdaPtDCAd", HistType::kTH2D, {v0PtAxis, v0DCAdAxis}); - - registry.add("data/PC/nV0sConePtEta", "nV0sConePtEta", HistType::kTH3D, {v0Count, jetPtAxis, etaAxis}); - registry.add("data/PC/ConePtEtaPhi", "ConePtEtaPhi", HistType::kTH3D, {jetPtAxis, etaAxis, phiAxis}); - registry.add("data/PC/JetPtEtaConePt", "JetPtEtaConePt", HistType::kTH3D, {jetPtAxis, etaAxis, jetPtAxis}); - } // doprocessDataV0PerpCone - - if (doprocessMcV0PerpCone) { registry.add("mcd/V0/nV0sEvent", "NV0s in event", HistType::kTH1D, {v0Count}); registry.add("mcd/V0/nV0sEventWeighted", "NV0s in event weighted", HistType::kTH1D, {v0Count}, true); @@ -1095,7 +1090,7 @@ struct JetFragmentation { // PerpCone - Matched AntiLambda in Matched jets registry.add("matching/PC/matchedJetPtAntiLambdaPtMass", "matchedJetPtAntiLambdaPtMass", HistType::kTH3D, {detJetPtAxis, v0PtAxis, lambdaMassAxis}, true); registry.add("matching/PC/matchedJetsPtAntiLambdaPtMass", "matchedJetsPtAntiLambdaPtMass", HistType::kTHnSparseD, {partJetPtAxis, detJetPtAxis, v0PtAxis, lambdaMassAxis}, true); - } // doprocessMcV0PerpCone + } // doprocessMcV0 } // init // --------------------------------------------------- @@ -1174,7 +1169,6 @@ struct JetFragmentation { // Convert state from uint32_t to std::vector containing the particle classes for that weight std::vector convertState(uint32_t state, int nParticles, int nClasses = 4) { - // TODO: Should we cast these to uint32_t? std::vector v(nParticles, nClasses); int nStates = std::round(std::pow(static_cast(nClasses), static_cast(nParticles))); int nBitsPerParticle = std::round(std::log2(nClasses)); @@ -1307,11 +1301,11 @@ struct JetFragmentation { return (jet.candidatesIds().size() > 0); } template - bool v0sAreMatched(T const& v0, U const& particle, V const& /*tracks*/) + // bool v0sAreMatched(T const& v0, U const& particle, V const& /*tracks*/) + bool v0sAreMatched(U const& v0, V const& particle) { - // FIXME: Should this check whether V0 is selected as correct species? Probably! - auto negId = v0.template negTrack_as().mcParticleId(); - auto posId = v0.template posTrack_as().mcParticleId(); + auto negId = v0.template negTrack_as().mcParticleId(); + auto posId = v0.template posTrack_as().mcParticleId(); auto daughters = particle.daughtersIds(); return ((negId == daughters[0] && posId == daughters[1]) || (posId == daughters[0] && negId == daughters[1])); } @@ -1353,7 +1347,7 @@ struct JetFragmentation { // --------------------------------------------------- // Data - Counts template - void fillDataV0Histograms(CollisionType const& collision, V0Type const& V0s) + void fillDataV0sInclusive(CollisionType const& coll, V0Type const& V0s) { // Fill histograms unweighted. Hists will be filled with V0 counts float nV0s = 0; @@ -1362,9 +1356,9 @@ struct JetFragmentation { continue; nV0s += 1; - double ctauLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0; - double ctauAntiLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0Bar; - double ctauK0s = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassK0Short; + double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; + double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; + double ctauK0s = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassK0Short; double massDiff = v0.mLambda() - v0.mAntiLambda(); double massRatio = v0.mAntiLambda() / v0.mLambda(); @@ -1426,17 +1420,23 @@ struct JetFragmentation { registry.fill(HIST("data/V0/nV0sEventAcc"), nV0s); } template - void fillDataJetHistograms(T const& jet) + void fillDataJet(T const& jet) + { + registry.fill(HIST("data/jets/inclJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi()); + if (jetContainsV0s(jet)) + registry.fill(HIST("data/jets/jetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi()); + } + void fillDataJetWeighted(double jetpt, double jeteta, double jetphi, double weight = 1.) { - registry.fill(HIST("data/jets/jetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi()); + registry.fill(HIST("data/jets/weighted/jetPtEtaPhi"), jetpt, jeteta, jetphi, weight); } template - void fillDataV0FragHistograms(CollisionType const& collision, JetType const& jet, V0Type const& v0) + void fillDataV0sInJet(CollisionType const& coll, JetType const& jet, V0Type const& v0) { double trackProj = getMomFrac(jet, v0); - double ctauLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0; - double ctauAntiLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0Bar; - double ctauK0s = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassK0Short; + double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; + double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; + double ctauK0s = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassK0Short; double massDiff = v0.mLambda() - v0.mAntiLambda(); double massRatio = v0.mAntiLambda() / v0.mLambda(); @@ -1553,7 +1553,7 @@ struct JetFragmentation { } } template - void fillDataPerpConeHists(T const& coll, U const& jet, V const& v0s) + void fillDataV0sInPerpCone(T const& coll, U const& jet, V const& v0s) { const int nCones = 2; double perpConeR = jet.r() * 1e-2; @@ -1633,7 +1633,7 @@ struct JetFragmentation { } // Data - Weighted template - void fillDataV0HistogramsWeighted(CollisionType const& collision, V0Type const& V0s) + void fillDataV0sInclusiveWeighted(CollisionType const& coll, V0Type const& V0s) { // Fill histograms with V0 signal weights float nV0s = 0; @@ -1644,9 +1644,9 @@ struct JetFragmentation { float weight = getV0SignalProb(v0); nV0s += weight; // Sum weights (purity) of V0s - double ctauLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0; - double ctauAntiLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0Bar; - double ctauK0s = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassK0Short; + double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; + double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; + double ctauK0s = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassK0Short; double massDiff = v0.mLambda() - v0.mAntiLambda(); double massRatio = v0.mAntiLambda() / v0.mLambda(); @@ -1705,25 +1705,21 @@ struct JetFragmentation { } // for v0 registry.fill(HIST("data/V0/nV0sEventAccWeighted"), nV0s); } - void fillDataJetHistogramsWeighted(double jetpt, double jeteta, double jetphi, double weight = 1.) - { - registry.fill(HIST("data/jets/weighted/jetPtEtaPhi"), jetpt, jeteta, jetphi, weight); - } template - void fillDataV0FragHistogramsWeighted(C const& collision, J const& jet, std::vector state, std::vector values, double weight) + void fillDataV0sInJetWeighted(C const& coll, J const& jet, std::vector state, std::vector values, double weight) { double jetpt = values[values.size() - 1]; int ip = 0; - for (const auto& v0 : jet.template candidates_as()) { + for (const auto& v0 : jet.template candidates_as()) { if (v0.isRejectedCandidate()) continue; double z = values[ip]; ip++; - double ctauK0s = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassK0Short; - double ctauLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0; - double ctauAntiLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0Bar; + double ctauK0s = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassK0Short; + double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; + double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; double massDiff = v0.mLambda() - v0.mAntiLambda(); double massRatio = v0.mAntiLambda() / v0.mLambda(); @@ -1818,7 +1814,7 @@ struct JetFragmentation { // MC - Counts (or event weights) template - void fillMCPV0Histograms(T const& pV0s, double weight = 1.) + void fillMcpV0sInclusive(T const& pV0s, double weight = 1.) { float nV0s = 0; for (const auto& pv0 : pV0s) { @@ -1838,12 +1834,15 @@ struct JetFragmentation { registry.fill(HIST("mcp/V0/nV0sEventAccWeighted"), nV0s, weight); } template - void fillMCPJetHistograms(T const& jet, double weight = 1.) + void fillMcpJet(T const& jet, double weight = 1.) { - registry.fill(HIST("mcp/jets/partJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + registry.fill(HIST("mcp/jets/inclPartJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + if (jetContainsV0s(jet)) { + registry.fill(HIST("mcp/jets/partJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + } } template - void fillMCDV0Histograms(T const& coll, U const& V0s, double weight = 1.) + void fillMcdV0sInclusive(T const& coll, U const& V0s, double weight = 1.) { float nV0s = 0; for (const auto& v0 : V0s) { @@ -1912,14 +1911,50 @@ struct JetFragmentation { registry.fill(HIST("mcd/V0/nV0sEventAcc"), nV0s); registry.fill(HIST("mcd/V0/nV0sEventAccWeighted"), nV0s, weight); } + template + void fillMatchingV0sInclusive(V const& coll, W const& V0s, X const& pV0s, double weight = 1.) + { + float nV0s = 0; + for (const auto& v0 : V0s) { + if (v0.isRejectedCandidate()) + continue; + + if (!v0.has_mcParticle()) { + fillFakeV0Inclusive(coll, v0, weight); + fillFakeV0DaughtersInclusive(v0, weight); + fillFakeV0DecayedInclusive(v0, weight); + continue; + } + for (const auto& pv0 : pV0s) { + if (v0sAreMatched(v0, pv0)) { + nV0s += 1; + fillMatchedV0Inclusive(coll, v0, pv0, weight); + fillMatchedV0DaughtersInclusive(v0, pv0, weight); + } + } + } // Reconstructed V0s + for (const auto& pv0 : pV0s) { + for (const auto& v0 : V0s) { + if (v0sAreMatched(v0, pv0)) + continue; + + fillMissV0Inclusive(pv0); + } + } + registry.fill(HIST("matching/V0/nV0sEventAcc"), nV0s); + registry.fill(HIST("matching/V0/nV0sEventAccWeighted"), nV0s, weight); + } template - void fillMCDJetHistograms(T const& jet, double weight = 1.) + void fillMcdJet(T const& jet, double weight = 1.) { - registry.fill(HIST("mcd/jets/detJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + registry.fill(HIST("mcd/jets/inclDetJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + if (jetContainsV0s(jet)) { + registry.fill(HIST("mcd/jets/detJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + } } // Reconstructed V0s in the cone of MCD jets template - void fillMcPerpConeHists(T const& coll, U const& mcdjet, V const& v0s, W const& /* V0 particles */, double weight = 1.) + void fillMcV0sInPerpCone(T const& coll, U const& mcdjet, V const& v0s, W const& /* V0 particles */, double weight = 1.) { const int nCones = 2; double perpConeR = mcdjet.r() * 1e-2; @@ -1989,7 +2024,7 @@ struct JetFragmentation { } // Reconstructed V0s in the cone of matched jets template - void fillMcPerpConeHists(T const& coll, U const& mcdjet, V const& mcpjet, W const& v0s, X const& /* V0 particles */, double weight = 1.) + void fillMcV0sInMatchedPerpCone(T const& coll, U const& mcdjet, V const& mcpjet, W const& v0s, X const& /* V0 particles */, double weight = 1.) { const int nCones = 2; double perpConeR = mcdjet.r() * 1e-2; @@ -2066,11 +2101,11 @@ struct JetFragmentation { } // Matched - Counts (or event weights) template // Reconstructed signal for inclusive V0s - void fillMatchingV0Histograms(CollisionType const& collision, V0Type const& v0, particleType const& particle, double weight = 1.) + void fillMatchedV0Inclusive(CollisionType const& coll, V0Type const& v0, particleType const& particle, double weight = 1.) { - double ctauLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0; - double ctauAntiLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0Bar; - double ctauK0s = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassK0Short; + double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; + double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; + double ctauK0s = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassK0Short; registry.fill(HIST("matching/V0/V0PartPtDetPt"), particle.pt(), v0.pt(), weight); registry.fill(HIST("matching/V0/V0PartPtRatioPtRelDiffPt"), particle.pt(), v0.pt() / particle.pt(), (v0.pt() - particle.pt()) / particle.pt(), weight); @@ -2107,7 +2142,7 @@ struct JetFragmentation { } } template // Reconstructed signal for inclusive V0s: daughters - void fillMatchingV0DauHistograms(V0Type const& v0, ParticleType const& pv0, double weight = 1.) + void fillMatchedV0DaughtersInclusive(V0Type const& v0, ParticleType const& pv0, double weight = 1.) { auto negTrack = v0.template negTrack_as(); auto posTrack = v0.template posTrack_as(); @@ -2125,12 +2160,27 @@ struct JetFragmentation { } } template // Reconstructed jets - void fillMatchingHistogramsJet(DetJet const& detJet, PartJet const& partJet, double weight = 1.) + void fillMatchedJet(DetJet const& detJet, PartJet const& partJet, double weight = 1.) { double deltaEta = detJet.eta() - partJet.eta(); double deltaPhi = RecoDecay::constrainAngle(detJet.phi() - partJet.phi(), -constants::math::PI); double dR = jetutilities::deltaR(detJet, partJet); + registry.fill(HIST("matching/jets/inclMatchDetJetPtEtaPhi"), detJet.pt(), detJet.eta(), detJet.phi(), weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtEtaPhi"), partJet.pt(), partJet.eta(), partJet.phi(), weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtEtaPhiMatchDist"), partJet.pt(), partJet.eta(), partJet.phi(), dR, weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtEnergyScale"), partJet.pt(), detJet.pt() / partJet.pt(), weight); + registry.fill(HIST("matching/jets/inclMatchDetJetPtPartJetPt"), detJet.pt(), partJet.pt(), weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtDetJetEtaPartJetEta"), partJet.pt(), detJet.eta(), partJet.eta(), weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtDetJetPhiPartJetPhi"), partJet.pt(), detJet.phi(), partJet.phi(), weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtResolutionPt"), partJet.pt(), (detJet.pt() - partJet.pt()), weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtResolutionEta"), partJet.pt(), partJet.eta(), deltaEta, weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtResolutionPhi"), partJet.pt(), partJet.phi(), deltaPhi, weight); + registry.fill(HIST("matching/jets/inclMatchPartJetPtRelDiffPt"), partJet.pt(), (detJet.pt() - partJet.pt()) / partJet.pt(), weight); + + if (!jetContainsV0s(detJet)) + return; + registry.fill(HIST("matching/jets/matchDetJetPtEtaPhi"), detJet.pt(), detJet.eta(), detJet.phi(), weight); registry.fill(HIST("matching/jets/matchPartJetPtEtaPhi"), partJet.pt(), partJet.eta(), partJet.phi(), weight); registry.fill(HIST("matching/jets/matchPartJetPtEtaPhiMatchDist"), partJet.pt(), partJet.eta(), partJet.phi(), dR, weight); @@ -2144,15 +2194,15 @@ struct JetFragmentation { registry.fill(HIST("matching/jets/matchPartJetPtRelDiffPt"), partJet.pt(), (detJet.pt() - partJet.pt()) / partJet.pt(), weight); } template // Reconstructed signal for in-jet V0s - void fillMatchingV0FragHistograms(CollisionType const& collision, DetJetType const& detJet, PartJetType const& partJet, V0Type const& v0, ParticleType const& particle, double weight = 1.) + void fillMatchedV0InJet(CollisionType const& coll, DetJetType const& detJet, PartJetType const& partJet, V0Type const& v0, ParticleType const& particle, double weight = 1.) { - bool correctCollision = (collision.mcCollisionId() == particle.mcCollisionId()); + bool correctCollision = (coll.mcCollisionId() == particle.mcCollisionId()); double detTrackProj = getMomFrac(detJet, v0); double partTrackProj = getMomFrac(partJet, particle); - double ctauLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0; - double ctauAntiLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0Bar; - double ctauK0s = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassK0Short; + double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; + double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; + double ctauK0s = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassK0Short; registry.fill(HIST("matching/jets/V0/matchDetJetPtV0TrackProjPartJetPtV0TrackProj"), detJet.pt(), detTrackProj, partJet.pt(), partTrackProj, weight); registry.fill(HIST("matching/jets/V0/partJetPtV0PtDetJetPtV0Pt"), partJet.pt(), particle.pt(), detJet.pt(), v0.pt(), weight); @@ -2296,7 +2346,7 @@ struct JetFragmentation { } // AntiLambda } template // Reconstructed signal for in-jet V0s: daughters - void fillMatchingV0DauJetHistograms(DetJetType const& detJet, PartJetType const& partJet, V0Type const& v0, ParticleType const& particle, double weight = 1.) + void fillMatchedV0DaughtersInJet(DetJetType const& detJet, PartJetType const& partJet, V0Type const& v0, ParticleType const& particle, double weight = 1.) { auto negTrack = v0.template negTrack_as(); auto posTrack = v0.template posTrack_as(); @@ -2306,8 +2356,8 @@ struct JetFragmentation { registry.fill(HIST("matching/jets/V0/partJetPtDetJetPtPartV0PtNegPtRatioPtRelDiffPt"), partJet.pt(), detJet.pt(), particle.pt(), negPart.pt(), negTrack.pt() / negPart.pt(), (negTrack.pt() - negPart.pt()) / negPart.pt(), weight); } // Misses - Counts (or event weights) - template // Missed inclusive V0s - void fillMatchingMissV0Histograms(T const& pv0, double weight = 1.) + template + void fillMissV0Inclusive(T const& pv0, double weight = 1.) { int pdg = pv0.pdgCode(); registry.fill(HIST("matching/V0/missV0PtEtaPhi"), pv0.pt(), pv0.eta(), pv0.phi(), weight); @@ -2319,8 +2369,17 @@ struct JetFragmentation { registry.fill(HIST("matching/V0/missAntiLambdaPtEtaPhi"), pv0.pt(), pv0.eta(), pv0.phi(), weight); } } - template // Missed V0s in jets - void fillMatchingV0Miss(JetType const& jet, V0Type const& v0, double weight = 1.) + template + void fillMissJet(T const& jet, double weight = 1.) + { + registry.fill(HIST("matching/jets/inclMissJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + if (!jetContainsV0s(jet)) + return; + + registry.fill(HIST("matching/jets/missJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + } + template + void fillMissV0InJet(JetType const& jet, V0Type const& v0, double weight = 1.) { double trackProj = getMomFrac(jet, v0); @@ -2338,8 +2397,8 @@ struct JetFragmentation { } } // Fakes - Counts (or event weights) - template // Fake inclusive V0s - void fillMatchingV0FakeHistograms(T const& coll, U const& v0, double weight = 1.) + template + void fillFakeV0Inclusive(T const& coll, U const& v0, double weight = 1.) { double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; @@ -2399,8 +2458,8 @@ struct JetFragmentation { registry.fill(HIST("matching/V0/fakeAntiLambdaPtDCAdMass"), v0.pt(), v0.dcaV0daughters(), v0.mAntiLambda(), weight); } } - template // Fake inclusive V0s daughters - void fillMatchingFakeV0DauHistograms(U const& v0, double weight = 1.) + template + void fillFakeV0DaughtersInclusive(U const& v0, double weight = 1.) { auto negTrack = v0.template negTrack_as(); auto posTrack = v0.template posTrack_as(); @@ -2427,7 +2486,7 @@ struct JetFragmentation { } } template // Check if inclusive V0 was missed because daughter decayed - void fillMatchingV0DecayedHistograms(V const& v0, double weight = 1.) + void fillFakeV0DecayedInclusive(V const& v0, double weight = 1.) { // Check if decayed daughter auto posTrack = v0.template posTrack_as(); @@ -2481,13 +2540,22 @@ struct JetFragmentation { } } } - template // Fake V0s in jets - void fillMatchingV0JetFakeHistograms(CollisionType const& collision, JetType const& jet, V0Type const& v0, double weight = 1.) + template + void fillFakeJet(T const& jet, double weight = 1.) + { + registry.fill(HIST("matching/jets/fakeJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + if (!jetContainsV0s(jet)) + return; + + registry.fill(HIST("matching/jets/inclFakeJetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi(), weight); + } + template + void fillFakeV0InJet(CollisionType const& coll, JetType const& jet, V0Type const& v0, double weight = 1.) { double trackProj = getMomFrac(jet, v0); - double ctauLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0; - double ctauAntiLambda = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassLambda0Bar; - double ctauK0s = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * constants::physics::MassK0Short; + double ctauLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0; + double ctauAntiLambda = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassLambda0Bar; + double ctauK0s = v0.distovertotmom(coll.posX(), coll.posY(), coll.posZ()) * constants::physics::MassK0Short; double massDiff = v0.mLambda() - v0.mAntiLambda(); double massRatio = v0.mAntiLambda() / v0.mLambda(); double massRelDiff = (v0.mLambda() - v0.mAntiLambda()) / v0.mLambda(); @@ -2566,7 +2634,7 @@ struct JetFragmentation { } } template // Check if V0 in jet was missed because daughter decayed - void fillMatchingV0JetDecayedHistograms(V const& partJet, W const& detJet, X const& v0, double weight = 1.) + void fillFakeV0DecayedInJet(V const& partJet, W const& detJet, X const& v0, double weight = 1.) { // Check if decayed daughter auto posTrack = v0.template posTrack_as(); @@ -2658,40 +2726,30 @@ struct JetFragmentation { void processDummy(aod::JetTracks const&) {} PROCESS_SWITCH(JetFragmentation, processDummy, "Dummy process function turned on by default", true); - // Data - void processDataV0(soa::Filtered::iterator const& jcoll, CandidatesV0DataWithFlags const& V0s) + void processDataV0(soa::Filtered::iterator const& coll, DataV0JetsWithConstituents const& jets, aod::CandidatesV0Data const& V0s, aod::JetTracks const&) { - if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) { + registry.fill(HIST("data/hEvents"), 0.5); + if (!jetderiveddatautilities::selectCollision(coll, eventSelectionBits)) return; - } - registry.fill(HIST("data/V0/nV0sEvent"), V0s.size()); - fillDataV0Histograms(jcoll, V0s); - fillDataV0HistogramsWeighted(jcoll, V0s); - } - PROCESS_SWITCH(JetFragmentation, processDataV0, "Data V0", false); - void processDataV0JetsFrag(soa::Filtered::iterator const& jcoll, DataV0JetsWithConstituents const& v0jets, CandidatesV0DataWithFlags const&, aod::JetTracks const&) - { - if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) { - return; - } - - for (const auto& jet : v0jets) { - if (!jetfindingutilities::isInEtaAcceptance(jet, -99., -99., v0EtaMin, v0EtaMax)) { - continue; - } + registry.fill(HIST("data/hEvents"), 1.5); + registry.fill(HIST("data/V0/nV0sEvent"), V0s.size()); + fillDataV0sInclusive(coll, V0s); + fillDataV0sInclusiveWeighted(coll, V0s); - if (!jetContainsV0s(jet)) // Double check if the jet contains V0s + for (const auto& jet : jets) { + if (!jetfindingutilities::isInEtaAcceptance(jet, -99., -99., v0EtaMin, v0EtaMax)) continue; - fillDataJetHistograms(jet); + fillDataJet(jet); + fillDataV0sInPerpCone(coll, jet, V0s); int nV0inJet = 0, nLambdainJet = 0, nAntiLambdainJet = 0, nK0SinJet = 0; // Counts float wV0inJet = 0, wLambdainJet = 0, wAntiLambdainJet = 0, wK0SinJet = 0; // Weights std::vector values; std::vector> weights; - for (const auto& v0 : jet.candidates_as()) { + for (const auto& v0 : jet.candidates_as()) { if (v0.isRejectedCandidate()) continue; @@ -2711,7 +2769,7 @@ struct JetFragmentation { wAntiLambdainJet += signalProb; } - fillDataV0FragHistograms(jcoll, jet, v0); + fillDataV0sInJet(coll, jet, v0); double z = getMomFrac(jet, v0); std::vector w; @@ -2729,144 +2787,94 @@ struct JetFragmentation { registry.fill(HIST("data/jets/V0/jetPtnV0nK0SnLambdanAntiLambda"), jet.pt(), nV0inJet, nK0SinJet, nLambdainJet, nAntiLambdainJet); registry.fill(HIST("data/jets/weighted/V0/jetPtnV0nK0SnLambdanAntiLambda"), jet.pt(), wV0inJet, wK0SinJet, wLambdainJet, wAntiLambdainJet); + if (nV0inJet == 0) + continue; + int nStates = std::round(std::pow(static_cast(nV0Classes), static_cast(nV0inJet))); for (int M = 0; M < nStates; M++) { std::vector state = convertState(M, nV0inJet, nV0Classes); std::vector corrected; if (doCorrectionWithTracks) - corrected = correctedValuesPlusTracks(state, jet); + corrected = correctedValuesPlusTracks(state, jet); else corrected = correctedValues(state, values); double ws = stateWeight(state, weights); double jetpt = corrected[nV0inJet]; - fillDataJetHistogramsWeighted(jetpt, jet.eta(), jet.phi(), ws); - fillDataV0FragHistogramsWeighted(jcoll, jet, state, corrected, ws); + fillDataJetWeighted(jetpt, jet.eta(), jet.phi(), ws); + fillDataV0sInJetWeighted(coll, jet, state, corrected, ws); } } } - PROCESS_SWITCH(JetFragmentation, processDataV0JetsFrag, "Data V0 jets fragmentation with weights", false); - - void processDataV0PerpCone(soa::Filtered::iterator const& jcoll, aod::V0ChargedJets const& v0jets, CandidatesV0DataWithFlags const& v0s) - { - if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) { - return; - } - - for (const auto& jet : v0jets) { - if (!jetfindingutilities::isInEtaAcceptance(jet, -99., -99., v0EtaMin, v0EtaMax)) { - continue; - } - fillDataPerpConeHists(jcoll, jet, v0s); - } - } - PROCESS_SWITCH(JetFragmentation, processDataV0PerpCone, "Perpendicular cone V0s in data", false); + PROCESS_SWITCH(JetFragmentation, processDataV0, "Data V0", false); - // Matching - void processMcMatchedV0(soa::Filtered::iterator const& jcoll, aod::JetMcCollisions const&, CandidatesV0MCDWithLabelsAndFlags const& v0s, aod::CandidatesV0MCP const& pv0s, aod::JetTracksMCD const& jTracks, aod::JetParticles const&) + void processMcV0(soa::Filtered::iterator const& coll, aod::JetMcCollisions const&, MatchedMCDV0JetsWithConstituents const& v0jetsMCD, MatchedMCPV0JetsWithConstituents const& v0jetsMCP, CandidatesV0MCDWithLabels const& V0s, aod::CandidatesV0MCP const& pV0s, aod::JetTracksMCD const&, aod::JetParticles const&, aod::McParticles const& particles) { - if (!jcoll.has_mcCollision()) { + registry.fill(HIST("matching/hEvents"), 0.5); + if (!coll.has_mcCollision()) return; - } - if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) { + if (!jetderiveddatautilities::selectCollision(coll, eventSelectionBits)) return; - } - double weight = jcoll.mcCollision().weight(); - registry.fill(HIST("matching/V0/nV0sEvent"), v0s.size()); - registry.fill(HIST("matching/V0/nV0sEventWeighted"), v0s.size(), weight); - // TODO: Fill mcd and mcp hists - fillMCDV0Histograms(jcoll, v0s, weight); - fillMCPV0Histograms(pv0s, weight); + registry.fill(HIST("matching/hEvents"), 1.5); + double weight = coll.mcCollision().weight(); + registry.fill(HIST("matching/hEvents"), 2.5, weight); + registry.fill(HIST("matching/V0/nV0sEvent"), V0s.size()); + registry.fill(HIST("matching/V0/nV0sEventWeighted"), V0s.size(), weight); - float nV0s = 0; - for (const auto& v0 : v0s) { - if (v0.isRejectedCandidate()) - continue; - - if (!v0.has_mcParticle()) { - fillMatchingV0FakeHistograms(jcoll, v0, weight); - fillMatchingFakeV0DauHistograms(v0, weight); - fillMatchingV0DecayedHistograms(v0, weight); - continue; - } - for (const auto& pv0 : pv0s) { - if (v0sAreMatched(v0, pv0, jTracks)) { - nV0s += 1; - fillMatchingV0Histograms(jcoll, v0, pv0, weight); - fillMatchingV0DauHistograms(v0, pv0, weight); - } - } - } // Reconstructed V0s - for (const auto& pv0 : pv0s) { - for (const auto& v0 : v0s) { - if (v0sAreMatched(v0, pv0, jTracks)) { - continue; - } - // Fill histograms for missed V0s - } - } - registry.fill(HIST("matching/V0/nV0sEventAcc"), nV0s); - registry.fill(HIST("matching/V0/nV0sEventAccWeighted"), nV0s, weight); - } - PROCESS_SWITCH(JetFragmentation, processMcMatchedV0, "Monte Carlo V0", false); - - void processMcMatchedV0JetsFrag(soa::Filtered::iterator const& jcoll, aod::JetMcCollisions const&, MatchedMCDV0JetsWithConstituents const& v0jetsMCD, MatchedMCPV0JetsWithConstituents const& v0jetsMCP, CandidatesV0MCDWithLabelsAndFlags const&, aod::CandidatesV0MCP const&, aod::JetTracksMCD const& jTracks, aod::JetParticles const&) - { - if (!jcoll.has_mcCollision()) - return; - if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) - return; - - double weight = jcoll.mcCollision().weight(); + fillMcdV0sInclusive(coll, V0s, weight); + fillMcpV0sInclusive(pV0s, weight); + fillMatchingV0sInclusive(coll, V0s, pV0s, weight); for (const auto& detJet : v0jetsMCD) { if (!jetfindingutilities::isInEtaAcceptance(detJet, -99., -99., v0EtaMin, v0EtaMax)) continue; - // Double check if the jet contains V0s - if (!jetContainsV0s(detJet)) - continue; - - fillMCDJetHistograms(detJet, weight); + fillMcdJet(detJet, weight); + fillMcV0sInPerpCone(coll, detJet, V0s, particles, weight); int nV0inJet = 0, nLambdainJet = 0, nAntiLambdainJet = 0, nK0SinJet = 0; if (!detJet.has_matchedJetGeo()) { - for (const auto& detV0 : detJet.candidates_as()) { - fillMatchingV0JetFakeHistograms(jcoll, detJet, detV0, weight); + fillFakeJet(detJet, weight); + for (const auto& v0 : detJet.candidates_as()) { + fillFakeV0InJet(coll, detJet, v0, weight); } continue; } // if jet not matched for (const auto& partJet : detJet.template matchedJetGeo_as()) { - fillMatchingHistogramsJet(detJet, partJet, weight); - for (const auto& detV0 : detJet.candidates_as()) { + fillMatchedJet(detJet, partJet, weight); + fillMcV0sInMatchedPerpCone(coll, detJet, partJet, V0s, particles, weight); + + for (const auto& detV0 : detJet.candidates_as()) { if (!detV0.has_mcParticle()) { - fillMatchingV0JetFakeHistograms(jcoll, detJet, detV0, weight); - fillMatchingV0JetDecayedHistograms(partJet, detJet, detV0, weight); + fillFakeV0InJet(coll, detJet, detV0, weight); + fillFakeV0DecayedInJet(partJet, detJet, detV0, weight); continue; } + bool isV0Matched = false; for (const auto& partV0 : partJet.template candidates_as()) { - if (v0sAreMatched(detV0, partV0, jTracks)) { - isV0Matched = true; - nV0inJet++; - fillMatchingV0FragHistograms(jcoll, detJet, partJet, detV0, partV0, weight); - fillMatchingV0DauJetHistograms(detJet, partJet, detV0, partV0, weight); - - if (std::abs(partV0.pdgCode()) == PDG_t::kK0Short) { - nK0SinJet++; - } else if (partV0.pdgCode() == PDG_t::kLambda0) { - nLambdainJet++; - } else if (partV0.pdgCode() == PDG_t::kLambda0Bar) { - nAntiLambdainJet++; - } - break; - } // if matched + if (!v0sAreMatched(detV0, partV0)) + continue; + + isV0Matched = true; + nV0inJet++; + fillMatchedV0InJet(coll, detJet, partJet, detV0, partV0, weight); + fillMatchedV0DaughtersInJet(detJet, partJet, detV0, partV0, weight); + + if (std::abs(partV0.pdgCode()) == PDG_t::kK0Short) { + nK0SinJet++; + } else if (partV0.pdgCode() == PDG_t::kLambda0) { + nLambdainJet++; + } else if (partV0.pdgCode() == PDG_t::kLambda0Bar) { + nAntiLambdainJet++; + } + break; } // partV0 loop if (!isV0Matched) { - fillMatchingV0JetFakeHistograms(jcoll, detJet, detV0, weight); + fillFakeV0InJet(coll, detJet, detV0, weight); } } // detV0 loop registry.fill(HIST("matching/jets/V0/jetPtnV0MatchednK0SnLambdanAntiLambda"), partJet.pt(), nV0inJet, nK0SinJet, nLambdainJet, nAntiLambdainJet, weight); @@ -2874,70 +2882,46 @@ struct JetFragmentation { } // detJet loop for (const auto& partJet : v0jetsMCP) { - if (!jetContainsV0s(partJet)) { - continue; - } - fillMCPJetHistograms(partJet, weight); + fillMcpJet(partJet, weight); if (!partJet.has_matchedJetGeo()) { + fillMissJet(partJet, weight); for (const auto& partV0 : partJet.candidates_as()) { - fillMatchingV0Miss(partJet, partV0, weight); + fillMissV0InJet(partJet, partV0, weight); } continue; } // if jet not matched bool isJetMatched = false; for (const auto& detJet : partJet.template matchedJetGeo_as()) { - if (!jetfindingutilities::isInEtaAcceptance(detJet, -99., -99., v0EtaMin, v0EtaMax)) { + if (!jetfindingutilities::isInEtaAcceptance(detJet, -99., -99., v0EtaMin, v0EtaMax)) continue; - } + isJetMatched = true; for (const auto& partV0 : partJet.candidates_as()) { bool isV0Matched = false; - for (const auto& detV0 : detJet.candidates_as()) { - if (v0sAreMatched(detV0, partV0, jTracks)) { + for (const auto& detV0 : detJet.candidates_as()) { + if (v0sAreMatched(detV0, partV0)) { isV0Matched = true; break; } } // detV0 loop - if (!isV0Matched) { - fillMatchingV0Miss(partJet, partV0, weight); - } + + // If V0 is matched, it has already been filled in the mcdjet loop + if (!isV0Matched) + fillMissV0InJet(partJet, partV0, weight); } // partV0 loop } // detJet loop // To account for matched jets where the detector level jet is outside of the eta range (cut applied within this task) if (!isJetMatched) { for (const auto& partV0 : partJet.candidates_as()) { - fillMatchingV0Miss(partJet, partV0, weight); + fillMissV0InJet(partJet, partV0, weight); } } } // partJet loop } - PROCESS_SWITCH(JetFragmentation, processMcMatchedV0JetsFrag, "Matched V0 jets fragmentation", false); - - void processMcV0PerpCone(soa::Filtered::iterator const& jcoll, aod::JetMcCollisions const&, MatchedMCDV0Jets const& v0jets, CandidatesV0MCDWithLabelsAndFlags const& v0s, aod::McParticles const& particles, MatchedMCPV0Jets const&) - { - if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) { - return; - } - double weight = jcoll.mcCollision().weight(); - registry.fill(HIST("mcd/V0/nV0sEvent"), v0s.size()); - registry.fill(HIST("mcd/V0/nV0sEventWeighted"), v0s.size(), weight); - - for (const auto& mcdjet : v0jets) { - if (!jetfindingutilities::isInEtaAcceptance(mcdjet, -99., -99., v0EtaMin, v0EtaMax)) - continue; - - fillMcPerpConeHists(jcoll, mcdjet, v0s, particles, weight); - - for (const auto& mcpjet : mcdjet.template matchedJetGeo_as()) { - fillMcPerpConeHists(jcoll, mcdjet, mcpjet, v0s, particles, weight); - break; // Make sure we only do this once - } - } - } - PROCESS_SWITCH(JetFragmentation, processMcV0PerpCone, "Perpendicular cone V0s in MC", false); + PROCESS_SWITCH(JetFragmentation, processMcV0, "MC V0", false); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) diff --git a/PWGJE/Tasks/v0QA.cxx b/PWGJE/Tasks/v0QA.cxx index c66bf164def..31b47a88a6e 100644 --- a/PWGJE/Tasks/v0QA.cxx +++ b/PWGJE/Tasks/v0QA.cxx @@ -45,6 +45,7 @@ using namespace o2::framework; using namespace o2::framework::expressions; // V0 jets +using V0ChargedJetsWithConstituents = soa::Join; using MCDV0Jets = aod::V0ChargedMCDetectorLevelJets; using MCDV0JetsWithConstituents = soa::Join; using MatchedMCDV0Jets = soa::Join; @@ -59,6 +60,10 @@ using MatchedMCPV0JetsWithConstituents = soa::Join; +// Tracks +using DaughterJTracks = soa::Join; +using DaughterTracks = soa::Join; + struct V0QA { HistogramRegistry registry{"registry"}; @@ -300,6 +305,44 @@ struct V0QA { registry.add("tests/sub/JetPtEtaAntiLambdaPt", "Jet Pt, Eta, AntiLambda Pt", HistType::kTH3D, {axisJetPt, axisEta, axisV0Pt}); registry.add("tests/sub/JetPtEtaAntiLambdaZ", "Jet Pt, Eta, AntiLambda Z", HistType::kTH3D, {axisJetPt, axisEta, axisV0Z}); } + if (doprocessTestV0DaughterSharing) { + registry.add("sharing/hEvents", "Events", {HistType::kTH1D, {{2, 0.0f, 2.0f}}}); + registry.add("sharing/V0PtEtaPhi", "V0 Pt, Eta, Phi", HistType::kTH3D, {axisV0Pt, axisEta, axisPhi}); + registry.add("sharing/K0SPtEtaPhi", "K0S Pt, Eta, Phi", HistType::kTH3D, {axisV0Pt, axisEta, axisPhi}); + registry.add("sharing/LambdaPtEtaPhi", "Lambda Pt, Eta, Phi", HistType::kTH3D, {axisV0Pt, axisEta, axisPhi}); + registry.add("sharing/AntiLambdaPtEtaPhi", "AntiLambda Pt, Eta, Phi", HistType::kTH3D, {axisV0Pt, axisEta, axisPhi}); + + registry.add("sharing/V0PtEtaPt", "V0s w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/V0PtEtaPtDaughterPt", "V0s w shared daughter and daughter pt", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/K0SK0S", "K0S-K0S w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/K0SLambda", "K0S-Lambda w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/K0SAntiLambda", "K0S-AntiLambda w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/LambdaK0S", "Lambda-K0S w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/LambdaLambda", "Lambda-Lambda w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/LambdaAntiLambda", "Lambda-AntiLambda w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/AntiLambdaK0S", "AntiLambda-K0S w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/AntiLambdaLambda", "AntiLambda-Lambda w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + registry.add("sharing/AntiLambdaAntiLambda", "AntiLambda-AntiLambda w shared daughter", HistType::kTHnSparseD, {axisV0Pt, axisEta, axisV0Pt}); + + registry.add("sharing/JetPtEtaPhi", "JetPtEtaPhi", HistType::kTH3D, {axisJetPt, axisEta, axisPhi}); + registry.add("sharing/JetPtEtaPhiSingle", "JetPtEtaPhiSingle", HistType::kTH3D, {axisJetPt, axisEta, axisPhi}); + registry.add("sharing/JetPtEtaPhiMultiple", "JetPtEtaPhiMultiple", HistType::kTH3D, {axisJetPt, axisEta, axisPhi}); + registry.add("sharing/JetPtEtaV0Pt", "JetPtEtaV0Pt", HistType::kTH3D, {axisJetPt, axisEta, axisPhi}); + registry.add("sharing/JetPtEtaK0SPt", "JetPtEtaK0SPt", HistType::kTH3D, {axisJetPt, axisEta, axisPhi}); + registry.add("sharing/JetPtEtaLambdaPt", "JetPtEtaLambdaPt", HistType::kTH3D, {axisJetPt, axisEta, axisPhi}); + registry.add("sharing/JetPtEtaAntiLambdaPt", "JetPtEtaAntiLambdaPt", HistType::kTH3D, {axisJetPt, axisEta, axisPhi}); + registry.add("sharing/JetPtEtaV0PtPt", "JetPtEtaV0PtPt", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetPtEtaV0PtPtDaughterPt", "JetPtEtaV0PtPtDaughterPt", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetK0SK0S", "JetK0SK0S", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetK0SLambda", "JetK0SLambda", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetK0SAntiLambda", "JetK0SAntiLambda", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetLambdaK0S", "JetLambdaK0S", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetLambdaLambda", "JetLambdaLambda", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetLambdaAntiLambda", "JetLambdaAntiLambda", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetAntiLambdaK0S", "JetAntiLambdaK0S", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetAntiLambdaLambda", "JetAntiLambdaLambda", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + registry.add("sharing/JetAntiLambdaAntiLambda", "JetAntiLambdaAntiLambda", HistType::kTHnSparseD, {axisJetPt, axisEta, axisV0Pt, axisV0Pt}); + } if (doprocessV0TrackQA) { registry.add("tracks/hEvents", "evts", {HistType::kTH1D, {{2, 0.0f, 2.0f}}}); registry.add("tracks/Pos", "pos", HistType::kTHnSparseD, {axisV0Pt, axisV0Pt, axisEta, axisPhi}); @@ -496,6 +539,16 @@ struct V0QA { auto daughters = particle.daughtersIds(); return ((negId == daughters[0] && posId == daughters[1]) || (posId == daughters[0] && negId == daughters[1])); } + template + bool v0sShareDaughter(U const& trigger, U const& associate) + { + // LOGF(info, "Checking if V0s share daughter"); + auto trigNeg = trigger.template negTrack_as(); + auto trigPos = trigger.template posTrack_as(); + auto assocNeg = associate.template negTrack_as(); + auto assocPos = associate.template posTrack_as(); + return (trigNeg == assocNeg || trigNeg == assocPos || trigPos == assocNeg || trigPos == assocPos); + } template bool genV0PassesEfficiencyCuts(T const& pv0) { @@ -784,6 +837,138 @@ struct V0QA { } } + template + void fillV0DaughterSharingInclusive(T const& v0) + { + registry.fill(HIST("sharing/V0PtEtaPhi"), v0.pt(), v0.eta(), v0.phi()); + if (v0.isK0SCandidate()) + registry.fill(HIST("sharing/K0SPtEtaPhi"), v0.pt(), v0.eta(), v0.phi()); + if (v0.isLambdaCandidate()) + registry.fill(HIST("sharing/LambdaPtEtaPhi"), v0.pt(), v0.eta(), v0.phi()); + if (v0.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/AntiLambdaPtEtaPhi"), v0.pt(), v0.eta(), v0.phi()); + } + + template + void fillV0DaughterSharingInclusive(V const& trigger, V const& associate) + { + double weight = 0.5; // To correct for double-counting + double pthard, etahard, ptsoft; + if (trigger.pt() > associate.pt()) { + pthard = trigger.pt(); + etahard = trigger.eta(); + ptsoft = associate.pt(); + } else { + pthard = associate.pt(); + etahard = associate.eta(); + ptsoft = trigger.pt(); + } + + registry.fill(HIST("sharing/V0PtEtaPt"), pthard, etahard, ptsoft, weight); + // Tried to get this in a function, but couldn't make it work + auto trigNeg = trigger.template negTrack_as().template track_as(); + auto trigPos = trigger.template posTrack_as().template track_as(); + auto assocNeg = associate.template negTrack_as().template track_as(); + auto assocPos = associate.template posTrack_as().template track_as(); + double sharedDaughterPt; + if (trigNeg == assocNeg || trigNeg == assocPos) + sharedDaughterPt = trigNeg.pt(); + else + sharedDaughterPt = trigPos.pt(); + registry.fill(HIST("sharing/V0PtEtaPtDaughterPt"), pthard, etahard, ptsoft, sharedDaughterPt, weight); + + if (trigger.isK0SCandidate() && associate.isK0SCandidate()) + registry.fill(HIST("sharing/K0SK0S"), pthard, etahard, ptsoft, weight); + if (trigger.isK0SCandidate() && associate.isLambdaCandidate()) + registry.fill(HIST("sharing/K0SLambda"), pthard, etahard, ptsoft, weight); + if (trigger.isK0SCandidate() && associate.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/K0SAntiLambda"), pthard, etahard, ptsoft, weight); + + if (trigger.isLambdaCandidate() && associate.isK0SCandidate()) + registry.fill(HIST("sharing/LambdaK0S"), pthard, etahard, ptsoft, weight); + if (trigger.isLambdaCandidate() && associate.isLambdaCandidate()) + registry.fill(HIST("sharing/LambdaLambda"), pthard, etahard, ptsoft, weight); + if (trigger.isLambdaCandidate() && associate.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/LambdaAntiLambda"), pthard, etahard, ptsoft, weight); + + if (trigger.isAntiLambdaCandidate() && associate.isK0SCandidate()) + registry.fill(HIST("sharing/AntiLambdaK0S"), pthard, etahard, ptsoft, weight); + if (trigger.isAntiLambdaCandidate() && associate.isLambdaCandidate()) + registry.fill(HIST("sharing/AntiLambdaLambda"), pthard, etahard, ptsoft, weight); + if (trigger.isAntiLambdaCandidate() && associate.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/AntiLambdaAntiLambda"), pthard, etahard, ptsoft, weight); + } + + template + void fillV0DaughterSharingJet(T const& jet, bool jetHasSingleV0) + { + registry.fill(HIST("sharing/JetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi()); + + if (jetHasSingleV0) + registry.fill(HIST("sharing/JetPtEtaPhiSingle"), jet.pt(), jet.eta(), jet.phi()); + else + registry.fill(HIST("sharing/JetPtEtaPhiMultiple"), jet.pt(), jet.eta(), jet.phi()); + } + + template + void fillV0DaughterSharingJet(T const& jet, U const& v0) + { + registry.fill(HIST("sharing/JetPtEtaV0Pt"), jet.pt(), jet.eta(), v0.pt()); + if (v0.isK0SCandidate()) + registry.fill(HIST("sharing/JetPtEtaK0SPt"), jet.pt(), jet.eta(), v0.pt()); + if (v0.isLambdaCandidate()) + registry.fill(HIST("sharing/JetPtEtaLambdaPt"), jet.pt(), jet.eta(), v0.pt()); + if (v0.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/JetPtEtaAntiLambdaPt"), jet.pt(), jet.eta(), v0.pt()); + } + + template + void fillV0DaughterSharingJet(V const& jet, W const& trigger, W const& associate) + { + double weight = 0.5; // To correct for double-counting + double pthard, ptsoft; + if (trigger.pt() > associate.pt()) { + pthard = trigger.pt(); + ptsoft = associate.pt(); + } else { + pthard = associate.pt(); + ptsoft = trigger.pt(); + } + + registry.fill(HIST("sharing/JetPtEtaV0PtPt"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + auto trigNeg = trigger.template negTrack_as().template track_as(); + auto trigPos = trigger.template posTrack_as().template track_as(); + auto assocNeg = associate.template negTrack_as().template track_as(); + auto assocPos = associate.template posTrack_as().template track_as(); + double sharedDaughterPt; + if (trigNeg == assocNeg || trigNeg == assocPos) + sharedDaughterPt = trigNeg.pt(); + else + sharedDaughterPt = trigPos.pt(); + registry.fill(HIST("sharing/JetPtEtaV0PtPtDaughterPt"), jet.pt(), jet.eta(), pthard, ptsoft, sharedDaughterPt, weight); + + if (trigger.isK0SCandidate() && associate.isK0SCandidate()) + registry.fill(HIST("sharing/JetK0SK0S"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + if (trigger.isK0SCandidate() && associate.isLambdaCandidate()) + registry.fill(HIST("sharing/JetK0SLambda"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + if (trigger.isK0SCandidate() && associate.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/JetK0SAntiLambda"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + + if (trigger.isLambdaCandidate() && associate.isK0SCandidate()) + registry.fill(HIST("sharing/JetLambdaK0S"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + if (trigger.isLambdaCandidate() && associate.isLambdaCandidate()) + registry.fill(HIST("sharing/JetLambdaLambda"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + if (trigger.isLambdaCandidate() && associate.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/JetLambdaAntiLambda"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + + if (trigger.isAntiLambdaCandidate() && associate.isK0SCandidate()) + registry.fill(HIST("sharing/JetAntiLambdaK0S"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + if (trigger.isAntiLambdaCandidate() && associate.isLambdaCandidate()) + registry.fill(HIST("sharing/JetAntiLambdaLambda"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + if (trigger.isAntiLambdaCandidate() && associate.isAntiLambdaCandidate()) + registry.fill(HIST("sharing/JetAntiLambdaAntiLambda"), jet.pt(), jet.eta(), pthard, ptsoft, weight); + } + template void fillTrackQa(V const& v0) { @@ -1477,7 +1662,7 @@ struct V0QA { PROCESS_SWITCH(V0QA, processFeeddownMatchedJets, "Jets feeddown", false); // Test the difference between excluding V0s from jet finding and subtracting V0s from jets afterwards - void processTestWeightedJetFinder(soa::Filtered::iterator const& jcoll, soa::Join const& jets, aod::CandidatesV0Data const& v0s) + void processTestWeightedJetFinder(soa::Filtered::iterator const& jcoll, V0ChargedJetsWithConstituents const& jets, aod::CandidatesV0Data const& v0s) { registry.fill(HIST("tests/weighted/hEvents"), 0.5); if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) @@ -1505,7 +1690,7 @@ struct V0QA { } PROCESS_SWITCH(V0QA, processTestWeightedJetFinder, "Test weighted jet finder", false); - void processTestSubtractedJetFinder(soa::Filtered::iterator const& jcoll, soa::Join const& jets, aod::CandidatesV0Data const& v0s) + void processTestSubtractedJetFinder(soa::Filtered::iterator const& jcoll, V0ChargedJetsWithConstituents const& jets, aod::CandidatesV0Data const& v0s) { registry.fill(HIST("tests/hEvents"), 0.5); if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) @@ -1552,8 +1737,63 @@ struct V0QA { } PROCESS_SWITCH(V0QA, processTestSubtractedJetFinder, "Test subtracted jet finder", false); - using DaughterJTracks = soa::Join; - using DaughterTracks = soa::Join; + void processTestV0DaughterSharing(soa::Filtered::iterator const& jcoll, V0ChargedJetsWithConstituents const& jets, aod::CandidatesV0Data const& v0s, DaughterJTracks const&, DaughterTracks const&) + { + registry.fill(HIST("sharing/hEvents"), 0.5); + if (!jetderiveddatautilities::selectCollision(jcoll, eventSelectionBits)) + return; + registry.fill(HIST("sharing/hEvents"), 1.5); + + // Check if V0s within the same event share daughters + for (const auto& trigger : v0s) { + if (trigger.isRejectedCandidate()) + continue; + + fillV0DaughterSharingInclusive(trigger); + + for (const auto& associate : v0s) { + if (associate.isRejectedCandidate()) + continue; + + if (trigger == associate) + continue; + + // Double-counting accounted for by filling histograms with weight 0.5 + if (v0sShareDaughter(trigger, associate)) { + fillV0DaughterSharingInclusive(trigger, associate); + } + } + } + + // Check if V0s within the same jet share daughters + for (const auto& jet : jets) { + bool jetHasSingleV0 = true; + + for (const auto& trigger : jet.template candidates_as()) { + if (trigger.isRejectedCandidate()) + continue; + + fillV0DaughterSharingJet(jet, trigger); + + for (const auto& associate : jet.template candidates_as()) { + if (associate.isRejectedCandidate()) + continue; + + if (trigger == associate) + continue; + + jetHasSingleV0 = false; + + // Double-counting accounted for by filling histograms with weight 0.5 + if (v0sShareDaughter(trigger, associate)) + fillV0DaughterSharingJet(jet, trigger, associate); + } + } + fillV0DaughterSharingJet(jet, jetHasSingleV0); + } + } + PROCESS_SWITCH(V0QA, processTestV0DaughterSharing, "Test V0s with shared daughters", false); + void processV0TrackQA(aod::JetCollision const& jcoll, aod::CandidatesV0Data const& v0s, DaughterJTracks const&, DaughterTracks const&) { registry.fill(HIST("tracks/hEvents"), 0.5);