PWG/FLOW/Tasks/AliFlowEventCuts.cxx

   1  /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   2  *                                                                        *
   3  * Author: The ALICE Off-line Project.                                    *
   4  * Contributors are mentioned in the code where appropriate.              *
   5  *                                                                        *
   6  * Permission to use, copy, modify and distribute this software and its   *
   7  * documentation strictly for non-commercial purposes is hereby granted   *
   8  * without fee, provided that the above copyright notice appears in all   *
   9  * copies and that both the copyright notice and this permission notice   *
  10  * appear in the supporting documentation. The authors make no claims     *
  11  * about the suitability of this software for any purpose. It is          *
  12  * provided "as is" without express or implied warranty.                  *
  13  **************************************************************************/
  14
  15 /* $Id$ */
  16
  17 // AliFlowEventCuts:
  18 // An event cut class for the flow framework
  19 //
  20 // origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
  21
  22 #include <limits.h>
  23 #include <float.h>
  24 #include <TList.h>
  25 #include <TH1F.h>
  26 #include <TH2F.h>
  27 #include <TBrowser.h>
  28 #include "TMath.h"
  29 #include "TNamed.h"
  30 #include "AliVVertex.h"
  31 #include "AliVEvent.h"
  32 #include "AliESDEvent.h"
  33 #include "AliAODEvent.h"
  34 #include "AliAODHeader.h"
  35 #include "AliCentrality.h"
  36 #include "AliESDVZERO.h"
  37 #include "AliMultiplicity.h"
  38 #include "AliMCEvent.h"
  39 #include "AliFlowEventCuts.h"
  40 #include "AliFlowTrackCuts.h"
  41 #include "AliTriggerAnalysis.h"
  42 #include "AliCollisionGeometry.h"
  43 #include "AliGenEventHeader.h"
  44 #include <iostream>
  45 using namespace std;
  46
  47 ClassImp(AliFlowEventCuts)
  48
  49 //-----------------------------------------------------------------------
  50 AliFlowEventCuts::AliFlowEventCuts():
  51   AliFlowEventSimpleCuts(),
  52   fQA(NULL),
  53   fCutNumberOfTracks(kFALSE),
  54   fNumberOfTracksMax(INT_MAX),
  55   fNumberOfTracksMin(INT_MIN),
  56   fCutRefMult(kFALSE),
  57   fRefMultMethod(kTPConly),
  58   fUseAliESDtrackCutsRefMult(kFALSE),
  59   fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
  60   fRefMultMax(INT_MAX),
  61   fRefMultMin(INT_MIN),
  62   fRefMultCuts(NULL),
  63   fMeanPtCuts(NULL),
  64   fStandardTPCcuts(NULL),
  65   fStandardGlobalCuts(NULL),
  66   fCutPrimaryVertexX(kFALSE),
  67   fPrimaryVertexXmax(INT_MAX),
  68   fPrimaryVertexXmin(INT_MIN),
  69   fCutPrimaryVertexY(kFALSE),
  70   fPrimaryVertexYmax(INT_MAX),
  71   fPrimaryVertexYmin(INT_MIN),
  72   fCutPrimaryVertexZ(kFALSE),
  73   fPrimaryVertexZmax(INT_MAX),
  74   fPrimaryVertexZmin(INT_MIN),
  75   fCutNContributors(kFALSE),
  76   fNContributorsMax(INT_MAX),
  77   fNContributorsMin(INT_MIN),
  78   fCutMeanPt(kFALSE),
  79   fMeanPtMax(-DBL_MAX),
  80   fMeanPtMin(DBL_MAX),
  81   fCutSPDvertexerAnomaly(kFALSE),
  82   fCutSPDTRKVtxZ(kFALSE),
  83   fCutTPCmultiplicityOutliers(kFALSE),
  84   fCutTPCmultiplicityOutliersAOD(kFALSE),
  85   fUseCentralityUnchecked(kFALSE),
  86   fCentralityPercentileMethod(kTPConly),
  87   fCutZDCtiming(kFALSE),
  88   fTrigAna(),
  89   fCutImpactParameter(kFALSE),
  90   fImpactParameterMin(0.0),
  91   fImpactParameterMax(100.0),
  92   fhistTPCvsGlobalMult(0),
  93   fData2011(kFALSE)
  94 {
  95   //constructor
  96 }
  97
  98 //-----------------------------------------------------------------------
  99 AliFlowEventCuts::AliFlowEventCuts(const char* name, const char* title):
 100   AliFlowEventSimpleCuts(name, title),
 101   fQA(NULL),
 102   fCutNumberOfTracks(kFALSE),
 103   fNumberOfTracksMax(INT_MAX),
 104   fNumberOfTracksMin(INT_MIN),
 105   fCutRefMult(kFALSE),
 106   fRefMultMethod(kTPConly),
 107   fUseAliESDtrackCutsRefMult(kFALSE),
 108   fRefMultMethodAliESDtrackCuts(AliESDtrackCuts::kTrackletsITSTPC),
 109   fRefMultMax(INT_MAX),
 110   fRefMultMin(INT_MIN),
 111   fRefMultCuts(NULL),
 112   fMeanPtCuts(NULL),
 113   fStandardTPCcuts(AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010()),
 114   fStandardGlobalCuts(AliFlowTrackCuts::GetStandardGlobalTrackCuts2010()),
 115   fCutPrimaryVertexX(kFALSE),
 116   fPrimaryVertexXmax(INT_MAX),
 117   fPrimaryVertexXmin(INT_MIN),
 118   fCutPrimaryVertexY(kFALSE),
 119   fPrimaryVertexYmax(INT_MAX),
 120   fPrimaryVertexYmin(INT_MIN),
 121   fCutPrimaryVertexZ(kFALSE),
 122   fPrimaryVertexZmax(INT_MAX),
 123   fPrimaryVertexZmin(INT_MIN),
 124   fCutNContributors(kFALSE),
 125   fNContributorsMax(INT_MAX),
 126   fNContributorsMin(INT_MIN),
 127   fCutMeanPt(kFALSE),
 128   fMeanPtMax(-DBL_MAX),
 129   fMeanPtMin(DBL_MAX),
 130   fCutSPDvertexerAnomaly(kFALSE),
 131   fCutSPDTRKVtxZ(kFALSE),
 132   fCutTPCmultiplicityOutliers(kFALSE),
 133   fCutTPCmultiplicityOutliersAOD(kFALSE),
 134   fUseCentralityUnchecked(kFALSE),
 135   fCentralityPercentileMethod(kTPConly),
 136   fCutZDCtiming(kFALSE),
 137   fTrigAna(),
 138   fCutImpactParameter(kFALSE),
 139   fImpactParameterMin(0.0),
 140   fImpactParameterMax(100.0),
 141   fhistTPCvsGlobalMult(0),
 142   fData2011(kFALSE)
 143 {
 144   //constructor
 145 }
 146
 147 ////-----------------------------------------------------------------------
 148 AliFlowEventCuts::AliFlowEventCuts(const AliFlowEventCuts& that):
 149   AliFlowEventSimpleCuts(that),
 150   fQA(NULL),
 151   fCutNumberOfTracks(that.fCutNumberOfTracks),
 152   fNumberOfTracksMax(that.fNumberOfTracksMax),
 153   fNumberOfTracksMin(that.fNumberOfTracksMin),
 154   fCutRefMult(that.fCutRefMult),
 155   fRefMultMethod(that.fRefMultMethod),
 156   fUseAliESDtrackCutsRefMult(that.fUseAliESDtrackCutsRefMult),
 157   fRefMultMethodAliESDtrackCuts(that.fRefMultMethodAliESDtrackCuts),
 158   fRefMultMax(that.fRefMultMax),
 159   fRefMultMin(that.fRefMultMin),
 160   fRefMultCuts(NULL),
 161   fMeanPtCuts(NULL),
 162   fStandardTPCcuts(NULL),
 163   fStandardGlobalCuts(NULL),
 164   fCutPrimaryVertexX(that.fCutPrimaryVertexX),
 165   fPrimaryVertexXmax(that.fPrimaryVertexXmax),
 166   fPrimaryVertexXmin(that.fPrimaryVertexXmin),
 167   fCutPrimaryVertexY(that.fCutPrimaryVertexX),
 168   fPrimaryVertexYmax(that.fPrimaryVertexYmax),
 169   fPrimaryVertexYmin(that.fPrimaryVertexYmin),
 170   fCutPrimaryVertexZ(that.fCutPrimaryVertexX),
 171   fPrimaryVertexZmax(that.fPrimaryVertexZmax),
 172   fPrimaryVertexZmin(that.fPrimaryVertexZmin),
 173   fCutNContributors(that.fCutNContributors),
 174   fNContributorsMax(that.fNContributorsMax),
 175   fNContributorsMin(that.fNContributorsMin),
 176   fCutMeanPt(that.fCutMeanPt),
 177   fMeanPtMax(that.fMeanPtMax),
 178   fMeanPtMin(that.fMeanPtMin),
 179   fCutSPDvertexerAnomaly(that.fCutSPDvertexerAnomaly),
 180   fCutSPDTRKVtxZ(that.fCutSPDTRKVtxZ),
 181   fCutTPCmultiplicityOutliers(that.fCutTPCmultiplicityOutliers),
 182   fCutTPCmultiplicityOutliersAOD(that.fCutTPCmultiplicityOutliersAOD),
 183   fUseCentralityUnchecked(that.fUseCentralityUnchecked),
 184   fCentralityPercentileMethod(that.fCentralityPercentileMethod),
 185   fCutZDCtiming(that.fCutZDCtiming),
 186   fTrigAna(),
 187   fCutImpactParameter(that.fCutImpactParameter),
 188   fImpactParameterMin(that.fImpactParameterMin),
 189   fImpactParameterMax(that.fImpactParameterMax),
 190   fhistTPCvsGlobalMult(that.fhistTPCvsGlobalMult),
 191   fData2011(that.fData2011)
 192 {
 193   if (that.fQA) DefineHistograms();
 194   //copy constructor
 195   if (that.fRefMultCuts)
 196     fRefMultCuts = new AliFlowTrackCuts(*(that.fRefMultCuts));
 197   if (that.fMeanPtCuts)
 198     fMeanPtCuts = new AliFlowTrackCuts(*(that.fMeanPtCuts));
 199   fStandardTPCcuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
 200   fStandardGlobalCuts = AliFlowTrackCuts::GetStandardGlobalTrackCuts2010();
 201 }
 202
 203 ////-----------------------------------------------------------------------
 204 AliFlowEventCuts::~AliFlowEventCuts()
 205 {
 206   //dtor
 207   delete fMeanPtCuts;
 208   delete fRefMultCuts;
 209   delete fStandardGlobalCuts;
 210   delete fStandardTPCcuts;
 211   if (fQA) { fQA->SetOwner(); fQA->Delete(); delete fQA; }
 212 }
 213
 214 ////-----------------------------------------------------------------------
 215 AliFlowEventCuts& AliFlowEventCuts::operator=(const AliFlowEventCuts& that)
 216 {
 217   //assignment
 218   if (this==&that) return *this;
 219
 220   if (that.fQA)
 221   {
 222     if (fQA)
 223     {
 224       fQA->Delete();
 225       delete fQA;
 226     }
 227     fQA = static_cast<TList*>(that.fQA->Clone());
 228   }
 229   else
 230   {
 231     fQA->Delete();
 232     delete fQA;
 233     fQA=NULL;
 234   }
 235
 236   fCutNumberOfTracks=that.fCutNumberOfTracks;
 237   fNumberOfTracksMax=that.fNumberOfTracksMax;
 238   fNumberOfTracksMin=that.fNumberOfTracksMin;
 239   fCutRefMult=that.fCutRefMult;
 240   fRefMultMethod=that.fRefMultMethod;
 241   fUseAliESDtrackCutsRefMult=that.fUseAliESDtrackCutsRefMult;
 242   fRefMultMethodAliESDtrackCuts=that.fRefMultMethodAliESDtrackCuts;
 243   fRefMultMax=that.fRefMultMax;
 244   fRefMultMin=that.fRefMultMin;
 245   if (that.fRefMultCuts) *fRefMultCuts=*(that.fRefMultCuts);
 246   if (that.fMeanPtCuts) *fMeanPtCuts=*(that.fMeanPtCuts);
 247   fStandardTPCcuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
 248   fStandardGlobalCuts = AliFlowTrackCuts::GetStandardGlobalTrackCuts2010();
 249   fCutPrimaryVertexX=that.fCutPrimaryVertexX;
 250   fPrimaryVertexXmax=that.fPrimaryVertexXmax;
 251   fPrimaryVertexXmin=that.fPrimaryVertexXmin;
 252   fCutPrimaryVertexY=that.fCutPrimaryVertexY;
 253   fPrimaryVertexYmax=that.fPrimaryVertexYmax;
 254   fPrimaryVertexYmin=that.fPrimaryVertexYmin;
 255   fCutPrimaryVertexZ=that.fCutPrimaryVertexZ;
 256   fPrimaryVertexZmax=that.fPrimaryVertexZmax;
 257   fPrimaryVertexZmin=that.fPrimaryVertexZmin;
 258   fCutNContributors=that.fCutNContributors;
 259   fNContributorsMax=that.fNContributorsMax;
 260   fNContributorsMin=that.fNContributorsMin;
 261   fCutMeanPt=that.fCutMeanPt;
 262   fMeanPtMax=that.fMeanPtMax;
 263   fMeanPtMin=that.fMeanPtMin;
 264   fCutSPDvertexerAnomaly=that.fCutSPDvertexerAnomaly;
 265   fCutSPDTRKVtxZ=that.fCutSPDTRKVtxZ;
 266   fCutTPCmultiplicityOutliers=that.fCutTPCmultiplicityOutliers;
 267   fCutTPCmultiplicityOutliersAOD=that.fCutTPCmultiplicityOutliersAOD;
 268   fUseCentralityUnchecked=that.fUseCentralityUnchecked;
 269   fCentralityPercentileMethod=that.fCentralityPercentileMethod;
 270   fCutZDCtiming=that.fCutZDCtiming;
 271   fhistTPCvsGlobalMult=that.fhistTPCvsGlobalMult;
 272   fData2011=that.fData2011;
 273   return *this;
 274 }
 275
 276 //-----------------------------------------------------------------------
 277 Bool_t AliFlowEventCuts::IsSelected(TObject* obj, TObject* objmc)
 278 {
 279   //check cuts
 280   AliVEvent* vevent = dynamic_cast<AliVEvent*>(obj);
 281   AliMCEvent* mcevent = dynamic_cast<AliMCEvent*>(objmc);
 282   if (vevent) return PassesCuts(vevent,mcevent);;
 283   return kFALSE;  //when passed wrong type of object
 284 }
 285 //-----------------------------------------------------------------------
 286 Bool_t AliFlowEventCuts::PassesCuts(AliVEvent *event, AliMCEvent *mcevent)
 287 {
 288   ///check if event passes cuts
 289   const AliVVertex* pvtx=event->GetPrimaryVertex();
 290   Double_t pvtxx = pvtx->GetX();
 291   Double_t pvtxy = pvtx->GetY();
 292   Double_t pvtxz = pvtx->GetZ();
 293   Int_t ncontrib = pvtx->GetNContributors();
 294   Bool_t pass=kTRUE;
 295   AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
 296   AliAODEvent* aodevent = dynamic_cast<AliAODEvent*>(event);
 297   Int_t multTPC = 0;
 298   Int_t multGlobal = 0;
 299
 300   // to remove multiplicity outliers, an explicit cut on the correlation
 301   // between global and tpc only tracks can be made by counting the two
 302   // sets. as counting is expensive, only count when qa is requested or cut is enabeled
 303   // the cut criteria are different for different data takign periods
 304   // and (of course) cut criteria, specific routines exist for 10h, 11h data
 305   // and esds (by calling AliFlowTrackCuts) or aods (evaluated here explicitely)
 306   if(esdevent && (fQA || fCutTPCmultiplicityOutliers))
 307   {
 308     //this is pretty slow as we check the event track by track twice
 309     //this cut will work for 2010 PbPb data and is dependent on
 310     //TPC and ITS reco efficiency (e.g. geometry, calibration etc)
 311     multTPC = fStandardTPCcuts->Count(event);
 312     multGlobal = fStandardGlobalCuts->Count(event);
 313     if(fCutTPCmultiplicityOutliers) {
 314       if (multTPC > ( 23+1.216*multGlobal)) pass = kFALSE;
 315       if (multTPC < (-20+1.087*multGlobal)) pass = kFALSE;
 316     }
 317   }
 318   if(aodevent && (fQA || fCutTPCmultiplicityOutliersAOD))
 319   {
 320     //similar (slow) cut for aod's. will work for both 2010 and 2010 pbpb data
 321     //but the user is responsible that this object is configured
 322     //correctly to select the dataset
 323     //FIXME data could dynamically be determined by this class via the
 324     //runnumber
 325     Int_t nTracks(aodevent->GetNumberOfTracks());
 326     for(Int_t iTracks = 0; iTracks < nTracks; iTracks++) {
 327       AliAODTrack* track = dynamic_cast<AliAODTrack*>(aodevent->GetTrack(iTracks));
 328       if(!track) continue;
 329       if (!track || track->Pt() < .2 || track->Pt() > 5.0 || TMath::Abs(track->Eta()) > .8 || track->GetTPCNcls() < 70 || !track->GetDetPid() || track->GetDetPid()->GetTPCsignal() < 10.0)  continue;  // general quality cut
 330       if (track->TestFilterBit(1) && track->Chi2perNDF() > 0.2) multTPC++;
 331       if (!track->TestFilterBit(16) || track->Chi2perNDF() < 0.1) continue;
 332       Double_t b[2] = {-99., -99.};
 333       Double_t bCov[3] = {-99., -99., -99.};
 334       AliAODTrack copy(*track);
 335       if (copy.PropagateToDCA(event->GetPrimaryVertex(), event->GetMagneticField(), 100., b, bCov) && TMath::Abs(b[0]) < 0.3 && TMath::Abs(b[1]) < 0.3) multGlobal++;
 336     }
 337     if(fCutTPCmultiplicityOutliersAOD)
 338     {
 339       if(!fData2011 && (multTPC < (-40.3+1.22*multGlobal) || multTPC > (32.1+1.59*multGlobal))) pass = kFALSE;
 340       else if(fData2011  && (multTPC < (-36.73 + 1.48*multGlobal) || multTPC > (62.87 + 1.78*multGlobal))) pass = kFALSE;
 341     }
 342   }
 343
 344   if (fQA)
 345   {
 346     QAbefore(0)->Fill(pvtxz);
 347     QAbefore(1)->Fill(multGlobal,multTPC);
 348   }
 349
 350   if (fCutNContributors)
 351   {
 352     if (ncontrib < fNContributorsMin || ncontrib >= fNContributorsMax) pass=kFALSE;
 353   }
 354   if (fCutPrimaryVertexX)
 355   {
 356     if (pvtxx < fPrimaryVertexXmin || pvtxx >= fPrimaryVertexXmax) pass=kFALSE;
 357   }
 358   if (fCutPrimaryVertexY)
 359   {
 360     if (pvtxy < fPrimaryVertexYmin || pvtxy >= fPrimaryVertexYmax) pass=kFALSE;
 361   }
 362   if (fCutPrimaryVertexZ)
 363   {
 364     if (pvtxz < fPrimaryVertexZmin || pvtxz >= fPrimaryVertexZmax)
 365       pass=kFALSE;
 366   }
 367   if (fCutCentralityPercentile&&esdevent)
 368   {
 369     AliCentrality* centr = esdevent->GetCentrality();
 370     if (fUseCentralityUnchecked)
 371     {
 372       if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
 373                                                      fCentralityPercentileMax,
 374                                                      CentrMethName(fCentralityPercentileMethod) ))
 375       {
 376         pass=kFALSE;
 377       }
 378     }
 379     else
 380     {
 381       if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
 382                                             fCentralityPercentileMax,
 383                                             CentrMethName(fCentralityPercentileMethod) ))
 384       {
 385         pass=kFALSE;
 386       }
 387     }
 388   }
 389   if (fCutSPDvertexerAnomaly&&esdevent)
 390   {
 391     const AliESDVertex* sdpvertex = esdevent->GetPrimaryVertexSPD();
 392     if (sdpvertex->GetNContributors()<1) pass=kFALSE;
 393     if (sdpvertex->GetDispersion()>0.04) pass=kFALSE;
 394     if (sdpvertex->GetZRes()>0.25) pass=kFALSE;
 395     const AliESDVertex* tpcvertex = esdevent->GetPrimaryVertexTPC();
 396     if (tpcvertex->GetNContributors()<1) pass=kFALSE;
 397     const AliMultiplicity* tracklets = esdevent->GetMultiplicity();
 398     if (tpcvertex->GetNContributors()<(-10.0+0.25*tracklets->GetNumberOfITSClusters(0)))
 399     {
 400       pass=kFALSE;
 401     }
 402   }
 403   if (fCutZDCtiming&&esdevent)
 404   {
 405     if (!fTrigAna.ZDCTimeTrigger(esdevent))
 406     {
 407       pass=kFALSE;
 408     }
 409   }
 410   if(fCutNumberOfTracks) {if ( event->GetNumberOfTracks() < fNumberOfTracksMin ||
 411                                event->GetNumberOfTracks() >= fNumberOfTracksMax ) pass=kFALSE;}
 412   if((fCutRefMult&&mcevent)||(fCutRefMult&&esdevent))
 413   {
 414     //reference multiplicity still to be defined
 415     Double_t refMult = RefMult(event,mcevent);
 416     if (refMult < fRefMultMin || refMult >= fRefMultMax )
 417     {
 418       pass=kFALSE;
 419     }
 420   }
 421
 422   // Handles AOD event
 423   if(aodevent) {
 424     if(fCutSPDTRKVtxZ) {
 425       Double_t tVtxZ = aodevent->GetPrimaryVertex()->GetZ();
 426       Double_t tSPDVtxZ = aodevent->GetPrimaryVertexSPD()->GetZ();
 427       if( TMath::Abs(tVtxZ-tSPDVtxZ) > 0.5 ) pass = kFALSE;
 428     }
 429     AliCentrality* centr = ((AliVAODHeader*)aodevent->GetHeader())->GetCentralityP();
 430     if(fCutTPCmultiplicityOutliers || fCutTPCmultiplicityOutliersAOD){
 431       Double_t v0Centr  = centr->GetCentralityPercentile("V0M");
 432       Double_t trkCentr = centr->GetCentralityPercentile("TRK");
 433       if(TMath::Abs(v0Centr-trkCentr) > 5) pass = kFALSE;
 434     }
 435     if (fCutCentralityPercentile) {
 436       if (fUseCentralityUnchecked) {
 437         if (!centr->IsEventInCentralityClassUnchecked( fCentralityPercentileMin,
 438                                                        fCentralityPercentileMax,
 439                                                        CentrMethName(fCentralityPercentileMethod) )) {
 440           pass = kFALSE;
 441         }
 442       } else {
 443         if (!centr->IsEventInCentralityClass( fCentralityPercentileMin,
 444                                               fCentralityPercentileMax,
 445                                               CentrMethName(fCentralityPercentileMethod) )) {
 446           pass = kFALSE;
 447         }
 448       }
 449     }
 450   }
 451
 452   if (fCutMeanPt)
 453   {
 454     Float_t meanpt=0.0;
 455     Int_t ntracks=event->GetNumberOfTracks();
 456     Int_t nselected=0;
 457     for (Int_t i=0; i<ntracks; i++)
 458     {
 459       AliVParticle* track = event->GetTrack(i);
 460       if (!track) continue;
 461       Bool_t localpass=kTRUE;
 462       if (fMeanPtCuts) localpass=fMeanPtCuts->IsSelected(track);
 463       if (localpass)
 464       {
 465         meanpt += track->Pt();
 466         nselected++;
 467       }
 468     }
 469     if (nselected) meanpt=meanpt/nselected;
 470     if (meanpt<fMeanPtMin || meanpt >= fMeanPtMax) pass=kFALSE;
 471   }
 472
 473   //impact parameter cut
 474   if(fCutImpactParameter) {
 475     Double_t gImpactParameter = 0.;
 476     if(mcevent) {
 477       AliCollisionGeometry* headerH = dynamic_cast<AliCollisionGeometry*>(dynamic_cast<AliMCEvent*>(mcevent)->GenEventHeader());
 478       if(headerH)
 479         gImpactParameter = headerH->ImpactParameter();
 480     }
 481     if ((gImpactParameter < fImpactParameterMin) || (gImpactParameter >= fImpactParameterMax ))
 482       pass=kFALSE;
 483   }
 484
 485   if (fQA&&pass)
 486   {
 487     QAafter(1)->Fill(multGlobal,multTPC);
 488     QAafter(0)->Fill(pvtxz);
 489   }
 490   return pass;
 491 }
 492
 493 //-----------------------------------------------------------------------
 494 Float_t AliFlowEventCuts::GetCentrality(AliVEvent* event, AliMCEvent* /*mcEvent*/)
 495 {
 496   //get the centrality percentile of the event
 497   AliESDEvent* esdEvent = dynamic_cast<AliESDEvent*>(event);
 498   AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(event);
 499
 500   Float_t centrality=-1.;
 501
 502   AliCentrality* centr = NULL;
 503   if (esdEvent)
 504     centr = esdEvent->GetCentrality();
 505   if (aodEvent)
 506     centr = ((AliVAODHeader*)aodEvent->GetHeader())->GetCentralityP();
 507
 508   if (!centr) return -1.;
 509
 510   if (fUseCentralityUnchecked)
 511     centrality=centr->GetCentralityPercentileUnchecked(CentrMethName(fCentralityPercentileMethod));
 512   else
 513     centrality=centr->GetCentralityPercentile(CentrMethName(fCentralityPercentileMethod));
 514
 515   return centrality;
 516 }
 517
 518 //-----------------------------------------------------------------------
 519 const char* AliFlowEventCuts::CentrMethName(refMultMethod method) const
 520 {
 521   //get the string for refmultmethod, for use with AliCentrality in
 522   //the cut on centrality percentile
 523   switch (method)
 524   {
 525     case kSPDtracklets:
 526       return "TKL";
 527     case kSPD1clusters:
 528       return "CL1";
 529     case kTPConly:
 530       return "TRK";
 531     case kVZERO:
 532       return "V0M";
 533     default:
 534       return "";
 535   }
 536 }
 537 //-----------------------------------------------------------------------
 538 AliFlowEventCuts* AliFlowEventCuts::StandardCuts()
 539 {
 540   //make a set of standard event cuts, caller becomes owner
 541   AliFlowEventCuts* cuts = new AliFlowEventCuts();
 542   return cuts;
 543 }
 544
 545 //-----------------------------------------------------------------------
 546 Int_t AliFlowEventCuts::RefMult(AliVEvent* event, AliMCEvent *mcEvent)
 547 {
 548   //calculate the reference multiplicity, if all fails return 0
 549   AliESDVZERO* vzero = NULL;
 550   AliESDEvent* esdevent = dynamic_cast<AliESDEvent*>(event);
 551
 552   if (fUseAliESDtrackCutsRefMult && esdevent)
 553   {
 554     //use the standard ALICE reference multiplicity with the default eta range
 555     return AliESDtrackCuts::GetReferenceMultiplicity(esdevent, fRefMultMethodAliESDtrackCuts);
 556   }
 557
 558   if (fRefMultMethod==kTPConly && !fRefMultCuts)
 559   {
 560     fRefMultCuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts();
 561     fRefMultCuts->SetEtaRange(-0.8,0.8);
 562     fRefMultCuts->SetPtMin(0.15);
 563   }
 564   else if (fRefMultMethod==kSPDtracklets && !fRefMultCuts)
 565   {
 566     fRefMultCuts = new AliFlowTrackCuts("tracklet refmult cuts");
 567     fRefMultCuts->SetParamType(AliFlowTrackCuts::kSPDtracklet);
 568     fRefMultCuts->SetEtaRange(-0.8,0.8);
 569   }
 570   else if (fRefMultMethod==kVZERO)
 571   {
 572     if (!esdevent) return 0;
 573     vzero=esdevent->GetVZEROData();
 574     if (!vzero) return 0;
 575     return TMath::Nint(vzero->GetMTotV0A()+vzero->GetMTotV0C());
 576   }
 577   else if (fRefMultMethod==kSPD1clusters)
 578   {
 579     if (!esdevent) return 0;
 580     const AliMultiplicity* mult = esdevent->GetMultiplicity();
 581     if (!mult) return 0;
 582     return mult->GetNumberOfITSClusters(1);
 583   }
 584
 585   Int_t refmult=0;
 586   fRefMultCuts->SetEvent(event,mcEvent);
 587   Int_t numberOfInputObjects = fRefMultCuts->GetNumberOfInputObjects();
 588   for (Int_t i=0; i<numberOfInputObjects; i++) {
 589     if (fRefMultCuts->IsSelected(fRefMultCuts->GetInputObject(i),i))
 590       refmult++;
 591   }
 592   return refmult;
 593 }
 594 //_____________________________________________________________________________
 595 void AliFlowEventCuts::DefineHistograms()
 596 {
 597   //define QA histos
 598   if (fQA) return;
 599
 600   Bool_t adddirstatus = TH1::AddDirectoryStatus();
 601   TH1::AddDirectory(kFALSE);
 602   fQA = new TList(); fQA->SetOwner();
 603   fQA->SetName(Form("%s QA",GetName()));
 604   TList* before = new TList(); before->SetOwner();
 605   before->SetName("before");
 606   TList* after = new TList(); after->SetOwner();
 607   after->SetName("after");
 608   fQA->Add(before);
 609   fQA->Add(after);
 610   before->Add(new TH1F("zvertex",";z;event cout",500,-15.,15.)); //0
 611   after->Add(new TH1F("zvertex",";z;event cout",500,-15.,15.)); //0
 612   before->Add(new TH2F("fTPCvsGlobalMult","TPC only vs Global track multiplicity;global;TPC only",500,0,2500,500,0,3500));//1
 613   after->Add(new TH2F("fTPCvsGlobalMult","TPC only vs Global track multiplicity;global;TPC only",500,0,2500,500,0,3500));//1
 614   TH1::AddDirectory(adddirstatus);
 615 }
 616
 617 //---------------------------------------------------------------//
 618 void AliFlowEventCuts::Browse(TBrowser* b)
 619 {
 620   //some browsing capabilities
 621   if (fQA) b->Add(fQA);
 622 }
 623
 624 //---------------------------------------------------------------//
 625 Long64_t AliFlowEventCuts::Merge(TCollection* list)
 626 {
 627   //merge
 628   Int_t number=0;
 629   AliFlowEventCuts* obj;
 630   if (!list) return 0;
 631   if (list->GetEntries()<1) return 0;
 632   TIter next(list);
 633   while ( (obj = dynamic_cast<AliFlowEventCuts*>(next())) )
 634   {
 635     if (obj==this) continue;
 636     TList listwrapper;
 637     listwrapper.Add(obj->GetQA());
 638     fQA->Merge(&listwrapper);
 639     number++;
 640   }
 641   return number;
 642 }
 643
 644