fc392675 |
1 | // Jet-Hadron Correlations PID |
2 | // Event plane dependence task. |
3 | // |
4 | // Author: J Mazer |
5 | |
6 | // task head include |
7 | #include "AliAnalysisTaskEmcalJetHadEPpid.h" |
8 | |
9 | // general ROOT includes |
10 | #include <TCanvas.h> |
11 | #include <TChain.h> |
12 | #include <TClonesArray.h> |
13 | #include <TH1F.h> |
14 | #include <TH2F.h> |
15 | #include <TH3F.h> |
16 | #include <THnSparse.h> |
17 | #include <TList.h> |
18 | #include <TLorentzVector.h> |
19 | #include <TParameter.h> |
20 | #include <TParticle.h> |
21 | #include <TTree.h> |
22 | #include <TVector3.h> |
23 | #include <TObjArray.h> |
24 | |
25 | // AliROOT includes |
26 | #include "AliAODEvent.h" |
27 | #include "AliESDEvent.h" |
28 | #include "AliAnalysisManager.h" |
29 | #include "AliAnalysisTask.h" |
30 | #include "AliCentrality.h" |
31 | #include "AliEmcalJet.h" |
32 | #include "AliAODJet.h" |
33 | #include "AliVCluster.h" |
34 | #include "AliVTrack.h" |
35 | #include <AliVEvent.h> |
36 | #include <AliVParticle.h> |
37 | #include "AliRhoParameter.h" |
38 | #include "AliEmcalParticle.h" |
39 | |
40 | // Localized Rho includes |
41 | #include "AliLocalRhoParameter.h" |
42 | #include "AliAnalysisTaskLocalRho.h" |
43 | |
44 | // event handler (and pico's) includes |
45 | #include <AliInputEventHandler.h> |
46 | #include <AliVEventHandler.h> |
47 | #include "AliESDInputHandler.h" |
48 | #include "AliPicoTrack.h" |
49 | #include "AliEventPoolManager.h" |
3263f852 |
50 | #include "AliESDtrackCuts.h" |
fc392675 |
51 | |
52 | // PID includes |
53 | #include "AliPIDResponse.h" |
54 | #include "AliTPCPIDResponse.h" |
55 | #include "AliESDpid.h" |
56 | |
57 | // magnetic field includes |
58 | #include "TGeoGlobalMagField.h" |
59 | #include "AliMagF.h" |
60 | |
8ee2d316 |
61 | // container includes |
62 | #include "AliJetContainer.h" |
63 | #include "AliParticleContainer.h" |
64 | #include "AliClusterContainer.h" |
65 | |
5e2539d9 |
66 | using std::cout; |
67 | using std::endl; |
68 | |
fc392675 |
69 | ClassImp(AliAnalysisTaskEmcalJetHadEPpid) |
70 | |
71 | //________________________________________________________________________ |
72 | AliAnalysisTaskEmcalJetHadEPpid::AliAnalysisTaskEmcalJetHadEPpid() : |
73 | AliAnalysisTaskEmcalJet("correlations",kFALSE), |
74 | fPhimin(-10), fPhimax(10), |
75 | fEtamin(-0.9), fEtamax(0.9), |
f1ce4cc3 |
76 | fAreacut(0.0), fTrkBias(5), fClusBias(5), fTrkEta(0.9), |
77 | fJetPtcut(15.0), fJetRad(0.4), fConstituentCut(0.15), |
3263f852 |
78 | fesdTrackCuts(0), |
a409c5ba |
79 | fDoEventMixing(0), fMixingTracks(50000), fNMIXtracks(5000), fNMIXevents(5), |
a4248490 |
80 | fCentBinSize(1), |
d60e724e |
81 | fTriggerEventType(AliVEvent::kAny), fMixingEventType(AliVEvent::kAny), |
a4248490 |
82 | fDoEffCorr(0), fEffFunctionCorr(0), |
8ee2d316 |
83 | doPlotGlobalRho(0), doVariableBinning(0), dovarbinTHnSparse(0), |
daa728b0 |
84 | makeQAhistos(0), makeBIAShistos(0), makeextraCORRhistos(0), makeoldJEThadhistos(0), |
139e13f4 |
85 | allpidAXIS(0), fcutType("EMCAL"), doPID(0), doPIDtrackBIAS(0), |
56523939 |
86 | doComments(0), |
87 | doFlavourJetAnalysis(0), fJetFlavTag(-99), |
88 | fBeam(kNA), |
fc392675 |
89 | fLocalRhoVal(0), |
56523939 |
90 | fTracksName(""), fTracksNameME(""), fJetsName(""), |
fc392675 |
91 | event(0), |
aa5f3717 |
92 | isPItpc(0), isKtpc(0), isPtpc(0), // pid TPC |
93 | isPIits(0), isKits(0), isPits(0), // pid ITS |
94 | isPItof(0), isKtof(0), isPtof(0), // pid TOF |
fc392675 |
95 | fPoolMgr(0x0), |
96 | fPIDResponse(0x0), fTPCResponse(), |
d60e724e |
97 | fESD(0), fAOD(0), fVevent(0), |
98 | fHistEventQA(0), fHistEventSelectionQA(0), |
a409c5ba |
99 | fHistCentZvertGA(0), fHistCentZvertJE(0), fHistCentZvertMB(0), fHistCentZvertAny(0), |
fc392675 |
100 | fHistTPCdEdX(0), fHistITSsignal(0), //fHistTOFsignal(0), |
101 | fHistRhovsCent(0), fHistNjetvsCent(0), fHistCentrality(0), |
102 | fHistZvtx(0), fHistMult(0), |
103 | fHistJetPhi(0), fHistTrackPhi(0), |
d60e724e |
104 | fHistLocalRhoJetpt(0), |
fc392675 |
105 | fHistJetHaddPhiIN(0), fHistJetHaddPhiOUT(0), fHistJetHaddPhiMID(0), |
106 | fHistJetHaddPhiBias(0), fHistJetHaddPhiINBias(0), fHistJetHaddPhiOUTBias(0), fHistJetHaddPhiMIDBias(0), |
107 | fHistMEdPHI(0), |
108 | fHistTrackPtallcent(0), |
109 | fHistJetEtaPhi(0), fHistJetHEtaPhi(0), |
110 | fHistSEphieta(0), fHistMEphieta(0), |
111 | fHistJetHaddPHI(0), |
112 | fHistPID(0), |
daa728b0 |
113 | fhnPID(0x0), fhnMixedEvents(0x0), fhnJH(0x0), fhnCorr(0x0), |
8ee2d316 |
114 | fJetsCont(0), fTracksCont(0), fCaloClustersCont(0), |
115 | fContainerAllJets(0), fContainerPIDJets(1) |
fc392675 |
116 | { |
117 | // Default Constructor |
118 | for(Int_t ilab=0; ilab<4; ilab++){ |
119 | for(Int_t ipta=0; ipta<7; ipta++){ |
72a4d320 |
120 | //fHistTrackEtaPhi[ilab][ipta]=0; // keep out for now |
fc392675 |
121 | } // end of pt-associated loop |
122 | } // end of lab loop |
123 | |
124 | for(Int_t itrackpt=0; itrackpt<9; itrackpt++){ |
125 | fHistJetHadbindPhi[itrackpt]=0; |
126 | fHistJetHadbindPhiIN[itrackpt]=0; |
127 | fHistJetHadbindPhiMID[itrackpt]=0; |
128 | fHistJetHadbindPhiOUT[itrackpt]=0; |
129 | } // end of trackpt bin loop |
130 | |
daa728b0 |
131 | for(Int_t icent = 0; icent<6; ++icent){ |
fc392675 |
132 | for(Int_t iptjet = 0; iptjet<5; ++iptjet){ |
133 | for(Int_t ieta = 0; ieta<3; ++ieta){ |
134 | fHistJetH[icent][iptjet][ieta]=0; |
135 | fHistJetHBias[icent][iptjet][ieta]=0; |
136 | fHistJetHTT[icent][iptjet][ieta]=0; |
137 | } // end of eta loop |
138 | } // end of pt-jet loop |
139 | } // end of centrality loop |
140 | |
141 | // centrality dependent histo's |
142 | for (Int_t i = 0;i<6;++i){ |
daa728b0 |
143 | fHistJetPt[i] = 0; |
144 | fHistJetPtBias[i] = 0; |
145 | fHistJetPtTT[i] = 0; |
146 | fHistAreavsRawPt[i] = 0; |
fc392675 |
147 | fHistJetPtvsTrackPt[i] = 0; |
148 | fHistRawJetPtvsTrackPt[i] = 0; |
149 | fHistTrackPt[i] = 0; |
150 | fHistEP0[i] = 0; |
151 | fHistEP0A[i] = 0; |
152 | fHistEP0C[i] = 0; |
153 | fHistEPAvsC[i] = 0; |
8ee2d316 |
154 | fHistJetPtcorrGlRho[i] = 0; |
fc392675 |
155 | fHistJetPtvsdEP[i] = 0; |
156 | fHistJetPtvsdEPBias[i] = 0; |
157 | fHistRhovsdEP[i] = 0; |
158 | fHistJetEtaPhiPt[i] = 0; |
159 | fHistJetEtaPhiPtBias[i] = 0; |
160 | fHistJetPtArea[i] = 0; |
161 | fHistJetPtAreaBias[i] = 0; |
162 | fHistJetPtNcon[i] = 0; |
163 | fHistJetPtNconBias[i] = 0; |
164 | fHistJetPtNconCh[i] = 0; |
165 | fHistJetPtNconBiasCh[i] = 0; |
166 | fHistJetPtNconEm[i] = 0; |
167 | fHistJetPtNconBiasEm[i] = 0; |
168 | fHistJetHaddPhiINcent[i] = 0; |
169 | fHistJetHaddPhiOUTcent[i] = 0; |
170 | fHistJetHaddPhiMIDcent[i] = 0; |
171 | } |
172 | |
173 | SetMakeGeneralHistograms(kTRUE); |
174 | |
fc392675 |
175 | } |
176 | |
177 | //________________________________________________________________________ |
178 | AliAnalysisTaskEmcalJetHadEPpid::AliAnalysisTaskEmcalJetHadEPpid(const char *name) : |
179 | AliAnalysisTaskEmcalJet(name,kTRUE), |
180 | fPhimin(-10), fPhimax(10), |
181 | fEtamin(-0.9), fEtamax(0.9), |
f1ce4cc3 |
182 | fAreacut(0.0), fTrkBias(5), fClusBias(5), fTrkEta(0.9), |
183 | fJetPtcut(15.0), fJetRad(0.4), fConstituentCut(0.15), |
3263f852 |
184 | fesdTrackCuts(0), |
a409c5ba |
185 | fDoEventMixing(0), fMixingTracks(50000), fNMIXtracks(5000), fNMIXevents(5), |
a4248490 |
186 | fCentBinSize(1), |
d60e724e |
187 | fTriggerEventType(AliVEvent::kAny), fMixingEventType(AliVEvent::kAny), |
a4248490 |
188 | fDoEffCorr(0), fEffFunctionCorr(0), |
8ee2d316 |
189 | doPlotGlobalRho(0), doVariableBinning(0), dovarbinTHnSparse(0), |
daa728b0 |
190 | makeQAhistos(0), makeBIAShistos(0), makeextraCORRhistos(0), makeoldJEThadhistos(0), |
139e13f4 |
191 | allpidAXIS(0), fcutType("EMCAL"), doPID(0), doPIDtrackBIAS(0), |
56523939 |
192 | doComments(0), |
193 | doFlavourJetAnalysis(0), fJetFlavTag(-99), |
194 | fBeam(kNA), |
fc392675 |
195 | fLocalRhoVal(0), |
56523939 |
196 | fTracksName(""), fTracksNameME(""), fJetsName(""), |
fc392675 |
197 | event(0), |
aa5f3717 |
198 | isPItpc(0), isKtpc(0), isPtpc(0), // pid TPC |
199 | isPIits(0), isKits(0), isPits(0), // pid ITS |
200 | isPItof(0), isKtof(0), isPtof(0), // pid TOF |
fc392675 |
201 | fPoolMgr(0x0), |
202 | fPIDResponse(0x0), fTPCResponse(), |
d60e724e |
203 | fESD(0), fAOD(0), fVevent(0), |
204 | fHistEventQA(0), fHistEventSelectionQA(0), |
a409c5ba |
205 | fHistCentZvertGA(0), fHistCentZvertJE(0), fHistCentZvertMB(0), fHistCentZvertAny(0), |
fc392675 |
206 | fHistTPCdEdX(0), fHistITSsignal(0), //fHistTOFsignal(0), |
207 | fHistRhovsCent(0), fHistNjetvsCent(0), fHistCentrality(0), |
208 | fHistZvtx(0), fHistMult(0), |
209 | fHistJetPhi(0), fHistTrackPhi(0), |
d60e724e |
210 | fHistLocalRhoJetpt(0), |
fc392675 |
211 | fHistJetHaddPhiIN(0), fHistJetHaddPhiOUT(0), fHistJetHaddPhiMID(0), |
212 | fHistJetHaddPhiBias(0), fHistJetHaddPhiINBias(0), fHistJetHaddPhiOUTBias(0), fHistJetHaddPhiMIDBias(0), |
213 | fHistMEdPHI(0), |
214 | fHistTrackPtallcent(0), |
215 | fHistJetEtaPhi(0), fHistJetHEtaPhi(0), |
216 | fHistSEphieta(0), fHistMEphieta(0), |
217 | fHistJetHaddPHI(0), |
218 | fHistPID(0), |
daa728b0 |
219 | fhnPID(0x0), fhnMixedEvents(0x0), fhnJH(0x0), fhnCorr(0x0), |
8ee2d316 |
220 | fJetsCont(0), fTracksCont(0), fCaloClustersCont(0), |
221 | fContainerAllJets(0), fContainerPIDJets(1) |
fc392675 |
222 | { |
223 | // Default Constructor |
224 | for(Int_t ilab=0; ilab<4; ilab++){ |
225 | for(Int_t ipta=0; ipta<7; ipta++){ |
72a4d320 |
226 | //fHistTrackEtaPhi[ilab][ipta]=0; //keep out for now |
fc392675 |
227 | } // end of pt-associated loop |
228 | } // end of lab loop |
229 | |
230 | for(Int_t itrackpt=0; itrackpt<9; itrackpt++){ |
231 | fHistJetHadbindPhi[itrackpt]=0; |
232 | fHistJetHadbindPhiIN[itrackpt]=0; |
233 | fHistJetHadbindPhiMID[itrackpt]=0; |
234 | fHistJetHadbindPhiOUT[itrackpt]=0; |
235 | } // end of trackpt bin loop |
fc392675 |
236 | |
daa728b0 |
237 | for(Int_t icent = 0; icent<6; ++icent){ |
fc392675 |
238 | for(Int_t iptjet = 0; iptjet<5; ++iptjet){ |
239 | for(Int_t ieta = 0; ieta<3; ++ieta){ |
240 | fHistJetH[icent][iptjet][ieta]=0; |
241 | fHistJetHBias[icent][iptjet][ieta]=0; |
242 | fHistJetHTT[icent][iptjet][ieta]=0; |
243 | } // end of eta loop |
244 | } // end of pt-jet loop |
245 | } // end of centrality loop |
246 | |
247 | // centrality dependent histo's |
248 | for (Int_t i = 0;i<6;++i){ |
daa728b0 |
249 | fHistJetPt[i] = 0; |
250 | fHistJetPtBias[i] = 0; |
251 | fHistJetPtTT[i] = 0; |
252 | fHistAreavsRawPt[i] = 0; |
fc392675 |
253 | fHistJetPtvsTrackPt[i] = 0; |
254 | fHistRawJetPtvsTrackPt[i] = 0; |
255 | fHistTrackPt[i] = 0; |
256 | fHistEP0[i] = 0; |
257 | fHistEP0A[i] = 0; |
258 | fHistEP0C[i] = 0; |
259 | fHistEPAvsC[i] = 0; |
8ee2d316 |
260 | fHistJetPtcorrGlRho[i] = 0; |
fc392675 |
261 | fHistJetPtvsdEP[i] = 0; |
262 | fHistJetPtvsdEPBias[i] = 0; |
263 | fHistRhovsdEP[i] = 0; |
264 | fHistJetEtaPhiPt[i] = 0; |
265 | fHistJetEtaPhiPtBias[i] = 0; |
266 | fHistJetPtArea[i] = 0; |
267 | fHistJetPtAreaBias[i] = 0; |
268 | fHistJetPtNcon[i] = 0; |
269 | fHistJetPtNconBias[i] = 0; |
270 | fHistJetPtNconCh[i] = 0; |
271 | fHistJetPtNconBiasCh[i] = 0; |
272 | fHistJetPtNconEm[i] = 0; |
273 | fHistJetPtNconBiasEm[i] = 0; |
274 | fHistJetHaddPhiINcent[i] = 0; |
275 | fHistJetHaddPhiOUTcent[i] = 0; |
276 | fHistJetHaddPhiMIDcent[i] = 0; |
277 | } |
278 | |
279 | SetMakeGeneralHistograms(kTRUE); |
280 | |
fc392675 |
281 | } |
282 | |
283 | //_______________________________________________________________________ |
284 | AliAnalysisTaskEmcalJetHadEPpid::~AliAnalysisTaskEmcalJetHadEPpid() |
285 | { |
286 | // destructor |
287 | if (fOutput) { |
288 | delete fOutput; |
289 | fOutput = 0; |
290 | } |
291 | } |
292 | |
293 | //________________________________________________________________________ |
294 | void AliAnalysisTaskEmcalJetHadEPpid::UserCreateOutputObjects() |
295 | { // This is called ONCE! |
296 | if (!fCreateHisto) return; |
297 | AliAnalysisTaskEmcalJet::UserCreateOutputObjects(); |
298 | OpenFile(1); // do I need the 1? |
299 | |
8ee2d316 |
300 | // char array for naming histograms |
301 | int nchar = 200; |
302 | char name[nchar]; |
303 | |
daa728b0 |
304 | // strings for titles |
305 | TString name1; |
306 | TString title1; |
307 | |
fc392675 |
308 | // create histograms that arn't array |
daa728b0 |
309 | fHistNjetvsCent = new TH2F("NjetvsCent", "NjetvsCent", 100, 0.0, 100.0, 100, 0, 100); |
310 | fHistJetHaddPHI = new TH1F("fHistJetHaddPHI", "Jet-Hadron #Delta#varphi", 128,-0.5*TMath::Pi(),1.5*TMath::Pi()); |
311 | fHistJetHaddPhiIN = new TH1F("fHistJetHaddPhiIN","Jet-Hadron #Delta#varphi IN PLANE", 128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
312 | fHistJetHaddPhiOUT = new TH1F("fHistJetHaddPhiOUT","Jet-Hadron #Delta#varphi OUT PLANE",128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
313 | fHistJetHaddPhiMID = new TH1F("fHistJetHaddPhiMID","Jet-Hadron #Delta#varphi MIDDLE of PLANE",128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
d60e724e |
314 | fHistLocalRhoJetpt = new TH1F("fHistLocalRhoJetpt","Local Rho corrected Jet p_{T}", 50, -50, 200); |
315 | |
a409c5ba |
316 | // Centrality and Zvertex distribution for various triggers - Event Mixing QA |
317 | fHistCentZvertGA = new TH2F("fHistCentZvertGA", "Centrality - Z-vertex distribution for GA trigger", 20, 0, 100, 10, -10, 10); |
318 | fOutput->Add(fHistCentZvertGA); |
319 | fHistCentZvertJE = new TH2F("fHistCentZvertJE", "Centrality - Z-vertex distribution for JE trigger", 20, 0, 100, 10, -10, 10); |
320 | fOutput->Add(fHistCentZvertJE); |
321 | fHistCentZvertMB = new TH2F("fHistCentZvertMB", "Centrality - Z-vertex distribution for MB trigger", 20, 0, 100, 10, -10, 10); |
322 | fOutput->Add(fHistCentZvertMB); |
323 | fHistCentZvertAny = new TH2F("fHistCentZvertAny", "Centrality - Z-vertex distribution for kAny trigger", 20, 0, 100, 10, -10, 10); |
324 | fOutput->Add(fHistCentZvertAny); |
325 | |
d60e724e |
326 | // Event QA histo |
327 | fHistEventQA = new TH1F("fHistEventQA", "Event Counter at checkpoints in code", 20, 0.5, 20.5); |
f9b4eee4 |
328 | SetfHistQAcounterLabels(fHistEventQA); |
329 | fOutput->Add(fHistEventQA); |
d60e724e |
330 | |
331 | // Event Selection QA histo |
332 | fHistEventSelectionQA = new TH1F("fHistEventSelectionQA", "Trigger Selection Counter", 20, 0.5, 20.5); |
333 | SetfHistEvtSelQALabels(fHistEventSelectionQA); |
334 | fOutput->Add(fHistEventSelectionQA); |
335 | |
daa728b0 |
336 | // add to output lists |
337 | fOutput->Add(fHistNjetvsCent); |
338 | fOutput->Add(fHistJetHaddPHI); |
339 | fOutput->Add(fHistJetHaddPhiIN); |
340 | fOutput->Add(fHistJetHaddPhiOUT); |
341 | fOutput->Add(fHistJetHaddPhiMID); |
d60e724e |
342 | fOutput->Add(fHistLocalRhoJetpt); |
daa728b0 |
343 | |
a4248490 |
344 | fHistTPCdEdX = new TH2F("TPCdEdX", "TPCdEdX", 400, 0.0, 20.0, 500, 0, 500); |
345 | fOutput->Add(fHistTPCdEdX); |
346 | |
daa728b0 |
347 | // create histo's used for general QA |
8ee2d316 |
348 | if (makeQAhistos) { |
a4248490 |
349 | //fHistTPCdEdX = new TH2F("TPCdEdX", "TPCdEdX", 2000, 0.0, 100.0, 500, 0, 500); |
8ee2d316 |
350 | fHistITSsignal = new TH2F("ITSsignal", "ITSsignal", 2000, 0.0, 100.0, 500, 0, 500); |
daa728b0 |
351 | // fHistTOFsignal = new TH2F("TOFsignal", "TOFsignal", 2000, 0.0, 100.0, 500, 0, 500); |
8ee2d316 |
352 | fHistCentrality = new TH1F("fHistCentrality","centrality",100,0,100); |
353 | fHistZvtx = new TH1F("fHistZvertex","z vertex",60,-30,30); |
354 | fHistJetPhi = new TH1F("fHistJetPhi", "Jet #phi Distribution", 128, -2.0*TMath::Pi(), 2.0*TMath::Pi()); |
355 | fHistTrackPhi = new TH1F("fHistTrackPhi", "Track #phi Distribution", 128, -2.0*TMath::Pi(), 2.0*TMath::Pi()); |
356 | fHistRhovsCent = new TH2F("RhovsCent", "RhovsCent", 100, 0.0, 100.0, 500, 0, 500); |
daa728b0 |
357 | fHistTrackPtallcent = new TH1F("fHistTrackPtallcent", "p_{T} distribution", 1000, 0.0, 100.0); |
358 | fHistJetEtaPhi = new TH2F("fHistJetEtaPhi","Jet #eta-#phi",512,-1.8,1.8,512,-3.2,3.2); |
359 | fHistJetHEtaPhi = new TH2F("fHistJetHEtaPhi","Jet-Hadron #Delta#eta-#Delta#phi", 72, -1.8, 1.8, 72, -1.6, 4.8); |
360 | fHistSEphieta = new TH2F("fHistSEphieta", "Single Event #phi-#eta distribution", 64,-0.5*TMath::Pi(), 1.5*TMath::Pi(), 64,-1.8,1.8); // was 64 bins |
361 | fHistMEphieta = new TH2F("fHistMEphieta", "Mixed Event #phi-#eta distribution", 64, -0.5*TMath::Pi(), 1.5*TMath::Pi(), 64,-1.8,1.8); // was 64 bins |
8ee2d316 |
362 | |
363 | // add to output list |
a4248490 |
364 | //fOutput->Add(fHistTPCdEdX); |
8ee2d316 |
365 | fOutput->Add(fHistITSsignal); |
366 | //fOutput->Add(fHistTOFsignal); |
367 | fOutput->Add(fHistCentrality); |
368 | fOutput->Add(fHistZvtx); |
369 | fOutput->Add(fHistJetPhi); |
370 | fOutput->Add(fHistTrackPhi); |
371 | //fOutput->Add(fHistTrackEtaPhi); |
daa728b0 |
372 | fOutput->Add(fHistTrackPtallcent); |
373 | fOutput->Add(fHistJetEtaPhi); |
374 | fOutput->Add(fHistJetHEtaPhi); |
375 | fOutput->Add(fHistSEphieta); |
376 | fOutput->Add(fHistMEphieta); |
8ee2d316 |
377 | |
72a4d320 |
378 | //for(Int_t ipta=0; ipta<7; ++ipta){ |
379 | // for(Int_t ilab=0; ilab<4; ++ilab){ |
380 | // snprintf(name, nchar, "fHistTrackEtaPhi_%i_%i", ilab,ipta); |
381 | // fHistTrackEtaPhi[ilab][ipta] = new TH2F(name,name,400,-1,1,640,0.0,2.*TMath::Pi()); |
382 | // fOutput->Add(fHistTrackEtaPhi[ilab][ipta]); |
383 | // } // end of lab loop |
384 | //} // end of pt-associated loop |
8ee2d316 |
385 | |
daa728b0 |
386 | for (Int_t i = 0;i<6;++i){ |
8ee2d316 |
387 | name1 = TString(Form("EP0_%i",i)); |
388 | title1 = TString(Form("EP VZero cent bin %i",i)); |
389 | fHistEP0[i] = new TH1F(name1,title1,144,-TMath::Pi(),TMath::Pi()); |
390 | fOutput->Add(fHistEP0[i]); |
391 | |
392 | name1 = TString(Form("EP0A_%i",i)); |
393 | title1 = TString(Form("EP VZero cent bin %i",i)); |
394 | fHistEP0A[i] = new TH1F(name1,title1,144,-TMath::Pi(),TMath::Pi()); |
395 | fOutput->Add(fHistEP0A[i]); |
396 | |
397 | name1 = TString(Form("EP0C_%i",i)); |
398 | title1 = TString(Form("EP VZero cent bin %i",i)); |
399 | fHistEP0C[i] = new TH1F(name1,title1,144,-TMath::Pi(),TMath::Pi()); |
400 | fOutput->Add(fHistEP0C[i]); |
401 | |
402 | name1 = TString(Form("EPAvsC_%i",i)); |
403 | title1 = TString(Form("EP VZero cent bin %i",i)); |
404 | fHistEPAvsC[i] = new TH2F(name1,title1,144,-TMath::Pi(),TMath::Pi(),144,-TMath::Pi(),TMath::Pi()); |
405 | fOutput->Add(fHistEPAvsC[i]); |
406 | |
407 | name1 = TString(Form("JetPtvsTrackPt_%i",i)); |
408 | title1 = TString(Form("Jet p_{T} vs Leading Track p_{T} cent bin %i",i)); |
409 | fHistJetPtvsTrackPt[i] = new TH2F(name1,title1,250,-50,200,250,0,50); |
410 | fOutput->Add(fHistJetPtvsTrackPt[i]); |
411 | |
412 | name1 = TString(Form("RawJetPtvsTrackPt_%i",i)); |
413 | title1 = TString(Form("Raw Jet p_{T} vs Leading Track p_{T} cent bin %i",i)); |
414 | fHistRawJetPtvsTrackPt[i] = new TH2F(name1,title1,250,-50,200,250,0,50); |
415 | fOutput->Add(fHistRawJetPtvsTrackPt[i]); |
416 | |
417 | name1 = TString(Form("TrackPt_%i",i)); |
418 | title1 = TString(Form("Track p_{T} cent bin %i",i)); |
419 | fHistTrackPt[i] = new TH1F(name1,title1,1000,0,100); // up to 200? |
420 | fOutput->Add(fHistTrackPt[i]); |
421 | |
422 | name1 = TString(Form("JetPtcorrGLrho_%i",i)); |
423 | title1 = TString(Form("Jet p_{T} corrected with Global #rho cent bin %i",i)); |
424 | fHistJetPtcorrGlRho[i] = new TH1F(name1,title1,300,-100,200); // up to 200? |
425 | fOutput->Add(fHistJetPtcorrGlRho[i]); |
daa728b0 |
426 | |
427 | name1 = TString(Form("JetPtvsdEP_%i",i)); |
428 | title1 = TString(Form("Jet p_{T} vs #DeltaEP cent bin %i",i)); |
429 | fHistJetPtvsdEP[i] = new TH2F(name1,title1,250,-50,200,288,-2*TMath::Pi(),2*TMath::Pi()); |
430 | fOutput->Add(fHistJetPtvsdEP[i]); |
431 | |
432 | name1 = TString(Form("RhovsdEP_%i",i)); |
433 | title1 = TString(Form("#rho vs #DeltaEP cent bin %i",i)); |
434 | fHistRhovsdEP[i] = new TH2F(name1,title1,500,0,500,288,-2*TMath::Pi(),2*TMath::Pi()); |
435 | fOutput->Add(fHistRhovsdEP[i]); |
436 | |
437 | name1 = TString(Form("JetEtaPhiPt_%i",i)); |
438 | title1 = TString(Form("Jet #eta-#phi p_{T} cent bin %i",i)); |
439 | fHistJetEtaPhiPt[i] = new TH3F(name1,title1,250,-50,200,100,-1,1,64,-3.2,3.2); |
440 | fOutput->Add(fHistJetEtaPhiPt[i]); |
441 | |
442 | name1 = TString(Form("JetPtArea_%i",i)); |
443 | title1 = TString(Form("Jet p_{T} Area cent bin %i",i)); |
444 | fHistJetPtArea[i] = new TH2F(name1,title1,250,-50,200,100,0,1); |
445 | fOutput->Add(fHistJetPtArea[i]); |
446 | |
447 | snprintf(name, nchar, "fHistAreavsRawPt_%i",i); |
448 | fHistAreavsRawPt[i] = new TH2F(name,name,100,0,1,200,0,200); |
449 | fOutput->Add(fHistAreavsRawPt[i]); |
450 | } // loop over centrality |
fc392675 |
451 | |
daa728b0 |
452 | } // QA histo switch |
453 | |
454 | if (makeBIAShistos) { |
455 | fHistJetHaddPhiBias = new TH1F("fHistJetHaddPhiBias","Jet-Hadron #Delta#varphi with bias",128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
456 | fHistJetHaddPhiINBias = new TH1F("fHistJetHaddPhiINBias","Jet-Hadron #Delta#varphi IN PLANE with bias", 128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
457 | fHistJetHaddPhiOUTBias = new TH1F("fHistJetHaddPhiOUTBias","Jet-Hadron #Delta#varphi OUT PLANE with bias",128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
458 | fHistJetHaddPhiMIDBias = new TH1F("fHistJetHaddPhiMIDBias","Jet-Hadron #Delta#varphi MIDDLE of PLANE with bias",128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
459 | |
460 | // add to output list |
461 | fOutput->Add(fHistJetHaddPhiBias); |
462 | fOutput->Add(fHistJetHaddPhiINBias); |
463 | fOutput->Add(fHistJetHaddPhiOUTBias); |
464 | fOutput->Add(fHistJetHaddPhiMIDBias); |
465 | |
466 | for (Int_t i = 0;i<6;++i){ |
467 | name1 = TString(Form("JetPtvsdEPBias_%i",i)); |
468 | title1 = TString(Form("Bias Jet p_{T} vs #DeltaEP cent bin %i",i)); |
469 | fHistJetPtvsdEPBias[i] = new TH2F(name1,title1,250,-50,200,288,-2*TMath::Pi(),2*TMath::Pi()); |
470 | fOutput->Add(fHistJetPtvsdEPBias[i]); |
fc392675 |
471 | |
daa728b0 |
472 | name1 = TString(Form("JetEtaPhiPtBias_%i",i)); |
473 | title1 = TString(Form("Jet #eta-#phi p_{T} Bias cent bin %i",i)); |
474 | fHistJetEtaPhiPtBias[i] = new TH3F(name1,title1,250,-50,200,100,-1,1,64,-3.2,3.2); |
475 | fOutput->Add(fHistJetEtaPhiPtBias[i]); |
fc392675 |
476 | |
8ee2d316 |
477 | name1 = TString(Form("JetPtAreaBias_%i",i)); |
478 | title1 = TString(Form("Jet p_{T} Area Bias cent bin %i",i)); |
479 | fHistJetPtAreaBias[i] = new TH2F(name1,title1,250,-50,200,100,0,1); |
480 | fOutput->Add(fHistJetPtAreaBias[i]); |
daa728b0 |
481 | } // end of centrality loop |
482 | } // bias histos |
483 | |
484 | if (makeoldJEThadhistos) { |
485 | for(Int_t icent = 0; icent<6; ++icent){ |
486 | snprintf(name, nchar, "fHistJetPtTT_%i",icent); |
487 | fHistJetPtTT[icent] = new TH1F(name,name,200,0,200); |
488 | fOutput->Add(fHistJetPtTT[icent]); |
489 | |
490 | snprintf(name, nchar, "fHistJetPt_%i",icent); |
491 | fHistJetPt[icent] = new TH1F(name,name,200,0,200); |
492 | fOutput->Add(fHistJetPt[icent]); |
493 | |
494 | snprintf(name, nchar, "fHistJetPtBias_%i",icent); |
495 | fHistJetPtBias[icent] = new TH1F(name,name,200,0,200); |
496 | fOutput->Add(fHistJetPtBias[icent]); |
497 | |
498 | for(Int_t iptjet = 0; iptjet<5; ++iptjet){ |
499 | for(Int_t ieta = 0; ieta<3; ++ieta){ |
500 | snprintf(name, nchar, "fHistJetH_%i_%i_%i",icent,iptjet,ieta); |
501 | fHistJetH[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30); |
502 | fOutput->Add(fHistJetH[icent][iptjet][ieta]); |
503 | |
504 | snprintf(name, nchar, "fHistJetHBias_%i_%i_%i",icent,iptjet,ieta); |
505 | fHistJetHBias[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30); |
506 | fOutput->Add(fHistJetHBias[icent][iptjet][ieta]); |
507 | |
508 | snprintf(name, nchar, "fHistJetHTT_%i_%i_%i",icent,iptjet,ieta); |
509 | fHistJetHTT[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30); |
510 | fOutput->Add(fHistJetHTT[icent][iptjet][ieta]); |
511 | } // end of eta loop |
512 | } // end of pt-jet loop |
513 | } // end of centrality loop |
514 | } // make JetHadhisto old |
515 | |
516 | if (makeextraCORRhistos) { |
517 | for(Int_t itrackpt=0; itrackpt<9; itrackpt++){ |
518 | snprintf(name, nchar, "fHistJetHadbindPhi_%i",itrackpt); |
519 | fHistJetHadbindPhi[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
520 | fOutput->Add(fHistJetHadbindPhi[itrackpt]); |
521 | |
522 | snprintf(name, nchar, "fHistJetHadbindPhiIN_%i",itrackpt); |
523 | fHistJetHadbindPhiIN[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
524 | fOutput->Add(fHistJetHadbindPhiIN[itrackpt]); |
525 | |
526 | snprintf(name, nchar, "fHistJetHadbindPhiMID_%i",itrackpt); |
527 | fHistJetHadbindPhiMID[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
528 | fOutput->Add(fHistJetHadbindPhiMID[itrackpt]); |
529 | |
530 | snprintf(name, nchar, "fHistJetHadbindPhiOUT_%i",itrackpt); |
531 | fHistJetHadbindPhiOUT[itrackpt] = new TH1F(name,name,128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
532 | fOutput->Add(fHistJetHadbindPhiOUT[itrackpt]); |
533 | } // end of trackpt bin loop |
534 | |
535 | for (Int_t i = 0;i<6;++i){ |
536 | name1 = TString(Form("JetHaddPhiINcent_%i",i)); |
537 | title1 = TString(Form("Jet Hadron #Delta#varphi Distribution IN PLANE cent bin %i",i)); |
538 | fHistJetHaddPhiINcent[i] = new TH1F(name1,title1,128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
539 | fOutput->Add(fHistJetHaddPhiINcent[i]); |
540 | |
541 | name1 = TString(Form("JetHaddPhiOUTcent_%i",i)); |
542 | title1 = TString(Form("Jet Hadron #Delta#varphi Distribution OUT PLANE cent bin %i",i)); |
543 | fHistJetHaddPhiOUTcent[i] = new TH1F(name1,title1,128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
544 | fOutput->Add(fHistJetHaddPhiOUTcent[i]); |
545 | |
546 | name1 = TString(Form("JetHaddPhiMIDcent_%i",i)); |
547 | title1 = TString(Form("Jet Hadron #Delta#varphi Distribution MIDDLE of PLANE cent bin %i",i)); |
548 | fHistJetHaddPhiMIDcent[i] = new TH1F(name1,title1,128,-0.5*TMath::Pi(), 1.5*TMath::Pi()); |
549 | fOutput->Add(fHistJetHaddPhiMIDcent[i]); |
8ee2d316 |
550 | |
551 | name1 = TString(Form("JetPtNcon_%i",i)); |
552 | title1 = TString(Form("Jet p_{T} Ncon cent bin %i",i)); |
553 | fHistJetPtNcon[i] = new TH2F(name1,title1,250,-50,200,100,0,2000); |
554 | fOutput->Add(fHistJetPtNcon[i]); |
555 | |
556 | name1 = TString(Form("JetPtNconBias_%i",i)); |
557 | title1 = TString(Form("Jet p_{T} NconBias cent bin %i",i)); |
558 | fHistJetPtNconBias[i] = new TH2F(name1,title1,250,-50,200,100,0,2000); |
559 | fOutput->Add(fHistJetPtNconBias[i]); |
560 | |
561 | name1 = TString(Form("JetPtNconCh_%i",i)); |
562 | title1 = TString(Form("Jet p_{T} NconCh cent bin %i",i)); |
563 | fHistJetPtNconCh[i] = new TH2F(name1,title1,250,-50,200,100,0,2000); |
564 | fOutput->Add(fHistJetPtNconCh[i]); |
565 | |
566 | name1 = TString(Form("JetPtNconBiasCh_%i",i)); |
567 | title1 = TString(Form("Jet p_{T} NconBiasCh cent bin %i",i)); |
568 | fHistJetPtNconBiasCh[i] = new TH2F(name1,title1,250,-50,200,100,0,2000); |
569 | fOutput->Add(fHistJetPtNconBiasCh[i]); |
570 | |
571 | name1 = TString(Form("JetPtNconEm_%i",i)); |
572 | title1 = TString(Form("Jet p_{T} NconEm cent bin %i",i)); |
573 | fHistJetPtNconEm[i] = new TH2F(name1,title1,250,-50,200,100,0,2000); |
574 | fOutput->Add(fHistJetPtNconEm[i]); |
575 | |
576 | name1 = TString(Form("JetPtNconBiasEm_%i",i)); |
577 | title1 = TString(Form("Jet p_{T} NconBiasEm cent bin %i",i)); |
578 | fHistJetPtNconBiasEm[i] = new TH2F(name1,title1,250,-50,200,100,0,2000); |
579 | fOutput->Add(fHistJetPtNconBiasEm[i]); |
580 | } // extra Correlations histos switch |
581 | } |
3263f852 |
582 | |
583 | // variable binned pt for THnSparse's |
e56c3f6e |
584 | Double_t xlowjetPT[] = {15, 16, 18, 20, 22, 24, 26, 28, 30, 32, 34, 36, 38, 40, 45, 50, 55, 60, 65, 70, 75, 80, 100, 150, 200, 300}; |
585 | Double_t xlowtrPT[] = {0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0,4.2,4.4,4.6,4.8,5.0,5.5,6.0,6.5,7.0,7.5,8.0,9.0,10.0,12.0,14.0,16.0,18.0,20.0,25.0,30.0,40.0,50.0,75.0}; |
3263f852 |
586 | |
e56c3f6e |
587 | // tracks: 51, jets: 26 |
3263f852 |
588 | // number of bins you tell histogram should be (# in array - 1) because the last bin |
589 | // is the right-most edge of the histogram |
590 | // i.e. this is for PT and there are 57 numbers (bins) thus we have 56 bins in our histo |
591 | Int_t nbinsjetPT = sizeof(xlowjetPT)/sizeof(Double_t) - 1; |
592 | Int_t nbinstrPT = sizeof(xlowtrPT)/sizeof(Double_t) - 1; |
8ee2d316 |
593 | |
fc392675 |
594 | // set up jet-hadron sparse |
3f135c41 |
595 | UInt_t bitcodeMESE = 0; // bit coded, see GetDimParams() below |
596 | UInt_t bitcodePID = 0; // bit coded, see GetDimParamsPID() below |
daa728b0 |
597 | UInt_t bitcodeCorr = 0; // bit coded, see GetDimparamsCorr() below |
3f135c41 |
598 | bitcodeMESE = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6; // | 1<<7 | 1<<8 | 1<<9; |
599 | if(fDoEventMixing) { |
600 | fhnJH = NewTHnSparseF("fhnJH", bitcodeMESE); |
3263f852 |
601 | |
602 | if(dovarbinTHnSparse){ |
603 | fhnJH->GetAxis(1)->Set(nbinsjetPT, xlowjetPT); |
604 | fhnJH->GetAxis(2)->Set(nbinstrPT, xlowtrPT); |
605 | } |
606 | |
607 | fOutput->Add(fhnJH); |
3f135c41 |
608 | } |
fc392675 |
609 | |
f9b4eee4 |
610 | bitcodeCorr = 1<<0 | 1<<1 | 1<<2 | 1<<3; // | 1<<4 | 1<<5; |
e56c3f6e |
611 | fhnCorr = NewTHnSparseFCorr("fhnCorr", bitcodeCorr); |
3263f852 |
612 | if(dovarbinTHnSparse) fhnCorr->GetAxis(1)->Set(nbinsjetPT, xlowjetPT); |
3badfe33 |
613 | fOutput->Add(fhnCorr); |
f1ce4cc3 |
614 | |
3badfe33 |
615 | /* |
616 | Double_t centralityBins[nCentralityBins+1]; |
617 | for(Int_t ic=0; ic<nCentralityBins+1; ic++){ |
618 | if(ic==nCentralityBins) centralityBins[ic]=500; |
619 | else centralityBins[ic]=10.0*ic; |
620 | } |
621 | */ |
622 | |
623 | // set up centrality bins for mixed events |
624 | // for pp we need mult bins for event mixing. Create binning here, to also make a histogram from it |
56523939 |
625 | Int_t nCentralityBinspp = 8; |
3badfe33 |
626 | //Double_t centralityBinspp[nCentralityBinspp+1]; |
56523939 |
627 | Double_t centralityBinspp[9] = {0.0, 4., 9, 15, 25, 35, 55, 100.0, 500.0}; |
3badfe33 |
628 | |
629 | // Setup for Pb-Pb collisions |
a4248490 |
630 | Int_t nCentralityBinsPbPb = 100; |
631 | Double_t mult = 1.0; |
632 | if(fCentBinSize==1) { |
633 | nCentralityBinsPbPb = 100; |
634 | mult = 1.0; |
635 | } else if(fCentBinSize==2){ |
636 | nCentralityBinsPbPb = 50; |
637 | mult = 2.0; |
638 | } else if(fCentBinSize==5){ |
639 | nCentralityBinsPbPb = 20; |
640 | mult = 5.0; |
641 | } else if(fCentBinSize==10){ |
642 | nCentralityBinsPbPb = 10; |
643 | mult = 10.0; |
644 | } |
645 | Double_t centralityBinsPbPb[nCentralityBinsPbPb]; // nCentralityBinsPbPb |
646 | for(Int_t ic=0; ic<nCentralityBinsPbPb; ic++){ |
647 | centralityBinsPbPb[ic]=mult*ic; |
648 | } |
649 | |
650 | /* |
a409c5ba |
651 | Int_t nCentralityBinsPbPb = 10; //100; |
3badfe33 |
652 | Double_t centralityBinsPbPb[nCentralityBinsPbPb+1]; |
653 | for(Int_t ic=0; ic<nCentralityBinsPbPb; ic++){ |
a409c5ba |
654 | centralityBinsPbPb[ic]=10.0*ic; //1.0*ic; |
fc392675 |
655 | } |
a4248490 |
656 | */ |
fc392675 |
657 | |
56523939 |
658 | if(fBeam == 0) fHistMult = new TH1F("fHistMult","multiplicity",nCentralityBinspp,centralityBinspp); |
659 | if(fBeam == 1) fHistMult = new TH1F("fHistMult","multiplicity",nCentralityBinsPbPb,centralityBinsPbPb); |
8ee2d316 |
660 | // fOutput->Add(fHistMult); |
fc392675 |
661 | |
662 | // Event Mixing |
663 | Int_t trackDepth = fMixingTracks; |
664 | Int_t poolsize = 1000; // Maximum number of events, ignored in the present implemented of AliEventPoolManager |
665 | Int_t nZvtxBins = 5+1+5; |
666 | Double_t vertexBins[] = { -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10}; |
667 | Double_t* zvtxbin = vertexBins; |
56523939 |
668 | if(fBeam == 0) fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBinspp, centralityBinspp, nZvtxBins, zvtxbin); |
669 | if(fBeam == 1) fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBinsPbPb, centralityBinsPbPb, nZvtxBins, zvtxbin); |
fc392675 |
670 | |
671 | // set up event mixing sparse |
672 | if(fDoEventMixing){ |
3f135c41 |
673 | bitcodeMESE = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6; // | 1<<7 | 1<<8 | 1<<9; |
674 | fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", bitcodeMESE); |
fc392675 |
675 | |
3f135c41 |
676 | // set up some variable binning of the sparse |
3263f852 |
677 | if(dovarbinTHnSparse){ |
678 | fhnMixedEvents->GetAxis(1)->Set(nbinsjetPT, xlowjetPT); |
679 | fhnMixedEvents->GetAxis(2)->Set(nbinstrPT, xlowtrPT); |
680 | } |
fc392675 |
681 | |
3263f852 |
682 | fOutput->Add(fhnMixedEvents); |
683 | } // end of do-eventmixing |
fc392675 |
684 | |
685 | // set up PID sparse |
686 | if(doPID){ |
8ee2d316 |
687 | // ****************************** PID ***************************************************** |
688 | // set up PID handler |
689 | AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager(); |
690 | AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler()); |
691 | if(!inputHandler) { |
692 | AliFatal("Input handler needed"); |
f1ce4cc3 |
693 | return; |
8ee2d316 |
694 | } |
695 | |
696 | // PID response object |
697 | fPIDResponse = inputHandler->GetPIDResponse(); |
698 | if (!fPIDResponse) { |
699 | AliError("PIDResponse object was not created"); |
f1ce4cc3 |
700 | return; |
8ee2d316 |
701 | } |
702 | // ***************************************************************************************** |
703 | |
704 | // PID counter |
d60e724e |
705 | fHistPID = new TH1F("fHistPID","PID Counter", 15, 0.5, 15.5); |
8ee2d316 |
706 | SetfHistPIDcounterLabels(fHistPID); |
707 | fOutput->Add(fHistPID); |
708 | |
139e13f4 |
709 | if(allpidAXIS) { |
3f135c41 |
710 | bitcodePID = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | |
a4248490 |
711 | 1<<10 | 1<<11 | 1<<12 | 1<<13 | 1<<14 | 1<<15 | 1<<16 | 1<<17 | 1<<18; // | 1<<19 | 1<<20; |
3f135c41 |
712 | fhnPID = NewTHnSparseFPID("fhnPID", bitcodePID); |
139e13f4 |
713 | } else { |
3f135c41 |
714 | bitcodePID = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | |
a4248490 |
715 | 1<<10 | 1<<11; // | 1<<12 | 1<<13; |
3f135c41 |
716 | fhnPID = NewTHnSparseFPID("fhnPID", bitcodePID); |
139e13f4 |
717 | } |
8ee2d316 |
718 | |
3f135c41 |
719 | // set some variable binning of sparse |
8ee2d316 |
720 | if(dovarbinTHnSparse){ |
139e13f4 |
721 | fhnPID->GetAxis(1)->Set(nbinstrPT, xlowtrPT); |
3263f852 |
722 | fhnPID->GetAxis(8)->Set(nbinsjetPT, xlowjetPT); |
8ee2d316 |
723 | } |
724 | |
fc392675 |
725 | fOutput->Add(fhnPID); |
726 | } // end of do-PID |
727 | |
fc392675 |
728 | // =========== Switch on Sumw2 for all histos =========== |
729 | for (Int_t i=0; i<fOutput->GetEntries(); ++i) { |
730 | TH1 *h1 = dynamic_cast<TH1*>(fOutput->At(i)); |
731 | if (h1){ |
732 | h1->Sumw2(); |
733 | continue; |
734 | } |
735 | TH2 *h2 = dynamic_cast<TH2*>(fOutput->At(i)); |
736 | if (h2){ |
737 | h2->Sumw2(); |
738 | continue; |
739 | } |
740 | TH3 *h3 = dynamic_cast<TH3*>(fOutput->At(i)); |
741 | if (h3){ |
742 | h3->Sumw2(); |
743 | continue; |
744 | } |
745 | THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i)); |
746 | if(hn)hn->Sumw2(); |
747 | } |
748 | |
749 | PostData(1, fOutput); |
750 | } |
751 | |
752 | //_________________________________________________________________________ |
753 | void AliAnalysisTaskEmcalJetHadEPpid::ExecOnce() |
754 | { |
755 | AliAnalysisTaskEmcalJet::ExecOnce(); |
8ee2d316 |
756 | |
757 | if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0; |
758 | if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0; |
759 | if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0; |
fc392675 |
760 | } |
761 | |
762 | //_________________________________________________________________________ |
763 | Bool_t AliAnalysisTaskEmcalJetHadEPpid::Run() |
764 | { // Main loop called for each event |
765 | // TEST TEST TEST TEST for OBJECTS! |
3badfe33 |
766 | |
f9b4eee4 |
767 | fHistEventQA->Fill(1); // All Events that get entered |
768 | |
3f135c41 |
769 | // check and fill a Event Selection QA histogram for different trigger selections |
d60e724e |
770 | UInt_t trig = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected(); |
771 | if(trig == 0) fHistEventSelectionQA->Fill(1); |
772 | if(trig & AliVEvent::kAny) fHistEventSelectionQA->Fill(2); |
773 | if(trig & AliVEvent::kAnyINT) fHistEventSelectionQA->Fill(3); |
774 | if(trig & AliVEvent::kMB) fHistEventSelectionQA->Fill(4); |
775 | if(trig & AliVEvent::kINT7) fHistEventSelectionQA->Fill(5); |
776 | if(trig & AliVEvent::kEMC1) fHistEventSelectionQA->Fill(6); |
777 | if(trig & AliVEvent::kEMC7) fHistEventSelectionQA->Fill(7); |
778 | if(trig & AliVEvent::kEMC8) fHistEventSelectionQA->Fill(8); |
779 | if(trig & AliVEvent::kEMCEJE) fHistEventSelectionQA->Fill(9); |
780 | if(trig & AliVEvent::kEMCEGA) fHistEventSelectionQA->Fill(10); |
781 | if(trig & AliVEvent::kCentral) fHistEventSelectionQA->Fill(11); |
782 | if(trig & AliVEvent::kSemiCentral) fHistEventSelectionQA->Fill(12); |
783 | if(trig & AliVEvent::kINT8) fHistEventSelectionQA->Fill(13); |
784 | |
785 | if(trig & (AliVEvent::kEMCEJE | AliVEvent::kMB)) fHistEventSelectionQA->Fill(14); |
786 | if(trig & (AliVEvent::kEMCEGA | AliVEvent::kMB)) fHistEventSelectionQA->Fill(15); |
787 | if(trig & (AliVEvent::kAnyINT | AliVEvent::kMB)) fHistEventSelectionQA->Fill(16); |
788 | |
3f135c41 |
789 | //if(trig & (AliVEvent::kEMCEJE && AliVEvent::kMB)) fHistEventSelectionQA->Fill(17); |
790 | //if(trig & (AliVEvent::kEMCEGA && AliVEvent::kMB)) fHistEventSelectionQA->Fill(18); |
791 | //if(trig & (AliVEvent::kAnyINT && AliVEvent::kMB)) fHistEventSelectionQA->Fill(19); |
d60e724e |
792 | |
3f135c41 |
793 | // see if we are running on PbPb and try and get LocalRho object |
d60e724e |
794 | if(GetBeamType() == 1) { |
795 | if(!fLocalRho){ |
796 | AliError(Form("Couldn't get fLocalRho object, try to get it from Event based on name\n")); |
797 | fLocalRho = GetLocalRhoFromEvent(fLocalRhoName); |
798 | if(!fLocalRho) return kTRUE; |
799 | } |
800 | } // check for LocalRho object if PbPb data |
801 | |
fc392675 |
802 | if(!fTracks){ |
803 | AliError(Form("No fTracks object!!\n")); |
804 | return kTRUE; |
805 | } |
806 | if(!fJets){ |
807 | AliError(Form("No fJets object!!\n")); |
808 | return kTRUE; |
809 | } |
810 | |
f9b4eee4 |
811 | fHistEventQA->Fill(2); // events after object check |
812 | |
fc392675 |
813 | // what kind of event do we have: AOD or ESD? |
139e13f4 |
814 | Bool_t useAOD; |
815 | if (dynamic_cast<AliAODEvent*>(InputEvent())) useAOD = kTRUE; |
816 | else useAOD = kFALSE; |
817 | |
fc392675 |
818 | // if we have ESD event, set up ESD object |
819 | if(!useAOD){ |
820 | fESD = dynamic_cast<AliESDEvent*>(InputEvent()); |
821 | if (!fESD) { |
822 | AliError(Form("ERROR: fESD not available\n")); |
823 | return kTRUE; |
824 | } |
825 | } |
826 | |
827 | // if we have AOD event, set up AOD object |
828 | if(useAOD){ |
829 | fAOD = dynamic_cast<AliAODEvent*>(InputEvent()); |
830 | if(!fAOD) { |
831 | AliError(Form("ERROR: fAOD not available\n")); |
832 | return kTRUE; |
833 | } |
834 | } |
835 | |
f9b4eee4 |
836 | fHistEventQA->Fill(3); // events after Aod/esd check |
837 | |
fc392675 |
838 | // get centrality |
839 | Int_t centbin = GetCentBin(fCent); |
8ee2d316 |
840 | if (makeQAhistos) fHistCentrality->Fill(fCent); // won't be filled in pp collision (Keep this in mind!) |
fc392675 |
841 | |
3badfe33 |
842 | // if we are on PbPb data do cut on centrality > 90%, else by default DON'T |
843 | if (GetBeamType() == 1) { |
844 | // apply cut to event on Centrality > 90% |
845 | if(fCent>90) return kTRUE; |
846 | } |
8ee2d316 |
847 | |
3f135c41 |
848 | // BEAM TYPE enumerator: kNA = -1, kpp = 0, kAA = 1, kpA = 2 |
849 | // for pp analyses we will just use the first centrality bin |
850 | if(GetBeamType() == 0) if (centbin == -1) centbin = 0; |
851 | if(GetBeamType() == 1) if (centbin == -1) return kTRUE; |
852 | |
f9b4eee4 |
853 | fHistEventQA->Fill(4); // events after centrality check |
854 | |
fc392675 |
855 | // get vertex information |
856 | Double_t fvertex[3]={0,0,0}; |
857 | InputEvent()->GetPrimaryVertex()->GetXYZ(fvertex); |
858 | Double_t zVtx=fvertex[2]; |
8ee2d316 |
859 | if (makeQAhistos) fHistZvtx->Fill(zVtx); |
fc392675 |
860 | |
861 | // get z-vertex bin |
862 | //Int_t zVbin = GetzVertexBin(zVtx); |
863 | |
864 | // apply zVtx cut |
8ee2d316 |
865 | if(fabs(zVtx)>10.0) return kTRUE; |
fc392675 |
866 | |
f9b4eee4 |
867 | fHistEventQA->Fill(5); // events after zvertex check |
868 | |
fc392675 |
869 | // create pointer to list of input event |
870 | TList *list = InputEvent()->GetList(); |
871 | if(!list) { |
872 | AliError(Form("ERROR: list not attached\n")); |
8ee2d316 |
873 | return kTRUE; |
fc392675 |
874 | } |
875 | |
f9b4eee4 |
876 | fHistEventQA->Fill(6); // events after list check |
877 | |
fc392675 |
878 | // initialize TClonesArray pointers to jets and tracks |
879 | TClonesArray *jets = 0; |
880 | TClonesArray *tracks = 0; |
3badfe33 |
881 | TClonesArray *tracksME = 0; |
fc392675 |
882 | |
883 | // get Tracks object |
884 | tracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracks)); |
885 | if (!tracks) { |
886 | AliError(Form("Pointer to tracks %s == 0", fTracks->GetName())); |
887 | return kTRUE; |
888 | } // verify existence of tracks |
889 | |
3badfe33 |
890 | // get ME Tracks object |
891 | tracksME = dynamic_cast<TClonesArray*>(list->FindObject(fTracksNameME)); |
892 | if (!tracksME) { |
893 | AliError(Form("Pointer to tracks %s == 0", fTracksNameME.Data())); |
894 | return kTRUE; |
895 | } // verify existence of tracks |
896 | |
fc392675 |
897 | // get Jets object |
898 | jets = dynamic_cast<TClonesArray*>(list->FindObject(fJets)); |
899 | if(!jets){ |
900 | AliError(Form("Pointer to jets %s == 0", fJets->GetName())); |
901 | return kTRUE; |
902 | } // verify existence of jets |
903 | |
f9b4eee4 |
904 | fHistEventQA->Fill(7); // events after track/jet pointer check |
905 | |
fc392675 |
906 | // get number of jets and tracks |
907 | const Int_t Njets = jets->GetEntries(); |
908 | const Int_t Ntracks = tracks->GetEntries(); |
909 | if(Ntracks<1) return kTRUE; |
910 | if(Njets<1) return kTRUE; |
911 | |
f9b4eee4 |
912 | fHistEventQA->Fill(8); // events after #track and jets < 1 check |
913 | |
8ee2d316 |
914 | if (makeQAhistos) fHistMult->Fill(Ntracks); // fill multiplicity distribution |
fc392675 |
915 | |
916 | // initialize track parameters |
917 | Int_t iTT=-1; |
918 | Double_t ptmax=-10; |
a409c5ba |
919 | Int_t NtrackAcc = 0; |
fc392675 |
920 | |
3badfe33 |
921 | fVevent = dynamic_cast<AliVEvent*>(InputEvent()); |
922 | if (!fVevent) { |
923 | printf("ERROR: fVEvent not available\n"); |
924 | return kTRUE; |
925 | } |
926 | |
a409c5ba |
927 | // fill event mixing QA |
928 | if(trig & AliVEvent::kEMCEGA) fHistCentZvertGA->Fill(fCent, zVtx); |
929 | if(trig & AliVEvent::kEMCEJE) fHistCentZvertJE->Fill(fCent, zVtx); |
930 | if(trig & AliVEvent::kMB) fHistCentZvertMB->Fill(fCent, zVtx); |
931 | if(trig & AliVEvent::kAny) fHistCentZvertAny->Fill(fCent, zVtx); |
3badfe33 |
932 | |
fc392675 |
933 | // loop over tracks - to get hardest track (highest pt) |
3f135c41 |
934 | for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++){ |
935 | AliVTrack *track = dynamic_cast<AliVTrack*>(tracks->At(iTracks)); |
936 | if(!track){ |
937 | AliError(Form("Couldn't get AliVtrack %d\n", iTracks)); |
e56c3f6e |
938 | continue; |
3263f852 |
939 | } |
3f135c41 |
940 | |
daa728b0 |
941 | // track cuts |
fc392675 |
942 | if(TMath::Abs(track->Eta())>0.9) continue; |
943 | if(track->Pt()<0.15) continue; |
944 | //iCount++; |
56523939 |
945 | NtrackAcc++; |
946 | |
fc392675 |
947 | if(track->Pt()>ptmax){ |
948 | ptmax=track->Pt(); // max pT track |
949 | iTT=iTracks; // trigger tracks |
950 | } // check if Pt>maxpt |
a4248490 |
951 | |
8ee2d316 |
952 | if (makeQAhistos) fHistTrackPhi->Fill(track->Phi()); |
953 | if (makeQAhistos) fHistTrackPt[centbin]->Fill(track->Pt()); |
daa728b0 |
954 | if (makeQAhistos) fHistTrackPtallcent->Fill(track->Pt()); |
fc392675 |
955 | } // end of loop over tracks |
956 | |
fc392675 |
957 | // get rho from event and fill relative histo's |
958 | fRho = GetRhoFromEvent(fRhoName); |
959 | fRhoVal = fRho->GetVal(); |
fc392675 |
960 | |
8ee2d316 |
961 | if (makeQAhistos) { |
daa728b0 |
962 | fHistRhovsdEP[centbin]->Fill(fRhoVal,fEPV0); // Global Rho vs delta Event Plane angle |
8ee2d316 |
963 | fHistRhovsCent->Fill(fCent,fRhoVal); // Global Rho vs Centrality |
964 | fHistEP0[centbin]->Fill(fEPV0); |
965 | fHistEP0A[centbin]->Fill(fEPV0A); |
966 | fHistEP0C[centbin]->Fill(fEPV0C); |
967 | fHistEPAvsC[centbin]->Fill(fEPV0A,fEPV0C); |
968 | } |
969 | |
fc392675 |
970 | // initialize jet parameters |
971 | Int_t ijethi=-1; |
972 | Double_t highestjetpt=0.0; |
973 | Int_t passedTTcut=0; |
974 | Int_t NjetAcc = 0; |
daa728b0 |
975 | Double_t leadhadronPT = 0; |
fc392675 |
976 | |
977 | // loop over jets in an event - to find highest jet pT and apply some cuts |
978 | for (Int_t ijet = 0; ijet < Njets; ijet++){ |
979 | // get our jets |
980 | AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet)); |
981 | if (!jet) continue; |
982 | |
f1ce4cc3 |
983 | // apply jet cuts |
8ee2d316 |
984 | if ((jet->Phi()<fPhimin)||(jet->Phi()>fPhimax)) continue; |
985 | if ((jet->Eta()<fEtamin)||(jet->Eta()>fEtamax)) continue; |
daa728b0 |
986 | if (makeQAhistos) fHistAreavsRawPt[centbin]->Fill(jet->Pt(),jet->Area()); |
8ee2d316 |
987 | if(!AcceptMyJet(jet)) continue; |
fc392675 |
988 | |
989 | NjetAcc++; // # of accepted jets |
56523939 |
990 | |
991 | // if FlavourJetAnalysis, get desired FlavTag and check against Jet |
992 | if(doFlavourJetAnalysis) { if(!AcceptFlavourJet(jet, fJetFlavTag)) continue;} |
993 | |
8ee2d316 |
994 | // use this to get total # of jets passing cuts in events!!!!!!!! |
995 | if (makeQAhistos) fHistJetPhi->Fill(jet->Phi()); // Jet Phi histogram (filled) |
fc392675 |
996 | |
997 | // get highest Pt jet in event |
998 | if(highestjetpt<jet->Pt()){ |
999 | ijethi=ijet; |
1000 | highestjetpt=jet->Pt(); |
1001 | } |
1002 | } // end of looping over jets |
1003 | |
1004 | // accepted jets |
1005 | fHistNjetvsCent->Fill(fCent,NjetAcc); |
1006 | Int_t NJETAcc = 0; |
f9b4eee4 |
1007 | fHistEventQA->Fill(9); // events after track/jet loop to get highest pt |
1008 | |
fc392675 |
1009 | // loop over jets in event and make appropriate cuts |
1010 | for (Int_t ijet = 0; ijet < Njets; ++ijet) { |
1011 | AliEmcalJet *jet = static_cast<AliEmcalJet*>(fJets->At(ijet)); |
8ee2d316 |
1012 | if (!jet) continue; |
1013 | |
3f135c41 |
1014 | // check our jet is in our selected trigger event |
d60e724e |
1015 | if (!(trig & fTriggerEventType)) continue; |
1016 | |
8ee2d316 |
1017 | // (should probably be higher..., but makes a cut on jet pT) |
1018 | if (jet->Pt()<0.1) continue; |
f1ce4cc3 |
1019 | // do we accept jet? apply jet cuts |
8ee2d316 |
1020 | if (!AcceptMyJet(jet)) continue; |
fc392675 |
1021 | |
56523939 |
1022 | // if FlavourJetAnalysis, get desired FlavTag and check against Jet |
1023 | if(doFlavourJetAnalysis) { if(!AcceptFlavourJet(jet, fJetFlavTag)) continue;} |
1024 | |
f9b4eee4 |
1025 | fHistEventQA->Fill(10); // accepted jets |
1026 | |
fc392675 |
1027 | // check on lead jet |
1028 | Double_t leadjet=0; |
1029 | if (ijet==ijethi) leadjet=1; |
1030 | |
daa728b0 |
1031 | // check on leading hadron pt |
139e13f4 |
1032 | if (ijet==ijethi) leadhadronPT = GetLeadingHadronPt(jet); |
daa728b0 |
1033 | |
fc392675 |
1034 | // initialize and calculate various parameters: pt, eta, phi, rho, etc... |
1035 | Double_t jetphi = jet->Phi(); // phi of jet |
1036 | NJETAcc++; // # accepted jets |
fc392675 |
1037 | Double_t jeteta = jet->Eta(); // ETA of jet |
aa5f3717 |
1038 | Double_t jetPt = -500; |
1039 | Double_t jetPtGlobal = -500; |
d60e724e |
1040 | Double_t jetPtLocal = -500; // initialize corr jet pt |
1041 | if(GetBeamType() == 1) { |
1042 | fLocalRhoVal = fLocalRho->GetLocalVal(jetphi, fJetRad); //GetJetRadius(0)); // get local rho value |
1043 | jetPtLocal = jet->Pt()-jet->Area()*fLocalRhoVal; // corrected pT of jet using Rho modulated for V2 and V3 |
1044 | } |
8ee2d316 |
1045 | jetPt = jet->Pt(); |
1046 | jetPtGlobal = jet->Pt()-jet->Area()*fRhoVal; // corrected pT of jet from rho value |
a8150e70 |
1047 | Double_t dEP = -500; // initialize angle between jet and event plane |
fc392675 |
1048 | dEP = RelativeEPJET(jetphi,fEPV0); // angle betweeen jet and event plane |
1049 | |
1050 | // make histo's |
8ee2d316 |
1051 | if(makeQAhistos) fHistJetPtvsTrackPt[centbin]->Fill(jetPt,jet->MaxTrackPt()); |
1052 | if(makeQAhistos) fHistRawJetPtvsTrackPt[centbin]->Fill(jetPt,jet->MaxTrackPt()); |
1053 | if(makeQAhistos) fHistJetPtcorrGlRho[centbin]->Fill(jetPtGlobal); |
daa728b0 |
1054 | if(makeQAhistos) fHistJetPtvsdEP[centbin]->Fill(jetPt, dEP); |
1055 | if(makeQAhistos) fHistJetEtaPhiPt[centbin]->Fill(jetPt,jet->Eta(),jet->Phi()); |
1056 | if(makeQAhistos) fHistJetPtArea[centbin]->Fill(jetPt,jet->Area()); |
1057 | if(makeQAhistos) fHistJetEtaPhi->Fill(jet->Eta(),jet->Phi()); // fill jet eta-phi distribution histo |
8ee2d316 |
1058 | if(makeextraCORRhistos) fHistJetPtNcon[centbin]->Fill(jetPt,jet->GetNumberOfConstituents()); |
1059 | if(makeextraCORRhistos) fHistJetPtNconCh[centbin]->Fill(jetPt,jet->GetNumberOfTracks()); |
1060 | if(makeextraCORRhistos) fHistJetPtNconEm[centbin]->Fill(jetPt,jet->GetNumberOfClusters()); |
daa728b0 |
1061 | if (makeoldJEThadhistos) fHistJetPt[centbin]->Fill(jet->Pt()); // fill #jets vs pT histo |
8ee2d316 |
1062 | //fHistDeltaPtvsArea->Fill(jetPt,jet->Area()); |
fc392675 |
1063 | |
1064 | // make histo's with BIAS applied |
1065 | if (jet->MaxTrackPt()>fTrkBias){ |
8ee2d316 |
1066 | if(makeBIAShistos) fHistJetPtvsdEPBias[centbin]->Fill(jetPt, dEP); |
1067 | if(makeBIAShistos) fHistJetEtaPhiPtBias[centbin]->Fill(jetPt,jet->Eta(),jet->Phi()); |
1068 | if(makeextraCORRhistos) fHistJetPtAreaBias[centbin]->Fill(jetPt,jet->Area()); |
1069 | if(makeextraCORRhistos) fHistJetPtNconBias[centbin]->Fill(jetPt,jet->GetNumberOfConstituents()); |
1070 | if(makeextraCORRhistos) fHistJetPtNconBiasCh[centbin]->Fill(jetPt,jet->GetNumberOfTracks()); |
1071 | if(makeextraCORRhistos) fHistJetPtNconBiasEm[centbin]->Fill(jetPt,jet->GetNumberOfClusters()); |
fc392675 |
1072 | } |
1073 | |
8ee2d316 |
1074 | //if(leadjet && centbin==0){ |
1075 | // if(makeextraCORRhistos) fHistJetPt[centbin+1]->Fill(jet->Pt()); |
1076 | //} |
fc392675 |
1077 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){ |
daa728b0 |
1078 | if (makeoldJEThadhistos){ |
8ee2d316 |
1079 | fHistJetPtBias[centbin]->Fill(jet->Pt()); |
1080 | //if(leadjet && centbin==0) fHistJetPtBias[centbin+1]->Fill(jet->Pt()); |
1081 | } |
fc392675 |
1082 | } // end of MaxTrackPt>ftrkBias or maxclusterPt>fclusBias |
1083 | |
1084 | // do we have trigger tracks |
1085 | if(iTT>0){ |
1086 | AliVTrack* TT = static_cast<AliVTrack*>(tracks->At(iTT)); |
1087 | if(TMath::Abs(jet->Phi()-TT->Phi()-TMath::Pi())<0.6) passedTTcut=1; |
1088 | else passedTTcut=0; |
1089 | } // end of check on iTT > 0 |
1090 | |
8ee2d316 |
1091 | if(passedTTcut) { |
daa728b0 |
1092 | if (makeoldJEThadhistos) fHistJetPtTT[centbin]->Fill(jet->Pt()); |
8ee2d316 |
1093 | } |
fc392675 |
1094 | |
1095 | // cut on HIGHEST jet pt in event (15 GeV default) |
139e13f4 |
1096 | //if (highestjetpt>fJetPtcut) { |
1097 | if (jet->Pt() > fJetPtcut) { |
f9b4eee4 |
1098 | fHistEventQA->Fill(11); // jets meeting pt threshold |
72a4d320 |
1099 | |
daa728b0 |
1100 | // does our max track or cluster pass the bias? |
1101 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){ |
1102 | // set up and fill Jet-Hadron Correction THnSparse |
f9b4eee4 |
1103 | Double_t CorrEntries[4] = {fCent, jet->Pt(), dEP, zVtx}; |
daa728b0 |
1104 | fhnCorr->Fill(CorrEntries); // fill Sparse Histo with Correction entries |
d60e724e |
1105 | |
1106 | if(GetBeamType() == 1) fHistLocalRhoJetpt->Fill(jetPtLocal); |
daa728b0 |
1107 | } |
1108 | |
fc392675 |
1109 | // loop over all track for an event containing a jet with a pT>fJetPtCut (15)GeV |
1110 | for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) { |
3f135c41 |
1111 | AliVTrack *track = dynamic_cast<AliVTrack*>(tracks->At(iTracks)); |
fc392675 |
1112 | if (!track) { |
8ee2d316 |
1113 | AliError(Form("Couldn't get AliVtrack %d\n", iTracks)); |
fc392675 |
1114 | continue; |
1115 | } |
1116 | |
8ee2d316 |
1117 | // apply track cuts |
daa728b0 |
1118 | if(TMath::Abs(track->Eta())>fTrkEta) continue; |
1119 | if (track->Pt()<0.15) continue; |
fc392675 |
1120 | |
f9b4eee4 |
1121 | fHistEventQA->Fill(12); // accepted tracks in events from trigger jets |
1122 | |
fc392675 |
1123 | // calculate and get some track parameters |
a8150e70 |
1124 | Double_t trCharge = -99; |
1125 | trCharge = track->Charge(); |
fc392675 |
1126 | Double_t tracketa=track->Eta(); // eta of track |
1127 | Double_t deta=tracketa-jeteta; // dETA between track and jet |
daa728b0 |
1128 | //Double_t dR=sqrt(deta*deta+dphijh*dphijh); // difference of R between jet and hadron track |
1129 | |
a4248490 |
1130 | // calculate single particle tracking efficiency for correlations |
1131 | //TF2 effCORR = fEffFunctionCorr; |
1132 | Double_t trefficiency = -999; |
1133 | trefficiency = EffCorrection(track->Eta(), track->Pt(), fDoEffCorr); |
1134 | |
daa728b0 |
1135 | Int_t ieta = -1; // initialize deta bin |
1136 | Int_t iptjet = -1; // initialize jet pT bin |
1137 | if (makeoldJEThadhistos) { |
1138 | ieta=GetEtaBin(deta); // bin of eta |
f1ce4cc3 |
1139 | if(ieta<0) continue; // double check we don't have a negative array index |
daa728b0 |
1140 | iptjet=GetpTjetBin(jet->Pt()); // bin of jet pt |
1141 | if(iptjet<0) continue; // double check we don't have a negative array index |
f1ce4cc3 |
1142 | } |
fc392675 |
1143 | |
1144 | // dPHI between jet and hadron |
1145 | Double_t dphijh = RelativePhi(jet->Phi(), track->Phi()); // angle between jet and hadron |
1146 | |
fc392675 |
1147 | // fill some jet-hadron histo's |
daa728b0 |
1148 | if (makeoldJEThadhistos) fHistJetH[centbin][iptjet][ieta]->Fill(dphijh,track->Pt()); // fill jet-hadron dPHI--track pT distribution |
1149 | if(makeQAhistos) fHistJetHEtaPhi->Fill(deta,dphijh); // fill jet-hadron eta--phi distribution |
8ee2d316 |
1150 | fHistJetHaddPHI->Fill(dphijh); |
1151 | if(passedTTcut){ |
daa728b0 |
1152 | if (makeoldJEThadhistos) fHistJetHTT[centbin][iptjet][ieta]->Fill(dphijh,track->Pt()); |
8ee2d316 |
1153 | } |
fc392675 |
1154 | |
1155 | // does our max track or cluster pass the bias? |
1156 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){ |
1157 | // set up and fill Jet-Hadron THnSparse |
daa728b0 |
1158 | Double_t triggerEntries[9] = {fCent, jet->Pt(), track->Pt(), deta, dphijh, dEP, zVtx, trCharge, leadjet}; |
a4248490 |
1159 | //cout<<"itracks#: "<<iTracks<<" tracking efficiency: "<<trefficiency<<endl; |
1160 | if(fDoEventMixing) fhnJH->Fill(triggerEntries, 1.0/trefficiency); // fill Sparse Histo with trigger entries |
d60e724e |
1161 | |
fc392675 |
1162 | // fill histo's |
daa728b0 |
1163 | if(makeQAhistos) fHistSEphieta->Fill(dphijh, deta); // single event distribution |
a4248490 |
1164 | if(makeoldJEThadhistos) fHistJetHBias[centbin][iptjet][ieta]->Fill(dphijh,track->Pt()); |
8ee2d316 |
1165 | |
1166 | if (makeBIAShistos) { |
1167 | fHistJetHaddPhiBias->Fill(dphijh); |
fc392675 |
1168 | |
8ee2d316 |
1169 | // in plane and out of plane histo's |
1170 | if( dEP>0 && dEP<=(TMath::Pi()/6) ){ |
1171 | // we are IN plane |
1172 | fHistJetHaddPhiINBias->Fill(dphijh); |
1173 | }else if( dEP>(TMath::Pi()/3) && dEP<=(TMath::Pi()/2) ){ |
1174 | // we are OUT of PLANE |
1175 | fHistJetHaddPhiOUTBias->Fill(dphijh); |
1176 | }else if( dEP>(TMath::Pi()/6) && dEP<=(TMath::Pi()/3) ){ |
1177 | // we are in middle of plane |
1178 | fHistJetHaddPhiMIDBias->Fill(dphijh); |
1179 | } |
1180 | } // BIAS histos switch |
a4248490 |
1181 | |
1182 | if (!fPIDResponse) break; // just return |
1183 | fHistTPCdEdX->Fill(track->Pt(), track->GetTPCsignal()); |
1184 | |
fc392675 |
1185 | } // end of check maxtrackpt>ftrackbias or maxclusterpt>fclustbias |
1186 | |
daa728b0 |
1187 | // ************************************************************************************************************** |
1188 | // *********************************** PID ********************************************************************** |
1189 | // ************************************************************************************************************** |
1190 | if(doPIDtrackBIAS){ |
1191 | //if(ptmax < fTrkBias) continue; // force PID to happen when max track pt > 5.0 GeV |
1192 | if(leadhadronPT < fTrkBias) continue; // force PID to happen when lead hadron pt > 5.0 GeV |
1193 | } |
8ee2d316 |
1194 | |
daa728b0 |
1195 | // some variables for PID |
a4248490 |
1196 | Double_t pt = -999, dEdx = -999, ITSsig = -999, TOFsig = -999, charge = -999; |
daa728b0 |
1197 | |
1198 | // nSigma of particles in TPC, TOF, and ITS |
1199 | Double_t nSigmaPion_TPC, nSigmaProton_TPC, nSigmaKaon_TPC; |
1200 | Double_t nSigmaPion_TOF, nSigmaProton_TOF, nSigmaKaon_TOF; |
1201 | Double_t nSigmaPion_ITS, nSigmaProton_ITS, nSigmaKaon_ITS; |
1202 | |
a4248490 |
1203 | if(doPID && ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)) ){ |
daa728b0 |
1204 | // get parameters of track |
1205 | charge = track->Charge(); // charge of track |
daa728b0 |
1206 | pt = track->Pt(); // pT of track |
1207 | |
1208 | // extra attempt |
1209 | AliVEvent *vevent=InputEvent(); |
1210 | if (!vevent||!fPIDResponse) return kTRUE; // just return, maybe put at beginning |
1211 | |
f9b4eee4 |
1212 | fHistEventQA->Fill(13); // check for AliVEvent and fPIDresponse objects |
1213 | |
3badfe33 |
1214 | // get detector signals |
1215 | dEdx = track->GetTPCsignal(); |
1216 | ITSsig = track->GetITSsignal(); |
1217 | TOFsig = track->GetTOFsignal(); |
daa728b0 |
1218 | |
3badfe33 |
1219 | // TPC nSigma's |
1220 | nSigmaPion_TPC = fPIDResponse->NumberOfSigmasTPC(track,AliPID::kPion); |
1221 | nSigmaKaon_TPC = fPIDResponse->NumberOfSigmasTPC(track,AliPID::kKaon); |
1222 | nSigmaProton_TPC = fPIDResponse->NumberOfSigmasTPC(track,AliPID::kProton); |
1223 | |
1224 | // TOF nSigma's |
1225 | nSigmaPion_TOF = fPIDResponse->NumberOfSigmasTOF(track,AliPID::kPion); |
1226 | nSigmaKaon_TOF = fPIDResponse->NumberOfSigmasTOF(track,AliPID::kKaon); |
1227 | nSigmaProton_TOF = fPIDResponse->NumberOfSigmasTOF(track,AliPID::kProton); |
1228 | |
1229 | // ITS nSigma's |
1230 | nSigmaPion_ITS = fPIDResponse->NumberOfSigmasITS(track,AliPID::kPion); |
1231 | nSigmaKaon_ITS = fPIDResponse->NumberOfSigmasITS(track,AliPID::kKaon); |
1232 | nSigmaProton_ITS = fPIDResponse->NumberOfSigmasITS(track,AliPID::kProton); |
e56c3f6e |
1233 | |
daa728b0 |
1234 | // fill detector signal histograms |
a4248490 |
1235 | //if (makeQAhistos) fHistTPCdEdX->Fill(pt, dEdx); // temp out |
daa728b0 |
1236 | if (makeQAhistos) fHistITSsignal->Fill(pt, ITSsig); |
1237 | //if (makeQAhistos) fHistTOFsignal->Fill(pt, TOFsig); |
1238 | |
1239 | // Tests to PID pions, kaons, and protons, (default is undentified tracks) |
a4248490 |
1240 | //Double_t nPIDtpc = 0, nPIDits = 0, nPIDtof = 0; |
daa728b0 |
1241 | Double_t nPID = -99; |
1242 | |
1243 | // check track has pT < 0.900 GeV - use TPC pid |
1244 | if (pt<0.900 && dEdx>0) { |
a4248490 |
1245 | //nPIDtpc = 4; |
d60e724e |
1246 | nPID = 1; |
daa728b0 |
1247 | |
1248 | // PION check - TPC |
1249 | if (TMath::Abs(nSigmaPion_TPC)<2 && TMath::Abs(nSigmaKaon_TPC)>2 && TMath::Abs(nSigmaProton_TPC)>2 ){ |
1250 | isPItpc = kTRUE; |
a4248490 |
1251 | //nPIDtpc = 1; |
d60e724e |
1252 | nPID=2; |
daa728b0 |
1253 | }else isPItpc = kFALSE; |
1254 | |
1255 | // KAON check - TPC |
1256 | if (TMath::Abs(nSigmaKaon_TPC)<2 && TMath::Abs(nSigmaPion_TPC)>3 && TMath::Abs(nSigmaProton_TPC)>2 ){ |
1257 | isKtpc = kTRUE; |
a4248490 |
1258 | //nPIDtpc = 2; |
d60e724e |
1259 | nPID=3; |
daa728b0 |
1260 | }else isKtpc = kFALSE; |
1261 | |
1262 | // PROTON check - TPC |
1263 | if (TMath::Abs(nSigmaProton_TPC)<2 && TMath::Abs(nSigmaPion_TPC)>3 && TMath::Abs(nSigmaKaon_TPC)>2 ){ |
1264 | isPtpc = kTRUE; |
a4248490 |
1265 | //nPIDtpc = 3; |
d60e724e |
1266 | nPID=4; |
daa728b0 |
1267 | }else isPtpc = kFALSE; |
1268 | } // cut on track pT for TPC |
1269 | |
1270 | // check track has pT < 0.500 GeV - use ITS pid |
1271 | if (pt<0.500 && ITSsig>0) { |
a4248490 |
1272 | //nPIDits = 4; |
d60e724e |
1273 | nPID = 5; |
daa728b0 |
1274 | |
1275 | // PION check - ITS |
1276 | if (TMath::Abs(nSigmaPion_ITS)<2 && TMath::Abs(nSigmaKaon_ITS)>2 && TMath::Abs(nSigmaProton_ITS)>2 ){ |
1277 | isPIits = kTRUE; |
a4248490 |
1278 | //nPIDits = 1; |
d60e724e |
1279 | nPID=6; |
daa728b0 |
1280 | }else isPIits = kFALSE; |
1281 | |
1282 | // KAON check - ITS |
1283 | if (TMath::Abs(nSigmaKaon_ITS)<2 && TMath::Abs(nSigmaPion_ITS)>3 && TMath::Abs(nSigmaProton_ITS)>2 ){ |
1284 | isKits = kTRUE; |
a4248490 |
1285 | //nPIDits = 2; |
d60e724e |
1286 | nPID=7; |
daa728b0 |
1287 | }else isKits = kFALSE; |
1288 | |
1289 | // PROTON check - ITS |
1290 | if (TMath::Abs(nSigmaProton_ITS)<2 && TMath::Abs(nSigmaPion_ITS)>3 && TMath::Abs(nSigmaKaon_ITS)>2 ){ |
1291 | isPits = kTRUE; |
a4248490 |
1292 | //nPIDits = 3; |
d60e724e |
1293 | nPID=8; |
daa728b0 |
1294 | }else isPits = kFALSE; |
1295 | } // cut on track pT for ITS |
1296 | |
1297 | // check track has 0.900 GeV < pT < 2.500 GeV - use TOF pid |
1298 | if (pt>0.900 && pt<2.500 && TOFsig>0) { |
a4248490 |
1299 | //nPIDtof = 4; |
d60e724e |
1300 | nPID = 9; |
daa728b0 |
1301 | |
1302 | // PION check - TOF |
1303 | if (TMath::Abs(nSigmaPion_TOF)<2 && TMath::Abs(nSigmaKaon_TOF)>2 && TMath::Abs(nSigmaProton_TOF)>2 ){ |
1304 | isPItof = kTRUE; |
a4248490 |
1305 | //nPIDtof = 1; |
d60e724e |
1306 | nPID=10; |
daa728b0 |
1307 | }else isPItof = kFALSE; |
1308 | |
1309 | // KAON check - TOF |
1310 | if (TMath::Abs(nSigmaKaon_TOF)<2 && TMath::Abs(nSigmaPion_TOF)>3 && TMath::Abs(nSigmaProton_TOF)>2 ){ |
1311 | isKtof = kTRUE; |
a4248490 |
1312 | //nPIDtof = 2; |
d60e724e |
1313 | nPID=11; |
daa728b0 |
1314 | }else isKtof = kFALSE; |
1315 | |
1316 | // PROTON check - TOF |
1317 | if (TMath::Abs(nSigmaProton_TOF)<2 && TMath::Abs(nSigmaPion_TOF)>3 && TMath::Abs(nSigmaKaon_TOF)>2 ){ |
1318 | isPtof = kTRUE; |
a4248490 |
1319 | //nPIDtof = 3; |
d60e724e |
1320 | nPID=12; |
daa728b0 |
1321 | }else isPtof = kFALSE; |
1322 | } // cut on track pT for TOF |
1323 | |
d60e724e |
1324 | if (nPID == -99) nPID = 14; |
daa728b0 |
1325 | fHistPID->Fill(nPID); |
fc392675 |
1326 | |
1327 | // PID sparse getting filled |
139e13f4 |
1328 | if (allpidAXIS) { // FILL ALL axis |
a4248490 |
1329 | Double_t pid_EntriesALL[19] = {fCent,pt,charge,deta,dphijh,leadjet,zVtx,dEP,jetPt, |
139e13f4 |
1330 | nSigmaPion_TPC, nSigmaPion_TOF, // pion nSig values in TPC/TOF |
a4248490 |
1331 | nPID, //nPIDtpc, nPIDits, nPIDtof, // PID label for each detector |
139e13f4 |
1332 | nSigmaProton_TPC, nSigmaKaon_TPC, // nSig in TPC |
1333 | nSigmaPion_ITS, nSigmaProton_ITS, nSigmaKaon_ITS, // nSig in ITS |
1334 | nSigmaProton_TOF, nSigmaKaon_TOF, // nSig in TOF |
1335 | }; //array for PID sparse |
a4248490 |
1336 | fhnPID->Fill(pid_EntriesALL, 1.0/trefficiency); |
139e13f4 |
1337 | } else { |
1338 | // PID sparse getting filled |
a4248490 |
1339 | Double_t pid_Entries[12] = {fCent,pt,charge,deta,dphijh,leadjet,zVtx,dEP,jetPt, |
139e13f4 |
1340 | nSigmaPion_TPC, nSigmaPion_TOF, // pion nSig values in TPC/TOF |
a4248490 |
1341 | nPID //nPIDtpc, nPIDits, nPIDtof // PID label for each detector |
fc392675 |
1342 | }; //array for PID sparse |
a4248490 |
1343 | fhnPID->Fill(pid_Entries, 1.0/trefficiency); // fill Sparse histo of PID tracks |
139e13f4 |
1344 | } // minimal pid sparse filling |
1345 | |
daa728b0 |
1346 | } // end of doPID check |
fc392675 |
1347 | |
8ee2d316 |
1348 | // get track pt bin |
fc392675 |
1349 | Int_t itrackpt = -500; // initialize track pT bin |
1350 | itrackpt = GetpTtrackBin(track->Pt()); |
1351 | |
8ee2d316 |
1352 | // all tracks: jet hadron correlations in hadron pt bins |
daa728b0 |
1353 | if(makeextraCORRhistos) fHistJetHadbindPhi[itrackpt]->Fill(dphijh); |
fc392675 |
1354 | |
1355 | // in plane and out of plane jet-hadron histo's |
1356 | if( dEP>0 && dEP<=(TMath::Pi()/6) ){ |
daa728b0 |
1357 | // we are IN plane |
1358 | if(makeextraCORRhistos) fHistJetHaddPhiINcent[centbin]->Fill(dphijh); |
1359 | fHistJetHaddPhiIN->Fill(dphijh); |
1360 | if(makeextraCORRhistos) fHistJetHadbindPhiIN[itrackpt]->Fill(dphijh); |
fc392675 |
1361 | }else if( dEP>(TMath::Pi()/3) && dEP<=(TMath::Pi()/2) ){ |
daa728b0 |
1362 | // we are OUT of PLANE |
1363 | if(makeextraCORRhistos) fHistJetHaddPhiOUTcent[centbin]->Fill(dphijh); |
1364 | fHistJetHaddPhiOUT->Fill(dphijh); |
1365 | if(makeextraCORRhistos) fHistJetHadbindPhiOUT[itrackpt]->Fill(dphijh); |
fc392675 |
1366 | }else if( dEP>(TMath::Pi()/6) && dEP<=(TMath::Pi()/3) ){ |
daa728b0 |
1367 | // we are in the middle of plane |
1368 | if(makeextraCORRhistos) fHistJetHaddPhiMIDcent[centbin]->Fill(dphijh); |
1369 | fHistJetHaddPhiMID->Fill(dphijh); |
1370 | if(makeextraCORRhistos) fHistJetHadbindPhiMID[itrackpt]->Fill(dphijh); |
fc392675 |
1371 | } |
1372 | } // loop over tracks found in event with highest JET pT > 10.0 GeV (change) |
1373 | } // jet pt cut |
1374 | } // jet loop |
1375 | |
f9b4eee4 |
1376 | fHistEventQA->Fill(14); // events right before event mixing |
1377 | |
fc392675 |
1378 | // *************************************************************************************************************** |
1379 | // ******************************** Event MIXING ***************************************************************** |
3badfe33 |
1380 | // *************************************************************************************************************** |
1381 | |
1382 | // initialize object array of cloned picotracks |
1383 | TObjArray* tracksClone = 0x0; |
1384 | |
1385 | // PbPb collisions - create cloned picotracks |
d60e724e |
1386 | //if(GetBeamType() == 1) tracksClone = CloneAndReduceTrackList(tracks); // TEST |
fc392675 |
1387 | |
1388 | //Prepare to do event mixing |
1389 | if(fDoEventMixing>0){ |
1390 | // event mixing |
1391 | |
1392 | // 1. First get an event pool corresponding in mult (cent) and |
1393 | // zvertex to the current event. Once initialized, the pool |
1394 | // should contain nMix (reduced) events. This routine does not |
1395 | // pre-scan the chain. The first several events of every chain |
1396 | // will be skipped until the needed pools are filled to the |
1397 | // specified depth. If the pool categories are not too rare, this |
1398 | // should not be a problem. If they are rare, you could lose |
1399 | // statistics. |
1400 | |
1401 | // 2. Collect the whole pool's content of tracks into one TObjArray |
1402 | // (bgTracks), which is effectively a single background super-event. |
1403 | |
1404 | // 3. The reduced and bgTracks arrays must both be passed into |
1405 | // FillCorrelations(). Also nMix should be passed in, so a weight |
1406 | // of 1./nMix can be applied. |
1407 | |
1408 | // mix jets from triggered events with tracks from MB events |
3badfe33 |
1409 | // get the trigger bit, need to change trigger bits between different runs |
1410 | UInt_t trigger = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected(); |
d60e724e |
1411 | // if event was not selected (triggered) for any reseason (should never happen) then return |
1412 | if (trigger==0) return kTRUE; |
fc392675 |
1413 | |
3badfe33 |
1414 | // initialize event pools |
1415 | AliEventPool* pool = 0x0; |
1416 | AliEventPool* poolpp = 0x0; |
1417 | Double_t Ntrks = -999; |
1418 | |
1419 | // pp collisions - get event pool |
1420 | if(GetBeamType() == 0) { |
1421 | Ntrks=(Double_t)Ntracks*1.0; |
1422 | //cout<<"Test.. Ntrks: "<<fPoolMgr->GetEventPool(Ntrks); |
3badfe33 |
1423 | poolpp = fPoolMgr->GetEventPool(Ntrks, zVtx); // for pp |
1424 | } |
fc392675 |
1425 | |
3badfe33 |
1426 | // PbPb collisions - get event pool |
1427 | if(GetBeamType() == 1) pool = fPoolMgr->GetEventPool(fCent, zVtx); // for PbPb? fcent |
fc392675 |
1428 | |
d60e724e |
1429 | // if we don't have a pool, return |
3badfe33 |
1430 | if (!pool && !poolpp){ |
1431 | if(GetBeamType() == 1) AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCent, zVtx)); |
1432 | if(GetBeamType() == 0) AliFatal(Form("No pool found for multiplicity = %f, zVtx = %f", Ntrks, zVtx)); |
fc392675 |
1433 | return kTRUE; |
1434 | } |
1435 | |
f9b4eee4 |
1436 | fHistEventQA->Fill(15); // mixed events cases that have pool |
1437 | |
3badfe33 |
1438 | // initialize background tracks array |
1439 | TObjArray* bgTracks; |
1440 | |
1441 | // next line might not apply for PbPb collisions |
1442 | // use only jets from EMCal-triggered events (for lhc11a use AliVEvent::kEMC1) |
1443 | //check for a trigger jet |
1444 | // fmixingtrack/10 ?? |
a409c5ba |
1445 | if(GetBeamType() == 1) if(trigger & fTriggerEventType) { //kEMCEJE)) { |
1446 | if (pool->IsReady() || pool->NTracksInPool() > fNMIXtracks || pool->GetCurrentNEvents() >= fNMIXevents) { |
3badfe33 |
1447 | |
1448 | // loop over jets (passing cuts?) |
1449 | for (Int_t ijet = 0; ijet < Njets; ijet++) { |
1450 | Double_t leadjet=0; |
1451 | if (ijet==ijethi) leadjet=1; |
1452 | |
1453 | // get jet object |
1454 | AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet)); |
1455 | if (!jet) continue; |
1456 | |
1457 | // (should probably be higher..., but makes a cut on jet pT) |
1458 | if (jet->Pt()<0.1) continue; |
1459 | if (!AcceptMyJet(jet)) continue; |
1460 | |
1461 | fHistEventQA->Fill(16); // event mixing jets |
1462 | |
1463 | // set cut to do event mixing only if we have a jet meeting our pt threshold (bias applied below) |
1464 | if (jet->Pt()<fJetPtcut) continue; |
1465 | |
1466 | // get number of current events in pool |
1467 | Int_t nMix = pool->GetCurrentNEvents(); // how many particles in pool to mix |
1468 | |
1469 | // Fill for biased jet triggers only |
1470 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)) { // && jet->Pt() > fJetPtcut) { |
1471 | // Fill mixed-event histos here |
1472 | for (Int_t jMix=0; jMix<nMix; jMix++) { |
1473 | fHistEventQA->Fill(17); // event mixing nMix |
1474 | |
1475 | // get jMix'th event |
1476 | bgTracks = pool->GetEvent(jMix); |
1477 | const Int_t Nbgtrks = bgTracks->GetEntries(); |
1478 | for(Int_t ibg=0; ibg<Nbgtrks; ibg++) { |
1479 | AliPicoTrack *part = static_cast<AliPicoTrack*>(bgTracks->At(ibg)); |
1480 | if(!part) continue; |
1481 | if(TMath::Abs(part->Eta())>0.9) continue; |
1482 | if(part->Pt()<0.15) continue; |
1483 | |
1484 | Double_t DEta = part->Eta()-jet->Eta(); // difference in eta |
1485 | Double_t DPhi = RelativePhi(jet->Phi(),part->Phi()); // difference in phi |
1486 | Double_t dEP = RelativeEPJET(jet->Phi(),fEPV0); // difference between jet and EP |
1487 | Double_t mixcharge = part->Charge(); |
1488 | //Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta); // difference in R |
a4248490 |
1489 | |
1490 | // calculate single particle tracking efficiency of mixed events for correlations |
1491 | Double_t mixefficiency = -999; |
1492 | mixefficiency = EffCorrection(part->Eta(), part->Pt(), fDoEffCorr); |
1493 | |
3badfe33 |
1494 | // create / fill mixed event sparse |
3f135c41 |
1495 | Double_t triggerEntries[9] = {fCent,jet->Pt(),part->Pt(),DEta,DPhi,dEP,zVtx, mixcharge, leadjet}; //array for ME sparse |
a4248490 |
1496 | fhnMixedEvents->Fill(triggerEntries,1./(nMix*mixefficiency)); // fill Sparse histo of mixed events |
3badfe33 |
1497 | |
d60e724e |
1498 | fHistEventQA->Fill(18); // event mixing - nbgtracks |
a4248490 |
1499 | if(makeextraCORRhistos) fHistMEphieta->Fill(DPhi,DEta, 1./(nMix*mixefficiency)); |
3badfe33 |
1500 | } // end of background track loop |
1501 | } // end of filling mixed-event histo's |
1502 | } // end of check for biased jet triggers |
1503 | } // end of jet loop |
d60e724e |
1504 | } // end of check for pool being ready |
1505 | } // end EMC triggered loop |
3badfe33 |
1506 | |
1507 | //============================================================================================================= |
1508 | |
fc392675 |
1509 | // use only jets from EMCal-triggered events (for lhc11a use AliVEvent::kEMC1) |
1510 | /// if (trigger & AliVEvent::kEMC1) { |
3badfe33 |
1511 | // pp collisions |
d60e724e |
1512 | if(GetBeamType() == 0) if(trigger & fTriggerEventType) { //kEMC1)) { |
a409c5ba |
1513 | if (poolpp->IsReady() || poolpp->NTracksInPool() > fNMIXtracks || poolpp->GetCurrentNEvents() >= fNMIXevents) { |
3badfe33 |
1514 | |
fc392675 |
1515 | // loop over jets (passing cuts?) |
1516 | for (Int_t ijet = 0; ijet < Njets; ijet++) { |
8ee2d316 |
1517 | Double_t leadjet=0; |
1518 | if (ijet==ijethi) leadjet=1; |
fc392675 |
1519 | |
1520 | // get jet object |
1521 | AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet)); |
8ee2d316 |
1522 | if (!jet) continue; |
fc392675 |
1523 | |
8ee2d316 |
1524 | // (should probably be higher..., but makes a cut on jet pT) |
1525 | if (jet->Pt()<0.1) continue; |
1526 | if (!AcceptMyJet(jet)) continue; |
fc392675 |
1527 | |
f9b4eee4 |
1528 | fHistEventQA->Fill(16); // event mixing jets |
1529 | |
3badfe33 |
1530 | // set cut to do event mixing only if we have a jet meeting our pt threshold (bias applied below) |
1531 | if (jet->Pt()<fJetPtcut) continue; |
fc392675 |
1532 | |
3badfe33 |
1533 | // get number of current events in pool |
1534 | Int_t nMix = poolpp->GetCurrentNEvents(); // how many particles in pool to mix |
e56c3f6e |
1535 | |
fc392675 |
1536 | // Fill for biased jet triggers only |
e56c3f6e |
1537 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)) { // && jet->Pt() > fJetPtcut) { |
fc392675 |
1538 | // Fill mixed-event histos here |
1539 | for (Int_t jMix=0; jMix<nMix; jMix++) { |
3badfe33 |
1540 | fHistEventQA->Fill(17); // event mixing nMix |
1541 | |
1542 | // get jMix'th event |
1543 | bgTracks = poolpp->GetEvent(jMix); |
1544 | const Int_t Nbgtrks = bgTracks->GetEntries(); |
1545 | for(Int_t ibg=0; ibg<Nbgtrks; ibg++) { |
1546 | AliPicoTrack *part = static_cast<AliPicoTrack*>(bgTracks->At(ibg)); |
1547 | if(!part) continue; |
1548 | if(TMath::Abs(part->Eta())>0.9) continue; |
1549 | if(part->Pt()<0.15) continue; |
1550 | |
1551 | Double_t DEta = part->Eta()-jet->Eta(); // difference in eta |
1552 | Double_t DPhi = RelativePhi(jet->Phi(),part->Phi()); // difference in phi |
3f135c41 |
1553 | Double_t dEP = RelativeEPJET(jet->Phi(),fEPV0); // difference between jet and EP |
3badfe33 |
1554 | Double_t mixcharge = part->Charge(); |
3f135c41 |
1555 | //Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta); // difference in R |
d60e724e |
1556 | |
3badfe33 |
1557 | // create / fill mixed event sparse |
3f135c41 |
1558 | Double_t triggerEntries[9] = {fCent,jet->Pt(),part->Pt(),DEta,DPhi,dEP,zVtx, mixcharge, leadjet}; //array for ME sparse |
3badfe33 |
1559 | fhnMixedEvents->Fill(triggerEntries,1./nMix); // fill Sparse histo of mixed events |
d60e724e |
1560 | |
1561 | fHistEventQA->Fill(18); // event mixing - nbgtracks |
3badfe33 |
1562 | if(makeextraCORRhistos) fHistMEphieta->Fill(DPhi,DEta, 1./nMix); |
1563 | } // end of background track loop |
1564 | } // end of filling mixed-event histo's |
fc392675 |
1565 | } // end of check for biased jet triggers |
1566 | } // end of jet loop |
d60e724e |
1567 | } // end of check for pool being ready |
3badfe33 |
1568 | } //end EMC triggered loop |
fc392675 |
1569 | |
3badfe33 |
1570 | // pp collisions |
3badfe33 |
1571 | if(GetBeamType() == 0) { |
1572 | |
3f135c41 |
1573 | // use only tracks from MB events (for lhc11a use AliVEvent::kMB) //kAnyINT as test |
d60e724e |
1574 | if(trigger & fMixingEventType) { //kMB) { |
3badfe33 |
1575 | |
1576 | // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool |
d60e724e |
1577 | tracksClone = CloneAndReduceTrackList(tracks); |
1578 | |
3badfe33 |
1579 | // update pool if jet in event or not |
1580 | poolpp->UpdatePool(tracksClone); |
1581 | |
1582 | } // check on track from MB events |
1583 | } |
fc392675 |
1584 | |
3badfe33 |
1585 | // PbPb collisions |
1586 | if(GetBeamType() == 1) { |
d60e724e |
1587 | |
1588 | // use only tracks from MB events |
1589 | if(trigger & fMixingEventType) { //kMB) { |
3badfe33 |
1590 | |
d60e724e |
1591 | // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool |
1592 | tracksClone = CloneAndReduceTrackList(tracks); |
1593 | |
1594 | // update pool if jet in event or not |
1595 | pool->UpdatePool(tracksClone); |
1596 | |
1597 | } // MB events |
1598 | } // PbPb collisions |
fc392675 |
1599 | } // end of event mixing |
1600 | |
1601 | // print some stats on the event |
1602 | event++; |
f9b4eee4 |
1603 | fHistEventQA->Fill(19); // events making it to end |
1604 | |
f1ce4cc3 |
1605 | if (doComments) { |
1606 | cout<<"Event #: "<<event<<" Jet Radius: "<<fJetRad<<" Constituent Pt Cut: "<<fConstituentCut<<endl; |
56523939 |
1607 | cout<<"# of jets: "<<Njets<<" NjetAcc: "<<NjetAcc<<" Highest jet pt: "<<highestjetpt<<" leading hadron pt: "<<leadhadronPT<<endl; |
1608 | cout<<"# tracks: "<<Ntracks<<" NtrackAcc: "<<NtrackAcc<<" Highest track pt: "<<ptmax<<endl; |
f1ce4cc3 |
1609 | cout<<" =============================================== "<<endl; |
1610 | } |
fc392675 |
1611 | |
1612 | return kTRUE; // used when the function is of type bool |
1613 | } // end of RUN |
1614 | |
1615 | //________________________________________________________________________ |
1616 | Int_t AliAnalysisTaskEmcalJetHadEPpid::GetCentBin(Double_t cent) const |
1617 | { // Get centrality bin. |
1618 | Int_t centbin = -1; |
1619 | if (cent>=0 && cent<10) centbin = 0; |
1620 | else if (cent>=10 && cent<20) centbin = 1; |
1621 | else if (cent>=20 && cent<30) centbin = 2; |
1622 | else if (cent>=30 && cent<40) centbin = 3; |
1623 | else if (cent>=40 && cent<50) centbin = 4; |
1624 | else if (cent>=50 && cent<90) centbin = 5; |
1625 | |
1626 | return centbin; |
1627 | } |
1628 | |
1629 | //________________________________________________________________________ |
a8150e70 |
1630 | Double_t AliAnalysisTaskEmcalJetHadEPpid::RelativePhi(Double_t mphi,Double_t vphi) const |
fc392675 |
1631 | { // function to calculate relative PHI |
1632 | double dphi = mphi-vphi; |
1633 | |
1634 | // set dphi to operate on adjusted scale |
1635 | if(dphi<-0.5*TMath::Pi()) dphi+=2.*TMath::Pi(); |
1636 | if(dphi>3./2.*TMath::Pi()) dphi-=2.*TMath::Pi(); |
1637 | |
1638 | // test |
1639 | if( dphi < -1.*TMath::Pi()/2 || dphi > 3.*TMath::Pi()/2 ) |
1640 | AliWarning(Form("%s: dPHI not in range [-0.5*Pi, 1.5*Pi]!", GetName())); |
1641 | |
1642 | return dphi; // dphi in [-0.5Pi, 1.5Pi] |
1643 | } |
1644 | |
fc392675 |
1645 | //_________________________________________________________________________ |
a4248490 |
1646 | Double_t AliAnalysisTaskEmcalJetHadEPpid::RelativeEPJET(Double_t jetAng, Double_t EPAng) const |
fc392675 |
1647 | { // function to calculate angle between jet and EP in the 1st quadrant (0,Pi/2) |
1648 | Double_t dphi = (EPAng - jetAng); |
1649 | |
1650 | // ran into trouble with a few dEP<-Pi so trying this... |
1651 | if( dphi<-1*TMath::Pi() ){ |
1652 | dphi = dphi + 1*TMath::Pi(); |
1653 | } // this assumes we are doing full jets currently |
1654 | |
1655 | if( (dphi>0) && (dphi<1*TMath::Pi()/2) ){ |
1656 | // Do nothing! we are in quadrant 1 |
1657 | }else if( (dphi>1*TMath::Pi()/2) && (dphi<1*TMath::Pi()) ){ |
1658 | dphi = 1*TMath::Pi() - dphi; |
1659 | }else if( (dphi<0) && (dphi>-1*TMath::Pi()/2) ){ |
1660 | dphi = fabs(dphi); |
1661 | }else if( (dphi<-1*TMath::Pi()/2) && (dphi>-1*TMath::Pi()) ){ |
1662 | dphi = dphi + 1*TMath::Pi(); |
1663 | } |
1664 | |
1665 | // test |
1666 | if( dphi < 0 || dphi > TMath::Pi()/2 ) |
1667 | AliWarning(Form("%s: dPHI not in range [0, 0.5*Pi]!", GetName())); |
1668 | |
1669 | return dphi; // dphi in [0, Pi/2] |
1670 | } |
1671 | |
1672 | //Int_t ieta=GetEtaBin(deta); |
1673 | //________________________________________________________________________ |
1674 | Int_t AliAnalysisTaskEmcalJetHadEPpid::GetEtaBin(Double_t eta) const |
1675 | { |
1676 | // Get eta bin for histos. |
1677 | Int_t etabin = -1; |
f1ce4cc3 |
1678 | if (TMath::Abs(eta)<=0.4) etabin = 0; |
fc392675 |
1679 | else if (TMath::Abs(eta)>0.4 && TMath::Abs(eta)<0.8) etabin = 1; |
f1ce4cc3 |
1680 | else if (TMath::Abs(eta)>=0.8) etabin = 2; |
fc392675 |
1681 | |
1682 | return etabin; |
1683 | } // end of get-eta-bin |
1684 | |
1685 | //________________________________________________________________________ |
1686 | Int_t AliAnalysisTaskEmcalJetHadEPpid::GetpTjetBin(Double_t pt) const |
1687 | { |
1688 | // Get jet pt bin for histos. |
1689 | Int_t ptbin = -1; |
1690 | if (pt>=15 && pt<20) ptbin = 0; |
1691 | else if (pt>=20 && pt<25) ptbin = 1; |
1692 | else if (pt>=25 && pt<40) ptbin = 2; |
1693 | else if (pt>=40 && pt<60) ptbin = 3; |
f1ce4cc3 |
1694 | else if (pt>=60) ptbin = 4; |
fc392675 |
1695 | |
1696 | return ptbin; |
1697 | } // end of get-jet-pt-bin |
1698 | |
1699 | //________________________________________________________________________ |
1700 | Int_t AliAnalysisTaskEmcalJetHadEPpid::GetpTtrackBin(Double_t pt) const |
1701 | { |
1702 | // May need to update bins for future runs... (testing locally) |
1703 | |
1704 | // Get track pt bin for histos. |
1705 | Int_t ptbin = -1; |
1706 | if (pt < 0.5) ptbin = 0; |
1707 | else if (pt>=0.5 && pt<1.0) ptbin = 1; |
1708 | else if (pt>=1.0 && pt<1.5) ptbin = 2; |
1709 | else if (pt>=1.5 && pt<2.0) ptbin = 3; |
1710 | else if (pt>=2.0 && pt<2.5) ptbin = 4; |
1711 | else if (pt>=2.5 && pt<3.0) ptbin = 5; |
1712 | else if (pt>=3.0 && pt<4.0) ptbin = 6; |
1713 | else if (pt>=4.0 && pt<5.0) ptbin = 7; |
1714 | else if (pt>=5.0) ptbin = 8; |
1715 | |
1716 | return ptbin; |
1717 | } // end of get-jet-pt-bin |
1718 | |
fc392675 |
1719 | //___________________________________________________________________________ |
1720 | Int_t AliAnalysisTaskEmcalJetHadEPpid::GetzVertexBin(Double_t zVtx) const |
1721 | { |
1722 | // get z-vertex bin for histo. |
1723 | int zVbin= -1; |
f1ce4cc3 |
1724 | if (zVtx>=-10 && zVtx<-8) zVbin = 0; |
fc392675 |
1725 | else if (zVtx>=-8 && zVtx<-6) zVbin = 1; |
1726 | else if (zVtx>=-6 && zVtx<-4) zVbin = 2; |
1727 | else if (zVtx>=-4 && zVtx<-2) zVbin = 3; |
1728 | else if (zVtx>=-2 && zVtx<0) zVbin = 4; |
1729 | else if (zVtx>=0 && zVtx<2) zVbin = 5; |
1730 | else if (zVtx>=2 && zVtx<4) zVbin = 6; |
1731 | else if (zVtx>=4 && zVtx<6) zVbin = 7; |
1732 | else if (zVtx>=6 && zVtx<8) zVbin = 8; |
1733 | else if (zVtx>=8 && zVtx<10) zVbin = 9; |
1734 | else zVbin = 10; |
1735 | |
1736 | return zVbin; |
1737 | } // end of get z-vertex bin |
1738 | |
1739 | //______________________________________________________________________ |
e56c3f6e |
1740 | THnSparse* AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseF(const char* name, UInt_t entries) |
fc392675 |
1741 | { |
e56c3f6e |
1742 | // generate new THnSparseF, axes are defined in GetDimParams() |
fc392675 |
1743 | Int_t count = 0; |
1744 | UInt_t tmp = entries; |
1745 | while(tmp!=0){ |
1746 | count++; |
1747 | tmp = tmp &~ -tmp; // clear lowest bit |
1748 | } |
1749 | |
1750 | TString hnTitle(name); |
1751 | const Int_t dim = count; |
1752 | Int_t nbins[dim]; |
1753 | Double_t xmin[dim]; |
1754 | Double_t xmax[dim]; |
1755 | |
1756 | Int_t i=0; |
1757 | Int_t c=0; |
1758 | while(c<dim && i<32){ |
1759 | if(entries&(1<<i)){ |
fc392675 |
1760 | TString label(""); |
1761 | GetDimParams(i, label, nbins[c], xmin[c], xmax[c]); |
1762 | hnTitle += Form(";%s",label.Data()); |
1763 | c++; |
1764 | } |
1765 | |
1766 | i++; |
1767 | } |
1768 | hnTitle += ";"; |
1769 | |
e56c3f6e |
1770 | return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax); |
1771 | } // end of NewTHnSparseF |
fc392675 |
1772 | |
1773 | void AliAnalysisTaskEmcalJetHadEPpid::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax) |
1774 | { |
1775 | // stores label and binning of axis for THnSparse |
1776 | const Double_t pi = TMath::Pi(); |
1777 | |
1778 | switch(iEntry){ |
1779 | |
1780 | case 0: |
1781 | label = "V0 centrality (%)"; |
1782 | nbins = 10; |
1783 | xmin = 0.; |
1784 | xmax = 100.; |
1785 | break; |
1786 | |
1787 | case 1: |
a8150e70 |
1788 | label = "Jet p_{T}"; |
e56c3f6e |
1789 | nbins = 50; // 216 |
a8150e70 |
1790 | xmin = 0.; |
e56c3f6e |
1791 | xmax = 250.; // 216 |
fc392675 |
1792 | break; |
1793 | |
1794 | case 2: |
1795 | label = "Track p_{T}"; |
e56c3f6e |
1796 | nbins = 80; //300; // 750 pid |
fc392675 |
1797 | xmin = 0.; |
e56c3f6e |
1798 | xmax = 20.; //75.; |
fc392675 |
1799 | break; |
1800 | |
1801 | case 3: |
1802 | label = "Relative Eta"; |
3f135c41 |
1803 | nbins = 56; // 42 |
1804 | xmin = -1.4; |
1805 | xmax = 1.4; |
fc392675 |
1806 | break; |
1807 | |
1808 | case 4: |
1809 | label = "Relative Phi"; |
daa728b0 |
1810 | nbins = 72; |
fc392675 |
1811 | xmin = -0.5*pi; |
1812 | xmax = 1.5*pi; |
1813 | break; |
1814 | |
1815 | case 5: |
1816 | label = "Relative angle of Jet and Reaction Plane"; |
3f135c41 |
1817 | nbins = 3; // (12) 72 |
e56c3f6e |
1818 | xmin = 0; |
1819 | xmax = 0.5*pi; |
fc392675 |
1820 | break; |
1821 | |
1822 | case 6: |
1823 | label = "z-vertex"; |
1824 | nbins = 10; |
1825 | xmin = -10; |
1826 | xmax = 10; |
1827 | break; |
1828 | |
1829 | case 7: |
1830 | label = "track charge"; |
a8150e70 |
1831 | nbins = 3; |
fc392675 |
1832 | xmin = -1.5; |
1833 | xmax = 1.5; |
1834 | break; |
1835 | |
daa728b0 |
1836 | case 8: |
fc392675 |
1837 | label = "leading jet"; |
1838 | nbins = 3; |
1839 | xmin = -0.5; |
1840 | xmax = 2.5; |
1841 | break; |
a8150e70 |
1842 | |
1843 | case 9: // need to update |
1844 | label = "leading track"; |
1845 | nbins = 10; |
1846 | xmin = 0; |
1847 | xmax = 50; |
1848 | break; |
1849 | |
fc392675 |
1850 | } // end of switch |
1851 | } // end of getting dim-params |
1852 | |
1853 | //_________________________________________________ |
1854 | // From CF event mixing code PhiCorrelations |
3badfe33 |
1855 | TObjArray* AliAnalysisTaskEmcalJetHadEPpid::CloneAndReduceTrackList(TObjArray* tracksME) |
fc392675 |
1856 | { |
1857 | // clones a track list by using AliPicoTrack which uses much less memory (used for event mixing) |
1858 | TObjArray* tracksClone = new TObjArray; |
1859 | tracksClone->SetOwner(kTRUE); |
1860 | |
3263f852 |
1861 | // =============================== |
1862 | |
1863 | // cout << "RM Hybrid track : " << i << " " << particle->Pt() << endl; |
fc392675 |
1864 | |
1865 | //cout << "nEntries " << tracks->GetEntriesFast() <<endl; |
3badfe33 |
1866 | for (Int_t i=0; i<tracksME->GetEntriesFast(); i++) { // AOD/general case |
1867 | AliVParticle* particle = (AliVParticle*) tracksME->At(i); // AOD/general case |
fc392675 |
1868 | if(TMath::Abs(particle->Eta())>fTrkEta) continue; |
1869 | if(particle->Pt()<0.15)continue; |
1870 | |
8ee2d316 |
1871 | /* |
1872 | // DON'T USE |
fc392675 |
1873 | Double_t trackpt=particle->Pt(); // track pT |
1874 | |
1875 | Int_t trklabel=-1; |
1876 | trklabel=particle->GetLabel(); |
1877 | //cout << "TRACK_LABEL: " << particle->GetLabel()<<endl; |
1878 | |
1879 | Int_t hadbin=-1; |
1880 | if(trackpt<0.5) hadbin=0; |
1881 | else if(trackpt<1) hadbin=1; |
1882 | else if(trackpt<2) hadbin=2; |
1883 | else if(trackpt<3) hadbin=3; |
1884 | else if(trackpt<5) hadbin=4; |
1885 | else if(trackpt<8) hadbin=5; |
1886 | else if(trackpt<20) hadbin=6; |
8ee2d316 |
1887 | // end of DON'T USE |
fc392675 |
1888 | |
1889 | //feb10 comment out |
1890 | if(hadbin>-1 && trklabel>-1 && trklabel <3) fHistTrackEtaPhi[trklabel][hadbin]->Fill(particle->Eta(),particle->Phi()); |
1891 | if(hadbin>-1) fHistTrackEtaPhi[3][hadbin]->Fill(particle->Eta(),particle->Phi()); |
1892 | |
8ee2d316 |
1893 | if(hadbin>-1) fHistTrackEtaPhi[hadbin]->Fill(particle->Eta(),particle->Phi()); // TEST |
1894 | */ |
fc392675 |
1895 | |
1896 | tracksClone->Add(new AliPicoTrack(particle->Pt(), particle->Eta(), particle->Phi(), particle->Charge(), 0, 0, 0, 0)); |
1897 | } // end of looping through tracks |
1898 | |
1899 | return tracksClone; |
1900 | } |
1901 | |
1902 | //____________________________________________________________________________________________ |
e56c3f6e |
1903 | THnSparse* AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseFPID(const char* name, UInt_t entries) |
fc392675 |
1904 | { |
e56c3f6e |
1905 | // generate new THnSparseF PID, axes are defined in GetDimParams() |
fc392675 |
1906 | Int_t count = 0; |
1907 | UInt_t tmp = entries; |
1908 | while(tmp!=0){ |
1909 | count++; |
1910 | tmp = tmp &~ -tmp; // clear lowest bit |
1911 | } |
1912 | |
1913 | TString hnTitle(name); |
1914 | const Int_t dim = count; |
1915 | Int_t nbins[dim]; |
1916 | Double_t xmin[dim]; |
1917 | Double_t xmax[dim]; |
1918 | |
1919 | Int_t i=0; |
1920 | Int_t c=0; |
1921 | while(c<dim && i<32){ |
1922 | if(entries&(1<<i)){ |
fc392675 |
1923 | TString label(""); |
1924 | GetDimParamsPID(i, label, nbins[c], xmin[c], xmax[c]); |
1925 | hnTitle += Form(";%s",label.Data()); |
1926 | c++; |
1927 | } |
1928 | |
1929 | i++; |
1930 | } |
1931 | hnTitle += ";"; |
1932 | |
e56c3f6e |
1933 | return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax); |
1934 | } // end of NewTHnSparseF PID |
fc392675 |
1935 | |
1936 | //________________________________________________________________________________ |
1937 | void AliAnalysisTaskEmcalJetHadEPpid::GetDimParamsPID(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax) |
1938 | { |
1939 | // stores label and binning of axis for THnSparse |
1940 | const Double_t pi = TMath::Pi(); |
1941 | |
1942 | switch(iEntry){ |
1943 | |
1944 | case 0: |
1945 | label = "V0 centrality (%)"; |
1946 | nbins = 10; |
1947 | xmin = 0.; |
1948 | xmax = 100.; |
1949 | break; |
1950 | |
1951 | case 1: |
fc392675 |
1952 | label = "Track p_{T}"; |
e56c3f6e |
1953 | nbins = 80; //300; // 750 |
fc392675 |
1954 | xmin = 0.; |
e56c3f6e |
1955 | xmax = 20.; //75.; |
fc392675 |
1956 | break; |
1957 | |
139e13f4 |
1958 | case 2: |
fc392675 |
1959 | label = "Charge of Track"; |
1960 | nbins = 3; |
1961 | xmin = -1.5; |
1962 | xmax = 1.5; |
1963 | break; |
1964 | |
139e13f4 |
1965 | case 3: |
fc392675 |
1966 | label = "Relative Eta of Track and Jet"; |
3f135c41 |
1967 | nbins = 56; // 48 |
1968 | xmin = -1.4; |
1969 | xmax = 1.4; |
fc392675 |
1970 | break; |
1971 | |
139e13f4 |
1972 | case 4: |
fc392675 |
1973 | label = "Relative Phi of Track and Jet"; |
daa728b0 |
1974 | nbins = 72; |
fc392675 |
1975 | xmin = -0.5*pi; |
1976 | xmax = 1.5*pi; |
1977 | break; |
1978 | |
139e13f4 |
1979 | case 5: |
fc392675 |
1980 | label = "leading jet"; |
1981 | nbins = 3; |
1982 | xmin = -.5; |
1983 | xmax = 2.5; |
1984 | break; |
1985 | |
139e13f4 |
1986 | case 6: |
fc392675 |
1987 | label = "Z-vertex"; |
1988 | nbins = 10; |
1989 | xmin = -10.; |
1990 | xmax = 10.; |
1991 | break; |
1992 | |
139e13f4 |
1993 | case 7: |
8ee2d316 |
1994 | label = "Relative angle: Jet and Reaction Plane"; |
3f135c41 |
1995 | nbins = 3; // (12) 48 |
fc392675 |
1996 | xmin = 0.; |
8ee2d316 |
1997 | xmax = 0.5*pi; |
fc392675 |
1998 | break; |
1999 | |
139e13f4 |
2000 | case 8: |
2001 | label = "Jet p_{T}"; |
e56c3f6e |
2002 | nbins = 50; // 216 |
fc392675 |
2003 | xmin = 0.; |
e56c3f6e |
2004 | xmax = 250.; // 216 |
fc392675 |
2005 | break; |
2006 | |
139e13f4 |
2007 | case 9: |
fc392675 |
2008 | label = "N-Sigma of pions in TPC"; |
8ee2d316 |
2009 | nbins = 200; |
daa728b0 |
2010 | xmin = -10.0; |
2011 | xmax = 10.0; |
fc392675 |
2012 | break; |
2013 | |
139e13f4 |
2014 | case 10: |
2015 | label = "N-Sigma of pions in TOF"; |
2016 | nbins = 200; |
2017 | xmin = -10.; |
2018 | xmax = 10.; |
2019 | break; |
2020 | |
2021 | case 11: |
2022 | label = "TPC PID determination"; |
a4248490 |
2023 | nbins = 15; //5; |
2024 | xmin = 0.5; //0.; |
2025 | xmax = 15.5; //5.; |
139e13f4 |
2026 | break; |
2027 | |
a4248490 |
2028 | /* |
fc392675 |
2029 | case 12: |
139e13f4 |
2030 | label = "ITS PID determination"; |
2031 | nbins = 5; |
2032 | xmin = 0.; |
2033 | xmax = 5.; |
2034 | break; |
2035 | |
2036 | case 13: |
2037 | label = "TOF PID determination"; |
2038 | nbins = 5; |
2039 | xmin = 0.; |
2040 | xmax = 5.; |
2041 | break; |
a4248490 |
2042 | */ |
139e13f4 |
2043 | |
a4248490 |
2044 | case 12: |
fc392675 |
2045 | label = "N-Sigma of protons in TPC"; |
daa728b0 |
2046 | nbins = 200; |
2047 | xmin = -10.; |
2048 | xmax = 10.; |
fc392675 |
2049 | break; |
2050 | |
a4248490 |
2051 | case 13: |
fc392675 |
2052 | label = "N-Sigma of kaons in TPC"; |
daa728b0 |
2053 | nbins = 200; |
2054 | xmin = -10.; |
2055 | xmax = 10.; |
fc392675 |
2056 | break; |
2057 | |
a4248490 |
2058 | case 14: |
fc392675 |
2059 | label = "N-Sigma of pions in ITS"; |
8ee2d316 |
2060 | nbins = 200; |
daa728b0 |
2061 | xmin = -10.0; |
2062 | xmax = 10.0; |
fc392675 |
2063 | break; |
2064 | |
a4248490 |
2065 | case 15: |
fc392675 |
2066 | label = "N-Sigma of protons in ITS"; |
daa728b0 |
2067 | nbins = 200; |
2068 | xmin = -10.; |
2069 | xmax = 10.; |
fc392675 |
2070 | break; |
2071 | |
a4248490 |
2072 | case 16: |
fc392675 |
2073 | label = "N-Sigma of kaons in ITS"; |
daa728b0 |
2074 | nbins = 200; |
2075 | xmin = -10.; |
2076 | xmax = 10.; |
fc392675 |
2077 | break; |
2078 | |
a4248490 |
2079 | case 17: |
fc392675 |
2080 | label = "N-Sigma of protons in TOF"; |
daa728b0 |
2081 | nbins = 200; |
2082 | xmin = -10.; |
2083 | xmax = 10.; |
fc392675 |
2084 | break; |
2085 | |
a4248490 |
2086 | case 18: |
fc392675 |
2087 | label = "N-Sigma of kaons in TOF"; |
daa728b0 |
2088 | nbins = 200; |
2089 | xmin = -10.; |
2090 | xmax = 10.; |
fc392675 |
2091 | break; |
2092 | |
fc392675 |
2093 | } // end of switch |
2094 | } // end of get dimension parameters PID |
2095 | |
8ee2d316 |
2096 | void AliAnalysisTaskEmcalJetHadEPpid::Terminate(Option_t *) { |
fc392675 |
2097 | cout<<"#########################"<<endl; |
2098 | cout<<"#### DONE RUNNING!!! ####"<<endl; |
2099 | cout<<"#########################"<<endl; |
fc392675 |
2100 | } // end of terminate |
2101 | |
8ee2d316 |
2102 | //________________________________________________________________________ |
2103 | Int_t AliAnalysisTaskEmcalJetHadEPpid::AcceptMyJet(AliEmcalJet *jet) { |
2104 | //applies all jet cuts except pt |
2105 | if ((jet->Phi()<fPhimin)||(jet->Phi()>fPhimax)) return 0; |
2106 | if ((jet->Eta()<fEtamin)||(jet->Eta()>fEtamax)) return 0; |
2107 | if (jet->Area()<fAreacut) return 0; |
2108 | // prevents 0 area jets from sneaking by when area cut == 0 |
2109 | if (jet->Area()==0) return 0; |
2110 | //exclude jets with extremely high pt tracks which are likely misreconstructed |
2111 | if(jet->MaxTrackPt()>100) return 0; |
8ee2d316 |
2112 | |
2113 | //passed all above cuts |
2114 | return 1; |
2115 | } |
fc392675 |
2116 | |
2117 | //void AliAnalysisTaskEmcalJetHadEPpid::FillAnalysisSummaryHistogram() const |
2118 | void AliAnalysisTaskEmcalJetHadEPpid::SetfHistPIDcounterLabels(TH1* h) const |
2119 | { |
2120 | // fill the analysis summary histrogram, saves all relevant analysis settigns |
d60e724e |
2121 | h->GetXaxis()->SetBinLabel(1, "TPC: Unidentified"); |
2122 | h->GetXaxis()->SetBinLabel(2, "TPC: Pion"); |
2123 | h->GetXaxis()->SetBinLabel(3, "TPC: Kaon"); |
2124 | h->GetXaxis()->SetBinLabel(4, "TPC: Proton"); |
2125 | h->GetXaxis()->SetBinLabel(5, "ITS: Unidentified"); |
2126 | h->GetXaxis()->SetBinLabel(6, "ITS: Pion"); |
2127 | h->GetXaxis()->SetBinLabel(7, "ITS: Kaon"); |
2128 | h->GetXaxis()->SetBinLabel(8, "ITS: Proton"); |
2129 | h->GetXaxis()->SetBinLabel(9, "TOF: Unidentified"); |
2130 | h->GetXaxis()->SetBinLabel(10, "TOF: Pion"); |
2131 | h->GetXaxis()->SetBinLabel(11, "TOF: Kaon"); |
2132 | h->GetXaxis()->SetBinLabel(12, "TOF: Proton"); |
2133 | h->GetXaxis()->SetBinLabel(14, "Unidentified tracks"); |
2134 | |
2135 | // set x-axis labels vertically |
2136 | h->LabelsOption("v"); |
fc392675 |
2137 | } |
daa728b0 |
2138 | |
f9b4eee4 |
2139 | //void AliAnalysisTaskEmcalJetHadEPpid::FillAnalysisSummaryHistogram() const |
2140 | void AliAnalysisTaskEmcalJetHadEPpid::SetfHistQAcounterLabels(TH1* h) const |
2141 | { |
d60e724e |
2142 | // label bins of the analysis event summary |
f9b4eee4 |
2143 | h->GetXaxis()->SetBinLabel(1, "All events started"); |
2144 | h->GetXaxis()->SetBinLabel(2, "object check"); |
2145 | h->GetXaxis()->SetBinLabel(3, "aod/esd check"); |
2146 | h->GetXaxis()->SetBinLabel(4, "centrality check"); |
2147 | h->GetXaxis()->SetBinLabel(5, "zvertex check"); |
2148 | h->GetXaxis()->SetBinLabel(6, "list check"); |
2149 | h->GetXaxis()->SetBinLabel(7, "track/jet pointer check"); |
e56c3f6e |
2150 | h->GetXaxis()->SetBinLabel(8, "tracks & jets < than 1 check"); |
f9b4eee4 |
2151 | h->GetXaxis()->SetBinLabel(9, "after track/jet loop to get highest pt"); |
2152 | h->GetXaxis()->SetBinLabel(10, "accepted jets"); |
2153 | h->GetXaxis()->SetBinLabel(11, "jets meeting pt threshold"); |
d60e724e |
2154 | h->GetXaxis()->SetBinLabel(12, "accepted tracks in events w/ trigger jet"); |
2155 | h->GetXaxis()->SetBinLabel(13, "after AliVEvent & fPIDResponse"); |
f9b4eee4 |
2156 | h->GetXaxis()->SetBinLabel(14, "events before event mixing"); |
d60e724e |
2157 | h->GetXaxis()->SetBinLabel(15, "mixed events w/ pool"); |
f9b4eee4 |
2158 | h->GetXaxis()->SetBinLabel(16, "event mixing: jets"); |
2159 | h->GetXaxis()->SetBinLabel(17, "event mixing: nMix"); |
2160 | h->GetXaxis()->SetBinLabel(18, "event mixing: nbackground tracks"); |
2161 | h->GetXaxis()->SetBinLabel(19, "event mixing: THE END"); |
d60e724e |
2162 | |
2163 | // set x-axis labels vertically |
2164 | h->LabelsOption("v"); |
2165 | } |
2166 | |
2167 | void AliAnalysisTaskEmcalJetHadEPpid::SetfHistEvtSelQALabels(TH1* h) const |
2168 | { |
2169 | // label bins of the analysis trigger selection summary |
2170 | h->GetXaxis()->SetBinLabel(1, "no trigger"); |
2171 | h->GetXaxis()->SetBinLabel(2, "kAny"); |
2172 | h->GetXaxis()->SetBinLabel(3, "kAnyINT"); |
2173 | h->GetXaxis()->SetBinLabel(4, "kMB"); |
2174 | h->GetXaxis()->SetBinLabel(5, "kINT7"); |
2175 | h->GetXaxis()->SetBinLabel(6, "kEMC1"); |
2176 | h->GetXaxis()->SetBinLabel(7, "kEMC7"); |
2177 | h->GetXaxis()->SetBinLabel(8, "kEMC8"); |
2178 | h->GetXaxis()->SetBinLabel(9, "kEMCEJE"); |
2179 | h->GetXaxis()->SetBinLabel(10, "kEMCEGA"); |
2180 | h->GetXaxis()->SetBinLabel(11, "kCentral"); |
2181 | h->GetXaxis()->SetBinLabel(12, "kSemiCentral"); |
2182 | h->GetXaxis()->SetBinLabel(13, "kINT8"); |
2183 | h->GetXaxis()->SetBinLabel(14, "kEMCEJE or kMB"); |
2184 | h->GetXaxis()->SetBinLabel(15, "kEMCEGA or kMB"); |
2185 | h->GetXaxis()->SetBinLabel(16, "kAnyINT or kMB"); |
2186 | h->GetXaxis()->SetBinLabel(17, "kEMCEJE & kMB"); |
2187 | h->GetXaxis()->SetBinLabel(18, "kEMCEGA & kMB"); |
2188 | h->GetXaxis()->SetBinLabel(19, "kAnyINT & kMB"); |
2189 | |
d60e724e |
2190 | // set x-axis labels vertically |
2191 | h->LabelsOption("v"); |
2192 | //h->LabelsDeflate("X"); |
f9b4eee4 |
2193 | } |
2194 | |
daa728b0 |
2195 | //______________________________________________________________________ |
e56c3f6e |
2196 | THnSparse* AliAnalysisTaskEmcalJetHadEPpid::NewTHnSparseFCorr(const char* name, UInt_t entries) { |
daa728b0 |
2197 | // generate new THnSparseD, axes are defined in GetDimParamsD() |
2198 | Int_t count = 0; |
2199 | UInt_t tmp = entries; |
2200 | while(tmp!=0){ |
2201 | count++; |
2202 | tmp = tmp &~ -tmp; // clear lowest bit |
2203 | } |
2204 | |
2205 | TString hnTitle(name); |
2206 | const Int_t dim = count; |
2207 | Int_t nbins[dim]; |
2208 | Double_t xmin[dim]; |
2209 | Double_t xmax[dim]; |
2210 | |
2211 | Int_t i=0; |
2212 | Int_t c=0; |
2213 | while(c<dim && i<32){ |
2214 | if(entries&(1<<i)){ |
daa728b0 |
2215 | TString label(""); |
2216 | GetDimParamsCorr(i, label, nbins[c], xmin[c], xmax[c]); |
2217 | hnTitle += Form(";%s",label.Data()); |
2218 | c++; |
2219 | } |
2220 | |
2221 | i++; |
2222 | } |
2223 | hnTitle += ";"; |
2224 | |
e56c3f6e |
2225 | return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax); |
2226 | } // end of NewTHnSparseF |
daa728b0 |
2227 | |
2228 | //______________________________________________________________________________________________ |
2229 | void AliAnalysisTaskEmcalJetHadEPpid::GetDimParamsCorr(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax) |
2230 | { |
2231 | //stores label and binning of axis for THnSparse |
2232 | const Double_t pi = TMath::Pi(); |
2233 | |
2234 | switch(iEntry){ |
2235 | |
2236 | case 0: |
2237 | label = "V0 centrality (%)"; |
2238 | nbins = 10; |
2239 | xmin = 0.; |
2240 | xmax = 100.; |
2241 | break; |
2242 | |
2243 | case 1: |
2244 | label = "Jet p_{T}"; |
e56c3f6e |
2245 | nbins = 50; // 216 |
daa728b0 |
2246 | xmin = 0.; |
e56c3f6e |
2247 | xmax = 250.; // 216 |
daa728b0 |
2248 | break; |
2249 | |
2250 | case 2: |
2251 | label = "Relative angle: Jet and Reaction Plane"; |
3f135c41 |
2252 | nbins = 3; // (12) 48 |
daa728b0 |
2253 | xmin = 0.; |
2254 | xmax = 0.5*pi; |
2255 | break; |
2256 | |
f9b4eee4 |
2257 | case 3: |
2258 | label = "Z-vertex"; |
2259 | nbins = 10; |
2260 | xmin = -10.; |
2261 | xmax = 10.; |
2262 | break; |
2263 | |
2264 | case 4: |
2265 | label = "Jet p_{T} corrected with Local Rho"; |
e56c3f6e |
2266 | nbins = 50; // 250 |
f9b4eee4 |
2267 | xmin = -50.; |
2268 | xmax = 200.; |
2269 | break; |
2270 | |
2271 | case 5: |
2272 | label = "Jet p_{T} corrected with Global Rho"; |
e56c3f6e |
2273 | nbins = 50; // 250 |
f9b4eee4 |
2274 | xmin = -50.; |
2275 | xmax = 200.; |
2276 | break; |
2277 | |
daa728b0 |
2278 | }// end of switch |
2279 | } // end of Correction (ME) sparse |
2280 | |
56523939 |
2281 | //________________________________________________________________________ |
2282 | //Int_t AliAnalysisTaskEmcalJetHadEPpid::AcceptFlavourJet(AliEmcalJet* fljet, Int_t NUM, Int_t NUM2, Int_t NUM3) { |
2283 | Int_t AliAnalysisTaskEmcalJetHadEPpid::AcceptFlavourJet(AliEmcalJet* fljet, Int_t NUM) { |
2284 | // Get jet if accepted for given flavour tag |
2285 | // If jet not accepted return 0 |
2286 | if(!fljet) { |
2287 | AliError(Form("%s:Jet not found",GetName())); |
2288 | return 0; |
2289 | } |
2290 | |
fdf5fcfe |
2291 | Int_t flavNUM = -99;//, flavNUM2 = -99, flavNUM3 = -99; FIXME commented out to avoid compiler warning |
56523939 |
2292 | flavNUM = NUM; |
2293 | //flavNUM2 = NUM2; |
2294 | //flavNUM3 = NUM3; |
daa728b0 |
2295 | |
56523939 |
2296 | /* |
2297 | // from the AliEmcalJet class, the tagging enumerator |
2298 | enum EFlavourTag{ |
2299 | kDStar = 1<<0, kD0 = 1<<1, |
2300 | kSig1 = 1<<2, kSig2 = 1<<3, |
2301 | kBckgrd1 = 1<<4, kBckgrd2 = 1<<5, kBckgrd3 = 1<<6 |
2302 | }; |
2303 | // bit 0 = no tag, bit 1 = Dstar, bit 2 = D0 and so forth... |
2304 | */ |
2305 | |
2306 | // get flavour of jet (if any) |
2307 | Int_t flav = -999; |
2308 | flav = fljet->GetFlavour(); |
2309 | |
2310 | // cases (for now..) |
2311 | // 3 = electron rich, 5 = hadron (bkgrd) rich |
2312 | // if flav < 1, its not tagged, so return kFALSE (0) |
2313 | if(flav < 1) return 0; |
2314 | |
2315 | // if flav is not equal to what we want then return kFALSE (0) |
2316 | //if(flav != flavNUM && flav != flavNUM2 && flav != flavNUM3) return 0; |
2317 | if(flav != flavNUM) return 0; |
2318 | |
2319 | // we have the flavour we want, so return kTRUE (1) |
2320 | //if(flav == flavNUM || flav == flavNUM2 || flav == flavNUM3) return 1; |
2321 | if(flav == flavNUM) return 1; |
2322 | |
2323 | // we by default have a flavour tagged jet |
2324 | return 1; |
2325 | } |
a4248490 |
2326 | |
2327 | //________________________________________________________________________ |
2328 | Double_t AliAnalysisTaskEmcalJetHadEPpid::EffCorrection(Double_t trackETA, Double_t trackPT, Int_t effSwitch) const { |
2329 | // default (current) parameters |
2330 | /* |
2331 | // set min/max range of x & y (track Eta & track pt) |
2332 | Double_t trETAmin = -0.9; Double_t trETAmax = 0.9; |
2333 | Double_t trPTmin = 0.0; Double_t trPTmax = 20.0; |
2334 | Int_t npar = 10; |
2335 | //TF2 TRefficiency = TRefficiency = new TF2("eff", EffFunctionCorr, trETAmin, trETAmax, trPTmin, trPTmax, npar); |
2336 | TRefficiency = new TF2("eff", "([0]*exp(-pow([1]/x,[2])) + [3]*x) * ( (y>-0.07)*([4] + [5]*y + [6]*y*y) + (y<=-0.07)*([7] + [8]*y + [9]*y*y) ) ", trETAmin, trETAmax, trPTmin, trPTmax, npar); |
2337 | */ |
2338 | |
bd67d847 |
2339 | Int_t runNUM = fCurrentRunNumber; |
2340 | |
2341 | if ((runNUM == 169975 || runNUM == 169981 || runNUM == 170038 || runNUM == 170040 || runNUM == 170083 || runNUM == 170084 || runNUM == 170085 || runNUM == 170088 || runNUM == 170089 || runNUM == 170091 || runNUM == 170155 || runNUM == 170159 || runNUM == 170163 || runNUM == 170193 || runNUM == 170195 || runNUM == 170203 || runNUM == 170204 || runNUM == 170228 || runNUM == 170230 || runNUM == 170268 || runNUM == 170269 || runNUM == 170270 || runNUM == 170306 || runNUM == 170308 || runNUM == 170309) && fCent < 30 && effSwitch!=1 && effSwitch!=6) effSwitch = 2; |
2342 | |
2343 | if ((runNUM == 169975 || runNUM == 169981 || runNUM == 170038 || runNUM == 170040 || runNUM == 170083 || runNUM == 170084 || runNUM == 170085 || runNUM == 170088 || runNUM == 170089 || runNUM == 170091 || runNUM == 170155 || runNUM == 170159 || runNUM == 170163 || runNUM == 170193 || runNUM == 170195 || runNUM == 170203 || runNUM == 170204 || runNUM == 170228 || runNUM == 170230 || runNUM == 170268 || runNUM == 170269 || runNUM == 170270 || runNUM == 170306 || runNUM == 170308 || runNUM == 170309) && fCent > 30 && effSwitch!=1 && effSwitch!=6) effSwitch = 3; |
2344 | |
2345 | if ((runNUM == 167903 || runNUM == 167915 || runNUM == 167987 || runNUM == 167988 || runNUM == 168066 || runNUM == 168068 || runNUM == 168069 || runNUM == 168076 || runNUM == 168104 || runNUM == 168107 || runNUM == 168108 || runNUM == 168115 || runNUM == 168212 || runNUM == 168310 || runNUM == 168311 || runNUM == 168322 || runNUM == 168325 || runNUM == 168341 || runNUM == 168342 || runNUM == 168361 || runNUM == 168362 || runNUM == 168458 || runNUM == 168460 || runNUM == 168461 || runNUM == 168464 || runNUM == 168467 || runNUM == 168511 || runNUM == 168512 || runNUM == 168777 || runNUM == 168826 || runNUM == 168984 || runNUM == 168988 || runNUM == 168992 || runNUM == 169035 || runNUM == 169091 || runNUM == 169094 || runNUM == 169138 || runNUM == 169143 || runNUM == 169144 || runNUM == 169145 || runNUM == 169148 || runNUM == 169156 || runNUM == 169160 || runNUM == 169167 || runNUM == 169238 || runNUM == 169411 || runNUM == 169415 || runNUM == 169417 || runNUM == 169835 || runNUM == 169837 || runNUM == 169838 || runNUM == 169846 || runNUM == 169855 || runNUM == 169858 || runNUM == 169859 || runNUM == 169923 || runNUM == 169956 || runNUM == 170027 || runNUM == 170036 || runNUM == 170081) && fCent < 30 && effSwitch!=1 && effSwitch!=6) effSwitch = 4; |
2346 | |
2347 | if ((runNUM == 167903 || runNUM == 167915 || runNUM == 167987 || runNUM == 167988 || runNUM == 168066 || runNUM == 168068 || runNUM == 168069 || runNUM == 168076 || runNUM == 168104 || runNUM == 168107 || runNUM == 168108 || runNUM == 168115 || runNUM == 168212 || runNUM == 168310 || runNUM == 168311 || runNUM == 168322 || runNUM == 168325 || runNUM == 168341 || runNUM == 168342 || runNUM == 168361 || runNUM == 168362 || runNUM == 168458 || runNUM == 168460 || runNUM == 168461 || runNUM == 168464 || runNUM == 168467 || runNUM == 168511 || runNUM == 168512 || runNUM == 168777 || runNUM == 168826 || runNUM == 168984 || runNUM == 168988 || runNUM == 168992 || runNUM == 169035 || runNUM == 169091 || runNUM == 169094 || runNUM == 169138 || runNUM == 169143 || runNUM == 169144 || runNUM == 169145 || runNUM == 169148 || runNUM == 169156 || runNUM == 169160 || runNUM == 169167 || runNUM == 169238 || runNUM == 169411 || runNUM == 169415 || runNUM == 169417 || runNUM == 169835 || runNUM == 169837 || runNUM == 169838 || runNUM == 169846 || runNUM == 169855 || runNUM == 169858 || runNUM == 169859 || runNUM == 169923 || runNUM == 169956 || runNUM == 170027 || runNUM == 170036 || runNUM == 170081) && fCent > 30 && effSwitch!=1 && effSwitch!=6) effSwitch = 5; |
2348 | |
a4248490 |
2349 | // set parameter values |
2350 | Double_t p[10] = {9.30132e-01, 9.78828e-02, 1.50626e+00, -1.30398e-02, 9.04935e-01, 2.93456e-01, -3.84654e-01, 8.64688e-01, -2.59244e-01, -3.50544e-01}; |
bd67d847 |
2351 | Double_t m[10] = {9.64508e-01, 9.53247e-01, 3.34584e+01, -4.44705e-03, 8.42526e-01, 2.64031e-01, -3.37155e-01, 8.00933e-01, -2.52274e-01, -3.32314e-01}; |
2352 | Double_t r[10] = {9.32462e-01, 9.65724e-02, 1.46878e+00, -1.58579e-02, 9.46241e-01, 1.98548e-01, -2.93280e-01, 9.05743e-01, -1.55519e-01, -2.54959e-01}; |
2353 | Double_t n[9] = {1.01595e+00, 4.46603e-05, 1.53587e+00, -1.14282e-03, 7.87525e-01, 3.35267e-01, -5.68939e-01, 5.19236e-01, -5.20195e-01}; |
2354 | Double_t q[9] = {19.72915e-01, 6.38626e-01, 8.47813e+00, -1.03864e-02, 9.01610e-01, 9.71765e-02, -4.94174e-01, 1.11112e+00, -9.64431e-01}; |
a4248490 |
2355 | |
2356 | // x-variable = track pt, y-variable = track eta |
2357 | Double_t x = trackPT; |
2358 | Double_t y = trackETA; |
2359 | Double_t TRefficiency = -999; |
2360 | |
a4248490 |
2361 | // set up a switch for different parameter values... |
2362 | switch(effSwitch) { |
2363 | case 1 : |
2364 | // first switch value, initial set of parameters on low statistics data for Semi-Good (LHC11h) Run |
2365 | //Double_t p[] = {9.30132e-01, 9.78828e-02, 1.50626e+00, -1.30398e-02, 9.04935e-01, 2.93456e-01, -3.84654e-01, 8.64688e-01, -2.59244e-01, -3.50544e-01}; |
2366 | |
2367 | // calculate tracking efficiency as function of track pt and track eta |
2368 | TRefficiency = (p[0]*exp(-pow(p[1]/x,p[2])) + p[3]*x) * ( (y>-0.07)*(p[4] + p[5]*y + p[6]*y*y) + (y<=-0.07)*(p[7] + p[8]*y + p[9]*y*y) ); |
2369 | break; |
2370 | |
2371 | case 2 : |
bd67d847 |
2372 | // Parameter values for Semi-GOOD TPC (LHC11h) runs (0-30%): |
a4248490 |
2373 | TRefficiency = (m[0]*exp(-pow(m[1]/x,m[2])) + m[3]*x) * ( (y>-0.07)*(m[4] + m[5]*y + m[6]*y*y) + (y<=-0.07)*(m[7] + m[8]*y + m[9]*y*y) ); |
a4248490 |
2374 | break; |
2375 | |
2376 | case 3 : |
bd67d847 |
2377 | // Parameter values for Semi-GOOD TPC (LHC11h) runs (30-90%): |
2378 | TRefficiency = (r[0]*exp(-pow(r[1]/x,r[2])) + r[3]*x) * ( (y>-0.07)*(r[4] + r[5]*y + r[6]*y*y) + (y<=-0.07)*(r[7] + r[8]*y + r[9]*y*y) ); |
a4248490 |
2379 | break; |
2380 | |
2381 | case 4 : |
bd67d847 |
2382 | // Parameter values for GOOD TPC (LHC11h) runs (0-30%): |
2383 | TRefficiency = (n[0]*exp(-pow(n[1]/x,n[2])) + n[3]*x) * (n[4] + n[5]*y + n[6]*y*y + n[7]*y*y*y + n[8]*y*y*y*y); |
2384 | break; |
2385 | |
2386 | case 5 : |
2387 | // Parameter values for GOOD TPC (LHC11h) runs (30-90%): |
2388 | TRefficiency = (q[0]*exp(-pow(q[1]/x,q[2])) + q[3]*x) * (q[4] + q[5]*y + q[6]*y*y + q[7]*y*y*y + q[8]*y*y*y*y); |
2389 | break; |
2390 | |
2391 | case 6 : |
a4248490 |
2392 | // Parameter values for tweaking efficiency function.. |
2393 | TRefficiency = (q[0]*exp(-pow(q[1]/x,q[2])) + q[3]*x) * ( (y>-0.07)*(q[4] + q[5]*y + q[6]*y*y) + (y<=-0.07)*(q[7] + q[8]*y + q[9]*y*y) ); |
bd67d847 |
2394 | TRefficiency = 1.0; // for now |
a4248490 |
2395 | break; |
2396 | |
2397 | default : |
2398 | // no Efficiency Switch option selected.. therefore don't correct, and set eff = 1 |
2399 | TRefficiency = 1.0; |
2400 | |
2401 | } |
2402 | |
2403 | return TRefficiency; |
2404 | } |
2405 | |