]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWG/FLOW/Tasks/AliFlowEventCuts.cxx
addition vzero ep for phi and extension to higher pt
[u/mrichter/AliRoot.git] / PWG / FLOW / Tasks / AliFlowEventCuts.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
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  **************************************************************************/
15
16 /* $Id$ */ 
17
18 // AliFlowEventCuts:
19 // An event cut class for the flow framework
20 //
21 // origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
22
23 #include <limits.h>
24 #include <float.h>
25 #include <TList.h>
26 #include <TH1F.h>
27 #include <TH2F.h>
28 #include <TBrowser.h>
29 #include "TMath.h"
30 #include "TNamed.h"
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"
43
44 ClassImp(AliFlowEventCuts)
45
46 //-----------------------------------------------------------------------
47 AliFlowEventCuts::AliFlowEventCuts():
48   TNamed(),
49   fQA(NULL),
50   fCutNumberOfTracks(kFALSE),
51   fNumberOfTracksMax(INT_MAX),
52   fNumberOfTracksMin(INT_MIN),
53   fCutRefMult(kFALSE),
54   fRefMultMethod(kTPConly),
55   fUseAliESDtrackCutsRefMult(kFALSE),
56   fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
57   fRefMultMax(INT_MAX),
58   fRefMultMin(INT_MIN),
59   fRefMultCuts(NULL),
60   fMeanPtCuts(NULL),
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),
75   fCutMeanPt(kFALSE),
76   fMeanPtMax(-DBL_MAX),
77   fMeanPtMin(DBL_MAX),
78   fCutSPDvertexerAnomaly(kFALSE),
79   fCutTPCmultiplicityOutliers(kFALSE),
80   fCutCentralityPercentile(kFALSE),
81   fUseCentralityUnchecked(kFALSE),
82   fCentralityPercentileMethod(kTPConly),
83   fCentralityPercentileMax(100.),
84   fCentralityPercentileMin(0.),
85   fCutZDCtiming(kFALSE),
86   fTrigAna()
87 {
88   //constructor 
89 }
90
91 //-----------------------------------------------------------------------
92 AliFlowEventCuts::AliFlowEventCuts(const char* name, const char* title):
93   TNamed(name, title),
94   fQA(NULL),
95   fCutNumberOfTracks(kFALSE),
96   fNumberOfTracksMax(INT_MAX),
97   fNumberOfTracksMin(INT_MIN),
98   fCutRefMult(kFALSE),
99   fRefMultMethod(kTPConly),
100   fUseAliESDtrackCutsRefMult(kFALSE),
101   fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
102   fRefMultMax(INT_MAX),
103   fRefMultMin(INT_MIN),
104   fRefMultCuts(NULL),
105   fMeanPtCuts(NULL),
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),
120   fCutMeanPt(kFALSE),
121   fMeanPtMax(-DBL_MAX),
122   fMeanPtMin(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),
131   fTrigAna()
132 {
133   //constructor 
134 }
135
136 ////-----------------------------------------------------------------------
137 AliFlowEventCuts::AliFlowEventCuts(const AliFlowEventCuts& that):
138   TNamed(that),
139   fQA(NULL),
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),
149   fRefMultCuts(NULL),
150   fMeanPtCuts(NULL),
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),
176   fTrigAna()
177 {
178   if (that.fQA) DefineHistograms();
179   //copy constructor 
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();
186 }
187
188 ////-----------------------------------------------------------------------
189 AliFlowEventCuts::~AliFlowEventCuts()
190 {
191   //dtor
192   delete fMeanPtCuts;
193   delete fRefMultCuts;
194   delete fStandardGlobalCuts;
195   delete fStandardTPCcuts;
196   if (fQA) { fQA->SetOwner(); fQA->Delete(); delete fQA; }
197 }
198
199 ////-----------------------------------------------------------------------
200 AliFlowEventCuts& AliFlowEventCuts::operator=(const AliFlowEventCuts& that)
201 {
202   //assignment
203   if (this==&that) return *this;
204
205   if (that.fQA)
206   {
207     if (fQA)
208     {
209       fQA->Delete();
210       delete fQA;
211     }
212     fQA = static_cast<TList*>(that.fQA->Clone());
213   }
214   else
215   {
216     fQA->Delete();
217     delete fQA;
218     fQA=NULL;
219   }
220
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;
257   return *this;
258 }
259
260 //----------------------------------------------------------------------- 
261 Bool_t AliFlowEventCuts::IsSelected(TObject* obj)
262 {
263   //check cuts
264   AliVEvent* vevent = dynamic_cast<AliVEvent*>(obj);
265   if (vevent) return PassesCuts(vevent);
266   return kFALSE;  //when passed wrong type of object
267 }
268 //----------------------------------------------------------------------- 
269 Bool_t AliFlowEventCuts::PassesCuts(AliVEvent *event)
270 {
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();
277   Bool_t pass=kTRUE;
278   AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
279   AliAODEvent* aodevent = dynamic_cast<AliAODEvent*>(event);
280   Int_t multTPC = 0;
281   Int_t multGlobal = 0; 
282   if (fQA)
283   {
284     multTPC = fStandardTPCcuts->Count(event);
285     multGlobal = fStandardGlobalCuts->Count(event);
286     QAbefore(0)->Fill(pvtxz);
287     QAbefore(1)->Fill(multGlobal,multTPC);
288   }
289   if (fCutTPCmultiplicityOutliers)
290   {
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)
294     if (!fQA)
295     {
296       multTPC = fStandardTPCcuts->Count(event);
297       multGlobal = fStandardGlobalCuts->Count(event);
298     }
299     if (multTPC > ( 23+1.216*multGlobal)) {pass=kFALSE;}
300     if (multTPC < (-20+1.087*multGlobal)) {pass=kFALSE;}
301   }
302   if (fCutNContributors)
303   {
304     if (ncontrib < fNContributorsMin || ncontrib >= fNContributorsMax) pass=kFALSE;
305   }
306   if (fCutPrimaryVertexX)
307   {
308     if (pvtxx < fPrimaryVertexXmin || pvtxx >= fPrimaryVertexXmax) pass=kFALSE;
309   }
310   if (fCutPrimaryVertexY)
311   {
312     if (pvtxy < fPrimaryVertexYmin || pvtxy >= fPrimaryVertexYmax) pass=kFALSE;
313   }
314   if (fCutPrimaryVertexZ)
315   {
316     if (pvtxz < fPrimaryVertexZmin || pvtxz >= fPrimaryVertexZmax)
317       pass=kFALSE;
318   }
319   if (fCutCentralityPercentile&&esdevent)
320   {
321     AliCentrality* centr = esdevent->GetCentrality();
322     if (fUseCentralityUnchecked)
323     {
324       if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
325                                                      fCentralityPercentileMax,
326                                                      CentrMethName(fCentralityPercentileMethod) ))
327       {
328         pass=kFALSE;
329       }
330     }
331     else
332     {
333       if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
334                                             fCentralityPercentileMax,
335                                             CentrMethName(fCentralityPercentileMethod) ))
336       {
337         pass=kFALSE;
338       }
339     }
340   }
341   if (fCutSPDvertexerAnomaly&&esdevent)
342   {
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)))
351     {
352       pass=kFALSE;
353     }
354   }
355   if (fCutZDCtiming&&esdevent)
356   {
357     if (!fTrigAna.ZDCTimeTrigger(esdevent))
358     {
359       pass=kFALSE;
360     }
361   }
362   if(fCutNumberOfTracks) {if ( event->GetNumberOfTracks() < fNumberOfTracksMin ||
363                                event->GetNumberOfTracks() >= fNumberOfTracksMax ) pass=kFALSE;}
364   if(fCutRefMult&&esdevent)
365   {
366     //reference multiplicity still to be defined
367     Double_t refMult = RefMult(event);
368     if (refMult < fRefMultMin || refMult >= fRefMultMax )
369     {
370       pass=kFALSE;
371     }
372   }
373
374   // Handles AOD event
375   if(aodevent) {
376     if (fCutCentralityPercentile) {
377       AliCentrality* centr = aodevent->GetHeader()->GetCentralityP();
378       if (fUseCentralityUnchecked) {
379         if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
380                                                        fCentralityPercentileMax,
381                                                        CentrMethName(fCentralityPercentileMethod) )) {
382           pass = kFALSE;
383         }
384       } else {
385         if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
386                                               fCentralityPercentileMax,
387                                               CentrMethName(fCentralityPercentileMethod) )) {
388           pass = kFALSE;
389         }
390       }
391     }
392   }
393
394   if (fCutMeanPt)
395   {
396     Float_t meanpt=0.0;
397     Int_t ntracks=event->GetNumberOfTracks();
398     Int_t nselected=0;
399     for (Int_t i=0; i<ntracks; i++)
400     {
401       AliVParticle* track = event->GetTrack(i);
402       if (!track) continue;
403       Bool_t localpass=kTRUE;
404       if (fMeanPtCuts) localpass=fMeanPtCuts->IsSelected(track);
405       if (localpass) 
406       {
407         meanpt += track->Pt();
408         nselected++;
409       }
410     }
411     if (nselected) meanpt=meanpt/nselected;
412     if (meanpt<fMeanPtMin || meanpt >= fMeanPtMax) pass=kFALSE;
413   }
414   if (fQA&&pass) 
415   {
416     QAafter(1)->Fill(multGlobal,multTPC);
417     QAafter(0)->Fill(pvtxz);
418   }
419   return pass;
420 }
421
422 //----------------------------------------------------------------------- 
423 const char* AliFlowEventCuts::CentrMethName(refMultMethod method) const
424 {
425   //get the string for refmultmethod, for use with AliCentrality in
426   //the cut on centrality percentile
427   switch (method)
428   {
429     case kSPDtracklets:
430       return "TKL";
431     case kSPD1clusters:
432       return "CL1";
433     case kTPConly:
434       return "TRK";
435     case kV0:
436       return "V0M";
437     default:
438       return "";
439   }
440 }
441 //----------------------------------------------------------------------- 
442 AliFlowEventCuts* AliFlowEventCuts::StandardCuts()
443 {
444   //make a set of standard event cuts, caller becomes owner
445   AliFlowEventCuts* cuts = new AliFlowEventCuts();
446   return cuts;
447 }
448
449 //----------------------------------------------------------------------- 
450 Int_t AliFlowEventCuts::RefMult(AliVEvent* event)
451 {
452   //calculate the reference multiplicity, if all fails return 0
453   AliESDVZERO* vzero = NULL;
454   AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
455
456   if (fUseAliESDtrackCutsRefMult && esdevent)
457   {
458     //use the standard ALICE reference multiplicity with the default eta range
459     return AliESDtrackCuts::GetReferenceMultiplicity(esdevent, fRefMultMethodAliESDtrackCuts);
460   }
461
462   if (fRefMultMethod==kTPConly && !fRefMultCuts)
463   {
464     fRefMultCuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts();
465     fRefMultCuts->SetEtaRange(-0.8,0.8);
466     fRefMultCuts->SetPtMin(0.15);
467   }
468   else if (fRefMultMethod==kSPDtracklets && !fRefMultCuts)
469   {
470     fRefMultCuts = new AliFlowTrackCuts("tracklet refmult cuts");
471     fRefMultCuts->SetParamType(AliFlowTrackCuts::kSPDtracklet);
472     fRefMultCuts->SetEtaRange(-0.8,0.8);
473   }
474   else if (fRefMultMethod==kV0)
475   {
476     if (!esdevent) return 0;
477     vzero=esdevent->GetVZEROData();
478     if (!vzero) return 0;
479     return TMath::Nint(vzero->GetMTotV0A()+vzero->GetMTotV0C());
480   }
481   else if (fRefMultMethod==kSPD1clusters)
482   {
483     if (!esdevent) return 0;
484     const AliMultiplicity* mult = esdevent->GetMultiplicity();
485     if (!mult) return 0;
486     return mult->GetNumberOfITSClusters(1);
487   }
488
489   Int_t refmult=0;
490   fRefMultCuts->SetEvent(event);
491   for (Int_t i=0; i<fRefMultCuts->GetNumberOfInputObjects(); i++)
492   {
493     if (fRefMultCuts->IsSelected(fRefMultCuts->GetInputObject(i),i))
494       refmult++;
495   }
496   return refmult;
497 }
498 //_____________________________________________________________________________
499 void AliFlowEventCuts::DefineHistograms()
500 {
501   //define QA histos
502   if (fQA) return;
503
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");
512   fQA->Add(before);
513   fQA->Add(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);
519 }
520
521 //---------------------------------------------------------------//
522 void AliFlowEventCuts::Browse(TBrowser* b)
523 {
524   //some browsing capabilities
525   if (fQA) b->Add(fQA);
526 }
527
528 //---------------------------------------------------------------//
529 Long64_t AliFlowEventCuts::Merge(TCollection* list)
530 {
531   //merge
532   Int_t number=0;
533   AliFlowEventCuts* obj;
534   if (!list) return 0;
535   if (list->GetEntries()<1) return 0;
536   TIter next(list);
537   while ( (obj = dynamic_cast<AliFlowEventCuts*>(next())) )
538   {
539     if (obj==this) continue;
540     TList listwrapper;
541     listwrapper.Add(obj->GetQA());
542     fQA->Merge(&listwrapper);
543     number++;
544   }
545   return number;
546 }
547