]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/Correlations/Base/AliAnalyseLeadingTrackUE.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
fixed reference run for pA in centrality OADB
[u/mrichter/AliRoot.git] / PWGCF / Correlations / Base / AliAnalyseLeadingTrackUE.cxx
CommitLineData
a75aacd6 1/*************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: A.Abrahantes, E.Lopez, S.Vallero *
5 * Version 1.0 *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
a75aacd6 15//#include <TBranch.h>
16//#include <TCanvas.h>
17//#include <TChain.h>
18//#include <TFile.h>
19//#include <TH1F.h>
20//#include <TH1I.h>
21//#include <TH2F.h>
22#include <TList.h>
23//#include <TLorentzVector.h>
24#include <TMath.h>
25#include <TObjArray.h>
26#include <TObject.h>
27//#include <TProfile.h>
28//#include <TRandom.h>
29//#include <TSystem.h>
30//#include <TTree.h>
31#include <TVector3.h>
32
33#include "AliAnalyseLeadingTrackUE.h"
34//#include "AliAnalysisTask.h"
35
36//#include "AliAnalysisHelperJetTasks.h"
37//#include "AliAnalysisManager.h"
38#include "AliAODEvent.h"
39//#include "AliAODHandler.h"
a75aacd6 40//#include "AliAODJet.h"
41#include "AliAODMCParticle.h"
42#include "AliAODTrack.h"
43#include "AliESDEvent.h"
44#include "AliESDtrack.h"
45#include "AliESDtrackCuts.h"
46//#include "AliGenPythiaEventHeader.h"
47#include "AliInputEventHandler.h"
48//#include "AliKFVertex.h"
144bd037 49//#include "AliLog.h"
a75aacd6 50#include "AliMCEvent.h"
a75aacd6 51#include "AliVParticle.h"
7dd4dec4 52#include "AliAODMCHeader.h"
a75aacd6 53
ea1919ac 54#include "AliAnalysisManager.h"
55#include "AliMCEventHandler.h"
56#include "AliStack.h"
57
58
a75aacd6 59////////////////////////////////////////////////
60//---------------------------------------------
61// Class for transverse regions analysis
62//---------------------------------------------
63////////////////////////////////////////////////
64
65
66using namespace std;
67
68ClassImp(AliAnalyseLeadingTrackUE)
69
70//-------------------------------------------------------------------
71AliAnalyseLeadingTrackUE::AliAnalyseLeadingTrackUE() :
72 TObject(),
73 fDebug(0),
74 fFilterBit(16),
75 fOnlyHadrons(kFALSE),
3712ba26 76 fTrackEtaCut(0.8),
e0331fd9 77 fTrackPtMin(0),
5cabd2cd 78 fEventSelection(AliVEvent::kMB|AliVEvent::kUserDefined),
3712ba26 79 fEsdTrackCuts(0x0),
80 fEsdTrackCutsSPD(0x0),
81 fEsdTrackCutsSDD(0x0)
a75aacd6 82{
83 // constructor
84}
85
86
87//-------------------------------------------------------------------
88AliAnalyseLeadingTrackUE & AliAnalyseLeadingTrackUE::operator = (const AliAnalyseLeadingTrackUE & /*source*/)
89{
90 // assignment operator
91 return *this;
92}
93
94
95//-------------------------------------------------------------------
96AliAnalyseLeadingTrackUE::~AliAnalyseLeadingTrackUE()
97{
98
99 //clear memory
100
101}
102
103
104//____________________________________________________________________
3712ba26 105Bool_t AliAnalyseLeadingTrackUE::ApplyCuts(TObject* track)
a75aacd6 106{
a75aacd6 107
108 // select track according to set of cuts
85bfac17 109 if (!fEsdTrackCuts->IsSelected(track) )return kFALSE;
110 if (fEsdTrackCutsSPD && fEsdTrackCutsSDD && !fEsdTrackCutsSPD->IsSelected(track) && fEsdTrackCutsSDD->IsSelected(track)) return kFALSE;
3712ba26 111
a75aacd6 112 return kTRUE;
113}
114
115
3712ba26 116//____________________________________________________________________
1be4eb33 117void AliAnalyseLeadingTrackUE::DefineESDCuts(Int_t filterbit) {
3712ba26 118
85bfac17 119 // Reproduces the cuts of the corresponding bit in the ESD->AOD filtering
120 // (see $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C)
1be4eb33 121
122 if (filterbit == -1)
123 filterbit = fFilterBit;
a75aacd6 124
1be4eb33 125 if (filterbit == 128)
85bfac17 126 {
127 fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
128 fEsdTrackCuts->SetMinNClustersTPC(70);
129 }
1be4eb33 130 else if (filterbit == 256)
131 {
132 // syst study
133 fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
134 fEsdTrackCuts->SetMinNClustersTPC(80);
135 fEsdTrackCuts->SetMaxChi2PerClusterTPC(3);
136 fEsdTrackCuts->SetMaxDCAToVertexZ(2.7);
137 fEsdTrackCuts->SetMaxDCAToVertexXY(1.9);
138 }
139 else if (filterbit == 512)
140 {
141 // syst study
142 fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
143 fEsdTrackCuts->SetMinNClustersTPC(60);
144 fEsdTrackCuts->SetMaxChi2PerClusterTPC(5);
145 fEsdTrackCuts->SetMaxDCAToVertexZ(3.7);
146 fEsdTrackCuts->SetMaxDCAToVertexXY(2.9);
147 }
4bbcfff9 148 else if (filterbit == 1024)
149 {
150 fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
151 fEsdTrackCuts->SetMinNClustersTPC(-1);
152 fEsdTrackCuts->SetMinNCrossedRowsTPC(70);
153 fEsdTrackCuts->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
154 }
85bfac17 155 else
156 {
157 fEsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
158 fEsdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kNone);
159
160 // Add SPD requirement
161 fEsdTrackCutsSPD = new AliESDtrackCuts("SPD", "Require 1 cluster in SPD");
162 fEsdTrackCutsSPD->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
163
164 // Add SDD requirement
165 fEsdTrackCutsSDD = new AliESDtrackCuts("SDD", "Require 1 cluster in first layer SDD");
166 fEsdTrackCutsSDD->SetClusterRequirementITS(AliESDtrackCuts::kSDD,AliESDtrackCuts::kFirst);
167 }
3712ba26 168}
a75aacd6 169
a75aacd6 170//____________________________________________________________________
171TObjArray* AliAnalyseLeadingTrackUE::FindLeadingObjects(TObject *obj)
172{
173
174 // Returns an array of charged particles (or jets) ordered according to their pT.
175
176 Int_t nTracks = NParticles(obj);
177
178
179 if( !nTracks ) return 0;
180
181 // Define array of AliVParticle objects
182 TObjArray* tracks = new TObjArray(nTracks);
183
184 // Loop over tracks or jets
185 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
186 AliVParticle* part = ParticleWithCuts( obj, ipart );
187 if (!part) continue;
188 // Accept leading-tracks in a limited pseudo-rapidity range
189 if( TMath::Abs(part->Eta()) > fTrackEtaCut ) continue;
190 tracks->AddLast( part );
191 }
192 // Order tracks by pT
193 QSortTracks( *tracks, 0, tracks->GetEntriesFast() );
194
195 nTracks = tracks->GetEntriesFast();
196 if( !nTracks ) return 0;
197
198 return tracks;
199 }
200
201
1ccd8a0a 202void AliAnalyseLeadingTrackUE::RemoveInjectedSignals(TObjArray* tracks, TObject* mcObj, Int_t maxLabel)
203{
204 // remove injected signals (primaries above <maxLabel>)
205 // <tracks> can be the following cases:
206 // a. tracks: in this case the label is taken and then case b.
207 // b. particles: the first stable mother is searched and checked if it is <= <maxLabel>
208 // <mcObj> can be AOD (TClonesArray) or ESD (AliMCEvent)
209
210 TClonesArray* arrayMC = 0;
211 AliMCEvent* mcEvent = 0;
212 if (mcObj->InheritsFrom("AliMCEvent"))
213 mcEvent = static_cast<AliMCEvent*>(mcObj);
214 else if (mcObj->InheritsFrom("TClonesArray"))
215 arrayMC = static_cast<TClonesArray*>(mcObj);
216 else
217 {
218 arrayMC->Dump();
219 AliFatal("Invalid object passed");
220 }
221
222 Int_t before = tracks->GetEntriesFast();
223
224 for (Int_t i=0; i<before; ++i)
225 {
226 AliVParticle* part = (AliVParticle*) tracks->At(i);
227
228 if (part->InheritsFrom("AliESDtrack") || part->InheritsFrom("AliAODTrack"))
229 part = ((mcEvent) ? mcEvent->GetTrack(TMath::Abs(part->GetLabel())) : (AliVParticle*)arrayMC->At(TMath::Abs(part->GetLabel())));
230
231 AliVParticle* mother = part;
232 if (mcEvent)
233 {
234 while (!mcEvent->IsPhysicalPrimary(mother->GetLabel()))
235 {
236 if (((AliMCParticle*)mother)->GetMother() < 0)
237 {
238 mother = 0;
239 break;
240 }
241
242 mother = (AliMCParticle*) mcEvent->GetTrack(((AliMCParticle*)mother)->GetMother());
243 if (!mother)
244 break;
245 }
246 }
247 else
248 {
249 // find the primary mother
250 while (!((AliAODMCParticle*)mother)->IsPhysicalPrimary())
251 {
252 if (((AliAODMCParticle*)mother)->GetMother() < 0)
253 {
254 mother = 0;
255 break;
256 }
257
258 mother = (AliVParticle*) arrayMC->At(((AliAODMCParticle*)mother)->GetMother());
259 if (!mother)
260 break;
261 }
262 }
263
264 if (!mother)
265 {
266 Printf("WARNING: No mother found for particle %d:", part->GetLabel());
267 continue;
268 }
269
270 if (mother->GetLabel() > maxLabel)
271 {
272// Printf("Removing %d with label %d", i, part->GetLabel()); part->Dump();
273 TObject* object = tracks->RemoveAt(i);
274 if (tracks->IsOwner())
275 delete object;
276 }
277 }
278
279 tracks->Compress();
280
281 AliInfo(Form("Reduced from %d to %d", before, tracks->GetEntriesFast()));
282}
283
b1831bcb 284//-------------------------------------------------------------------
e0331fd9 285TObjArray* AliAnalyseLeadingTrackUE::GetAcceptedParticles(TObject* obj, TObject* arrayMC, Bool_t onlyprimaries, Int_t particleSpecies, Bool_t useEtaPtCuts)
b1831bcb 286{
287 // Returns an array of particles that pass the cuts, if arrayMC is given each reconstructed particle is replaced by its corresponding MC particles, depending on the parameter onlyprimaries only for primaries
288 // particleSpecies: -1 all particles are returned
289 // 0 (pions) 1 (kaons) 2 (protons) 3 (others) particles
290
291 Int_t nTracks = NParticles(obj);
292 TObjArray* tracks = new TObjArray;
85bfac17 293
294 // for TPC only tracks
eed401dc 295 Bool_t hasOwnership = kFALSE;
1bba939a 296 if ((fFilterBit == 128 || fFilterBit == 256 || fFilterBit == 512 || fFilterBit == 1024) && obj->InheritsFrom("AliESDEvent"))
eed401dc 297 hasOwnership = kTRUE;
298
299 if (!arrayMC)
300 tracks->SetOwner(hasOwnership);
b1831bcb 301
302 // Loop over tracks or jets
303 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
304 AliVParticle* part = ParticleWithCuts( obj, ipart, onlyprimaries, particleSpecies );
305 if (!part) continue;
306
e0331fd9 307 if (useEtaPtCuts)
308 if (TMath::Abs(part->Eta()) > fTrackEtaCut || part->Pt() < fTrackPtMin)
92834c86 309 {
310 if (hasOwnership)
311 delete part;
e0331fd9 312 continue;
92834c86 313 }
e0331fd9 314
eed401dc 315 if (arrayMC) {
316 Int_t label = part->GetLabel();
317 if (hasOwnership)
318 delete part;
b1831bcb 319 // re-define part as the matched MC particle
eed401dc 320 part = ParticleWithCuts(arrayMC, TMath::Abs(label),onlyprimaries, particleSpecies);
b1831bcb 321 if (!part)continue;
322 }
323
324 tracks->Add(part);
325 }
326
327 return tracks;
328}
329
a75aacd6 330//-------------------------------------------------------------------
331TObjArray* AliAnalyseLeadingTrackUE::GetMinMaxRegion(TList *transv1, TList *transv2)
332{
333
334 // Returns two lists of particles, one for MIN and one for MAX region
335 Double_t sumpT1 = 0.;
336 Double_t sumpT2 = 0.;
337
338 Int_t particles1 = transv1->GetEntries();
339 Int_t particles2 = transv2->GetEntries();
340
341 // Loop on transverse region 1
342 for (Int_t i=0; i<particles1; i++){
343 AliVParticle *part = (AliVParticle*)transv1->At(i);
344 sumpT1 += part->Pt();
345 }
346
347 // Loop on transverse region 2
348 for (Int_t i=0; i<particles2; i++){
349 AliVParticle *part = (AliVParticle*)transv2->At(i);
350 sumpT2 += part->Pt();
351 }
352
353 TObjArray *regionParticles = new TObjArray;
354 if ( sumpT2 >= sumpT1 ){
355 regionParticles->AddLast(transv1); // MIN
356 regionParticles->AddLast(transv2); // MAX
357 }else {
358 regionParticles->AddLast(transv2); // MIN
359 regionParticles->AddLast(transv1); // MAX
360 }
361
362 return regionParticles;
363}
364
365//-------------------------------------------------------------------
366Int_t AliAnalyseLeadingTrackUE::NParticles(TObject* obj)
367{
368
369 //Returns the number of particles in AliAODMCParticle array or AliAODTracks or AliESDTracks
370
371 Int_t nTracks;
372
373 if (obj->InheritsFrom("TClonesArray")){ // MC particles
ea1919ac 374 TClonesArray *arrayMC = static_cast<TClonesArray*>(obj);
a75aacd6 375 nTracks = arrayMC->GetEntriesFast();
376 }else if (obj->InheritsFrom("TObjArray")){ // list of AliVParticle
ea1919ac 377 TObjArray *array = static_cast<TObjArray*>(obj);
a75aacd6 378 nTracks = array->GetEntriesFast();
379 }else if (obj->InheritsFrom("AliAODEvent")){ // RECO AOD tracks
ea1919ac 380 AliAODEvent *aodEvent = static_cast<AliAODEvent*>(obj);
a75aacd6 381 nTracks = aodEvent->GetNTracks();
382 }else if (obj->InheritsFrom("AliESDEvent")){ // RECO ESD tracks
ea1919ac 383 AliESDEvent *esdEvent = static_cast<AliESDEvent*>(obj);
a75aacd6 384 nTracks = esdEvent->GetNumberOfTracks();
3712ba26 385 }else if (obj->InheritsFrom("AliMCEvent")){ // RECO ESD tracks
ea1919ac 386 AliMCEvent *mcEvent = static_cast<AliMCEvent*>(obj);
3712ba26 387 nTracks = mcEvent->GetNumberOfTracks();
a75aacd6 388 }else {
389 if (fDebug > 1) AliFatal(" Analysis type not defined !!! ");
390 return 0;
391 }
392
393 return nTracks;
394}
395
396
397//-------------------------------------------------------------------
b1831bcb 398AliVParticle* AliAnalyseLeadingTrackUE::ParticleWithCuts(TObject* obj, Int_t ipart, Bool_t onlyprimaries, Int_t particleSpecies)
a75aacd6 399{
400 // Returns track or MC particle at position "ipart" if passes selection criteria
b1831bcb 401 // particleSpecies: -1 all particles are returned
402 // 0 (pions) 1 (kaons) 2 (protons) 3 (others) particles
a75aacd6 403 AliVParticle *part=0;
404
405 if (obj->InheritsFrom("TClonesArray")){ // AOD-MC PARTICLE
ea1919ac 406 TClonesArray *arrayMC = static_cast<TClonesArray*>(obj);
a75aacd6 407 part = (AliVParticle*)arrayMC->At( ipart );
408 if (!part)return 0;
409 // eventually only primaries
410 if (onlyprimaries && !( ((AliAODMCParticle*)part)->IsPhysicalPrimary()) )return 0;
411 // eventually only hadrons
412 if (fOnlyHadrons){
413 Int_t pdgCode = ((AliAODMCParticle*)part)->GetPdgCode();
414 Bool_t isHadron = TMath::Abs(pdgCode)==211 || // Pion
415 TMath::Abs(pdgCode)==2212 || // Proton
416 TMath::Abs(pdgCode)==321; // Kaon
417 if (!isHadron) return 0;
418 }
b1831bcb 419 if (particleSpecies != -1) {
420 // find the primary mother
421 AliVParticle* mother = part;
422 while (!((AliAODMCParticle*)mother)->IsPhysicalPrimary())
423 {
424 if (((AliAODMCParticle*)mother)->GetMother() < 0)
425 {
426 mother = 0;
427 break;
428 }
429
430 mother = (AliVParticle*) arrayMC->At(((AliAODMCParticle*)mother)->GetMother());
431 if (!mother)
432 break;
433 }
434
435 if (mother)
436 {
437 Int_t pdgCode = ((AliAODMCParticle*)mother)->GetPdgCode();
438 if (particleSpecies == 0 && TMath::Abs(pdgCode)!=211)
439 return 0;
440 if (particleSpecies == 1 && TMath::Abs(pdgCode)!=321)
441 return 0;
442 if (particleSpecies == 2 && TMath::Abs(pdgCode)!=2212)
443 return 0;
444 if (particleSpecies == 3 && (TMath::Abs(pdgCode)==211 || TMath::Abs(pdgCode)==321 || TMath::Abs(pdgCode)==2212))
445 return 0;
446 }
1c6b6c97 447 else
448 {
449 // if mother not found, accept particle only in case of particleSpecies == 3. To include it in all or no sample is no solution
ea1919ac 450 Printf("WARNING: No mother found for particle %d:", part->GetLabel());
451 part->Print();
452
453 /*
454 // this code prints the details of the mother that is missing in the AOD
455 AliMCEventHandler* fMcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
456
457 AliMCEvent* fMcEvent = fMcHandler->MCEvent();
458
459 fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->Print();
460 fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Print();
461 Printf("eta = %f", fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Eta());
462 */
463
1c6b6c97 464 if (particleSpecies != 3)
465 return 0;
466 }
b1831bcb 467 }
a75aacd6 468
469 }else if (obj->InheritsFrom("TObjArray")){ // list of AliVParticle
ea1919ac 470 TObjArray *array = static_cast<TObjArray*>(obj);
a75aacd6 471 part = (AliVParticle*)array->At( ipart );
472 if (!part)return 0;
473 }else if (obj->InheritsFrom("AliMCEvent")){ // MC PARTICLE
ea1919ac 474 AliMCEvent* mcEvent = static_cast<AliMCEvent*>(obj);
a75aacd6 475 part = mcEvent->GetTrack( ipart );
476 if (!part) return 0;
477 // eventually only primaries
478 if (onlyprimaries && !( mcEvent->IsPhysicalPrimary(ipart)) )return 0;
479 // eventually only hadrons
480 //TO-DO
481 /*if (fOnlyHadrons){
482 Int_t pdgCode = part->GetPdgCode();
483 Bool_t isHadron = TMath::Abs(pdgCode)==211 || // Pion
484 TMath::Abs(pdgCode)==2212 || // Proton
485 TMath::Abs(pdgCode)==321; // Kaon
486 if (!isHadron) return 0;
487 }
488 */
eed401dc 489 if (particleSpecies != -1) {
490 // find the primary mother
491 AliMCParticle* mother = (AliMCParticle*) part;
492// Printf("");
493 while (!mcEvent->IsPhysicalPrimary(mother->GetLabel()))
494 {
495// Printf("pdg = %d; mother = %d", mother->PdgCode(), mother->GetMother());
496 if (mother->GetMother() < 0)
497 {
498 mother = 0;
499 break;
500 }
501
502 mother = (AliMCParticle*) mcEvent->GetTrack(mother->GetMother());
503 if (!mother)
504 break;
505 }
506
507 if (mother)
508 {
509 Int_t pdgCode = mother->PdgCode();
510 if (particleSpecies == 0 && TMath::Abs(pdgCode)!=211)
511 return 0;
512 if (particleSpecies == 1 && TMath::Abs(pdgCode)!=321)
513 return 0;
514 if (particleSpecies == 2 && TMath::Abs(pdgCode)!=2212)
515 return 0;
516 if (particleSpecies == 3 && (TMath::Abs(pdgCode)==211 || TMath::Abs(pdgCode)==321 || TMath::Abs(pdgCode)==2212))
517 return 0;
518 }
519 else
520 {
521 // if mother not found, accept particle only in case of particleSpecies == 3. To include it in all or no sample is no solution
522 Printf("WARNING: No mother found for particle %d:", part->GetLabel());
523 //part->Dump();
524 //part->Print();
525
526 /*
527 // this code prints the details of the mother that is missing in the AOD
528 AliMCEventHandler* fMcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
529
530 AliMCEvent* fMcEvent = fMcHandler->MCEvent();
531
532 fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->Print();
533 fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Print();
534 Printf("eta = %f", fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Eta());
535 */
536
537 if (particleSpecies != 3)
538 return 0;
539 }
540 }
a75aacd6 541 }else if (obj->InheritsFrom("AliAODEvent")){ // RECO AOD TRACKS
ea1919ac 542 AliAODEvent *aodEvent = static_cast<AliAODEvent*>(obj);
a75aacd6 543 part = aodEvent->GetTrack(ipart);
544 // track selection cuts
144bd037 545 if ( !(((AliAODTrack*)part)->TestFilterBit(fFilterBit)) ) return 0;
546 //if ( !(((AliAODTrack*)part)->TestFilterBit(fFilterBit)) && !(((AliAODTrack*)part)->TestFilterBit(32)) ) return 0;
a75aacd6 547 // only primary candidates
548 //if ( ((AliAODTrack*)part)->IsPrimaryCandidate() )return 0;
549 // eventually only hadrons
550 if (fOnlyHadrons){
551 Bool_t isHadron = ((AliAODTrack*)part)->GetMostProbablePID()==AliAODTrack::kPion ||
552 ((AliAODTrack*)part)->GetMostProbablePID()==AliAODTrack::kKaon ||
553 ((AliAODTrack*)part)->GetMostProbablePID()==AliAODTrack::kProton;
554 if (!isHadron) return 0;
555 }
556
557 }else if (obj->InheritsFrom("AliESDEvent")){ // RECO ESD TRACKS
ea1919ac 558 AliESDEvent *esdEvent = static_cast<AliESDEvent*>(obj);
a75aacd6 559 part = esdEvent->GetTrack(ipart);
560 if (!part)return 0;
561 // track selection cuts
e0331fd9 562
563 if (!( ApplyCuts(part)) )
564 return 0;
a75aacd6 565
1bba939a 566 if (fFilterBit == 128 || fFilterBit == 256 || fFilterBit == 512 || fFilterBit == 1024)
85bfac17 567 {
568 // create TPC only tracks constrained to the SPD vertex
569
570 const AliESDVertex *vtxSPD = esdEvent->GetPrimaryVertexSPD();
571
572 AliESDtrack* track = AliESDtrackCuts::GetTPCOnlyTrack(esdEvent, ipart);
573 if(!track) return 0;
574
85bfac17 575 if(track->Pt()>0.){
576 // only constrain tracks above threshold
577 AliExternalTrackParam exParam;
578 // take the B-feild from the ESD, no 3D fieldMap available at this point
579 Bool_t relate = kFALSE;
caa8c7b2 580 relate = track->RelateToVertexTPC(vtxSPD,esdEvent->GetMagneticField(),kVeryBig,&exParam);
85bfac17 581 if(!relate)
582 {
583// Printf("relating failed");
584 delete track;
585 return 0;
586 }
587 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
588 }
589
590 part = track;
591 }
a75aacd6 592
593 // eventually only hadrons
594 //TO-DO
595 /*if (fOnlyHadrons){
596 Bool_t isHadron = ((AliAODTrack*)part)->GetMostProbablePID()==AliAODTrack::kPion ||
597 ((AliAODTrack*)part)->GetMostProbablePID()==AliAODTrack::kKaon ||
598 ((AliAODTrack*)part)->GetMostProbablePID()==AliAODTrack::kProton;
599 if (!isHadron) return 0;
600 }
601 */
602 }else {
603 if (fDebug > 1) AliFatal(" Analysis type not defined !!! ");
604 return 0;
605 }
606
607 // only charged
608 if (!part->Charge())return 0;
609
610 return part;
611}
612
613
614//-------------------------------------------------------------------
615void AliAnalyseLeadingTrackUE::QSortTracks(TObjArray &a, Int_t first, Int_t last)
616{
617 // Sort array of TObjArray of tracks by Pt using a quicksort algorithm.
618
619 static TObject *tmp;
620 static int i; // "static" to save stack space
621 int j;
622
623 while (last - first > 1) {
624 i = first;
625 j = last;
626 for (;;) {
627 while (++i < last && ((AliVParticle*)a[i])->Pt() > ((AliVParticle*)a[first])->Pt() )
628 ;
629 while (--j > first && ((AliVParticle*)a[j])->Pt() < ((AliVParticle*)a[first])->Pt() )
630 ;
631 if (i >= j)
632 break;
633
634 tmp = a[i];
635 a[i] = a[j];
636 a[j] = tmp;
637 }
638 if (j == first) {
639 ++first;
640 continue;
641 }
642 tmp = a[first];
643 a[first] = a[j];
644 a[j] = tmp;
645 if (j - first < last - (j + 1)) {
646 QSortTracks(a, first, j);
647 first = j + 1; // QSortTracks(j + 1, last);
648 } else {
649 QSortTracks(a, j + 1, last);
650 last = j; // QSortTracks(first, j);
651 }
652 }
653}
654
655//____________________________________________________________________
656TObjArray* AliAnalyseLeadingTrackUE::SortRegions(const AliVParticle* leading, TObject* obj, TObject* arrayMC, Bool_t onlyprimaries)
657{
658
659 // Assign particles to towards, away or transverse regions.
660 // Returns a lists of particles for each region.
661
662 static const Double_t k60rad = 60.*TMath::Pi()/180.;
663 static const Double_t k120rad = 120.*TMath::Pi()/180.;
664
665 // Define output lists of particles
666 TList *toward = new TList();
667 TList *away = new TList();
668 // Two transverse regions, for the moment those are not yet MIN and MAX!!!
669 // MIN and MAX can be sorted in GetMinMaxRegion function
670 TList *transverse1 = new TList();
671 TList *transverse2 = new TList();
672
673 TObjArray *regionParticles = new TObjArray;
05c47aff 674 regionParticles->SetOwner(kTRUE);
675
a75aacd6 676 regionParticles->AddLast(toward);
677 regionParticles->AddLast(away);
678 regionParticles->AddLast(transverse1);
679 regionParticles->AddLast(transverse2);
680
681 if (!leading)
682 return regionParticles;
683
684 // Switch to vector for leading particle
685 TVector3 leadVect(leading->Px(),leading->Py(),leading->Pz());
686
687 Int_t nTracks = NParticles(obj);
688 if( !nTracks ) return 0;
689 // Loop over tracks
690 for (Int_t ipart=0; ipart<nTracks; ++ipart) {
691 AliVParticle* part = ParticleWithCuts(obj, ipart);
692 if (!part)continue;
693 //Switch to vectors for particles
694 TVector3 partVect(part->Px(), part->Py(), part->Pz());
695
696 Int_t region = 0;
697 if( TMath::Abs(partVect.Eta()) > fTrackEtaCut ) continue;
698 // transverse regions
699 if (leadVect.DeltaPhi(partVect) < -k60rad && leadVect.DeltaPhi(partVect) > -k120rad )region = -1; //left
700 if (leadVect.DeltaPhi(partVect) > k60rad && leadVect.DeltaPhi(partVect) < k120rad ) region = 1; //right
701
702 if (TMath::Abs(leadVect.DeltaPhi(partVect)) < k60rad ) region = 2; //forward
703 if (TMath::Abs(leadVect.DeltaPhi(partVect)) > k120rad ) region = -2; //backward
704
705 // skip leading particle
b1831bcb 706 if (leading == part)
a75aacd6 707 continue;
708
709 if (!region)continue;
710 if (arrayMC && arrayMC->InheritsFrom("TClonesArray") && obj->InheritsFrom("AliAODEvent")){
711 Int_t label = ((AliAODTrack*)part)->GetLabel();
712 // re-define part as the matched MC particle
713 part = (AliAODMCParticle*)ParticleWithCuts(arrayMC, TMath::Abs(label),onlyprimaries);
714 if (!part)continue;
715 // skip leading particle
716 if (leading == part)
717 continue;
718 }
719 if (arrayMC && arrayMC->InheritsFrom("AliMCEvent") && obj->InheritsFrom("AliESDEvent")){
720 Int_t label = ((AliESDtrack*)part)->GetLabel();
721 // look for the matched MC particle (but do not re-define part)
722 if (!ParticleWithCuts(arrayMC, TMath::Abs(label),onlyprimaries)) continue;
723 }
724
725 if ( region == 1 ) transverse1->Add(part);
726 if ( region == -1 ) transverse2->Add(part);
727 if ( region == 2 ) toward->Add(part);
728 if ( region == -2 ) away->Add(part);
729
730 }//end loop on tracks
731
732 return regionParticles;
733
734}
735
736
737//____________________________________________________________________
738Bool_t AliAnalyseLeadingTrackUE::TriggerSelection(const TObject* obj)
739{
1c6b6c97 740 if (!obj) // MC
741 return kFALSE;
a75aacd6 742
85bfac17 743 // Use AliPhysicsSelection to select good events, works for ESD and AOD
5cabd2cd 744 if (!(((AliInputEventHandler*)obj)->IsEventSelected()&(fEventSelection)))
85bfac17 745 return kFALSE;
05c47aff 746
a75aacd6 747 return kTRUE;
a75aacd6 748}
749
750//____________________________________________________________________
751Bool_t AliAnalyseLeadingTrackUE::VertexSelection(const TObject* obj, Int_t ntracks, Double_t zed)
752{
753
754 //Require 1 vertex (no TPC stand-alone) with a minimum number of tracks and z-coordinate in a limited range
755
756 if (obj->InheritsFrom("AliAODEvent")){
757 Int_t nVertex = ((AliAODEvent*)obj)->GetNumberOfVertices();
758 if( nVertex > 0 ) {
759 AliAODVertex* vertex = (AliAODVertex*)((AliAODEvent*)obj)->GetPrimaryVertex();
144bd037 760 Int_t nTracksPrim = vertex->GetNContributors();
a75aacd6 761 Double_t zVertex = vertex->GetZ();
762 if (fDebug > 1)AliInfo(Form(" Vertex in = %f with %d particles by %s data ...",zVertex,nTracksPrim,vertex->GetName()));
144bd037 763 // Reject TPC only vertex
764 TString name(vertex->GetName());
765 if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex"))return kFALSE;
766
767 // Select a quality vertex by number of tracks?
a75aacd6 768 if( nTracksPrim < ntracks || TMath::Abs(zVertex) > zed ) {
769 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event REJECTED ...");
770 return kFALSE;
771 }
772 // TODO remove vertexer Z events with dispersion > 0.02: Doesn't work for AOD at present
773 //if (strcmp(vertex->GetTitle(), "AliVertexerZ") == 0 && vertex->GetDispersion() > 0.02)
774 // return kFALSE;
775 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event ACCEPTED...");
776 } else {
777 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event REJECTED ...");
778 return kFALSE;
779 }
780 }
781
782 if (obj->InheritsFrom("AliMCEvent"))
783 {
784 if (TMath::Abs(((AliMCEvent*) obj)->GetPrimaryVertex()->GetZ()) > zed)
785 {
786 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event (based on MC) REJECTED ...");
787 return kFALSE;
788 }
789 }
144bd037 790
7dd4dec4 791 if (obj->InheritsFrom("AliAODMCHeader"))
792 {
793 if (TMath::Abs(((AliAODMCHeader*) obj)->GetVtxZ()) > zed)
794 {
795 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event (based on MC) REJECTED ...");
796 return kFALSE;
797 }
798 }
799
144bd037 800 // ESD case for DCA studies
801 if (obj->InheritsFrom("AliESDEvent")){
802 AliESDVertex* vertex = (AliESDVertex*)((AliESDEvent*)obj)->GetPrimaryVertex();
803 if ( vertex){
804 Int_t nTracksPrim = vertex->GetNContributors();
805 Double_t zVertex = vertex->GetZ();
806 if (fDebug > 1)AliInfo(Form(" Vertex in = %f with %d particles by %s data ...",zVertex,nTracksPrim,vertex->GetName()));
807 // Reject SPD or TPC only vertex
808 TString name(vertex->GetName());
809 if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex"))return kFALSE;
810
811 // Select a quality vertex by number of tracks?
812 if( nTracksPrim < ntracks || TMath::Abs(zVertex) > zed ) {
813 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event REJECTED ...");
814 return kFALSE;
815 }
816 // TODO remove vertexer Z events with dispersion > 0.02: Doesn't work for AOD at present
817 //if (strcmp(vertex->GetTitle(), "AliVertexerZ") == 0 && vertex->GetDispersion() > 0.02)
818 // return kFALSE;
819 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event ACCEPTED...");
820 } else {
821 if (fDebug > 1) AliInfo(" Primary-vertex Selection: event REJECTED ...");
822 return kFALSE;
823 }
824 }
a75aacd6 825
826 return kTRUE;
827}