[u/mrichter/AliRoot.git] / PWG4 / totEt / AliAnalysisEt.cxx

//_________________________________________________________________________
//  Utility Class for transverse energy studies
//  Base class for ESD & MC analysis
//  - reconstruction and MonteCarlo output
// implementation file
//
//*-- Authors: Oystein Djuvsland (Bergen), David Silvermyr (ORNL)
//_________________________________________________________________________

#include "AliAnalysisEt.h"
#include "TMath.h"
#include "TList.h"
#include "TH1F.h"
#include "TH2F.h"
#include "TTree.h"
#include <iostream>
#include "AliAnalysisEtCuts.h"
#include "AliVEvent.h"
#include "TDatabasePDG.h"
#include "Rtypes.h"

using namespace std;
ClassImp(AliAnalysisEt);


AliAnalysisEt::AliAnalysisEt() :
        fHistogramNameSuffix("")
        ,fCuts(0)
        ,fPdgDB(0)
        ,fPiPlusCode(0)
        ,fPiMinusCode(0)
        ,fKPlusCode(0)
        ,fKMinusCode(0)
        ,fProtonCode(0)
        ,fAntiProtonCode(0)
        ,fLambdaCode(0)
        ,fAntiLambdaCode(0)
        ,fK0SCode(0)
        ,fOmegaCode(0)
        ,fAntiOmegaCode(0)
        ,fXi0Code(0)
        ,fAntiXi0Code(0)
        ,fXiCode(0)
        ,fAntiXiCode(0)
        ,fSigmaCode(0)
        ,fAntiSigmaCode(0)
        ,fK0LCode(0)
        ,fNeutronCode(0)
        ,fAntiNeutronCode(0)
        ,fEPlusCode(0)
        ,fEMinusCode(0)
        ,fPionMass(0)
        ,fTotEt(0)
        ,fTotEtAcc(0)
        ,fTotNeutralEt(0)
        ,fTotNeutralEtAcc(0)
        ,fTotChargedEt(0)
        ,fTotChargedEtAcc(0)
        ,fMultiplicity(0)
        ,fChargedMultiplicity(0)
        ,fNeutralMultiplicity(0)
        ,fBaryonEt(0)
        ,fAntiBaryonEt(0)
        ,fMesonEt(0)
        ,fBaryonEtAcc(0)
        ,fAntiBaryonEtAcc(0)
        ,fMesonEtAcc(0)
        ,fProtonEt(0)
        ,fChargedKaonEt(0)
        ,fMuonEt(0)
        ,fElectronEt(0)
        ,fProtonEtAcc(0)
        ,fChargedKaonEtAcc(0)
        ,fMuonEtAcc(0)
        ,fElectronEtAcc(0)
        ,fEtaCut(0)
	,fEtaCutAcc(0)
	,fPhiCutAccMin(0)
	,fPhiCutAccMax(0)
	,fDetectorRadius(0)
	,fClusterEnergyCut(0) 
	,fSingleCellEnergyCut(0)
        ,fHistEt(0)
        ,fHistChargedEt(0)
        ,fHistNeutralEt(0)
        ,fHistEtAcc(0)
        ,fHistChargedEtAcc(0)
        ,fHistNeutralEtAcc(0)
        ,fHistMult(0)
        ,fHistChargedMult(0)
        ,fHistNeutralMult(0)
        ,fHistPhivsPtPos(0)
        ,fHistPhivsPtNeg(0)
        ,fHistBaryonEt(0)
        ,fHistAntiBaryonEt(0)
        ,fHistMesonEt(0)
        ,fHistBaryonEtAcc(0)
        ,fHistAntiBaryonEtAcc(0)
        ,fHistMesonEtAcc(0)
        ,fHistProtonEt(0)
        ,fHistChargedKaonEt(0)
        ,fHistMuonEt(0)
        ,fHistElectronEt(0)
        ,fHistProtonEtAcc(0)
        ,fHistChargedKaonEtAcc(0)
        ,fHistMuonEtAcc(0)
        ,fHistElectronEtAcc(0)
        ,fHistEtRecvsEtMC(0)
        ,fHistTMDeltaR(0)
        ,fTree(0)
{
}

AliAnalysisEt::~AliAnalysisEt()
{

}

void AliAnalysisEt::FillOutputList(TList *list)
{ // histograms to be added to output
    list->Add(fHistEt);
    list->Add(fHistChargedEt);
    list->Add(fHistNeutralEt);

    list->Add(fHistEtAcc);
    list->Add(fHistChargedEtAcc);
    list->Add(fHistNeutralEtAcc);

    list->Add(fHistMult);
    list->Add(fHistChargedMult);
    list->Add(fHistNeutralMult);

    list->Add(fHistPhivsPtPos);
    list->Add(fHistPhivsPtNeg);

    list->Add(fHistBaryonEt);
    list->Add(fHistAntiBaryonEt);
    list->Add(fHistMesonEt);

    list->Add(fHistBaryonEtAcc);
    list->Add(fHistAntiBaryonEtAcc);
    list->Add(fHistMesonEtAcc);

    list->Add(fHistProtonEtAcc);
    list->Add(fHistChargedKaonEtAcc);
    list->Add(fHistMuonEtAcc);
    list->Add(fHistElectronEtAcc);

    list->Add(fHistEtRecvsEtMC);

    list->Add(fHistTMDeltaR);

    if (fCuts) {
      if (fCuts->GetHistMakeTree()) {
	list->Add(fTree);
      }
    }

}

void AliAnalysisEt::Init()
{// clear variables, set up cuts and PDG info
  ResetEventValues();
}

