[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliAnalysisTaskScale.cxx

// $Id$
//
// Scale task.
//
// Author: R.Reed, M.Connors

#include "AliAnalysisTaskScale.h"

#include <TClonesArray.h>
#include <TF1.h>
#include <TH2F.h>
#include <TLorentzVector.h>
#include <TMath.h>

#include "AliEMCALGeometry.h"
#include "AliLog.h"
#include "AliVCluster.h"
#include "AliVParticle.h"
#include "AliParticleContainer.h"

ClassImp(AliAnalysisTaskScale)

//________________________________________________________________________
AliAnalysisTaskScale::AliAnalysisTaskScale() : 
  AliAnalysisTaskEmcal("AliAnalysisTaskScale", kTRUE), 
  fScaleFunction(0),
  fEmcalArea(1),
  fTpcArea(1),
  fHistPtTPCvsCent(0), 
  fHistPtEMCALvsCent(0), 
  fHistEtvsCent(0),  
  fHistScalevsCent(0),  
  fHistDeltaScalevsCent(0), 
  fHistScaleEmcalvsCent(0),      
  fHistScale2EmcalvsCent(0),     
  fHistChScalevsCent(0),          
  fHistChScale2EmcalvsCent(0),   
  fHistPtTPCvsNtrack(0), 
  fHistPtEMCALvsNtrack(0), 
  fHistEtvsNtrack(0),  
  fHistScalevsNtrack(0),  
  fHistDeltaScalevsNtrack(0),
  fHistScaleEmcalvsNtrack(0),      
  fHistScale2EmcalvsNtrack(0),     
  fHistChScalevsNtrack(0),          
  fHistChScale2EmcalvsNtrack(0),   
  fHistTrackPtvsCent(0), 
  fHistClusterPtvsCent(0),
  fHistTrackEtaPhi(0), 
  fHistClusterEtaPhi(0),
  fHistScalevsScale2Emcal(0),      
  fHistScalevsScaleEmcal(0),       
  fHistScaleEmcalvsScale2Emcal(0)
{
  // Default constructor.

  SetMakeGeneralHistograms(kTRUE);
}

//________________________________________________________________________
AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) :
  AliAnalysisTaskEmcal(name, kTRUE), 
  fScaleFunction(0),
  fEmcalArea(1),
  fTpcArea(1),
  fHistPtTPCvsCent(0), 
  fHistPtEMCALvsCent(0), 
  fHistEtvsCent(0),  
  fHistScalevsCent(0),  
  fHistDeltaScalevsCent(0), 
  fHistScaleEmcalvsCent(0),      
  fHistScale2EmcalvsCent(0),     
  fHistChScalevsCent(0),          
  fHistChScale2EmcalvsCent(0),   
  fHistPtTPCvsNtrack(0), 
  fHistPtEMCALvsNtrack(0), 
  fHistEtvsNtrack(0),  
  fHistScalevsNtrack(0),  
  fHistDeltaScalevsNtrack(0),
  fHistScaleEmcalvsNtrack(0),      
  fHistScale2EmcalvsNtrack(0),     
  fHistChScalevsNtrack(0),          
  fHistChScale2EmcalvsNtrack(0),   
  fHistTrackPtvsCent(0), 
  fHistClusterPtvsCent(0),
  fHistTrackEtaPhi(0), 
  fHistClusterEtaPhi(0),
  fHistScalevsScale2Emcal(0),      
  fHistScalevsScaleEmcal(0),       
  fHistScaleEmcalvsScale2Emcal(0)
{
  // Constructor.

  SetMakeGeneralHistograms(kTRUE);
}

