[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTaskAODCentralityMaker.cxx

/**************************************************************************
 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/////////////////////////////////////////////////////////////
//
// AliAnalysisTaskSE to make AOD centrality
// Author: Alberica Toia, CERN, Alberica.Toia@cern.ch
//
/////////////////////////////////////////////////////////////

#include <TROOT.h>
#include <TSystem.h>

#include "AliVEvent.h"
#include "AliESDEvent.h"
#include "AliAODEvent.h"
#include "AliAODCentrality.h"
#include "AliAnalysisManager.h"
#include "AliESDInputHandler.h"
#include "AliESDZDC.h"
#include "AliESDFMD.h"
#include "AliESDVZERO.h"
#include "AliMultiplicity.h"
#include "AliAODHandler.h"
#include "AliAODHeader.h"
#include "AliAODEvent.h"
#include "AliAODVertex.h"
#include "AliMCEvent.h"
#include "AliMCEventHandler.h"
#include "AliMCParticle.h"
#include "AliStack.h"
#include "AliHeader.h"
#include "AliGenEventHeader.h"
#include "AliGenHijingEventHeader.h"
#include "AliAnalysisTaskAODCentralityMaker.h"

ClassImp(AliAnalysisTaskAODCentralityMaker)


//________________________________________________________________________
AliAnalysisTaskAODCentralityMaker::AliAnalysisTaskAODCentralityMaker():
AliAnalysisTaskSE(),
fAODCentrality(0),
fDeltaAODFileName("AliAOD.Centrality.root"),
fAODHeader        (0),
fIsMCInput        (0)
{
  // Default constructor
}

//________________________________________________________________________
AliAnalysisTaskAODCentralityMaker::AliAnalysisTaskAODCentralityMaker(const char *name):
AliAnalysisTaskSE(name),
fAODCentrality(0),
fDeltaAODFileName("AliAOD.Centrality.root"),
fAODHeader        (0),
fIsMCInput        (0)
{
  // Standard constructor
}


//________________________________________________________________________
AliAnalysisTaskAODCentralityMaker::~AliAnalysisTaskAODCentralityMaker()
{
  // Destructor
}  

//________________________________________________________________________
void AliAnalysisTaskAODCentralityMaker::Init()
{
  // Initialization
  if(fDebug > 1) printf("AnalysisTaskAODCentralityMaker::Init() \n");
  AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(fDeltaAODFileName.Data());

  return;
}

//________________________________________________________________________
void AliAnalysisTaskAODCentralityMaker::UserCreateOutputObjects()
{
  
  // Create the output container
  //
  if(fDebug > 1) printf("AnalysisTaskAODCentralityMaker::UserCreateOutPutData() \n");
  // Support both the case when the AOD + deltaAOD are produced in an ESD
  // analysis or if the deltaAOD is produced on an analysis on AOD's. (A.G. 27/04/09)
  if(!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) {
    Fatal("UserCreateOutputObjects", "This task needs an AOD handler");
    return;
  }   
  TString filename = fDeltaAODFileName;
  // When running on standard AOD to produce deltas, IsStandardAOD is never set,
  // If AODEvent is NULL, new branches have to be added to the new file(s) (A.G. 15/01/10)
  if(!IsStandardAOD() && AODEvent()) filename = "";

  fAODCentrality = new AliAODCentrality();
  fAODCentrality->SetName("AODCentrality");
  AddAODBranch("AliAODCentrality", &fAODCentrality, filename);
  

  fAODHeader = new AliAODHeader();
  AddAODBranch("AliAODHeader", &fAODHeader, filename);
  return;
}


//________________________________________________________________________
void AliAnalysisTaskAODCentralityMaker::UserExec(Option_t */*option*/)
{
  AliVEvent*   event = InputEvent();
  AliESDEvent* esd   = dynamic_cast<AliESDEvent*>(event);

  Float_t beamEnergy = esd->GetBeamEnergy();
  Int_t   nTracks    = event->GetNumberOfTracks();     
  Int_t   nPmdTracks = esd->GetNumberOfPmdTracks();     
    
  // ***** V0 info
  AliESDVZERO* esdV0 = esd->GetVZEROData();
  Double_t multV0A = esdV0->GetMTotV0A();
  Double_t multV0C = esdV0->GetMTotV0C();
    
  // ***** Trigger selection
  char     trigClass[100];	//  fired trigger classes
  TString triggerClass = esd->GetFiredTriggerClasses();
  sprintf(trigClass,"%s",triggerClass.Data());
  
  // ***** vertex info
  const AliESDVertex *vertex = esd->GetPrimaryVertexSPD();
  Double_t xVertex = vertex->GetX();
  Double_t yVertex = vertex->GetY();
  Double_t zVertex = vertex->GetZ();
  Bool_t vertexer3d;
  
  if(vertex->IsFromVertexer3D()) vertexer3d = kTRUE;
  else vertexer3d = kFALSE;
  Double_t vertex3[3];
  vertex->GetXYZ(vertex3);
  
  // ***** CB info (tracklets, clusters, chips)
  const AliMultiplicity *mult = esd->GetMultiplicity();
  Int_t nTracklets = mult->GetNumberOfTracklets();
  Int_t nSingleClusters;
  Int_t nClusters[6];
  
  for(Int_t ilay = 0; ilay < 6; ilay++){
    nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
  }
  nSingleClusters = mult->GetNumberOfSingleClusters();

  Int_t nChips[2];
  for(Int_t ilay = 0; ilay < 2; ilay++){
    nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
  }
  
  // ***** FMD info
  AliESDFMD *fmd = esd->GetFMDData();
  Float_t totalMultA = 0;
  Float_t totalMultC = 0;
  const Float_t fmdLowCut = 0.4;
  
  for(UShort_t det = 1;det <= 3; det++) {
      Int_t nRings = (det==1 ? 1 : 2);
      for (UShort_t ir = 0; ir < nRings; ir++) {	  
	  Char_t   ring = (ir == 0 ? 'I' : 'O');
	  UShort_t nsec = (ir == 0 ? 20  : 40);
	  UShort_t nstr = (ir == 0 ? 512 : 256);
	  for(UShort_t sec =0; sec < nsec;  sec++)  {
	      for(UShort_t strip = 0; strip < nstr; strip++) {
		  Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
		  if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
		  Float_t nParticles=0;
		  if(fmdMult > fmdLowCut) {
		      nParticles = 1.;
		  }
	  
		  if (det<3) totalMultA = totalMultA + nParticles;
		  else totalMultC = totalMultC + nParticles;
		  
	      }
	  }
      }
  }
  Float_t multFMDA = totalMultA;
  Float_t multFMDC = totalMultC;
  
  // ***** ZDC info
  AliESDZDC *esdZDC = esd->GetESDZDC();
  UInt_t esdFlag =  esdZDC->GetESDQuality();   
  
  Float_t znCEnergy  = (Float_t) (esdZDC->GetZDCN1Energy());
  Float_t zpCEnergy  = (Float_t) (esdZDC->GetZDCP1Energy());
  Float_t znAEnergy  = (Float_t) (esdZDC->GetZDCN2Energy());
  Float_t zpAEnergy  = (Float_t) (esdZDC->GetZDCP2Energy());
  Float_t zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0));
  Float_t zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1));
  
  Double_t bZDC      = esdZDC->GetImpactParameter();
  Int_t    nPartZDC  = esdZDC->GetZDCParticipants();
  Double_t bZDCA     = esdZDC->GetImpactParamSideA();
  Int_t    nPartZDCA = esdZDC->GetZDCPartSideA();
  Double_t bZDCC     = esdZDC->GetImpactParamSideC();
  Int_t    nPartZDCC = esdZDC->GetZDCPartSideC();
  
  const Double_t * towZNC = esdZDC->GetZN1TowerEnergy();
  const Double_t * towZPC = esdZDC->GetZP1TowerEnergy();
  const Double_t * towZNA = esdZDC->GetZN2TowerEnergy();
  const Double_t * towZPA = esdZDC->GetZP2TowerEnergy();
  //
  Float_t  znCtower[5];	//  ZNC 5 tower signals
  Float_t  zpCtower[5];	//  ZPC 5 tower signals
  Float_t  znAtower[5];	//  ZNA 5 tower signals
  Float_t  zpAtower[5];	//  ZPA 5 tower signals
  Float_t  centrZNC[2];	//  centroid over ZNC
  Float_t  centrZNA[2];	//  centroid over ZNA

  for(Int_t it = 0; it < 5; it++){
    znCtower[it] = (Float_t) (towZNC[it]);
    zpCtower[it] = (Float_t) (towZPC[it]);
    znAtower[it] = (Float_t) (towZNA[it]); 
    zpAtower[it] = (Float_t) (towZPA[it]);  
  }
  
  Double_t xyZNC[2] = {-99.,-99.};
  Double_t xyZNA[2] = {-99.,-99.};

  esdZDC->GetZNCentroidInPbPb(beamEnergy, xyZNC, xyZNA);
  for(Int_t it = 0; it < 2; it++){
      centrZNC[it] = xyZNC[it];
      centrZNA[it] = xyZNA[it];
  }

  // ***** MC info
  Double_t bMC          = 0.;
  Int_t specNeutronProj = 0;
  Int_t specProtonProj  = 0;
  Int_t specNeutronTarg = 0;
  Int_t specProtonTarg  = 0;
  Int_t nPartTargMC     = 0;
  Int_t nPartProjMC     = 0;
  Int_t nnColl          = 0;
  Int_t nnwColl         = 0;
  Int_t nwNColl         = 0;
  Int_t nwNwColl        = 0;

  if(fIsMCInput){
    
    AliMCEvent* mcEvent = MCEvent();
    if (!mcEvent) {
      printf("   Could not retrieve MC event!!!\n");
      return;
    }
    
    Int_t nMyTracks_gen = 0;
    AliStack *stack = 0x0; // needed for MC studies
    stack = MCEvent()->Stack();
    for (Int_t iTrack = 0; iTrack < MCEvent()->GetNumberOfTracks(); iTrack++) {
      //get properties of mc particle
      AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(iTrack);
      // Primaries only
      if (!(stack->IsPhysicalPrimary(mcP->Label()))) continue;
      //charged tracks only
      if (mcP->Particle()->GetPDG()->Charge() == 0) continue;
      //same cuts as on ESDtracks
      // 	  if(TMath::Abs(mcP->Eta())>0.9)continue;
      // 	  if(mcP->Pt()<0.2)continue;
      // 	  if(mcP->Pt()>200)continue;
      
      nMyTracks_gen ++;
    } 
    
    AliGenEventHeader* genHeader = mcEvent->GenEventHeader();
    if(!genHeader){
      printf("  Event generator header not available!!!\n");
      return;
    }
	

    if(genHeader->InheritsFrom(AliGenHijingEventHeader::Class())){
	bMC = ((AliGenHijingEventHeader*) genHeader)->ImpactParameter();
	specNeutronProj = ((AliGenHijingEventHeader*) genHeader)->ProjSpectatorsn();
	specProtonProj  = ((AliGenHijingEventHeader*) genHeader)->ProjSpectatorsp();
	specNeutronTarg = ((AliGenHijingEventHeader*) genHeader)->TargSpectatorsn();
	specProtonTarg  = ((AliGenHijingEventHeader*) genHeader)->TargSpectatorsp();
	nPartTargMC = Int_t (208.-(specNeutronTarg+specProtonTarg));
	nPartProjMC = Int_t (208.-(specNeutronProj+specProtonProj));
	nnColl   = ((AliGenHijingEventHeader*) genHeader)->NN();
	nnwColl  = ((AliGenHijingEventHeader*) genHeader)->NNw();
	nwNColl  = ((AliGenHijingEventHeader*) genHeader)->NwN();
	nwNwColl = ((AliGenHijingEventHeader*) genHeader)->NwNw();
    }  
  }
  
  fAODCentrality->SetTrigClass         (trigClass[100]);
  fAODCentrality->SetxVertex	       (xVertex       );
  fAODCentrality->SetyVertex	       (yVertex       );
  fAODCentrality->SetzVertex	       (zVertex       );
  fAODCentrality->SetVertexer3d        (vertexer3d    );
  fAODCentrality->SetbMC               (bMC);
  fAODCentrality->SetNpartTargMC       (nPartTargMC);
  fAODCentrality->SetNpartProjMC       (nPartProjMC);
  fAODCentrality->SetNNColl            (nnColl);
  fAODCentrality->SetNNwColl           (nnwColl);
  fAODCentrality->SetNwNColl           (nwNColl);
  fAODCentrality->SetNwNwColl          (nwNwColl);
  fAODCentrality->SetNTracklets        (nTracklets);
  fAODCentrality->SetNSingleClusters   (nSingleClusters);
  fAODCentrality->SetNClusters         (
					   nClusters[0],
					   nClusters[1],
					   nClusters[2],
					   nClusters[3],
					   nClusters[4],
					   nClusters[5]);
  fAODCentrality->SetNChips            (
					   nChips[0],
					   nChips[1]);
  fAODCentrality->SetbZDC              (bZDC);
  fAODCentrality->SetNpartZDC          (nPartZDC);
  fAODCentrality->SetbZDCA             (bZDCA);
  fAODCentrality->SetNpartZDCA         (nPartZDCA);
  fAODCentrality->SetbZDCC             (bZDCC);
  fAODCentrality->SetNpartZDCC         (nPartZDCC);
  fAODCentrality->SetESDFlag  	  (esdFlag);
  fAODCentrality->SetZNCEnergy 	  (znCEnergy);
  fAODCentrality->SetZPCEnergy	  (zpCEnergy);
  fAODCentrality->SetZNAEnergy	  (znAEnergy);
  fAODCentrality->SetZPAEnergy	  (zpAEnergy);
  fAODCentrality->SetZEM1Energy	  (zem1Energy);
  fAODCentrality->SetZEM2Energy	  (zem2Energy);
  fAODCentrality->SetZNCtower          (
					   znCtower[0],
					   znCtower[1],
					   znCtower[2],
					   znCtower[3],
					   znCtower[4]);
  fAODCentrality->SetZPCtower          (
					 zpCtower[0],
					 zpCtower[1],
					 zpCtower[2],
					 zpCtower[3],
					 zpCtower[4]);
 
  fAODCentrality-> SetZNAtower          (
			     znAtower[0],  
			     znAtower[1], 
			     znAtower[2],  
			     znAtower[3],  
			     znAtower[4]); 
  fAODCentrality-> SetZPAtower          (
					    zpAtower[0], 
					    zpAtower[1], 
					    zpAtower[2], 
					    zpAtower[3],
					    zpAtower[4]);
  fAODCentrality-> SetCentrZNC          (
					    centrZNC[0],
					    centrZNC[1]);
  fAODCentrality-> SetCentrZNA          (
					    centrZNA[0],
					    centrZNA[1]);
  fAODCentrality-> SetNTracks           (nTracks);
  fAODCentrality-> SetNPmdTracks        (nPmdTracks);
  fAODCentrality-> SetMultV0A           (multV0A);
  fAODCentrality-> SetMultV0C           (multV0C);
  fAODCentrality-> SetMultFMDA          (multFMDA);
  fAODCentrality-> SetMultFMDC          (multFMDC);

