-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathDiJetAnalysis.h
227 lines (200 loc) · 8.7 KB
/
DiJetAnalysis.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
/**
* @file DiJetAnalysis.h
* @author Grigory Nigmatkulov ([email protected])
* @brief Dijet analysis
* @version 1.1
* @date 2025-01-09
*
* @copyright Copyright (c) 2025
*
*/
#ifndef DiJetAnalysis_h
#define DiJetAnalysis_h
// Load ROOT libraries
#include "TObject.h"
#include "TString.h"
#include "Rtypes.h"
#include "TChain.h"
#include "TF1.h"
// Jet analysis headers
#include "BaseAnalysis.h"
#include "HistoManagerDiJet.h"
#include "Event.h"
// C++ headers
#include <vector>
//________________
class DiJetAnalysis : public BaseAnalysis {
public:
/// @brief Default constructor
DiJetAnalysis();
/// @brief Destructor
virtual ~DiJetAnalysis();
/// @brief Initialize variables and functions
void init();
/// @brief Process event
void processEvent(const Event* ev);
/// @brief Finish analysis
void finish();
/// @brief Returns reports of all cuts applied and correlation functions being done
virtual void report();
/// @brief Return a TList of objects to be written as output
virtual TList* getOutputList();
/// @brief Add histogram manager to the analysis
void addHistoManager(HistoManagerDiJet *hm) { fHM = hm; }
/// @brief Add lorentz shift
void setEtaShift(const double& shift) { fEtaShift = shift; }
/// @brief Set dataset to be MC
void setIsMc(const bool& isMc) { fIsMc = isMc; }
/// @brief Set collision system: 0 - pp, 1 - pPb, 2 - PbPb
void setCollisionSystem(const int& syst) { fCollisionSystem = syst; }
/// @brief Set collision energy in GeV (default: 8160 GeV)
void setCollisionEnergyInGeV(const int& en) { fCollisionEnergy = en; }
/// @brief Set cut on the ptHat of the event (for MC in pPb only due to the xsection matching)
void setPtHatRange(const double& lo, const double& hi) { fPtHatRange[0] = lo; fPtHatRange[1] = hi; }
/// @brief Cut on the lowest momentum of leading jet
void setLeadJetPtLow(const double& lo) { fLeadJetPtLow = lo; }
/// @brief Cut on the lowest momentum of subleading jet
void setSubLeadJetPtLow(const double& lo) { fSubleadJetPtLow = lo; }
/// @brief Cut on angle between leading and subleading jet
void setDijetPhiCut(const double& cut) { fDijetPhiCut = cut; }
/// @brief Set the direction of Pb-going ion
void setPbGoing() { fIsPbGoingDir = {true}; }
/// @brief Set the direction of p-going ion
void setpGoing() { fIsPbGoingDir = {false}; }
/// @brief Set verbose mode
void setVerbose() { fVerbose = {true}; }
/// @brief Set number of events in the embedding sample (for the given ptHat)
void setNEventsInSample(const int& n) { fNEventsInSample = n; }
/// @brief Set loose jetId cut
void setLooseJetIdCut() { fIsLooseJetIdCut = {true}; }
/// @brief Select inclusive jets in the center-of-mass frame
void selectJetsInCMFrame() { fSelectJetsInCMFrame = {true}; }
/// @brief Set eta range to select jets in the lab frame
void setJetEtaLabRange(const double& lo, const double& hi) { fJetEtaLab[0]=lo; fJetEtaLab[1]=hi; }
/// @brief Set eta range to select jets in the center-of-mass frame
void setJetEtaCMRange(const double& lo, const double& hi) { fJetEtaCM[0]=lo; fJetEtaCM[1]=hi; }
/// @brief Reweight MC to data (trigger-dependent):
/// 0 - do not reweight (default)
/// 1 - MB
/// 2 - Jet60
/// 3 - Jet80
/// 4 - Jet100
void setUseMcReweighting(const int& w = 0) { fUseMcReweighting = (short)w; }
/// @brief Set jetId selection of the jets (default: trkMax)
void useJetIdSelection() { fUseJetIdSelection = {true}; }
void findMcWeight(const double& ptLead, const double& ptSublead);
/// @brief Print DiJetAnalysis setup
void print();
/// @brief Return collision system name
TString collisionSystem() const;
private:
/// @brief Calculate event weight
double eventWeight(const double& ptHat, const double& vz, const double& centWeight, const double& ptHatW);
/// @brief Process gen jets
void processGenJets(const Event* event, double ptHatW);
/// @brief Process reco jets
void processRecoJets(const Event* event, double ptHatW);
/// @brief Process ref jets
void processRefJets(const Event* event, double ptHatW);
/// @brief Dijet selection
bool isGoodDijet(const double& ptLead, const double& ptSublead, const double& dphi);
/// @brief Dijet selection
bool isGoodDijet(const double& ptLead, const double& etaLead, const double& ptSubLead,
const double& etaSubLead, const double& dphi, const bool& isCM = false);
/// @brief Calculate delta phi between two jets in the range [-pi, pi]
double deltaPhi(const double& phi1, const double &phi2);
/// @brief Single gen/ref jet selection criteria
bool isGoodGenJet(const GenJet* jet);
/// @brief Single reco jet selection criteria
bool isGoodRecoJet(const RecoJet* jet);
/// @brief Check if jet passes jetId requirements
bool isGoodJetId(const RecoJet* jet);
/// @brief Check if good track max cut
bool isGoodTrkMax(const RecoJet* jet);
/// @brief Boost eta to the center-of-mass frame
double boostEta2CM(const double &etaLab);
/// @brief Get proper eta in the lab frame depending on beam direction
double etaLab(const double &eta);
/// @brief Dijet eta calculation
double dijetEtaInFrame(const double& eta1, const double& eta2, bool isCM = false);
/// @brief Pass pt of the jet and check if it is leading or subleading jet
void findLeadSubleadJets(const double &pt, const int &counter, double &ptLead, double &ptSublead,
int &idLead, int &idSubLead);
/// @brief Find dijet ptAve bin
int findDijetPtAveBin(const double& pt);
/// @brief Find dijet ptAve bin (old binning)
int findDijetPtAveOldBin(const double& pt);
/// @brief Initialize vz weight function
void initVzWeightFunction();
/// @brief Vz weight to match MC to data
TF1 *fVzWeight;
/// @brief Dijet ptAve weight (to match PYTHIA 2 pp data)
TF1 *fDijetPtAveWeight;
/// @brief Centrality weight
bool fUseCentralityWeight;
/// @brief Histogram manager
HistoManagerDiJet *fHM;
/// @brief Pseudorapidity shift for asymmetric collisions (pPb)
double fEtaShift;
/// @brief Is MC sample (needed for event weight corrections)
bool fIsMc;
/// @brief Type of collision system: 0 - pp, 1 - pPb, 2 - PbPb. Default is PbPb
int fCollisionSystem;
/// @brief Collision energy in GeV
int fCollisionEnergy;
/// @brief ptHat range for the generated events (must cut events on this one)
double fPtHatRange[2];
/// @brief Momentum selection of the leading jet
double fLeadJetPtLow;
/// @brief Momentum selection of the subleading jet
double fSubleadJetPtLow;
/// @brief Angular selection of dijet
double fDijetPhiCut;
/// @brief Lead going direction for pPb collisions
bool fIsPbGoingDir;
/// @brief Verbose mode
bool fVerbose;
/// @brief Number of events in the embedding sample
int fNEventsInSample;
/// @brief Use jetId selection (default - false, i.e. trkMax)
bool fUseJetIdSelection;
/// @brief Is loose/tight jetId cut (default: false = tight)
bool fIsLooseJetIdCut;
/// @brief Gen dijet in the lab frame found (default: false)
bool fIsGenDijetLabFound;
/// @brief Gen dijet in the center-of-mass frame found
bool fIsGenDijetCMFound;
/// @brief Reco dijet in the lab frame found
bool fIsRecoDijetLabFound;
/// @brief Reco dijet in the center-of-mass frame found
bool fIsRecoDijetCMFound;
/// @brief Ref-selected dijet in the lab frame found
bool fIsRefSelDijetLabFound;
/// @brief Ref-selected dijet in the center-of-mass frame found
bool fIsRefSelDijetCMFound;
/// @brief Select jets in the center-of-mass frame (default: false)
bool fSelectJetsInCMFrame;
/// @brief Reweight MC to data (trigger-dependent):
/// 0 - do not reweight (default)
/// 1 - MB
/// 2 - Jet60
/// 3 - Jet80
/// 4 - Jet100
short fUseMcReweighting;
int fJetPtBins;
double fJetPtLow;
double fJetPtHi;
double fJetPtStep;
double fJetPtLeadPtSubleadReweightMatrix[75][75];
double fMcReweight;
/// Range of eta selection in the lab frame
double fJetEtaLab[2];
/// Range of eta selection in the center-of-mass frame
double fJetEtaCM[2];
/// @brief Values for new dijet ptAve binning
std::vector<double> fPtAveBins;
/// @brief Values for old dijet ptAve binning
std::vector<double> fPtAveOldBins;
ClassDef(DiJetAnalysis, 0)
};
#endif // #define DiJetAnalysis_h