Skip to content

Commit

Permalink
Final functionality required for hyperloop usage
Browse files Browse the repository at this point in the history
  • Loading branch information
@smaff92
smaff92 committed May 9, 2024
1 parent 7fbfcab commit da169da
Showing 1 changed file with 63 additions and 31 deletions.
94 changes: 63 additions & 31 deletions PWGJE/Tasks/phiInJets.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ struct phiInJets {
HistogramRegistry JEhistos{"JEhistos", {}, OutputObjHandlingPolicy::AnalysisObject};

HistogramRegistry registry{"registry",
{{"h_jet_pt", "jet pT;#it{p}_{T,jet} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0., 200.}}}},
{{"h_jet_pt", "jet pT;#it{p}_{T,jet} (GeV/#it{c});entries", {HistType::kTH1F, {{4000, 0., 200.}}}},
{"h_jet_eta", "jet #eta;#eta_{jet};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}}},
{"h_jet_phi", "jet #phi;#phi_{jet};entries", {HistType::kTH1F, {{80, -1.0, 7.}}}},
{"h_matched_REC_jet_pt", "matched_REC level jet pT;#it{p}_{T,jet part} (GeV/#it{c});Delta", {HistType::kTH2F, {{200, 0., 200.}, {400, -20., 20.}}}},
Expand All @@ -59,14 +59,15 @@ struct phiInJets {
{"h_matched_GEN_jet_pt", "matched_GEN level jet pT;#it{p}_{T,jet part} (GeV/#it{c});Delta", {HistType::kTH2F, {{200, 0., 200.}, {400, -20., 20.}}}},
{"h_matched_GEN_jet_eta", "matched_GEN level jet #eta;#eta_{jet part};Delta", {HistType::kTH2F, {{100, -1.0, 1.0}, {400, -20., 20.}}}},
{"h_matched_GEN_jet_phi", "matched_GEN level jet #phi;#phi_{jet part};Delta", {HistType::kTH2F, {{80, -1.0, 7.}, {400, -20., 20.}}}},
{"h_part_jet_pt", "particle level jet pT;#it{p}_{T,jet part} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0., 200.}}}},
{"h_part_jet_pt", "particle level jet pT;#it{p}_{T,jet part} (GeV/#it{c});entries", {HistType::kTH1F, {{4000, 0., 200.}}}},
{"h_part_jet_eta", "particle level jet #eta;#eta_{jet part};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}}},
{"h_part_jet_phi", "particle level jet #phi;#phi_{jet part};entries", {HistType::kTH1F, {{80, -1.0, 7.}}}}}};

Configurable<std::string> cfgeventSelections{"cfgeventSelections", "sel8", "choose event selection"};
Configurable<std::string> cfgtrackSelections{"cfgtrackSelections", "globalTracks", "set track selections"};

Configurable<double> cfgtrkMinPt{"cfgtrkMinPt", 0.15, "set track min pT"};
Configurable<double> cfgtrkMaxEta{"cfgtrkMaxEta", 0.9, "set track max Eta"};
Configurable<double> cfgMaxDCArToPVcut{"cfgMaxDCArToPVcut", 0.5, "Track DCAr cut to PV Maximum"};
Configurable<double> cfgMaxDCAzToPVcut{"cfgMaxDCAzToPVcut", 2.0, "Track DCAz cut to PV Maximum"};
Configurable<bool> cfgPrimaryTrack{"cfgPrimaryTrack", true, "Primary track selection"}; // kGoldenChi2 | kDCAxy | kDCAz
Expand Down Expand Up @@ -131,7 +132,10 @@ struct phiInJets {
JEhistos.add("nEvents_MCGen", "nEvents_MCGen", kTH1F, {{4, 0.0, 4.0}});
JEhistos.add("nEvents_MCRec_MATCHED", "nEvents_MCRec_MATCHED", kTH1F, {{4, 0.0, 4.0}});
JEhistos.add("nEvents_MCGen_MATCHED", "nEvents_MCGen_MATCHED", kTH1F, {{4, 0.0, 4.0}});


JEhistos.add("JetVsPhi_GEN","JetVsPhi_GEN",kTH2F,{{4000, 0., 200.},{200, 0, 20.0}});
JEhistos.add("JetVsPhi_REC","JetVsPhi_REC",kTH2F,{{4000, 0., 200.},{200, 0, 20.0}});

JEhistos.add("nJetsPerEvent", "nJetsPerEvent", kTH1F, {{10, 0.0, 10.0}});

JEhistos.add("hDCArToPv", "DCArToPv", kTH1F, {{300, 0.0, 3.0}});
Expand Down Expand Up @@ -234,6 +238,9 @@ struct phiInJets {
if (track.pt() < cfgtrkMinPt)
return false;

if (std::abs(track.eta()) > cfgtrkMaxEta)
return false;

if (std::abs(track.dcaXY()) > cfgMaxDCArToPVcut)
return false;

Expand Down Expand Up @@ -306,7 +313,7 @@ struct phiInJets {
lDecayDaughter2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massKa);
lResonance = lDecayDaughter1 + lDecayDaughter2;

if (std::abs(lResonance.Rapidity()) > 0.5)
if (std::abs(lResonance.Eta()) > cfgtrkMaxEta)
return;

/////////////////////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -455,7 +462,7 @@ struct phiInJets {

if (!trackSelection(originalTrack))
continue;

JEhistos.fill(HIST("etaHistogram"), trackC.eta());
JEhistos.fill(HIST("phiHistogram"), trackC.phi());
} // JTrack Loop
Expand Down Expand Up @@ -489,7 +496,7 @@ struct phiInJets {

bool INELgt0 = false;
for (const auto& track : tracks) {
if (TMath::Abs(track.eta()) < 0.8) {
if (fabs(track.eta()) < cfgtrkMaxEta) {
INELgt0 = true;
break;
}
Expand Down Expand Up @@ -523,7 +530,7 @@ struct phiInJets {
continue;
if (track.has_mcParticle()) {
auto mcParticle = track.mcParticle();
if (mcParticle.isPhysicalPrimary() && fabs(mcParticle.y()) <= 0.5) {
if (mcParticle.isPhysicalPrimary() && fabs(mcParticle.eta()) <= cfgtrkMaxEta) {
if (abs(mcParticle.pdgCode()) == 211)
JEhistos.fill(HIST("ptJEHistogramPion"), mcParticle.pt());
if (abs(mcParticle.pdgCode()) == 321)
Expand All @@ -540,7 +547,7 @@ struct phiInJets {
if (originalTrack.index() >= originalTrack2.index())
continue;

if (fabs(originalTrack.eta()) > 0.8 || fabs(originalTrack2.eta()) > 0.8)
if (fabs(originalTrack.eta()) > cfgtrkMaxEta || fabs(originalTrack2.eta()) > cfgtrkMaxEta)
continue;

// check PID
Expand Down Expand Up @@ -582,13 +589,13 @@ struct phiInJets {
lDecayDaughter1.SetXYZM(originalTrack.px(), originalTrack.py(), originalTrack.pz(), massKa);
lDecayDaughter2.SetXYZM(originalTrack2.px(), originalTrack2.py(), originalTrack2.pz(), massKa);
lResonance = lDecayDaughter1 + lDecayDaughter2;
if (lResonance.Rapidity() > 0.5)
if (fabs(lResonance.Eta()) > cfgtrkMaxEta)
continue;
JEhistos.fill(HIST("ptJEHistogramPhi"), lResonance.Pt());

bool jetFlag = false;
for (int i = 0; i < mcd_pt.size(); i++) {
double phidiff = TVector2::Phi_mpi_pi(mcd_pt[i] - lResonance.Phi());
double phidiff = TVector2::Phi_mpi_pi(mcd_phi[i] - lResonance.Phi());
double etadiff = mcd_eta[i] - lResonance.Eta();
double R = TMath::Sqrt((etadiff * etadiff) + (phidiff * phidiff));
if (R < cfgjetR)
Expand Down Expand Up @@ -618,8 +625,13 @@ struct phiInJets {
JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_OUTSIDE"), 1.0, lResonance.Pt(), lResonance.M());
} //! jetflag

if (hasJets)
if (hasJets){
JEhistos.fill(HIST("ptJEHistogramPhi_JetTrigger"), lResonance.Pt());
auto triggerjet = std::min_element(mcd_pt.begin(), mcd_pt.end());
double triggerjet_pt = *triggerjet;
JEhistos.fill(HIST("JetVsPhi_REC"), triggerjet_pt, lResonance.Pt());

}
JEhistos.fill(HIST("minvJEHistogramPhi"), lResonance.M());
} // mcpart check
} // tracks2
Expand All @@ -631,8 +643,9 @@ struct phiInJets {
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
int nprocessSimEvents = 0;
void processSim(o2::aod::JMcCollision const& collision, aod::JMcParticles const& mcParticles, soa::Filtered<aod::ChargedMCParticleLevelJets> const& mcpjets)
{
// Preslice<aod::JCollisions> slice = o2::aod::JCollision::collisionId;
void processSim(o2::aod::JMcCollision const& collision, soa::SmallGroups<soa::Join<aod::JMcCollisionLbs,aod::JCollisions>> const& recocolls, aod::JMcParticles const& mcParticles, soa::Filtered<aod::ChargedMCParticleLevelJets> const& mcpjets)
{
if (cDebugLevel > 0) {
nprocessSimEvents++;
if ((nprocessSimEvents + 1) % 10000 == 0)
Expand All @@ -641,18 +654,29 @@ struct phiInJets {

JEhistos.fill(HIST("nEvents_MCGen"), 0.5);

if (fabs(collision.posZ()) > 10)
if(recocolls.size()<=0)//not reconstructed
return;

for(auto& recocoll : recocolls){ //poorly reconstructed
if (!jetderiveddatautilities::selectCollision(recocoll, jetderiveddatautilities::JCollisionSel::sel8))
return;
}



if (fabs(collision.posZ()) > 10)//bad vertex
return;
bool INELgt0 = false;
for (const auto& mcParticle : mcParticles) {
if (TMath::Abs(mcParticle.eta()) < 0.8) {
if (fabs(mcParticle.eta()) < cfgtrkMaxEta) {
INELgt0 = true;
break;
}
}
if (!INELgt0)
if (!INELgt0) //not INEL
return;


std::vector<double> mcp_pt{};
std::vector<double> mcp_phi{};
std::vector<double> mcp_eta{};
Expand All @@ -675,15 +699,15 @@ struct phiInJets {

// Check pikp and phi
for (const auto& mcParticle : mcParticles) {
if (mcParticle.isPhysicalPrimary() && fabs(mcParticle.y()) <= 0.5) { // watch out for context!!!
if (mcParticle.isPhysicalPrimary() && fabs(mcParticle.eta()) <= cfgtrkMaxEta) { // watch out for context!!!
if (abs(mcParticle.pdgCode()) == 211)
JEhistos.fill(HIST("ptGeneratedPion"), mcParticle.pt());
if (abs(mcParticle.pdgCode()) == 321)
JEhistos.fill(HIST("ptGeneratedKaon"), mcParticle.pt());
if (abs(mcParticle.pdgCode()) == 2212)
JEhistos.fill(HIST("ptGeneratedProton"), mcParticle.pt());
}
if (fabs(mcParticle.y()) <= 0.5) { // watch out for context!!!
if (fabs(mcParticle.eta()) <= cfgtrkMaxEta) { // watch out for context!!!
TLorentzVector lResonance;
lResonance.SetPxPyPzE(mcParticle.px(), mcParticle.py(), mcParticle.pz(), mcParticle.e());
if (abs(mcParticle.pdgCode()) == 333) {
Expand All @@ -695,7 +719,7 @@ struct phiInJets {
lResonance.SetPxPyPzE(mcParticle.px(), mcParticle.py(), mcParticle.pz(), mcParticle.e());
bool jetFlag = false;
for (int i = 0; i < mcp_pt.size(); i++) {
double phidiff = TVector2::Phi_mpi_pi(mcp_pt[i] - lResonance.Phi());
double phidiff = TVector2::Phi_mpi_pi(mcp_phi[i] - lResonance.Phi());
double etadiff = mcp_eta[i] - lResonance.Eta();
double R = TMath::Sqrt((etadiff * etadiff) + (phidiff * phidiff));
if (R < cfgjetR)
Expand Down Expand Up @@ -727,9 +751,11 @@ struct phiInJets {
} //! jetflag

////////////////////////////Phi found

if (hasJets) {
JEhistos.fill(HIST("ptGeneratedPhi_JetTrigger"), mcParticle.pt());
auto triggerjet = std::min_element(mcp_pt.begin(), mcp_pt.end());
double triggerjet_pt = *triggerjet;
JEhistos.fill(HIST("JetVsPhi_GEN"), triggerjet_pt, mcParticle.pt());
} // check for jets

} // check for phi
Expand All @@ -747,6 +773,7 @@ struct phiInJets {
// void processMatchedGen(o2::aod::JMcCollision const& collision, aod::JMcParticles const& mcParticles, soa::Filtered<aod::ChargedMCParticleLevelJets> const& mcpjets)
int nprocessSimJEEvents = 0;
void processMatchedGen(aod::JMcCollision const& collision,
soa::SmallGroups<soa::Join<aod::JMcCollisionLbs,aod::JCollisions>> const& recocolls,
JetMCDTable const& mcdjets,
JetMCPTable const& mcpjets,
aod::JMcParticles const& mcParticles)
Expand All @@ -763,11 +790,17 @@ struct phiInJets {
if (fabs(collision.posZ()) > 10)
return;

if (fabs(collision.posZ()) > 10)
if(recocolls.size()<=0)//not reconstructed
return;
for(auto& recocoll : recocolls){ //poorly reconstructed
if (!jetderiveddatautilities::selectCollision(recocoll, jetderiveddatautilities::JCollisionSel::sel8))
return;
}


bool INELgt0 = false;
for (const auto& mcParticle : mcParticles) {
if (TMath::Abs(mcParticle.eta()) < 0.8) {
if (TMath::Abs(mcParticle.eta()) < cfgtrkMaxEta) {
INELgt0 = true;
break;
}
Expand Down Expand Up @@ -812,17 +845,16 @@ struct phiInJets {

// First we do GEN part
for (const auto& mcParticle : mcParticles) {
if (fabs(mcParticle.y()) > 0.5)
continue;
if (fabs(mcParticle.eta() > 0.8))
if (fabs(mcParticle.eta()) > cfgtrkMaxEta)
continue;
if (abs(mcParticle.pdgCode()) == 333) {

if (fabs(mcParticle.pdgCode()) == 333) {

TLorentzVector lResonance;
lResonance.SetPxPyPzE(mcParticle.px(), mcParticle.py(), mcParticle.pz(), mcParticle.e());
bool jetFlag = false;
for (int i = 0; i < mcp_pt.size(); i++) {
double phidiff = TVector2::Phi_mpi_pi(mcp_pt[i] - lResonance.Phi());
double phidiff = TVector2::Phi_mpi_pi(mcp_phi[i] - lResonance.Phi());
double etadiff = mcp_eta[i] - lResonance.Eta();
double R = TMath::Sqrt((etadiff * etadiff) + (phidiff * phidiff));
if (R < cfgjetR)
Expand Down Expand Up @@ -883,7 +915,7 @@ struct phiInJets {

bool INELgt0 = false;
for (const auto& track : tracks) {
if (TMath::Abs(track.eta()) < 0.8) {
if (fabs(track.eta()) < cfgtrkMaxEta) {
INELgt0 = true;
break;
}
Expand Down Expand Up @@ -931,7 +963,7 @@ struct phiInJets {
auto trk2 = track2.track_as<myCompleteTracks>();
if (trk1.index() >= trk2.index())
continue;
if (fabs(trk1.eta()) > 0.8 || fabs(trk2.eta()) > 0.8)
if (fabs(trk1.eta()) > cfgtrkMaxEta || fabs(trk2.eta()) > cfgtrkMaxEta)
continue;
if ((trk1.sign() * trk2.sign()) > 0)
continue; // Not K+K-
Expand Down Expand Up @@ -972,12 +1004,12 @@ struct phiInJets {
lDecayDaughter1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massKa);
lDecayDaughter2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massKa);
lResonance = lDecayDaughter1 + lDecayDaughter2;
if (lResonance.Rapidity() > 0.5)
if (fabs(lResonance.Eta()) > cfgtrkMaxEta)
continue;

bool jetFlag = false;
for (int i = 0; i < mcd_pt.size(); i++) {
double phidiff = TVector2::Phi_mpi_pi(mcd_pt[i] - lResonance.Phi());
double phidiff = TVector2::Phi_mpi_pi(mcd_phi[i] - lResonance.Phi());
double etadiff = mcd_eta[i] - lResonance.Eta();
double R = TMath::Sqrt((etadiff * etadiff) + (phidiff * phidiff));
if (R < cfgjetR)
Expand Down

0 comments on commit da169da

Please sign in to comment.