New AliFlowTrackCuts and AliFlowEventCuts, allow running of centrality train
authorsnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 11 Oct 2010 20:23:01 +0000 (20:23 +0000)
committersnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 11 Oct 2010 20:23:01 +0000 (20:23 +0000)
18 files changed:
PWG2/FLOW/AliFlowCommon/AliFlowTrackSimpleCuts.cxx
PWG2/FLOW/AliFlowCommon/AliFlowTrackSimpleCuts.h
PWG2/FLOW/AliFlowTasks/AliAnalysisTaskFlowEvent.cxx
PWG2/FLOW/AliFlowTasks/AliAnalysisTaskFlowEvent.h
PWG2/FLOW/AliFlowTasks/AliFlowEvent.cxx
PWG2/FLOW/AliFlowTasks/AliFlowEvent.h
PWG2/FLOW/AliFlowTasks/AliFlowEventCuts.cxx [new file with mode: 0644]
PWG2/FLOW/AliFlowTasks/AliFlowEventCuts.h [new file with mode: 0644]
PWG2/FLOW/AliFlowTasks/AliFlowTrack.cxx
PWG2/FLOW/AliFlowTasks/AliFlowTrack.h
PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.cxx [new file with mode: 0644]
PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.h [new file with mode: 0644]
PWG2/FLOW/Documentation/FLOW.pdf
PWG2/FLOW/macros/AddTaskFlowCentrality.C [new file with mode: 0644]
PWG2/FLOW/macros/runFlowTask.C
PWG2/FLOW/macros/runFlowTaskCentralityTrain.C [new file with mode: 0644]
PWG2/PWG2flowTasksLinkDef.h
PWG2/libPWG2flowTasks.pkg

index 3c77ed5..2128385 100644 (file)
@@ -94,14 +94,25 @@ AliFlowTrackSimpleCuts::AliFlowTrackSimpleCuts():
 //}
 
 //----------------------------------------------------------------------- 
+Bool_t AliFlowTrackSimpleCuts::IsSelected(TObject* obj)
+{
+  //check cuts
+  TParticle* p = dynamic_cast<TParticle*>(obj);
+  if (p) return PassesCuts(p);
+  AliFlowTrackSimple* ts = dynamic_cast<AliFlowTrackSimple*>(obj);
+  if (ts) return PassesCuts(ts);
+  return kFALSE; //default when passed a wrong type of object
+}
+
+//----------------------------------------------------------------------- 
 Bool_t AliFlowTrackSimpleCuts::PassesCuts(const AliFlowTrackSimple *track) const
 {
   //simple method to check if the simple track passes the simple cuts
   if(fCutPt) {if (track->Pt() < fPtMin || track->Pt() >= fPtMax ) return kFALSE;}
   if(fCutEta) {if (track->Eta() < fEtaMin || track->Eta() >= fEtaMax ) return kFALSE;}
   if(fCutPhi) {if (track->Phi() < fPhiMin || track->Phi() >= fPhiMax ) return kFALSE;}
-  //if(fCutPID) {if (track->PID() != fPID) return kFALSE;}
   if(fCutCharge) {if (track->Charge() != fCharge) return kFALSE;}
+  //if(fCutPID) {if (track->PID() != fPID) return kFALSE;}
   return kTRUE;
 }
 
@@ -119,7 +130,7 @@ Bool_t AliFlowTrackSimpleCuts::PassesCuts(TParticle* track) const
   if (fCutCharge) 
   {
     TParticlePDG* ppdg = track->GetPDG();
-    Int_t charge = TMath::Nint(ppdg->Charge()/3.0);
+    Int_t charge = TMath::Nint(ppdg->Charge()/3.0); //mc particles have charge in units of 1/3e
     return (charge==fCharge);
   }
   return kTRUE;
index b022a1e..699c1d1 100644 (file)
@@ -48,7 +48,9 @@ class AliFlowTrackSimpleCuts : public TNamed {
   Bool_t PassesCuts(const AliFlowTrackSimple *track) const;
   Bool_t PassesCuts(TParticle* p) const;
 
- private:
+  virtual Bool_t IsSelected(TObject* obj);
+
+ protected:
   Bool_t   fCutPt;
   Double_t fPtMax;
   Double_t fPtMin;
index eb38e97..b9a4a3a 100644 (file)
@@ -62,6 +62,8 @@
 
 // Interface to make the Flow Event Simple used in the flow analysis methods
 #include "AliFlowEvent.h"
+#include "AliFlowTrackCuts.h"
+#include "AliFlowEventCuts.h"
 #include "AliFlowCommonConstants.h"
 #include "AliAnalysisTaskFlowEvent.h"
 
@@ -77,6 +79,9 @@ AliAnalysisTaskFlowEvent::AliAnalysisTaskFlowEvent() :
   fRPType("Global"),
   fCFManager1(NULL),
   fCFManager2(NULL),
+  fCutsEvent(NULL),
+  fCutsRP(NULL),
+  fCutsPOI(NULL),
   fQAInt(NULL),
   fQADiff(NULL),
   fMinMult(0),
@@ -125,6 +130,9 @@ AliAnalysisTaskFlowEvent::AliAnalysisTaskFlowEvent(const char *name, TString RPt
   fRPType(RPtype),
   fCFManager1(NULL),
   fCFManager2(NULL),
+  fCutsEvent(NULL),
+  fCutsRP(NULL),
+  fCutsPOI(NULL),
   fQAInt(NULL),
   fQADiff(NULL),
   fMinMult(0),
@@ -253,6 +261,20 @@ void AliAnalysisTaskFlowEvent::UserExec(Option_t *)
     }
   }
   
+  //use the new and temporarily inclomplete way of doing things
+  if (fAnalysisType == "MK")
+  {
+    if (!(fCutsRP&&fCutsPOI))
+    {
+      AliError("cuts not set");
+      return;
+    }
+    fCutsRP->SetMCevent( MCEvent() );
+    fCutsPOI->SetMCevent( MCEvent() );
+    flowEvent = new AliFlowEvent( InputEvent(), fCutsRP, fCutsPOI );
+    if (myESD)
+      flowEvent->SetReferenceMultiplicity(AliESDtrackCuts::GetReferenceMultiplicity(myESD,kTRUE));
+  }
 
   // Make the FlowEvent for MC input
   if (fAnalysisType == "MC")
index c3055ee..c03c9ac 100644 (file)
@@ -12,6 +12,8 @@
 #define AliAnalysisTaskFlowEvent_H
 
 class AliCFManager;
+class AliFlowEventCuts;
+class AliFlowTrackCuts;
 class AliFlowEventSimpleMaker;
 class TList;
 class TRandom3;
@@ -50,6 +52,13 @@ class AliAnalysisTaskFlowEvent : public AliAnalysisTaskSE {
   {this->fExcludedEtaMin = etaMin; this->fExcludedEtaMax = etaMax; 
     this->fExcludedPhiMin = phiMin; this->fExcludedPhiMax = phiMax; }
 
+  void          SetCutsEvent(AliFlowEventCuts* cutsEvent) {fCutsEvent=cutsEvent;}
+  AliFlowEventCuts* GetCutsEvent() const {return fCutsEvent;}
+  void          SetCutsRP(AliFlowTrackCuts* cutsRP) {fCutsRP=cutsRP;}
+  AliFlowTrackCuts* GetCutsRP() const {return fCutsRP;}
+  void          SetCutsPOI(AliFlowTrackCuts* cutsPOI) {fCutsPOI=cutsPOI;}
+  AliFlowTrackCuts* GetCutsPOI() const {return fCutsPOI;}
+
   void          SetCFManager1(AliCFManager* cfmgr) {this->fCFManager1 = cfmgr; } 
   AliCFManager* GetCFManager1()           {return this->fCFManager1; }
   void          SetCFManager2(AliCFManager* cfmgr) {this->fCFManager2 = cfmgr; } 
@@ -98,6 +107,9 @@ class AliAnalysisTaskFlowEvent : public AliAnalysisTaskSE {
   TString       fRPType;            // can be Global or Tracklet or FMD
   AliCFManager* fCFManager1;        // correction framework manager
   AliCFManager* fCFManager2;        // correction framework manager
+  AliFlowEventCuts* fCutsEvent;     //event cuts
+  AliFlowTrackCuts* fCutsRP;        //cuts for RPs
+  AliFlowTrackCuts* fCutsPOI;       //cuts for POIs
   TList*        fQAInt;             // QA histogram list
   TList*        fQADiff;            // QA histogram list
   Int_t         fMinMult;           // Minimum multiplicity from tracks selected using CORRFW
index 517f067..537ddf8 100644 (file)
@@ -33,7 +33,7 @@
 #include "AliGenHijingEventHeader.h"
 #include "AliGenGeVSimEventHeader.h"
 #include "AliMultiplicity.h"
-#include "AliFlowTrackSimpleCuts.h"
+#include "AliFlowTrackCuts.h"
 #include "AliFlowEventSimple.h"
 #include "AliFlowTrack.h"
 #include "AliFlowEvent.h"
@@ -580,3 +580,59 @@ AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
 
 }
 
+//-----------------------------------------------------------------------
+AliFlowEvent::AliFlowEvent( AliVEvent* event,
+                            AliFlowTrackCuts* rpCuts,
+                            AliFlowTrackCuts* poiCuts ):
+  AliFlowEventSimple(20)
+{
+  //Fills the event from a vevent: AliESDEvent,AliAODEvent,AliMCEvent
+
+  if (!rpCuts || !poiCuts) return;
+
+  //if input event empty try to do MC analysis
+  if (!event) event = rpCuts->GetMCevent();
+  if (!event) return;
+
+  Int_t numberOfTracks = event->GetNumberOfTracks() ;
+
+  //loop over tracks
+  for (Int_t i=0; i<numberOfTracks; i++)
+  {
+    AliVParticle* particle = event->GetTrack(i);   //get input particle
+
+    Bool_t rp = rpCuts->IsSelected(particle);
+    Bool_t poi = poiCuts->IsSelected(particle);
+    
+    if (!(rp||poi)) continue;
+
+    //make new AliFLowTrack
+    //here we need to be careful: if selected particle passes both rp and poi cuts
+    //then both cuts should have been done on the same set of parameters, e.g. global
+    //or TPConly. Otherwise we would have to introduce the same particle twice.
+    //this means that in a sane scenario when we pass both rp and poi cuts we get our
+    //parameters from any one of them (here rp).
+    AliFlowTrack* pTrack = NULL;
+    if (rp&&poi)
+    {
+      pTrack = rpCuts->MakeFlowTrack();
+      pTrack->TagRP();
+      pTrack->TagPOI();
+    }
+    else
+    if (rp)
+    {
+      pTrack = rpCuts->MakeFlowTrack();
+      pTrack->TagRP();
+    }
+    else
+    if (poi)
+    {
+      pTrack = poiCuts->MakeFlowTrack();
+      pTrack->TagPOI();
+    }
+
+    AddTrack(pTrack);
+  }//end of while (i < numberOfTracks)
+}
+
index 35e6b23..b6fd397 100644 (file)
 #ifndef ALIFLOWEVENT_H
 #define ALIFLOWEVENT_H
 
-class TTree;
-class AliFlowTrackSimpleCuts;
+class AliFlowTrackCuts;
 class AliFlowTrack;
 class AliCFManager;
+class AliVEvent;
 class AliMCEvent;
 class AliESDEvent;
-class AliMCEvent;
 class AliAODEvent;
 class AliMultiplicity;
 class TH2F;
@@ -57,6 +56,9 @@ public:
   AliFlowEvent( const AliESDEvent* anInput,
                 const TH2F* anInputFMDhist,
                 const AliCFManager* poiCFManager );
+  AliFlowEvent( AliVEvent* anInput,
+                AliFlowTrackCuts* rpCuts,
+                AliFlowTrackCuts* poiCuts );
 
 
   void SetMCReactionPlaneAngle(const AliMCEvent* mcEvent);
diff --git a/PWG2/FLOW/AliFlowTasks/AliFlowEventCuts.cxx b/PWG2/FLOW/AliFlowTasks/AliFlowEventCuts.cxx
new file mode 100644 (file)
index 0000000..8e0b562
--- /dev/null
@@ -0,0 +1,84 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id$ */ 
+
+// AliFlowEventCuts:
+// An event cut class for the flow framework
+//
+// origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
+
+#include <limits.h>
+#include <float.h>
+#include "TNamed.h"
+#include "AliVEvent.h"
+#include "AliFlowEventCuts.h"
+
+ClassImp(AliFlowEventCuts)
+
+//-----------------------------------------------------------------------
+AliFlowEventCuts::AliFlowEventCuts():
+  TNamed(),
+  fCutNumberOfTracks(kFALSE),
+  fNumberOfTracksMax(INT_MAX),
+  fNumberOfTracksMin(INT_MIN)
+{
+  //constructor 
+}
+
+////-----------------------------------------------------------------------
+//AliFlowEventCuts::AliFlowEventCuts(const AliFlowEventCuts& someCuts):
+//  TNamed(),
+//  fCutNumberOfTracks(that.fCutNumberOfTracks),
+//  fNumberOfTracksMax(that.fNumberOfTracksMax),
+//  fNumberOfTracksMin(that.fNumberOfTracksMin),
+//{
+//  //copy constructor 
+//}
+//
+////-----------------------------------------------------------------------
+//AliFlowEventCuts& AliFlowEventCuts::operator=(const AliFlowEventCuts& someCuts)
+//{
+//  //assignment
+//  fCutNumberOfTracks=that.fCutNumberOfTracks;
+//  fNumberOfTracksMax=that.fNumberOfTracksMax;
+//  fNumberOfTracksMin=that.fNumberOfTracksMin;
+//
+//  return *this;
+//}
+
+//----------------------------------------------------------------------- 
+Bool_t AliFlowEventCuts::IsSelected(const TObject* obj)
+{
+  //check cuts
+  const AliVEvent* vevent = dynamic_cast<const AliVEvent*>(obj);
+  if (vevent) return PassesCuts(vevent);
+  return kFALSE;  //when passed wrong type of object
+}
+//----------------------------------------------------------------------- 
+Bool_t AliFlowEventCuts::PassesCuts(const AliVEvent *event)
+{
+  ///check if event passes cuts
+  if(fCutNumberOfTracks) {if (event->GetNumberOfTracks() < fNumberOfTracksMin || event->GetNumberOfTracks() > fNumberOfTracksMax ) return kFALSE;}
+  return kTRUE;
+}
+
+//----------------------------------------------------------------------- 
+AliFlowEventCuts* AliFlowEventCuts::StandardCuts()
+{
+  //make a set of standard event cuts, caller becomes owner
+  AliFlowEventCuts* cuts = new AliFlowEventCuts();
+  return cuts;
+}
diff --git a/PWG2/FLOW/AliFlowTasks/AliFlowEventCuts.h b/PWG2/FLOW/AliFlowTasks/AliFlowEventCuts.h
new file mode 100644 (file)
index 0000000..87b5404
--- /dev/null
@@ -0,0 +1,47 @@
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. */
+/* See cxx source for full Copyright notice */
+/* $Id$ */
+
+// AliFlowEventCuts:
+// An event cut class
+// origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
+
+#ifndef ALIFLOWEVENTCUTS_H
+#define ALIFLOWEVENTCUTS_H
+
+#include <float.h>
+#include "TNamed.h"
+
+class AliVEvent;;
+
+class AliFlowEventCuts : public TNamed {
+
+ public:
+  AliFlowEventCuts();
+  //AliFlowEventCuts(const AliFlowEventCuts& someCuts);
+  //AliFlowEventCuts& operator=(const AliFlowEventCuts& someCuts);
+  virtual  ~AliFlowEventCuts() {}
+  
+  virtual Bool_t IsSelected(const TObject* obj);
+
+  Bool_t PassesCuts(const AliVEvent* event);
+  
+  static AliFlowEventCuts* StandardCuts();
+  
+  void SetNumberOfTracksMax(const Int_t value) {fNumberOfTracksMax=value;fCutNumberOfTracks=kTRUE;}
+  void SetNumberOfTracksMin(const Int_t value) {fNumberOfTracksMin=value;fCutNumberOfTracks=kTRUE;}
+
+  Int_t GetNumberOfTracksMax() const {return fNumberOfTracksMax;}
+  Int_t GetNumberOfTracksMin() const {return fNumberOfTracksMin;}
+
+ private:
+  Bool_t   fCutNumberOfTracks;//cut on # of tracks
+  Int_t fNumberOfTracksMax;  //limits
+  Int_t fNumberOfTracksMin;  //limits
+
+  ClassDef(AliFlowEventCuts,1)
+};
+
+#endif
+
+
index c78f7c3..ee2b8cc 100644 (file)
@@ -34,7 +34,7 @@ AliFlowTrack::AliFlowTrack():
 }
 
 //-----------------------------------------------------------------------