void AliAnalysisEt::CreateHistograms()
{ // create histograms..
  // histogram binning for E_T, p_T and Multiplicity: defaults for p+p
  Int_t nbinsEt = 1000;
  Double_t minEt = 0.0001;
  Double_t maxEt = 100;
  Int_t nbinsPt = 200;
  Double_t minPt = 0;
  Double_t maxPt = 20;
  Int_t nbinsMult = 200;
  Double_t minMult = -0.5; // offset -0.5 to have integer bins centered around 0
  Double_t maxMult = nbinsMult + minMult; // 1 bin per integer value

  // see if we should change histogram limits etc, and possibly create a tree
  if (fCuts) {
    if (fCuts->GetHistMakeTree()) {
      CreateTree();
    }
  }

    TString histname = "fHistEt" + fHistogramNameSuffix;
    fHistEt = new TH1F(histname.Data(), "Total E_{T} Distribution", nbinsEt, minEt, maxEt);
    fHistEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
    fHistEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");

    histname = "fHistChargedEt" + fHistogramNameSuffix;
    fHistChargedEt = new TH1F(histname.Data(), "Total Charged E_{T} Distribution", nbinsEt, minEt, maxEt);
    fHistChargedEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
    fHistChargedEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");

    histname = "fHistNeutralEt" + fHistogramNameSuffix;
    fHistNeutralEt = new TH1F(histname.Data(), "Total Neutral E_{T} Distribution", nbinsEt, minEt, maxEt);
    fHistNeutralEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
    fHistNeutralEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");

    histname = "fHistEtAcc" + fHistogramNameSuffix;
    fHistEtAcc = new TH1F(histname.Data(), "Total E_{T} Distribution in Acceptance", nbinsEt, minEt, maxEt);
    fHistEtAcc->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
    fHistEtAcc->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");

    histname = "fHistChargedEtAcc" + fHistogramNameSuffix;
    fHistChargedEtAcc = new TH1F(histname.Data(), "Total Charged E_{T} Distribution in Acceptance", nbinsEt, minEt, maxEt);
    fHistChargedEtAcc->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
    fHistChargedEtAcc->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");

    histname = "fHistNeutralEtAcc" + fHistogramNameSuffix;
    fHistNeutralEtAcc = new TH1F(histname.Data(), "Total Neutral E_{T} Distribution in Acceptance", nbinsEt, minEt, maxEt);
    fHistNeutralEtAcc->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
    fHistNeutralEtAcc->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
    std::cout << histname << std::endl;
    histname = "fHistMult" + fHistogramNameSuffix;
    fHistMult = new TH1F(histname.Data(), "Total Multiplicity", nbinsMult, minMult, maxMult);
    fHistMult->GetXaxis()->SetTitle("N");
    fHistMult->GetYaxis()->SetTitle("Multiplicity");

    histname = "fHistChargedMult" + fHistogramNameSuffix;
    fHistChargedMult = new TH1F(histname.Data(), "Charged Multiplicity", nbinsMult, minMult, maxMult);
    fHistChargedMult->GetXaxis()->SetTitle("N");
    fHistChargedMult->GetYaxis()->SetTitle("Multiplicity");

    histname = "fHistNeutralMult" + fHistogramNameSuffix;
    fHistNeutralMult = new TH1F(histname.Data(), "Neutral Multiplicity", nbinsMult, minMult, maxMult);
    fHistNeutralMult->GetXaxis()->SetTitle("N");
    fHistNeutralMult->GetYaxis()->SetTitle("Multiplicity");

    histname = "fHistPhivsPtPos" + fHistogramNameSuffix;
    fHistPhivsPtPos = new TH2F(histname.Data(), "Phi vs pT of positively charged tracks hitting the calorimeter", 	200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);

    histname = "fHistPhivsPtNeg" + fHistogramNameSuffix;
    fHistPhivsPtNeg = new TH2F(histname.Data(), "Phi vs pT of negatively charged tracks hitting the calorimeter", 	200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);

    histname = "fHistBaryonEt" + fHistogramNameSuffix;
    fHistBaryonEt = new TH1F(histname.Data(), "E_{T} for baryons",  nbinsEt, minEt, maxEt);

    histname = "fHistAntiBaryonEt" + fHistogramNameSuffix;
    fHistAntiBaryonEt = new TH1F(histname.Data(), "E_{T} for anti baryons",  nbinsEt, minEt, maxEt);

    histname = "fHistMesonEt" + fHistogramNameSuffix;
    fHistMesonEt = new TH1F(histname.Data(), "E_{T} for mesons",  nbinsEt, minEt, maxEt);

    histname = "fHistBaryonEtAcc" + fHistogramNameSuffix;
    fHistBaryonEtAcc = new TH1F(histname.Data(), "E_{T} for baryons in calorimeter acceptance",  nbinsEt, minEt, maxEt);

    histname = "fHistAntiBaryonEtAcc" + fHistogramNameSuffix;
    fHistAntiBaryonEtAcc = new TH1F(histname.Data(), "E_{T} for anti baryons in calorimeter acceptance",  nbinsEt, minEt, maxEt);

    histname = "fHistMesonEtAcc" + fHistogramNameSuffix;
    fHistMesonEtAcc = new TH1F(histname.Data(), "E_{T} for mesons in calorimeter acceptance",  nbinsEt, minEt, maxEt);

    histname = "fHistProtonEt" + fHistogramNameSuffix;
    fHistProtonEt = new TH1F(histname.Data(), "E_{T} for (anti-)protons", nbinsEt, minEt, maxEt);

    histname = "fHistKaonEt" + fHistogramNameSuffix;
    fHistChargedKaonEt = new TH1F(histname.Data(), "E_{T} for charged kaons", nbinsEt, minEt, maxEt);

    histname = "fHistMuonEt" + fHistogramNameSuffix;
    fHistMuonEt = new TH1F(histname.Data(), "E_{T} for muons", nbinsEt, minEt, maxEt);

    histname = "fHistElectronEt" + fHistogramNameSuffix;
    fHistElectronEt = new TH1F(histname.Data(), "E_{T} for electrons/positrons", nbinsEt, minEt, maxEt);

    histname = "fHistProtonEtAcc" + fHistogramNameSuffix;
    fHistProtonEtAcc = new TH1F(histname.Data(), "E_{T} for (anti-)protons in calorimeter acceptance", nbinsEt, minEt, maxEt);

    histname = "fHistKaonEtAcc" + fHistogramNameSuffix;
    fHistChargedKaonEtAcc = new TH1F(histname.Data(), "E_{T} for charged kaons in calorimeter acceptance", nbinsEt, minEt, maxEt);

    histname = "fHistMuonEtAcc" + fHistogramNameSuffix;
    fHistMuonEtAcc = new TH1F(histname.Data(), "E_{T} for muons in calorimeter acceptance", nbinsEt, minEt, maxEt);

    histname = "fHistElectronEtAcc" + fHistogramNameSuffix;
    fHistElectronEtAcc = new TH1F(histname.Data(), "E_{T} for electrons/positrons in calorimeter acceptance", nbinsEt, minEt, maxEt);

    histname = "fHistEtRecvsEtMC" + fHistogramNameSuffix;
    fHistEtRecvsEtMC = new TH2F(histname.Data(), "Reconstructed E_{t} vs MC E_{t}", nbinsEt, minEt, maxEt, nbinsEt, minEt, maxEt);

    //
    histname = "fHistTMDeltaR" + fHistogramNameSuffix;
    fHistTMDeltaR = new TH1F(histname.Data(), "#Delta R for calorimeter clusters", 200, 0, 50);

}

