Skip to content

Commit acf7757

Browse files
smaff92smaff92
committed
[PWGJE]: Adding MC weights for JEJE injections, and adding mask for software triggers
1 parent c20e1db commit acf7757

File tree

1 file changed

+82
-87
lines changed

1 file changed

+82
-87
lines changed

‎PWGJE/Tasks/statPromptPhoton.cxx‎

Lines changed: 82 additions & 87 deletions
Original file line numberDiff line numberDiff line change
@@ -13,40 +13,58 @@
1313
/// \brief Reconstruction of Phi yield through track-track Minv correlations for resonance hadrochemistry analysis.
1414
///
1515
///
16-
/// \author Adrian Fereydon Nassirpour <adrian.fereydon.nassirpour@cern.ch>
16+
/// \author Adrian Fereydon Nassirpour <adrian.fereydon.nassirpour@cern.ch>
1717

18-
#include "PWGJE/Core/JetDerivedDataUtilities.h"
19-
#include "PWGJE/DataModel/EMCALClusters.h"
20-
#include "PWGJE/DataModel/JetReducedData.h"
18+
#include <CCDB/BasicCCDBManager.h>
2119

22-
#include "Common/DataModel/EventSelection.h"
23-
#include "Common/DataModel/TrackSelectionTables.h"
20+
#include <iostream>
21+
#include <vector>
22+
#include <string>
23+
#include <optional>
24+
25+
#include <TLorentzVector.h>
26+
#include <TVector2.h>
2427

2528
#include "Framework/ASoA.h"
2629
#include "Framework/AnalysisDataModel.h"
2730
#include "Framework/AnalysisTask.h"
2831
#include "Framework/HistogramRegistry.h"
29-
#include <Framework/Configurable.h>
30-
#include <Framework/HistogramSpec.h>
31-
#include <Framework/InitContext.h>
32-
#include <Framework/OutputObjHeader.h>
33-
#include <Framework/runDataProcessing.h>
32+
#include "Framework/runDataProcessing.h"
33+
#include "ReconstructionDataFormats/Track.h"
3434

35-
#include <TLorentzVector.h>
36-
#include <TMath.h>
37-
#include <TVector2.h>
35+
#include "Common/Core/RecoDecay.h"
36+
#include "Common/Core/TrackSelection.h"
37+
#include "Common/Core/TrackSelectionDefaults.h"
38+
#include "Common/Core/trackUtilities.h"
39+
#include "Common/DataModel/EventSelection.h"
40+
#include "Common/DataModel/TrackSelectionTables.h"
41+
#include "Common/DataModel/Multiplicity.h"
42+
#include "Common/DataModel/PIDResponse.h"
43+
#include "CommonConstants/PhysicsConstants.h"
3844

39-
#include <cmath>
40-
#include <iostream>
41-
#include <optional>
42-
#include <string>
43-
#include <vector>
45+
#include "PWGJE/Core/FastJetUtilities.h"
46+
#include "PWGJE/Core/JetDerivedDataUtilities.h"
47+
#include "PWGJE/DataModel/Jet.h"
48+
#include "PWGJE/DataModel/EMCALClusters.h"
49+
#include "EMCALBase/Geometry.h"
50+
#include "EMCALCalib/BadChannelMap.h"
51+
52+
#include "DataFormatsEMCAL/Cell.h"
53+
#include "DataFormatsEMCAL/Constants.h"
54+
#include "DataFormatsEMCAL/AnalysisCluster.h"
55+
#include "DataFormatsParameters/GRPObject.h"
56+
#include "DataFormatsParameters/GRPMagField.h"
57+
58+
#include "DetectorsBase/Propagator.h"
59+
60+
#include "CommonDataFormat/InteractionRecord.h"
4461

4562
using namespace o2;
4663
using namespace o2::framework;
4764
using namespace o2::framework::expressions;
4865