-AliFlowTrack::AliFlowTrack(AliVParticle* p):
+AliFlowTrack::AliFlowTrack(const AliVParticle* p):
   AliFlowTrackSimple(),
   fTrackSourceBits()
 {
index 87633c7..ee9f444 100644 (file)
@@ -21,7 +21,7 @@ public:
                      kFromTracklet=3,
                      kFromFMD=4 };
   AliFlowTrack();
-  AliFlowTrack(AliVParticle* p);
+  AliFlowTrack(const AliVParticle* p);
   AliFlowTrack& operator=(const AliFlowTrack& aTrack);
   //virtual AliFlowTrackSimple& operator=(const AliFlowTrackSimple& aTrack);
   AliFlowTrack(const AliFlowTrack& aTrack);
diff --git a/PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.cxx b/PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.cxx
new file mode 100644 (file)
index 0000000..e1a7417
--- /dev/null
@@ -0,0 +1,265 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id$ */ 
+
+// AliFlowTrackCuts:
+// ESD track cuts for flow framework 
+//
+// origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
+
+#include <limits.h>
+#include <float.h>
+#include "AliMCEvent.h"
+#include "AliVParticle.h"
+#include "AliMCParticle.h"
+#include "AliESDtrack.h"
+#include "AliAODTrack.h"
+#include "AliFlowTrack.h"
+#include "AliFlowTrackCuts.h"
+#include "AliLog.h"
+
+ClassImp(AliFlowTrackCuts)
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts::AliFlowTrackCuts():
+  AliFlowTrackSimpleCuts(),
+  fAliESDtrackCuts(new AliESDtrackCuts()),
+  fCutMCprocessType(kFALSE),
+  fMCprocessType(kPNoProcess),
+  fCutMCPID(kFALSE),
+  fMCPID(0),
+  fCutMCisPrimary(kFALSE),
+  fMCisPrimary(kFALSE),
+  fParamType(kESD_Global),
+  fParamMix(kPure),
+  fMCevent(NULL),
+  fCleanupTrack(kFALSE),
+  fTrack(NULL),
+  fMCparticle(NULL)
+{
+  //constructor 
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& someCuts):
+  AliFlowTrackSimpleCuts(someCuts),
+  fAliESDtrackCuts(new AliESDtrackCuts(*(someCuts.fAliESDtrackCuts))),
+  fCutMCprocessType(someCuts.fCutMCprocessType),
+  fMCprocessType(someCuts.fMCprocessType),
+  fCutMCPID(someCuts.fCutMCPID),
+  fMCPID(someCuts.fMCPID),
+  fCutMCisPrimary(someCuts.fCutMCisPrimary),
+  fMCisPrimary(someCuts.fMCisPrimary),
+  fParamType(someCuts.fParamType),
+  fParamMix(someCuts.fParamMix),
+  fMCevent(NULL),
+  fCleanupTrack(kFALSE),
+  fTrack(NULL),
+  fMCparticle(NULL)
+{
+  //copy constructor
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& someCuts)
+{
+  //assignment
+  AliFlowTrackSimpleCuts::operator=(someCuts);
+  *fAliESDtrackCuts=*(someCuts.fAliESDtrackCuts);
+  fCutMCprocessType=someCuts.fCutMCprocessType;
+  fMCprocessType=someCuts.fMCprocessType;
+  fCutMCPID=someCuts.fCutMCPID;
+  fMCPID=someCuts.fMCPID;
+  fCutMCisPrimary=someCuts.fCutMCisPrimary;
+  fMCisPrimary=someCuts.fMCisPrimary;
+  fParamType=someCuts.fParamType;
+  fParamMix=someCuts.fParamMix;
+
+  fMCevent=NULL;
+  fCleanupTrack=kFALSE;
+  fTrack=NULL;
+  fMCparticle=NULL;
+
+  return *this;
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts::~AliFlowTrackCuts()
+{
+  //dtor
+  if (fCleanupTrack) delete fTrack;
+  delete fAliESDtrackCuts;
+}
+
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::IsSelected(TObject* obj)
+{
+  //check cuts
+  AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj);
+  if (vparticle) return PassesCuts(vparticle);
+  AliFlowTrackSimple* flowtrack = dynamic_cast<AliFlowTrackSimple*>(obj);
+  if (flowtrack) return PassesCuts(flowtrack);
+  return kFALSE;  //default when passed wrong type of object
+}
+
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesCuts(AliFlowTrackSimple* track)
+{
+  //check cuts on a flowtracksimple
+  if (fCleanupTrack) delete fTrack;
+  fTrack = NULL;
+  return AliFlowTrackSimpleCuts::PassesCuts(track);
+}
+
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesCuts(AliVParticle* vparticle)
+{
+  //check cuts for an ESD vparticle
+
+  Int_t mcLabel = vparticle->GetLabel();
+  if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(mcLabel));
+  else fMCparticle=NULL;
+
+  AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*>(vparticle);
+  if (esdTrack)
+    HandleESDtrack(esdTrack);
+  else
+    HandleVParticle(vparticle);
+
+  //check the common cuts for the current particle (MC,AOD,ESD)
+  if (fCutPt) {if (fTrack->Pt() < fPtMin || fTrack->Pt() >= fPtMax ) return kFALSE;}
+  if (fCutEta) {if (fTrack->Eta() < fEtaMin || fTrack->Eta() >= fEtaMax ) return kFALSE;}
+  if (fCutPhi) {if (fTrack->Phi() < fPhiMin || fTrack->Phi() >= fPhiMax ) return kFALSE;}
+  if (fCutCharge) {if (fTrack->Charge() != fCharge) return kFALSE;}
+  //if(fCutPID) {if (fTrack->PID() != fPID) return kFALSE;}
+
+  //when MC info is required
+  if (fCutMCisPrimary)
+  {
+    if (!fMCevent) {AliError("no MC info"); return kFALSE;}
+    if (fMCevent->IsPhysicalPrimary(mcLabel) != fMCisPrimary) return kFALSE;
+  }
+  if (fCutMCPID)
+  {
+    if (!fMCparticle) {AliError("no MC info"); return kFALSE;}
+    Int_t pdgCode = fMCparticle->PdgCode();
+    if (fMCPID != pdgCode) return kFALSE;
+  }
+  if ( fCutMCprocessType )
+  {
+    if (!fMCparticle) {AliError("no MC info"); return kFALSE;}
+    TParticle* particle = fMCparticle->Particle();
+    Int_t processID = particle->GetUniqueID();
+    if (processID != fMCprocessType ) return kFALSE;
+  }
+
+  //check all else for ESDs using aliesdtrackcuts
+  if (esdTrack && (fParamType!=kMC) ) return fAliESDtrackCuts->IsSelected(static_cast<AliESDtrack*>(fTrack));
+
+  return kTRUE; //true by default, if we didn't set any cuts
+}
+
+//-----------------------------------------------------------------------
+void AliFlowTrackCuts::HandleVParticle(AliVParticle* track)
+{
+  //handle the general case
+  if (fCleanupTrack) delete fTrack;
+  switch (fParamType)
+  {
+    case kMC:
+      fCleanupTrack = kFALSE;
+      fTrack = fMCparticle;
+      break;
+    default:
+      fCleanupTrack = kFALSE;
+      fTrack = track;
+  }
+}
+
+//-----------------------------------------------------------------------
+void AliFlowTrackCuts::HandleESDtrack(AliESDtrack* track)
+{
+  //handle esd track
+  if (fCleanupTrack) delete fTrack;
+  switch (fParamType)
+  {
+    case kESD_Global:
+      fTrack = track;
+      fCleanupTrack = kFALSE;
+      break;
+    case kESD_TPConly:
+      fTrack = new AliESDtrack();
+      track->FillTPCOnlyTrack(*(static_cast<AliESDtrack*>(fTrack)));
+      fCleanupTrack = kTRUE;
+      break;
+    case kMC:
+      fCleanupTrack = kFALSE;
+      fTrack = fMCparticle;
+      break;
+    default:
+      fTrack = track;
+      fCleanupTrack = kFALSE;
+  }
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCOnlyTrackCuts()
+{
+  //get standard cuts
+  AliFlowTrackCuts* cuts = new AliFlowTrackCuts();
+  cuts->SetName("standard global track cuts 2009");
+  delete cuts->fAliESDtrackCuts;
+  cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+  return cuts;
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts* AliFlowTrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
+{
+  //get standard cuts
+  AliFlowTrackCuts* cuts = new AliFlowTrackCuts();
+  cuts->SetName("standard global track cuts 2009");
+  delete cuts->fAliESDtrackCuts;
+  cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selPrimaries);
+  return cuts;
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrack* AliFlowTrackCuts::MakeFlowTrack() const
+{
+  //get a flow track constructed from whatever we applied cuts on
+  //caller is resposible for deletion
+  AliFlowTrack* flowtrack=NULL;
+  switch(fParamMix)
+  {
+    case kPure:
+      flowtrack = new AliFlowTrack(fTrack);
+      break;
+    case kTrackWithMCkine:
+      flowtrack = new AliFlowTrack(fMCparticle);
+      break;
+    case kTrackWithMCPID:
+      flowtrack = new AliFlowTrack(fTrack);
+      break;
+    default:
+      flowtrack = new AliFlowTrack(fTrack);
+  }
+  if (fParamType==kMC) flowtrack->SetSource(AliFlowTrack::kFromMC);
+  else if (dynamic_cast<AliMCParticle*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromMC);
+  else if (dynamic_cast<AliESDtrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromESD);
+  else if (dynamic_cast<AliAODTrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromAOD);
+  return flowtrack;
+}
diff --git a/PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.h b/PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.h
new file mode 100644 (file)
index 0000000..5c4bc13
--- /dev/null
@@ -0,0 +1,122 @@
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. */
+/* See cxx source for full Copyright notice */
+/* $Id$ */
+
+// AliFlowTrackESDCuts:
+// A cut class for ESD, AOD and MC particles for the flow framework
+// author: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
+
+#ifndef ALIFLOWTRACKCUTS_H
+#define ALIFLOWTRACKCUTS_H
+
+#include "AliFlowTrackSimpleCuts.h"
+#include "AliESDtrackCuts.h"
+#include "TMCProcess.h"
+
+class AliVParticle;
+class AliMCParticle;
+class AliFlowTrack;
+class AliMCEvent;
+class AliVEvent;
+
+class AliFlowTrackCuts : public AliFlowTrackSimpleCuts {
+
+ public:
+  AliFlowTrackCuts();
+  AliFlowTrackCuts(const AliFlowTrackCuts& someCuts);
+  AliFlowTrackCuts& operator=(const AliFlowTrackCuts& someCuts);
+  virtual ~AliFlowTrackCuts();
+
+  static AliFlowTrackCuts* GetStandardTPCOnlyTrackCuts();
+  static AliFlowTrackCuts* GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries=kTRUE);
+
+  enum trackParameterType { kMC, kESD_Global, kESD_TPConly };
+  enum trackParameterMix  { kPure, kTrackWithMCkine, kTrackWithMCPID };
+
+  //setters (interface to AliESDtrackCuts)
+  void SetMinNClustersTPC( Int_t a ) {fAliESDtrackCuts->SetMinNClustersTPC(a);}
+  void SetMinNClustersITS( Int_t a ) {fAliESDtrackCuts->SetMinNClustersITS(a);}
+  void SetClusterRequirementITS( AliESDtrackCuts::Detector det,
+                                 AliESDtrackCuts::ITSClusterRequirement req = AliESDtrackCuts::kOff )
+                                 { fAliESDtrackCuts->SetClusterRequirementITS(det,req); } 
+  void SetMaxChi2PerClusterTPC( Float_t a ) {fAliESDtrackCuts->SetMaxChi2PerClusterTPC(a);}
+  void SetMaxChi2PerClusterITS( Float_t a ) {fAliESDtrackCuts->SetMaxChi2PerClusterITS(a);}
+  void SetRequireTPCRefit( Bool_t a ) {fAliESDtrackCuts->SetRequireTPCRefit(a);}
+  void SetRequireTPCStandAlone( Bool_t a) {fAliESDtrackCuts->SetRequireTPCStandAlone(a);}
+  void SetRequireITSRefit( Bool_t a ) {fAliESDtrackCuts->SetRequireITSRefit(a);}
+  void SetRequireITSStandAlone( Bool_t a) {fAliESDtrackCuts->SetRequireITSStandAlone(a);}
+  void SetAcceptKinkDaughters( Bool_t a ) {fAliESDtrackCuts->SetAcceptKinkDaughters(a);}
+  void SetMaxDCAToVertexZ( Float_t a ) {fAliESDtrackCuts->SetMaxDCAToVertexZ(a);}
+  void SetMaxDCAToVertexXY( Float_t a ) {fAliESDtrackCuts->SetMaxDCAToVertexXY(a);}
+  void SetMaxDCAToVertexXYPtDep( const char* a ) {fAliESDtrackCuts->SetMaxDCAToVertexXYPtDep(a);}
+  void SetRequireSigmaToVertex(Bool_t a) {fAliESDtrackCuts->SetRequireSigmaToVertex(a);}
+  void SetMaxNsigmaToVertex(Float_t sigma=1e10) { fAliESDtrackCuts->SetMaxNsigmaToVertex(sigma); }
+  void SetDCAToVertex2D( Bool_t a ) {fAliESDtrackCuts->SetDCAToVertex2D(a);}
+  void SetEtaRange( Float_t r1, Float_t r2 ) { SetEtaMin(r1); SetEtaMax(r2); }
+  void SetPtRange( Float_t r1, Float_t r2 ) { SetPtMin(r1); SetPtMax(r2); }
+
+  Int_t GetMinNClustersTPC() const {return fAliESDtrackCuts->GetMinNClusterTPC();}
+  Int_t GetMinNClustersITS() const {return fAliESDtrackCuts->GetMinNClustersITS();}
+  AliESDtrackCuts::ITSClusterRequirement GetClusterRequirementITS( AliESDtrackCuts::Detector det ) const
+                                 { return fAliESDtrackCuts->GetClusterRequirementITS(det); } 
+  Float_t GetMaxChi2PerClusterTPC() const {return fAliESDtrackCuts->GetMaxChi2PerClusterTPC();}
+  Float_t GetMaxChi2PerClusterITS() const {return fAliESDtrackCuts->GetMaxChi2PerClusterITS();}
+  Bool_t GetRequireTPCRefit() const {return fAliESDtrackCuts->GetRequireTPCRefit();}
+  Bool_t GetRequireTPCStandAlone() const {return fAliESDtrackCuts->GetRequireTPCStandAlone();}
+  Bool_t GetRequireITSRefit() const {return fAliESDtrackCuts->GetRequireITSRefit();}
+  Bool_t GetRequireITSStandAlone() const {return fAliESDtrackCuts->GetRequireITSStandAlone();}
+  Bool_t GetAcceptKinkDaughters() const {return fAliESDtrackCuts->GetAcceptKinkDaughters();}
+  Float_t GetMaxDCAToVertexZ() const {return fAliESDtrackCuts->GetMaxDCAToVertexZ();}
+  Float_t GetMaxDCAToVertexXY() const {return fAliESDtrackCuts->GetMaxDCAToVertexXY();}
+  const char* GetMaxDCAToVertexXYPtDep() const {return fAliESDtrackCuts->GetMaxDCAToVertexXYPtDep();}
+  Bool_t GetRequireSigmaToVertex() const {return fAliESDtrackCuts->GetRequireSigmaToVertex();}
+  Float_t GetMaxNsigmaToVertex() const {return fAliESDtrackCuts->GetMaxNsigmaToVertex(); }
+  Bool_t GetDCAToVertex2D() const {return fAliESDtrackCuts->GetDCAToVertex2D();}
+  void GetEtaRange( Float_t& r1, Float_t& r2 ) const { r1=GetEtaMin(); r2=GetEtaMax(); }
+  void GetPtRange( Float_t& r1, Float_t& r2 ) const { r1=GetPtMin(); r2=GetPtMax(); }
+
+  //MC stuff
+  void SetMCprocessType( TMCProcess t ) { fMCprocessType = t; fCutMCprocessType=kTRUE; }
+  TMCProcess GetMCprocessType() const { return fMCprocessType; }
+  void SetMCisPrimary( Bool_t b ) { fMCisPrimary=b; fCutMCisPrimary=kTRUE; }
+  Bool_t GetMCisPrimary() const {return fMCisPrimary;}
+
+  void SetParamType(trackParameterType paramType) {fParamType=paramType;}
+  trackParameterType GetParamType() const {return fParamType;}
+  void SetParamMix(trackParameterMix paramMix) {fParamMix=paramMix;}
+  trackParameterMix GetParamMix() const {return fParamMix;}
+
+  virtual Bool_t IsSelected(TObject* obj);
+  const AliVParticle* GetTrack() const {return fTrack;}
+  AliFlowTrack* MakeFlowTrack() const;
+  
+  void SetMCevent(AliMCEvent* mcEvent) {fMCevent=mcEvent;}
+  AliMCEvent* GetMCevent() const {return fMCevent;}
+
+ protected:
+  Bool_t PassesCuts(AliVParticle* track);
+  Bool_t PassesCuts(AliFlowTrackSimple* track);
+  void HandleESDtrack(AliESDtrack* track);
+  void HandleVParticle(AliVParticle* track);
+
+  AliESDtrackCuts* fAliESDtrackCuts; //alianalysis cuts
+  Bool_t fCutMCprocessType;          //do we cut on mc process type?
+  TMCProcess fMCprocessType;         //mc process type
+  Bool_t fCutMCPID;                  //cut on MC pid?
+  Int_t fMCPID;                      //MC PID
+  Bool_t fCutMCisPrimary;            //do we cut on primaryness?
+  Bool_t fMCisPrimary;               //is MC primary
+
+  trackParameterType fParamType;     //parameter type tu cut on
+  trackParameterMix fParamMix;       //parameter mixing
+  AliMCEvent* fMCevent;              //!mc event
+  Bool_t fCleanupTrack;              //check if we need to delete
+  AliVParticle* fTrack;              //!the track to apply cuts on
+  AliMCParticle* fMCparticle;        //!mc particle
+
+  ClassDef(AliFlowTrackCuts,1)
+};
+
+#endif
+
+
index 1e70008..65c4c8f 100644 (file)
Binary files a/PWG2/FLOW/Documentation/FLOW.pdf and b/PWG2/FLOW/Documentation/FLOW.pdf differ
diff --git a/PWG2/FLOW/macros/AddTaskFlowCentrality.C b/PWG2/FLOW/macros/AddTaskFlowCentrality.C
new file mode 100644 (file)
index 0000000..1b02a71
--- /dev/null
@@ -0,0 +1,1003 @@
+/////////////////////////////////////////////////////////////////////////////////////////////
+//
+// AddTask* macro for flow analysis
+// Creates a Flow Event task and adds it to the analysis manager.
+// Sets the cuts using the correction framework (CORRFW) classes.
+// Also creates Flow Analysis tasks and connects them to the output of the flow event task.
+//
+/////////////////////////////////////////////////////////////////////////////////////////////
+
+// Define the range for eta subevents (for SP method)
+//-----(FMD 1.7 - 5.0)-----
+//Double_t minA = -5.0;
+//Double_t maxA = -1.7;
+//Double_t minB = 1.7;
+//Double_t maxB = 5.0;
+//-----(Tracklets 0.9 - 2.0)-----
+//Double_t minA = -2.0;
+//Double_t maxA = -0.9;
+//Double_t minB = 0.9;
+//Double_t maxB = 2.0;
+//-----(Global 0.5 - 0.9)-----
+Double_t minA = -0.9;
+Double_t maxA = -0.5;
+Double_t minB = 0.5;
+Double_t maxB = 0.9;
+
+// AFTERBURNER
+Bool_t useAfterBurner=kFALSE;
+Double_t v1=0.0;
+Double_t v2=0.0;
+Double_t v3=0.0;
+Double_t v4=0.0;
+Int_t numberOfTrackClones=0; //non-flow
+
+// Define a range of the detector to exclude
+Bool_t ExcludeRegion = kFALSE;
+Double_t excludeEtaMin = -0.;
+Double_t excludeEtaMax = 0.;
+Double_t excludePhiMin = 0.;
+Double_t excludePhiMax = 0.;
+
+// use physics selection class
+Bool_t  UsePhysicsSelection = kTRUE;
+
+// SETTING THE CUTS
+
+//----------Event selection----------
+Bool_t UseMultCutforESD = kTRUE;
+//Bool_t UseMultCutforESD = kFALSE;
+const Int_t multminESD = 1;  //used for CORRFW cuts 
+const Int_t multmaxESD = 1000000; //used for CORRFW cuts 
+
+Bool_t requireVtxCuts = kTRUE;
+const Double_t vertexXmin = -1.e99; 
+const Double_t vertexXmax = 1.e99;
+const Double_t vertexYmin = -1.e99;
+const Double_t vertexYmax = 1.e99;
+const Double_t vertexZmin = -10.; 
+const Double_t vertexZmax = 10.; 
+const Int_t vertexNContributorsmin = 1;
+const Int_t vertexNContributorsmax = 10000;
+
+//Bool_t UseMultCut = kFALSE;
+Bool_t UseMultCut = kTRUE;
+const Int_t multmin = 1;     //used for AliFlowEventSimple (to set the centrality)
+const Int_t multmax = 10000;     //used for AliFlowEventSimple (to set the centrality)
+//const Int_t multmin = 10;     //used for AliFlowEventSimple (to set the centrality)
+//const Int_t multmax = 1000000;     //used for AliFlowEventSimple (to set the centrality)
+
+
+//----------For RP selection----------
+// Use Global tracks ("Global"), or SPD tracklets ("Tracklet") 
+// or FMD hits ("FMD") for the RP selection
+const TString rptype = "Global";
+//const TString rptype = "Tracklet";
+//const TString rptype = "FMD";
+
+//KINEMATICS (on generated and reconstructed tracks)
+Bool_t UseKineforRP =  kFALSE;
+const Double_t ptminRP = 0.0;
+const Double_t ptmaxRP = 10.0;
+const Double_t etaminRP  = -0.9;
+const Double_t etamaxRP  = 0.9;
+const Int_t    chargeRP = 1;  //not used
+const Bool_t   isChargedRP = kTRUE;
+
+//PID (on generated and reconstructed tracks)
+Bool_t UsePIDforRP = kFALSE;
+const Int_t PdgRP = 211;
+
+//TRACK QUALITY (on reconstructed tracks only)
+//see /CORRFW/AliCFTrackQualityCuts class
+Bool_t UseTrackQualityforRP =  kFALSE;
+const Int_t    minClustersTpcRP = 80;           //default = -1; 
+const Double_t maxChi2PerClusterTpcRP = 4.0;    //default = 1.e+09;
+const UShort_t minDedxClusterTpcRP = 0;
+const Int_t    minClustersItsRP = 2;            //panos
+const Double_t maxChi2PerClusterItsRP = 1.e+09; 
+const Int_t    minClustersTrdRP = -1;
+const Int_t    minTrackletTrdRP = -1;
+const Int_t    minTrackletTrdPidRP = -1;
+const Double_t maxChi2PerClusterTrdRP = 1.e+09;
+const ULong_t  statusRP = AliESDtrack::kTPCrefit;   //AliESDtrack::kTPCrefit &  AliESDtrack::kITSrefit 
+
+//PRIMARY (on reconstructed tracks only)
+//see /CORRFW/AliCFTrackIsPrimaryCuts class
+Bool_t UsePrimariesforRP = kFALSE;
+const Bool_t   spdVertexRP = kFALSE;
+const Bool_t   tpcVertexRP = kFALSE;
+const Float_t  minDcaToVertexXyRP = 0.;
+const Float_t  minDcaToVertexZRP = 0.;
+const Float_t  maxDcaToVertexXyRP = 2.4;         //default = 1.e+10;  //2.4;
+const Float_t  maxDcaToVertexZRP = 3.2;          //default = 1.e+10;  //3.2;
+const Bool_t   dcaToVertex2dRP = kFALSE;         //default = kFALSE;
+const Bool_t   absDcaToVertexRP = kTRUE;         //default = kTRUE;
+const Double_t minNSigmaToVertexRP = 0.;
+const Double_t maxNSigmaToVertexRP = 1.e+10; //3.; //1.e+10
+const Double_t maxSigmaDcaXySP = 1.e+10;
+const Double_t maxSigmaDcaZSP = 1.e+10;
+const Bool_t   requireSigmaToVertexSP = kFALSE;
+const Bool_t   acceptKinkDaughtersSP = kFALSE;  //default = kTRUE;
+
+//ACCEPTANCE (on generated tracks only : AliMCParticle)
+//see /CORRFW/AliCFAcceptanceCuts class
+Bool_t UseAcceptanceforRP =  kFALSE; 
+const Int_t  minTrackrefsItsRP = 0;//3;
+const Int_t  minTrackrefsTpcRP = 0;//2;
+const Int_t  minTrackrefsTrdRP = 0; 
+const Int_t  minTrackrefsTofRP = 0; 
+const Int_t  minTrackrefsMuonRP = 0; 
+//default for all is 0
+
+//----------For POI selection----------
+//KINEMATICS (on generated and reconstructed tracks)
+Bool_t UseKineforPOI = kTRUE;
+const Double_t ptminPOI = 0.0;
+const Double_t ptmaxPOI = 10.0;
+const Double_t etaminPOI  = -0.5;
+const Double_t etamaxPOI  = 0.5;
+const Int_t    chargePOI = 1;  //not used
+const Bool_t   isChargedPOI = kTRUE;
+
+//PID (on generated and reconstructed tracks)
+Bool_t UsePIDforPOI = kFALSE;
+const Int_t PdgPOI = 321;
+
+//TRACK QUALITY (on reconstructed tracks only)
+//see /CORRFW/AliCFTrackQualityCuts class
+Bool_t UseTrackQualityforPOI = kTRUE;
+const Int_t    minClustersTpcPOI = 80;
+const Double_t maxChi2PerClusterTpcPOI = 4.0;    
+const UShort_t minDedxClusterTpcPOI = 0;
+const Int_t    minClustersItsPOI = 2;
+const Double_t maxChi2PerClusterItsPOI = 1.e+09;
+const Int_t    minClustersTrdPOI = -1;
+const Int_t    minTrackletTrdPOI = -1;
+const Int_t    minTrackletTrdPidPOI = -1;
+const Double_t maxChi2PerClusterTrdPOI = 1.e+09;
+const ULong_t  statusPOI = AliESDtrack::kTPCrefit;   
+
+//PRIMARY (on reconstructed tracks only)
+//see /CORRFW/AliCFTrackIsPrimaryCuts class
+Bool_t UsePrimariesforPOI = kTRUE;
+const Bool_t   spdVertexPOI = kFALSE;
+const Bool_t   tpcVertexPOI = kFALSE;
+const Float_t  minDcaToVertexXyPOI = 0.;
+const Float_t  minDcaToVertexZPOI = 0.;
+const Float_t  maxDcaToVertexXyPOI = 2.4;
+const Float_t  maxDcaToVertexZPOI = 3.2;
+const Bool_t   dcaToVertex2dPOI =  kFALSE;
+const Bool_t   absDcaToVertexPOI = kTRUE;
+const Double_t minNSigmaToVertexPOI = 0.;
+const Double_t maxNSigmaToVertexPOI = 1.e+10;  
+const Double_t maxSigmaDcaXyPOI = 1.e+10;
+const Double_t maxSigmaDcaZPOI = 1.e+10;
+const Bool_t   requireSigmaToVertexPOI = kFALSE;
+const Bool_t   acceptKinkDaughtersPOI = kFALSE;
+
+//ACCEPTANCE (on generated tracks only : AliMCParticle)
+//see /CORRFW/AliCFAcceptanceCuts class
+Bool_t UseAcceptanceforPOI = kFALSE;
+const Int_t minTrackrefsItsPOI = 3;
+const Int_t minTrackrefsTpcPOI = 2;
+const Int_t minTrackrefsTrdPOI = 0; 
+const Int_t minTrackrefsTofPOI = 0; 
+const Int_t minTrackrefsMuonPOI = 0; 
+
+void AddTaskFlowCentrality( TString type,
+                                                 Bool_t* METHODS,
+                                                 Bool_t QA,
+                                                 Bool_t* WEIGHTS,
+                                                 Int_t refMultMin=0,
+                                                 Int_t refMultMax=1e10,
+                                                 TString fileName="AnalysisResults.root" )
+{
+  //boleans for the methods
+  Bool_t SP       = METHODS[0];
+  Bool_t LYZ1SUM  = METHODS[1];
+  Bool_t LYZ1PROD = METHODS[2];
+  Bool_t LYZ2SUM  = METHODS[3];
+  Bool_t LYZ2PROD = METHODS[4];
+  Bool_t LYZEP    = METHODS[5];
+  Bool_t GFC      = METHODS[6];
+  Bool_t QC       = METHODS[7];
+  Bool_t FQD      = METHODS[8];
+  Bool_t MCEP     = METHODS[9];      
+  Bool_t MH       = METHODS[10];
+  Bool_t NL       = METHODS[11];  
+  //for using weights
+  Bool_t useWeights  = WEIGHTS[0] || WEIGHTS[1] || WEIGHTS[2];
+  if (useWeights) cout<<"Weights are used"<<endl;
+  else cout<<"Weights are not used"<<endl;
+
+
+  // Get the pointer to the existing analysis manager via the static access method.
+  //==============================================================================
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr) {
+    Error("AddTaskFlowEvent", "No analysis manager to connect to.");
+    return NULL;
+  }   
+  
+  // Check the analysis type using the event handlers connected to the analysis
+  // manager. The availability of MC handler cann also be checked here.
+  //==============================================================================
+  if (!mgr->GetInputEventHandler()) {
+    ::Error("AddTaskFlowEvent", "This task requires an input event handler");
+    return NULL;
+  }  
+    
+  // Open external input files
+  //===========================================================================
+  //weights: 
+  TFile *weightsFile = NULL;
+  TList *weightsList = NULL;
+
+  if(useWeights) {
+    //open the file with the weights:
+    weightsFile = TFile::Open("weights.root","READ");
+    if(weightsFile) {
+      //access the list which holds the histos with weigths:
+      weightsList = (TList*)weightsFile->Get("weights");
+    }
+    else {
+      cout<<" WARNING: the file <weights.root> with weights from the previous run was not available."<<endl;
+      break;
+    } 
+  }
+    
+  //LYZ2
+  if (LYZ2SUM || LYZ2PROD) {
+    //read the outputfile of the first run
+    TString outputFileName = "AnalysisResults1.root";
+    TString pwd(gSystem->pwd());
+    pwd+="/";
+    pwd+=outputFileName.Data();
+    TFile *outputFile = NULL;
+    if(gSystem->AccessPathName(pwd.Data(),kFileExists)) {
+      cout<<"WARNING: You do not have an output file:"<<endl;
+      cout<<"         "<<pwd.Data()<<endl;
+      exit(0);
+    } else {
+      outputFile = TFile::Open(pwd.Data(),"READ");
+    }
+    
+    if (LYZ2SUM){  
+      // read the output directory from LYZ1SUM 
+      TString inputFileNameLYZ2SUM = "outputLYZ1SUManalysis" ;
+      inputFileNameLYZ2SUM += type;
+      cout<<"The input directory is "<<inputFileNameLYZ2SUM.Data()<<endl;
+      TFile* fInputFileLYZ2SUM = (TFile*)outputFile->FindObjectAny(inputFileNameLYZ2SUM.Data());
+      if(!fInputFileLYZ2SUM || fInputFileLYZ2SUM->IsZombie()) { 
+       cerr << " ERROR: To run LYZ2SUM you need the output file from LYZ1SUM. This file is not there! Please run LYZ1SUM first." << endl ; 
+       break;
+      }
+      else {
+       TList* fInputListLYZ2SUM = (TList*)fInputFileLYZ2SUM->Get("cobjLYZ1SUM");
+       if (!fInputListLYZ2SUM) {cout<<"list is NULL pointer!"<<endl;}
+      }
+      cout<<"LYZ2SUM input file/list read..."<<endl;
+    }
+
+    if (LYZ2PROD){  
+      // read the output directory from LYZ1PROD 
+      TString inputFileNameLYZ2PROD = "outputLYZ1PRODanalysis" ;
+      inputFileNameLYZ2PROD += type;
+      cout<<"The input directory is "<<inputFileNameLYZ2PROD.Data()<<endl;
+      TFile* fInputFileLYZ2PROD = (TFile*)outputFile->FindObjectAny(inputFileNameLYZ2PROD.Data());
+      if(!fInputFileLYZ2PROD || fInputFileLYZ2PROD->IsZombie()) { 
+       cerr << " ERROR: To run LYZ2PROD you need the output file from LYZ1PROD. This file is not there! Please run LYZ1PROD first." << endl ; 
+       break;
+      }
+      else {
+       TList* fInputListLYZ2PROD = (TList*)fInputFileLYZ2PROD->Get("cobjLYZ1PROD");
+       if (!fInputListLYZ2PROD) {cout<<"list is NULL pointer!"<<endl;}
+      }
+      cout<<"LYZ2PROD input file/list read..."<<endl;
+    }
+  }
+
+
+  if (LYZEP) {
+    //read the outputfile of the second run
+    TString outputFileName = "AnalysisResults2.root";
+    TString pwd(gSystem->pwd());
+    pwd+="/";
+    pwd+=outputFileName.Data();
+    TFile *outputFile = NULL;
+    if(gSystem->AccessPathName(pwd.Data(),kFileExists)) {
+      cout<<"WARNING: You do not have an output file:"<<endl;
+      cout<<"         "<<pwd.Data()<<endl;
+      exit(0);
+    } else {
+      outputFile = TFile::Open(pwd.Data(),"READ");
+    }
+
+    // read the output file from LYZ2SUM
+    TString inputFileNameLYZEP = "outputLYZ2SUManalysis" ;
+    inputFileNameLYZEP += type;
+    cout<<"The input file is "<<inputFileNameLYZEP.Data()<<endl;
+    TFile* fInputFileLYZEP = (TFile*)outputFile->FindObjectAny(inputFileNameLYZEP.Data());
+    if(!fInputFileLYZEP || fInputFileLYZEP->IsZombie()) { 
+      cerr << " ERROR: To run LYZEP you need the output file from LYZ2SUM. This file is not there! Please run LYZ2SUM first." << endl ; 
+      break;
+    }
+    else {
+      TList* fInputListLYZEP = (TList*)fInputFileLYZEP->Get("cobjLYZ2SUM");
+      if (!fInputListLYZEP) {cout<<"list is NULL pointer!"<<endl;}
+    }
+    cout<<"LYZEP input file/list read..."<<endl;
+  }
+  
+  
+  // Create the FMD task and add it to the manager
+  //===========================================================================
+  if (rptype == "FMD") {
+    AliFMDAnalysisTaskSE *taskfmd = NULL;
+    if (rptype == "FMD") {
+      taskfmd = new AliFMDAnalysisTaskSE("TaskFMD");
+      mgr->AddTask(taskfmd);
+  
+      AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
+      pars->Init();
+      pars->SetProcessPrimary(kTRUE); //for MC only
+      pars->SetProcessHits(kFALSE);
+
+      //pars->SetRealData(kTRUE); //for real data
+      //pars->SetProcessPrimary(kFALSE); //for real data
+
+    }
+  }
+  
+
+  // Create the task, add it to the manager.
+  //===========================================================================
+  AliAnalysisTaskFlowEvent *taskFE = NULL;
+  if (QA) { 
+    if(useAfterBurner)
+    { 
+      taskFE = new AliAnalysisTaskFlowEvent("TaskFlowEvent",rptype,kTRUE,1);
+      taskFE->SetFlow(v1,v2,v3,v4); 
+      taskFE->SetNonFlowNumberOfTrackClones(numberOfTrackClones);
+      taskFE->SetAfterburnerOn();
+    }
+    else {taskFE = new AliAnalysisTaskFlowEvent("TaskFlowEvent",rptype,kTRUE); }
+    taskFE->SetAnalysisType(type);
+    taskFE->SetRPType(rptype); //only for ESD
+    if (UseMultCut) {
+      taskFE->SetMinMult(multmin);
+      taskFE->SetMaxMult(multmax);
+    }
+    if (ExcludeRegion) {
+      taskFE->DefineDeadZone(excludeEtaMin, excludeEtaMax, excludePhiMin, excludePhiMax); 
+    }
+
+    taskFE->SetSubeventEtaRange(minA, maxA, minB, maxB);
+    if (UsePhysicsSelection) {
+      taskFE->SelectCollisionCandidates();
+      cout<<"Using Physics Selection"<<endl;
+    }
+    mgr->AddTask(taskFE);
+  }
+  else { 
+    if(useAfterBurner)
+    { 
+      taskFE = new AliAnalysisTaskFlowEvent("TaskFlowEvent",rptype,kFALSE,1);
+      taskFE->SetFlow(v1,v2,v3,v4); 
+      taskFE->SetNonFlowNumberOfTrackClones(numberOfTrackClones);
+      taskFE->SetAfterburnerOn();
+    }
+    else {taskFE = new AliAnalysisTaskFlowEvent("TaskFlowEvent",rptype,kFALSE); }
+    taskFE->SetAnalysisType(type);
+    taskFE->SetRPType(rptype); //only for ESD
+    if (UseMultCut) {
+      taskFE->SetMinMult(multmin);
+      taskFE->SetMaxMult(multmax);
+    }
+    if (ExcludeRegion) {
+      taskFE->DefineDeadZone(excludeEtaMin, excludeEtaMax, excludePhiMin, excludePhiMax); 
+    }
+    taskFE->SetSubeventEtaRange(minA, maxA, minB, maxB);
+    if (UsePhysicsSelection) {
+      taskFE->SelectCollisionCandidates();
+      cout<<"Using Physics Selection"<<endl;
+    }
+    mgr->AddTask(taskFE);
+  }
+  // Create cuts using the correction framework (CORRFW)
+  //===========================================================================
+  if (QA){
+    //Set TList for the QA histograms
+    TList* qaRP  = new TList(); 
+    TList* qaPOI = new TList();
+  }
+
+  //----------Event cuts----------
+  AliCFEventGenCuts* mcEventCuts = new AliCFEventGenCuts("mcEventCuts","MC-level event cuts");
+  mcEventCuts->SetNTracksCut(refMultMin,refMultMax); 
+  mcEventCuts->SetRequireVtxCuts(requireVtxCuts);
+  mcEventCuts->SetVertexXCut(vertexXmin, vertexXmax);
+  mcEventCuts->SetVertexYCut(vertexYmin, vertexYmax);
+  mcEventCuts->SetVertexZCut(vertexZmin, vertexZmax);
+  if (QA) { 
+    mcEventCuts->SetQAOn(qaRP);
+  }
+  AliCFEventRecCuts* recEventCuts = new AliCFEventRecCuts("recEventCuts","rec-level event cuts");
+  recEventCuts->SetNTracksCut(refMultMin,refMultMax); 
+  recEventCuts->SetRequireVtxCuts(requireVtxCuts);
+  recEventCuts->SetVertexXCut(vertexXmin, vertexXmax);
+  recEventCuts->SetVertexYCut(vertexYmin, vertexYmax);
+  recEventCuts->SetVertexZCut(vertexZmin, vertexZmax);
+  recEventCuts->SetVertexNContributors(vertexNContributorsmin,vertexNContributorsmax);
+  if (QA) { 
+    recEventCuts->SetQAOn(qaRP);
+  }
+  
+  //----------Cuts for RP----------
+  //KINEMATICS (MC and reconstructed)
+  AliCFTrackKineCuts* mcKineCutsRP = new AliCFTrackKineCuts("mcKineCutsRP","MC-level kinematic cuts");
+  mcKineCutsRP->SetPtRange(ptminRP,ptmaxRP);
+  mcKineCutsRP->SetEtaRange(etaminRP,etamaxRP);
+  //mcKineCutsRP->SetChargeMC(chargeRP);
+  mcKineCutsRP->SetRequireIsCharged(isChargedRP);
+  if (QA) { 
+    mcKineCutsRP->SetQAOn(qaRP);
+  }
+
+  AliCFTrackKineCuts *recKineCutsRP = new AliCFTrackKineCuts("recKineCutsRP","rec-level kine cuts");
+  recKineCutsRP->SetPtRange(ptminRP,ptmaxRP);
+  recKineCutsRP->SetEtaRange(etaminRP,etamaxRP);
+  //recKineCutsRP->SetChargeRec(chargeRP);
+  recKineCutsRP->SetRequireIsCharged(isChargedRP);
+  if (QA) { 
+    recKineCutsRP->SetQAOn(qaRP);
+  }
+
+  //PID (MC and reconstructed)
+  AliCFParticleGenCuts* mcGenCutsRP = new AliCFParticleGenCuts("mcGenCutsRP","MC particle generation cuts for RP");
+  mcGenCutsRP->SetRequireIsPrimary();
+  if (UsePIDforRP) {mcGenCutsRP->SetRequirePdgCode(PdgRP);}
+  if (QA) { 
+    mcGenCutsRP->SetQAOn(qaRP);
+  }
+
+  int n_species = AliPID::kSPECIES ;
+  Double_t* prior = new Double_t[n_species];
+  
+  prior[0] = 0.0244519 ;
+  prior[1] = 0.0143988 ;
+  prior[2] = 0.805747  ;
+  prior[3] = 0.0928785 ;
+  prior[4] = 0.0625243 ;
+  
+  AliCFTrackCutPid* cutPidRP = NULL;
+  if(UsePIDforRP) {
+    cutPidRP = new AliCFTrackCutPid("cutPidRP","ESD_PID for RP") ;
+    cutPidRP->SetPriors(prior);
+    cutPidRP->SetProbabilityCut(0.0);
+    cutPidRP->SetDetectors("TPC TOF");
+    switch(TMath::Abs(PDG1)) {
+    case 11   : cutPidRP->SetParticleType(AliPID::kElectron, kTRUE); break;
+    case 13   : cutPidRP->SetParticleType(AliPID::kMuon    , kTRUE); break;
+    case 211  : cutPidRP->SetParticleType(AliPID::kPion    , kTRUE); break;
+    case 321  : cutPidRP->SetParticleType(AliPID::kKaon    , kTRUE); break;
+    case 2212 : cutPidRP->SetParticleType(AliPID::kProton  , kTRUE); break;
+    default   : printf("UNDEFINED PID\n"); break;
+    }
+    if (QA) { 
+      cutPidRP->SetQAOn(qaRP); 
+    }
+  }
+  
+  //TRACK QUALITY
+  AliCFTrackQualityCuts *recQualityCutsRP = new AliCFTrackQualityCuts("recQualityCutsRP","rec-level quality cuts");
+  recQualityCutsRP->SetMinNClusterTPC(minClustersTpcRP);
+  //recQualityCutsRP->SetMinFoundClusterTPC(minFoundClustersTpcRP); //only for internal TPC QA
+  recQualityCutsRP->SetMaxChi2PerClusterTPC(maxChi2PerClusterTpcRP);
+  recQualityCutsRP->SetMinNdEdxClusterTPC(minDedxClusterTpcRP);     //to reject secondaries
+
+  recQualityCutsRP->SetMinNClusterITS(minClustersItsRP);
+  recQualityCutsRP->SetMaxChi2PerClusterITS(maxChi2PerClusterItsRP);
+
+  recQualityCutsRP->SetMinNClusterTRD(minClustersTrdRP);
+  recQualityCutsRP->SetMinNTrackletTRD(minTrackletTrdRP);
+  recQualityCutsRP->SetMinNTrackletTRDpid(minTrackletTrdPidRP);
+  recQualityCutsRP->SetMaxChi2PerTrackletTRD(maxChi2PerClusterTrdRP);
+  recQualityCutsRP->SetStatus(statusRP);  
+  if (QA) { 
+    recQualityCutsRP->SetQAOn(qaRP);
+  }
+
+  /* 
+  //How to set this?
+  void SetMaxCovDiagonalElements(Float_t c1=1.e+09, Float_t c2=1.e+09, Float_t c3=1.e+09, Float_t c4=1.e+09, Float_t c5=1.e+09)
+    {fCovariance11Max=c1;fCovariance22Max=c2;fCovariance33Max=c3;fCovariance44Max=c4;fCovariance55Max=c5;}
+  */
+
+  //PRIMARIES
+  AliCFTrackIsPrimaryCuts *recIsPrimaryCutsRP = new AliCFTrackIsPrimaryCuts("recIsPrimaryCutsRP","rec-level isPrimary cuts");
+  recIsPrimaryCutsRP->UseSPDvertex(spdVertexRP);
+  recIsPrimaryCutsRP->UseTPCvertex(tpcVertexRP);
+  recIsPrimaryCutsRP->SetMinDCAToVertexXY(minDcaToVertexXyRP); 
+  recIsPrimaryCutsRP->SetMinDCAToVertexZ(minDcaToVertexZRP);
+  recIsPrimaryCutsRP->SetMaxDCAToVertexXY(maxDcaToVertexXyRP);
+  recIsPrimaryCutsRP->SetMaxDCAToVertexZ(maxDcaToVertexZRP); 
+  recIsPrimaryCutsRP->SetDCAToVertex2D(dcaToVertex2dRP);
+  recIsPrimaryCutsRP->SetAbsDCAToVertex(absDcaToVertexRP);
+  recIsPrimaryCutsRP->SetMinNSigmaToVertex(minNSigmaToVertexRP); 
+  recIsPrimaryCutsRP->SetMaxNSigmaToVertex(maxNSigmaToVertexRP); 
+  recIsPrimaryCutsRP->SetMaxSigmaDCAxy(maxSigmaDcaXySP);
+  recIsPrimaryCutsRP->SetMaxSigmaDCAz(maxSigmaDcaZSP);
+  recIsPrimaryCutsRP->SetRequireSigmaToVertex(requireSigmaToVertexSP);
+  recIsPrimaryCutsRP->SetAcceptKinkDaughters(acceptKinkDaughtersSP);
+  if (QA) { 
+    recIsPrimaryCutsRP->SetQAOn(qaRP);
+  }
+  
+  //ACCEPTANCE
+  AliCFAcceptanceCuts *mcAccCutsRP = new AliCFAcceptanceCuts("mcAccCutsRP","MC acceptance cuts");
+  mcAccCutsRP->SetMinNHitITS(minTrackrefsItsRP);
+  mcAccCutsRP->SetMinNHitTPC(minTrackrefsTpcRP);
+  mcAccCutsRP->SetMinNHitTRD(minTrackrefsTrdRP); 
+  mcAccCutsRP->SetMinNHitTOF(minTrackrefsTofRP);
+  mcAccCutsRP->SetMinNHitMUON(minTrackrefsMuonRP);
+  if (QA) { 
+    mcAccCutsRP->SetQAOn(qaRP);
+  }
+
+  
+  //----------Cuts for POI----------
+  //KINEMATICS (MC and reconstructed)
+  AliCFTrackKineCuts* mcKineCutsPOI = new AliCFTrackKineCuts("mcKineCutsPOI","MC-level kinematic cuts");
+  mcKineCutsPOI->SetPtRange(ptminPOI,ptmaxPOI);
+  mcKineCutsPOI->SetEtaRange(etaminPOI,etamaxPOI);
+  //mcKineCutsPOI->SetChargeMC(chargePOI);
+  mcKineCutsPOI->SetRequireIsCharged(isChargedPOI);
+  if (QA) { 
+    mcKineCutsPOI->SetQAOn(qaPOI);
+  }
+  
+  AliCFTrackKineCuts *recKineCutsPOI = new AliCFTrackKineCuts("recKineCutsPOI","rec-level kine cuts");
+  recKineCutsPOI->SetPtRange(ptminPOI,ptmaxPOI);
+  recKineCutsPOI->SetEtaRange(etaminPOI,etamaxPOI);
+  //recKineCutsPOI->SetChargeRec(chargePOI);
+  recKineCutsPOI->SetRequireIsCharged(isChargedPOI);
+  if (QA) { 
+    recKineCutsPOI->SetQAOn(qaPOI);
+  }
+  
+  //PID (MC and reconstructed)
+  AliCFParticleGenCuts* mcGenCutsPOI = new AliCFParticleGenCuts("mcGenCutsPOI","MC particle generation cuts for POI");
+  mcGenCutsPOI->SetRequireIsPrimary();
+  if (UsePIDforPOI) {mcGenCutsPOI->SetRequirePdgCode(PdgPOI);}
+  if (QA) { 
+    mcGenCutsPOI->SetQAOn(qaPOI);
+  }
+
+  AliCFTrackCutPid* cutPidPOI = NULL;
+  if (UsePIDforPOI) {
+    cutPidPOI = new AliCFTrackCutPid("cutPidPOI","ESD_PID for POI") ;
+    cutPidPOI->SetPriors(prior);
+    cutPidPOI->SetProbabilityCut(0.0);
+    cutPidPOI->SetDetectors("TPC TOF");
+    switch(TMath::Abs(PDG2)) {
+    case 11   : cutPidPOI->SetParticleType(AliPID::kElectron, kTRUE); break;
+    case 13   : cutPidPOI->SetParticleType(AliPID::kMuon    , kTRUE); break;
+    case 211  : cutPidPOI->SetParticleType(AliPID::kPion    , kTRUE); break;
+    case 321  : cutPidPOI->SetParticleType(AliPID::kKaon    , kTRUE); break;
+    case 2212 : cutPidPOI->SetParticleType(AliPID::kProton  , kTRUE); break;
+    default   : printf("UNDEFINED PID\n"); break;
+    }
+    if (QA) { 
+      cutPidPOI->SetQAOn(qaPOI);
+    }
+  }
+
+  //TRACK QUALITY
+  AliCFTrackQualityCuts *recQualityCutsPOI = new AliCFTrackQualityCuts("recQualityCutsPOI","rec-level quality cuts");
+  recQualityCutsPOI->SetMinNClusterTPC(minClustersTpcPOI);
+  //recQualityCutsPOI->SetMinFoundClusterTPC(minFoundClustersTpcPOI); //only for internal TPC QA
+  recQualityCutsPOI->SetMaxChi2PerClusterTPC(maxChi2PerClusterTpcPOI);
+  recQualityCutsPOI->SetMinNdEdxClusterTPC(minDedxClusterTpcPOI);     //to reject secondaries
+
+  recQualityCutsPOI->SetMinNClusterITS(minClustersItsPOI);
+  recQualityCutsPOI->SetMaxChi2PerClusterITS(maxChi2PerClusterItsPOI);
+
+  recQualityCutsPOI->SetMinNClusterTRD(minClustersTrdPOI);
+  recQualityCutsPOI->SetMinNTrackletTRD(minTrackletTrdPOI);
+  recQualityCutsPOI->SetMinNTrackletTRDpid(minTrackletTrdPidPOI);
+  recQualityCutsPOI->SetMaxChi2PerTrackletTRD(maxChi2PerClusterTrdPOI);
+  recQualityCutsPOI->SetStatus(statusPOI); 
+  if (QA) { 
+    recQualityCutsPOI->SetQAOn(qaPOI);
+  }
+
+  //PRIMARIES
+  AliCFTrackIsPrimaryCuts *recIsPrimaryCutsPOI = new AliCFTrackIsPrimaryCuts("recIsPrimaryCutsPOI","rec-level isPrimary cuts");
+  recIsPrimaryCutsPOI->UseSPDvertex(spdVertexPOI);
+  recIsPrimaryCutsPOI->UseTPCvertex(tpcVertexPOI);
+  recIsPrimaryCutsPOI->SetMinDCAToVertexXY(minDcaToVertexXyPOI); 
+  recIsPrimaryCutsPOI->SetMinDCAToVertexZ(minDcaToVertexZPOI);
+  recIsPrimaryCutsPOI->SetMaxDCAToVertexXY(maxDcaToVertexXyPOI);
+  recIsPrimaryCutsPOI->SetMaxDCAToVertexZ(maxDcaToVertexZPOI); 
+  recIsPrimaryCutsPOI->SetDCAToVertex2D(dcaToVertex2dPOI);
+  recIsPrimaryCutsPOI->SetAbsDCAToVertex(absDcaToVertexPOI);
+  recIsPrimaryCutsPOI->SetMinNSigmaToVertex(minNSigmaToVertexPOI); 
+  recIsPrimaryCutsPOI->SetMaxNSigmaToVertex(maxNSigmaToVertexPOI); 
+  recIsPrimaryCutsPOI->SetMaxSigmaDCAxy(maxSigmaDcaXyPOI);
+  recIsPrimaryCutsPOI->SetMaxSigmaDCAz(maxSigmaDcaZPOI);
+  recIsPrimaryCutsPOI->SetRequireSigmaToVertex(requireSigmaToVertexPOI);
+  recIsPrimaryCutsPOI->SetAcceptKinkDaughters(acceptKinkDaughtersPOI);
+  if (QA) { 
+    recIsPrimaryCutsPOI->SetQAOn(qaPOI);
+  }
+
+  //ACCEPTANCE
+  AliCFAcceptanceCuts *mcAccCutsPOI = new AliCFAcceptanceCuts("mcAccCutsPOI","MC acceptance cuts");
+  mcAccCutsPOI->SetMinNHitITS(minTrackrefsItsPOI);
+  mcAccCutsPOI->SetMinNHitTPC(minTrackrefsTpcPOI);
+  mcAccCutsPOI->SetMinNHitTRD(minTrackrefsTrdPOI); 
+  mcAccCutsPOI->SetMinNHitTOF(minTrackrefsTofPOI);
+  mcAccCutsPOI->SetMinNHitMUON(minTrackrefsMuonPOI);
+  if (QA) { 
+    mcAccCutsPOI->SetQAOn(qaPOI);
+  }
+
+     
+  
+  //----------Create Cut Lists----------
+  printf("CREATE EVENT CUTS\n");
+  TObjArray* mcEventList = new TObjArray(0);  
+  if (UseMultCutforESD) mcEventList->AddLast(mcEventCuts);//cut on mult and vertex
+    
+  TObjArray* recEventList = new TObjArray(0);
+  if (UseMultCutforESD) recEventList->AddLast(recEventCuts);//cut on mult and vertex
+
+  printf("CREATE MC KINE CUTS\n");
+  TObjArray* mcListRP = new TObjArray(0);
+  if (UseKineforRP) mcListRP->AddLast(mcKineCutsRP); //cut on pt/eta/phi
+  mcListRP->AddLast(mcGenCutsRP); //cut on primary and if (UsePIDforRP) MC PID
+  
+  TObjArray* mcListPOI = new TObjArray(0);
+  if (UseKineforPOI) mcListPOI->AddLast(mcKineCutsPOI); //cut on pt/eta/phi
+  mcListPOI->AddLast(mcGenCutsPOI); //cut on primary and if (UsePIDforPOI) MC PID
+  
+  printf("CREATE MC ACCEPTANCE CUTS\n");
+  TObjArray* accListRP = new TObjArray(0) ;
+  if (UseAcceptanceforRP) accListRP->AddLast(mcAccCutsRP); //cut on number of track references
+  
+  TObjArray* accListPOI = new TObjArray(0) ;
+  if (UseAcceptanceforPOI) accListPOI->AddLast(mcAccCutsPOI); //cut on number of track references
+  
+  printf("CREATE ESD RECONSTRUCTION CUTS\n");
+  TObjArray* recListRP = new TObjArray(0) ;
+  if (UseTrackQualityforRP) recListRP->AddLast(recQualityCutsRP);   //track quality
+  if (UsePrimariesforRP)    recListRP->AddLast(recIsPrimaryCutsRP); //cut if it is a primary
+  if (UseKineforRP)         recListRP->AddLast(recKineCutsRP);      //cut on pt/eta/phi  
+
+  TObjArray* recListPOI = new TObjArray(0) ;
+  if (UseTrackQualityforPOI) recListPOI->AddLast(recQualityCutsPOI);   //track quality
+  if (UsePrimariesforPOI)    recListPOI->AddLast(recIsPrimaryCutsPOI); //cut if it is a primary
+  if (UseKineforPOI)         recListPOI->AddLast(recKineCutsPOI);      //cut on pt/eta/phi
+
+  printf("CREATE ESD PID CUTS\n");
+  TObjArray* fPIDCutListRP = new TObjArray(0) ;
+  if(UsePIDforRP) {fPIDCutListRP->AddLast(cutPidRP);} //cut on ESD PID
+  
+  TObjArray* fPIDCutListPOI = new TObjArray(0) ;
+  if (UsePIDforPOI)  {fPIDCutListPOI->AddLast(cutPidPOI);} //cut on ESD PID
+  
+
+  //----------Add Cut Lists to the CF Manager----------
+  printf("CREATE INTERFACE AND CUTS\n");
+  AliCFManager* cfmgrRP = new AliCFManager();
+  cfmgrRP->SetNStepEvent(3);
+  cfmgrRP->SetEventCutsList(AliCFManager::kEvtGenCuts,mcEventList); 
+  cfmgrRP->SetEventCutsList(AliCFManager::kEvtRecCuts,recEventList); 
+  cfmgrRP->SetNStepParticle(4); 
+  cfmgrRP->SetParticleCutsList(AliCFManager::kPartGenCuts,mcListRP);
+  cfmgrRP->SetParticleCutsList(AliCFManager::kPartAccCuts,accListRP);
+  cfmgrRP->SetParticleCutsList(AliCFManager::kPartRecCuts,recListRP);
+  cfmgrRP->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutListRP);
+  
+  AliCFManager* cfmgrPOI = new AliCFManager();
+  cfmgrPOI->SetNStepEvent(3);
+  cfmgrPOI->SetEventCutsList(AliCFManager::kEvtGenCuts,mcEventList); 
+  cfmgrPOI->SetEventCutsList(AliCFManager::kEvtRecCuts,recEventList); 
+  cfmgrPOI->SetNStepParticle(4); 
+  cfmgrPOI->SetParticleCutsList(AliCFManager::kPartGenCuts,mcListPOI);
+  cfmgrPOI->SetParticleCutsList(AliCFManager::kPartAccCuts,accListPOI);
+  cfmgrPOI->SetParticleCutsList(AliCFManager::kPartRecCuts,recListPOI);
+  cfmgrPOI->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutListPOI);
+  
+  if (QA) {
+    taskFE->SetQAList1(qaRP);
+    taskFE->SetQAList2(qaPOI);
+  }
+  taskFE->SetCFManager1(cfmgrRP);
+  taskFE->SetCFManager2(cfmgrPOI);
+
+
+
+  // Create the analysis tasks, add them to the manager.
+  //===========================================================================
+  if (SP){
+    AliAnalysisTaskScalarProduct *taskSP = new AliAnalysisTaskScalarProduct("TaskScalarProduct",WEIGHTS[0]);
+    taskSP->SetRelDiffMsub(1.0);
+    taskSP->SetApplyCorrectionForNUA(kFALSE);
+    mgr->AddTask(taskSP);
+  }
+  if (LYZ1SUM){
+    AliAnalysisTaskLeeYangZeros *taskLYZ1SUM = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZerosSUM",kTRUE);
+    taskLYZ1SUM->SetFirstRunLYZ(kTRUE);
+    taskLYZ1SUM->SetUseSumLYZ(kTRUE);
+    mgr->AddTask(taskLYZ1SUM);
+  }
+  if (LYZ1PROD){
+    AliAnalysisTaskLeeYangZeros *taskLYZ1PROD = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZerosPROD",kTRUE);
+    taskLYZ1PROD->SetFirstRunLYZ(kTRUE);
+    taskLYZ1PROD->SetUseSumLYZ(kFALSE);
+    mgr->AddTask(taskLYZ1PROD);
+  }
+  if (LYZ2SUM){
+    AliAnalysisTaskLeeYangZeros *taskLYZ2SUM = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZerosSUM",kFALSE);
+    taskLYZ2SUM->SetFirstRunLYZ(kFALSE);
+    taskLYZ2SUM->SetUseSumLYZ(kTRUE);
+    mgr->AddTask(taskLYZ2SUM);
+  }
+  if (LYZ2PROD){
+    AliAnalysisTaskLeeYangZeros *taskLYZ2PROD = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZerosPROD",kFALSE);
+    taskLYZ2PROD->SetFirstRunLYZ(kFALSE);
+    taskLYZ2PROD->SetUseSumLYZ(kFALSE);
+    mgr->AddTask(taskLYZ2PROD);
+  }
+  if (LYZEP){
+    AliAnalysisTaskLYZEventPlane *taskLYZEP = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane");
+    mgr->AddTask(taskLYZEP);
+  }
+  if (GFC){
+    AliAnalysisTaskCumulants *taskGFC = new AliAnalysisTaskCumulants("TaskCumulants",useWeights);
+    taskGFC->SetUsePhiWeights(WEIGHTS[0]); 
+    taskGFC->SetUsePtWeights(WEIGHTS[1]);
+    taskGFC->SetUseEtaWeights(WEIGHTS[2]); 
+    mgr->AddTask(taskGFC);
+  }
+  if (QC){
+    AliAnalysisTaskQCumulants *taskQC = new AliAnalysisTaskQCumulants("TaskQCumulants",useWeights);
+    taskQC->SetUsePhiWeights(WEIGHTS[0]); 
+    taskQC->SetUsePtWeights(WEIGHTS[1]);
+    taskQC->SetUseEtaWeights(WEIGHTS[2]); 
+    taskQC->SetnBinsMult(10000);
+    taskQC->SetMinMult(0.);
+    taskQC->SetMaxMult(10000.);
+    mgr->AddTask(taskQC);
+  }
+  if (FQD){
+    AliAnalysisTaskFittingQDistribution *taskFQD = new AliAnalysisTaskFittingQDistribution("TaskFittingQDistribution",kFALSE);
+    taskFQD->SetUsePhiWeights(WEIGHTS[0]); 
+    taskFQD->SetqMin(0.);
+    taskFQD->SetqMax(1000.);
+    taskFQD->SetqNbins(10000);
+    mgr->AddTask(taskFQD);
+  }
+  if (MCEP){
+    AliAnalysisTaskMCEventPlane *taskMCEP = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane");
+    mgr->AddTask(taskMCEP);
+  }
+  if (MH){
+    AliAnalysisTaskMixedHarmonics *taskMH = new AliAnalysisTaskMixedHarmonics("TaskMixedHarmonics",useWeights);
+    taskMH->SetHarmonic(1); // n in cos[n(phi1+phi2-2phi3)] and cos[n(psi1+psi2-2phi3)]
+    taskMH->SetNoOfMultipicityBins(10);
+    taskMH->SetMultipicityBinWidth(2);
+    taskMH->SetMinMultiplicity(3);
+    taskMH->SetCorrectForDetectorEffects(kTRUE);
+    taskMH->SetEvaluateDifferential3pCorrelator(kFALSE); // evaluate <<cos[n(psi1+psi2-2phi3)]>> (Remark: two nested loops)    
+    taskMH->SetOppositeChargesPOI(kFALSE); // POIs psi1 and psi2 in cos[n(psi1+psi2-2phi3)] will have opposite charges  
+    mgr->AddTask(taskMH);
+  }  
+  if (NL){
+    AliAnalysisTaskNestedLoops *taskNL = new AliAnalysisTaskNestedLoops("TaskNestedLoops",useWeights);
+    taskNL->SetHarmonic(1); // n in cos[n(phi1+phi2-2phi3)] and cos[n(psi1+psi2-2phi3)]
+    taskNL->SetEvaluateNestedLoopsForRAD(kTRUE); // RAD = Relative Angle Distribution
+    taskNL->SetEvaluateNestedLoopsForMH(kTRUE); // evalaute <<cos[n(phi1+phi2-2phi3)]>> (Remark: three nested loops)   
+    taskNL->SetEvaluateDifferential3pCorrelator(kFALSE); // evaluate <<cos[n(psi1+psi2-2phi3)]>>  (Remark: three nested loops)   
+    taskNL->SetOppositeChargesPOI(kFALSE); // POIs psi1 and psi2 in cos[n(psi1+psi2-2phi3)] will have opposite charges  
+    mgr->AddTask(taskNL);
+  }
+  
+  // Create the output container for the data produced by the task
+  // Connect to the input and output containers
+  //===========================================================================
+  AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer();
+  
+  if (rptype == "FMD") {
+    AliAnalysisDataContainer *coutputFMD = 
+      mgr->CreateContainer(Form("BackgroundCorrected_%s",fileName.Data()), TList::Class(), AliAnalysisManager::kExchangeContainer);                        
+    //input and output taskFMD     
+    mgr->ConnectInput(taskfmd, 0, cinput1);
+    mgr->ConnectOutput(taskfmd, 1, coutputFMD);
+    //input into taskFE
+    mgr->ConnectInput(taskFE,1,coutputFMD);
+  }
+  
+  AliAnalysisDataContainer *coutputFE = 
+    mgr->CreateContainer(Form("cobjFlowEventSimple_%s",fileName.Data()),  AliFlowEventSimple::Class(),AliAnalysisManager::kExchangeContainer);
+  mgr->ConnectInput(taskFE,0,cinput1); 
+  mgr->ConnectOutput(taskFE,1,coutputFE);
+
+  if (QA) { 
+    TString qaNameRPFE = fileName;
+    qaNameRPFE += ":QAforRP_FE_";
+    qaNameRPFE += type;
+
+    AliAnalysisDataContainer *coutputQA1FE =
+      mgr->CreateContainer(Form("QARPFE_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,qaNameRPFE); 
+    
+    TString qaNamePOIFE = fileName;
+    qaNamePOIFE += ":QAforPOI_FE_";
+    qaNamePOIFE += type;
+        
+    AliAnalysisDataContainer *coutputQA2FE =
+      mgr->CreateContainer(Form("QAPOIFE_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,qaNamePOIFE); 
+
+    mgr->ConnectOutput(taskFE,2,coutputQA1FE); 
+    mgr->ConnectOutput(taskFE,3,coutputQA2FE); 
+  }
+
+  // Create the output containers for the data produced by the analysis tasks
+  // Connect to the input and output containers
+  //===========================================================================
+  if (useWeights) {    
+    AliAnalysisDataContainer *cinputWeights = mgr->CreateContainer(Form("cobjWeights_%s",fileName.Data()),TList::Class(),AliAnalysisManager::kInputContainer); 
+  }
+
+  if(SP) {
+    TString outputSP = fileName;
+    outputSP += ":outputSPanalysis";
+    outputSP+= type;
+    
+    AliAnalysisDataContainer *coutputSP = mgr->CreateContainer(Form("cobjSP_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputSP); 
+    mgr->ConnectInput(taskSP,0,coutputFE); 
+    mgr->ConnectOutput(taskSP,1,coutputSP); 
+    if (WEIGHTS[0]) {
+      mgr->ConnectInput(taskSP,1,cinputWeights);
+      cinputWeights->SetData(weightsList);
+    } 
+  }
+  if(LYZ1SUM) {
+    TString outputLYZ1SUM = fileName;
+    outputLYZ1SUM += ":outputLYZ1SUManalysis";
+    outputLYZ1SUM+= type;
+    
+    AliAnalysisDataContainer *coutputLYZ1SUM = mgr->CreateContainer(Form("cobjLYZ1SUM_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZ1SUM); 
+    mgr->ConnectInput(taskLYZ1SUM,0,coutputFE); 
+    mgr->ConnectOutput(taskLYZ1SUM,1,coutputLYZ1SUM); 
+  }
+  if(LYZ1PROD) {
+    TString outputLYZ1PROD = fileName;
+    outputLYZ1PROD += ":outputLYZ1PRODanalysis";
+    outputLYZ1PROD+= type;
+    
+    AliAnalysisDataContainer *coutputLYZ1PROD = mgr->CreateContainer(Form("cobjLYZ1PROD_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZ1PROD); 
+    mgr->ConnectInput(taskLYZ1PROD,0,coutputFE); 
+    mgr->ConnectOutput(taskLYZ1PROD,1,coutputLYZ1PROD);
+  }
+  if(LYZ2SUM) {
+    AliAnalysisDataContainer *cinputLYZ2SUM = mgr->CreateContainer(Form("cobjLYZ2SUMin_%s",fileName.Data()),TList::Class(),AliAnalysisManager::kInputContainer);
+    TString outputLYZ2SUM = fileName;
+    outputLYZ2SUM += ":outputLYZ2SUManalysis";
+    outputLYZ2SUM+= type;
+    
+    AliAnalysisDataContainer *coutputLYZ2SUM = mgr->CreateContainer(Form("cobjLYZ2SUM_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZ2SUM); 
+    mgr->ConnectInput(taskLYZ2SUM,0,coutputFE); 
+    mgr->ConnectInput(taskLYZ2SUM,1,cinputLYZ2SUM);
+    mgr->ConnectOutput(taskLYZ2SUM,1,coutputLYZ2SUM); 
+    cinputLYZ2SUM->SetData(fInputListLYZ2SUM);
+  }
+  if(LYZ2PROD) {
+    AliAnalysisDataContainer *cinputLYZ2PROD = mgr->CreateContainer(Form("cobjLYZ2PRODin_%s",fileName.Data()),TList::Class(),AliAnalysisManager::kInputContainer);
+    TString outputLYZ2PROD = fileName;
+    outputLYZ2PROD += ":outputLYZ2PRODanalysis";
+    outputLYZ2PROD+= type;
+    
+    AliAnalysisDataContainer *coutputLYZ2PROD = mgr->CreateContainer(Form("cobjLYZ2PROD_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZ2PROD); 
+    mgr->ConnectInput(taskLYZ2PROD,0,coutputFE); 
+    mgr->ConnectInput(taskLYZ2PROD,1,cinputLYZ2PROD);
+    mgr->ConnectOutput(taskLYZ2PROD,1,coutputLYZ2PROD); 
+    cinputLYZ2PROD->SetData(fInputListLYZ2PROD);
+  }
+  if(LYZEP) {
+    AliAnalysisDataContainer *cinputLYZEP = mgr->CreateContainer(Form("cobjLYZEPin_%s",fileName.Data()),TList::Class(),AliAnalysisManager::kInputContainer);
+    TString outputLYZEP = fileName;
+    outputLYZEP += ":outputLYZEPanalysis";
+    outputLYZEP+= type;
+    
+    AliAnalysisDataContainer *coutputLYZEP = mgr->CreateContainer(Form("cobjLYZEP_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZEP); 
+    mgr->ConnectInput(taskLYZEP,0,coutputFE); 
+    mgr->ConnectInput(taskLYZEP,1,cinputLYZEP);
+    mgr->ConnectOutput(taskLYZEP,1,coutputLYZEP); 
+    cinputLYZEP->SetData(fInputListLYZEP);
+  }
+  if(GFC) {
+    TString outputGFC = fileName;
+    outputGFC += ":outputGFCanalysis";
+    outputGFC+= type;
+    
+    AliAnalysisDataContainer *coutputGFC = mgr->CreateContainer(Form("cobjGFC_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputGFC); 
+    mgr->ConnectInput(taskGFC,0,coutputFE); 
+    mgr->ConnectOutput(taskGFC,1,coutputGFC);
+    if (useWeights) {
+      mgr->ConnectInput(taskGFC,1,cinputWeights);
+      cinputWeights->SetData(weightsList);
+    } 
+  }
+  if(QC) {
+    TString outputQC = fileName;
+    outputQC += ":outputQCanalysis";
+    outputQC+= type;
+
+    AliAnalysisDataContainer *coutputQC = mgr->CreateContainer(Form("cobjQC_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputQC); 
+    mgr->ConnectInput(taskQC,0,coutputFE); 
+    mgr->ConnectOutput(taskQC,1,coutputQC);
+    if (useWeights) {
+      mgr->ConnectInput(taskQC,1,cinputWeights);
+      cinputWeights->SetData(weightsList);
+    }    
+  }
+  if(FQD) {
+    TString outputFQD = fileName;
+    outputFQD += ":outputFQDanalysis";
+    outputFQD+= type;
+    
+    AliAnalysisDataContainer *coutputFQD = mgr->CreateContainer(Form("cobjFQD_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputFQD); 
+    mgr->ConnectInput(taskFQD,0,coutputFE); 
+    mgr->ConnectOutput(taskFQD,1,coutputFQD);
+    if(useWeights) {
+      mgr->ConnectInput(taskFQD,1,cinputWeights);
+      cinputWeights->SetData(weightsList);
+    } 
+  }
+  if(MCEP) {
+    TString outputMCEP = fileName;
+    outputMCEP += ":outputMCEPanalysis";
+    outputMCEP+= type;
+    
+    AliAnalysisDataContainer *coutputMCEP = mgr->CreateContainer(Form("cobjMCEP_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputMCEP); 
+    mgr->ConnectInput(taskMCEP,0,coutputFE); 
+    mgr->ConnectOutput(taskMCEP,1,coutputMCEP); 
+  }
+  if(MH) {
+    TString outputMH = fileName;
+    outputMH += ":outputMHanalysis";
+    outputMH += type;
+        
+    AliAnalysisDataContainer *coutputMH = mgr->CreateContainer(Form("cobjMH_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputMH); 
+    mgr->ConnectInput(taskMH,0,coutputFE); 
+    mgr->ConnectOutput(taskMH,1,coutputMH); 
+    //if (useWeights) {
+    //  mgr->ConnectInput(taskMH,1,cinputWeights);
+    //  cinputWeights->SetData(weightsList);
+    //} 
+  }
+  if(NL) {
+    TString outputNL = fileName;
+    outputNL += ":outputNLanalysis";
+    outputNL += type;
+        
+    AliAnalysisDataContainer *coutputNL = mgr->CreateContainer(Form("cobjNL_%s",fileName.Data()), TList::Class(),AliAnalysisManager::kOutputContainer,outputNL); 
+    mgr->ConnectInput(taskNL,0,coutputFE); 
+    mgr->ConnectOutput(taskNL,1,coutputNL); 
+    //if (useWeights) {
+    //  mgr->ConnectInput(taskNL,1,cinputWeights);
+    //  cinputWeights->SetData(weightsList);
+    //} 
+  }
+
+}
+
+
+
+
+
index e6d6114..fb5ebc7 100644 (file)
@@ -21,7 +21,7 @@ Bool_t NL       = kTRUE;  // nested loops (for instance distribution of phi1-phi
 
 Bool_t METHODS[] = {SP,LYZ1SUM,LYZ1PROD,LYZ2SUM,LYZ2PROD,LYZEP,GFC,QC,FQD,MCEP,MH,NL};
 
-// Analysis type can be ESD, AOD, MC, ESDMCkineESD, ESDMCkineMC
+// Analysis type can be ESD, AOD, MC, ESDMCkineESD, ESDMCkineMC, MK
 const TString type = "ESD";
 
 // Boolean to fill/not fill the QA histograms
diff --git a/PWG2/FLOW/macros/runFlowTaskCentralityTrain.C b/PWG2/FLOW/macros/runFlowTaskCentralityTrain.C
new file mode 100644 (file)
index 0000000..491ac75
--- /dev/null
@@ -0,0 +1,529 @@
+enum anaModes {mLocal,mLocalPAR,mPROOF,mGRID};
+//mLocal: Analyze locally files in your computer using aliroot
+//mLocalPAR: Analyze locally files in your computer using root + PAR files
+//mPROOF: Analyze CAF files with PROOF
+
+//CENTRALITY DEFINITION
+const Int_t numberOfCentralityBins = 3;
+Int_t centralityArray[numberOfCentralityBins+1] = {0,4000,6000,10000};
+
+// RUN SETTINGS
+
+// Flow analysis method can be:(set to kTRUE or kFALSE)
+Bool_t MCEP     = kTRUE;  // correlation with Monte Carlo reaction plane
+Bool_t SP       = kTRUE;  // scalar product method (similar to eventplane method)
+Bool_t GFC      = kTRUE;  // cumulants based on generating function
+Bool_t QC       = kTRUE;  // cumulants using Q vectors
+Bool_t FQD      = kTRUE;  // fit of the distribution of the Q vector (only integrated v)
+Bool_t LYZ1SUM  = kTRUE;  // Lee Yang Zeroes using sum generating function (integrated v)
+Bool_t LYZ1PROD = kTRUE;  // Lee Yang Zeroes using product generating function (integrated v)
+Bool_t LYZ2SUM  = kFALSE; // Lee Yang Zeroes using sum generating function (second pass differential v)
+Bool_t LYZ2PROD = kFALSE; // Lee Yang Zeroes using product generating function (second pass differential v)
+Bool_t LYZEP    = kFALSE; // Lee Yang Zeroes Event plane using sum generating function (gives eventplane + weight)
+Bool_t MH       = kTRUE;  // azimuthal correlators in mixed harmonics  
+Bool_t NL       = kTRUE;  // nested loops (for instance distribution of phi1-phi2 for all distinct pairs)
+
+Bool_t METHODS[] = {SP,LYZ1SUM,LYZ1PROD,LYZ2SUM,LYZ2PROD,LYZEP,GFC,QC,FQD,MCEP,MH,NL};
+
+// Analysis type can be ESD, AOD, MC, ESDMCkineESD, ESDMCkineMC
+const TString type = "ESD";
+
+// Boolean to fill/not fill the QA histograms
+Bool_t QA = kTRUE;   
+
+// Boolean to use/not use weights for the Q vector
+Bool_t WEIGHTS[] = {kFALSE,kFALSE,kFALSE}; //Phi, v'(pt), v'(eta)
+
+
+//void runFlowTask(Int_t mode=mLocal, Int_t nRuns = 2, 
+//Bool_t DATA = kTRUE, const Char_t* dataDir="/data/alice2/kolk/PP/data/LHC09d/104892/test", Int_t offset = 0)
+//              Bool_t DATA = kFALSE, const Char_t* dataDir="/data/alice2/kolk/PP/LHC09d10/104873", Int_t offset = 0)
+
+void runFlowTaskCentralityTrain(Int_t mode = mPROOF, Int_t nRuns = 50000000, 
+                //Bool_t DATA = kFALSE, const Char_t* dataDir="/PWG2/akisiel/Therminator_midcentral_ESD", Int_t offset=0)
+                //Bool_t DATA = kFALSE, const Char_t* dataDir="/PWG2/akisiel/LHC10d6_0.9TeV_EPOS_12502X", Int_t offset=0)
+                //Bool_t DATA = kFALSE, const Char_t* dataDir="/alice/sim/LHC10d2_117048", Int_t offset=0) //phojet 7 TeV               
+                //Bool_t DATA = kTRUE, const Char_t* dataDir="/alice/data/LHC09d_000104792_p6", Int_t offset=0) //data 0.9 TeV
+                Bool_t DATA = kFALSE, const Char_t* dataDir="/PWG4/morsch/HIJING_CENT_4EV", Int_t offset=0) //hijing Pb Pb pilot
+
+//void runFlowTask(Int_t mode = mGRID, Bool_t DATA = kTRUE)
+{
+  TStopwatch timer;
+  timer.Start();
+  
+  CrossCheckUserSettings(DATA);
+
+  LoadLibraries(mode);
+
+  if (mode == mGRID) {
+    // Create and configure the alien handler plugin
+    gROOT->LoadMacro("CreateAlienHandler.C");
+    AliAnalysisGrid *alienHandler = CreateAlienHandler();  
+    if (!alienHandler) return;
+  }
+  
+  if (mode==mLocal || mode == mLocalPAR) {
+    if (type!="AOD") { TChain* chain = CreateESDChain(dataDir, nRuns, offset);}
+    else { TChain* chain = CreateAODChain(dataDir, nRuns, offset);}
+  }
+  //____________________________________________//
+  // Make the analysis manager
+  AliAnalysisManager *mgr = new AliAnalysisManager("FlowAnalysisManager");
+  if (mode == mGRID) { 
+    // Connect plug-in to the analysis manager
+    mgr->SetGridHandler(alienHandler);
+  }
+
+  if (type == "ESD"){
+    AliVEventHandler* esdH = new AliESDInputHandler;
+    mgr->SetInputEventHandler(esdH);
+    if (MCEP) { 
+      AliMCEventHandler *mc = new AliMCEventHandler();
+      mgr->SetMCtruthEventHandler(mc); 
+    }
+  }
+  
+  if (type == "AOD"){
+    AliVEventHandler* aodH = new AliAODInputHandler;
+    mgr->SetInputEventHandler(aodH); 
+    if (MCEP) { 
+      AliMCEventHandler *mc = new AliMCEventHandler();
+      mgr->SetMCtruthEventHandler(mc);
+    } 
+  }
+  
+  if (type == "MC" || type == "ESDMCkineESD" || type == "ESDMCkineMC"){
+    AliVEventHandler* esdH = new AliESDInputHandler;
+    mgr->SetInputEventHandler(esdH);
+    
+    AliMCEventHandler *mc = new AliMCEventHandler();
+    mgr->SetMCtruthEventHandler(mc); 
+  }
+    
+  //____________________________________________//
+  // Load the analysis task
+  gROOT->LoadMacro("AddTaskFlow.C");
+
+  for (Int_t i=0; i<numberOfCentralityBins; i++)
+  {
+    Int_t lowCentralityBinEdge = centralityArray[i];
+    Int_t highCentralityBinEdge = centralityArray[i+1];
+    TString filename("outputCentrality");
+    filename += i;
+    //TDirectory* dir = new TDirectory(filename.Data(),"");
+    filename += ".root";
+    AddTaskFlowCentrality( type,
+                           METHODS,
+                           QA,
+                           WEIGHTS,
+                           lowCentralityBinEdge,
+                           highCentralityBinEdge,
+                           filename );
+  }
+
+  
+  // Task to check the offline trigger
+  if (mode == mLocal || mode == mGRID) {
+    gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C"); }
+  else if (mode == mPROOF || mode == mLocalPAR) {
+    gROOT->LoadMacro("AddTaskPhysicsSelection.C"); }
+  AliPhysicsSelectionTask* physicsSelTask = AddTaskPhysicsSelection();
+  if (!DATA) {physicsSelTask->GetPhysicsSelection()->SetAnalyzeMC();}
+  
+  // Enable debug printouts
+  mgr->SetDebugLevel(2);
+
+
+  //____________________________________________//
+  // Run the analysis
+  if (!mgr->InitAnalysis()) return;
+  mgr->PrintStatus();
+  
+  if (mode == mLocal || mode == mLocalPAR) {
+    mgr->StartAnalysis("local",chain);
+  }
+  else if (mode == mPROOF) {
+    mgr->StartAnalysis("proof",dataDir,nRuns,offset);
+  }
+  else if (mode == mGRID) { 
+    mgr->StartAnalysis("grid");
+  }
+  
+  timer.Stop();
+  timer.Print();
+  
+}
+
+void CrossCheckUserSettings(Bool_t bData) 
+{
+ // Check in this method if the user settings make sense.
+ if(MCEP==kTRUE && bData==kTRUE)
+ {
+  cout<<endl;
+  cout<<"WARNING: In real datasets there is no Monte Carlo information available !!!!"<<endl;
+  cout<<"         Set for real data analysis DATA = kTRUE and MCEP = kFALSE in the macro."<<endl;
+  cout<<endl;
+  exit(0);
+ }
+
+} // end of void CrossCheckUserSettings()
+
+void LoadLibraries(const anaModes mode) {
+  
+  //--------------------------------------
+  // Load the needed libraries most of them already loaded by aliroot
+  //--------------------------------------
+  gSystem->Load("libTree");
+  gSystem->Load("libGeom");
+  gSystem->Load("libVMC");
+  gSystem->Load("libXMLIO");
+  gSystem->Load("libPhysics");
+  
+  //----------------------------------------------------------
+  // >>>>>>>>>>> Local mode <<<<<<<<<<<<<< 
+  //----------------------------------------------------------
+  if (mode==mLocal || mode==mGRID) {
+    //--------------------------------------------------------
+    // If you want to use already compiled libraries 
+    // in the aliroot distribution
+    //--------------------------------------------------------
+    gSystem->Load("libSTEERBase");
+    gSystem->Load("libESD");
+    gSystem->Load("libAOD");
+    gSystem->Load("libANALYSIS");
+    gSystem->Load("libANALYSISalice");
+    gSystem->Load("libCORRFW");
+    gSystem->Load("libPWG2forward");
+    if (mode==mLocal) {
+      gSystem->Load("libPWG2flowCommon");
+      cerr<<"libPWG2flowCommon loaded..."<<endl;
+      gSystem->Load("libPWG2flowTasks");
+      cerr<<"libPWG2flowTasks loaded..."<<endl;
+    }
+    if (mode==mGRID) {
+      SetupPar("PWG2flowCommon");
+      cerr<<"PWG2flowCommon.par loaded..."<<endl;
+      SetupPar("PWG2flowTasks");
+      cerr<<"PWG2flowTasks.par loaded..."<<endl;
+    }
+  }
+  
+  else if (mode == mLocalPAR) {
+    //--------------------------------------------------------
+    //If you want to use root and par files from aliroot
+    //--------------------------------------------------------  
+    SetupPar("STEERBase");
+    SetupPar("ESD");
+    SetupPar("AOD");
+    SetupPar("ANALYSIS");
+    SetupPar("ANALYSISalice");
+        SetupPar("CORRFW");
+    SetupPar("PWG2flowCommon");
+    cerr<<"PWG2flowCommon.par loaded..."<<endl;
+    SetupPar("PWG2flowTasks");
+    cerr<<"PWG2flowTasks.par loaded..."<<endl;
+  }
+  
+  //---------------------------------------------------------
+  // <<<<<<<<<< PROOF mode >>>>>>>>>>>>
+  //---------------------------------------------------------
+  else if (mode==mPROOF) {
+    //
+    //gEnv->SetValue("XSec.GSI.DelegProxy","2");    
+    //  set to debug root versus if needed
+    //TProof::Mgr("alicecaf")->SetROOTVersion("v5-24-00a_dbg");
+    //TProof::Mgr("alicecaf")->SetROOTVersion("v5-24-00a");
+    //TProof::Reset("proof://snelling@alicecaf.cern.ch");     
+    // Connect to proof
+    printf("*** Connect to PROOF ***\n");
+    gEnv->SetValue("XSec.GSI.DelegProxy","2");
+    // Put appropriate username here
+    //TProof::Open("abilandz@alicecaf.cern.ch");
+    //TProof::Open("nkolk@alicecaf.cern.ch");
+    //TProof::Open("snelling@localhost");
+    TProof::Open("alice-caf.cern.ch");
+    //TProof::Open("skaf.saske.sk");
+    //TProof::Open("prf000-iep-grid.saske.sk");
+    //Info("runSKAF.C","Loading libs on proof (may take while, around 1 min) ...");
+    // list the data available
+    //gProof->ShowDataSets("/*/*"); 
+    //gProof->ShowDataSets("/alice/sim/"); //for MC Data
+    //gProof->ShowDataSets("/alice/data/"); //for REAL Data 
+    // Clear the Packages
+    
+    gProof->ClearPackage("STEERBase.par");
+    gProof->ClearPackage("ESD.par");
+    gProof->ClearPackage("AOD.par");
+    gProof->ClearPackage("ANALYSIS.par");
+    gProof->ClearPackage("ANALYSISalice.par");
+    gProof->ClearPackage("CORRFW.par");
+    
+    gProof->ClearPackage("PWG2flowCommon");
+    gProof->ClearPackage("PWG2flowTasks");
+    
+    // Upload the Packages
+    gProof->UploadPackage("STEERBase.par");
+    gProof->UploadPackage("ESD.par");    
+    gProof->UploadPackage("AOD.par");
+       
+    gProof->UploadPackage("ANALYSIS.par"); 
+    gProof->UploadPackage("ANALYSISalice.par");
+    gProof->UploadPackage("CORRFW.par");
+    gProof->UploadPackage("PWG2flowCommon.par");
+    gProof->UploadPackage("PWG2flowTasks.par");
+
+    // Enable the Packages 
+    // The global package
+    //gProof->EnablePackage("aliroot_v4-19-05-AN",kTRUE);
+    // Or separate
+    
+    gProof->EnablePackage("STEERBase");
+    gProof->EnablePackage("ESD");
+    gProof->EnablePackage("AOD");
+    
+    // Always needed
+    gProof->EnablePackage("ANALYSIS");
+    gProof->EnablePackage("ANALYSISalice");
+    gProof->EnablePackage("CORRFW");
+    gProof->EnablePackage("PWG2flowCommon");
+    gProof->EnablePackage("PWG2flowTasks");
+
+    // Show enables Packages
+    gProof->ShowEnabledPackages();
+  }  
+  
+}
+
+void SetupPar(char* pararchivename) {
+  //Load par files, create analysis libraries
+  //For testing, if par file already decompressed and modified
+  //classes then do not decompress.
+  
+  TString cdir(Form("%s", gSystem->WorkingDirectory() )) ; 
+  TString parpar(Form("%s.par", pararchivename)) ; 
+  if ( gSystem->AccessPathName(parpar.Data()) ) {
+    gSystem->ChangeDirectory(gSystem->Getenv("ALICE_ROOT")) ;
+    TString processline(Form(".! make %s", parpar.Data())) ; 
+    gROOT->ProcessLine(processline.Data()) ;
+    gSystem->ChangeDirectory(cdir) ; 
+    processline = Form(".! mv /tmp/%s .", parpar.Data()) ;
+    gROOT->ProcessLine(processline.Data()) ;
+  } 
+  if ( gSystem->AccessPathName(pararchivename) ) {  
+    TString processline = Form(".! tar xvzf %s",parpar.Data()) ;
+    gROOT->ProcessLine(processline.Data());
+  }
+  
+  TString ocwd = gSystem->WorkingDirectory();
+  gSystem->ChangeDirectory(pararchivename);
+  
+  // check for BUILD.sh and execute
+  if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
+    printf("*******************************\n");
+    printf("*** Building PAR archive    ***\n");
+    cout<<pararchivename<<endl;
+    printf("*******************************\n");
+    if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
+      Error("runProcess","Cannot Build the PAR Archive! - Abort!");
+      return -1;
+    }
+  }
+  // check for SETUP.C and execute
+  if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
+    printf("*******************************\n");
+    printf("*** Setup PAR archive       ***\n");
+    cout<<pararchivename<<endl;
+    printf("*******************************\n");
+    gROOT->Macro("PROOF-INF/SETUP.C");
+  }
+  
+  gSystem->ChangeDirectory(ocwd.Data());
+  printf("Current dir: %s\n", ocwd.Data());
+}
+
+
+// Helper macros for creating chains
+// from: CreateESDChain.C,v 1.10 jgrosseo Exp
+
+TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
+{
+  // creates chain of files in a given directory or file containing a list.
+  // In case of directory the structure is expected as:
+  // <aDataDir>/<dir0>/AliESDs.root
+  // <aDataDir>/<dir1>/AliESDs.root
+  // ...
+  
+  if (!aDataDir)
+    return 0;
+  
+  Long_t id, size, flags, modtime;
+  if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
+    {
+      printf("%s not found.\n", aDataDir);
+      return 0;
+    }
+  
+  TChain* chain = new TChain("esdTree");
+  TChain* chaingAlice = 0;
+  
+  if (flags & 2)
+    {
+      TString execDir(gSystem->pwd());
+      TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
+      TList* dirList            = baseDir->GetListOfFiles();
+      Int_t nDirs               = dirList->GetEntries();
+      gSystem->cd(execDir);
+      
+      Int_t count = 0;
+      
+      for (Int_t iDir=0; iDir<nDirs; ++iDir)
+       {
+         TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
+         if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
+           continue;
+         
+         if (offset > 0)
+           {
+             --offset;
+             continue;
+           }
+         
+         if (count++ == aRuns)
+           break;
+         
+         TString presentDirName(aDataDir);
+         presentDirName += "/";
+         presentDirName += presentDir->GetName();        
+         chain->Add(presentDirName + "/AliESDs.root/esdTree");
+         //  cerr<<presentDirName<<endl;
+       }
+      
+    }
+  else
+    {
+      // Open the input stream
+      ifstream in;
+      in.open(aDataDir);
+      
+      Int_t count = 0;
+      
+      // Read the input list of files and add them to the chain
+      TString esdfile;
+      while(in.good()) {
+       in >> esdfile;
+       if (!esdfile.Contains("root")) continue; // protection
+       
+       if (offset > 0)
+         {
+           --offset;
+           continue;
+         }
+       
+       if (count++ == aRuns)
+         break;
+       
+       // add esd file
+       chain->Add(esdfile);
+      }
+      
+      in.close();
+    }
+  
+  return chain;
+}
+
+
+// Helper macros for creating chains
+// from: CreateESDChain.C,v 1.10 jgrosseo Exp
+
+TChain* CreateAODChain(const char* aDataDir, Int_t aRuns, Int_t offset)
+{
+  // creates chain of files in a given directory or file containing a list.
+  // In case of directory the structure is expected as:
+  // <aDataDir>/<dir0>/AliAOD.root
+  // <aDataDir>/<dir1>/AliAOD.root
+  // ...
+  
+  if (!aDataDir)
+    return 0;
+  
+  Long_t id, size, flags, modtime;
+  if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
+    {
+      printf("%s not found.\n", aDataDir);
+      return 0;
+    }
+  
+  TChain* chain = new TChain("aodTree");
+  TChain* chaingAlice = 0;
+  
+  if (flags & 2)
+    {
+      TString execDir(gSystem->pwd());
+      TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
+      TList* dirList            = baseDir->GetListOfFiles();
+      Int_t nDirs               = dirList->GetEntries();
+      gSystem->cd(execDir);
+      
+      Int_t count = 0;
+      
+      for (Int_t iDir=0; iDir<nDirs; ++iDir)
+       {
+         TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
+         if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
+           continue;
+         
+         if (offset > 0)
+           {
+             --offset;
+             continue;
+           }
+         
+         if (count++ == aRuns)
+           break;
+         
+         TString presentDirName(aDataDir);
+         presentDirName += "/";
+         presentDirName += presentDir->GetName();        
+         chain->Add(presentDirName + "/AliAOD.root/aodTree");
+         // cerr<<presentDirName<<endl;
+       }
+      
+    }
+  else
+    {
+      // Open the input stream
+      ifstream in;
+      in.open(aDataDir);
+      
+      Int_t count = 0;
+      
+      // Read the input list of files and add them to the chain
+      TString aodfile;
+      while(in.good()) {
+       in >> aodfile;
+       if (!aodfile.Contains("root")) continue; // protection
+       
+       if (offset > 0)
+         {
+           --offset;
+           continue;
+         }
+       
+       if (count++ == aRuns)
+         break;
+       
+       // add aod file
+       chain->Add(aodfile);
+      }
+      
+      in.close();
+    }
+  
+  return chain;
+}
+
index 638ee56..19e973e 100644 (file)
@@ -5,7 +5,9 @@
 #pragma link off all functions;
 
 #pragma link C++ class AliFlowEvent+;
+#pragma link C++ class AliFlowEventCuts+;
 #pragma link C++ class AliFlowTrack+;
+#pragma link C++ class AliFlowTrackCuts+;
 #pragma link C++ class AliFlowEventSimpleMaker+;
 
 #pragma link C++ class AliAnalysisTaskScalarProduct+;
index 1cfd730..85b0a27 100644 (file)
@@ -2,7 +2,9 @@
 
 SRCS= FLOW/AliFlowTasks/AliFlowEventSimpleMaker.cxx \
       FLOW/AliFlowTasks/AliFlowEvent.cxx \
+      FLOW/AliFlowTasks/AliFlowEventCuts.cxx \
       FLOW/AliFlowTasks/AliFlowTrack.cxx \
+      FLOW/AliFlowTasks/AliFlowTrackCuts.cxx \
       FLOW/AliFlowTasks/AliAnalysisTaskScalarProduct.cxx \
       FLOW/AliFlowTasks/AliAnalysisTaskMCEventPlane.cxx \
       FLOW/AliFlowTasks/AliAnalysisTaskLYZEventPlane.cxx \