void AliAnalysisEt::CreateTree()
{ // create tree..
  TString treename = "fTree" + fHistogramNameSuffix;
  fTree = new TTree(treename, treename);
  fTree->Branch("fTotEt",&fTotEt,"fTotEt/D");
  fTree->Branch("fTotEtAcc",&fTotEtAcc,"fTotEtAcc/D");
  fTree->Branch("fTotNeutralEt",&fTotNeutralEt,"fTotNeutralEt/D");
  fTree->Branch("fTotNeutralEtAcc",&fTotNeutralEtAcc,"fTotNeutralEtAcc/D");
  fTree->Branch("fTotChargedEt",&fTotChargedEt,"fTotChargedEt/D");
  fTree->Branch("fTotChargedEtAcc",&fTotChargedEtAcc,"fTotChargedEtAcc/D");
  fTree->Branch("fMultiplicity",&fMultiplicity,"fMultiplicity/I");
  fTree->Branch("fChargedMultiplicity",&fChargedMultiplicity,"fChargedMultiplicity/I");
  fTree->Branch("fNeutralMultiplicity",&fNeutralMultiplicity,"fNeutralMultiplicity/I");
  fTree->Branch("fBaryonEt",&fBaryonEt,"fBaryonEt/D");
  fTree->Branch("fAntiBaryonEt",&fAntiBaryonEt,"fAntiBaryonEt/D");
  fTree->Branch("fMesonEt",&fMesonEt,"fMesonEt/D");
  fTree->Branch("fBaryonEtAcc",&fBaryonEtAcc,"fBaryonEtAcc/D");
  fTree->Branch("fAntiBaryonEtAcc",&fAntiBaryonEtAcc,"fAntiBaryonEtAcc/D");
  fTree->Branch("fMesonEtAcc",&fMesonEtAcc,"fMesonEtAcc/D");
  fTree->Branch("fProtonEt",&fProtonEt,"fProtonEt/D");
  fTree->Branch("fChargedKaonEt",&fChargedKaonEt,"fChargedKaonEt/D");
  fTree->Branch("fMuonEt",&fMuonEt,"fMuonEt/D");
  fTree->Branch("fElectronEt",&fElectronEt,"fElectronEt/D");
  fTree->Branch("fProtonEtAcc",&fProtonEtAcc,"fProtonEtAcc/D");
  fTree->Branch("fChargedKaonEtAcc",&fChargedKaonEtAcc,"fChargedKaonEtAcc/D");
  fTree->Branch("fMuonEtAcc",&fMuonEtAcc,"fMuonEtAcc/D");
  fTree->Branch("fElectronEtAcc",&fElectronEtAcc,"fElectronEtAcc/D");

  return;
}

void AliAnalysisEt::FillHistograms()
{ // fill histograms..
    fHistEt->Fill(fTotEt);
    fHistChargedEt->Fill(fTotChargedEt);
    fHistNeutralEt->Fill(fTotNeutralEt);

    fHistEtAcc->Fill(fTotEtAcc);
    fHistChargedEtAcc->Fill(fTotChargedEtAcc);
    fHistNeutralEtAcc->Fill(fTotNeutralEtAcc);

    fHistMult->Fill(fMultiplicity);
    fHistChargedMult->Fill(fChargedMultiplicity);
    fHistNeutralMult->Fill(fNeutralMultiplicity);

    fHistBaryonEt->Fill(fBaryonEt);
    fHistAntiBaryonEt->Fill(fAntiBaryonEt);
    fHistMesonEt->Fill(fMesonEt);

    fHistBaryonEtAcc->Fill(fBaryonEtAcc);
    fHistAntiBaryonEtAcc->Fill(fAntiBaryonEtAcc);
    fHistMesonEtAcc->Fill(fMesonEtAcc);
 
    fHistProtonEt->Fill(fProtonEt);
    fHistChargedKaonEt->Fill(fChargedKaonEt);
    fHistMuonEt->Fill(fMuonEt);
    fHistElectronEt->Fill(fElectronEt);
    
    fHistProtonEtAcc->Fill(fProtonEtAcc);
    fHistChargedKaonEtAcc->Fill(fChargedKaonEtAcc);
    fHistMuonEtAcc->Fill(fMuonEtAcc);
    fHistElectronEtAcc->Fill(fElectronEtAcc);

    if (fCuts) {
      if (fCuts->GetHistMakeTree()) {
	fTree->Fill();
      }
    }

}

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;
  ResetEventValues();
  return 0;
}

