]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaParticleIsolation.cxx
add pile-up related histograms plus other minor fixes
[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//////////////////////////////////////////////////////////////////////////////
db6fb352 27
28
1a31a9ab 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"
2ad19c3d 45#include "AliESDEvent.h"
46#include "AliAODEvent.h"
1a31a9ab 47
48ClassImp(AliAnaParticleIsolation)
db6fb352 49
803d06a8 50//______________________________________________________________________________
db6fb352 51AliAnaParticleIsolation::AliAnaParticleIsolation() :
52AliAnaCaloTrackCorrBaseClass(), fCalorimeter(""),
53fReMakeIC(0), fMakeSeveralIC(0),
2ad19c3d 54fFillPileUpHistograms(0),
db6fb352 55fFillTMHisto(0), fFillSSHisto(0),
56// Several IC
57fNCones(0), fNPtThresFrac(0),
58fConeSizes(), fPtThresholds(),
59fPtFractions(), fSumPtThresholds(),
60// Histograms
8736d400 61fhEIso(0), fhPtIso(0), fhPtNLocMaxIso(0),
db6fb352 62fhPhiIso(0), fhEtaIso(0), fhEtaPhiIso(0),
63fhEtaPhiNoIso(0),
b1f720a7 64fhENoIso(0), fhPtNoIso(0), fhPtNLocMaxNoIso(0),
8736d400 65fhPtDecayIso(0), fhPtDecayNoIso(0),
db6fb352 66fhEtaPhiDecayIso(0), fhEtaPhiDecayNoIso(0),
6c80c1bf 67fhConeSumPt(0), fhPtInCone(0),
68fhPtInConePileUp(0), fhPtInConeCent(0),
b7ce43b4 69fhFRConeSumPt(0), fhPtInFRCone(0), fhPhiUEConeSumPt(0),
70fhEtaUEConeSumPt(0), fhEtaBand(0), fhPhiBand(0),
71fhConeSumPtEtaUESub(0), fhConeSumPtPhiUESub(0),
db6fb352 72// MC histograms
73fhPtIsoPrompt(0), fhPhiIsoPrompt(0), fhEtaIsoPrompt(0),
74fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(),
75fhPtIsoFragmentation(0), fhPhiIsoFragmentation(0), fhEtaIsoFragmentation(0),
76fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(),
764ab1f4 77fhPtIsoPi0(0), fhPhiIsoPi0(0), fhEtaIsoPi0(0),
78fhPtThresIsolatedPi0(), fhPtFracIsolatedPi0(), fhPtSumIsolatedPi0(),
db6fb352 79fhPtIsoPi0Decay(0), fhPhiIsoPi0Decay(0), fhEtaIsoPi0Decay(0),
80fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(),
81fhPtIsoEtaDecay(0), fhPhiIsoEtaDecay(0), fhEtaIsoEtaDecay(0),
82fhPtThresIsolatedEtaDecay(), fhPtFracIsolatedEtaDecay(), fhPtSumIsolatedEtaDecay(),
83fhPtIsoOtherDecay(0), fhPhiIsoOtherDecay(0), fhEtaIsoOtherDecay(0),
84fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(),
764ab1f4 85//fhPtIsoConversion(0), fhPhiIsoConversion(0), fhEtaIsoConversion(0),
86//fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(),
87fhPtIsoHadron(0), fhPhiIsoHadron(0), fhEtaIsoHadron(0),
88fhPtThresIsolatedHadron(), fhPtFracIsolatedHadron(), fhPtSumIsolatedHadron(),
89fhPtNoIsoPi0(0), fhPtNoIsoPi0Decay(0),
90fhPtNoIsoEtaDecay(0), fhPtNoIsoOtherDecay(0),
db6fb352 91fhPtNoIsoPrompt(0), fhPtIsoMCPhoton(0), fhPtNoIsoMCPhoton(0),
764ab1f4 92//fhPtNoIsoConversion(0),
93fhPtNoIsoFragmentation(0), fhPtNoIsoHadron(0),
db6fb352 94// Hist several IC
764ab1f4 95fhSumPtLeadingPt(), fhPtLeadingPt(),
96fhFRSumPtLeadingPt(), fhFRPtLeadingPt(),
db6fb352 97fhPtThresIsolated(), fhPtFracIsolated(), fhPtSumIsolated(),
98fhEtaPhiPtThresIso(), fhEtaPhiPtThresDecayIso(), fhPtPtThresDecayIso(),
99fhEtaPhiPtFracIso(), fhEtaPhiPtFracDecayIso(), fhPtPtFracDecayIso(),
b0a31c92 100fhPtPtSumDecayIso(), fhEtaPhiSumDensityIso(), fhEtaPhiSumDensityDecayIso(),
101fhPtSumDensityIso(), fhPtSumDensityDecayIso(),
102fhPtFracPtSumIso(), fhPtFracPtSumDecayIso(),
103fhEtaPhiFracPtSumIso(), fhEtaPhiFracPtSumDecayIso(),
db6fb352 104// Cluster control histograms
ca134929 105fhTrackMatchedDEta(), fhTrackMatchedDPhi(), fhTrackMatchedDEtaDPhi(),
106fhdEdx(), fhEOverP(), fhTrackMatchedMCParticle(),
db7b861a 107fhELambda0() , fhELambda1(), fhELambda0SSBkg(),
ca134929 108fhELambda0TRD(), fhELambda1TRD(),
764ab1f4 109fhELambda0MCPhoton(), fhELambda0MCPi0(), fhELambda0MCPi0Decay(),
110fhELambda0MCEtaDecay(), fhELambda0MCOtherDecay(), fhELambda0MCHadron(),
db6fb352 111// Number of local maxima in cluster
ca134929 112fhNLocMax(),
113fhELambda0LocMax1(), fhELambda1LocMax1(),
114fhELambda0LocMax2(), fhELambda1LocMax2(),
115fhELambda0LocMaxN(), fhELambda1LocMaxN(),
2ad19c3d 116// PileUp
b1f720a7 117fhEIsoPileUp(0), fhPtIsoPileUp(0),
118fhENoIsoPileUp(0), fhPtNoIsoPileUp(0),
2ad19c3d 119fhTimeENoCut(0), fhTimeESPD(0), fhTimeESPDMulti(0),
120fhTimeNPileUpVertSPD(0), fhTimeNPileUpVertTrack(0),
121fhTimeNPileUpVertContributors(0),
122fhTimePileUpMainVertexZDistance(0), fhTimePileUpMainVertexZDiamond(0),
db6fb352 123// Histograms settings
124fHistoNPtSumBins(0), fHistoPtSumMax(0.), fHistoPtSumMin(0.),
125fHistoNPtInConeBins(0), fHistoPtInConeMax(0.), fHistoPtInConeMin(0.)
1a31a9ab 126{
127 //default ctor
128
129 //Initialize parameters
130 InitParameters();
db6fb352 131
b5dbb99b 132 for(Int_t i = 0; i < 5 ; i++)
133 {
803d06a8 134 fConeSizes[i] = 0 ;
1a31a9ab 135
db6fb352 136 fhPtSumIsolatedPrompt [i] = 0 ;
1a31a9ab 137 fhPtSumIsolatedFragmentation[i] = 0 ;
db6fb352 138 fhPtSumIsolatedPi0Decay [i] = 0 ;
764ab1f4 139 fhPtSumIsolatedPi0 [i] = 0 ;
db6fb352 140 fhPtSumIsolatedEtaDecay [i] = 0 ;
141 fhPtSumIsolatedOtherDecay [i] = 0 ;
764ab1f4 142// fhPtSumIsolatedConversion [i] = 0 ;
143 fhPtSumIsolatedHadron [i] = 0 ;
1a31a9ab 144
b5dbb99b 145 for(Int_t j = 0; j < 5 ; j++)
146 {
db6fb352 147 fhPtThresIsolated [i][j] = 0 ;
148 fhPtFracIsolated [i][j] = 0 ;
149 fhPtSumIsolated [i][j] = 0 ;
150
151 fhEtaPhiPtThresIso [i][j] = 0 ;
152 fhEtaPhiPtThresDecayIso[i][j] = 0 ;
153 fhPtPtThresDecayIso [i][j] = 0 ;
154
155 fhEtaPhiPtFracIso [i][j] = 0 ;
156 fhEtaPhiPtFracDecayIso [i][j] = 0 ;
157 fhPtPtFracDecayIso [i][j] = 0 ;
158 fhPtPtSumDecayIso [i][j] = 0 ;
159 fhPtSumDensityIso [i][j] = 0 ;
160 fhPtSumDensityDecayIso [i][j] = 0 ;
b0a31c92 161 fhEtaPhiSumDensityIso [i][j] = 0 ;
162 fhEtaPhiSumDensityDecayIso [i][j] = 0 ;
ca134929 163 fhPtFracPtSumIso [i][j] = 0 ;
164 fhPtFracPtSumDecayIso [i][j] = 0 ;
b0a31c92 165 fhEtaPhiFracPtSumIso [i][j] = 0 ;
166 fhEtaPhiFracPtSumDecayIso [i][j] = 0 ;
db6fb352 167
168 fhPtThresIsolatedPrompt [i][j] = 0 ;
169 fhPtThresIsolatedFragmentation[i][j] = 0 ;
170 fhPtThresIsolatedPi0Decay [i][j] = 0 ;
764ab1f4 171 fhPtThresIsolatedPi0 [i][j] = 0 ;
172 fhPtThresIsolatedEtaDecay [i][j] = 0 ;
db6fb352 173 fhPtThresIsolatedOtherDecay [i][j] = 0 ;
764ab1f4 174// fhPtThresIsolatedConversion [i][j] = 0 ;
175 fhPtThresIsolatedHadron [i][j] = 0 ;
db6fb352 176
177 fhPtFracIsolatedPrompt [i][j] = 0 ;
178 fhPtFracIsolatedFragmentation [i][j] = 0 ;
764ab1f4 179 fhPtFracIsolatedPi0 [i][j] = 0 ;
db6fb352 180 fhPtFracIsolatedPi0Decay [i][j] = 0 ;
181 fhPtFracIsolatedEtaDecay [i][j] = 0 ;
182 fhPtFracIsolatedOtherDecay [i][j] = 0 ;
764ab1f4 183// fhPtFracIsolatedConversion [i][j] = 0 ;
184 fhPtFracIsolatedHadron [i][j] = 0 ;
db6fb352 185
1a31a9ab 186 }
187 }
188
db6fb352 189 for(Int_t i = 0; i < 5 ; i++)
190 {
191 fPtFractions [i] = 0 ;
192 fPtThresholds [i] = 0 ;
193 fSumPtThresholds[i] = 0 ;
1a31a9ab 194 }
ca134929 195
196
197 for(Int_t i = 0; i < 2 ; i++)
198 {
199 fhTrackMatchedDEta[i] = 0 ; fhTrackMatchedDPhi[i] = 0 ; fhTrackMatchedDEtaDPhi [i] = 0 ;
200 fhdEdx [i] = 0 ; fhEOverP [i] = 0 ; fhTrackMatchedMCParticle[i] = 0 ;
201 fhELambda0 [i] = 0 ; fhELambda1 [i] = 0 ;
202 fhELambda0TRD [i] = 0 ; fhELambda1TRD [i] = 0 ;
203
764ab1f4 204 fhELambda0MCPhoton [i] = 0 ; fhELambda0MCPi0 [i] = 0 ; fhELambda0MCPi0Decay[i] = 0 ;
205 fhELambda0MCEtaDecay[i] = 0 ; fhELambda0MCOtherDecay[i] = 0 ; fhELambda0MCHadron [i] = 0 ;
206
207
ca134929 208 // Number of local maxima in cluster
209 fhNLocMax [i] = 0 ;
210 fhELambda0LocMax1[i] = 0 ; fhELambda1LocMax1[i] = 0 ;
211 fhELambda0LocMax2[i] = 0 ; fhELambda1LocMax2[i] = 0 ;
212 fhELambda0LocMaxN[i] = 0 ; fhELambda1LocMaxN[i] = 0 ;
213
214 }
db6fb352 215
1a31a9ab 216}
217
2ad19c3d 218//_________________________________________________________________
219void AliAnaParticleIsolation::FillPileUpHistograms(Int_t clusterID)
220{
221 // Fill some histograms to understand pile-up
222 if(!fFillPileUpHistograms) return;
223
224 if(clusterID < 0 )
225 {
226 printf("AliAnaParticleIsolation::FillPileUpHistograms(), ID of cluster = %d, not possible! ", clusterID);
227 return;
228 }
229
230 Int_t iclus = -1;
231 TObjArray* clusters = 0x0;
232 if (fCalorimeter == "EMCAL") clusters = GetEMCALClusters();
233 else if(fCalorimeter == "PHOS" ) clusters = GetPHOSClusters();
234
235 Float_t energy = 0;
236 Float_t time = -1000;
237
238 if(clusters)
239 {
240 AliVCluster *cluster = FindCluster(clusters,clusterID,iclus);
241 energy = cluster->E();
994051fa 242 time = cluster->GetTOF()*1e9;
2ad19c3d 243 }
244
245 //printf("E %f, time %f\n",energy,time);
246 AliVEvent * event = GetReader()->GetInputEvent();
247
248 fhTimeENoCut->Fill(energy,time);
249 if(GetReader()->IsPileUpFromSPD()) fhTimeESPD ->Fill(energy,time);
250 if(event->IsPileupFromSPDInMultBins()) fhTimeESPDMulti->Fill(energy,time);
251
de101942 252 if(energy < 8) return; // Fill time figures for high energy clusters not too close to trigger threshold
2ad19c3d 253
254 AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
255 AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
256
257 // N pile up vertices
258 Int_t nVerticesSPD = -1;
259 Int_t nVerticesTracks = -1;
260
261 if (esdEv)
262 {
263 nVerticesSPD = esdEv->GetNumberOfPileupVerticesSPD();
264 nVerticesTracks = esdEv->GetNumberOfPileupVerticesTracks();
265
266 }//ESD
267 else if (aodEv)
268 {
269 nVerticesSPD = aodEv->GetNumberOfPileupVerticesSPD();
270 nVerticesTracks = aodEv->GetNumberOfPileupVerticesTracks();
271 }//AOD
272
273 fhTimeNPileUpVertSPD ->Fill(time,nVerticesSPD);
274 fhTimeNPileUpVertTrack->Fill(time,nVerticesTracks);
275
276 //printf("Is SPD %d, Is SPD Multi %d, n spd %d, n track %d\n",
277 // GetReader()->IsPileUpFromSPD(),event->IsPileupFromSPDInMultBins(),nVerticesSPD,nVerticesTracks);
278
279 Int_t ncont = -1;
5559f30a 280 Float_t z1 = -1, z2 = -1;
2ad19c3d 281 Float_t diamZ = -1;
282 for(Int_t iVert=0; iVert<nVerticesSPD;iVert++)
283 {
284 if (esdEv)
285 {
286 const AliESDVertex* pv=esdEv->GetPileupVertexSPD(iVert);
287 ncont=pv->GetNContributors();
288 z1 = esdEv->GetPrimaryVertexSPD()->GetZ();
289 z2 = pv->GetZ();
290 diamZ = esdEv->GetDiamondZ();
291 }//ESD
292 else if (aodEv)
293 {
294 AliAODVertex *pv=aodEv->GetVertex(iVert);
295 if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
296 ncont=pv->GetNContributors();
297 z1=aodEv->GetPrimaryVertexSPD()->GetZ();
298 z2=pv->GetZ();
299 diamZ = aodEv->GetDiamondZ();
300 }// AOD
301
302 Double_t distZ = TMath::Abs(z2-z1);
303 diamZ = TMath::Abs(z2-diamZ);
304
305 fhTimeNPileUpVertContributors ->Fill(time,ncont);
306 fhTimePileUpMainVertexZDistance->Fill(time,distZ);
307 fhTimePileUpMainVertexZDiamond ->Fill(time,diamZ);
308
309 }// loop
310}
311
b5dbb99b 312//________________________________________________________________________________________________
ca134929 313void AliAnaParticleIsolation::FillTrackMatchingShowerShapeControlHistograms(const Bool_t isolated,
314 const Int_t clusterID,
315 const Int_t nMaxima,
db7b861a 316 const Int_t mcTag,
764ab1f4 317 const TObjArray * plCTS,
318 const TObjArray * plNe,
319 AliAODPWG4ParticleCorrelation *pCandidate,
320 const AliCaloTrackReader * reader,
db7b861a 321 const AliCaloPID * pid)
b5dbb99b 322{
db7b861a 323 // Fill Track matching and Shower Shape control histograms
b5dbb99b 324 if(!fFillTMHisto && !fFillSSHisto) return;
325
547c2f01 326 if(clusterID < 0 )
327 {
81a5e27b 328 printf("AliAnaParticleIsolation::FillTrackMatchingShowerShapeControlHistograms(), ID of cluster = %d, not possible! \n", clusterID);
547c2f01 329 return;
330 }
331
b5dbb99b 332 Int_t iclus = -1;
333 TObjArray* clusters = 0x0;
334 if (fCalorimeter == "EMCAL") clusters = GetEMCALClusters();
335 else if(fCalorimeter == "PHOS" ) clusters = GetPHOSClusters();
336
337 if(clusters)
338 {
339
340 AliVCluster *cluster = FindCluster(clusters,clusterID,iclus);
341 Float_t energy = cluster->E();
342
343 if(fFillSSHisto)
344 {
ca134929 345 fhELambda0[isolated]->Fill(energy, cluster->GetM02() );
346 fhELambda1[isolated]->Fill(energy, cluster->GetM20() );
b5dbb99b 347
764ab1f4 348 if(IsDataMC())
349 {
350 if (GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt) ||
351 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation)) fhELambda0MCPhoton [isolated]->Fill(energy, cluster->GetM02());
352 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0)) fhELambda0MCPi0 [isolated]->Fill(energy, cluster->GetM02());
353 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay)) fhELambda0MCPi0Decay [isolated]->Fill(energy, cluster->GetM02());
354 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay)) fhELambda0MCEtaDecay [isolated]->Fill(energy, cluster->GetM02());
355 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay)) fhELambda0MCOtherDecay[isolated]->Fill(energy, cluster->GetM02());
356
2dac6a45 357 // else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion)) fhPtNoIsoConversion ->Fill(energy, cluster->GetM02());
358 else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron)) fhELambda0MCHadron [isolated]->Fill(energy, cluster->GetM02());
764ab1f4 359
360 }
361
ca134929 362 if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5) // TO DO: CHANGE FOR 2012
b5dbb99b 363 {
ca134929 364 fhELambda0TRD[isolated]->Fill(energy, cluster->GetM02() );
365 fhELambda1TRD[isolated]->Fill(energy, cluster->GetM20() );
b5dbb99b 366 }
5c46c992 367
ca134929 368 fhNLocMax[isolated]->Fill(energy,nMaxima);
369 if (nMaxima==1) { fhELambda0LocMax1[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMax1[isolated]->Fill(energy,cluster->GetM20()); }
370 else if(nMaxima==2) { fhELambda0LocMax2[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMax2[isolated]->Fill(energy,cluster->GetM20()); }
371 else { fhELambda0LocMaxN[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMaxN[isolated]->Fill(energy,cluster->GetM20()); }
db7b861a 372
373 if(isolated==0)
374 {
375 //Analyse non-isolated events
376 Int_t n = 0;
377 Int_t nfrac = 0;
378 Bool_t iso = kFALSE ;
379 Float_t coneptsum = 0 ;
380 GetIsolationCut()->SetPtThresholdMax(1.);
381 GetIsolationCut()->MakeIsolationCut(plCTS, plNe,
382 reader, pid,
383 kFALSE, pCandidate, "",
384 n,nfrac,coneptsum, iso);
385 if (!iso) fhELambda0SSBkg->Fill(energy, cluster->GetM02());
386
5c46c992 387
db7b861a 388 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Energy Sum in Isolation Cone %2.2f\n", coneptsum);
389 }
390 GetIsolationCut()->SetPtThresholdMax(10000.);
391
b5dbb99b 392 } // SS histo fill
393
394
395 if(fFillTMHisto)
396 {
397 Float_t dZ = cluster->GetTrackDz();
398 Float_t dR = cluster->GetTrackDx();
399
400 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
401 {
402 dR = 2000., dZ = 2000.;
403 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
404 }
405
406 //printf("ParticleIsolation: dPhi %f, dEta %f\n",dR,dZ);
ec58c056 407 if(fhTrackMatchedDEta[isolated] && TMath::Abs(dR) < 999)
b5dbb99b 408 {
ca134929 409 fhTrackMatchedDEta[isolated]->Fill(energy,dZ);
410 fhTrackMatchedDPhi[isolated]->Fill(energy,dR);
411 if(energy > 0.5) fhTrackMatchedDEtaDPhi[isolated]->Fill(dZ,dR);
b5dbb99b 412 }
413
414 // Check dEdx and E/p of matched clusters
415
416 if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
417 {
db6fb352 418
4bfeae64 419 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
b5dbb99b 420
421 if(track)
422 {
423 Float_t dEdx = track->GetTPCsignal();
ca134929 424 fhdEdx[isolated]->Fill(cluster->E(), dEdx);
b5dbb99b 425
426 Float_t eOverp = cluster->E()/track->P();
ca134929 427 fhEOverP[isolated]->Fill(cluster->E(), eOverp);
b5dbb99b 428 }
4bfeae64 429 //else
430 // printf("AliAnaParticleIsolation::FillTrackMatchingShowerShapeHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
b5dbb99b 431
db6fb352 432
433 if(IsDataMC())
434 {
b5dbb99b 435 if ( !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) )
436 {
437 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
ca134929 438 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 2.5 );
439 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 0.5 );
440 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 1.5 );
441 else fhTrackMatchedMCParticle[isolated]->Fill(energy, 3.5 );
b5dbb99b 442
443 }
444 else
445 {
446 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
ca134929 447 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 6.5 );
448 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 4.5 );
449 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 5.5 );
450 else fhTrackMatchedMCParticle[isolated]->Fill(energy, 7.5 );
b5dbb99b 451 }
452
453 } // MC
454
455 } // match window
456
457 }// TM histos fill
458
459 } // clusters array available
460
461}
462
803d06a8 463//______________________________________________________
1a31a9ab 464TObjString * AliAnaParticleIsolation::GetAnalysisCuts()
465{
b0a31c92 466 //Save parameters used for analysis
1a31a9ab 467 TString parList ; //this will be list of parameters used for this analysis.
468 const Int_t buffersize = 255;
469 char onePar[buffersize] ;
470
471 snprintf(onePar, buffersize,"--- AliAnaParticleIsolation ---\n") ;
472 parList+=onePar ;
473 snprintf(onePar, buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
474 parList+=onePar ;
475 snprintf(onePar, buffersize,"fReMakeIC =%d (Flag for reisolation during histogram filling) \n",fReMakeIC) ;
476 parList+=onePar ;
477 snprintf(onePar, buffersize,"fMakeSeveralIC=%d (Flag for isolation with several cuts at the same time ) \n",fMakeSeveralIC) ;
09273901 478 parList+=onePar ;
479 snprintf(onePar, buffersize,"fFillTMHisto=%d (Flag for track matching histograms) \n",fFillTMHisto) ;
1a31a9ab 480 parList+=onePar ;
09273901 481 snprintf(onePar, buffersize,"fFillSSHisto=%d (Flag for shower shape histograms) \n",fFillSSHisto) ;
482 parList+=onePar ;
db6fb352 483
b5dbb99b 484 if(fMakeSeveralIC)
485 {
1a31a9ab 486 snprintf(onePar, buffersize,"fNCones =%d (Number of cone sizes) \n",fNCones) ;
487 parList+=onePar ;
488 snprintf(onePar, buffersize,"fNPtThresFrac=%d (Flag for isolation with several cuts at the same time ) \n",fNPtThresFrac) ;
489 parList+=onePar ;
490
b5dbb99b 491 for(Int_t icone = 0; icone < fNCones ; icone++)
492 {
1a31a9ab 493 snprintf(onePar, buffersize,"fConeSizes[%d]=%1.2f (isolation cone size) \n",icone, fConeSizes[icone]) ;
494 parList+=onePar ;
495 }
b5dbb99b 496 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
497 {
1a31a9ab 498 snprintf(onePar, buffersize,"fPtThresholds[%d]=%1.2f (isolation pt threshold) \n",ipt, fPtThresholds[ipt]) ;
499 parList+=onePar ;
500 }
b5dbb99b 501 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
502 {
1a31a9ab 503 snprintf(onePar, buffersize,"fPtFractions[%d]=%1.2f (isolation pt fraction threshold) \n",ipt, fPtFractions[ipt]) ;
504 parList+=onePar ;
db6fb352 505 }
506 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
507 {
508 snprintf(onePar, buffersize,"fSumPtThresholds[%d]=%1.2f (isolation sum pt threshold) \n",ipt, fSumPtThresholds[ipt]) ;
509 parList+=onePar ;
1a31a9ab 510 }
511 }
512
513 //Get parameters set in base class.
514 parList += GetBaseParametersList() ;
515
516 //Get parameters set in IC class.
517 if(!fMakeSeveralIC)parList += GetIsolationCut()->GetICParametersList() ;
518
519 return new TObjString(parList) ;
b0a31c92 520
1a31a9ab 521}
522
803d06a8 523//________________________________________________________
1a31a9ab 524TList * AliAnaParticleIsolation::GetCreateOutputObjects()
525{
526 // Create histograms to be saved in output file and
527 // store them in outputContainer
528 TList * outputContainer = new TList() ;
529 outputContainer->SetName("IsolatedParticleHistos") ;
530
745913ae 531 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins();
532 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins();
533 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
534 Float_t ptmax = GetHistogramRanges()->GetHistoPtMax();
535 Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();
536 Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
537 Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
538 Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
539 Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
09273901 540 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins();
541 Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax();
542 Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
2ad19c3d 543 Int_t ntimebins= GetHistogramRanges()->GetHistoTimeBins();
544 Float_t timemax = GetHistogramRanges()->GetHistoTimeMax();
545 Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
546
09273901 547 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
548 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
549 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
550 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
551 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
552 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
803d06a8 553
31ae6d59 554 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins();
555 Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax();
556 Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
557 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
558 Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax();
559 Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
560
803d06a8 561 Int_t nptsumbins = fHistoNPtSumBins;
562 Float_t ptsummax = fHistoPtSumMax;
563 Float_t ptsummin = fHistoPtSumMin;
564 Int_t nptinconebins = fHistoNPtInConeBins;
565 Float_t ptinconemax = fHistoPtInConeMax;
566 Float_t ptinconemin = fHistoPtInConeMin;
1a31a9ab 567
db6fb352 568 Float_t ptthre = GetIsolationCut()->GetPtThreshold();
569 Float_t ptfrac = GetIsolationCut()->GetPtFraction();
570 Float_t r = GetIsolationCut()->GetConeSize();
571
b5dbb99b 572 if(!fMakeSeveralIC)
573 {
ca134929 574 TString hName [] = {"NoIso",""};
575 TString hTitle[] = {"Not isolated" ,"isolated"};
764ab1f4 576 if(fFillSSHisto)
577 {
db7b861a 578 fhELambda0SSBkg = new TH2F
764ab1f4 579 ("hELambda0SSBkg","Non isolated clusters : E vs #lambda_{0}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
580 fhELambda0SSBkg->SetYTitle("#lambda_{0}^{2}");
581 fhELambda0SSBkg->SetXTitle("E (GeV)");
582 outputContainer->Add(fhELambda0SSBkg) ;
583 }
584
ca134929 585 for(Int_t iso = 0; iso < 2; iso++)
b5dbb99b 586 {
ca134929 587 if(fFillTMHisto)
31ae6d59 588 {
ca134929 589 fhTrackMatchedDEta[iso] = new TH2F
590 (Form("hTrackMatchedDEta%s",hName[iso].Data()),
591 Form("%s - d#eta of cluster-track vs cluster energy for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
592 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
593 fhTrackMatchedDEta[iso]->SetYTitle("d#eta");
594 fhTrackMatchedDEta[iso]->SetXTitle("E_{cluster} (GeV)");
595
596 fhTrackMatchedDPhi[iso] = new TH2F
597 (Form("hTrackMatchedDPhi%s",hName[iso].Data()),
598 Form("%s - d#phi of cluster-track vs cluster energy for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
599 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
600 fhTrackMatchedDPhi[iso]->SetYTitle("d#phi (rad)");
601 fhTrackMatchedDPhi[iso]->SetXTitle("E_{cluster} (GeV)");
602
603 fhTrackMatchedDEtaDPhi[iso] = new TH2F
604 (Form("hTrackMatchedDEtaDPhi%s",hName[iso].Data()),
605 Form("%s - d#eta vs d#phi of cluster-track for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
606 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
607 fhTrackMatchedDEtaDPhi[iso]->SetYTitle("d#phi (rad)");
608 fhTrackMatchedDEtaDPhi[iso]->SetXTitle("d#eta");
609
610 outputContainer->Add(fhTrackMatchedDEta[iso]) ;
611 outputContainer->Add(fhTrackMatchedDPhi[iso]) ;
612 outputContainer->Add(fhTrackMatchedDEtaDPhi[iso]) ;
31ae6d59 613
ca134929 614 fhdEdx[iso] = new TH2F
615 (Form("hdEdx%s",hName[iso].Data()),
616 Form("%s - Matched track <dE/dx> vs cluster E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
617 nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
618 fhdEdx[iso]->SetXTitle("E (GeV)");
619 fhdEdx[iso]->SetYTitle("<dE/dx>");
620 outputContainer->Add(fhdEdx[iso]);
31ae6d59 621
ca134929 622 fhEOverP[iso] = new TH2F
623 (Form("hEOverP%s",hName[iso].Data()),
624 Form("%s - Matched track E/p vs cluster E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
625 nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
626 fhEOverP[iso]->SetXTitle("E (GeV)");
627 fhEOverP[iso]->SetYTitle("E/p");
628 outputContainer->Add(fhEOverP[iso]);
629
630 if(IsDataMC())
631 {
632 fhTrackMatchedMCParticle[iso] = new TH2F
633 (Form("hTrackMatchedMCParticle%s",hName[iso].Data()),
634 Form("%s - Origin of particle vs energy vs cluster E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
635 nptbins,ptmin,ptmax,8,0,8);
636 fhTrackMatchedMCParticle[iso]->SetXTitle("E (GeV)");
637 //fhTrackMatchedMCParticle[iso]->SetYTitle("Particle type");
638
639 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(1 ,"Photon");
640 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(2 ,"Electron");
641 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
642 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(4 ,"Rest");
643 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
644 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
645 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
646 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
647
648 outputContainer->Add(fhTrackMatchedMCParticle[iso]);
649 }
31ae6d59 650 }
b5dbb99b 651
ca134929 652 if(fFillSSHisto)
b5dbb99b 653 {
ca134929 654 fhELambda0[iso] = new TH2F
655 (Form("hELambda0%s",hName[iso].Data()),
656 Form("%s cluster : E vs #lambda_{0}",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
657 fhELambda0[iso]->SetYTitle("#lambda_{0}^{2}");
658 fhELambda0[iso]->SetXTitle("E (GeV)");
659 outputContainer->Add(fhELambda0[iso]) ;
660
764ab1f4 661 if(IsDataMC())
662 {
663 fhELambda0MCPhoton[iso] = new TH2F
664 (Form("hELambda0%s_MCPhoton",hName[iso].Data()),
665 Form("%s cluster : E vs #lambda_{0}: Origin is final state photon",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
666 fhELambda0MCPhoton[iso]->SetYTitle("#lambda_{0}^{2}");
667 fhELambda0MCPhoton[iso]->SetXTitle("E (GeV)");
668 outputContainer->Add(fhELambda0MCPhoton[iso]) ;
669
670 fhELambda0MCPi0[iso] = new TH2F
671 (Form("hELambda0%s_MCPi0",hName[iso].Data()),
672 Form("%s cluster : E vs #lambda_{0}: Origin is pi0 (2 #gamma)",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
673 fhELambda0MCPi0[iso]->SetYTitle("#lambda_{0}^{2}");
674 fhELambda0MCPi0[iso]->SetXTitle("E (GeV)");
675 outputContainer->Add(fhELambda0MCPi0[iso]) ;
676
677 fhELambda0MCPi0Decay[iso] = new TH2F
678 (Form("hELambda0%s_MCPi0Decay",hName[iso].Data()),
679 Form("%s cluster : E vs #lambda_{0}: Origin is pi0 decay",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
680 fhELambda0MCPi0Decay[iso]->SetYTitle("#lambda_{0}^{2}");
681 fhELambda0MCPi0Decay[iso]->SetXTitle("E (GeV)");
682 outputContainer->Add(fhELambda0MCPi0Decay[iso]) ;
683
684 fhELambda0MCEtaDecay[iso] = new TH2F
685 (Form("hELambda0%s_MCEtaDecay",hName[iso].Data()),
686 Form("%s cluster : E vs #lambda_{0}: Origin is eta decay",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
687 fhELambda0MCEtaDecay[iso]->SetYTitle("#lambda_{0}^{2}");
688 fhELambda0MCEtaDecay[iso]->SetXTitle("E (GeV)");
689 outputContainer->Add(fhELambda0MCEtaDecay[iso]) ;
690
691 fhELambda0MCOtherDecay[iso] = new TH2F
692 (Form("hELambda0%s_MCOtherDecay",hName[iso].Data()),
693 Form("%s cluster : E vs #lambda_{0}: Origin is other decay",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
694 fhELambda0MCOtherDecay[iso]->SetYTitle("#lambda_{0}^{2}");
695 fhELambda0MCOtherDecay[iso]->SetXTitle("E (GeV)");
696 outputContainer->Add(fhELambda0MCOtherDecay[iso]) ;
697
698 fhELambda0MCHadron[iso] = new TH2F
699 (Form("hELambda0%s_MCHadron",hName[iso].Data()),
700 Form("%s cluster : E vs #lambda_{0}: Origin is hadron",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
701 fhELambda0MCHadron[iso]->SetYTitle("#lambda_{0}^{2}");
702 fhELambda0MCHadron[iso]->SetXTitle("E (GeV)");
703 outputContainer->Add(fhELambda0MCHadron[iso]) ;
704 }
705
ca134929 706 fhELambda1[iso] = new TH2F
707 (Form("hELambda1%s",hName[iso].Data()),
708 Form("%s cluster: E vs #lambda_{1}",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
709 fhELambda1[iso]->SetYTitle("#lambda_{1}^{2}");
710 fhELambda1[iso]->SetXTitle("E (GeV)");
db7b861a 711 outputContainer->Add(fhELambda1[iso]) ;
ca134929 712
713 if(fCalorimeter=="EMCAL")
714 {
715 fhELambda0TRD[iso] = new TH2F
716 (Form("hELambda0TRD%s",hName[iso].Data()),
717 Form("%s cluster: E vs #lambda_{0}, SM behind TRD",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
718 fhELambda0TRD[iso]->SetYTitle("#lambda_{0}^{2}");
719 fhELambda0TRD[iso]->SetXTitle("E (GeV)");
720 outputContainer->Add(fhELambda0TRD[iso]) ;
721
722 fhELambda1TRD[iso] = new TH2F
723 (Form("hELambda1TRD%s",hName[iso].Data()),
724 Form("%s cluster: E vs #lambda_{1}, SM behind TRD",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
725 fhELambda1TRD[iso]->SetYTitle("#lambda_{1}^{2}");
726 fhELambda1TRD[iso]->SetXTitle("E (GeV)");
727 outputContainer->Add(fhELambda1TRD[iso]) ;
728 }
729
730 fhNLocMax[iso] = new TH2F
731 (Form("hNLocMax%s",hName[iso].Data()),
732 Form("%s - Number of local maxima in cluster",hTitle[iso].Data()),
733 nptbins,ptmin,ptmax,10,0,10);
734 fhNLocMax[iso]->SetYTitle("N maxima");
735 fhNLocMax[iso]->SetXTitle("E (GeV)");
736 outputContainer->Add(fhNLocMax[iso]) ;
737
738 fhELambda0LocMax1[iso] = new TH2F
739 (Form("hELambda0LocMax1%s",hName[iso].Data()),
740 Form("%s cluster (#eta) pairs: E vs #lambda_{0}, 1 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
741 fhELambda0LocMax1[iso]->SetYTitle("#lambda_{0}^{2}");
742 fhELambda0LocMax1[iso]->SetXTitle("E (GeV)");
743 outputContainer->Add(fhELambda0LocMax1[iso]) ;
744
745 fhELambda1LocMax1[iso] = new TH2F
746 (Form("hELambda1LocMax1%s",hName[iso].Data()),
747 Form("%s cluster (#eta) pairs: E vs #lambda_{1}, 1 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
748 fhELambda1LocMax1[iso]->SetYTitle("#lambda_{1}^{2}");
749 fhELambda1LocMax1[iso]->SetXTitle("E (GeV)");
750 outputContainer->Add(fhELambda1LocMax1[iso]) ;
751
752 fhELambda0LocMax2[iso] = new TH2F
753 (Form("hELambda0LocMax2%s",hName[iso].Data()),
754 Form("%s cluster (#eta) pairs: E vs #lambda_{0}, 2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
755 fhELambda0LocMax2[iso]->SetYTitle("#lambda_{0}^{2}");
756 fhELambda0LocMax2[iso]->SetXTitle("E (GeV)");
757 outputContainer->Add(fhELambda0LocMax2[iso]) ;
758
759 fhELambda1LocMax2[iso] = new TH2F
760 (Form("hELambda1LocMax2%s",hName[iso].Data()),
761 Form("%s cluster (#eta) pairs: E vs #lambda_{1}, 2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
762 fhELambda1LocMax2[iso]->SetYTitle("#lambda_{1}^{2}");
763 fhELambda1LocMax2[iso]->SetXTitle("E (GeV)");
764 outputContainer->Add(fhELambda1LocMax2[iso]) ;
765
766 fhELambda0LocMaxN[iso] = new TH2F
767 ( Form("hELambda0LocMaxN%s",hName[iso].Data()),
768 Form("%s cluster (#eta) pairs: E vs #lambda_{0}, N>2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
769 fhELambda0LocMaxN[iso]->SetYTitle("#lambda_{0}^{2}");
770 fhELambda0LocMaxN[iso]->SetXTitle("E (GeV)");
771 outputContainer->Add(fhELambda0LocMaxN[iso]) ;
772
773 fhELambda1LocMaxN[iso] = new TH2F
774 (Form("hELambda1LocMaxN%s",hName[iso].Data()),
775 Form("%s cluster (#eta) pairs: E vs #lambda_{1}, N>2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
776 fhELambda1LocMaxN[iso]->SetYTitle("#lambda_{1}^{2}");
777 fhELambda1LocMaxN[iso]->SetXTitle("E (GeV)");
778 outputContainer->Add(fhELambda1LocMaxN[iso]) ;
b5dbb99b 779
b5dbb99b 780 }
ca134929 781 } // control histograms for isolated and non isolated objects
764ab1f4 782
db6fb352 783 fhConeSumPt = new TH2F("hConePtSum",
784 Form("#Sigma p_{T} in isolation cone for R = %2.2f",r),
785 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1a31a9ab 786 fhConeSumPt->SetYTitle("#Sigma p_{T}");
787 fhConeSumPt->SetXTitle("p_{T} (GeV/c)");
788 outputContainer->Add(fhConeSumPt) ;
789
db6fb352 790 fhPtInCone = new TH2F("hPtInCone",
791 Form("p_{T} in isolation cone for R = %2.2f",r),
792 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1a31a9ab 793 fhPtInCone->SetYTitle("p_{T in cone} (GeV/c)");
794 fhPtInCone->SetXTitle("p_{T} (GeV/c)");
795 outputContainer->Add(fhPtInCone) ;
764ab1f4 796
6c80c1bf 797 if(fFillPileUpHistograms)
798 {
799 fhPtInConePileUp = new TH2F("hPtInConePileUp",
800 Form("p_{T} in isolation cone for R = %2.2f, from pile-up (SPD)",r),
801 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
802 fhPtInConePileUp->SetYTitle("p_{T in cone} (GeV/c)");
803 fhPtInConePileUp->SetXTitle("p_{T} (GeV/c)");
804 outputContainer->Add(fhPtInConePileUp) ;
805 }
806
b7ce43b4 807 fhPtInConeCent = new TH2F("hPtInConeCent",
764ab1f4 808 Form("p_{T} in isolation cone for R = %2.2f",r),
809 100,0,100,nptinconebins,ptinconemin,ptinconemax);
b7ce43b4 810 fhPtInConeCent->SetYTitle("p_{T in cone} (GeV/c)");
811 fhPtInConeCent->SetXTitle("centrality");
812 outputContainer->Add(fhPtInConeCent) ;
1a31a9ab 813
db6fb352 814 fhFRConeSumPt = new TH2F("hFRConePtSum",
815 Form("#Sigma p_{T} in the forward region isolation cone for R = %2.2f",r),
816 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1a31a9ab 817 fhFRConeSumPt->SetYTitle("#Sigma p_{T}");
818 fhFRConeSumPt->SetXTitle("p_{T} (GeV/c)");
819 outputContainer->Add(fhFRConeSumPt) ;
820
db6fb352 821 fhPtInFRCone = new TH2F("hPtInFRCone",
822 Form("p_{T} in forward region isolation cone for R = %2.2f",r),
823 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1a31a9ab 824 fhPtInFRCone->SetYTitle("p_{T in cone} (GeV/c)");
825 fhPtInFRCone->SetXTitle("p_{T} (GeV/c)");
b7ce43b4 826 outputContainer->Add(fhPtInFRCone) ;
764ab1f4 827
b7ce43b4 828 fhPhiUEConeSumPt = new TH2F("hPhiUEConeSumPt",
764ab1f4 829 Form("p_{T} in phi band around isolation cone for R = %2.2f",r),
830 nptbins,ptmin,ptmax,4*nptinconebins,ptinconemin,4*ptinconemax);
b7ce43b4 831 fhPhiUEConeSumPt->SetYTitle("p_{T in band} (GeV/c)");
832 fhPhiUEConeSumPt->SetXTitle("p_{T} (GeV/c)");
833 outputContainer->Add(fhPhiUEConeSumPt) ;
764ab1f4 834
b7ce43b4 835 fhEtaUEConeSumPt = new TH2F("hEtaUEConeSumPt",
764ab1f4 836 Form("p_{T} in eta band around isolation cone for R = %2.2f",r),
837 nptbins,ptmin,ptmax,4*nptinconebins,ptinconemin,4*ptinconemax);
b7ce43b4 838 fhEtaUEConeSumPt->SetYTitle("p_{T in band} (GeV/c)");
839 fhEtaUEConeSumPt->SetXTitle("p_{T} (GeV/c)");
840 outputContainer->Add(fhEtaUEConeSumPt) ;
764ab1f4 841
b7ce43b4 842 fhEtaBand = new TH2F("fhEtaBand",
764ab1f4 843 Form("Eta/Phi of particle in Eta band isolation cone for R = %2.2f",r),
844 netabins,etamin,etamax,nphibins,phimin,phimax);
b7ce43b4 845 fhEtaBand->SetXTitle("#eta");
846 fhEtaBand->SetYTitle("#phi");
847 outputContainer->Add(fhEtaBand) ;
764ab1f4 848
b7ce43b4 849 fhPhiBand = new TH2F("fhPhiBand",
764ab1f4 850 Form("Eta/Phi of particle in Phi band isolation cone for R = %2.2f",r),
851 netabins,etamin,etamax,nphibins,phimin,phimax);
b7ce43b4 852 fhPhiBand->SetXTitle("#eta");
853 fhPhiBand->SetYTitle("#phi");
854 outputContainer->Add(fhPhiBand) ;
764ab1f4 855
b7ce43b4 856 fhConeSumPtEtaUESub = new TH2F("hConeSumPtEtaUESub",
764ab1f4 857 Form("#Sigma p_{T} after bkg subtraction from eta band in the isolation cone for R = %2.2f",r),
858 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
b7ce43b4 859 fhConeSumPtEtaUESub->SetYTitle("#Sigma p_{T}");
860 fhConeSumPtEtaUESub->SetXTitle("p_{T} (GeV/c)");
861 outputContainer->Add(fhConeSumPtEtaUESub) ;
764ab1f4 862
b7ce43b4 863 fhConeSumPtPhiUESub = new TH2F("hConeSumPtPhiUESub",
764ab1f4 864 Form("#Sigma p_{T} after bkg subtraction from phi band in the isolation cone for R = %2.2f",r),
865 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
b7ce43b4 866 fhConeSumPtPhiUESub->SetYTitle("#Sigma p_{T}");
867 fhConeSumPtPhiUESub->SetXTitle("p_{T} (GeV/c)");
868 outputContainer->Add(fhConeSumPtPhiUESub) ;
764ab1f4 869
db6fb352 870 fhEIso = new TH1F("hE",
871 Form("Number of isolated particles vs E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
872 nptbins,ptmin,ptmax);
0fb69ade 873 fhEIso->SetYTitle("dN / dE");
874 fhEIso->SetXTitle("E (GeV/c)");
875 outputContainer->Add(fhEIso) ;
876
db6fb352 877 fhPtIso = new TH1F("hPt",
878 Form("Number of isolated particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
879 nptbins,ptmin,ptmax);
0fb69ade 880 fhPtIso->SetYTitle("dN / p_{T}");
881 fhPtIso->SetXTitle("p_{T} (GeV/c)");
1a31a9ab 882 outputContainer->Add(fhPtIso) ;
883
8736d400 884 fhPtNLocMaxIso = new TH2F("hPtNLocMax",
885 Form("Number of isolated particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f vs NLM, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
886 nptbins,ptmin,ptmax,10,0,10);
887 fhPtNLocMaxIso->SetYTitle("NLM");
888 fhPtNLocMaxIso->SetXTitle("p_{T} (GeV/c)");
889 outputContainer->Add(fhPtNLocMaxIso) ;
890
db6fb352 891 fhPhiIso = new TH2F("hPhi",
892 Form("Number of isolated particles vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
893 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 894 fhPhiIso->SetYTitle("#phi");
895 fhPhiIso->SetXTitle("p_{T} (GeV/c)");
896 outputContainer->Add(fhPhiIso) ;
897
db6fb352 898 fhEtaIso = new TH2F("hEta",
899 Form("Number of isolated particles vs #eta for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
900 nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 901 fhEtaIso->SetYTitle("#eta");
902 fhEtaIso->SetXTitle("p_{T} (GeV/c)");
903 outputContainer->Add(fhEtaIso) ;
904
db6fb352 905 fhEtaPhiIso = new TH2F("hEtaPhiIso",
906 Form("Number of isolated particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
907 netabins,etamin,etamax,nphibins,phimin,phimax);
0fb69ade 908 fhEtaPhiIso->SetXTitle("#eta");
909 fhEtaPhiIso->SetYTitle("#phi");
910 outputContainer->Add(fhEtaPhiIso) ;
b0a31c92 911
db6fb352 912 fhPtDecayIso = new TH1F("hPtDecayIso",
913 Form("Number of isolated #pi^{0} decay particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
914 nptbins,ptmin,ptmax);
66e64043 915 fhPtDecayIso->SetYTitle("N");
916 fhPtDecayIso->SetXTitle("p_{T}(GeV/c)");
803d06a8 917 outputContainer->Add(fhPtDecayIso) ;
1a31a9ab 918
db6fb352 919 fhEtaPhiDecayIso = new TH2F("hEtaPhiDecayIso",
920 Form("Number of isolated Pi0 decay particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
921 netabins,etamin,etamax,nphibins,phimin,phimax);
d0a4f937 922 fhEtaPhiDecayIso->SetXTitle("#eta");
923 fhEtaPhiDecayIso->SetYTitle("#phi");
924 outputContainer->Add(fhEtaPhiDecayIso) ;
b0a31c92 925
b5dbb99b 926 if(IsDataMC())
927 {
732895a6 928 fhPtIsoPrompt = new TH1F("hPtMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax);
1a31a9ab 929 fhPtIsoPrompt->SetYTitle("N");
930 fhPtIsoPrompt->SetXTitle("p_{T #gamma}(GeV/c)");
931 outputContainer->Add(fhPtIsoPrompt) ;
932
933 fhPhiIsoPrompt = new TH2F
732895a6 934 ("hPhiMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 935 fhPhiIsoPrompt->SetYTitle("#phi");
936 fhPhiIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
937 outputContainer->Add(fhPhiIsoPrompt) ;
938
939 fhEtaIsoPrompt = new TH2F
732895a6 940 ("hEtaMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 941 fhEtaIsoPrompt->SetYTitle("#eta");
942 fhEtaIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
943 outputContainer->Add(fhEtaIsoPrompt) ;
944
732895a6 945 fhPtIsoFragmentation = new TH1F("hPtMCFragmentation","Number of isolated #gamma",nptbins,ptmin,ptmax);
1a31a9ab 946 fhPtIsoFragmentation->SetYTitle("N");
947 fhPtIsoFragmentation->SetXTitle("p_{T #gamma}(GeV/c)");
948 outputContainer->Add(fhPtIsoFragmentation) ;
949
950 fhPhiIsoFragmentation = new TH2F
732895a6 951 ("hPhiMCFragmentation","Number of isolated fragmentation #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 952 fhPhiIsoFragmentation->SetYTitle("#phi");
953 fhPhiIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
954 outputContainer->Add(fhPhiIsoFragmentation) ;
955
956 fhEtaIsoFragmentation = new TH2F
732895a6 957 ("hEtaMCFragmentation","Number of isolated fragmentation #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 958 fhEtaIsoFragmentation->SetYTitle("#eta");
959 fhEtaIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
960 outputContainer->Add(fhEtaIsoFragmentation) ;
961
764ab1f4 962 fhPtIsoPi0 = new TH1F("hPtMCPi0","Number of isolated #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax);
963 fhPtIsoPi0->SetYTitle("N");
964 fhPtIsoPi0->SetXTitle("p_{T #gamma}(GeV/c)");
965 outputContainer->Add(fhPtIsoPi0) ;
966
967 fhPhiIsoPi0 = new TH2F
968 ("hPhiMCPi0","Number of isolated #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
969 fhPhiIsoPi0->SetYTitle("#phi");
970 fhPhiIsoPi0->SetXTitle("p_{T #gamma} (GeV/c)");
971 outputContainer->Add(fhPhiIsoPi0) ;
972
973 fhEtaIsoPi0 = new TH2F
974 ("hEtaMCPi0","Number of isolated #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax,netabins,etamin,etamax);
975 fhEtaIsoPi0->SetYTitle("#eta");
976 fhEtaIsoPi0->SetXTitle("p_{T #gamma} (GeV/c)");
977 outputContainer->Add(fhEtaIsoPi0) ;
978
732895a6 979 fhPtIsoPi0Decay = new TH1F("hPtMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax);
1a31a9ab 980 fhPtIsoPi0Decay->SetYTitle("N");
981 fhPtIsoPi0Decay->SetXTitle("p_{T #gamma}(GeV/c)");
982 outputContainer->Add(fhPtIsoPi0Decay) ;
983
984 fhPhiIsoPi0Decay = new TH2F
732895a6 985 ("hPhiMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 986 fhPhiIsoPi0Decay->SetYTitle("#phi");
987 fhPhiIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
988 outputContainer->Add(fhPhiIsoPi0Decay) ;
989
990 fhEtaIsoPi0Decay = new TH2F
732895a6 991 ("hEtaMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 992 fhEtaIsoPi0Decay->SetYTitle("#eta");
993 fhEtaIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
994 outputContainer->Add(fhEtaIsoPi0Decay) ;
995
803d06a8 996 fhPtIsoEtaDecay = new TH1F("hPtMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax);
997 fhPtIsoEtaDecay->SetYTitle("N");
998 fhPtIsoEtaDecay->SetXTitle("p_{T #gamma}(GeV/c)");
999 outputContainer->Add(fhPtIsoEtaDecay) ;
1000
1001 fhPhiIsoEtaDecay = new TH2F
1002 ("hPhiMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1003 fhPhiIsoEtaDecay->SetYTitle("#phi");
1004 fhPhiIsoEtaDecay->SetXTitle("p_{T #gamma} (GeV/c)");
1005 outputContainer->Add(fhPhiIsoEtaDecay) ;
1006
1007 fhEtaIsoEtaDecay = new TH2F
1008 ("hEtaMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1009 fhEtaIsoEtaDecay->SetYTitle("#eta");
1010 fhEtaIsoEtaDecay->SetXTitle("p_{T #gamma} (GeV/c)");
1011 outputContainer->Add(fhEtaIsoEtaDecay) ;
1012
732895a6 1013 fhPtIsoOtherDecay = new TH1F("hPtMCOtherDecay","Number of isolated #gamma from non-#pi^{0} decay",nptbins,ptmin,ptmax);
1a31a9ab 1014 fhPtIsoOtherDecay->SetYTitle("N");
1015 fhPtIsoOtherDecay->SetXTitle("p_{T #gamma}(GeV/c)");
1016 outputContainer->Add(fhPtIsoOtherDecay) ;
1017
1018 fhPhiIsoOtherDecay = new TH2F
732895a6 1019 ("hPhiMCOtherDecay","Number of isolated #gamma from non-#pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1a31a9ab 1020 fhPhiIsoOtherDecay->SetYTitle("#phi");
1021 fhPhiIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
1022 outputContainer->Add(fhPhiIsoOtherDecay) ;
1023
1024 fhEtaIsoOtherDecay = new TH2F
732895a6 1025 ("hEtaMCOtherDecay","Number of isolated #gamma non-#pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1a31a9ab 1026 fhEtaIsoOtherDecay->SetYTitle("#eta");
1027 fhEtaIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
1028 outputContainer->Add(fhEtaIsoOtherDecay) ;
1029
764ab1f4 1030 // fhPtIsoConversion = new TH1F("hPtMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax);
1031 // fhPtIsoConversion->SetYTitle("N");
1032 // fhPtIsoConversion->SetXTitle("p_{T #gamma}(GeV/c)");
1033 // outputContainer->Add(fhPtIsoConversion) ;
1034 //
1035 // fhPhiIsoConversion = new TH2F
1036 // ("hPhiMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1037 // fhPhiIsoConversion->SetYTitle("#phi");
1038 // fhPhiIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
1039 // outputContainer->Add(fhPhiIsoConversion) ;
1040 //
1041 // fhEtaIsoConversion = new TH2F
1042 // ("hEtaMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1043 // fhEtaIsoConversion->SetYTitle("#eta");
1044 // fhEtaIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
1045 // outputContainer->Add(fhEtaIsoConversion) ;
1a31a9ab 1046
764ab1f4 1047 fhPtIsoHadron = new TH1F("hPtMCHadron","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax);
1048 fhPtIsoHadron->SetYTitle("N");
1049 fhPtIsoHadron->SetXTitle("p_{T}(GeV/c)");
1050 outputContainer->Add(fhPtIsoHadron) ;
1a31a9ab 1051
8736d400 1052
764ab1f4 1053 fhPhiIsoHadron = new TH2F
1054 ("hPhiMCHadron","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1055 fhPhiIsoHadron->SetYTitle("#phi");
1056 fhPhiIsoHadron->SetXTitle("p_{T} (GeV/c)");
1057 outputContainer->Add(fhPhiIsoHadron) ;
1a31a9ab 1058
764ab1f4 1059 fhEtaIsoHadron = new TH2F
1060 ("hEtaMCHadron","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1061 fhEtaIsoHadron->SetYTitle("#eta");
1062 fhEtaIsoHadron->SetXTitle("p_{T} (GeV/c)");
1063 outputContainer->Add(fhEtaIsoHadron) ;
1a31a9ab 1064
1a31a9ab 1065 }//Histos with MC
1066
1067 }
1068
db6fb352 1069 // Not Isolated histograms, reference histograms
1070
b1f720a7 1071 fhENoIso = new TH1F("hENoIso",
1072 Form("Number of not isolated leading particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1073 nptbins,ptmin,ptmax);
1074 fhENoIso->SetYTitle("N");
1075 fhENoIso->SetXTitle("p_{T}(GeV/c)");
1076 outputContainer->Add(fhENoIso) ;
1077
db6fb352 1078 fhPtNoIso = new TH1F("hPtNoIso",
1079 Form("Number of not isolated leading particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1080 nptbins,ptmin,ptmax);
1081 fhPtNoIso->SetYTitle("N");
1082 fhPtNoIso->SetXTitle("p_{T}(GeV/c)");
1083 outputContainer->Add(fhPtNoIso) ;
1084
8736d400 1085 fhPtNLocMaxNoIso = new TH2F("hPtNLocMaxNoIso",
1086 Form("Number of not isolated particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f vs NLM, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1087 nptbins,ptmin,ptmax,10,0,10);
1088 fhPtNLocMaxNoIso->SetYTitle("NLM");
1089 fhPtNLocMaxNoIso->SetXTitle("p_{T} (GeV/c)");
1090 outputContainer->Add(fhPtNLocMaxNoIso) ;
db6fb352 1091
1092 fhEtaPhiNoIso = new TH2F("hEtaPhiNoIso",
1093 Form("Number of not isolated leading particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1094 netabins,etamin,etamax,nphibins,phimin,phimax);
1095 fhEtaPhiNoIso->SetXTitle("#eta");
1096 fhEtaPhiNoIso->SetYTitle("#phi");
1097 outputContainer->Add(fhEtaPhiNoIso) ;
1098
1099 fhPtDecayNoIso = new TH1F("hPtDecayNoIso",
1100 Form("Number of not isolated leading pi0 decay particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1101 nptbins,ptmin,ptmax);
1102 fhPtDecayNoIso->SetYTitle("N");
1103 fhPtDecayNoIso->SetXTitle("p_{T}(GeV/c)");
1104 outputContainer->Add(fhPtDecayNoIso) ;
1105
1106 fhEtaPhiDecayNoIso = new TH2F("hEtaPhiDecayNoIso",
1107 Form("Number of not isolated leading Pi0 decay particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1108 netabins,etamin,etamax,nphibins,phimin,phimax);
1109 fhEtaPhiDecayNoIso->SetXTitle("#eta");
1110 fhEtaPhiDecayNoIso->SetYTitle("#phi");
1111 outputContainer->Add(fhEtaPhiDecayNoIso) ;
1112
44e48e82 1113
1114
db6fb352 1115 if(IsDataMC())
1116 {
764ab1f4 1117 fhPtNoIsoPi0 = new TH1F
1118 ("hPtNoIsoPi0","Number of not isolated leading #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax);
1119 fhPtNoIsoPi0->SetYTitle("N");
1120 fhPtNoIsoPi0->SetXTitle("p_{T} (GeV/c)");
1121 outputContainer->Add(fhPtNoIsoPi0) ;
1122
db6fb352 1123 fhPtNoIsoPi0Decay = new TH1F
1124 ("hPtNoIsoPi0Decay","Number of not isolated leading #gamma from #pi^{0} decay",nptbins,ptmin,ptmax);
1125 fhPtNoIsoPi0Decay->SetYTitle("N");
1126 fhPtNoIsoPi0Decay->SetXTitle("p_{T} (GeV/c)");
1127 outputContainer->Add(fhPtNoIsoPi0Decay) ;
1128
1129 fhPtNoIsoEtaDecay = new TH1F
1130 ("hPtNoIsoEtaDecay","Number of not isolated leading #gamma from eta decay",nptbins,ptmin,ptmax);
1131 fhPtNoIsoEtaDecay->SetYTitle("N");
1132 fhPtNoIsoEtaDecay->SetXTitle("p_{T} (GeV/c)");
1133 outputContainer->Add(fhPtNoIsoEtaDecay) ;
1134
1135 fhPtNoIsoOtherDecay = new TH1F
1136 ("hPtNoIsoOtherDecay","Number of not isolated leading #gamma from other decay",nptbins,ptmin,ptmax);
1137 fhPtNoIsoOtherDecay->SetYTitle("N");
1138 fhPtNoIsoOtherDecay->SetXTitle("p_{T} (GeV/c)");
1139 outputContainer->Add(fhPtNoIsoOtherDecay) ;
1140
1141 fhPtNoIsoPrompt = new TH1F
1142 ("hPtNoIsoPrompt","Number of not isolated leading prompt #gamma",nptbins,ptmin,ptmax);
1143 fhPtNoIsoPrompt->SetYTitle("N");
1144 fhPtNoIsoPrompt->SetXTitle("p_{T} (GeV/c)");
1145 outputContainer->Add(fhPtNoIsoPrompt) ;
1146
1147 fhPtIsoMCPhoton = new TH1F
1148 ("hPtIsoMCPhoton","Number of isolated leading #gamma",nptbins,ptmin,ptmax);
1149 fhPtIsoMCPhoton->SetYTitle("N");
1150 fhPtIsoMCPhoton->SetXTitle("p_{T} (GeV/c)");
1151 outputContainer->Add(fhPtIsoMCPhoton) ;
1152
1153 fhPtNoIsoMCPhoton = new TH1F
1154 ("hPtNoIsoMCPhoton","Number of not isolated leading #gamma",nptbins,ptmin,ptmax);
1155 fhPtNoIsoMCPhoton->SetYTitle("N");
1156 fhPtNoIsoMCPhoton->SetXTitle("p_{T} (GeV/c)");
1157 outputContainer->Add(fhPtNoIsoMCPhoton) ;
1158
764ab1f4 1159// fhPtNoIsoConversion = new TH1F
1160// ("hPtNoIsoConversion","Number of not isolated leading conversion #gamma",nptbins,ptmin,ptmax);
1161// fhPtNoIsoConversion->SetYTitle("N");
1162// fhPtNoIsoConversion->SetXTitle("p_{T} (GeV/c)");
1163// outputContainer->Add(fhPtNoIsoConversion) ;
db6fb352 1164
1165 fhPtNoIsoFragmentation = new TH1F
1166 ("hPtNoIsoFragmentation","Number of not isolated leading fragmentation #gamma",nptbins,ptmin,ptmax);
1167 fhPtNoIsoFragmentation->SetYTitle("N");
1168 fhPtNoIsoFragmentation->SetXTitle("p_{T} (GeV/c)");
1169 outputContainer->Add(fhPtNoIsoFragmentation) ;
1170
764ab1f4 1171 fhPtNoIsoHadron = new TH1F
1172 ("hPtNoIsoHadron","Number of not isolated leading hadrons",nptbins,ptmin,ptmax);
1173 fhPtNoIsoHadron->SetYTitle("N");
1174 fhPtNoIsoHadron->SetXTitle("p_{T} (GeV/c)");
1175 outputContainer->Add(fhPtNoIsoHadron) ;
db6fb352 1176
1177 }//Histos with MC
1178
1179
b5dbb99b 1180 if(fMakeSeveralIC)
1181 {
1a31a9ab 1182 const Int_t buffersize = 255;
e4ef72be 1183 char name[buffersize];
1184 char title[buffersize];
1185 for(Int_t icone = 0; icone<fNCones; icone++)
44e48e82 1186 {
1187 // sum pt in cone vs. pt leading
1188 snprintf(name, buffersize,"hSumPtLeadingPt_Cone_%d",icone);
1189 snprintf(title, buffersize,"#Sigma p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1190 fhSumPtLeadingPt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1191 fhSumPtLeadingPt[icone] ->SetYTitle("#sum_{cone}p_{T} (GeV/c)");//#Sigma p_{T}
1192 fhSumPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1193 outputContainer->Add(fhSumPtLeadingPt[icone]) ;
1194
1195 // pt in cone vs. pt leading
1196 snprintf(name, buffersize,"hPtLeadingPt_Cone_%d",icone);
1197 snprintf(title, buffersize,"p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1198 fhPtLeadingPt[icone] = new TH2F(name, title, nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1199 fhPtLeadingPt[icone] ->SetYTitle("p_{T}^{cone} (GeV/c)");
1200 fhPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1201 outputContainer->Add(fhPtLeadingPt[icone]) ;
1202
1203 // sum pt in cone vs. pt leading in the forward region (for background subtraction studies)
1204 snprintf(name, buffersize,"hFRSumPtLeadingPt_Cone_%d",icone);
1205 snprintf(title, buffersize,"#Sigma p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1206 fhFRSumPtLeadingPt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1207 fhFRSumPtLeadingPt[icone] ->SetYTitle("#sum_{cone}p_{T} (GeV/c)");//#Sigma p_{T}
1208 fhFRSumPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1209 outputContainer->Add(fhFRSumPtLeadingPt[icone]) ;
1210
1211 // pt in cone vs. pt leading in the forward region (for background subtraction studies)
1212 snprintf(name, buffersize,"hFRPtLeadingPt_Cone_%d",icone);
1213 snprintf(title, buffersize,"p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1214 fhFRPtLeadingPt[icone] = new TH2F(name, title, nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1215 fhFRPtLeadingPt[icone] ->SetYTitle("p_{T}^{cone} (GeV/c)");
1216 fhFRPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1217 outputContainer->Add(fhFRPtLeadingPt[icone]) ;
1218
1219
e4ef72be 1220 if(IsDataMC())
db6fb352 1221 {
b0a31c92 1222 snprintf(name, buffersize,"hPtSumPrompt_Cone_%d",icone);
1223 snprintf(title, buffersize,"Candidate Prompt cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1224 fhPtSumIsolatedPrompt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1225 fhPtSumIsolatedPrompt[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1226 fhPtSumIsolatedPrompt[icone]->SetXTitle("p_{T} (GeV/c)");
1227 outputContainer->Add(fhPtSumIsolatedPrompt[icone]) ;
1228
1229 snprintf(name, buffersize,"hPtSumFragmentation_Cone_%d",icone);
1230 snprintf(title, buffersize,"Candidate Fragmentation cone sum p_{T} for R = %2.2fvs candidate p_{T}",fConeSizes[icone]);
1231 fhPtSumIsolatedFragmentation[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1232 fhPtSumIsolatedFragmentation[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1233 fhPtSumIsolatedFragmentation[icone]->SetXTitle("p_{T} (GeV/c)");
1234 outputContainer->Add(fhPtSumIsolatedFragmentation[icone]) ;
1235
764ab1f4 1236 snprintf(name, buffersize,"hPtSumPi0_Cone_%d",icone);
1237 snprintf(title, buffersize,"Candidate Pi0 cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1238 fhPtSumIsolatedPi0[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1239 fhPtSumIsolatedPi0[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1240 fhPtSumIsolatedPi0[icone]->SetXTitle("p_{T} (GeV/c)");
1241 outputContainer->Add(fhPtSumIsolatedPi0[icone]) ;
1242
b0a31c92 1243 snprintf(name, buffersize,"hPtSumPi0Decay_Cone_%d",icone);
1244 snprintf(title, buffersize,"Candidate Pi0Decay cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1245 fhPtSumIsolatedPi0Decay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1246 fhPtSumIsolatedPi0Decay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1247 fhPtSumIsolatedPi0Decay[icone]->SetXTitle("p_{T} (GeV/c)");
1248 outputContainer->Add(fhPtSumIsolatedPi0Decay[icone]) ;
1249
1250 snprintf(name, buffersize,"hPtSumEtaDecay_Cone_%d",icone);
1251 snprintf(title, buffersize,"Candidate EtaDecay cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1252 fhPtSumIsolatedEtaDecay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1253 fhPtSumIsolatedEtaDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1254 fhPtSumIsolatedEtaDecay[icone]->SetXTitle("p_{T} (GeV/c)");
1255 outputContainer->Add(fhPtSumIsolatedEtaDecay[icone]) ;
1256
1257 snprintf(name, buffersize,"hPtSumOtherDecay_Cone_%d",icone);
1258 snprintf(title, buffersize,"Candidate OtherDecay cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1259 fhPtSumIsolatedOtherDecay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1260 fhPtSumIsolatedOtherDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1261 fhPtSumIsolatedOtherDecay[icone]->SetXTitle("p_{T} (GeV/c)");
1262 outputContainer->Add(fhPtSumIsolatedOtherDecay[icone]) ;
1263
764ab1f4 1264// snprintf(name, buffersize,"hPtSumConversion_Cone_%d",icone);
1265// snprintf(title, buffersize,"Candidate Conversion cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1266// fhPtSumIsolatedConversion[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1267// fhPtSumIsolatedConversion[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1268// fhPtSumIsolatedConversion[icone]->SetXTitle("p_{T} (GeV/c)");
1269// outputContainer->Add(fhPtSumIsolatedConversion[icone]) ;
b0a31c92 1270
764ab1f4 1271 snprintf(name, buffersize,"hPtSumHadron_Cone_%d",icone);
1272 snprintf(title, buffersize,"Candidate Hadron cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1273 fhPtSumIsolatedHadron[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1274 fhPtSumIsolatedHadron[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1275 fhPtSumIsolatedHadron[icone]->SetXTitle("p_{T} (GeV/c)");
1276 outputContainer->Add(fhPtSumIsolatedHadron[icone]) ;
b0a31c92 1277
e4ef72be 1278 }//Histos with MC
1279
1280 for(Int_t ipt = 0; ipt<fNPtThresFrac;ipt++)
1281 {
44e48e82 1282
b0a31c92 1283 snprintf(name, buffersize,"hPtThres_Cone_%d_Pt%d",icone,ipt);
1284 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1285 fhPtThresIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1286 fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1287 outputContainer->Add(fhPtThresIsolated[icone][ipt]) ;
1288
1289 snprintf(name, buffersize,"hPtFrac_Cone_%d_Pt%d",icone,ipt);
1290 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1291 fhPtFracIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1292 fhPtFracIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1293 outputContainer->Add(fhPtFracIsolated[icone][ipt]) ;
1294
1295
1296 snprintf(name, buffersize,"hPtSum_Cone_%d_Pt%d",icone,ipt);
1297 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1298 fhPtSumIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1299 // fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma p_{T} (GeV/c)");
1300 fhPtSumIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1301 outputContainer->Add(fhPtSumIsolated[icone][ipt]) ;
1302
1303 snprintf(name, buffersize,"hPtSumDensity_Cone_%d_Pt%d",icone,ipt);
1304 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for density in R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1305 fhPtSumDensityIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1306 //fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma p_{T} (GeV/c)");
1307 fhPtSumDensityIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1308 outputContainer->Add(fhPtSumDensityIso[icone][ipt]) ;
1309
1310 snprintf(name, buffersize,"hPtFracPtSum_Cone_%d_Pt%d",icone,ipt);
1311 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for PtFracPtSum in R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1312 fhPtFracPtSumIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1313 //fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma p_{T} (GeV/c)");
1314 fhPtFracPtSumIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1315 outputContainer->Add(fhPtFracPtSumIso[icone][ipt]) ;
1316
1317 // pt decays isolated
1318 snprintf(name, buffersize,"hPtThres_Decay_Cone_%d_Pt%d",icone,ipt);
1319 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1320 fhPtPtThresDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1321 fhPtPtThresDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1322 outputContainer->Add(fhPtPtThresDecayIso[icone][ipt]) ;
1323
1324 snprintf(name, buffersize,"hPtFrac_Decay_Cone_%d_Pt%d",icone,ipt);
1325 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1326 fhPtPtFracDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1327 fhPtPtFracDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1328 outputContainer->Add(fhPtPtFracDecayIso[icone][ipt]) ;
1329
1330 snprintf(name, buffersize,"hPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1331 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1332 fhPtPtSumDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1333 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1334 fhPtPtSumDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1335 outputContainer->Add(fhPtPtSumDecayIso[icone][ipt]) ;
1336
1337 snprintf(name, buffersize,"hPtSumDensity_Decay_Cone_%d_Pt%d",icone,ipt);
1338 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for density in R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1339 fhPtSumDensityDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1340 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1341 fhPtSumDensityDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1342 outputContainer->Add(fhPtSumDensityDecayIso[icone][ipt]) ;
1343
1344 snprintf(name, buffersize,"hPtFracPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1345 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for PtFracPtSum in R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1346 fhPtFracPtSumDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1347 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1348 fhPtFracPtSumDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1349 outputContainer->Add(fhPtFracPtSumDecayIso[icone][ipt]) ;
1350
1351
1352 // eta:phi
1353 snprintf(name, buffersize,"hEtaPhiPtThres_Cone_%d_Pt%d",icone,ipt);
1354 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1355 fhEtaPhiPtThresIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1356 fhEtaPhiPtThresIso[icone][ipt]->SetXTitle("#eta");
1357 fhEtaPhiPtThresIso[icone][ipt]->SetYTitle("#phi");
1358 outputContainer->Add(fhEtaPhiPtThresIso[icone][ipt]) ;
1359
1360 snprintf(name, buffersize,"hEtaPhiPtFrac_Cone_%d_Pt%d",icone,ipt);
1361 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1362 fhEtaPhiPtFracIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1363 fhEtaPhiPtFracIso[icone][ipt]->SetXTitle("#eta");
1364 fhEtaPhiPtFracIso[icone][ipt]->SetYTitle("#phi");
1365 outputContainer->Add(fhEtaPhiPtFracIso[icone][ipt]) ;
1366
1367 snprintf(name, buffersize,"hEtaPhiPtSum_Cone_%d_Pt%d",icone,ipt);
1368 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1369 fhEtaPhiPtSumIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1370 fhEtaPhiPtSumIso[icone][ipt]->SetXTitle("#eta");
1371 fhEtaPhiPtSumIso[icone][ipt]->SetYTitle("#phi");
1372 outputContainer->Add(fhEtaPhiPtSumIso[icone][ipt]) ;
1373
1374 snprintf(name, buffersize,"hEtaPhiSumDensity_Cone_%d_Pt%d",icone,ipt);
1375 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for density R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1376 fhEtaPhiSumDensityIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1377 fhEtaPhiSumDensityIso[icone][ipt]->SetXTitle("#eta");
1378 fhEtaPhiSumDensityIso[icone][ipt]->SetYTitle("#phi");
1379 outputContainer->Add(fhEtaPhiSumDensityIso[icone][ipt]) ;
1380
1381 snprintf(name, buffersize,"hEtaPhiFracPtSum_Cone_%d_Pt%d",icone,ipt);
1382 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for FracPtSum R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1383 fhEtaPhiFracPtSumIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1384 fhEtaPhiFracPtSumIso[icone][ipt]->SetXTitle("#eta");
1385 fhEtaPhiFracPtSumIso[icone][ipt]->SetYTitle("#phi");
1386 outputContainer->Add(fhEtaPhiFracPtSumIso[icone][ipt]) ;
1387
1388 // eta:phi decays
1389 snprintf(name, buffersize,"hEtaPhiPtThres_Decay_Cone_%d_Pt%d",icone,ipt);
1390 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1391 fhEtaPhiPtThresDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1392 fhEtaPhiPtThresDecayIso[icone][ipt]->SetXTitle("#eta");
1393 fhEtaPhiPtThresDecayIso[icone][ipt]->SetYTitle("#phi");
1394 outputContainer->Add(fhEtaPhiPtThresDecayIso[icone][ipt]) ;
1395
1396 snprintf(name, buffersize,"hEtaPhiPtFrac_Decay_Cone_%d_Pt%d",icone,ipt);
1397 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1398 fhEtaPhiPtFracDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1399 fhEtaPhiPtFracDecayIso[icone][ipt]->SetXTitle("#eta");
1400 fhEtaPhiPtFracDecayIso[icone][ipt]->SetYTitle("#phi");
1401 outputContainer->Add(fhEtaPhiPtFracDecayIso[icone][ipt]) ;
1402
1403
1404 snprintf(name, buffersize,"hEtaPhiPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1405 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1406 fhEtaPhiPtSumDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1407 fhEtaPhiPtSumDecayIso[icone][ipt]->SetXTitle("#eta");
1408 fhEtaPhiPtSumDecayIso[icone][ipt]->SetYTitle("#phi");
1409 outputContainer->Add(fhEtaPhiPtSumDecayIso[icone][ipt]) ;
1410
1411 snprintf(name, buffersize,"hEtaPhiSumDensity_Decay_Cone_%d_Pt%d",icone,ipt);
1412 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for density R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1413 fhEtaPhiSumDensityDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1414 fhEtaPhiSumDensityDecayIso[icone][ipt]->SetXTitle("#eta");
1415 fhEtaPhiSumDensityDecayIso[icone][ipt]->SetYTitle("#phi");
1416 outputContainer->Add(fhEtaPhiSumDensityDecayIso[icone][ipt]) ;
1417
1418 snprintf(name, buffersize,"hEtaPhiFracPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1419 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for FracPtSum R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1420 fhEtaPhiFracPtSumDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1421 fhEtaPhiFracPtSumDecayIso[icone][ipt]->SetXTitle("#eta");
1422 fhEtaPhiFracPtSumDecayIso[icone][ipt]->SetYTitle("#phi");
1423 outputContainer->Add(fhEtaPhiFracPtSumDecayIso[icone][ipt]) ;
1424
1425
1426 if(IsDataMC())
1427 {
1428 snprintf(name, buffersize,"hPtThresMCPrompt_Cone_%d_Pt%d",icone,ipt);
1429 snprintf(title, buffersize,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1430 fhPtThresIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1431 fhPtThresIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1432 outputContainer->Add(fhPtThresIsolatedPrompt[icone][ipt]) ;
1433
1434 snprintf(name, buffersize,"hPtFracMCPrompt_Cone_%d_Pt%d",icone,ipt);
1435 snprintf(title, buffersize,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1436 fhPtFracIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1437 fhPtFracIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1438 outputContainer->Add(fhPtFracIsolatedPrompt[icone][ipt]) ;
1439
1440 snprintf(name, buffersize,"hPtThresMCFragmentation_Cone_%d_Pt%d",icone,ipt);
1441 snprintf(title, buffersize,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1442 fhPtThresIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1443 fhPtThresIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1444 outputContainer->Add(fhPtThresIsolatedFragmentation[icone][ipt]) ;
1445
1446 snprintf(name, buffersize,"hPtFracMCFragmentation_Cone_%d_Pt%d",icone,ipt);
1447 snprintf(title, buffersize,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1448 fhPtFracIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1449 fhPtFracIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1450 outputContainer->Add(fhPtFracIsolatedFragmentation[icone][ipt]) ;
1451
764ab1f4 1452 snprintf(name, buffersize,"hPtThresMCPi0_Cone_%d_Pt%d",icone,ipt);
1453 snprintf(title, buffersize,"Isolated candidate Pi0 p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1454 fhPtThresIsolatedPi0[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1455 fhPtThresIsolatedPi0[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1456 outputContainer->Add(fhPtThresIsolatedPi0[icone][ipt]) ;
1457
1458 snprintf(name, buffersize,"hPtFracMCPi0_Cone_%d_Pt%d",icone,ipt);
1459 snprintf(title, buffersize,"Isolated candidate Pi0 p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1460 fhPtFracIsolatedPi0[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1461 fhPtFracIsolatedPi0[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1462 outputContainer->Add(fhPtFracIsolatedPi0[icone][ipt]) ;
1463
b0a31c92 1464 snprintf(name, buffersize,"hPtThresMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
1465 snprintf(title, buffersize,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1466 fhPtThresIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1467 fhPtThresIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1468 outputContainer->Add(fhPtThresIsolatedPi0Decay[icone][ipt]) ;
1469
1470 snprintf(name, buffersize,"hPtFracMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
1471 snprintf(title, buffersize,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1472 fhPtFracIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1473 fhPtFracIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1474 outputContainer->Add(fhPtFracIsolatedPi0Decay[icone][ipt]) ;
1475
1476 snprintf(name, buffersize,"hPtThresMCEtaDecay_Cone_%d_Pt%d",icone,ipt);
1477 snprintf(title, buffersize,"Isolated candidate EtaDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1478 fhPtThresIsolatedEtaDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1479 fhPtThresIsolatedEtaDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1480 outputContainer->Add(fhPtThresIsolatedEtaDecay[icone][ipt]) ;
1481
1482 snprintf(name, buffersize,"hPtFracMCEtaDecay_Cone_%d_Pt%d",icone,ipt);
1483 snprintf(title, buffersize,"Isolated candidate EtaDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1484 fhPtFracIsolatedEtaDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1485 fhPtFracIsolatedEtaDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1486 outputContainer->Add(fhPtFracIsolatedEtaDecay[icone][ipt]) ;
1487
1488
1489 snprintf(name, buffersize,"hPtThresMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
1490 snprintf(title, buffersize,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1491 fhPtThresIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1492 fhPtThresIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1493 outputContainer->Add(fhPtThresIsolatedOtherDecay[icone][ipt]) ;
1494
1495 snprintf(name, buffersize,"hPtFracMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
1496 snprintf(title, buffersize,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1497 fhPtFracIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1498 fhPtFracIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1499 outputContainer->Add(fhPtFracIsolatedOtherDecay[icone][ipt]) ;
1500
764ab1f4 1501// snprintf(name, buffersize,"hPtThresMCConversion_Cone_%d_Pt%d",icone,ipt);
1502// snprintf(title, buffersize,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1503// fhPtThresIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1504// fhPtThresIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1505// outputContainer->Add(fhPtThresIsolatedConversion[icone][ipt]) ;
b0a31c92 1506
764ab1f4 1507// snprintf(name, buffersize,"hPtFracMCConversion_Cone_%d_Pt%d",icone,ipt);
1508// snprintf(title, buffersize,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1509// fhPtFracIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1510// fhPtFracIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1511// outputContainer->Add(fhPtFracIsolatedConversion[icone][ipt]) ;
b0a31c92 1512
764ab1f4 1513 snprintf(name, buffersize,"hPtThresMCHadron_Cone_%d_Pt%d",icone,ipt);
1514 snprintf(title, buffersize,"Isolated candidate Hadron p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1515 fhPtThresIsolatedHadron[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1516 fhPtThresIsolatedHadron[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1517 outputContainer->Add(fhPtThresIsolatedHadron[icone][ipt]) ;
b0a31c92 1518
764ab1f4 1519 snprintf(name, buffersize,"hPtFracMCHadron_Cone_%d_Pt%d",icone,ipt);
1520 snprintf(title, buffersize,"Isolated candidate Hadron p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1521 fhPtFracIsolatedHadron[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1522 fhPtFracIsolatedHadron[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1523 outputContainer->Add(fhPtFracIsolatedHadron[icone][ipt]) ;
b0a31c92 1524
e4ef72be 1525 }//Histos with MC
1526 }//icone loop
1527 }//ipt loop
1a31a9ab 1528 }
1529
2ad19c3d 1530 if(fFillPileUpHistograms)
1531 {
b1f720a7 1532
1533 fhEIsoPileUp = new TH1F("hEPileUp",
1534 Form("Number of isolated particles vs E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f, pile-up event by SDD",r,ptthre,ptfrac),
1535 nptbins,ptmin,ptmax);
1536 fhEIsoPileUp->SetYTitle("dN / dE");
1537 fhEIsoPileUp->SetXTitle("E (GeV/c)");
1538 outputContainer->Add(fhEIsoPileUp) ;
1539
1540 fhPtIsoPileUp = new TH1F("hPtPileUp",
1541 Form("Number of isolated particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr}, pile-up event by SDD = %2.2f",r,ptthre,ptfrac),
1542 nptbins,ptmin,ptmax);
1543 fhPtIsoPileUp->SetYTitle("dN / p_{T}");
1544 fhPtIsoPileUp->SetXTitle("p_{T} (GeV/c)");
1545 outputContainer->Add(fhPtIsoPileUp) ;
1546
1547 fhENoIsoPileUp = new TH1F("hENoIsoPileUp",
1548 Form("Number of not isolated particles vs E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f, pile-up event by SDD",r,ptthre,ptfrac),
1549 nptbins,ptmin,ptmax);
1550 fhENoIsoPileUp->SetYTitle("dN / dE");
1551 fhENoIsoPileUp->SetXTitle("E (GeV/c)");
1552 outputContainer->Add(fhENoIsoPileUp) ;
1553
1554 fhPtNoIsoPileUp = new TH1F("hPtNoIsoPileUp",
1555 Form("Number of not isolated particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr}, pile-up event by SDD = %2.2f",r,ptthre,ptfrac),
1556 nptbins,ptmin,ptmax);
1557 fhPtNoIsoPileUp->SetYTitle("dN / p_{T}");
1558 fhPtNoIsoPileUp->SetXTitle("p_{T} (GeV/c)");
1559 outputContainer->Add(fhPtNoIsoPileUp) ;
1560
2ad19c3d 1561 fhTimeENoCut = new TH2F ("hTimeE_NoCut","time of cluster vs E of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1562 fhTimeENoCut->SetXTitle("E (GeV)");
1563 fhTimeENoCut->SetYTitle("time (ns)");
1564 outputContainer->Add(fhTimeENoCut);
1565
1566 fhTimeESPD = new TH2F ("hTimeE_SPD","time of cluster vs E of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1567 fhTimeESPD->SetXTitle("E (GeV)");
1568 fhTimeESPD->SetYTitle("time (ns)");
1569 outputContainer->Add(fhTimeESPD);
1570
1571 fhTimeESPDMulti = new TH2F ("hTimeE_SPDMulti","time of cluster vs E of clusters, SPD multi cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1572 fhTimeESPDMulti->SetXTitle("E (GeV)");
1573 fhTimeESPDMulti->SetYTitle("time (ns)");
1574 outputContainer->Add(fhTimeESPDMulti);
1575
1576 fhTimeNPileUpVertSPD = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
1577 fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
1578 fhTimeNPileUpVertSPD->SetXTitle("time (ns)");
1579 outputContainer->Add(fhTimeNPileUpVertSPD);
1580
1581 fhTimeNPileUpVertTrack = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 50,0,50 );
1582 fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
1583 fhTimeNPileUpVertTrack->SetXTitle("time (ns)");
1584 outputContainer->Add(fhTimeNPileUpVertTrack);
1585
1586 fhTimeNPileUpVertContributors = new TH2F ("hTime_NPileUpVertContributors","time of cluster vs N constributors to pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
1587 fhTimeNPileUpVertContributors->SetYTitle("# vertex ");
1588 fhTimeNPileUpVertContributors->SetXTitle("time (ns)");
1589 outputContainer->Add(fhTimeNPileUpVertContributors);
1590
1591 fhTimePileUpMainVertexZDistance = new TH2F ("hTime_PileUpMainVertexZDistance","time of cluster vs distance in Z pile-up SPD vertex - main SPD vertex",ntimebins,timemin,timemax,100,0,50);
1592 fhTimePileUpMainVertexZDistance->SetYTitle("distance Z (cm) ");
1593 fhTimePileUpMainVertexZDistance->SetXTitle("time (ns)");
1594 outputContainer->Add(fhTimePileUpMainVertexZDistance);
1595
1596 fhTimePileUpMainVertexZDiamond = new TH2F ("hTime_PileUpMainVertexZDiamond","time of cluster vs distance in Z pile-up SPD vertex - z diamond",ntimebins,timemin,timemax,100,0,50);
1597 fhTimePileUpMainVertexZDiamond->SetYTitle("diamond distance Z (cm) ");
1598 fhTimePileUpMainVertexZDiamond->SetXTitle("time (ns)");
1599 outputContainer->Add(fhTimePileUpMainVertexZDiamond);
1600
1601 }
1602
1a31a9ab 1603 return outputContainer ;
1604
1605}
1606
03bae431 1607//__________________________________
1608void AliAnaParticleIsolation::Init()
1609{
1610 // Do some checks and init stuff
1611
1612 // In case of several cone and thresholds analysis, open the cuts for the filling of the
1613 // track and cluster reference arrays in cone when done in the MakeAnalysisFillAOD().
1614 // The different cones, thresholds are tested for this list of tracks, clusters.
1615 if(fMakeSeveralIC)
1616 {
1617 printf("AliAnaParticleIsolation::Init() - Open default isolation cuts for multiple Isolation analysis\n");
1618 GetIsolationCut()->SetPtThreshold(100);
1619 GetIsolationCut()->SetPtFraction(100);
1620 GetIsolationCut()->SetConeSize(1);
1621 }
1622}
1623
803d06a8 1624//____________________________________________
1625void AliAnaParticleIsolation::InitParameters()
1626{
1627
1628 //Initialize the parameters of the analysis.
1629 SetInputAODName("PWG4Particle");
1630 SetAODObjArrayName("IsolationCone");
1631 AddToHistogramsName("AnaIsolation_");
1632
1633 fCalorimeter = "PHOS" ;
1634 fReMakeIC = kFALSE ;
1635 fMakeSeveralIC = kFALSE ;
1636
1637 //----------- Several IC-----------------
db6fb352 1638 fNCones = 5 ;
1639 fNPtThresFrac = 5 ;
1640 fConeSizes [0] = 0.1; fConeSizes [1] = 0.2; fConeSizes [2] = 0.3; fConeSizes [3] = 0.4; fConeSizes [4] = 0.5;
1641 fPtThresholds [0] = 1.; fPtThresholds [1] = 2.; fPtThresholds [2] = 3.; fPtThresholds [3] = 4.; fPtThresholds [4] = 5.;
1642 fPtFractions [0] = 0.05; fPtFractions [1] = 0.075; fPtFractions [2] = 0.1; fPtFractions [3] = 1.25; fPtFractions [4] = 1.5;
1643 fSumPtThresholds[0] = 1.; fSumPtThresholds[1] = 2.; fSumPtThresholds[2] = 3.; fSumPtThresholds[3] = 4.; fSumPtThresholds[4] = 5.;
803d06a8 1644
1645 //------------- Histograms settings -------
1646 fHistoNPtSumBins = 100 ;
1647 fHistoPtSumMax = 50 ;
1648 fHistoPtSumMin = 0. ;
1649
1650 fHistoNPtInConeBins = 100 ;
1651 fHistoPtInConeMax = 50 ;
1652 fHistoPtInConeMin = 0. ;
1653
1654}
1655
1656//__________________________________________________
1a31a9ab 1657void AliAnaParticleIsolation::MakeAnalysisFillAOD()
1658{
1659 //Do analysis and fill aods
1660 //Search for the isolated photon in fCalorimeter with pt > GetMinPt()
1661
b5dbb99b 1662 if(!GetInputAODBranch())
1663 {
1a31a9ab 1664 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
1665 abort();
1666 }
1667
b5dbb99b 1668 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation"))
1669 {
1a31a9ab 1670 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());
1671 abort();
1672 }
b0a31c92 1673
1a31a9ab 1674 Int_t n = 0, nfrac = 0;
1675 Bool_t isolated = kFALSE ;
1a31a9ab 1676 Float_t coneptsum = 0 ;
1677 TObjArray * pl = 0x0; ;
1678
1679 //Select the calorimeter for candidate isolation with neutral particles
b5dbb99b 1680 if (fCalorimeter == "PHOS" )
1a31a9ab 1681 pl = GetPHOSClusters();
1682 else if (fCalorimeter == "EMCAL")
1683 pl = GetEMCALClusters();
1684
1685 //Loop on AOD branch, filled previously in AliAnaPhoton, find leading particle to do isolation only with it
1686 Double_t ptLeading = 0. ;
1687 Int_t idLeading = -1 ;
1688 TLorentzVector mom ;
1689 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1690 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Input aod branch entries %d\n", naod);
1691
b5dbb99b 1692 for(Int_t iaod = 0; iaod < naod; iaod++)
1693 {
1a31a9ab 1694 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
81a5e27b 1695
1a31a9ab 1696 //If too small or too large pt, skip
1697 if(aodinput->Pt() < GetMinPt() || aodinput->Pt() > GetMaxPt() ) continue ;
1698
03bae431 1699 //check if it is low pt trigger particle
1700 if((aodinput->Pt() < GetIsolationCut()->GetPtThreshold() ||
1701 aodinput->Pt() < GetIsolationCut()->GetSumPtThreshold()) &&
1702 !fMakeSeveralIC)
b5dbb99b 1703 {
1a31a9ab 1704 continue ; //trigger should not come from underlying event
b5dbb99b 1705 }
1a31a9ab 1706
1707 //vertex cut in case of mixing
04f7a616 1708 Int_t check = CheckMixedEventVertex(aodinput->GetCaloLabel(0), aodinput->GetTrackLabel(0));
1709 if(check == 0) continue;
1710 if(check == -1) return;
1a31a9ab 1711
1712 //find the leading particles with highest momentum
547c2f01 1713 if ( aodinput->Pt() > ptLeading )
b5dbb99b 1714 {
1a31a9ab 1715 ptLeading = aodinput->Pt() ;
226b95ba 1716 idLeading = iaod ;
1a31a9ab 1717 }
b5dbb99b 1718
226b95ba 1719 aodinput->SetLeadingParticle(kFALSE);
b5dbb99b 1720
1a31a9ab 1721 }//finish searching for leading trigger particle
1722
1723 // Check isolation of leading particle
1724 if(idLeading < 0) return;
226b95ba 1725
1a31a9ab 1726 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(idLeading));
db6fb352 1727 aodinput->SetLeadingParticle(kTRUE);
547c2f01 1728
1729 // Check isolation only of clusters in fiducial region
1730 if(IsFiducialCutOn())
1731 {
050ad675 1732 Bool_t in = GetFiducialCut()->IsInFiducialCut(*aodinput->Momentum(),aodinput->GetDetector()) ;
547c2f01 1733 if(! in ) return ;
1734 }
1735
1a31a9ab 1736 //After cuts, study isolation
1737 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
ac5111f9 1738 GetIsolationCut()->MakeIsolationCut(GetCTSTracks(),pl,
1739 GetReader(), GetCaloPID(),
b5dbb99b 1740 kTRUE, aodinput, GetAODObjArrayName(),
1741 n,nfrac,coneptsum, isolated);
3d187b6c 1742
1743 if(!fMakeSeveralIC) aodinput->SetIsolated(isolated);
2ad19c3d 1744
b5dbb99b 1745 if(GetDebug() > 1)
1746 {
1a31a9ab 1747 if(isolated)printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() : Particle %d IS ISOLATED \n",idLeading);
1748 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - End fill AODs \n");
1749 }
1750
1751}
1752
803d06a8 1753//_________________________________________________________
1a31a9ab 1754void AliAnaParticleIsolation::MakeAnalysisFillHistograms()
1755{
1756 //Do analysis and fill histograms
db6fb352 1757
803d06a8 1758 Int_t n = 0, nfrac = 0;
1759 Bool_t isolated = kFALSE ;
1a31a9ab 1760 Float_t coneptsum = 0 ;
b7ce43b4 1761 Float_t etaUEptsum = 0 ;
1762 Float_t phiUEptsum = 0 ;
803d06a8 1763
1a31a9ab 1764 //Loop on stored AOD
1765 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1766 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Histo aod branch entries %d\n", naod);
1767
1768 //Get vertex for photon momentum calculation
3d187b6c 1769 Double_t vertex[] = {0,0,0} ; //vertex ;
b5dbb99b 1770 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1771 {
b0a31c92 1772 GetReader()->GetVertex(vertex);
1a31a9ab 1773 }
b0a31c92 1774
b5dbb99b 1775 for(Int_t iaod = 0; iaod < naod ; iaod++)
1776 {
1a31a9ab 1777 AliAODPWG4ParticleCorrelation* aod = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1778
1779 if(!aod->IsLeadingParticle()) continue; // Try to isolate only leading cluster or track
1780
547c2f01 1781 // Check isolation only of clusters in fiducial region
1782 if(IsFiducialCutOn())
1783 {
050ad675 1784 Bool_t in = GetFiducialCut()->IsInFiducialCut(*aod->Momentum(),aod->GetDetector()) ;
547c2f01 1785 if(! in ) continue ;
1786 }
1787
1a31a9ab 1788 Bool_t isolation = aod->IsIsolated();
803d06a8 1789 Bool_t decay = aod->IsTagged();
0fb69ade 1790 Float_t energy = aod->E();
1a31a9ab 1791 Float_t pt = aod->Pt();
1792 Float_t phi = aod->Phi();
1793 Float_t eta = aod->Eta();
1794 Float_t conesize = GetIsolationCut()->GetConeSize();
1795
1796 //Recover reference arrays with clusters and tracks
1797 TObjArray * refclusters = aod->GetObjArray(GetAODObjArrayName()+"Clusters");
1798 TObjArray * reftracks = aod->GetObjArray(GetAODObjArrayName()+"Tracks");
0fb69ade 1799
1a31a9ab 1800 //If too small or too large pt, skip
1801 if(pt < GetMinPt() || pt > GetMaxPt() ) continue ;
b0a31c92 1802
1a31a9ab 1803 // --- In case of redoing isolation from delta AOD ----
db6fb352 1804
03bae431 1805 if(fMakeSeveralIC)
1806 {
1a31a9ab 1807 //Analysis of multiple IC at same time
1808 MakeSeveralICAnalysis(aod);
03bae431 1809 continue;
1a31a9ab 1810 }
b5dbb99b 1811 else if(fReMakeIC)
1812 {
1a31a9ab 1813 //In case a more strict IC is needed in the produced AOD
1814 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
ac5111f9 1815 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
1816 GetReader(), GetCaloPID(),
b5dbb99b 1817 kFALSE, aod, "",
1818 n,nfrac,coneptsum, isolated);
1a31a9ab 1819 fhConeSumPt->Fill(pt,coneptsum);
1820 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Energy Sum in Isolation Cone %2.2f\n", coneptsum);
1821 }
db6fb352 1822
547c2f01 1823 // ---------------------------------------------------
1a31a9ab 1824
1825 //Fill pt distribution of particles in cone
1826 //Tracks
1827 coneptsum = 0;
1828 Double_t sumptFR = 0. ;
1829 TObjArray * trackList = GetCTSTracks() ;
b5dbb99b 1830 for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++)
1831 {
1a31a9ab 1832 AliVTrack* track = (AliVTrack *) trackList->At(itrack);
b7ce43b4 1833
3f150b4b 1834 if(!track)
1835 {
1a31a9ab 1836 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Track not available?");
1837 continue;
1838 }
0fb69ade 1839
b7ce43b4 1840 //fill histogram for UE in phi band
1841 if(track->Eta() > (eta-conesize) && track->Eta() < (eta+conesize))
6c80c1bf 1842 {
1843 phiUEptsum+=track->Pt();
1844 fhPhiBand->Fill(track->Eta(),track->Phi());
1845 }
b7ce43b4 1846 //fill histogram for UE in eta band in EMCal acceptance
1847 if(track->Phi() > (phi-conesize) && track->Phi() < (phi+conesize) && track->Eta() > -0.6 && track->Eta() < 0.6)
6c80c1bf 1848 {
1849 etaUEptsum+=track->Pt();
1850 fhEtaBand->Fill(track->Eta(),track->Phi());
1851 }
1852
b7ce43b4 1853 //fill the histograms at forward range
1a31a9ab 1854 Double_t dPhi = phi - track->Phi() + TMath::PiOver2();
1855 Double_t dEta = eta - track->Eta();
1856 Double_t arg = dPhi*dPhi + dEta*dEta;
b5dbb99b 1857 if(TMath::Sqrt(arg) < conesize)
1858 {
1a31a9ab 1859 fhPtInFRCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1860 sumptFR+=track->Pt();
0fb69ade 1861 }
1862
1a31a9ab 1863 dPhi = phi - track->Phi() - TMath::PiOver2();
1864 arg = dPhi*dPhi + dEta*dEta;
b5dbb99b 1865 if(TMath::Sqrt(arg) < conesize)
1866 {
1a31a9ab 1867 fhPtInFRCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1868 sumptFR+=track->Pt();
1869 }
1870 }
0fb69ade 1871
1a31a9ab 1872 fhFRConeSumPt->Fill(pt,sumptFR);
b5dbb99b 1873 if(reftracks)
1874 {
1875 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++)
1876 {
1a31a9ab 1877 AliVTrack* track = (AliVTrack *) reftracks->At(itrack);
6c80c1bf 1878 Float_t pTtrack = track->Pt();
1879 fhPtInCone->Fill(pt,pTtrack);
1880 if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD()) fhPtInConePileUp->Fill(pt,pTtrack);
1881 if (GetEventCentrality()) fhPtInConeCent->Fill(GetEventCentrality(),pTtrack);
1882 coneptsum+=pTtrack;
1a31a9ab 1883 }
1884 }
1885
1886 //CaloClusters
b5dbb99b 1887 if(refclusters)
1888 {
1a31a9ab 1889 TLorentzVector mom ;
b5dbb99b 1890 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++)
1891 {
1a31a9ab 1892 AliVCluster* calo = (AliVCluster *) refclusters->At(icalo);
1893 calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
6c80c1bf 1894
1895 //fill histogram for UE in phi band
1896 if(mom.Eta() > (eta-conesize) && mom.Eta() < (eta+conesize))
1897 {
1898 phiUEptsum+=mom.Pt();
1899 fhPhiBand->Fill(mom.Eta(),mom.Phi());
1900 }
1901 //fill histogram for UE in eta band in EMCal acceptance
1902 if(mom.Phi() > (phi-conesize) && mom.Phi() < (phi+conesize))
1903 {
1904 etaUEptsum+=mom.Pt();
1905 fhEtaBand->Fill(mom.Eta(),mom.Phi());
1906 }
1a31a9ab 1907
1908 fhPtInCone->Fill(pt, mom.Pt());
6c80c1bf 1909 if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD()) fhPtInConePileUp->Fill(pt,mom.Pt());
1910 if (GetEventCentrality()) fhPtInConeCent->Fill(GetEventCentrality(),mom.Pt());
1a31a9ab 1911 coneptsum+=mom.Pt();
1912 }
1913 }
6c80c1bf 1914
b7ce43b4 1915 //normalize phi/eta band per area unit
1916 fhPhiUEConeSumPt->Fill(pt, phiUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*TMath::Pi()));
1917 fhEtaUEConeSumPt->Fill(pt, etaUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*0.6));
6c80c1bf 1918
ec58c056 1919 Double_t sumPhiUESub = coneptsum-(phiUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*TMath::Pi()));
1920 Double_t sumEtaUESub = coneptsum-(etaUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*0.6));
6c80c1bf 1921
ec58c056 1922 fhConeSumPtPhiUESub->Fill(pt,sumPhiUESub);
1923 fhConeSumPtEtaUESub->Fill(pt,sumEtaUESub);
6c80c1bf 1924
1925
1a31a9ab 1926 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d Energy Sum in Isolation Cone %2.2f\n", iaod, coneptsum);
1927
1928 if(!fReMakeIC) fhConeSumPt->Fill(pt,coneptsum);
1929
b5dbb99b 1930 Int_t mcTag = aod->GetTag() ;
1931 Int_t clID = aod->GetCaloLabel(0) ;
8736d400 1932 Int_t nlm = aod->GetFiducialArea();
db6fb352 1933 if(GetDebug() > 0) printf(" AliAnaParticleIsolation::MakeAnalysisFillHistograms() - pt %1.1f, eta %1.1f, phi %1.1f\n",pt, eta, phi);
1934
8736d400 1935 FillTrackMatchingShowerShapeControlHistograms(isolation, clID,nlm,mcTag,reftracks,refclusters,aod,GetReader(), GetCaloPID());
2ad19c3d 1936
b5dbb99b 1937 if(isolation)
1938 {
db6fb352 1939 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d ISOLATED: fill histograms\n", iaod);
6c80c1bf 1940
2ad19c3d 1941 // Fill histograms to undertand pile-up before other cuts applied
1942 // Remember to relax time cuts in the reader
1943 FillPileUpHistograms(clID);
1944
b5dbb99b 1945 fhEIso ->Fill(energy);
1946 fhPtIso ->Fill(pt);
1947 fhPhiIso ->Fill(pt,phi);
1948 fhEtaIso ->Fill(pt,eta);
0fb69ade 1949 fhEtaPhiIso ->Fill(eta,phi);
8736d400 1950 fhPtNLocMaxIso->Fill(pt,nlm);
db6fb352 1951
1952 if(decay)
3f150b4b 1953 {
1954 fhPtDecayIso->Fill(pt);
1955 fhEtaPhiDecayIso->Fill(eta,phi);
1956 }
1a31a9ab 1957
b1f720a7 1958 if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD())
1959 {
1960 fhEIsoPileUp ->Fill(energy);
1961 fhPtIsoPileUp->Fill(pt);
1962 }
1963
b5dbb99b 1964 if(IsDataMC())
1965 {
1a31a9ab 1966
b5dbb99b 1967 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
803d06a8 1968 {
1969 fhPtIsoMCPhoton ->Fill(pt);
1970 }
1971
db6fb352 1972 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt))
1973 {
1a31a9ab 1974 fhPtIsoPrompt ->Fill(pt);
1975 fhPhiIsoPrompt ->Fill(pt,phi);
1976 fhEtaIsoPrompt ->Fill(pt,eta);
1977 }
b5dbb99b 1978 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation))
1a31a9ab 1979 {
1980 fhPtIsoFragmentation ->Fill(pt);
1981 fhPhiIsoFragmentation ->Fill(pt,phi);
1982 fhEtaIsoFragmentation ->Fill(pt,eta);
1983 }
764ab1f4 1984 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0))
1985 {
1986 fhPtIsoPi0 ->Fill(pt);
1987 fhPhiIsoPi0 ->Fill(pt,phi);
1988 fhEtaIsoPi0 ->Fill(pt,eta);
1989 }
b5dbb99b 1990 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay))
1a31a9ab 1991 {
1992 fhPtIsoPi0Decay ->Fill(pt);
1993 fhPhiIsoPi0Decay ->Fill(pt,phi);
1994 fhEtaIsoPi0Decay ->Fill(pt,eta);
1995 }
b5dbb99b 1996 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay))
803d06a8 1997 {
1998 fhPtIsoEtaDecay ->Fill(pt);
1999 fhPhiIsoEtaDecay ->Fill(pt,phi);
2000 fhEtaIsoEtaDecay ->Fill(pt,eta);
2001 }
b5dbb99b 2002 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay))
1a31a9ab 2003 {
2004 fhPtIsoOtherDecay ->Fill(pt);
2005 fhPhiIsoOtherDecay ->Fill(pt,phi);
2006 fhEtaIsoOtherDecay ->Fill(pt,eta);
2007 }
6c80c1bf 2008 // else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion))
2009 // {
2010 // fhPtIsoConversion ->Fill(pt);
2011 // fhPhiIsoConversion ->Fill(pt,phi);
2012 // fhEtaIsoConversion ->Fill(pt,eta);
2013 // }
764ab1f4 2014 else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron))// anything else but electrons
1a31a9ab 2015 {
764ab1f4 2016 fhPtIsoHadron ->Fill(pt);
2017 fhPhiIsoHadron ->Fill(pt,phi);
2018 fhEtaIsoHadron ->Fill(pt,eta);
1a31a9ab 2019 }
2020 }//Histograms with MC
2021
2022 }//Isolated histograms
ca134929 2023 else // NON isolated
1a31a9ab 2024 {
db6fb352 2025 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d NOT ISOLATED, fill histograms\n", iaod);
2026
b1f720a7 2027 fhENoIso ->Fill(energy);
2028 fhPtNoIso ->Fill(pt);
2029 fhEtaPhiNoIso ->Fill(eta,phi);
8736d400 2030 fhPtNLocMaxNoIso->Fill(pt,nlm);
2031
b1f720a7 2032 if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD())
2033 {
2034 fhENoIsoPileUp ->Fill(energy);
2035 fhPtNoIsoPileUp->Fill(pt);
2036 }
2037
db6fb352 2038 if(decay)
3f150b4b 2039 {
db6fb352 2040 fhPtDecayNoIso ->Fill(pt);
3f150b4b 2041 fhEtaPhiDecayNoIso->Fill(eta,phi);
2042 }
1a31a9ab 2043
b5dbb99b 2044 if(IsDataMC())
2045 {
db6fb352 2046 if (GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton)) fhPtNoIsoMCPhoton ->Fill(pt);
764ab1f4 2047 if (GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0)) fhPtNoIsoPi0 ->Fill(pt);
2048 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtNoIsoPi0Decay ->Fill(pt);
db6fb352 2049 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtNoIsoEtaDecay ->Fill(pt);
2050 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtNoIsoOtherDecay ->Fill(pt);
2051 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt)) fhPtNoIsoPrompt ->Fill(pt);
2052 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation)) fhPtNoIsoFragmentation->Fill(pt);
6c80c1bf 2053 // else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion)) fhPtNoIsoConversion ->Fill(pt);
764ab1f4 2054 else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron)) fhPtNoIsoHadron ->Fill(pt);
1a31a9ab 2055 }
2056 }
1a31a9ab 2057 }// aod loop
2058
2059}
2060
1a31a9ab 2061
803d06a8 2062//_____________________________________________________________________________________
1a31a9ab 2063void AliAnaParticleIsolation::MakeSeveralICAnalysis(AliAODPWG4ParticleCorrelation* ph)
2064{
1a31a9ab 2065
db6fb352 2066 //Isolation Cut Analysis for both methods and different pt cuts and cones
2067 Float_t ptC = ph->Pt();
2068 Float_t etaC = ph->Eta();
2069 Float_t phiC = ph->Phi();
2070 Int_t tag = ph->GetTag();
2071 Bool_t decay = ph->IsTagged();
b0a31c92 2072
03bae431 2073 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeSeveralICAnalysis() - Isolate pT %2.2f\n",ptC);
2074
1a31a9ab 2075 //Keep original setting used when filling AODs, reset at end of analysis
2076 Float_t ptthresorg = GetIsolationCut()->GetPtThreshold();
2077 Float_t ptfracorg = GetIsolationCut()->GetPtFraction();
2078 Float_t rorg = GetIsolationCut()->GetConeSize();
2079
b7ce43b4 2080 Float_t coneptsum = 0 ;
db6fb352 2081 Int_t n [10][10];//[fNCones][fNPtThresFrac];
2082 Int_t nfrac[10][10];//[fNCones][fNPtThresFrac];
ca134929 2083 Bool_t isolated = kFALSE;
2084 Int_t nCone = 0;
2085 Int_t nFracCone = 0;
44e48e82 2086
db6fb352 2087 // fill hist with all particles before isolation criteria
2088 fhPtNoIso ->Fill(ptC);
2089 fhEtaPhiNoIso->Fill(etaC,phiC);
2090
2091 if(IsDataMC())
2092 {
2093 if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)) fhPtNoIsoMCPhoton ->Fill(ptC);
764ab1f4 2094 if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtNoIsoPi0 ->Fill(ptC);
2095 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtNoIsoPi0Decay ->Fill(ptC);
db6fb352 2096 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtNoIsoEtaDecay ->Fill(ptC);
2097 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtNoIsoOtherDecay ->Fill(ptC);
2098 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtNoIsoPrompt ->Fill(ptC);
2099 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtNoIsoFragmentation->Fill(ptC);
764ab1f4 2100// else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtNoIsoConversion ->Fill(ptC);
2101 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtNoIsoHadron ->Fill(ptC);
db6fb352 2102 }
2103
2104 if(decay)
2105 {
2106 fhPtDecayNoIso ->Fill(ptC);
2107 fhEtaPhiDecayNoIso->Fill(etaC,phiC);
2108 }
44e48e82 2109 //Get vertex for photon momentum calculation
2110 Double_t vertex[] = {0,0,0} ; //vertex ;
2111 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
2112 {
2113 GetReader()->GetVertex(vertex);
2114 }
2115
1a31a9ab 2116 //Loop on cone sizes
b5dbb99b 2117 for(Int_t icone = 0; icone<fNCones; icone++)
2118 {
44e48e82 2119 //Recover reference arrays with clusters and tracks
2120 TObjArray * refclusters = ph->GetObjArray(GetAODObjArrayName()+"Clusters");
2121 TObjArray * reftracks = ph->GetObjArray(GetAODObjArrayName()+"Tracks");
2122
2123 //If too small or too large pt, skip
2124 if(ptC < GetMinPt() || ptC > GetMaxPt() ) continue ;
2125
2126 //In case a more strict IC is needed in the produced AOD
2127
b7ce43b4 2128 nCone=0; nFracCone = 0; isolated = kFALSE; coneptsum = 0;
44e48e82 2129
2130 GetIsolationCut()->SetSumPtThreshold(100);
2131 GetIsolationCut()->SetPtThreshold(100);
2132 GetIsolationCut()->SetPtFraction(100);
2133 GetIsolationCut()->SetConeSize(fConeSizes[icone]);
2134 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
2135 GetReader(), GetCaloPID(),
2136 kFALSE, ph, "",
b7ce43b4 2137 nCone,nFracCone,coneptsum, isolated);
44e48e82 2138
2139
2140 fhSumPtLeadingPt[icone]->Fill(ptC,coneptsum);
db6fb352 2141
44e48e82 2142 // retreive pt tracks to fill histo vs. pt leading
2143 //Fill pt distribution of particles in cone
2144 //fhPtLeadingPt(),fhFRSumPtLeadingPt(),fhFRPtLeadingPt(),
2145
2146 //Tracks
2147 coneptsum = 0;
2148 Double_t sumptFR = 0. ;
2149 TObjArray * trackList = GetCTSTracks() ;
2150 for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++)
2151 {
2152 AliVTrack* track = (AliVTrack *) trackList->At(itrack);
2153 //fill the histograms at forward range
2154 if(!track)
2155 {
2156 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Track not available?");
2157 continue;
2158 }
2159
2160 Double_t dPhi = phiC - track->Phi() + TMath::PiOver2();
2161 Double_t dEta = etaC - track->Eta();
2162 Double_t arg = dPhi*dPhi + dEta*dEta;
2163 if(TMath::Sqrt(arg) < fConeSizes[icone])
2164 {
2165 fhFRPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
2166 sumptFR+=track->Pt();
2167 }
2168
2169 dPhi = phiC - track->Phi() - TMath::PiOver2();
2170 arg = dPhi*dPhi + dEta*dEta;
2171 if(TMath::Sqrt(arg) < fConeSizes[icone])
2172 {
2173 fhFRPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
2174 sumptFR+=track->Pt();
2175 }
2176 }
2177 fhFRSumPtLeadingPt[icone]->Fill(ptC,sumptFR);
2178 if(reftracks)
2179 {
2180 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++)
2181 {
2182 AliVTrack* track = (AliVTrack *) reftracks->At(itrack);
2183 fhPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
2184 coneptsum+=track->Pt();
2185 }
2186 }
2187 //CaloClusters
2188 if(refclusters)
2189 {
2190 TLorentzVector mom ;
2191 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++)
2192 {
2193 AliVCluster* calo = (AliVCluster *) refclusters->At(icalo);
2194 calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
2195
2196 fhPtLeadingPt[icone]->Fill(ptC, mom.Pt());
2197 coneptsum+=mom.Pt();
2198 }
2199 }
2200 ///////////////////
2201
2202
1a31a9ab 2203 //Loop on ptthresholds
b5dbb99b 2204 for(Int_t ipt = 0; ipt<fNPtThresFrac ;ipt++)
2205 {
db6fb352 2206 n [icone][ipt]=0;
1a31a9ab 2207 nfrac[icone][ipt]=0;
2208 GetIsolationCut()->SetPtThreshold(fPtThresholds[ipt]);
db6fb352 2209 GetIsolationCut()->SetPtFraction(fPtFractions[ipt]) ;
2210 GetIsolationCut()->SetSumPtThreshold(fSumPtThresholds[ipt]);
2211
44e48e82 2212 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
ac5111f9 2213 GetReader(), GetCaloPID(),
b5dbb99b 2214 kFALSE, ph, "",
2215 n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);
1a31a9ab 2216
db6fb352 2217 if(!isolated) continue;
1a31a9ab 2218 //Normal ptThreshold cut
db6fb352 2219
2220 if(GetDebug() > 0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - cone size %1.1f, ptThres %1.1f, sumptThresh %1.1f, n %d, nfrac %d, coneptsum %2.2f, isolated %d\n",
2221 fConeSizes[icone],fPtThresholds[ipt],fSumPtThresholds[ipt],n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);
2222 if(GetDebug() > 0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - pt %1.1f, eta %1.1f, phi %1.1f\n",ptC, etaC, phiC);
2223
b5dbb99b 2224 if(n[icone][ipt] == 0)
2225 {
db6fb352 2226 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling pt threshold loop\n");
1a31a9ab 2227 fhPtThresIsolated[icone][ipt]->Fill(ptC);
db6fb352 2228 fhEtaPhiPtThresIso[icone][ipt]->Fill(etaC,phiC);
2229
2230 if(decay)
2231 {
2232 fhPtPtThresDecayIso[icone][ipt]->Fill(ptC);
2233 // fhEtaPhiPtThresDecayIso[icone][ipt]->Fill(etaC,phiC);
2234 }
2235
b5dbb99b 2236 if(IsDataMC())
2237 {
803d06a8 2238 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtThresIsolatedPrompt[icone][ipt] ->Fill(ptC) ;
764ab1f4 2239// else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtThresIsolatedConversion[icone][ipt] ->Fill(ptC) ;
803d06a8 2240 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtThresIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
764ab1f4 2241 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtThresIsolatedPi0[icone][ipt] ->Fill(ptC) ;
803d06a8 2242 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtThresIsolatedPi0Decay[icone][ipt] ->Fill(ptC) ;
2243 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtThresIsolatedEtaDecay[icone][ipt] ->Fill(ptC) ;
2244 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtThresIsolatedOtherDecay[icone][ipt] ->Fill(ptC) ;
764ab1f4 2245 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtThresIsolatedHadron[icone][ipt] ->Fill(ptC) ;
1a31a9ab 2246 }
2247 }
2248
db6fb352 2249 // pt in cone fraction
b5dbb99b 2250 if(nfrac[icone][ipt] == 0)
2251 {
db6fb352 2252 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling frac loop\n");
1a31a9ab 2253 fhPtFracIsolated[icone][ipt]->Fill(ptC);
db6fb352 2254 fhEtaPhiPtFracIso[icone][ipt]->Fill(etaC,phiC);
2255
2256 if(decay)
2257 {
2258 fhPtPtFracDecayIso[icone][ipt]->Fill(ptC);
b0a31c92 2259 fhEtaPhiPtFracDecayIso[icone][ipt]->Fill(etaC,phiC);
db6fb352 2260 }
2261
b5dbb99b 2262 if(IsDataMC())
2263 {
803d06a8 2264 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtFracIsolatedPrompt[icone][ipt] ->Fill(ptC) ;
764ab1f4 2265// else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtFracIsolatedConversion[icone][ipt] ->Fill(ptC) ;
803d06a8 2266 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtFracIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
764ab1f4 2267 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtFracIsolatedPi0[icone][ipt] ->Fill(ptC) ;
803d06a8 2268 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtFracIsolatedPi0Decay[icone][ipt] ->Fill(ptC) ;
2269 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtFracIsolatedEtaDecay[icone][ipt] ->Fill(ptC) ;
2270 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtFracIsolatedOtherDecay[icone][ipt] ->Fill(ptC) ;
764ab1f4 2271 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtFracIsolatedHadron[icone][ipt]->Fill(ptC) ;
1a31a9ab 2272 }
2273 }
db6fb352 2274
2275 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - checking IC method : %i\n",GetIsolationCut()->GetICMethod());
2276
2277 //Pt threshold on pt cand/ sum in cone histograms
2278 if(coneptsum<fSumPtThresholds[ipt])
2279 {// if((GetIsolationCut()->GetICMethod())==1){//kSumPtIC){
2280 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling sum loop\n");
2281 fhPtSumIsolated[icone][ipt]->Fill(ptC) ;
b0a31c92 2282 fhEtaPhiPtSumIso[icone][ipt]->Fill(etaC, phiC) ;
db6fb352 2283 if(decay)
2284 {
2285 fhPtPtSumDecayIso[icone][ipt]->Fill(ptC);
b0a31c92 2286 fhEtaPhiPtSumDecayIso[icone][ipt]->Fill(etaC, phiC) ;
db6fb352 2287 }
2288 }
2289
44e48e82 2290 // pt sum pt frac method
bb5fc123 2291// if( ((fPtFractions[ipt]*ptC < fSumPtThresholds[ipt]) && (coneptsum < fSumPtThresholds[ipt])) || ((fPtFractions[ipt]*ptC > fSumPtThresholds[ipt]) && (coneptsum < fPtFractions[ipt]*ptC)) )
2292
2293 if(coneptsum < fPtFractions[ipt]*ptC)
2294 {
b0a31c92 2295 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling PtFrac PtSum loop\n");
2296 fhPtFracPtSumIso[icone][ipt]->Fill(ptC) ;
2297 fhEtaPhiFracPtSumIso[icone][ipt]->Fill(etaC,phiC) ;
2298
2299 if(decay)
2300 {
2301 fhPtFracPtSumDecayIso[icone][ipt]->Fill(ptC);
2302 fhEtaPhiFracPtSumDecayIso[icone][ipt]->Fill(etaC,phiC);
2303 }
2304 }
2305
2306 // density method
db6fb352 2307 Float_t cellDensity = GetIsolationCut()->GetCellDensity( ph, GetReader());
2308 if(coneptsum<fSumPtThresholds[ipt]*cellDensity)
2309 {//(GetIsolationCut()->GetICMethod())==4){//kSumDensityIC) {
2310 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling density loop\n");
2311 fhPtSumDensityIso[icone][ipt]->Fill(ptC) ;
b0a31c92 2312 fhEtaPhiSumDensityIso[icone][ipt]->Fill(etaC,phiC) ;
2313
db6fb352 2314 if(decay)
2315 {
2316 fhPtSumDensityDecayIso[icone][ipt]->Fill(ptC);
b0a31c92 2317 fhEtaPhiSumDensityDecayIso[icone][ipt]->Fill(etaC, phiC);
db6fb352 2318 }
2319
2320 }
1a31a9ab 2321 }//pt thresh loop
2322
b5dbb99b 2323 if(IsDataMC())
2324 {
803d06a8 2325 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtSumIsolatedPrompt[icone] ->Fill(ptC,coneptsum) ;
764ab1f4 2326// else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtSumIsolatedConversion[icone] ->Fill(ptC,coneptsum) ;
803d06a8 2327 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtSumIsolatedFragmentation[icone]->Fill(ptC,coneptsum) ;
764ab1f4 2328 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtSumIsolatedPi0[icone] ->Fill(ptC,coneptsum) ;
803d06a8 2329 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtSumIsolatedPi0Decay[icone] ->Fill(ptC,coneptsum) ;
2330 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtSumIsolatedEtaDecay[icone] ->Fill(ptC,coneptsum) ;
2331 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtSumIsolatedOtherDecay[icone] ->Fill(ptC,coneptsum) ;
764ab1f4 2332 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtSumIsolatedHadron[icone]->Fill(ptC,coneptsum) ;
1a31a9ab 2333 }
2334
2335 }//cone size loop
2336
2337 //Reset original parameters for AOD analysis
2338 GetIsolationCut()->SetPtThreshold(ptthresorg);
2339 GetIsolationCut()->SetPtFraction(ptfracorg);
2340 GetIsolationCut()->SetConeSize(rorg);
2341
2342}
2343
803d06a8 2344//_____________________________________________________________
1a31a9ab 2345void AliAnaParticleIsolation::Print(const Option_t * opt) const
2346{
2347
2348 //Print some relevant parameters set for the analysis
2349 if(! opt)
2350 return;
2351
2352 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
745913ae 2353 AliAnaCaloTrackCorrBaseClass::Print(" ");
1a31a9ab 2354
2355 printf("ReMake Isolation = %d \n", fReMakeIC) ;
2356 printf("Make Several Isolation = %d \n", fMakeSeveralIC) ;
2357 printf("Calorimeter for isolation = %s \n", fCalorimeter.Data()) ;
2358
b5dbb99b 2359 if(fMakeSeveralIC)
2360 {
1a31a9ab 2361 printf("N Cone Sizes = %d\n", fNCones) ;
2362 printf("Cone Sizes = \n") ;
2363 for(Int_t i = 0; i < fNCones; i++)
2364 printf(" %1.2f;", fConeSizes[i]) ;
2365 printf(" \n") ;
2366
2367 printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ;
2368 printf(" pT thresholds = \n") ;
2369 for(Int_t i = 0; i < fNPtThresFrac; i++)
2370 printf(" %2.2f;", fPtThresholds[i]) ;
2371
2372 printf(" \n") ;
2373
2374 printf(" pT fractions = \n") ;
2375 for(Int_t i = 0; i < fNPtThresFrac; i++)
2376 printf(" %2.2f;", fPtFractions[i]) ;
2377
db6fb352 2378 printf(" \n") ;
2379
2380 printf("sum pT thresholds = \n") ;
2381 for(Int_t i = 0; i < fNPtThresFrac; i++)
2382 printf(" %2.2f;", fSumPtThresholds[i]) ;
2383
2384
1a31a9ab 2385 }
2386
b5dbb99b 2387 printf("Histograms: %3.1f < pT sum < %3.1f, Nbin = %d\n", fHistoPtSumMin, fHistoPtSumMax, fHistoNPtSumBins );
1a31a9ab 2388 printf("Histograms: %3.1f < pT in cone < %3.1f, Nbin = %d\n", fHistoPtInConeMin, fHistoPtInConeMax, fHistoNPtInConeBins);
2389
2390 printf(" \n") ;
2391
2392}
2393