[u/mrichter/AliRoot.git] / PWGGA / PHOSTasks / UserTasks / AliAnalysisTaskSEPHOSpPbPi0.cxx

/**************************************************************************
 * Copyright(c) 1998-2006, 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 for the gamma and pi0 from pPb collision analysis
//
// Author: H-S. Zhu, hongsheng.zhu@cern.ch
//                   hszhu@iopp.ccnu.edu.cn
///////////////////////////////////////////////////////////////////////////

#include <TH2I.h>
#include <TList.h>
#include <TMath.h>
#include <TArray.h>
#include <TClonesArray.h>

#include "AliAnalysisManager.h"
#include "AliMCEventHandler.h"
#include "AliMCEvent.h"
#include "AliStack.h"
#include "AliVEvent.h"
#include "AliESDEvent.h"
#include "AliAODEvent.h"
#include "AliESDHeader.h"
#include "AliAODHeader.h"
#include "AliVCluster.h"
#include "AliVCaloCells.h"
#include "AliPHOSGeoUtils.h"
#include "AliPHOSGeometry.h"
#include "AliOADBContainer.h"
#include "AliInputEventHandler.h"
#include "AliAnalysisManager.h"
#include "AliPHOSpPbPi0Header.h"
#include "AliCaloClusterInfo.h"
#include "AliAnalysisTaskSEPHOSpPbPi0.h"

ClassImp(AliAnalysisTaskSEPHOSpPbPi0)

Int_t    AliAnalysisTaskSEPHOSpPbPi0::fgMinNCells        = 2;
Double_t AliAnalysisTaskSEPHOSpPbPi0::fgMinClusterEnergy = 0.3;
Double_t AliAnalysisTaskSEPHOSpPbPi0::fgMinM02           = 0.2;
Double_t AliAnalysisTaskSEPHOSpPbPi0::fgMinDistToBad     = 2.5;

//________________________________________________________________________
AliAnalysisTaskSEPHOSpPbPi0::AliAnalysisTaskSEPHOSpPbPi0():
  AliAnalysisTaskSE(), fIsMC(kFALSE), fCentralityBin(10), fBufferSize(10), fRunNumber(-1), fPHOSGeo(0),
  fListEvent(0), fListCaloCl(0), fListPi0(0), fListMC(0), fHeader(0), fCaloClArr(0)
{
  //
  // Default constructor
  //
  for (Int_t i=0; i<10; i++) { for (Int_t j=0; j<10; j++) fEventList[i][j] = 0; }

  // Set bad channel map
  for (Int_t i=0; i<5; i++)
    fPHOSBadMap[i] = new TH2I(Form("PHOS_BadMap_mod%d", i), Form("PHOS_BadMap_mod%d", i), 64, 0., 64., 56, 0., 56.);
}
//________________________________________________________________________
AliAnalysisTaskSEPHOSpPbPi0::AliAnalysisTaskSEPHOSpPbPi0(const char *name):
  AliAnalysisTaskSE(name), fIsMC(kFALSE), fCentralityBin(10), fBufferSize(10), fRunNumber(-1), fPHOSGeo(0),
  fListEvent(0), fListCaloCl(0), fListPi0(0), fListMC(0), fHeader(0), fCaloClArr(0)
{
  // Constructor
  for (Int_t i=0; i<10; i++) { for (Int_t j=0; j<10; j++) fEventList[i][j] = 0; }

  // Set bad channel map
  for (Int_t i=0; i<5; i++)
    fPHOSBadMap[i] = new TH2I(Form("PHOS_BadMap_mod%d", i), Form("PHOS_BadMap_mod%d", i), 64, 0., 64., 56, 0., 56.);

  DefineOutput(1, TList::Class());
  DefineOutput(2, TList::Class());
  DefineOutput(3, TList::Class());
  DefineOutput(4, TList::Class());
}

//_____________________________________________________________________________
AliAnalysisTaskSEPHOSpPbPi0::~AliAnalysisTaskSEPHOSpPbPi0()
{
  //
  // Default destructor
  //
  if (fListEvent)   { delete fListEvent;   fListEvent  = NULL; }
  if (fListCaloCl)  { delete fListCaloCl;  fListCaloCl = NULL; }
  if (fListPi0)     { delete fListPi0;     fListPi0    = NULL; }
  if (fListMC)      { delete fListMC;      fListMC     = NULL; }
  if (fHeader)      { delete fHeader;      fHeader     = NULL; }
  if (fCaloClArr)   { delete fCaloClArr;   fCaloClArr  = NULL; }
}

//________________________________________________________________________
void AliAnalysisTaskSEPHOSpPbPi0::UserCreateOutputObjects()
{
  // Create the output container
  // Initialize the PHOS geometry
//fPHOSGeo = new AliPHOSGeoUtils("PHOSGeo");
//fPHOSGeo =  AliPHOSGeometry::GetInstance("IHEP");

  if (!fHeader)     fHeader     = new AliPHOSpPbPi0Header();
  if (!fCaloClArr)  fCaloClArr  = new TClonesArray("AliCaloClusterInfo", 0);
  if (!fListEvent)  fListEvent  = new TList();
  if (!fListCaloCl) fListCaloCl = new TList();
  if (!fListPi0)    fListPi0    = new TList();
  if (!fListMC)     fListMC     = new TList();

  fHeader->SetIsMC(fIsMC);
  fHeader->SetNCent(fCentralityBin.GetSize()-1);
  fHeader->CreateHistograms(fListEvent, fListCaloCl, fListPi0, fListMC);

  // Post output data.
  PostData(1, fListEvent);
  PostData(2, fListCaloCl);
  PostData(3, fListPi0);
  PostData(4, fListMC);

  return;
}

//________________________________________________________________________
void AliAnalysisTaskSEPHOSpPbPi0::UserExec(Option_t *) 
{
  // Main loop
  // Called for each event

  if (fIsMC) {
    if (MCEvent()) {
      if (MCEvent()->GetNumberOfTracks()<=0)
           { AliError("MC event not found. Nothing done!"); return; }
    } else { AliError("MC event not found. Nothing done!"); return; }
  }

  AliAODEvent *aod   = 0;
  AliESDEvent *esd   = 0;

  if (((TString)fInputEvent->IsA()->GetName())=="AliAODEvent") {
    aod = dynamic_cast<AliAODEvent*>(fInputEvent);
    if (!aod) { AliError("AOD event not found. Nothing done!");   return; }
    if (!fIsMC && (aod->GetHeader()->GetEventType()!=7))          return; // check event type; should be PHYSICS = 7 for data and 0 for MC
  } else {
    esd = dynamic_cast<AliESDEvent*>(fInputEvent);
    if (!esd) { AliError("ESD event not found. Nothing done!");   return; }
    if (!fIsMC && (esd->GetHeader()->GetEventType()!=7))          return; // check event type; should be PHYSICS = 7 for data and 0 for MC
  }

  // Fill Event info
  fHeader->SetEventInfo(fInputHandler);
  fHeader->FillHistosEvent(fListEvent);

  if (!fHeader->IsSelected()) return; // event selection

  // PHOS Geometry and Misalignment initialization at the first time it runs
  if(fRunNumber != fInputEvent->GetRunNumber()) {
    fRunNumber = fInputEvent->GetRunNumber();
    PHOSInitialize(esd);
  }

  // Fill PHOS cells QA histograms
  fHeader->FillHistosCaloCellsQA(fListCaloCl, fInputEvent->GetPHOSCells(), fPHOSGeo);
  aod = 0;

  // Fill PHOS cluster Clones Array
  FillCaloClusterInfo(esd);

  if (!fCaloClArr->GetEntriesFast()) return;

  // vertex bining and centrality bining
  Int_t zvtx = (Int_t)((fHeader->Vz() + 10.)/2.);   if (zvtx<0) zvtx = 0;   if (zvtx>9) zvtx = 9;
  Int_t cent = TMath::BinarySearch<Float_t>(fCentralityBin.GetSize()-1, fCentralityBin.GetArray(), fHeader->Centrality());
  if (!fEventList[zvtx][cent]) fEventList[zvtx][cent] = new TList(); 
  TList *eventList = fEventList[zvtx][cent];

  // Fill cluster histograms
  fHeader->FillHistosCaloCluster(fListCaloCl, fCaloClArr, cent);

  // Fill pi0 histograms
  fHeader->FillHistosPi0(fListPi0, fCaloClArr, cent);

  // Fill mixed pi0 histograms
  fHeader->FillHistosMixPi0(fListPi0, fCaloClArr, eventList, cent);

  // Fill MC info
  AliStack *stack = 0x0;
  if (fIsMC) {
    if (esd) {
      if (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()) {
      if (static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent())
        stack = static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent()->Stack();
      }
      fHeader->FillHistosMC(fListMC, stack, fPHOSGeo, cent);
    } else
      fHeader->FillHistosMC(fListMC, MCEvent(), fPHOSGeo, cent);
  }
  esd = 0;

  // Fill event list for mixing
  if (fCaloClArr->GetEntriesFast()>0) {
    eventList->AddFirst(fCaloClArr);   fCaloClArr = 0;
    if (eventList->GetSize()>fBufferSize[cent]) { // Remove redundant events
      TClonesArray *tmp = static_cast<TClonesArray*>(eventList->Last());
      eventList->RemoveLast();
      delete tmp;
    }
  }

  return;
}

//________________________________________________________________________
void AliAnalysisTaskSEPHOSpPbPi0::Terminate(Option_t *) 
{
  // Terminate analysis

  // add the correction matrix

  return;
}

//________________________________________________________________________
void AliAnalysisTaskSEPHOSpPbPi0::PHOSInitialize(AliESDEvent* const esd)
{
  // Initialize PHOS Geometry and misalignment
  // PHOS Geometry
  AliOADBContainer geomContainer("phosGeo");
  geomContainer.InitFromFile("$ALICE_ROOT/OADB/PHOS/PHOSGeometry.root","PHOSRotationMatrixes");
  TObjArray *matrixes = (TObjArray*)geomContainer.GetObject(fRunNumber,"PHOSRotationMatrixes");
  fPHOSGeo =  AliPHOSGeometry::GetInstance("IHEP");
  for (Int_t mod=0; mod<5; mod++) {
    if (!matrixes->At(mod)) continue;
    else fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod);
  }

  // sets the PHOS Misalignment vertex if ESD
  if (esd) {
    for (Int_t mod=0; mod<5; mod++) {
      const TGeoHMatrix* modMatrix = fInputEvent->GetPHOSMatrix(mod);
      if (!modMatrix) continue;
      else fPHOSGeo->SetMisalMatrix(modMatrix, mod);
    }
  }
}

//________________________________________________________________________
void AliAnalysisTaskSEPHOSpPbPi0::FillCaloClusterInfo(AliESDEvent* const esd)
{
  // Fill calo cluster info
  if (fCaloClArr) fCaloClArr->Clear();
  else fCaloClArr = new TClonesArray("AliCaloClusterInfo", 0);

  Int_t    nclsts       = fInputEvent->GetNumberOfCaloClusters();
  Int_t    countN       = 0;
  Int_t    relId[4]     = {0,0,0,0}; // module = relId[0]; cellX = relId[2]; cellZ = relId[3];
  Float_t  position[3]  = {0,0,0};
  Double_t vtx[3]       = {0,0,0};   fHeader->GetXYZ(vtx);
  TClonesArray &caloRef = *fCaloClArr;
  TLorentzVector momentum;
  AliVCluster        *clust       = 0;
  AliCaloClusterInfo *caloCluster = 0;
  for (Int_t iclst=0; iclst<nclsts; iclst++) {  // loop over all clusters
    clust = fInputEvent->GetCaloCluster(iclst);
    if (!(clust && clust->IsPHOS() && clust->E()>fgMinClusterEnergy))     { clust=0; continue; }
    if (!(clust->GetNCells()>fgMinNCells && clust->GetM02()>fgMinM02))    { clust=0; continue; } // To remove exotic clusters
    if (!(clust->GetDistanceToBadChannel()>fgMinDistToBad))               { clust=0; continue; }
    clust->GetPosition(position); TVector3 global(position);
    fPHOSGeo->GlobalPos2RelId(global,relId);
    if (!IsGoodCaloCluster(relId[0], relId[2], relId[3]))                 { clust=0; continue; } // calo cluster selection 
    if (relId[0] == 2)                                                    { clust=0; continue; } // !remove module 2

    caloCluster = new AliCaloClusterInfo(clust, esd, fPHOSGeo, vtx);
    if (esd && !(caloCluster->LorentzVector().E()>fgMinClusterEnergy)) { delete caloCluster; caloCluster=0; clust=0; continue; } // check again for ESD
    if (caloCluster->GetTrackPt()==-1. || caloCluster->TestCPV(fHeader->MagneticField()))    caloCluster->SetPIDBit(BIT(0));     // set CPV bit

    clust = 0;

    new(caloRef[countN++]) AliCaloClusterInfo(*caloCluster);

    delete caloCluster; caloCluster=0;
  }  // end loop of all clusters
}

//________________________________________________________________________
Bool_t AliAnalysisTaskSEPHOSpPbPi0::IsGoodCaloCluster(Int_t iMod, Int_t cellX, Int_t cellZ)
{
  // !Check whether this cluster is not in bad channel

  if (!(iMod>0 && iMod<5 && fPHOSBadMap[iMod]))          return kTRUE;   //No bad maps for this Module
  if (fPHOSBadMap[iMod]->GetBinContent(cellX,cellZ)>0)   return kFALSE;

  return kTRUE;
}
Commit	Line	Data
a01cd2c9	1	/**************************************************************************
	2	* Copyright(c) 1998-2006, 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 for the gamma and pi0 from pPb collision analysis
	19	//
	20	// Author: H-S. Zhu, hongsheng.zhu@cern.ch
	21	// hszhu@iopp.ccnu.edu.cn
	22	///////////////////////////////////////////////////////////////////////////
	23
	24	#include <TH2I.h>
	25	#include <TList.h>
	26	#include <TMath.h>
	27	#include <TArray.h>
	28	#include <TClonesArray.h>
	29
	30	#include "AliAnalysisManager.h"
	31	#include "AliMCEventHandler.h"
	32	#include "AliMCEvent.h"
	33	#include "AliStack.h"
	34	#include "AliVEvent.h"
	35	#include "AliESDEvent.h"
	36	#include "AliAODEvent.h"
	37	#include "AliESDHeader.h"
	38	#include "AliAODHeader.h"
	39	#include "AliVCluster.h"
	40	#include "AliVCaloCells.h"
	41	#include "AliPHOSGeoUtils.h"
	42	#include "AliPHOSGeometry.h"
	43	#include "AliOADBContainer.h"
	44	#include "AliInputEventHandler.h"
	45	#include "AliAnalysisManager.h"
	46	#include "AliPHOSpPbPi0Header.h"
	47	#include "AliCaloClusterInfo.h"
	48	#include "AliAnalysisTaskSEPHOSpPbPi0.h"
	49
	50	ClassImp(AliAnalysisTaskSEPHOSpPbPi0)
	51
	52	Int_t AliAnalysisTaskSEPHOSpPbPi0::fgMinNCells = 2;
	53	Double_t AliAnalysisTaskSEPHOSpPbPi0::fgMinClusterEnergy = 0.3;
	54	Double_t AliAnalysisTaskSEPHOSpPbPi0::fgMinM02 = 0.2;
	55	Double_t AliAnalysisTaskSEPHOSpPbPi0::fgMinDistToBad = 2.5;
	56
	57	//________________________________________________________________________
	58	AliAnalysisTaskSEPHOSpPbPi0::AliAnalysisTaskSEPHOSpPbPi0():
0f8c33c1	59	AliAnalysisTaskSE(), fIsMC(kFALSE), fCentralityBin(10), fBufferSize(10), fRunNumber(-1), fPHOSGeo(0),
0f8c33c1	60	fListEvent(0), fListCaloCl(0), fListPi0(0), fListMC(0), fHeader(0), fCaloClArr(0)
a01cd2c9	61	{
	62	//
	63	// Default constructor
	64	//
	65	for (Int_t i=0; i<10; i++) { for (Int_t j=0; j<10; j++) fEventList[i][j] = 0; }
	66
	67	// Set bad channel map
	68	for (Int_t i=0; i<5; i++)
	69	fPHOSBadMap[i] = new TH2I(Form("PHOS_BadMap_mod%d", i), Form("PHOS_BadMap_mod%d", i), 64, 0., 64., 56, 0., 56.);
	70	}
	71	//________________________________________________________________________
	72	AliAnalysisTaskSEPHOSpPbPi0::AliAnalysisTaskSEPHOSpPbPi0(const char *name):
0f8c33c1	73	AliAnalysisTaskSE(name), fIsMC(kFALSE), fCentralityBin(10), fBufferSize(10), fRunNumber(-1), fPHOSGeo(0),
0f8c33c1	74	fListEvent(0), fListCaloCl(0), fListPi0(0), fListMC(0), fHeader(0), fCaloClArr(0)
a01cd2c9	75	{
	76	// Constructor
	77	for (Int_t i=0; i<10; i++) { for (Int_t j=0; j<10; j++) fEventList[i][j] = 0; }
	78
	79	// Set bad channel map
	80	for (Int_t i=0; i<5; i++)
	81	fPHOSBadMap[i] = new TH2I(Form("PHOS_BadMap_mod%d", i), Form("PHOS_BadMap_mod%d", i), 64, 0., 64., 56, 0., 56.);
	82
	83	DefineOutput(1, TList::Class());
0f8c33c1	84	DefineOutput(2, TList::Class());
	85	DefineOutput(3, TList::Class());
	86	DefineOutput(4, TList::Class());
a01cd2c9	87	}
	88
	89	//_____________________________________________________________________________
	90	AliAnalysisTaskSEPHOSpPbPi0::~AliAnalysisTaskSEPHOSpPbPi0()
	91	{
	92	//
	93	// Default destructor
	94	//
0f8c33c1	95	if (fListEvent) { delete fListEvent; fListEvent = NULL; }
	96	if (fListCaloCl) { delete fListCaloCl; fListCaloCl = NULL; }
	97	if (fListPi0) { delete fListPi0; fListPi0 = NULL; }
	98	if (fListMC) { delete fListMC; fListMC = NULL; }
	99	if (fHeader) { delete fHeader; fHeader = NULL; }
	100	if (fCaloClArr) { delete fCaloClArr; fCaloClArr = NULL; }
a01cd2c9	101	}
	102
	103	//________________________________________________________________________
	104	void AliAnalysisTaskSEPHOSpPbPi0::UserCreateOutputObjects()
	105	{
	106	// Create the output container
	107	// Initialize the PHOS geometry
	108	//fPHOSGeo = new AliPHOSGeoUtils("PHOSGeo");
	109	//fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP");
	110
0f8c33c1	111	if (!fHeader) fHeader = new AliPHOSpPbPi0Header();
	112	if (!fCaloClArr) fCaloClArr = new TClonesArray("AliCaloClusterInfo", 0);
	113	if (!fListEvent) fListEvent = new TList();
	114	if (!fListCaloCl) fListCaloCl = new TList();
	115	if (!fListPi0) fListPi0 = new TList();
	116	if (!fListMC) fListMC = new TList();
a01cd2c9	117
a01cd2c9	118	fHeader->SetIsMC(fIsMC);
0f8c33c1	119	fHeader->SetNCent(fCentralityBin.GetSize()-1);
0f8c33c1	120	fHeader->CreateHistograms(fListEvent, fListCaloCl, fListPi0, fListMC);
a01cd2c9	121
a01cd2c9	122	// Post output data.
0f8c33c1	123	PostData(1, fListEvent);
	124	PostData(2, fListCaloCl);
	125	PostData(3, fListPi0);
	126	PostData(4, fListMC);
	127
a01cd2c9	128	return;
	129	}
	130
	131	//________________________________________________________________________
	132	void AliAnalysisTaskSEPHOSpPbPi0::UserExec(Option_t *)
	133	{
	134	// Main loop
	135	// Called for each event
	136
	137	if (fIsMC) {
	138	if (MCEvent()) {
	139	if (MCEvent()->GetNumberOfTracks()<=0)
	140	{ AliError("MC event not found. Nothing done!"); return; }
	141	} else { AliError("MC event not found. Nothing done!"); return; }
	142	}
	143
a01cd2c9	144	AliAODEvent *aod = 0;
	145	AliESDEvent *esd = 0;
	146
	147	if (((TString)fInputEvent->IsA()->GetName())=="AliAODEvent") {
	148	aod = dynamic_cast<AliAODEvent*>(fInputEvent);
	149	if (!aod) { AliError("AOD event not found. Nothing done!"); return; }
	150	if (!fIsMC && (aod->GetHeader()->GetEventType()!=7)) return; // check event type; should be PHYSICS = 7 for data and 0 for MC
a01cd2c9	151	} else {
	152	esd = dynamic_cast<AliESDEvent*>(fInputEvent);
	153	if (!esd) { AliError("ESD event not found. Nothing done!"); return; }
	154	if (!fIsMC && (esd->GetHeader()->GetEventType()!=7)) return; // check event type; should be PHYSICS = 7 for data and 0 for MC
a01cd2c9	155	}
	156
	157	// Fill Event info
	158	fHeader->SetEventInfo(fInputHandler);
0f8c33c1	159	fHeader->FillHistosEvent(fListEvent);
	160
	161	if (!fHeader->IsSelected()) return; // event selection
a01cd2c9	162
	163	// PHOS Geometry and Misalignment initialization at the first time it runs
	164	if(fRunNumber != fInputEvent->GetRunNumber()) {
	165	fRunNumber = fInputEvent->GetRunNumber();
	166	PHOSInitialize(esd);
	167	}
	168
	169	// Fill PHOS cells QA histograms
0f8c33c1	170	fHeader->FillHistosCaloCellsQA(fListCaloCl, fInputEvent->GetPHOSCells(), fPHOSGeo);
0f8c33c1	171	aod = 0;
a01cd2c9	172
a01cd2c9	173	// Fill PHOS cluster Clones Array
0f8c33c1	174	FillCaloClusterInfo(esd);
a01cd2c9	175
a01cd2c9	176	if (!fCaloClArr->GetEntriesFast()) return;
0f8c33c1	177
	178	// vertex bining and centrality bining
	179	Int_t zvtx = (Int_t)((fHeader->Vz() + 10.)/2.); if (zvtx<0) zvtx = 0; if (zvtx>9) zvtx = 9;
	180	Int_t cent = TMath::BinarySearch<Float_t>(fCentralityBin.GetSize()-1, fCentralityBin.GetArray(), fHeader->Centrality());
a01cd2c9	181	if (!fEventList[zvtx][cent]) fEventList[zvtx][cent] = new TList();
	182	TList *eventList = fEventList[zvtx][cent];
	183
	184	// Fill cluster histograms
0f8c33c1	185	fHeader->FillHistosCaloCluster(fListCaloCl, fCaloClArr, cent);
a01cd2c9	186
a01cd2c9	187	// Fill pi0 histograms
0f8c33c1	188	fHeader->FillHistosPi0(fListPi0, fCaloClArr, cent);
a01cd2c9	189
a01cd2c9	190	// Fill mixed pi0 histograms
0f8c33c1	191	fHeader->FillHistosMixPi0(fListPi0, fCaloClArr, eventList, cent);
a01cd2c9	192
a01cd2c9	193	// Fill MC info
0f8c33c1	194	AliStack *stack = 0x0;
a01cd2c9	195	if (fIsMC) {
a01cd2c9	196	if (esd) {
0f8c33c1	197	if (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()) {
0f8c33c1	198	if (static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent())
a01cd2c9	199	stack = static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent()->Stack();
a01cd2c9	200	}
0f8c33c1	201	fHeader->FillHistosMC(fListMC, stack, fPHOSGeo, cent);
a01cd2c9	202	} else
0f8c33c1	203	fHeader->FillHistosMC(fListMC, MCEvent(), fPHOSGeo, cent);
a01cd2c9	204	}
0f8c33c1	205	esd = 0;
a01cd2c9	206
a01cd2c9	207	// Fill event list for mixing
0f8c33c1	208	if (fCaloClArr->GetEntriesFast()>0) {
a01cd2c9	209	eventList->AddFirst(fCaloClArr); fCaloClArr = 0;
0f8c33c1	210	if (eventList->GetSize()>fBufferSize[cent]) { // Remove redundant events
0f8c33c1	211	TClonesArray tmp = static_cast<TClonesArray>(eventList->Last());
a01cd2c9	212	eventList->RemoveLast();
0f8c33c1	213	delete tmp;
a01cd2c9	214	}
	215	}
	216
	217	return;
0f8c33c1	218	}
a01cd2c9	219
	220	//________________________________________________________________________
	221	void AliAnalysisTaskSEPHOSpPbPi0::Terminate(Option_t *)
	222	{
	223	// Terminate analysis
	224
	225	// add the correction matrix
	226
	227	return;
	228	}
	229
	230	//________________________________________________________________________
	231	void AliAnalysisTaskSEPHOSpPbPi0::PHOSInitialize(AliESDEvent* const esd)
	232	{
	233	// Initialize PHOS Geometry and misalignment
	234	// PHOS Geometry
	235	AliOADBContainer geomContainer("phosGeo");
	236	geomContainer.InitFromFile("$ALICE_ROOT/OADB/PHOS/PHOSGeometry.root","PHOSRotationMatrixes");
	237	TObjArray matrixes = (TObjArray)geomContainer.GetObject(fRunNumber,"PHOSRotationMatrixes");
	238	fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP");
	239	for (Int_t mod=0; mod<5; mod++) {
	240	if (!matrixes->At(mod)) continue;
	241	else fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod);
	242	}
	243
	244	// sets the PHOS Misalignment vertex if ESD
	245	if (esd) {
	246	for (Int_t mod=0; mod<5; mod++) {
	247	const TGeoHMatrix* modMatrix = fInputEvent->GetPHOSMatrix(mod);
	248	if (!modMatrix) continue;
	249	else fPHOSGeo->SetMisalMatrix(modMatrix, mod);
	250	}
	251	}
	252	}
	253
	254	//________________________________________________________________________
0f8c33c1	255	void AliAnalysisTaskSEPHOSpPbPi0::FillCaloClusterInfo(AliESDEvent* const esd)
a01cd2c9	256	{
	257	// Fill calo cluster info
	258	if (fCaloClArr) fCaloClArr->Clear();
	259	else fCaloClArr = new TClonesArray("AliCaloClusterInfo", 0);
	260
0f8c33c1	261	Int_t nclsts = fInputEvent->GetNumberOfCaloClusters();
	262	Int_t countN = 0;
	263	Int_t relId[4] = {0,0,0,0}; // module = relId[0]; cellX = relId[2]; cellZ = relId[3];
	264	Float_t position[3] = {0,0,0};
	265	Double_t vtx[3] = {0,0,0}; fHeader->GetXYZ(vtx);
	266	TClonesArray &caloRef = *fCaloClArr;
	267	TLorentzVector momentum;
a01cd2c9	268	AliVCluster *clust = 0;
	269	AliCaloClusterInfo *caloCluster = 0;
	270	for (Int_t iclst=0; iclst<nclsts; iclst++) { // loop over all clusters
	271	clust = fInputEvent->GetCaloCluster(iclst);
	272	if (!(clust && clust->IsPHOS() && clust->E()>fgMinClusterEnergy)) { clust=0; continue; }
	273	if (!(clust->GetNCells()>fgMinNCells && clust->GetM02()>fgMinM02)) { clust=0; continue; } // To remove exotic clusters
	274	if (!(clust->GetDistanceToBadChannel()>fgMinDistToBad)) { clust=0; continue; }
	275	clust->GetPosition(position); TVector3 global(position);
	276	fPHOSGeo->GlobalPos2RelId(global,relId);
	277	if (!IsGoodCaloCluster(relId[0], relId[2], relId[3])) { clust=0; continue; } // calo cluster selection
	278	if (relId[0] == 2) { clust=0; continue; } // !remove module 2
	279
	280	caloCluster = new AliCaloClusterInfo(clust, esd, fPHOSGeo, vtx);
0f8c33c1	281	if (esd && !(caloCluster->LorentzVector().E()>fgMinClusterEnergy)) { delete caloCluster; caloCluster=0; clust=0; continue; } // check again for ESD
0f8c33c1	282	if (caloCluster->GetTrackPt()==-1. \|\| caloCluster->TestCPV(fHeader->MagneticField())) caloCluster->SetPIDBit(BIT(0)); // set CPV bit
a01cd2c9	283
	284	clust = 0;
	285
	286	new(caloRef[countN++]) AliCaloClusterInfo(*caloCluster);
	287
	288	delete caloCluster; caloCluster=0;
	289	} // end loop of all clusters
	290	}
	291
	292	//________________________________________________________________________
	293	Bool_t AliAnalysisTaskSEPHOSpPbPi0::IsGoodCaloCluster(Int_t iMod, Int_t cellX, Int_t cellZ)
	294	{
	295	// !Check whether this cluster is not in bad channel
	296
	297	if (!(iMod>0 && iMod<5 && fPHOSBadMap[iMod])) return kTRUE; //No bad maps for this Module
	298	if (fPHOSBadMap[iMod]->GetBinContent(cellX,cellZ)>0) return kFALSE;
	299
	300	return kTRUE;
	301	}