6da8866d20905cf10d91f893572c53b1574588ca
[u/mrichter/AliRoot.git] / PWGCF / EBYE / BalanceFunctions / AliAnalysisTaskTriggeredBF.cxx
1 #include "TChain.h"\r
2 #include "TList.h"\r
3 #include "TCanvas.h"\r
4 #include "TLorentzVector.h"\r
5 #include "TGraphErrors.h"\r
6 #include "TH1F.h"\r
7 #include "TH2F.h"\r
8 #include "TArrayF.h"\r
9 #include "TF1.h"\r
10 #include "TRandom.h"\r
11 \r
12 #include "AliLog.h"\r
13 \r
14 #include "AliAnalysisTaskSE.h"\r
15 #include "AliAnalysisManager.h"\r
16 \r
17 #include "AliESDVertex.h"\r
18 #include "AliESDEvent.h"\r
19 #include "AliESDInputHandler.h"\r
20 #include "AliAODEvent.h"\r
21 #include "AliAODTrack.h"\r
22 #include "AliAODInputHandler.h"\r
23 #include "AliGenEventHeader.h"\r
24 #include "AliGenHijingEventHeader.h"\r
25 #include "AliMCEventHandler.h"\r
26 #include "AliMCEvent.h"\r
27 #include "AliMixInputEventHandler.h"\r
28 #include "AliStack.h"\r
29 \r
30 #include "TH2D.h"    \r
31 #include "AliTHn.h"             \r
32 \r
33 #include "AliEventPoolManager.h" \r
34 \r
35 #include "AliAnalysisTaskTriggeredBF.h"\r
36 #include "AliBalanceTriggered.h"\r
37 \r
38 \r
39 // Analysis task for the TriggeredBF code\r
40 // Authors: Panos.Christakoglou@nikhef.nl, m.weber@cern.ch\r
41 \r
42 ClassImp(AliAnalysisTaskTriggeredBF)\r
43 \r
44 //________________________________________________________________________\r
45 AliAnalysisTaskTriggeredBF::AliAnalysisTaskTriggeredBF(const char *name) \r
46 : AliAnalysisTaskSE(name), \r
47   fBalance(0),\r
48   fRunShuffling(kFALSE),\r
49   fShuffledBalance(0),\r
50   fRunMixing(kFALSE),\r
51   fMixingTracks(50000),\r
52   fMixedBalance(0),\r
53   fList(0),\r
54   fListTriggeredBF(0),\r
55   fListTriggeredBFS(0),\r
56   fHistListPIDQA(0),\r
57   fHistEventStats(0),\r
58   fHistCentStats(0),\r
59   fHistTriggerStats(0),\r
60   fHistTrackStats(0),\r
61   fHistVx(0),\r
62   fHistVy(0),\r
63   fHistVz(0),\r
64   fHistClus(0),\r
65   fHistDCA(0),\r
66   fHistChi2(0),\r
67   fHistPt(0),\r
68   fHistEta(0),\r
69   fHistPhi(0),\r
70   fHistPhiBefore(0),\r
71   fHistPhiAfter(0),\r
72   fHistV0M(0),\r
73   fHistRefTracks(0),\r
74   fCentralityEstimator("V0M"),\r
75   fUseCentrality(kFALSE),\r
76   fCentralityPercentileMin(0.), \r
77   fCentralityPercentileMax(5.),\r
78   fImpactParameterMin(0.),\r
79   fImpactParameterMax(20.),\r
80   fUseMultiplicity(kFALSE),\r
81   fNumberOfAcceptedTracksMin(0),\r
82   fNumberOfAcceptedTracksMax(10000),\r
83   fHistNumberOfAcceptedTracks(0),\r
84   fUseOfflineTrigger(kFALSE),\r
85   fVxMax(0.3),\r
86   fVyMax(0.3),\r
87   fVzMax(10.),\r
88   nAODtrackCutBit(128),\r
89   fPtMin(0.3),\r
90   fPtMax(1.5),\r
91   fEtaMin(-0.8),\r
92   fEtaMax(-0.8),\r
93   fDCAxyCut(-1),\r
94   fDCAzCut(-1),\r
95   fTPCchi2Cut(-1),\r
96   fNClustersTPCCut(-1)\r
97 {\r
98   // Constructor\r
99   // Define input and output slots here\r
100   // Input slot #0 works with a TChain\r
101   DefineInput(0, TChain::Class());\r
102   // Output slot #0 writes into a TH1 container\r
103   DefineOutput(1, TList::Class());\r
104   DefineOutput(2, TList::Class());\r
105   DefineOutput(3, TList::Class());\r
106 }\r
107 \r
108 //________________________________________________________________________\r
109 AliAnalysisTaskTriggeredBF::~AliAnalysisTaskTriggeredBF() {\r
110 \r
111   // Destructor\r
112 \r
113 }\r
114 \r
115 //________________________________________________________________________\r
116 void AliAnalysisTaskTriggeredBF::UserCreateOutputObjects() {\r
117   // Create histograms\r
118   // Called once\r
119   if(!fBalance) {\r
120     fBalance = new AliBalanceTriggered();\r
121     fBalance->SetAnalysisLevel("AOD");\r
122   }\r
123   if(fRunShuffling) {\r
124     if(!fShuffledBalance) {\r
125       fShuffledBalance = new AliBalanceTriggered();\r
126       fShuffledBalance->SetAnalysisLevel("AOD");\r
127     }\r
128   }\r
129   if(fRunMixing) {\r
130     if(!fMixedBalance) {\r
131       fMixedBalance = new AliBalanceTriggered();\r
132       fMixedBalance->SetAnalysisLevel("AOD");\r
133     }\r
134   }\r
135 \r
136   //QA list\r
137   fList = new TList();\r
138   fList->SetName("listQA");\r
139   fList->SetOwner();\r
140 \r
141   //Balance Function list\r
142   fListTriggeredBF = new TList();\r
143   fListTriggeredBF->SetName("listTriggeredBF");\r
144   fListTriggeredBF->SetOwner();\r
145 \r
146   if(fRunShuffling) {\r
147     fListTriggeredBFS = new TList();\r
148     fListTriggeredBFS->SetName("listTriggeredBFShuffled");\r
149     fListTriggeredBFS->SetOwner();\r
150   }\r
151   if(fRunMixing) {\r
152     fListTriggeredBFM = new TList();\r
153     fListTriggeredBFM->SetName("listTriggeredBFMixed");\r
154     fListTriggeredBFM->SetOwner();\r
155   }\r
156   \r
157   \r
158   //Event stats.\r
159   TString gCutName[4] = {"Total","Offline trigger",\r
160                          "Vertex","Analyzed"};\r
161   fHistEventStats = new TH1F("fHistEventStats",\r
162                              "Event statistics;;N_{events}",\r
163                              4,0.5,4.5);\r
164   for(Int_t i = 1; i <= 4; i++)\r
165     fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());\r
166   fList->Add(fHistEventStats);\r
167   \r
168   TString gCentName[9] = {"V0M","FMD","TRK","TKL","CL0","CL1","V0MvsFMD","TKLvsV0M","ZEMvsZDC"};\r
169   fHistCentStats = new TH2F("fHistCentStats",\r
170                             "Centrality statistics;;Cent percentile",\r
171                             9,-0.5,8.5,220,-5,105);\r
172   for(Int_t i = 1; i <= 9; i++)\r
173     fHistCentStats->GetXaxis()->SetBinLabel(i,gCentName[i-1].Data());\r
174   fList->Add(fHistCentStats);\r
175   \r
176   fHistTriggerStats = new TH1F("fHistTriggerStats","Trigger statistics;TriggerBit;N_{events}",130,0,130);\r
177   fList->Add(fHistTriggerStats);\r
178   \r
179   fHistTrackStats = new TH1F("fHistTrackStats","Event statistics;TrackFilterBit;N_{events}",130,0,130);\r
180   fList->Add(fHistTrackStats);\r
181 \r
182   fHistNumberOfAcceptedTracks = new TH1F("fHistNumberOfAcceptedTracks",";N_{acc.};Entries",4001,-0.5,4000.5);\r
183   fList->Add(fHistNumberOfAcceptedTracks);\r
184 \r
185   // Vertex distributions\r
186   fHistVx = new TH1F("fHistVx","Primary vertex distribution - x coordinate;V_{x} (cm);Entries",100,-0.5,0.5);\r
187   fList->Add(fHistVx);\r
188   fHistVy = new TH1F("fHistVy","Primary vertex distribution - y coordinate;V_{y} (cm);Entries",100,-0.5,0.5);\r
189   fList->Add(fHistVy);\r
190   fHistVz = new TH1F("fHistVz","Primary vertex distribution - z coordinate;V_{z} (cm);Entries",100,-20.,20.);\r
191   fList->Add(fHistVz);\r
192 \r
193   // QA histograms\r
194   fHistClus = new TH2F("fHistClus","# Cluster (TPC vs. ITS)",10,0,10,200,0,200);\r
195   fList->Add(fHistClus);\r
196   fHistChi2 = new TH1F("fHistChi2","Chi2/NDF distribution",200,0,10);\r
197   fList->Add(fHistChi2);\r
198   fHistDCA  = new TH2F("fHistDCA","DCA (xy vs. z)",400,-5,5,400,-5,5); \r
199   fList->Add(fHistDCA);\r
200   fHistPt   = new TH1F("fHistPt","p_{T} distribution",200,0,10);\r
201   fList->Add(fHistPt);\r
202   fHistEta  = new TH1F("fHistEta","#eta distribution",200,-2,2);\r
203   fList->Add(fHistEta);\r
204   fHistPhi  = new TH1F("fHistPhi","#phi distribution",200,-20,380);\r
205   fList->Add(fHistPhi);\r
206   fHistPhiBefore  = new TH1F("fHistPhiBefore","#phi distribution",200,0.,2*TMath::Pi());\r
207   fList->Add(fHistPhiBefore);\r
208   fHistPhiAfter  = new TH1F("fHistPhiAfter","#phi distribution",200,0.,2*TMath::Pi());\r
209   fList->Add(fHistPhiAfter);\r
210   fHistV0M  = new TH2F("fHistV0M","V0 Multiplicity C vs. A",500, 0, 20000, 500, 0, 20000);\r
211   fList->Add(fHistV0M);\r
212   TString gRefTrackName[6] = {"tracks","tracksPos","tracksNeg","tracksTPConly","clusITS0","clusITS1"};\r
213   fHistRefTracks  = new TH2F("fHistRefTracks","Nr of Ref tracks/event vs. ref track estimator;;Nr of tracks",6, 0, 6, 400, 0, 20000);\r
214   for(Int_t i = 1; i <= 6; i++)\r
215     fHistRefTracks->GetXaxis()->SetBinLabel(i,gRefTrackName[i-1].Data());\r
216   fList->Add(fHistRefTracks);\r
217 \r
218 \r
219 \r
220   // Balance function histograms\r
221   // Initialize histograms if not done yet\r
222   if(!fBalance->GetHistNp()){\r
223     AliWarning("Histograms not yet initialized! --> Will be done now");\r
224     AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
225     fBalance->InitHistograms();\r
226   }\r
227 \r
228   if(fRunShuffling) {\r
229     if(!fShuffledBalance->GetHistNp()) {\r
230       AliWarning("Histograms (shuffling) not yet initialized! --> Will be done now");\r
231       AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
232       fShuffledBalance->InitHistograms();\r
233     }\r
234   }\r
235 \r
236   if(fRunMixing) {\r
237     if(!fMixedBalance->GetHistNp()) {\r
238       AliWarning("Histograms (mixing) not yet initialized! --> Will be done now");\r
239       AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
240       fMixedBalance->InitHistograms();\r
241     }\r
242   }\r
243 \r
244   fListTriggeredBF->Add(fBalance->GetHistNp());\r
245   fListTriggeredBF->Add(fBalance->GetHistNn());\r
246   fListTriggeredBF->Add(fBalance->GetHistNpn());\r
247   fListTriggeredBF->Add(fBalance->GetHistNnn());\r
248   fListTriggeredBF->Add(fBalance->GetHistNpp());\r
249   fListTriggeredBF->Add(fBalance->GetHistNnp());\r
250   \r
251   if(fRunShuffling) {\r
252     fListTriggeredBFS->Add(fShuffledBalance->GetHistNp());\r
253     fListTriggeredBFS->Add(fShuffledBalance->GetHistNn());\r
254     fListTriggeredBFS->Add(fShuffledBalance->GetHistNpn());\r
255     fListTriggeredBFS->Add(fShuffledBalance->GetHistNnn());\r
256     fListTriggeredBFS->Add(fShuffledBalance->GetHistNpp());\r
257     fListTriggeredBFS->Add(fShuffledBalance->GetHistNnp());\r
258   }  \r
259 \r
260   if(fRunMixing) {\r
261     fListTriggeredBFM->Add(fMixedBalance->GetHistNp());\r
262     fListTriggeredBFM->Add(fMixedBalance->GetHistNn());\r
263     fListTriggeredBFM->Add(fMixedBalance->GetHistNpn());\r
264     fListTriggeredBFM->Add(fMixedBalance->GetHistNnn());\r
265     fListTriggeredBFM->Add(fMixedBalance->GetHistNpp());\r
266     fListTriggeredBFM->Add(fMixedBalance->GetHistNnp());\r
267   }  \r
268 \r
269 \r
270   // Event Mixing\r
271   Int_t trackDepth = fMixingTracks; \r
272   Int_t poolsize   = 1000;  // Maximum number of events, ignored in the present implemented of AliEventPoolManager\r
273    \r
274   Double_t centralityBins[] = {0,1,2,3,4,5,6,7,8,9,10,15,20,25,30,35,40,45,50,55,60,65,70,75,80,90,100}; // SHOULD BE DEDUCED FROM CREATED ALITHN!!!\r
275   Double_t* centbins = centralityBins;\r
276   Int_t nCentralityBins  = 26;\r
277 \r
278   \r
279   // bins for second buffer are shifted by 100 cm\r
280   Double_t vertexBins[] = {-10, -7, -5, -3, -1, 1, 3, 5, 7, 10}; // SHOULD BE DEDUCED FROM CREATED ALITHN!!!\r
281   Double_t* vtxbins = vertexBins;\r
282   Int_t nVertexBins  = 9;\r
283 \r
284   fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, centbins, nVertexBins, vtxbins);\r
285 \r
286 \r
287   // Post output data.\r
288   PostData(1, fList);\r
289   PostData(2, fListTriggeredBF);\r
290   if(fRunShuffling) PostData(3, fListTriggeredBFS);\r
291   if(fRunMixing) PostData(4, fListTriggeredBFM);\r
292 }\r
293 \r
294 //________________________________________________________________________\r
295 void AliAnalysisTaskTriggeredBF::UserExec(Option_t *) {\r
296   // Main loop\r
297   // Called for each event\r
298 \r
299   TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
300   Float_t fCentrality = 0.;  \r
301 \r
302   // -------------------------------------------------------------                   \r
303   // AOD analysis (vertex and track cuts also here!!!!)\r
304   if(gAnalysisLevel == "AOD") {\r
305     AliVEvent* eventMain = dynamic_cast<AliVEvent*>(InputEvent()); \r
306     if(!eventMain) {\r
307       AliError("eventMain not available");\r
308       return;\r
309     }\r
310 \r
311     // check event cuts and fill event histograms\r
312     if((fCentrality = IsEventAccepted(eventMain)) < 0){\r
313       return;\r
314     }\r
315     \r
316     // get the accepted tracks in main event\r
317     TObjArray *tracksMain = GetAcceptedTracks(eventMain);\r
318 \r
319     // store charges of all accepted tracks, shuffle and reassign (two extra loops!)\r
320     TObjArray* tracksShuffled = NULL;\r
321     if(fRunShuffling){\r
322       tracksShuffled = GetShuffledTracks(tracksMain);\r
323     }\r
324     \r
325     // Event mixing --> UPDATE POOL IS MISSING!!!\r
326     if (fRunMixing)\r
327       {\r
328         // 1. First get an event pool corresponding in mult (cent) and\r
329         //    zvertex to the current event. Once initialized, the pool\r
330         //    should contain nMix (reduced) events. This routine does not\r
331         //    pre-scan the chain. The first several events of every chain\r
332         //    will be skipped until the needed pools are filled to the\r
333         //    specified depth. If the pool categories are not too rare, this\r
334         //    should not be a problem. If they are rare, you could lose`\r
335         //    statistics.\r
336         \r
337         // 2. Collect the whole pool's content of tracks into one TObjArray\r
338         //    (bgTracks), which is effectively a single background super-event.\r
339         \r
340         // 3. The reduced and bgTracks arrays must both be passed into\r
341         //    FillCorrelations(). Also nMix should be passed in, so a weight\r
342         //    of 1./nMix can be applied.\r
343         \r
344         AliEventPool* pool = fPoolMgr->GetEventPool(fCentrality, eventMain->GetPrimaryVertex()->GetZ());\r
345         \r
346         if (!pool)\r
347           AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCentrality, eventMain->GetPrimaryVertex()->GetZ()));\r
348         \r
349         //pool->SetDebug(1);\r
350         \r
351         if (pool->IsReady() || pool->NTracksInPool() > fMixingTracks / 10 || pool->GetCurrentNEvents() >= 5) \r
352           {\r
353             \r
354             Int_t nMix = pool->GetCurrentNEvents();\r
355             cout << "nMix = " << nMix << " tracks in pool = " << pool->NTracksInPool() << endl;\r
356             \r
357             //((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(2);\r
358             //((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());\r
359             //if (pool->IsReady())\r
360             //((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(3);\r
361             \r
362             // Fill mixed-event histos here  \r
363             for (Int_t jMix=0; jMix<nMix; jMix++) \r
364               {\r
365                 TObjArray* tracksMixed = pool->GetEvent(jMix);\r
366                 fMixedBalance->FillBalance(fCentrality,tracksMixed); //NOW ONLY THE MIXED EVENT ITSELF IS FILLED --> DO ONE TRACK OF MAIN AND ONE OF MIXED (LIKE UEHISTOGRAMS!!!!)\r
367               }\r
368           }\r
369       }\r
370     \r
371     // calculate balance function\r
372     fBalance->FillBalance(fCentrality,tracksMain);//,chargeVectorMixed); // here comes the mixing... in some time\r
373     \r
374     // calculate shuffled balance function\r
375     if(fRunShuffling && tracksShuffled != NULL) {\r
376        fShuffledBalance->FillBalance(fCentrality,tracksShuffled);\r
377     }\r
378     \r
379   }//AOD analysis\r
380   else{\r
381     AliError("Triggered Balance Function analysis only for AODs!");\r
382   }\r
383 }     \r
384 \r
385 //________________________________________________________________________\r
386 Float_t AliAnalysisTaskTriggeredBF::IsEventAccepted(AliVEvent *event){\r
387   // Checks the Event cuts\r
388   // Fills Event statistics histograms\r
389   \r
390   // event selection done in AliAnalysisTaskSE::Exec() --> this is not used\r
391   fHistEventStats->Fill(1); //all events\r
392 \r
393   Bool_t isSelectedMain = kTRUE;\r
394   Float_t fCentrality = -1.;\r
395   TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
396   \r
397   if(fUseOfflineTrigger)\r
398     isSelectedMain = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();\r
399   \r
400   if(isSelectedMain) {\r
401     fHistEventStats->Fill(2); //triggered events\r
402     \r
403     //Centrality stuff \r
404     if(fUseCentrality) {\r
405       if(gAnalysisLevel == "AOD") { //centrality in AOD header\r
406         AliAODHeader *header = (AliAODHeader*) event->GetHeader();\r
407         fCentrality = header->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());\r
408 \r
409         // QA for centrality estimators\r
410         fHistCentStats->Fill(0.,header->GetCentralityP()->GetCentralityPercentile("V0M"));\r
411         fHistCentStats->Fill(1.,header->GetCentralityP()->GetCentralityPercentile("FMD"));\r
412         fHistCentStats->Fill(2.,header->GetCentralityP()->GetCentralityPercentile("TRK"));\r
413         fHistCentStats->Fill(3.,header->GetCentralityP()->GetCentralityPercentile("TKL"));\r
414         fHistCentStats->Fill(4.,header->GetCentralityP()->GetCentralityPercentile("CL0"));\r
415         fHistCentStats->Fill(5.,header->GetCentralityP()->GetCentralityPercentile("CL1"));\r
416         fHistCentStats->Fill(6.,header->GetCentralityP()->GetCentralityPercentile("V0MvsFMD"));\r
417         fHistCentStats->Fill(7.,header->GetCentralityP()->GetCentralityPercentile("TKLvsV0M"));\r
418         fHistCentStats->Fill(8.,header->GetCentralityP()->GetCentralityPercentile("ZEMvsZDC"));\r
419         \r
420         // centrality QA (V0M)\r
421         fHistV0M->Fill(event->GetVZEROData()->GetMTotV0A(), event->GetVZEROData()->GetMTotV0C());\r
422         \r
423         // centrality QA (reference tracks)\r
424         fHistRefTracks->Fill(0.,header->GetRefMultiplicity());\r
425         fHistRefTracks->Fill(1.,header->GetRefMultiplicityPos());\r
426         fHistRefTracks->Fill(2.,header->GetRefMultiplicityNeg());\r
427         fHistRefTracks->Fill(3.,header->GetTPConlyRefMultiplicity());\r
428         fHistRefTracks->Fill(4.,header->GetNumberOfITSClusters(0));\r
429         fHistRefTracks->Fill(5.,header->GetNumberOfITSClusters(1));\r
430         fHistRefTracks->Fill(6.,header->GetNumberOfITSClusters(2));\r
431         fHistRefTracks->Fill(7.,header->GetNumberOfITSClusters(3));\r
432         fHistRefTracks->Fill(8.,header->GetNumberOfITSClusters(4));\r
433       }\r
434     }\r
435     \r
436     \r
437     const AliVVertex *vertex = event->GetPrimaryVertex();\r
438     \r
439     if(vertex) {\r
440       Double32_t fCov[6];\r
441       vertex->GetCovarianceMatrix(fCov);\r
442       if(vertex->GetNContributors() > 0) {\r
443         if(fCov[5] != 0) {\r
444           fHistEventStats->Fill(3); //events with a proper vertex\r
445           if(TMath::Abs(vertex->GetX()) < fVxMax) {\r
446             if(TMath::Abs(vertex->GetY()) < fVyMax) {\r
447               if(TMath::Abs(vertex->GetZ()) < fVzMax) {\r
448                 fHistEventStats->Fill(4); //analyzed events\r
449                 fHistVx->Fill(vertex->GetX());\r
450                 fHistVy->Fill(vertex->GetY());\r
451                 fHistVz->Fill(vertex->GetZ());\r
452 \r
453                 // take only events inside centrality class\r
454                 if((fCentrality > fCentralityPercentileMin) && (fCentrality < fCentralityPercentileMax)){\r
455                   return fCentrality;           \r
456                 }//centrality class\r
457               }//Vz cut\r
458             }//Vy cut\r
459           }//Vx cut\r
460         }//proper vertex resolution\r
461       }//proper number of contributors\r
462     }//vertex object valid\r
463   }//triggered event \r
464   \r
465   // in all other cases return -1 (event not accepted)\r
466   return -1;\r
467 }\r
468 \r
469 //________________________________________________________________________\r
470 TObjArray* AliAnalysisTaskTriggeredBF::GetAcceptedTracks(AliVEvent *event){\r
471   // Returns TObjArray with tracks after all track cuts (only for AOD!)\r
472   // Fills QA histograms\r
473 \r
474   //output TObjArray holding all good tracks\r
475   TObjArray* tracksAccepted = new TObjArray;\r
476   tracksAccepted->SetOwner(kTRUE);\r
477 \r
478   Double_t v_charge;\r
479   Double_t v_eta;\r
480   Double_t v_phi;\r
481   Double_t v_pt;\r
482   \r
483   // Loop over tracks in event\r
484   for (Int_t iTracks = 0; iTracks < event->GetNumberOfTracks(); iTracks++) {\r
485     AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(event->GetTrack(iTracks));\r
486     if (!aodTrack) {\r
487       AliError(Form("Could not receive track %d", iTracks));\r
488       continue;\r
489     }\r
490     \r
491     // AOD track cuts\r
492     \r
493     // For ESD Filter Information: ANALYSIS/macros/AddTaskESDfilter.C\r
494     // take only TPC only tracks \r
495     fHistTrackStats->Fill(aodTrack->GetFilterMap());\r
496     if(!aodTrack->TestFilterBit(nAODtrackCutBit)) continue;\r
497     \r
498     v_charge = aodTrack->Charge();\r
499     v_eta    = aodTrack->Eta();\r
500     v_phi    = aodTrack->Phi() * TMath::RadToDeg();\r
501     v_pt     = aodTrack->Pt();\r
502     \r
503     Float_t DCAxy = aodTrack->DCA();      // this is the DCA from global track (not exactly what is cut on)\r
504     Float_t DCAz  = aodTrack->ZAtDCA();   // this is the DCA from global track (not exactly what is cut on)\r
505     \r
506     \r
507     // Kinematics cuts from ESD track cuts\r
508     if( v_pt < fPtMin || v_pt > fPtMax)      continue;\r
509     if( v_eta < fEtaMin || v_eta > fEtaMax)  continue;\r
510     \r
511     // Extra DCA cuts (for systematic studies [!= -1])\r
512     if( fDCAxyCut != -1 && fDCAzCut != -1){\r
513       if(TMath::Sqrt((DCAxy*DCAxy)/(fDCAxyCut*fDCAxyCut)+(DCAz*DCAz)/(fDCAzCut*fDCAzCut)) > 1 ){\r
514         continue;  // 2D cut\r
515       }\r
516     }\r
517     \r
518     // Extra TPC cuts (for systematic studies [!= -1])\r
519     if( fTPCchi2Cut != -1 && aodTrack->Chi2perNDF() > fTPCchi2Cut){\r
520       continue;\r
521     }\r
522     if( fNClustersTPCCut != -1 && aodTrack->GetTPCNcls() < fNClustersTPCCut){\r
523       continue;\r
524     }\r
525     \r
526     // fill QA histograms\r
527     fHistClus->Fill(aodTrack->GetITSNcls(),aodTrack->GetTPCNcls());\r
528     fHistDCA->Fill(DCAz,DCAxy);\r
529     fHistChi2->Fill(aodTrack->Chi2perNDF());\r
530     fHistPt->Fill(v_pt);\r
531     fHistEta->Fill(v_eta);\r
532     fHistPhi->Fill(v_phi);\r
533     \r
534     // add the track to the TObjArray\r
535     tracksAccepted->Add(new AliBFBasicParticle(v_eta, v_phi, v_pt, v_charge));\r
536   }\r
537 \r
538   return tracksAccepted;\r
539 }\r
540 \r
541 //________________________________________________________________________\r
542 TObjArray* AliAnalysisTaskTriggeredBF::GetShuffledTracks(TObjArray *tracks){\r
543   // Clones TObjArray and returns it with tracks after shuffling the charges\r
544 \r
545   TObjArray* tracksShuffled = new TObjArray;\r
546   tracksShuffled->SetOwner(kTRUE);\r
547 \r
548   vector<Short_t> *chargeVector = new vector<Short_t>;   //original charge of accepted tracks \r
549 \r
550   for (Int_t i=0; i<tracks->GetEntriesFast(); i++)\r
551   {\r
552     AliVParticle* track = (AliVParticle*) tracks->At(i);\r
553     chargeVector->push_back(track->Charge());\r
554   }  \r
555  \r
556   random_shuffle(chargeVector->begin(), chargeVector->end());\r
557   \r
558   for(Int_t i = 0; i < tracks->GetEntriesFast(); i++){\r
559     AliVParticle* track = (AliVParticle*) tracks->At(i);\r
560     tracksShuffled->Add(new AliBFBasicParticle(track->Eta(), track->Phi(), track->Pt(),chargeVector->at(i)));\r
561   }\r
562    \r
563   return tracksShuffled;\r
564 }\r
565 \r
566 //________________________________________________________________________\r
567 void  AliAnalysisTaskTriggeredBF::FinishTaskOutput(){\r
568 \r
569   if (!fBalance) {\r
570     AliError("fBalance not available");\r
571     return;\r
572   }  \r
573   if(fRunShuffling) {\r
574     if (!fShuffledBalance) {\r
575       AliError("fShuffledBalance not available");\r
576       return;\r
577     }\r
578   }\r
579 \r
580 }\r
581 \r
582 //________________________________________________________________________\r
583 void AliAnalysisTaskTriggeredBF::Terminate(Option_t *) {\r
584   // Called once at the end of the query\r
585 \r
586   // not implemented ...\r
587 \r
588 }\r
589 \r
590 void AliAnalysisTaskTriggeredBF::UserExecMix(Option_t *)\r
591 {\r
592 \r
593   // not yet done for event mixing!\r
594   return;\r
595 \r
596 }\r
597 \r