//________________________________________________________________________
void AliAnalysisTaskScale::UserCreateOutputObjects()
{
  // Create my user objects.

  AliAnalysisTaskEmcal::UserCreateOutputObjects();

  fHistPtTPCvsCent             = new TH2F("PtTPCvsCent","rho vs cent",            101, -1, 100,   500,   0, 1000);
  fHistPtEMCALvsCent           = new TH2F("PtEMCALvsCent","rho vs cent",          101, -1, 100,   500,   0, 1000);
  fHistEtvsCent                = new TH2F("EtvsCent","rho vs cent",               101, -1, 100,   500,   0, 1000);
  fHistScalevsCent             = new TH2F("ScalevsCent","rho vs cent",            101, -1, 100,   500,   0, 5);
  fHistDeltaScalevsCent        = new TH2F("DeltaScalevsCent","rho vs cent",       101, -1, 100,   500,  -2.5, 2.5);
  fHistScaleEmcalvsCent        = new TH2F("ScaleEmcalvsCent","",                  101, -1, 100,   500,   0, 5);
  fHistScale2EmcalvsCent       = new TH2F("Scale2EmcalvsCent","",                 101, -1, 100,   500,   0, 5);
  fHistChScalevsCent           = new TH2F("ChScalevsCent","",                     101, -1, 100,   500,   0, 5);
  fHistChScale2EmcalvsCent     = new TH2F("ChScale2EmcalvsCent","",               101, -1, 100,   500,   0, 5);
  fHistPtTPCvsNtrack           = new TH2F("PtTPCvsNtrack","rho vs cent",          500,  0, 2500,  500,   0, 1000);
  fHistPtEMCALvsNtrack         = new TH2F("PtEMCALvsNtrack","rho vs cent",        500,  0, 2500,  500,   0, 1000);
  fHistEtvsNtrack              = new TH2F("EtvsNtrack","rho vs cent",             500,  0, 2500,  500,   0, 1000);
  fHistScalevsNtrack           = new TH2F("ScalevsNtrack","rho vs cent",          500,  0, 2500,  500,   0, 5);
  fHistDeltaScalevsNtrack      = new TH2F("DeltaScalevsNtrack","rho vs cent",     500,  0, 2500,  500,  -2.5, 2.5);
  fHistScaleEmcalvsNtrack      = new TH2F("ScaleEmcalvsNtrack","",                500,  0, 2500,  500,   0, 5);
  fHistScale2EmcalvsNtrack     = new TH2F("Scale2EmcalvsNtrack","",               500,  0, 2500,  500,   0, 5);
  fHistChScalevsNtrack         = new TH2F("ChScalevsNtrack","",                   500,  0, 2500,  500,   0, 5);
  fHistChScale2EmcalvsNtrack   = new TH2F("ChScale2EmcalvsNtrack","",             500,  0, 2500,  500,   0, 5);
  fHistTrackPtvsCent           = new TH2F("TrackPtvsCent","Track pt vs cent",     101, -1, 100,   500,   0, 100);
  fHistClusterPtvsCent         = new TH2F("ClusterPtvsCent","Cluster pt vs cent", 101, -1, 100,   500,   0, 100);
  fHistTrackEtaPhi             = new TH2F("TrackEtaPhi","Track eta phi",          100, -1.0, 1.0, 101,   0, 2.02*TMath::Pi());
  fHistClusterEtaPhi           = new TH2F("ClusterEtaPhi","Cluster eta phi",      100, -1.0, 1.0, 101,   0, 2.02*TMath::Pi());
  fHistScalevsScale2Emcal      = new TH2F("ScalevsScale2Emcal","",                500,  0, 5,     500,   0, 5);
  fHistScalevsScaleEmcal       = new TH2F("ScalevsScaleEmcal","",                 500,  0, 5,     500,   0, 5);
  fHistScaleEmcalvsScale2Emcal = new TH2F("ScaleEmcalvsScale2Emcal","",           500,  0, 5,     500,   0, 5);

  fOutput->Add(fHistPtTPCvsCent);
  fOutput->Add(fHistPtEMCALvsCent);
  fOutput->Add(fHistEtvsCent);
  fOutput->Add(fHistScalevsCent);
  fOutput->Add(fHistDeltaScalevsCent);
  fOutput->Add(fHistScaleEmcalvsCent);      
  fOutput->Add(fHistScale2EmcalvsCent);     
  fOutput->Add(fHistChScalevsCent);    
  fOutput->Add(fHistChScale2EmcalvsCent);   
  fOutput->Add(fHistPtTPCvsNtrack);
  fOutput->Add(fHistPtEMCALvsNtrack);
  fOutput->Add(fHistEtvsNtrack);
  fOutput->Add(fHistScalevsNtrack);
  fOutput->Add(fHistDeltaScalevsNtrack);
  fOutput->Add(fHistScaleEmcalvsNtrack);      
  fOutput->Add(fHistScale2EmcalvsNtrack);     
  fOutput->Add(fHistChScalevsNtrack);    
  fOutput->Add(fHistChScale2EmcalvsNtrack);   
  fOutput->Add(fHistTrackPtvsCent);
  fOutput->Add(fHistClusterPtvsCent);
  fOutput->Add(fHistTrackEtaPhi);
  fOutput->Add(fHistClusterEtaPhi);
  fOutput->Add(fHistScalevsScale2Emcal);      
  fOutput->Add(fHistScalevsScaleEmcal);       
  fOutput->Add(fHistScaleEmcalvsScale2Emcal);

  PostData(1, fOutput);
}

