]>
Commit | Line | Data |
---|---|---|
ad869500 | 1 | |
2 | // ****************************************** | |
3 | // This task computes several jet observables like | |
4 | // the fraction of energy in inner and outer coronnas, | |
5 | // jet-track correlations,triggered jet shapes and | |
6 | // correlation strength distribution of particles inside jets. | |
7 | // Author: lcunquei@cern.ch | |
8 | // ******************************************* | |
9 | ||
10 | ||
11 | /************************************************************************** | |
12 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
13 | * * | |
14 | * Author: The ALICE Off-line Project. * | |
15 | * Contributors are mentioned in the code where appropriate. * | |
16 | * * | |
17 | * Permission to use, copy, modify and distribute this software and its * | |
18 | * documentation strictly for non-commercial purposes is hereby granted * | |
19 | * without fee, provided that the above copyright notice appears in all * | |
20 | * copies and that both the copyright notice and this permission notice * | |
21 | * appear in the supporting documentation. The authors make no claims * | |
22 | * about the suitability of this software for any purpose. It is * | |
23 | * provided "as is" without express or implied warranty. * | |
24 | **************************************************************************/ | |
25 | ||
26 | ||
27 | #include "TChain.h" | |
28 | #include "TTree.h" | |
29 | #include "TMath.h" | |
30 | #include "TH1F.h" | |
31 | #include "TH1D.h" | |
32 | #include "TH2F.h" | |
33 | #include "TH3F.h" | |
34 | #include "THnSparse.h" | |
35 | #include "TCanvas.h" | |
36 | ||
37 | #include "AliLog.h" | |
38 | ||
39 | #include "AliAnalysisTask.h" | |
40 | #include "AliAnalysisManager.h" | |
41 | ||
42 | #include "AliVEvent.h" | |
43 | #include "AliESDEvent.h" | |
44 | #include "AliESDInputHandler.h" | |
45 | #include "AliCentrality.h" | |
46 | #include "AliAnalysisHelperJetTasks.h" | |
47 | #include "AliInputEventHandler.h" | |
48 | #include "AliAODJetEventBackground.h" | |
49 | #include "AliAODMCParticle.h" | |
50 | //#include "AliAnalysisTaskFastEmbedding.h" | |
51 | #include "AliAODEvent.h" | |
52 | #include "AliAODHandler.h" | |
53 | #include "AliAODJet.h" | |
54 | ||
55 | #include "AliAnalysisTaskJetCorePP.h" | |
56 | ||
57 | using std::cout; | |
58 | using std::endl; | |
59 | ||
60 | ClassImp(AliAnalysisTaskJetCorePP) | |
61 | ||
62 | //Filip Krizek 1st March 2013 | |
63 | ||
64 | //--------------------------------------------------------------------- | |
65 | AliAnalysisTaskJetCorePP::AliAnalysisTaskJetCorePP() : | |
66 | AliAnalysisTaskSE(), | |
67 | fESD(0x0), | |
68 | fAODIn(0x0), | |
69 | fAODOut(0x0), | |
70 | fAODExtension(0x0), | |
71 | fJetBranchName(""), | |
72 | //fListJets(0x0), | |
73 | fNonStdFile(""), | |
74 | fSystem(0), //pp=0 pPb=1 | |
75 | fJetParamR(0.4), | |
76 | fOfflineTrgMask(AliVEvent::kAny), | |
77 | fMinContribVtx(1), | |
78 | fVtxZMin(-10.0), | |
79 | fVtxZMax(10.0), | |
80 | fFilterMask(0), | |
81 | fCentMin(0.0), | |
82 | fCentMax(100.0), | |
83 | fJetEtaMin(-0.5), | |
84 | fJetEtaMax(0.5), | |
85 | fTriggerEtaCut(0.9), | |
86 | fTrackEtaCut(0.9), | |
87 | fTrackLowPtCut(0.15), | |
88 | fOutputList(0x0), | |
89 | fHistEvtSelection(0x0), | |
90 | fh2Ntriggers(0x0), | |
91 | fHJetSpec(0x0), | |
92 | fHJetDensity(0x0), | |
93 | fHJetDensityA4(0x0), | |
94 | fhJetPhi(0x0), | |
95 | fhTriggerPhi(0x0), | |
96 | fhJetEta(0x0), | |
97 | fhTriggerEta(0x0), | |
98 | fhVertexZ(0x0), | |
99 | fhVertexZAccept(0x0), | |
100 | fhContribVtx(0x0), | |
101 | fhContribVtxAccept(0x0), | |
102 | fhDphiTriggerJet(0x0), | |
103 | fhDphiTriggerJetAccept(0x0), | |
104 | fhCentrality(0x0), | |
105 | fhCentralityAccept(0x0), | |
106 | fkAcceptance(2.0*TMath::Pi()*1.8) | |
107 | { | |
108 | // default Constructor | |
109 | fListJets = new TList(); | |
110 | } | |
111 | ||
112 | //--------------------------------------------------------------------- | |
113 | ||
114 | AliAnalysisTaskJetCorePP::AliAnalysisTaskJetCorePP(const char *name) : | |
115 | AliAnalysisTaskSE(name), | |
116 | fESD(0x0), | |
117 | fAODIn(0x0), | |
118 | fAODOut(0x0), | |
119 | fAODExtension(0x0), | |
120 | fJetBranchName(""), | |
121 | //fListJets(0x0), | |
122 | fNonStdFile(""), | |
123 | fSystem(0), //pp=0 pPb=1 | |
124 | fJetParamR(0.4), | |
125 | fOfflineTrgMask(AliVEvent::kAny), | |
126 | fMinContribVtx(1), | |
127 | fVtxZMin(-10.0), | |
128 | fVtxZMax(10.0), | |
129 | fFilterMask(0), | |
130 | fCentMin(0.0), | |
131 | fCentMax(100.0), | |
132 | fJetEtaMin(-0.5), | |
133 | fJetEtaMax(0.5), | |
134 | fTriggerEtaCut(0.9), | |
135 | fTrackEtaCut(0.9), | |
136 | fTrackLowPtCut(0.15), | |
137 | fOutputList(0x0), | |
138 | fHistEvtSelection(0x0), | |
139 | fh2Ntriggers(0x0), | |
140 | fHJetSpec(0x0), | |
141 | fHJetDensity(0x0), | |
142 | fHJetDensityA4(0x0), | |
143 | fhJetPhi(0x0), | |
144 | fhTriggerPhi(0x0), | |
145 | fhJetEta(0x0), | |
146 | fhTriggerEta(0x0), | |
147 | fhVertexZ(0x0), | |
148 | fhVertexZAccept(0x0), | |
149 | fhContribVtx(0x0), | |
150 | fhContribVtxAccept(0x0), | |
151 | fhDphiTriggerJet(0x0), | |
152 | fhDphiTriggerJetAccept(0x0), | |
153 | fhCentrality(0x0), | |
154 | fhCentralityAccept(0x0), | |
155 | fkAcceptance(2.0*TMath::Pi()*1.8) | |
156 | { | |
157 | // Constructor | |
158 | fListJets = new TList(); | |
159 | ||
160 | DefineOutput(1, TList::Class()); | |
161 | } | |
162 | ||
163 | //-------------------------------------------------------------- | |
164 | AliAnalysisTaskJetCorePP::AliAnalysisTaskJetCorePP(const AliAnalysisTaskJetCorePP& a): | |
165 | AliAnalysisTaskSE(a.GetName()), | |
166 | fESD(a.fESD), | |
167 | fAODIn(a.fAODIn), | |
168 | fAODOut(a.fAODOut), | |
169 | fAODExtension(a.fAODExtension), | |
170 | fJetBranchName(a.fJetBranchName), | |
171 | fListJets(a.fListJets), | |
172 | fNonStdFile(a.fNonStdFile), | |
173 | fSystem(a.fSystem), | |
174 | fJetParamR(a.fJetParamR), | |
175 | fOfflineTrgMask(a.fOfflineTrgMask), | |
176 | fMinContribVtx(a.fMinContribVtx), | |
177 | fVtxZMin(a.fVtxZMin), | |
178 | fVtxZMax(a.fVtxZMax), | |
179 | fFilterMask(a.fFilterMask), | |
180 | fCentMin(a.fCentMin), | |
181 | fCentMax(a.fCentMax), | |
182 | fJetEtaMin(a.fJetEtaMin), | |
183 | fJetEtaMax(a.fJetEtaMax), | |
184 | fTriggerEtaCut(a.fTriggerEtaCut), | |
185 | fTrackEtaCut(a.fTrackEtaCut), | |
186 | fTrackLowPtCut(a.fTrackLowPtCut), | |
187 | fOutputList(a.fOutputList), | |
188 | fHistEvtSelection(a.fHistEvtSelection), | |
189 | fh2Ntriggers(a.fh2Ntriggers), | |
190 | fHJetSpec(a.fHJetSpec), | |
191 | fHJetDensity(a.fHJetDensity), | |
192 | fHJetDensityA4(a.fHJetDensityA4), | |
193 | fhJetPhi(a.fhJetPhi), | |
194 | fhTriggerPhi(a.fhTriggerPhi), | |
195 | fhJetEta(a.fhJetEta), | |
196 | fhTriggerEta(a.fhTriggerEta), | |
197 | fhVertexZ(a.fhVertexZ), | |
198 | fhVertexZAccept(a.fhVertexZAccept), | |
199 | fhContribVtx(a.fhContribVtx), | |
200 | fhContribVtxAccept(a.fhContribVtxAccept), | |
201 | fhDphiTriggerJet(a.fhDphiTriggerJet), | |
202 | fhDphiTriggerJetAccept(a.fhDphiTriggerJetAccept), | |
203 | fhCentrality(a.fhCentrality), | |
204 | fhCentralityAccept(a.fhCentralityAccept), | |
205 | fkAcceptance(a.fkAcceptance) | |
206 | { | |
207 | //Copy Constructor | |
208 | } | |
209 | //-------------------------------------------------------------- | |
210 | ||
211 | AliAnalysisTaskJetCorePP& AliAnalysisTaskJetCorePP::operator = (const AliAnalysisTaskJetCorePP& a){ | |
212 | // assignment operator | |
213 | this->~AliAnalysisTaskJetCorePP(); | |
214 | new(this) AliAnalysisTaskJetCorePP(a); | |
215 | return *this; | |
216 | } | |
217 | //-------------------------------------------------------------- | |
218 | ||
219 | AliAnalysisTaskJetCorePP::~AliAnalysisTaskJetCorePP() | |
220 | { | |
221 | //Destructor | |
222 | delete fListJets; | |
223 | delete fOutputList; // ???? | |
224 | } | |
225 | ||
226 | //-------------------------------------------------------------- | |
227 | ||
228 | void AliAnalysisTaskJetCorePP::Init() | |
229 | { | |
230 | // check for jet branches | |
231 | if(!strlen(fJetBranchName.Data())){ | |
232 | AliError("Jet branch name not set."); | |
233 | } | |
234 | ||
235 | } | |
236 | ||
237 | //-------------------------------------------------------------- | |
238 | ||
239 | void AliAnalysisTaskJetCorePP::UserCreateOutputObjects() | |
240 | { | |
241 | ||
242 | ||
243 | // Create histograms | |
244 | // Called once | |
245 | OpenFile(1); | |
246 | if(!fOutputList) fOutputList = new TList(); | |
247 | fOutputList->SetOwner(kTRUE); | |
248 | ||
249 | Bool_t oldStatus = TH1::AddDirectoryStatus(); | |
250 | TH1::AddDirectory(kFALSE); | |
251 | ||
252 | fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5); | |
253 | fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED"); | |
254 | fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN"); | |
255 | fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)"); | |
256 | fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)"); | |
257 | fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)"); | |
258 | fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)"); | |
259 | ||
260 | fOutputList->Add(fHistEvtSelection); | |
261 | ||
262 | Int_t nBinsCentrality = (fSystem==0) ? 1 : 10; // pp=1 else 10 | |
263 | ||
264 | fh2Ntriggers = new TH2F("fh2Ntriggers","# of triggers", | |
265 | nBinsCentrality,0.0,100.0,50,0.0,50.0); | |
266 | ||
267 | //Centrality, A, pT - rho*A, pTtrigg, pT, rho*A | |
268 | const Int_t dimSpec = 6; | |
269 | const Int_t nBinsSpec[dimSpec] = {nBinsCentrality, 100, 140, 50, 100, 50}; | |
270 | const Double_t lowBinSpec[dimSpec] = {0.0, 0.0,-80.0, 0.0, 0.0, 0.0}; | |
271 | const Double_t hiBinSpec[dimSpec] = {100.0, 1.0,200.0,50.0, 200, 20.0}; | |
272 | fHJetSpec = new THnSparseF("fHJetSpec","Recoil jet spectrum [cent,A,pT-rho*A,pTtrig]", | |
273 | dimSpec,nBinsSpec,lowBinSpec,hiBinSpec); | |
274 | fOutputList->Add(fHJetSpec); | |
275 | ||
276 | //------------------- HISTOS FOR DIAGNOSTIC ---------------------- | |
277 | //Jet number density histos [Trk Mult, jet density, pT trigger] | |
278 | const Int_t dimJetDens = 3; | |
279 | const Int_t nBinsJetDens[dimJetDens] = {100, 100, 10}; | |
280 | const Double_t lowBinJetDens[dimJetDens] = {0.0, 0.0, 0.0}; | |
281 | const Double_t hiBinJetDens[dimJetDens] = {500.0, 5.0, 50.0 }; | |
282 | ||
283 | fHJetDensity = new THnSparseF("fHJetDensity","Jet dens vs trk mult A>0.07", | |
284 | dimJetDens,nBinsJetDens,lowBinJetDens,hiBinJetDens); | |
285 | ||
286 | fHJetDensityA4 =new THnSparseF("fHJetDensityA4","Jet dens vs trk mult A>0.4", | |
287 | dimJetDens,nBinsJetDens,lowBinJetDens,hiBinJetDens); | |
288 | ||
289 | fOutputList->Add(fh2Ntriggers); | |
290 | fOutputList->Add(fHJetDensity); | |
291 | fOutputList->Add(fHJetDensityA4); | |
292 | ||
293 | ||
294 | //inclusive azimuthal and pseudorapidity histograms | |
295 | fhJetPhi = new TH1D("fhJetPhi","Azim dist jets",50,-TMath::Pi(),TMath::Pi()); | |
296 | fhTriggerPhi= new TH1D("fhTriggerPhi","azim dist trig had",50,-TMath::Pi(),TMath::Pi()); | |
297 | fhJetEta = new TH1D("fhJetEta","Eta dist jets",40,-0.9,0.9); | |
298 | fhTriggerEta = new TH1D("fhTriggerEta","Eta dist trig had",40,-0.9,0.9); | |
299 | fhVertexZ = new TH1D("fhVertexZ","z vertex",40,-20,20); | |
300 | fhVertexZAccept = new TH1D("fhVertexZAccept","z vertex after cut",40,-20,20); | |
301 | fhContribVtx = new TH1D("fhContribVtx","contrib to vtx",200,0,200); | |
302 | fhContribVtxAccept = new TH1D("fhContribVtxAccept","contrib to vtx after cut",200,0,200); | |
303 | fhDphiTriggerJet = new TH1D("fhDphiTriggerJet","Deltaphi trig-jet",50, -TMath::Pi(),TMath::Pi()); | |
304 | fhDphiTriggerJetAccept = new TH1D("fhDphiTriggerJetAccept","Deltaphi trig-jet after cut",50, -TMath::Pi(),TMath::Pi()); | |
305 | fhCentrality = new TH1D("fhCentrality","Centrality",20,0,100); | |
306 | fhCentralityAccept = new TH1D("fhCentralityAccept","Centrality after cut",20,0,100); | |
307 | ||
308 | fOutputList->Add(fhJetPhi); | |
309 | fOutputList->Add(fhTriggerPhi); | |
310 | fOutputList->Add(fhJetEta); | |
311 | fOutputList->Add(fhTriggerEta); | |
312 | fOutputList->Add(fhVertexZ); | |
313 | fOutputList->Add(fhVertexZAccept); | |
314 | fOutputList->Add(fhContribVtx); | |
315 | fOutputList->Add(fhContribVtxAccept); | |
316 | fOutputList->Add(fhDphiTriggerJet); | |
317 | fOutputList->Add(fhDphiTriggerJetAccept); | |
318 | fOutputList->Add(fhCentrality); | |
319 | fOutputList->Add(fhCentralityAccept); | |
320 | ||
321 | ||
322 | // =========== Switch on Sumw2 for all histos =========== | |
323 | for(Int_t i=0; i<fOutputList->GetEntries(); i++){ | |
324 | TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i)); | |
325 | if(h1){ | |
326 | h1->Sumw2(); | |
327 | continue; | |
328 | } | |
329 | THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i)); | |
330 | if(hn){ | |
331 | hn->Sumw2(); | |
332 | } | |
333 | } | |
334 | TH1::AddDirectory(oldStatus); | |
335 | ||
336 | PostData(1, fOutputList); | |
337 | } | |
338 | ||
339 | //-------------------------------------------------------------------- | |
340 | ||
341 | void AliAnalysisTaskJetCorePP::UserExec(Option_t *) | |
342 | { | |
343 | ||
344 | //Event loop | |
345 | ||
346 | if(TMath::Abs((Float_t) fJetParamR)<0.00001){ | |
347 | AliError("Cone radius is set to zero."); | |
348 | return; | |
349 | } | |
350 | if(!strlen(fJetBranchName.Data())){ | |
351 | AliError("Jet branch name not set."); | |
352 | return; | |
353 | } | |
354 | ||
355 | fESD = dynamic_cast<AliESDEvent*>(InputEvent()); | |
356 | if(!fESD){ | |
357 | AliError("ESD not available"); | |
358 | fAODIn = dynamic_cast<AliAODEvent*>(InputEvent()); | |
359 | } | |
360 | ||
361 | fAODOut = dynamic_cast<AliAODEvent*>(AODEvent()); | |
362 | AliAODEvent* aod = NULL; | |
363 | // take all other information from the aod we take the tracks from | |
364 | if(!aod){ | |
365 | if(!fESD) aod = fAODIn; | |
366 | else aod = fAODOut; | |
367 | } | |
368 | ||
369 | if(fNonStdFile.Length()!=0){ | |
370 | // case that we have an AOD extension we can fetch the jets from the extended output | |
371 | AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()); | |
372 | fAODExtension = aodH ? aodH->GetExtension(fNonStdFile.Data()) : 0; | |
373 | if(!fAODExtension){ | |
374 | if(fDebug>1) Printf("AODExtension found for %s",fNonStdFile.Data()); | |
375 | } | |
376 | } | |
377 | ||
378 | // ----------------- event selection -------------------------- | |
379 | fHistEvtSelection->Fill(1); // number of events before event selection | |
380 | ||
381 | // physics selection | |
382 | AliInputEventHandler* inputHandler = (AliInputEventHandler*) | |
383 | ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler()); | |
384 | ||
385 | if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){ | |
386 | if(fDebug) Printf(" Trigger Selection: event REJECTED ... "); | |
387 | fHistEvtSelection->Fill(2); | |
388 | PostData(1, fOutputList); | |
389 | return; | |
390 | } | |
391 | ||
392 | //check AOD pointer | |
393 | if(!aod){ | |
394 | if(fDebug) Printf("%s:%d No AOD",(char*)__FILE__,__LINE__); | |
395 | fHistEvtSelection->Fill(3); | |
396 | PostData(1, fOutputList); | |
397 | return; | |
398 | } | |
399 | ||
400 | // vertex selection | |
401 | AliAODVertex* primVtx = aod->GetPrimaryVertex(); | |
402 | ||
403 | if(!primVtx){ | |
404 | if(fDebug) Printf("%s:%d No primVtx",(char*)__FILE__,__LINE__); | |
405 | fHistEvtSelection->Fill(3); | |
406 | PostData(1, fOutputList); | |
407 | return; | |
408 | } | |
409 | ||
410 | Int_t nTracksPrim = primVtx->GetNContributors(); | |
411 | Float_t vtxz = primVtx->GetZ(); | |
412 | //Input events | |
413 | fhContribVtx->Fill(nTracksPrim); | |
414 | if( nTracksPrim > 0 ) fhVertexZ->Fill(vtxz); | |
415 | ||
416 | if((nTracksPrim < fMinContribVtx) || | |
417 | (vtxz < fVtxZMin) || | |
418 | (vtxz > fVtxZMax)){ | |
419 | if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...", | |
420 | (char*)__FILE__,__LINE__,vtxz); | |
421 | fHistEvtSelection->Fill(3); | |
422 | PostData(1, fOutputList); | |
423 | return; | |
424 | }else{ | |
425 | //Accepted events | |
426 | fhContribVtxAccept->Fill(nTracksPrim); | |
427 | fhVertexZAccept->Fill(vtxz); | |
428 | } | |
429 | //FK// No event class selection imposed | |
430 | // event class selection (from jet helper task) | |
431 | //Int_t eventClass = AliAnalysisHelperJetTasks::EventClass(); | |
432 | //if(fDebug) Printf("Event class %d", eventClass); | |
433 | //if(eventClass < fEvtClassMin || eventClass > fEvtClassMax){ | |
434 | // fHistEvtSelection->Fill(4); | |
435 | // PostData(1, fOutputList); | |
436 | // return; | |
437 | //} | |
438 | ||
439 | // centrality selection | |
440 | AliCentrality *cent = 0x0; | |
441 | Double_t centValue = 0.0; | |
442 | if(fSystem){ //fSystem=0 for pp, fSystem=1 for pPb | |
443 | if(fESD){ | |
444 | cent = fESD->GetCentrality(); | |
445 | if(cent) centValue = cent->GetCentralityPercentile("V0M"); | |
446 | }else{ | |
447 | centValue = aod->GetHeader()->GetCentrality(); | |
448 | } | |
449 | if(fDebug) printf("centrality: %f\n", centValue); | |
450 | //Input events | |
451 | fhCentrality->Fill(centValue); | |
452 | ||
453 | if(centValue < fCentMin || centValue > fCentMax){ | |
454 | fHistEvtSelection->Fill(4); | |
455 | PostData(1, fOutputList); | |
456 | return; | |
457 | }else{ | |
458 | //Accepted events | |
459 | fhCentralityAccept->Fill( centValue ); | |
460 | } | |
461 | } | |
462 | ||
463 | if(fDebug) std::cout<<" ACCEPTED EVENT "<<endl; | |
464 | ||
465 | fHistEvtSelection->Fill(0); // accepted events | |
466 | ||
467 | // ------------------- end event selection -------------------- | |
468 | ||
469 | // fetch jets | |
470 | TClonesArray *aodJets = 0x0; | |
471 | ||
472 | if(fAODOut && !aodJets){ | |
473 | aodJets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName.Data())); | |
474 | } | |
475 | if(fAODExtension && !aodJets){ | |
476 | aodJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName.Data())); | |
477 | } | |
478 | if(fAODIn && !aodJets){ | |
479 | aodJets = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName.Data())); | |
480 | } | |
481 | ||
482 | // ------------- Hadron trigger -------------- | |
483 | TList particleList; //list of tracks | |
484 | Int_t indexTrigg = GetListOfTracks(&particleList); //index of trigger hadron in Particle list | |
485 | ||
486 | if(indexTrigg<0) return; // no trigger track found | |
487 | ||
488 | AliVParticle *triggerHadron = (AliVParticle*) particleList.At(indexTrigg); | |
489 | if(!triggerHadron){ | |
490 | PostData(1, fOutputList); | |
491 | return; | |
492 | } | |
493 | ||
494 | ||
495 | fh2Ntriggers->Fill(centValue,triggerHadron->Pt()); //trigger pT | |
496 | ||
497 | // if(triggerHadron->Pt() > 10.0){} | |
498 | if(triggerHadron->Pt() > 3.0){ | |
499 | //Trigger Diagnostics--------------------------------- | |
500 | fhTriggerPhi->Fill(RelativePhi(triggerHadron->Phi(),0.0)); //phi -pi,pi | |
501 | fhTriggerEta->Fill(triggerHadron->Eta()); | |
502 | } | |
503 | ||
504 | //--------- Fill list of jets -------------- | |
505 | fListJets->Clear(); | |
506 | if(aodJets){ | |
507 | if(fDebug) Printf("########## %s: %d jets",fJetBranchName.Data(), | |
508 | aodJets->GetEntriesFast()); | |
509 | for(Int_t iJet = 0; iJet < aodJets->GetEntriesFast(); iJet++) { | |
510 | AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets)[iJet]); | |
511 | if (jet) fListJets->Add(jet); | |
512 | } | |
513 | } | |
514 | ||
515 | Double_t etaJet = 0.0; | |
516 | Double_t pTJet = 0.0; | |
517 | Double_t areaJet = 0.0; | |
518 | Double_t phiJet = 0.0; | |
519 | Int_t injet4 = 0; | |
520 | Int_t injet = 0; | |
521 | //---------- jet loop --------- | |
522 | for(Int_t ij=0; ij<fListJets->GetEntries(); ij++){ | |
523 | AliAODJet* jet = (AliAODJet*)(fListJets->At(ij)); | |
524 | if(!jet) continue; | |
525 | etaJet = jet->Eta(); | |
526 | phiJet = jet->Phi(); | |
527 | pTJet = jet->Pt(); | |
528 | if(pTJet==0) continue; | |
529 | ||
530 | if((etaJet<fJetEtaMin) || (etaJet>fJetEtaMax)) continue; | |
531 | areaJet = jet->EffectiveAreaCharged(); | |
532 | ||
533 | //Jet Diagnostics--------------------------------- | |
534 | fhJetPhi->Fill(RelativePhi(phiJet,0.0)); //phi -pi,pi | |
535 | fhJetEta->Fill(etaJet); | |
536 | if(areaJet >= 0.07) injet++; | |
537 | if(areaJet >= 0.4) injet4++; | |
538 | //-------------------------------------------------- | |
539 | ||
540 | Double_t dphi = RelativePhi(triggerHadron->Phi(), phiJet); | |
541 | ||
542 | fhDphiTriggerJet->Fill(dphi); //Input | |
543 | if(TMath::Abs((Double_t) dphi) < TMath::Pi()-0.6) continue; | |
544 | fhDphiTriggerJetAccept->Fill(dphi); //Accepted | |
545 | ||
546 | //Background w.r.t. jet axis | |
547 | Double_t pTBckWrtJet = | |
548 | GetBackgroundInPerpCone(fJetParamR, phiJet, etaJet, &particleList); | |
549 | ||
550 | Double_t ratioOfAreas = areaJet/(TMath::Pi()*fJetParamR*fJetParamR); | |
551 | Double_t rhoA = ratioOfAreas*pTBckWrtJet; //bg activity in a cone of similar size | |
552 | Double_t ptcorr = pTJet - rhoA; //Pt Jet UE subtr | |
553 | ||
554 | ||
555 | //Centrality, A, pT - rho*A, pTtrigg, pT, rho*A | |
556 | Double_t fillspec[] = { centValue, | |
557 | areaJet, | |
558 | ptcorr, | |
559 | triggerHadron->Pt(), | |
560 | pTJet, | |
561 | rhoA | |
562 | }; | |
563 | fHJetSpec->Fill(fillspec); | |
564 | ||
565 | }//jet loop | |
566 | ||
567 | //Fill Jet Density In the Event A>0.07 | |
568 | if(injet>0){ | |
569 | Double_t filldens[]={ (Double_t) particleList.GetEntries(), | |
570 | injet/fkAcceptance, | |
571 | triggerHadron->Pt() | |
572 | }; | |
573 | fHJetDensity->Fill(filldens); | |
574 | } | |
575 | ||
576 | //Fill Jet Density In the Event A>0.4 | |
577 | if(injet4>0){ | |
578 | Double_t filldens4[]={ (Double_t) particleList.GetEntries(), | |
579 | injet4/fkAcceptance, | |
580 | triggerHadron->Pt() | |
581 | }; | |
582 | fHJetDensityA4->Fill(filldens4); | |
583 | } | |
584 | ||
585 | PostData(1, fOutputList); | |
586 | } | |
587 | ||
588 | //---------------------------------------------------------------------------- | |
589 | void AliAnalysisTaskJetCorePP::Terminate(const Option_t *) | |
590 | { | |
591 | // Draw result to the screen | |
592 | // Called once at the end of the query | |
593 | ||
594 | if(fDebug) printf("AliAnalysisTaskJetCorePP DONE\n"); | |
595 | if(!GetOutputData(1)) return; | |
596 | } | |
597 | ||
598 | ||
599 | //---------------------------------------------------------------------------- | |
600 | Int_t AliAnalysisTaskJetCorePP::GetListOfTracks(TList *list){ | |
601 | //Fill the list of accepted tracks (passed track cut) | |
602 | //return consecutive index of the hardest ch hadron in the list | |
603 | Int_t iCount = 0; | |
604 | AliAODEvent *aodevt = NULL; | |
605 | ||
606 | if(!fESD) aodevt = fAODIn; | |
607 | else aodevt = fAODOut; | |
608 | ||
609 | if(!aodevt) return -1; | |
610 | ||
611 | Int_t index = -1; //index of the highest particle in the list | |
612 | Double_t ptmax = -10; | |
613 | ||
614 | for(int it = 0; it < aodevt->GetNumberOfTracks(); it++){ | |
615 | AliAODTrack *tr = aodevt->GetTrack(it); | |
616 | ||
617 | //if((fFilterMask > 0) && !(tr->TestFilterBit(fFilterMask))) continue; | |
618 | if((fFilterMask > 0) && !(tr->IsHybridGlobalConstrainedGlobal())) continue; | |
619 | if(TMath::Abs((Float_t) tr->Eta()) > fTrackEtaCut) continue; | |
620 | if(tr->Pt() < fTrackLowPtCut) continue; | |
621 | list->Add(tr); | |
622 | if(tr->Pt()>ptmax){ | |
623 | ptmax = tr->Pt(); | |
624 | index = iCount; | |
625 | } | |
626 | iCount++; | |
627 | } | |
628 | ||
629 | if(index>-1){ //check pseudorapidity cut on trigger | |
630 | AliAODTrack *trigger = (AliAODTrack*) list->At(index); | |
631 | if(trigger && TMath::Abs((Float_t) trigger->Eta())< fTriggerEtaCut){ return index;} | |
632 | return -1; | |
633 | }else{ | |
634 | return -1; | |
635 | } | |
636 | } | |
637 | ||
638 | //---------------------------------------------------------------------------- | |
639 | ||
640 | Double_t AliAnalysisTaskJetCorePP::GetBackgroundInPerpCone(Float_t jetR, Double_t jetPhi, Double_t jetEta, TList* trkList){ | |
641 | //calculate sum of track pT in the cone perpendicular in phi to the jet | |
642 | //jetR = cone radius | |
643 | // jetPhi, jetEta = direction of the jet | |
644 | Int_t numberOfTrks = trkList->GetEntries(); | |
645 | Double_t pTsum = 0.0; | |
646 | Double_t perpConePhi = jetPhi + TMath::Pi()/2;//perp cone w.r.t. jet in phi | |
647 | for(Int_t it=0; it<numberOfTrks; it++){ | |
648 | AliVParticle *trk = (AliVParticle*) trkList->At(it); | |
649 | Double_t dphi = RelativePhi(perpConePhi,trk->Phi()); | |
650 | Double_t deta = trk->Eta()-jetEta; | |
651 | ||
652 | if( (dphi*dphi + deta*deta)< (jetR*jetR)){ | |
653 | pTsum += trk->Pt(); | |
654 | } | |
655 | } | |
656 | ||
657 | return pTsum; | |
658 | } | |
659 | ||
660 | //---------------------------------------------------------------------------- | |
661 | ||
662 | Double_t AliAnalysisTaskJetCorePP::RelativePhi(Double_t mphi,Double_t vphi){ | |
663 | //Get relative azimuthal angle of two particles -pi to pi | |
664 | if (vphi < -TMath::Pi()) vphi += TMath::TwoPi(); | |
665 | else if (vphi > TMath::Pi()) vphi -= TMath::TwoPi(); | |
666 | ||
667 | if (mphi < -TMath::Pi()) mphi += TMath::TwoPi(); | |
668 | else if (mphi > TMath::Pi()) mphi -= TMath::TwoPi(); | |
669 | ||
670 | Double_t dphi = mphi - vphi; | |
671 | if (dphi < -TMath::Pi()) dphi += TMath::TwoPi(); | |
672 | else if (dphi > TMath::Pi()) dphi -= TMath::TwoPi(); | |
673 | ||
674 | return dphi;//dphi in [-Pi, Pi] | |
675 | } | |
676 | ||
677 |