Skip to content

Commit

Permalink
PWGJE: add weighted track and number of collisions option (#4179)
Browse files Browse the repository at this point in the history 
* add weighted track and number of collisions options to QA tasks

* fix exponent in weighting
  • Loading branch information
@jaimenorman
jaimenorman authored Dec 24, 2023
1 parent b7ea22e commit 1b8043c
Show file tree
Hide file tree
Showing 3 changed files with 120 additions and 49 deletions.
52 changes: 38 additions & 14 deletions PWGJE/Tasks/jetfinderQA.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,7 @@ struct JetFinderQATask {
Configurable<std::string> trackSelections{"trackSelections", "globalTracks", "set track selections"};
Configurable<float> pTHatMaxMCD{"pTHatMaxMCD", 999.0, "maximum fraction of hard scattering for jet acceptance in detector MC"};
Configurable<float> pTHatMaxMCP{"pTHatMaxMCP", 999.0, "maximum fraction of hard scattering for jet acceptance in particle MC"};
Configurable<float> pTHatExponent{"pTHatExponent", 0.1666, "exponent of the event weight for the calculation of pTHat"};
Configurable<float> pTHatExponent{"pTHatExponent", 6.0, "exponent of the event weight for the calculation of pTHat"};

std::vector<bool> filledJetR;
std::vector<double> jetRadiiValues;
Expand Down Expand Up @@ -158,11 +158,15 @@ struct JetFinderQATask {
registry.add("h3_jet_r_jet_pt_track_eta_Triggered_Both", "#it{R}_{jet};#it{p}_{T,jet} (GeV/#it{c});#eta_{jet tracks}", {HistType::kTH3F, {{jetRadiiBins, ""}, {200, 0., 200.}, {100, -1.0, 1.0}}});
registry.add("h3_jet_r_jet_pt_track_phi_Triggered_Both", "#it{R}_{jet};#it{p}_{T,jet} (GeV/#it{c});#varphi_{jet tracks}", {HistType::kTH3F, {{jetRadiiBins, ""}, {200, 0., 200.}, {160, -1.0, 7.}}});
}
if (doprocessTracks) {

if (doprocessTracks || doprocessTracksWeighted) {
registry.add("h_collisions", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
registry.add("h_track_pt", "track pT;#it{p}_{T,track} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0., 200.}}});
registry.add("h_track_eta", "track #eta;#eta_{track};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}});
registry.add("h_track_phi", "track #varphi;#varphi_{track};entries", {HistType::kTH1F, {{160, -1.0, 7.}}});
if (doprocessTracksWeighted) {
registry.add("h_collisions_weighted", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
}
}

if (doprocessMCCollisionsWeighted) {
Expand All @@ -180,7 +184,7 @@ struct JetFinderQATask {
void fillHistograms(T const& jet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (jet.pt() > pTHatMaxMCD * pTHat) {
return;
}
Expand Down Expand Up @@ -210,7 +214,7 @@ struct JetFinderQATask {
void fillMCPHistograms(T const& jet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (jet.pt() > pTHatMaxMCP * pTHat) {
return;
}
Expand Down Expand Up @@ -238,7 +242,7 @@ struct JetFinderQATask {
template <typename T, typename U>
void fillMCMatchedHistograms(T const& mcdjet, float weight = 1.0)
{
float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (mcdjet.pt() > pTHatMaxMCD * pTHat) {
return;
}
Expand All @@ -265,6 +269,18 @@ struct JetFinderQATask {
}
}

void fillTrackHistograms(soa::Filtered<JetTracks> const& tracks, float weight = 1.0)
{
for (auto const& track : tracks) {
if (!JetDerivedDataUtilities::selectTrack(track, trackSelection)) {
continue;
}
registry.fill(HIST("h_track_pt"), track.pt(), weight);
registry.fill(HIST("h_track_eta"), track.eta(), weight);
registry.fill(HIST("h_track_phi"), track.phi(), weight);
}
}

void processJetsData(soa::Join<aod::ChargedJets, aod::ChargedJetConstituents>::iterator const& jet, JetTracks const& tracks)
{
fillHistograms(jet);
Expand Down Expand Up @@ -452,22 +468,30 @@ struct JetFinderQATask {
void processTracks(aod::JCollision const& collision,
soa::Filtered<JetTracks> const& tracks)
{

registry.fill(HIST("h_collisions"), 0.5);
if (!JetDerivedDataUtilities::selectCollision(collision, eventSelection)) {
return;
}
registry.fill(HIST("h_collisions"), 1.5);
for (auto const& track : tracks) {
if (!JetDerivedDataUtilities::selectTrack(track, trackSelection)) {
continue;
}
registry.fill(HIST("h_track_pt"), track.pt());
registry.fill(HIST("h_track_eta"), track.eta());
registry.fill(HIST("h_track_phi"), track.phi());
}
fillTrackHistograms(tracks);
}
PROCESS_SWITCH(JetFinderQATask, processTracks, "QA for charged tracks", false);

void processTracksWeighted(soa::Join<aod::JCollisions, aod::JMcCollisionLbs>::iterator const& collision,
aod::JMcCollisions const& mcCollisions,
soa::Filtered<JetTracks> const& tracks)
{
float eventWeight = collision.mcCollision().weight();
registry.fill(HIST("h_collisions"), 0.5);
registry.fill(HIST("h_collisions_weighted"), 0.5, eventWeight);
if (!JetDerivedDataUtilities::selectCollision(collision, eventSelection)) {
return;
}
registry.fill(HIST("h_collisions"), 1.5);
registry.fill(HIST("h_collisions_weighted"), 1.5, eventWeight);
fillTrackHistograms(tracks, eventWeight);
}
PROCESS_SWITCH(JetFinderQATask, processTracksWeighted, "QA for charged tracks weighted", false);
};

WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { return WorkflowSpec{adaptAnalysisTask<JetFinderQATask>(cfgc, TaskName{"jet-finder-charged-qa"})}; }
66 changes: 45 additions & 21 deletions PWGJE/Tasks/jetfinderfullQA.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ struct JetFinderFullQATask {

Configurable<float> pTHatMaxMCD{"pTHatMaxMCD", 999.0, "maximum fraction of hard scattering for jet acceptance in detector MC"};
Configurable<float> pTHatMaxMCP{"pTHatMaxMCP", 999.0, "maximum fraction of hard scattering for jet acceptance in particle MC"};
Configurable<float> pTHatExponent{"pTHatExponent", 0.1666, "exponent of the event weight for the calculation of pTHat"};
Configurable<float> pTHatExponent{"pTHatExponent", 6.0, "exponent of the event weight for the calculation of pTHat"};

std::vector<bool> filledJetR;
std::vector<double> jetRadiiValues;
Expand Down Expand Up @@ -138,7 +138,7 @@ struct JetFinderFullQATask {
registry.add("h3_jet_pt_part_jet_ntracks_part_jet_ntracks", ";#it{p}_{T,jet}^{part} (GeV/#it{c}); N_{jet tracks}^{part}; N_{jet tracks}", {HistType::kTH3F, {{200, 0.0, 200}, {100, -0.5, 99.5}, {100, -0.5, 99.5}}});
}

if (doprocessTracks) {
if (doprocessTracks || doprocessTracksWeighted) {
registry.add("h_collisions", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
registry.add("h_track_pt", "track pT;#it{p}_{T,track} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0., 200.}}});
registry.add("h_track_eta", "track #eta;#eta_{track};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}});
Expand All @@ -147,6 +147,9 @@ struct JetFinderFullQATask {
registry.add("h_cluster_eta", "cluster #eta;#eta_{cluster};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}});
registry.add("h_cluster_phi", "cluster #varphi;#varphi_{cluster};entries", {HistType::kTH1F, {{160, -1.0, 7.}}});
registry.add("h_cluster_energy", "cluster E;E_{cluster} (GeV);entries", {HistType::kTH1F, {{200, 0., 200.}}});
if (doprocessTracksWeighted) {
registry.add("h_collisions_weighted", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
}
}

if (doprocessMCCollisionsWeighted) {
Expand Down Expand Up @@ -176,7 +179,7 @@ struct JetFinderFullQATask {
void fillHistograms(T const& jet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (jet.pt() > pTHatMaxMCD * pTHat) {
return;
}
Expand Down Expand Up @@ -216,7 +219,7 @@ struct JetFinderFullQATask {
void fillMCPHistograms(T const& jet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (jet.pt() > pTHatMaxMCP * pTHat) {
return;
}
Expand Down Expand Up @@ -245,7 +248,7 @@ struct JetFinderFullQATask {
void fillMCMatchedHistograms(T const& mcdjet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (mcdjet.pt() > pTHatMaxMCD * pTHat) {
return;
}
Expand All @@ -272,6 +275,25 @@ struct JetFinderFullQATask {
}
}

void fillTrackHistograms(soa::Filtered<JetTracks> const& tracks, soa::Filtered<JetClusters> const& clusters, float weight = 1.0)
{
for (auto const& track : tracks) {
if (!JetDerivedDataUtilities::selectTrack(track, trackSelection) || !JetDerivedDataUtilities::applyTrackKinematics(track, trackPtMin, trackPtMax, trackEtaMin, trackEtaMax)) {
continue;
}
registry.fill(HIST("h_track_pt"), track.pt(), weight);
registry.fill(HIST("h_track_eta"), track.eta(), weight);
registry.fill(HIST("h_track_phi"), track.phi(), weight);
}
for (auto const& cluster : clusters) {
double clusterpt = cluster.energy() / std::cosh(cluster.eta());
registry.fill(HIST("h_cluster_pt"), clusterpt, weight);
registry.fill(HIST("h_cluster_eta"), cluster.eta(), weight);
registry.fill(HIST("h_cluster_phi"), cluster.phi(), weight);
registry.fill(HIST("h_cluster_energy"), cluster.energy(), weight);
}
}

void processDummy(aod::JCollisions const& collision)
{
}
Expand Down Expand Up @@ -348,29 +370,31 @@ struct JetFinderFullQATask {
soa::Filtered<JetTracks> const& tracks,
soa::Filtered<JetClusters> const& clusters)
{

registry.fill(HIST("h_collisions"), 0.5);
if (!JetDerivedDataUtilities::eventEMCAL(collision)) {
return;
}
registry.fill(HIST("h_collisions"), 1.5);
for (auto const& track : tracks) {
if (!JetDerivedDataUtilities::selectTrack(track, trackSelection) || !JetDerivedDataUtilities::applyTrackKinematics(track, trackPtMin, trackPtMax, trackEtaMin, trackEtaMax)) {
continue;
}
registry.fill(HIST("h_track_pt"), track.pt());
registry.fill(HIST("h_track_eta"), track.eta());
registry.fill(HIST("h_track_phi"), track.phi());
}
for (auto const& cluster : clusters) {
double clusterpt = cluster.energy() / std::cosh(cluster.eta());
registry.fill(HIST("h_cluster_pt"), clusterpt);
registry.fill(HIST("h_cluster_eta"), cluster.eta());
registry.fill(HIST("h_cluster_phi"), cluster.phi());
registry.fill(HIST("h_cluster_energy"), cluster.energy());
}
fillTrackHistograms(tracks, clusters);
}
PROCESS_SWITCH(JetFinderFullQATask, processTracks, "QA for charged tracks", false);

void processTracksWeighted(soa::Join<aod::JCollisions, aod::JMcCollisionLbs>::iterator const& collision,
aod::JMcCollisions const& mcCollisions,
soa::Filtered<JetTracks> const& tracks,
soa::Filtered<JetClusters> const& clusters)
{
float eventWeight = collision.mcCollision().weight();
registry.fill(HIST("h_collisions"), 0.5);
registry.fill(HIST("h_collisions_weighted"), 0.5, eventWeight);
if (!JetDerivedDataUtilities::eventEMCAL(collision)) {
return;
}
registry.fill(HIST("h_collisions"), 1.5);
registry.fill(HIST("h_collisions_weighted"), 1.5, eventWeight);
fillTrackHistograms(tracks, clusters, eventWeight);
}
PROCESS_SWITCH(JetFinderFullQATask, processTracksWeighted, "QA for charged tracks weighted", false);
};

using JetFinderFullJetsQATask = JetFinderFullQATask<aod::FullJets, aod::FullJetConstituents, aod::FullMCDetectorLevelJets, aod::FullMCDetectorLevelJetConstituents, aod::FullMCDetectorLevelJetsMatchedToFullMCParticleLevelJets, aod::FullMCDetectorLevelJetEventWeights, aod::FullMCParticleLevelJets, aod::FullMCParticleLevelJetConstituents, aod::FullMCParticleLevelJetsMatchedToFullMCDetectorLevelJets, aod::FullMCParticleLevelJetEventWeights>;
Expand Down
51 changes: 37 additions & 14 deletions PWGJE/Tasks/jetfinderhfQA.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -60,7 +60,7 @@ struct JetFinderHFQATask {
Configurable<int> selectionFlagBplus{"selectionFlagBplus", 1, "Selection Flag for B+"};
Configurable<float> pTHatMaxMCD{"pTHatMaxMCD", 999.0, "maximum fraction of hard scattering for jet acceptance in detector MC"};
Configurable<float> pTHatMaxMCP{"pTHatMaxMCP", 999.0, "maximum fraction of hard scattering for jet acceptance in particle MC"};
Configurable<float> pTHatExponent{"pTHatExponent", 0.1666, "exponent of the event weight for the calculation of pTHat"};
Configurable<float> pTHatExponent{"pTHatExponent", 6.0, "exponent of the event weight for the calculation of pTHat"};

HfHelper hfHelper;
std::vector<bool> filledJetR;
Expand Down Expand Up @@ -204,11 +204,14 @@ struct JetFinderHFQATask {
registry.add("h3_jet_r_jet_pt_candidate_y_Triggered_Both", "#it{R}_{jet};#it{p}_{T,jet} (GeV/#it{c});y_{candidate}", {HistType::kTH3F, {{jetRadiiBins, ""}, {200, 0., 200.}, {100, -1.0, 1.0}}});
}

if (doprocessTracks) {
if (doprocessTracks || doprocessTracksWeighted) {
registry.add("h_collisions", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
registry.add("h_track_pt", "track pT;#it{p}_{T,track} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0., 200.}}});
registry.add("h_track_eta", "track #eta;#eta_{track};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}});
registry.add("h_track_phi", "track #varphi;#varphi_{track};entries", {HistType::kTH1F, {{160, -1.0, 7.}}});
if (doprocessTracksWeighted) {
registry.add("h_collisions_weighted", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}});
}
}

if (doprocessMCCollisionsWeighted) {
Expand All @@ -234,7 +237,7 @@ struct JetFinderHFQATask {
void fillHistograms(T const& jet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (jet.pt() > pTHatMaxMCD * pTHat) {
return;
}
Expand Down Expand Up @@ -303,7 +306,7 @@ struct JetFinderHFQATask {
void fillMCPHistograms(T const& jet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (jet.pt() > pTHatMaxMCP * pTHat) {
return;
}
Expand Down Expand Up @@ -344,7 +347,7 @@ struct JetFinderHFQATask {
void fillMCMatchedHistograms(T const& mcdjet, float weight = 1.0)
{

float pTHat = 10. / (std::pow(weight, pTHatExponent));
float pTHat = 10. / (std::pow(weight, 1.0 / pTHatExponent));
if (mcdjet.pt() > pTHatMaxMCD * pTHat) {
return;
}
Expand Down Expand Up @@ -399,6 +402,18 @@ struct JetFinderHFQATask {
}
}

void fillTrackHistograms(soa::Filtered<JetTracks> const& tracks, float weight = 1.0)
{
for (auto const& track : tracks) {
if (!JetDerivedDataUtilities::selectTrack(track, trackSelection)) {
continue;
}
registry.fill(HIST("h_track_pt"), track.pt(), weight);
registry.fill(HIST("h_track_eta"), track.eta(), weight);
registry.fill(HIST("h_track_phi"), track.phi(), weight);
}
}

void processDummy(aod::Collisions const& collision)
{
}
Expand Down Expand Up @@ -611,22 +626,30 @@ struct JetFinderHFQATask {
void processTracks(aod::JCollision const& collision,
soa::Filtered<JetTracks> const& tracks)
{

registry.fill(HIST("h_collisions"), 0.5);
if (!JetDerivedDataUtilities::selectCollision(collision, eventSelection)) {
return;
}
registry.fill(HIST("h_collisions"), 1.5);
for (auto const& track : tracks) {
if (!JetDerivedDataUtilities::selectTrack(track, trackSelection)) {
continue;
}
registry.fill(HIST("h_track_pt"), track.pt());
registry.fill(HIST("h_track_eta"), track.eta());
registry.fill(HIST("h_track_phi"), track.phi());
}
fillTrackHistograms(tracks);
}
PROCESS_SWITCH(JetFinderHFQATask, processTracks, "QA for charged tracks", false);

void processTracksWeighted(soa::Join<aod::JCollisions, aod::JMcCollisionLbs>::iterator const& collision,
aod::JMcCollisions const& mcCollisions,
soa::Filtered<JetTracks> const& tracks)
{
float eventWeight = collision.mcCollision().weight();
registry.fill(HIST("h_collisions"), 0.5);
registry.fill(HIST("h_collisions_weighted"), 0.5, eventWeight);
if (!JetDerivedDataUtilities::selectCollision(collision, eventSelection)) {
return;
}
registry.fill(HIST("h_collisions"), 1.5);
registry.fill(HIST("h_collisions_weighted"), 1.5, eventWeight);
fillTrackHistograms(tracks, eventWeight);
}
PROCESS_SWITCH(JetFinderHFQATask, processTracksWeighted, "QA for charged tracks weighted", false);
};

using JetFinderD0QATask = JetFinderHFQATask<aod::D0ChargedJets, aod::D0ChargedJetConstituents, soa::Join<aod::HfCand2Prong, aod::HfSelD0>, aod::D0ChargedMCDetectorLevelJets, aod::D0ChargedMCDetectorLevelJetConstituents, aod::D0ChargedMCDetectorLevelJetsMatchedToD0ChargedMCParticleLevelJets, aod::D0ChargedMCDetectorLevelJetEventWeights, soa::Join<aod::HfCand2Prong, aod::HfSelD0, aod::HfCand2ProngMcRec>, aod::D0ChargedMCParticleLevelJets, aod::D0ChargedMCParticleLevelJetConstituents, aod::D0ChargedMCParticleLevelJetsMatchedToD0ChargedMCDetectorLevelJets, aod::D0ChargedMCParticleLevelJetEventWeights, soa::Join<aod::McParticles, aod::HfCand2ProngMcGen>>;
Expand Down

0 comments on commit 1b8043c

Please sign in to comment.