4966
struct statPromptPhoton {
67+
5068
SliceCache cache;
5169
HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject};
5270
Configurable<double> cfgMaxDCArToPVcut{"cfgMaxDCArToPVcut", 0.5, "Track DCAr cut to PV Maximum"};
@@ -85,15 +103,20 @@ struct statPromptPhoton {
85103
Configurable<bool> cfgGenHistograms{"cfgGenHistograms", false, "Enables Generated histograms"};
86104
Configurable<bool> cfgRecHistograms{"cfgRecHistograms", false, "Enables Reconstructed histograms"};
87105
Configurable<bool> cfgDataHistograms{"cfgDataHistograms", false, "Enables Data histograms"};
106+
Configurable<bool> cfgSkimmedTrigger{"cfgSkimmedTrigger", false, "Enables trigger for skimmied datasets (2023 onwards)"};
107+
Configurable<std::string> cfgTriggerMasks{"cfgTriggerMasks", "", "possible JE Trigger masks: fJetChLowPt,fJetChHighPt,fTrackLowPt,fTrackHighPt,fJetD0ChLowPt,fJetD0ChHighPt,fJetLcChLowPt,fJetLcChHighPt,fEMCALReadout,fJetFullHighPt,fJetFullLowPt,fJetNeutralHighPt,fJetNeutralLowPt,fGammaVeryHighPtEMCAL,fGammaVeryHighPtDCAL,fGammaHighPtEMCAL,fGammaHighPtDCAL,fGammaLowPtEMCAL,fGammaLowPtDCAL,fGammaVeryLowPtEMCAL,fGammaVeryLowPtDCAL"};
88108
Configurable<bool> cfgDebug{"cfgDebug", false, "Enables debug information for local running"};
109+
89110
int trackFilter = -1;
111+
std::vector<int> triggerMaskBits;
90112

91113
// INIT
92114
void init(InitContext const&)
93115
{
94116
std::vector<double> ptBinning = {0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.5, 3.0, 4.0, 5.0, 6.0, 8.0, 12.0, 16.0, 20.0, 25.0, 30.0, 40.0, 50.0, 75.0, 100.0, 150.0, 200.0, 300.0, 500.0};
95117
AxisSpec pthadAxis = {ptBinning, "#it{p}_{T}^{had sum} [GeV/c]"};
96118

119+
triggerMaskBits = jetderiveddatautilities::initialiseTriggerMaskBits(cfgTriggerMasks);
97120
if (cfgJETracks) {
98121
trackFilter = jetderiveddatautilities::initialiseTrackSelection(static_cast<std::string>(cfgTrackFilter));
99122
}
@@ -235,8 +258,11 @@ struct statPromptPhoton {
235258

236259
using jMCClusters = o2::soa::Join<o2::aod::JMcClusterLbs, o2::aod::JClusters, o2::aod::JClusterTracks>;
237260
using jClusters = o2::soa::Join<o2::aod::JClusters, o2::aod::JClusterTracks>;
238-
using jselectedCollisions = soa::Join<aod::JCollisions, aod::JCollisionBCs, aod::JCollisionPIs, aod::EvSels, aod::JEMCCollisionLbs>;
261+
using jselectedCollisions = soa::Join<aod::JCollisions, aod::JCollisionBCs, aod::JCollisionPIs, aod::EvSels, aod::JEMCCollisionLbs, aod::JMcCollisionLbs>;
262+
using jselectedDataCollisions = soa::Join<aod::JCollisions, aod::JCollisionBCs, aod::JCollisionPIs, aod::EvSels, aod::JEMCCollisionLbs>;
263+
// using jselectedDataCollisions = soa::Join<aod::JCollisions, aod::JCollisionBCs, aod::JCollisionPIs, aod::JCollisionMcInfos, aod::EvSels, aod::JEMCCollisionLbs>;
239264
using jfilteredCollisions = soa::Filtered<jselectedCollisions>;
265+
using jfilteredDataCollisions = soa::Filtered<jselectedDataCollisions>;
240266
using jfilteredMCClusters = soa::Filtered<jMCClusters>;
241267
using jfilteredClusters = soa::Filtered<jClusters>;
242268

@@ -588,7 +614,7 @@ struct statPromptPhoton {
588614

589615
PresliceUnsorted<jEMCtracks> EMCTrackPerTrack = aod::jemctrack::trackId;
590616
int nEventsRecMC_JE = 0;
591-
void processMCRec_JE(jfilteredCollisions::iterator const& collision, jfilteredMCClusters const& mcclusters, jTrackCandidates const& tracks, soa::Join<aod::JTracks, aod::JTrackExtras, aod::JTrackPIs> const&, TrackCandidates const&, aod::JMcParticles const&, BcCandidates const&, jEMCtracks const& emctracks)
617+
void processMCRec_JE(jfilteredCollisions::iterator const& collision, jfilteredMCClusters const& mcclusters, jTrackCandidates const& tracks, soa::Join<aod::JTracks, aod::JTrackExtras, aod::JTrackPIs> const&, TrackCandidates const&, aod::JMcParticles const&, BcCandidates const&, jEMCtracks const& emctracks, aod::JetMcCollisions const&)
592618
{
593619

594620
nEventsRecMC_JE++;
@@ -613,6 +639,20 @@ struct statPromptPhoton {
613639
}
614640
histos.fill(HIST("REC_nEvents"), 2.5);
615641

642+
643+
if(cfgSkimmedTrigger){
644+
if (!jetderiveddatautilities::selectTrigger(collision, triggerMaskBits)) {
645+
return;
646+
}
647+
}//JE Software Triggers
648+
649+
histos.fill(HIST("REC_nEvents"), 3.5);
650+
651+
double weight = 1;
652+
if (collision.has_mcCollision()) {
653+
weight = collision.mcCollision().weight();
654+
}
655+
616656
bool noTrk = true;
617657
for (auto& track : tracks) {
618658
if (cfgJETracks) {
@@ -828,9 +868,9 @@ struct statPromptPhoton {
828868
histos.fill(HIST("REC_Cluster_QA"), 4.5);
829869
clustertrigger = true;
830870
double pthadsum = GetPtHadSum(tracks, mccluster, cfgMinR, cfgMaxR, false, false, true);
831-
histos.fill(HIST("REC_Trigger_V_PtHadSum_Photon"), mccluster.energy(), pthadsum);
832-
histos.fill(HIST("REC_PtHadSum_Photon"), pthadsum);
833-
histos.fill(HIST("REC_Trigger_Energy"), mccluster.energy());
871+
histos.fill(HIST("REC_Trigger_V_PtHadSum_Photon"), mccluster.energy(), pthadsum, weight);
872+
histos.fill(HIST("REC_PtHadSum_Photon"), pthadsum, weight);
873+
histos.fill(HIST("REC_Trigger_Energy"), mccluster.energy(), weight);
834874
}
835875

836876
auto ClusterParticles = mccluster.mcParticles_as<aod::JMcParticles>();
@@ -839,11 +879,6 @@ struct statPromptPhoton {
839879
bool goodgentrigger = true;
840880
double chPe = 0;
841881
for (auto& clusterparticle : ClusterParticles) {
842-
// double etaP = clusterparticle.eta();
843-
// double etaC = mccluster.eta();
844-
// double phiP = clusterparticle.phi();
845-
// double phiC = mccluster.phi();
846-
// double ptP = clusterparticle.pt();
847882
int cindex = clusterparticle.globalIndex();
848883
double pdgcode = fabs(clusterparticle.pdgCode());
849884
if (!clusterparticle.isPhysicalPrimary()) {
@@ -884,8 +919,6 @@ struct statPromptPhoton {
884919
histos.fill(HIST("REC_Cluster_ParticleWITHtrack_Phi"), clusterparticle.phi());
885920
histos.fill(HIST("REC_Cluster_ParticleWITHtrack_Eta"), clusterparticle.eta());
886921
histos.fill(HIST("REC_Cluster_ParticleWITHtrack_Pt_Phi"), clusterparticle.pt(), clusterparticle.phi());
887-
// if (phiPrimeP > (0.12/ptP + TMath::Pi()/18. + 0.035) ||
888-
// phiPrimeP < (0.1/ptP/ptP + TMath::Pi()/18. - 0.025) ) {
889922
histos.fill(HIST("REC_Cluster_ParticleWITHtrack_Pt_PhiPrime"), ptP, phiPrimeP);
890923
if (photontrigger) {
891924
histos.fill(HIST("REC_Impurity_ParticleWITHtrack_Pt_PhiPrime"), ptP, phiPrimeP);
@@ -904,8 +937,6 @@ struct statPromptPhoton {
904937
histos.fill(HIST("REC_Cluster_ParticleWITHOUTtrack_Phi"), clusterparticle.phi());
905938
histos.fill(HIST("REC_Cluster_ParticleWITHOUTtrack_Eta"), clusterparticle.eta());
906939
histos.fill(HIST("REC_Cluster_ParticleWITHOUTtrack_Pt_Phi"), clusterparticle.pt(), clusterparticle.phi());
907-
// if (phiPrimeP > (0.12/ptP + TMath::Pi()/18. + 0.035) ||
908-
// phiPrimeP < (0.1/ptP/ptP + TMath::Pi()/18. - 0.025) ) {
909940
histos.fill(HIST("REC_Cluster_ParticleWITHOUTtrack_Pt_PhiPrime"), ptP, phiPrimeP);
910941
if (photontrigger) {
911942
histos.fill(HIST("REC_Impurity_ParticleWITHOUTtrack_Pt_PhiPrime"), ptP, phiPrimeP);
@@ -964,12 +995,12 @@ struct statPromptPhoton {
964995
std::cout << "Photon mom 2: " << mom2 << std::endl;
965996
}
966997
if (std::abs(clusterparticle.getGenStatusCode()) > 19 && std::abs(clusterparticle.getGenStatusCode()) < 90) {
967-
histos.fill(HIST("REC_True_Prompt_Trigger_Energy"), clusterparticle.e());
998+
histos.fill(HIST("REC_True_Prompt_Trigger_Energy"), clusterparticle.e(), weight);
968999
TLorentzVector lRealPhoton;
9691000
lRealPhoton.SetPxPyPzE(clusterparticle.px(), clusterparticle.py(), clusterparticle.pz(), clusterparticle.e());
9701001
double truepthadsum = GetPtHadSum(tracks, lRealPhoton, cfgMinR, cfgMaxR, false, false, false);
9711002
truephotonPt = clusterparticle.e();
972-
histos.fill(HIST("REC_TrueTrigger_V_PtHadSum_Photon"), truephotonPt, truepthadsum);
1003+
histos.fill(HIST("REC_TrueTrigger_V_PtHadSum_Photon"), truephotonPt, truepthadsum, weight);
9731004
}
9741005
} // photon check
9751006
} // clusterparticle loop
@@ -1000,39 +1031,6 @@ struct statPromptPhoton {
10001031
}
10011032
} // cluster loop
10021033

1003-
// auto bc = collision.bc_as<BcCandidates>();
1004-
// int rnr = bc.runNumber();
1005-
1006-
// std::string rnrstring = std::to_string(rnr);
1007-
// if (runs.find(rnrstring) == std::string::npos) {
1008-
// std::cout<<"++++++++++++++++++++++++++++++++"<<std::endl;
1009-
// std::cout<<"FETCHING NEW RUN NUMBER FOR RUN: "<<rnr<<std::endl;
1010-
// std::cout<<"++++++++++++++++++++++++++++++++"<<std::endl;
1011-
1012-
// std::string ccdbpath="GLO/Config/GRPMagField";
1013-
// static o2::parameters::GRPMagField* grpmag = ccdb->getForTimeStamp<o2::parameters::GRPMagField>(ccdbpath, bc.timestamp());
1014-
// if(grpmag) {
1015-
// bfield = std::lround(5.f * grpmag->getL3Current() / 30000.f);
1016-
// std::cout<<"++++++++++++++++++++++++++++++++"<<std::endl;
1017-
// std::cout<<"MAG FIELD IS: "<<bfield<<std::endl;
1018-
// std::cout<<"++++++++++++++++++++++++++++++++"<<std::endl;
1019-
// }
1020-
// else {
1021-
// ccdbpath="GLO/GRP/GRP";
1022-
// static o2::parameters::GRPObject* grpo = ccdb->getForTimeStamp<o2::parameters::GRPObject>(ccdbpath, bc.timestamp());
1023-
// if(!grpo) {
1024-
// std::cout<<"WE CAN NEITHER FETCH GRPMAG OR GRPO!!! SHIT IS SCREWED"<<std::endl;
1025-
// }
1026-
// bfield = grpo->getNominalL3Field();
1027-
// }
1028-
// bfield = 5;
1029-
// runs += rnrstring;
1030-
// std::cout << "++++++++++++++++++++++++++++++++" << std::endl;
1031-
// std::cout << "Run is now appended to string: " << runs << std::endl;
1032-
// std::cout << "++++++++++++++++++++++++++++++++" << std::endl;
1033-
1034-
// } // check mag field for current run number: done!
1035-
10361034
// clusters done, now we do the sternheimer tracks
10371035
for (auto& track : tracks) {
10381036
bool sterntrigger = false;
@@ -1057,19 +1055,9 @@ struct statPromptPhoton {
10571055
phiPrime = 2 * TMath::Pi() - phiPrime;
10581056
}
10591057

1060-
// if (bfield < 0) {
1061-
// phiPrime = 2 * TMath::Pi() - phiPrime;
1062-
// }
1063-
10641058
phiPrime = phiPrime + TMath::Pi() / 18.;
10651059
phiPrime = fmod(phiPrime, 2 * TMath::Pi() / 18.);
1066-
// double pt = track.pt();
1067-
// if (phiPrime > (0.12/pt + TMath::Pi()/18. + 0.035) ||
1068-
// phiPrime < (0.1/pt/pt + TMath::Pi()/18. - 0.025) ) {
10691060
histos.fill(HIST("REC_Track_PhiPrime_Pt"), phiPrime, track.pt());
1070-
// }//geo cut
1071-
// Done with geometric cuts
1072-
10731061
histos.fill(HIST("REC_Track_Pt"), track.pt());
10741062
histos.fill(HIST("REC_Track_Phi"), track.phi());
10751063
if (clustertrigger) {
@@ -1083,12 +1071,12 @@ struct statPromptPhoton {
10831071
}
10841072
}
10851073
double pthadsum = GetPtHadSum(tracks, track, cfgMinR, cfgMaxR, true, false, true);
1086-
histos.fill(HIST("REC_Trigger_V_PtHadSum_Nch"), sternPt, pthadsum);
1074+
histos.fill(HIST("REC_Trigger_V_PtHadSum_Nch"), sternPt, pthadsum,weight);
10871075
if (sterntrigger) {
10881076
bool doStern = true;
10891077
double sterncount = 1.0;
10901078
while (doStern) {
1091-
histos.fill(HIST("REC_Trigger_V_PtHadSum_Stern"), sterncount, pthadsum, 2.0 / sternPt);
1079+
histos.fill(HIST("REC_Trigger_V_PtHadSum_Stern"), sterncount, pthadsum, (2.0 / sternPt)*weight);
10921080
if (sterncount < sternPt) {
10931081
sterncount++;
10941082
} else {
@@ -1103,18 +1091,21 @@ struct statPromptPhoton {
11031091
PROCESS_SWITCH(statPromptPhoton, processMCRec_JE, "processJE MC data", false);
11041092

11051093
int nEventsData = 0;
1106-
void processData(jfilteredCollisions::iterator const& collision, jfilteredClusters const& clusters, jDataTrackCandidates const& tracks, soa::Join<aod::JTracks, aod::JTrackExtras, aod::JTrackPIs> const&, TrackCandidates const&, BcCandidates const&, jEMCtracks const& emctracks)
1094+
void processData(jfilteredDataCollisions::iterator const& collision, jfilteredClusters const& clusters, jDataTrackCandidates const& tracks, soa::Join<aod::JTracks, aod::JTrackExtras, aod::JTrackPIs> const&, TrackCandidates const&, BcCandidates const&, jEMCtracks const& emctracks)
11071095
{
1108-
11091096
nEventsData++;
11101097
if (cfgDebug) {
11111098
if (nEventsData == 1) {
11121099
std::cout << "Starting Data Processing: " << nEventsData << std::endl;
11131100
}
11141101
if ((nEventsData + 1) % 10000 == 0) {
11151102
std::cout << "Processed Data Events: " << nEventsData << std::endl;
1103+
std::cout << "Events Trigger Bit: " << collision.triggerSel() << std::endl;
1104+
std::cout << "Trigger Mask Bit: " << triggerMaskBits[0] << std::endl;
1105+
std::cout << "Trigger Mask Cfg Line: " << cfgTriggerMasks << std::endl;
11161106
}
11171107
}
1108+
11181109
histos.fill(HIST("DATA_nEvents"), 0.5);
11191110

11201111
// required cuts
@@ -1124,13 +1115,21 @@ struct statPromptPhoton {
11241115
return;
11251116

11261117
histos.fill(HIST("DATA_nEvents"), 1.5);
1127-
11281118
if (cfgEmcTrigger) {
11291119
if (!collision.isEmcalReadout())
11301120
return;
11311121
}
1122+
11321123
histos.fill(HIST("DATA_nEvents"), 2.5);
11331124

1125+
if(cfgSkimmedTrigger){
1126+
if (!jetderiveddatautilities::selectTrigger(collision, triggerMaskBits)) {
1127+
return;
1128+
}
1129+
}//JE Software Triggers
1130+
1131+
histos.fill(HIST("DATA_nEvents"), 3.5);
1132+
11341133
bool noTrk = true;
11351134
for (auto& track : tracks) {
11361135

@@ -1331,10 +1330,6 @@ struct statPromptPhoton {
13311330
phiPrime = 2 * TMath::Pi() - phiPrime;
13321331
}
13331332

1334-
// if (bfield < 0) {
1335-
// phiPrime = 2 * TMath::Pi() - phiPrime;
1336-
// }
1337-
13381333
phiPrime = phiPrime + TMath::Pi() / 18.;
13391334
phiPrime = fmod(phiPrime, 2 * TMath::Pi() / 18.);
13401335
double pt = track.pt();

0 commit comments

Comments
 (0)