//________________________________________________________________________
Double_t AliAnalysisTaskScale::GetScaleFactor(Double_t cent)
{
  // Get scale function.

  Double_t scale = -1;
  if (fScaleFunction)
    scale = fScaleFunction->Eval(cent);
  return scale;
}

//________________________________________________________________________
Bool_t AliAnalysisTaskScale::FillHistograms() 
{
  // Execute on each event.

  const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
  const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
  const Double_t EmcalWidth = (EmcalMaxPhi-EmcalMinPhi)/2.0;

  Double_t ptTPC   = 0;
  Double_t ptEMCAL = 0;
  Double_t ptEMCAL2 = 0;

  const Int_t Ntracks = fTracks->GetEntries();
  for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
    AliVParticle *track = static_cast<AliVParticle*>(fTracks->At(iTracks));

    if (!track)
      continue;

    if (!AcceptTrack(track))
      continue;

    if (TMath::Abs(track->Eta()) > 0.7)   // only accept tracks in the EMCal eta range
      continue;

    fHistTrackPtvsCent->Fill(fCent,track->Pt());
    fHistTrackEtaPhi->Fill(track->Eta(),track->Phi());
    ptTPC += track->Pt();
    if ((track->Phi() > (EmcalMaxPhi+EmcalWidth)) || (track->Phi() < (EmcalMinPhi-EmcalWidth))) continue;
    ptEMCAL2 += track->Pt();
    if ((track->Phi() > EmcalMaxPhi) || (track->Phi() < EmcalMinPhi)) continue;
    ptEMCAL += track->Pt();
  }

  if (ptTPC == 0) 
    return kFALSE;
  
  Double_t Et = 0;
  const Int_t Nclus = fCaloClusters->GetEntries();
  for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
    AliVCluster *c = static_cast<AliVCluster*>(fCaloClusters->At(iClus));
    if (!c)
      continue;

    if (!AcceptCluster(c))
      continue;

    TLorentzVector nPart;
    c->GetMomentum(nPart, fVertex);

    fHistClusterPtvsCent->Fill(fCent, nPart.Pt());
    fHistClusterEtaPhi->Fill(nPart.Eta(), nPart.Phi());

    Et += nPart.Pt();
  }
 
  Double_t scalecalc         = -1;
  if (ptEMCAL > 0 && Et > 0 && ptTPC > 0)
    scalecalc         =  ((Et + ptEMCAL) / fEmcalArea) * (fTpcArea / ptTPC);
  const Double_t scale             = GetScaleFactor(fCent);
  Double_t scalecalcemcal          = -1;
  if (ptEMCAL > 0)
    scalecalcemcal                 = (Et+ptEMCAL)/ptEMCAL;
  Double_t scalecalcemcal2         = -1;
  Double_t Chscalecalcemcal2       = -1;
  if (ptEMCAL2 > 0){
    scalecalcemcal2                = 2*(Et+ptEMCAL)/ptEMCAL2;
    Chscalecalcemcal2              = 2*ptEMCAL/ptEMCAL2;}
  const Double_t Chscalecalcemcal  = ((ptEMCAL) / fEmcalArea) * (fTpcArea / ptTPC);

  fHistScaleEmcalvsCent->Fill(fCent,scalecalcemcal);      
  fHistScale2EmcalvsCent->Fill(fCent,scalecalcemcal2);     
  fHistChScalevsCent->Fill(fCent,Chscalecalcemcal);    
  fHistChScale2EmcalvsCent->Fill(fCent,Chscalecalcemcal2);   
  fHistScaleEmcalvsNtrack->Fill(Ntracks,scalecalcemcal);      
  fHistScale2EmcalvsNtrack->Fill(Ntracks,scalecalcemcal2);     
  fHistChScalevsNtrack->Fill(Ntracks,Chscalecalcemcal);    
  fHistChScale2EmcalvsNtrack->Fill(Ntracks,Chscalecalcemcal2);   
  fHistPtTPCvsCent->Fill(fCent, ptTPC);
  fHistPtEMCALvsCent->Fill(fCent, ptEMCAL);
  fHistEtvsCent->Fill(fCent, Et);
  fHistScalevsCent->Fill(fCent, scalecalc);
  fHistDeltaScalevsCent->Fill(fCent, scalecalc - scale);
  fHistPtTPCvsNtrack->Fill(Ntracks, ptTPC);
  fHistPtEMCALvsNtrack->Fill(Ntracks, ptEMCAL);
  fHistEtvsNtrack->Fill(Ntracks, Et);
  fHistScalevsNtrack->Fill(Ntracks, scalecalc);
  fHistDeltaScalevsNtrack->Fill(Ntracks, scalecalc - scale);
  fHistScalevsScale2Emcal->Fill(scalecalc,scalecalcemcal2);      
  fHistScalevsScaleEmcal->Fill(scalecalc,scalecalcemcal); 
  fHistScaleEmcalvsScale2Emcal->Fill(scalecalcemcal,scalecalcemcal2);

  return kTRUE;
}


