totEt updates from Christine
authordsilverm <dsilverm@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 2 Sep 2010 10:07:56 +0000 (10:07 +0000)
committerdsilverm <dsilverm@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 2 Sep 2010 10:07:56 +0000 (10:07 +0000)
14 files changed:
PWG4/PWG4totEtLinkDef.h
PWG4/libPWG4totEt.pkg
PWG4/totEt/AliAnalysisEt.cxx
PWG4/totEt/AliAnalysisEt.h
PWG4/totEt/AliAnalysisEtCuts.h
PWG4/totEt/AliAnalysisEtMonteCarlo.cxx
PWG4/totEt/AliAnalysisHadEt.cxx [new file with mode: 0644]
PWG4/totEt/AliAnalysisHadEt.h [new file with mode: 0644]
PWG4/totEt/AliAnalysisHadEtMonteCarlo.cxx [new file with mode: 0644]
PWG4/totEt/AliAnalysisHadEtMonteCarlo.h [new file with mode: 0644]
PWG4/totEt/AliAnalysisHadEtReconstructed.cxx [new file with mode: 0644]
PWG4/totEt/AliAnalysisHadEtReconstructed.h [new file with mode: 0644]
PWG4/totEt/AliAnalysisTaskHadEt.cxx [new file with mode: 0644]
PWG4/totEt/AliAnalysisTaskHadEt.h [new file with mode: 0644]

index b02f06e..9b7bd56 100644 (file)
 #pragma link C++ class AliAnalysisEtReconstructedEmcal+;
 #pragma link C++ class AliAnalysisEtReconstructedPhos+;
 #pragma link C++ class AliAnalysisTaskTotEt+;
+
+#pragma link C++ class AliAnalysisHadEt+;
+#pragma link C++ class AliAnalysisHadEtMonteCarlo+;
+#pragma link C++ class AliAnalysisHadEtReconstructed+;
+#pragma link C++ class AliAnalysisTaskHadEt+;
+
       
 #endif
index 261a957..e4ff165 100644 (file)
@@ -8,6 +8,10 @@ totEt/AliAnalysisEtReconstructed.cxx \
 totEt/AliAnalysisEtReconstructedEmcal.cxx \
 totEt/AliAnalysisEtReconstructedPhos.cxx \
 totEt/AliAnalysisTaskTotEt.cxx \
+totEt/AliAnalysisHadEt.cxx \
+totEt/AliAnalysisHadEtMonteCarlo.cxx \
+totEt/AliAnalysisHadEtReconstructed.cxx \
+totEt/AliAnalysisTaskHadEt.cxx \
 
 HDRS:= $(SRCS:.cxx=.h) 
 
index 8072c72..cc17d08 100644 (file)
@@ -5,12 +5,38 @@
 #include "TH2F.h"
 #include <iostream>
 #include "AliAnalysisEtCuts.h"
+#include "AliVEvent.h"
+
+using namespace std;
 ClassImp(AliAnalysisEt);
 
 
 AliAnalysisEt::AliAnalysisEt() :
         fHistogramNameSuffix("")
         ,fPdgDB(0)
+        ,PiPlusCode(0)
+        ,PiMinusCode(0)
+        ,KPlusCode(0)
+        ,KMinusCode(0)
+        ,ProtonCode(0)
+        ,AntiProtonCode(0)
+        ,LambdaCode(0)
+        ,AntiLambdaCode(0)
+        ,K0SCode(0)
+        ,OmegaCode(0)
+        ,AntiOmegaCode(0)
+        ,Xi0Code(0)
+        ,AntiXi0Code(0)
+        ,XiCode(0)
+        ,AntiXiCode(0)
+        ,SigmaCode(0)
+        ,AntiSigmaCode(0)
+        ,K0LCode(0)
+        ,NeutronCode(0)
+        ,AntiNeutronCode(0)
+        ,EPlusCode(0)
+        ,EMinusCode(0)
+        ,PionMass(0)
         ,fSumEt(0)
         ,fSumEtAcc(0)
         ,fTotEt(0)
@@ -95,7 +121,32 @@ void AliAnalysisEt::FillOutputList(TList *list)
 
 void AliAnalysisEt::Init()
 {
-    fPdgDB = new TDatabasePDG();
+
+  if(!fPdgDB) fPdgDB = new TDatabasePDG();
+  PionMass = fPdgDB->GetParticle("pi+")->Mass();
+  PiPlusCode = fPdgDB->GetParticle("pi+")->PdgCode();
+    PiMinusCode = fPdgDB->GetParticle("pi-")->PdgCode();
+    KPlusCode = fPdgDB->GetParticle("K+")->PdgCode();
+    KMinusCode = fPdgDB->GetParticle("K-")->PdgCode();
+    ProtonCode = fPdgDB->GetParticle("proton")->PdgCode();
+    AntiProtonCode = fPdgDB->GetParticle("antiproton")->PdgCode();
+    LambdaCode = fPdgDB->GetParticle("Lambda0")->PdgCode();
+    AntiLambdaCode = fPdgDB->GetParticle("Lambda0_bar")->PdgCode();
+    K0SCode = fPdgDB->GetParticle("K_S0")->PdgCode();
+    OmegaCode = fPdgDB->GetParticle("Omega-")->PdgCode();
+    AntiOmegaCode = fPdgDB->GetParticle("Omega+")->PdgCode();
+    Xi0Code = fPdgDB->GetParticle("Xi0")->PdgCode();
+    AntiXi0Code = fPdgDB->GetParticle("Xi0_bar")->PdgCode();
+    XiCode = fPdgDB->GetParticle("Xi-")->PdgCode();
+    AntiXiCode = fPdgDB->GetParticle("Xi-_bar")->PdgCode();
+    SigmaCode = fPdgDB->GetParticle("Sigma-")->PdgCode();
+    AntiSigmaCode = fPdgDB->GetParticle("Sigma+")->PdgCode();
+    K0LCode = fPdgDB->GetParticle("K_L0")->PdgCode();
+    NeutronCode = fPdgDB->GetParticle("neutron")->PdgCode();
+    AntiNeutronCode = fPdgDB->GetParticle("antineutron")->PdgCode();
+    EPlusCode = fPdgDB->GetParticle("e+")->PdgCode();
+    EMinusCode = fPdgDB->GetParticle("e-")->PdgCode();
+
 }
 
 void AliAnalysisEt::CreateHistograms()
@@ -208,6 +259,12 @@ void AliAnalysisEt::FillHistograms()
     fHistTMDeltaR;
     */
 }
