1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 // An event cut class for the flow framework
21 // origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
31 #include "AliVVertex.h"
32 #include "AliVEvent.h"
33 #include "AliESDEvent.h"
34 #include "AliAODEvent.h"
35 #include "AliAODHeader.h"
36 #include "AliCentrality.h"
37 #include "AliESDVZERO.h"
38 #include "AliMultiplicity.h"
39 #include "AliMCEvent.h"
40 #include "AliFlowEventCuts.h"
41 #include "AliFlowTrackCuts.h"
42 #include "AliTriggerAnalysis.h"
44 ClassImp(AliFlowEventCuts)
46 //-----------------------------------------------------------------------
47 AliFlowEventCuts::AliFlowEventCuts():
50 fCutNumberOfTracks(kFALSE),
51 fNumberOfTracksMax(INT_MAX),
52 fNumberOfTracksMin(INT_MIN),
54 fRefMultMethod(kTPConly),
55 fUseAliESDtrackCutsRefMult(kFALSE),
56 fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
61 fStandardTPCcuts(NULL),
62 fStandardGlobalCuts(NULL),
63 fCutPrimaryVertexX(kFALSE),
64 fPrimaryVertexXmax(INT_MAX),
65 fPrimaryVertexXmin(INT_MIN),
66 fCutPrimaryVertexY(kFALSE),
67 fPrimaryVertexYmax(INT_MAX),
68 fPrimaryVertexYmin(INT_MIN),
69 fCutPrimaryVertexZ(kFALSE),
70 fPrimaryVertexZmax(INT_MAX),
71 fPrimaryVertexZmin(INT_MIN),
72 fCutNContributors(kFALSE),
73 fNContributorsMax(INT_MAX),
74 fNContributorsMin(INT_MIN),
78 fCutSPDvertexerAnomaly(kFALSE),
79 fCutTPCmultiplicityOutliers(kFALSE),
80 fCutCentralityPercentile(kFALSE),
81 fUseCentralityUnchecked(kFALSE),
82 fCentralityPercentileMethod(kTPConly),
83 fCentralityPercentileMax(100.),
84 fCentralityPercentileMin(0.),
85 fCutZDCtiming(kFALSE),
91 //-----------------------------------------------------------------------
92 AliFlowEventCuts::AliFlowEventCuts(const char* name, const char* title):
95 fCutNumberOfTracks(kFALSE),
96 fNumberOfTracksMax(INT_MAX),
97 fNumberOfTracksMin(INT_MIN),
99 fRefMultMethod(kTPConly),
100 fUseAliESDtrackCutsRefMult(kFALSE),
101 fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
102 fRefMultMax(INT_MAX),
103 fRefMultMin(INT_MIN),
106 fStandardTPCcuts(AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010()),
107 fStandardGlobalCuts(AliFlowTrackCuts::GetStandardGlobalTrackCuts2010()),
108 fCutPrimaryVertexX(kFALSE),
109 fPrimaryVertexXmax(INT_MAX),
110 fPrimaryVertexXmin(INT_MIN),
111 fCutPrimaryVertexY(kFALSE),
112 fPrimaryVertexYmax(INT_MAX),
113 fPrimaryVertexYmin(INT_MIN),
114 fCutPrimaryVertexZ(kFALSE),
115 fPrimaryVertexZmax(INT_MAX),
116 fPrimaryVertexZmin(INT_MIN),
117 fCutNContributors(kFALSE),
118 fNContributorsMax(INT_MAX),
119 fNContributorsMin(INT_MIN),
121 fMeanPtMax(-DBL_MAX),
123 fCutSPDvertexerAnomaly(kFALSE),
124 fCutTPCmultiplicityOutliers(kFALSE),
125 fCutCentralityPercentile(kFALSE),
126 fUseCentralityUnchecked(kFALSE),
127 fCentralityPercentileMethod(kTPConly),
128 fCentralityPercentileMax(100.),
129 fCentralityPercentileMin(0.),
130 fCutZDCtiming(kFALSE),
136 ////-----------------------------------------------------------------------
137 AliFlowEventCuts::AliFlowEventCuts(const AliFlowEventCuts& that):
140 fCutNumberOfTracks(that.fCutNumberOfTracks),
141 fNumberOfTracksMax(that.fNumberOfTracksMax),
142 fNumberOfTracksMin(that.fNumberOfTracksMin),
143 fCutRefMult(that.fCutRefMult),
144 fRefMultMethod(that.fRefMultMethod),
145 fUseAliESDtrackCutsRefMult(that.fUseAliESDtrackCutsRefMult),
146 fRefMultMethodAliESDtrackCuts(that.fRefMultMethodAliESDtrackCuts),
147 fRefMultMax(that.fRefMultMax),
148 fRefMultMin(that.fRefMultMin),
151 fStandardTPCcuts(NULL),
152 fStandardGlobalCuts(NULL),
153 fCutPrimaryVertexX(that.fCutPrimaryVertexX),
154 fPrimaryVertexXmax(that.fPrimaryVertexXmax),
155 fPrimaryVertexXmin(that.fPrimaryVertexXmin),
156 fCutPrimaryVertexY(that.fCutPrimaryVertexX),
157 fPrimaryVertexYmax(that.fPrimaryVertexYmax),
158 fPrimaryVertexYmin(that.fPrimaryVertexYmin),
159 fCutPrimaryVertexZ(that.fCutPrimaryVertexX),
160 fPrimaryVertexZmax(that.fPrimaryVertexZmax),
161 fPrimaryVertexZmin(that.fPrimaryVertexZmin),
162 fCutNContributors(that.fCutNContributors),
163 fNContributorsMax(that.fNContributorsMax),
164 fNContributorsMin(that.fNContributorsMin),
165 fCutMeanPt(that.fCutMeanPt),
166 fMeanPtMax(that.fMeanPtMax),
167 fMeanPtMin(that.fMeanPtMin),
168 fCutSPDvertexerAnomaly(that.fCutSPDvertexerAnomaly),
169 fCutTPCmultiplicityOutliers(that.fCutTPCmultiplicityOutliers),
170 fCutCentralityPercentile(that.fCutCentralityPercentile),
171 fUseCentralityUnchecked(that.fUseCentralityUnchecked),
172 fCentralityPercentileMethod(that.fCentralityPercentileMethod),
173 fCentralityPercentileMax(that.fCentralityPercentileMax),
174 fCentralityPercentileMin(that.fCentralityPercentileMin),
175 fCutZDCtiming(that.fCutZDCtiming),
178 if (that.fQA) DefineHistograms();
180 if (that.fRefMultCuts)
181 fRefMultCuts = new AliFlowTrackCuts(*(that.fRefMultCuts));
182 if (that.fMeanPtCuts)
183 fMeanPtCuts = new AliFlowTrackCuts(*(that.fMeanPtCuts));
184 fStandardTPCcuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
185 fStandardGlobalCuts = AliFlowTrackCuts::GetStandardGlobalTrackCuts2010();
188 ////-----------------------------------------------------------------------
189 AliFlowEventCuts::~AliFlowEventCuts()
194 delete fStandardGlobalCuts;
195 delete fStandardTPCcuts;
196 if (fQA) { fQA->SetOwner(); fQA->Delete(); delete fQA; }
199 ////-----------------------------------------------------------------------
200 AliFlowEventCuts& AliFlowEventCuts::operator=(const AliFlowEventCuts& that)
203 if (this==&that) return *this;
212 fQA = static_cast<TList*>(that.fQA->Clone());
221 fCutNumberOfTracks=that.fCutNumberOfTracks;
222 fNumberOfTracksMax=that.fNumberOfTracksMax;
223 fNumberOfTracksMin=that.fNumberOfTracksMin;
224 fCutRefMult=that.fCutRefMult;
225 fRefMultMethod=that.fRefMultMethod;
226 fUseAliESDtrackCutsRefMult=that.fUseAliESDtrackCutsRefMult;
227 fRefMultMethodAliESDtrackCuts=that.fRefMultMethodAliESDtrackCuts;
228 fRefMultMax=that.fRefMultMax;
229 fRefMultMin=that.fRefMultMin;
230 if (that.fRefMultCuts) *fRefMultCuts=*(that.fRefMultCuts);
231 if (that.fMeanPtCuts) *fMeanPtCuts=*(that.fMeanPtCuts);
232 fStandardTPCcuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
233 fStandardGlobalCuts = AliFlowTrackCuts::GetStandardGlobalTrackCuts2010();
234 fCutPrimaryVertexX=that.fCutPrimaryVertexX;
235 fPrimaryVertexXmax=that.fPrimaryVertexXmax;
236 fPrimaryVertexXmin=that.fPrimaryVertexXmin;
237 fCutPrimaryVertexY=that.fCutPrimaryVertexY;
238 fPrimaryVertexYmax=that.fPrimaryVertexYmax;
239 fPrimaryVertexYmin=that.fPrimaryVertexYmin;
240 fCutPrimaryVertexZ=that.fCutPrimaryVertexZ;
241 fPrimaryVertexZmax=that.fPrimaryVertexZmax;
242 fPrimaryVertexZmin=that.fPrimaryVertexZmin;
243 fCutNContributors=that.fCutNContributors;
244 fNContributorsMax=that.fNContributorsMax;
245 fNContributorsMin=that.fNContributorsMin;
246 fCutMeanPt=that.fCutMeanPt;
247 fMeanPtMax=that.fMeanPtMax;
248 fMeanPtMin=that.fMeanPtMin;
249 fCutSPDvertexerAnomaly=that.fCutSPDvertexerAnomaly;
250 fCutTPCmultiplicityOutliers=that.fCutTPCmultiplicityOutliers;
251 fCutCentralityPercentile=that.fCutCentralityPercentile;
252 fUseCentralityUnchecked=that.fUseCentralityUnchecked;
253 fCentralityPercentileMethod=that.fCentralityPercentileMethod;
254 fCentralityPercentileMax=that.fCentralityPercentileMax;
255 fCentralityPercentileMin=that.fCentralityPercentileMin;
256 fCutZDCtiming=that.fCutZDCtiming;
260 //-----------------------------------------------------------------------
261 Bool_t AliFlowEventCuts::IsSelected(TObject* obj)
264 AliVEvent* vevent = dynamic_cast<AliVEvent*>(obj);
265 if (vevent) return PassesCuts(vevent);
266 return kFALSE; //when passed wrong type of object
268 //-----------------------------------------------------------------------
269 Bool_t AliFlowEventCuts::PassesCuts(AliVEvent *event)
271 ///check if event passes cuts
272 const AliVVertex* pvtx=event->GetPrimaryVertex();
273 Double_t pvtxx = pvtx->GetX();
274 Double_t pvtxy = pvtx->GetY();
275 Double_t pvtxz = pvtx->GetZ();
276 Int_t ncontrib = pvtx->GetNContributors();
278 AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
279 AliAODEvent* aodevent = dynamic_cast<AliAODEvent*>(event);
281 Int_t multGlobal = 0;
284 multTPC = fStandardTPCcuts->Count(event);
285 multGlobal = fStandardGlobalCuts->Count(event);
286 QAbefore(0)->Fill(pvtxz);
287 QAbefore(1)->Fill(multGlobal,multTPC);
289 if (fCutTPCmultiplicityOutliers)
291 //this is pretty slow as we check the event track by track twice
292 //this cut will work for 2010 PbPb data and is dependent on
293 //TPC and ITS reco efficiency (e.g. geometry, calibration etc)
296 multTPC = fStandardTPCcuts->Count(event);
297 multGlobal = fStandardGlobalCuts->Count(event);
299 if (multTPC > ( 23+1.216*multGlobal)) {pass=kFALSE;}
300 if (multTPC < (-20+1.087*multGlobal)) {pass=kFALSE;}
302 if (fCutNContributors)
304 if (ncontrib < fNContributorsMin || ncontrib >= fNContributorsMax) pass=kFALSE;
306 if (fCutPrimaryVertexX)
308 if (pvtxx < fPrimaryVertexXmin || pvtxx >= fPrimaryVertexXmax) pass=kFALSE;
310 if (fCutPrimaryVertexY)
312 if (pvtxy < fPrimaryVertexYmin || pvtxy >= fPrimaryVertexYmax) pass=kFALSE;
314 if (fCutPrimaryVertexZ)
316 if (pvtxz < fPrimaryVertexZmin || pvtxz >= fPrimaryVertexZmax)
319 if (fCutCentralityPercentile&&esdevent)
321 AliCentrality* centr = esdevent->GetCentrality();
322 if (fUseCentralityUnchecked)
324 if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
325 fCentralityPercentileMax,
326 CentrMethName(fCentralityPercentileMethod) ))
333 if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
334 fCentralityPercentileMax,
335 CentrMethName(fCentralityPercentileMethod) ))
341 if (fCutSPDvertexerAnomaly&&esdevent)
343 const AliESDVertex* sdpvertex = esdevent->GetPrimaryVertexSPD();
344 if (sdpvertex->GetNContributors()<1) pass=kFALSE;
345 if (sdpvertex->GetDispersion()>0.04) pass=kFALSE;
346 if (sdpvertex->GetZRes()>0.25) pass=kFALSE;
347 const AliESDVertex* tpcvertex = esdevent->GetPrimaryVertexTPC();
348 if (tpcvertex->GetNContributors()<1) pass=kFALSE;
349 const AliMultiplicity* tracklets = esdevent->GetMultiplicity();
350 if (tpcvertex->GetNContributors()<(-10.0+0.25*tracklets->GetNumberOfITSClusters(0)))
355 if (fCutZDCtiming&&esdevent)
357 if (!fTrigAna.ZDCTimeTrigger(esdevent))
362 if(fCutNumberOfTracks) {if ( event->GetNumberOfTracks() < fNumberOfTracksMin ||
363 event->GetNumberOfTracks() >= fNumberOfTracksMax ) pass=kFALSE;}
364 if(fCutRefMult&&esdevent)
366 //reference multiplicity still to be defined
367 Double_t refMult = RefMult(event);
368 if (refMult < fRefMultMin || refMult >= fRefMultMax )
376 if (fCutCentralityPercentile) {
377 AliCentrality* centr = aodevent->GetHeader()->GetCentralityP();
378 if (fUseCentralityUnchecked) {
379 if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
380 fCentralityPercentileMax,
381 CentrMethName(fCentralityPercentileMethod) )) {
385 if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
386 fCentralityPercentileMax,
387 CentrMethName(fCentralityPercentileMethod) )) {
397 Int_t ntracks=event->GetNumberOfTracks();
399 for (Int_t i=0; i<ntracks; i++)
401 AliVParticle* track = event->GetTrack(i);
402 if (!track) continue;
403 Bool_t localpass=kTRUE;
404 if (fMeanPtCuts) localpass=fMeanPtCuts->IsSelected(track);
407 meanpt += track->Pt();
411 if (nselected) meanpt=meanpt/nselected;
412 if (meanpt<fMeanPtMin || meanpt >= fMeanPtMax) pass=kFALSE;
416 QAafter(1)->Fill(multGlobal,multTPC);
417 QAafter(0)->Fill(pvtxz);
422 //-----------------------------------------------------------------------
423 const char* AliFlowEventCuts::CentrMethName(refMultMethod method) const
425 //get the string for refmultmethod, for use with AliCentrality in
426 //the cut on centrality percentile
441 //-----------------------------------------------------------------------
442 AliFlowEventCuts* AliFlowEventCuts::StandardCuts()
444 //make a set of standard event cuts, caller becomes owner
445 AliFlowEventCuts* cuts = new AliFlowEventCuts();
449 //-----------------------------------------------------------------------
450 Int_t AliFlowEventCuts::RefMult(AliVEvent* event)
452 //calculate the reference multiplicity, if all fails return 0
453 AliESDVZERO* vzero = NULL;
454 AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
456 if (fUseAliESDtrackCutsRefMult && esdevent)
458 //use the standard ALICE reference multiplicity with the default eta range
459 return AliESDtrackCuts::GetReferenceMultiplicity(esdevent, fRefMultMethodAliESDtrackCuts);
462 if (fRefMultMethod==kTPConly && !fRefMultCuts)
464 fRefMultCuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts();
465 fRefMultCuts->SetEtaRange(-0.8,0.8);
466 fRefMultCuts->SetPtMin(0.15);
468 else if (fRefMultMethod==kSPDtracklets && !fRefMultCuts)
470 fRefMultCuts = new AliFlowTrackCuts("tracklet refmult cuts");
471 fRefMultCuts->SetParamType(AliFlowTrackCuts::kSPDtracklet);
472 fRefMultCuts->SetEtaRange(-0.8,0.8);
474 else if (fRefMultMethod==kV0)
476 if (!esdevent) return 0;
477 vzero=esdevent->GetVZEROData();
478 if (!vzero) return 0;
479 return TMath::Nint(vzero->GetMTotV0A()+vzero->GetMTotV0C());
481 else if (fRefMultMethod==kSPD1clusters)
483 if (!esdevent) return 0;
484 const AliMultiplicity* mult = esdevent->GetMultiplicity();
486 return mult->GetNumberOfITSClusters(1);
490 fRefMultCuts->SetEvent(event);
491 for (Int_t i=0; i<fRefMultCuts->GetNumberOfInputObjects(); i++)
493 if (fRefMultCuts->IsSelected(fRefMultCuts->GetInputObject(i),i))
498 //_____________________________________________________________________________
499 void AliFlowEventCuts::DefineHistograms()
504 Bool_t adddirstatus = TH1::AddDirectoryStatus();
505 TH1::AddDirectory(kFALSE);
506 fQA = new TList(); fQA->SetOwner();
507 fQA->SetName(Form("%s QA",GetName()));
508 TList* before = new TList(); before->SetOwner();
509 before->SetName("before");
510 TList* after = new TList(); after->SetOwner();
511 after->SetName("after");
514 before->Add(new TH1F("zvertex",";z;event cout",500,-15.,15.)); //0
515 after->Add(new TH1F("zvertex",";z;event cout",500,-15.,15.)); //0
516 before->Add(new TH2F("fTPCvsGlobalMult","TPC only vs Global track multiplicity;global;TPC only",500,0,2500,500,0,3500));//1
517 after->Add(new TH2F("fTPCvsGlobalMult","TPC only vs Global track multiplicity;global;TPC only",500,0,2500,500,0,3500));//1
518 TH1::AddDirectory(adddirstatus);
521 //---------------------------------------------------------------//
522 void AliFlowEventCuts::Browse(TBrowser* b)
524 //some browsing capabilities
525 if (fQA) b->Add(fQA);
528 //---------------------------------------------------------------//
529 Long64_t AliFlowEventCuts::Merge(TCollection* list)
533 AliFlowEventCuts* obj;
535 if (list->GetEntries()<1) return 0;
537 while ( (obj = dynamic_cast<AliFlowEventCuts*>(next())) )
539 if (obj==this) continue;
541 listwrapper.Add(obj->GetQA());
542 fQA->Merge(&listwrapper);