Extension for AOD
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 6 Sep 2011 18:54:11 +0000 (18:54 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 6 Sep 2011 18:54:11 +0000 (18:54 +0000)
Robert Grajcarek <grajcarek@physi.uni-heidelberg.de>

ANALYSIS/AliEPSelectionTask.cxx
ANALYSIS/AliEPSelectionTask.h
ANALYSIS/macros/AddTaskEventplane.C
STEER/STEERBase/AliEventplane.cxx
STEER/STEERBase/AliEventplane.h

index 76d379d..2764442 100644 (file)
@@ -56,7 +56,9 @@
 #include "AliBackgroundSelection.h"
 #include "AliESDUtils.h"
 #include "AliOADBContainer.h"
-
+#include "AliAODMCHeader.h"
+#include "AliAODTrack.h"
+#include "AliVTrack.h"
 #include "AliEventplane.h"
 
 ClassImp(AliEPSelectionTask)
@@ -73,6 +75,9 @@ AliAnalysisTaskSE(),
   fUserphidist(kFALSE),
   fUsercuts(kFALSE),
   fRunNumber(-15),
+  fAODfilterbit(1),
+  fEtaGap(0.),
+  fSplitMethod(0),
   fESDtrackCuts(0),
   fEPContainer(0),
   fPhiDist(0),
@@ -109,6 +114,9 @@ AliEPSelectionTask::AliEPSelectionTask(const char *name):
   fUserphidist(kFALSE),
   fUsercuts(kFALSE),
   fRunNumber(-15),
+  fAODfilterbit(1),
+  fEtaGap(0.),
+  fSplitMethod(0),
   fESDtrackCuts(0),
   fEPContainer(0),
   fPhiDist(0),
@@ -189,8 +197,15 @@ void AliEPSelectionTask::UserCreateOutputObjects()
   
     if(!fUserphidist) { // if it's already set and custom class is required, we use the one provided by the user
 
-    TString oadbfilename = (Form("%s/COMMON/EVENTPLANE/data/epphidist.root", AliAnalysisManager::GetOADBPath()));
+    
+    TString oadbfilename; 
 
+    if (fAnalysisInput.CompareTo("AOD")==0){
+      oadbfilename = (Form("%s/COMMON/EVENTPLANE/data/epphidist.aod.root", AliAnalysisManager::GetOADBPath()));
+    } else if (fAnalysisInput.CompareTo("ESD")==0){
+      oadbfilename = (Form("%s/COMMON/EVENTPLANE/data/epphidist.root", AliAnalysisManager::GetOADBPath()));
+    }
     TFile foadb(oadbfilename); 
     if(!foadb.IsOpen()) AliFatal(Form("Cannot open OADB file %s", oadbfilename.Data()));
 
@@ -199,6 +214,7 @@ void AliEPSelectionTask::UserCreateOutputObjects()
     if (!fEPContainer) AliFatal("Cannot fetch OADB container for EP selection");
     foadb.Close();
     }
+
 }
 
 //________________________________________________________________________
@@ -209,7 +225,7 @@ void AliEPSelectionTask::UserExec(Option_t */*option*/)
   
 //   fRunNumber = -15;
  
-  AliEventplane* esdEP = 0;
+  AliEventplane *esdEP;
   TVector2 qq1;
   TVector2 qq2;
   Double_t fRP = 0.; // the monte carlo reaction plane angle
@@ -237,6 +253,10 @@ void AliEPSelectionTask::UserExec(Option_t */*option*/)
       if (fSaveTrackContribution) {
        esdEP->GetQContributionXArray()->Set(esd->GetNumberOfTracks());
        esdEP->GetQContributionYArray()->Set(esd->GetNumberOfTracks());
+        esdEP->GetQContributionXArraysub1()->Set(esd->GetNumberOfTracks());
+       esdEP->GetQContributionYArraysub1()->Set(esd->GetNumberOfTracks());
+        esdEP->GetQContributionXArraysub2()->Set(esd->GetNumberOfTracks());
+       esdEP->GetQContributionYArraysub2()->Set(esd->GetNumberOfTracks());
       }
       
       TObjArray* tracklist = new TObjArray;
