[u/mrichter/AliRoot.git] / PWG4 / AliAnaGammaIsolCut.cxx

/**************************************************************************
 * 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$ */

/* History of cvs commits:
 *
 * $Log$
 * Revision 1.1  2007/01/23 17:17:29  schutz
 * New Gamma package
 *
 *
 */

//_________________________________________________________________________
// Class for the analysis of gamma correlations (gamma-jet, 
// gamma-hadron and isolation cut.
// This class makes isolation cut analysis for 2 IC methods 
//(cone pt sum and particle pt threshold), for different cone sizes 
//and pt thresholds
//
//  Class created from old AliPHOSGammaJet 
//  (see AliRoot versions previous Release 4-09)
//
//*-- Author: Gustavo Conesa (LNF-INFN) 
//////////////////////////////////////////////////////////////////////////////


// --- ROOT system ---

#include <TFile.h>
#include <TParticle.h>
#include <TH2.h>

#include "AliAnaGammaIsolCut.h" 
#include "AliESD.h"
#include "AliESDtrack.h"
#include "AliESDCaloCluster.h"
#include "Riostream.h"
#include "AliLog.h"

ClassImp(AliAnaGammaIsolCut)

//____________________________________________________________________________
  AliAnaGammaIsolCut::AliAnaGammaIsolCut(const char *name) : 
    AliAnaGammaDirect(name),
    fOutputContainer(new TObjArray(100)),  
    fNCones(0),fNPtThres(0)
{
  //Ctor        
  TList * list = gDirectory->GetListOfKeys() ; 
  TIter next(list) ; 
  TH2F * h = 0 ;
  Int_t index ; 
  for (index = 0 ; index < list->GetSize()-1 ; index++) { 
    //-1 to avoid GammaJet Task
    h = dynamic_cast<TH2F*>(gDirectory->Get(list->At(index)->GetName())) ; 
    fOutputContainer->Add(h) ; 
  }
  
  for(Int_t i = 0; i < 10 ; i++){
    fConeSizes[i]=0;
    fPtThresholds[i]=0;
  }
  
  // Input slot #0 works with an Ntuple
  DefineInput(0, TChain::Class());
  // Output slot #0 writes into a TH1 container
  DefineOutput(0,  TObjArray::Class()) ; 
  
}


//____________________________________________________________________________
AliAnaGammaIsolCut::AliAnaGammaIsolCut(const AliAnaGammaIsolCut & ic) : 
  AliAnaGammaDirect(ic),
  fOutputContainer(ic. fOutputContainer), 
  fNCones(ic.fNCones),fNPtThres(ic.fNPtThres)
{
  // cpy ctor
  SetName (ic.GetName()) ; 
  SetTitle(ic.GetTitle()) ; 

  for(Int_t i = 0; i < 10 ; i++){
    fConeSizes[i]=  ic.fConeSizes[i];
    fPtThresholds[i]=   ic.fPtThresholds[i];
  }
}

//____________________________________________________________________________
AliAnaGammaIsolCut::~AliAnaGammaIsolCut() 
{
  // Remove all pointers
  fOutputContainer->Clear() ; 
  delete fOutputContainer ;

  delete fhPtCandidate ;
  delete [] fhPtThresIsolated ;
  delete [] fhPtSumIsolated ;

}


