AOD analysis for protons - analysis train
authorpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 8 Feb 2008 14:29:00 +0000 (14:29 +0000)
committerpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 8 Feb 2008 14:29:00 +0000 (14:29 +0000)
PWG2/SPECTRA/AliProtonAnalysis.cxx
PWG2/SPECTRA/AliProtonAnalysis.h

index a127d53..275df00 100644 (file)
 
 #include "AliProtonAnalysis.h"
 
+#include <AliAODEvent.h>
 #include <AliESDEvent.h>
 #include <AliLog.h>
+#include <AliPID.h>
 
 ClassImp(AliProtonAnalysis)
 
@@ -299,6 +301,31 @@ void AliProtonAnalysis::Analyze(AliESDEvent* fESD) {
 }
 
 //____________________________________________________________________//
+void AliProtonAnalysis::Analyze(AliAODEvent* fAOD) {
+  //Main analysis part
+  Int_t nGoodTracks = fAOD->GetNumberOfTracks();
+  for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) {
+    AliAODTrack* track = fAOD->GetTrack(iTracks);
+    Double_t Pt = track->Pt();
+    Double_t P = track->P();
+    
+    //pid
+    Double_t probability[10];
+    track->GetPID(probability);
+    Double_t rcc = 0.0;
+    for(Int_t i = 0; i < AliPID::kSPECIESN; i++) rcc += probability[i]*GetParticleFraction(i,P);
+    if(rcc == 0.0) continue;
+    Double_t w[5];
+    for(Int_t i = 0; i < AliPID::kSPECIESN; i++) w[i] = probability[i]*GetParticleFraction(i,P)/rcc;
+    Long64_t fParticleType = TMath::LocMax(AliPID::kSPECIESN,w);
+    if(fParticleType == 4) {
+      if(track->Charge() > 0) fHistYPtProtons->Fill(track->Y(fParticleType),Pt);
+      else if(track->Charge() < 0) fHistYPtAntiProtons->Fill(track->Y(fParticleType),Pt);
+    }//proton check
+  }//track loop 
+}
+
+//____________________________________________________________________//
 Bool_t AliProtonAnalysis::IsAccepted(AliESDtrack* track) {
   // Checks if the track is excluded from the cuts
   Int_t  fIdxInt[200];
index fbaca56..f910e20 100644 (file)
 //-------------------------------------------------------------------------
 
 #include <TObject.h>
+#include "AliPID.h"
 
 class TF1;
 class TH2F;
 class TH1D;
+class AliAODEvent;
+class AliAODtrack;
 class AliESDEvent;
 class AliESDtrack;
 
@@ -32,6 +35,7 @@ class AliProtonAnalysis : public TObject {
                      Int_t nbinsPt, Float_t fLowPt, Float_t fHighPt);
   void ReadFromFile(const char* filename);
   void Analyze(AliESDEvent* fESD);
+  void Analyze(AliAODEvent* fAOD);
   
   TH2F *GetProtonYPtHistogram() {return fHistYPtProtons;}
   TH2F *GetAntiProtonYPtHistogram() {return fHistYPtAntiProtons;}
@@ -74,7 +78,7 @@ class AliProtonAnalysis : public TObject {
   void    SetTPCRefit() {fTPCRefitFlag = kTRUE;}
   
   //Prior probabilities
-  void    SetPriorProbabilities(Double_t *partFrac) {for(Int_t i = 0; i < 5; i++) fPartFrac[i] = partFrac[i];} 
+  void    SetPriorProbabilities(Double_t *partFrac) {for(Int_t i = 0; i < AliPID::kSPECIESN; i++) fPartFrac[i] = partFrac[i];} 
   void    SetPriorProbabilityFunctions(TF1 *felectron, TF1 *fmuon, TF1 *fpion, TF1 *fkaon, TF1 *fproton) {
     fFunctionProbabilityFlag = kTRUE;
     fElectronFunction = felectron; 
@@ -108,7 +112,7 @@ class AliProtonAnalysis : public TObject {
   
   //pid
   Bool_t fFunctionProbabilityFlag; //flag: kTRUE if functions used
-  Double_t fPartFrac[5]; //prior probabilities
+  Double_t fPartFrac[10]; //prior probabilities
   TF1  *fElectronFunction; //momentum dependence of the prior probs
   TF1  *fMuonFunction; //momentum dependence of the prior probs
   TF1  *fPionFunction; //momentum dependence of the prior probs