@@ -247,7 +267,7 @@ void AliEPSelectionTask::UserExec(Option_t */*option*/)
       if (nt>4){
        fQVector = new TVector2(GetQ(esdEP,tracklist));
        fEventplaneQ = fQVector->Phi()/2; 
-       GetQsub(qq1, qq2, tracklist);
+       GetQsub(qq1, qq2, tracklist, esdEP);
        fQsub1 = new TVector2(qq1);
        fQsub2 = new TVector2(qq2);
        fQsubRes = (fQsub1->Phi()/2 - fQsub2->Phi()/2);
@@ -287,12 +307,83 @@ void AliEPSelectionTask::UserExec(Option_t */*option*/)
     }
   }
   
-  else if (fAnalysisInput.CompareTo("AOD")==0){
-    //AliAODEvent *aod =  dynamic_cast<AliAODEvent*> (InputEvent());
-    // to be implemented
-    printf("  AOD analysis not yet implemented!!!\n\n");
-    return;
+    else if (fAnalysisInput.CompareTo("AOD")==0){
+    AliAODEvent *aod =  dynamic_cast<AliAODEvent*> (InputEvent());
+
+    if (!(fRunNumber == aod->GetRunNumber())) {
+      fRunNumber = aod->GetRunNumber();
+      SetPhiDist();      
+    }
+  
+    if (fUseMCRP) {
+      AliAODMCHeader *headerH = dynamic_cast<AliAODMCHeader*>(aod->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
+      if (headerH) fRP = headerH->GetReactionPlaneAngle();
+    }
+  
+    if (aod){
+      esdEP = aod->GetHeader()->GetEventplaneP();
+      if(esdEP) {esdEP->Reset();} // reset eventplane if not NULL       
+     
+    Int_t maxID = 0;
+    TObjArray* tracklist = GetAODTracksAndMaxID(aod,maxID);
+       
+    if (fSaveTrackContribution) {
+      esdEP->GetQContributionXArray()->Set(maxID+1);
+      esdEP->GetQContributionYArray()->Set(maxID+1);
+      esdEP->GetQContributionXArraysub1()->Set(maxID+1);
+      esdEP->GetQContributionYArraysub1()->Set(maxID+1);
+      esdEP->GetQContributionXArraysub2()->Set(maxID+1);
+      esdEP->GetQContributionYArraysub2()->Set(maxID+1);
+    }
+       
+    const int NT = tracklist->GetEntries();
+      
+  if (NT>4){
+    fQVector = new TVector2(GetQ(esdEP,tracklist));
+    fEventplaneQ = fQVector->Phi()/2.; 
+    GetQsub(qq1, qq2, tracklist, esdEP);
+    fQsub1 = new TVector2(qq1);
+    fQsub2 = new TVector2(qq2);
+    fQsubRes = (fQsub1->Phi()/2. - fQsub2->Phi()/2.);
+       
+    esdEP->SetQVector(fQVector);
+    esdEP->SetEventplaneQ(fEventplaneQ);
+    esdEP->SetQsub(fQsub1,fQsub2);
+    esdEP->SetQsubRes(fQsubRes);
+       
+    fHOutEventplaneQ->Fill(fEventplaneQ);
+    fHOutsub1sub2->Fill(fQsub1->Phi()/2.,fQsub2->Phi()/2.);
+    fHOutNTEPRes->Fill(NT,fQsubRes);
+
+       if (fUseMCRP) fHOutDiff->Fill(fEventplaneQ, fRP);
+       
+       for (int iter = 0; iter<NT;iter++){
+         AliAODTrack* track = dynamic_cast<AliAODTrack*> (tracklist->At(iter));
+         if (track) {
+           float delta = track->Phi()-fEventplaneQ;
+           while (delta < 0) delta += TMath::Pi();
+           while (delta > TMath::Pi()) delta -= TMath::Pi();
+           fHOutPTPsi->Fill(track->Pt(),delta);
+           fHOutPhi->Fill(track->Phi());
+           fHOutPhiCorr->Fill(track->Phi(),GetPhiWeight(track));
+         }
+       }
+       
+       AliAODTrack* trmax = aod->GetTrack(0);
+       for (int iter = 1; iter<NT;iter++){
+         AliAODTrack* track = dynamic_cast<AliAODTrack*> (tracklist->At(iter));
+         if (track && (track->Pt() > trmax->Pt())) trmax = track;
+       }
+       fHOutleadPTPsi->Fill(trmax->Phi(),fEventplaneQ);      
+      }     
+      delete tracklist;
+      tracklist = 0;
+    }  
+       
+    
   }  
+
+  
   else {
     printf(" Analysis Input not known!\n\n ");
     return;
@@ -312,22 +403,25 @@ TVector2 AliEPSelectionTask::GetQ(AliEventplane* EP, TObjArray* tracklist)
 // Get the Q vector
   TVector2 mQ;
   float mQx=0, mQy=0;
-  AliESDtrack* track;
+  AliVTrack* track;
   Double_t weight;
+  Int_t idtemp = -1;
   
   int nt = tracklist->GetEntries();
 
   for (int i=0; i<nt; i++){
     weight = 1;
-    track = dynamic_cast<AliESDtrack*> (tracklist->At(i));
+    track = dynamic_cast<AliVTrack*> (tracklist->At(i));
     if (track) {
       weight = GetWeight(track);
-      if (fSaveTrackContribution){
-       EP->GetQContributionXArray()->AddAt(weight*cos(2*track->Phi()),track->GetID());
-       EP->GetQContributionYArray()->AddAt(weight*sin(2*track->Phi()),track->GetID());
-      }
-      mQx += (weight*cos(2*track->Phi()));
-      mQy += (weight*sin(2*track->Phi()));
+    if (fSaveTrackContribution){
+      idtemp = track->GetID(); 
+      if ((fAnalysisInput.CompareTo("AOD")==0) && (fAODfilterbit == 128)) idtemp = idtemp*(-1) - 1;
+      EP->GetQContributionXArray()->AddAt(weight*cos(2*track->Phi()),idtemp);
+      EP->GetQContributionYArray()->AddAt(weight*sin(2*track->Phi()),idtemp);
+     }
+     mQx += (weight*cos(2*track->Phi()));
+     mQy += (weight*sin(2*track->Phi()));
     }
   }
   mQ.Set(mQx,mQy);
@@ -335,50 +429,103 @@ TVector2 AliEPSelectionTask::GetQ(AliEventplane* EP, TObjArray* tracklist)
 }
   
   //________________________________________________________________________
-void AliEPSelectionTask::GetQsub(TVector2 &Q1, TVector2 &Q2, TObjArray* tracklist)
+void AliEPSelectionTask::GetQsub(TVector2 &Q1, TVector2 &Q2, TObjArray* tracklist,AliEventplane* EP)
 {
 // Get Qsub
   TVector2 mQ[2];
   float mQx1=0, mQy1=0, mQx2=0, mQy2=0;
   Double_t weight;
 
-  AliESDtrack* track;
+  AliVTrack* track;
   TRandom2 rn = 0;
   
   int nt = tracklist->GetEntries();
   int trackcounter1=0, trackcounter2=0;
-  
-  for (Int_t i = 0; i < nt; i++) {
-    weight = 1;
-    track = dynamic_cast<AliESDtrack*> (tracklist->At(i));
-    if (!track) continue;
-    weight = GetWeight(track);
+  int idtemp = 0;
+
+  if (fSplitMethod == AliEPSelectionTask::kRandom){
+    
+    for (Int_t i = 0; i < nt; i++) {
+      weight = 1;
+      track = dynamic_cast<AliVTrack*> (tracklist->At(i));
+      if (!track) continue;
+      weight = GetWeight(track);
+      idtemp = track->GetID(); 
+      if ((fAnalysisInput.CompareTo("AOD")==0) && (fAODfilterbit == 128)) idtemp = idtemp*(-1) - 1;
     
-    // This loop splits the track set into 2 random subsets
-    if( trackcounter1 < int(nt/2.) && trackcounter2 < int(nt/2.)){
-      float random = rn.Rndm();
-      if(random < .5){
+      // This loop splits the track set into 2 random subsets
+      if( trackcounter1 < int(nt/2.) && trackcounter2 < int(nt/2.)){
+        float random = rn.Rndm();
+        if(random < .5){
+          mQx1 += (weight*cos(2*track->Phi()));
+          mQy1 += (weight*sin(2*track->Phi()));
+          if (fSaveTrackContribution){
+            EP->GetQContributionXArraysub1()->AddAt(weight*cos(2*track->Phi()),idtemp);
+            EP->GetQContributionYArraysub1()->AddAt(weight*sin(2*track->Phi()),idtemp);
+          }
+          trackcounter1++;
+        }
+        else {
+          mQx2 += (weight*cos(2*track->Phi()));
+          mQy2 += (weight*sin(2*track->Phi()));
+          if (fSaveTrackContribution){
+            EP->GetQContributionXArraysub2()->AddAt(weight*cos(2*track->Phi()),idtemp);
+            EP->GetQContributionYArraysub2()->AddAt(weight*sin(2*track->Phi()),idtemp);
+          }
+          trackcounter2++;
+        }
+      }
+      else if( trackcounter1 >= int(nt/2.)){
+        mQx2 += (weight*cos(2*track->Phi()));
+        mQy2 += (weight*sin(2*track->Phi()));
+        if (fSaveTrackContribution){
+          EP->GetQContributionXArraysub2()->AddAt(weight*cos(2*track->Phi()),idtemp);
+          EP->GetQContributionYArraysub2()->AddAt(weight*sin(2*track->Phi()),idtemp);
+        }
+        trackcounter2++;
+      }
+      else {
         mQx1 += (weight*cos(2*track->Phi()));
         mQy1 += (weight*sin(2*track->Phi()));
+        if (fSaveTrackContribution){
+          EP->GetQContributionXArraysub1()->AddAt(weight*cos(2*track->Phi()),idtemp);
+          EP->GetQContributionYArraysub1()->AddAt(weight*sin(2*track->Phi()),idtemp);
+        }
         trackcounter1++;
       }
-      else {
+    }
+  } else if (fSplitMethod == AliEPSelectionTask::kEta) {
+     
+    for (Int_t i = 0; i < nt; i++) {
+      weight = 1;
+      track = dynamic_cast<AliVTrack*> (tracklist->At(i));
+      if (!track) continue;
+      weight = GetWeight(track);
+      Double_t eta = track->Eta();
+      idtemp = track->GetID(); 
+      if ((fAnalysisInput.CompareTo("AOD")==0) && (fAODfilterbit == 128)) idtemp = idtemp*(-1) - 1;
+
+      if (eta > fEtaGap/2.) {  
+        mQx1 += (weight*cos(2*track->Phi()));
+        mQy1 += (weight*sin(2*track->Phi()));
+        if (fSaveTrackContribution){
+          EP->GetQContributionXArraysub1()->AddAt(weight*cos(2*track->Phi()),idtemp);
+          EP->GetQContributionYArraysub1()->AddAt(weight*sin(2*track->Phi()),idtemp);
+        }
+      } else if (eta < -1.*fEtaGap/2.) {
         mQx2 += (weight*cos(2*track->Phi()));
         mQy2 += (weight*sin(2*track->Phi()));
-        trackcounter2++;
+        if (fSaveTrackContribution){
+          EP->GetQContributionXArraysub2()->AddAt(weight*cos(2*track->Phi()),idtemp);
+          EP->GetQContributionYArraysub2()->AddAt(weight*sin(2*track->Phi()),idtemp);
+        }
       }
     }
-    else if( trackcounter1 >= int(nt/2.)){
-      mQx2 += (weight*cos(2*track->Phi()));
-      mQy2 += (weight*sin(2*track->Phi()));
-      trackcounter2++;
-    }
-    else {
-      mQx1 += (weight*cos(2*track->Phi()));
-      mQy1 += (weight*sin(2*track->Phi()));
-      trackcounter1++;
-    }
+  } else {
+    printf("plane resolution determination method not available!\n\n ");
+    return;
   }
+     
   mQ[0].Set(mQx1,mQy1);
   mQ[1].Set(mQx2,mQy2);
   Q1 = mQ[0];
@@ -392,29 +539,59 @@ void AliEPSelectionTask::SetPersonalESDtrackCuts(AliESDtrackCuts* trackcuts){
     delete fESDtrackCuts;
     fESDtrackCuts = 0;
   }
-    
+  if (fAnalysisInput.CompareTo("AOD")==0){
+    AliInfo("ESD track cuts not possible for AOD analysis; please use SetPersonalAODtrackCuts(); using TPC only track cuts");  
+    fUsercuts = kFALSE;
+    SetTrackType("TPC");
+    return;
+  } 
   fUsercuts = kTRUE;
   fESDtrackCuts = trackcuts;
 }
 
+//________________________________________________________________________
+void AliEPSelectionTask::SetPersonalAODtrackCuts(UInt_t filterbit, Float_t etalow, Float_t etaup, Float_t ptlow, Float_t ptup){
+  
+  if(fESDtrackCuts){ 
+    delete fESDtrackCuts;
+    fESDtrackCuts = 0;
+  }
+  if (fAnalysisInput.CompareTo("ESD")==0){
+    AliInfo("AOD track cuts not possible for ESD analysis; please use SetPersonalESDtrackCuts(); using TPC only track cuts");  
+    fUsercuts = kFALSE;
+    SetTrackType("TPC");
+    return;
+  }
+  fUsercuts = kTRUE;
+  fESDtrackCuts = new AliESDtrackCuts();
+  fESDtrackCuts->SetPtRange(ptlow,ptup);
+  fESDtrackCuts->SetEtaRange(etalow,etaup);
+  fAODfilterbit = filterbit;
+}
+
 //_____________________________________________________________________________
+
 void AliEPSelectionTask::SetTrackType(TString tracktype){
-// Set the track type
   fTrackType = tracktype;
   if (!fUsercuts) {
-  if (fTrackType.CompareTo("GLOBAL")==0) fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
-  if (fTrackType.CompareTo("TPC")==0)    fESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
-  fESDtrackCuts->SetPtRange(0.15,20);
+  if (fTrackType.CompareTo("GLOBAL")==0){ 
+    fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
+    fAODfilterbit = 32;
+  }    
+  if (fTrackType.CompareTo("TPC")==0){  
+    fESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+    fAODfilterbit = 128;
+  }
+  fESDtrackCuts->SetPtRange(0.15,20.);
   fESDtrackCuts->SetEtaRange(-0.8,0.8);
   }
 }
 
 //________________________________________________________________________
-Double_t AliEPSelectionTask::GetWeight(AliESDtrack* track)
+Double_t AliEPSelectionTask::GetWeight(TObject* track1)
 {
-// Get weight for track
   Double_t ptweight=1;
-
+  AliVTrack* track = dynamic_cast<AliVTrack*>(track1);
   if (fUsePtWeight) {      
     if (track->Pt()<2) ptweight=track->Pt();
     else ptweight=2;
@@ -423,23 +600,23 @@ Double_t AliEPSelectionTask::GetWeight(AliESDtrack* track)
 }
 
 //________________________________________________________________________
-Double_t AliEPSelectionTask::GetPhiWeight(AliESDtrack* track)
+Double_t AliEPSelectionTask::GetPhiWeight(TObject* track1)
 {
-// Get phi weight for track
   Double_t phiweight=1;
+  AliVTrack* track = dynamic_cast<AliVTrack*>(track1);
   
   if (fUsePhiWeight && fPhiDist && track) {
     Double_t nParticles = fPhiDist->Integral();
     Double_t nPhibins = fPhiDist->GetNbinsX();
   
-    Double_t phi = track->Phi();
+    Double_t Phi = track->Phi();
     
-    while (phi<0) phi += TMath::TwoPi();
-    while (phi>TMath::TwoPi()) phi -= TMath::TwoPi();
+    while (Phi<0) Phi += TMath::TwoPi();
+    while (Phi>TMath::TwoPi()) Phi -= TMath::TwoPi();
       
-    Double_t phiDistValue = fPhiDist->GetBinContent(1+TMath::FloorNint((track->Phi())*nPhibins/TMath::TwoPi()));
+    Double_t PhiDistValue = fPhiDist->GetBinContent(1+TMath::FloorNint((track->Phi())*nPhibins/TMath::TwoPi()));
     
-    if (phiDistValue > 0) phiweight = nParticles/nPhibins/phiDistValue;
+    if (PhiDistValue > 0) phiweight = nParticles/nPhibins/PhiDistValue;
   }
   return phiweight;
 }
@@ -447,7 +624,6 @@ Double_t AliEPSelectionTask::GetPhiWeight(AliESDtrack* track)
 //__________________________________________________________________________
 void AliEPSelectionTask::SetPhiDist() 
 {
-// Set the phi distribution
   if(!fUserphidist) { // if it's already set and custom class is required, we use the one provided by the user
 
     fPhiDist = (TH1F*) fEPContainer->GetObject(fRunNumber, "Default");
@@ -485,7 +661,7 @@ void AliEPSelectionTask::SetPhiDist()
 //__________________________________________________________________________
 void AliEPSelectionTask::SetPersonalPhiDistribution(const char* infilename, char* listname)
 {
-    // Set a personal phi distribution
+  
   fUserphidist = kTRUE;
   
   TFile f(infilename);
@@ -495,3 +671,37 @@ void AliEPSelectionTask::SetPersonalPhiDistribution(const char* infilename, char
 
   f.Close();
 } 
+
+
+//_________________________________________________________________________
+TObjArray* AliEPSelectionTask::GetAODTracksAndMaxID(AliAODEvent* aod, Int_t& maxid)
+{
+  TObjArray *acctracks = new TObjArray();
+  
+  AliAODTrack *tr = 0;
+  Int_t maxid1 = 0;
+  Int_t maxidtemp = -1;
+  Float_t ptlow = 0;
+  Float_t ptup = 0;
+  Float_t etalow = 0;
+  Float_t etaup = 0;
+  fESDtrackCuts->GetPtRange(ptlow,ptup);
+  fESDtrackCuts->GetEtaRange(etalow,etaup);
+  
+  for (Int_t i = 0; i < aod->GetNumberOfTracks() ; i++){
+     tr = aod->GetTrack(i);
+     maxidtemp = tr->GetID(); 
+     if(maxidtemp < 0 && fAODfilterbit != 128) continue;
+     if(maxidtemp > -1 && fAODfilterbit == 128) continue;
+     if (fAODfilterbit == 128) maxidtemp = maxidtemp*(-1) - 1;
+     if (maxidtemp > maxid1) maxid1 = maxidtemp;
+     if(tr->TestFilterBit(fAODfilterbit) && tr->Pt() < ptup && tr->Pt() > ptlow && tr->Eta() < etaup && tr->Eta() > etalow){
+     acctracks->Add(tr);
+     }
+  }
+  
+  maxid = maxid1;
+  
+  return acctracks;
+  
+}
index 7f56495..2eb4c61 100644 (file)
@@ -27,6 +27,8 @@ class AliOADBContainer;
 class AliEPSelectionTask : public AliAnalysisTaskSE {
 
  public:
+  
+  enum ResoMethod{kRandom,kEta};
 
   AliEPSelectionTask();
   AliEPSelectionTask(const char *name);
@@ -38,26 +40,31 @@ class AliEPSelectionTask : public AliAnalysisTaskSE {
   virtual void Terminate(Option_t *option);
   
   TVector2 GetQ(AliEventplane* EP, TObjArray* event);
-  void GetQsub(TVector2& Qsub1, TVector2& Qsub2, TObjArray* event);
-  Double_t GetWeight(AliESDtrack* track);
-  Double_t GetPhiWeight(AliESDtrack* track);
-
-  virtual void  SetDebugLevel(Int_t level) {fDebug = level;}
-  void SetInput(const char* input)         {fAnalysisInput = input;}
-  void SetUseMCRP()                       {fUseMCRP = kTRUE;}
-  void SetUsePhiWeight()                  {fUsePhiWeight = kTRUE;}
-  void SetUsePtWeight()                           {fUsePtWeight = kTRUE;}
-  void SetSaveTrackContribution()         {fSaveTrackContribution = kTRUE;}
+  void GetQsub(TVector2& Qsub1, TVector2& Qsub2, TObjArray* event,AliEventplane* EP);
+  Double_t GetWeight(TObject* track1);
+  Double_t GetPhiWeight(TObject* track1);
+
+  virtual void  SetDebugLevel(Int_t level)   {fDebug = level;}
+  void SetInput(const char* input)           {fAnalysisInput = input;}
+  void SetUseMCRP()                         {fUseMCRP = kTRUE;}
+  void SetUsePhiWeight()                    {fUsePhiWeight = kTRUE;}
+  void SetUsePtWeight()                             {fUsePtWeight = kTRUE;}
+  void SetSaveTrackContribution()           {fSaveTrackContribution = kTRUE;}
   void SetTrackType(TString tracktype);
   void SetPhiDist();
   void SetPersonalESDtrackCuts(AliESDtrackCuts* trackcuts);
+  void SetPersonalAODtrackCuts(UInt_t filterbit = 1, Float_t etalow = -0.8, Float_t etaup = 0.8, Float_t ptlow = 0.5, Float_t ptup = 20.);
   void SetPersonalPhiDistribution(const char* filename, char* listname);
+  void SetEtaGap(Float_t etagap)             {fEtaGap = etagap;}
+  void SetSubeventsSplitMethod(Int_t method) {fSplitMethod = method;}
   
  private:
    
   AliEPSelectionTask(const AliEPSelectionTask& ep);
   AliEPSelectionTask& operator= (const AliEPSelectionTask& ep); 
 
+  TObjArray* GetAODTracksAndMaxID(AliAODEvent* aod, Int_t& maxid);
+
   TString  fAnalysisInput;             // "ESD", "AOD"
   TString  fTrackType;                 // "GLOBAL", "TPC"
   Bool_t   fUseMCRP;                   // i.e. usable for Therminator, when MC RP is provided
@@ -67,7 +74,11 @@ class AliEPSelectionTask : public AliAnalysisTaskSE {
   Bool_t   fUserphidist;               // bool, if personal phi distribution should be used
   Bool_t   fUsercuts;                  // bool, if personal cuts should be used
   Int_t    fRunNumber;                 // runnumber
+  UInt_t   fAODfilterbit;               // AOD filter bit for AOD track selection  
+  Float_t  fEtaGap;                     // Eta Gap between Subevent A and B
+  Int_t    fSplitMethod;                // Splitting Method for subevents
   
+
   AliESDtrackCuts* fESDtrackCuts;       // track cuts
   
   AliOADBContainer* fEPContainer;      //! OADB Container
index fe35e6a..3d33454 100644 (file)
@@ -13,13 +13,13 @@ AliEPSelectionTask *AddTaskEventplane()
     return NULL;
   }
   TString inputDataType = mgr->GetInputEventHandler()->GetDataType(); // can be "ESD" or "AOD"
-  if (inputDataType != "ESD") {
-    ::Error("AddTaskEventplane", "This task works only on ESD analysis");
-    return NULL;
-  }
   
   AliEPSelectionTask *eventplaneTask = new AliEPSelectionTask("EventplaneSelection");
   eventplaneTask->SelectCollisionCandidates(AliVEvent::kMB);
+  if (inputDataType == "AOD"){
+    eventplaneTask->SetInput("AOD");
+  }
   eventplaneTask->SetTrackType("TPC");
   eventplaneTask->SetUsePtWeight();
   eventplaneTask->SetUsePhiWeight();
index fa24f14..c1875d1 100644 (file)
@@ -31,6 +31,10 @@ AliEventplane::AliEventplane() : TNamed("Eventplane", "Eventplane"),
   fQVector(0),
   fQContributionX(0),
   fQContributionY(0),
+  fQContributionXsub1(0),
+  fQContributionYsub1(0),
+  fQContributionXsub2(0),
+  fQContributionYsub2(0),
   fEventplaneQ(-1),
   fQsub1(0),
   fQsub2(0),
@@ -39,6 +43,10 @@ AliEventplane::AliEventplane() : TNamed("Eventplane", "Eventplane"),
   /// constructor
   fQContributionX = new TArrayF(0);
   fQContributionY = new TArrayF(0);
+  fQContributionXsub1 = new TArrayF(0);
+  fQContributionYsub1 = new TArrayF(0);
+  fQContributionXsub2 = new TArrayF(0);
+  fQContributionYsub2 = new TArrayF(0);
 }
 
 AliEventplane::AliEventplane(const AliEventplane& ep) : 
@@ -46,6 +54,10 @@ AliEventplane::AliEventplane(const AliEventplane& ep) :
   fQVector(0),
   fQContributionX(0),
   fQContributionY(0),
+  fQContributionXsub1(0),
+  fQContributionYsub1(0),
+  fQContributionXsub2(0),
+  fQContributionYsub2(0),
   fEventplaneQ(0),
   fQsub1(0),
   fQsub2(0),
@@ -72,6 +84,14 @@ void AliEventplane::CopyEP(AliEventplane& ep) const
       target.fQContributionX = fQContributionX;
   if (fQContributionY)
       target.fQContributionY = fQContributionY;
+  if (fQContributionXsub1)
+      target.fQContributionXsub1 = fQContributionXsub1;
+  if (fQContributionYsub1)
+      target.fQContributionYsub1 = fQContributionYsub1;
+  if (fQContributionXsub2)
+      target.fQContributionXsub2 = fQContributionXsub2;
+  if (fQContributionYsub2)
+      target.fQContributionYsub2 = fQContributionYsub2;
   if (fEventplaneQ)
       target.fEventplaneQ = fEventplaneQ;
   if (fQVector)
@@ -95,6 +115,22 @@ AliEventplane::~AliEventplane()
       delete fQContributionY;
       fQContributionY = 0;
   }
+  if (fQContributionXsub1){
+      delete fQContributionXsub1;
+      fQContributionXsub1 = 0;
+  }
+  if (fQContributionYsub1){
+      delete fQContributionYsub1;
+      fQContributionYsub1 = 0;
+  }
+  if (fQContributionXsub2){
+      delete fQContributionXsub2;
+      fQContributionXsub2 = 0;
+  }
+  if (fQContributionYsub2){
+      delete fQContributionYsub2;
+      fQContributionYsub2 = 0;
+  }
   if (fQVector){
       delete fQVector;
       fQVector = 0;
@@ -153,11 +189,35 @@ Double_t AliEventplane::GetQContributionY(AliVTrack* track)
   return fQContributionY->GetAt(track->GetID());
 }
 
+Double_t AliEventplane::GetQContributionXsub1(AliVTrack* track)
+{ 
+  return fQContributionXsub1->GetAt(track->GetID());
+}
+
+Double_t AliEventplane::GetQContributionYsub1(AliVTrack* track)
+{ 
+  return fQContributionYsub1->GetAt(track->GetID());
+}
+
+Double_t AliEventplane::GetQContributionXsub2(AliVTrack* track)
+{ 
+  return fQContributionXsub2->GetAt(track->GetID());
+}
+
+Double_t AliEventplane::GetQContributionYsub2(AliVTrack* track)
+{ 
+  return fQContributionYsub2->GetAt(track->GetID());
+}
+
 void AliEventplane::Reset()
 { 
   delete fQVector; fQVector=0;
   fQContributionX->Reset();
   fQContributionY->Reset();
+  fQContributionXsub1->Reset();
+  fQContributionYsub1->Reset();
+  fQContributionXsub2->Reset();
+  fQContributionYsub2->Reset();
   fEventplaneQ = -1;
   delete fQsub1; fQsub1=0;
   delete fQsub2; fQsub2=0;
index edb21a6..276a105 100644 (file)
@@ -1,5 +1,5 @@
-#ifndef ALIEventplane_H
-#define ALIEventplane_H
+#ifndef ALIEVENTPLANE_H
+#define ALIEVENTPLANE_H
 
 /* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
@@ -36,8 +36,16 @@ class AliEventplane : public TNamed
   TVector2* GetQVector(); 
   Double_t  GetQContributionX(AliVTrack* track);
   Double_t  GetQContributionY(AliVTrack* track);
+  Double_t  GetQContributionXsub1(AliVTrack* track);
+  Double_t  GetQContributionYsub1(AliVTrack* track);
+  Double_t  GetQContributionXsub2(AliVTrack* track);
+  Double_t  GetQContributionYsub2(AliVTrack* track);
   TArrayF*  GetQContributionXArray() { return fQContributionX; }
   TArrayF*  GetQContributionYArray() { return fQContributionY; }
+  TArrayF*  GetQContributionXArraysub1() { return fQContributionXsub1; }
+  TArrayF*  GetQContributionYArraysub1() { return fQContributionYsub1; }
+  TArrayF*  GetQContributionXArraysub2() { return fQContributionXsub2; }
+  TArrayF*  GetQContributionYArraysub2() { return fQContributionYsub2; }
   Double_t  GetEventplane(const char *method);
   TVector2* GetQsub1();
   TVector2* GetQsub2();
@@ -47,13 +55,17 @@ class AliEventplane : public TNamed
   void Reset();
 
  private:
-   TVector2* fQVector;         // Q-Vector of event
-   TArrayF* fQContributionX;   // array of the tracks' contributions to X component of Q-Vector - index = track ID
-   TArrayF* fQContributionY;   // array of the tracks' contributions to Y component of Q-Vector - index = track ID
-   Double_t fEventplaneQ;      // Event plane angle from Q-Vector
-   TVector2* fQsub1;           // Q-Vector of subevent 1
-   TVector2* fQsub2;           // Q-Vector of subevent 2
-   Double_t fQsubRes;          // Difference of EP angles of subevents
+   TVector2* fQVector;          // Q-Vector of event
+   TArrayF* fQContributionX;    // array of the tracks' contributions to X component of Q-Vector - index = track ID
+   TArrayF* fQContributionY;    // array of the tracks' contributions to Y component of Q-Vector - index = track ID
+   TArrayF* fQContributionXsub1; // array of the tracks' contributions to X component of Q-Vectorsub1 - index = track ID
+   TArrayF* fQContributionYsub1; // array of the tracks' contributions to Y component of Q-Vectorsub1 - index = track ID
+   TArrayF* fQContributionXsub2; // array of the tracks' contributions to X component of Q-Vectorsub2 - index = track ID
+   TArrayF* fQContributionYsub2; // array of the tracks' contributions to Y component of Q-Vectorsub2 - index = track ID
+   Double_t fEventplaneQ;       // Event plane angle from Q-Vector
+   TVector2* fQsub1;            // Q-Vector of subevent 1
+   TVector2* fQsub2;            // Q-Vector of subevent 2
+   Double_t fQsubRes;           // Difference of EP angles of subevents
  
   ClassDef(AliEventplane, 1)
 };