]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/PiKaPr/TestAOD/AliSpectraBothEventCuts.cxx
AliAODEvent::GetHeader() returns AliVHeader
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / PiKaPr / TestAOD / AliSpectraBothEventCuts.cxx
CommitLineData
239a080a 1
2/**************************************************************************
3 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
17//-----------------------------------------------------------------
18// AliSpectraBothEventCuts class
19//-----------------------------------------------------------------
20
21#include "TChain.h"
22#include "TTree.h"
23#include "TLegend.h"
24#include "TH1F.h"
25#include "TH1I.h"
26#include "TH2F.h"
27#include "TCanvas.h"
28#include "AliAnalysisTask.h"
29#include "AliAnalysisManager.h"
30#include "AliAODTrack.h"
31#include "AliAODMCParticle.h"
32#include "AliAODEvent.h"
33#include "AliESDEvent.h"
34#include "AliAODInputHandler.h"
35#include "AliAnalysisTaskESDfilter.h"
36#include "AliAnalysisDataContainer.h"
37#include "AliSpectraBothEventCuts.h"
38#include "AliSpectraBothTrackCuts.h"
49ed9c10 39#include "AliAnalysisUtils.h"
40#include "AliPPVsMultUtils.h"
239a080a 41//#include "AliSpectraBothHistoManager.h"
42#include <iostream>
43
44using namespace std;
45
46ClassImp(AliSpectraBothEventCuts)
47
b3ea73e1 48AliSpectraBothEventCuts::AliSpectraBothEventCuts(const char *name) : TNamed(name, "AOD Event Cuts"), fAOD(0),fAODEvent(AliSpectraBothTrackCuts::kAODobject), fTrackBits(0),fIsMC(0),fCentEstimator(""), fUseCentPatchAOD049(0), fUseSDDPatchforLHC11a(kDoNotCheckforSDD),fTriggerSettings(AliVEvent::kMB),fTrackCuts(0),
239a080a 49fIsSelected(0), fCentralityCutMin(0), fCentralityCutMax(0), fQVectorCutMin(0), fQVectorCutMax(0), fVertexCutMin(0), fVertexCutMax(0), fMultiplicityCutMin(0), fMultiplicityCutMax(0),fMaxChi2perNDFforVertex(0),
49ed9c10 50fMinRun(0),fMaxRun(0),fetarangeofmultiplicitycut(0.0),fUseAliPPVsMultUtils(false),
239a080a 51fHistoCuts(0),fHistoVtxBefSel(0),fHistoVtxAftSel(0),fHistoEtaBefSel(0),fHistoEtaAftSel(0),fHistoNChAftSel(0),fHistoQVector(0)
0ab8c127 52,fHistoEP(0),fHistoVtxAftSelwithoutZvertexCut(0),fHistoVtxalltriggerEventswithMCz(0),fHistoVtxAftSelwithoutZvertexCutusingMCz(0),fHistoRunNumbers(0),
49ed9c10 53fHistoCentrality(0),fHistoMultiplicty(0),fAnalysisUtils(0),fAliPPVsMultUtils(0)
0ab8c127 54
239a080a 55{
8fda510f 56 // Constructori
2d98dd91 57 // Bool_t oldStatus = TH1::AddDirectoryStatus();
58 //TH1::AddDirectory(kFALSE);
59 // fHistoCuts = new TH1I("fEventCuts", "Event Cuts", kNVtxCuts, -0.5, kNVtxCuts - 0.5);
60 // fHistoVtxBefSel = new TH1F("fHistoVtxBefSel", "Vtx distr before event selection",300,-15,15);
61// fHistoVtxAftSel = new TH1F("fHistoVtxAftSel", "Vtx distr after event selection",300,-15,15);
62// fHistoVtxAftSelwithoutZvertexCut=new TH1F("fHistoVtxAftSelwithoutZvertexcut", "Vtx distr after event selection without Z vertex cut",300,-15,15);
63// fHistoVtxalltriggerEventswithMCz=new TH1F("fHistoVtxalltriggerEventswithMCz", "generated z vertex position",300,-15,15);
64// fHistoVtxAftSelwithoutZvertexCutusingMCz=new TH1F("fHistoVtxAftSelwithoutZvertexCutusingMCz", "Vtx distr after event selection without Z vertex cut using MC z",300,-15,15);
65// fHistoEtaBefSel = new TH1F("fHistoEtaBefSel", "Eta distr before event selection",500,-2,2);
66// fHistoEtaAftSel = new TH1F("fHistoEtaAftSel", "Eta distr after event selection",500,-2,2);
67 // fHistoNChAftSel = new TH1F("fHistoNChAftSel", "NCh distr after event selection",2000,-0.5,1999.5);
239a080a 68 //fHistoQVectorPos = new TH1F("fHistoQVectorPos", "QVectorPos distribution",100,0,10);
69 //fHistoQVectorNeg = new TH1F("fHistoQVectorNeg", "QVectorNeg distribution",100,0,10);
2d98dd91 70 // fHistoQVector = new TH1F("fHistoQVector", "QVector with VZERO distribution",100,0,10);
71 // fHistoEP = new TH1F("fHistoEP", "EP with VZERO distribution",100,-10,10);
239a080a 72 fCentralityCutMin = 0.0; // default value of centrality cut minimum, 0 ~ no cut
73 fCentralityCutMax = 10000.0; // default value of centrality cut maximum, ~ no cut
74 // fQVectorPosCutMin=0.0;
75 // fQVectorPosCutMax=10000.0;
76 // fQVectorNegCutMin=0.0;
77 // fQVectorNegCutMax=10000.0;
78 fQVectorCutMin=0.0;
79 fQVectorCutMax=10000.0;
80 fVertexCutMin=-10.0;
81 fVertexCutMax=10.0;
82 fMultiplicityCutMin=-1.0;
83 fMultiplicityCutMax=-1.0;
84 fTrackBits=1;
b3ea73e1 85 fCentEstimator="V0M";
239a080a 86 fMaxChi2perNDFforVertex=-1;
2d98dd91 87 // TH1::AddDirectory(oldStatus);
8fda510f 88}
89//______________________________________________________
90
91AliSpectraBothEventCuts::~AliSpectraBothEventCuts()
92{
93 if(fHistoCuts)
94 delete fHistoCuts;
95 if(fHistoVtxBefSel)
96 delete fHistoVtxBefSel;
97 if(fHistoVtxAftSel)
98 delete fHistoVtxAftSel;
99 if(fHistoVtxAftSelwithoutZvertexCut)
100 delete fHistoVtxAftSelwithoutZvertexCut;
101 if(fHistoVtxalltriggerEventswithMCz)
102 delete fHistoVtxalltriggerEventswithMCz;
103 if(fHistoVtxAftSelwithoutZvertexCutusingMCz)
104 delete fHistoVtxAftSelwithoutZvertexCutusingMCz;
105 if(fHistoEtaBefSel)
106 delete fHistoEtaBefSel;
107 if(fHistoEtaAftSel)
108 delete fHistoEtaAftSel ;
109 if(fHistoNChAftSel)
110 delete fHistoNChAftSel;
111 if(fHistoQVector)
112 delete fHistoQVector;
113 if(fHistoEP)
114 delete fHistoEP;
2d98dd91 115 if(fHistoRunNumbers)
116 delete fHistoRunNumbers;
0ab8c127 117 if(fHistoCentrality)
118 delete fHistoCentrality;
119 if(fHistoMultiplicty)
120 delete fHistoMultiplicty;
4954491e 121 if(fAnalysisUtils)
122 delete fAnalysisUtils;
49ed9c10 123 if(fAliPPVsMultUtils)
124 delete fAliPPVsMultUtils;
0ab8c127 125
239a080a 126}
2d98dd91 127//______________________________________________________
128void AliSpectraBothEventCuts::InitHisto()
129{
130 Bool_t oldStatus = TH1::AddDirectoryStatus();
131 TH1::AddDirectory(kFALSE);
132 if(!fHistoCuts)
133 fHistoCuts = new TH1I("fEventCuts", "Event Cuts", kNVtxCuts, -0.5, kNVtxCuts - 0.5);
134 if(!fHistoVtxBefSel )
135 fHistoVtxBefSel = new TH1F("fHistoVtxBefSel", "Vtx distr before event selection",300,-15,15);
136 if(!fHistoVtxAftSel)
137 fHistoVtxAftSel = new TH1F("fHistoVtxAftSel", "Vtx distr after event selection",300,-15,15);
138 if(!fHistoVtxAftSelwithoutZvertexCut)
139 fHistoVtxAftSelwithoutZvertexCut=new TH1F("fHistoVtxAftSelwithoutZvertexcut", "Vtx distr after event selection without Z vertex cut",300,-15,15);
140 if(!fHistoVtxalltriggerEventswithMCz)
141 fHistoVtxalltriggerEventswithMCz=new TH1F("fHistoVtxalltriggerEventswithMCz", "generated z vertex position",300,-15,15);
142 if(!fHistoVtxAftSelwithoutZvertexCutusingMCz)
143 fHistoVtxAftSelwithoutZvertexCutusingMCz=new TH1F("fHistoVtxAftSelwithoutZvertexCutusingMCz", "Vtx distr after event selection without Z vertex cut using MC z",300,-15,15);
144 if(!fHistoEtaBefSel)
145 fHistoEtaBefSel = new TH1F("fHistoEtaBefSel", "Eta distr before event selection",500,-2,2);
146 if(!fHistoEtaAftSel)
147 fHistoEtaAftSel = new TH1F("fHistoEtaAftSel", "Eta distr after event selection",500,-2,2);
148 if(!fHistoNChAftSel)
149 fHistoNChAftSel = new TH1F("fHistoNChAftSel", "NCh distr after event selection",2000,-0.5,1999.5);
150 //fHistoQVectorPos = new TH1F("fHistoQVectorPos", "QVectorPos distribution",100,0,10);
151 //fHistoQVectorNeg = new TH1F("fHistoQVectorNeg", "QVectorNeg distribution",100,0,10);
152 if(!fHistoQVector)
153 fHistoQVector = new TH1F("fHistoQVector", "QVector with VZERO distribution",100,0,10);
154 if(!fHistoEP)
155 fHistoEP = new TH1F("fHistoEP", "EP with VZERO distribution",100,-10,10);
156 if(!fHistoRunNumbers)
157 {
158 if(fMaxRun>fMinRun&&fMinRun>=0)
159 fHistoRunNumbers=new TH1F("fHistoRunNumbers","Run numbers",fMaxRun-fMinRun+1,fMinRun-0.5,fMaxRun+0.5);
160 else
161 fHistoRunNumbers=new TH1F("fHistoRunNumbers","Run numbers",1001,120000-.5,121000+0.5);
239a080a 162
2d98dd91 163 }
0ab8c127 164 if(!fHistoCentrality)
165 fHistoCentrality = new TH2F("fHistoCentrality", "centrality",2,0,2,100,0.0,100);
166
167 if(!fHistoMultiplicty)
49ed9c10 168 fHistoMultiplicty= new TH2F("fHistoMultiplicty", "multiplicty estimator",2,0,2,155,-4.5,150.5);
0ab8c127 169
2d98dd91 170 TH1::AddDirectory(oldStatus);
171}
239a080a 172//______________________________________________________
34347952 173Bool_t AliSpectraBothEventCuts::IsSelected(AliVEvent * aod,AliSpectraBothTrackCuts* trackcuts,Bool_t isMC,Double_t mcZ,TH1F* managerhisteventcuts)
239a080a 174{
175 // Returns true if Event Cuts are selected and applied
176 fAOD = aod;
177 fTrackCuts = trackcuts;
178 fHistoCuts->Fill(kProcessedEvents);
34347952 179 if(managerhisteventcuts)
180 managerhisteventcuts->Fill(0);
2d98dd91 181 fHistoRunNumbers->Fill(aod->GetRunNumber());
1c23d41e 182 Bool_t IsPhysSel = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fTriggerSettings);//FIXME we can add the trigger mask here
d0592deb 183 if(!IsPhysSel)
184 return IsPhysSel;
185 if(fAnalysisUtils) // we check for pile-up
4954491e 186 IsPhysSel = (!fAnalysisUtils->IsPileUpEvent(fAOD));
187 if(!IsPhysSel)
188 return IsPhysSel;
189
8bb435b0 190 if(isMC)
191 fHistoVtxalltriggerEventswithMCz->Fill(mcZ);
239a080a 192 //loop on tracks, before event selection, filling QA histos
4954491e 193 AliESDEvent* esdevent=0x0;
239a080a 194 AliAODEvent* aodevent=0x0;
195 Bool_t isSDD=kFALSE;
4954491e 196 TString nameoftrack(fAOD->ClassName());
239a080a 197 if(!nameoftrack.CompareTo("AliESDEvent"))
198 {
199 fAODEvent=AliSpectraBothTrackCuts::kESDobject;
200
201 if(fUseSDDPatchforLHC11a!=kDoNotCheckforSDD)
202 {
203 esdevent=dynamic_cast<AliESDEvent*>(fAOD);
737ad8ac 204 if(!esdevent)
94a8ed71 205 return kFALSE;
239a080a 206 if(esdevent->GetFiredTriggerClasses().Contains("ALLNOTRD"))
207 isSDD=kTRUE;
208 }
209 }
210 else if(!nameoftrack.CompareTo("AliAODEvent"))
211 {
212 fAODEvent=AliSpectraBothTrackCuts::kAODobject;
213 if(fUseSDDPatchforLHC11a!=kDoNotCheckforSDD)
214 {
215 aodevent=dynamic_cast<AliAODEvent*>(fAOD);
737ad8ac 216 if(!aodevent)
94a8ed71 217 return kFALSE;
239a080a 218 if(aodevent->GetFiredTriggerClasses().Contains("ALLNOTRD"))
219 isSDD=kTRUE;
220 }
221 }
222 else
223 return false;
224 if(fUseSDDPatchforLHC11a==kwithSDD&&isSDD==kFALSE)
225 return false;
226 if(fUseSDDPatchforLHC11a==kwithoutSDD&&isSDD==kTRUE)
227 return false;
228
239a080a 229 fHistoCuts->Fill(kPhysSelEvents);
34347952 230 if(managerhisteventcuts)
231 managerhisteventcuts->Fill(1);
239a080a 232
233
234
235
bca0c28e 236 const AliVVertex * vertex = fAOD->GetPrimaryVertex();//FIXME vertex is recreated
237
4954491e 238 if(vertex)
239 fHistoVtxBefSel->Fill(vertex->GetZ());
239a080a 240 fIsSelected =kFALSE;
49ed9c10 241 if(CheckVtx())
242 { //selection on vertex
547e52a0 243
4954491e 244 fIsSelected=kTRUE;
239a080a 245 }
8bb435b0 246 if(fIsSelected&&vertex)
49ed9c10 247 {
248 fHistoVtxAftSelwithoutZvertexCut->Fill(vertex->GetZ());
249 if(isMC)
250 fHistoVtxAftSelwithoutZvertexCutusingMCz->Fill(mcZ);
251 if (vertex->GetZ() > fVertexCutMin && vertex->GetZ() < fVertexCutMax)
252 {
547e52a0 253 fIsSelected=kTRUE;
254 fHistoVtxAftSel->Fill(vertex->GetZ());
49ed9c10 255 }
256 else
257 {
547e52a0 258 fIsSelected=kFALSE;
49ed9c10 259 }
239a080a 260 }
4954491e 261 if(fIsSelected)
49ed9c10 262 {
263 if( CheckCentralityCut() && CheckMultiplicityCut() && CheckQVectorCut())
264 fIsSelected=kTRUE;
77898649 265 else
266 fIsSelected=kFALSE;
49ed9c10 267 }
4954491e 268
269
270
239a080a 271 Int_t Nch=0;
49ed9c10 272 for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++)
273 {
274 AliVTrack* track =dynamic_cast<AliVTrack*>(fAOD->GetTrack(iTracks));
239a080a 275 /* if(fAODEvent==AliSpectraBothTrackCuts::kESDobject)
276 track=dynamic_cast<AliVTrack*>(esdevent->GetTrack(iTracks));
277 else if (fAODEvent==AliSpectraBothTrackCuts::kAODobject)
278 track=dynamic_cast<AliVTrack*>(aodevent->GetTrack(iTracks));
279 else return false;*/
280
49ed9c10 281 if (!fTrackCuts->IsSelected(track,kFALSE))
282 continue;
283 fHistoEtaBefSel->Fill(track->Eta());
284 if(fIsSelected)
285 {
286 fHistoEtaAftSel->Fill(track->Eta());
287 Nch++;
288 }
239a080a 289 }
4954491e 290 if(fIsSelected)
77898649 291 {
292 fHistoNChAftSel->Fill(Nch);
293 fHistoCuts->Fill(kAcceptedEvents);
34347952 294 if(managerhisteventcuts)
295 managerhisteventcuts->Fill(2);
296
77898649 297 }
239a080a 298 return fIsSelected;
299}
300
301//______________________________________________________
547e52a0 302Bool_t AliSpectraBothEventCuts::CheckVtx()
239a080a 303{
304 // reject events outside of range
305 const AliVVertex * vertex = fAOD->GetPrimaryVertex();
306 //when moving to 2011 wìone has to add a cut using SPD vertex.
307 //The point is that for events with |z|>20 the vertexer tracks is not working (only 2011!). One has to put a safety cut using SPD vertex large e.g. 15cm
308 if (!vertex)
309 {
310 fHistoCuts->Fill(kVtxNoEvent);
311 return kFALSE;
312 }
bca0c28e 313 if(vertex->GetNContributors()<1)
239a080a 314 {
bca0c28e 315
316 fHistoCuts->Fill(kZeroCont);
239a080a 317 return kFALSE;
318
319 }
bca0c28e 320
321 TString tmp(vertex->GetTitle());
322 if(tmp.Contains("NoConstraint"))
239a080a 323 {
bca0c28e 324 fHistoCuts->Fill(kTPCasPV);
325 return kFALSE;
326 }
327
328
547e52a0 329 // if (vertex->GetZ() > fVertexCutMin && vertex->GetZ() < fVertexCutMax)
330 // {
331 // return kTRUE;
332 // }
4954491e 333 if(!CheckVtxChi2perNDF())
334 return kFALSE;
335
336 fHistoCuts->Fill(kGoodVtx);
547e52a0 337 //return kFALSE;
338 return kTRUE;
239a080a 339}
340
341//______________________________________________________
342Bool_t AliSpectraBothEventCuts::CheckCentralityCut()
343{
344 // Check centrality cut
49ed9c10 345 if ( fCentralityCutMax<0.0 && fCentralityCutMin<0.0 )
346 return kTRUE;
347 Double_t cent=0;
348 if(fUseAliPPVsMultUtils)
349 {
350 if(!fAliPPVsMultUtils)
351 fAliPPVsMultUtils=new AliPPVsMultUtils();
352 cent=fAliPPVsMultUtils->GetMultiplicityPercentile(fAOD,fCentEstimator.Data());
353 }
354 else
355 {
356 if(!fUseCentPatchAOD049)
357 cent=fAOD->GetCentrality()->GetCentralityPercentile(fCentEstimator.Data());
358 else
359 cent=ApplyCentralityPatchAOD049();
360 }
0ab8c127 361 fHistoCentrality->Fill(0.5,cent);
77898649 362 if ( (cent < fCentralityCutMax) && (cent >= fCentralityCutMin) )
0ab8c127 363 {
364 fHistoCentrality->Fill(1.5,cent);
365 return kTRUE;
366 }
239a080a 367 fHistoCuts->Fill(kVtxCentral);
368
369 return kFALSE;
370}
371
372//______________________________________________________
373Bool_t AliSpectraBothEventCuts::CheckMultiplicityCut()
374{
375 // Check multiplicity cut
0ab8c127 376 if(fMultiplicityCutMin<0 && fMultiplicityCutMax<0)
377 return kTRUE;
378 Int_t Ncharged=-1;
379 if(fAODEvent==AliSpectraBothTrackCuts::kESDobject)
380 {
381 AliESDEvent* esdevent=dynamic_cast<AliESDEvent*>(fAOD);
382 AliESDtrackCuts::MultEstTrackType estType = esdevent->GetPrimaryVertexTracks()->GetStatus() ? AliESDtrackCuts::kTrackletsITSTPC : AliESDtrackCuts::kTracklets;
383 Ncharged=AliESDtrackCuts::GetReferenceMultiplicity(esdevent,estType,fetarangeofmultiplicitycut);
384 }
385 else if(fAODEvent==AliSpectraBothTrackCuts::kAODobject)
386 {
387 AliAODEvent* aodevent=0x0;
388 aodevent=dynamic_cast<AliAODEvent*>(fAOD);
0a918d8d 389 AliAODHeader * header = dynamic_cast<AliAODHeader*>(aodevent->GetHeader());
390 if(!header) AliFatal("Not a standard AOD");
391
0ab8c127 392 if(TMath::Abs(0.8-fetarangeofmultiplicitycut)<0.1)
0a918d8d 393 Ncharged=header->GetRefMultiplicityComb08();
0ab8c127 394 else if (TMath::Abs(0.5-fetarangeofmultiplicitycut)<0.1)
0a918d8d 395 Ncharged=header->GetRefMultiplicityComb05();
0ab8c127 396 else
397 Ncharged=-1;
398 }
399 else
400 return kFALSE;
239a080a 401
0ab8c127 402 fHistoMultiplicty->Fill(0.5,Ncharged);
403 if(Ncharged>=fMultiplicityCutMin && Ncharged<fMultiplicityCutMax&& Ncharged>0)
404 {
405 fHistoMultiplicty->Fill(1.5,Ncharged);
406 return kTRUE;
407 }
408 fHistoCuts->Fill(kVtxCentral);
239a080a 409 return kFALSE;
410}
411
412//______________________________________________________
413Bool_t AliSpectraBothEventCuts::CheckQVectorCut()
414{
415 if(fQVectorCutMin<0.0 && fQVectorCutMax<0.0)return kTRUE;
416 // Check qvector cut
417 /// FIXME: Q vector
418 // //Selection on QVector, before ANY other selection on the event
419 // //Spectra MUST be normalized wrt events AFTER the selection on Qvector
420 // Double_t Qx2EtaPos = 0, Qy2EtaPos = 0;
421 // Double_t Qx2EtaNeg = 0, Qy2EtaNeg = 0;
422
423 // Int_t multPos = 0;
424 // Int_t multNeg = 0;
425 // for(Int_t iT = 0; iT < fAOD->GetNumberOfTracks(); iT++) {
426 // AliAODTrack* aodTrack = fAOD->GetTrack(iT);
427 // if (!aodTrack->TestFilterBit(fTrackBits)) continue;
428 // if (aodTrack->Eta() >= 0){
429 // multPos++;
430 // Qx2EtaPos += TMath::Cos(2*aodTrack->Phi());
431 // Qy2EtaPos += TMath::Sin(2*aodTrack->Phi());
432 // } else {
433 // multNeg++;
434 // Qx2EtaNeg += TMath::Cos(2*aodTrack->Phi());
435 // Qy2EtaNeg += TMath::Sin(2*aodTrack->Phi());
436 // }
437 // }
438 // Double_t qPos=-999;
439 // if(multPos!=0)qPos= TMath::Sqrt((Qx2EtaPos*Qx2EtaPos + Qy2EtaPos*Qy2EtaPos)/multPos);
440 // Double_t qNeg=-999;
441 // if(multNeg!=0)qNeg= TMath::Sqrt((Qx2EtaNeg*Qx2EtaNeg + Qy2EtaNeg*Qy2EtaNeg)/multNeg);
442 //if(qPos<fQVectorPosCutMin || qPos>fQVectorPosCutMax || qNeg<fQVectorNegCutMin || qNeg>fQVectorNegCutMax)return kFALSE;
443
444 Double_t qxEPVZERO = 0, qyEPVZERO = 0;
445 Double_t qVZERO = -999;
446 Double_t psi=fAOD->GetEventplane()->CalculateVZEROEventPlane(fAOD,10,2,qxEPVZERO,qyEPVZERO);
447
448 qVZERO= TMath::Sqrt(qxEPVZERO*qxEPVZERO + qyEPVZERO*qyEPVZERO);
449 if(qVZERO<fQVectorCutMin || qVZERO>fQVectorCutMax)return kFALSE;
450 fHistoQVector->Fill(qVZERO);
451 fHistoEP->Fill(psi);
452
453 fHistoCuts->Fill(kQVector);
454 // fHistoQVectorPos->Fill(qPos);
455 // fHistoQVectorNeg->Fill(qNeg);
456 return kTRUE;
457}
458//____________________________________________________________
459 Bool_t AliSpectraBothEventCuts::CheckVtxChi2perNDF()
460 {
461 if(fMaxChi2perNDFforVertex<0)
462 return kTRUE;
463 const AliVVertex * vertex = fAOD->GetPrimaryVertex();
464 if(TMath::Abs(vertex->GetChi2perNDF())>fMaxChi2perNDFforVertex)
465 return kFALSE;
466 return kTRUE;
467 }
468
469
470
471//______________________________________________________
472void AliSpectraBothEventCuts::PrintCuts()
473{
474 // print info about event cuts
475 cout << "Event Stats" << endl;
2d98dd91 476 if(fHistoCuts)
477 {
478 cout << " > Number of accepted events: " << fHistoCuts->GetBinContent(kAcceptedEvents + 1) << endl;
479 cout << " > Number of processed events: " << fHistoCuts->GetBinContent(kProcessedEvents + 1) << endl;
480 cout << " > Number of PhysSel events: " << fHistoCuts->GetBinContent(kPhysSelEvents + 1) << endl;
4954491e 481 cout << " > With good veretx: " << fHistoCuts->GetBinContent(kGoodVtx + 1) << endl;
2d98dd91 482 cout << " > Events cut by centrality: " << fHistoCuts->GetBinContent(kVtxCentral + 1) << endl;
483 cout << " > Events without vertex: " << fHistoCuts->GetBinContent(kVtxNoEvent + 1) << endl;
484 cout << " > QVector cut: " << fHistoCuts->GetBinContent(kQVector + 1) << endl;
485 }
239a080a 486 cout << " > Track type used for the QVector calculation: " << fTrackBits << endl;
487 // cout << " > QPosRange: [" << fQVectorPosCutMin <<"," <<fQVectorPosCutMax<<"]"<< endl;
488 // cout << " > QNegRange: [" << fQVectorNegCutMin <<"," <<fQVectorNegCutMax<<"]"<< endl;
489 cout << " > QRange: [" << fQVectorCutMin <<"," <<fQVectorCutMax<<"]"<< endl;
490 cout << " > Vertex: [" << fVertexCutMin <<"," <<fVertexCutMax<<"]"<< endl;
491 cout << " > Multiplicity: [" << fMultiplicityCutMin <<"," <<fMultiplicityCutMax<<"]"<< endl;
492 cout << " > Centrality: [" << fCentralityCutMin <<"," <<fCentralityCutMax<<"]"<< endl;
493}
494//______________________________________________________
495
496Double_t AliSpectraBothEventCuts::ApplyCentralityPatchAOD049()
497{
498 //
499 //Apply centrality patch for AOD049 outliers
500 //
501 // if (fCentralityType!="V0M") {
502 // AliWarning("Requested patch forAOD049 for wrong value (not centrality from V0).");
503 // return -999.0;
504 // }
505 AliCentrality *centrality = fAOD->GetCentrality();
506 if (!centrality) {
507 AliWarning("Cannot get centrality from AOD event.");
508 return -999.0;
509 }
510
511 Float_t cent = (Float_t)(centrality->GetCentralityPercentile("V0M"));
512 /*
513 Bool_t isSelRun = kFALSE;
514 Int_t selRun[5] = {138364, 138826, 138828, 138836, 138871};
515 if(cent<0){
516 Int_t quality = centrality->GetQuality();
517 if(quality<=1){
518 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
519 } else {
520 Int_t runnum=aodEvent->GetRunNumber();
521 for(Int_t ir=0;ir<5;ir++){
522 if(runnum==selRun[ir]){
523 isSelRun=kTRUE;
524 break;
525 }
526 }
527 if((quality==8||quality==9)&&isSelRun) cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
528 }
529 }
530 */
531 if(cent>=0.0) {
532 Float_t v0 = 0.0;
533 AliAODEvent *aodEvent = (AliAODEvent *)fAOD;
534 AliAODVZERO *aodV0 = (AliAODVZERO *) aodEvent->GetVZEROData();
535 v0+=aodV0->GetMTotV0A();
536 v0+=aodV0->GetMTotV0C();
537 if ( (cent==0) && (v0<19500) ) {
538 AliDebug(3, Form("Filtering issue in centrality -> cent = %5.2f",cent));
539 return -999.0;
540 }
541 Float_t tkl = (Float_t)(aodEvent->GetTracklets()->GetNumberOfTracklets());
542 Float_t val = 1.30552 + 0.147931 * v0;
543
544 Float_t tklSigma[101] = {176.644, 156.401, 153.789, 153.015, 142.476, 137.951, 136.127, 129.852, 127.436, 124.86,
545 120.788, 115.611, 113.172, 110.496, 109.127, 104.421, 102.479, 99.9766, 97.5152, 94.0654,
546 92.4602, 89.3364, 87.1342, 83.3497, 82.6216, 81.1084, 78.0793, 76.1234, 72.9434, 72.1334,
547 68.0056, 68.2755, 66.0376, 62.9666, 62.4274, 59.65, 58.3776, 56.6361, 54.5184, 53.4224,
548 51.932, 50.8922, 48.2848, 47.912, 46.5717, 43.4114, 43.2083, 41.3065, 40.1863, 38.5255,
549 37.2851, 37.5396, 34.4949, 33.8366, 31.8043, 31.7412, 30.8392, 30.0274, 28.8793, 27.6398,
550 26.6488, 25.0183, 25.1489, 24.4185, 22.9107, 21.2002, 21.6977, 20.1242, 20.4963, 19.0235,
551 19.298, 17.4103, 16.868, 15.2939, 15.2939, 16.0295, 14.186, 14.186, 15.2173, 12.9504, 12.9504,
552 12.9504, 15.264, 12.3674, 12.3674, 12.3674, 12.3674, 12.3674, 18.3811, 13.7544, 13.7544,
553 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544
554 };
555
556 if ( TMath::Abs(tkl-val) > 6.*tklSigma[(Int_t)cent] ) {
557 AliDebug(3, Form("Outlier event in centrality -> cent = %5.2f",cent));
558 return -999.0;
559 }
560 } else {
561 //force it to be -999. whatever the negative value was
562 cent = -999.;
563 }
564 return cent;
565}
566
567//______________________________________________________
568
569
570Long64_t AliSpectraBothEventCuts::Merge(TCollection* list)
571{
572 // Merge a list of AliSpectraBothEventCuts objects with this.
573 // Returns the number of merged objects (including this).
574
575 // AliInfo("Merging");
576
577 if (!list)
578 return 0;
579
580 if (list->IsEmpty())
581 return 1;
582
583 TIterator* iter = list->MakeIterator();
584 TObject* obj;
585
586 // collections of all histograms
587 TList collections;//FIXME we should use only 1 collection
588 TList collections_histoVtxBefSel;
589 TList collections_histoVtxAftSel;
590 TList collections_histoEtaBefSel;
591 TList collections_histoEtaAftSel;
592 TList collections_histoNChAftSel;
593 // TList collections_histoQVectorPos;
594 // TList collections_histoQVectorNeg;
595 TList collections_histoQVector;
596 TList collections_histoEP;
3a71f081 597 TList collections_histoVtxAftSelwithoutZvertexCut;
598 TList collections_histoVtxalltriggerEventswithMCz;
599 TList collections_histoVtxAftSelwithoutZvertexCutusingMCz;
2d98dd91 600 TList collections_histoRunNumbers;
0ab8c127 601 TList collections_histoCentrality;
602 TList collections_histoMultiplicty;
603
239a080a 604 Int_t count = 0;
605
606 while ((obj = iter->Next())) {
607 AliSpectraBothEventCuts* entry = dynamic_cast<AliSpectraBothEventCuts*> (obj);
608 if (entry == 0)
609 continue;
610
611 TH1I * histo = entry->GetHistoCuts();
612 collections.Add(histo);
613 TH1F * histo_histoVtxBefSel = entry->GetHistoVtxBefSel();
614 collections_histoVtxBefSel.Add(histo_histoVtxBefSel);
615 TH1F * histo_histoVtxAftSel = entry->GetHistoVtxAftSel();
616 collections_histoVtxAftSel.Add(histo_histoVtxAftSel);
617 TH1F * histo_histoEtaBefSel = entry->GetHistoEtaBefSel();
618 collections_histoEtaBefSel.Add(histo_histoEtaBefSel);
619 TH1F * histo_histoEtaAftSel = entry->GetHistoEtaAftSel();
620 collections_histoEtaAftSel.Add(histo_histoEtaAftSel);
621 TH1F * histo_histoNChAftSel = entry->GetHistoNChAftSel();
622 collections_histoNChAftSel.Add(histo_histoNChAftSel);
623 // TH1F * histo_histoQVectorPos = entry->GetHistoQVectorPos();
624 // collections_histoQVectorPos.Add(histo_histoQVectorPos);
625 // TH1F * histo_histoQVectorNeg = entry->GetHistoQVectorNeg();
626 // collections_histoQVectorNeg.Add(histo_histoQVectorNeg);
627 TH1F * histo_histoQVector = entry->GetHistoQVector();
628 collections_histoQVector.Add(histo_histoQVector);
629 TH1F * histo_histoEP = entry->GetHistoEP();
630 collections_histoEP.Add(histo_histoEP);
3a71f081 631 TH1F* histo_histoVtxAftSelwithoutZvertexCut=entry->GetHistoVtxAftSelwithoutZvertexCut();
632 collections_histoVtxAftSelwithoutZvertexCut.Add(histo_histoVtxAftSelwithoutZvertexCut);
633 TH1F* histo_histoVtxalltriggerEventswithMCz=entry->GetHistoVtxGenerated();
634 collections_histoVtxalltriggerEventswithMCz.Add(histo_histoVtxalltriggerEventswithMCz);
2d98dd91 635
636 TH1F* histo_histoVtxAftSelwithoutZvertexCutusingMCz=entry->GetHistoVtxAftSelwithoutZvertexCutusingMCz();
3a71f081 637 collections_histoVtxAftSelwithoutZvertexCutusingMCz.Add(histo_histoVtxAftSelwithoutZvertexCutusingMCz);
2d98dd91 638
639 TH1F* histo_histoRunNumbers=entry->GetHistoRunNumbers();
640 collections_histoRunNumbers.Add(histo_histoRunNumbers);
0ab8c127 641
642 TH2F* histo_histoCentrality=entry->GetHistoCentrality();
643 collections_histoCentrality.Add(histo_histoCentrality);
644
645TH2F* histo_histoMultiplicty=entry->GetHistoMultiplicty();
646 collections_histoMultiplicty.Add(histo_histoMultiplicty);
647
648
239a080a 649 count++;
650 }
651
652 fHistoCuts->Merge(&collections);
653 fHistoVtxBefSel->Merge(&collections_histoVtxBefSel);
654 fHistoVtxAftSel->Merge(&collections_histoVtxAftSel);
655 fHistoEtaBefSel->Merge(&collections_histoEtaBefSel);
656 fHistoEtaAftSel->Merge(&collections_histoEtaAftSel);
657 fHistoNChAftSel->Merge(&collections_histoNChAftSel);
658 // fHistoQVectorPos->Merge(&collections_histoQVectorPos);
659 // fHistoQVectorNeg->Merge(&collections_histoQVectorNeg);
660 fHistoQVector->Merge(&collections_histoQVector);
661 fHistoEP->Merge(&collections_histoEP);
3a71f081 662
663 fHistoVtxAftSelwithoutZvertexCut->Merge(&collections_histoVtxAftSelwithoutZvertexCut);
664 fHistoVtxalltriggerEventswithMCz->Merge(&collections_histoVtxalltriggerEventswithMCz);
665 fHistoVtxAftSelwithoutZvertexCutusingMCz->Merge(&collections_histoVtxAftSelwithoutZvertexCutusingMCz);
2d98dd91 666 fHistoRunNumbers->Merge(&collections_histoRunNumbers);
0ab8c127 667 fHistoCentrality->Merge(&collections_histoCentrality);
668 fHistoMultiplicty->Merge(&collections_histoMultiplicty);
669
239a080a 670
671 delete iter;
672
673 return count+1;
674}
2d98dd91 675//__________________________________________________________________________________________________________
676void AliSpectraBothEventCuts::SetRunNumberRange(Int_t min, Int_t max)
677{
678 if(max>min&&min>=0)
679 {
680 fMinRun=min;
681 fMaxRun=max;
682 }
683}