//________________________________________________________________________
void AliAnalysisTaskScale::ExecOnce() 
{
  AliAnalysisTaskEmcal::ExecOnce();

  const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
  const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
  const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
  const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();

  fEmcalArea  = (EmcalMaxPhi - EmcalMinPhi) * (EmcalMinEta - EmcalMaxEta);

  AliParticleContainer *partCont = GetParticleContainer(0);
  if (!partCont) {
    AliError(Form("%s: No particle container found! Assuming tpc area = 1...",GetName()));
    fTpcArea = 1;
    return;
  }

  Float_t TpcMaxPhi = partCont->GetParticlePhiMax();
  Float_t TpcMinPhi = partCont->GetParticlePhiMin();
  
  if (TpcMaxPhi > TMath::Pi()*2) TpcMaxPhi = TMath::Pi()*2;
  if (TpcMinPhi < 0) TpcMinPhi = 0;

  fTpcArea = (TpcMaxPhi - TpcMinPhi) * (EmcalMinEta - EmcalMaxEta);
}
Commit	Line	Data
0b777a09	1	// $Id$
4a6b693d	2	//
	3	// Scale task.
	4	//
68263112	5	// Author: R.Reed, M.Connors
0b777a09	6
c25dc25a	7	#include "AliAnalysisTaskScale.h"
c25dc25a	8
0b777a09	9	#include <TClonesArray.h>
c25dc25a	10	#include <TF1.h>
8154249f	11	#include <TH2F.h>
0b777a09	12	#include <TLorentzVector.h>
73e2fd59	13	#include <TMath.h>
b1dfb3f4	14
40e2fd10	15	#include "AliEMCALGeometry.h"
40e2fd10	16	#include "AliLog.h"
c25dc25a	17	#include "AliVCluster.h"
7cd832c7	18	#include "AliVParticle.h"
7cd832c7	19	#include "AliParticleContainer.h"
0b777a09	20
	21	ClassImp(AliAnalysisTaskScale)
	22
c25dc25a	23	//________________________________________________________________________
c25dc25a	24	AliAnalysisTaskScale::AliAnalysisTaskScale() :
9239b066	25	AliAnalysisTaskEmcal("AliAnalysisTaskScale", kTRUE),
c25dc25a	26	fScaleFunction(0),
7cd832c7	27	fEmcalArea(1),
7cd832c7	28	fTpcArea(1),
c25dc25a	29	fHistPtTPCvsCent(0),
	30	fHistPtEMCALvsCent(0),
	31	fHistEtvsCent(0),
	32	fHistScalevsCent(0),
	33	fHistDeltaScalevsCent(0),