//
// Header Replication
//
  AliAODHeader* hdr = AODEvent()->GetHeader();
  *fAODHeader =  *hdr;
//
  return;
}

//________________________________________________________________________
void AliAnalysisTaskAODCentralityMaker::Terminate(Option_t */*option*/)
{
  // Terminate analysis
  //
  if(fDebug > 1) printf("AnalysisTaskAODCentralityMaker: Terminate() \n");
}
Commit	Line	Data
c94a2509	1	/**************************************************************************
	2	* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16	/////////////////////////////////////////////////////////////
	17	//
	18	// AliAnalysisTaskSE to make AOD centrality
	19	// Author: Alberica Toia, CERN, Alberica.Toia@cern.ch
	20	//
	21	/////////////////////////////////////////////////////////////
	22
	23	#include <TROOT.h>
	24	#include <TSystem.h>
	25
	26	#include "AliVEvent.h"
	27	#include "AliESDEvent.h"
	28	#include "AliAODEvent.h"
	29	#include "AliAODCentrality.h"
c94a2509	30	#include "AliAnalysisManager.h"
c94a2509	31	#include "AliESDInputHandler.h"
	32	#include "AliESDZDC.h"
	33	#include "AliESDFMD.h"
	34	#include "AliESDVZERO.h"
	35	#include "AliMultiplicity.h"
	36	#include "AliAODHandler.h"
42f0d071	37	#include "AliAODHeader.h"
c94a2509	38	#include "AliAODEvent.h"
c94a2509	39	#include "AliAODVertex.h"
c94a2509	40	#include "AliMCEvent.h"
	41	#include "AliMCEventHandler.h"
	42	#include "AliMCParticle.h"
	43	#include "AliStack.h"
	44	#include "AliHeader.h"
c94a2509	45	#include "AliGenEventHeader.h"
c94a2509	46	#include "AliGenHijingEventHeader.h"
c94a2509	47	#include "AliAnalysisTaskAODCentralityMaker.h"
	48
	49	ClassImp(AliAnalysisTaskAODCentralityMaker)
	50
	51
	52	//________________________________________________________________________
	53	AliAnalysisTaskAODCentralityMaker::AliAnalysisTaskAODCentralityMaker():
	54	AliAnalysisTaskSE(),
	55	fAODCentrality(0),
	56	fDeltaAODFileName("AliAOD.Centrality.root"),
42f0d071	57	fAODHeader (0),
42f0d071	58	fIsMCInput (0)
c94a2509	59	{
c94a2509	60	// Default constructor
c94a2509	61	}
	62
	63	//________________________________________________________________________
	64	AliAnalysisTaskAODCentralityMaker::AliAnalysisTaskAODCentralityMaker(const char *name):
	65	AliAnalysisTaskSE(name),
	66	fAODCentrality(0),
	67	fDeltaAODFileName("AliAOD.Centrality.root"),
42f0d071	68	fAODHeader (0),
42f0d071	69	fIsMCInput (0)
c94a2509	70	{
c94a2509	71	// Standard constructor
c94a2509	72	}
	73
	74
	75	//________________________________________________________________________
	76	AliAnalysisTaskAODCentralityMaker::~AliAnalysisTaskAODCentralityMaker()
	77	{
	78	// Destructor
	79	}
	80
	81	//________________________________________________________________________
	82	void AliAnalysisTaskAODCentralityMaker::Init()
	83	{
	84	// Initialization
	85	if(fDebug > 1) printf("AnalysisTaskAODCentralityMaker::Init() \n");
	86	AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(fDeltaAODFileName.Data());
	87
	88	return;
	89	}
	90
	91	//________________________________________________________________________
	92	void AliAnalysisTaskAODCentralityMaker::UserCreateOutputObjects()
	93	{
	94
	95	// Create the output container
	96	//
	97	if(fDebug > 1) printf("AnalysisTaskAODCentralityMaker::UserCreateOutPutData() \n");
	98	// Support both the case when the AOD + deltaAOD are produced in an ESD
	99	// analysis or if the deltaAOD is produced on an analysis on AOD's. (A.G. 27/04/09)
	100	if(!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) {
	101	Fatal("UserCreateOutputObjects", "This task needs an AOD handler");
	102	return;
	103	}
	104	TString filename = fDeltaAODFileName;
	105	// When running on standard AOD to produce deltas, IsStandardAOD is never set,
	106	// If AODEvent is NULL, new branches have to be added to the new file(s) (A.G. 15/01/10)
	107	if(!IsStandardAOD() && AODEvent()) filename = "";
	108
	109	fAODCentrality = new AliAODCentrality();
	110	fAODCentrality->SetName("AODCentrality");
	111	AddAODBranch("AliAODCentrality", &fAODCentrality, filename);
	112
42f0d071	113
	114	fAODHeader = new AliAODHeader();
	115	AddAODBranch("AliAODHeader", &fAODHeader, filename);
c94a2509	116	return;
	117	}
	118
42f0d071	119
c94a2509	120	//________________________________________________________________________
	121	void AliAnalysisTaskAODCentralityMaker::UserExec(Option_t /option*/)
	122	{
42f0d071	123	AliVEvent* event = InputEvent();
42f0d071	124	AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
c94a2509	125
42f0d071	126	Float_t beamEnergy = esd->GetBeamEnergy();
	127	Int_t nTracks = event->GetNumberOfTracks();
	128	Int_t nPmdTracks = esd->GetNumberOfPmdTracks();
c94a2509	129
	130	// ***** V0 info
	131	AliESDVZERO* esdV0 = esd->GetVZEROData();
42f0d071	132	Double_t multV0A = esdV0->GetMTotV0A();
42f0d071	133	Double_t multV0C = esdV0->GetMTotV0C();
c94a2509	134
c94a2509	135	// ***** Trigger selection
42f0d071	136	char trigClass[100]; // fired trigger classes
c94a2509	137	TString triggerClass = esd->GetFiredTriggerClasses();
42f0d071	138	sprintf(trigClass,"%s",triggerClass.Data());
c94a2509	139
	140	// ***** vertex info
	141	const AliESDVertex *vertex = esd->GetPrimaryVertexSPD();
42f0d071	142	Double_t xVertex = vertex->GetX();
	143	Double_t yVertex = vertex->GetY();
	144	Double_t zVertex = vertex->GetZ();
	145	Bool_t vertexer3d;
	146
	147	if(vertex->IsFromVertexer3D()) vertexer3d = kTRUE;
	148	else vertexer3d = kFALSE;
c94a2509	149	Double_t vertex3[3];
	150	vertex->GetXYZ(vertex3);
	151
	152	// ***** CB info (tracklets, clusters, chips)
	153	const AliMultiplicity *mult = esd->GetMultiplicity();
42f0d071	154	Int_t nTracklets = mult->GetNumberOfTracklets();
	155	Int_t nSingleClusters;
	156	Int_t nClusters[6];
c94a2509	157
42f0d071	158	for(Int_t ilay = 0; ilay < 6; ilay++){
42f0d071	159	nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
c94a2509	160	}
42f0d071	161	nSingleClusters = mult->GetNumberOfSingleClusters();
	162
	163	Int_t nChips[2];
	164	for(Int_t ilay = 0; ilay < 2; ilay++){
	165	nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
c94a2509	166	}
	167
	168	// ***** FMD info
	169	AliESDFMD *fmd = esd->GetFMDData();
	170	Float_t totalMultA = 0;
	171	Float_t totalMultC = 0;
42f0d071	172	const Float_t fmdLowCut = 0.4;
c94a2509	173
42f0d071	174	for(UShort_t det = 1;det <= 3; det++) {
	175	Int_t nRings = (det==1 ? 1 : 2);
	176	for (UShort_t ir = 0; ir < nRings; ir++) {
	177	Char_t ring = (ir == 0 ? 'I' : 'O');
	178	UShort_t nsec = (ir == 0 ? 20 : 40);
	179	UShort_t nstr = (ir == 0 ? 512 : 256);
	180	for(UShort_t sec =0; sec < nsec; sec++) {
	181	for(UShort_t strip = 0; strip < nstr; strip++) {
	182	Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
	183	if(fmdMult == 0 \|\| fmdMult == AliESDFMD::kInvalidMult) continue;
	184	Float_t nParticles=0;
	185	if(fmdMult > fmdLowCut) {
	186	nParticles = 1.;
	187	}
c94a2509	188
42f0d071	189	if (det<3) totalMultA = totalMultA + nParticles;
	190	else totalMultC = totalMultC + nParticles;
	191
	192	}
c94a2509	193	}
c94a2509	194	}
c94a2509	195	}
42f0d071	196	Float_t multFMDA = totalMultA;
42f0d071	197	Float_t multFMDC = totalMultC;
c94a2509	198
	199	// ***** ZDC info
	200	AliESDZDC *esdZDC = esd->GetESDZDC();
42f0d071	201	UInt_t esdFlag = esdZDC->GetESDQuality();
c94a2509	202
42f0d071	203	Float_t znCEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
	204	Float_t zpCEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
	205	Float_t znAEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
	206	Float_t zpAEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
	207	Float_t zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0));
	208	Float_t zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1));