void AliAnalysisEt::ResetEventValues()
{ // clear
  fTotEt = 0;
  fTotEtAcc = 0;
  fTotNeutralEt = 0;
  fTotNeutralEtAcc = 0;
  fTotChargedEt  = 0;
  fTotChargedEtAcc = 0;
  fMultiplicity = 0;
  fChargedMultiplicity = 0;
  fNeutralMultiplicity = 0;
  fBaryonEt = 0;
  fAntiBaryonEt = 0;
  fMesonEt = 0;
  fBaryonEtAcc = 0;
  fAntiBaryonEtAcc = 0;
  fMesonEtAcc = 0;
  fProtonEt = 0;
  fChargedKaonEt = 0;
  fMuonEt = 0;
  fElectronEt = 0;
  fProtonEtAcc = 0;
  fChargedKaonEtAcc = 0;
  fMuonEtAcc = 0;
  fElectronEtAcc = 0;
    
  if (!fCuts || !fPdgDB || fPiPlusCode==0) { // some Init's needed
    cout << __FILE__ << ":" << __LINE__ << " : Init " << endl;
    if (!fCuts) {
      cout << " setting up Cuts " << endl;
      fCuts = new AliAnalysisEtCuts();
    }
    if(!fPdgDB) {
      cout << " setting up PdgDB " << endl;
      fPdgDB = new TDatabasePDG();
    }
    
    if (fPiPlusCode==0) {
      SetParticleCodes();
    }
  }
  return;
}


void AliAnalysisEt::SetParticleCodes()
{ // set PDG info    
  fPionMass = fPdgDB->GetParticle("pi+")->Mass();
  fPiPlusCode = fPdgDB->GetParticle("pi+")->PdgCode();
  fPiMinusCode = fPdgDB->GetParticle("pi-")->PdgCode();
  fKPlusCode = fPdgDB->GetParticle("K+")->PdgCode();
  fKMinusCode = fPdgDB->GetParticle("K-")->PdgCode();
  fProtonCode = fPdgDB->GetParticle("proton")->PdgCode();
  fAntiProtonCode = fPdgDB->GetParticle("antiproton")->PdgCode();
  fLambdaCode = fPdgDB->GetParticle("Lambda0")->PdgCode();
  fAntiLambdaCode = fPdgDB->GetParticle("Lambda0_bar")->PdgCode();
  fK0SCode = fPdgDB->GetParticle("K_S0")->PdgCode();
  fOmegaCode = fPdgDB->GetParticle("Omega-")->PdgCode();
  fAntiOmegaCode = fPdgDB->GetParticle("Omega+")->PdgCode();
  fXi0Code = fPdgDB->GetParticle("Xi0")->PdgCode();
  fAntiXi0Code = fPdgDB->GetParticle("Xi0_bar")->PdgCode();
  fXiCode = fPdgDB->GetParticle("Xi-")->PdgCode();
  fAntiXiCode = fPdgDB->GetParticle("Xi-_bar")->PdgCode();
  fSigmaCode = fPdgDB->GetParticle("Sigma-")->PdgCode();
  fAntiSigmaCode = fPdgDB->GetParticle("Sigma+")->PdgCode();
  fK0LCode = fPdgDB->GetParticle("K_L0")->PdgCode();
  fNeutronCode = fPdgDB->GetParticle("neutron")->PdgCode();
  fAntiNeutronCode = fPdgDB->GetParticle("antineutron")->PdgCode();
  fEPlusCode = fPdgDB->GetParticle("e+")->PdgCode();
  fEMinusCode = fPdgDB->GetParticle("e-")->PdgCode();
  cout << "Resetting Codes: Pion " << fPiPlusCode
       << "," << fPiMinusCode 
       << " Kaon " << fKPlusCode 
       << "," << fKMinusCode << endl;
}
Commit	Line	Data
cf6522d1	1	//_________________________________________________________________________
	2	// Utility Class for transverse energy studies
	3	// Base class for ESD & MC analysis
	4	// - reconstruction and MonteCarlo output
	5	// implementation file
	6	//
	7	//*-- Authors: Oystein Djuvsland (Bergen), David Silvermyr (ORNL)
	8	//_________________________________________________________________________
	9
2fbf38ac	10	#include "AliAnalysisEt.h"
	11	#include "TMath.h"
	12	#include "TList.h"
	13	#include "TH1F.h"
	14	#include "TH2F.h"
ce546038	15	#include "TTree.h"
2fbf38ac	16	#include <iostream>
2fbf38ac	17	#include "AliAnalysisEtCuts.h"
641e1e0c	18	#include "AliVEvent.h"
cf6522d1	19	#include "TDatabasePDG.h"
cf6522d1	20	#include "Rtypes.h"
641e1e0c	21
641e1e0c	22	using namespace std;
2fbf38ac	23	ClassImp(AliAnalysisEt);
	24
	25
	26	AliAnalysisEt::AliAnalysisEt() :
	27	fHistogramNameSuffix("")
4998becf	28	,fCuts(0)
2fbf38ac	29	,fPdgDB(0)
8985c506	30	,fPiPlusCode(0)
	31	,fPiMinusCode(0)
	32	,fKPlusCode(0)
	33	,fKMinusCode(0)
	34	,fProtonCode(0)
	35	,fAntiProtonCode(0)
	36	,fLambdaCode(0)
	37	,fAntiLambdaCode(0)
	38	,fK0SCode(0)
	39	,fOmegaCode(0)
	40	,fAntiOmegaCode(0)
	41	,fXi0Code(0)
	42	,fAntiXi0Code(0)
	43	,fXiCode(0)
	44	,fAntiXiCode(0)
	45	,fSigmaCode(0)
	46	,fAntiSigmaCode(0)
	47	,fK0LCode(0)
	48	,fNeutronCode(0)
	49	,fAntiNeutronCode(0)
	50	,fEPlusCode(0)
	51	,fEMinusCode(0)
	52	,fPionMass(0)