a828ea11	34	fHistScaleEmcalvsCent(0),
	35	fHistScale2EmcalvsCent(0),
	36	fHistChScalevsCent(0),
	37	fHistChScale2EmcalvsCent(0),
c25dc25a	38	fHistPtTPCvsNtrack(0),
	39	fHistPtEMCALvsNtrack(0),
	40	fHistEtvsNtrack(0),
	41	fHistScalevsNtrack(0),
	42	fHistDeltaScalevsNtrack(0),
a828ea11	43	fHistScaleEmcalvsNtrack(0),
	44	fHistScale2EmcalvsNtrack(0),
	45	fHistChScalevsNtrack(0),
	46	fHistChScale2EmcalvsNtrack(0),
c25dc25a	47	fHistTrackPtvsCent(0),
	48	fHistClusterPtvsCent(0),
	49	fHistTrackEtaPhi(0),
a828ea11	50	fHistClusterEtaPhi(0),
	51	fHistScalevsScale2Emcal(0),
	52	fHistScalevsScaleEmcal(0),
	53	fHistScaleEmcalvsScale2Emcal(0)
c25dc25a	54	{
c25dc25a	55	// Default constructor.
73e2fd59	56
73e2fd59	57	SetMakeGeneralHistograms(kTRUE);
c25dc25a	58	}
c25dc25a	59
0b777a09	60	//________________________________________________________________________
b1dfb3f4	61	AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) :
9239b066	62	AliAnalysisTaskEmcal(name, kTRUE),
b1dfb3f4	63	fScaleFunction(0),
7cd832c7	64	fEmcalArea(1),
7cd832c7	65	fTpcArea(1),
b1dfb3f4	66	fHistPtTPCvsCent(0),
	67	fHistPtEMCALvsCent(0),
	68	fHistEtvsCent(0),
	69	fHistScalevsCent(0),
	70	fHistDeltaScalevsCent(0),
a828ea11	71	fHistScaleEmcalvsCent(0),
	72	fHistScale2EmcalvsCent(0),
	73	fHistChScalevsCent(0),
	74	fHistChScale2EmcalvsCent(0),
b1dfb3f4	75	fHistPtTPCvsNtrack(0),
	76	fHistPtEMCALvsNtrack(0),
	77	fHistEtvsNtrack(0),
	78	fHistScalevsNtrack(0),
a5e84976	79	fHistDeltaScalevsNtrack(0),
a828ea11	80	fHistScaleEmcalvsNtrack(0),
	81	fHistScale2EmcalvsNtrack(0),
	82	fHistChScalevsNtrack(0),
	83	fHistChScale2EmcalvsNtrack(0),
a5e84976	84	fHistTrackPtvsCent(0),
	85	fHistClusterPtvsCent(0),
	86	fHistTrackEtaPhi(0),
a828ea11	87	fHistClusterEtaPhi(0),
	88	fHistScalevsScale2Emcal(0),
	89	fHistScalevsScaleEmcal(0),
	90	fHistScaleEmcalvsScale2Emcal(0)