c94a2509	209
42f0d071	210	Double_t bZDC = esdZDC->GetImpactParameter();
	211	Int_t nPartZDC = esdZDC->GetZDCParticipants();
	212	Double_t bZDCA = esdZDC->GetImpactParamSideA();
	213	Int_t nPartZDCA = esdZDC->GetZDCPartSideA();
	214	Double_t bZDCC = esdZDC->GetImpactParamSideC();
	215	Int_t nPartZDCC = esdZDC->GetZDCPartSideC();
c94a2509	216
	217	const Double_t * towZNC = esdZDC->GetZN1TowerEnergy();
	218	const Double_t * towZPC = esdZDC->GetZP1TowerEnergy();
	219	const Double_t * towZNA = esdZDC->GetZN2TowerEnergy();
	220	const Double_t * towZPA = esdZDC->GetZP2TowerEnergy();
	221	//
42f0d071	222	Float_t znCtower[5]; // ZNC 5 tower signals
	223	Float_t zpCtower[5]; // ZPC 5 tower signals
	224	Float_t znAtower[5]; // ZNA 5 tower signals
	225	Float_t zpAtower[5]; // ZPA 5 tower signals
	226	Float_t centrZNC[2]; // centroid over ZNC
	227	Float_t centrZNA[2]; // centroid over ZNA
	228
	229	for(Int_t it = 0; it < 5; it++){
	230	znCtower[it] = (Float_t) (towZNC[it]);
	231	zpCtower[it] = (Float_t) (towZPC[it]);
	232	znAtower[it] = (Float_t) (towZNA[it]);
	233	zpAtower[it] = (Float_t) (towZPA[it]);
c94a2509	234	}
c94a2509	235
42f0d071	236	Double_t xyZNC[2] = {-99.,-99.};
	237	Double_t xyZNA[2] = {-99.,-99.};
	238
	239	esdZDC->GetZNCentroidInPbPb(beamEnergy, xyZNC, xyZNA);
	240	for(Int_t it = 0; it < 2; it++){
	241	centrZNC[it] = xyZNC[it];
	242	centrZNA[it] = xyZNA[it];
c94a2509	243	}
	244
	245	// ***** MC info
