]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaParticleIsolation.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
AliAnaInsideClusterInvariantMass : Improved calculation of splitted clusters position...
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaParticleIsolation.cxx
CommitLineData
1a31a9ab 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 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 **************************************************************************/
1a31a9ab 15
16//_________________________________________________________________________
17// Class for analysis of particle isolation
18// Input is selected particles put in AOD branch (AliAODPWG4ParticleCorrelation)
19//
20// Class created from old AliPHOSGammaJet
21// (see AliRoot versions previous Release 4-09)
22//
23// -- Author: Gustavo Conesa (LNF-INFN)
24
25//-Yaxian Mao (add the possibility for different IC method with different pt range, 01/10/2010)
26//////////////////////////////////////////////////////////////////////////////
27
28
29// --- ROOT system ---
30#include <TClonesArray.h>
31#include <TList.h>
32#include <TObjString.h>
33#include <TH2F.h>
1a31a9ab 34#include <TClass.h>
35
36// --- Analysis system ---
37#include "AliAnaParticleIsolation.h"
38#include "AliCaloTrackReader.h"
39#include "AliIsolationCut.h"
40#include "AliNeutralMesonSelection.h"
41#include "AliAODPWG4ParticleCorrelation.h"
42#include "AliMCAnalysisUtils.h"
43#include "AliVTrack.h"
44#include "AliVCluster.h"
45
46ClassImp(AliAnaParticleIsolation)
47
803d06a8 48//______________________________________________________________________________
1a31a9ab 49 AliAnaParticleIsolation::AliAnaParticleIsolation() :
09273901 50 AliAnaCaloTrackCorrBaseClass(), fCalorimeter(""),
51 fReMakeIC(0), fMakeSeveralIC(0),
52 fFillTMHisto(0), fFillSSHisto(0),
803d06a8 53 // Several IC
54 fNCones(0), fNPtThresFrac(0),
55 fConeSizes(), fPtThresholds(), fPtFractions(),
56 // Histograms
0fb69ade 57 fhEIso(0), fhPtIso(0),
58 fhPhiIso(0), fhEtaIso(0), fhEtaPhiIso(0),
803d06a8 59 fhPtNoIso(0), fhPtDecayIso(0), fhPtDecayNoIso(0),
60 fhConeSumPt(0), fhPtInCone(0),
61 fhFRConeSumPt(0), fhPtInFRCone(0),
62 // MC histograms
63 fhPtIsoPrompt(0), fhPhiIsoPrompt(0), fhEtaIsoPrompt(0),
64 fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(),
65 fhPtIsoFragmentation(0), fhPhiIsoFragmentation(0), fhEtaIsoFragmentation(0),
66 fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(),
67 fhPtIsoPi0Decay(0), fhPhiIsoPi0Decay(0), fhEtaIsoPi0Decay(0),
68 fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(),
69 fhPtIsoEtaDecay(0), fhPhiIsoEtaDecay(0), fhEtaIsoEtaDecay(0),
70 fhPtThresIsolatedEtaDecay(), fhPtFracIsolatedEtaDecay(), fhPtSumIsolatedEtaDecay(),
71 fhPtIsoOtherDecay(0), fhPhiIsoOtherDecay(0), fhEtaIsoOtherDecay(0),
72 fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(),
73 fhPtIsoConversion(0), fhPhiIsoConversion(0), fhEtaIsoConversion(0),
74 fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(),
75 fhPtIsoUnknown(0), fhPhiIsoUnknown(0), fhEtaIsoUnknown(0),
76 fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(), fhPtSumIsolatedUnknown(),
77 fhPtNoIsoPi0Decay(0), fhPtNoIsoEtaDecay(0), fhPtNoIsoOtherDecay(0),
78 fhPtNoIsoPrompt(0), fhPtIsoMCPhoton(0), fhPtNoIsoMCPhoton(0),
0fb69ade 79 fhPtNoIsoConversion(0), fhPtNoIsoFragmentation(0), fhPtNoIsoUnknown(0),
5c46c992 80 // Cluster control histograms
09273901 81 fhTrackMatchedDEta(0x0), fhTrackMatchedDPhi(0x0), fhTrackMatchedDEtaDPhi(0x0),
31ae6d59 82 fhdEdx(0), fhEOverP(0), fhTrackMatchedMCParticle(0),
09273901 83 fhELambda0(0), fhELambda1(0),
5c46c992 84 fhELambda0TRD(0), fhELambda1TRD(0),
85 // Number of local maxima in cluster
86 fhNLocMax(0),
87 fhELambda0LocMax1(0), fhELambda1LocMax1(0),
88 fhELambda0LocMax2(0), fhELambda1LocMax2(0),
89 fhELambda0LocMaxN(0), fhELambda1LocMaxN(0),
90 // Histograms settings
803d06a8 91 fHistoNPtSumBins(0), fHistoPtSumMax(0.), fHistoPtSumMin(0.),
92 fHistoNPtInConeBins(0), fHistoPtInConeMax(0.), fHistoPtInConeMin(0.)
1a31a9ab 93{
94 //default ctor
95
96 //Initialize parameters
97 InitParameters();
98
b5dbb99b 99 for(Int_t i = 0; i < 5 ; i++)
100 {
803d06a8 101 fConeSizes[i] = 0 ;
1a31a9ab 102 fhPtSumIsolated[i] = 0 ;
103
803d06a8 104 fhPtSumIsolatedPrompt[i] = 0 ;
1a31a9ab 105 fhPtSumIsolatedFragmentation[i] = 0 ;
803d06a8 106 fhPtSumIsolatedPi0Decay[i] = 0 ;
107 fhPtSumIsolatedEtaDecay[i] = 0 ;
108 fhPtSumIsolatedOtherDecay[i] = 0 ;
109 fhPtSumIsolatedConversion[i] = 0 ;
110 fhPtSumIsolatedUnknown[i] = 0 ;
1a31a9ab 111
b5dbb99b 112 for(Int_t j = 0; j < 5 ; j++)
113 {
1a31a9ab 114 fhPtThresIsolated[i][j] = 0 ;
803d06a8 115 fhPtFracIsolated[i][j] = 0 ;
1a31a9ab 116
803d06a8 117 fhPtThresIsolatedPrompt[i][j] = 0 ;
118 fhPtThresIsolatedFragmentation[i][j]= 0 ;
119 fhPtThresIsolatedPi0Decay[i][j] = 0 ;
120 fhPtThresIsolatedEtaDecay[i][j] = 0 ;
121 fhPtThresIsolatedOtherDecay[i][j] = 0 ;
122 fhPtThresIsolatedConversion[i][j] = 0 ;
123 fhPtThresIsolatedUnknown[i][j] = 0 ;
124
125 fhPtFracIsolatedPrompt[i][j] = 0 ;
1a31a9ab 126 fhPtFracIsolatedFragmentation[i][j] = 0 ;
803d06a8 127 fhPtFracIsolatedPi0Decay[i][j] = 0 ;
128 fhPtFracIsolatedEtaDecay[i][j] = 0 ;
129 fhPtFracIsolatedOtherDecay[i][j] = 0 ;
130 fhPtFracIsolatedConversion[i][j] = 0 ;
131 fhPtFracIsolatedUnknown[i][j] = 0 ;
1a31a9ab 132
133 }
134 }
135
136 for(Int_t i = 0; i < 5 ; i++){
b5dbb99b 137 fPtFractions [i] = 0 ;
138 fPtThresholds[i] = 0 ;
1a31a9ab 139 }
140
1a31a9ab 141}
142
b5dbb99b 143//________________________________________________________________________________________________
144void AliAnaParticleIsolation::FillTrackMatchingShowerShapeControlHistograms(
145 const Int_t clusterID,
5c46c992 146 const Int_t nMaxima,
b5dbb99b 147 const Int_t mcTag
148 )
149{
150 // Fill Track matching and Shower Shape control histograms
151
152 if(!fFillTMHisto && !fFillSSHisto) return;
153
154 Int_t iclus = -1;
155 TObjArray* clusters = 0x0;
156 if (fCalorimeter == "EMCAL") clusters = GetEMCALClusters();
157 else if(fCalorimeter == "PHOS" ) clusters = GetPHOSClusters();
158
159 if(clusters)
160 {
161
162 AliVCluster *cluster = FindCluster(clusters,clusterID,iclus);
163 Float_t energy = cluster->E();
164
165 if(fFillSSHisto)
166 {
167 fhELambda0 ->Fill(energy, cluster->GetM02() );
168 fhELambda1 ->Fill(energy, cluster->GetM20() );
169
170 if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5)
171 {
172 fhELambda0TRD ->Fill(energy, cluster->GetM02() );
173 fhELambda1TRD ->Fill(energy, cluster->GetM20() );
174 }
5c46c992 175
176 fhNLocMax->Fill(energy,nMaxima);
177 if (nMaxima==1) { fhELambda0LocMax1->Fill(energy,cluster->GetM02()); fhELambda1LocMax1->Fill(energy,cluster->GetM20()); }
178 else if(nMaxima==2) { fhELambda0LocMax2->Fill(energy,cluster->GetM02()); fhELambda1LocMax2->Fill(energy,cluster->GetM20()); }
179 else { fhELambda0LocMaxN->Fill(energy,cluster->GetM02()); fhELambda1LocMaxN->Fill(energy,cluster->GetM20()); }
180
b5dbb99b 181 } // SS histo fill
182
183
184 if(fFillTMHisto)
185 {
186 Float_t dZ = cluster->GetTrackDz();
187 Float_t dR = cluster->GetTrackDx();
188
189 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
190 {
191 dR = 2000., dZ = 2000.;
192 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
193 }
194
195 //printf("ParticleIsolation: dPhi %f, dEta %f\n",dR,dZ);
196 if(fhTrackMatchedDEta && TMath::Abs(dR) < 999)
197 {
198 fhTrackMatchedDEta->Fill(energy,dZ);
199 fhTrackMatchedDPhi->Fill(energy,dR);
200 if(energy > 0.5) fhTrackMatchedDEtaDPhi->Fill(dZ,dR);
201 }
202
203 // Check dEdx and E/p of matched clusters
204
205 if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
206 {
4bfeae64 207
208 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
b5dbb99b 209
210 if(track)
211 {
212 Float_t dEdx = track->GetTPCsignal();
213 fhdEdx->Fill(cluster->E(), dEdx);
214
215 Float_t eOverp = cluster->E()/track->P();
216 fhEOverP->Fill(cluster->E(), eOverp);
217 }
4bfeae64 218 //else
219 // printf("AliAnaParticleIsolation::FillTrackMatchingShowerShapeHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
220
b5dbb99b 221
222 if(IsDataMC()){
223
224 if ( !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) )
225 {
226 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
227 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle->Fill(energy, 2.5 );
228 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle->Fill(energy, 0.5 );
229 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle->Fill(energy, 1.5 );
230 else fhTrackMatchedMCParticle->Fill(energy, 3.5 );
231
232 }
233 else
234 {
235 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
236 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle->Fill(energy, 6.5 );
237 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle->Fill(energy, 4.5 );
238 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle->Fill(energy, 5.5 );
239 else fhTrackMatchedMCParticle->Fill(energy, 7.5 );
240 }
241
242 } // MC
243
244 } // match window
245
246 }// TM histos fill
247
248 } // clusters array available
249
250}
251
803d06a8 252//______________________________________________________
1a31a9ab 253TObjString * AliAnaParticleIsolation::GetAnalysisCuts()
254{
255 //Save parameters used for analysis
256 TString parList ; //this will be list of parameters used for this analysis.
257 const Int_t buffersize = 255;
258 char onePar[buffersize] ;
259
260 snprintf(onePar, buffersize,"--- AliAnaParticleIsolation ---\n") ;
261 parList+=onePar ;
262 snprintf(onePar, buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
263 parList+=onePar ;
264 snprintf(onePar, buffersize,"fReMakeIC =%d (Flag for reisolation during histogram filling) \n",fReMakeIC) ;
265 parList+=onePar ;
266 snprintf(onePar, buffersize,"fMakeSeveralIC=%d (Flag for isolation with several cuts at the same time ) \n",fMakeSeveralIC) ;
09273901 267 parList+=onePar ;
268 snprintf(onePar, buffersize,"fFillTMHisto=%d (Flag for track matching histograms) \n",fFillTMHisto) ;
1a31a9ab 269 parList+=onePar ;
09273901 270 snprintf(onePar, buffersize,"fFillSSHisto=%d (Flag for shower shape histograms) \n",fFillSSHisto) ;
271 parList+=onePar ;
272
b5dbb99b 273 if(fMakeSeveralIC)
274 {
1a31a9ab 275 snprintf(onePar, buffersize,"fNCones =%d (Number of cone sizes) \n",fNCones) ;
276 parList+=onePar ;
277 snprintf(onePar, buffersize,"fNPtThresFrac=%d (Flag for isolation with several cuts at the same time ) \n",fNPtThresFrac) ;
278 parList+=onePar ;
279
b5dbb99b 280 for(Int_t icone = 0; icone < fNCones ; icone++)
281 {
1a31a9ab 282 snprintf(onePar, buffersize,"fConeSizes[%d]=%1.2f (isolation cone size) \n",icone, fConeSizes[icone]) ;
283 parList+=onePar ;
284 }
b5dbb99b 285 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
286 {
1a31a9ab 287 snprintf(onePar, buffersize,"fPtThresholds[%d]=%1.2f (isolation pt threshold) \n",ipt, fPtThresholds[ipt]) ;
288 parList+=onePar ;
289 }
b5dbb99b 290 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
291 {
1a31a9ab 292 snprintf(onePar, buffersize,"fPtFractions[%d]=%1.2f (isolation pt fraction threshold) \n",ipt, fPtFractions[ipt]) ;
293 parList+=onePar ;
294 }
295 }
296
297 //Get parameters set in base class.
298 parList += GetBaseParametersList() ;
299
300 //Get parameters set in IC class.
301 if(!fMakeSeveralIC)parList += GetIsolationCut()->GetICParametersList() ;
302
303 return new TObjString(parList) ;
304
305}
306
803d06a8 307//________________________________________________________
1a31a9ab 308TList * AliAnaParticleIsolation::GetCreateOutputObjects()
309{
310 // Create histograms to be saved in output file and
311 // store them in outputContainer
312 TList * outputContainer = new TList() ;
313 outputContainer->SetName("IsolatedParticleHistos") ;
314
745913ae 315 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins();
316 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins();
317 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
318 Float_t ptmax = GetHistogramRanges()->GetHistoPtMax();
319 Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();
320 Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
321 Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
322 Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
323 Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
09273901 324 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins();
325 Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax();
326 Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
327
328 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
329 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
330 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
331 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
332 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
333 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
803d06a8 334
31ae6d59 335 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins();
336 Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax();
337 Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
338 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
339 Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax();
340 Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
341
342
803d06a8 343 Int_t nptsumbins = fHistoNPtSumBins;
344 Float_t ptsummax = fHistoPtSumMax;
345 Float_t ptsummin = fHistoPtSumMin;
346 Int_t nptinconebins = fHistoNPtInConeBins;
347 Float_t ptinconemax = fHistoPtInConeMax;
348 Float_t ptinconemin = fHistoPtInConeMin;
1a31a9ab 349
b5dbb99b 350 if(!fMakeSeveralIC)
351 {
1a31a9ab 352
b5dbb99b 353 if(fFillTMHisto)
354 {
09273901 355 fhTrackMatchedDEta = new TH2F
31ae6d59 356 ("hTrackMatchedDEta",
09273901 357 "d#eta of cluster-track vs cluster energy",
358 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
359 fhTrackMatchedDEta->SetYTitle("d#eta");
360 fhTrackMatchedDEta->SetXTitle("E_{cluster} (GeV)");
361
362 fhTrackMatchedDPhi = new TH2F
31ae6d59 363 ("hTrackMatchedDPhi",
09273901 364 "d#phi of cluster-track vs cluster energy",
365 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
366 fhTrackMatchedDPhi->SetYTitle("d#phi (rad)");
367 fhTrackMatchedDPhi->SetXTitle("E_{cluster} (GeV)");
368
369 fhTrackMatchedDEtaDPhi = new TH2F
31ae6d59 370 ("hTrackMatchedDEtaDPhi",
09273901 371 "d#eta vs d#phi of cluster-track vs cluster energy",
372 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
373 fhTrackMatchedDEtaDPhi->SetYTitle("d#phi (rad)");
374 fhTrackMatchedDEtaDPhi->SetXTitle("d#eta");
375
376 outputContainer->Add(fhTrackMatchedDEta) ;
377 outputContainer->Add(fhTrackMatchedDPhi) ;
378 outputContainer->Add(fhTrackMatchedDEtaDPhi) ;
31ae6d59 379
380 fhdEdx = new TH2F ("hdEdx","matched track <dE/dx> vs cluster E ", nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
381 fhdEdx->SetXTitle("E (GeV)");
382 fhdEdx->SetYTitle("<dE/dx>");
383 outputContainer->Add(fhdEdx);
384
385 fhEOverP = new TH2F ("hEOverP","matched track E/p vs cluster E ", nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
386 fhEOverP->SetXTitle("E (GeV)");
387 fhEOverP->SetYTitle("E/p");
388 outputContainer->Add(fhEOverP);
389
390 if(IsDataMC())
391 {
392 fhTrackMatchedMCParticle = new TH2F
393 ("hTrackMatchedMCParticle",
394 "Origin of particle vs energy",
395 nptbins,ptmin,ptmax,8,0,8);
396 fhTrackMatchedMCParticle->SetXTitle("E (GeV)");
397 //fhTrackMatchedMCParticle->SetYTitle("Particle type");
398
399 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(1 ,"Photon");
400 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(2 ,"Electron");
401 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
402 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(4 ,"Rest");
403 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
404 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
405 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
406 fhTrackMatchedMCParticle->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
407
408 outputContainer->Add(fhTrackMatchedMCParticle);
409 }
09273901 410 }
411
b5dbb99b 412 if(fFillSSHisto)
413 {
09273901 414 fhELambda0 = new TH2F
415 ("hELambda0","Selected #pi^{0} pairs: E vs #lambda_{0}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
416 fhELambda0->SetYTitle("#lambda_{0}^{2}");
417 fhELambda0->SetXTitle("E (GeV)");
418 outputContainer->Add(fhELambda0) ;
419
420 fhELambda1 = new TH2F
421 ("hELambda1","Selected #pi^{0} pairs: E vs #lambda_{1}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
422 fhELambda1->SetYTitle("#lambda_{1}^{2}");
423 fhELambda1->SetXTitle("E (GeV)");
b5dbb99b 424 outputContainer->Add(fhELambda1) ;
425
426 if(fCalorimeter=="EMCAL")
427 {
428 fhELambda0TRD = new TH2F
429 ("hELambda0TRD","Selected #pi^{0} pairs: E vs #lambda_{0}, SM behind TRD",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
430 fhELambda0TRD->SetYTitle("#lambda_{0}^{2}");
431 fhELambda0TRD->SetXTitle("E (GeV)");
432 outputContainer->Add(fhELambda0TRD) ;
433
434 fhELambda1TRD = new TH2F
435 ("hELambda1TRD","Selected #pi^{0} pairs: E vs #lambda_{1}, SM behind TRD",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
436 fhELambda1TRD->SetYTitle("#lambda_{1}^{2}");
437 fhELambda1TRD->SetXTitle("E (GeV)");
438 outputContainer->Add(fhELambda1TRD) ;
439 }
5c46c992 440
441 fhNLocMax = new TH2F("hNLocMax","Number of local maxima in cluster",
442 nptbins,ptmin,ptmax,10,0,10);
443 fhNLocMax ->SetYTitle("N maxima");
444 fhNLocMax ->SetXTitle("E (GeV)");
445 outputContainer->Add(fhNLocMax) ;
446
447 fhELambda0LocMax1 = new TH2F
448 ("hELambda0LocMax1","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, 1 Local maxima",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
449 fhELambda0LocMax1->SetYTitle("#lambda_{0}^{2}");
450 fhELambda0LocMax1->SetXTitle("E (GeV)");
451 outputContainer->Add(fhELambda0LocMax1) ;
452
453 fhELambda1LocMax1 = new TH2F
454 ("hELambda1LocMax1","Selected #pi^{0} (#eta) pairs: E vs #lambda_{1}, 1 Local maxima",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
455 fhELambda1LocMax1->SetYTitle("#lambda_{1}^{2}");
456 fhELambda1LocMax1->SetXTitle("E (GeV)");
457 outputContainer->Add(fhELambda1LocMax1) ;
458
459 fhELambda0LocMax2 = new TH2F
460 ("hELambda0LocMax2","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, 2 Local maxima",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
461 fhELambda0LocMax2->SetYTitle("#lambda_{0}^{2}");
462 fhELambda0LocMax2->SetXTitle("E (GeV)");
463 outputContainer->Add(fhELambda0LocMax2) ;
464
465 fhELambda1LocMax2 = new TH2F
466 ("hELambda1LocMax2","Selected #pi^{0} (#eta) pairs: E vs #lambda_{1}, 2 Local maxima",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
467 fhELambda1LocMax2->SetYTitle("#lambda_{1}^{2}");
468 fhELambda1LocMax2->SetXTitle("E (GeV)");
469 outputContainer->Add(fhELambda1LocMax2) ;
470
471 fhELambda0LocMaxN = new TH2F
472 ("hELambda0LocMaxN","Selected #pi^{0} (#eta) pairs: E vs #lambda_{0}, N>2 Local maxima",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
473 fhELambda0LocMaxN->SetYTitle("#lambda_{0}^{2}");
474 fhELambda0LocMaxN->SetXTitle("E (GeV)");
475 outputContainer->Add(fhELambda0LocMaxN) ;
476
477 fhELambda1LocMaxN = new TH2F
478 ("hELambda1LocMaxN","Selected #pi^{0} (#eta) pairs: E vs #lambda_{1}, N>2 Local maxima",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
479 fhELambda1LocMaxN->SetYTitle("#lambda_{1}^{2}");
480 fhELambda1LocMaxN->SetXTitle("E (GeV)");
481 outputContainer->Add(fhELambda1LocMaxN) ;
482
09273901 483 }
484
1a31a9ab 485 fhConeSumPt = new TH2F
486 ("hConePtSum","#Sigma p_{T} in isolation cone ",nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
487 fhConeSumPt->SetYTitle("#Sigma p_{T}");
488 fhConeSumPt->SetXTitle("p_{T} (GeV/c)");
489 outputContainer->Add(fhConeSumPt) ;
490
491 fhPtInCone = new TH2F
492 ("hPtInCone","p_{T} in isolation cone ",nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
493 fhPtInCone->SetYTitle("p_{T in cone} (GeV/c)");
494 fhPtInCone->SetXTitle("p_{T} (GeV/c)");
495 outputContainer->Add(fhPtInCone) ;
496
497 fhFRConeSumPt = new TH2F
498 ("hFRConePtSum","#Sigma p_{T} in the froward region isolation cone ",nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
499 fhFRConeSumPt->SetYTitle("#Sigma p_{T}");
500 fhFRConeSumPt->SetXTitle("p_{T} (GeV/c)");
501 outputContainer->Add(fhFRConeSumPt) ;
502
503 fhPtInFRCone = new TH2F
504 ("hPtInFRCone","p_{T} in froward region isolation cone ",nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
505 fhPtInFRCone->SetYTitle("p_{T in cone} (GeV/c)");
506 fhPtInFRCone->SetXTitle("p_{T} (GeV/c)");
507 outputContainer->Add(fhPtInFRCone) ;
508
0fb69ade 509 fhEIso = new TH1F("hE","Number of isolated particles vs E",nptbins,ptmin,ptmax);
510 fhEIso->SetYTitle("dN / dE");
511 fhEIso->SetXTitle("E (GeV/c)");
512 outputContainer->Add(fhEIso) ;
513
514 fhPtIso = new TH1F("hPt","Number of isolated particles vs p_{T}",nptbins,ptmin,ptmax);
515 fhPtIso->SetYTitle("dN / p_{T}");
516 fhPtIso->SetXTitle("p_{T} (GeV/c)");
1a31a9ab 517 outputContainer->Add(fhPtIso) ;
518
519 fhPhiIso = new TH2F
732895a6 520 ("hPhi","Number of isolated particles",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 521 fhPhiIso->SetYTitle("#phi");
522 fhPhiIso->SetXTitle("p_{T} (GeV/c)");
523 outputContainer->Add(fhPhiIso) ;
524
525 fhEtaIso = new TH2F
732895a6 526 ("hEta","Number of isolated particles",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 527 fhEtaIso->SetYTitle("#eta");
528 fhEtaIso->SetXTitle("p_{T} (GeV/c)");
529 outputContainer->Add(fhEtaIso) ;
530
0fb69ade 531 fhEtaPhiIso = new TH2F
532 ("hEtaPhi","Number of isolated particlesm #eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
533 fhEtaPhiIso->SetXTitle("#eta");
534 fhEtaPhiIso->SetYTitle("#phi");
535 outputContainer->Add(fhEtaPhiIso) ;
536
1a31a9ab 537 fhPtNoIso = new TH1F("hPtNoIso","Number of not isolated leading particles",nptbins,ptmin,ptmax);
538 fhPtNoIso->SetYTitle("N");
539 fhPtNoIso->SetXTitle("p_{T}(GeV/c)");
540 outputContainer->Add(fhPtNoIso) ;
541
803d06a8 542 fhPtDecayIso = new TH1F("hPtDecayIso","Number of isolated #pi^{0} decay particles",nptbins,ptmin,ptmax);
66e64043 543 fhPtDecayIso->SetYTitle("N");
544 fhPtDecayIso->SetXTitle("p_{T}(GeV/c)");
803d06a8 545 outputContainer->Add(fhPtDecayIso) ;
1a31a9ab 546
803d06a8 547 fhPtDecayNoIso = new TH1F("hPtDecayNoIso","Number of not isolated leading pi0 decay particles",nptbins,ptmin,ptmax);
66e64043 548 fhPtDecayNoIso->SetYTitle("N");
549 fhPtDecayNoIso->SetXTitle("p_{T}(GeV/c)");
803d06a8 550 outputContainer->Add(fhPtDecayNoIso) ;
1a31a9ab 551
b5dbb99b 552 if(IsDataMC())
553 {
732895a6 554 fhPtIsoPrompt = new TH1F("hPtMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax);
1a31a9ab 555 fhPtIsoPrompt->SetYTitle("N");
556 fhPtIsoPrompt->SetXTitle("p_{T #gamma}(GeV/c)");
557 outputContainer->Add(fhPtIsoPrompt) ;
558
559 fhPhiIsoPrompt = new TH2F
732895a6 560 ("hPhiMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 561 fhPhiIsoPrompt->SetYTitle("#phi");
562 fhPhiIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
563 outputContainer->Add(fhPhiIsoPrompt) ;
564
565 fhEtaIsoPrompt = new TH2F
732895a6 566 ("hEtaMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 567 fhEtaIsoPrompt->SetYTitle("#eta");
568 fhEtaIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
569 outputContainer->Add(fhEtaIsoPrompt) ;
570
732895a6 571 fhPtIsoFragmentation = new TH1F("hPtMCFragmentation","Number of isolated #gamma",nptbins,ptmin,ptmax);
1a31a9ab 572 fhPtIsoFragmentation->SetYTitle("N");
573 fhPtIsoFragmentation->SetXTitle("p_{T #gamma}(GeV/c)");
574 outputContainer->Add(fhPtIsoFragmentation) ;
575
576 fhPhiIsoFragmentation = new TH2F
732895a6 577 ("hPhiMCFragmentation","Number of isolated fragmentation #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 578 fhPhiIsoFragmentation->SetYTitle("#phi");
579 fhPhiIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
580 outputContainer->Add(fhPhiIsoFragmentation) ;
581
582 fhEtaIsoFragmentation = new TH2F
732895a6 583 ("hEtaMCFragmentation","Number of isolated fragmentation #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 584 fhEtaIsoFragmentation->SetYTitle("#eta");
585 fhEtaIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
586 outputContainer->Add(fhEtaIsoFragmentation) ;
587
732895a6 588 fhPtIsoPi0Decay = new TH1F("hPtMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax);
1a31a9ab 589 fhPtIsoPi0Decay->SetYTitle("N");
590 fhPtIsoPi0Decay->SetXTitle("p_{T #gamma}(GeV/c)");
591 outputContainer->Add(fhPtIsoPi0Decay) ;
592
593 fhPhiIsoPi0Decay = new TH2F
732895a6 594 ("hPhiMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 595 fhPhiIsoPi0Decay->SetYTitle("#phi");
596 fhPhiIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
597 outputContainer->Add(fhPhiIsoPi0Decay) ;
598
599 fhEtaIsoPi0Decay = new TH2F
732895a6 600 ("hEtaMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 601 fhEtaIsoPi0Decay->SetYTitle("#eta");
602 fhEtaIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
603 outputContainer->Add(fhEtaIsoPi0Decay) ;
604
803d06a8 605 fhPtIsoEtaDecay = new TH1F("hPtMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax);
606 fhPtIsoEtaDecay->SetYTitle("N");
607 fhPtIsoEtaDecay->SetXTitle("p_{T #gamma}(GeV/c)");
608 outputContainer->Add(fhPtIsoEtaDecay) ;
609
610 fhPhiIsoEtaDecay = new TH2F
611 ("hPhiMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
612 fhPhiIsoEtaDecay->SetYTitle("#phi");
613 fhPhiIsoEtaDecay->SetXTitle("p_{T #gamma} (GeV/c)");
614 outputContainer->Add(fhPhiIsoEtaDecay) ;
615
616 fhEtaIsoEtaDecay = new TH2F
617 ("hEtaMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
618 fhEtaIsoEtaDecay->SetYTitle("#eta");
619 fhEtaIsoEtaDecay->SetXTitle("p_{T #gamma} (GeV/c)");
620 outputContainer->Add(fhEtaIsoEtaDecay) ;
621
732895a6 622 fhPtIsoOtherDecay = new TH1F("hPtMCOtherDecay","Number of isolated #gamma from non-#pi^{0} decay",nptbins,ptmin,ptmax);
1a31a9ab 623 fhPtIsoOtherDecay->SetYTitle("N");
624 fhPtIsoOtherDecay->SetXTitle("p_{T #gamma}(GeV/c)");
625 outputContainer->Add(fhPtIsoOtherDecay) ;
626
627 fhPhiIsoOtherDecay = new TH2F
732895a6 628 ("hPhiMCOtherDecay","Number of isolated #gamma from non-#pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 629 fhPhiIsoOtherDecay->SetYTitle("#phi");
630 fhPhiIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
631 outputContainer->Add(fhPhiIsoOtherDecay) ;
632
633 fhEtaIsoOtherDecay = new TH2F
732895a6 634 ("hEtaMCOtherDecay","Number of isolated #gamma non-#pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 635 fhEtaIsoOtherDecay->SetYTitle("#eta");
636 fhEtaIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
637 outputContainer->Add(fhEtaIsoOtherDecay) ;
638
732895a6 639 fhPtIsoConversion = new TH1F("hPtMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax);
1a31a9ab 640 fhPtIsoConversion->SetYTitle("N");
641 fhPtIsoConversion->SetXTitle("p_{T #gamma}(GeV/c)");
642 outputContainer->Add(fhPtIsoConversion) ;
643
644 fhPhiIsoConversion = new TH2F
732895a6 645 ("hPhiMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 646 fhPhiIsoConversion->SetYTitle("#phi");
647 fhPhiIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
648 outputContainer->Add(fhPhiIsoConversion) ;
649
650 fhEtaIsoConversion = new TH2F
732895a6 651 ("hEtaMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 652 fhEtaIsoConversion->SetYTitle("#eta");
653 fhEtaIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
654 outputContainer->Add(fhEtaIsoConversion) ;
655
732895a6 656 fhPtIsoUnknown = new TH1F("hPtMCUnknown","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax);
1a31a9ab 657 fhPtIsoUnknown->SetYTitle("N");
658 fhPtIsoUnknown->SetXTitle("p_{T}(GeV/c)");
659 outputContainer->Add(fhPtIsoUnknown) ;
660
661 fhPhiIsoUnknown = new TH2F
732895a6 662 ("hPhiMCUnknown","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 663 fhPhiIsoUnknown->SetYTitle("#phi");
664 fhPhiIsoUnknown->SetXTitle("p_{T} (GeV/c)");
665 outputContainer->Add(fhPhiIsoUnknown) ;
666
667 fhEtaIsoUnknown = new TH2F
732895a6 668 ("hEtaMCUnknown","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 669 fhEtaIsoUnknown->SetYTitle("#eta");
670 fhEtaIsoUnknown->SetXTitle("p_{T} (GeV/c)");
671 outputContainer->Add(fhEtaIsoUnknown) ;
672
673 fhPtNoIsoPi0Decay = new TH1F
732895a6 674 ("hPtNoIsoPi0Decay","Number of not isolated leading #gamma from #pi^{0} decay",nptbins,ptmin,ptmax);
1a31a9ab 675 fhPtNoIsoPi0Decay->SetYTitle("N");
676 fhPtNoIsoPi0Decay->SetXTitle("p_{T} (GeV/c)");
677 outputContainer->Add(fhPtNoIsoPi0Decay) ;
678
803d06a8 679 fhPtNoIsoEtaDecay = new TH1F
680 ("hPtNoIsoEtaDecay","Number of not isolated leading #gamma from eta decay",nptbins,ptmin,ptmax);
681 fhPtNoIsoEtaDecay->SetYTitle("N");
682 fhPtNoIsoEtaDecay->SetXTitle("p_{T} (GeV/c)");
683 outputContainer->Add(fhPtNoIsoEtaDecay) ;
684
685 fhPtNoIsoOtherDecay = new TH1F
686 ("hPtNoIsoOtherDecay","Number of not isolated leading #gamma from other decay",nptbins,ptmin,ptmax);
687 fhPtNoIsoOtherDecay->SetYTitle("N");
688 fhPtNoIsoOtherDecay->SetXTitle("p_{T} (GeV/c)");
689 outputContainer->Add(fhPtNoIsoOtherDecay) ;
690
1a31a9ab 691 fhPtNoIsoPrompt = new TH1F
692 ("hPtNoIsoPrompt","Number of not isolated leading prompt #gamma",nptbins,ptmin,ptmax);
693 fhPtNoIsoPrompt->SetYTitle("N");
694 fhPtNoIsoPrompt->SetXTitle("p_{T} (GeV/c)");
695 outputContainer->Add(fhPtNoIsoPrompt) ;
696
697 fhPtIsoMCPhoton = new TH1F
698 ("hPtIsoMCPhoton","Number of isolated leading #gamma",nptbins,ptmin,ptmax);
699 fhPtIsoMCPhoton->SetYTitle("N");
700 fhPtIsoMCPhoton->SetXTitle("p_{T} (GeV/c)");
701 outputContainer->Add(fhPtIsoMCPhoton) ;
702
703 fhPtNoIsoMCPhoton = new TH1F
732895a6 704 ("hPtNoIsoMCPhoton","Number of not isolated leading #gamma",nptbins,ptmin,ptmax);
1a31a9ab 705 fhPtNoIsoMCPhoton->SetYTitle("N");
706 fhPtNoIsoMCPhoton->SetXTitle("p_{T} (GeV/c)");
707 outputContainer->Add(fhPtNoIsoMCPhoton) ;
0fb69ade 708
709 fhPtNoIsoConversion = new TH1F
710 ("hPtNoIsoConversion","Number of not isolated leading conversion #gamma",nptbins,ptmin,ptmax);
711 fhPtNoIsoConversion->SetYTitle("N");
712 fhPtNoIsoConversion->SetXTitle("p_{T} (GeV/c)");
713 outputContainer->Add(fhPtNoIsoConversion) ;
714
715 fhPtNoIsoFragmentation = new TH1F
716 ("hPtNoIsoFragmentation","Number of not isolated leading fragmentation #gamma",nptbins,ptmin,ptmax);
717 fhPtNoIsoFragmentation->SetYTitle("N");
718 fhPtNoIsoFragmentation->SetXTitle("p_{T} (GeV/c)");
719 outputContainer->Add(fhPtNoIsoFragmentation) ;
720
721 fhPtNoIsoUnknown = new TH1F
722 ("hPtNoIsoUnknown","Number of not isolated leading hadrons",nptbins,ptmin,ptmax);
723 fhPtNoIsoUnknown->SetYTitle("N");
724 fhPtNoIsoUnknown->SetXTitle("p_{T} (GeV/c)");
725 outputContainer->Add(fhPtNoIsoUnknown) ;
726
1a31a9ab 727 }//Histos with MC
728
729 }
730
b5dbb99b 731 if(fMakeSeveralIC)
732 {
1a31a9ab 733 const Int_t buffersize = 255;
734 char name[buffersize];
735 char title[buffersize];
736 for(Int_t icone = 0; icone<fNCones; icone++){
737 snprintf(name, buffersize,"hPtSum_Cone_%d",icone);
738 snprintf(title, buffersize,"Candidate cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
739 fhPtSumIsolated[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
740 fhPtSumIsolated[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
741 fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)");
742 outputContainer->Add(fhPtSumIsolated[icone]) ;
743
744 if(IsDataMC()){
745 snprintf(name, buffersize,"hPtSumPrompt_Cone_%d",icone);
746 snprintf(title, buffersize,"Candidate Prompt cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
747 fhPtSumIsolatedPrompt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
748 fhPtSumIsolatedPrompt[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
749 fhPtSumIsolatedPrompt[icone]->SetXTitle("p_{T} (GeV/c)");
750 outputContainer->Add(fhPtSumIsolatedPrompt[icone]) ;
751
752 snprintf(name, buffersize,"hPtSumFragmentation_Cone_%d",icone);
753 snprintf(title, buffersize,"Candidate Fragmentation cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
754 fhPtSumIsolatedFragmentation[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
755 fhPtSumIsolatedFragmentation[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
756 fhPtSumIsolatedFragmentation[icone]->SetXTitle("p_{T} (GeV/c)");
757 outputContainer->Add(fhPtSumIsolatedFragmentation[icone]) ;
758
759 snprintf(name, buffersize,"hPtSumPi0Decay_Cone_%d",icone);
760 snprintf(title, buffersize,"Candidate Pi0Decay cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
761 fhPtSumIsolatedPi0Decay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
762 fhPtSumIsolatedPi0Decay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
763 fhPtSumIsolatedPi0Decay[icone]->SetXTitle("p_{T} (GeV/c)");
764 outputContainer->Add(fhPtSumIsolatedPi0Decay[icone]) ;
765
803d06a8 766 snprintf(name, buffersize,"hPtSumEtaDecay_Cone_%d",icone);
767 snprintf(title, buffersize,"Candidate EtaDecay cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
768 fhPtSumIsolatedEtaDecay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
769 fhPtSumIsolatedEtaDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
770 fhPtSumIsolatedEtaDecay[icone]->SetXTitle("p_{T} (GeV/c)");
771 outputContainer->Add(fhPtSumIsolatedEtaDecay[icone]) ;
772
1a31a9ab 773 snprintf(name, buffersize,"hPtSumOtherDecay_Cone_%d",icone);
774 snprintf(title, buffersize,"Candidate OtherDecay cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
775 fhPtSumIsolatedOtherDecay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
776 fhPtSumIsolatedOtherDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
777 fhPtSumIsolatedOtherDecay[icone]->SetXTitle("p_{T} (GeV/c)");
778 outputContainer->Add(fhPtSumIsolatedOtherDecay[icone]) ;
779
780 snprintf(name, buffersize,"hPtSumConversion_Cone_%d",icone);
781 snprintf(title, buffersize,"Candidate Conversion cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
782 fhPtSumIsolatedConversion[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
783 fhPtSumIsolatedConversion[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
784 fhPtSumIsolatedConversion[icone]->SetXTitle("p_{T} (GeV/c)");
785 outputContainer->Add(fhPtSumIsolatedConversion[icone]) ;
786
787 snprintf(name, buffersize,"hPtSumUnknown_Cone_%d",icone);
788 snprintf(title, buffersize,"Candidate Unknown cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
789 fhPtSumIsolatedUnknown[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
790 fhPtSumIsolatedUnknown[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
791 fhPtSumIsolatedUnknown[icone]->SetXTitle("p_{T} (GeV/c)");
792 outputContainer->Add(fhPtSumIsolatedUnknown[icone]) ;
793
794 }//Histos with MC
795
b5dbb99b 796 for(Int_t ipt = 0; ipt<fNPtThresFrac;ipt++)
797 {
1a31a9ab 798 snprintf(name, buffersize,"hPtThres_Cone_%d_Pt%d",icone,ipt);
799 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
800 fhPtThresIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
801 fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
802 outputContainer->Add(fhPtThresIsolated[icone][ipt]) ;
803
804 snprintf(name, buffersize,"hPtFrac_Cone_%d_Pt%d",icone,ipt);
805 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
806 fhPtFracIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
807 fhPtFracIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
808 outputContainer->Add(fhPtFracIsolated[icone][ipt]) ;
809
b5dbb99b 810 if(IsDataMC())
811 {
1a31a9ab 812 snprintf(name, buffersize,"hPtThresMCPrompt_Cone_%d_Pt%d",icone,ipt);
813 snprintf(title, buffersize,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
814 fhPtThresIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
815 fhPtThresIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
816 outputContainer->Add(fhPtThresIsolatedPrompt[icone][ipt]) ;
817
818 snprintf(name, buffersize,"hPtFracMCPrompt_Cone_%d_Pt%d",icone,ipt);
819 snprintf(title, buffersize,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
820 fhPtFracIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
821 fhPtFracIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
822 outputContainer->Add(fhPtFracIsolatedPrompt[icone][ipt]) ;
823
824 snprintf(name, buffersize,"hPtThresMCFragmentation_Cone_%d_Pt%d",icone,ipt);
825 snprintf(title, buffersize,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
826 fhPtThresIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
827 fhPtThresIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
828 outputContainer->Add(fhPtThresIsolatedFragmentation[icone][ipt]) ;
829
830 snprintf(name, buffersize,"hPtFracMCFragmentation_Cone_%d_Pt%d",icone,ipt);
831 snprintf(title, buffersize,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
832 fhPtFracIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
833 fhPtFracIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
834 outputContainer->Add(fhPtFracIsolatedFragmentation[icone][ipt]) ;
835
836 snprintf(name, buffersize,"hPtThresMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
837 snprintf(title, buffersize,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
838 fhPtThresIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
839 fhPtThresIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
840 outputContainer->Add(fhPtThresIsolatedPi0Decay[icone][ipt]) ;
841
842 snprintf(name, buffersize,"hPtFracMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
843 snprintf(title, buffersize,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
844 fhPtFracIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
845 fhPtFracIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
846 outputContainer->Add(fhPtFracIsolatedPi0Decay[icone][ipt]) ;
847
803d06a8 848 snprintf(name, buffersize,"hPtThresMCEtaDecay_Cone_%d_Pt%d",icone,ipt);
849 snprintf(title, buffersize,"Isolated candidate EtaDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
850 fhPtThresIsolatedEtaDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
851 fhPtThresIsolatedEtaDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
852 outputContainer->Add(fhPtThresIsolatedEtaDecay[icone][ipt]) ;
853
854 snprintf(name, buffersize,"hPtFracMCEtaDecay_Cone_%d_Pt%d",icone,ipt);
855 snprintf(title, buffersize,"Isolated candidate EtaDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
856 fhPtFracIsolatedEtaDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
857 fhPtFracIsolatedEtaDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
858 outputContainer->Add(fhPtFracIsolatedEtaDecay[icone][ipt]) ;
859
860
1a31a9ab 861 snprintf(name, buffersize,"hPtThresMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
862 snprintf(title, buffersize,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
863 fhPtThresIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
864 fhPtThresIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
865 outputContainer->Add(fhPtThresIsolatedOtherDecay[icone][ipt]) ;
866
867 snprintf(name, buffersize,"hPtFracMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
868 snprintf(title, buffersize,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
869 fhPtFracIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
870 fhPtFracIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
871 outputContainer->Add(fhPtFracIsolatedOtherDecay[icone][ipt]) ;
872
873 snprintf(name, buffersize,"hPtThresMCConversion_Cone_%d_Pt%d",icone,ipt);
874 snprintf(title, buffersize,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
875 fhPtThresIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
876 fhPtThresIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
877 outputContainer->Add(fhPtThresIsolatedConversion[icone][ipt]) ;
878
879 snprintf(name, buffersize,"hPtFracMCConversion_Cone_%d_Pt%d",icone,ipt);
880 snprintf(title, buffersize,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
881 fhPtFracIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
882 fhPtFracIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
883 outputContainer->Add(fhPtFracIsolatedConversion[icone][ipt]) ;
884
885 snprintf(name, buffersize,"hPtThresMCUnknown_Cone_%d_Pt%d",icone,ipt);
886 snprintf(title, buffersize,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
887 fhPtThresIsolatedUnknown[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
888 fhPtThresIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
889 outputContainer->Add(fhPtThresIsolatedUnknown[icone][ipt]) ;
890
891 snprintf(name, buffersize,"hPtFracMCUnknown_Cone_%d_Pt%d",icone,ipt);
892 snprintf(title, buffersize,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
893 fhPtFracIsolatedUnknown[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
894 fhPtFracIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
895 outputContainer->Add(fhPtFracIsolatedUnknown[icone][ipt]) ;
896
897 }//Histos with MC
898
899 }//icone loop
900 }//ipt loop
901 }
902
1a31a9ab 903
904 return outputContainer ;
905
906}
907
803d06a8 908//____________________________________________
909void AliAnaParticleIsolation::InitParameters()
910{
911
912 //Initialize the parameters of the analysis.
913 SetInputAODName("PWG4Particle");
914 SetAODObjArrayName("IsolationCone");
915 AddToHistogramsName("AnaIsolation_");
916
917 fCalorimeter = "PHOS" ;
918 fReMakeIC = kFALSE ;
919 fMakeSeveralIC = kFALSE ;
920
921 //----------- Several IC-----------------
922 fNCones = 5 ;
923 fNPtThresFrac = 5 ;
924 fConeSizes[0] = 0.1; fConeSizes[1] = 0.2; fConeSizes[2] = 0.3; fConeSizes[3] = 0.4; fConeSizes[4] = 0.5;
925 fPtThresholds[0] = 1.; fPtThresholds[1] = 2.; fPtThresholds[2] = 3.; fPtThresholds[3] = 4.; fPtThresholds[4] = 5.;
926 fPtFractions[0] = 0.05; fPtFractions[1] = 0.075; fPtFractions[2] = 0.1; fPtFractions[3] = 1.25; fPtFractions[4] = 1.5;
927
928 //------------- Histograms settings -------
929 fHistoNPtSumBins = 100 ;
930 fHistoPtSumMax = 50 ;
931 fHistoPtSumMin = 0. ;
932
933 fHistoNPtInConeBins = 100 ;
934 fHistoPtInConeMax = 50 ;
935 fHistoPtInConeMin = 0. ;
936
937}
938
939//__________________________________________________
1a31a9ab 940void AliAnaParticleIsolation::MakeAnalysisFillAOD()
941{
942 //Do analysis and fill aods
943 //Search for the isolated photon in fCalorimeter with pt > GetMinPt()
944
b5dbb99b 945 if(!GetInputAODBranch())
946 {
1a31a9ab 947 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
948 abort();
949 }
950
b5dbb99b 951 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation"))
952 {
1a31a9ab 953 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
954 abort();
955 }
956
957 Int_t n = 0, nfrac = 0;
958 Bool_t isolated = kFALSE ;
1a31a9ab 959 Float_t coneptsum = 0 ;
960 TObjArray * pl = 0x0; ;
961
962 //Select the calorimeter for candidate isolation with neutral particles
b5dbb99b 963 if (fCalorimeter == "PHOS" )
1a31a9ab 964 pl = GetPHOSClusters();
965 else if (fCalorimeter == "EMCAL")
966 pl = GetEMCALClusters();
967
968 //Loop on AOD branch, filled previously in AliAnaPhoton, find leading particle to do isolation only with it
969 Double_t ptLeading = 0. ;
970 Int_t idLeading = -1 ;
971 TLorentzVector mom ;
972 Int_t naod = GetInputAODBranch()->GetEntriesFast();
973 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Input aod branch entries %d\n", naod);
974
b5dbb99b 975 for(Int_t iaod = 0; iaod < naod; iaod++)
976 {
1a31a9ab 977 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
978
979 //If too small or too large pt, skip
980 if(aodinput->Pt() < GetMinPt() || aodinput->Pt() > GetMaxPt() ) continue ;
981
982 //check if it is low pt trigger particle, then adjust the isolation method
b5dbb99b 983 if(aodinput->Pt() < GetIsolationCut()->GetPtThreshold() ||
984 aodinput->Pt() < GetIsolationCut()->GetSumPtThreshold())
985 {
1a31a9ab 986 continue ; //trigger should not come from underlying event
b5dbb99b 987 }
1a31a9ab 988
989 //vertex cut in case of mixing
04f7a616 990 Int_t check = CheckMixedEventVertex(aodinput->GetCaloLabel(0), aodinput->GetTrackLabel(0));
991 if(check == 0) continue;
992 if(check == -1) return;
1a31a9ab 993
994 //find the leading particles with highest momentum
b5dbb99b 995 if ((aodinput->Pt())>ptLeading)
996 {
1a31a9ab 997 ptLeading = aodinput->Pt() ;
226b95ba 998 idLeading = iaod ;
1a31a9ab 999 }
b5dbb99b 1000
226b95ba 1001 aodinput->SetLeadingParticle(kFALSE);
b5dbb99b 1002
1a31a9ab 1003 }//finish searching for leading trigger particle
1004
1005 // Check isolation of leading particle
1006 if(idLeading < 0) return;
226b95ba 1007
1a31a9ab 1008 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(idLeading));
226b95ba 1009 aodinput->SetLeadingParticle(kTRUE);
803d06a8 1010
1a31a9ab 1011 //After cuts, study isolation
1012 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
ac5111f9 1013 GetIsolationCut()->MakeIsolationCut(GetCTSTracks(),pl,
1014 GetReader(), GetCaloPID(),
b5dbb99b 1015 kTRUE, aodinput, GetAODObjArrayName(),
1016 n,nfrac,coneptsum, isolated);
1a31a9ab 1017 aodinput->SetIsolated(isolated);
1a31a9ab 1018
b5dbb99b 1019 if(GetDebug() > 1)
1020 {
1a31a9ab 1021 if(isolated)printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() : Particle %d IS ISOLATED \n",idLeading);
1022 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - End fill AODs \n");
1023 }
1024
1025}
1026
803d06a8 1027//_________________________________________________________
1a31a9ab 1028void AliAnaParticleIsolation::MakeAnalysisFillHistograms()
1029{
1030 //Do analysis and fill histograms
803d06a8 1031 Int_t n = 0, nfrac = 0;
1032 Bool_t isolated = kFALSE ;
1a31a9ab 1033 Float_t coneptsum = 0 ;
803d06a8 1034
1a31a9ab 1035 //Loop on stored AOD
1036 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1037 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Histo aod branch entries %d\n", naod);
1038
1039 //Get vertex for photon momentum calculation
1040 Double_t vertex[]={0,0,0} ; //vertex ;
b5dbb99b 1041 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1042 {
1a31a9ab 1043 GetReader()->GetVertex(vertex);
1a31a9ab 1044 }
1045
b5dbb99b 1046 for(Int_t iaod = 0; iaod < naod ; iaod++)
1047 {
1a31a9ab 1048 AliAODPWG4ParticleCorrelation* aod = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1049
1050 if(!aod->IsLeadingParticle()) continue; // Try to isolate only leading cluster or track
1051
1052 Bool_t isolation = aod->IsIsolated();
803d06a8 1053 Bool_t decay = aod->IsTagged();
0fb69ade 1054 Float_t energy = aod->E();
1a31a9ab 1055 Float_t pt = aod->Pt();
1056 Float_t phi = aod->Phi();
1057 Float_t eta = aod->Eta();
1058 Float_t conesize = GetIsolationCut()->GetConeSize();
1059
1060 //Recover reference arrays with clusters and tracks
1061 TObjArray * refclusters = aod->GetObjArray(GetAODObjArrayName()+"Clusters");
1062 TObjArray * reftracks = aod->GetObjArray(GetAODObjArrayName()+"Tracks");
0fb69ade 1063
1a31a9ab 1064 //If too small or too large pt, skip
1065 if(pt < GetMinPt() || pt > GetMaxPt() ) continue ;
1066
1067 // --- In case of redoing isolation from delta AOD ----
1068 if(fMakeSeveralIC) {
1069 //Analysis of multiple IC at same time
1070 MakeSeveralICAnalysis(aod);
1071 }
b5dbb99b 1072 else if(fReMakeIC)
1073 {
1a31a9ab 1074 //In case a more strict IC is needed in the produced AOD
1075 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
ac5111f9 1076 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
1077 GetReader(), GetCaloPID(),
b5dbb99b 1078 kFALSE, aod, "",
1079 n,nfrac,coneptsum, isolated);
1a31a9ab 1080 fhConeSumPt->Fill(pt,coneptsum);
1081 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Energy Sum in Isolation Cone %2.2f\n", coneptsum);
1082 }
1083 // --- -------------------------------------------- ----
1084
1085 //Fill pt distribution of particles in cone
1086 //Tracks
1087 coneptsum = 0;
1088 Double_t sumptFR = 0. ;
1089 TObjArray * trackList = GetCTSTracks() ;
b5dbb99b 1090 for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++)
1091 {
1a31a9ab 1092 AliVTrack* track = (AliVTrack *) trackList->At(itrack);
1093 //fill the histograms at forward range
1094 if(!track){
1095 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Track not available?");
1096 continue;
1097 }
0fb69ade 1098
1a31a9ab 1099 Double_t dPhi = phi - track->Phi() + TMath::PiOver2();
1100 Double_t dEta = eta - track->Eta();
1101 Double_t arg = dPhi*dPhi + dEta*dEta;
b5dbb99b 1102 if(TMath::Sqrt(arg) < conesize)
1103 {
1a31a9ab 1104 fhPtInFRCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1105 sumptFR+=track->Pt();
0fb69ade 1106 }
1107
1a31a9ab 1108 dPhi = phi - track->Phi() - TMath::PiOver2();
1109 arg = dPhi*dPhi + dEta*dEta;
b5dbb99b 1110 if(TMath::Sqrt(arg) < conesize)
1111 {
1a31a9ab 1112 fhPtInFRCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1113 sumptFR+=track->Pt();
1114 }
1115 }
0fb69ade 1116
1a31a9ab 1117 fhFRConeSumPt->Fill(pt,sumptFR);
b5dbb99b 1118 if(reftracks)
1119 {
1120 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++)
1121 {
1a31a9ab 1122 AliVTrack* track = (AliVTrack *) reftracks->At(itrack);
1123 fhPtInCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1124 coneptsum+=track->Pt();
1125 }
1126 }
1127
1128 //CaloClusters
b5dbb99b 1129 if(refclusters)
1130 {
1a31a9ab 1131 TLorentzVector mom ;
b5dbb99b 1132 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++)
1133 {
1a31a9ab 1134 AliVCluster* calo = (AliVCluster *) refclusters->At(icalo);
1135 calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
1136
1137 fhPtInCone->Fill(pt, mom.Pt());
1138 coneptsum+=mom.Pt();
1139 }
1140 }
1141
1142 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d Energy Sum in Isolation Cone %2.2f\n", iaod, coneptsum);
1143
1144 if(!fReMakeIC) fhConeSumPt->Fill(pt,coneptsum);
1145
b5dbb99b 1146 Int_t mcTag = aod->GetTag() ;
1147 Int_t clID = aod->GetCaloLabel(0) ;
1148
1149 if(isolation)
1150 {
1a31a9ab 1151 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d ISOLATED, fill histograms\n", iaod);
0fb69ade 1152
5c46c992 1153 FillTrackMatchingShowerShapeControlHistograms(clID,aod->GetFiducialArea(),mcTag);
b5dbb99b 1154
1155 fhEIso ->Fill(energy);
1156 fhPtIso ->Fill(pt);
1157 fhPhiIso ->Fill(pt,phi);
1158 fhEtaIso ->Fill(pt,eta);
0fb69ade 1159 fhEtaPhiIso ->Fill(eta,phi);
1160
803d06a8 1161 if (decay) fhPtDecayIso->Fill(pt);
1a31a9ab 1162
b5dbb99b 1163 if(IsDataMC())
1164 {
1a31a9ab 1165
b5dbb99b 1166 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
803d06a8 1167 {
1168 fhPtIsoMCPhoton ->Fill(pt);
1169 }
1170
b5dbb99b 1171 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt)){
1a31a9ab 1172 fhPtIsoPrompt ->Fill(pt);
1173 fhPhiIsoPrompt ->Fill(pt,phi);
1174 fhEtaIsoPrompt ->Fill(pt,eta);
1175 }
b5dbb99b 1176 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation))
1a31a9ab 1177 {
1178 fhPtIsoFragmentation ->Fill(pt);
1179 fhPhiIsoFragmentation ->Fill(pt,phi);
1180 fhEtaIsoFragmentation ->Fill(pt,eta);
1181 }
b5dbb99b 1182 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay))
1a31a9ab 1183 {
1184 fhPtIsoPi0Decay ->Fill(pt);
1185 fhPhiIsoPi0Decay ->Fill(pt,phi);
1186 fhEtaIsoPi0Decay ->Fill(pt,eta);
1187 }
b5dbb99b 1188 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay))
803d06a8 1189 {
1190 fhPtIsoEtaDecay ->Fill(pt);
1191 fhPhiIsoEtaDecay ->Fill(pt,phi);
1192 fhEtaIsoEtaDecay ->Fill(pt,eta);
1193 }
b5dbb99b 1194 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay))
1a31a9ab 1195 {
1196 fhPtIsoOtherDecay ->Fill(pt);
1197 fhPhiIsoOtherDecay ->Fill(pt,phi);
1198 fhEtaIsoOtherDecay ->Fill(pt,eta);
1199 }
b5dbb99b 1200 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion))
1a31a9ab 1201 {
1202 fhPtIsoConversion ->Fill(pt);
1203 fhPhiIsoConversion ->Fill(pt,phi);
1204 fhEtaIsoConversion ->Fill(pt,eta);
1205 }
803d06a8 1206 else // anything else
1a31a9ab 1207 {
1208 fhPtIsoUnknown ->Fill(pt);
1209 fhPhiIsoUnknown ->Fill(pt,phi);
1210 fhEtaIsoUnknown ->Fill(pt,eta);
1211 }
1212 }//Histograms with MC
1213
1214 }//Isolated histograms
1215
1216 if(!isolation)
1217 {
1218 fhPtNoIso ->Fill(pt);
803d06a8 1219 if (decay) fhPtDecayNoIso->Fill(pt);
1a31a9ab 1220
b5dbb99b 1221 if(IsDataMC())
1222 {
803d06a8 1223
b5dbb99b 1224 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
803d06a8 1225 {
1226 fhPtNoIsoMCPhoton->Fill(pt);
1227 }
1228
b5dbb99b 1229 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay))
1a31a9ab 1230 {
1231 fhPtNoIsoPi0Decay->Fill(pt);
1232 }
b5dbb99b 1233 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay))
1a31a9ab 1234 {
803d06a8 1235 fhPtNoIsoEtaDecay->Fill(pt);
1a31a9ab 1236 }
b5dbb99b 1237 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay))
1a31a9ab 1238 {
803d06a8 1239 fhPtNoIsoOtherDecay->Fill(pt);
1240 }
b5dbb99b 1241 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt))
803d06a8 1242 {
1243 fhPtNoIsoPrompt->Fill(pt);
1a31a9ab 1244 }
b5dbb99b 1245 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation))
0fb69ade 1246 {
1247 fhPtNoIsoFragmentation->Fill(pt);
1248 }
b5dbb99b 1249 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion))
0fb69ade 1250 {
1251 fhPtNoIsoConversion->Fill(pt);
1252 }
1253 else
1254 {
1255 fhPtNoIsoUnknown->Fill(pt);
1256 }
1a31a9ab 1257
1258 }
1259 }
1260
1261 }// aod loop
1262
1263}
1264
1a31a9ab 1265
803d06a8 1266//_____________________________________________________________________________________
1a31a9ab 1267void AliAnaParticleIsolation::MakeSeveralICAnalysis(AliAODPWG4ParticleCorrelation* ph)
1268{
1269 //Isolation Cut Analysis for both methods and different pt cuts and cones
1270 Float_t ptC = ph->Pt();
1271 Int_t tag = ph->GetTag();
1272
1273 //Keep original setting used when filling AODs, reset at end of analysis
1274 Float_t ptthresorg = GetIsolationCut()->GetPtThreshold();
1275 Float_t ptfracorg = GetIsolationCut()->GetPtFraction();
1276 Float_t rorg = GetIsolationCut()->GetConeSize();
1277
1278 Float_t coneptsum = 0 ;
1279 Int_t n[10][10];//[fNCones][fNPtThresFrac];
1280 Int_t nfrac[10][10];//[fNCones][fNPtThresFrac];
1281 Bool_t isolated = kFALSE;
1282
1283 //Loop on cone sizes
b5dbb99b 1284 for(Int_t icone = 0; icone<fNCones; icone++)
1285 {
1a31a9ab 1286 GetIsolationCut()->SetConeSize(fConeSizes[icone]);
1287 coneptsum = 0 ;
1288
1289 //Loop on ptthresholds
b5dbb99b 1290 for(Int_t ipt = 0; ipt<fNPtThresFrac ;ipt++)
1291 {
1a31a9ab 1292 n[icone][ipt]=0;
1293 nfrac[icone][ipt]=0;
1294 GetIsolationCut()->SetPtThreshold(fPtThresholds[ipt]);
1295 GetIsolationCut()->MakeIsolationCut(ph->GetObjArray(GetAODObjArrayName()+"Tracks"),
1296 ph->GetObjArray(GetAODObjArrayName()+"Clusters"),
ac5111f9 1297 GetReader(), GetCaloPID(),
b5dbb99b 1298 kFALSE, ph, "",
1299 n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);
1a31a9ab 1300
1301 //Normal ptThreshold cut
b5dbb99b 1302 if(n[icone][ipt] == 0)
1303 {
1a31a9ab 1304 fhPtThresIsolated[icone][ipt]->Fill(ptC);
b5dbb99b 1305 if(IsDataMC())
1306 {
803d06a8 1307 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtThresIsolatedPrompt[icone][ipt] ->Fill(ptC) ;
1308 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtThresIsolatedConversion[icone][ipt] ->Fill(ptC) ;
1309 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtThresIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
1310 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtThresIsolatedPi0Decay[icone][ipt] ->Fill(ptC) ;
1311 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtThresIsolatedEtaDecay[icone][ipt] ->Fill(ptC) ;
1312 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtThresIsolatedOtherDecay[icone][ipt] ->Fill(ptC) ;
1a31a9ab 1313 else fhPtThresIsolatedUnknown[icone][ipt]->Fill(ptC) ;
1314 }
1315 }
1316
1317 //Pt threshold on pt cand/ pt in cone fraction
b5dbb99b 1318 if(nfrac[icone][ipt] == 0)
1319 {
1a31a9ab 1320 fhPtFracIsolated[icone][ipt]->Fill(ptC);
b5dbb99b 1321 if(IsDataMC())
1322 {
803d06a8 1323 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtFracIsolatedPrompt[icone][ipt] ->Fill(ptC) ;
1324 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtFracIsolatedConversion[icone][ipt] ->Fill(ptC) ;
1325 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtFracIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
1326 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtFracIsolatedPi0Decay[icone][ipt] ->Fill(ptC) ;
1327 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtFracIsolatedEtaDecay[icone][ipt] ->Fill(ptC) ;
1328 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtFracIsolatedOtherDecay[icone][ipt] ->Fill(ptC) ;
1a31a9ab 1329 else fhPtFracIsolatedUnknown[icone][ipt]->Fill(ptC) ;
1330 }
1331 }
1332 }//pt thresh loop
1333
1334 //Sum in cone histograms
1335 fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ;
b5dbb99b 1336 if(IsDataMC())
1337 {
803d06a8 1338 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtSumIsolatedPrompt[icone] ->Fill(ptC,coneptsum) ;
1339 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtSumIsolatedConversion[icone] ->Fill(ptC,coneptsum) ;
1340 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtSumIsolatedFragmentation[icone]->Fill(ptC,coneptsum) ;
1341 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtSumIsolatedPi0Decay[icone] ->Fill(ptC,coneptsum) ;
1342 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtSumIsolatedEtaDecay[icone] ->Fill(ptC,coneptsum) ;
1343 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtSumIsolatedOtherDecay[icone] ->Fill(ptC,coneptsum) ;
1a31a9ab 1344 else fhPtSumIsolatedUnknown[icone]->Fill(ptC,coneptsum) ;
1345 }
1346
1347 }//cone size loop
1348
1349 //Reset original parameters for AOD analysis
1350 GetIsolationCut()->SetPtThreshold(ptthresorg);
1351 GetIsolationCut()->SetPtFraction(ptfracorg);
1352 GetIsolationCut()->SetConeSize(rorg);
1353
1354}
1355
803d06a8 1356//_____________________________________________________________
1a31a9ab 1357void AliAnaParticleIsolation::Print(const Option_t * opt) const
1358{
1359
1360 //Print some relevant parameters set for the analysis
1361 if(! opt)
1362 return;
1363
1364 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
745913ae 1365 AliAnaCaloTrackCorrBaseClass::Print(" ");
1a31a9ab 1366
1367 printf("ReMake Isolation = %d \n", fReMakeIC) ;
1368 printf("Make Several Isolation = %d \n", fMakeSeveralIC) ;
1369 printf("Calorimeter for isolation = %s \n", fCalorimeter.Data()) ;
1370
b5dbb99b 1371 if(fMakeSeveralIC)
1372 {
1a31a9ab 1373 printf("N Cone Sizes = %d\n", fNCones) ;
1374 printf("Cone Sizes = \n") ;
1375 for(Int_t i = 0; i < fNCones; i++)
1376 printf(" %1.2f;", fConeSizes[i]) ;
1377 printf(" \n") ;
1378
1379 printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ;
1380 printf(" pT thresholds = \n") ;
1381 for(Int_t i = 0; i < fNPtThresFrac; i++)
1382 printf(" %2.2f;", fPtThresholds[i]) ;
1383
1384 printf(" \n") ;
1385
1386 printf(" pT fractions = \n") ;
1387 for(Int_t i = 0; i < fNPtThresFrac; i++)
1388 printf(" %2.2f;", fPtFractions[i]) ;
1389
1390 }
1391
b5dbb99b 1392 printf("Histograms: %3.1f < pT sum < %3.1f, Nbin = %d\n", fHistoPtSumMin, fHistoPtSumMax, fHistoNPtSumBins );
1a31a9ab 1393 printf("Histograms: %3.1f < pT in cone < %3.1f, Nbin = %d\n", fHistoPtInConeMin, fHistoPtInConeMax, fHistoNPtInConeBins);
1394
1395 printf(" \n") ;
1396
1397}
1398