0b777a09	91	{
4a6b693d	92	// Constructor.
0b777a09	93
73e2fd59	94	SetMakeGeneralHistograms(kTRUE);
0b777a09	95	}
	96
	97	//________________________________________________________________________
	98	void AliAnalysisTaskScale::UserCreateOutputObjects()
	99	{
4a6b693d	100	// Create my user objects.
4a6b693d	101
9239b066	102	AliAnalysisTaskEmcal::UserCreateOutputObjects();
0b777a09	103
a828ea11	104	fHistPtTPCvsCent = new TH2F("PtTPCvsCent","rho vs cent", 101, -1, 100, 500, 0, 1000);
	105	fHistPtEMCALvsCent = new TH2F("PtEMCALvsCent","rho vs cent", 101, -1, 100, 500, 0, 1000);
	106	fHistEtvsCent = new TH2F("EtvsCent","rho vs cent", 101, -1, 100, 500, 0, 1000);
	107	fHistScalevsCent = new TH2F("ScalevsCent","rho vs cent", 101, -1, 100, 500, 0, 5);
	108	fHistDeltaScalevsCent = new TH2F("DeltaScalevsCent","rho vs cent", 101, -1, 100, 500, -2.5, 2.5);
	109	fHistScaleEmcalvsCent = new TH2F("ScaleEmcalvsCent","", 101, -1, 100, 500, 0, 5);
	110	fHistScale2EmcalvsCent = new TH2F("Scale2EmcalvsCent","", 101, -1, 100, 500, 0, 5);
	111	fHistChScalevsCent = new TH2F("ChScalevsCent","", 101, -1, 100, 500, 0, 5);
	112	fHistChScale2EmcalvsCent = new TH2F("ChScale2EmcalvsCent","", 101, -1, 100, 500, 0, 5);
	113	fHistPtTPCvsNtrack = new TH2F("PtTPCvsNtrack","rho vs cent", 500, 0, 2500, 500, 0, 1000);
	114	fHistPtEMCALvsNtrack = new TH2F("PtEMCALvsNtrack","rho vs cent", 500, 0, 2500, 500, 0, 1000);
	115	fHistEtvsNtrack = new TH2F("EtvsNtrack","rho vs cent", 500, 0, 2500, 500, 0, 1000);
	116	fHistScalevsNtrack = new TH2F("ScalevsNtrack","rho vs cent", 500, 0, 2500, 500, 0, 5);
	117	fHistDeltaScalevsNtrack = new TH2F("DeltaScalevsNtrack","rho vs cent", 500, 0, 2500, 500, -2.5, 2.5);
	118	fHistScaleEmcalvsNtrack = new TH2F("ScaleEmcalvsNtrack","", 500, 0, 2500, 500, 0, 5);
	119	fHistScale2EmcalvsNtrack = new TH2F("Scale2EmcalvsNtrack","", 500, 0, 2500, 500, 0, 5);
	120	fHistChScalevsNtrack = new TH2F("ChScalevsNtrack","", 500, 0, 2500, 500, 0, 5);
	121	fHistChScale2EmcalvsNtrack = new TH2F("ChScale2EmcalvsNtrack","", 500, 0, 2500, 500, 0, 5);
	122	fHistTrackPtvsCent = new TH2F("TrackPtvsCent","Track pt vs cent", 101, -1, 100, 500, 0, 100);
	123	fHistClusterPtvsCent = new TH2F("ClusterPtvsCent","Cluster pt vs cent", 101, -1, 100, 500, 0, 100);
	124	fHistTrackEtaPhi = new TH2F("TrackEtaPhi","Track eta phi", 100, -1.0, 1.0, 101, 0, 2.02*TMath::Pi());
	125	fHistClusterEtaPhi = new TH2F("ClusterEtaPhi","Cluster eta phi", 100, -1.0, 1.0, 101, 0, 2.02*TMath::Pi());
	126	fHistScalevsScale2Emcal = new TH2F("ScalevsScale2Emcal","", 500, 0, 5, 500, 0, 5);
	127	fHistScalevsScaleEmcal = new TH2F("ScalevsScaleEmcal","", 500, 0, 5, 500, 0, 5);
	128	fHistScaleEmcalvsScale2Emcal = new TH2F("ScaleEmcalvsScale2Emcal","", 500, 0, 5, 500, 0, 5);
	129
e44e8726	130	fOutput->Add(fHistPtTPCvsCent);
	131	fOutput->Add(fHistPtEMCALvsCent);
	132	fOutput->Add(fHistEtvsCent);
	133	fOutput->Add(fHistScalevsCent);
	134	fOutput->Add(fHistDeltaScalevsCent);
a828ea11	135	fOutput->Add(fHistScaleEmcalvsCent);
	136	fOutput->Add(fHistScale2EmcalvsCent);
	137	fOutput->Add(fHistChScalevsCent);
	138	fOutput->Add(fHistChScale2EmcalvsCent);
e44e8726	139	fOutput->Add(fHistPtTPCvsNtrack);
	140	fOutput->Add(fHistPtEMCALvsNtrack);
	141	fOutput->Add(fHistEtvsNtrack);
	142	fOutput->Add(fHistScalevsNtrack);
	143	fOutput->Add(fHistDeltaScalevsNtrack);
a828ea11	144	fOutput->Add(fHistScaleEmcalvsNtrack);
	145	fOutput->Add(fHistScale2EmcalvsNtrack);
	146	fOutput->Add(fHistChScalevsNtrack);
	147	fOutput->Add(fHistChScale2EmcalvsNtrack);
e44e8726	148	fOutput->Add(fHistTrackPtvsCent);
	149	fOutput->Add(fHistClusterPtvsCent);
	150	fOutput->Add(fHistTrackEtaPhi);
	151	fOutput->Add(fHistClusterEtaPhi);
a828ea11	152	fOutput->Add(fHistScalevsScale2Emcal);
	153	fOutput->Add(fHistScalevsScaleEmcal);
	154	fOutput->Add(fHistScaleEmcalvsScale2Emcal);
e44e8726	155
e44e8726	156	PostData(1, fOutput);
0b777a09	157	}
0b777a09	158
b1dfb3f4	159	//________________________________________________________________________
	160	Double_t AliAnalysisTaskScale::GetScaleFactor(Double_t cent)
	161	{
c25dc25a	162	// Get scale function.
c25dc25a	163
b1dfb3f4	164	Double_t scale = -1;
	165	if (fScaleFunction)
	166	scale = fScaleFunction->Eval(cent);
	167	return scale;
	168	}
	169