42f0d071	246	Double_t bMC = 0.;
	247	Int_t specNeutronProj = 0;
	248	Int_t specProtonProj = 0;
	249	Int_t specNeutronTarg = 0;
	250	Int_t specProtonTarg = 0;
	251	Int_t nPartTargMC = 0;
	252	Int_t nPartProjMC = 0;
	253	Int_t nnColl = 0;
	254	Int_t nnwColl = 0;
	255	Int_t nwNColl = 0;
	256	Int_t nwNwColl = 0;
	257
c94a2509	258	if(fIsMCInput){
	259
	260	AliMCEvent* mcEvent = MCEvent();
	261	if (!mcEvent) {
	262	printf(" Could not retrieve MC event!!!\n");
	263	return;
	264	}
	265
42f0d071	266	Int_t nMyTracks_gen = 0;
c94a2509	267	AliStack *stack = 0x0; // needed for MC studies
	268	stack = MCEvent()->Stack();
	269	for (Int_t iTrack = 0; iTrack < MCEvent()->GetNumberOfTracks(); iTrack++) {
	270	//get properties of mc particle
	271	AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(iTrack);
	272	// Primaries only
	273	if (!(stack->IsPhysicalPrimary(mcP->Label()))) continue;
	274	//charged tracks only
	275	if (mcP->Particle()->GetPDG()->Charge() == 0) continue;
	276	//same cuts as on ESDtracks
	277	// if(TMath::Abs(mcP->Eta())>0.9)continue;
	278	// if(mcP->Pt()<0.2)continue;
	279	// if(mcP->Pt()>200)continue;
	280
42f0d071	281	nMyTracks_gen ++;
c94a2509	282	}
	283
	284	AliGenEventHeader* genHeader = mcEvent->GenEventHeader();
	285	if(!genHeader){
	286	printf(" Event generator header not available!!!\n");
	287	return;
	288	}
	289
42f0d071	290
c94a2509	291	if(genHeader->InheritsFrom(AliGenHijingEventHeader::Class())){
42f0d071	292	bMC = ((AliGenHijingEventHeader*) genHeader)->ImpactParameter();
	293	specNeutronProj = ((AliGenHijingEventHeader*) genHeader)->ProjSpectatorsn();
	294	specProtonProj = ((AliGenHijingEventHeader*) genHeader)->ProjSpectatorsp();
	295	specNeutronTarg = ((AliGenHijingEventHeader*) genHeader)->TargSpectatorsn();
	296	specProtonTarg = ((AliGenHijingEventHeader*) genHeader)->TargSpectatorsp();
	297	nPartTargMC = Int_t (208.-(specNeutronTarg+specProtonTarg));
	298	nPartProjMC = Int_t (208.-(specNeutronProj+specProtonProj));
	299	nnColl = ((AliGenHijingEventHeader*) genHeader)->NN();
	300	nnwColl = ((AliGenHijingEventHeader*) genHeader)->NNw();
	301	nwNColl = ((AliGenHijingEventHeader*) genHeader)->NwN();
	302	nwNwColl = ((AliGenHijingEventHeader*) genHeader)->NwNw();
c94a2509	303	}
c94a2509	304	}
c94a2509	305
42f0d071	306	fAODCentrality->SetTrigClass (trigClass[100]);
	307	fAODCentrality->SetxVertex (xVertex );
	308	fAODCentrality->SetyVertex (yVertex );
	309	fAODCentrality->SetzVertex (zVertex );
	310	fAODCentrality->SetVertexer3d (vertexer3d );
	311	fAODCentrality->SetbMC (bMC);
	312	fAODCentrality->SetNpartTargMC (nPartTargMC);
	313	fAODCentrality->SetNpartProjMC (nPartProjMC);
	314	fAODCentrality->SetNNColl (nnColl);
	315	fAODCentrality->SetNNwColl (nnwColl);
	316	fAODCentrality->SetNwNColl (nwNColl);
	317	fAODCentrality->SetNwNwColl (nwNwColl);
	318	fAODCentrality->SetNTracklets (nTracklets);
	319	fAODCentrality->SetNSingleClusters (nSingleClusters);