2fbf38ac	53	,fTotEt(0)
	54	,fTotEtAcc(0)
	55	,fTotNeutralEt(0)
	56	,fTotNeutralEtAcc(0)
	57	,fTotChargedEt(0)
	58	,fTotChargedEtAcc(0)
	59	,fMultiplicity(0)
	60	,fChargedMultiplicity(0)
	61	,fNeutralMultiplicity(0)
87efb15c	62	,fBaryonEt(0)
	63	,fAntiBaryonEt(0)
	64	,fMesonEt(0)
	65	,fBaryonEtAcc(0)
	66	,fAntiBaryonEtAcc(0)
	67	,fMesonEtAcc(0)
	68	,fProtonEt(0)
	69	,fChargedKaonEt(0)
	70	,fMuonEt(0)
	71	,fElectronEt(0)
	72	,fProtonEtAcc(0)
	73	,fChargedKaonEtAcc(0)
	74	,fMuonEtAcc(0)
	75	,fElectronEtAcc(0)
	76	,fEtaCut(0)
4998becf	77	,fEtaCutAcc(0)
	78	,fPhiCutAccMin(0)
	79	,fPhiCutAccMax(0)
	80	,fDetectorRadius(0)
	81	,fClusterEnergyCut(0)
	82	,fSingleCellEnergyCut(0)
2fbf38ac	83	,fHistEt(0)
	84	,fHistChargedEt(0)
	85	,fHistNeutralEt(0)
	86	,fHistEtAcc(0)
	87	,fHistChargedEtAcc(0)
	88	,fHistNeutralEtAcc(0)
	89	,fHistMult(0)
	90	,fHistChargedMult(0)
	91	,fHistNeutralMult(0)
	92	,fHistPhivsPtPos(0)
	93	,fHistPhivsPtNeg(0)
	94	,fHistBaryonEt(0)
	95	,fHistAntiBaryonEt(0)
	96	,fHistMesonEt(0)
	97	,fHistBaryonEtAcc(0)
	98	,fHistAntiBaryonEtAcc(0)
	99	,fHistMesonEtAcc(0)
87efb15c	100	,fHistProtonEt(0)
	101	,fHistChargedKaonEt(0)
	102	,fHistMuonEt(0)
	103	,fHistElectronEt(0)
	104	,fHistProtonEtAcc(0)
	105	,fHistChargedKaonEtAcc(0)
	106	,fHistMuonEtAcc(0)
	107	,fHistElectronEtAcc(0)
	108	,fHistEtRecvsEtMC(0)
2fbf38ac	109	,fHistTMDeltaR(0)
ce546038	110	,fTree(0)
2fbf38ac	111	{
2fbf38ac	112	}
	113
	114	AliAnalysisEt::~AliAnalysisEt()
	115	{
	116
	117	}
	118
	119	void AliAnalysisEt::FillOutputList(TList *list)
cf6522d1	120	{ // histograms to be added to output
2fbf38ac	121	list->Add(fHistEt);
	122	list->Add(fHistChargedEt);
	123	list->Add(fHistNeutralEt);
	124
	125	list->Add(fHistEtAcc);
	126	list->Add(fHistChargedEtAcc);
	127	list->Add(fHistNeutralEtAcc);
	128
	129	list->Add(fHistMult);
	130	list->Add(fHistChargedMult);
	131	list->Add(fHistNeutralMult);
	132
	133	list->Add(fHistPhivsPtPos);
	134	list->Add(fHistPhivsPtNeg);
	135
	136	list->Add(fHistBaryonEt);
	137	list->Add(fHistAntiBaryonEt);
	138	list->Add(fHistMesonEt);
	139
	140	list->Add(fHistBaryonEtAcc);
	141	list->Add(fHistAntiBaryonEtAcc);
	142	list->Add(fHistMesonEtAcc);
	143
87efb15c	144	list->Add(fHistProtonEtAcc);
	145	list->Add(fHistChargedKaonEtAcc);
	146	list->Add(fHistMuonEtAcc);
	147	list->Add(fHistElectronEtAcc);
	148
	149	list->Add(fHistEtRecvsEtMC);
	150
2fbf38ac	151	list->Add(fHistTMDeltaR);
ce546038	152
	153	if (fCuts) {
	154	if (fCuts->GetHistMakeTree()) {
	155	list->Add(fTree);
	156	}
	157	}
	158
2fbf38ac	159	}
	160
	161	void AliAnalysisEt::Init()
0065531a	162	{// clear variables, set up cuts and PDG info
0065531a	163	ResetEventValues();
2fbf38ac	164	}
	165
	166	void AliAnalysisEt::CreateHistograms()