e44e8726	170	//________________________________________________________________________
	171	Bool_t AliAnalysisTaskScale::FillHistograms()
	172	{
	173	// Execute on each event.
7cd832c7	174
40e2fd10	175	const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
40e2fd10	176	const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
a828ea11	177	const Double_t EmcalWidth = (EmcalMaxPhi-EmcalMinPhi)/2.0;
0b777a09	178
85211d8e	179	Double_t ptTPC = 0;
85211d8e	180	Double_t ptEMCAL = 0;
a828ea11	181	Double_t ptEMCAL2 = 0;
b1dfb3f4	182
e44e8726	183	const Int_t Ntracks = fTracks->GetEntries();
0b777a09	184	for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
7cd832c7	185	AliVParticle track = static_cast<AliVParticle>(fTracks->At(iTracks));
b1dfb3f4	186
0b777a09	187	if (!track)
0b777a09	188	continue;
b1dfb3f4	189
e44e8726	190	if (!AcceptTrack(track))
	191	continue;
	192
b1dfb3f4	193	if (TMath::Abs(track->Eta()) > 0.7) // only accept tracks in the EMCal eta range
0b777a09	194	continue;
b1dfb3f4	195
e44e8726	196	fHistTrackPtvsCent->Fill(fCent,track->Pt());
a5e84976	197	fHistTrackEtaPhi->Fill(track->Eta(),track->Phi());
b1dfb3f4	198	ptTPC += track->Pt();
7cd832c7	199	if ((track->Phi() > (EmcalMaxPhi+EmcalWidth)) \|\| (track->Phi() < (EmcalMinPhi-EmcalWidth))) continue;
a828ea11	200	ptEMCAL2 += track->Pt();
7cd832c7	201	if ((track->Phi() > EmcalMaxPhi) \|\| (track->Phi() < EmcalMinPhi)) continue;
b1dfb3f4	202	ptEMCAL += track->Pt();
85211d8e	203	}
e44e8726	204
	205	if (ptTPC == 0)
	206	return kFALSE;
0b777a09	207
85211d8e	208	Double_t Et = 0;
e44e8726	209	const Int_t Nclus = fCaloClusters->GetEntries();
4a6b693d	210	for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
e44e8726	211	AliVCluster c = static_cast<AliVCluster>(fCaloClusters->At(iClus));
d03084cd	212	if (!c)
d03084cd	213	continue;
e44e8726	214
e44e8726	215	if (!AcceptCluster(c))
0b777a09	216	continue;
a5e84976	217
e44e8726	218	TLorentzVector nPart;
e44e8726	219	c->GetMomentum(nPart, fVertex);
a5e84976	220
e44e8726	221	fHistClusterPtvsCent->Fill(fCent, nPart.Pt());
e44e8726	222	fHistClusterEtaPhi->Fill(nPart.Eta(), nPart.Phi());
a5e84976	223
85211d8e	224	Et += nPart.Pt();
0b777a09	225	}
a828ea11	226
a38a6344	227	Double_t scalecalc = -1;
a38a6344	228	if (ptEMCAL > 0 && Et > 0 && ptTPC > 0)
7cd832c7	229	scalecalc = ((Et + ptEMCAL) / fEmcalArea) * (fTpcArea / ptTPC);
a828ea11	230	const Double_t scale = GetScaleFactor(fCent);
	231	Double_t scalecalcemcal = -1;
	232	if (ptEMCAL > 0)
	233	scalecalcemcal = (Et+ptEMCAL)/ptEMCAL;
	234	Double_t scalecalcemcal2 = -1;
	235	Double_t Chscalecalcemcal2 = -1;
	236	if (ptEMCAL2 > 0){
	237	scalecalcemcal2 = 2*(Et+ptEMCAL)/ptEMCAL2;
	238	Chscalecalcemcal2 = 2*ptEMCAL/ptEMCAL2;}
