]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Adding the trigger and vertex checking in both the QA and the actual analysis
authorpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 29 Dec 2008 18:05:30 +0000 (18:05 +0000)
committerpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 29 Dec 2008 18:05:30 +0000 (18:05 +0000)
PWG2/AliAnalysisTaskProtons.cxx
PWG2/AliAnalysisTaskProtons.h
PWG2/AliAnalysisTaskProtonsQA.cxx
PWG2/SPECTRA/AliProtonAnalysis.cxx
PWG2/SPECTRA/AliProtonAnalysis.h
PWG2/SPECTRA/AliProtonQAAnalysis.cxx
PWG2/SPECTRA/AliProtonQAAnalysis.h
PWG2/runProtonAnalysis.C

index 55a2ed1feaf7c3cba7351231fc994921e7ac5e68..de5fc253e1f9cec50138c1cc80539dcc2252b1a7 100644 (file)
@@ -19,6 +19,8 @@
 #include "AliMCEvent.h"
 #include "AliStack.h"
 #include "AliCFContainer.h"
+#include "AliVertexerTracks.h"
+#include "AliESDVertex.h"
 
 #include "AliProtonAnalysis.h"
 #include "AliAnalysisTaskProtons.h"
@@ -31,7 +33,7 @@ ClassImp(AliAnalysisTaskProtons)
 //________________________________________________________________________ 
 AliAnalysisTaskProtons::AliAnalysisTaskProtons()
   : AliAnalysisTask(), fESD(0), fAOD(0), fMC(0), fAnalysisType("ESD"),
