AliAODEvent::GetHeader() returns AliVHeader
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / AliAnalysisTaskTwoPlusOne.cxx
1 #include <TROOT.h>
2 #include <TChain.h>
3 #include <TFile.h>
4 #include <TList.h>
5 #include <TMath.h>
6 #include <TTree.h>
7 //#include <TParameter.h>
8
9 #include "AliAnalysisTaskTwoPlusOne.h"
10 #include "AliCFParticle.h"
11 #include "AliAnalyseLeadingTrackUE.h"
12 #include "AliTwoPlusOneContainer.h"
13 #include "AliUEHist.h"
14
15 #include "AliAODHandler.h"
16 #include "AliAODInputHandler.h"
17 #include "AliVParticle.h"
18 #include "AliCFContainer.h"
19
20 #include "AliEventPoolManager.h"
21
22
23 ClassImp( AliAnalysisTaskTwoPlusOne )
24
25 //________________________________________________________________________
26 AliAnalysisTaskTwoPlusOne::AliAnalysisTaskTwoPlusOne(const char *name)
27 : AliAnalysisTaskSE(name),
28   fMixingTracks(10000),
29   fAnalyseUE(0x0),
30 // pointers to UE classes
31   fHistos(0x0),
32 // handlers and events
33   fAOD(0x0),
34   fPoolMgr(0x0),
35 // histogram settings
36   fListOfHistos(0x0), 
37 // event QA
38   fnTracksVertex(1),  // QA tracks pointing to principal vertex (= 3 default) 
39   fZVertex(7.),
40   fCentralityMethod("V0M"),
41 // track cuts
42   fTrackEtaCut(0.9),
43   fTrackEtaCutMin(-1.),
44   fPtMin(0.5),
45   fDCAXYCut(0),
46   fSharedClusterCut(-1),
47   fCrossedRowsCut(-1),
48   fFoundFractionCut(-1),
49   fFilterBit(0xFF),
50   fTrackStatus(0),
51   fThreeParticleMixed(0),
52   fCustomBinning(),
53   fAlpha(0.2)
54 {
55
56   DefineOutput(1, TList::Class());
57 }
58
59 //________________________________________________________________________
60 AliAnalysisTaskTwoPlusOne::~AliAnalysisTaskTwoPlusOne()
61 {
62     // Destructor
63
64 }
65
66 //________________________________________________________________________
67 void AliAnalysisTaskTwoPlusOne::UserCreateOutputObjects()
68 {
69   // Initialize class with main algorithms, event and track selection. 
70   fAnalyseUE = new AliAnalyseLeadingTrackUE();
71   fAnalyseUE->SetParticleSelectionCriteria(fFilterBit, kFALSE, fTrackEtaCut, fTrackEtaCutMin, fPtMin);
72   fAnalyseUE->SetDCAXYCut(fDCAXYCut);
73   fAnalyseUE->SetSharedClusterCut(fSharedClusterCut);
74   fAnalyseUE->SetCrossedRowsCut(fCrossedRowsCut);
75   fAnalyseUE->SetFoundFractionCut(fFoundFractionCut);
76   fAnalyseUE->SetTrackStatus(fTrackStatus);
77   fAnalyseUE->SetDebug(fDebug); 
78   fAnalyseUE->DefineESDCuts(fFilterBit);
79
80   fListOfHistos = new TList();
81   fListOfHistos->SetOwner(kTRUE); 
82
83   fHistos = new AliTwoPlusOneContainer("AliTwoPlusOneContainer", fCustomBinning, fAlpha);
84   fHistos->GetData()->SetTrackEtaCut(fTrackEtaCut);
85
86   fListOfHistos->Add(fHistos);
87
88   fListOfHistos->Add(new TH1F("eventStat", ";;events", 4, -0.5, 3.5));
89   fListOfHistos->Add(new TH2F("mixedDist", ";centrality;tracks;events", 101, 0, 101, 200, 0, fMixingTracks * 1.5));
90
91   PostData(1,fListOfHistos);
92
93   // Add task configuration to output list 
94   AddSettingsTree();
95   
96   // event mixing
97   Int_t poolsize   = 1000;  // Maximum number of events, ignored in the present implemention of AliEventPoolManager
98    
99   Int_t nCentralityBins  = fHistos->GetData()->GetTrackHist(AliUEHist::kToward)->GetNBins(3);
100
101   
102   Double_t* centralityBins = (Double_t*) fHistos->GetData()->GetTrackHist(AliUEHist::kToward)->GetAxis(3, 0)->GetXbins()->GetArray();
103   
104   Int_t nZvtxBins = fHistos->GetData()->GetTrackHist(AliUEHist::kToward)->GetNBins(5);
105   
106   Double_t* zvtxbin = (Double_t*) fHistos->GetData()->GetTrackHist(AliUEHist::kToward)->GetAxis(5, 0)->GetXbins()->GetArray();
107
108   fPoolMgr = new AliEventPoolManager(poolsize, fMixingTracks, nCentralityBins, centralityBins, nZvtxBins, zvtxbin);
109   fPoolMgr->SetTargetValues(fMixingTracks, 0.1, 5);
110 }
111
112 //________________________________________________________________________
113 void AliAnalysisTaskTwoPlusOne::UserExec(Option_t *)
114 {
115   // exec (per event)
116   fAnalyseUE->NextEvent();
117
118   ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(0);
119
120   // Support for AOD based analysis
121   fAOD = dynamic_cast<AliAODEvent*> (InputEvent());
122
123   Double_t centrality = 0;
124   AliCentrality *centralityObj = 0;
125
126   centralityObj = ((AliVAODHeader*)fAOD->GetHeader())->GetCentralityP();
127
128
129  if (centralityObj)
130    centrality = centralityObj->GetCentralityPercentile(fCentralityMethod);
131  else
132    centrality = -1;
133       
134   if (fAOD){
135     // remove outliers
136     if (centrality == 0)
137       {
138         if (fAOD->GetVZEROData())
139           {
140             Float_t multV0 = 0;
141             for (Int_t i=0; i<64; i++)
142               multV0 += fAOD->GetVZEROData()->GetMultiplicity(i);
143             if (multV0 < 19500)
144               {
145                 centrality = -1;
146                 AliInfo("Rejecting event due to too small V0 multiplicity");
147               }
148           }
149       }
150   }
151
152   
153   AliInfo(Form("Centrality is %f", centrality));
154
155
156   // Vertex selection *************************************************
157   if(!fAnalyseUE->VertexSelection(InputEvent(), fnTracksVertex, fZVertex)) return;
158
159   // optimization
160   if (centrality < 0)
161     return;
162
163   TObjArray* tracks = fAnalyseUE->GetAcceptedParticles(InputEvent(), 0, kTRUE, -1, kTRUE);
164   // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool
165   TObjArray* tracksClone = CloneAndReduceTrackList(tracks);
166   delete tracks;
167
168   const AliVVertex* vertex = InputEvent()->GetPrimaryVertex();
169   Double_t zVtx = vertex->GetZ();
170
171   fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kSameNS, tracksClone, tracksClone, tracksClone, tracksClone, 1.0, kFALSE, kFALSE);//same event for near and away side
172
173   fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::k1plus1, tracksClone, tracksClone, tracksClone, tracksClone, 1.0, kTRUE, kFALSE);//get number of possible away side triggers in the trigger area and outside of it
174
175   fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kBackgroundSameNS, tracksClone, tracksClone, tracksClone, tracksClone, 1.0, kFALSE, kTRUE);//background estimation for the same event
176
177   ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(1);
178   
179   // event mixing
180
181   // 1. First get an event pool corresponding in mult (cent) and
182   //    zvertex to the current event. Once initialized, the pool
183   //    should contain nMix (reduced) events. This routine does not
184   //    pre-scan the chain. The first several events of every chain
185   //    will be skipped until the needed pools are filled to the
186   //    specified depth. If the pool categories are not too rare, this
187   //    should not be a problem. If they are rare, you could lose
188   //    statistics.
189   
190   // 2. Collect the whole pool's content of tracks into one TObjArray
191   //    (bgTracks), which is effectively a single background super-event.
192   
193   // 3. The reduced and bgTracks arrays must both be passed into
194   //    FillCorrelations(). Also nMix should be passed in, so a weight
195   //    of 1./nMix can be applied.
196   
197   AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx);
198
199   if (!pool)
200     AliFatal(Form("No pool found for centrality = %f, zVtx = %f", centrality, zVtx));
201   if (pool->IsReady()){    
202     ((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());
203     Int_t nMix = pool->GetCurrentNEvents();
204
205     ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(2);
206     
207     // Fill mixed-event histos here  
208     for (Int_t jMix=0; jMix<nMix; jMix++){
209       TObjArray* bgTracks = pool->GetEvent(jMix);
210       
211       //standard mixed event
212       if(!fThreeParticleMixed){
213         fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kMixedNS, tracksClone, tracksClone, bgTracks, bgTracks, 1.0 / (2*nMix), kFALSE, kFALSE);
214         fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kMixedNS, bgTracks, bgTracks, tracksClone, tracksClone, 1.0 / (2*nMix), kFALSE, kFALSE);
215       }
216
217       //1plus1 mixed event
218       fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kMixed1plus1, tracksClone, tracksClone, bgTracks, bgTracks, 1.0 / (2*nMix), kTRUE, kFALSE);
219       fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kMixed1plus1, bgTracks, bgTracks, tracksClone, tracksClone, 1.0 / (2*nMix), kTRUE, kFALSE);
220
221       //mixed combinatorics
222       fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kMixedCombNS, tracksClone, bgTracks, tracksClone, bgTracks, 1.0 / (2*nMix), kFALSE, kFALSE);
223       fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kMixedCombNS, bgTracks, tracksClone, bgTracks, tracksClone, 1.0 / (2*nMix), kFALSE, kFALSE);
224     }
225     
226     
227     // use 3 particle mixed event
228     if(fThreeParticleMixed && nMix>1){
229       TObjArray* tracks_t2 = pool->GetEvent(0);
230       for (Int_t jMix=1; jMix<nMix; jMix++){
231
232         TObjArray* bgTracks = pool->GetEvent(jMix);
233
234         fHistos->FillCorrelations(centrality, zVtx, AliTwoPlusOneContainer::kMixedNS, tracksClone, tracks_t2, bgTracks, bgTracks, 1.0 / (nMix-1), kFALSE, kFALSE);
235       }
236     }
237     
238     
239   }
240
241   // ownership is with the pool now
242   pool->UpdatePool(tracksClone);
243 }
244
245 //________________________________________________________________________
246 void AliAnalysisTaskTwoPlusOne::Terminate(Option_t*)
247 {
248     //terminate function is called at the end
249     //can be used to draw histograms etc.
250     
251     
252 }
253
254
255 TObjArray* AliAnalysisTaskTwoPlusOne::CloneAndReduceTrackList(TObjArray* tracks)
256 {
257   // clones a track list by using AliCFParticle which uses much less memory (used for event mixing)
258   
259   TObjArray* tracksClone = new TObjArray;
260   tracksClone->SetOwner(kTRUE);
261   
262   for (Int_t i=0; i<tracks->GetEntriesFast(); i++)
263   {
264     AliVParticle* particle = (AliVParticle*) tracks->UncheckedAt(i);
265     AliCFParticle* copy = new AliCFParticle(particle->Pt(), particle->Eta(), particle->Phi(), particle->Charge(), 0);
266     copy->SetUniqueID(particle->GetUniqueID());
267     tracksClone->Add(copy);
268   }
269   
270   return tracksClone;
271 }
272
273
274 //____________________________________________________________________
275 void  AliAnalysisTaskTwoPlusOne::AddSettingsTree()
276 {
277   //Write settings to output list
278   TTree *settingsTree   = new TTree("UEAnalysisSettings","Analysis Settings in UE estimation");
279   settingsTree->Branch("fnTracksVertex", &fnTracksVertex,"nTracksVertex/I");
280   settingsTree->Branch("fZVertex", &fZVertex,"ZVertex/D");
281   //settingsTree->Branch("fCentralityMethod", fCentralityMethod.Data(),"CentralityMethod/C");
282   settingsTree->Branch("fTrackEtaCut", &fTrackEtaCut, "TrackEtaCut/D");
283   settingsTree->Branch("fTrackEtaCutMin", &fTrackEtaCutMin, "TrackEtaCutMin/D");
284   settingsTree->Branch("fPtMin", &fPtMin, "PtMin/D");
285   settingsTree->Branch("fFilterBit", &fFilterBit,"FilterBit/I");
286   settingsTree->Branch("fSharedClusterCut", &fSharedClusterCut,"SharedClusterCut/D");
287   settingsTree->Branch("fCrossedRowsCut", &fCrossedRowsCut,"CrossedRowsCut/I");
288   settingsTree->Branch("fFoundFractionCut", &fFoundFractionCut,"FoundFractionCut/D");
289   settingsTree->Branch("fTrackStatus", &fTrackStatus,"TrackStatus/I");
290   settingsTree->Branch("fThreeParticleMixed", &fThreeParticleMixed,"fThreeParticleMixed/I");
291   settingsTree->Branch("fMixingTracks", &fMixingTracks,"MixingTracks/I");
292   
293   settingsTree->Fill();
294   fListOfHistos->Add(settingsTree);
295 }