cf6522d1	167	{ // create histograms..
13b0d3c1	168	// histogram binning for E_T, p_T and Multiplicity: defaults for p+p
	169	Int_t nbinsEt = 1000;
	170	Double_t minEt = 0.0001;
	171	Double_t maxEt = 100;
	172	Int_t nbinsPt = 200;
	173	Double_t minPt = 0;
	174	Double_t maxPt = 20;
	175	Int_t nbinsMult = 200;
	176	Double_t minMult = -0.5; // offset -0.5 to have integer bins centered around 0
	177	Double_t maxMult = nbinsMult + minMult; // 1 bin per integer value
2fbf38ac	178
ce546038	179	// see if we should change histogram limits etc, and possibly create a tree
	180	if (fCuts) {
	181	if (fCuts->GetHistMakeTree()) {
	182	CreateTree();
	183	}
	184	}
	185
2fbf38ac	186	TString histname = "fHistEt" + fHistogramNameSuffix;
13b0d3c1	187	fHistEt = new TH1F(histname.Data(), "Total E_{T} Distribution", nbinsEt, minEt, maxEt);
2fbf38ac	188	fHistEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
	189	fHistEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
	190
	191	histname = "fHistChargedEt" + fHistogramNameSuffix;
13b0d3c1	192	fHistChargedEt = new TH1F(histname.Data(), "Total Charged E_{T} Distribution", nbinsEt, minEt, maxEt);
2fbf38ac	193	fHistChargedEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
	194	fHistChargedEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
	195
	196	histname = "fHistNeutralEt" + fHistogramNameSuffix;
13b0d3c1	197	fHistNeutralEt = new TH1F(histname.Data(), "Total Neutral E_{T} Distribution", nbinsEt, minEt, maxEt);
2fbf38ac	198	fHistNeutralEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
	199	fHistNeutralEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
	200
	201	histname = "fHistEtAcc" + fHistogramNameSuffix;
13b0d3c1	202	fHistEtAcc = new TH1F(histname.Data(), "Total E_{T} Distribution in Acceptance", nbinsEt, minEt, maxEt);
2fbf38ac	203	fHistEtAcc->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
	204	fHistEtAcc->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
	205
	206	histname = "fHistChargedEtAcc" + fHistogramNameSuffix;
13b0d3c1	207	fHistChargedEtAcc = new TH1F(histname.Data(), "Total Charged E_{T} Distribution in Acceptance", nbinsEt, minEt, maxEt);
2fbf38ac	208	fHistChargedEtAcc->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
	209	fHistChargedEtAcc->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
	210
	211	histname = "fHistNeutralEtAcc" + fHistogramNameSuffix;
13b0d3c1	212	fHistNeutralEtAcc = new TH1F(histname.Data(), "Total Neutral E_{T} Distribution in Acceptance", nbinsEt, minEt, maxEt);
2fbf38ac	213	fHistNeutralEtAcc->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
	214	fHistNeutralEtAcc->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
	215	std::cout << histname << std::endl;
	216	histname = "fHistMult" + fHistogramNameSuffix;
13b0d3c1	217	fHistMult = new TH1F(histname.Data(), "Total Multiplicity", nbinsMult, minMult, maxMult);
2fbf38ac	218	fHistMult->GetXaxis()->SetTitle("N");
	219	fHistMult->GetYaxis()->SetTitle("Multiplicity");
	220
	221	histname = "fHistChargedMult" + fHistogramNameSuffix;
13b0d3c1	222	fHistChargedMult = new TH1F(histname.Data(), "Charged Multiplicity", nbinsMult, minMult, maxMult);
2fbf38ac	223	fHistChargedMult->GetXaxis()->SetTitle("N");
	224	fHistChargedMult->GetYaxis()->SetTitle("Multiplicity");
	225
	226	histname = "fHistNeutralMult" + fHistogramNameSuffix;
13b0d3c1	227	fHistNeutralMult = new TH1F(histname.Data(), "Neutral Multiplicity", nbinsMult, minMult, maxMult);
2fbf38ac	228	fHistNeutralMult->GetXaxis()->SetTitle("N");
	229	fHistNeutralMult->GetYaxis()->SetTitle("Multiplicity");
	230
	231	histname = "fHistPhivsPtPos" + fHistogramNameSuffix;
13b0d3c1	232	fHistPhivsPtPos = new TH2F(histname.Data(), "Phi vs pT of positively charged tracks hitting the calorimeter", 200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
2fbf38ac	233
2fbf38ac	234	histname = "fHistPhivsPtNeg" + fHistogramNameSuffix;
13b0d3c1	235	fHistPhivsPtNeg = new TH2F(histname.Data(), "Phi vs pT of negatively charged tracks hitting the calorimeter", 200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
2fbf38ac	236
2fbf38ac	237	histname = "fHistBaryonEt" + fHistogramNameSuffix;
13b0d3c1	238	fHistBaryonEt = new TH1F(histname.Data(), "E_{T} for baryons", nbinsEt, minEt, maxEt);
2fbf38ac	239
2fbf38ac	240	histname = "fHistAntiBaryonEt" + fHistogramNameSuffix;
13b0d3c1	241	fHistAntiBaryonEt = new TH1F(histname.Data(), "E_{T} for anti baryons", nbinsEt, minEt, maxEt);
2fbf38ac	242
2fbf38ac	243	histname = "fHistMesonEt" + fHistogramNameSuffix;
13b0d3c1	244	fHistMesonEt = new TH1F(histname.Data(), "E_{T} for mesons", nbinsEt, minEt, maxEt);
2fbf38ac	245
2fbf38ac	246	histname = "fHistBaryonEtAcc" + fHistogramNameSuffix;
13b0d3c1	247	fHistBaryonEtAcc = new TH1F(histname.Data(), "E_{T} for baryons in calorimeter acceptance", nbinsEt, minEt, maxEt);
2fbf38ac	248
2fbf38ac	249	histname = "fHistAntiBaryonEtAcc" + fHistogramNameSuffix;
13b0d3c1	250	fHistAntiBaryonEtAcc = new TH1F(histname.Data(), "E_{T} for anti baryons in calorimeter acceptance", nbinsEt, minEt, maxEt);
2fbf38ac	251
2fbf38ac	252	histname = "fHistMesonEtAcc" + fHistogramNameSuffix;
13b0d3c1	253	fHistMesonEtAcc = new TH1F(histname.Data(), "E_{T} for mesons in calorimeter acceptance", nbinsEt, minEt, maxEt);
2fbf38ac	254
87efb15c	255	histname = "fHistProtonEt" + fHistogramNameSuffix;
	256	fHistProtonEt = new TH1F(histname.Data(), "E_{T} for (anti-)protons", nbinsEt, minEt, maxEt);
	257
	258	histname = "fHistKaonEt" + fHistogramNameSuffix;
	259	fHistChargedKaonEt = new TH1F(histname.Data(), "E_{T} for charged kaons", nbinsEt, minEt, maxEt);
	260
	261	histname = "fHistMuonEt" + fHistogramNameSuffix;
	262	fHistMuonEt = new TH1F(histname.Data(), "E_{T} for muons", nbinsEt, minEt, maxEt);
	263
	264	histname = "fHistElectronEt" + fHistogramNameSuffix;
	265	fHistElectronEt = new TH1F(histname.Data(), "E_{T} for electrons/positrons", nbinsEt, minEt, maxEt);
	266
	267	histname = "fHistProtonEtAcc" + fHistogramNameSuffix;
	268	fHistProtonEtAcc = new TH1F(histname.Data(), "E_{T} for (anti-)protons in calorimeter acceptance", nbinsEt, minEt, maxEt);
	269
	270	histname = "fHistKaonEtAcc" + fHistogramNameSuffix;
	271	fHistChargedKaonEtAcc = new TH1F(histname.Data(), "E_{T} for charged kaons in calorimeter acceptance", nbinsEt, minEt, maxEt);
	272
	273	histname = "fHistMuonEtAcc" + fHistogramNameSuffix;
	274	fHistMuonEtAcc = new TH1F(histname.Data(), "E_{T} for muons in calorimeter acceptance", nbinsEt, minEt, maxEt);
	275
	276	histname = "fHistElectronEtAcc" + fHistogramNameSuffix;
	277	fHistElectronEtAcc = new TH1F(histname.Data(), "E_{T} for electrons/positrons in calorimeter acceptance", nbinsEt, minEt, maxEt);
	278
	279	histname = "fHistEtRecvsEtMC" + fHistogramNameSuffix;
	280	fHistEtRecvsEtMC = new TH2F(histname.Data(), "Reconstructed E_{t} vs MC E_{t}", nbinsEt, minEt, maxEt, nbinsEt, minEt, maxEt);
	281
cf6522d1	282	//
2fbf38ac	283	histname = "fHistTMDeltaR" + fHistogramNameSuffix;
	284	fHistTMDeltaR = new TH1F(histname.Data(), "#Delta R for calorimeter clusters", 200, 0, 50);
	285
	286	}
	287
ce546038	288	void AliAnalysisEt::CreateTree()
	289	{ // create tree..
	290	TString treename = "fTree" + fHistogramNameSuffix;
	291	fTree = new TTree(treename, treename);
	292	fTree->Branch("fTotEt",&fTotEt,"fTotEt/D");
	293	fTree->Branch("fTotEtAcc",&fTotEtAcc,"fTotEtAcc/D");
	294	fTree->Branch("fTotNeutralEt",&fTotNeutralEt,"fTotNeutralEt/D");
	295	fTree->Branch("fTotNeutralEtAcc",&fTotNeutralEtAcc,"fTotNeutralEtAcc/D");
	296	fTree->Branch("fTotChargedEt",&fTotChargedEt,"fTotChargedEt/D");
	297	fTree->Branch("fTotChargedEtAcc",&fTotChargedEtAcc,"fTotChargedEtAcc/D");
	298	fTree->Branch("fMultiplicity",&fMultiplicity,"fMultiplicity/I");
	299	fTree->Branch("fChargedMultiplicity",&fChargedMultiplicity,"fChargedMultiplicity/I");
	300	fTree->Branch("fNeutralMultiplicity",&fNeutralMultiplicity,"fNeutralMultiplicity/I");
	301	fTree->Branch("fBaryonEt",&fBaryonEt,"fBaryonEt/D");
	302	fTree->Branch("fAntiBaryonEt",&fAntiBaryonEt,"fAntiBaryonEt/D");
	303	fTree->Branch("fMesonEt",&fMesonEt,"fMesonEt/D");
	304	fTree->Branch("fBaryonEtAcc",&fBaryonEtAcc,"fBaryonEtAcc/D");
	305	fTree->Branch("fAntiBaryonEtAcc",&fAntiBaryonEtAcc,"fAntiBaryonEtAcc/D");
	306	fTree->Branch("fMesonEtAcc",&fMesonEtAcc,"fMesonEtAcc/D");
	307	fTree->Branch("fProtonEt",&fProtonEt,"fProtonEt/D");
	308	fTree->Branch("fChargedKaonEt",&fChargedKaonEt,"fChargedKaonEt/D");
	309	fTree->Branch("fMuonEt",&fMuonEt,"fMuonEt/D");
	310	fTree->Branch("fElectronEt",&fElectronEt,"fElectronEt/D");
	311	fTree->Branch("fProtonEtAcc",&fProtonEtAcc,"fProtonEtAcc/D");
	312	fTree->Branch("fChargedKaonEtAcc",&fChargedKaonEtAcc,"fChargedKaonEtAcc/D");
	313	fTree->Branch("fMuonEtAcc",&fMuonEtAcc,"fMuonEtAcc/D");
	314	fTree->Branch("fElectronEtAcc",&fElectronEtAcc,"fElectronEtAcc/D");
	315
	316	return;
	317	}
	318
2fbf38ac	319	void AliAnalysisEt::FillHistograms()
cf6522d1	320	{ // fill histograms..
2fbf38ac	321	fHistEt->Fill(fTotEt);
	322	fHistChargedEt->Fill(fTotChargedEt);
	323	fHistNeutralEt->Fill(fTotNeutralEt);
	324
	325	fHistEtAcc->Fill(fTotEtAcc);
	326	fHistChargedEtAcc->Fill(fTotChargedEtAcc);
	327	fHistNeutralEtAcc->Fill(fTotNeutralEtAcc);
	328
	329	fHistMult->Fill(fMultiplicity);
	330	fHistChargedMult->Fill(fChargedMultiplicity);
	331	fHistNeutralMult->Fill(fNeutralMultiplicity);
	332
ce546038	333	fHistBaryonEt->Fill(fBaryonEt);
	334	fHistAntiBaryonEt->Fill(fAntiBaryonEt);
	335	fHistMesonEt->Fill(fMesonEt);
2fbf38ac	336
ce546038	337	fHistBaryonEtAcc->Fill(fBaryonEtAcc);
	338	fHistAntiBaryonEtAcc->Fill(fAntiBaryonEtAcc);
	339	fHistMesonEtAcc->Fill(fMesonEtAcc);
	340
87efb15c	341	fHistProtonEt->Fill(fProtonEt);
	342	fHistChargedKaonEt->Fill(fChargedKaonEt);
	343	fHistMuonEt->Fill(fMuonEt);
	344	fHistElectronEt->Fill(fElectronEt);
	345
	346	fHistProtonEtAcc->Fill(fProtonEtAcc);
	347	fHistChargedKaonEtAcc->Fill(fChargedKaonEtAcc);
	348	fHistMuonEtAcc->Fill(fMuonEtAcc);
	349	fHistElectronEtAcc->Fill(fElectronEtAcc);
ce546038	350
	351	if (fCuts) {
	352	if (fCuts->GetHistMakeTree()) {
	353	fTree->Fill();
	354	}
	355	}
	356
2fbf38ac	357	}
0065531a	358
641e1e0c	359	Int_t AliAnalysisEt::AnalyseEvent(AliVEvent *event)
0065531a	360	{ //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.
	361	cout << "This event has " << event->GetNumberOfTracks() << " tracks" << endl;
	362	ResetEventValues();
641e1e0c	363	return 0;
641e1e0c	364	}
2fbf38ac	365
2fbf38ac	366	void AliAnalysisEt::ResetEventValues()
cf6522d1	367	{ // clear
0065531a	368	fTotEt = 0;
	369	fTotEtAcc = 0;
	370	fTotNeutralEt = 0;
	371	fTotNeutralEtAcc = 0;
	372	fTotChargedEt = 0;
	373	fTotChargedEtAcc = 0;
	374	fMultiplicity = 0;
	375	fChargedMultiplicity = 0;
	376	fNeutralMultiplicity = 0;
87efb15c	377	fBaryonEt = 0;
	378	fAntiBaryonEt = 0;
	379	fMesonEt = 0;
	380	fBaryonEtAcc = 0;
	381	fAntiBaryonEtAcc = 0;
	382	fMesonEtAcc = 0;
	383	fProtonEt = 0;
	384	fChargedKaonEt = 0;
	385	fMuonEt = 0;
	386	fElectronEt = 0;
	387	fProtonEtAcc = 0;
	388	fChargedKaonEtAcc = 0;
	389	fMuonEtAcc = 0;
ce546038	390	fElectronEtAcc = 0;
87efb15c	391
0065531a	392	if (!fCuts \|\| !fPdgDB \|\| fPiPlusCode==0) { // some Init's needed
	393	cout << __FILE__ << ":" << __LINE__ << " : Init " << endl;
	394	if (!fCuts) {
	395	cout << " setting up Cuts " << endl;
	396	fCuts = new AliAnalysisEtCuts();
	397	}
	398	if(!fPdgDB) {
	399	cout << " setting up PdgDB " << endl;
	400	fPdgDB = new TDatabasePDG();
	401	}
	402
	403	if (fPiPlusCode==0) {
99a6613d	404	SetParticleCodes();
99a6613d	405	}
0065531a	406	}
0065531a	407	return;
2fbf38ac	408	}
99a6613d	409
0065531a	410
99a6613d	411	void AliAnalysisEt::SetParticleCodes()
cf6522d1	412	{ // set PDG info
8985c506	413	fPionMass = fPdgDB->GetParticle("pi+")->Mass();
	414	fPiPlusCode = fPdgDB->GetParticle("pi+")->PdgCode();
	415	fPiMinusCode = fPdgDB->GetParticle("pi-")->PdgCode();
	416	fKPlusCode = fPdgDB->GetParticle("K+")->PdgCode();
	417	fKMinusCode = fPdgDB->GetParticle("K-")->PdgCode();
	418	fProtonCode = fPdgDB->GetParticle("proton")->PdgCode();
	419	fAntiProtonCode = fPdgDB->GetParticle("antiproton")->PdgCode();
	420	fLambdaCode = fPdgDB->GetParticle("Lambda0")->PdgCode();
	421	fAntiLambdaCode = fPdgDB->GetParticle("Lambda0_bar")->PdgCode();
	422	fK0SCode = fPdgDB->GetParticle("K_S0")->PdgCode();
	423	fOmegaCode = fPdgDB->GetParticle("Omega-")->PdgCode();
	424	fAntiOmegaCode = fPdgDB->GetParticle("Omega+")->PdgCode();
	425	fXi0Code = fPdgDB->GetParticle("Xi0")->PdgCode();
	426	fAntiXi0Code = fPdgDB->GetParticle("Xi0_bar")->PdgCode();
	427	fXiCode = fPdgDB->GetParticle("Xi-")->PdgCode();
	428	fAntiXiCode = fPdgDB->GetParticle("Xi-_bar")->PdgCode();
	429	fSigmaCode = fPdgDB->GetParticle("Sigma-")->PdgCode();
	430	fAntiSigmaCode = fPdgDB->GetParticle("Sigma+")->PdgCode();
	431	fK0LCode = fPdgDB->GetParticle("K_L0")->PdgCode();
	432	fNeutronCode = fPdgDB->GetParticle("neutron")->PdgCode();
	433	fAntiNeutronCode = fPdgDB->GetParticle("antineutron")->PdgCode();
	434	fEPlusCode = fPdgDB->GetParticle("e+")->PdgCode();
	435	fEMinusCode = fPdgDB->GetParticle("e-")->PdgCode();
	436	cout << "Resetting Codes: Pion " << fPiPlusCode
	437	<< "," << fPiMinusCode
	438	<< " Kaon " << fKPlusCode
	439	<< "," << fKMinusCode << endl;
99a6613d	440	}
99a6613d	441