c94a2509	320	fAODCentrality->SetNClusters (
42f0d071	321	nClusters[0],
	322	nClusters[1],
	323	nClusters[2],
	324	nClusters[3],
	325	nClusters[4],
	326	nClusters[5]);
c94a2509	327	fAODCentrality->SetNChips (
42f0d071	328	nChips[0],
	329	nChips[1]);
	330	fAODCentrality->SetbZDC (bZDC);
	331	fAODCentrality->SetNpartZDC (nPartZDC);
	332	fAODCentrality->SetbZDCA (bZDCA);
	333	fAODCentrality->SetNpartZDCA (nPartZDCA);
	334	fAODCentrality->SetbZDCC (bZDCC);
	335	fAODCentrality->SetNpartZDCC (nPartZDCC);
	336	fAODCentrality->SetESDFlag (esdFlag);
	337	fAODCentrality->SetZNCEnergy (znCEnergy);
	338	fAODCentrality->SetZPCEnergy (zpCEnergy);
	339	fAODCentrality->SetZNAEnergy (znAEnergy);
	340	fAODCentrality->SetZPAEnergy (zpAEnergy);
	341	fAODCentrality->SetZEM1Energy (zem1Energy);
	342	fAODCentrality->SetZEM2Energy (zem2Energy);
c94a2509	343	fAODCentrality->SetZNCtower (
42f0d071	344	znCtower[0],
	345	znCtower[1],
	346	znCtower[2],
	347	znCtower[3],
	348	znCtower[4]);