//____________________________________________________________________________
void AliAnaGammaIsolCut::Exec(Option_t *) 
{
  
  // Processing of one event
    
  //Get ESDs
  Long64_t entry = GetChain()->GetReadEntry() ;
  
  if (!GetESD()) {
    AliError("fESD is not connected to the input!") ; 
    return ; 
  }
  
  if (GetPrintInfo()) 
    AliInfo(Form("%s ----> Processing event # %lld",  (dynamic_cast<TChain *>(GetChain()))->GetFile()->GetName(), entry)) ; 

  //CreateTLists with arrays of TParticles. Filled with particles only relevant for the analysis.

  TClonesArray * particleList = new TClonesArray("TParticle",1000); // All particles refitted in CTS and detected in EMCAL (jet)
  TClonesArray * plCTS         = new TClonesArray("TParticle",1000); // All particles refitted in Central Tracking System (ITS+TPC)
  TClonesArray * plNe         = new TClonesArray("TParticle",1000);   // All particles measured in Jet Calorimeter (EMCAL)
  TClonesArray * plPHOS     = new TClonesArray("TParticle",1000);  // All particles measured in PHOS as Gamma calorimeter
  TClonesArray * plEMCAL   = new TClonesArray("TParticle",1000);  // All particles measured in EMCAL as Gamma calorimeter
  
  
  AliDebug(2, "Fill particle lists");
  
  //Fill particle lists 
  CreateParticleList(particleList, plCTS,plEMCAL,plPHOS); 
  
  if(GetCalorimeter() == "PHOS")
    plNe = plPHOS;
  if(GetCalorimeter() == "EMCAL")
    plNe = plEMCAL;
  
  //Isolation Cut Analysis for both methods and different pt cuts and cones
  
  for(Int_t ipr = 0; ipr < plNe->GetEntries() ; ipr ++ ){
    TParticle * pCandidate = dynamic_cast<TParticle *>(plNe->At(ipr)) ;
    
    if(pCandidate->Pt() > GetMinGammaPt()){
      
      Bool_t  icPtThres   = kFALSE;
      Bool_t  icPtSum     = kFALSE;
      
      Float_t ptC             = pCandidate->Pt() ;
      
      fhPtCandidate->Fill(ptC);
      
      for(Int_t icone = 0; icone<fNCones; icone++){
	SetConeSize(fConeSizes[icone]) ;
	Float_t coneptsum = 0 ;
	for(Int_t ipt = 0; ipt<fNPtThres;ipt++){ 
	  SetPtThreshold(fPtThresholds[ipt]) ;
	  MakeIsolationCut(plCTS,plNe, pCandidate, ipr, icPtThres, icPtSum,coneptsum);
	  AliDebug(4,Form("Candidate pt %f, pt in cone %f, Isolated? ICPt %d, ICSum %d",
			  pCandidate->Pt(), coneptsum, icPtThres, icPtSum));
	  
	  fhPtThresIsolated[icone][ipt]->Fill(ptC); 
	}//pt thresh loop
	fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ;
      }//cone size loop
    }//min pt candidate
  }//candidate loop
  
  AliDebug(2, "End of analysis, delete pointers");
  
  particleList->Delete() ; 
  plCTS->Delete() ;
  plNe->Delete() ;
  plPHOS->Delete() ;
  plEMCAL->Delete() ;

  delete plNe ;
  delete plCTS ;
  //delete plPHOS ;
  //delete plEMCAL ;
  delete particleList ;

  PostData(0, fOutputContainer);
}    

  //____________________________________________________________________________
void AliAnaGammaIsolCut::Init(const Option_t * )
{
  // Initialisation of branch container 
  AliAnaGammaDirect::Init();

  fNCones           = 4 ; 
  fNPtThres         = 4 ; 
  fConeSizes[0] = 0.1; fConeSizes[0] = 0.2; fConeSizes[2] = 0.3; fConeSizes[3] = 0.4;
  fPtThresholds[0]=1.; fPtThresholds[0]=2.; fPtThresholds[0]=3.; fPtThresholds[0]=4.;

  //Initialization of histograms 
  MakeHistos() ;
}