-    fList(0), fAnalysis(0),
+    fList(0), fProtonAnalysis(0),
     fElectronFunction(0), fMuonFunction(0),
     fPionFunction(0), fKaonFunction(0), fProtonFunction(0),
     fFunctionUsed(kFALSE) {
@@ -42,10 +44,12 @@ AliAnalysisTaskProtons::AliAnalysisTaskProtons()
 //________________________________________________________________________
 AliAnalysisTaskProtons::AliAnalysisTaskProtons(const char *name) 
 : AliAnalysisTask(name, ""), fESD(0), fAOD(0), fMC(0), fAnalysisType("ESD"),
-  fList(0), fAnalysis(0), 
+  fList(0), fProtonAnalysis(0), 
   fElectronFunction(0), fMuonFunction(0),
   fPionFunction(0), fKaonFunction(0), fProtonFunction(0),
-  fFunctionUsed(kFALSE) { 
+  fFunctionUsed(kFALSE),
+  fTriggerMode(kMB2), fProtonAnalysisMode(kTPC),
+  fVxMax(0), fVyMax(0), fVzMax(0) { 
   // Constructor
 
   // Define input and output slots here
@@ -65,9 +69,6 @@ void AliAnalysisTaskProtons::ConnectInputData(Option_t *) {
     Printf("ERROR: Could not read chain from input slot 0");
   } else {
     if(fAnalysisType == "ESD") {
-      tree->SetBranchStatus("*", kFALSE);
-      tree->SetBranchStatus("Tracks.*", kTRUE);
-      
       AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
      
       if (!esdH) {
@@ -104,44 +105,89 @@ void AliAnalysisTaskProtons::CreateOutputObjects() {
   Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
   
   //proton analysis object
-  fAnalysis = new AliProtonAnalysis();
-  fAnalysis->InitAnalysisHistograms(10,-1.0,1.0,26,0.5,3.1);
+  fProtonAnalysis = new AliProtonAnalysis();
 
   if(fAnalysisType == "ESD") {
     //Use of TPConly tracks
-    fAnalysis->UseTPCOnly();
-
-    //TPC related cuts       
-    fAnalysis->SetMinTPCClusters(50);
-    fAnalysis->SetMaxChi2PerTPCCluster(3.5);
-    fAnalysis->SetMaxCov11(2.0);
-    fAnalysis->SetMaxCov22(2.0);
-    fAnalysis->SetMaxCov33(0.5);
-    fAnalysis->SetMaxCov44(0.5);
-    fAnalysis->SetMaxCov55(2.0);
-    fAnalysis->SetMaxSigmaToVertex(2.5);
-    fAnalysis->SetTPCRefit();
-
-    //ITS related cuts - to be used for the analysis of global tracking 
-    //fAnalysis->SetMinITSClusters(5);
-    //fAnalysis->SetITSRefit();
-  }
+    if(fProtonAnalysisMode == kTPC) {
+      fProtonAnalysis->InitAnalysisHistograms(10, -0.5, 0.5, 16, 0.5, 0.9);
+      fProtonAnalysis->UseTPCOnly();
+      fProtonAnalysis->SetTPCpid();
+      fProtonAnalysis->SetMinTPCClusters(100);
+      fProtonAnalysis->SetMaxChi2PerTPCCluster(2.2);
+      fProtonAnalysis->SetMaxCov11(0.5);
+      fProtonAnalysis->SetMaxCov22(0.5);
+      fProtonAnalysis->SetMaxCov33(0.5);
+      fProtonAnalysis->SetMaxCov44(0.5);
+      fProtonAnalysis->SetMaxCov55(0.5);
+      fProtonAnalysis->SetMaxSigmaToVertexTPC(2.0);
+      //fProtonAnalysis->SetMaxDCAXYTPC(1.5);
+      //fProtonAnalysis->SetMaxDCAZTPC(1.5);
+    }
+    //Use of HybridTPC tracks
+    else if(fProtonAnalysisMode == kHybrid) {
+      fProtonAnalysis->InitAnalysisHistograms(10, -0.5, 0.5, 16, 0.5, 0.9);
+      fProtonAnalysis->UseHybridTPC();
+      fProtonAnalysis->SetTPCpid();
+      fProtonAnalysis->SetMinTPCClusters(110);
+      fProtonAnalysis->SetMaxChi2PerTPCCluster(2.2);
+      fProtonAnalysis->SetMaxCov11(0.5);
+      fProtonAnalysis->SetMaxCov22(0.5);
+      fProtonAnalysis->SetMaxCov33(0.5);
+      fProtonAnalysis->SetMaxCov44(0.5);
+      fProtonAnalysis->SetMaxCov55(0.5);
+      fProtonAnalysis->SetMaxSigmaToVertex(2.0);
+      /*fProtonAnalysis->SetMaxDCAXY(1.5);
+       fProtonAnalysis->SetMaxDCAZ(1.5);*/
+      fProtonAnalysis->SetPointOnITSLayer6();
+      fProtonAnalysis->SetPointOnITSLayer5();
+      //fProtonAnalysis->SetPointOnITSLayer4();
+      //fProtonAnalysis->SetPointOnITSLayer3();
+      fProtonAnalysis->SetPointOnITSLayer2();
+      fProtonAnalysis->SetPointOnITSLayer1();
+      fProtonAnalysis->SetMinITSClusters(5);
+    }
+    //Combined tracking
+    else if(fProtonAnalysisMode == kGlobal) {
+      fProtonAnalysis->InitAnalysisHistograms(10, -0.5, 0.5, 16, 0.5, 0.9);
+      fProtonAnalysis->SetMinTPCClusters(110);
+      fProtonAnalysis->SetMaxChi2PerTPCCluster(2.2);
+      fProtonAnalysis->SetMaxCov11(0.5);
+      fProtonAnalysis->SetMaxCov22(0.5);
+      fProtonAnalysis->SetMaxCov33(0.5);
+      fProtonAnalysis->SetMaxCov44(0.5);
+      fProtonAnalysis->SetMaxCov55(0.5);
+      fProtonAnalysis->SetMaxSigmaToVertex(2.0);
+      //fProtonAnalysis->SetMaxDCAXY(2.0);
+      //fProtonAnalysis->SetMaxDCAZ(2.0);
+      fProtonAnalysis->SetTPCRefit();
+      fProtonAnalysis->SetPointOnITSLayer1();
+      fProtonAnalysis->SetPointOnITSLayer2();
+      //fProtonAnalysis->SetPointOnITSLayer3();
+      //fProtonAnalysis->SetPointOnITSLayer4();
+      fProtonAnalysis->SetPointOnITSLayer5();
+      fProtonAnalysis->SetPointOnITSLayer6();
+      fProtonAnalysis->SetMinITSClusters(5);
+      fProtonAnalysis->SetITSRefit();
+      fProtonAnalysis->SetESDpid();
+    }
 
-  if(fFunctionUsed)
-    fAnalysis->SetPriorProbabilityFunctions(fElectronFunction,
-                                           fMuonFunction,
-                                           fPionFunction,
-                                           fKaonFunction,
-                                           fProtonFunction);
-  else
-    fAnalysis->SetPriorProbabilities(partFrac);
+    if(fFunctionUsed)
+      fProtonAnalysis->SetPriorProbabilityFunctions(fElectronFunction,
+                                             fMuonFunction,
+                                             fPionFunction,
+                                             fKaonFunction,
+                                             fProtonFunction);
+    else
+      fProtonAnalysis->SetPriorProbabilities(partFrac);
+  }//ESD analysis
 
   fList = new TList();
-  fList->Add(fAnalysis->GetProtonYPtHistogram());
-  fList->Add(fAnalysis->GetAntiProtonYPtHistogram());
-  fList->Add(fAnalysis->GetEventHistogram());
-  fList->Add(fAnalysis->GetProtonContainer());
-  fList->Add(fAnalysis->GetAntiProtonContainer());
+  fList->Add(fProtonAnalysis->GetProtonYPtHistogram());
+  fList->Add(fProtonAnalysis->GetAntiProtonYPtHistogram());
+  fList->Add(fProtonAnalysis->GetEventHistogram());
+  fList->Add(fProtonAnalysis->GetProtonContainer());
+  fList->Add(fProtonAnalysis->GetAntiProtonContainer());
 }
 
 //________________________________________________________________________
@@ -155,8 +201,14 @@ void AliAnalysisTaskProtons::Exec(Option_t *) {
       return;
     }
 
-    Printf("Proton ESD analysis task: There are %d tracks in this event", fESD->GetNumberOfTracks());
-    fAnalysis->Analyze(fESD);
+    if(IsEventTriggered(fESD,fTriggerMode)) {
+      const AliESDVertex *vertex = GetVertex(fESD,fProtonAnalysisMode,
+                                            fVxMax,fVyMax,fVzMax);
+      if(vertex) {
+       Printf("Proton ESD analysis task: There are %d tracks in this event", fESD->GetNumberOfTracks());
+       fProtonAnalysis->Analyze(fESD,vertex);
+      }//reconstructed vertex
+    }//triggered event
   }//ESD analysis              
   
   else if(fAnalysisType == "AOD") {
@@ -166,7 +218,7 @@ void AliAnalysisTaskProtons::Exec(Option_t *) {
     }
     
     Printf("Proton AOD analysis task: There are %d tracks in this event", fAOD->GetNumberOfTracks());
-    fAnalysis->Analyze(fAOD);
+    fProtonAnalysis->Analyze(fAOD);
   }//AOD analysis                                                                                                                                                                
   
   else if(fAnalysisType == "MC") {
@@ -181,7 +233,7 @@ void AliAnalysisTaskProtons::Exec(Option_t *) {
       return;
     }
     Printf("Proton MC analysis task: There are %d primaries in this event", stack->GetNprimary());
-    fAnalysis->Analyze(stack);
+    fProtonAnalysis->Analyze(stack);
   }//MC analysis                      
 
   // Post output data.
@@ -211,4 +263,83 @@ void AliAnalysisTaskProtons::Terminate(Option_t *) {
   if (fHistYPtAntiProtons) fHistYPtAntiProtons->DrawCopy("colz");
 }
 
+//________________________________________________________________________
+Bool_t AliAnalysisTaskProtons::IsEventTriggered(const AliESDEvent *esd, 
+                                               TriggerMode trigger) {
+  // check if the event was triggered
+  ULong64_t triggerMask = esd->GetTriggerMask();
+  
+  // definitions from p-p.cfg
+  ULong64_t spdFO = (1 << 14);
+  ULong64_t v0left = (1 << 11);
+  ULong64_t v0right = (1 << 12);
+
+  switch (trigger) {
+  case kMB1: {
+    if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
+      return kTRUE;
+    break;
+  }
+  case kMB2: {
+    if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
+      return kTRUE;
+    break;
+  }
+  case kSPDFASTOR: {
+    if (triggerMask & spdFO)
+      return kTRUE;
+    break;
+  }
+  }//switch
+
+  return kFALSE;
+}
+
+//________________________________________________________________________
+const AliESDVertex* AliAnalysisTaskProtons::GetVertex(AliESDEvent* esd, 
+                                                     AnalysisMode mode,
+                                                     Double_t gVxMax,
+                                                     Double_t gVyMax,
+                                                     Double_t gVzMax) {
+  // Get the vertex from the ESD and returns it if the vertex is valid
+  // Second argument decides which vertex is used (this selects
+  // also the quality criteria that are applied)
+  const AliESDVertex* vertex = 0;
+  if(mode == kHybrid)
+    vertex = esd->GetPrimaryVertexSPD();
+  else if(mode == kTPC){
+    Double_t kBz = esd->GetMagneticField();
+    AliVertexerTracks vertexer(kBz);
+    vertexer.SetTPCMode();
+    AliESDVertex *vTPC = vertexer.FindPrimaryVertex(esd);
+    esd->SetPrimaryVertexTPC(vTPC);
+    for (Int_t i=0; i<esd->GetNumberOfTracks(); i++) {
+      AliESDtrack *t = esd->GetTrack(i);
+      t->RelateToVertexTPC(vTPC, kBz, kVeryBig);
+    }
+    delete vTPC;
+    vertex = esd->GetPrimaryVertexTPC();
+  }
+  else if(mode == kGlobal)
+    vertex = esd->GetPrimaryVertex();
+  else
+    Printf("GetVertex: ERROR: Invalid second argument %d", mode);
+  
+  if(!vertex) return 0;
+  
+  // check Ncontributors
+  if(vertex->GetNContributors() <= 0) return 0;
+  
+  // check resolution
+  Double_t zRes = vertex->GetZRes();
+  if(zRes == 0) return 0;
+  
+  //check position
+  if(TMath::Abs(vertex->GetXv()) > gVxMax) return 0;
+  if(TMath::Abs(vertex->GetYv()) > gVyMax) return 0;
+  if(TMath::Abs(vertex->GetZv()) > gVzMax) return 0;
+  
+  return vertex;
+}
+
 
index 73a4baacf0844c3de17a6c9fd5936e7a3f404183..5fe4c8cbd615520606b07a15c73a4abe69524c83 100644 (file)
@@ -15,6 +15,9 @@ class TF1;
 
 class AliAnalysisTaskProtons : public AliAnalysisTask {
  public:
+  enum TriggerMode { kMB1 = 0, kMB2, kSPDFASTOR }; 
+  enum AnalysisMode { kInvalid = -1, kTPC = 0, kHybrid, kGlobal };
+  
   AliAnalysisTaskProtons();
   AliAnalysisTaskProtons(const char *name);
   virtual ~AliAnalysisTaskProtons() {}
@@ -37,28 +40,45 @@ class AliAnalysisTaskProtons : public AliAnalysisTask {
     fKaonFunction = fkaons;
     fProtonFunction = fprotons;
   }
+
+  void SetTriggerMode(TriggerMode triggermode) {fTriggerMode = triggermode;}
+  void SetAnalysisMode(AnalysisMode analysismode) {fProtonAnalysisMode = analysismode;}
+  void SetAcceptedVertexDiamond(Double_t gVx, Double_t gVy, Double_t gVz) {
+    fVxMax = gVx; fVyMax = gVy; fVzMax = gVz;}
+
+  static Bool_t IsEventTriggered(const AliESDEvent *esd, 
+                                TriggerMode trigger = kMB2);
+  static const  AliESDVertex *GetVertex(AliESDEvent *esd, 
+                                       AnalysisMode mode = kTPC,
+                                       Double_t gVx = 100.,
+                                       Double_t gVy = 100.,
+                                       Double_t gVz = 100.);
   
  private:
-  AliESDEvent *fESD;    //ESD object                                                                                 
-  AliAODEvent *fAOD;    //AOD object                                                                                  
-  AliMCEvent  *fMC;     //MC object                                                                                   
+  AliESDEvent *fESD;    //ESD object 
+  AliAODEvent *fAOD;    //AOD object
+  AliMCEvent  *fMC;     //MC object 
   
-  TString fAnalysisType;//"ESD", "AOD" or "MC"                                                                        
+  TString fAnalysisType;//"ESD", "AOD" or "MC"
 
-  TList  *fList; //TList output object                                                                                
+  TList  *fList; //TList output object 
   
-  AliProtonAnalysis *fAnalysis; //analysis object                                                                     
+  AliProtonAnalysis *fProtonAnalysis; //analysis object 
   
-  TF1 *fElectronFunction; //TF1 for e                                                                                 
-  TF1 *fMuonFunction; //TF1 for mu                                                                                    
-  TF1 *fPionFunction; //TF1 for pi                                                                                    
-  TF1 *fKaonFunction; //TF1 for K                                                                                     
-  TF1 *fProtonFunction; //TF1 for p                                                                                   
+  TF1 *fElectronFunction; //TF1 for e
+  TF1 *fMuonFunction; //TF1 for mu
+  TF1 *fPionFunction; //TF1 for pi
+  TF1 *fKaonFunction; //TF1 for K
+  TF1 *fProtonFunction; //TF1 for p
   
-  Bool_t fFunctionUsed; //kTRUE if Functions are used                                                                 
+  Bool_t fFunctionUsed; //kTRUE if Functions are used
   
-  AliAnalysisTaskProtons(const AliAnalysisTaskProtons&); // not implemented                                           
-  AliAnalysisTaskProtons& operator=(const AliAnalysisTaskProtons&); // not implemented                                 
+  TriggerMode   fTriggerMode; //Trigger mode
+  AnalysisMode  fProtonAnalysisMode; //Analysis mode
+  Double_t      fVxMax, fVyMax, fVzMax; //vertex diamond constrain
+
+  AliAnalysisTaskProtons(const AliAnalysisTaskProtons&); // not implemented
+  AliAnalysisTaskProtons& operator=(const AliAnalysisTaskProtons&); // not implemented
   
   ClassDef(AliAnalysisTaskProtons, 1); // example of analysis
 };
index f7d61ecc103814c1e6a9c9b99afa01285e9fe9c7..69a8c42f96a7223a49d59d581ece78f2ea6e08b6 100644 (file)
@@ -99,7 +99,7 @@ void AliAnalysisTaskProtonsQA::CreateOutputObjects() {
 
   //Use of TPConly tracks
   if(fProtonAnalysisMode == kTPC) {
-    fProtonQAAnalysis->SetQAYPtBins(10, -0.5, 0.5, 12, 0.5, 0.9); //TPC only
+    fProtonQAAnalysis->SetQAYPtBins(10, -0.5, 0.5, 16, 0.5, 0.9); //TPC only
     fProtonQAAnalysis->UseTPCOnly();
     //fProtonQAAnalysis->SetTPCpid();
     fProtonQAAnalysis->SetMinTPCClusters(100);
@@ -115,7 +115,7 @@ void AliAnalysisTaskProtonsQA::CreateOutputObjects() {
   }
   //Use of HybridTPC tracks
   else if(fProtonAnalysisMode == kHybrid) {
-    fProtonQAAnalysis->SetQAYPtBins(10, -0.5, 0.5, 12, 0.5, 0.9); //HybridTPC
+    fProtonQAAnalysis->SetQAYPtBins(10, -0.5, 0.5, 16, 0.5, 0.9); //HybridTPC
     fProtonQAAnalysis->UseHybridTPC();
     fProtonQAAnalysis->SetTPCpid();
     fProtonQAAnalysis->SetMinTPCClusters(110);
@@ -207,9 +207,9 @@ void AliAnalysisTaskProtonsQA::Exec(Option_t *) {
     const AliESDVertex *vertex = GetVertex(fESD,fProtonAnalysisMode,
                                           fVxMax,fVyMax,fVzMax);
     if(vertex) {
-      fProtonQAAnalysis->RunQAAnalysis(stack, fESD);
+      fProtonQAAnalysis->RunQAAnalysis(stack, fESD, vertex);
       fProtonQAAnalysis->RunMCAnalysis(stack);
-      fProtonQAAnalysis->RunEfficiencyAnalysis(stack, fESD);
+      fProtonQAAnalysis->RunEfficiencyAnalysis(stack, fESD, vertex);
     }//accepted vertex
   }//triggered event
   
index ce0c90ebfbd317ee9735b4f5384efaac0ae86143..e3887b80e4cfea226a43942bfe327fbd4218736f 100644 (file)
@@ -38,6 +38,7 @@
 #include <AliCFContainer.h>
 #include <AliCFEffGrid.h>
 #include <AliCFDataGrid.h>
+#include <AliESDVertex.h>
 
 ClassImp(AliProtonAnalysis)
 
@@ -514,14 +515,15 @@ Double_t AliProtonAnalysis::GetParticleFraction(Int_t i, Double_t p) {
 }
 
 //____________________________________________________________________//
-void AliProtonAnalysis::Analyze(AliESDEvent* fESD) {
+void AliProtonAnalysis::Analyze(AliESDEvent* esd,
+                               const AliESDVertex *vertex) {
   //Main analysis part - ESD
   fHistEvents->Fill(0); //number of analyzed events
   Double_t containerInput[2] ;
   Double_t Pt = 0.0, P = 0.0;
-  Int_t nGoodTracks = fESD->GetNumberOfTracks();
+  Int_t nGoodTracks = esd->GetNumberOfTracks();
   for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) {
-    AliESDtrack* track = fESD->GetTrack(iTracks);
+    AliESDtrack* track = esd->GetTrack(iTracks);
     Double_t probability[5];
     AliESDtrack trackTPC;
 
@@ -530,14 +532,14 @@ void AliProtonAnalysis::Analyze(AliESDEvent* fESD) {
       Float_t p[2],cov[3];
       track->GetImpactParametersTPC(p,cov);
       if (p[0]==0 && p[1]==0)  
-       track->RelateToVertexTPC(((AliESDEvent*)fESD)->GetPrimaryVertexTPC(),fESD->GetMagneticField(),kVeryBig);
+       track->RelateToVertexTPC(((AliESDEvent*)esd)->GetPrimaryVertexTPC(),esd->GetMagneticField(),kVeryBig);
       if (!track->FillTPCOnlyTrack(trackTPC)) {
        continue;
       }
       track = &trackTPC ;
     }
 
-    if(IsAccepted(track)) {    
+    if(IsAccepted(esd,vertex,track)) { 
       if(fUseTPCOnly) {
        AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
        if(!tpcTrack) continue;
@@ -680,23 +682,45 @@ void AliProtonAnalysis::Analyze(AliStack* stack) {
 }
 
 //____________________________________________________________________//
-Bool_t AliProtonAnalysis::IsAccepted(AliESDtrack* track) {
+Bool_t AliProtonAnalysis::IsAccepted(AliESDEvent *esd,
+                                    const AliESDVertex *vertex, 
+                                    AliESDtrack* track) {
   // Checks if the track is excluded from the cuts
   Double_t Pt = 0.0, Px = 0.0, Py = 0.0, Pz = 0.0;
-  Float_t dcaXY = 0.0, dcaZ = 0.0;
-
-  if(fUseHybridTPC) {
+  Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};  //The impact parameters and their covariance.
+  
+  if((fUseTPCOnly)&&(!fUseHybridTPC)) {
+    AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
+    if(!tpcTrack) {
+      Pt = 0.0; Px = 0.0; Py = 0.0; Pz = 0.0;
+      dca[0] = -100.; dca[1] = -100.;
+      cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
+    }
+    else {
+      Pt = tpcTrack->Pt();
+      Px = tpcTrack->Px();
+      Py = tpcTrack->Py();
+      Pz = tpcTrack->Pz();
+      tpcTrack->PropagateToDCA(vertex,
+                              esd->GetMagneticField(),
+                              100.,dca,cov);
+    }
+  }
+  else if(fUseHybridTPC) {
      AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
     if(!tpcTrack) {
       Pt = 0.0; Px = 0.0; Py = 0.0; Pz = 0.0;
-      dcaXY = -100.0, dcaZ = -100.0;
+      dca[0] = -100.; dca[1] = -100.;
+      cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
     }
     else {
       Pt = tpcTrack->Pt();
       Px = tpcTrack->Px();
       Py = tpcTrack->Py();
       Pz = tpcTrack->Pz();
-      track->GetImpactParameters(dcaXY,dcaZ);
+      tpcTrack->PropagateToDCA(vertex,
+                              esd->GetMagneticField(),
+                              100.,dca,cov);
     }
   }
   else{
@@ -704,7 +728,9 @@ Bool_t AliProtonAnalysis::IsAccepted(AliESDtrack* track) {
     Px = track->Px();
     Py = track->Py();
     Pz = track->Pz();
-    track->GetImpactParameters(dcaXY,dcaZ);
+    track->PropagateToDCA(vertex,
+                         esd->GetMagneticField(),
+                         100.,dca,cov);
   }
      
   Int_t  fIdxInt[200];
@@ -756,15 +782,13 @@ Bool_t AliProtonAnalysis::IsAccepted(AliESDtrack* track) {
   if(fMaxSigmaToVertexTPCFlag)
     if(GetSigmaToVertex(track) > fMaxSigmaToVertexTPC) return kFALSE;
   if(fMaxDCAXYFlag) 
-    if(dcaXY > fMaxDCAXY) return kFALSE;
+    if(TMath::Abs(dca[0]) > fMaxDCAXY) return kFALSE;
   if(fMaxDCAXYTPCFlag) 
-    if(dcaXY > fMaxDCAXYTPC) return kFALSE;
+    if(TMath::Abs(dca[0]) > fMaxDCAXYTPC) return kFALSE;
     if(fMaxDCAZFlag) 
-    if(dcaZ > fMaxDCAZ) return kFALSE;
+    if(TMath::Abs(dca[1]) > fMaxDCAZ) return kFALSE;
   if(fMaxDCAZTPCFlag) 
-    if(dcaZ > fMaxDCAZTPC) return kFALSE;
-  if(fMaxDCAXYFlag) 
-    if(dcaXY > fMaxDCAXY) return kFALSE;
+    if(TMath::Abs(dca[1]) > fMaxDCAZTPC) return kFALSE;
   if(fMaxConstrainChi2Flag) {
     if(track->GetConstrainedChi2() > 0) 
       if(TMath::Log(track->GetConstrainedChi2()) > fMaxConstrainChi2) return kFALSE;
index ef6de10cdf4916cd7570cb9f43874b4c4f9bcbc3..440092756d2a0b689ba464a954899de09d402ddf 100644 (file)
@@ -31,6 +31,7 @@ class AliESDEvent;
 class AliESDtrack;
 class AliExternalTrackParam;
 class AliStack;
+class AliESDVertex;
 
 class AliProtonAnalysis : public TObject {
  public:
@@ -45,7 +46,8 @@ class AliProtonAnalysis : public TObject {
   void InitAnalysisHistograms(Int_t nbinsY, Float_t fLowY, Float_t fHighY,
                              Int_t nbinsPt, Float_t fLowPt, Float_t fHighPt);
   Bool_t ReadFromFile(const char* filename);
-  void Analyze(AliESDEvent *fESD);
+  void Analyze(AliESDEvent *fESD, 
+              const AliESDVertex *vertex);
   void Analyze(AliAODEvent *fAOD);
   void Analyze(AliStack *stack);
   
@@ -173,7 +175,9 @@ class AliProtonAnalysis : public TObject {
   AliProtonAnalysis(const AliProtonAnalysis&); // Not implemented
   AliProtonAnalysis& operator=(const AliProtonAnalysis&); // Not implemented
 
-  Bool_t   IsAccepted(AliESDtrack *track);
+  Bool_t   IsAccepted(AliESDEvent *esd,
+                     const AliESDVertex *vertex, 
+                     AliESDtrack *track);
   Float_t  GetSigmaToVertex(AliESDtrack* esdTrack); 
   Double_t Rapidity(Double_t Px, Double_t Py, Double_t Pz);
   
index 37f5ca152a2bba3cee2ee674fb85b88ca3816feb..d3d97abe3e53a602430910f215242c237ad24040 100644 (file)
@@ -35,6 +35,7 @@
 #include <AliLog.h>
 #include <AliPID.h>
 #include <AliStack.h>
+#include <AliESDVertex.h>
 
 ClassImp(AliProtonQAAnalysis)
 
@@ -134,37 +135,45 @@ Double_t AliProtonQAAnalysis::GetParticleFraction(Int_t i, Double_t p) {
 }
 
 //____________________________________________________________________//
-Bool_t AliProtonQAAnalysis::IsAccepted(AliESDtrack* track) {
+Bool_t AliProtonQAAnalysis::IsAccepted(AliESDEvent *esd,
+                                      const AliESDVertex *vertex, 
+                                      AliESDtrack* track) {
   // Checks if the track is excluded from the cuts
   Double_t Pt = 0.0, Px = 0.0, Py = 0.0, Pz = 0.0;
-  Float_t dcaXY = 0.0, dcaZ = 0.0;
-
+  Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};  //The impact parameters and their covariance.
+  
   if((fUseTPCOnly)&&(!fUseHybridTPC)) {
     AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
     if(!tpcTrack) {
       Pt = 0.0; Px = 0.0; Py = 0.0; Pz = 0.0;
-      dcaXY = -100.0, dcaZ = -100.0;
+      dca[0] = -100.; dca[1] = -100.;
+      cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
     }
     else {
       Pt = tpcTrack->Pt();
       Px = tpcTrack->Px();
       Py = tpcTrack->Py();
       Pz = tpcTrack->Pz();
-      track->GetImpactParametersTPC(dcaXY,dcaZ);
+      tpcTrack->PropagateToDCA(vertex,
+                              esd->GetMagneticField(),
+                              100.,dca,cov);
     }
   }
   else if(fUseHybridTPC) {
      AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
     if(!tpcTrack) {
       Pt = 0.0; Px = 0.0; Py = 0.0; Pz = 0.0;
-      dcaXY = -100.0, dcaZ = -100.0;
+      dca[0] = -100.; dca[1] = -100.;
+      cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
     }
     else {
       Pt = tpcTrack->Pt();
       Px = tpcTrack->Px();
       Py = tpcTrack->Py();
       Pz = tpcTrack->Pz();
-      track->GetImpactParameters(dcaXY,dcaZ);
+      tpcTrack->PropagateToDCA(vertex,
+                              esd->GetMagneticField(),
+                              100.,dca,cov);
     }
   }
   else{
@@ -172,7 +181,9 @@ Bool_t AliProtonQAAnalysis::IsAccepted(AliESDtrack* track) {
     Px = track->Px();
     Py = track->Py();
     Pz = track->Pz();
-    track->GetImpactParameters(dcaXY,dcaZ);
+    track->PropagateToDCA(vertex,
+                         esd->GetMagneticField(),
+                         100.,dca,cov);
   }
      
   Int_t  fIdxInt[200];
@@ -224,13 +235,13 @@ Bool_t AliProtonQAAnalysis::IsAccepted(AliESDtrack* track) {
   if(fMaxSigmaToVertexTPCFlag)
     if(GetSigmaToVertex(track) > fMaxSigmaToVertexTPC) return kFALSE;
   if(fMaxDCAXYFlag) 
-    if(TMath::Abs(dcaXY) > fMaxDCAXY) return kFALSE;
+    if(TMath::Abs(dca[0]) > fMaxDCAXY) return kFALSE;
   if(fMaxDCAXYTPCFlag) 
-    if(TMath::Abs(dcaXY) > fMaxDCAXYTPC) return kFALSE;
+    if(TMath::Abs(dca[0]) > fMaxDCAXYTPC) return kFALSE;
     if(fMaxDCAZFlag) 
-    if(TMath::Abs(dcaZ) > fMaxDCAZ) return kFALSE;
+    if(TMath::Abs(dca[1]) > fMaxDCAZ) return kFALSE;
   if(fMaxDCAZTPCFlag) 
-    if(TMath::Abs(dcaZ) > fMaxDCAZTPC) return kFALSE;
+    if(TMath::Abs(dca[1]) > fMaxDCAZTPC) return kFALSE;
   if(fMaxConstrainChi2Flag) {
     if(track->GetConstrainedChi2() > 0) 
       if(TMath::Log(track->GetConstrainedChi2()) > fMaxConstrainChi2) return kFALSE;
@@ -1166,13 +1177,12 @@ void AliProtonQAAnalysis::FillQA(AliESDtrack* track, AliStack *stack) {
 //____________________________________________________________________//
 Float_t AliProtonQAAnalysis::GetSigmaToVertex(AliESDtrack* esdTrack) {
   // Calculates the number of sigma to the vertex.
-  
   Float_t b[2];
   Float_t bRes[2];
   Float_t bCov[3];
   if((fUseTPCOnly)&&(!fUseHybridTPC))
     esdTrack->GetImpactParametersTPC(b,bCov);
-  else
+  else 
     esdTrack->GetImpactParameters(b,bCov);
   
   if (bCov[0]<=0 || bCov[2]<=0) {
@@ -2630,7 +2640,8 @@ void AliProtonQAAnalysis::InitQA() {
 
 //____________________________________________________________________//
 void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack, 
-                                               AliESDEvent *esd) {
+                                               AliESDEvent *esd,
+                                               const AliESDVertex *vertex) {
   //Runs the efficiency code
   //MC loop
   Int_t nMCProtons = 0, nESDProtons = 0;
@@ -2728,19 +2739,19 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
   
   }//MC loop
 
-  //ESD loop
+  //ESD track loop
   Int_t nGoodTracks = esd->GetNumberOfTracks();
   TArrayI labelArray(nGoodTracks);
   Int_t labelCounter = 0;
   for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) {
     AliESDtrack* track = esd->GetTrack(iTracks);
     if(!track) continue;
-
+    
     Int_t label = TMath::Abs(track->GetLabel());
     if(IsLabelUsed(labelArray,label)) continue;
     labelArray.AddAt(label,labelCounter);
     labelCounter += 1;
-
+    
     TParticle *particle = stack->Particle(label);
     if(!particle) continue;
     Int_t pdgcode = particle->GetPdgCode();
@@ -2755,7 +2766,7 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
       if(!tpcTrack) continue;
       Pt = tpcTrack->Pt();
       P = tpcTrack->P();
-
+      
       if((particle->Pt() > fMaxPt)||(particle->Pt() < fMinPt)) continue;
       if(fRunEfficiencyAnalysisEtaMode) {
        if((particle->Eta() > fMaxY)|| (particle->Eta() < fMinY)) continue;
@@ -2764,8 +2775,8 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
        if((Rapidity(particle->Px(),particle->Py(),particle->Pz()) > fMaxY)||(Rapidity(particle->Px(),particle->Py(),particle->Pz()) < fMinY)) continue;
       
       if(fUseCutsInEfficiency) 
-       if(!IsAccepted(track)) continue;
-
+       if(!IsAccepted(esd,vertex,track)) continue;
+      
       //reconstructed primary (anti)protons
       if(pdgcode == 2212) {
        if(fRunEfficiencyAnalysisEtaMode)
@@ -2793,7 +2804,7 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
          lPartMother = particle->GetFirstMother();
          TParticle *motherParticle = stack->Particle(lPartMother);
          if(motherParticle) motherPDGCode = motherParticle->GetPdgCode();
-
+         
          if((particle->GetUniqueID() == 4)&&(TMath::Abs(motherPDGCode) == 3122)) {
            if(fRunEfficiencyAnalysisEtaMode)
              ((TH2D *)(fEfficiencyList->At(8)))->Fill(particle->Eta(),
@@ -2841,7 +2852,7 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
          lPartMother = particle->GetFirstMother();
          TParticle *motherParticle = stack->Particle(lPartMother);
          if(motherParticle) motherPDGCode = motherParticle->GetPdgCode();
-
+         
          if((particle->GetUniqueID() == 4)&&(TMath::Abs(motherPDGCode) == 3122)) {
            if(fRunEfficiencyAnalysisEtaMode)
              ((TH2D *)(fEfficiencyList->At(9)))->Fill(particle->Eta(),
@@ -2916,18 +2927,18 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
       }
       
       if(fUseCutsInEfficiency) 
-       if(!IsAccepted(track)) continue;
-
+       if(!IsAccepted(esd,vertex,track)) continue;
+      
       //reconstructed primary (anti)protons
       if(pdgcode == 2212) {
        if(fRunEfficiencyAnalysisEtaMode)
          ((TH2D *)(fEfficiencyList->At(12)))->Fill(particle->Eta(),
-                                                  particle->Pt());
+                                                   particle->Pt());
        else
          ((TH2D *)(fEfficiencyList->At(12)))->Fill(Rapidity(particle->Px(),
-                                                           particle->Py(),
-                                                           particle->Pz()),
-                                                  particle->Pt());
+                                                            particle->Py(),
+                                                            particle->Pz()),
+                                                   particle->Pt());
        if(label <= stack->GetNprimary()) {
          nESDProtons += 1;
          if(fRunEfficiencyAnalysisEtaMode)
@@ -2945,7 +2956,7 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
          lPartMother = particle->GetFirstMother();
          TParticle *motherParticle = stack->Particle(lPartMother);
          if(motherParticle) motherPDGCode = motherParticle->GetPdgCode();
-
+         
          if((particle->GetUniqueID() == 4)&&(TMath::Abs(motherPDGCode) == 3122)) {
            if(fRunEfficiencyAnalysisEtaMode)
              ((TH2D *)(fEfficiencyList->At(8)))->Fill(particle->Eta(),
@@ -2971,12 +2982,12 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
       if(pdgcode == -2212) {   
        if(fRunEfficiencyAnalysisEtaMode)
          ((TH2D *)(fEfficiencyList->At(12)))->Fill(particle->Eta(),
-                                                  particle->Pt());
+                                                   particle->Pt());
        else
          ((TH2D *)(fEfficiencyList->At(12)))->Fill(Rapidity(particle->Px(),
-                                                           particle->Py(),
-                                                           particle->Pz()),
-                                                  particle->Pt());
+                                                            particle->Py(),
+                                                            particle->Pz()),
+                                                   particle->Pt());
        if(label <= stack->GetNprimary()) {
          if(fRunEfficiencyAnalysisEtaMode)
            ((TH2D *)(fEfficiencyList->At(7)))->Fill(particle->Eta(),
@@ -2993,7 +3004,7 @@ void AliProtonQAAnalysis::RunEfficiencyAnalysis(AliStack *stack,
          lPartMother = particle->GetFirstMother();
          TParticle *motherParticle = stack->Particle(lPartMother);
          if(motherParticle) motherPDGCode = motherParticle->GetPdgCode();
-
+         
          if((particle->GetUniqueID() == 4)&&(TMath::Abs(motherPDGCode) == 3122)) {
            if(fRunEfficiencyAnalysisEtaMode)
              ((TH2D *)(fEfficiencyList->At(9)))->Fill(particle->Eta(),
@@ -3077,7 +3088,8 @@ Bool_t AliProtonQAAnalysis::IsLabelUsed(TArrayI labelArray,
 
 //____________________________________________________________________//
 void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack, 
-                                       AliESDEvent *esd) {
+                                       AliESDEvent *esd,
+                                       const AliESDVertex *vertex) {
   //Runs the QA code
   //MC loop
   for(Int_t iParticle = 0; iParticle < stack->GetNprimary(); iParticle++) {
@@ -3101,7 +3113,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
                                         particle->Pt());
   }//MC loop
 
-  //ESD loop
+  //ESD track loop
   Int_t nGoodTracks = esd->GetNumberOfTracks();
   TArrayI labelArray(nGoodTracks);
   Int_t labelCounter = 0;
@@ -3113,8 +3125,9 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
     if(IsLabelUsed(labelArray,label)) continue;
     labelArray.AddAt(label,labelCounter);
     labelCounter += 1;
-
+    
     TParticle *particle = stack->Particle(label);
+    if(!particle) continue;
     if(TMath::Abs(particle->Eta()) > 1.0) continue;//acceptance
     
     AliESDtrack trackTPC;
@@ -3130,7 +3143,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
       }
       track = &trackTPC ;
     }
-
+    
     Double_t Pt = 0.0, P = 0.0;
     Double_t probability[5];
     Float_t dcaXY = 0.0, dcaZ = 0.0;
@@ -3138,7 +3151,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
     Int_t  fIdxInt[200];
     Int_t nClustersITS = track->GetITSclusters(fIdxInt);
     Int_t nClustersTPC = track->GetTPCclusters(fIdxInt);
-
+    
     Float_t chi2PerClusterITS = -1;
     if (nClustersITS!=0)
       chi2PerClusterITS = track->GetITSchi2()/Float_t(nClustersITS);
@@ -3148,7 +3161,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
     Double_t chi2ConstrainVertex = TMath::Log(track->GetConstrainedChi2());    
     Double_t extCov[15];
     track->GetExternalCovariance(extCov);
-
+    
     //TPC only
     if(fUseTPCOnly) {
       AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
@@ -3170,9 +3183,9 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
       Long64_t fParticleType = TMath::LocMax(AliPID::kSPECIES,w);
       if(fParticleType == 4) {
        FillQA(track, stack);
-       if(IsAccepted(track)) {
+       if(IsAccepted(esd,vertex,track)) {
          if(label <= stack->GetNprimary()) {
-            if(track->Charge() > 0) {
+           if(track->Charge() > 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
                  ((TH1F *)(fAcceptedCutList->At(0)))->Fill(iLayer+1);
@@ -3187,16 +3200,16 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(32)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(36)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(40)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(0)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(4)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(8)))->Fill(nSigmaToVertex);
-              ((TH2D *)(fQA2DList->At(0)))->Fill(Rapidity(tpcTrack->Px(),
+             ((TH2D *)(fQA2DList->At(0)))->Fill(Rapidity(tpcTrack->Px(),
                                                          tpcTrack->Py(),
                                                          tpcTrack->Pz()),
                                                 Pt);
            }
-            else if(track->Charge() < 0) {
+           else if(track->Charge() < 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
                  ((TH1F *)(fAcceptedCutList->At(1)))->Fill(iLayer+1);
@@ -3211,20 +3224,17 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(33)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(37)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(41)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(1)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(5)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(9)))->Fill(nSigmaToVertex);
-              ((TH2D *)(fQA2DList->At(4)))->Fill(Rapidity(tpcTrack->Px(),
+             ((TH2D *)(fQA2DList->At(4)))->Fill(Rapidity(tpcTrack->Px(),
                                                          tpcTrack->Py(),
                                                          tpcTrack->Pz()),
                                                 Pt);
            }
          }//primary particles
          else if(label > stack->GetNprimary()) {
-           TParticle *particle = stack->Particle(label);
-           if(!particle) continue;
-
            Int_t lPartMother = -1;
            Int_t motherPDGCode = -1;
            if(particle) {
@@ -3232,12 +3242,12 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              TParticle *motherParticle = stack->Particle(lPartMother);
              if(motherParticle) motherPDGCode = motherParticle->GetPdgCode();
            }
-
+           
            if(fMCProcessIdFlag)
              if(particle->GetUniqueID() != fMCProcessId) continue;
            if(fMotherParticlePDGCodeFlag)
              if(TMath::Abs(motherPDGCode) != fMotherParticlePDGCode) continue;
-
+           
            if(track->Charge() > 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
@@ -3253,7 +3263,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(34)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(38)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(42)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(2)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(6)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(10)))->Fill(nSigmaToVertex);
@@ -3267,7 +3277,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
                                                  Pt,
                                                  ConvertPDGToInt(motherPDGCode));
            }
-            else if(track->Charge() < 0) {
+           else if(track->Charge() < 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
                  ((TH1F *)(fAcceptedCutList->At(3)))->Fill(iLayer+1);
@@ -3282,11 +3292,11 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(35)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(39)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(43)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(3)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(7)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(11)))->Fill(nSigmaToVertex);
-              ((TH2D *)(fQA2DList->At(6)))->Fill(Rapidity(tpcTrack->Px(),
+             ((TH2D *)(fQA2DList->At(6)))->Fill(Rapidity(tpcTrack->Px(),
                                                          tpcTrack->Py(),
                                                          tpcTrack->Pz()),
                                                 Pt);
@@ -3298,27 +3308,27 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
            }
          }//secondary particles
        }//accepted - track cuts
-       else if(!IsAccepted(track)) {
+       else if(!IsAccepted(esd,vertex,track)) {
          if(label <= stack->GetNprimary()) {
-            if(track->Charge() > 0)
-              ((TH2D *)(fQA2DList->At(1)))->Fill(Rapidity(tpcTrack->Px(),
+           if(track->Charge() > 0)
+             ((TH2D *)(fQA2DList->At(1)))->Fill(Rapidity(tpcTrack->Px(),
                                                          tpcTrack->Py(),
                                                          tpcTrack->Pz()),
                                                 Pt);
-            else if(track->Charge() < 0)
-              ((TH2D *)(fQA2DList->At(5)))->Fill(Rapidity(tpcTrack->Px(),
+           else if(track->Charge() < 0)
+             ((TH2D *)(fQA2DList->At(5)))->Fill(Rapidity(tpcTrack->Px(),
                                                          tpcTrack->Py(),
                                                          tpcTrack->Pz()),
                                                 Pt);
          }//primary particles
          else if(label > stack->GetNprimary()) {
            if(track->Charge() > 0)
-              ((TH2D *)(fQA2DList->At(3)))->Fill(Rapidity(tpcTrack->Px(),
+             ((TH2D *)(fQA2DList->At(3)))->Fill(Rapidity(tpcTrack->Px(),
                                                          tpcTrack->Py(),
                                                          tpcTrack->Pz()),
                                                 Pt);
-            else if(track->Charge() < 0)
-              ((TH2D *)(fQA2DList->At(7)))->Fill(Rapidity(tpcTrack->Px(),
+           else if(track->Charge() < 0)
+             ((TH2D *)(fQA2DList->At(7)))->Fill(Rapidity(tpcTrack->Px(),
                                                          tpcTrack->Py(),
                                                          tpcTrack->Pz()),
                                                 Pt);
@@ -3331,7 +3341,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
       Pt = track->Pt();
       P = track->P();
       track->GetImpactParameters(dcaXY,dcaZ);
-
+      
       //pid
       track->GetESDpid(probability);
       Double_t rcc = 0.0;
@@ -3344,9 +3354,9 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
       Long64_t fParticleType = TMath::LocMax(AliPID::kSPECIES,w);
       if(fParticleType == 4) {
        FillQA(track, stack);
-       if(IsAccepted(track)) {
+       if(IsAccepted(esd,vertex,track)) {
          if(label <= stack->GetNprimary()) {
-            if(track->Charge() > 0) {
+           if(track->Charge() > 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
                  ((TH1F *)(fAcceptedCutList->At(0)))->Fill(iLayer+1);
@@ -3361,16 +3371,16 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(32)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(36)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(40)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(0)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(4)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(8)))->Fill(nSigmaToVertex);
-              ((TH2D *)(fQA2DList->At(0)))->Fill(Rapidity(track->Px(),
+             ((TH2D *)(fQA2DList->At(0)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
            }
-            else if(track->Charge() < 0) {
+           else if(track->Charge() < 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
                  ((TH1F *)(fAcceptedCutList->At(1)))->Fill(iLayer+1);
@@ -3385,20 +3395,17 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(33)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(37)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(41)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(1)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(5)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(9)))->Fill(nSigmaToVertex);
-              ((TH2D *)(fQA2DList->At(4)))->Fill(Rapidity(track->Px(),
+             ((TH2D *)(fQA2DList->At(4)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
            }
          }//primary particles
          else if(label > stack->GetNprimary()) {
-           TParticle *particle = stack->Particle(label);
-           if(!particle) continue;
-
            Int_t lPartMother = -1;
            Int_t motherPDGCode = -1;
            if(particle) {
@@ -3406,12 +3413,12 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              TParticle *motherParticle = stack->Particle(lPartMother);
              if(motherParticle) motherPDGCode = motherParticle->GetPdgCode();
            }
-
+           
            if(fMCProcessIdFlag)
              if(particle->GetUniqueID() != fMCProcessId) continue;
            if(fMotherParticlePDGCodeFlag)
              if(TMath::Abs(motherPDGCode) != fMotherParticlePDGCode) continue;
-
+           
            if(track->Charge() > 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
@@ -3427,11 +3434,11 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(34)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(38)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(42)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(2)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(6)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(10)))->Fill(nSigmaToVertex);
-              ((TH2D *)(fQA2DList->At(2)))->Fill(Rapidity(track->Px(),
+             ((TH2D *)(fQA2DList->At(2)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
@@ -3441,7 +3448,7 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
                                                  Pt,
                                                  ConvertPDGToInt(motherPDGCode));
            }
-            else if(track->Charge() < 0) {
+           else if(track->Charge() < 0) {
              for(Int_t iLayer = 0; iLayer < 6; iLayer++) {
                if(track->HasPointOnITSLayer(iLayer))
                  ((TH1F *)(fAcceptedCutList->At(3)))->Fill(iLayer+1);
@@ -3456,44 +3463,44 @@ void AliProtonQAAnalysis::RunQAAnalysis(AliStack *stack,
              ((TH1F *)(fAcceptedCutList->At(35)))->Fill(extCov[5]);
              ((TH1F *)(fAcceptedCutList->At(39)))->Fill(extCov[9]);
              ((TH1F *)(fAcceptedCutList->At(43)))->Fill(extCov[14]);
-
+             
              ((TH1F *)(fAcceptedDCAList->At(3)))->Fill(TMath::Abs(dcaXY));
              ((TH1F *)(fAcceptedDCAList->At(7)))->Fill(TMath::Abs(dcaZ));
              ((TH1F *)(fAcceptedDCAList->At(11)))->Fill(nSigmaToVertex);
-              ((TH2D *)(fQA2DList->At(6)))->Fill(Rapidity(track->Px(),
+             ((TH2D *)(fQA2DList->At(6)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
-
+             
              ((TH3F *)(fQA2DList->At(11)))->Fill(Rapidity(track->Px(),
-                                                         track->Py(),
-                                                         track->Pz()),
-                                                Pt,
-                                                ConvertPDGToInt(motherPDGCode));
+                                                          track->Py(),
+                                                          track->Pz()),
+                                                 Pt,
+                                                 ConvertPDGToInt(motherPDGCode));
            }
          }//secondary particles
        }//accepted - track cuts
-       else if(!IsAccepted(track)) {
+       else if(!IsAccepted(esd,vertex,track)) {
          if(label <= stack->GetNprimary()) {
-            if(track->Charge() > 0)
-              ((TH2D *)(fQA2DList->At(1)))->Fill(Rapidity(track->Px(),
+           if(track->Charge() > 0)
+             ((TH2D *)(fQA2DList->At(1)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
-            else if(track->Charge() < 0)
-              ((TH2D *)(fQA2DList->At(5)))->Fill(Rapidity(track->Px(),
+           else if(track->Charge() < 0)
+             ((TH2D *)(fQA2DList->At(5)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
          }//primary particles
          else if(label > stack->GetNprimary()) {
            if(track->Charge() > 0)
-              ((TH2D *)(fQA2DList->At(3)))->Fill(Rapidity(track->Px(),
+             ((TH2D *)(fQA2DList->At(3)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
-            else if(track->Charge() < 0)
-              ((TH2D *)(fQA2DList->At(7)))->Fill(Rapidity(track->Px(),
+           else if(track->Charge() < 0)
+             ((TH2D *)(fQA2DList->At(7)))->Fill(Rapidity(track->Px(),
                                                          track->Py(),
                                                          track->Pz()),
                                                 Pt);
@@ -3747,6 +3754,50 @@ Int_t AliProtonQAAnalysis::ConvertPDGToInt(Int_t pdgCode) {
   return code;
 }
 
+//________________________________________________________________________
+/*const AliESDVertex* AliProtonQAAnalysis::GetVertex(AliESDEvent* esd,
+                                                  Double_t gVxMax,
+                                                  Double_t gVyMax,
+                                                  Double_t gVzMax) {
+  // Get the vertex from the ESD and returns it if the vertex is valid
+  // depending on the analysis mode: TPC - Hybrid - Global
+  const AliESDVertex* vertex = 0;
+  if((fUseTPCOnly)&&(!fUseHybridTPC)) {
+    Double_t kBz = esd->GetMagneticField();
+    AliVertexerTracks vertexer(kBz);
+    vertexer.SetTPCMode();
+    AliESDVertex *vTPC = vertexer.FindPrimaryVertex(esd);
+    esd->SetPrimaryVertexTPC(vTPC);
+    for (Int_t i=0; i<esd->GetNumberOfTracks(); i++) {
+      AliESDtrack *t = esd->GetTrack(i);
+      t->RelateToVertexTPC(vTPC, kBz, kVeryBig);
+    }
+    delete vTPC;
+    vertex = esd->GetPrimaryVertexTPC();
+  }
+  else if(fUseHybridTPC)
+    vertex = esd->GetPrimaryVertexSPD();
+  else if(!fUseTPCOnly)
+    vertex = esd->GetPrimaryVertex();
+  else
+    Printf("GetVertex: ERROR: Invalid analysis mode");
+  
+  if(!vertex) return 0;
+  
+  // check Ncontributors
+  if(vertex->GetNContributors() <= 0) return 0;
+  
+  // check resolution
+  Double_t zRes = vertex->GetZRes();
+  if(zRes == 0) return 0;
+  
+  //check position
+  if(TMath::Abs(vertex->GetXv()) > gVxMax) return 0;
+  if(TMath::Abs(vertex->GetYv()) > gVyMax) return 0;
+  if(TMath::Abs(vertex->GetZv()) > gVzMax) return 0;
+
+  return vertex;
+  }*/
 
 
 
index 33c9aa2836cb8bfdd9cc8e1a7d5958f88709f371..a978146d9b5297fddbbc104bfe195523694c0499 100644 (file)
@@ -27,6 +27,7 @@ class TH3F;
 class AliESDEvent;
 class AliESDtrack;
 class AliStack;
+class AliESDVertex;
 
 class AliProtonQAAnalysis : public TObject {
  public:
@@ -35,7 +36,7 @@ class AliProtonQAAnalysis : public TObject {
 
   void UseTPCOnly() {fUseTPCOnly = kTRUE;}
   void UseHybridTPC() {fUseTPCOnly = kTRUE; fUseHybridTPC = kTRUE;}
-  
+
   //Cut functions
   void    SetPointOnITSLayer1() {fPointOnITSLayer1Flag = kTRUE;}
   void    SetPointOnITSLayer2() {fPointOnITSLayer2Flag = kTRUE;}
@@ -118,24 +119,22 @@ class AliProtonQAAnalysis : public TObject {
   //QA histograms
   void SetQAYPtBins(Int_t nbinsY, Double_t minY, Double_t maxY,
                    Int_t nbinsPt, Double_t minPt, Double_t maxPt);
-  void RunQAAnalysis(AliStack *stack, AliESDEvent *esd);
+  void RunQAAnalysis(AliStack *stack, 
+                    AliESDEvent *esd,
+                    const AliESDVertex *vertex);
   void SetRunQAAnalysis();
   TList *GetGlobalQAList() {return fGlobalQAList;}
 
   //Efficiency plots (reconstruction & PID)
-  void RunEfficiencyAnalysis(AliStack *stack, AliESDEvent *esd);
+  void RunEfficiencyAnalysis(AliStack *stack, 
+                            AliESDEvent *esd,
+                            const AliESDVertex *vertex);
   void SetRunEfficiencyAnalysis(Bool_t gEtaMode, Bool_t gUseCuts) {
     fRunEfficiencyAnalysis = kTRUE;
     fRunEfficiencyAnalysisEtaMode = gEtaMode;
     fUseCutsInEfficiency = gUseCuts;
   }
   TList *GetEfficiencyQAList() {return fEfficiencyList;}
-  /*TH1F *GetReconstructionEfficiency(const char *variable, 
-                                   const char *particle);
-  TH1F *GetPIDEfficiencyY();
-  TH1F *GetPIDEfficiencyPt();
-  TH1F *GetPIDContaminationY();
-  TH1F *GetPIDContaminationPt();*/
 
   //MC analysis
   void RunMCAnalysis(AliStack* stack);
@@ -160,7 +159,9 @@ class AliProtonQAAnalysis : public TObject {
   AliProtonQAAnalysis(const AliProtonQAAnalysis&); // Not implemented
   AliProtonQAAnalysis& operator=(const AliProtonQAAnalysis&);// Not implemented
 
-  Bool_t   IsAccepted(AliESDtrack *track);
+  Bool_t   IsAccepted(AliESDEvent *esd,
+                     const AliESDVertex *vertex, 
+                     AliESDtrack *track);
 
   void     FillQA(AliESDtrack *track, AliStack *stack);
 
index b592b42bed7ca5cf14d7a3554db621126ac5dc0e..ab0c6eda8cacc151b4d5db6898f08d5f16ee79de 100644 (file)
@@ -2,9 +2,9 @@ void runProtonAnalysis() {
   TStopwatch timer;
   timer.Start();
   
-  //runLocal();
-  //runInteractive();
-  //runBatch();
+  //runLocal("ESD");
+  //runInteractive("ESD");
+  //runBatch("ESD");
   
   runProof("ESD",200000,"/COMMON/COMMON/LHC08c11_10TeV_0.5T"); //use data sets
   //runProof("ESD",200); //use ascii files
@@ -14,7 +14,7 @@ void runProtonAnalysis() {
 }
 
 //_________________________________________________//
-void runLocal() {
+void runLocal(const char* mode = "ESD") {
   TStopwatch timer;
   timer.Start();
   gSystem->Load("libTree.so");
@@ -81,9 +81,19 @@ void runLocal() {
   AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
   AliVEventHandler* esdH = new AliESDInputHandler;
   mgr->SetInputEventHandler(esdH);  
+  if(smode == "MC") {
+    AliMCEventHandler *mc = new AliMCEventHandler();
+    mgr->SetMCtruthEventHandler(mc);
+  }
+
   //____________________________________________//
   // 1st Proton task
   AliAnalysisTaskProtons *taskProtons = new AliAnalysisTaskProtons("TaskProtons");
+  taskProtons->SetType(mode);
+  taskProtons->SetTriggerMode(AliAnalysisTaskProtonsQA::kMB2);
+  taskProtons->SetAnalysisMode(AliAnalysisTaskProtonsQA::kTPC);
+  taskProtons->SetAcceptedVertexDiamond(5.,5.,15.);
+  //Momentum dependent priors
   /*TFile *f = TFile::Open("PriorProb/PriorProbabilities.root ");
   TF1 *fitElectrons = (TF1 *)f->Get("fitElectrons");
   TF1 *fitMuons = (TF1 *)f->Get("fitMuons");
@@ -118,7 +128,7 @@ void runLocal() {
 }
 
 //_________________________________________________//
-void runInteractive() {
+void runInteractive(const char* mode = "ESD") {
   TStopwatch timer;
   timer.Start();
   gSystem->Load("libProofPlayer.so");
@@ -192,9 +202,19 @@ void runInteractive() {
   AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
   AliVEventHandler* esdH = new AliESDInputHandler;
   mgr->SetInputEventHandler(esdH);  
+  if(smode == "MC") {
+    AliMCEventHandler *mc = new AliMCEventHandler();
+    mgr->SetMCtruthEventHandler(mc);
+  }
+
   //____________________________________________//
   // 1st Proton task
   AliAnalysisTaskProtons *taskProtons = new AliAnalysisTaskProtons("TaskProtons");
+  taskProtons->SetType(mode);
+  taskProtons->SetTriggerMode(AliAnalysisTaskProtonsQA::kMB2);
+  taskProtons->SetAnalysisMode(AliAnalysisTaskProtonsQA::kTPC);
+  taskProtons->SetAcceptedVertexDiamond(5.,5.,15.);
+  //Momentum dependent priors
   /*TFile *f = TFile::Open("PriorProb/PriorProbabilities.root ");
   TF1 *fitElectrons = (TF1 *)f->Get("fitElectrons");
   TF1 *fitMuons = (TF1 *)f->Get("fitMuons");
@@ -229,7 +249,7 @@ void runInteractive() {
 }
 
 //_________________________________________________//
-void runBatch() {
+void runBatch(const char* mode = "ESD") {
   TStopwatch timer;
   timer.Start();
 
@@ -296,9 +316,19 @@ void runBatch() {
   AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
   AliVEventHandler* esdH = new AliESDInputHandler;
   mgr->SetInputEventHandler(esdH);  
+  if(smode == "MC") {
+    AliMCEventHandler *mc = new AliMCEventHandler();
+    mgr->SetMCtruthEventHandler(mc);
+  }
+  
   //____________________________________________//
   // 1st Proton task
   AliAnalysisTaskProtons *taskProtons = new AliAnalysisTaskProtons("TaskProtons");
+  taskProtons->SetType(mode);
+  taskProtons->SetTriggerMode(AliAnalysisTaskProtonsQA::kMB2);
+  taskProtons->SetAnalysisMode(AliAnalysisTaskProtonsQA::kTPC);
+  taskProtons->SetAcceptedVertexDiamond(5.,5.,15.);
+  //Momentum dependent priors
   /*TFile *f = TFile::Open("PriorProb/PriorProbabilities.root ");
   TF1 *fitElectrons = (TF1 *)f->Get("fitElectrons");
   TF1 *fitMuons = (TF1 *)f->Get("fitMuons");
@@ -376,6 +406,10 @@ void runProof(const char* mode = "ESD",
   // 1st Proton task
   AliAnalysisTaskProtons *taskProtons = new AliAnalysisTaskProtons("TaskProtons");
   taskProtons->SetType(mode);
+  taskProtons->SetTriggerMode(AliAnalysisTaskProtonsQA::kMB2);
+  taskProtons->SetAnalysisMode(AliAnalysisTaskProtonsQA::kTPC);
+  taskProtons->SetAcceptedVertexDiamond(5.,5.,15.);
+  //Momentum dependent priors
   /*TFile *f = TFile::Open("PriorProb/PriorProbabilities.root ");
   TF1 *fitElectrons = (TF1 *)f->Get("fitElectrons");
   TF1 *fitMuons = (TF1 *)f->Get("fitMuons");