c94a2509	349	fAODCentrality->SetZPCtower (
42f0d071	350	zpCtower[0],
	351	zpCtower[1],
	352	zpCtower[2],
	353	zpCtower[3],
	354	zpCtower[4]);
c94a2509	355
c94a2509	356	fAODCentrality-> SetZNAtower (
42f0d071	357	znAtower[0],
	358	znAtower[1],
	359	znAtower[2],
	360	znAtower[3],
	361	znAtower[4]);
c94a2509	362	fAODCentrality-> SetZPAtower (
42f0d071	363	zpAtower[0],
	364	zpAtower[1],
	365	zpAtower[2],
	366	zpAtower[3],
	367	zpAtower[4]);
c94a2509	368	fAODCentrality-> SetCentrZNC (
42f0d071	369	centrZNC[0],
42f0d071	370	centrZNC[1]);
c94a2509	371	fAODCentrality-> SetCentrZNA (
42f0d071	372	centrZNA[0],
	373	centrZNA[1]);
	374	fAODCentrality-> SetNTracks (nTracks);
	375	fAODCentrality-> SetNPmdTracks (nPmdTracks);
	376	fAODCentrality-> SetMultV0A (multV0A);
	377	fAODCentrality-> SetMultV0C (multV0C);
	378	fAODCentrality-> SetMultFMDA (multFMDA);
	379	fAODCentrality-> SetMultFMDC (multFMDC);
c94a2509	380
42f0d071	381	//
	382	// Header Replication
	383	//
	384	AliAODHeader* hdr = AODEvent()->GetHeader();
	385	fAODHeader = hdr;
	386	//
c94a2509	387	return;
	388	}
	389
	390	//________________________________________________________________________
	391	void AliAnalysisTaskAODCentralityMaker::Terminate(Option_t /option*/)
	392	{
	393	// Terminate analysis
	394	//
	395	if(fDebug > 1) printf("AnalysisTaskAODCentralityMaker: Terminate() \n");
	396	}