+Int_t AliAnalysisEt::AnalyseEvent(AliVEvent *event)
+{
+  //this line is basically here to eliminate a compiler warning that event is not used.  Making it a virtual function did not work with the plugin.
+  cout<<"This event has "<<event->GetNumberOfTracks()<<" tracks"<<endl;
+  return 0;
+}
 
 void AliAnalysisEt::ResetEventValues()
 {
index 3aa941c..6087334 100644 (file)
@@ -25,13 +25,13 @@ private:
 public:
   
     /** Analyse the event! */
-    virtual Int_t AnalyseEvent(AliVEvent *event) = 0;
+    virtual Int_t AnalyseEvent(AliVEvent *event);
 
     /** Fill the objects you want to output, classes which add new histograms should overload this. */
     virtual void FillOutputList(TList* list);
 
     /** Initialise the analysis, must be overloaded. */
-    virtual void Init() = 0;
+    virtual void Init();
 
     /** 
     * Creates the histograms, must be overloaded if you want to add your own. 
@@ -78,6 +78,30 @@ protected:
     /** PDG Database */
     TDatabasePDG *fPdgDB;
 
+    Int_t PiPlusCode;
+    Int_t PiMinusCode;
+    Int_t KPlusCode;
+    Int_t KMinusCode;
+    Int_t ProtonCode;
+    Int_t AntiProtonCode;
+    Int_t LambdaCode;
+    Int_t AntiLambdaCode;
+    Int_t K0SCode;
+    Int_t OmegaCode;
+    Int_t AntiOmegaCode;
+    Int_t Xi0Code;
+    Int_t AntiXi0Code;
+    Int_t XiCode;
+    Int_t AntiXiCode;
+    Int_t SigmaCode;
+    Int_t AntiSigmaCode;
+    Int_t K0LCode;
+    Int_t NeutronCode;
+    Int_t AntiNeutronCode;
+    Int_t EPlusCode;
+    Int_t EMinusCode;
+    Float_t PionMass;
+
     /** Sum of the total Et for all events */
     Double_t fSumEt;
 
index d579b0f..f2c8720 100644 (file)
@@ -6,7 +6,7 @@
 namespace EtCommonCuts
 {
    
-   const Double_t kEtaCut = 0.5;
+   const Double_t kEtaCut = 0.8;
    const Double_t kClusterEnergyCut = 0.0;
    const Double_t kTrackPtCut = 0.0;
    
index d7be236..50c479e 100644 (file)
@@ -45,12 +45,12 @@ Int_t AliAnalysisEtMonteCarlo::AnalyseEvent(AliVEvent* ev)
 
             TParticlePDG *pdgCode =  part->GetPDG(0);
             if (
-                TMath::Abs(pdgCode->PdgCode()) == fPdgDB->GetParticle("proton")->PdgCode() ||
-                TMath::Abs(pdgCode->PdgCode()) == fPdgDB->GetParticle("neutron")->PdgCode() ||
-                TMath::Abs(pdgCode->PdgCode()) == fPdgDB->GetParticle("Lambda0")->PdgCode() ||
-                TMath::Abs(pdgCode->PdgCode()) == fPdgDB->GetParticle("Xi-")->PdgCode() ||
-                TMath::Abs(pdgCode->PdgCode()) == fPdgDB->GetParticle("Xi0")->PdgCode() ||
-                TMath::Abs(pdgCode->PdgCode()) == fPdgDB->GetParticle("Omega-")->PdgCode()
+                TMath::Abs(pdgCode->PdgCode()) == ProtonCode ||
+                TMath::Abs(pdgCode->PdgCode()) == NeutronCode ||
+                TMath::Abs(pdgCode->PdgCode()) == LambdaCode ||
+                TMath::Abs(pdgCode->PdgCode()) == XiCode ||
+                TMath::Abs(pdgCode->PdgCode()) == Xi0Code ||
+                TMath::Abs(pdgCode->PdgCode()) == OmegaCode
               )
             {
                 particleMassPart = -TMath::Sign(pdgCode->PdgCode(), pdgCode->PdgCode())*pdgCode->Mass();
diff --git a/PWG4/totEt/AliAnalysisHadEt.cxx b/PWG4/totEt/AliAnalysisHadEt.cxx
new file mode 100644 (file)
index 0000000..22212c3
--- /dev/null
@@ -0,0 +1,314 @@
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
+//University of Tennessee at Knoxville
+#include "AliAnalysisHadEt.h"
+#include "TMath.h"
+#include "TList.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include <iostream>
+#include "AliAnalysisEtCuts.h"
+#include "AliVEvent.h"
+
+using namespace std;
+
+ClassImp(AliAnalysisHadEt);
+
+
+Int_t AliAnalysisHadEt::numOfEtaBins = 46;
+Float_t AliAnalysisHadEt::etaAxis[47]={-0.78, -0.74, -0.7, -0.66, -0.62, -0.58, -0.54, -0.5, -0.46, -0.42, -0.38, -0.34, -0.3, -0.26, -0.22, -0.18, -0.14, -0.12, -0.1, -0.08, -0.06, -0.04, -0.02, -0.0, 0.02, 0.04, 0.06, 0.08, 0.1, 0.12, 0.14, 0.18, 0.22, 0.26, 0.3, 0.34, 0.38, 0.42, 0.46, 0.5, 0.54, 0.58, 0.62, 0.66, 0.7, 0.74, 0.78};
+Int_t AliAnalysisHadEt::numOfPtBins = 111;
+Float_t AliAnalysisHadEt::ptAxis[117]=
+  {0.0,0.01,0.02,0.03,0.04, 0.05, 0.06,0.07,0.08,0.09, 0.10,0.11, .12,0.13, .14,0.15, .16,0.17, .18,0.19,
+   0.2, .22, .24, .26, .28, 0.30, 0.32, .34, .36, .38, 0.40, .42, .44, .46, .48,
+   0.5, .52, .54, .56, .58, 0.60, 0.62, .64, .66, .68, 0.70, .72, .74, .76, .78,
+   .80, .82, .84, .86, .88, 0.90, 0.92, .94, .96, .98, 1.00,1.05, 1.1,1.15, 1.2,
+  1.25, 1.3,1.35,1.40,1.45, 1.50, 1.55, 1.6,1.65, 1.7, 1.75, 1.8,1.85, 1.9,1.95,
+   2.0, 2.2, 2.4, 2.6, 2.8, 3.00, 3.20, 3.4, 3.6, 3.8, 4.00, 4.2, 4.4, 4.6, 4.8,
+   5.0, 5.5, 6.0, 6.5, 7.0, 7.50, 8.00, 8.5, 9.0, 9.5, 10.0,12.0,14.0,16.0,18.0,
+  20.0,25.0,30.0,35.0,40.0, 45.0, 50.0}; 
+
+AliAnalysisHadEt::AliAnalysisHadEt() :
+        fHistogramNameSuffix("")
+        ,fPdgDB(0)
+        ,PiPlusCode(0)
+        ,PiMinusCode(0)
+        ,KPlusCode(0)
+        ,KMinusCode(0)
+        ,ProtonCode(0)
+        ,AntiProtonCode(0)
+        ,LambdaCode(0)
+        ,AntiLambdaCode(0)
+        ,K0SCode(0)
+        ,OmegaCode(0)
+        ,AntiOmegaCode(0)
+        ,Xi0Code(0)
+        ,AntiXi0Code(0)
+        ,XiCode(0)
+        ,AntiXiCode(0)
+        ,SigmaCode(0)
+        ,AntiSigmaCode(0)
+        ,K0LCode(0)
+        ,NeutronCode(0)
+        ,AntiNeutronCode(0)
+        ,EPlusCode(0)
+        ,EMinusCode(0)
+        ,PionMass(0)
+       ,fSumEt(0)
+       ,fSumEtAcc(0)
+       ,fTotEt(0)
+       ,fTotEtAcc(0)
+       ,fTotNeutralEt(0)
+       ,fTotNeutralEtAcc(0)
+       ,fTotChargedEt(0)
+       ,fTotChargedEtAcc(0)
+       ,fMultiplicity(0)
+       ,fChargedMultiplicity(0)
+       ,fNeutralMultiplicity(0)
+        ,fEtaCut(EtCommonCuts::kEtaCut)
+        ,fEtaCutAcc(0)
+                                 //,fPhiCutAccMin(0)
+                                 //,fPhiCutAccMax(360.)
+        ,fVertexXCut(0)
+        ,fVertexYCut(0)
+        ,fVertexZCut(0)
+                                 ,fIPxyCut(0)
+                                 ,fIPzCut(0)
+                                 //,fSingleCellEnergyCut(0)
+                                 //,fClusterEnergyCut(EtCommonCuts::kClusterEnergyCut)
+                                 //,fTrackPtCut(EtCommonCuts::kTrackPtCut)
+        ,esdtrackCutsITSTPC(0)
+        ,esdtrackCutsTPC(0)
+       ,esdtrackCutsITS(0)
+        ,histoList(0)
+{
+
+}
+
+AliAnalysisHadEt::~AliAnalysisHadEt()
+{
+
+}
+
+Int_t AliAnalysisHadEt::AnalyseEvent(AliVEvent *event)
+{
+  //this line is basically here to eliminate a compiler warning that event is not used.  Making it a virtual function did not work with the plugin.
+  cout<<"This event has "<<event->GetNumberOfTracks()<<" tracks"<<endl;
+  return 0;
+}
+
+void AliAnalysisHadEt::FillOutputList()
+{
+}
+
+void AliAnalysisHadEt::Init()
+{
+  if(!fPdgDB) fPdgDB = new TDatabasePDG();
+  //the codes are defined in $ROOTSYS/etc/pdg_table.txt
+  PionMass = fPdgDB->GetParticle("pi+")->Mass();
+  PiPlusCode = fPdgDB->GetParticle("pi+")->PdgCode();
+    PiMinusCode = fPdgDB->GetParticle("pi-")->PdgCode();
+    KPlusCode = fPdgDB->GetParticle("K+")->PdgCode();
+    KMinusCode = fPdgDB->GetParticle("K-")->PdgCode();
+    ProtonCode = fPdgDB->GetParticle("proton")->PdgCode();
+    AntiProtonCode = fPdgDB->GetParticle("antiproton")->PdgCode();
+    LambdaCode = fPdgDB->GetParticle("Lambda0")->PdgCode();
+    AntiLambdaCode = fPdgDB->GetParticle("Lambda0_bar")->PdgCode();
+    K0SCode = fPdgDB->GetParticle("K_S0")->PdgCode();
+    OmegaCode = fPdgDB->GetParticle("Omega-")->PdgCode();
+    AntiOmegaCode = fPdgDB->GetParticle("Omega+")->PdgCode();
+    Xi0Code = fPdgDB->GetParticle("Xi0")->PdgCode();
+    AntiXi0Code = fPdgDB->GetParticle("Xi0_bar")->PdgCode();
+    XiCode = fPdgDB->GetParticle("Xi-")->PdgCode();
+    AntiXiCode = fPdgDB->GetParticle("Xi-_bar")->PdgCode();
+    SigmaCode = fPdgDB->GetParticle("Sigma-")->PdgCode();
+    AntiSigmaCode = fPdgDB->GetParticle("Sigma+")->PdgCode();
+    K0LCode = fPdgDB->GetParticle("K_L0")->PdgCode();
+    NeutronCode = fPdgDB->GetParticle("neutron")->PdgCode();
+    AntiNeutronCode = fPdgDB->GetParticle("antineutron")->PdgCode();
+    EPlusCode = fPdgDB->GetParticle("e+")->PdgCode();
+    EMinusCode = fPdgDB->GetParticle("e-")->PdgCode();
+}
+
+void AliAnalysisHadEt::CreateHistograms()
+{
+}
+
+void AliAnalysisHadEt::FillHistograms()
+{
+}
+
+void AliAnalysisHadEt::ResetEventValues()
+{
+    fTotEt = 0;
+    fTotEtAcc = 0;
+    fTotNeutralEt = 0;
+    fTotNeutralEtAcc = 0;
+    fTotChargedEt  = 0;
+    fTotChargedEtAcc = 0;
+    fMultiplicity = 0;
+    fChargedMultiplicity = 0;
+    fNeutralMultiplicity = 0;
+}
+void AliAnalysisHadEt::CreateEtaPtHisto2D(TString name, TString title)
+{     
+  TString *histoname   = new TString();
+  TString *histotitle   = new TString();
+
+  histoname->Append(name);
+  histotitle->Append(title);
+  //TH2F *h1 = new TH2F("h1", "Histogram with Gaussian random distribution", numOfPtBins, ptBinsArray, numOfEtaBins, etaBinsArray);
+
+  TH2F *histo = new TH2F(histoname->Data(),histotitle->Data(),numOfPtBins, ptAxis, numOfEtaBins, etaAxis);
+  histo->SetYTitle("#eta");
+  histo->SetXTitle("p_{T}");
+  histo->SetZTitle("E_{T}");
+  histo->Sumw2();
+  histoList->Add(histo);
+  delete histoname;
+  delete histotitle;
+    
+}
+
+void AliAnalysisHadEt::CreateHisto1D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Float_t xlow,Float_t xhigh)
+{     
+  TString *histoname   = new TString();
+  TString *histotitle   = new TString();
+  
+  //cout<<"creating "<<name<<endl;
+
+  histoname->Append(name);
+  histotitle->Append(title);
+  // printf("%s \n ",histoname->Data());
+  TH1F *histo = new TH1F(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh);
+  histo->SetYTitle(ytitle);
+  histo->SetXTitle(xtitle);
+  histo->Sumw2();
+  histoList->Add(histo);
+  delete histoname;
+  delete histotitle;
+    
+}
+void AliAnalysisHadEt::CreateIntHisto1D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Int_t xlow,Int_t xhigh)
+{     
+  TString *histoname   = new TString();
+  TString *histotitle   = new TString();
+  
+  //cout<<"creating "<<name<<endl;
+
+  histoname->Append(name);
+  histotitle->Append(title);
+  // printf("%s \n ",histoname->Data());
+  TH1I *histo = new TH1I(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh);
+  histo->SetYTitle(ytitle);
+  histo->SetXTitle(xtitle);
+  histo->Sumw2();
+  histoList->Add(histo);
+  delete histoname;
+  delete histotitle;
+    
+}
+void AliAnalysisHadEt::CreateHisto2D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Float_t xlow,Float_t xhigh,Int_t ybins,Float_t ylow,Float_t yhigh)
+{     
+  TString *histoname   = new TString();
+  TString *histotitle   = new TString();
+  
+  //cout<<"creating "<<name<<endl;
+
+  histoname->Append(name);
+  histotitle->Append(title);
+  // printf("%s \n ",histoname->Data());
+  TH2F *histo = new TH2F(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh,ybins,ylow,yhigh);
+  histo->SetYTitle(ytitle);
+  histo->SetXTitle(xtitle);
+  histo->Sumw2();
+  histoList->Add(histo);
+  delete histoname;
+  delete histotitle;
+    
+}
+void AliAnalysisHadEt::CreateIntHisto2D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Int_t xlow,Int_t xhigh,Int_t ybins,Int_t ylow,Int_t yhigh)
+{     
+  TString *histoname   = new TString();
+  TString *histotitle   = new TString();
+  
+  //cout<<"creating "<<name<<endl;
+
+  histoname->Append(name);
+  histotitle->Append(title);
+  // printf("%s \n ",histoname->Data());
+  TH2I *histo = new TH2I(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh,ybins,ylow,yhigh);
+  histo->SetYTitle(ytitle);
+  histo->SetXTitle(xtitle);
+  histo->Sumw2();
+  histoList->Add(histo);
+  delete histoname;
+  delete histotitle;
+    
+}
+
+void AliAnalysisHadEt::CreateEtaHisto1D(TString name, TString title)
+{     
+  TString *histoname   = new TString();
+  TString *histotitle   = new TString();
+  
+
+  histoname->Append(name);
+  histotitle->Append(title);
+  TH1F *histo = new TH1F(histoname->Data(),histotitle->Data(),numOfEtaBins, etaAxis);
+  histo->SetYTitle("E_{T}");
+  histo->SetXTitle("#eta");
+  histo->Sumw2();
+  histoList->Add(histo);
+  delete histoname;
+  delete histotitle;
+    
+}
+void AliAnalysisHadEt::FillHisto1D(TString histname, Float_t x, Float_t weight)
+{
+  TH1F     *histo; 
+  TString  *name   = new TString();
+
+  name->Append(histname);       
+  histo = (TH1F *)histoList->FindObject(name->Data()); 
+  if(histo){
+    histo->Fill((Double_t)x, weight);
+  }
+  else{cerr<<"CorrelationMaker::FillHisto1D: no histogram "<< name->Data()<<endl;}
+  delete name;
+}
+void AliAnalysisHadEt::FillHisto2D(TString histname, Float_t x, Float_t y, Float_t weight)
+{
+  TH2F     *histo; 
+  TString  *name   = new TString();
+  
+  name->Append(histname);       
+  histo = (TH2F *)histoList->FindObject(name->Data()); 
+  if(histo){
+    histo->Fill((Double_t)x,(Double_t)y, weight);
+  }
+  else{cerr<<"CorrelationMaker::FillHisto2D: no histogram "<< name->Data()<<endl;}
+  delete name;
+}
+
+
+Float_t AliAnalysisHadEt::Et(TParticle *part, float mass){
+  if(mass == -1000){//if no mass given return default 
+    if(TMath::Abs(part->GetPDG(0)->PdgCode())==2212 || TMath::Abs(part->GetPDG(0)->PdgCode())==2112){
+      if(part->GetPDG(0)->PdgCode()==-2212 || part->GetPDG(0)->PdgCode()==-2112){//antiproton or antineutron
+       //for antinucleons we specifically want to return the kinetic energy plus twice the rest mass
+       return (part->Energy()+part->GetMass())*TMath::Sin(part->Theta());
+      }
+      if(part->GetPDG(0)->PdgCode()==2212 || part->GetPDG(0)->PdgCode()==2112){//antiproton or antineutron
+       //for nucleons we specifically want to return the kinetic energy only
+       return (part->Energy()-part->GetMass())*TMath::Sin(part->Theta());
+      }
+    }
+    else{//otherwise go to the default
+      return part->Energy()*TMath::Sin(part->Theta());
+    }
+  }
+  else{//otherwise use the mass that was given
+    return (TMath::Sqrt(TMath::Power(part->P(),2.0)+TMath::Power(mass,2.0)))*TMath::Sin(part->Theta());
+  }
+  return 0.0;
+}
diff --git a/PWG4/totEt/AliAnalysisHadEt.h b/PWG4/totEt/AliAnalysisHadEt.h
new file mode 100644 (file)
index 0000000..3199462
--- /dev/null
@@ -0,0 +1,188 @@
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
+//University of Tennessee at Knoxville
+#ifndef ALIANALYSISHADET_H
+#define ALIANALYSISHADET_H
+
+#include "TString.h"
+#include "TDatabasePDG.h"
+#include "TParticle.h"
+#include "Rtypes.h"
+#include "TString.h"
+#include "AliESDtrackCuts.h"
+
+class TH2F;
+class TH1F;
+class AliVEvent;
+class TList;
+
+class AliAnalysisHadEt
+{
+public:
+   
+  AliAnalysisHadEt();
+    virtual ~AliAnalysisHadEt();
+
+    /** Analyse the event! */
+    virtual Int_t AnalyseEvent(AliVEvent *event);
+
+    /** Fill the objects you want to output, classes which add new histograms should overload this. */
+    virtual void FillOutputList();
+
+    /** Initialise the analysis, must be overloaded. */
+    virtual void Init();
+
+    /** 
+    * Creates the histograms, must be overloaded if you want to add your own. 
+    * Uses the fHistogramNameSuffix to create proper histogram names
+    */
+    virtual void CreateHistograms();
+    
+    /** Fills the histograms, must be overloaded if you want to add your own */
+    virtual void FillHistograms();
+
+    /** Reset event specific values (Et etc.) */
+    virtual void ResetEventValues();
+    
+    /** Sum of the total Et for all events */
+    Double_t GetSumEt() { return fSumEt; }
+
+    /** Sum of the total Et within our acceptance for all events */
+    Double_t GetSumEtAcc() { return fSumEtAcc; }
+
+    /** Total Et in the event (without acceptance cuts) */
+    Double_t GetTotEt() { return fTotEt; }
+
+    /** Total Et in the event within the acceptance cuts */
+    Double_t GetTotEtAcc() { return fTotEtAcc; }
+
+   /** Total neutral Et in the event (without acceptance cuts) */
+    Double_t GetTotNeutralEt() { return fTotNeutralEt; }
+
+    /** Total neutral Et in the event within the acceptance cuts */
+    Double_t GetTotNeutralEtAcc() { return fTotNeutralEtAcc; }
+    
+    /** Total charged Et in the event (without acceptance cuts) */
+    Double_t GetTotChargedEt() { return fTotChargedEt; }
+
+    /** Total charged Et in the event within the acceptance cuts */
+    Double_t GetTotChargedEtAcc() { return fTotChargedEtAcc; }
+
+
+    void SetHistoList(TList *mylist){histoList = mylist;}
+
+    void SetTPCITSTrackCuts(AliESDtrackCuts *cuts){ esdtrackCutsITSTPC = cuts;}
+    void SetTPCOnlyTrackCuts(AliESDtrackCuts *cuts){ esdtrackCutsTPC = cuts;}
+    void SetITSTrackCuts(AliESDtrackCuts *cuts){ esdtrackCutsITS = cuts;}
+
+protected:
+   
+    /** The suffix for the histogram names */
+    TString fHistogramNameSuffix;
+
+    /** PDG Database */
+    TDatabasePDG *fPdgDB;
+    Int_t PiPlusCode;
+    Int_t PiMinusCode;
+    Int_t KPlusCode;
+    Int_t KMinusCode;
+    Int_t ProtonCode;
+    Int_t AntiProtonCode;
+    Int_t LambdaCode;
+    Int_t AntiLambdaCode;
+    Int_t K0SCode;
+    Int_t OmegaCode;
+    Int_t AntiOmegaCode;
+    Int_t Xi0Code;
+    Int_t AntiXi0Code;
+    Int_t XiCode;
+    Int_t AntiXiCode;
+    Int_t SigmaCode;
+    Int_t AntiSigmaCode;
+    Int_t K0LCode;
+    Int_t NeutronCode;
+    Int_t AntiNeutronCode;
+    Int_t EPlusCode;
+    Int_t EMinusCode;
+    Float_t PionMass;
+
+    /** Sum of the total Et for all events */
+    Double_t fSumEt;
+
+    /** Sum of the total Et within our acceptance for all events */
+    Double_t fSumEtAcc;
+
+    /** Total Et in the event (without acceptance cuts) */
+    Double_t fTotEt;
+
+    /** Total Et in the event within the acceptance cuts */
+    Double_t fTotEtAcc;
+
+    /** Total neutral Et in the event */
+    Double_t fTotNeutralEt;
+
+    /** Total neutral Et in the event within the acceptance cuts */
+    Double_t fTotNeutralEtAcc;
+
+    /** Total charged Et in the event */
+    Double_t fTotChargedEt;
+
+    /** Total charged Et in the event within the acceptance cuts */
+    Double_t fTotChargedEtAcc;
+
+    /** Multiplicity of particles in the event */
+    Int_t fMultiplicity;
+    
+    /** Multiplicity of charged particles in the event */
+    Int_t fChargedMultiplicity;
+    
+    /** Multiplicity of neutral particles in the event */
+    Int_t fNeutralMultiplicity; 
+    
+    /** Cut in eta ( normally |eta| < 0.5 */
+    Double_t fEtaCut;
+
+    /** Eta cut for our acceptance */
+    Double_t fEtaCutAcc;
+
+    /** Vertex cuts */
+    Double_t fVertexXCut;
+    Double_t fVertexYCut;
+    Double_t fVertexZCut;
+
+    /** Impact parameter cuts */
+    Double_t fIPxyCut;
+    Double_t fIPzCut;
+
+
+    void CreateEtaPtHisto2D(TString name, TString title);
+    void CreateEtaHisto1D(TString name, TString title);
+    void CreateHisto2D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Float_t xlow,Float_t xhigh,Int_t ybins,Float_t ylow,Float_t yhigh);
+    void CreateHisto1D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Float_t xlow,Float_t xhigh);
+    void CreateIntHisto1D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Int_t xlow,Int_t xhigh);
+    void CreateIntHisto2D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Int_t xlow,Int_t xhigh,Int_t ybins,Int_t ylow,Int_t yhigh);
+    void FillHisto1D(TString histname, Float_t x, Float_t weight);
+    void FillHisto2D(TString histname, Float_t x, Float_t y, Float_t weight);
+
+    Float_t Et(TParticle *part, float mass = -1000);
+    AliESDtrackCuts* esdtrackCutsITSTPC;
+    AliESDtrackCuts* esdtrackCutsTPC;
+    AliESDtrackCuts* esdtrackCutsITS;
+
+    TList *histoList;
+    static Float_t etaAxis[47];
+    static Int_t numOfEtaBins;
+    static Float_t ptAxis[117];
+    static Int_t numOfPtBins;
+    
+
+ private:
+
+private:
+  //Declare it private to avoid compilation warning
+  AliAnalysisHadEt & operator = (const AliAnalysisHadEt & g) ;//cpy assignment
+  AliAnalysisHadEt(const AliAnalysisHadEt & g) ; // cpy ctor
+
+    ClassDef(AliAnalysisHadEt, 0);
+};
+
+#endif // ALIANALYSISHADET_H
diff --git a/PWG4/totEt/AliAnalysisHadEtMonteCarlo.cxx b/PWG4/totEt/AliAnalysisHadEtMonteCarlo.cxx
new file mode 100644 (file)
index 0000000..cbe871c
--- /dev/null
@@ -0,0 +1,896 @@
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
+//University of Tennessee at Knoxville
+#include "AliAnalysisHadEtMonteCarlo.h"
+#include "AliAnalysisEtCuts.h"
+
+#include "AliStack.h"
+#include "AliMCEvent.h"
+#include "AliESDEvent.h"
+#include "AliESDtrackCuts.h"
+#include "AliESDpid.h"
+
+#include <iostream>
+using namespace std;
+
+
+Int_t AliAnalysisHadEtMonteCarlo::AnalyseEvent(AliVEvent* ev,AliVEvent* ev2){
+  FillHisto1D("NEvents",0.5,1);
+
+  AnalyseEvent(ev);
+  AliMCEvent *mcEvent = dynamic_cast<AliMCEvent*>(ev);
+  AliESDEvent *realEvent = dynamic_cast<AliESDEvent*>(ev2);
+  AliStack *stack = mcEvent->Stack();
+
+  //for PID
+  AliESDpid *PID = new AliESDpid();
+  PID->MakePID(realEvent);
+
+  //This code taken from https://twiki.cern.ch/twiki/bin/view/ALICE/SelectionOfPrimaryTracksForPp2009DataAnalysis
+  //Gets good tracks
+  //=============================================
+  // Primary vertex
+  const AliESDVertex *vertex = realEvent->GetPrimaryVertexTracks();
+  if(vertex->GetNContributors()<1) {
+    // SPD vertex
+    vertex = realEvent->GetPrimaryVertexSPD();
+    if(vertex->GetNContributors()<1) {
+      // NO GOOD VERTEX, SKIP EVENT 
+    }
+  }
+  // apply a cut |zVertex| < CUT, if needed
+
+  
+
+  //esdtrackCutsITSTPC->SetEtaRange(-0.8,0.8); // normally, |eta|<0.8
+  //=============================================
+
+  //Roughly following $ALICE_ROOT/PWG0/dNdEta/AlidNdEtaCorrectionTask
+
+  //=============================================TPC&&ITS=============================================
+  TString *TPC = new TString("TPC");
+  TString *ITS = new TString("ITS");
+  TString *TPCITS = new TString("TPCITS");
+  for(Int_t cutset=0;cutset<3;cutset++){
+    TString *CutName;
+    TObjArray* list;
+    switch(cutset){
+    case 0:
+      CutName = TPC;
+      list = esdtrackCutsTPC->GetAcceptedTracks(realEvent);
+      break;
+    case 1:
+      CutName = ITS;
+      list = esdtrackCutsITS->GetAcceptedTracks(realEvent);
+      break;
+    case 2:
+      CutName = TPCITS;
+      list = esdtrackCutsITSTPC->GetAcceptedTracks(realEvent);
+      break;
+    default:
+      cerr<<"Error:  cannot fill histograms!"<<endl;
+      return -1;
+    }
+    Int_t nGoodTracks = list->GetEntries();
+    for (Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++)
+      {
+       AliESDtrack *track = dynamic_cast<AliESDtrack*> (list->At(iTrack));
+       if (!track)
+         {
+           Printf("ERROR: Could not get track %d", iTrack);
+           continue;
+         }
+       else{
+         bool IsGood = true;
+         if(cutset==1){//if these are ITS stand alone tracks, apply some specific cuts
+           ULong_t trStatus=track->GetStatus();
+           if(trStatus&AliESDtrack::kTPCin) IsGood=false; // reject tracks found in TPC
+           if(trStatus&AliESDtrack::kITSpureSA) IsGood=false; // reject "pure standalone" ITS tracks
+           if(!(trStatus&AliESDtrack::kITSrefit)) IsGood = false; // require proper refit in ITS 
+         }
+         if(!IsGood) continue;
+         Float_t nSigmaPion,nSigmaProton,nSigmaKaon,nSigmaElectron;
+         if(cutset!=1){
+           nSigmaPion = TMath::Abs(PID->NumberOfSigmasTPC(track,AliPID::kPion));
+           nSigmaProton = TMath::Abs(PID->NumberOfSigmasTPC(track,AliPID::kProton));
+           nSigmaKaon = TMath::Abs(PID->NumberOfSigmasTPC(track,AliPID::kKaon));
+           nSigmaElectron = TMath::Abs(PID->NumberOfSigmasTPC(track,AliPID::kElectron));
+         }
+         else{
+           nSigmaPion = TMath::Abs(PID->NumberOfSigmasITS(track,AliPID::kPion));
+           nSigmaProton = TMath::Abs(PID->NumberOfSigmasITS(track,AliPID::kProton));
+           nSigmaKaon = TMath::Abs(PID->NumberOfSigmasITS(track,AliPID::kKaon));
+           nSigmaElectron = TMath::Abs(PID->NumberOfSigmasITS(track,AliPID::kElectron));
+         }
+         bool IsPion = (nSigmaPion<3.0 && nSigmaProton>2.0 && nSigmaKaon>2.0);
+         bool IsElectron = (nSigmaElectron<2.0 && nSigmaPion>4.0 && nSigmaProton>3.0 && nSigmaKaon>3.0);
+         bool IsKaon = (nSigmaPion>3.0 && nSigmaProton>2.0 && nSigmaKaon<2.0);
+         bool IsProton = (nSigmaPion>3.0 && nSigmaProton<2.0 && nSigmaKaon>2.0);
+
+         //bool IsElectron = false;
+         bool Unidentified = (!IsProton && !IsKaon && !IsElectron);
+         Float_t dEdx = track->GetTPCsignal();
+         if(cutset==1) dEdx = track->GetITSsignal();
+
+         FillHisto2D(Form("dEdxAll%s",CutName->Data()),track->P(),dEdx,1.0);
+         //if(cutset==1) cout<<"filling "<<track->P()<<" "<<dEdx<<endl;
+
+         UInt_t label = (UInt_t)TMath::Abs(track->GetLabel());
+         TParticle  *simPart  = stack->Particle(label);
+         if(!simPart) {
+           Printf("no MC particle\n");                 
+           continue;           
+         }
+         else{//analysis
+           if(stack->IsPhysicalPrimary(label)){
+             if (TMath::Abs(simPart->Eta()) < fEtaCut)     {
+
+               Int_t pdgCode =  simPart->GetPDG(0)->PdgCode();
+               Int_t mypid = 0;
+               if(pdgCode==PiPlusCode) mypid = 1;
+               if(pdgCode==ProtonCode) mypid = 2;
+               if(pdgCode==KPlusCode) mypid = 3;
+               if(pdgCode==EPlusCode) mypid = 4;
+               if(pdgCode==PiMinusCode) mypid = 1;
+               if(pdgCode==AntiProtonCode) mypid = 2;
+               if(pdgCode==KMinusCode) mypid = 3;
+               if(pdgCode==EMinusCode) mypid = 4;
+               //cout<<pdgCode->PdgCode()<<" ";
+               //fPdgDB->GetSimParticle("pi+")->PdgCode();
+               bool filled = false;      
+               //============Charged hadrons===================================
+               //identified...
+               if(IsPion){
+                 if(pdgCode!=PiPlusCode && pdgCode!=PiMinusCode){
+                   FillHisto2D(Form("MisidentifiedPIDs%s",CutName->Data()),1,mypid,1);
+                   //if(mypid==0)cerr<<"I was misidentified! I'm not a pion! I am a "<<simPart->GetName()<<endl;
+                 }
+                 float myEt = Et(simPart);
+                 if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedPiPlus",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 else{ FillHisto2D(Form("EtReconstructed%sIdentifiedPiMinus",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 FillHisto2D(Form("dEdxPion%s",CutName->Data()),track->P(),dEdx,1.0);
+               }
+               if(IsProton){
+                 if(pdgCode!=ProtonCode && pdgCode!=AntiProtonCode){
+                   FillHisto2D(Form("MisidentifiedPIDs%s",CutName->Data()),2,mypid,1);
+                   // if(mypid==0)cerr<<"I was misidentified!  I'm not a proton! I am a "<<simPart->GetName()<<endl;
+                 }
+                 float myEt = Et(simPart);
+                 if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedProton",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 else{ FillHisto2D(Form("EtReconstructed%sIdentifiedAntiProton",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 FillHisto2D(Form("dEdxProton%s",CutName->Data()),track->P(),dEdx,1.0);
+               }
+               if(IsKaon){
+                 if(pdgCode!=KMinusCode && pdgCode!=KPlusCode){
+                   FillHisto2D(Form("MisidentifiedPIDs%s",CutName->Data()),3,mypid,1);
+                   //if(mypid==0)cerr<<"I was misidentified!  I'm not a kaon! I am a "<<simPart->GetName()<<" p "<<track->P()<<" nSigmaProton "<<nSigmaProton<<" nSigmaPion "<<nSigmaPion<<" nSigmaKaon "<<nSigmaKaon<<endl;
+                 }
+                 float myEt = Et(simPart);
+                 if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedKPlus",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 else{ FillHisto2D(Form("EtReconstructed%sIdentifiedKMinus",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 FillHisto2D(Form("dEdxKaon%s",CutName->Data()),track->P(),dEdx,1.0);
+               }
+               if(IsElectron){
+                 if(pdgCode!=EMinusCode && pdgCode!=EPlusCode){
+                   FillHisto2D(Form("MisidentifiedPIDs%s",CutName->Data()),4,mypid,1);
+                   //cerr<<"I was misidentified!  I'm not an electron! I am a "<<simPart->GetName()<<endl;
+                 }
+                 float myEt = Et(simPart);
+                 if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedEPlus",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 else{ FillHisto2D(Form("EtReconstructed%sIdentifiedEMinus",CutName->Data()),track->Pt(),track->Eta(),myEt);}
+                 FillHisto2D(Form("dEdxElectron%s",CutName->Data()),track->P(),dEdx,1.0);
+               }
+               if(Unidentified){
+                 if(pdgCode!=EMinusCode && pdgCode!=EPlusCode){
+                   float myEtPi = Et(simPart,PionMass);
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sUnidentifiedAssumingPion",CutName->Data()),track->Pt(),track->Eta(),myEtPi);
+                   FillHisto2D(Form("EtReconstructed%sUnidentified",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                   FillHisto2D(Form("EtNReconstructed%sUnidentified",CutName->Data()),track->Pt(),track->Eta(),1.0);
+                 }
+                 FillHisto2D(Form("dEdxUnidentified%s",CutName->Data()),track->P(),dEdx,1.0);
+                 //cout<<"I was not identified.  I am a "<<simPart->GetName()<<" PID "<<pdgCode<<endl;
+                 //track what was not identified successfully
+                 FillHisto1D(Form("UnidentifiedPIDs%s",CutName->Data()),mypid,1);
+               }
+               //...simulated
+               if(pdgCode == PiPlusCode){
+                 //cout<<"I'm a real primary "<<simPart->GetName()<<"! "<<"my label is "<<simPart->GetFirstMother()<<" track no "<<iTrack<<"/"<<realEvent->GetNumberOfTracks()<<endl;//<<" "<<label<<" "<<pdgCode<<endl;
+               
+                 float myEt = Et(simPart);
+                 FillHisto2D(Form("EtReconstructed%sPiPlus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sPiPlus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 filled = true;
+               }
+               if(pdgCode == PiMinusCode){
+                 float myEt = Et(simPart);
+                 FillHisto2D(Form("EtReconstructed%sPiMinus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sPiMinus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 filled = true;
+               }
+               if(pdgCode == KPlusCode){
+                 float myEt = Et(simPart);
+                 float myEtPi = Et(simPart,PionMass);
+                 FillHisto2D(Form("EtReconstructed%sKPlus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sKPlus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 FillHisto2D(Form("EtReconstructed%sKPlusAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 filled = true;
+               }
+               if(pdgCode == KMinusCode){
+                 float myEt = Et(simPart);
+                 float myEtPi = Et(simPart,PionMass);
+                 FillHisto2D(Form("EtReconstructed%sKMinus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sKMinus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 FillHisto2D(Form("EtReconstructed%sKMinusAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 filled = true;
+               }
+               if(pdgCode == ProtonCode){
+                 float myEt = Et(simPart);
+                 float myEtPi = Et(simPart,PionMass);
+                 FillHisto2D(Form("EtReconstructed%sProton",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sProton",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 FillHisto2D(Form("EtReconstructed%sProtonAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 filled = true;
+               }
+               if(pdgCode == AntiProtonCode){
+                 float myEt = Et(simPart);
+                 float myEtPi = Et(simPart,PionMass);
+                 FillHisto2D(Form("EtReconstructed%sAntiProton",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sAntiProton",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtNReconstructed%sChargedHadron",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 FillHisto2D(Form("EtReconstructed%sAntiProtonAssumingPion",CutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
+                 filled = true;
+               }
+               if(pdgCode == EPlusCode){
+                 float myEt = Et(simPart);
+                 FillHisto2D(Form("EtReconstructed%sEPlus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                 if(!IsElectron || Unidentified){
+                   float myEtPi = Et(simPart,PionMass);
+                   FillHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",CutName->Data()),track->Pt(),track->Eta(),myEtPi);
+                 }
+                 filled = true;
+               }
+               if(pdgCode == EMinusCode){
+                 if(!IsElectron || Unidentified){
+                   float myEtPi = Et(simPart,PionMass);
+                   FillHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",CutName->Data()),track->Pt(),track->Eta(),myEtPi);
+                 }
+                 float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sEMinus",CutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
+                   filled = true;
+               }
+               //if(!filled){
+                 //TParticlePDG *pc = simPart->GetPDG(0);
+                 //if( strcmp(pc->ParticleClass(),"Baryon")==0 || strcmp(pc->ParticleClass(),"Meson")==0 ){
+                 //cout<<"Did not find a place for "<<simPart->GetName()<<" "<<pdgCode<<" which is a "<<pc->ParticleClass()<<endl;
+                 //}
+                 //}
+             }
+             
+           }
+           else{//not a primary - we're after V0 daughters!
+             //cout<<"I'm a secondary "<<simPart->GetName()<<"!";//<<endl;
+             TParticle *mom = stack->Particle(simPart->GetFirstMother());
+             if(mom){
+               TParticlePDG *pc = mom->GetPDG(0);
+               if(pc){
+                 Int_t pdgCode =  mom->GetPDG(0)->PdgCode();
+                 if(pdgCode == LambdaCode){
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sLambdaDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                 }
+                 if(pdgCode == AntiLambdaCode){
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sAntiLambdaDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                 }
+                 if(pdgCode == K0SCode){
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sK0SDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                 }
+                 if(pdgCode == XiCode){
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sXiDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                 }
+                 if(pdgCode == AntiXiCode){
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sAntiXiDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                 }
+                 if(pdgCode == OmegaCode){
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sOmegaDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                 }
+                 if(pdgCode == XiCode){
+                   float myEt = Et(simPart);
+                   FillHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                 }
+
+                 if(mom->GetFirstMother()>0){
+                   TParticle *grandma = stack->Particle(mom->GetFirstMother());
+                   if(grandma){
+                     Int_t pdgCodeMom =  mom->GetPDG(0)->PdgCode();
+                     if(pdgCodeMom==PiPlusCode || pdgCodeMom==PiMinusCode || pdgCodeMom==ProtonCode ||pdgCodeMom==AntiProtonCode || pdgCodeMom==KPlusCode || pdgCode==KMinusCode){
+                       //cout<<" my grandmother is "<<grandma->GetName()<<" "<<endl;
+                       Int_t pdgCodeGrandma =  grandma->GetPDG(0)->PdgCode();
+                     
+                       if(pdgCodeGrandma == XiCode){
+                         float myEt = Et(simPart);
+                         FillHisto2D(Form("EtReconstructed%sXiDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                       }
+                       if(pdgCodeGrandma == AntiXiCode){
+                         float myEt = Et(simPart);
+                         FillHisto2D(Form("EtReconstructed%sAntiXiDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                       }
+                       if(pdgCodeGrandma == OmegaCode){
+                         float myEt = Et(simPart);
+                         FillHisto2D(Form("EtReconstructed%sOmegaDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                       }
+                       if(pdgCodeGrandma == XiCode){
+                         float myEt = Et(simPart);
+                         FillHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",CutName->Data()),track->Pt(),track->Eta(),myEt);
+                       }
+
+                     }
+                   }
+                 }
+               }
+             }
+           }
+         }
+
+       }
+      }
+    delete list;
+  }
+  //delete AliESDpid;
+  return 1;
+}
+Int_t AliAnalysisHadEtMonteCarlo::AnalyseEvent(AliVEvent* ev)
+{
+     ResetEventValues();
+     
+    // Get us an mc event
+    AliMCEvent *mcEvent = dynamic_cast<AliMCEvent*>(ev);
+
+    // Let's play with the stack!
+    AliStack *stack = mcEvent->Stack();
+
+    Int_t nPrim = stack->GetNtrack();
+
+    //=================Tracks which may or may not have been reconstructed=================
+
+    for (Int_t iPart = 0; iPart < nPrim; iPart++)
+    {
+
+        TParticle *part = stack->Particle(iPart);
+
+        if (!part)
+         {
+            Printf("ERROR: Could not get particle %d", iPart);
+            continue;
+         }
+
+        TParticlePDG *pc = part->GetPDG(0);
+
+        // Check if it is a primary particle
+        if (!stack->IsPhysicalPrimary(iPart)){//secondaries...
+       }
+       else{//primaries
+         // Check for reasonable (for now neutral and singly charged) charge on the particle
+         //note that the charge is stored in units of e/3
+         //if (TMath::Abs(pc->Charge()) != EtMonteCarloCuts::kSingleChargedParticle && pc->Charge() != EtMonteCarloCuts::kNeutralParticle) continue;
+
+
+         if (TMath::Abs(part->Eta()) < fEtaCut)            {
+
+           Int_t pdgCode =  part->GetPDG(0)->PdgCode();
+           //cout<<pdgCode->PdgCode()<<" ";
+           //fPdgDB->GetParticle("pi+")->PdgCode();
+           bool filled = false;
+           //============Charged hadrons===================================
+           if(pdgCode == PiPlusCode){
+             //cout<<"I'm a simulated primary "<<part->GetName()<<"! "<<"my label is "<<part->GetFirstMother()<<" pt "<<part->Pt()<<endl;
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedPiPlus",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedPiPlus",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == PiMinusCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedPiMinus",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedPiMinus",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == KPlusCode){
+             float myEt = Et(part);
+             float myEtPi = Et(part,PionMass);
+             FillHisto2D("EtSimulatedKPlus",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedKPlus",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedKPlusAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == KMinusCode){
+             float myEt = Et(part);
+             float myEtPi = Et(part,PionMass);
+             FillHisto2D("EtSimulatedKMinus",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedKMinus",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedKMinusAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == ProtonCode){
+             float myEt = Et(part);
+             float myEtPi = Et(part,PionMass);
+             FillHisto2D("EtSimulatedProton",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedProton",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedProtonAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == AntiProtonCode){
+             float myEt = Et(part);
+             float myEtPi = Et(part,PionMass);
+             FillHisto2D("EtSimulatedAntiProton",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedAntiProton",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
+             FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedAntiProtonAssumingPion",part->Pt(),part->Eta(),myEtPi);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           //============Other hadrons===================================
+
+           if(pdgCode == NeutronCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedNeutron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == AntiNeutronCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedAntiNeutron",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == LambdaCode){
+             float myEt = Et(part);
+             //cout<<"I am a simulated lambda! pt "<<part->Pt()<<" eta "<<part->Eta()<<endl;
+             FillHisto2D("EtSimulatedLambda",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             Int_t ndaughters = part->GetNDaughters();
+             for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
+               Int_t daughterindex = part->GetDaughter(idaughter);
+               if(daughterindex<0 || daughterindex>1e5) continue;
+               TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
+               if(daughter){
+                 if(daughter->GetPDG(0)){
+                   Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
+                   if(daughtercode==PiMinusCode || daughtercode==ProtonCode){
+                     myEt = Et(daughter);
+                     FillHisto2D("EtSimulatedLambdaDaughters",daughter->Pt(),daughter->Eta(),myEt);
+                     //cout<<"Lambda daughter is a "<<daughter->GetName()<<endl;
+                   }
+                 }
+                 else{
+                   //cout<<"Lambda daughter is a "<<daughter->GetName()<<endl;
+                 }
+               }
+             }
+             filled = true;
+           }
+           if(pdgCode == AntiLambdaCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedAntiLambda",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             Int_t ndaughters = part->GetNDaughters();
+             for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
+               Int_t daughterindex = part->GetDaughter(idaughter);
+               if(daughterindex<0 || daughterindex>1e5) continue;
+               TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
+               if(daughter){
+                 if(daughter->GetPDG(0)){
+                   Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
+                   if(daughtercode==PiPlusCode || daughtercode==AntiProtonCode){
+                     myEt = Et(daughter);
+                     FillHisto2D("EtSimulatedAntiLambdaDaughters",daughter->Pt(),daughter->Eta(),myEt);
+                     //cout<<"AntiLambda daughter is a "<<daughter->GetName()<<endl;
+                   }
+                 }
+                 else{
+                   //cout<<"AntiLambda daughter is a "<<daughter->GetName()<<endl;
+                 }
+               }
+             }
+             filled = true;
+           }
+           if(pdgCode == K0SCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedK0S",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             Int_t ndaughters = part->GetNDaughters();
+             for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
+               Int_t daughterindex = part->GetDaughter(idaughter);
+               if(daughterindex<0 || daughterindex>1e5) continue;
+               TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
+               if(daughter){
+                 if(daughter->GetPDG(0)){
+
+                   Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
+                   if(daughtercode==PiMinusCode || daughtercode==PiPlusCode){
+                     myEt = Et(daughter);
+                     FillHisto2D("EtSimulatedK0SDaughters",daughter->Pt(),daughter->Eta(),myEt);
+                     //cout<<"K0S daughter is a "<<daughter->GetName()<<endl;
+                   }
+                 }
+                 else{
+                   //cout<<"K0S daughter is a "<<daughter->GetName()<<endl;
+                 }
+               }
+             }
+             filled = true;
+           }
+           if(pdgCode == K0LCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedK0L",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == OmegaCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedOmega",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             Int_t ndaughters = part->GetNDaughters();
+             for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
+               Int_t daughterindex = part->GetDaughter(idaughter);
+               if(daughterindex<0 || daughterindex>1e5) continue;
+               TParticle *daughter = stack->Particle(daughterindex);
+               if(daughter){
+                 if(daughter->GetPDG(0)){
+
+                   Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
+                   if(daughtercode==PiPlusCode || daughtercode==ProtonCode || daughtercode==KMinusCode){
+                     myEt = Et(daughter);
+                     FillHisto2D("EtSimulatedOmegaDaughters",daughter->Pt(),daughter->Eta(),myEt);
+                   //cout<<"Omega daughter is a "<<daughter->GetName()<<endl;
+                   }
+                 }
+                 else{
+                   //cout<<"Omega daughter is a "<<daughter->GetName()<<endl;
+                 }
+               }
+             }
+             filled = true;
+           }
+           if(pdgCode == AntiOmegaCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedOmega",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             Int_t ndaughters = part->GetNDaughters();
+             for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
+               Int_t daughterindex = part->GetDaughter(idaughter);
+               if(daughterindex<0 || daughterindex>1e5) continue;
+               TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
+               if(daughter){
+                 if(daughter->GetPDG(0)){
+                   Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
+                   if(daughtercode==PiMinusCode || daughtercode==AntiProtonCode || daughtercode==KPlusCode){
+                     myEt = Et(daughter);
+                     FillHisto2D("EtSimulatedAntiOmegaDaughters",daughter->Pt(),daughter->Eta(),myEt);
+                     //cout<<"AntiOmega daughter is a "<<daughter->GetName()<<endl;
+                   }
+                 }
+                 else{
+                   //cout<<"AntiOmega daughter is a "<<daughter->GetName()<<endl;
+                 }
+               }
+             }
+             filled = true;
+           }
+           //There are two codes for Sigmas
+           if(pdgCode == SigmaCode || pdgCode == -3222){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedSigma",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == AntiSigmaCode || pdgCode == 3222){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedAntiSigma",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == XiCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedXi",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             Int_t ndaughters = part->GetNDaughters();
+             for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
+               Int_t daughterindex = part->GetDaughter(idaughter);
+               if(daughterindex<0 || daughterindex>1e5 || daughterindex>1e5) continue;
+               //cerr<<"Daughter index "<<daughterindex<<" npart "<<nPrim<<endl;
+               TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
+               if(daughter){
+                 if(daughter->GetPDG(0)){
+
+                   Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
+                   if(daughtercode==PiPlusCode || daughtercode==ProtonCode || daughtercode==PiMinusCode){
+                     myEt = Et(daughter);
+                     FillHisto2D("EtSimulatedXiDaughters",daughter->Pt(),daughter->Eta(),myEt);
+                   //cout<<"Xi daughter is a "<<daughter->GetName()<<endl;
+                   }
+                 }
+                 else{
+                   //cout<<"Xi daughter is a "<<daughter->GetName()<<endl;
+                 }
+               }
+             }
+             filled = true;
+           }
+           if(pdgCode == AntiXiCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedAntiXi",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             Int_t ndaughters = part->GetNDaughters();
+             for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
+               Int_t daughterindex = part->GetDaughter(idaughter);
+               if(daughterindex<0 || daughterindex>1e5) continue;
+               TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
+               if(daughter){
+                 if(daughter->GetPDG(0)){
+                   Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
+                   if(daughtercode==PiPlusCode || daughtercode==AntiProtonCode || daughtercode==PiMinusCode){
+                     myEt = Et(daughter);
+                     FillHisto2D("EtSimulatedAntiXiDaughters",daughter->Pt(),daughter->Eta(),myEt);
+                     //cout<<"AntiXi daughter is a "<<daughter->GetName()<<endl;
+                   }
+                 }
+                 else{
+                   //cout<<"AntiXi daughter is a "<<daughter->GetName()<<endl;
+                 }
+               }
+             }
+             filled = true;
+           }
+           if(pdgCode == Xi0Code){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedXi0",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == AntiXi0Code){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedAntiXi0",part->Pt(),part->Eta(),myEt);
+             FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           //============electrons===================================
+
+           if(pdgCode == EPlusCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedEPlus",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(pdgCode == EMinusCode){
+             float myEt = Et(part);
+             FillHisto2D("EtSimulatedEMinus",part->Pt(),part->Eta(),myEt);
+             filled = true;
+           }
+           if(!filled){
+             if( strcmp(pc->ParticleClass(),"Baryon")==0 || strcmp(pc->ParticleClass(),"Meson")==0 ){
+               //cout<<"Did not find a place for "<<part->GetName()<<" "<<pdgCode<<" which is a "<<pc->ParticleClass()<<endl;
+             }
+           }
+         }
+       }
+    }
+
+
+
+
+//     fTotNeutralEtAcc = fTotNeutralEt;
+//     fTotEt = fTotChargedEt + fTotNeutralEt;
+//     fTotEtAcc = fTotChargedEtAcc + fTotNeutralEtAcc;
+    
+//     FillHistograms();
+
+    return 1;
+    
+}
+
+void AliAnalysisHadEtMonteCarlo::Init()
+{
+
+    AliAnalysisHadEt::Init();
+
+    fVertexXCut = EtReconstructedCuts::kVertexXCut;
+    fVertexYCut = EtReconstructedCuts::kVertexYCut;
+    fVertexZCut = EtReconstructedCuts::kVertexZCut;
+    fIPxyCut = EtReconstructedCuts::kIPxyCut;
+    fIPzCut = EtReconstructedCuts::kIPzCut;
+    // Track cuts
+    //Bool_t selectPrimaries=kTRUE;
+    //esdtrackCutsITSTPC = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selectPrimaries);
+    //esdtrackCutsITSTPC = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+
+}
+void AliAnalysisHadEtMonteCarlo::CreateHistograms(){
+  //for simulated Et only (no reconstruction)
+  CreateEtaPtHisto2D(TString("EtSimulatedPiPlus"),TString("Simulated E_{T} from #pi^{+}"));
+  CreateEtaPtHisto2D("EtSimulatedPiMinus","Simulated E_{T} from #pi^{-}");
+  CreateEtaPtHisto2D("EtSimulatedKPlus","Simulated E_{T} from K^{+}");
+  CreateEtaPtHisto2D("EtSimulatedKMinus","Simulated E_{T} from K^{-}");
+  CreateEtaPtHisto2D("EtSimulatedProton","Simulated E_{T} from p");
+  CreateEtaPtHisto2D("EtSimulatedAntiProton","Simulated E_{T} from #bar{p}");
+  CreateEtaPtHisto2D("EtSimulatedChargedHadron","Simulated E_{T} from charged hadrons");
+  CreateEtaPtHisto2D(TString("EtNSimulatedPiPlus"),TString("Number of Simulated #pi^{+}"));
+  CreateEtaPtHisto2D("EtNSimulatedPiMinus","Number of simulated #pi^{-}");
+  CreateEtaPtHisto2D("EtNSimulatedKPlus","Number of simulated K^{+}");
+  CreateEtaPtHisto2D("EtNSimulatedKMinus","Number of simulated K^{-}");
+  CreateEtaPtHisto2D("EtNSimulatedProton","Number of simulated p");
+  CreateEtaPtHisto2D("EtNSimulatedAntiProton","Number of simulated #bar{p}");
+  CreateEtaPtHisto2D("EtNSimulatedChargedHadron","Number of simulated charged hadrons");
+
+  CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingPion","Simulated E_{T} from charged hadrons assuming they are all pions");
+  CreateEtaPtHisto2D("EtSimulatedKPlusAssumingPion","Simulated E_{T} from K^{+} assuming #pi mass");
+  CreateEtaPtHisto2D("EtSimulatedKMinusAssumingPion","Simulated E_{T} from K^{-} assuming #pi mass");
+  CreateEtaPtHisto2D("EtSimulatedProtonAssumingPion","Simulated E_{T} from p assuming #pi mass");
+  CreateEtaPtHisto2D("EtSimulatedAntiProtonAssumingPion","Simulated E_{T} from #bar{p} assuming #pi mass");
+
+  CreateEtaPtHisto2D("EtSimulatedLambda","Simulated E_{T} from #Lambda");
+  CreateEtaPtHisto2D("EtSimulatedAntiLambda","Simulated E_{T} from #bar{#Lambda}");
+  CreateEtaPtHisto2D("EtSimulatedK0S","Simulated E_{T} from K^{0}_{S}");
+  CreateEtaPtHisto2D("EtSimulatedK0L","Simulated E_{T} from K^{0}_{L}");
+  CreateEtaPtHisto2D("EtSimulatedNeutron","Simulated E_{T} from neutrons");
+  CreateEtaPtHisto2D("EtSimulatedAntiNeutron","Simulated E_{T} from #bar{n}");
+  CreateEtaPtHisto2D("EtSimulatedEPlus","Simulated E_{T} from e^{+}");
+  CreateEtaPtHisto2D("EtSimulatedEMinus","Simulated E_{T} from e^{-}");
+  CreateEtaPtHisto2D("EtSimulatedOmega","Simulated E_{T} from #Omega^{-}");
+  CreateEtaPtHisto2D("EtSimulatedAntiOmega","Simulated E_{T} from #Omega^{+}");
+  CreateEtaPtHisto2D("EtSimulatedXi","Simulated E_{T} from #Xi^{-}");
+  CreateEtaPtHisto2D("EtSimulatedAntiXi","Simulated E_{T} from #Xi^{+}");
+  CreateEtaPtHisto2D("EtSimulatedSigma","Simulated E_{T} from #Xi^{-}");
+  CreateEtaPtHisto2D("EtSimulatedAntiSigma","Simulated E_{T} from #Xi^{+}");
+  CreateEtaPtHisto2D("EtSimulatedXi0","Simulated E_{T} from #Xi^{0}");
+  CreateEtaPtHisto2D("EtSimulatedAntiXi0","Simulated E_{T} from #Xi^{0}");
+  CreateEtaPtHisto2D("EtSimulatedAllHadron","Simulated E_{T} from all hadrons");
+
+
+  CreateEtaPtHisto2D("EtSimulatedLambdaDaughters","Simulated E_{T} from #Lambda Daughters");
+  CreateEtaPtHisto2D("EtSimulatedAntiLambdaDaughters","Simulated E_{T} from #bar{#Lambda} Daughters");
+  CreateEtaPtHisto2D("EtSimulatedK0SDaughters","Simulated E_{T} from K^{0}_{S} Daughters");
+  CreateEtaPtHisto2D("EtSimulatedOmegaDaughters","Simulated E_{T} from #Omega^{-} Daughters");
+  CreateEtaPtHisto2D("EtSimulatedAntiOmegaDaughters","Simulated E_{T} from #Omega^{+} Daughters");
+  CreateEtaPtHisto2D("EtSimulatedXiDaughters","Simulated E_{T} from #Xi^{-} Daughters");
+  CreateEtaPtHisto2D("EtSimulatedAntiXiDaughters","Simulated E_{T} from #Xi^{+} Daughters");
+
+  TString *TPC = new TString("TPC");
+  TString *ITS = new TString("ITS");
+  TString *TPCITS = new TString("TPCITS");
+  for(Int_t i=0;i<3;i++){
+    TString *CutName;
+    Float_t maxPtdEdx = 10;
+    Float_t mindEdx = 35;
+    Float_t maxdEdx = 150.0;
+    switch(i){
+    case 0:
+      CutName = TPC;
+      break;
+    case 1:
+      CutName = ITS;
+      maxPtdEdx = 5;
+      maxdEdx = 500.0;
+      break;
+    case 2:
+      CutName = TPCITS;
+      break;
+    default:
+      cerr<<"Error:  cannot make histograms!"<<endl;
+      return;
+    }
+
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedPiPlus",CutName->Data()),"Reconstructed E_{T} from identified #pi^{+}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedPiMinus",CutName->Data()),"Reconstructed E_{T} from identified #pi^{-}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedKPlus",CutName->Data()),"Reconstructed E_{T} from identified K^{+}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedEMinus",CutName->Data()),"Reconstructed E_{T} from identified e^{-}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedEPlus",CutName->Data()),"Reconstructed E_{T} from identified e^{+}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedKMinus",CutName->Data()),"Reconstructed E_{T} from identified K^{-}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedProton",CutName->Data()),"Reconstructed E_{T} from identified p");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedAntiProton",CutName->Data()),"Reconstructed E_{T} from identified #bar{p}");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sUnidentified",CutName->Data()),"Number of Reconstructed unidentified particles");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedAssumingPion",CutName->Data()),"Reconstructed E_{T} from unidentified particles assuming pion mass");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentified",CutName->Data()),"Reconstructed E_{T} from unidentified particles using real mass");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",CutName->Data()),"Reconstructed E_{T} from misidentified electrons");
+
+
+    CreateEtaPtHisto2D(Form("EtReconstructed%sPiPlus",CutName->Data()),"Reconstructed E_{T} from #pi^{+}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sPiMinus",CutName->Data()),"Reconstructed E_{T} from #pi^{-}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sKPlus",CutName->Data()),"Reconstructed E_{T} from K^{+}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sKMinus",CutName->Data()),"Reconstructed E_{T} from K^{-}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sProton",CutName->Data()),"Reconstructed E_{T} from p");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProton",CutName->Data()),"Reconstructed E_{T} from #bar{p}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadron",CutName->Data()),"Reconstructed E_{T} from charged hadrons");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sPiPlus",CutName->Data()),"Reconstructed E_{T} from #pi^{+}");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sPiMinus",CutName->Data()),"Reconstructed E_{T} from #pi^{-}");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sKPlus",CutName->Data()),"Reconstructed E_{T} from K^{+}");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sKMinus",CutName->Data()),"Reconstructed E_{T} from K^{-}");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sProton",CutName->Data()),"Reconstructed E_{T} from p");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sAntiProton",CutName->Data()),"Reconstructed E_{T} from #bar{p}");
+    CreateEtaPtHisto2D(Form("EtNReconstructed%sChargedHadron",CutName->Data()),"Reconstructed E_{T} from charged hadrons");
+
+    CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",CutName->Data()),"Reconstructed E_{T} from charged hadrons assuming they are all pions");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sKPlusAssumingPion",CutName->Data()),"Reconstructed E_{T} from K^{+} assuming #pi mass");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sKMinusAssumingPion",CutName->Data()),"Reconstructed E_{T} from K^{-} assuming #pi mass");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sProtonAssumingPion",CutName->Data()),"Reconstructed E_{T} from p assuming #pi mass");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProtonAssumingPion",CutName->Data()),"Reconstructed E_{T} from #bar{p} assuming #pi mass");
+
+    CreateEtaPtHisto2D(Form("EtReconstructed%sEPlus",CutName->Data()),"Reconstructed E_{T} from e^{+}");
+    CreateEtaPtHisto2D(Form("EtReconstructed%sEMinus",CutName->Data()),"Reconstructed E_{T} from e^{-}");
+
+
+
+    CreateEtaPtHisto2D(Form("EtReconstructed%sLambdaDaughters",CutName->Data()),"Reconstructed E_{T} from #Lambda Daughters");
+  CreateEtaPtHisto2D(Form("EtReconstructed%sAntiLambdaDaughters",CutName->Data()),"Reconstructed E_{T} from #bar{#Lambda} Daughters");
+  CreateEtaPtHisto2D(Form("EtReconstructed%sK0SDaughters",CutName->Data()),"Reconstructed E_{T} from K^{0}_{S} Daughters");
+  CreateEtaPtHisto2D(Form("EtReconstructed%sOmegaDaughters",CutName->Data()),"Reconstructed E_{T} from #Omega^{-} Daughters");
+  CreateEtaPtHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",CutName->Data()),"Reconstructed E_{T} from #Omega^{+} Daughters");
+  CreateEtaPtHisto2D(Form("EtReconstructed%sXiDaughters",CutName->Data()),"Reconstructed E_{T} from #Xi^{-} Daughters");
+  CreateEtaPtHisto2D(Form("EtReconstructed%sAntiXiDaughters",CutName->Data()),"Reconstructed E_{T} from #Xi^{+} Daughters");
+
+    CreateIntHisto1D(Form("UnidentifiedPIDs%s",CutName->Data()),"PIDs of unidentified particles", "PID", "Number of particles",9, -4,4);
+    CreateHisto2D(Form("MisidentifiedPIDs%s",CutName->Data()),"PIDs of misidentified particles", "PID real","PID identified",5, -.5,4.5,5, -.5,4.5);
+    CreateHisto2D(Form("dEdxAll%s",CutName->Data()),"dE/dx for all particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxPion%s",CutName->Data()),"dE/dx for #pi^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxKaon%s",CutName->Data()),"dE/dx for K^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxProton%s",CutName->Data()),"dE/dx for p(#bar{p})","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxElectron%s",CutName->Data()),"dE/dx for e^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxUnidentified%s",CutName->Data()),"dE/dx for unidentified particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+  }
+
+  CreateIntHisto1D("NEvents","Number of events","number of events","Number of events",1,0,1);
+
+  //CreateHisto1D("MisidentifiedPIDs","PIDs for particles misidentified that are not a #pi,K,p","PID","number of entries",3000,0.5,3000.5);
+
+
+
+//     list->Add(fHistEt);
+//     TString histname = "fHistEt" + fHistogramNameSuffix;
+
+//     fHistEt = new TH1F(histname.Data(), "Total E_{T} Distribution", 1000, 0.00, 99);
+//     fHistEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
+//     fHistEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
+}
+
diff --git a/PWG4/totEt/AliAnalysisHadEtMonteCarlo.h b/PWG4/totEt/AliAnalysisHadEtMonteCarlo.h
new file mode 100644 (file)
index 0000000..4043568
--- /dev/null
@@ -0,0 +1,28 @@
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
+//University of Tennessee at Knoxville
+#ifndef ALIANALYSISHADETMONTECARLO_H
+#define ALIANALYSISHADETMONTECARLO_H
+
+#include "AliAnalysisHadEt.h"
+#include "AliESDtrackCuts.h"
+#include <iostream>
+
+class AliAnalysisHadEtMonteCarlo : public AliAnalysisHadEt
+{
+
+public:
+   
+    AliAnalysisHadEtMonteCarlo() {}
+   
+    virtual Int_t AnalyseEvent(AliVEvent* event);
+    virtual Int_t AnalyseEvent(AliVEvent* event,AliVEvent* event2);
+
+    //void FillHistograms();
+    void CreateHistograms();
+    virtual void Init();
+    
+ private:
+    //AliESDtrackCuts* esdtrackCutsITSTPC;
+};
+
+#endif // ALIANALYSISHADETMONTECARLO_H
diff --git a/PWG4/totEt/AliAnalysisHadEtReconstructed.cxx b/PWG4/totEt/AliAnalysisHadEtReconstructed.cxx
new file mode 100644 (file)
index 0000000..feaa7d4
--- /dev/null
@@ -0,0 +1,97 @@
+
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
+//University of Tennessee at Knoxville
+#include "AliAnalysisHadEtReconstructed.h"
+#include "AliAnalysisEtCuts.h"
+#include "AliESDtrack.h"
+#include "AliESDCaloCluster.h"
+#include "TVector3.h"
+#include "AliVEvent.h"
+#include "AliESDEvent.h"
+#include "AliVParticle.h"
+#include <iostream>
+#include "TH2F.h"
+
+AliAnalysisHadEtReconstructed::AliAnalysisHadEtReconstructed() :
+        AliAnalysisHadEt()
+        ,fNTpcClustersCut(EtReconstructedCuts::kNTpcClustersCut)
+        ,fNItsClustersCut(EtReconstructedCuts::knItsClustersCut)
+        ,fTrackDistanceCut(0)
+{
+
+}
+
+Int_t AliAnalysisHadEtReconstructed::AnalyseEvent(AliVEvent* ev)
+{
+    ResetEventValues();
+    AliESDEvent *event = dynamic_cast<AliESDEvent*>(ev);
+
+    for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++)
+    {
+        AliVParticle *track = event->GetTrack(iTrack);
+        if (!track)
+        {
+            Printf("ERROR: Could not get track %d", iTrack);
+            continue;
+        }
+
+        fMultiplicity++;
+
+        const Double_t *pidWeights = track->PID();
+        if (pidWeights)
+        {
+            Int_t maxpid = -1;
+            Float_t maxpidweight = 0;
+            for (Int_t p =0; p < AliPID::kSPECIES; p++)
+            {
+                if (pidWeights[p] > maxpidweight)
+                {
+                    maxpidweight = pidWeights[p];
+                    maxpid = p;
+                }
+            }
+            if (maxpid == AliPID::kProton)
+            {
+                //     massPart = -0.938*track->Charge();
+            }
+
+        }
+
+    }
+
+    fTotNeutralEtAcc = fTotNeutralEt;
+    fTotEt = fTotChargedEt + fTotNeutralEt;
+    fTotEtAcc = fTotChargedEtAcc + fTotNeutralEtAcc;
+
+    std::cout << fTotChargedEtAcc << std::endl;
+    // Fill the histograms...
+    FillHistograms();
+
+    return 1;
+}
+
+bool AliAnalysisHadEtReconstructed::CheckGoodVertex(AliVParticle* track)
+{
+
+    Float_t bxy = 999.;
+    Float_t bz = 999.;
+    dynamic_cast<AliESDtrack*>(track)->GetImpactParametersTPC(bxy,bz);
+
+    return TMath::Abs(track->Xv()) < fVertexXCut && TMath::Abs(track->Yv()) < fVertexYCut && TMath::Abs(track->Zv()) < fVertexZCut && TMath::Abs(bxy) < fIPxyCut && TMath::Abs(bz) < fIPzCut;;
+
+}
+
+void AliAnalysisHadEtReconstructed::Init()
+{
+
+    AliAnalysisHadEt::Init();
+
+    fVertexXCut = EtReconstructedCuts::kVertexXCut;
+    fVertexYCut = EtReconstructedCuts::kVertexYCut;
+    fVertexZCut = EtReconstructedCuts::kVertexZCut;
+    fIPxyCut = EtReconstructedCuts::kIPxyCut;
+    fIPzCut = EtReconstructedCuts::kIPzCut;
+
+}
+
+
diff --git a/PWG4/totEt/AliAnalysisHadEtReconstructed.h b/PWG4/totEt/AliAnalysisHadEtReconstructed.h
new file mode 100644 (file)
index 0000000..d71bf1f
--- /dev/null
@@ -0,0 +1,33 @@
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
+//University of Tennessee at Knoxville
+#ifndef ALIANALYSISHADETRECONSTRUCTED_H
+#define ALIANALYSISHADETRECONSTRUCTED_H
+
+#include "AliAnalysisHadEt.h"
+
+class AliVParticle;
+
+class AliAnalysisHadEtReconstructed : public AliAnalysisHadEt
+{
+
+public:
+   
+    AliAnalysisHadEtReconstructed();
+   
+    virtual Int_t AnalyseEvent(AliVEvent* event);
+
+    virtual void Init();
+    
+protected:
+
+    bool CheckGoodVertex(AliVParticle *track);
+    //virtual bool TrackHitsCalorimeter(AliVParticle *track, Double_t magField) = 0;
+
+    Int_t fNTpcClustersCut;
+    Int_t fNItsClustersCut;
+   
+    Double_t fTrackDistanceCut;
+    
+};
+
+#endif // ALIANALYSISHADETRECONSTRUCTED_H
diff --git a/PWG4/totEt/AliAnalysisTaskHadEt.cxx b/PWG4/totEt/AliAnalysisTaskHadEt.cxx
new file mode 100644 (file)
index 0000000..937a074
--- /dev/null
@@ -0,0 +1,160 @@
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
+//University of Tennessee at Knoxville
+#include "TChain.h"
+#include "TTree.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TNtuple.h"
+#include "TCanvas.h"
+#include "TMath.h"
+#include "TObjArray.h"
+
+#include "AliAnalysisTask.h"
+#include "AliAnalysisManager.h"
+
+#include "AliESDEvent.h"
+#include "AliAODEvent.h"
+#include "AliVEvent.h"
+#include "AliVTrack.h"
+#include "AliVParticle.h"
+#include "AliESDtrack.h"
+#include "AliESDInputHandler.h"
+#include "AliMCEvent.h"
+#include "AliMCParticle.h"
+#include "TDatabasePDG.h"
+#include "AliAnalysisTaskHadEt.h"
+
+#include <iostream>
+#include "AliStack.h"
+
+using namespace std;
+
+ClassImp(AliAnalysisTaskHadEt)
+
+
+
+//________________________________________________________________________
+AliAnalysisTaskHadEt::AliAnalysisTaskHadEt(const char *name) :
+        AliAnalysisTaskSE(name)
+        ,fESD(0)
+        ,fOutputList(0)
+        ,fRecAnalysis(0)
+        ,fMCAnalysis(0)
+        ,fHistEtRecvsEtMC(0)
+        ,fTriggerSelection(false)
+        ,fCount(0)
+        ,fkPhotonPdg(22)
+        ,fkProtonMass(.938)
+        ,fPdgDB(0)
+        ,fRecEventVars(0)
+        ,fSimEventVars(0)
+        ,esdtrackCutsITSTPC(0)
+        ,esdtrackCutsTPC(0)
+        ,esdtrackCutsITS(0)
+{
+    // Constructor
+
+    fRecAnalysis = new AliAnalysisHadEtReconstructed();
+    fRecAnalysis->Init();
+    fMCAnalysis = new AliAnalysisHadEtMonteCarlo();
+    fMCAnalysis->Init();
+    if(!fPdgDB) fPdgDB = new TDatabasePDG();
+
+    fTriggerSelection = false;
+    fCount = 0;
+
+    // Define input and output slots here
+    // Input slot #0 works with a TChain
+    DefineInput(0, TChain::Class());
+    // Output slot #1 writes into a TH1 container
+
+    DefineOutput(1, TList::Class());
+
+}
+
+
+//________________________________________________________________________
+void AliAnalysisTaskHadEt::UserCreateOutputObjects()
+{
+    // Create histograms
+    // Called once
+  fOutputList = new TList;
+  fMCAnalysis->SetHistoList(fOutputList);
+  fRecAnalysis->SetHistoList(fOutputList);
+  fMCAnalysis->CreateHistograms();
+  fRecAnalysis->CreateHistograms();
+  fRecAnalysis->FillOutputList();
+  fMCAnalysis->FillOutputList();
+
+
+    Bool_t selectPrimaries=kTRUE;
+    esdtrackCutsITSTPC = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selectPrimaries);
+    esdtrackCutsITSTPC->SetName("fEsdTrackCuts");
+    esdtrackCutsTPC = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+    esdtrackCutsTPC->SetName("fEsdTrackCutsTPCOnly");
+    //ITS stand alone cuts - similar to 2009 cuts but with only ITS hits required
+    esdtrackCutsITS =  new AliESDtrackCuts;
+    esdtrackCutsITS->SetName("fEsdTrackCutsITS");
+    esdtrackCutsITS->SetRequireITSRefit(kTRUE);
+    esdtrackCutsITS->SetRequireITSStandAlone(kTRUE);
+    esdtrackCutsITS->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
+    esdtrackCutsITS->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/pt^0.9");
+    esdtrackCutsITS->SetMaxDCAToVertexZ(1.e6);
+    esdtrackCutsITS->SetDCAToVertex2D(kFALSE);
+    esdtrackCutsITS->SetRequireSigmaToVertex(kFALSE);
+    esdtrackCutsITS->SetAcceptKinkDaughters(kFALSE);
+
+    fOutputList->Add(esdtrackCutsITSTPC);
+    fOutputList->Add(esdtrackCutsTPC);
+    fOutputList->Add(esdtrackCutsITS);
+    if(esdtrackCutsITSTPC && esdtrackCutsTPC){
+      fRecAnalysis->SetITSTrackCuts( GetITSTrackCuts());
+      fMCAnalysis->SetITSTrackCuts( GetITSTrackCuts());
+      fRecAnalysis->SetTPCITSTrackCuts( GetTPCITSTrackCuts());
+      fMCAnalysis->SetTPCITSTrackCuts( GetTPCITSTrackCuts());
+      fRecAnalysis->SetTPCOnlyTrackCuts( GetTPCOnlyTrackCuts());
+      fMCAnalysis->SetTPCOnlyTrackCuts( GetTPCOnlyTrackCuts());
+      //add ITS stuff!
+    }
+    else{
+      Printf("Error: no track cuts!");
+    }
+
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskHadEt::UserExec(Option_t *)
+{
+  AliESDEvent *event = dynamic_cast<AliESDEvent*>(InputEvent());
+if (!event) {
+  Printf("ERROR: Could not retrieve event");
+  return;
+ }
+
+fRecAnalysis->AnalyseEvent(event);
+
+AliMCEvent* mcEvent = MCEvent();
+if (mcEvent)
+  {
+    ((AliAnalysisHadEtMonteCarlo*)fMCAnalysis)->AnalyseEvent((AliVEvent*)mcEvent,(AliVEvent*)event);
+  }
+
+// Post output data.
+PostData(1, fOutputList);
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskHadEt::Terminate(Option_t *)
+{
+    // Draw result to the screen
+    // Called once at the end of the query
+
+    fOutputList = dynamic_cast<TList*> (GetOutputData(1));
+    if (!fOutputList) {
+        printf("ERROR: Output list not available\n");
+        return;
+    }
+}
+
+
diff --git a/PWG4/totEt/AliAnalysisTaskHadEt.h b/PWG4/totEt/AliAnalysisTaskHadEt.h
new file mode 100644 (file)
index 0000000..73dcc07
--- /dev/null
@@ -0,0 +1,113 @@
+//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili\r
+//University of Tennessee at Knoxville\r
+#ifndef AliAnalysisTaskHadEt_cxx\r
+#define AliAnalysisTaskHadEt_cxx\r
+\r
+class AliAnalysisHadEt;\r
+class TTree;\r
+class AliVParticle;\r
+class TH1F;\r
+class TH2F;\r
+class TNtuple;\r
+class TObjArray;\r
+class AliESDEvent;\r
+class AliMCParticle;\r
+class TDatabasePDG;\r
+\r
+#include "AliAnalysisTaskSE.h"\r
+#include "TObject.h"\r
+#include "AliAnalysisHadEtReconstructed.h"\r
+#include "AliAnalysisHadEtMonteCarlo.h"\r
+\r
+/* class ParticleVars : public TObject        // Inherit from TObject to put in TClonesArray */\r
+/*     { */\r
+/*        public: */\r
+         \r
+/*       ParticleVars() : TObject(){} */\r
+/*        Int_t fPdgCode; // from MC */\r
+/*        Int_t fPid; //from ESDs */\r
+/*        Int_t fMass; */\r
+/*        Int_t fCharge; */\r
+/*        Double_t fEt; */\r
+/*        Double_t fPhi; */\r
+/*        Double_t fEta; */\r
+\r
+/*        ClassDef(ParticleVars, 1); */\r
+       \r
+/*     }; */\r
+\r
+class AliAnalysisTaskHadEt : public AliAnalysisTaskSE {\r
+public:\r
+    AliAnalysisTaskHadEt(const char *name = "AliAnalysisTaskHadEt");\r
+    virtual ~AliAnalysisTaskHadEt() {}\r
+\r
+    //  virtual void   ConnectInputData(Option_t *);\r
+    virtual void   UserCreateOutputObjects();\r
+    virtual void   UserExec(Option_t *option);\r
+    virtual void   Terminate(Option_t *);\r
+\r
+    virtual void SetTriggerSelection(Bool_t v) {\r
+        fTriggerSelection = v;\r
+    }\r
+\r
+    AliESDtrackCuts* GetTPCITSTrackCuts(){return (AliESDtrackCuts*) fOutputList->FindObject("fEsdTrackCuts");}\r
+    AliESDtrackCuts* GetTPCOnlyTrackCuts(){return (AliESDtrackCuts*) fOutputList->FindObject("fEsdTrackCutsTPCOnly");}\r
+    AliESDtrackCuts* GetITSTrackCuts(){return (AliESDtrackCuts*) fOutputList->FindObject("fEsdTrackCutsITS");}\r
+\r
+private:\r
+\r
+  //Declare it private to avoid compilation warning\r
+    AliAnalysisTaskHadEt & operator = (const AliAnalysisTaskHadEt & g) ;//cpy assignment\r
+    AliAnalysisTaskHadEt(const AliAnalysisTaskHadEt & g) ; // cpy ctor\r
+\r
+    bool CheckGoodVertex(AliVParticle *track);\r
+    bool TrackHits(AliVParticle *track, Double_t magField);\r
+\r
+\r
+    AliESDEvent *fESD;    //ESD object\r
+\r
+    TList *fOutputList;\r
+\r
+    AliAnalysisHadEtReconstructed *fRecAnalysis;\r
+    AliAnalysisHadEtMonteCarlo *fMCAnalysis;\r
+\r
+    TH2F *fHistEtRecvsEtMC;\r
+    \r
+    Bool_t fTriggerSelection;\r
+\r
+    Int_t fCount;\r
+\r
+    const int fkPhotonPdg;\r
+\r
+    const Float_t fkProtonMass;\r
+\r
+    TDatabasePDG *fPdgDB;\r
+\r
+    class EventVars\r
+    {\r
+       public:\r
+        Double_t fTotEt;\r
+       Double_t fTotEtAcc;\r
+        Double_t fTotEnergy;\r
+\r
+        Double_t fTotNeutralEt;\r
+       Double_t fTotNeutralEtAcc;\r
+\r
+        Double_t fTotChargedEt;\r
+        Double_t fTotChargedEtAcc;\r
+\r
+        Int_t fChargedMultiplicity;\r
+        Int_t fNeutralMultiplicity;\r
+\r
+    };\r
+    \r
+    EventVars *fRecEventVars;\r
+    EventVars *fSimEventVars;\r
+    AliESDtrackCuts* esdtrackCutsITSTPC;\r
+    AliESDtrackCuts* esdtrackCutsTPC;\r
+    AliESDtrackCuts* esdtrackCutsITS;\r
+    \r
+    ClassDef(AliAnalysisTaskHadEt, 1); // example of analysis\r
+};\r
+\r
+#endif\r