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