7cd832c7	239	const Double_t Chscalecalcemcal = ((ptEMCAL) / fEmcalArea) * (fTpcArea / ptTPC);
a828ea11	240
	241	fHistScaleEmcalvsCent->Fill(fCent,scalecalcemcal);
	242	fHistScale2EmcalvsCent->Fill(fCent,scalecalcemcal2);
	243	fHistChScalevsCent->Fill(fCent,Chscalecalcemcal);
	244	fHistChScale2EmcalvsCent->Fill(fCent,Chscalecalcemcal2);
	245	fHistScaleEmcalvsNtrack->Fill(Ntracks,scalecalcemcal);
	246	fHistScale2EmcalvsNtrack->Fill(Ntracks,scalecalcemcal2);
	247	fHistChScalevsNtrack->Fill(Ntracks,Chscalecalcemcal);
	248	fHistChScale2EmcalvsNtrack->Fill(Ntracks,Chscalecalcemcal2);
e44e8726	249	fHistPtTPCvsCent->Fill(fCent, ptTPC);
	250	fHistPtEMCALvsCent->Fill(fCent, ptEMCAL);
	251	fHistEtvsCent->Fill(fCent, Et);
	252	fHistScalevsCent->Fill(fCent, scalecalc);
	253	fHistDeltaScalevsCent->Fill(fCent, scalecalc - scale);
b1dfb3f4	254	fHistPtTPCvsNtrack->Fill(Ntracks, ptTPC);
	255	fHistPtEMCALvsNtrack->Fill(Ntracks, ptEMCAL);
	256	fHistEtvsNtrack->Fill(Ntracks, Et);
	257	fHistScalevsNtrack->Fill(Ntracks, scalecalc);
	258	fHistDeltaScalevsNtrack->Fill(Ntracks, scalecalc - scale);
a828ea11	259	fHistScalevsScale2Emcal->Fill(scalecalc,scalecalcemcal2);
	260	fHistScalevsScaleEmcal->Fill(scalecalc,scalecalcemcal);
	261	fHistScaleEmcalvsScale2Emcal->Fill(scalecalcemcal,scalecalcemcal2);
0b777a09	262
e44e8726	263	return kTRUE;
4a6b693d	264	}
7cd832c7	265
	266
	267	//________________________________________________________________________
	268	void AliAnalysisTaskScale::ExecOnce()
	269	{
9239b066	270	AliAnalysisTaskEmcal::ExecOnce();
7cd832c7	271
	272	const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
	273	const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
	274	const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
	275	const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
	276
	277	fEmcalArea = (EmcalMaxPhi - EmcalMinPhi) * (EmcalMinEta - EmcalMaxEta);
	278
	279	AliParticleContainer *partCont = GetParticleContainer(0);
	280	if (!partCont) {
	281	AliError(Form("%s: No particle container found! Assuming tpc area = 1...",GetName()));
	282	fTpcArea = 1;
	283	return;
	284	}
	285
	286	Float_t TpcMaxPhi = partCont->GetParticlePhiMax();
	287	Float_t TpcMinPhi = partCont->GetParticlePhiMin();
	288
	289	if (TpcMaxPhi > TMath::Pi()2) TpcMaxPhi = TMath::Pi()2;
	290	if (TpcMinPhi < 0) TpcMinPhi = 0;
	291
	292	fTpcArea = (TpcMaxPhi - TpcMinPhi) * (EmcalMinEta - EmcalMaxEta);
	293	}