]>
Commit | Line | Data |
---|---|---|
1 | // $Id$ | |
2 | ||
3 | ////////// | |
4 | //Measure Jet-hadron correlations | |
5 | //Does event Mixing using AliEventPoolManager | |
6 | ///////// | |
7 | ||
8 | #include "AliAnalysisTaskEmcalJetHMEC.h" | |
9 | ||
10 | #include "TChain.h" | |
11 | #include "TTree.h" | |
12 | #include "TList.h" | |
13 | #include "TH1F.h" | |
14 | #include "TH2F.h" | |
15 | #include "THnSparse.h" | |
16 | #include "TCanvas.h" | |
17 | #include <TClonesArray.h> | |
18 | #include <TParticle.h> | |
19 | #include "AliVTrack.h" | |
20 | #include "TParameter.h" | |
21 | ||
22 | #include "AliAODEvent.h" | |
23 | #include "AliAnalysisTask.h" | |
24 | #include "AliAnalysisManager.h" | |
25 | ||
26 | #include "AliESDEvent.h" | |
27 | #include "AliESDInputHandler.h" | |
28 | #include "AliESDCaloCluster.h" | |
29 | #include "AliESDVertex.h" | |
30 | #include "AliCentrality.h" | |
31 | #include "AliAODJet.h" | |
32 | #include "AliEmcalJet.h" | |
33 | #include "AliESDtrackCuts.h" | |
34 | ||
35 | #include "TVector3.h" | |
36 | #include "AliPicoTrack.h" | |
37 | #include "AliEventPoolManager.h" | |
38 | ||
39 | ClassImp(AliAnalysisTaskEmcalJetHMEC) | |
40 | ||
41 | //________________________________________________________________________ | |
42 | AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC() : | |
43 | AliAnalysisTaskEmcalJet("HMEC",kFALSE), | |
44 | fTracksName(""), | |
45 | fJetsName(""), | |
46 | fPhimin(-10), | |
47 | fPhimax(10), | |
48 | fEtamin(-0.9), | |
49 | fEtamax(0.9), | |
50 | fAreacut(0.0), | |
51 | fTrkBias(5), | |
52 | fClusBias(5), | |
53 | fTrkEta(0.9), | |
54 | fDoEventMixing(0), | |
55 | fMixingTracks(50000), | |
56 | fESD(0), | |
57 | fAOD(0), | |
58 | fPoolMgr(0x0), | |
59 | fHistTrackPt(0), | |
60 | fHistCentrality(0), | |
61 | fHistJetEtaPhi(0), | |
62 | fHistJetHEtaPhi(0), | |
63 | fhnMixedEvents(0x0), | |
64 | fhnJH(0x0) | |
65 | { | |
66 | // Default Constructor | |
67 | ||
68 | for(Int_t ipta=0; ipta<7; ipta++){ | |
69 | fHistTrackEtaPhi[ipta]=0; | |
70 | } | |
71 | ||
72 | for(Int_t icent = 0; icent<6; ++icent){ | |
73 | fHistJetPt[icent]=0; | |
74 | fHistJetPtBias[icent]=0; | |
75 | fHistLeadJetPt[icent]=0; | |
76 | fHistLeadJetPtBias[icent]=0; | |
77 | fHistJetPtTT[icent]=0; | |
78 | for(Int_t iptjet = 0; iptjet<5; ++iptjet){ | |
79 | for(Int_t ieta = 0; ieta<3; ++ieta){ | |
80 | fHistJetH[icent][iptjet][ieta]=0; | |
81 | fHistJetHBias[icent][iptjet][ieta]=0; | |
82 | fHistJetHTT[icent][iptjet][ieta]=0; | |
83 | } | |
84 | } | |
85 | } | |
86 | ||
87 | } | |
88 | //________________________________________________________________________ | |
89 | AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC(const char *name) : | |
90 | AliAnalysisTaskEmcalJet(name,kTRUE), | |
91 | fTracksName(""), | |
92 | fJetsName(""), | |
93 | fPhimin(-10), | |
94 | fPhimax(10), | |
95 | fEtamin(-0.9), | |
96 | fEtamax(0.9), | |
97 | fAreacut(0.0), | |
98 | fTrkBias(5), | |
99 | fClusBias(5), | |
100 | fTrkEta(0.9), | |
101 | fDoEventMixing(0), | |
102 | fMixingTracks(50000), | |
103 | fESD(0), | |
104 | fAOD(0), | |
105 | fPoolMgr(0x0), | |
106 | fHistTrackPt(0), | |
107 | fHistCentrality(0), | |
108 | fHistJetEtaPhi(0), | |
109 | fHistJetHEtaPhi(0), | |
110 | fhnMixedEvents(0x0), | |
111 | fhnJH(0x0) | |
112 | { | |
113 | // Constructor | |
114 | for(Int_t ipta=0; ipta<7; ipta++){ | |
115 | fHistTrackEtaPhi[ipta]=0; | |
116 | } | |
117 | for(Int_t icent = 0; icent<6; ++icent){ | |
118 | fHistJetPt[icent]=0; | |
119 | fHistJetPtBias[icent]=0; | |
120 | fHistLeadJetPt[icent]=0; | |
121 | fHistLeadJetPtBias[icent]=0; | |
122 | fHistJetPtTT[icent]=0; | |
123 | for(Int_t iptjet = 0; iptjet<5; ++iptjet){ | |
124 | for(Int_t ieta = 0; ieta<3; ++ieta){ | |
125 | fHistJetH[icent][iptjet][ieta]=0; | |
126 | fHistJetHBias[icent][iptjet][ieta]=0; | |
127 | fHistJetHTT[icent][iptjet][ieta]=0; | |
128 | } | |
129 | } | |
130 | } | |
131 | ||
132 | } | |
133 | ||
134 | //________________________________________________________________________ | |
135 | void AliAnalysisTaskEmcalJetHMEC::UserCreateOutputObjects() { | |
136 | // Called once | |
137 | AliAnalysisTaskEmcalJet::UserCreateOutputObjects(); | |
138 | OpenFile(1); | |
139 | ||
140 | // Create histograms | |
141 | fHistTrackPt = new TH1F("fHistTrackPt", "P_{T} distribution", 1000, 0.0, 100.0); | |
142 | fHistCentrality = new TH1F("fHistCentrality","centrality",100,0,100); | |
143 | fHistJetEtaPhi = new TH2F("fHistJetEtaPhi","Jet eta-phi",900,-1.8,1.8,720,-3.2,3.2); | |
144 | fHistJetHEtaPhi = new TH2F("fHistJetHEtaPhi","Jet-Hadron deta-dphi",900,-1.8,1.8,720,-1.6,4.8); | |
145 | ||
146 | TString name; | |
147 | ||
148 | for(Int_t ipta=0; ipta<7; ++ipta){ | |
149 | name = Form("fHistTrackEtaPhi_%i", ipta); | |
150 | fHistTrackEtaPhi[ipta] = new TH2F(name,name,400,-1,1,720,0.0,2.0*TMath::Pi()); | |
151 | fOutput->Add(fHistTrackEtaPhi[ipta]); | |
152 | } | |
153 | ||
154 | for(Int_t icent = 0; icent<6; ++icent){ | |
155 | name = Form("fHistJetPt_%i",icent); | |
156 | fHistJetPt[icent] = new TH1F(name,name,200,0,200); | |
157 | fOutput->Add(fHistJetPt[icent]); | |
158 | ||
159 | name = Form("fHistJetPtBias_%i",icent); | |
160 | fHistJetPtBias[icent] = new TH1F(name,name,200,0,200); | |
161 | fOutput->Add(fHistJetPtBias[icent]); | |
162 | ||
163 | name = Form("fHistLeadJetPt_%i",icent); | |
164 | fHistLeadJetPt[icent] = new TH1F(name,name,200,0,200); | |
165 | fOutput->Add(fHistLeadJetPt[icent]); | |
166 | ||
167 | name = Form("fHistLeadJetPtBias_%i",icent); | |
168 | fHistLeadJetPtBias[icent] = new TH1F(name,name,200,0,200); | |
169 | fOutput->Add(fHistLeadJetPtBias[icent]); | |
170 | ||
171 | name = Form("fHistJetPtTT_%i",icent); | |
172 | fHistJetPtTT[icent] = new TH1F(name,name,200,0,200); | |
173 | fOutput->Add(fHistJetPtTT[icent]); | |
174 | ||
175 | for(Int_t iptjet = 0; iptjet<5; ++iptjet){ | |
176 | for(Int_t ieta = 0; ieta<3; ++ieta){ | |
177 | name = Form("fHistJetH_%i_%i_%i",icent,iptjet,ieta); | |
178 | fHistJetH[icent][iptjet][ieta]=new TH2F(name,name,72,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30); | |
179 | fOutput->Add(fHistJetH[icent][iptjet][ieta]); | |
180 | ||
181 | name = Form("fHistJetHBias_%i_%i_%i",icent,iptjet,ieta); | |
182 | fHistJetHBias[icent][iptjet][ieta]=new TH2F(name,name,72,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30); | |
183 | fOutput->Add(fHistJetHBias[icent][iptjet][ieta]); | |
184 | ||
185 | name = Form("fHistJetHTT_%i_%i_%i",icent,iptjet,ieta); | |
186 | fHistJetHTT[icent][iptjet][ieta]=new TH2F(name,name,72,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30); | |
187 | fOutput->Add(fHistJetHTT[icent][iptjet][ieta]); | |
188 | ||
189 | } | |
190 | } | |
191 | } | |
192 | ||
193 | UInt_t cifras = 0; // bit coded, see GetDimParams() below | |
194 | cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8; | |
195 | fhnJH = NewTHnSparseF("fhnJH", cifras); | |
196 | fhnJH->Sumw2(); | |
197 | fOutput->Add(fhnJH); | |
198 | ||
199 | if(fDoEventMixing){ | |
200 | cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8; | |
201 | fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras); | |
202 | fhnMixedEvents->Sumw2(); | |
203 | fOutput->Add(fhnMixedEvents); | |
204 | } | |
205 | ||
206 | fOutput->Add(fHistTrackPt); | |
207 | fOutput->Add(fHistCentrality); | |
208 | fOutput->Add(fHistJetEtaPhi); | |
209 | fOutput->Add(fHistJetHEtaPhi); | |
210 | ||
211 | PostData(1, fOutput); | |
212 | ||
213 | //Event Mixing | |
214 | Int_t trackDepth = fMixingTracks; | |
215 | Int_t poolsize = 1000; // Maximum number of events, ignored in the present implemented of AliEventPoolManager | |
216 | ||
217 | Int_t nZvtxBins = 5+1+5; | |
218 | // bins for second buffer are shifted by 100 cm | |
219 | Double_t vertexBins[] = { -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10, }; | |
220 | Double_t* zvtxbin = vertexBins; | |
221 | ||
222 | Int_t nCentralityBins = 100; | |
223 | Double_t centralityBins[nCentralityBins]; | |
224 | for(Int_t ic=0; ic<nCentralityBins; ic++){ | |
225 | centralityBins[ic]=1.0*ic; | |
226 | } | |
227 | ||
228 | fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, centralityBins, nZvtxBins, zvtxbin); | |
229 | ||
230 | } | |
231 | ||
232 | //________________________________________________________________________ | |
233 | ||
234 | Double_t AliAnalysisTaskEmcalJetHMEC:: RelativePhi(Double_t mphi,Double_t vphi) { | |
235 | ||
236 | if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi()); | |
237 | else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi()); | |
238 | if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi()); | |
239 | else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi()); | |
240 | Double_t dphi = vphi-mphi; | |
241 | if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi()); | |
242 | else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi()); | |
243 | ||
244 | return dphi;//dphi in [-Pi, Pi] | |
245 | } | |
246 | ||
247 | //________________________________________________________________________ | |
248 | Int_t AliAnalysisTaskEmcalJetHMEC::GetCentBin(Double_t cent) const { | |
249 | // Get centrality bin. | |
250 | ||
251 | Int_t centbin = -1; | |
252 | if (cent>=0 && cent<10) centbin = 0; | |
253 | else if (cent>=10 && cent<20) centbin = 1; | |
254 | else if (cent>=20 && cent<30) centbin = 2; | |
255 | else if (cent>=30 && cent<40) centbin = 3; | |
256 | else if (cent>=40 && cent<50) centbin = 4; | |
257 | else if (cent>=50 && cent<90) centbin = 5; | |
258 | return centbin; | |
259 | } | |
260 | ||
261 | //________________________________________________________________________ | |
262 | Int_t AliAnalysisTaskEmcalJetHMEC::GetEtaBin(Double_t eta) const { | |
263 | // Get eta bin for histos. | |
264 | ||
265 | Int_t etabin = -1; | |
266 | if (TMath::Abs(eta)<=0.4) etabin = 0; | |
267 | else if (TMath::Abs(eta)>0.4 && TMath::Abs(eta)<0.8) etabin = 1; | |
268 | else if (TMath::Abs(eta)>=0.8) etabin = 2; | |
269 | return etabin; | |
270 | } | |
271 | //________________________________________________________________________ | |
272 | Int_t AliAnalysisTaskEmcalJetHMEC::GetpTjetBin(Double_t pt) const | |
273 | { | |
274 | // Get jet pt bin for histos. | |
275 | ||
276 | Int_t ptbin = -1; | |
277 | if (pt>=15 && pt<20) ptbin = 0; | |
278 | else if (pt>=20 && pt<25) ptbin = 1; | |
279 | else if (pt>=25 && pt<30) ptbin = 2; | |
280 | else if (pt>=30 && pt<60) ptbin = 3; | |
281 | else if (pt>=60) ptbin = 4; | |
282 | ||
283 | return ptbin; | |
284 | } | |
285 | ||
286 | //________________________________________________________________________ | |
287 | void AliAnalysisTaskEmcalJetHMEC::ExecOnce() { | |
288 | AliAnalysisTaskEmcalJet::ExecOnce(); | |
289 | ||
290 | } | |
291 | ||
292 | //________________________________________________________________________ | |
293 | Bool_t AliAnalysisTaskEmcalJetHMEC::Run() { | |
294 | // Main loop called for each event | |
295 | if(!fTracks){ | |
296 | AliError(Form("No fTracks object!!\n")); | |
297 | return kTRUE; | |
298 | } | |
299 | if(!fJets){ | |
300 | AliError(Form("No fJets object!!\n")); | |
301 | return kTRUE; | |
302 | } | |
303 | ||
304 | // what kind of event do we have: AOD or ESD? | |
305 | Bool_t esdMode = kTRUE; | |
306 | if (dynamic_cast<AliAODEvent*>(InputEvent())) esdMode = kFALSE; | |
307 | ||
308 | // if we have ESD event, set up ESD object | |
309 | if(esdMode){ | |
310 | fESD = dynamic_cast<AliESDEvent*>(InputEvent()); | |
311 | if (!fESD) { | |
312 | AliError(Form("ERROR: fESD not available\n")); | |
313 | return kTRUE; | |
314 | } | |
315 | } | |
316 | ||
317 | // if we have AOD event, set up AOD object | |
318 | if(!esdMode){ | |
319 | fAOD = dynamic_cast<AliAODEvent*>(InputEvent()); | |
320 | if(!fAOD) { | |
321 | AliError(Form("ERROR: fAOD not available\n")); | |
322 | return kTRUE; | |
323 | } | |
324 | } | |
325 | ||
326 | TList *list = InputEvent()->GetList(); | |
327 | if(!list) { | |
328 | AliError(Form("ERROR: list not attached\n")); | |
329 | return kTRUE; | |
330 | } | |
331 | ||
332 | // get centrality | |
333 | if (fCent<0) { | |
334 | AliError(Form("Centrality negative: %f", fCent)); | |
335 | return kTRUE; | |
336 | } | |
337 | ||
338 | Int_t centbin = GetCentBin(fCent); | |
339 | if(centbin<0) return kTRUE; | |
340 | ||
341 | Double_t fvertex[3]={0,0,0}; | |
342 | InputEvent()->GetPrimaryVertex()->GetXYZ(fvertex); | |
343 | Double_t zVtx=fvertex[2]; | |
344 | ||
345 | if(fabs(zVtx)>10.0) return kTRUE; | |
346 | ||
347 | fHistCentrality->Fill(fCent); | |
348 | ||
349 | TClonesArray *jets = 0; | |
350 | TClonesArray *tracks = 0; | |
351 | ||
352 | tracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracks)); | |
353 | if (!tracks) { | |
354 | AliError(Form("Pointer to tracks %s == 0", fTracks->GetName() )); | |
355 | return kTRUE; | |
356 | } | |
357 | const Int_t Ntracks=tracks->GetEntries(); | |
358 | ||
359 | jets= dynamic_cast<TClonesArray*>(list->FindObject(fJets)); | |
360 | if (!jets) { | |
361 | AliError(Form("Pointer to tracks %s == 0", fJets->GetName() )); | |
362 | return kTRUE; | |
363 | } | |
364 | const Int_t Njets = jets->GetEntries(); | |
365 | ||
366 | //Leticia's loop to find hardest track | |
367 | Int_t iTT=-1; | |
368 | Double_t ptmax=-10; | |
369 | ||
370 | for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++) { | |
371 | AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks)); | |
372 | if (!track) { | |
373 | printf("ERROR: Could not receive track %d\n", iTracks); | |
374 | continue; | |
375 | } | |
376 | ||
377 | if(TMath::Abs(track->Eta())>0.9) continue; | |
378 | if(track->Pt()<0.15) continue; | |
379 | //iCount++; | |
380 | if(track->Pt()>ptmax){ | |
381 | ptmax=track->Pt(); | |
382 | iTT=iTracks; | |
383 | } | |
384 | } | |
385 | ||
386 | Int_t ijethi=-1; | |
387 | Double_t highestjetpt=0.0; | |
388 | Int_t passedTTcut=0; | |
389 | ||
390 | for (Int_t ijets = 0; ijets < Njets; ijets++){ | |
391 | AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijets)); | |
392 | if (!jet) continue; | |
393 | if(!AcceptthisJet(jet)) continue; | |
394 | ||
395 | Double_t jetPt = jet->Pt(); | |
396 | ||
397 | if(highestjetpt<jetPt){ | |
398 | ijethi=ijets; | |
399 | highestjetpt=jetPt; | |
400 | } | |
401 | } | |
402 | ||
403 | for (Int_t ijet = 0; ijet < Njets; ijet++){ | |
404 | AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet)); | |
405 | if (!jet) continue; | |
406 | if(!AcceptthisJet(jet)) continue; | |
407 | ||
408 | Double_t jetphi = jet->Phi(); | |
409 | Double_t jetPt = jet->Pt(); | |
410 | Double_t jeteta=jet->Eta(); | |
411 | ||
412 | Double_t leadjet=0; | |
413 | if (ijet==ijethi) leadjet=1; | |
414 | ||
415 | fHistJetPt[centbin]->Fill(jet->Pt()); | |
416 | fHistLeadJetPt[centbin]->Fill(jet->Pt()); | |
417 | ||
418 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){ | |
419 | fHistJetPtBias[centbin]->Fill(jet->Pt()); | |
420 | fHistLeadJetPtBias[centbin]->Fill(jet->Pt()); | |
421 | } | |
422 | ||
423 | fHistJetEtaPhi->Fill(jet->Eta(),jetphi); | |
424 | ||
425 | if(iTT>0){ | |
426 | AliVTrack* TT = static_cast<AliVTrack*>(tracks->At(iTT)); | |
427 | if(TMath::Abs(jetphi-TT->Phi()-TMath::Pi())<0.6) passedTTcut=1; | |
428 | else passedTTcut=0; | |
429 | } | |
430 | ||
431 | if(passedTTcut) | |
432 | fHistJetPtTT[centbin]->Fill(jet->Pt()); | |
433 | ||
434 | Int_t iptjet=-1; | |
435 | iptjet=GetpTjetBin(jetPt); | |
436 | if(iptjet<0) continue; | |
437 | ||
438 | //if (highestjetpt>15) { | |
439 | if (highestjetpt>8) { | |
440 | ||
441 | for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++) { | |
442 | AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks)); | |
443 | if (!track) { | |
444 | printf("ERROR: Could not receive track %d\n", iTracks); | |
445 | continue; | |
446 | } | |
447 | ||
448 | if(TMath::Abs(track->Eta())>fTrkEta) continue; | |
449 | ||
450 | fHistTrackPt->Fill(track->Pt()); | |
451 | ||
452 | if (track->Pt()<0.15) continue; | |
453 | ||
454 | Double_t trackphi = track->Phi(); | |
455 | if (trackphi > TMath::Pi()) | |
456 | trackphi = trackphi-2*TMath::Pi(); | |
457 | ||
458 | Double_t tracketa=track->Eta(); | |
459 | Double_t trackpt=track->Pt(); | |
460 | Double_t deta=tracketa-jeteta; | |
461 | Int_t ieta=GetEtaBin(deta); | |
462 | if (ieta<0) { | |
463 | AliError(Form("Eta Bin negative: %f", deta)); | |
464 | continue; | |
465 | } | |
466 | ||
467 | //Jet pt, track pt, dPhi,deta,fCent | |
468 | Double_t dphijh = RelativePhi(jetphi,trackphi); | |
469 | if (dphijh < -0.5*TMath::Pi()) | |
470 | dphijh+= 2*TMath::Pi(); | |
471 | if (dphijh > 1.5*TMath::Pi()) | |
472 | dphijh-=2.*TMath::Pi(); | |
473 | ||
474 | fHistJetH[centbin][iptjet][ieta]->Fill(dphijh,track->Pt()); | |
475 | fHistJetHEtaPhi->Fill(deta,dphijh); | |
476 | ||
477 | Double_t dR=sqrt(deta*deta+dphijh*dphijh); | |
478 | ||
479 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){ | |
480 | fHistJetHBias[centbin][iptjet][ieta]->Fill(dphijh,trackpt); | |
481 | ||
482 | Double_t triggerEntries[8] = {fCent,jetPt,trackpt,dR,deta,dphijh,0.0,leadjet}; | |
483 | fhnJH->Fill(triggerEntries); | |
484 | } | |
485 | ||
486 | if(passedTTcut) | |
487 | fHistJetHTT[centbin][iptjet][ieta]->Fill(dphijh,trackpt); | |
488 | ||
489 | } //track loop | |
490 | }//jet pt cut | |
491 | }//jet loop | |
492 | ||
493 | //Prepare to do event mixing | |
494 | ||
495 | // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool | |
496 | TObjArray* tracksClone = CloneAndReduceTrackList(tracks); | |
497 | //delete tracks; | |
498 | ||
499 | if(fDoEventMixing>0){ | |
500 | ||
501 | // event mixing | |
502 | ||
503 | // 1. First get an event pool corresponding in mult (cent) and | |
504 | // zvertex to the current event. Once initialized, the pool | |
505 | // should contain nMix (reduced) events. This routine does not | |
506 | // pre-scan the chain. The first several events of every chain | |
507 | // will be skipped until the needed pools are filled to the | |
508 | // specified depth. If the pool categories are not too rare, this | |
509 | // should not be a problem. If they are rare, you could lose | |
510 | // statistics. | |
511 | ||
512 | // 2. Collect the whole pool's content of tracks into one TObjArray | |
513 | // (bgTracks), which is effectively a single background super-event. | |
514 | ||
515 | // 3. The reduced and bgTracks arrays must both be passed into | |
516 | // FillCorrelations(). Also nMix should be passed in, so a weight | |
517 | // of 1./nMix can be applied. | |
518 | ||
519 | AliEventPool* pool = fPoolMgr->GetEventPool(fCent, zVtx); | |
520 | ||
521 | if (!pool){ | |
522 | AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCent, zVtx)); | |
523 | return kTRUE; | |
524 | } | |
525 | ||
526 | //check for a trigger jet | |
527 | if (pool->IsReady() || pool->NTracksInPool() > fMixingTracks / 10 || pool->GetCurrentNEvents() >= 5) { | |
528 | ||
529 | for (Int_t ijet = 0; ijet < Njets; ijet++){ | |
530 | Double_t leadjet=0; | |
531 | if (ijet==ijethi) leadjet=1; | |
532 | ||
533 | AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet)); | |
534 | if(!AcceptthisJet(jet)) continue; | |
535 | ||
536 | Double_t jetPt = jet->Pt(); | |
537 | Double_t jetphi = jet->Phi(); | |
538 | Double_t jeteta=jet->Eta(); | |
539 | ||
540 | Int_t nMix = pool->GetCurrentNEvents(); | |
541 | ||
542 | //Fill for biased jet triggers only | |
543 | if ((jet->MaxTrackPt()>fTrkBias) || (jet->MaxClusterPt()>fClusBias)){ | |
544 | ||
545 | // Fill mixed-event histos here | |
546 | for (Int_t jMix=0; jMix<nMix; jMix++) { | |
547 | TObjArray* bgTracks = pool->GetEvent(jMix); | |
548 | const Int_t Nbgtrks = bgTracks->GetEntries(); | |
549 | ||
550 | for(Int_t ibg=0; ibg<Nbgtrks; ibg++){ | |
551 | AliPicoTrack *part = static_cast<AliPicoTrack*>(bgTracks->At(ibg)); | |
552 | if(!part) continue; | |
553 | ||
554 | Double_t DEta = part->Eta()-jeteta; | |
555 | Double_t DPhi = RelativePhi(jetphi,part->Phi()); | |
556 | ||
557 | Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta); | |
558 | if(DPhi<-0.5*TMath::Pi()) DPhi+=2.*TMath::Pi(); | |
559 | if(DPhi>3./2.*TMath::Pi()) DPhi-=2.*TMath::Pi(); | |
560 | Double_t triggerEntries[8] = {fCent,jetPt,part->Pt(),DR,DEta,DPhi,0.0,leadjet}; | |
561 | fhnMixedEvents->Fill(triggerEntries,1./nMix); | |
562 | } | |
563 | } | |
564 | } | |
565 | } | |
566 | } | |
567 | ||
568 | //update pool if jet in event or not | |
569 | pool->UpdatePool(tracksClone); | |
570 | } | |
571 | ||
572 | return kTRUE; | |
573 | } | |
574 | ||
575 | //________________________________________________________________________ | |
576 | void AliAnalysisTaskEmcalJetHMEC::Terminate(Option_t *) | |
577 | { | |
578 | //just terminate | |
579 | } | |
580 | ||
581 | //________________________________________________________________________ | |
582 | Int_t AliAnalysisTaskEmcalJetHMEC::AcceptthisJet(AliEmcalJet *jet) | |
583 | { | |
584 | //applies all jet cuts except pt | |
585 | float jetphi = jet->Phi(); | |
586 | if (jetphi>TMath::Pi()) | |
587 | jetphi = jetphi-2*TMath::Pi(); | |
588 | ||
589 | if ((jet->Phi()<fPhimin)||(jet->Phi()>fPhimax)) | |
590 | return 0; | |
591 | if ((jet->Eta()<fEtamin)||(jet->Eta()>fEtamax)) | |
592 | return 0; | |
593 | if (jet->Area()<fAreacut) | |
594 | return 0; | |
595 | //prevents 0 area jets from sneaking by when area cut == 0 | |
596 | if (jet->Area()==0) | |
597 | return 0; | |
598 | //exclude jets with extremely high pt tracks which are likely misreconstructed | |
599 | if(jet->MaxTrackPt()>100) | |
600 | return 0; | |
601 | ||
602 | //passed all above cuts | |
603 | return 1; | |
604 | } | |
605 | ||
606 | //________________________________________________________________________ | |
607 | ||
608 | THnSparse* AliAnalysisTaskEmcalJetHMEC::NewTHnSparseF(const char* name, UInt_t entries){ | |
609 | // generate new THnSparseF, axes are defined in GetDimParams() | |
610 | ||
611 | Int_t count = 0; | |
612 | UInt_t tmp = entries; | |
613 | while(tmp!=0){ | |
614 | count++; | |
615 | tmp = tmp &~ -tmp; // clear lowest bit | |
616 | } | |
617 | ||
618 | TString hnTitle(name); | |
619 | const Int_t dim = count; | |
620 | Int_t nbins[dim]; | |
621 | Double_t xmin[dim]; | |
622 | Double_t xmax[dim]; | |
623 | ||
624 | Int_t i=0; | |
625 | Int_t c=0; | |
626 | while(c<dim && i<32){ | |
627 | if(entries&(1<<i)){ | |
628 | ||
629 | TString label(""); | |
630 | GetDimParams(i, label, nbins[c], xmin[c], xmax[c]); | |
631 | hnTitle += Form(";%s",label.Data()); | |
632 | c++; | |
633 | } | |
634 | ||
635 | i++; | |
636 | } | |
637 | hnTitle += ";"; | |
638 | ||
639 | return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax); | |
640 | } | |
641 | ||
642 | void AliAnalysisTaskEmcalJetHMEC::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax) | |
643 | { | |
644 | // stores label and binning of axis for THnSparse | |
645 | ||
646 | const Double_t pi = TMath::Pi(); | |
647 | ||
648 | switch(iEntry){ | |
649 | ||
650 | case 0: | |
651 | label = "V0 centrality (%)"; | |
652 | nbins = 10; | |
653 | xmin = 0.; | |
654 | xmax = 100.; | |
655 | break; | |
656 | ||
657 | case 1: | |
658 | label = "corrected jet pt"; | |
659 | nbins = 20; | |
660 | xmin = 0.; | |
661 | xmax = 200.; | |
662 | break; | |
663 | ||
664 | case 2: | |
665 | label = "track pT"; | |
666 | nbins = 100; | |
667 | xmin = 0.; | |
668 | xmax = 10; | |
669 | break; | |
670 | ||
671 | case 3: | |
672 | label = "deltaR"; | |
673 | nbins = 10; | |
674 | xmin = 0.; | |
675 | xmax = 5.0; | |
676 | break; | |
677 | ||
678 | case 4: | |
679 | label = "deltaEta"; | |
680 | nbins = 24; | |
681 | xmin = -1.2; | |
682 | xmax = 1.2; | |
683 | break; | |
684 | ||
685 | case 5: | |
686 | label = "deltaPhi"; | |
687 | nbins = 72; | |
688 | xmin = -0.5*pi; | |
689 | xmax = 1.5*pi; | |
690 | break; | |
691 | ||
692 | case 6: | |
693 | label = "leading track"; | |
694 | nbins = 13; | |
695 | xmin = 0; | |
696 | xmax = 50; | |
697 | break; | |
698 | ||
699 | case 7: | |
700 | ||
701 | label = "trigger track"; | |
702 | nbins =10; | |
703 | xmin = 0; | |
704 | xmax = 50; | |
705 | break; | |
706 | ||
707 | case 8: | |
708 | label = "leading jet"; | |
709 | nbins = 3; | |
710 | xmin = -0.5; | |
711 | xmax = 2.5; | |
712 | break; | |
713 | ||
714 | } | |
715 | } | |
716 | ||
717 | //_________________________________________________ | |
718 | // From CF event mixing code PhiCorrelations | |
719 | TObjArray* AliAnalysisTaskEmcalJetHMEC::CloneAndReduceTrackList(TObjArray* tracks) | |
720 | { | |
721 | // clones a track list by using AliPicoTrack which uses much less memory (used for event mixing) | |
722 | ||
723 | TObjArray* tracksClone = new TObjArray; | |
724 | tracksClone->SetOwner(kTRUE); | |
725 | ||
726 | for (Int_t i=0; i<tracks->GetEntriesFast(); i++) | |
727 | { | |
728 | AliVParticle* particle = (AliVParticle*) tracks->At(i); | |
729 | if(TMath::Abs(particle->Eta())>fTrkEta) continue; | |
730 | if(particle->Pt()<0.15) continue; | |
731 | ||
732 | Double_t trackpt=particle->Pt(); | |
733 | ||
734 | Int_t hadbin=-1; | |
735 | if(trackpt<0.5) hadbin=0; | |
736 | else if(trackpt<1) hadbin=1; | |
737 | else if(trackpt<2) hadbin=2; | |
738 | else if(trackpt<3) hadbin=3; | |
739 | else if(trackpt<5) hadbin=4; | |
740 | else if(trackpt<8) hadbin=5; | |
741 | else if(trackpt<20) hadbin=6; | |
742 | ||
743 | if(hadbin>-1) fHistTrackEtaPhi[hadbin]->Fill(particle->Eta(),particle->Phi()); | |
744 | ||
745 | tracksClone->Add(new AliPicoTrack(particle->Pt(), particle->Eta(), particle->Phi(), particle->Charge(), 0, 0, 0, 0)); | |
746 | } | |
747 | ||
748 | return tracksClone; | |
749 | } |