//___________________________________________________________________
void AliAnaGammaIsolCut::MakeHistos()
{
  // Create histograms to be saved in output file and 
  // stores them in fOutputContainer
  
  fOutputContainer = new TObjArray(10000) ; 

  //Isolation cut histograms
  fhPtCandidate  = new TH1F
    ("PtCandidate","p_{T} of candidate particles for isolation",240,0,120); 
  fhPtCandidate->SetXTitle("p_{T} (GeV/c)");
  fOutputContainer->Add(fhPtCandidate) ;
  
  char name[128];
  char title[128];
  for(Int_t icone = 0; icone<fNCones; icone++){
    sprintf(name,"PtSumIsolated_Cone_%d",icone);
    sprintf(title,"Candidate cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
    fhPtSumIsolated[icone]  = new TH2F(name, title,240,0,120,120,0,10);
    fhPtSumIsolated[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
    fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)");
    fOutputContainer->Add(fhPtSumIsolated[icone]) ; 
    
    for(Int_t ipt = 0; ipt<fNPtThres;ipt++){ 
      sprintf(name,"PtThresIsol_Cone_%d_Pt%d",icone,ipt);
      sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
      fhPtThresIsolated[icone][ipt]  = new TH1F(name, title,240,0,120);
      fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
      fOutputContainer->Add(fhPtThresIsolated[icone][ipt]) ; 
    }//icone loop
  }//ipt loop

}


void AliAnaGammaIsolCut::Print(const Option_t * opt) const
{

  //Print some relevant parameters set for the analysis
  if(! opt)
    return;

  Info("Print", "%s %s", GetName(), GetTitle() ) ;
  printf("N Cone Sizes               =     %d\n", fNCones) ; 
  printf("N pT thresholds           =     %d\n", fNPtThres) ;
  printf("Cone Sizes                  =    \n") ;
  for(Int_t i = 0; i < fNCones; i++)
    printf("   %f;",  fConeSizes[i]) ;
  printf("    \n") ;
  for(Int_t i = 0; i < fNPtThres; i++)
    printf("   %f;",  fPtThresholds[i]) ;
  printf("    \n") ;
  
} 

void AliAnaGammaIsolCut::Terminate(Option_t *)
{
  // The Terminate() function is the last function to be called during
  // a query. It always runs on the client, it can be used to present
  // the results graphically or save the results to file.
 
}
Commit	Line	Data
2d11ded3	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, 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	/* $Id$ */
	16
	17	/* History of cvs commits:
	18	*
	19	* $Log$
	20	* Revision 1.1 2007/01/23 17:17:29 schutz
	21	* New Gamma package
	22	*
	23	*
	24	*/
	25
	26	//_________________________________________________________________________
	27	// Class for the analysis of gamma correlations (gamma-jet,
	28	// gamma-hadron and isolation cut.
	29	// This class makes isolation cut analysis for 2 IC methods
	30	//(cone pt sum and particle pt threshold), for different cone sizes
	31	//and pt thresholds
	32	//
	33	// Class created from old AliPHOSGammaJet
	34	// (see AliRoot versions previous Release 4-09)
	35	//
	36	//*-- Author: Gustavo Conesa (LNF-INFN)
	37	//////////////////////////////////////////////////////////////////////////////
	38
	39
	40	// --- ROOT system ---
	41
	42	#include <TFile.h>
	43	#include <TParticle.h>
	44	#include <TH2.h>
	45
	46	#include "AliAnaGammaIsolCut.h"
	47	#include "AliESD.h"
	48	#include "AliESDtrack.h"
	49	#include "AliESDCaloCluster.h"
	50	#include "Riostream.h"
	51	#include "AliLog.h"
	52
	53	ClassImp(AliAnaGammaIsolCut)
	54
	55	//____________________________________________________________________________
	56	AliAnaGammaIsolCut::AliAnaGammaIsolCut(const char *name) :
	57	AliAnaGammaDirect(name),
	58	fOutputContainer(new TObjArray(100)),
	59	fNCones(0),fNPtThres(0)
	60	{
	61	//Ctor
	62	TList * list = gDirectory->GetListOfKeys() ;
	63	TIter next(list) ;
	64	TH2F * h = 0 ;
65	Int_t index ;
66	for (index = 0 ; index < list->GetSize()-1 ; index++) {
67	//-1 to avoid GammaJet Task
68	h = dynamic_cast<TH2F*>(gDirectory->Get(list->At(index)->GetName())) ;
69	fOutputContainer->Add(h) ;
70	}
71
72	for(Int_t i = 0; i < 10 ; i++){
73	fConeSizes[i]=0;
74	fPtThresholds[i]=0;
75	}
76
77	// Input slot #0 works with an Ntuple
78	DefineInput(0, TChain::Class());
79	// Output slot #0 writes into a TH1 container
80	DefineOutput(0, TObjArray::Class()) ;
81
82	}
83
84
85	//____________________________________________________________________________
86	AliAnaGammaIsolCut::AliAnaGammaIsolCut(const AliAnaGammaIsolCut & ic) :
87	AliAnaGammaDirect(ic),
88	fOutputContainer(ic. fOutputContainer),
89	fNCones(ic.fNCones),fNPtThres(ic.fNPtThres)
90	{
91	// cpy ctor
92	SetName (ic.GetName()) ;
93	SetTitle(ic.GetTitle()) ;
94
95	for(Int_t i = 0; i < 10 ; i++){
96	fConeSizes[i]= ic.fConeSizes[i];
97	fPtThresholds[i]= ic.fPtThresholds[i];
98	}
99	}
100
101	//____________________________________________________________________________
102	AliAnaGammaIsolCut::~AliAnaGammaIsolCut()
103	{
104	// Remove all pointers
105	fOutputContainer->Clear() ;
106	delete fOutputContainer ;
107
108	delete fhPtCandidate ;
109	delete [] fhPtThresIsolated ;
110	delete [] fhPtSumIsolated ;
111
112	}
113
114
115
116	//____________________________________________________________________________
117	void AliAnaGammaIsolCut::Exec(Option_t *)
118	{
119
120	// Processing of one event
121
122	//Get ESDs
123	Long64_t entry = GetChain()->GetReadEntry() ;
124
125	if (!GetESD()) {
126	AliError("fESD is not connected to the input!") ;
127	return ;
128	}
129
130	if (GetPrintInfo())
131	AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(GetChain()))->GetFile()->GetName(), entry)) ;
132
133	//CreateTLists with arrays of TParticles. Filled with particles only relevant for the analysis.
134
135	TClonesArray * particleList = new TClonesArray("TParticle",1000); // All particles refitted in CTS and detected in EMCAL (jet)
136	TClonesArray * plCTS = new TClonesArray("TParticle",1000); // All particles refitted in Central Tracking System (ITS+TPC)
137	TClonesArray * plNe = new TClonesArray("TParticle",1000); // All particles measured in Jet Calorimeter (EMCAL)
138	TClonesArray * plPHOS = new TClonesArray("TParticle",1000); // All particles measured in PHOS as Gamma calorimeter
139	TClonesArray * plEMCAL = new TClonesArray("TParticle",1000); // All particles measured in EMCAL as Gamma calorimeter
140
141
142	AliDebug(2, "Fill particle lists");
143
144	//Fill particle lists
145	CreateParticleList(particleList, plCTS,plEMCAL,plPHOS);
146
147	if(GetCalorimeter() == "PHOS")
148	plNe = plPHOS;
149	if(GetCalorimeter() == "EMCAL")
150	plNe = plEMCAL;
151
152	//Isolation Cut Analysis for both methods and different pt cuts and cones
153
154	for(Int_t ipr = 0; ipr < plNe->GetEntries() ; ipr ++ ){
155	TParticle * pCandidate = dynamic_cast<TParticle *>(plNe->At(ipr)) ;
156
157	if(pCandidate->Pt() > GetMinGammaPt()){
158
159	Bool_t icPtThres = kFALSE;
160	Bool_t icPtSum = kFALSE;
161
162	Float_t ptC = pCandidate->Pt() ;
163
164	fhPtCandidate->Fill(ptC);
165
166	for(Int_t icone = 0; icone<fNCones; icone++){
167	SetConeSize(fConeSizes[icone]) ;
168	Float_t coneptsum = 0 ;
169	for(Int_t ipt = 0; ipt<fNPtThres;ipt++){
170	SetPtThreshold(fPtThresholds[ipt]) ;
171	MakeIsolationCut(plCTS,plNe, pCandidate, ipr, icPtThres, icPtSum,coneptsum);
172	AliDebug(4,Form("Candidate pt %f, pt in cone %f, Isolated? ICPt %d, ICSum %d",
173	pCandidate->Pt(), coneptsum, icPtThres, icPtSum));
174
175	fhPtThresIsolated[icone][ipt]->Fill(ptC);
176	}//pt thresh loop
177	fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ;
178	}//cone size loop
179	}//min pt candidate
180	}//candidate loop
181
182	AliDebug(2, "End of analysis, delete pointers");
183
184	particleList->Delete() ;
185	plCTS->Delete() ;
186	plNe->Delete() ;
187	plPHOS->Delete() ;
188	plEMCAL->Delete() ;
189
190	delete plNe ;
191	delete plCTS ;
192	//delete plPHOS ;
193	//delete plEMCAL ;
194	delete particleList ;
195
196	PostData(0, fOutputContainer);
197	}
198
199	//____________________________________________________________________________
200	void AliAnaGammaIsolCut::Init(const Option_t * )
201	{
202	// Initialisation of branch container
203	AliAnaGammaDirect::Init();
204
205	fNCones = 4 ;
206	fNPtThres = 4 ;
207	fConeSizes[0] = 0.1; fConeSizes[0] = 0.2; fConeSizes[2] = 0.3; fConeSizes[3] = 0.4;
208	fPtThresholds[0]=1.; fPtThresholds[0]=2.; fPtThresholds[0]=3.; fPtThresholds[0]=4.;
209
210	//Initialization of histograms
211	MakeHistos() ;
212	}
213
214	//___________________________________________________________________
215	void AliAnaGammaIsolCut::MakeHistos()
216	{
217	// Create histograms to be saved in output file and
218	// stores them in fOutputContainer
219
220	fOutputContainer = new TObjArray(10000) ;
221
222	//Isolation cut histograms
223	fhPtCandidate = new TH1F
224	("PtCandidate","p_{T} of candidate particles for isolation",240,0,120);
225	fhPtCandidate->SetXTitle("p_{T} (GeV/c)");
226	fOutputContainer->Add(fhPtCandidate) ;
227
228	char name[128];
229	char title[128];
230	for(Int_t icone = 0; icone<fNCones; icone++){
231	sprintf(name,"PtSumIsolated_Cone_%d",icone);
232	sprintf(title,"Candidate cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
233	fhPtSumIsolated[icone] = new TH2F(name, title,240,0,120,120,0,10);
234	fhPtSumIsolated[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
235	fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)");
236	fOutputContainer->Add(fhPtSumIsolated[icone]) ;
237
238	for(Int_t ipt = 0; ipt<fNPtThres;ipt++){
239	sprintf(name,"PtThresIsol_Cone_%d_Pt%d",icone,ipt);
240	sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
241	fhPtThresIsolated[icone][ipt] = new TH1F(name, title,240,0,120);
242	fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
243	fOutputContainer->Add(fhPtThresIsolated[icone][ipt]) ;
244	}//icone loop
245	}//ipt loop
246
247	}
248
249
250	void AliAnaGammaIsolCut::Print(const Option_t * opt) const
251	{
252
253	//Print some relevant parameters set for the analysis
254	if(! opt)
255	return;
256
257	Info("Print", "%s %s", GetName(), GetTitle() ) ;
258	printf("N Cone Sizes = %d\n", fNCones) ;
259	printf("N pT thresholds = %d\n", fNPtThres) ;
260	printf("Cone Sizes = \n") ;
261	for(Int_t i = 0; i < fNCones; i++)
262	printf(" %f;", fConeSizes[i]) ;
263	printf(" \n") ;
264	for(Int_t i = 0; i < fNPtThres; i++)
265	printf(" %f;", fPtThresholds[i]) ;
266	printf(" \n") ;
267
268	}
269
270	void AliAnaGammaIsolCut::Terminate(Option_t *)
271	{
272	// The Terminate() function is the last function to be called during
273	// a query. It always runs on the client, it can be used to present
274	// the results graphically or save the results to file.
275
276	}