else if(clus->E() < 6 ) fhDeltaIEtaDeltaIPhiE2[matched]->Fill(dIeta,dIphi);
else fhDeltaIEtaDeltaIPhiE6[matched]->Fill(dIeta,dIphi);
- Float_t dIA = 1.*(dIphi-dIeta)/(dIeta+dIphi);
+ Float_t dIA = 0;
+ if( dIeta+dIphi > 0 ) dIA = 1.*(dIphi-dIeta)/(dIeta+dIphi);
fhDeltaIA[matched]->Fill(clus->E(),dIA);
if(clus->E() > 0.5){
fhNCellsE [index] = 0;
fhTimeE [index] = 0;
fhMaxCellDiffClusterE[index] = 0;
- fhE [index] = 0;
+ fhE [index] = 0;
+ fhPt [index] = 0;
fhPhi [index] = 0;
fhEta [index] = 0;
fhEtaPhi [index] = 0;
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, 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. *
- **************************************************************************/
-/* $Id: $ */
-
-//_________________________________________________________________________
-// Example class on how to read CaloClusters, and AODTracks and how
-// fill AODs with PWG4PartCorr analysis frame
-// Select the type of detector information that you want to analyze, CTS (tracking), PHOS or EMCAL
-// Select the PID custer type of the calorimeters
-// Set min momentum of the cluster/tracks
-// Fill few histograms
-//
-//-- Author: Gustavo Conesa (INFN-LNF)
-//_________________________________________________________________________
-
-
-// --- ROOT system ---
-//#include "Riostream.h"
-#include "TObjArray.h"
-#include "TParticle.h"
-#include "TCanvas.h"
-#include "TPad.h"
-#include "TROOT.h"
-#include "TH2F.h"
-
-//---- AliRoot system ----
-#include "AliAnaExample.h"
-#include "AliCaloTrackReader.h"
-#include "AliAODPWG4Particle.h"
-#include "AliStack.h"
-#include "AliCaloPID.h"
-#include "AliFiducialCut.h"
-#include "AliVCluster.h"
-#include "AliVTrack.h"
-
-ClassImp(AliAnaExample)
-
-//____________________________________________________________________________
- AliAnaExample::AliAnaExample() :
- AliAnaPartCorrBaseClass(),fPdg(0), fDetector(""), fhPt(0),fhPhi(0),fhEta(0), fh2Pt(0),fh2Phi(0),fh2Eta(0)
-{
- //Default Ctor
-
- //Initialize parameters
- InitParameters();
-}
-
-//________________________________________________________________________
-TList * AliAnaExample::GetCreateOutputObjects()
-{
- // Create histograms to be saved in output file and
- // store them in fOutputContainer
-
- TList * outputContainer = new TList() ;
- outputContainer->SetName("ExampleHistos") ;
-
- Int_t nptbins = GetHistoPtBins();
- Int_t nphibins = GetHistoPhiBins();
- Int_t netabins = GetHistoEtaBins();
- Float_t ptmax = GetHistoPtMax();
- Float_t phimax = GetHistoPhiMax();
- Float_t etamax = GetHistoEtaMax();
- Float_t ptmin = GetHistoPtMin();
- Float_t phimin = GetHistoPhiMin();
- Float_t etamin = GetHistoEtaMin();
-
- fhPt = new TH1F ("hPt","p_T distribution", nptbins,ptmin,ptmax);
- fhPt->SetXTitle("p_{T} (GeV/c)");
- outputContainer->Add(fhPt);
-
- fhPhi = new TH1F ("hPhi","#phi distribution",nphibins,phimin,phimax);
- fhPhi->SetXTitle("#phi (rad)");
- outputContainer->Add(fhPhi);
-
- fhEta = new TH1F ("hEta","#eta distribution",netabins,etamin,etamax);
- fhEta->SetXTitle("#eta ");
- outputContainer->Add(fhEta);
-
- if(IsDataMC()){
- fh2Pt = new TH2F ("h2Pt","p_T distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
- fh2Pt->SetXTitle("p_{T,rec} (GeV/c)");
- fh2Pt->SetYTitle("p_{T,gen} (GeV/c)");
- outputContainer->Add(fh2Pt);
-
- fh2Phi = new TH2F ("h2Phi","#phi distribution, reconstructed vs generated", nphibins,phimin,phimax, nphibins,phimin,phimax);
- fh2Phi->SetXTitle("#phi_{rec} (rad)");
- fh2Phi->SetYTitle("#phi_{gen} (rad)");
- outputContainer->Add(fh2Phi);
-
- fh2Eta = new TH2F ("h2Eta","#eta distribution, reconstructed vs generated", netabins,etamin,etamax,netabins,etamin,etamax);
- fh2Eta->SetXTitle("#eta_{rec} ");
- fh2Eta->SetYTitle("#eta_{gen} ");
- outputContainer->Add(fh2Eta);
-
- }
-
- return outputContainer;
-}
-
-//__________________________________________________
-void AliAnaExample::InitParameters()
-{
- //Initialize the parameters of the analysis.
- SetOutputAODClassName("AliAODPWG4Particle");
- SetOutputAODName("PWG4Particle");
- AddToHistogramsName("AnaExample_");
-
- fPdg = 22; //Keep photons
- fDetector = "PHOS";
-
-}
-
-//__________________________________________________________________
-void AliAnaExample::Print(const Option_t * opt) const
-{
- //Print some relevant parameters set for the analysis
- if(! opt)
- return;
-
- printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
- AliAnaPartCorrBaseClass::Print(" ");
-
- printf("Select clusters with pdg %d \n",fPdg);
- printf("Select detector %s \n",fDetector.Data());
-
-}
-
-//__________________________________________________________________
-void AliAnaExample::MakeAnalysisFillAOD()
-{
- //Do analysis and fill aods
-
- //Some prints
- if(GetDebug() > 0){
- if(fDetector == "EMCAL" && GetEMCALClusters())printf("AliAnaExample::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetEMCALClusters()->GetEntriesFast());
- if(fDetector == "CTS" && GetCTSTracks()) printf("AliAnaExample::MakeAnalysisFillAOD() - In CTS aod entries %d\n", GetCTSTracks() ->GetEntriesFast());
- if(fDetector == "PHOS" && GetPHOSClusters()) printf("AliAnaExample::MakeAnalysisFillAOD() - In PHOS aod entries %d\n", GetPHOSClusters() ->GetEntriesFast());
- }
-
- //Get List with tracks or clusters
- TObjArray * partList = 0x0;
- if(fDetector == "CTS") partList = GetCTSTracks();
- else if(fDetector == "EMCAL") partList = GetEMCALClusters();
- else if(fDetector == "PHOS") partList = GetPHOSClusters();
-
- if(!partList || partList->GetEntriesFast() == 0) return ;
-
- //Fill AODParticle with PHOS/EMCAL aods
- if(fDetector == "EMCAL" || fDetector == "PHOS"){
-
- //Get vertex for photon momentum calculation
-
-
- TLorentzVector mom ;
- for(Int_t i = 0; i < partList->GetEntriesFast(); i++){
-
- AliVCluster * calo = (AliVCluster*) (partList->At(i));
-
- //Get the index where the cluster comes, to retrieve the corresponding vertex
- Int_t evtIndex = 0 ;
- if (GetMixedEvent()) {
- evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
- }
-
- //Fill AODParticle after some selection
- if(GetReader()->GetDataType() != AliCaloTrackReader::kMC){
- calo->GetMomentum(mom,GetVertex(evtIndex)) ;}//Assume that come from vertex in straight line
- else{
- Double_t vertex[]={0,0,0};
- calo->GetMomentum(mom,vertex) ;
- }
-
- Int_t pdg = fPdg;
-
- if(IsCaloPIDOn()){
- const Double_t *pid = calo->GetPID();
- pdg = GetCaloPID()->GetIdentifiedParticleType(fDetector,pid,mom.E());
- //pdg = GetCaloPID()->GetIdentifiedParticleType(fDetector,mom,
- // calo->GetM02(), calo->GetM02(),
- // calo->GetDispersion(), 0, 0);
- }
-
- //Acceptance selection
- Bool_t in = kTRUE;
- if(IsFiducialCutOn())
- in = GetFiducialCut()->IsInFiducialCut(mom,fDetector) ;
-
- if(GetDebug() > 1) printf("AliAnaExample::MakeAnalysisFillAOD() - Cluster pt %2.2f, phi %2.2f, pdg %d, in fiducial cut %d \n",mom.Pt(), mom.Phi(), pdg, in);
-
- //Select cluster if momentum, pdg and acceptance are good
- if(mom.Pt() > GetMinPt() && pdg ==fPdg && in) {
- AliAODPWG4Particle ph = AliAODPWG4Particle(mom);
- //AddAODParticleCorrelation(AliAODPWG4ParticleCorrelation(mom));
- ph.SetLabel(calo->GetLabel());
- ph.SetIdentifiedParticleType(pdg);
- ph.SetDetector(fDetector);
- AddAODParticle(ph);
- }//selection
- }//loop
- }//Calorimeters
- else if(fDetector == "CTS"){ //Track analysis
- //Fill AODParticle with CTS aods
- TVector3 p3;
- for(Int_t i = 0; i < GetCTSTracks()->GetEntriesFast(); i++){
-
- AliVTrack * track = (AliVTrack*) (GetCTSTracks()->At(i));
-
- //Fill AODParticle after some selection
- Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
- p3.SetXYZ(mom[0],mom[1],mom[2]);
-
- //Acceptance selection
- Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,"CTS") ;
- if(GetDebug() > 1) printf("AliAnaExample::MakeAnalysisFillAOD() - Track pt %2.2f, phi %2.2f, in fiducial cut %d\n",p3.Pt(), p3.Phi(), in);
- if(p3.Pt() > GetMinPt() && in) {
- AliAODPWG4Particle tr = AliAODPWG4Particle(mom[0],mom[1],mom[2],0);
- tr.SetDetector("CTS");
- AddAODParticle(tr);
- }//selection
- }//loop
- }//CTS analysis
-
- if(GetDebug() > 0) {
- if(fDetector!="CTS" && GetReader()->GetDataType()!= AliCaloTrackReader::kMC)
- printf("AliAnaExample::MakeAnalysisFillAOD() - Final aod branch entries %d\n", GetOutputAODBranch()->GetEntriesFast());
- // if(fDetector!="CTS" && GetReader()->GetDataType()!= AliCaloTrackReader::kMC)
- //printf("Example: final aod cell entries %d\n", GetAODCaloCells()->GetNumberOfCells());
- }
-}
-
-//__________________________________________________________________
-void AliAnaExample::MakeAnalysisFillHistograms()
-{
- //Do analysis and fill histograms
-
- //Loop on stored AODParticles
- Int_t naod = GetOutputAODBranch()->GetEntriesFast();
- if(GetDebug() > 0) printf("AliAnaExample::MakeAnalysisFillHistograms() - Histo aod branch entries %d\n", naod);
- for(Int_t iaod = 0; iaod < naod ; iaod++){
- AliAODPWG4Particle* ph = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
-
- fhPt->Fill(ph->Pt());
- fhPhi->Fill(ph->Phi());
- fhEta->Fill(ph->Eta());
-
- if(IsDataMC()){
- //Play with the MC stack if available
- Float_t ptprim = 0;
- Float_t phiprim = 0;
- Float_t etaprim = 0;
- if(GetReader()->ReadStack()){
- AliStack * stack = GetMCStack() ;
-
- if(ph->GetLabel() < 0 || !stack) {
- printf("AliAnaExample::MakeAnalysisFillHistograms() *** bad label or no stack ***: label %d \n", ph->GetLabel());
- continue;
- }
-
- if(ph->GetLabel() >= stack->GetNtrack()) {
- printf("AliAnaExample::MakeAnalysisFillHistograms() *** large label ***: label %d, n tracks %d \n", ph->GetLabel(), stack->GetNtrack());
- continue ;
- }
-
- TParticle * mom = GetMCStack()->Particle(ph->GetLabel());
- ptprim = mom->Pt();
- phiprim = mom->Phi();
- etaprim = mom->Eta();
- }
- else if(GetReader()->ReadAODMCParticles()){
- printf("AliAnaExample::MakeAnalysisFillHistograms() - Acceptance calculation with MCParticles not implemented yet\n");
- }
- fh2Pt->Fill(ph->Pt(), ptprim);
- fh2Phi->Fill(ph->Phi(), phiprim);
- fh2Eta->Fill(ph->Eta(), etaprim);
- }//Work with stack also
- }// aod branch loop
-
- }
-
-//________________________________________________________________________
-void AliAnaExample::ReadHistograms(TList* outputList)
-{
- // Needed when Terminate is executed in distributed environment
- // Refill analysis histograms of this class with corresponding histograms in output list.
-
- // Histograms of this analsys are kept in the same list as other analysis, recover the position of
- // the first one and then add the next
- Int_t index = outputList->IndexOf(outputList->FindObject(GetAddedHistogramsStringToName()+"hPt"));
-
- //Read histograms, must be in the same order as in GetCreateOutputObject.
- fhPt = (TH1F *) outputList->At(index++);
- fhPhi = (TH1F *) outputList->At(index++);
- fhEta = (TH1F *) outputList->At(index++);
-
- if(IsDataMC()){
- fh2Pt = (TH2F *) outputList->At(index++);
- fh2Phi = (TH2F *) outputList->At(index++);
- fh2Eta = (TH2F *) outputList->At(index);
- }
-
-}
-
-//__________________________________________________________________
-void AliAnaExample::Terminate(TList* outputList)
-{
-
- //Do some plots to end
-
- //Recover histograms from output histograms list, needed for distributed analysis.
- ReadHistograms(outputList);
-
- printf(" AliAnaExample::Terminate() *** %s Report:", GetName()) ;
- printf(" AliAnaExample::Terminate() pt : %5.3f , RMS : %5.3f \n", fhPt->GetMean(), fhPt->GetRMS() ) ;
-
- TCanvas * c = new TCanvas("c", "PHOS ESD Test", 400, 10, 600, 700) ;
- c->Divide(1, 3);
-
- c->cd(1) ;
- gPad->SetLogy();
- fhPt->SetLineColor(2);
- fhPt->Draw();
-
- c->cd(2) ;
- fhPhi->SetLineColor(2);
- fhPhi->Draw();
-
- c->cd(3) ;
- fhEta->SetLineColor(2);
- fhEta->Draw();
-
- c->Print("Example.eps");
-
- const Int_t buffersize = 1024;
- char line[buffersize] ;
- snprintf(line,buffersize, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
- gROOT->ProcessLine(line);
- snprintf(line,buffersize, ".!rm -fR *.eps");
- gROOT->ProcessLine(line);
-
- printf("AliAnaExample::Terminate() - !! All the eps files are in %s.tar.gz !!!", GetName());
-
-}
+++ /dev/null
-#ifndef ALIANAEXAMPLE_H
-#define ALIANAEXAMPLE_H
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice */
-/* $Id: $ */
-
-//_________________________________________________________________________
-// Example class on how to read AODCaloClusters, ESDCaloCells and AODTracks and how
-// fill AODs with PWG4PartCorr analysis frame
-// Select the type of detector information that you want to analyze, CTS (tracking), PHOS or EMCAL
-// Select the PID custer type of the calorimeters
-// Set min momentum of the cluster/tracks
-// Fill few histograms
-//
-//-- Author: Gustavo Conesa (INFN-LNF)
-
-// --- Root system ---
-class TH2F;
-class TH1F;
-
-// --- Analysis system ---
-#include "AliAnaPartCorrBaseClass.h"
-
-class AliAnaExample : public AliAnaPartCorrBaseClass {
-
- public:
- AliAnaExample() ; // default ctor
- virtual ~AliAnaExample() {;} //virtual dtor
-
- private :
- AliAnaExample(const AliAnaExample & g) ; // cpy ctor
- AliAnaExample & operator = (const AliAnaExample & g) ;//cpy assignment
-
- public:
-
- TList * GetCreateOutputObjects();
-
- void InitParameters();
-
- void Print(const Option_t * opt) const;
-
- void MakeAnalysisFillAOD() ;
-
- void MakeAnalysisFillHistograms() ;
-
- Int_t GetPdg() const {return fPdg ;}
- void SetPdg( Int_t pdg ) {fPdg = pdg; }
-
- TString GetDetector() const {return fDetector ;}
- void SetDetector( TString calo ) {fDetector = calo; }
-
- void Terminate(TList * outputList);
- void ReadHistograms(TList * outputList); //Fill histograms with histograms in ouput list, needed in Terminate.
-
- private:
-
- Int_t fPdg ; //identified particle id
- TString fDetector ; //detector selection
- //Histograms
- //CaloClusters
- TH1F * fhPt; //! pT distribution
- TH1F * fhPhi; //! phi distribution
- TH1F * fhEta; //! eta distribution
- TH2F * fh2Pt; //!pT distribution, reconstructed vs generated
- TH2F * fh2Phi; //! phi distribution, reconstructed vs generated
- TH2F * fh2Eta; //! eta distribution, reconstructed vs generated
-
- ClassDef(AliAnaExample,2)
-} ;
-
-
-#endif //ALIANAEXAMPLE_H
-
-
-
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-1999, 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. *
+ **************************************************************************/
+
+//_________________________________________________________________________
+//
+// Split clusters with some criteria and calculate invariant mass
+// to identify them as pi0 or conversion
+//
+//
+//-- Author: Gustavo Conesa (LPSC-Grenoble)
+//_________________________________________________________________________
+
+//////////////////////////////////////////////////////////////////////////////
+
+
+// --- ROOT system ---
+#include <TList.h>
+#include <TClonesArray.h>
+#include <TObjString.h>
+#include <TH3F.h>
+//#include "Riostream.h"
+
+// --- Analysis system ---
+#include "AliAnaInsideClusterInvariantMass.h"
+#include "AliCaloTrackReader.h"
+#include "AliMCAnalysisUtils.h"
+#include "AliStack.h"
+#include "AliFiducialCut.h"
+#include "TParticle.h"
+#include "AliVCluster.h"
+#include "AliAODEvent.h"
+#include "AliAODMCParticle.h"
+#include "AliEMCALGeoParams.h"
+
+ClassImp(AliAnaInsideClusterInvariantMass)
+
+//__________________________________________________________________
+AliAnaInsideClusterInvariantMass::AliAnaInsideClusterInvariantMass() :
+ AliAnaPartCorrBaseClass(),
+ fCalorimeter(""),
+ fM02Cut(0),
+ fMinNCells(0)
+{
+ //default ctor
+
+ // Init array of histograms
+ for(Int_t i = 0; i < 7; i++){
+ fhMassNLocMax1[i] = 0;
+ fhMassNLocMax2[i] = 0;
+ fhMassNLocMaxN[i] = 0;
+ fhNLocMax[i] = 0;
+ fhNLocMaxM02Cut[i] = 0;
+ fhM02NLocMax1[i] = 0;
+ fhM02NLocMax2[i] = 0;
+ fhM02NLocMaxN[i] = 0;
+ fhNCellNLocMax1[i] = 0;
+ fhNCellNLocMax2[i] = 0;
+ fhNCellNLocMaxN[i] = 0;
+ fhM02Pi0[i] = 0;
+ fhM02Eta[i] = 0;
+ fhM02Con[i] = 0;
+ fhInvMassAllCells[i]=0;
+ }
+
+ InitParameters();
+
+}
+
+//_______________________________________________________________
+TObjString * AliAnaInsideClusterInvariantMass::GetAnalysisCuts()
+{
+ //Save parameters used for analysis
+ TString parList ; //this will be list of parameters used for this analysis.
+ const Int_t buffersize = 255;
+ char onePar[buffersize] ;
+
+ snprintf(onePar,buffersize,"--- AliAnaInsideClusterInvariantMass ---\n") ;
+ parList+=onePar ;
+
+ snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"fM02Cut =%f \n" ,fM02Cut) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"fMinNCells =%f \n",fMinNCells) ;
+ parList+=onePar ;
+
+ return new TObjString(parList) ;
+
+}
+
+//____________________________________________________________________________________________________
+Bool_t AliAnaInsideClusterInvariantMass::AreNeighbours( const Int_t absId1, const Int_t absId2 ) const
+{
+ // Tells if (true) or not (false) two digits are neighbours
+ // A neighbour is defined as being two digits which share a corner
+
+ Bool_t areNeighbours = kFALSE ;
+ Int_t nSupMod =0, nModule =0, nIphi =0, nIeta =0;
+ Int_t nSupMod1=0, nModule1=0, nIphi1=0, nIeta1=0;
+ Int_t relid1[2], relid2[2] ;
+ Int_t rowdiff=0, coldiff=0;
+
+ areNeighbours = kFALSE ;
+
+ GetEMCALGeometry()->GetCellIndex(absId1, nSupMod,nModule,nIphi,nIeta);
+ GetEMCALGeometry()->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, relid1[0],relid1[1]);
+
+ GetEMCALGeometry()->GetCellIndex(absId2, nSupMod1,nModule1,nIphi1,nIeta1);
+ GetEMCALGeometry()->GetCellPhiEtaIndexInSModule(nSupMod1,nModule1,nIphi1,nIeta1, relid2[0],relid2[1]);
+
+ // In case of a shared cluster, index of SM in C side, columns start at 48 and ends at 48*2-1
+ // C Side impair SM, nSupMod%2=1; A side pair SM nSupMod%2=0
+ if(nSupMod1!=nSupMod){
+ if(nSupMod1%2) relid1[1]+=AliEMCALGeoParams::fgkEMCALCols;
+ else relid2[1]+=AliEMCALGeoParams::fgkEMCALCols;
+ }
+
+ rowdiff = TMath::Abs( relid1[0] - relid2[0] ) ;
+ coldiff = TMath::Abs( relid1[1] - relid2[1] ) ;
+
+ if (( coldiff <= 1 ) && ( rowdiff <= 1 ) && (coldiff + rowdiff > 0))
+ areNeighbours = kTRUE ;
+
+ return areNeighbours;
+}
+
+//_____________________________________________________________________________________
+TLorentzVector AliAnaInsideClusterInvariantMass::GetCellMomentum(const Int_t absId,
+ AliVCaloCells * cells)
+{
+
+ // Cell momentum calculation for invariant mass
+
+ Double_t cellpos[] = {0, 0, 0};
+ GetEMCALGeometry()->GetGlobal(absId, cellpos);
+
+ if(GetVertex(0)){//calculate direction from vertex
+ cellpos[0]-=GetVertex(0)[0];
+ cellpos[1]-=GetVertex(0)[1];
+ cellpos[2]-=GetVertex(0)[2];
+ }
+
+ Double_t r = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]+cellpos[2]*cellpos[2] ) ;
+
+ Float_t en = cells->GetCellAmplitude(absId);
+ RecalibrateCellAmplitude(en,absId);
+ TLorentzVector cellMom ;
+ cellMom.SetPxPyPzE( en*cellpos[0]/r, en*cellpos[1]/r, en*cellpos[2]/r, en) ;
+
+ return cellMom;
+
+}
+
+//________________________________________________________________
+TList * AliAnaInsideClusterInvariantMass::GetCreateOutputObjects()
+{
+ // Create histograms to be saved in output file and
+ // store them in outputContainer
+ TList * outputContainer = new TList() ;
+ outputContainer->SetName("InsideClusterHistos") ;
+
+ Int_t nptbins = GetHistoPtBins(); Float_t ptmax = GetHistoPtMax(); Float_t ptmin = GetHistoPtMin();
+ Int_t ssbins = GetHistoShowerShapeBins(); Float_t ssmax = GetHistoShowerShapeMax(); Float_t ssmin = GetHistoShowerShapeMin();
+ Int_t mbins = GetHistoMassBins(); Float_t mmax = GetHistoMassMax(); Float_t mmin = GetHistoMassMin();
+ Int_t ncbins = GetHistoNClusterCellBins(); Int_t ncmax = GetHistoNClusterCellMax(); Int_t ncmin = GetHistoNClusterCellMin();
+
+ TString ptype[] ={"","#gamma","#gamma->e^{#pm}","#pi^{0}","#eta","e^{#pm}", "hadron"};
+ TString pname[] ={"","Photon","Conversion", "Pi0", "Eta", "Electron","Hadron"};
+
+ Int_t n = 1;
+
+ if(IsDataMC()) n = 7;
+
+ for(Int_t i = 0; i < n; i++){
+
+ fhInvMassAllCells[i] = new TH2F(Form("hInvMassAllCells%s",pname[i].Data()),
+ Form("Invariant mass of all cells %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,mbins,mmin,mmax);
+ fhInvMassAllCells[i]->SetYTitle("M (MeV/c^2)");
+ fhInvMassAllCells[i]->SetXTitle("E (GeV)");
+ outputContainer->Add(fhInvMassAllCells[i]) ;
+
+
+ fhMassNLocMax1[i] = new TH2F(Form("hMassNLocMax1%s",pname[i].Data()),
+ Form("Invariant mass of 2 highest energy cells %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,mbins,mmin,mmax);
+ fhMassNLocMax1[i]->SetYTitle("M (MeV/c^2)");
+ fhMassNLocMax1[i]->SetXTitle("E (GeV)");
+ outputContainer->Add(fhMassNLocMax1[i]) ;
+
+ fhMassNLocMax2[i] = new TH2F(Form("hMassNLocMax2%s",pname[i].Data()),
+ Form("Invariant mass of 2 local maxima cells %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,mbins,mmin,mmax);
+ fhMassNLocMax2[i]->SetYTitle("M (MeV/c^2)");
+ fhMassNLocMax2[i]->SetXTitle("E (GeV)");
+ outputContainer->Add(fhMassNLocMax2[i]) ;
+
+ fhMassNLocMaxN[i] = new TH2F(Form("hMassNLocMaxN%s",pname[i].Data()),
+ Form("Invariant mass of 2 local maxima cells %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,mbins,mmin,mmax);
+ fhMassNLocMaxN[i]->SetYTitle("M (MeV/c^2)");
+ fhMassNLocMax2[i]->SetXTitle("E (GeV)");
+ outputContainer->Add(fhMassNLocMaxN[i]) ;
+
+
+ fhNLocMax[i] = new TH2F(Form("hNLocMax%s",pname[i].Data()),
+ Form("Number of local maxima in cluster %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,5,0,5);
+ fhNLocMax[i] ->SetYTitle("N maxima");
+ fhNLocMax[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhNLocMax[i]) ;
+
+ fhNLocMaxM02Cut[i] = new TH2F(Form("hNLocMaxM02Cut%s",pname[i].Data()),
+ Form("Number of local maxima in cluster %s for M02 > %2.2f",ptype[i].Data(),fM02Cut),
+ nptbins,ptmin,ptmax,5,0,5);
+ fhNLocMaxM02Cut[i]->SetYTitle("N maxima");
+ fhNLocMaxM02Cut[i]->SetXTitle("E (GeV)");
+ outputContainer->Add(fhNLocMaxM02Cut[i]) ;
+
+
+ fhM02NLocMax1[i] = new TH2F(Form("hM02NLocMax1%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for N max = 1 %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
+ fhM02NLocMax1[i] ->SetYTitle("#lambda_{0}^{2}");
+ fhM02NLocMax1[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhM02NLocMax1[i]) ;
+
+ fhM02NLocMax2[i] = new TH2F(Form("hM02NLocMax2%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for N max = 2 %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
+ fhM02NLocMax2[i] ->SetYTitle("#lambda_{0}^{2}");
+ fhM02NLocMax2[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhM02NLocMax2[i]) ;
+
+
+ fhM02NLocMaxN[i] = new TH2F(Form("hM02NLocMaxN%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for N max > 2 %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
+ fhM02NLocMaxN[i] ->SetYTitle("#lambda_{0}^{2}");
+ fhM02NLocMaxN[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhM02NLocMaxN[i]) ;
+
+
+ fhNCellNLocMax1[i] = new TH2F(Form("hNCellNLocMax1%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for N max = 1 %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ncbins,ncmin,ncmax);
+ fhNCellNLocMax1[i] ->SetYTitle("N cells");
+ fhNCellNLocMax1[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhNCellNLocMax1[i]) ;
+
+ fhNCellNLocMax2[i] = new TH2F(Form("hNCellNLocMax2%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for N max = 2 %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ncbins,ncmin,ncmax);
+ fhNCellNLocMax2[i] ->SetYTitle("N cells");
+ fhNCellNLocMax2[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhNCellNLocMax2[i]) ;
+
+
+ fhNCellNLocMaxN[i] = new TH2F(Form("hNCellNLocMaxN%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for N max > 2 %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ncbins,ncmin,ncmax);
+ fhNCellNLocMaxN[i] ->SetYTitle("N cells");
+ fhNCellNLocMaxN[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhNCellNLocMaxN[i]) ;
+
+
+ fhM02Pi0[i] = new TH2F(Form("hM02Pi0%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E ffor mass [0.1-0.2] MeV/c2 %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
+ fhM02Pi0[i] ->SetYTitle("#lambda_{0}^{2}");
+ fhM02Pi0[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhM02Pi0[i]) ;
+
+ fhM02Eta[i] = new TH2F(Form("hM02Eta%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for mass [0.4-0.7] MeV/c2, %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
+ fhM02Eta[i] ->SetYTitle("#lambda_{0}^{2}");
+ fhM02Eta[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhM02Eta[i]) ;
+
+
+ fhM02Con[i] = new TH2F(Form("hM02Con%s",pname[i].Data()),
+ Form("#lambda_{0}^{2} vs E for mass < 0.5 MeV/c2, %s",ptype[i].Data()),
+ nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
+ fhM02Con[i] ->SetYTitle("#lambda_{0}^{2}");
+ fhM02Con[i] ->SetXTitle("E (GeV)");
+ outputContainer->Add(fhM02Con[i]) ;
+
+ }
+
+ return outputContainer ;
+
+}
+
+//________________________________________________________________________________________________________
+Int_t AliAnaInsideClusterInvariantMass::GetNumberOfLocalMaxima(AliVCluster* cluster, AliVCaloCells* cells,
+ Int_t *absIdList, Float_t *maxEList)
+{
+ // Find local maxima in cluster
+
+ Float_t locMaxCut = 0; // not used for the moment
+
+ Int_t iDigitN = 0 ;
+ Int_t iDigit = 0 ;
+ Int_t absId1 = -1 ;
+ Int_t absId2 = -1 ;
+ const Int_t nCells = cluster->GetNCells();
+ //printf("cluster :");
+ for(iDigit = 0; iDigit < nCells ; iDigit++){
+ absIdList[iDigit] = cluster->GetCellsAbsId()[iDigit] ;
+
+ //Float_t en = cells->GetCellAmplitude(absIdList[iDigit]);
+ //RecalibrateCellAmplitude(en,absIdList[iDigit]);
+ //Int_t icol = -1, irow = -1, iRCU = -1;
+ //Int_t sm = GetCaloUtils()->GetModuleNumberCellIndexes(absIdList[iDigit], fCalorimeter, icol, irow, iRCU) ;
+
+ //printf("\t cell %d, id %d, sm %d, col %d, row %d, e %f\n", iDigit, absIdList[iDigit], sm, icol, irow, en );
+ }
+
+ for(iDigit = 0 ; iDigit < nCells; iDigit++) {
+ if(absIdList[iDigit]>=0) {
+
+ absId1 = absIdList[iDigit] ;
+ Float_t en1 = cells->GetCellAmplitude(absId1);
+ RecalibrateCellAmplitude(en1,absId1);
+
+ for(iDigitN = 0; iDigitN < nCells; iDigitN++) {
+ absId2 = absIdList[iDigitN] ;
+
+ Float_t en2 = cells->GetCellAmplitude(absId2);
+ RecalibrateCellAmplitude(en2,absId2);
+
+ if ( AreNeighbours(absId1, absId2) ) {
+
+ if (en1 > en2 ) {
+ absIdList[iDigitN] = -1 ;
+ // but may be digit too is not local max ?
+ if(en1 < en2 + locMaxCut)
+ absIdList[iDigit] = -1 ;
+ }
+ else {
+ absIdList[iDigit] = -1 ;
+ // but may be digitN too is not local max ?
+ if(en1 > en2 - locMaxCut)
+ absIdList[iDigitN] = -1 ;
+ }
+ } // if Areneighbours
+ } // while digitN
+ } // slot not empty
+ } // while digit
+
+ iDigitN = 0 ;
+ for(iDigit = 0; iDigit < nCells; iDigit++) {
+ if(absIdList[iDigit]>=0 ){
+ absIdList[iDigitN] = absIdList[iDigit] ;
+ Float_t en = cells->GetCellAmplitude(absIdList[iDigit]);
+ RecalibrateCellAmplitude(en,absIdList[iDigit]);
+ if(en < 0.1) continue; // Maxima only with seed energy at least
+ maxEList[iDigitN] = en ;
+ //printf("Local max %d, id %d, en %f\n", iDigit,absIdList[iDigitN],en);
+ iDigitN++ ;
+ }
+ }
+
+ return iDigitN ;
+
+}
+
+
+//___________________________________________
+void AliAnaInsideClusterInvariantMass::Init()
+{
+ //Init
+ //Do some checks
+ if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
+ printf("AliAnaInsideClusterInvariantMass::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
+ abort();
+ }
+ else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){
+ printf("AliAnaInsideClusterInvariantMass::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
+ abort();
+ }
+
+ if( GetReader()->GetDataType() == AliCaloTrackReader::kMC ){
+ printf("AliAnaInsideClusterInvariantMass::Init() - !!STOP: You want to use pure MC data!!\n");
+ abort();
+
+ }
+
+}
+
+//_____________________________________________________
+void AliAnaInsideClusterInvariantMass::InitParameters()
+{
+ //Initialize the parameters of the analysis.
+ AddToHistogramsName("AnaPi0InsideClusterInvariantMass_");
+
+ fCalorimeter = "EMCAL" ;
+ fM02Cut = 0.26 ;
+ fMinNCells = 4 ;
+}
+
+
+//__________________________________________________________________
+void AliAnaInsideClusterInvariantMass::MakeAnalysisFillHistograms()
+{
+ //Search for pi0 in fCalorimeter with shower shape analysis
+
+ TObjArray * pl = 0x0;
+ AliVCaloCells* cells = 0x0;
+
+ //Select the Calorimeter of the photon
+ if(fCalorimeter == "PHOS"){
+ pl = GetPHOSClusters();
+ cells = GetPHOSCells();
+ }
+ else if (fCalorimeter == "EMCAL"){
+ pl = GetEMCALClusters();
+ cells = GetEMCALCells();
+ }
+
+ if(!pl || !cells) {
+ Info("MakeAnalysisFillHistograms","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
+ return;
+ }
+
+ if(fCalorimeter == "PHOS") return; // Not implemented for PHOS yet
+
+ for(Int_t icluster = 0; icluster < pl->GetEntriesFast(); icluster++){
+ AliVCluster * cluster = (AliVCluster*) (pl->At(icluster));
+
+ // Study clusters with large shape parameter
+ Float_t en = cluster->E();
+ Float_t l0 = cluster->GetM02();
+ Int_t nc = cluster->GetNCells();
+
+ //If too small or big E or low number of cells, skip it
+ if( ( en < GetMinEnergy() || en > GetMaxEnergy() ) && nc < fMinNCells) continue ;
+
+ Int_t absId1 = -1; Int_t absId2 = -1;
+ Int_t *absIdList = new Int_t [nc];
+ Float_t *maxEList = new Float_t[nc];
+ Int_t nMax = GetNumberOfLocalMaxima(cluster, cells, absIdList, maxEList) ;
+
+ if (nMax <= 0) {
+ printf("AliAnaInsideClusterInvariantMass::MakeAnalysisFillHistograms() - No local maximum found!");
+ delete absIdList ;
+ delete maxEList ;
+ return;
+ }
+
+ fhNLocMax[0]->Fill(en,nMax);
+
+ if ( nMax == 1 ) { fhM02NLocMax1[0]->Fill(en,l0) ; fhNCellNLocMax1[0]->Fill(en,nc) ; }
+ else if( nMax == 2 ) { fhM02NLocMax2[0]->Fill(en,l0) ; fhNCellNLocMax2[0]->Fill(en,nc) ; }
+ else if( nMax >= 3 ) { fhM02NLocMaxN[0]->Fill(en,l0) ; fhNCellNLocMaxN[0]->Fill(en,nc) ; }
+ else printf("N max smaller than 1 -> %d \n",nMax);
+
+ // Play with the MC stack if available
+ // Check origin of the candidates
+ Int_t mcindex = -1;
+ if(IsDataMC()){
+
+ Int_t tag = GetMCAnalysisUtils()->CheckOrigin(cluster->GetLabels(),cluster->GetNLabels(), GetReader(), 0);
+
+ if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ) mcindex = mcPi0;
+ else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) ) mcindex = mcEta;
+ else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) &&
+ !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) mcindex = mcPhoton;
+ else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) &&
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) mcindex = mcConversion;
+ else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) mcindex = mcElectron;
+ else mcindex = mcHadron;
+
+/* printf("AliAnaInsideClusterInvariantMass::FillAnalysisMakeHistograms() - tag %d, photon %d, prompt %d, frag %d, isr %d, pi0 decay %d, eta decay %d, other decay %d \n conv %d, pi0 %d, hadron %d, electron %d, unk %d, muon %d,pion %d, proton %d, neutron %d, kaon %d, antiproton %d, antineutron %d, bad %d\n",tag,
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay),
+
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOther),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCUnknown),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCMuon),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCProton),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCKaon),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron),
+ GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCBadLabel)
+);
+*/
+
+ fhNLocMax[mcindex]->Fill(en,nMax);
+
+ if (nMax == 1 ) { fhM02NLocMax1[mcindex]->Fill(en,l0) ; fhNCellNLocMax1[mcindex]->Fill(en,nc) ; }
+ else if(nMax == 2 ) { fhM02NLocMax2[mcindex]->Fill(en,l0) ; fhNCellNLocMax2[mcindex]->Fill(en,nc) ; }
+ else if(nMax >= 3 ) { fhM02NLocMaxN[mcindex]->Fill(en,l0) ; fhNCellNLocMaxN[mcindex]->Fill(en,nc) ; }
+
+ }
+
+ if( l0 < fM02Cut) continue;
+
+ // Get the 2 max indeces and do inv mass
+
+ absId1 = absIdList[0];
+ TLorentzVector cellMomi = GetCellMomentum(absId1, cells);
+
+ if ( nMax == 2 ) absId2 = absIdList[1];
+ else if( nMax == 1 ){
+ //Find second highest energy cell
+ Int_t enmax = 0 ;
+ for(Int_t iDigit = 0 ; iDigit < cluster->GetNCells() ; iDigit++){
+ Int_t absId = cluster->GetCellsAbsId()[iDigit];
+ if( absId == absId1 ) continue ;
+ Float_t endig = cells->GetCellAmplitude(absId);
+ RecalibrateCellAmplitude(endig,absId);
+ if(endig > enmax) {
+ enmax = endig ;
+ absId2 = absId ;
+ }
+
+ //Get mass of all cell in cluster combinations
+
+
+ Float_t enj = cells->GetCellAmplitude(absId);
+ RecalibrateCellAmplitude(enj,absId);
+
+ if(enj<0.3) continue;
+
+ TLorentzVector cellMomj = GetCellMomentum(absId, cells);
+
+ fhInvMassAllCells[0]->Fill(en,(cellMomj+cellMomi).M());
+
+ if(IsDataMC()) fhInvMassAllCells[mcindex]->Fill(en,(cellMomj+cellMomi).M());
+
+ }// cell loop
+
+ }
+ else { // loop on maxima, find 2 highest
+
+ Int_t enmax = 0 ;
+ for(Int_t iDigit = 0 ; iDigit < nMax ; iDigit++){
+ Float_t endig = maxEList[iDigit];
+ if(endig > enmax) {
+ endig = enmax ;
+ absId2 = absIdList[iDigit];
+ }
+ }// maxima loop
+
+ // First max is not highest, check if there is other higher
+ if(maxEList[0] < enmax){
+
+ for(Int_t iDigit = 0 ; iDigit < nMax ; iDigit++){
+ if(absId2 == absIdList[iDigit]) continue;
+ Float_t endig = maxEList[iDigit];
+ if(endig > enmax) {
+ endig = enmax ;
+ absId1 = absIdList[iDigit];
+ }
+ }// maxima loop
+
+ }
+
+ }
+
+ fhNLocMaxM02Cut[0]->Fill(en,nMax);
+
+ TLorentzVector cellMom1 = GetCellMomentum(absId1, cells);
+ TLorentzVector cellMom2 = GetCellMomentum(absId2, cells);
+
+ Float_t mass = (cellMom1+cellMom2).M();
+
+ if (nMax==1) fhMassNLocMax1[0]->Fill(en,mass);
+ else if(nMax==2) fhMassNLocMax2[0]->Fill(en,mass);
+ else if(nMax >2) fhMassNLocMaxN[0]->Fill(en,mass);
+
+ if (mass < 0.1) fhM02Con[0]->Fill(en,l0);
+ else if(mass < 0.2) fhM02Pi0[0]->Fill(en,l0);
+ else if(mass < 0.6 && mass > 0.4) fhM02Eta[0]->Fill(en,l0);
+
+ if(IsDataMC()){
+
+ fhNLocMaxM02Cut[mcindex]->Fill(en,nMax);
+
+ if (nMax==1) fhMassNLocMax1[mcindex]->Fill(en,mass);
+ else if(nMax==2) fhMassNLocMax2[mcindex]->Fill(en,mass);
+ else if(nMax >2) fhMassNLocMaxN[mcindex]->Fill(en,mass);
+
+ if (mass < 0.1) fhM02Con[mcindex]->Fill(en,l0);
+ else if(mass < 0.2) fhM02Pi0[mcindex]->Fill(en,l0);
+ else if(mass < 0.6 && mass > 0.4) fhM02Eta[mcindex]->Fill(en,l0);
+
+
+ }//Work with MC truth first
+
+ delete absIdList ;
+ delete maxEList ;
+
+ }//loop
+
+ if(GetDebug() > 1) printf("AliAnaInsideClusterInvariantMass::MakeAnalysisFillHistograms() - END \n");
+
+}
+
+//______________________________________________________________________
+void AliAnaInsideClusterInvariantMass::Print(const Option_t * opt) const
+{
+ //Print some relevant parameters set for the analysis
+ if(! opt)
+ return;
+
+ printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
+ AliAnaPartCorrBaseClass::Print("");
+ printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
+ printf("lambda 0 sqared > %2.1f\n", fM02Cut);
+ printf(" \n") ;
+
+}
+
+//____________________________________________________________________________________________
+void AliAnaInsideClusterInvariantMass::RecalibrateCellAmplitude(Float_t & amp, const Int_t id)
+{
+ //Recaculate cell energy if recalibration factor
+
+ Int_t icol = -1; Int_t irow = -1; Int_t iRCU = -1;
+ Int_t nModule = GetModuleNumberCellIndexes(id,fCalorimeter, icol, irow, iRCU);
+
+ if (GetCaloUtils()->IsRecalibrationOn()) {
+ if(fCalorimeter == "PHOS") {
+ amp *= GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModule,icol,irow);
+ }
+ else {
+ amp *= GetCaloUtils()->GetEMCALChannelRecalibrationFactor(nModule,icol,irow);
+ }
+ }
+}
+
+//____________________________________________________________________________________________
+void AliAnaInsideClusterInvariantMass::SplitEnergy(const Int_t absId1, const Int_t absId2, AliVCaloCells* cells,
+ Float_t & e1, Float_t & e2 )
+{
+
+ // Split energy of cluster between the 2 local maxima.
+
+ Int_t icol1 = -1, irow1 = -1, iRCU1 = -1;
+ Int_t sm1 = GetCaloUtils()->GetModuleNumberCellIndexes(absId1, fCalorimeter, icol1, irow1, iRCU1) ;
+ Int_t icol2 = -1, irow2 = -1, iRCU2 = -1;
+ Int_t sm2 = GetCaloUtils()->GetModuleNumberCellIndexes(absId2, fCalorimeter, icol2, irow2, iRCU2) ;
+
+ if(sm1!=sm2) {
+ if(sm1%2) icol1+=AliEMCALGeoParams::fgkEMCALCols;
+ else icol2+=AliEMCALGeoParams::fgkEMCALCols;
+ }
+
+/// continue here
+
+}
+
+
--- /dev/null
+#ifndef ALIANAINSIDECLUSTERINVARIANTMASS_H
+#define ALIANAINSIDECLUSTERINVARIANTMASS_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+/* $Id: AliAnaInsideClusterInvariantMass.h 27413 2008-07-18 13:28:12Z gconesab $ */
+
+//_________________________________________________________________________
+//
+// Split clusters with some criteria and calculate invariant mass
+// to identify them as pi0 or conversion
+//
+//
+//-- Author: Gustavo Conesa (LPSC-Grenoble)
+//_________________________________________________________________________
+
+
+// --- ROOT system ---
+class TList ;
+class TObjString;
+class TLorentzVector;
+
+// --- ANALYSIS system ---
+#include "AliAnaPartCorrBaseClass.h"
+
+class AliAnaInsideClusterInvariantMass : public AliAnaPartCorrBaseClass {
+
+ public:
+ AliAnaInsideClusterInvariantMass() ; // default ctor
+ virtual ~AliAnaInsideClusterInvariantMass() { ; } //virtual dtor
+ private:
+ AliAnaInsideClusterInvariantMass(const AliAnaInsideClusterInvariantMass & g) ; // cpy ctor
+ AliAnaInsideClusterInvariantMass & operator = (const AliAnaInsideClusterInvariantMass & g) ;//cpy assignment
+
+ public:
+
+ Bool_t AreNeighbours(const Int_t absId1, const Int_t absId2) const ;
+
+ TObjString * GetAnalysisCuts();
+
+ TList * GetCreateOutputObjects();
+
+ TLorentzVector GetCellMomentum(const Int_t absId, AliVCaloCells* cells);
+
+ Int_t GetNumberOfLocalMaxima(AliVCluster* cluster, AliVCaloCells* cells,
+ Int_t *absIdList, Float_t *maxEList) ;
+
+ void Init();
+
+ void InitParameters();
+
+ void MakeAnalysisFillHistograms() ;
+
+ void RecalibrateCellAmplitude(Float_t & amp, const Int_t absId);
+
+ void SetCalorimeter(TString & det) { fCalorimeter = det ; }
+
+ void SetM02Cut(Float_t cut) { fM02Cut = cut ; }
+
+ void SetMinNCells(Float_t cut) { fMinNCells = cut ; }
+
+ void SplitEnergy(const Int_t absId1, const Int_t absId2, AliVCaloCells* cells,
+ Float_t & e1, Float_t & e2 );
+
+ void Print(const Option_t * opt) const;
+
+ //For histograms
+ enum mcTypes { mcPhoton = 1, mcConversion = 2, mcPi0 = 3,
+ mcEta = 4, mcElectron = 5, mcHadron = 6 };
+
+ private:
+
+ TString fCalorimeter ; // Calorimeter where the gamma is searched
+ Float_t fM02Cut ; // Study clusters with l0 larger than cut
+ Float_t fMinNCells ; // Study clusters with ncells larger than cut
+
+ //Histograms
+
+ TH2F * fhMassNLocMax1[7] ; //! Mass of 2 highest energy cells when 1 local max, 1-6 for different MC particle types
+ TH2F * fhMassNLocMax2[7] ; //! Mass of 2 cells local maxima vs E, 1-6 for different MC particle types
+ TH2F * fhMassNLocMaxN[7] ; //! Mass of >2 cells local maxima vs E, 1-6 for different MC particle types
+
+ TH2F * fhNLocMax[7] ; //! Number of maxima in cluster vs E, 1-6 for different MC particle types
+ TH2F * fhNLocMaxM02Cut[7]; //! Number of maxima in cluster vs E, 1-6 for different MC particle types, after SS cut
+
+ TH2F * fhM02NLocMax1[7] ; //! M02 vs E for N max in cluster = 1, 1-6 for different MC particle types
+ TH2F * fhM02NLocMax2[7] ; //! M02 vs E for N max in cluster = 2, 1-6 for different MC particle types
+ TH2F * fhM02NLocMaxN[7] ; //! M02 vs E for N max in cluster > 2, 1-6 for different MC particle types
+
+ TH2F * fhNCellNLocMax1[7] ; //! n cells in cluster vs E for N max in cluster = 1, 1-6 for different MC particle types
+ TH2F * fhNCellNLocMax2[7] ; //! n cells in cluster vs E for N max in cluster = 2, 1-6 for different MC particle types
+ TH2F * fhNCellNLocMaxN[7] ; //! n cells in cluster vs E for N max in cluster > 2, 1-6 for different MC particle types
+
+ TH2F * fhM02Pi0[7] ; //! M02 for Mass around pi0
+ TH2F * fhM02Eta[7] ; //! M02 for Mass around eta
+ TH2F * fhM02Con[7] ; //! M02 for Mass around close to 0
+
+ TH2F * fhInvMassAllCells[7] ; //! Inv mass of all cells
+
+
+ ClassDef(AliAnaInsideClusterInvariantMass,1)
+
+} ;
+
+#endif //ALIANAINSIDECLUSTERINVARIANTMASS_H
+
+
+
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, 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. *
- **************************************************************************/
-
-/* $Id: */
-
-//________________________________________________________________________
-// A basic analysis task to analyse photon detected by PHOS
-// A basic analysis task to analyse photon detected by PHOS
-// A basic analysis task to analyse photon detected by PHOS
-// Adapted for AliAnalysisTaskSE and AOD production
-// by Gustavo Conesa
-//
-//*-- Yves Schutz
-//////////////////////////////////////////////////////////////////////////////
-//Root system
-#include <TCanvas.h>
-#include <TH1D.h>
-#include <TROOT.h>
-#include <TNtuple.h>
-#include <TVector3.h>
-
-//Analysis system
-#include "AliAnalysisTaskPHOSExample.h"
-#include "AliESDEvent.h"
-#include "AliESDCaloCluster.h"
-#include "AliAODEvent.h"
-#include "AliAODPhoton.h"
-#include "AliLog.h"
-
-
-//______________________________________________________________________________
-AliAnalysisTaskPHOSExample::AliAnalysisTaskPHOSExample() :
- fDebug(0),
- fAODPhotons(0x0),
- fPhotonsInPhos(0),
- fPhotonId(1.0),
- fOutputList(0x0),
- fhPHOSPos(0),
- fhPHOS(0),
- fhPHOSEnergy(0),
- fhPHOSDigits(0),
- fhPHOSRecParticles(0),
- fhPHOSPhotons(0),
- fhPHOSInvariantMass(0),
- fhPHOSDigitsEvent(0)
-{
- //Default constructor
-}
-//______________________________________________________________________________
-AliAnalysisTaskPHOSExample::AliAnalysisTaskPHOSExample(const char *name) :
- AliAnalysisTaskSE(name),
- fDebug(0),
- fAODPhotons(0x0),
- fPhotonsInPhos(0),
- fPhotonId(1.0),
- fOutputList(0x0),
- fhPHOSPos(0),
- fhPHOS(0),
- fhPHOSEnergy(0),
- fhPHOSDigits(0),
- fhPHOSRecParticles(0),
- fhPHOSPhotons(0),
- fhPHOSInvariantMass(0),
- fhPHOSDigitsEvent(0)
-{
- // Constructor.
-
- // Output slots
- DefineOutput(1, TList::Class()) ;
-}
-
-//____________________________________________________________________________
-AliAnalysisTaskPHOSExample::AliAnalysisTaskPHOSExample(const AliAnalysisTaskPHOSExample& ap) :
- AliAnalysisTaskSE(ap.GetName()),
- fDebug(ap.fDebug),
- fAODPhotons(ap.fAODPhotons),
- fPhotonsInPhos(ap.fPhotonsInPhos),
- fPhotonId(ap.fPhotonId),
- fOutputList(ap.fOutputList),
- fhPHOSPos(ap.fhPHOSPos),
- fhPHOS(ap.fhPHOS),
- fhPHOSEnergy(ap.fhPHOSEnergy),
- fhPHOSDigits(ap.fhPHOSDigits),
- fhPHOSRecParticles(ap.fhPHOSRecParticles),
- fhPHOSPhotons(ap.fhPHOSPhotons),
- fhPHOSInvariantMass(ap.fhPHOSInvariantMass),
- fhPHOSDigitsEvent(ap.fhPHOSDigitsEvent)
-{
- // cpy ctor
-}
-
-//_____________________________________________________________________________
-AliAnalysisTaskPHOSExample& AliAnalysisTaskPHOSExample::operator = (const AliAnalysisTaskPHOSExample& ap)
-{
-// assignment operator
-
- this->~AliAnalysisTaskPHOSExample();
- new(this) AliAnalysisTaskPHOSExample(ap);
-
- fDebug = ap.fDebug;
- fAODPhotons = ap.fAODPhotons;
- fPhotonsInPhos = ap.fPhotonsInPhos;
- fPhotonId = ap.fPhotonId;
- fOutputList = ap.fOutputList;
- fhPHOSPos = ap.fhPHOSPos;
- fhPHOS = ap.fhPHOS;
- fhPHOSEnergy = ap.fhPHOSEnergy;
- fhPHOSDigits = ap.fhPHOSDigits;
- fhPHOSRecParticles = ap.fhPHOSRecParticles;
- fhPHOSPhotons = ap.fhPHOSPhotons;
- fhPHOSInvariantMass = ap.fhPHOSInvariantMass;
- fhPHOSDigitsEvent = ap.fhPHOSDigitsEvent;
-
- return *this;
-}
-
-//______________________________________________________________________________
-AliAnalysisTaskPHOSExample::~AliAnalysisTaskPHOSExample()
-{
- // dtor
- if(fOutputList) {
- fOutputList->Clear() ;
- delete fOutputList ;
- }
-}
-
-
-//________________________________________________________________________
-void AliAnalysisTaskPHOSExample::UserCreateOutputObjects()
-{
- // Create the outputs containers
- //AODs
- fAODPhotons = new TClonesArray("AliAODPhoton", 0);
- fAODPhotons->SetName("Photons");
- AddAODBranch("TClonesArray", &fAODPhotons);
-
- OpenFile(1) ;
-
- fhPHOSPos = new TNtuple("PHOSPos" , "Position in PHOS" , "x:y:z");
- fhPHOS = new TNtuple("PHOS" , "PHOS" , "event:digits:clusters:photons");
- fhPHOSEnergy = new TH1D("PHOSEnergy" , "PHOSEnergy" , 100, 0., 100. ) ;
- fhPHOSDigits = new TH1I("PHOSDigitsCluster" , "PHOSDigits" , 20 , 0 , 20 ) ;
- fhPHOSRecParticles = new TH1D("PHOSRecParticles" , "PHOSRecParticles" , 20 , 0., 20. ) ;
- fhPHOSPhotons = new TH1I("PHOSPhotons" , "PHOSPhotons" , 20 , 0 , 20 ) ;
- fhPHOSInvariantMass = new TH1D("PHOSInvariantMass" , "PHOSInvariantMass" , 400, 0., 400.) ;
- fhPHOSDigitsEvent = new TH1I("PHOSDigitsEvent" , "PHOSDigitsEvent" , 30 , 0 , 30 ) ;
-
- // create output container
-
- fOutputList = new TList() ;
- fOutputList->SetName(GetName()) ;
-
- fOutputList->AddAt(fhPHOSPos, 0) ;
- fOutputList->AddAt(fhPHOS, 1) ;
- fOutputList->AddAt(fhPHOSEnergy, 2) ;
- fOutputList->AddAt(fhPHOSDigits, 3) ;
- fOutputList->AddAt(fhPHOSRecParticles, 4) ;
- fOutputList->AddAt(fhPHOSPhotons, 5) ;
- fOutputList->AddAt(fhPHOSInvariantMass, 6) ;
- fOutputList->AddAt(fhPHOSDigitsEvent, 7) ;
-
-
-}
-
-//______________________________________________________________________________
-void AliAnalysisTaskPHOSExample::UserExec(Option_t *)
-{
- // Processing of one event
-
- // if ( !((Entry()-1)%100) )
- AliInfo(Form(" Processing event # %lld", Entry())) ;
- AliESDEvent* esd = (AliESDEvent*)InputEvent();
-
- //************************ PHOS *************************************
- TRefArray * caloClustersArr = new TRefArray();
- esd->GetPHOSClusters(caloClustersArr);
-
- const Int_t kNumberOfPhosClusters = caloClustersArr->GetEntries() ;
-
- TVector3 ** phosVector = new TVector3*[kNumberOfPhosClusters] ;
- Float_t * phosPhotonsEnergy = new Float_t[kNumberOfPhosClusters] ;
- Int_t phosCluster ;
- Int_t numberOfDigitsInPhos = 0 ;
- Double_t v[3] ; //vertex ;
- esd->GetVertex()->GetXYZ(v) ;
-
- fPhotonsInPhos = 0 ;
- // loop over the PHOS Cluster
- for(phosCluster = 0 ; phosCluster < kNumberOfPhosClusters ; phosCluster++) {
- AliESDCaloCluster * caloCluster = (AliESDCaloCluster *) caloClustersArr->At(phosCluster) ;
-
- //AliESDCaloCluster * caloCluster = fESD->GetCaloCluster(phosCluster) ;
- if (caloCluster) {
- Float_t pos[3] ;
- caloCluster->GetPosition( pos ) ;
- fhPHOSEnergy->Fill( caloCluster->E() ) ;
- fhPHOSPos->Fill( pos[0], pos[1], pos[2] ) ;
- fhPHOSDigits->Fill(Entry(), caloCluster->GetNCells() ) ;
- numberOfDigitsInPhos += caloCluster->GetNCells() ;
- const Double_t * pid = caloCluster->GetPID() ;
- if(pid[AliPID::kPhoton] > GetPhotonId() ) {
- phosVector[fPhotonsInPhos] = new TVector3(pos[0],pos[1],pos[2]) ;
- phosPhotonsEnergy[fPhotonsInPhos]=caloCluster->E() ;
- TLorentzVector momentum ;
- caloCluster->GetMomentum(momentum, v);
- new ((*fAODPhotons)[fPhotonsInPhos++]) AliAODPhoton (momentum);
- }
- }
- } //PHOS clusters
-
- fhPHOSRecParticles->Fill(kNumberOfPhosClusters);
- fhPHOSPhotons->Fill(fPhotonsInPhos);
- fhPHOSDigitsEvent->Fill(numberOfDigitsInPhos);
- fhPHOS->Fill(Entry(), numberOfDigitsInPhos, kNumberOfPhosClusters, fPhotonsInPhos) ;
-
- // invariant Mass
- if (fPhotonsInPhos > 1 ) {
- Int_t phosPhoton1, phosPhoton2 ;
- for(phosPhoton1 = 0 ; phosPhoton1 < fPhotonsInPhos ; phosPhoton1++) {
- for(phosPhoton2 = phosPhoton1 + 1 ; phosPhoton2 < fPhotonsInPhos ; phosPhoton2++) {
- Float_t tempMass = TMath::Sqrt( 2 * phosPhotonsEnergy[phosPhoton1] * phosPhotonsEnergy[phosPhoton2] *
- ( 1 - TMath::Cos(phosVector[phosPhoton1]->Angle(*phosVector[phosPhoton2])) )
- );
- fhPHOSInvariantMass->Fill(tempMass*1000.);
- }
- }
- }
-
- PostData(1, fOutputList);
-
- delete [] phosVector ;
- delete [] phosPhotonsEnergy ;
-
-}
-
-
-//______________________________________________________________________________
-void AliAnalysisTaskPHOSExample::Init()
-{
- // Intialisation of parameters
- AliInfo("Doing initialisation") ;
- fPhotonId = 0.9 ;
-}
-
-//______________________________________________________________________________
-void AliAnalysisTaskPHOSExample::Terminate(Option_t *)
-{
- // Processing when the event loop is ended
-
-// Bool_t problem = kFALSE ;
- AliInfo(Form(" *** %s Report:", GetName())) ;
- printf(" PHOSEnergy Mean : %5.3f , RMS : %5.3f \n", fhPHOSEnergy->GetMean(), fhPHOSEnergy->GetRMS() ) ;
- printf(" PHOSDigits Mean : %5.3f , RMS : %5.3f \n", fhPHOSDigits->GetMean(), fhPHOSDigits->GetRMS() ) ;
- printf(" PHOSRecParticles Mean : %5.3f , RMS : %5.3f \n", fhPHOSRecParticles->GetMean(), fhPHOSRecParticles->GetRMS() ) ;
- printf(" PHOSPhotons Mean : %5.3f , RMS : %5.3f \n", fhPHOSPhotons->GetMean(), fhPHOSPhotons->GetRMS() ) ;
- printf(" PHOSInvariantMass Mean : %5.3f , RMS : %5.3f \n", fhPHOSInvariantMass->GetMean(), fhPHOSInvariantMass->GetRMS() ) ;
- printf(" PHOSDigitsEvent Mean : %5.3f , RMS : %5.3f \n", fhPHOSDigitsEvent->GetMean(), fhPHOSDigitsEvent->GetRMS() ) ;
-
- TCanvas * cPHOS = new TCanvas("cPHOS", "PHOS ESD Test", 400, 10, 600, 700) ;
- cPHOS->Divide(3, 2);
-
- cPHOS->cd(1) ;
- if ( fhPHOSEnergy->GetMaximum() > 0. )
- gPad->SetLogy();
- fhPHOSEnergy->SetAxisRange(0, 25.);
- fhPHOSEnergy->SetLineColor(2);
- fhPHOSEnergy->Draw();
-
- cPHOS->cd(2) ;
- fhPHOSDigits->SetAxisRange(0,25.);
- fhPHOSDigits->SetLineColor(2);
- fhPHOSDigits->Draw();
-
- cPHOS->cd(3) ;
- if ( fhPHOSRecParticles->GetMaximum() > 0. )
- gPad->SetLogy();
- fhPHOSRecParticles->SetAxisRange(0, 25.);
- fhPHOSRecParticles->SetLineColor(2);
- fhPHOSRecParticles->Draw();
-
- cPHOS->cd(4) ;
- if ( fhPHOSPhotons->GetMaximum() > 0. )
- gPad->SetLogy();
- fhPHOSPhotons->SetAxisRange(0,25.);
- fhPHOSPhotons->SetLineColor(2);
- fhPHOSPhotons->Draw();
-
- cPHOS->cd(5) ;
- fhPHOSInvariantMass->SetLineColor(2);
- fhPHOSInvariantMass->Draw();
-
- cPHOS->cd(6) ;
- if ( fhPHOSDigitsEvent->GetMaximum() > 0. )
- gPad->SetLogy();
- fhPHOSDigitsEvent->SetAxisRange(0,40.);
- fhPHOSDigitsEvent->SetLineColor(2);
- fhPHOSDigitsEvent->Draw();
-
- cPHOS->Print("PHOS.eps");
-
- char line[1024] ;
- snprintf(line,1024, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
- gROOT->ProcessLine(line);
- snprintf(line,1024, ".!rm -fR *.eps");
- gROOT->ProcessLine(line);
-
- AliInfo(Form("!!! All the eps files are in %s.tar.gz !!!", GetName())) ;
-
-// char * report = 0x0 ;
-// if(problem)
-// sprintf(report,"Problems found, please check!!!");
-// else
-// sprintf(report,"OK");
-
-// AliInfo(Form("*** %s Summary Report: %s \n",GetName(), report)) ;
-}
+++ /dev/null
-#ifndef ALIANALYSISTASKPHOSEXAMPLE_H
-#define ALIANALYSISTASKPHOSEXAMPLE_H
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice */
-//______________________________________________________________________________
-// A basic analysis task to analyse photon detected by PHOS
-// A basic analysis task to analyse photon detected by PHOS
-// Adapted for AliAnalysisTaskSE and AOD production
-// by Gustavo Conesa
-//
-//*-- Yves Schutz
-//////////////////////////////////////////////////////////////////////////////
-
-class TTree ;
-#include "AliAnalysisTaskSE.h"
-
-class AliESDEvent ;
-class AliAODEvent ;
-class TNtuple ;
-class TH1D ;
-class TH1I ;
-
-class AliAnalysisTaskPHOSExample : public AliAnalysisTaskSE {
-
-public:
- AliAnalysisTaskPHOSExample() ;
- AliAnalysisTaskPHOSExample(const char *name) ;
- AliAnalysisTaskPHOSExample(const AliAnalysisTaskPHOSExample& ap) ;
- AliAnalysisTaskPHOSExample& operator = (const AliAnalysisTaskPHOSExample& ap) ;
- virtual ~AliAnalysisTaskPHOSExample() ;
-
- virtual void UserCreateOutputObjects();
- virtual void Init() ;
- virtual void LocalInit() { Init() ; }
- virtual void UserExec(Option_t * opt = "") ;
- Float_t GetPhotonId() const { return fPhotonId ; }
- void SetDebugLevel(Int_t level) { fDebug = level ; }
- void SetPhotonId(Float_t threshold) { fPhotonId = threshold ; }
- virtual void Terminate(Option_t * opt = "") ;
-
-private:
- // input and output
- Int_t fDebug ; // Debug flag
- TClonesArray * fAODPhotons ; //! reconstructed photons
- Int_t fPhotonsInPhos ; //! number of photons found
- // task parameters
- Float_t fPhotonId ; // threshold for photon identification
-
- // Histograms
- TList * fOutputList ; //! output data list
- TNtuple * fhPHOSPos ; //! PHOS (x,y)
- TNtuple * fhPHOS ; //! all PHOS parameters
- TH1D * fhPHOSEnergy ; //! PHOS energy
- TH1I * fhPHOSDigits ; //! PHOS numer of SDigits
- TH1D * fhPHOSRecParticles ;//! PHOS number of RecParticles
- TH1I * fhPHOSPhotons ; //! PHOS number of photons
- TH1D * fhPHOSInvariantMass ; //! PHOS invariant mass
- TH1I * fhPHOSDigitsEvent ; //! PHOS numbet of Sdigits per event
-
- ClassDef(AliAnalysisTaskPHOSExample, 1); // a PHOS photon analysis task
-};
-#endif // ALIANALYSISTASKPHOSEXAMPLE_H