1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //_________________________________________________________________________
17 // Class for analysis of particle isolation
18 // Input is selected particles put in AOD branch (AliAODPWG4ParticleCorrelation)
20 // Class created from old AliPHOSGammaJet
21 // (see AliRoot versions previous Release 4-09)
23 // -- Author: Gustavo Conesa (LNF-INFN)
25 //-Yaxian Mao (add the possibility for different IC method with different pt range, 01/10/2010)
26 //////////////////////////////////////////////////////////////////////////////
29 // --- ROOT system ---
30 #include <TClonesArray.h>
32 #include <TObjString.h>
36 #include "TParticle.h"
37 #include "TDatabasePDG.h"
41 // --- Analysis system ---
42 #include "AliAnaParticleIsolation.h"
43 #include "AliCaloTrackReader.h"
45 #include "AliIsolationCut.h"
46 #include "AliFiducialCut.h"
47 #include "AliMCAnalysisUtils.h"
48 #include "AliNeutralMesonSelection.h"
49 #include "AliAODMCParticle.h"
50 #include "AliAODPWG4ParticleCorrelation.h"
51 #include "AliMCAnalysisUtils.h"
52 #include "AliVTrack.h"
53 #include "AliVCluster.h"
54 #include "AliESDEvent.h"
55 #include "AliAODEvent.h"
57 #include "AliEMCALGeometry.h"
58 #include "AliPHOSGeoUtils.h"
60 ClassImp(AliAnaParticleIsolation)
62 //______________________________________________________________________________
63 AliAnaParticleIsolation::AliAnaParticleIsolation() :
64 AliAnaCaloTrackCorrBaseClass(),
65 fCalorimeter(""), fIsoDetector(""),
66 fReMakeIC(0), fMakeSeveralIC(0),
67 fFillPileUpHistograms(0),
68 fFillTMHisto(0), fFillSSHisto(1),
69 fFillUEBandSubtractHistograms(1), fFillCellHistograms(0),
70 fFillHighMultHistograms(0), fFillTaggedDecayHistograms(0),
71 fFillNLMHistograms(0),
73 fNCones(0), fNPtThresFrac(0),
74 fConeSizes(), fPtThresholds(),
75 fPtFractions(), fSumPtThresholds(),
77 fhEIso(0), fhPtIso(0),
78 fhPtCentralityIso(0), fhPtEventPlaneIso(0),
80 fhPhiIso(0), fhEtaIso(0), fhEtaPhiIso(0),
82 fhENoIso(0), fhPtNoIso(0), fhPtNLocMaxNoIso(0),
83 fhPtDecayIso(0), fhPtDecayNoIso(0),
84 fhEtaPhiDecayIso(0), fhEtaPhiDecayNoIso(0),
86 fhPtClusterInCone(0), fhPtCellInCone(0), fhPtTrackInCone(0),
87 fhPtTrackInConeOtherBC(0), fhPtTrackInConeOtherBCPileUpSPD(0),
88 fhPtTrackInConeBC0(0), fhPtTrackInConeVtxBC0(0),
89 fhPtTrackInConeBC0PileUpSPD(0),
90 fhPtInConePileUp(), fhPtInConeCent(0),
91 fhPerpConeSumPt(0), fhPtInPerpCone(0),
92 fhEtaPhiInConeCluster(0), fhEtaPhiCluster(0),
93 fhEtaPhiInConeTrack(0), fhEtaPhiTrack(0),
94 fhEtaBandCluster(0), fhPhiBandCluster(0),
95 fhEtaBandTrack(0), fhPhiBandTrack(0),
96 fhEtaBandCell(0), fhPhiBandCell(0),
97 fhConeSumPt(0), fhConeSumPtCellTrack(0),
98 fhConeSumPtCell(0), fhConeSumPtCluster(0), fhConeSumPtTrack(0),
99 fhConeSumPtEtaBandUECluster(0), fhConeSumPtPhiBandUECluster(0),
100 fhConeSumPtEtaBandUETrack(0), fhConeSumPtPhiBandUETrack(0),
101 fhConeSumPtEtaBandUECell(0), fhConeSumPtPhiBandUECell(0),
102 fhConeSumPtTrigEtaPhi(0),
103 fhConeSumPtCellTrackTrigEtaPhi(0),
104 fhConeSumPtEtaBandUEClusterTrigEtaPhi(0), fhConeSumPtPhiBandUEClusterTrigEtaPhi(0),
105 fhConeSumPtEtaBandUETrackTrigEtaPhi(0), fhConeSumPtPhiBandUETrackTrigEtaPhi(0),
106 fhConeSumPtEtaBandUECellTrigEtaPhi(0), fhConeSumPtPhiBandUECellTrigEtaPhi(0),
107 fhConeSumPtEtaUESub(0), fhConeSumPtPhiUESub(0),
108 fhConeSumPtEtaUESubTrigEtaPhi(0), fhConeSumPtPhiUESubTrigEtaPhi(0),
109 fhConeSumPtEtaUESubTrackCell(0), fhConeSumPtPhiUESubTrackCell(0),
110 fhConeSumPtEtaUESubTrackCellTrigEtaPhi(0), fhConeSumPtPhiUESubTrackCellTrigEtaPhi(0),
111 fhConeSumPtEtaUESubCluster(0), fhConeSumPtPhiUESubCluster(0),
112 fhConeSumPtEtaUESubClusterTrigEtaPhi(0), fhConeSumPtPhiUESubClusterTrigEtaPhi(0),
113 fhConeSumPtEtaUESubCell(0), fhConeSumPtPhiUESubCell(0),
114 fhConeSumPtEtaUESubCellTrigEtaPhi(0), fhConeSumPtPhiUESubCellTrigEtaPhi(0),
115 fhConeSumPtEtaUESubTrack(0), fhConeSumPtPhiUESubTrack(0),
116 fhConeSumPtEtaUESubTrackTrigEtaPhi(0), fhConeSumPtPhiUESubTrackTrigEtaPhi(0),
117 fhFractionTrackOutConeEta(0), fhFractionTrackOutConeEtaTrigEtaPhi(0),
118 fhFractionClusterOutConeEta(0), fhFractionClusterOutConeEtaTrigEtaPhi(0),
119 fhFractionClusterOutConePhi(0), fhFractionClusterOutConePhiTrigEtaPhi(0),
120 fhFractionCellOutConeEta(0), fhFractionCellOutConeEtaTrigEtaPhi(0),
121 fhFractionCellOutConePhi(0), fhFractionCellOutConePhiTrigEtaPhi(0),
122 fhConeSumPtClustervsTrack(0),
123 fhConeSumPtEtaUESubClustervsTrack(0), fhConeSumPtPhiUESubClustervsTrack(0),
124 fhConeSumPtCellvsTrack(0),
125 fhConeSumPtEtaUESubCellvsTrack(0), fhConeSumPtPhiUESubCellvsTrack(0),
126 fhEtaBandClustervsTrack(0), fhPhiBandClustervsTrack(0),
127 fhEtaBandNormClustervsTrack(0), fhPhiBandNormClustervsTrack(0),
128 fhEtaBandCellvsTrack(0), fhPhiBandCellvsTrack(0),
129 fhEtaBandNormCellvsTrack(0), fhPhiBandNormCellvsTrack(0),
130 fhConeSumPtSubvsConeSumPtTotPhiTrack(0), fhConeSumPtSubNormvsConeSumPtTotPhiTrack(0),
131 fhConeSumPtSubvsConeSumPtTotEtaTrack(0), fhConeSumPtSubNormvsConeSumPtTotEtaTrack(0),
132 fhConeSumPtSubvsConeSumPtTotPhiCluster(0), fhConeSumPtSubNormvsConeSumPtTotPhiCluster(0),
133 fhConeSumPtSubvsConeSumPtTotEtaCluster(0), fhConeSumPtSubNormvsConeSumPtTotEtaCluster(0),
134 fhConeSumPtSubvsConeSumPtTotPhiCell(0), fhConeSumPtSubNormvsConeSumPtTotPhiCell(0),
135 fhConeSumPtSubvsConeSumPtTotEtaCell(0), fhConeSumPtSubNormvsConeSumPtTotEtaCell(0),
136 fhConeSumPtVSUETracksEtaBand(0), fhConeSumPtVSUETracksPhiBand(0),
137 fhConeSumPtVSUEClusterEtaBand(0), fhConeSumPtVSUEClusterPhiBand(0),
139 // Number of local maxima in cluster
141 fhELambda0LocMax1(), fhELambda1LocMax1(),
142 fhELambda0LocMax2(), fhELambda1LocMax2(),
143 fhELambda0LocMaxN(), fhELambda1LocMaxN(),
145 fhEIsoPileUp(), fhPtIsoPileUp(),
146 fhENoIsoPileUp(), fhPtNoIsoPileUp(),
147 fhTimeENoCut(0), fhTimeESPD(0), fhTimeESPDMulti(0),
148 fhTimeNPileUpVertSPD(0), fhTimeNPileUpVertTrack(0),
149 fhTimeNPileUpVertContributors(0),
150 fhTimePileUpMainVertexZDistance(0), fhTimePileUpMainVertexZDiamond(0)
154 //Initialize parameters
157 for(Int_t i = 0; i < 5 ; i++)
161 for(Int_t imc = 0; imc < 9; imc++)
162 fhSumPtLeadingPtMC[imc][i] = 0 ;
164 for(Int_t j = 0; j < 5 ; j++)
166 fhPtThresIsolated [i][j] = 0 ;
167 fhPtFracIsolated [i][j] = 0 ;
168 fhSumPtIsolated [i][j] = 0 ;
170 fhEtaPhiPtThresIso [i][j] = 0 ;
171 fhEtaPhiPtThresDecayIso [i][j] = 0 ;
172 fhPtPtThresDecayIso [i][j] = 0 ;
174 fhEtaPhiPtFracIso [i][j] = 0 ;
175 fhEtaPhiPtFracDecayIso [i][j] = 0 ;
176 fhPtPtFracDecayIso [i][j] = 0 ;
177 fhPtPtSumDecayIso [i][j] = 0 ;
178 fhPtSumDensityIso [i][j] = 0 ;
179 fhPtSumDensityDecayIso [i][j] = 0 ;
180 fhEtaPhiSumDensityIso [i][j] = 0 ;
181 fhEtaPhiSumDensityDecayIso [i][j] = 0 ;
182 fhPtFracPtSumIso [i][j] = 0 ;
183 fhPtFracPtSumDecayIso [i][j] = 0 ;
184 fhEtaPhiFracPtSumIso [i][j] = 0 ;
185 fhEtaPhiFracPtSumDecayIso [i][j] = 0 ;
187 for(Int_t imc = 0; imc < 9; imc++)
189 fhPtThresIsolatedMC[imc][i][j] = 0 ;
190 fhPtFracIsolatedMC [imc][i][j] = 0 ;
191 fhSumPtIsolatedMC [imc][i][j] = 0 ;
197 for(Int_t i = 0; i < 5 ; i++)
199 fPtFractions [i] = 0 ;
200 fPtThresholds [i] = 0 ;
201 fSumPtThresholds[i] = 0 ;
203 fhSumPtLeadingPt [i] = 0 ;
204 fhPtLeadingPt [i] = 0 ;
205 fhPerpSumPtLeadingPt[i] = 0 ;
206 fhPerpPtLeadingPt [i] = 0 ;
209 for(Int_t imc = 0; imc < 9; imc++)
211 fhPtNoIsoMC [imc] = 0;
213 fhPhiIsoMC [imc] = 0;
214 fhEtaIsoMC [imc] = 0;
215 fhPtLambda0MC[imc][0] = 0;
216 fhPtLambda0MC[imc][1] = 0;
218 fhPtDecayIsoMC [imc] = 0;
219 fhPtDecayNoIsoMC[imc] = 0;
222 for(Int_t i = 0; i < 2 ; i++)
224 fhTrackMatchedDEta[i] = 0 ; fhTrackMatchedDPhi[i] = 0 ; fhTrackMatchedDEtaDPhi [i] = 0 ;
225 fhdEdx [i] = 0 ; fhEOverP [i] = 0 ; fhTrackMatchedMCParticle[i] = 0 ;
226 fhELambda0 [i] = 0 ; fhELambda1 [i] = 0 ; fhPtLambda0 [i] = 0 ;
227 fhELambda0TRD [i] = 0 ; fhELambda1TRD [i] = 0 ; fhPtLambda0TRD [i] = 0 ;
229 // Number of local maxima in cluster
231 fhELambda0LocMax1[i] = 0 ; fhELambda1LocMax1[i] = 0 ;
232 fhELambda0LocMax2[i] = 0 ; fhELambda1LocMax2[i] = 0 ;
233 fhELambda0LocMaxN[i] = 0 ; fhELambda1LocMaxN[i] = 0 ;
237 for(Int_t i = 0; i < 6; i++)
239 fhPtPrimMCiso[i] = 0;
247 for(Int_t i = 0 ; i < 7 ; i++)
249 fhPtInConePileUp[i] = 0 ;
250 fhEIsoPileUp [i] = 0 ;
251 fhPtIsoPileUp [i] = 0 ;
252 fhENoIsoPileUp [i] = 0 ;
253 fhPtNoIsoPileUp [i] = 0 ;
258 //_______________________________________________________________________________________________
259 void AliAnaParticleIsolation::CalculateCaloUEBand(AliAODPWG4ParticleCorrelation * pCandidate,
260 Float_t & etaBandPtSum, Float_t & phiBandPtSum)
262 // Get the clusters pT or sum of pT in phi/eta bands or at 45 degrees from trigger
264 if( GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kOnlyCharged ) return ;
266 Float_t conesize = GetIsolationCut()->GetConeSize();
269 //Select the Calorimeter
270 TObjArray * pl = 0x0;
271 if (fCalorimeter == "PHOS" )
272 pl = GetPHOSClusters();
273 else if (fCalorimeter == "EMCAL")
274 pl = GetEMCALClusters();
278 //Get vertex for cluster momentum calculation
279 Double_t vertex[] = {0,0,0} ; //vertex ;
280 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
281 GetReader()->GetVertex(vertex);
283 Float_t ptTrig = pCandidate->Pt() ;
284 Float_t phiTrig = pCandidate->Phi();
285 Float_t etaTrig = pCandidate->Eta();
287 for(Int_t icluster=0; icluster < pl->GetEntriesFast(); icluster++)
289 AliVCluster* cluster = (AliVCluster *) pl->At(icluster);
293 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Cluster not available?");
297 //Do not count the candidate (photon or pi0) or the daughters of the candidate
298 if(cluster->GetID() == pCandidate->GetCaloLabel(0) ||
299 cluster->GetID() == pCandidate->GetCaloLabel(1) ) continue ;
301 //Remove matched clusters to tracks if Neutral and Track info is used
302 if( GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kNeutralAndCharged &&
303 IsTrackMatched(cluster,GetReader()->GetInputEvent())) continue ;
305 cluster->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
307 //exclude particles in cone
308 Float_t rad = GetIsolationCut()->Radius(etaTrig, phiTrig, mom.Eta(), mom.Phi());
310 // histo of eta and phi for all clusters
311 fhEtaPhiCluster->Fill(mom.Eta(), mom.Phi());
313 // histos for all clusters in cone
314 fhEtaPhiInConeCluster->Fill(mom.Eta(), mom.Phi());
317 //fill histogram for UE in phi band in EMCal acceptance
318 if(mom.Eta() > (etaTrig-conesize) && mom.Eta() < (etaTrig+conesize))
320 phiBandPtSum+=mom.Pt();
321 fhPhiBandCluster->Fill(mom.Eta(),mom.Phi());
325 //fill histogram for UE in eta band in EMCal acceptance
326 if(mom.Phi() > (phiTrig-conesize) && mom.Phi() < (phiTrig+conesize))
328 etaBandPtSum+=mom.Pt();
329 fhEtaBandCluster->Fill(mom.Eta(),mom.Phi());
333 fhConeSumPtEtaBandUECluster ->Fill(ptTrig , etaBandPtSum);
334 fhConeSumPtPhiBandUECluster ->Fill(ptTrig , phiBandPtSum);
335 fhConeSumPtEtaBandUEClusterTrigEtaPhi->Fill(etaTrig,phiTrig,etaBandPtSum);
336 fhConeSumPtPhiBandUEClusterTrigEtaPhi->Fill(etaTrig,phiTrig,phiBandPtSum);
340 //________________________________________________________________________________________________
341 void AliAnaParticleIsolation::CalculateCaloCellUEBand(AliAODPWG4ParticleCorrelation * pCandidate,
342 Float_t & etaBandPtSumCells, Float_t & phiBandPtSumCells)
344 // Get the cells amplitude or sum of amplitude in phi/eta bands or at 45 degrees from trigger
346 if( GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kOnlyCharged ) return ;
348 Float_t conesize = GetIsolationCut()->GetConeSize();
350 Float_t phiTrig = pCandidate->Phi();
351 if(phiTrig<0) phiTrig += TMath::TwoPi();
352 Float_t etaTrig = pCandidate->Eta();
354 if(pCandidate->GetDetector()=="EMCAL")
356 AliEMCALGeometry* eGeom = AliEMCALGeometry::GetInstance();
359 if (eGeom->GetAbsCellIdFromEtaPhi(etaTrig,phiTrig,absId))
361 if(!eGeom->CheckAbsCellId(absId)) return ;
363 // Get absolute (col,row) of trigger particle
364 Int_t nSupMod = eGeom->GetSuperModuleNumber(absId);
366 Int_t imEta=-1, imPhi=-1;
367 Int_t ieta =-1, iphi =-1;
369 if (eGeom->GetCellIndex(absId,nSupMod,nModule,imPhi,imEta))
371 eGeom->GetCellPhiEtaIndexInSModule(nSupMod,nModule,imPhi,imEta,iphi,ieta);
373 Int_t colTrig = ieta;
374 if (nSupMod % 2) colTrig = AliEMCALGeoParams::fgkEMCALCols + ieta ;
375 Int_t rowTrig = iphi + AliEMCALGeoParams::fgkEMCALRows*int(nSupMod/2);
377 Int_t sqrSize = int(conesize/0.0143);
379 AliVCaloCells * cells = GetEMCALCells();
381 Int_t nTotalRows = AliEMCALGeoParams::fgkEMCALRows*16/3 ; // 24*(16/3) 5 full-size Sectors (2 SM) + 1 one-third Sector (2 SM)
382 Int_t nTotalCols = 2*AliEMCALGeoParams::fgkEMCALCols;
383 // printf("nTotalRows %i, nTotalCols %i\n",nTotalRows,nTotalCols);
384 // Loop on cells in eta band
386 Int_t irowmin = rowTrig-sqrSize;
387 if(irowmin<0) irowmin=0;
388 Int_t irowmax = rowTrig+sqrSize;
389 if(irowmax>AliEMCALGeoParams::fgkEMCALRows) irowmax=AliEMCALGeoParams::fgkEMCALRows;
392 for(Int_t irow = irowmin; irow <irowmax; irow++)
394 for(Int_t icol = 0; icol < nTotalCols; icol++)
396 Int_t inSector = int(irow/AliEMCALGeoParams::fgkEMCALRows);
397 if(inSector==5) continue;
400 if(icol < AliEMCALGeoParams::fgkEMCALCols)
402 inSupMod = 2*inSector + 1;
405 else if(icol > AliEMCALGeoParams::fgkEMCALCols - 1)
407 inSupMod = 2*inSector;
408 icolLoc = icol-AliEMCALGeoParams::fgkEMCALCols;
411 Int_t irowLoc = irow - AliEMCALGeoParams::fgkEMCALRows*inSector ;
413 // Exclude cells in cone
414 if(TMath::Abs(icol-colTrig) < sqrSize || TMath::Abs(irow-rowTrig) < sqrSize){
417 Int_t iabsId = eGeom->GetAbsCellIdFromCellIndexes(inSupMod,irowLoc,icolLoc);
418 if(!eGeom->CheckAbsCellId(iabsId)) continue;
419 etaBandPtSumCells += cells->GetCellAmplitude(iabsId);
420 fhEtaBandCell->Fill(colTrig,rowTrig);
422 // printf("ETA inSupMod %i,irowLoc %i,icolLoc %i, iabsId %i, etaBandPtSumCells %f\n",nSupMod,irowLoc,icolLoc,iabsId,etaBandPtSumCells);
425 Int_t icolmin = colTrig-sqrSize;
426 if(icolmin<0) icolmin=0;
427 Int_t icolmax = colTrig+sqrSize;
428 if(icolmax>AliEMCALGeoParams::fgkEMCALCols) icolmax=AliEMCALGeoParams::fgkEMCALCols;
430 // Loop on cells in phi band
431 for(Int_t icol = icolmin; icol < icolmax; icol++)
433 for(Int_t irow = 0; irow < nTotalRows; irow++)
435 Int_t inSector = int(irow/AliEMCALGeoParams::fgkEMCALRows);
436 if(inSector==5) continue;
439 // printf("icol %i, irow %i, inSector %i\n",icol,irow ,inSector);
440 if(icol < AliEMCALGeoParams::fgkEMCALCols)
442 // printf("icol < AliEMCALGeoParams::fgkEMCALCols %i\n",AliEMCALGeoParams::fgkEMCALCols );
443 inSupMod = 2*inSector + 1;
446 else if(icol > AliEMCALGeoParams::fgkEMCALCols - 1)
448 // printf("icol > AliEMCALGeoParams::fgkEMCALCols -1 %i\n",AliEMCALGeoParams::fgkEMCALCols -1 );
449 inSupMod = 2*inSector;
450 icolLoc = icol-AliEMCALGeoParams::fgkEMCALCols;
453 Int_t irowLoc = irow - AliEMCALGeoParams::fgkEMCALRows*inSector ; // Stesso problema di sopra //
455 // Exclude cells in cone
456 if(TMath::Abs(icol-colTrig) < sqrSize) {
457 //printf("TMath::Abs(icol-colTrig) %i < sqrSize %i\n",TMath::Abs(icol-colTrig) ,sqrSize);continue ;
459 if(TMath::Abs(irow-rowTrig) < sqrSize) {
460 //printf("TMath::Abs(irow-rowTrig) %i < sqrSize %i\n",TMath::Abs(irow-rowTrig) ,sqrSize);continue ;
463 Int_t iabsId = eGeom->GetAbsCellIdFromCellIndexes(inSupMod,irowLoc,icolLoc);
464 if(!eGeom->CheckAbsCellId(iabsId)) {printf("!eGeom->CheckAbsCellId(iabsId=%i) inSupMod %i irowLoc %i icolLoc %i \n",iabsId,inSupMod, irowLoc, icolLoc);continue;}
465 phiBandPtSumCells += cells->GetCellAmplitude(iabsId);
466 fhPhiBandCell->Fill(colTrig,rowTrig);
467 //printf("inSupMod %i,irowLoc %i,icolLoc %i, iabsId %i, phiBandPtSumCells %f\n",nSupMod,irowLoc,icolLoc,iabsId,phiBandPtSumCells);
474 Float_t ptTrig = pCandidate->Pt();
476 fhConeSumPtEtaBandUECell ->Fill(ptTrig , etaBandPtSumCells);
477 fhConeSumPtPhiBandUECell ->Fill(ptTrig , phiBandPtSumCells);
478 fhConeSumPtEtaBandUECellTrigEtaPhi->Fill(etaTrig,phiTrig,etaBandPtSumCells);
479 fhConeSumPtPhiBandUECellTrigEtaPhi->Fill(etaTrig,phiTrig,phiBandPtSumCells);
483 //________________________________________________________________________________________________
484 void AliAnaParticleIsolation::CalculateTrackUEBand(AliAODPWG4ParticleCorrelation * pCandidate,
485 Float_t & etaBandPtSum, Float_t & phiBandPtSum)
487 // Get the track pT or sum of pT in phi/eta bands or at 45 degrees from trigger
489 if( GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kOnlyNeutral ) return ;
491 Float_t conesize = GetIsolationCut()->GetConeSize();
493 Double_t sumptPerp= 0. ;
494 Float_t ptTrig = pCandidate->Pt() ;
495 Float_t phiTrig = pCandidate->Phi();
496 Float_t etaTrig = pCandidate->Eta();
498 TObjArray * trackList = GetCTSTracks() ;
499 for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++)
501 AliVTrack* track = (AliVTrack *) trackList->At(itrack);
505 printf("AliAnaParticleIsolation::CalculateTrackUEBand() - Track not available?");
509 //Do not count the candidate (pion, conversion photon) or the daughters of the candidate
510 if(track->GetID() == pCandidate->GetTrackLabel(0) || track->GetID() == pCandidate->GetTrackLabel(1) ||
511 track->GetID() == pCandidate->GetTrackLabel(2) || track->GetID() == pCandidate->GetTrackLabel(3) ) continue ;
513 // histo of eta:phi for all tracks
514 fhEtaPhiTrack->Fill(track->Eta(),track->Phi());
516 //exclude particles in cone
517 Float_t rad = GetIsolationCut()->Radius(etaTrig, phiTrig, track->Eta(), track->Phi());
519 // histo of eta:phi for all tracks in cone
520 fhEtaPhiInConeTrack->Fill(track->Eta(),track->Phi());
524 //fill histogram for UE in phi band
525 if(track->Eta() > (etaTrig-conesize) && track->Eta() < (etaTrig+conesize))
527 phiBandPtSum+=track->Pt();
528 fhPhiBandTrack->Fill(track->Eta(),track->Phi());
531 //fill histogram for UE in eta band in EMCal acceptance
532 if(track->Phi() > (phiTrig-conesize) && track->Phi() < (phiTrig+conesize))
534 etaBandPtSum+=track->Pt();
535 fhEtaBandTrack->Fill(track->Eta(),track->Phi());
538 //fill the histograms at +-45 degrees in phi from trigger particle, perpedicular to trigger axis in phi
539 Double_t dPhi = phiTrig - track->Phi() + TMath::PiOver2();
540 Double_t dEta = etaTrig - track->Eta();
541 Double_t arg = dPhi*dPhi + dEta*dEta;
542 if(TMath::Sqrt(arg) < conesize)
544 fhPtInPerpCone->Fill(ptTrig,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
545 sumptPerp+=track->Pt();
548 dPhi = phiTrig - track->Phi() - TMath::PiOver2();
549 arg = dPhi*dPhi + dEta*dEta;
550 if(TMath::Sqrt(arg) < conesize)
552 fhPtInPerpCone->Fill(ptTrig,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
553 sumptPerp+=track->Pt();
557 fhPerpConeSumPt ->Fill(ptTrig , sumptPerp );
558 fhConeSumPtEtaBandUETrack ->Fill(ptTrig , etaBandPtSum);
559 fhConeSumPtPhiBandUETrack ->Fill(ptTrig , phiBandPtSum);
560 fhConeSumPtEtaBandUETrackTrigEtaPhi->Fill(etaTrig,phiTrig,etaBandPtSum);
561 fhConeSumPtPhiBandUETrackTrigEtaPhi->Fill(etaTrig,phiTrig,phiBandPtSum);
567 //_____________________________________________________________________________________________________________________________________
568 void AliAnaParticleIsolation::CalculateNormalizeUEBandPerUnitArea(AliAODPWG4ParticleCorrelation * pCandidate, Float_t coneptsumCluster,
569 Float_t coneptsumCell, Float_t coneptsumTrack,
570 Float_t &etaBandptsumTrackNorm, Float_t &etaBandptsumClusterNorm)
572 //normalize phi/eta band per area unit
574 Float_t etaUEptsumTrack = 0 ;
575 Float_t phiUEptsumTrack = 0 ;
576 Float_t etaUEptsumCluster = 0 ;
577 Float_t phiUEptsumCluster = 0 ;
578 Float_t etaUEptsumCell = 0 ;
579 Float_t phiUEptsumCell = 0 ;
581 Int_t partTypeInCone = GetIsolationCut()->GetParticleTypeInCone();
583 // Do the normalization
585 Float_t conesize = GetIsolationCut()->GetConeSize();
586 Float_t coneA = conesize*conesize*TMath::Pi(); // A = pi R^2, isolation cone area
587 Float_t ptTrig = pCandidate->Pt() ;
588 Float_t phiTrig = pCandidate->Phi();
589 Float_t etaTrig = pCandidate->Eta();
595 Float_t phiUEptsumTrackNorm = 0 ;
596 Float_t etaUEptsumTrackNorm = 0 ;
597 Float_t coneptsumTrackSubPhi = 0 ;
598 Float_t coneptsumTrackSubEta = 0 ;
599 Float_t coneptsumTrackSubPhiNorm = 0 ;
600 Float_t coneptsumTrackSubEtaNorm = 0 ;
601 etaBandptsumTrackNorm = 0 ;
603 if( partTypeInCone!=AliIsolationCut::kOnlyNeutral )
605 // Sum the pT in the phi or eta band for clusters or tracks
606 CalculateTrackUEBand (pCandidate,etaUEptsumTrack ,phiUEptsumTrack );// rajouter ici l'histo eta phi
609 fhConeSumPtVSUETracksEtaBand->Fill(coneptsumTrack,etaUEptsumTrack);
610 fhConeSumPtVSUETracksPhiBand->Fill(coneptsumTrack,phiUEptsumTrack);
613 Float_t correctConeSumTrack = 1;
614 Float_t correctConeSumTrackPhi = 1;
616 GetIsolationCut()->CalculateUEBandTrackNormalization(GetReader(),etaTrig, phiTrig,
617 phiUEptsumTrack,etaUEptsumTrack,
618 phiUEptsumTrackNorm,etaUEptsumTrackNorm,
619 correctConeSumTrack,correctConeSumTrackPhi);
621 coneptsumTrackSubPhi = coneptsumTrack - phiUEptsumTrackNorm;
622 coneptsumTrackSubEta = coneptsumTrack - etaUEptsumTrackNorm;
624 etaBandptsumTrackNorm = etaUEptsumTrackNorm;
626 fhConeSumPtPhiUESubTrack ->Fill(ptTrig , coneptsumTrackSubPhi);
627 fhConeSumPtPhiUESubTrackTrigEtaPhi ->Fill(etaTrig, phiTrig, coneptsumTrackSubPhi);
628 fhConeSumPtEtaUESubTrack ->Fill(ptTrig , coneptsumTrackSubEta);
629 fhConeSumPtEtaUESubTrackTrigEtaPhi ->Fill(etaTrig, phiTrig, coneptsumTrackSubEta);
631 fhFractionTrackOutConeEta ->Fill(ptTrig , correctConeSumTrack-1);
632 fhFractionTrackOutConeEtaTrigEtaPhi->Fill(etaTrig, phiTrig,correctConeSumTrack-1);
634 if(coneptsumTrack > 0)
636 coneptsumTrackSubPhiNorm = coneptsumTrackSubPhi/coneptsumTrack;
637 coneptsumTrackSubEtaNorm = coneptsumTrackSubEta/coneptsumTrack;
640 fhConeSumPtSubvsConeSumPtTotPhiTrack ->Fill(coneptsumTrack,coneptsumTrackSubPhi);
641 fhConeSumPtSubNormvsConeSumPtTotPhiTrack->Fill(coneptsumTrack,coneptsumTrackSubPhiNorm);
642 fhConeSumPtSubvsConeSumPtTotEtaTrack ->Fill(coneptsumTrack,coneptsumTrackSubEta);
643 fhConeSumPtSubNormvsConeSumPtTotEtaTrack->Fill(coneptsumTrack,coneptsumTrackSubEtaNorm);
647 // ------------------------ //
648 // EMCal Clusters and cells //
649 // ------------------------ //
650 Float_t phiUEptsumClusterNorm = 0 ;
651 Float_t etaUEptsumClusterNorm = 0 ;
652 Float_t coneptsumClusterSubPhi = 0 ;
653 Float_t coneptsumClusterSubEta = 0 ;
654 Float_t coneptsumClusterSubPhiNorm = 0 ;
655 Float_t coneptsumClusterSubEtaNorm = 0 ;
656 Float_t phiUEptsumCellNorm = 0 ;
657 Float_t etaUEptsumCellNorm = 0 ;
658 Float_t coneptsumCellSubPhi = 0 ;
659 Float_t coneptsumCellSubEta = 0 ;
660 Float_t coneptsumCellSubPhiNorm = 0 ;
661 Float_t coneptsumCellSubEtaNorm = 0 ;
662 etaBandptsumClusterNorm = 0;
664 if( partTypeInCone!=AliIsolationCut::kOnlyCharged )
671 // Sum the pT in the phi or eta band for clusters or tracks
672 CalculateCaloUEBand (pCandidate,etaUEptsumCluster,phiUEptsumCluster);// rajouter ici l'histo eta phi
675 fhConeSumPtVSUEClusterEtaBand->Fill(coneptsumCluster,etaUEptsumCluster);
676 fhConeSumPtVSUEClusterPhiBand->Fill(coneptsumCluster,phiUEptsumCluster);
679 Float_t correctConeSumClusterEta = 1;
680 Float_t correctConeSumClusterPhi = 1;
682 GetIsolationCut()->CalculateUEBandClusterNormalization(GetReader(),etaTrig, phiTrig,
683 phiUEptsumCluster,etaUEptsumCluster,
684 phiUEptsumClusterNorm,etaUEptsumClusterNorm,
685 correctConeSumClusterEta,correctConeSumClusterPhi);
687 // In case that cone is out of eta and phi side, we are over correcting, not too often with the current cuts ...
688 // Comment if not used
689 // Float_t coneBadCellsCoeff =1;
690 // Float_t etaBandBadCellsCoeff=1;
691 // Float_t phiBandBadCellsCoeff=1;
692 // GetIsolationCut()->GetCoeffNormBadCell(pCandidate, GetReader(),coneBadCellsCoeff,etaBandBadCellsCoeff,phiBandBadCellsCoeff) ;
694 //coneptsumCluster=coneptsumCluster*coneBadCellsCoeff*correctConeSumClusterEta*correctConeSumClusterPhi;
696 coneptsumClusterSubPhi = coneptsumCluster - phiUEptsumClusterNorm;
697 coneptsumClusterSubEta = coneptsumCluster - etaUEptsumClusterNorm;
699 etaBandptsumClusterNorm = etaUEptsumClusterNorm;
701 fhConeSumPtPhiUESubCluster ->Fill(ptTrig , coneptsumClusterSubPhi);
702 fhConeSumPtPhiUESubClusterTrigEtaPhi ->Fill(etaTrig, phiTrig, coneptsumClusterSubPhi);
703 fhConeSumPtEtaUESubCluster ->Fill(ptTrig , coneptsumClusterSubEta);
704 fhConeSumPtEtaUESubClusterTrigEtaPhi ->Fill(etaTrig, phiTrig, coneptsumClusterSubEta);
706 fhFractionClusterOutConeEta ->Fill(ptTrig , correctConeSumClusterEta-1);
707 fhFractionClusterOutConeEtaTrigEtaPhi->Fill(etaTrig, phiTrig, correctConeSumClusterEta-1);
708 fhFractionClusterOutConePhi ->Fill(ptTrig , correctConeSumClusterPhi-1);
709 fhFractionClusterOutConePhiTrigEtaPhi->Fill(etaTrig, phiTrig, correctConeSumClusterPhi-1);
711 if(coneptsumCluster!=0)
713 coneptsumClusterSubPhiNorm = coneptsumClusterSubPhi/coneptsumCluster;
714 coneptsumClusterSubEtaNorm = coneptsumClusterSubEta/coneptsumCluster;
717 fhConeSumPtSubvsConeSumPtTotPhiCluster ->Fill(coneptsumCluster,coneptsumClusterSubPhi);
718 fhConeSumPtSubNormvsConeSumPtTotPhiCluster->Fill(coneptsumCluster,coneptsumClusterSubPhiNorm);
719 fhConeSumPtSubvsConeSumPtTotEtaCluster ->Fill(coneptsumCluster,coneptsumClusterSubEta);
720 fhConeSumPtSubNormvsConeSumPtTotEtaCluster->Fill(coneptsumCluster,coneptsumClusterSubEtaNorm);
726 if(fFillCellHistograms)
728 // Sum the pT in the phi or eta band for clusters or tracks
729 CalculateCaloCellUEBand(pCandidate,etaUEptsumCell ,phiUEptsumCell );
731 // Move to AliIsolationCut the calculation not the histograms??
733 //Careful here if EMCal limits changed .. 2010 (4 SM) to 2011-12 (10 SM), for the moment consider 100 deg in phi
734 Float_t emcEtaSize = 0.7*2; // TO FIX
735 Float_t emcPhiSize = TMath::DegToRad()*100.; // TO FIX
737 if(((2*conesize*emcPhiSize)-coneA)!=0)phiUEptsumCellNorm = phiUEptsumCell*(coneA / ((2*conesize*emcPhiSize)-coneA));
738 if(((2*conesize*emcEtaSize)-coneA)!=0)etaUEptsumCellNorm = etaUEptsumCell*(coneA / ((2*conesize*emcEtaSize)-coneA));
740 // Need to correct coneptsumCluster by the fraction of the cone out of the calorimeter cut acceptance!
742 Float_t correctConeSumCellEta = 1;
743 if(TMath::Abs(etaTrig)+conesize > emcEtaSize/2.)
745 Float_t excess = TMath::Abs(etaTrig) + conesize - emcEtaSize/2.;
746 correctConeSumCellEta = GetIsolationCut()->CalculateExcessAreaFraction(excess);
747 //printf("Excess EMC-Eta %2.3f, coneA %2.2f, excessA %2.2f, angle %2.2f,factor %2.2f\n",excess,coneA, excessA, angle*TMath::RadToDeg(), correctConeSumClusterEta);
748 // Need to correct phi band surface if part of the cone falls out of track cut acceptance!
749 if(((2*(conesize-excess)*emcPhiSize)-(coneA-correctConeSumCellEta))!=0)phiUEptsumCellNorm = phiUEptsumCell*(coneA / ((2*(conesize-excess)*emcPhiSize)-(coneA-correctConeSumCellEta)));
752 Float_t correctConeSumCellPhi = 1;
753 //printf("EMCPhiTrig %2.2f, conesize %2.2f, sum %2.2f, rest %2.2f \n",phiTrig*TMath::RadToDeg(),conesize*TMath::RadToDeg(),(phiTrig+conesize)*TMath::RadToDeg(),(phiTrig-conesize)*TMath::RadToDeg() );
754 if((phiTrig+conesize > 180*TMath::DegToRad()) ||
755 (phiTrig-conesize < 80*TMath::DegToRad()))
758 if( phiTrig+conesize > 180*TMath::DegToRad() ) excess = conesize + phiTrig - 180*TMath::DegToRad() ;
759 else excess = conesize - phiTrig + 80*TMath::DegToRad() ;
761 correctConeSumCellPhi = GetIsolationCut()->CalculateExcessAreaFraction(excess);
762 //printf("Excess EMC-Phi %2.3f, coneA %2.2f, excessA %2.2f, angle %2.2f,factor %2.2f\n",excess,coneA, excessA, angle*TMath::RadToDeg(), correctConeSumClusterPhi);
764 // Need to correct eta band surface if part of the cone falls out of track cut acceptance!
765 if(((2*(conesize-excess)*emcEtaSize)-(coneA-correctConeSumCellPhi))!=0)etaUEptsumCellNorm = etaUEptsumCell*(coneA / ((2*(conesize-excess)*emcEtaSize)-(coneA-correctConeSumCellPhi)));
769 // In case that cone is out of eta and phi side, we are over correcting, not too often with the current cuts ...
770 coneptsumCellSubPhi = coneptsumCell*correctConeSumCellEta*correctConeSumCellPhi - phiUEptsumCellNorm;
771 coneptsumCellSubEta = coneptsumCell*correctConeSumCellEta*correctConeSumCellPhi - etaUEptsumCellNorm;
773 fhConeSumPtPhiUESubCell ->Fill(ptTrig , coneptsumCellSubPhi);
774 fhConeSumPtPhiUESubCellTrigEtaPhi ->Fill(etaTrig, phiTrig, coneptsumCellSubPhi);
775 fhConeSumPtEtaUESubCell ->Fill(ptTrig , coneptsumCellSubEta);
776 fhConeSumPtEtaUESubCellTrigEtaPhi ->Fill(etaTrig, phiTrig, coneptsumCellSubEta);
778 fhFractionCellOutConeEta ->Fill(ptTrig , correctConeSumCellEta-1);
779 fhFractionCellOutConeEtaTrigEtaPhi->Fill(etaTrig, phiTrig, correctConeSumCellEta-1);
780 fhFractionCellOutConePhi ->Fill(ptTrig , correctConeSumCellPhi-1);
781 fhFractionCellOutConePhiTrigEtaPhi->Fill(etaTrig, phiTrig, correctConeSumCellPhi-1);
784 coneptsumCellSubPhiNorm = coneptsumCellSubPhi/coneptsumCell;
785 coneptsumCellSubEtaNorm = coneptsumCellSubEta/coneptsumCell;
788 fhConeSumPtSubvsConeSumPtTotPhiCell ->Fill(coneptsumCell,coneptsumCellSubPhi);
789 fhConeSumPtSubNormvsConeSumPtTotPhiCell->Fill(coneptsumCell,coneptsumCellSubPhiNorm);
790 fhConeSumPtSubvsConeSumPtTotEtaCell ->Fill(coneptsumCell,coneptsumCellSubEta);
791 fhConeSumPtSubNormvsConeSumPtTotEtaCell->Fill(coneptsumCell,coneptsumCellSubEtaNorm);
795 if( partTypeInCone==AliIsolationCut::kNeutralAndCharged )
797 // --------------------------- //
798 // Tracks and clusters in cone //
799 // --------------------------- //
801 Double_t sumPhiUESub = coneptsumClusterSubPhi + coneptsumTrackSubPhi;
802 Double_t sumEtaUESub = coneptsumClusterSubEta + coneptsumTrackSubEta;
804 fhConeSumPtPhiUESub ->Fill(ptTrig , sumPhiUESub);
805 fhConeSumPtPhiUESubTrigEtaPhi->Fill(etaTrig, phiTrig, sumPhiUESub);
806 fhConeSumPtEtaUESub ->Fill(ptTrig , sumEtaUESub);
807 fhConeSumPtEtaUESubTrigEtaPhi->Fill(etaTrig, phiTrig, sumEtaUESub);
809 fhEtaBandClustervsTrack ->Fill(etaUEptsumCluster ,etaUEptsumTrack );
810 fhPhiBandClustervsTrack ->Fill(phiUEptsumCluster ,phiUEptsumTrack );
811 fhEtaBandNormClustervsTrack->Fill(etaUEptsumClusterNorm,etaUEptsumTrackNorm);
812 fhPhiBandNormClustervsTrack->Fill(phiUEptsumClusterNorm,phiUEptsumTrackNorm);
814 fhConeSumPtEtaUESubClustervsTrack->Fill(coneptsumClusterSubEta,coneptsumTrackSubEta);
815 fhConeSumPtPhiUESubClustervsTrack->Fill(coneptsumClusterSubPhi,coneptsumTrackSubPhi);
817 // ------------------------ //
818 // Tracks and cells in cone //
819 // ------------------------ //
821 if(fFillCellHistograms)
823 Double_t sumPhiUESubTrackCell = coneptsumCellSubPhi + coneptsumTrackSubPhi;
824 Double_t sumEtaUESubTrackCell = coneptsumCellSubEta + coneptsumTrackSubEta;
826 fhConeSumPtPhiUESubTrackCell ->Fill(ptTrig , sumPhiUESubTrackCell);
827 fhConeSumPtPhiUESubTrackCellTrigEtaPhi->Fill(etaTrig, phiTrig, sumPhiUESubTrackCell);
828 fhConeSumPtEtaUESubTrackCell ->Fill(ptTrig , sumEtaUESubTrackCell);
829 fhConeSumPtEtaUESubTrackCellTrigEtaPhi->Fill(etaTrig, phiTrig, sumEtaUESubTrackCell);
831 fhEtaBandCellvsTrack ->Fill(etaUEptsumCell ,etaUEptsumTrack );
832 fhPhiBandCellvsTrack ->Fill(phiUEptsumCell ,phiUEptsumTrack );
833 fhEtaBandNormCellvsTrack->Fill(etaUEptsumCellNorm,etaUEptsumTrackNorm);
834 fhPhiBandNormCellvsTrack->Fill(phiUEptsumCellNorm,phiUEptsumTrackNorm);
836 fhConeSumPtEtaUESubCellvsTrack->Fill(coneptsumCellSubEta,coneptsumTrackSubEta);
837 fhConeSumPtPhiUESubCellvsTrack->Fill(coneptsumCellSubPhi,coneptsumTrackSubPhi);
844 //__________________________________________________________________________________________________
845 void AliAnaParticleIsolation::CalculateCaloSignalInCone(AliAODPWG4ParticleCorrelation * aodParticle,
846 Float_t & coneptsumCluster)
848 // Get the cluster pT or sum of pT in isolation cone
850 if( GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kOnlyCharged ) return ;
852 //Recover reference arrays with clusters and tracks
853 TObjArray * refclusters = aodParticle->GetObjArray(GetAODObjArrayName()+"Clusters");
854 if(!refclusters) return ;
856 Float_t ptTrig = aodParticle->Pt();
858 //Get vertex for cluster momentum calculation
859 Double_t vertex[] = {0,0,0} ; //vertex ;
860 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
861 GetReader()->GetVertex(vertex);
864 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++)
866 AliVCluster* calo = (AliVCluster *) refclusters->At(icalo);
867 calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
869 fhPtInCone ->Fill(ptTrig, mom.Pt());
870 fhPtClusterInCone->Fill(ptTrig, mom.Pt());
872 if(fFillPileUpHistograms)
874 if(GetReader()->IsPileUpFromSPD()) fhPtInConePileUp[0]->Fill(ptTrig,mom.Pt());
875 if(GetReader()->IsPileUpFromEMCal()) fhPtInConePileUp[1]->Fill(ptTrig,mom.Pt());
876 if(GetReader()->IsPileUpFromSPDOrEMCal()) fhPtInConePileUp[2]->Fill(ptTrig,mom.Pt());
877 if(GetReader()->IsPileUpFromSPDAndEMCal()) fhPtInConePileUp[3]->Fill(ptTrig,mom.Pt());
878 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhPtInConePileUp[4]->Fill(ptTrig,mom.Pt());
879 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhPtInConePileUp[5]->Fill(ptTrig,mom.Pt());
880 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtInConePileUp[6]->Fill(ptTrig,mom.Pt());
883 if(fFillHighMultHistograms) fhPtInConeCent->Fill(GetEventCentrality(),mom.Pt());
885 coneptsumCluster+=mom.Pt();
888 fhConeSumPtCluster ->Fill(ptTrig, coneptsumCluster);
891 //______________________________________________________________________________________________________
892 void AliAnaParticleIsolation::CalculateCaloCellSignalInCone(AliAODPWG4ParticleCorrelation * aodParticle,
893 Float_t & coneptsumCell)
895 // Get the cell amplityde or sum of amplitudes in isolation cone
896 // Mising: Remove signal cells in cone in case the trigger is a cluster!
898 if( GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kOnlyCharged ) return ;
900 Float_t conesize = GetIsolationCut()->GetConeSize();
902 Float_t ptTrig = aodParticle->Pt();
903 Float_t phiTrig = aodParticle->Phi();
904 if(phiTrig<0) phiTrig += TMath::TwoPi();
905 Float_t etaTrig = aodParticle->Eta();
907 if(aodParticle->GetDetector()=="EMCAL")
909 AliEMCALGeometry* eGeom = AliEMCALGeometry::GetInstance();
912 if (eGeom->GetAbsCellIdFromEtaPhi(etaTrig,phiTrig,absId))
914 if(!eGeom->CheckAbsCellId(absId)) return ;
916 // Get absolute (col,row) of trigger particle
917 Int_t nSupMod = eGeom->GetSuperModuleNumber(absId);
919 Int_t imEta=-1, imPhi=-1;
920 Int_t ieta =-1, iphi =-1;
922 if (eGeom->GetCellIndex(absId,nSupMod,nModule,imPhi,imEta))
924 Int_t iEta=-1, iPhi=-1;
925 eGeom->GetCellPhiEtaIndexInSModule(nSupMod,nModule,imPhi,imEta,iphi,ieta);
927 Int_t colTrig = iEta;
928 if (nSupMod % 2) colTrig = AliEMCALGeoParams::fgkEMCALCols + iEta ;
929 Int_t rowTrig = iPhi + AliEMCALGeoParams::fgkEMCALRows*int(nSupMod/2);
931 Int_t sqrSize = int(conesize/0.0143);
933 AliVCaloCells * cells = GetEMCALCells();
935 // Loop on cells in cone
936 for(Int_t irow = rowTrig-sqrSize; irow < rowTrig+sqrSize; irow++)
938 for(Int_t icol = colTrig-sqrSize; icol < colTrig+sqrSize; icol++)
940 Int_t inSector = int(irow/AliEMCALGeoParams::fgkEMCALRows);
941 if(inSector==5) continue;
945 if(icol < AliEMCALGeoParams::fgkEMCALCols)
947 inSupMod = 2*inSector + 1;
950 else if(icol > AliEMCALGeoParams::fgkEMCALCols - 1)
952 inSupMod = 2*inSector;
953 icolLoc = icol-AliEMCALGeoParams::fgkEMCALCols;
956 Int_t irowLoc = irow - AliEMCALGeoParams::fgkEMCALRows*inSector ;
958 Int_t iabsId = eGeom->GetAbsCellIdFromCellIndexes(inSupMod,irowLoc,icolLoc);
959 if(!eGeom->CheckAbsCellId(iabsId)) continue;
961 fhPtCellInCone->Fill(ptTrig, cells->GetCellAmplitude(iabsId));
962 coneptsumCell += cells->GetCellAmplitude(iabsId);
969 fhConeSumPtCell->Fill(ptTrig,coneptsumCell);
973 //___________________________________________________________________________________________________
974 void AliAnaParticleIsolation::CalculateTrackSignalInCone(AliAODPWG4ParticleCorrelation * aodParticle,
975 Float_t & coneptsumTrack)
977 // Get the track pT or sum of pT in isolation cone
979 if( GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kOnlyNeutral ) return ;
981 //Recover reference arrays with clusters and tracks
982 TObjArray * reftracks = aodParticle->GetObjArray(GetAODObjArrayName()+"Tracks");
983 if(!reftracks) return ;
985 Float_t ptTrig = aodParticle->Pt();
986 Double_t bz = GetReader()->GetInputEvent()->GetMagneticField();
988 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++)
990 AliVTrack* track = (AliVTrack *) reftracks->At(itrack);
991 Float_t pTtrack = track->Pt();
993 fhPtInCone ->Fill(ptTrig,pTtrack);
994 fhPtTrackInCone->Fill(ptTrig,pTtrack);
996 if(fFillPileUpHistograms)
998 ULong_t status = track->GetStatus();
999 Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
1000 //Double32_t tof = track->GetTOFsignal()*1e-3;
1001 Int_t trackBC = track->GetTOFBunchCrossing(bz);
1003 if ( okTOF && trackBC!=0 ) fhPtTrackInConeOtherBC->Fill(ptTrig,pTtrack);
1004 else if( okTOF && trackBC==0 ) fhPtTrackInConeBC0 ->Fill(ptTrig,pTtrack);
1006 Int_t vtxBC = GetReader()->GetVertexBC();
1007 if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA) fhPtTrackInConeVtxBC0->Fill(ptTrig,pTtrack);
1009 if(GetReader()->IsPileUpFromSPD()) { fhPtInConePileUp[0]->Fill(ptTrig,pTtrack);
1010 if(okTOF && trackBC!=0 ) fhPtTrackInConeOtherBCPileUpSPD->Fill(ptTrig,pTtrack);
1011 if(okTOF && trackBC==0 ) fhPtTrackInConeBC0PileUpSPD ->Fill(ptTrig,pTtrack); }
1012 if(GetReader()->IsPileUpFromEMCal()) fhPtInConePileUp[1]->Fill(ptTrig,pTtrack);
1013 if(GetReader()->IsPileUpFromSPDOrEMCal()) fhPtInConePileUp[2]->Fill(ptTrig,pTtrack);
1014 if(GetReader()->IsPileUpFromSPDAndEMCal()) fhPtInConePileUp[3]->Fill(ptTrig,pTtrack);
1015 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhPtInConePileUp[4]->Fill(ptTrig,pTtrack);
1016 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhPtInConePileUp[5]->Fill(ptTrig,pTtrack);
1017 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtInConePileUp[6]->Fill(ptTrig,pTtrack);
1020 if(fFillHighMultHistograms) fhPtInConeCent->Fill(GetEventCentrality(),pTtrack);
1022 coneptsumTrack+=pTtrack;
1025 fhConeSumPtTrack->Fill(ptTrig, coneptsumTrack);
1029 //_________________________________________________________________
1030 void AliAnaParticleIsolation::FillPileUpHistograms(Int_t clusterID)
1032 // Fill some histograms to understand pile-up
1036 printf("AliAnaParticleIsolation::FillPileUpHistograms(), ID of cluster = %d, not possible! ", clusterID);
1041 TObjArray* clusters = 0x0;
1042 if (fCalorimeter == "EMCAL") clusters = GetEMCALClusters();
1043 else if(fCalorimeter == "PHOS" ) clusters = GetPHOSClusters();
1046 Float_t time = -1000;
1050 AliVCluster *cluster = FindCluster(clusters,clusterID,iclus);
1051 energy = cluster->E();
1052 time = cluster->GetTOF()*1e9;
1055 //printf("E %f, time %f\n",energy,time);
1056 AliVEvent * event = GetReader()->GetInputEvent();
1058 fhTimeENoCut->Fill(energy,time);
1059 if(GetReader()->IsPileUpFromSPD()) fhTimeESPD ->Fill(energy,time);
1060 if(event->IsPileupFromSPDInMultBins()) fhTimeESPDMulti->Fill(energy,time);
1062 if(energy < 8) return; // Fill time figures for high energy clusters not too close to trigger threshold
1064 AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
1065 AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
1067 // N pile up vertices
1068 Int_t nVerticesSPD = -1;
1069 Int_t nVerticesTracks = -1;
1073 nVerticesSPD = esdEv->GetNumberOfPileupVerticesSPD();
1074 nVerticesTracks = esdEv->GetNumberOfPileupVerticesTracks();
1079 nVerticesSPD = aodEv->GetNumberOfPileupVerticesSPD();
1080 nVerticesTracks = aodEv->GetNumberOfPileupVerticesTracks();
1083 fhTimeNPileUpVertSPD ->Fill(time,nVerticesSPD);
1084 fhTimeNPileUpVertTrack->Fill(time,nVerticesTracks);
1086 //printf("Is SPD %d, Is SPD Multi %d, n spd %d, n track %d\n",
1087 // GetReader()->IsPileUpFromSPD(),event->IsPileupFromSPDInMultBins(),nVerticesSPD,nVerticesTracks);
1090 Float_t z1 = -1, z2 = -1;
1092 for(Int_t iVert=0; iVert<nVerticesSPD;iVert++)
1096 const AliESDVertex* pv=esdEv->GetPileupVertexSPD(iVert);
1097 ncont=pv->GetNContributors();
1098 z1 = esdEv->GetPrimaryVertexSPD()->GetZ();
1100 diamZ = esdEv->GetDiamondZ();
1104 AliAODVertex *pv=aodEv->GetVertex(iVert);
1105 if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
1106 ncont=pv->GetNContributors();
1107 z1=aodEv->GetPrimaryVertexSPD()->GetZ();
1109 diamZ = aodEv->GetDiamondZ();
1112 Double_t distZ = TMath::Abs(z2-z1);
1113 diamZ = TMath::Abs(z2-diamZ);
1115 fhTimeNPileUpVertContributors ->Fill(time,ncont);
1116 fhTimePileUpMainVertexZDistance->Fill(time,distZ);
1117 fhTimePileUpMainVertexZDiamond ->Fill(time,diamZ);
1122 //_____________________________________________________________________________________________________________________
1123 void AliAnaParticleIsolation::FillTrackMatchingShowerShapeControlHistograms(AliAODPWG4ParticleCorrelation *pCandidate,
1126 // Fill Track matching and Shower Shape control histograms
1127 if(!fFillTMHisto && !fFillSSHisto) return;
1129 Int_t clusterID = pCandidate->GetCaloLabel(0) ;
1130 Int_t nMaxima = pCandidate->GetFiducialArea(); // bad name, just place holder for the moment
1131 Int_t mcTag = pCandidate->GetTag() ;
1132 Bool_t isolated = pCandidate->IsIsolated();
1136 printf("AliAnaParticleIsolation::FillTrackMatchingShowerShapeControlHistograms(), ID of cluster = %d, not possible! \n", clusterID);
1141 TObjArray* clusters = 0x0;
1142 if (fCalorimeter == "EMCAL") clusters = GetEMCALClusters();
1143 else if(fCalorimeter == "PHOS" ) clusters = GetPHOSClusters();
1145 Float_t energy = pCandidate->E();
1146 Float_t pt = pCandidate->Pt();
1150 AliVCluster *cluster = FindCluster(clusters,clusterID,iclus);
1154 fhELambda0 [isolated]->Fill(energy, cluster->GetM02() );
1155 fhPtLambda0[isolated]->Fill(pt, cluster->GetM02() );
1156 fhELambda1 [isolated]->Fill(energy, cluster->GetM20() );
1160 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
1161 fhPtLambda0MC[kmcPhoton][isolated]->Fill(pt, cluster->GetM02());
1163 fhPtLambda0MC[mcIndex][isolated]->Fill(pt, cluster->GetM02());
1166 if(fCalorimeter == "EMCAL" && GetFirstSMCoveredByTRD() >= 0 &&
1167 GetModuleNumber(cluster) >= GetFirstSMCoveredByTRD() )
1169 fhELambda0TRD [isolated]->Fill(energy, cluster->GetM02() );
1170 fhPtLambda0TRD[isolated]->Fill(pt , cluster->GetM02() );
1171 fhELambda1TRD [isolated]->Fill(energy, cluster->GetM20() );
1174 if(fFillNLMHistograms)
1176 fhNLocMax[isolated]->Fill(energy,nMaxima);
1177 if (nMaxima==1) { fhELambda0LocMax1[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMax1[isolated]->Fill(energy,cluster->GetM20()); }
1178 else if(nMaxima==2) { fhELambda0LocMax2[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMax2[isolated]->Fill(energy,cluster->GetM20()); }
1179 else { fhELambda0LocMaxN[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMaxN[isolated]->Fill(energy,cluster->GetM20()); }
1186 Float_t dZ = cluster->GetTrackDz();
1187 Float_t dR = cluster->GetTrackDx();
1189 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
1191 dR = 2000., dZ = 2000.;
1192 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
1195 //printf("ParticleIsolation: dPhi %f, dEta %f\n",dR,dZ);
1196 if(fhTrackMatchedDEta[isolated] && TMath::Abs(dR) < 999)
1198 fhTrackMatchedDEta[isolated]->Fill(energy,dZ);
1199 fhTrackMatchedDPhi[isolated]->Fill(energy,dR);
1200 if(energy > 0.5) fhTrackMatchedDEtaDPhi[isolated]->Fill(dZ,dR);
1203 // Check dEdx and E/p of matched clusters
1205 if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
1208 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
1212 Float_t dEdx = track->GetTPCsignal();
1213 fhdEdx[isolated]->Fill(cluster->E(), dEdx);
1215 Float_t eOverp = cluster->E()/track->P();
1216 fhEOverP[isolated]->Fill(cluster->E(), eOverp);
1219 // printf("AliAnaParticleIsolation::FillTrackMatchingShowerShapeHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
1224 if ( !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) )
1226 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
1227 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 2.5 );
1228 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 0.5 );
1229 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 1.5 );
1230 else fhTrackMatchedMCParticle[isolated]->Fill(energy, 3.5 );
1235 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
1236 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 6.5 );
1237 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 4.5 );
1238 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 5.5 );
1239 else fhTrackMatchedMCParticle[isolated]->Fill(energy, 7.5 );
1248 } // clusters array available
1252 //______________________________________________________
1253 TObjString * AliAnaParticleIsolation::GetAnalysisCuts()
1255 //Save parameters used for analysis
1256 TString parList ; //this will be list of parameters used for this analysis.
1257 const Int_t buffersize = 255;
1258 char onePar[buffersize] ;
1260 snprintf(onePar, buffersize,"--- AliAnaParticleIsolation ---\n") ;
1262 snprintf(onePar, buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
1264 snprintf(onePar, buffersize,"Isolation Cand Detector: %s\n",fIsoDetector.Data()) ;
1266 snprintf(onePar, buffersize,"fReMakeIC =%d (Flag for reisolation during histogram filling) \n",fReMakeIC) ;
1268 snprintf(onePar, buffersize,"fMakeSeveralIC=%d (Flag for isolation with several cuts at the same time ) \n",fMakeSeveralIC) ;
1270 snprintf(onePar, buffersize,"fFillTMHisto=%d (Flag for track matching histograms) \n",fFillTMHisto) ;
1272 snprintf(onePar, buffersize,"fFillSSHisto=%d (Flag for shower shape histograms) \n",fFillSSHisto) ;
1277 snprintf(onePar, buffersize,"fNCones =%d (Number of cone sizes) \n",fNCones) ;
1279 snprintf(onePar, buffersize,"fNPtThresFrac=%d (Flag for isolation with several cuts at the same time ) \n",fNPtThresFrac) ;
1282 for(Int_t icone = 0; icone < fNCones ; icone++)
1284 snprintf(onePar, buffersize,"fConeSizes[%d]=%1.2f (isolation cone size) \n",icone, fConeSizes[icone]) ;
1287 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
1289 snprintf(onePar, buffersize,"fPtThresholds[%d]=%1.2f (isolation pt threshold) \n",ipt, fPtThresholds[ipt]) ;
1292 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
1294 snprintf(onePar, buffersize,"fPtFractions[%d]=%1.2f (isolation pt fraction threshold) \n",ipt, fPtFractions[ipt]) ;
1297 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
1299 snprintf(onePar, buffersize,"fSumPtThresholds[%d]=%1.2f (isolation sum pt threshold) \n",ipt, fSumPtThresholds[ipt]) ;
1304 //Get parameters set in base class.
1305 parList += GetBaseParametersList() ;
1307 //Get parameters set in IC class.
1308 if(!fMakeSeveralIC)parList += GetIsolationCut()->GetICParametersList() ;
1310 return new TObjString(parList) ;
1314 //________________________________________________________
1315 TList * AliAnaParticleIsolation::GetCreateOutputObjects()
1317 // Create histograms to be saved in output file and
1318 // store them in outputContainer
1319 TList * outputContainer = new TList() ;
1320 outputContainer->SetName("IsolatedParticleHistos") ;
1322 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins();
1323 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins();
1324 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
1325 Float_t ptmax = GetHistogramRanges()->GetHistoPtMax();
1326 Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();
1327 Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
1328 Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
1329 Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
1330 Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
1331 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins();
1332 Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax();
1333 Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
1334 Int_t ntimebins= GetHistogramRanges()->GetHistoTimeBins();
1335 Float_t timemax = GetHistogramRanges()->GetHistoTimeMax();
1336 Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
1338 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
1339 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
1340 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
1341 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
1342 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
1343 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
1345 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins();
1346 Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax();
1347 Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
1348 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
1349 Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax();
1350 Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
1352 Int_t nptsumbins = GetHistogramRanges()->GetHistoNPtSumBins();
1353 Float_t ptsummax = GetHistogramRanges()->GetHistoPtSumMax();
1354 Float_t ptsummin = GetHistogramRanges()->GetHistoPtSumMin();
1355 Int_t nptinconebins = GetHistogramRanges()->GetHistoNPtInConeBins();
1356 Float_t ptinconemax = GetHistogramRanges()->GetHistoPtInConeMax();
1357 Float_t ptinconemin = GetHistogramRanges()->GetHistoPtInConeMin();
1359 //Float_t ptthre = GetIsolationCut()->GetPtThreshold();
1360 //Float_t ptsumthre = GetIsolationCut()->GetSumPtThreshold();
1361 //Float_t ptfrac = GetIsolationCut()->GetPtFraction();
1362 Float_t r = GetIsolationCut()->GetConeSize();
1363 Int_t method = GetIsolationCut()->GetICMethod() ;
1364 Int_t particle = GetIsolationCut()->GetParticleTypeInCone() ;
1366 TString sThreshold = "";
1367 if ( method == AliIsolationCut::kSumPtIC )
1369 sThreshold = Form(", %2.2f < #Sigma #it{p}_{T}^{in cone} < %2.2f GeV/#it{c}",
1370 GetIsolationCut()->GetSumPtThreshold(), GetIsolationCut()->GetSumPtThresholdMax());
1371 if(GetIsolationCut()->GetSumPtThresholdMax() > 200)
1372 sThreshold = Form(", #Sigma #it{p}_{T}^{in cone} = %2.2f GeV/#it{c}",
1373 GetIsolationCut()->GetSumPtThreshold());
1375 else if ( method == AliIsolationCut::kPtThresIC)
1377 sThreshold = Form(", %2.2f < #it{p}_{T}^{th} < %2.2f GeV/#it{c}",
1378 GetIsolationCut()->GetPtThreshold(),GetIsolationCut()->GetPtThresholdMax());
1379 if(GetIsolationCut()->GetSumPtThreshold() > 200)
1380 sThreshold = Form(", #it{p}_{T}^{th} = %2.2f GeV/#it{c}",
1381 GetIsolationCut()->GetPtThreshold());
1383 else if ( method == AliIsolationCut::kPtFracIC)
1384 sThreshold = Form(", #Sigma #it{p}_{T}^{in cone}/#it{p}_{T}^{trig} = %2.2f" ,
1385 GetIsolationCut()->GetPtFraction());
1387 TString sParticle = ", x^{0,#pm}";
1388 if ( particle == AliIsolationCut::kOnlyNeutral ) sParticle = ", x^{0}";
1389 else if ( particle == AliIsolationCut::kOnlyCharged ) sParticle = ", x^{#pm}";
1391 TString parTitle = Form("#it{R} = %2.2f%s%s",GetIsolationCut()->GetConeSize(), sThreshold.Data(),sParticle.Data());
1393 TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
1395 // MC histograms title and name
1396 TString mcPartType[] = { "#gamma", "#gamma_{prompt}", "#gamma_{fragmentation}",
1397 "#pi^{0} (merged #gamma)","#gamma_{#pi decay}",
1398 "#gamma_{#eta decay}","#gamma_{other decay}",
1399 "e^{#pm}","hadrons?"} ;
1401 TString mcPartName[] = { "Photon","PhotonPrompt","PhotonFrag",
1402 "Pi0","Pi0Decay","EtaDecay","OtherDecay",
1403 "Electron","Hadron"} ;
1405 // Primary MC histograms title and name
1406 TString pptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}",
1407 "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}"} ;
1409 TString ppname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay",
1410 "PhotonPrompt","PhotonFrag","PhotonISR"} ;
1412 // Not Isolated histograms, reference histograms
1414 fhENoIso = new TH1F("hENoIso",
1415 Form("Number of not isolated leading particles vs #it{p}_{T}, %s",parTitle.Data()),
1416 nptbins,ptmin,ptmax);
1417 fhENoIso->SetYTitle("#it{counts}");
1418 fhENoIso->SetXTitle("E (GeV/#it{c})");
1419 outputContainer->Add(fhENoIso) ;
1421 fhPtNoIso = new TH1F("hPtNoIso",
1422 Form("Number of not isolated leading particles vs #it{p}_{T}, %s",parTitle.Data()),
1423 nptbins,ptmin,ptmax);
1424 fhPtNoIso->SetYTitle("#it{counts}");
1425 fhPtNoIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1426 outputContainer->Add(fhPtNoIso) ;
1428 fhEtaPhiNoIso = new TH2F("hEtaPhiNoIso",
1429 Form("Number of not isolated leading particles #eta vs #phi, %s",parTitle.Data()),
1430 netabins,etamin,etamax,nphibins,phimin,phimax);
1431 fhEtaPhiNoIso->SetXTitle("#eta");
1432 fhEtaPhiNoIso->SetYTitle("#phi");
1433 outputContainer->Add(fhEtaPhiNoIso) ;
1437 // For histograms in arrays, index in the array, corresponding to any particle origin
1439 for(Int_t imc = 0; imc < 9; imc++)
1442 fhPtNoIsoMC[imc] = new TH1F(Form("hPtNoIsoMC%s",mcPartName[imc].Data()),
1443 Form("#it{p}_{T} of NOT isolated %s, %s",mcPartType[imc].Data(),parTitle.Data()),
1444 nptbins,ptmin,ptmax);
1445 fhPtNoIsoMC[imc]->SetYTitle("#it{counts}");
1446 fhPtNoIsoMC[imc]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1447 outputContainer->Add(fhPtNoIsoMC[imc]) ;
1449 fhPtIsoMC[imc] = new TH1F(Form("hPtMC%s",mcPartName[imc].Data()),
1450 Form("#it{p}_{T} of isolated %s, %s",mcPartType[imc].Data(),parTitle.Data()),
1451 nptbins,ptmin,ptmax);
1452 fhPtIsoMC[imc]->SetYTitle("#it{counts}");
1453 fhPtIsoMC[imc]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1454 outputContainer->Add(fhPtIsoMC[imc]) ;
1456 fhPhiIsoMC[imc] = new TH2F(Form("hPhiMC%s",mcPartName[imc].Data()),
1457 Form("#phi vs #it{p}_{T} of isolated %s, %s",mcPartType[imc].Data(),parTitle.Data()),
1458 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1459 fhPhiIsoMC[imc]->SetYTitle("#phi");
1460 fhPhiIsoMC[imc]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1461 outputContainer->Add(fhPhiIsoMC[imc]) ;
1463 fhEtaIsoMC[imc] = new TH2F(Form("hEtaMC%s",mcPartName[imc].Data()),
1464 Form("#phi vs #it{p}_{T} of isolated %s, %s",mcPartType[imc].Data(),parTitle.Data()),
1465 nptbins,ptmin,ptmax,netabins,etamin,etamax);
1466 fhEtaIsoMC[imc]->SetYTitle("#eta");
1467 fhEtaIsoMC[imc]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1468 outputContainer->Add(fhEtaIsoMC[imc]) ;
1472 // Histograms for tagged candidates as decay
1473 if(fFillTaggedDecayHistograms)
1475 fhPtDecayNoIso = new TH1F("hPtDecayNoIso",
1476 Form("Number of not isolated leading pi0 decay particles vs #it{p}_{T}, %s",parTitle.Data()),
1477 nptbins,ptmin,ptmax);
1478 fhPtDecayNoIso->SetYTitle("#it{counts}");
1479 fhPtDecayNoIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1480 outputContainer->Add(fhPtDecayNoIso) ;
1482 fhEtaPhiDecayNoIso = new TH2F("hEtaPhiDecayNoIso",
1483 Form("Number of not isolated leading Pi0 decay particles #eta vs #phi, %s",parTitle.Data()),
1484 netabins,etamin,etamax,nphibins,phimin,phimax);
1485 fhEtaPhiDecayNoIso->SetXTitle("#eta");
1486 fhEtaPhiDecayNoIso->SetYTitle("#phi");
1487 outputContainer->Add(fhEtaPhiDecayNoIso) ;
1491 fhPtDecayIso = new TH1F("hPtDecayIso",
1492 Form("Number of isolated #pi^{0} decay particles vs #it{p}_{T}, %s",parTitle.Data()),
1493 nptbins,ptmin,ptmax);
1494 fhPtDecayIso->SetYTitle("#it{counts}");
1495 fhPtDecayIso->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1496 outputContainer->Add(fhPtDecayIso) ;
1498 fhEtaPhiDecayIso = new TH2F("hEtaPhiDecayIso",
1499 Form("Number of isolated Pi0 decay particles #eta vs #phi, %s",parTitle.Data()),
1500 netabins,etamin,etamax,nphibins,phimin,phimax);
1501 fhEtaPhiDecayIso->SetXTitle("#eta");
1502 fhEtaPhiDecayIso->SetYTitle("#phi");
1503 outputContainer->Add(fhEtaPhiDecayIso) ;
1509 for(Int_t imc = 0; imc < 9; imc++)
1512 fhPtDecayNoIsoMC[imc] = new TH1F(Form("hPtDecayNoIsoMC%s",mcPartName[imc].Data()),
1513 Form("#it{p}_{T} of NOT isolated, decay tag, %s, %s",mcPartType[imc].Data(),parTitle.Data()),
1514 nptbins,ptmin,ptmax);
1515 fhPtDecayNoIsoMC[imc]->SetYTitle("#it{counts}");
1516 fhPtDecayNoIsoMC[imc]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1517 outputContainer->Add(fhPtDecayNoIsoMC[imc]) ;
1521 fhPtDecayIsoMC[imc] = new TH1F(Form("hPtDecayMC%s",mcPartName[imc].Data()),
1522 Form("#it{p}_{T} of isolated %s, decay tag, %s",mcPartType[imc].Data(),parTitle.Data()),
1523 nptbins,ptmin,ptmax);
1524 fhPtDecayIsoMC[imc]->SetYTitle("#it{counts}");
1525 fhPtDecayIsoMC[imc]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1526 outputContainer->Add(fhPtDecayIsoMC[imc]) ;
1534 TString isoName [] = {"NoIso",""};
1535 TString isoTitle[] = {"Not isolated" ,"isolated"};
1537 fhEIso = new TH1F("hE",
1538 Form("Number of isolated particles vs E, %s",parTitle.Data()),
1539 nptbins,ptmin,ptmax);
1540 fhEIso->SetYTitle("d#it{N} / d#it{E}");
1541 fhEIso->SetXTitle("#it{E} (GeV/#it{c})");
1542 outputContainer->Add(fhEIso) ;
1544 fhPtIso = new TH1F("hPt",
1545 Form("Number of isolated particles vs #it{p}_{T}, %s",parTitle.Data()),
1546 nptbins,ptmin,ptmax);
1547 fhPtIso->SetYTitle("d#it{N} / #it{p}_{T}");
1548 fhPtIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1549 outputContainer->Add(fhPtIso) ;
1551 fhPhiIso = new TH2F("hPhi",
1552 Form("Number of isolated particles vs #phi, %s",parTitle.Data()),
1553 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1554 fhPhiIso->SetYTitle("#phi");
1555 fhPhiIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1556 outputContainer->Add(fhPhiIso) ;
1558 fhEtaIso = new TH2F("hEta",
1559 Form("Number of isolated particles vs #eta, %s",parTitle.Data()),
1560 nptbins,ptmin,ptmax,netabins,etamin,etamax);
1561 fhEtaIso->SetYTitle("#eta");
1562 fhEtaIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1563 outputContainer->Add(fhEtaIso) ;
1565 fhEtaPhiIso = new TH2F("hEtaPhiIso",
1566 Form("Number of isolated particles #eta vs #phi, %s",parTitle.Data()),
1567 netabins,etamin,etamax,nphibins,phimin,phimax);
1568 fhEtaPhiIso->SetXTitle("#eta");
1569 fhEtaPhiIso->SetYTitle("#phi");
1570 outputContainer->Add(fhEtaPhiIso) ;
1572 if(fFillHighMultHistograms)
1574 fhPtCentralityIso = new TH2F("hPtCentrality",
1575 Form("centrality vs #it{p}_{T} for isolated particles, %s",parTitle.Data()),
1576 nptbins,ptmin,ptmax, 100,0,100);
1577 fhPtCentralityIso->SetYTitle("centrality");
1578 fhPtCentralityIso->SetXTitle("#it{p}_{T}(GeV/#it{c})");
1579 outputContainer->Add(fhPtCentralityIso) ;
1581 fhPtEventPlaneIso = new TH2F("hPtEventPlane",
1582 Form("event plane angle vs #it{p}_{T} for isolated particles, %s",parTitle.Data()),
1583 nptbins,ptmin,ptmax, 100,0,TMath::Pi());
1584 fhPtEventPlaneIso->SetYTitle("Event plane angle (rad)");
1585 fhPtEventPlaneIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1586 outputContainer->Add(fhPtEventPlaneIso) ;
1589 if(fFillNLMHistograms)
1591 fhPtNLocMaxIso = new TH2F("hPtNLocMax",
1592 Form("Number of isolated particles vs #it{p}_{T}, %s",parTitle.Data()),
1593 nptbins,ptmin,ptmax,10,0,10);
1594 fhPtNLocMaxIso->SetYTitle("#it{NLM}");
1595 fhPtNLocMaxIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1597 fhPtNLocMaxNoIso = new TH2F("hPtNLocMaxNoIso",
1598 Form("Number of not isolated particles vs #it{p}_{T}, %s",parTitle.Data()),
1599 nptbins,ptmin,ptmax,10,0,10);
1600 fhPtNLocMaxNoIso->SetYTitle("#it{NLM}");
1601 fhPtNLocMaxNoIso->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1602 outputContainer->Add(fhPtNLocMaxNoIso) ;
1605 fhConeSumPt = new TH2F("hConePtSum",
1606 Form("Track and Cluster #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
1607 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1608 fhConeSumPt->SetYTitle("#Sigma #it{p}_{T}");
1609 fhConeSumPt->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1610 outputContainer->Add(fhConeSumPt) ;
1612 fhConeSumPtTrigEtaPhi = new TH2F("hConePtSumTrigEtaPhi",
1613 Form("Trigger #eta vs #phi, #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
1614 netabins,etamin,etamax,nphibins,phimin,phimax);
1615 fhConeSumPtTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1616 fhConeSumPtTrigEtaPhi->SetXTitle("#eta_{trigger}");
1617 fhConeSumPtTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1618 outputContainer->Add(fhConeSumPtTrigEtaPhi) ;
1620 fhPtInCone = new TH2F("hPtInCone",
1621 Form("#it{p}_{T} of clusters and tracks in isolation cone for #it{R} = %2.2f",r),
1622 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1623 fhPtInCone->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
1624 fhPtInCone->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1625 outputContainer->Add(fhPtInCone) ;
1627 if(fFillHighMultHistograms)
1629 fhPtInConeCent = new TH2F("hPtInConeCent",
1630 Form("#it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
1631 100,0,100,nptinconebins,ptinconemin,ptinconemax);
1632 fhPtInConeCent->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
1633 fhPtInConeCent->SetXTitle("centrality");
1634 outputContainer->Add(fhPtInConeCent) ;
1637 // Cluster only histograms
1638 if(GetIsolationCut()->GetParticleTypeInCone()!=AliIsolationCut::kOnlyCharged)
1640 fhConeSumPtCluster = new TH2F("hConePtSumCluster",
1641 Form("Cluster #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
1642 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1643 fhConeSumPtCluster->SetYTitle("#Sigma #it{p}_{T}");
1644 fhConeSumPtCluster->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1645 outputContainer->Add(fhConeSumPtCluster) ;
1647 if(fFillCellHistograms)
1649 fhConeSumPtCell = new TH2F("hConePtSumCell",
1650 Form("Cell #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
1651 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1652 fhConeSumPtCell->SetYTitle("#Sigma #it{p}_{T}");
1653 fhConeSumPtCell->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1654 outputContainer->Add(fhConeSumPtCell) ;
1657 if(fFillUEBandSubtractHistograms)
1659 fhConeSumPtEtaBandUECluster = new TH2F("hConePtSumEtaBandUECluster",
1660 "#Sigma cluster #it{p}_{T} in UE Eta Band",
1661 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1662 fhConeSumPtEtaBandUECluster->SetYTitle("#Sigma #it{p}_{T}");
1663 fhConeSumPtEtaBandUECluster->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1664 outputContainer->Add(fhConeSumPtEtaBandUECluster) ;
1666 fhConeSumPtPhiBandUECluster = new TH2F("hConePtSumPhiBandUECluster",
1667 "#Sigma cluster #it{p}_{T} UE Phi Band",
1668 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1669 fhConeSumPtPhiBandUECluster->SetYTitle("#Sigma #it{p}_{T}");
1670 fhConeSumPtPhiBandUECluster->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1671 outputContainer->Add(fhConeSumPtPhiBandUECluster) ;
1673 fhConeSumPtEtaBandUEClusterTrigEtaPhi = new TH2F("hConePtSumEtaBandUEClusterTrigEtaPhi",
1674 "Trigger #eta vs #phi, #Sigma cluster #it{p}_{T} in UE Eta Band",
1675 netabins,etamin,etamax,nphibins,phimin,phimax);
1676 fhConeSumPtEtaBandUEClusterTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1677 fhConeSumPtEtaBandUEClusterTrigEtaPhi->SetXTitle("#eta_{trigger}");
1678 fhConeSumPtEtaBandUEClusterTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1679 outputContainer->Add(fhConeSumPtEtaBandUEClusterTrigEtaPhi) ;
1681 fhConeSumPtPhiBandUEClusterTrigEtaPhi = new TH2F("hConePtSumPhiBandUEClusterTrigEtaPhi",
1682 "Trigger #eta vs #phi, #Sigma cluster #it{p}_{T} UE Phi Band",
1683 netabins,etamin,etamax,nphibins,phimin,phimax);
1684 fhConeSumPtPhiBandUEClusterTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1685 fhConeSumPtPhiBandUEClusterTrigEtaPhi->SetXTitle("#eta_{trigger}");
1686 fhConeSumPtPhiBandUEClusterTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1687 outputContainer->Add(fhConeSumPtPhiBandUEClusterTrigEtaPhi) ;
1688 if(fFillCellHistograms)
1691 fhConeSumPtEtaBandUECell = new TH2F("hConePtSumEtaBandUECell",
1692 "#Sigma cell #it{p}_{T} in UE Eta Band",
1693 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1694 fhConeSumPtEtaBandUECell->SetYTitle("#Sigma #it{p}_{T}");
1695 fhConeSumPtEtaBandUECell->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1696 outputContainer->Add(fhConeSumPtEtaBandUECell) ;
1698 fhConeSumPtPhiBandUECell = new TH2F("hConePtSumPhiBandUECell",
1699 "#Sigma cell #it{p}_{T} UE Phi Band",
1700 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1701 fhConeSumPtPhiBandUECell->SetYTitle("#Sigma #it{p}_{T}");
1702 fhConeSumPtPhiBandUECell->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1703 outputContainer->Add(fhConeSumPtPhiBandUECell) ;
1705 fhConeSumPtEtaBandUECellTrigEtaPhi = new TH2F("hConePtSumEtaBandUECellTrigEtaPhi",
1706 "Trigger #eta vs #phi, #Sigma cell #it{p}_{T} in UE Eta Band",
1707 netabins,etamin,etamax,nphibins,phimin,phimax);
1708 fhConeSumPtEtaBandUECellTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1709 fhConeSumPtEtaBandUECellTrigEtaPhi->SetXTitle("#eta_{trigger}");
1710 fhConeSumPtEtaBandUECellTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1711 outputContainer->Add(fhConeSumPtEtaBandUECellTrigEtaPhi) ;
1713 fhConeSumPtPhiBandUECellTrigEtaPhi = new TH2F("hConePtSumPhiBandUECellTrigEtaPhi",
1714 "Trigger #eta vs #phi, #Sigma cell #it{p}_{T} UE Phi Band",
1715 netabins,etamin,etamax,nphibins,phimin,phimax);
1716 fhConeSumPtPhiBandUECellTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1717 fhConeSumPtPhiBandUECellTrigEtaPhi->SetXTitle("#eta_{trigger}");
1718 fhConeSumPtPhiBandUECellTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1719 outputContainer->Add(fhConeSumPtPhiBandUECellTrigEtaPhi) ;
1722 fhEtaBandCluster = new TH2F("hEtaBandCluster",
1723 Form("#eta vs #phi of clusters in #eta band isolation cone for #it{R} = %2.2f",r),
1724 netabins,-1,1,nphibins,0,TMath::TwoPi());
1725 fhEtaBandCluster->SetXTitle("#eta");
1726 fhEtaBandCluster->SetYTitle("#phi");
1727 outputContainer->Add(fhEtaBandCluster) ;
1729 fhPhiBandCluster = new TH2F("hPhiBandCluster",
1730 Form("#eta vs #phi of clusters in #phi band isolation cone for #it{R} = %2.2f",r),
1731 netabins,-1,1,nphibins,0,TMath::TwoPi());
1732 fhPhiBandCluster->SetXTitle("#eta");
1733 fhPhiBandCluster->SetYTitle("#phi");
1734 outputContainer->Add(fhPhiBandCluster) ;
1736 fhEtaPhiInConeCluster= new TH2F("hEtaPhiInConeCluster",
1737 Form("#eta vs #phi of clusters in cone for #it{R} = %2.2f",r),
1738 netabins,-1,1,nphibins,0,TMath::TwoPi());
1739 fhEtaPhiInConeCluster->SetXTitle("#eta");
1740 fhEtaPhiInConeCluster->SetYTitle("#phi");
1741 outputContainer->Add(fhEtaPhiInConeCluster) ;
1743 fhEtaPhiCluster= new TH2F("hEtaPhiCluster",
1744 Form("#eta vs #phi of all clusters"),
1745 netabins,-1,1,nphibins,0,TMath::TwoPi());
1746 fhEtaPhiCluster->SetXTitle("#eta");
1747 fhEtaPhiCluster->SetYTitle("#phi");
1748 outputContainer->Add(fhEtaPhiCluster) ;
1752 fhPtClusterInCone = new TH2F("hPtClusterInCone",
1753 Form("#it{p}_{T} of clusters in isolation cone for #it{R} = %2.2f",r),
1754 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1755 fhPtClusterInCone->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
1756 fhPtClusterInCone->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1757 outputContainer->Add(fhPtClusterInCone) ;
1759 if(fFillCellHistograms)
1761 fhPtCellInCone = new TH2F("hPtCellInCone",
1762 Form("#it{p}_{T} of cells in isolation cone for #it{R} = %2.2f",r),
1763 nptbins,ptmin,ptmax,1000,0,50);
1764 fhPtCellInCone->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
1765 fhPtCellInCone->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1766 outputContainer->Add(fhPtCellInCone) ;
1768 fhEtaBandCell = new TH2F("hEtaBandCell",
1769 Form("#col vs #row of cells in #eta band isolation cone for #it{R} = %2.2f",r),
1771 fhEtaBandCell->SetXTitle("#col");
1772 fhEtaBandCell->SetYTitle("#row");
1773 outputContainer->Add(fhEtaBandCell) ;
1775 fhPhiBandCell = new TH2F("hPhiBandCell",
1776 Form("#col vs #row of cells in #phi band isolation cone for #it{R} = %2.2f",r),
1778 fhPhiBandCell->SetXTitle("#col");
1779 fhPhiBandCell->SetYTitle("#row");
1780 outputContainer->Add(fhPhiBandCell) ;
1783 if(fFillUEBandSubtractHistograms)
1785 fhConeSumPtEtaUESubCluster = new TH2F("hConeSumPtEtaUESubCluster",
1786 Form("Clusters #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
1787 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
1788 fhConeSumPtEtaUESubCluster->SetYTitle("#Sigma #it{p}_{T}");
1789 fhConeSumPtEtaUESubCluster->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1790 outputContainer->Add(fhConeSumPtEtaUESubCluster) ;
1792 fhConeSumPtPhiUESubCluster = new TH2F("hConeSumPtPhiUESubCluster",
1793 Form("Clusters #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
1794 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
1795 fhConeSumPtPhiUESubCluster->SetYTitle("#Sigma #it{p}_{T}");
1796 fhConeSumPtPhiUESubCluster->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1797 outputContainer->Add(fhConeSumPtPhiUESubCluster) ;
1799 fhConeSumPtEtaUESubClusterTrigEtaPhi = new TH2F("hConeSumPtEtaUESubClusterTrigEtaPhi",
1800 Form("Trigger #eta vs #phi, Clusters #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
1801 netabins,etamin,etamax,nphibins,phimin,phimax);
1802 fhConeSumPtEtaUESubClusterTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1803 fhConeSumPtEtaUESubClusterTrigEtaPhi->SetXTitle("#eta_{trigger}");
1804 fhConeSumPtEtaUESubClusterTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1805 outputContainer->Add(fhConeSumPtEtaUESubClusterTrigEtaPhi) ;
1807 fhConeSumPtPhiUESubClusterTrigEtaPhi = new TH2F("hConeSumPtPhiUESubClusterTrigEtaPhi",
1808 Form("Trigger #eta vs #phi, Clusters #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
1809 netabins,etamin,etamax,nphibins,phimin,phimax);
1810 fhConeSumPtPhiUESubClusterTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1811 fhConeSumPtPhiUESubClusterTrigEtaPhi->SetXTitle("#eta_{trigger}");
1812 fhConeSumPtPhiUESubClusterTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1813 outputContainer->Add(fhConeSumPtPhiUESubClusterTrigEtaPhi) ;
1815 if(fFillCellHistograms)
1817 fhConeSumPtEtaUESubCell = new TH2F("hConeSumPtEtaUESubCell",
1818 Form("Cells #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
1819 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
1820 fhConeSumPtEtaUESubCell->SetYTitle("#Sigma #it{p}_{T}");
1821 fhConeSumPtEtaUESubCell->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1822 outputContainer->Add(fhConeSumPtEtaUESubCell) ;
1824 fhConeSumPtPhiUESubCell = new TH2F("hConeSumPtPhiUESubCell",
1825 Form("Cells #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
1826 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
1827 fhConeSumPtPhiUESubCell->SetYTitle("#Sigma #it{p}_{T}");
1828 fhConeSumPtPhiUESubCell->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1829 outputContainer->Add(fhConeSumPtPhiUESubCell) ;
1831 fhConeSumPtEtaUESubCellTrigEtaPhi = new TH2F("hConeSumPtEtaUESubCellTrigEtaPhi",
1832 Form("Trigger #eta vs #phi, Cells #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
1833 netabins,etamin,etamax,nphibins,phimin,phimax);
1834 fhConeSumPtEtaUESubCellTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1835 fhConeSumPtEtaUESubCellTrigEtaPhi->SetXTitle("#eta_{trigger}");
1836 fhConeSumPtEtaUESubCellTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1837 outputContainer->Add(fhConeSumPtEtaUESubCellTrigEtaPhi) ;
1839 fhConeSumPtPhiUESubCellTrigEtaPhi = new TH2F("hConeSumPtPhiUESubCellTrigEtaPhi",
1840 Form("Trigger #eta vs #phi, Cells #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
1841 netabins,etamin,etamax,nphibins,phimin,phimax);
1842 fhConeSumPtPhiUESubCellTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
1843 fhConeSumPtPhiUESubCellTrigEtaPhi->SetXTitle("#eta_{trigger}");
1844 fhConeSumPtPhiUESubCellTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1845 outputContainer->Add(fhConeSumPtPhiUESubCellTrigEtaPhi) ;
1848 fhFractionClusterOutConeEta = new TH2F("hFractionClusterOutConeEta",
1849 Form("Fraction of the isolation cone #it{R} = %2.2f, out of clusters #eta acceptance",r),
1850 nptbins,ptmin,ptmax,100,0,1);
1851 fhFractionClusterOutConeEta->SetYTitle("#it{fraction}");
1852 fhFractionClusterOutConeEta->SetXTitle("#it{p}_{T,trigger} (GeV/#it{c})");
1853 outputContainer->Add(fhFractionClusterOutConeEta) ;
1855 fhFractionClusterOutConeEtaTrigEtaPhi = new TH2F("hFractionClusterOutConeEtaTrigEtaPhi",
1856 Form("Fraction of the isolation cone #it{R} = %2.2f, out of clusters #eta acceptance, in trigger #eta-#phi ",r),
1857 netabins,etamin,etamax,nphibins,phimin,phimax);
1858 fhFractionClusterOutConeEtaTrigEtaPhi->SetZTitle("#it{fraction}");
1859 fhFractionClusterOutConeEtaTrigEtaPhi->SetXTitle("#eta_{trigger}");
1860 fhFractionClusterOutConeEtaTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1861 outputContainer->Add(fhFractionClusterOutConeEtaTrigEtaPhi) ;
1863 fhFractionClusterOutConePhi = new TH2F("hFractionClusterOutConePhi",
1864 Form("Fraction of the isolation cone #it{R} = %2.2f, out of clusters #phi acceptance",r),
1865 nptbins,ptmin,ptmax,100,0,1);
1866 fhFractionClusterOutConePhi->SetYTitle("#it{fraction}");
1867 fhFractionClusterOutConePhi->SetXTitle("#it{p}_{T,trigger} (GeV/#it{c})");
1868 outputContainer->Add(fhFractionClusterOutConePhi) ;
1870 fhFractionClusterOutConePhiTrigEtaPhi = new TH2F("hFractionClusterOutConePhiTrigEtaPhi",
1871 Form("Fraction of the isolation cone #it{R} = %2.2f, out of clusters #phi acceptance, in trigger #eta-#phi ",r),
1872 netabins,etamin,etamax,nphibins,phimin,phimax);
1873 fhFractionClusterOutConePhiTrigEtaPhi->SetZTitle("#it{fraction}");
1874 fhFractionClusterOutConePhiTrigEtaPhi->SetXTitle("#eta_{trigger}");
1875 fhFractionClusterOutConePhiTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1876 outputContainer->Add(fhFractionClusterOutConePhiTrigEtaPhi) ;
1878 fhConeSumPtSubvsConeSumPtTotPhiCluster = new TH2F("hConeSumPtSubvsConeSumPtTotPhiCluster",
1879 Form("#Sigma #it{p}_{T} in cone after bkg sub from #phi band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1880 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1881 fhConeSumPtSubvsConeSumPtTotPhiCluster->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1882 fhConeSumPtSubvsConeSumPtTotPhiCluster->SetYTitle("#Sigma #it{p}_{T, sub} (GeV/#it{c})");
1883 outputContainer->Add(fhConeSumPtSubvsConeSumPtTotPhiCluster);
1885 fhConeSumPtSubNormvsConeSumPtTotPhiCluster = new TH2F("hConeSumPtSubNormvsConeSumPtTotPhiCluster",
1886 Form("#Sigma #it{p}_{T, norm} in cone after bkg sub from #phi band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1887 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1888 fhConeSumPtSubNormvsConeSumPtTotPhiCluster->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1889 fhConeSumPtSubNormvsConeSumPtTotPhiCluster->SetYTitle("#Sigma #it{p}_{T, sub norm} (GeV/#it{c})");
1890 outputContainer->Add(fhConeSumPtSubNormvsConeSumPtTotPhiCluster);
1892 fhConeSumPtSubvsConeSumPtTotEtaCluster = new TH2F("hConeSumPtSubvsConeSumPtTotEtaCluster",
1893 Form("#Sigma #it{p}_{T} in cone after bkg sub from #eta band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1894 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1895 fhConeSumPtSubvsConeSumPtTotEtaCluster->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1896 fhConeSumPtSubvsConeSumPtTotEtaCluster->SetYTitle("#Sigma #it{p}_{T, sub} (GeV/#it{c})");
1897 outputContainer->Add(fhConeSumPtSubvsConeSumPtTotEtaCluster);
1899 fhConeSumPtSubNormvsConeSumPtTotEtaCluster = new TH2F("hConeSumPtSubNormvsConeSumPtTotEtaCluster",
1900 Form("#Sigma #it{p}_{T, norm} in cone after bkg sub from #eta band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1901 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1902 fhConeSumPtSubNormvsConeSumPtTotEtaCluster->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1903 fhConeSumPtSubNormvsConeSumPtTotEtaCluster->SetYTitle("#Sigma #it{p}_{T, sub norm} (GeV/#it{c})");
1904 outputContainer->Add(fhConeSumPtSubNormvsConeSumPtTotEtaCluster);
1906 fhConeSumPtVSUEClusterEtaBand = new TH2F("hConeSumPtVSUEClusterEtaBand",
1907 Form("#Sigma #it{p}_{T} in cone versus #Sigma #it{p}_{T} in eta band for cluster (before normalization), R=%2.2f",r),
1908 nptsumbins,ptsummin,ptsummax,2*nptsumbins,ptsummin,2*ptsummax);
1909 fhConeSumPtVSUEClusterEtaBand->SetXTitle("#Sigma #it{p}_{T} cone (GeV/#it{c})");
1910 fhConeSumPtVSUEClusterEtaBand->SetYTitle("#Sigma #it{p}_{T} UE (GeV/#it{c})");
1911 outputContainer->Add(fhConeSumPtVSUEClusterEtaBand);
1913 fhConeSumPtVSUEClusterPhiBand = new TH2F("hConeSumPtVSUEClusterPhiBand",
1914 Form("#Sigma #it{p}_{T} in cone versus #Sigma #it{p}_{T} in phi band for cluster (before normalization), R=%2.2f",r),
1915 nptsumbins,ptsummin,ptsummax,8*nptsumbins,ptsummin,8*ptsummax);
1916 fhConeSumPtVSUEClusterPhiBand->SetXTitle("#Sigma #it{p}_{T} cone (GeV/#it{c})");
1917 fhConeSumPtVSUEClusterPhiBand->SetYTitle("#Sigma #it{p}_{T} UE (GeV/#it{c})");
1918 outputContainer->Add(fhConeSumPtVSUEClusterPhiBand);
1920 if(fFillCellHistograms)
1922 fhFractionCellOutConeEta = new TH2F("hFractionCellOutConeEta",
1923 Form("Fraction of the isolation cone #it{R} = %2.2f, out of cells #eta acceptance",r),
1924 nptbins,ptmin,ptmax,100,0,1);
1925 fhFractionCellOutConeEta->SetYTitle("#it{fraction}");
1926 fhFractionCellOutConeEta->SetXTitle("#it{p}_{T,trigger} (GeV/#it{c})");
1927 outputContainer->Add(fhFractionCellOutConeEta) ;
1929 fhFractionCellOutConeEtaTrigEtaPhi = new TH2F("hFractionCellOutConeEtaTrigEtaPhi",
1930 Form("Fraction of the isolation cone #it{R} = %2.2f, out of cells #eta acceptance, in trigger #eta-#phi ",r),
1931 netabins,etamin,etamax,nphibins,phimin,phimax);
1932 fhFractionCellOutConeEtaTrigEtaPhi->SetZTitle("#it{fraction}");
1933 fhFractionCellOutConeEtaTrigEtaPhi->SetXTitle("#eta_{trigger}");
1934 fhFractionCellOutConeEtaTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1935 outputContainer->Add(fhFractionCellOutConeEtaTrigEtaPhi) ;
1937 fhFractionCellOutConePhi = new TH2F("hFractionCellOutConePhi",
1938 Form("Fraction of the isolation cone #it{R} = %2.2f, out of cells #phi acceptance",r),
1939 nptbins,ptmin,ptmax,100,0,1);
1940 fhFractionCellOutConePhi->SetYTitle("#it{fraction}");
1941 fhFractionCellOutConePhi->SetXTitle("#it{p}_{T,trigger} (GeV/#it{c})");
1942 outputContainer->Add(fhFractionCellOutConePhi) ;
1944 fhFractionCellOutConePhiTrigEtaPhi = new TH2F("hFractionCellOutConePhiTrigEtaPhi",
1945 Form("Fraction of the isolation cone #it{R} = %2.2f, out of cells #phi acceptance, in trigger #eta-#phi ",r),
1946 netabins,etamin,etamax,nphibins,phimin,phimax);
1947 fhFractionCellOutConePhiTrigEtaPhi->SetZTitle("#it{fraction}");
1948 fhFractionCellOutConePhiTrigEtaPhi->SetXTitle("#eta_{trigger}");
1949 fhFractionCellOutConePhiTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
1950 outputContainer->Add(fhFractionCellOutConePhiTrigEtaPhi) ;
1953 fhConeSumPtSubvsConeSumPtTotPhiCell = new TH2F("hConeSumPtSubvsConeSumPtTotPhiCell",
1954 Form("#Sigma #it{p}_{T} in cone after bkg sub from #phi band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1955 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1956 fhConeSumPtSubvsConeSumPtTotPhiCell->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1957 fhConeSumPtSubvsConeSumPtTotPhiCell->SetYTitle("#Sigma #it{p}_{T, sub} (GeV/#it{c})");
1958 outputContainer->Add(fhConeSumPtSubvsConeSumPtTotPhiCell);
1960 fhConeSumPtSubNormvsConeSumPtTotPhiCell = new TH2F("hConeSumPtSubNormvsConeSumPtTotPhiCell",
1961 Form("#Sigma #it{p}_{T, norm} in cone after bkg sub from #phi band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1962 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1963 fhConeSumPtSubNormvsConeSumPtTotPhiCell->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1964 fhConeSumPtSubNormvsConeSumPtTotPhiCell->SetYTitle("#Sigma #it{p}_{T, sub norm} (GeV/#it{c})");
1965 outputContainer->Add(fhConeSumPtSubNormvsConeSumPtTotPhiCell);
1967 fhConeSumPtSubvsConeSumPtTotEtaCell = new TH2F("hConeSumPtSubvsConeSumPtTotEtaCell",
1968 Form("#Sigma #it{p}_{T} in cone after bkg sub from #eta band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1969 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1970 fhConeSumPtSubvsConeSumPtTotEtaCell->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1971 fhConeSumPtSubvsConeSumPtTotEtaCell->SetYTitle("#Sigma #it{p}_{T, sub} (GeV/#it{c})");
1972 outputContainer->Add(fhConeSumPtSubvsConeSumPtTotEtaCell);
1974 fhConeSumPtSubNormvsConeSumPtTotEtaCell = new TH2F("hConeSumPtSubNormvsConeSumPtTotEtaCell",
1975 Form("#Sigma #it{p}_{T, norm} in cone after bkg sub from #eta band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
1976 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
1977 fhConeSumPtSubNormvsConeSumPtTotEtaCell->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
1978 fhConeSumPtSubNormvsConeSumPtTotEtaCell->SetYTitle("#Sigma #it{p}_{T, sub norm} (GeV/#it{c})");
1979 outputContainer->Add(fhConeSumPtSubNormvsConeSumPtTotEtaCell);
1984 // Track only histograms
1985 if(GetIsolationCut()->GetParticleTypeInCone()!=AliIsolationCut::kOnlyNeutral)
1987 fhConeSumPtTrack = new TH2F("hConePtSumTrack",
1988 Form("Track #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
1989 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1990 fhConeSumPtTrack->SetYTitle("#Sigma #it{p}_{T}");
1991 fhConeSumPtTrack->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
1992 outputContainer->Add(fhConeSumPtTrack) ;
1994 fhPtTrackInCone = new TH2F("hPtTrackInCone",
1995 Form("#it{p}_{T} of tracks in isolation cone for #it{R} = %2.2f",r),
1996 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1997 fhPtTrackInCone->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
1998 fhPtTrackInCone->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1999 outputContainer->Add(fhPtTrackInCone) ;
2002 if(fFillUEBandSubtractHistograms)
2004 fhConeSumPtEtaBandUETrack = new TH2F("hConePtSumEtaBandUETrack",
2005 "#Sigma track #it{p}_{T} in UE Eta Band",
2006 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
2007 fhConeSumPtEtaBandUETrack->SetYTitle("#Sigma #it{p}_{T}");
2008 fhConeSumPtEtaBandUETrack->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
2009 outputContainer->Add(fhConeSumPtEtaBandUETrack) ;
2011 fhConeSumPtPhiBandUETrack = new TH2F("hConePtSumPhiBandUETrack",
2012 "#Sigma track #it{p}_{T} in UE Phi Band",
2013 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax*8);
2014 fhConeSumPtPhiBandUETrack->SetYTitle("#Sigma #it{p}_{T}");
2015 fhConeSumPtPhiBandUETrack->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
2016 outputContainer->Add(fhConeSumPtPhiBandUETrack) ;
2019 fhConeSumPtEtaBandUETrackTrigEtaPhi = new TH2F("hConePtSumEtaBandUETrackTrigEtaPhi",
2020 "Trigger #eta vs #phi, #Sigma track #it{p}_{T} in UE Eta Band",
2021 netabins,etamin,etamax,nphibins,phimin,phimax);
2022 fhConeSumPtEtaBandUETrackTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2023 fhConeSumPtEtaBandUETrackTrigEtaPhi->SetXTitle("#eta_{trigger}");
2024 fhConeSumPtEtaBandUETrackTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2025 outputContainer->Add(fhConeSumPtEtaBandUETrackTrigEtaPhi) ;
2027 fhConeSumPtPhiBandUETrackTrigEtaPhi = new TH2F("hConePtSumPhiBandUETrackTrigEtaPhi",
2028 "Trigger #eta vs #phi, #Sigma track #it{p}_{T} in UE Phi Band",
2029 netabins,etamin,etamax,nphibins,phimin,phimax);
2030 fhConeSumPtPhiBandUETrackTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2031 fhConeSumPtPhiBandUETrackTrigEtaPhi->SetXTitle("#eta_{trigger}");
2032 fhConeSumPtPhiBandUETrackTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2033 outputContainer->Add(fhConeSumPtPhiBandUETrackTrigEtaPhi) ;
2035 fhEtaBandTrack = new TH2F("hEtaBandTrack",
2036 Form("#eta vs #phi of tracks in #eta band isolation cone for #it{R} = %2.2f",r),
2037 netabins,-1,1,nphibins,0,TMath::TwoPi());
2038 fhEtaBandTrack->SetXTitle("#eta");
2039 fhEtaBandTrack->SetYTitle("#phi");
2040 outputContainer->Add(fhEtaBandTrack) ;
2042 fhPhiBandTrack = new TH2F("hPhiBandTrack",
2043 Form("#eta vs #phi of tracks in #phi band isolation cone for #it{R} = %2.2f",r),
2044 netabins,-1,1,nphibins,0,TMath::TwoPi());
2045 fhPhiBandTrack->SetXTitle("#eta");
2046 fhPhiBandTrack->SetYTitle("#phi");
2047 outputContainer->Add(fhPhiBandTrack) ;
2049 fhConeSumPtEtaUESubTrack = new TH2F("hConeSumPtEtaUESubTrack",
2050 Form("Tracks #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
2051 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
2052 fhConeSumPtEtaUESubTrack->SetYTitle("#Sigma #it{p}_{T}");
2053 fhConeSumPtEtaUESubTrack->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2054 outputContainer->Add(fhConeSumPtEtaUESubTrack) ;
2056 fhConeSumPtPhiUESubTrack = new TH2F("hConeSumPtPhiUESubTrack",
2057 Form("Tracks #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
2058 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
2059 fhConeSumPtPhiUESubTrack->SetYTitle("#Sigma #it{p}_{T}");
2060 fhConeSumPtPhiUESubTrack->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2061 outputContainer->Add(fhConeSumPtPhiUESubTrack) ;
2063 fhConeSumPtEtaUESubTrackTrigEtaPhi = new TH2F("hConeSumPtEtaUESubTrackTrigEtaPhi",
2064 Form("Trigger #eta vs #phi, Tracks #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
2065 netabins,etamin,etamax,nphibins,phimin,phimax);
2066 fhConeSumPtEtaUESubTrackTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2067 fhConeSumPtEtaUESubTrackTrigEtaPhi->SetXTitle("#eta_{trigger}");
2068 fhConeSumPtEtaUESubTrackTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2069 outputContainer->Add(fhConeSumPtEtaUESubTrackTrigEtaPhi) ;
2071 fhConeSumPtPhiUESubTrackTrigEtaPhi = new TH2F("hConeSumPtPhiUESubTrackTrigEtaPhi",
2072 Form("Trigger #eta vs #phi, Tracks #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
2073 netabins,etamin,etamax,nphibins,phimin,phimax);
2074 fhConeSumPtPhiUESubTrackTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2075 fhConeSumPtPhiUESubTrackTrigEtaPhi->SetXTitle("#eta_{trigger}");
2076 fhConeSumPtPhiUESubTrackTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2077 outputContainer->Add(fhConeSumPtPhiUESubTrackTrigEtaPhi) ;
2079 fhFractionTrackOutConeEta = new TH2F("hFractionTrackOutConeEta",
2080 Form("Fraction of the isolation cone #it{R} = %2.2f, out of tracks #eta acceptance",r),
2081 nptbins,ptmin,ptmax,100,0,1);
2082 fhFractionTrackOutConeEta->SetYTitle("#it{fraction}");
2083 fhFractionTrackOutConeEta->SetXTitle("#it{p}_{T,trigger} (GeV/#it{c})");
2084 outputContainer->Add(fhFractionTrackOutConeEta) ;
2086 fhFractionTrackOutConeEtaTrigEtaPhi = new TH2F("hFractionTrackOutConeEtaTrigEtaPhi",
2087 Form("Fraction of the isolation cone #it{R} = %2.2f, out of tracks #eta acceptance, in trigger #eta-#phi ",r),
2088 netabins,etamin,etamax,nphibins,phimin,phimax);
2089 fhFractionTrackOutConeEtaTrigEtaPhi->SetZTitle("#it{fraction}");
2090 fhFractionTrackOutConeEtaTrigEtaPhi->SetXTitle("#eta_{trigger}");
2091 fhFractionTrackOutConeEtaTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2092 outputContainer->Add(fhFractionTrackOutConeEtaTrigEtaPhi) ;
2094 fhConeSumPtSubvsConeSumPtTotPhiTrack = new TH2F("hConeSumPtSubvsConeSumPtTotPhiTrack",
2095 Form("#Sigma #it{p}_{T} in cone after bkg sub from #phi band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
2096 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2097 fhConeSumPtSubvsConeSumPtTotPhiTrack->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
2098 fhConeSumPtSubvsConeSumPtTotPhiTrack->SetYTitle("#Sigma #it{p}_{T, sub} (GeV/#it{c})");
2099 outputContainer->Add(fhConeSumPtSubvsConeSumPtTotPhiTrack);
2101 fhConeSumPtSubNormvsConeSumPtTotPhiTrack = new TH2F("hConeSumPtSubNormvsConeSumPtTotPhiTrack",
2102 Form("#Sigma #it{p}_{T, norm} in cone after bkg sub from #phi band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
2103 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2104 fhConeSumPtSubNormvsConeSumPtTotPhiTrack->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
2105 fhConeSumPtSubNormvsConeSumPtTotPhiTrack->SetYTitle("#Sigma #it{p}_{T, sub norm} (GeV/#it{c})");
2106 outputContainer->Add(fhConeSumPtSubNormvsConeSumPtTotPhiTrack);
2108 fhConeSumPtSubvsConeSumPtTotEtaTrack = new TH2F("hConeSumPtSubvsConeSumPtTotEtaTrack",
2109 Form("#Sigma #it{p}_{T} in cone after bkg sub from #eta band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
2110 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2111 fhConeSumPtSubvsConeSumPtTotEtaTrack->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
2112 fhConeSumPtSubvsConeSumPtTotEtaTrack->SetYTitle("#Sigma #it{p}_{T, sub} (GeV/#it{c})");
2113 outputContainer->Add(fhConeSumPtSubvsConeSumPtTotEtaTrack);
2115 fhConeSumPtSubNormvsConeSumPtTotEtaTrack = new TH2F("hConeSumPtSubNormvsConeSumPtTotEtaTrack",
2116 Form("#Sigma #it{p}_{T, norm} in cone after bkg sub from #eta band vs #Sigma #it{p}_{T} in cone before bkg sub, R=%2.2f",r),
2117 nptsumbins,ptsummin,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2118 fhConeSumPtSubNormvsConeSumPtTotEtaTrack->SetXTitle("#Sigma #it{p}_{T, tot} (GeV/#it{c})");
2119 fhConeSumPtSubNormvsConeSumPtTotEtaTrack->SetYTitle("#Sigma #it{p}_{T, sub norm} (GeV/#it{c})");
2120 outputContainer->Add(fhConeSumPtSubNormvsConeSumPtTotEtaTrack);
2123 // UE in perpendicular cone
2124 fhPerpConeSumPt = new TH2F("hPerpConePtSum",
2125 Form("#Sigma #it{p}_{T} in isolation cone at #pm 45 degree phi from trigger particle, #it{R} = %2.2f",r),
2126 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
2127 fhPerpConeSumPt->SetYTitle("#Sigma #it{p}_{T}");
2128 fhPerpConeSumPt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2129 outputContainer->Add(fhPerpConeSumPt) ;
2131 fhPtInPerpCone = new TH2F("hPtInPerpCone",
2132 Form("#it{p}_{T} in isolation cone at #pm 45 degree phi from trigger particle, #it{R} = %2.2f",r),
2133 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2134 fhPtInPerpCone->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
2135 fhPtInPerpCone->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2136 outputContainer->Add(fhPtInPerpCone) ;
2138 fhEtaPhiTrack= new TH2F("hEtaPhiTrack",
2139 Form("#eta vs #phi of all Tracks"),
2140 netabins,-1,1,nphibins,0,TMath::TwoPi());
2141 fhEtaPhiTrack->SetXTitle("#eta");
2142 fhEtaPhiTrack->SetYTitle("#phi");
2143 outputContainer->Add(fhEtaPhiTrack) ;
2145 fhEtaPhiInConeTrack= new TH2F("hEtaPhiInConeTrack",
2146 Form("#eta vs #phi of Tracks in cone for #it{R} = %2.2f",r),
2147 netabins,-1,1,nphibins,0,TMath::TwoPi());
2148 fhEtaPhiInConeTrack->SetXTitle("#eta");
2149 fhEtaPhiInConeTrack->SetYTitle("#phi");
2150 outputContainer->Add(fhEtaPhiInConeTrack) ;
2152 fhConeSumPtVSUETracksEtaBand = new TH2F("hConeSumPtVSUETracksEtaBand",
2153 Form("#Sigma #it{p}_{T} in cone versus #Sigma #it{p}_{T} in eta band for tracks (before normalization), R=%2.2f",r),
2154 nptsumbins,ptsummin,ptsummax,2*nptsumbins,ptsummin,2*ptsummax);
2155 fhConeSumPtVSUETracksEtaBand->SetXTitle("#Sigma #it{p}_{T} cone (GeV/#it{c})");
2156 fhConeSumPtVSUETracksEtaBand->SetYTitle("#Sigma #it{p}_{T} UE (GeV/#it{c})");
2157 outputContainer->Add(fhConeSumPtVSUETracksEtaBand);
2159 fhConeSumPtVSUETracksPhiBand = new TH2F("hConeSumPtVSUETracksPhiBand",
2160 Form("#Sigma #it{p}_{T} in cone versus #Sigma #it{p}_{T} in phi band for tracks (before normalization), R=%2.2f",r),
2161 nptsumbins,ptsummin,ptsummax,8*nptsumbins,ptsummin,8*ptsummax);
2162 fhConeSumPtVSUETracksPhiBand->SetXTitle("#Sigma #it{p}_{T} cone (GeV/#it{c})");
2163 fhConeSumPtVSUETracksPhiBand->SetYTitle("#Sigma #it{p}_{T} UE (GeV/#it{c})");
2164 outputContainer->Add(fhConeSumPtVSUETracksPhiBand);
2168 if(GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kNeutralAndCharged )
2170 fhConeSumPtClustervsTrack = new TH2F("hConePtSumClustervsTrack",
2171 Form("Track vs Cluster #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
2172 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2173 fhConeSumPtClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2174 fhConeSumPtClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2175 outputContainer->Add(fhConeSumPtClustervsTrack) ;
2177 if(fFillCellHistograms)
2179 fhConeSumPtCellvsTrack = new TH2F("hConePtSumCellvsTrack",
2180 Form("Track vs cell #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
2181 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2182 fhConeSumPtCellvsTrack->SetXTitle("#Sigma #it{p}_{T} cell");
2183 fhConeSumPtCellvsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2184 outputContainer->Add(fhConeSumPtCellvsTrack) ;
2186 fhConeSumPtCellTrack = new TH2F("hConePtSumCellTrack",
2187 Form("Track and Cell #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
2188 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
2189 fhConeSumPtCellTrack->SetYTitle("#Sigma #it{p}_{T}");
2190 fhConeSumPtCellTrack->SetXTitle("#it{p}_{T, trigger} (GeV/#it{c})");
2191 outputContainer->Add(fhConeSumPtCellTrack) ;
2193 fhConeSumPtCellTrackTrigEtaPhi = new TH2F("hConePtSumCellTrackTrigEtaPhi",
2194 Form("Trigger #eta vs #phi, #Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",r),
2195 netabins,etamin,etamax,nphibins,phimin,phimax);
2196 fhConeSumPtCellTrackTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2197 fhConeSumPtCellTrackTrigEtaPhi->SetXTitle("#eta_{trigger}");
2198 fhConeSumPtCellTrackTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2199 outputContainer->Add(fhConeSumPtCellTrackTrigEtaPhi) ;
2203 if(fFillUEBandSubtractHistograms)
2205 fhConeSumPtEtaUESub = new TH2F("hConeSumPtEtaUESub",
2206 Form("#Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
2207 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
2208 fhConeSumPtEtaUESub->SetYTitle("#Sigma #it{p}_{T}");
2209 fhConeSumPtEtaUESub->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2210 outputContainer->Add(fhConeSumPtEtaUESub) ;
2212 fhConeSumPtPhiUESub = new TH2F("hConeSumPtPhiUESub",
2213 Form("#Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
2214 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
2215 fhConeSumPtPhiUESub->SetYTitle("#Sigma #it{p}_{T}");
2216 fhConeSumPtPhiUESub->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2217 outputContainer->Add(fhConeSumPtPhiUESub) ;
2219 fhConeSumPtEtaUESubTrigEtaPhi = new TH2F("hConeSumPtEtaUESubTrigEtaPhi",
2220 Form("Trigger #eta vs #phi, #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
2221 netabins,etamin,etamax,nphibins,phimin,phimax);
2222 fhConeSumPtEtaUESubTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2223 fhConeSumPtEtaUESubTrigEtaPhi->SetXTitle("#eta_{trigger}");
2224 fhConeSumPtEtaUESubTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2225 outputContainer->Add(fhConeSumPtEtaUESubTrigEtaPhi) ;
2227 fhConeSumPtPhiUESubTrigEtaPhi = new TH2F("hConeSumPtPhiUESubTrigEtaPhi",
2228 Form("Trigger #eta vs #phi, #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
2229 netabins,etamin,etamax,nphibins,phimin,phimax);
2230 fhConeSumPtPhiUESubTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2231 fhConeSumPtPhiUESubTrigEtaPhi->SetXTitle("#eta_{trigger}");
2232 fhConeSumPtPhiUESubTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2233 outputContainer->Add(fhConeSumPtPhiUESubTrigEtaPhi) ;
2235 fhConeSumPtEtaUESubClustervsTrack = new TH2F("hConePtSumEtaUESubClustervsTrack",
2236 Form("Track vs Cluster #Sigma #it{p}_{T} UE sub eta band in isolation cone for #it{R} = %2.2f",r),
2237 2*nptsumbins,-ptsummax,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2238 fhConeSumPtEtaUESubClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2239 fhConeSumPtEtaUESubClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2240 outputContainer->Add(fhConeSumPtEtaUESubClustervsTrack) ;
2242 fhConeSumPtPhiUESubClustervsTrack = new TH2F("hConePhiUESubPtSumClustervsTrack",
2243 Form("Track vs Cluster #Sigma #it{p}_{T} UE sub phi band in isolation cone for #it{R} = %2.2f",r),
2244 2*nptsumbins,-ptsummax,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2245 fhConeSumPtPhiUESubClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2246 fhConeSumPtPhiUESubClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2247 outputContainer->Add(fhConeSumPtPhiUESubClustervsTrack) ;
2249 fhEtaBandClustervsTrack = new TH2F("hEtaBandClustervsTrack",
2250 Form("Track vs Cluster #Sigma #it{p}_{T} in Eta band in isolation cone for #it{R} = %2.2f",r),
2251 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2252 fhEtaBandClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2253 fhEtaBandClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2254 outputContainer->Add(fhEtaBandClustervsTrack) ;
2256 fhPhiBandClustervsTrack = new TH2F("hPhiBandClustervsTrack",
2257 Form("Track vs Cluster #Sigma #it{p}_{T} in Phi band in isolation cone for #it{R} = %2.2f",r),
2258 nptsumbins,ptsummin,ptsummax*4,nptsumbins,ptsummin,ptsummax*8);
2259 fhPhiBandClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2260 fhPhiBandClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2261 outputContainer->Add(fhPhiBandClustervsTrack) ;
2263 fhEtaBandNormClustervsTrack = new TH2F("hEtaBandNormClustervsTrack",
2264 Form("Track vs Cluster #Sigma #it{p}_{T} in Eta band in isolation cone for #it{R} = %2.2f",r),
2265 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2266 fhEtaBandNormClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2267 fhEtaBandNormClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2268 outputContainer->Add(fhEtaBandNormClustervsTrack) ;
2270 fhPhiBandNormClustervsTrack = new TH2F("hPhiBandNormClustervsTrack",
2271 Form("Track vs Cluster #Sigma #it{p}_{T} in Phi band in isolation cone for #it{R} = %2.2f",r),
2272 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2273 fhPhiBandNormClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2274 fhPhiBandNormClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2275 outputContainer->Add(fhPhiBandNormClustervsTrack) ;
2277 fhConeSumPtEtaUESubClustervsTrack = new TH2F("hConePtSumEtaUESubClustervsTrack",
2278 Form("Track vs Cluster #Sigma #it{p}_{T} UE sub eta band in isolation cone for #it{R} = %2.2f",r),
2279 2*nptsumbins,-ptsummax,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2280 fhConeSumPtEtaUESubClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2281 fhConeSumPtEtaUESubClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2282 outputContainer->Add(fhConeSumPtEtaUESubClustervsTrack) ;
2284 fhConeSumPtPhiUESubClustervsTrack = new TH2F("hConePhiUESubPtSumClustervsTrack",
2285 Form("Track vs Cluster #Sigma #it{p}_{T} UE sub phi band in isolation cone for #it{R} = %2.2f",r),
2286 2*nptsumbins,-ptsummax,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2287 fhConeSumPtPhiUESubClustervsTrack->SetXTitle("#Sigma #it{p}_{T} cluster");
2288 fhConeSumPtPhiUESubClustervsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2289 outputContainer->Add(fhConeSumPtPhiUESubClustervsTrack) ;
2291 if(fFillCellHistograms)
2294 fhConeSumPtEtaUESubCellvsTrack = new TH2F("hConePtSumEtaUESubCellvsTrack",
2295 Form("Track vs Cell #Sigma #it{p}_{T} UE sub eta band in isolation cone for #it{R} = %2.2f",r),
2296 2*nptsumbins,-ptsummax,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2297 fhConeSumPtEtaUESubCellvsTrack->SetXTitle("#Sigma #it{p}_{T} cell");
2298 fhConeSumPtEtaUESubCellvsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2299 outputContainer->Add(fhConeSumPtEtaUESubCellvsTrack) ;
2301 fhConeSumPtPhiUESubCellvsTrack = new TH2F("hConePhiUESubPtSumCellvsTrack",
2302 Form("Track vs Cell #Sigma #it{p}_{T} UE sub phi band in isolation cone for #it{R} = %2.2f",r),
2303 2*nptsumbins,-ptsummax,ptsummax,2*nptsumbins,-ptsummax,ptsummax);
2304 fhConeSumPtPhiUESubCellvsTrack->SetXTitle("#Sigma #it{p}_{T} cell");
2305 fhConeSumPtPhiUESubCellvsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2306 outputContainer->Add(fhConeSumPtPhiUESubCellvsTrack) ;
2308 fhEtaBandCellvsTrack = new TH2F("hEtaBandCellvsTrack",
2309 Form("Track vs Cell #Sigma #it{p}_{T} in Eta band in isolation cone for #it{R} = %2.2f",r),
2310 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2311 fhEtaBandCellvsTrack->SetXTitle("#Sigma #it{p}_{T} cell");
2312 fhEtaBandCellvsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2313 outputContainer->Add(fhEtaBandCellvsTrack) ;
2315 fhPhiBandCellvsTrack = new TH2F("hPhiBandCellvsTrack",
2316 Form("Track vs Cell #Sigma #it{p}_{T} in Phi band in isolation cone for #it{R} = %2.2f",r),
2317 nptsumbins,ptsummin,ptsummax*4,nptsumbins,ptsummin,ptsummax*8);
2318 fhPhiBandCellvsTrack->SetXTitle("#Sigma #it{p}_{T} cell");
2319 fhPhiBandCellvsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2320 outputContainer->Add(fhPhiBandCellvsTrack) ;
2322 fhEtaBandNormCellvsTrack = new TH2F("hEtaBandNormCellvsTrack",
2323 Form("Track vs Cell #Sigma #it{p}_{T} in Eta band in isolation cone for #it{R} = %2.2f",r),
2324 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2325 fhEtaBandNormCellvsTrack->SetXTitle("#Sigma #it{p}_{T} cell");
2326 fhEtaBandNormCellvsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2327 outputContainer->Add(fhEtaBandNormCellvsTrack) ;
2329 fhPhiBandNormCellvsTrack = new TH2F("hPhiBandNormCellvsTrack",
2330 Form("Track vs Cell #Sigma #it{p}_{T} in Phi band in isolation cone for #it{R} = %2.2f",r),
2331 nptsumbins,ptsummin,ptsummax,nptsumbins,ptsummin,ptsummax);
2332 fhPhiBandNormCellvsTrack->SetXTitle("#Sigma #it{p}_{T} cell");
2333 fhPhiBandNormCellvsTrack->SetYTitle("#Sigma #it{p}_{T} track");
2334 outputContainer->Add(fhPhiBandNormCellvsTrack) ;
2336 fhConeSumPtEtaUESubTrackCell = new TH2F("hConeSumPtEtaUESubTrackCell",
2337 Form("Tracks #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
2338 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
2339 fhConeSumPtEtaUESubTrackCell->SetYTitle("#Sigma #it{p}_{T}");
2340 fhConeSumPtEtaUESubTrackCell->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2341 outputContainer->Add(fhConeSumPtEtaUESubTrackCell) ;
2343 fhConeSumPtPhiUESubTrackCell = new TH2F("hConeSumPtPhiUESubTrackCell",
2344 Form("Tracks #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
2345 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
2346 fhConeSumPtPhiUESubTrackCell->SetYTitle("#Sigma #it{p}_{T}");
2347 fhConeSumPtPhiUESubTrackCell->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2348 outputContainer->Add(fhConeSumPtPhiUESubTrackCell) ;
2350 fhConeSumPtEtaUESubTrackCellTrigEtaPhi = new TH2F("hConeSumPtEtaUESubTrackCellTrigEtaPhi",
2351 Form("Trigger #eta vs #phi, Tracks #Sigma #it{p}_{T} after bkg subtraction from eta band in the isolation cone for #it{R} = %2.2f",r),
2352 netabins,etamin,etamax,nphibins,phimin,phimax);
2353 fhConeSumPtEtaUESubTrackCellTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2354 fhConeSumPtEtaUESubTrackCellTrigEtaPhi->SetXTitle("#eta_{trigger}");
2355 fhConeSumPtEtaUESubTrackCellTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2356 outputContainer->Add(fhConeSumPtEtaUESubTrackCellTrigEtaPhi) ;
2358 fhConeSumPtPhiUESubTrackCellTrigEtaPhi = new TH2F("hConeSumPtPhiUESubTrackCellTrigEtaPhi",
2359 Form("Trigger #eta vs #phi, Tracks #Sigma #it{p}_{T} after bkg subtraction from phi band in the isolation cone for #it{R} = %2.2f",r),
2360 netabins,etamin,etamax,nphibins,phimin,phimax);
2361 fhConeSumPtPhiUESubTrackCellTrigEtaPhi->SetZTitle("#Sigma #it{p}_{T}");
2362 fhConeSumPtPhiUESubTrackCellTrigEtaPhi->SetXTitle("#eta_{trigger}");
2363 fhConeSumPtPhiUESubTrackCellTrigEtaPhi->SetYTitle("#phi_{trigger} (rad)");
2364 outputContainer->Add(fhConeSumPtPhiUESubTrackCellTrigEtaPhi) ;
2369 for(Int_t iso = 0; iso < 2; iso++)
2373 fhTrackMatchedDEta[iso] = new TH2F
2374 (Form("hTrackMatchedDEta%s",isoName[iso].Data()),
2375 Form("%s - d#eta of cluster-track vs cluster energy, %s",isoTitle[iso].Data(),parTitle.Data()),
2376 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
2377 fhTrackMatchedDEta[iso]->SetYTitle("d#eta");
2378 fhTrackMatchedDEta[iso]->SetXTitle("E_{cluster} (GeV)");
2380 fhTrackMatchedDPhi[iso] = new TH2F
2381 (Form("hTrackMatchedDPhi%s",isoName[iso].Data()),
2382 Form("%s - d#phi of cluster-track vs cluster energy, %s",isoTitle[iso].Data(),parTitle.Data()),
2383 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
2384 fhTrackMatchedDPhi[iso]->SetYTitle("d#phi (rad)");
2385 fhTrackMatchedDPhi[iso]->SetXTitle("E_{cluster} (GeV)");
2387 fhTrackMatchedDEtaDPhi[iso] = new TH2F
2388 (Form("hTrackMatchedDEtaDPhi%s",isoName[iso].Data()),
2389 Form("%s - d#eta vs d#phi of cluster-track, %s",isoTitle[iso].Data(),parTitle.Data()),
2390 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
2391 fhTrackMatchedDEtaDPhi[iso]->SetYTitle("d#phi (rad)");
2392 fhTrackMatchedDEtaDPhi[iso]->SetXTitle("d#eta");
2394 outputContainer->Add(fhTrackMatchedDEta[iso]) ;
2395 outputContainer->Add(fhTrackMatchedDPhi[iso]) ;
2396 outputContainer->Add(fhTrackMatchedDEtaDPhi[iso]) ;
2398 fhdEdx[iso] = new TH2F
2399 (Form("hdEdx%s",isoName[iso].Data()),
2400 Form("%s - Matched track <d#it{E}/d#it{x}> vs cluster #it{E}, %s",isoTitle[iso].Data(),parTitle.Data()),
2401 nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
2402 fhdEdx[iso]->SetXTitle("#it{E} (GeV)");
2403 fhdEdx[iso]->SetYTitle("<d#it{E}/d#it{x}>");
2404 outputContainer->Add(fhdEdx[iso]);
2406 fhEOverP[iso] = new TH2F
2407 (Form("hEOverP%s",isoName[iso].Data()),
2408 Form("%s - Matched track #it{E}/#it{p} vs cluster, %s",isoTitle[iso].Data(),parTitle.Data()),
2409 nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
2410 fhEOverP[iso]->SetXTitle("#it{E} (GeV)");
2411 fhEOverP[iso]->SetYTitle("#it{E}/#it{p}");
2412 outputContainer->Add(fhEOverP[iso]);
2416 fhTrackMatchedMCParticle[iso] = new TH2F
2417 (Form("hTrackMatchedMCParticle%s",isoName[iso].Data()),
2418 Form("%s - Origin of particle vs cluster #it{E}, %s",isoTitle[iso].Data(),parTitle.Data()),
2419 nptbins,ptmin,ptmax,8,0,8);
2420 fhTrackMatchedMCParticle[iso]->SetXTitle("#it{E} (GeV)");
2421 //fhTrackMatchedMCParticle[iso]->SetYTitle("Particle type");
2423 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(1 ,"Photon");
2424 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(2 ,"Electron");
2425 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
2426 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(4 ,"Rest");
2427 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
2428 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
2429 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
2430 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
2432 outputContainer->Add(fhTrackMatchedMCParticle[iso]);
2438 fhELambda0[iso] = new TH2F
2439 (Form("hELambda0%s",isoName[iso].Data()),
2440 Form("%s cluster : #it{E} vs #lambda_{0}, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2441 fhELambda0[iso]->SetYTitle("#lambda_{0}^{2}");
2442 fhELambda0[iso]->SetXTitle("#it{E} (GeV)");
2443 outputContainer->Add(fhELambda0[iso]) ;
2445 fhELambda1[iso] = new TH2F
2446 (Form("hELambda1%s",isoName[iso].Data()),
2447 Form("%s cluster: #it{E} vs #lambda_{1}, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2448 fhELambda1[iso]->SetYTitle("#lambda_{1}^{2}");
2449 fhELambda1[iso]->SetXTitle("#it{E} (GeV)");
2450 outputContainer->Add(fhELambda1[iso]) ;
2452 fhPtLambda0[iso] = new TH2F
2453 (Form("hPtLambda0%s",isoName[iso].Data()),
2454 Form("%s cluster : #it{p}_{T} vs #lambda_{0}, %s",isoTitle[iso].Data(), parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2455 fhPtLambda0[iso]->SetYTitle("#lambda_{0}^{2}");
2456 fhPtLambda0[iso]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2457 outputContainer->Add(fhPtLambda0[iso]) ;
2461 for(Int_t imc = 0; imc < 9; imc++)
2463 fhPtLambda0MC[imc][iso] = new TH2F(Form("hPtLambda0%s_MC%s",isoName[iso].Data(),mcPartName[imc].Data()),
2464 Form("%s cluster : #it{p}_{T} vs #lambda_{0}: %s %s",isoTitle[iso].Data(),mcPartType[imc].Data(),parTitle.Data()),
2465 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2466 fhPtLambda0MC[imc][iso]->SetYTitle("#lambda_{0}^{2}");
2467 fhPtLambda0MC[imc][iso]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
2468 outputContainer->Add( fhPtLambda0MC[imc][iso]) ;
2472 if(fIsoDetector=="EMCAL" && GetFirstSMCoveredByTRD() >= 0)
2474 fhPtLambda0TRD[iso] = new TH2F
2475 (Form("hPtLambda0TRD%s",isoName[iso].Data()),
2476 Form("%s cluster: #it{p}_{T} vs #lambda_{0}, SM behind TRD, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2477 fhPtLambda0TRD[iso]->SetYTitle("#lambda_{0}^{2}");
2478 fhPtLambda0TRD[iso]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2479 outputContainer->Add(fhPtLambda0TRD[iso]) ;
2481 fhELambda0TRD[iso] = new TH2F
2482 (Form("hELambda0TRD%s",isoName[iso].Data()),
2483 Form("%s cluster: #it{E} vs #lambda_{0}, SM behind TRD, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2484 fhELambda0TRD[iso]->SetYTitle("#lambda_{0}^{2}");
2485 fhELambda0TRD[iso]->SetXTitle("#it{E} (GeV)");
2486 outputContainer->Add(fhELambda0TRD[iso]) ;
2488 fhELambda1TRD[iso] = new TH2F
2489 (Form("hELambda1TRD%s",isoName[iso].Data()),
2490 Form("%s cluster: #it{E} vs #lambda_{1}, SM behind TRD, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2491 fhELambda1TRD[iso]->SetYTitle("#lambda_{1}^{2}");
2492 fhELambda1TRD[iso]->SetXTitle("#it{E} (GeV)");
2493 outputContainer->Add(fhELambda1TRD[iso]) ;
2496 if(fFillNLMHistograms)
2498 fhNLocMax[iso] = new TH2F
2499 (Form("hNLocMax%s",isoName[iso].Data()),
2500 Form("%s - Number of local maxima in cluster, %s",isoTitle[iso].Data(),parTitle.Data()),
2501 nptbins,ptmin,ptmax,10,0,10);
2502 fhNLocMax[iso]->SetYTitle("#it{NLM}");
2503 fhNLocMax[iso]->SetXTitle("#it{E} (GeV)");
2504 outputContainer->Add(fhNLocMax[iso]) ;
2506 fhELambda0LocMax1[iso] = new TH2F
2507 (Form("hELambda0LocMax1%s",isoName[iso].Data()),
2508 Form("%s cluster (#eta) pairs: #it{E} vs #lambda_{0}, #it{NLM}=1, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2509 fhELambda0LocMax1[iso]->SetYTitle("#lambda_{0}^{2}");
2510 fhELambda0LocMax1[iso]->SetXTitle("#it{E} (GeV)");
2511 outputContainer->Add(fhELambda0LocMax1[iso]) ;
2513 fhELambda1LocMax1[iso] = new TH2F
2514 (Form("hELambda1LocMax1%s",isoName[iso].Data()),
2515 Form("%s cluster (#eta) pairs: #it{E} vs #lambda_{1}, #it{NLM}=1, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2516 fhELambda1LocMax1[iso]->SetYTitle("#lambda_{1}^{2}");
2517 fhELambda1LocMax1[iso]->SetXTitle("#it{E} (GeV)");
2518 outputContainer->Add(fhELambda1LocMax1[iso]) ;
2520 fhELambda0LocMax2[iso] = new TH2F
2521 (Form("hELambda0LocMax2%s",isoName[iso].Data()),
2522 Form("%s cluster (#eta) pairs: #it{E} vs #lambda_{0}, #it{NLM}=2, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2523 fhELambda0LocMax2[iso]->SetYTitle("#lambda_{0}^{2}");
2524 fhELambda0LocMax2[iso]->SetXTitle("#it{E} (GeV)");
2525 outputContainer->Add(fhELambda0LocMax2[iso]) ;
2527 fhELambda1LocMax2[iso] = new TH2F
2528 (Form("hELambda1LocMax2%s",isoName[iso].Data()),
2529 Form("%s cluster (#eta) pairs: #it{E} vs #lambda_{1}, #it{NLM}=2, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2530 fhELambda1LocMax2[iso]->SetYTitle("#lambda_{1}^{2}");
2531 fhELambda1LocMax2[iso]->SetXTitle("#it{E} (GeV)");
2532 outputContainer->Add(fhELambda1LocMax2[iso]) ;
2534 fhELambda0LocMaxN[iso] = new TH2F
2535 ( Form("hELambda0LocMaxN%s",isoName[iso].Data()),
2536 Form("%s cluster (#eta) pairs: #it{E} vs #lambda_{0}, #it{NLM}>2, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2537 fhELambda0LocMaxN[iso]->SetYTitle("#lambda_{0}^{2}");
2538 fhELambda0LocMaxN[iso]->SetXTitle("#it{E} (GeV)");
2539 outputContainer->Add(fhELambda0LocMaxN[iso]) ;
2541 fhELambda1LocMaxN[iso] = new TH2F
2542 (Form("hELambda1LocMaxN%s",isoName[iso].Data()),
2543 Form("%s cluster (#eta) pairs: #it{E} vs #lambda_{1}, #it{NLM}>2, %s",isoTitle[iso].Data(),parTitle.Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
2544 fhELambda1LocMaxN[iso]->SetYTitle("#lambda_{1}^{2}");
2545 fhELambda1LocMaxN[iso]->SetXTitle("#it{E} (GeV)");
2546 outputContainer->Add(fhELambda1LocMaxN[iso]) ;
2549 } // control histograms for isolated and non isolated objects
2552 if(fFillPileUpHistograms)
2554 fhPtTrackInConeOtherBC = new TH2F("hPtTrackInConeOtherBC",
2555 Form("#it{p}_{T} of tracks in isolation cone for #it{R} = %2.2f, TOF from BC!=0",r),
2556 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2557 fhPtTrackInConeOtherBC->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
2558 fhPtTrackInConeOtherBC->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2559 outputContainer->Add(fhPtTrackInConeOtherBC) ;
2561 fhPtTrackInConeOtherBCPileUpSPD = new TH2F("hPtTrackInConeOtherBCPileUpSPD",
2562 Form("#it{p}_{T} of tracks in isolation cone for #it{R} = %2.2f, TOF from BC!=0, pile-up from SPD",r),
2563 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2564 fhPtTrackInConeOtherBCPileUpSPD->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
2565 fhPtTrackInConeOtherBCPileUpSPD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2566 outputContainer->Add(fhPtTrackInConeOtherBCPileUpSPD) ;
2568 fhPtTrackInConeBC0 = new TH2F("hPtTrackInConeBC0",
2569 Form("#it{p}_{T} of tracks in isolation cone for #it{R} = %2.2f, TOF from BC==0",r),
2570 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2571 fhPtTrackInConeBC0->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
2572 fhPtTrackInConeBC0->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2573 outputContainer->Add(fhPtTrackInConeBC0) ;
2575 fhPtTrackInConeVtxBC0 = new TH2F("hPtTrackInConeVtxBC0",
2576 Form("#it{p}_{T} of tracks in isolation cone for #it{R} = %2.2f, TOF from BC==0",r),
2577 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2578 fhPtTrackInConeVtxBC0->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
2579 fhPtTrackInConeVtxBC0->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2580 outputContainer->Add(fhPtTrackInConeVtxBC0) ;
2583 fhPtTrackInConeBC0PileUpSPD = new TH2F("hPtTrackInConeBC0PileUpSPD",
2584 Form("#it{p}_{T} of tracks in isolation cone for #it{R} = %2.2f, TOF from BC==0, pile-up from SPD",r),
2585 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2586 fhPtTrackInConeBC0PileUpSPD->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
2587 fhPtTrackInConeBC0PileUpSPD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2588 outputContainer->Add(fhPtTrackInConeBC0PileUpSPD) ;
2591 for (Int_t i = 0; i < 7 ; i++)
2593 fhPtInConePileUp[i] = new TH2F(Form("hPtInConePileUp%s",pileUpName[i].Data()),
2594 Form("#it{p}_{T} in isolation cone for #it{R} = %2.2f, from pile-up (%s)",r,pileUpName[i].Data()),
2595 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2596 fhPtInConePileUp[i]->SetYTitle("#it{p}_{T in cone} (GeV/#it{c})");
2597 fhPtInConePileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2598 outputContainer->Add(fhPtInConePileUp[i]) ;
2604 // For histograms in arrays, index in the array, corresponding to any particle origin
2606 for(Int_t i = 0; i < 6; i++)
2608 fhEPrimMC[i] = new TH1F(Form("hEPrim_MC%s",ppname[i].Data()),
2609 Form("primary photon %s : #it{E}, %s",pptype[i].Data(),parTitle.Data()),
2610 nptbins,ptmin,ptmax);
2611 fhEPrimMC[i]->SetXTitle("#it{E} (GeV)");
2612 outputContainer->Add(fhEPrimMC[i]) ;
2614 fhPtPrimMCiso[i] = new TH1F(Form("hPtPrim_MCiso%s",ppname[i].Data()),
2615 Form("primary isolated photon %s : #it{p}_{T}, %s",pptype[i].Data(),parTitle.Data()),
2616 nptbins,ptmin,ptmax);
2617 fhPtPrimMCiso[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2618 outputContainer->Add(fhPtPrimMCiso[i]) ;
2620 fhEtaPrimMC[i] = new TH2F(Form("hEtaPrim_MC%s",ppname[i].Data()),
2621 Form("primary photon %s : #eta vs #it{p}_{T}, %s",pptype[i].Data(),parTitle.Data()),
2622 nptbins,ptmin,ptmax,200,-2,2);
2623 fhEtaPrimMC[i]->SetYTitle("#eta");
2624 fhEtaPrimMC[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2625 outputContainer->Add(fhEtaPrimMC[i]) ;
2627 fhPhiPrimMC[i] = new TH2F(Form("hPhiPrim_MC%s",ppname[i].Data()),
2628 Form("primary photon %s : #phi vs #it{p}_{T}, %s",pptype[i].Data(),parTitle.Data()),
2629 nptbins,ptmin,ptmax,200,0.,TMath::TwoPi());
2630 fhPhiPrimMC[i]->SetYTitle("#phi");
2631 fhPhiPrimMC[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2632 outputContainer->Add(fhPhiPrimMC[i]) ;
2641 const Int_t buffersize = 255;
2642 char name[buffersize];
2643 char title[buffersize];
2644 for(Int_t icone = 0; icone<fNCones; icone++)
2646 // sum pt in cone vs. pt leading
2647 snprintf(name, buffersize,"hSumPtLeadingPt_Cone_%d",icone);
2648 snprintf(title, buffersize,"#Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",fConeSizes[icone]);
2649 fhSumPtLeadingPt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
2650 fhSumPtLeadingPt[icone] ->SetYTitle("#sum_{cone}#it{p}_{T} (GeV/#it{c})");//#Sigma #it{p}_{T}
2651 fhSumPtLeadingPt[icone] ->SetXTitle("#it{p}_{T}^{leading} (GeV/#it{c})");
2652 outputContainer->Add(fhSumPtLeadingPt[icone]) ;
2654 // pt in cone vs. pt leading
2655 snprintf(name, buffersize,"hPtLeadingPt_Cone_%d",icone);
2656 snprintf(title, buffersize,"#it{p}_{T} in isolation cone for #it{R} = %2.2f",fConeSizes[icone]);
2657 fhPtLeadingPt[icone] = new TH2F(name, title, nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2658 fhPtLeadingPt[icone] ->SetYTitle("#it{p}_{T}^{cone} (GeV/#it{c})");
2659 fhPtLeadingPt[icone] ->SetXTitle("#it{p}_{T}^{leading} (GeV/#it{c})");
2660 outputContainer->Add(fhPtLeadingPt[icone]) ;
2662 // sum pt in cone vs. pt leading in the forward region (for background subtraction studies)
2663 snprintf(name, buffersize,"hPerpSumPtLeadingPt_Cone_%d",icone);
2664 snprintf(title, buffersize,"#Sigma #it{p}_{T} in isolation cone for #it{R} = %2.2f",fConeSizes[icone]);
2665 fhPerpSumPtLeadingPt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
2666 fhPerpSumPtLeadingPt[icone] ->SetYTitle("#sum_{cone}#it{p}_{T} (GeV/#it{c})");//#Sigma #it{p}_{T}
2667 fhPerpSumPtLeadingPt[icone] ->SetXTitle("#it{p}_{T}^{leading} (GeV/#it{c})");
2668 outputContainer->Add(fhPerpSumPtLeadingPt[icone]) ;
2670 // pt in cone vs. pt leading in the forward region (for background subtraction studies)
2671 snprintf(name, buffersize,"hPerpPtLeadingPt_Cone_%d",icone);
2672 snprintf(title, buffersize,"#it{p}_{T} in isolation cone for #it{R} = %2.2f",fConeSizes[icone]);
2673 fhPerpPtLeadingPt[icone] = new TH2F(name, title, nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
2674 fhPerpPtLeadingPt[icone] ->SetYTitle("#it{p}_{T}^{cone} (GeV/#it{c})");
2675 fhPerpPtLeadingPt[icone] ->SetXTitle("#it{p}_{T}^{leading} (GeV/#it{c})");
2676 outputContainer->Add(fhPerpPtLeadingPt[icone]) ;
2680 for(Int_t imc = 0; imc < 9; imc++)
2682 snprintf(name , buffersize,"hSumPtLeadingPt_MC%s_Cone_%d",mcPartName[imc].Data(),icone);
2683 snprintf(title, buffersize,"Candidate %s #it{p}_{T} vs cone #Sigma #it{p}_{T} for #it{R}=%2.2f",mcPartType[imc].Data(),fConeSizes[icone]);
2684 fhSumPtLeadingPtMC[imc][icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
2685 fhSumPtLeadingPtMC[imc][icone]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
2686 fhSumPtLeadingPtMC[imc][icone]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
2687 outputContainer->Add(fhSumPtLeadingPtMC[imc][icone]) ;
2691 for(Int_t ipt = 0; ipt<fNPtThresFrac;ipt++)
2693 snprintf(name, buffersize,"hPtThres_Cone_%d_Pt%d",icone,ipt);
2694 snprintf(title, buffersize,"Isolated candidate #it{p}_{T} distribution for #it{R} = %2.2f and #it{p}_{T}^{th} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtThresholds[ipt]);
2695 fhPtThresIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2696 fhPtThresIsolated[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2697 outputContainer->Add(fhPtThresIsolated[icone][ipt]) ;
2699 snprintf(name, buffersize,"hPtFrac_Cone_%d_Pt%d",icone,ipt);
2700 snprintf(title, buffersize,"Isolated candidate #it{p}_{T} distribution for #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2701 fhPtFracIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2702 fhPtFracIsolated[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2703 outputContainer->Add(fhPtFracIsolated[icone][ipt]) ;
2705 snprintf(name, buffersize,"hSumPt_Cone_%d_Pt%d",icone,ipt);
2706 snprintf(title, buffersize,"Isolated candidate #it{p}_{T} distribution for #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2707 fhSumPtIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2708 // fhSumPtIsolated[icone][ipt]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
2709 fhSumPtIsolated[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2710 outputContainer->Add(fhSumPtIsolated[icone][ipt]) ;
2712 snprintf(name, buffersize,"hPtSumDensity_Cone_%d_Pt%d",icone,ipt);
2713 snprintf(title, buffersize,"Isolated candidate #it{p}_{T} distribution for density in #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2714 fhPtSumDensityIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
2715 //fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
2716 fhPtSumDensityIso[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2717 outputContainer->Add(fhPtSumDensityIso[icone][ipt]) ;
2719 snprintf(name, buffersize,"hPtFracPtSum_Cone_%d_Pt%d",icone,ipt);
2720 snprintf(title, buffersize,"Isolated candidate #it{p}_{T} distribution for PtFracPtSum in #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2721 fhPtFracPtSumIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
2722 //fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
2723 fhPtFracPtSumIso[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2724 outputContainer->Add(fhPtFracPtSumIso[icone][ipt]) ;
2727 snprintf(name, buffersize,"hEtaPhiPtThres_Cone_%d_Pt%d",icone,ipt);
2728 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for #it{R} = %2.2f and #it{p}_{T}^{th} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtThresholds[ipt]);
2729 fhEtaPhiPtThresIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2730 fhEtaPhiPtThresIso[icone][ipt]->SetXTitle("#eta");
2731 fhEtaPhiPtThresIso[icone][ipt]->SetYTitle("#phi");
2732 outputContainer->Add(fhEtaPhiPtThresIso[icone][ipt]) ;
2734 snprintf(name, buffersize,"hEtaPhiPtFrac_Cone_%d_Pt%d",icone,ipt);
2735 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2736 fhEtaPhiPtFracIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2737 fhEtaPhiPtFracIso[icone][ipt]->SetXTitle("#eta");
2738 fhEtaPhiPtFracIso[icone][ipt]->SetYTitle("#phi");
2739 outputContainer->Add(fhEtaPhiPtFracIso[icone][ipt]) ;
2741 snprintf(name, buffersize,"hEtaPhiPtSum_Cone_%d_Pt%d",icone,ipt);
2742 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2743 fhEtaPhiPtSumIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2744 fhEtaPhiPtSumIso[icone][ipt]->SetXTitle("#eta");
2745 fhEtaPhiPtSumIso[icone][ipt]->SetYTitle("#phi");
2746 outputContainer->Add(fhEtaPhiPtSumIso[icone][ipt]) ;
2748 snprintf(name, buffersize,"hEtaPhiSumDensity_Cone_%d_Pt%d",icone,ipt);
2749 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for density #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2750 fhEtaPhiSumDensityIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2751 fhEtaPhiSumDensityIso[icone][ipt]->SetXTitle("#eta");
2752 fhEtaPhiSumDensityIso[icone][ipt]->SetYTitle("#phi");
2753 outputContainer->Add(fhEtaPhiSumDensityIso[icone][ipt]) ;
2755 snprintf(name, buffersize,"hEtaPhiFracPtSum_Cone_%d_Pt%d",icone,ipt);
2756 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for FracPtSum #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2757 fhEtaPhiFracPtSumIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2758 fhEtaPhiFracPtSumIso[icone][ipt]->SetXTitle("#eta");
2759 fhEtaPhiFracPtSumIso[icone][ipt]->SetYTitle("#phi");
2760 outputContainer->Add(fhEtaPhiFracPtSumIso[icone][ipt]) ;
2762 if(fFillTaggedDecayHistograms)
2764 // pt decays isolated
2765 snprintf(name, buffersize,"hPtThres_Decay_Cone_%d_Pt%d",icone,ipt);
2766 snprintf(title, buffersize,"Isolated decay candidate #it{p}_{T} distribution for #it{R} = %2.2f and #it{p}_{T}^{th} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtThresholds[ipt]);
2767 fhPtPtThresDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2768 fhPtPtThresDecayIso[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2769 outputContainer->Add(fhPtPtThresDecayIso[icone][ipt]) ;
2771 snprintf(name, buffersize,"hPtFrac_Decay_Cone_%d_Pt%d",icone,ipt);
2772 snprintf(title, buffersize,"Isolated decay candidate #it{p}_{T} distribution for #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2773 fhPtPtFracDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2774 fhPtPtFracDecayIso[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2775 outputContainer->Add(fhPtPtFracDecayIso[icone][ipt]) ;
2777 snprintf(name, buffersize,"hPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
2778 snprintf(title, buffersize,"Isolated decay candidate #it{p}_{T} distribution for #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2779 fhPtPtSumDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
2780 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
2781 fhPtPtSumDecayIso[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2782 outputContainer->Add(fhPtPtSumDecayIso[icone][ipt]) ;
2784 snprintf(name, buffersize,"hPtSumDensity_Decay_Cone_%d_Pt%d",icone,ipt);
2785 snprintf(title, buffersize,"Isolated decay candidate #it{p}_{T} distribution for density in #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2786 fhPtSumDensityDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
2787 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
2788 fhPtSumDensityDecayIso[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2789 outputContainer->Add(fhPtSumDensityDecayIso[icone][ipt]) ;
2791 snprintf(name, buffersize,"hPtFracPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
2792 snprintf(title, buffersize,"Isolated decay candidate #it{p}_{T} distribution for PtFracPtSum in #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2793 fhPtFracPtSumDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
2794 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
2795 fhPtFracPtSumDecayIso[icone][ipt]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2796 outputContainer->Add(fhPtFracPtSumDecayIso[icone][ipt]) ;
2799 snprintf(name, buffersize,"hEtaPhiPtThres_Decay_Cone_%d_Pt%d",icone,ipt);
2800 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for #it{R} = %2.2f and #it{p}_{T}^{th} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtThresholds[ipt]);
2801 fhEtaPhiPtThresDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2802 fhEtaPhiPtThresDecayIso[icone][ipt]->SetXTitle("#eta");
2803 fhEtaPhiPtThresDecayIso[icone][ipt]->SetYTitle("#phi");
2804 outputContainer->Add(fhEtaPhiPtThresDecayIso[icone][ipt]) ;
2806 snprintf(name, buffersize,"hEtaPhiPtFrac_Decay_Cone_%d_Pt%d",icone,ipt);
2807 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2808 fhEtaPhiPtFracDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2809 fhEtaPhiPtFracDecayIso[icone][ipt]->SetXTitle("#eta");
2810 fhEtaPhiPtFracDecayIso[icone][ipt]->SetYTitle("#phi");
2811 outputContainer->Add(fhEtaPhiPtFracDecayIso[icone][ipt]) ;
2814 snprintf(name, buffersize,"hEtaPhiPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
2815 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2816 fhEtaPhiPtSumDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2817 fhEtaPhiPtSumDecayIso[icone][ipt]->SetXTitle("#eta");
2818 fhEtaPhiPtSumDecayIso[icone][ipt]->SetYTitle("#phi");
2819 outputContainer->Add(fhEtaPhiPtSumDecayIso[icone][ipt]) ;
2821 snprintf(name, buffersize,"hEtaPhiSumDensity_Decay_Cone_%d_Pt%d",icone,ipt);
2822 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for density #it{R} = %2.2f and #it{p}_{T}^{sum} = %2.2f GeV/#it{c}",fConeSizes[icone],fSumPtThresholds[ipt]);
2823 fhEtaPhiSumDensityDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2824 fhEtaPhiSumDensityDecayIso[icone][ipt]->SetXTitle("#eta");
2825 fhEtaPhiSumDensityDecayIso[icone][ipt]->SetYTitle("#phi");
2826 outputContainer->Add(fhEtaPhiSumDensityDecayIso[icone][ipt]) ;
2828 snprintf(name, buffersize,"hEtaPhiFracPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
2829 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for FracPtSum #it{R} = %2.2f and #it{p}_{T}^{fr} = %2.2f GeV/#it{c}",fConeSizes[icone],fPtFractions[ipt]);
2830 fhEtaPhiFracPtSumDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
2831 fhEtaPhiFracPtSumDecayIso[icone][ipt]->SetXTitle("#eta");
2832 fhEtaPhiFracPtSumDecayIso[icone][ipt]->SetYTitle("#phi");
2833 outputContainer->Add(fhEtaPhiFracPtSumDecayIso[icone][ipt]) ;
2839 for(Int_t imc = 0; imc < 9; imc++)
2841 snprintf(name , buffersize,"hPtThreshMC%s_Cone_%d_Pt%d",mcPartName[imc].Data(),icone,ipt);
2842 snprintf(title, buffersize,"Isolated %s #it{p}_{T} for #it{R}=%2.2f and #it{p}_{T}^{th}=%2.2f",
2843 mcPartType[imc].Data(),fConeSizes[icone], fPtThresholds[ipt]);
2844 fhPtThresIsolatedMC[imc][icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2845 fhPtThresIsolatedMC[imc][icone][ipt]->SetYTitle("#it{counts}");
2846 fhPtThresIsolatedMC[imc][icone][ipt]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
2847 outputContainer->Add(fhPtThresIsolatedMC[imc][icone][ipt]) ;
2850 snprintf(name , buffersize,"hPtFracMC%s_Cone_%d_Pt%d",mcPartName[imc].Data(),icone,ipt);
2851 snprintf(title, buffersize,"Isolated %s #it{p}_{T} for #it{R}=%2.2f and #Sigma #it{p}_{T}^{in cone}/#it{p}_{T}^{trig}=%2.2f",
2852 mcPartType[imc].Data(),fConeSizes[icone], fPtFractions[ipt]);
2853 fhPtFracIsolatedMC[imc][icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2854 fhPtFracIsolatedMC[imc][icone][ipt]->SetYTitle("#it{counts}");
2855 fhPtFracIsolatedMC[imc][icone][ipt]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
2856 outputContainer->Add(fhPtFracIsolatedMC[imc][icone][ipt]) ;
2858 snprintf(name , buffersize,"hSumPtMC%s_Cone_%d_Pt%d",mcPartName[imc].Data(),icone,ipt);
2859 snprintf(title, buffersize,"Isolated %s #it{p}_{T} for #it{R}=%2.2f and #Sigma #it{p}_{T}^{in cone}=%2.2f",
2860 mcPartType[imc].Data(),fConeSizes[icone], fSumPtThresholds[ipt]);
2861 fhSumPtIsolatedMC[imc][icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
2862 fhSumPtIsolatedMC[imc][icone][ipt]->SetYTitle("#it{counts}");
2863 fhSumPtIsolatedMC[imc][icone][ipt]->SetXTitle("#it{p}_{T}(GeV/#it{c})");
2864 outputContainer->Add(fhSumPtIsolatedMC[imc][icone][ipt]) ;
2871 if(fFillPileUpHistograms)
2873 for (Int_t i = 0; i < 7 ; i++)
2875 fhEIsoPileUp[i] = new TH1F(Form("hEPileUp%s",pileUpName[i].Data()),
2876 Form("Number of isolated particles vs E, %s, pile-up event by %s",parTitle.Data(),pileUpName[i].Data()),
2877 nptbins,ptmin,ptmax);
2878 fhEIsoPileUp[i]->SetYTitle("d#it{N} / d#it{E}");
2879 fhEIsoPileUp[i]->SetXTitle("#it{E} (GeV)");
2880 outputContainer->Add(fhEIsoPileUp[i]) ;
2882 fhPtIsoPileUp[i] = new TH1F(Form("hPtPileUp%s",pileUpName[i].Data()),
2883 Form("Number of isolated particles vs #it{p}_{T}, %s, pile-up event by %s",parTitle.Data(),pileUpName[i].Data()),
2884 nptbins,ptmin,ptmax);
2885 fhPtIsoPileUp[i]->SetYTitle("d#it{N} / #it{p}_{T}");
2886 fhPtIsoPileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2887 outputContainer->Add(fhPtIsoPileUp[i]) ;
2889 fhENoIsoPileUp[i] = new TH1F(Form("hENoIsoPileUp%s",pileUpName[i].Data()),
2890 Form("Number of not isolated particles vs E, %s, pile-up event by %s",parTitle.Data(),pileUpName[i].Data()),
2891 nptbins,ptmin,ptmax);
2892 fhENoIsoPileUp[i]->SetYTitle("d#it{N} / dE");
2893 fhENoIsoPileUp[i]->SetXTitle("#it{E} (GeV)");
2894 outputContainer->Add(fhENoIsoPileUp[i]) ;
2896 fhPtNoIsoPileUp[i] = new TH1F(Form("hPtNoIsoPileUp%s",pileUpName[i].Data()),
2897 Form("Number of not isolated particles vs #it{p}_{T}, %s, pile-up event by %s",parTitle.Data(),pileUpName[i].Data()),
2898 nptbins,ptmin,ptmax);
2899 fhPtNoIsoPileUp[i]->SetYTitle("d#it{N} / #it{p}_{T}");
2900 fhPtNoIsoPileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
2901 outputContainer->Add(fhPtNoIsoPileUp[i]) ;
2904 fhTimeENoCut = new TH2F ("hTimeE_NoCut","time of cluster vs E of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2905 fhTimeENoCut->SetXTitle("#it{E} (GeV)");
2906 fhTimeENoCut->SetYTitle("#it{time} (ns)");
2907 outputContainer->Add(fhTimeENoCut);
2909 fhTimeESPD = new TH2F ("hTimeE_SPD","time of cluster vs E of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2910 fhTimeESPD->SetXTitle("#it{E} (GeV)");
2911 fhTimeESPD->SetYTitle("#it{time} (ns)");
2912 outputContainer->Add(fhTimeESPD);
2914 fhTimeESPDMulti = new TH2F ("hTimeE_SPDMulti","time of cluster vs E of clusters, SPD multi cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
2915 fhTimeESPDMulti->SetXTitle("#it{E} (GeV)");
2916 fhTimeESPDMulti->SetYTitle("#it{time} (ns)");
2917 outputContainer->Add(fhTimeESPDMulti);
2919 fhTimeNPileUpVertSPD = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
2920 fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
2921 fhTimeNPileUpVertSPD->SetXTitle("#it{time} (ns)");
2922 outputContainer->Add(fhTimeNPileUpVertSPD);
2924 fhTimeNPileUpVertTrack = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 50,0,50 );
2925 fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
2926 fhTimeNPileUpVertTrack->SetXTitle("#it{time} (ns)");
2927 outputContainer->Add(fhTimeNPileUpVertTrack);
2929 fhTimeNPileUpVertContributors = new TH2F ("hTime_NPileUpVertContributors","time of cluster vs N constributors to pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
2930 fhTimeNPileUpVertContributors->SetYTitle("# vertex ");
2931 fhTimeNPileUpVertContributors->SetXTitle("#it{time} (ns)");
2932 outputContainer->Add(fhTimeNPileUpVertContributors);
2934 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);
2935 fhTimePileUpMainVertexZDistance->SetYTitle("distance #it{z} (cm) ");
2936 fhTimePileUpMainVertexZDistance->SetXTitle("#it{time} (ns)");
2937 outputContainer->Add(fhTimePileUpMainVertexZDistance);
2939 fhTimePileUpMainVertexZDiamond = new TH2F ("hTime_PileUpMainVertexZDiamond","time of cluster vs distance in Z pile-up SPD vertex - z diamond",ntimebins,timemin,timemax,100,0,50);
2940 fhTimePileUpMainVertexZDiamond->SetYTitle("diamond distance #it{z} (cm) ");
2941 fhTimePileUpMainVertexZDiamond->SetXTitle("#it{time} (ns)");
2942 outputContainer->Add(fhTimePileUpMainVertexZDiamond);
2945 return outputContainer ;
2949 //____________________________________________________
2950 Int_t AliAnaParticleIsolation::GetMCIndex(Int_t mcTag)
2952 // Histogram index depending on origin of candidate
2954 if(!IsDataMC()) return -1;
2956 if (GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt))
2960 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation))
2964 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0))
2968 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay))
2972 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay))
2976 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay))
2978 return kmcOtherDecay;
2980 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron))
2984 else // anything else
2986 // careful can contain also other decays, to be checked.
2991 //__________________________________
2992 void AliAnaParticleIsolation::Init()
2994 // Do some checks and init stuff
2996 // In case of several cone and thresholds analysis, open the cuts for the filling of the
2997 // track and cluster reference arrays in cone when done in the MakeAnalysisFillAOD().
2998 // The different cones, thresholds are tested for this list of tracks, clusters.
3001 printf("AliAnaParticleIsolation::Init() - Open default isolation cuts for multiple Isolation analysis\n");
3002 GetIsolationCut()->SetPtThreshold(100);
3003 GetIsolationCut()->SetPtFraction(100);
3004 GetIsolationCut()->SetConeSize(1);
3007 if(!GetReader()->IsCTSSwitchedOn() && GetIsolationCut()->GetParticleTypeInCone()!=AliIsolationCut::kOnlyNeutral)
3008 AliFatal("STOP!: You want to use CTS tracks in analysis but not read!! \n!!Check the configuration file!!\n");
3013 //____________________________________________
3014 void AliAnaParticleIsolation::InitParameters()
3017 //Initialize the parameters of the analysis.
3018 SetInputAODName("PWG4Particle");
3019 SetAODObjArrayName("IsolationCone");
3020 AddToHistogramsName("AnaIsolation_");
3022 fCalorimeter = "EMCAL" ;
3023 fIsoDetector = "EMCAL" ;
3025 fReMakeIC = kFALSE ;
3026 fMakeSeveralIC = kFALSE ;
3028 //----------- Several IC-----------------
3031 fConeSizes [0] = 0.1; fConeSizes [1] = 0.2; fConeSizes [2] = 0.3; fConeSizes [3] = 0.4; fConeSizes [4] = 0.5;
3032 fPtThresholds [0] = 1.; fPtThresholds [1] = 2.; fPtThresholds [2] = 3.; fPtThresholds [3] = 4.; fPtThresholds [4] = 5.;
3033 fPtFractions [0] = 0.05; fPtFractions [1] = 0.075; fPtFractions [2] = 0.1; fPtFractions [3] = 1.25; fPtFractions [4] = 1.5;
3034 fSumPtThresholds[0] = 1.; fSumPtThresholds[1] = 2.; fSumPtThresholds[2] = 3.; fSumPtThresholds[3] = 4.; fSumPtThresholds[4] = 5.;
3038 //__________________________________________________
3039 void AliAnaParticleIsolation::MakeAnalysisFillAOD()
3041 //Do analysis and fill aods
3042 //Search for the isolated photon in fCalorimeter with pt > GetMinPt()
3044 if(!GetInputAODBranch())
3045 AliFatal(Form("No input particles in AOD with name branch < %s >, STOP",GetInputAODName().Data()));
3047 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation"))
3048 AliFatal(Form("Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName()));
3050 Int_t n = 0, nfrac = 0;
3051 Bool_t isolated = kFALSE ;
3052 Float_t coneptsum = 0 ;
3053 TObjArray * pl = 0x0; ;
3055 //Select the calorimeter for candidate isolation with neutral particles
3056 if (fCalorimeter == "PHOS" )
3057 pl = GetPHOSClusters();
3058 else if (fCalorimeter == "EMCAL")
3059 pl = GetEMCALClusters();
3061 //Loop on AOD branch, filled previously in AliAnaPhoton, find leading particle to do isolation only with it
3062 Double_t ptLeading = 0. ;
3063 Int_t idLeading = -1 ;
3064 TLorentzVector mom ;
3065 Int_t naod = GetInputAODBranch()->GetEntriesFast();
3067 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Input aod branch entries %d\n", naod);
3069 for(Int_t iaod = 0; iaod < naod; iaod++)
3071 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
3073 //If too small or too large pt, skip
3074 if(aodinput->Pt() < GetMinPt() || aodinput->Pt() > GetMaxPt() ) continue ;
3076 //check if it is low pt trigger particle
3077 if((aodinput->Pt() < GetIsolationCut()->GetPtThreshold() ||
3078 aodinput->Pt() < GetIsolationCut()->GetSumPtThreshold()) &&
3081 continue ; //trigger should not come from underlying event
3084 //vertex cut in case of mixing
3085 Int_t check = CheckMixedEventVertex(aodinput->GetCaloLabel(0), aodinput->GetTrackLabel(0));
3086 if(check == 0) continue;
3087 if(check == -1) return;
3089 //find the leading particles with highest momentum
3090 if ( aodinput->Pt() > ptLeading )
3092 ptLeading = aodinput->Pt() ;
3096 aodinput->SetLeadingParticle(kFALSE);
3098 }//finish searching for leading trigger particle
3100 // Check isolation of leading particle
3101 if(idLeading < 0) return;
3103 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(idLeading));
3104 aodinput->SetLeadingParticle(kTRUE);
3106 // Check isolation only of clusters in fiducial region
3107 if(IsFiducialCutOn())
3109 Bool_t in = GetFiducialCut()->IsInFiducialCut(*aodinput->Momentum(),aodinput->GetDetector()) ;
3113 //After cuts, study isolation
3114 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
3115 GetIsolationCut()->MakeIsolationCut(GetCTSTracks(),pl,
3116 GetReader(), GetCaloPID(),
3117 kTRUE, aodinput, GetAODObjArrayName(),
3118 n,nfrac,coneptsum, isolated);
3120 if(!fMakeSeveralIC) aodinput->SetIsolated(isolated);
3124 if(isolated)printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() : Particle %d IS ISOLATED \n",idLeading);
3125 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - End fill AODs \n");
3130 //_________________________________________________________
3131 void AliAnaParticleIsolation::MakeAnalysisFillHistograms()
3133 //Do analysis and fill histograms
3135 // In case of simulated data, fill acceptance histograms
3136 // Not ready for multiple case analysis.
3137 if(IsDataMC() && !fMakeSeveralIC) FillAcceptanceHistograms();
3139 //Loop on stored AOD
3140 Int_t naod = GetInputAODBranch()->GetEntriesFast();
3142 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Histo aod branch entries %d\n", naod);
3144 for(Int_t iaod = 0; iaod < naod ; iaod++)
3146 AliAODPWG4ParticleCorrelation* aod = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
3148 if(!aod->IsLeadingParticle()) continue; // Try to isolate only leading cluster or track
3150 // Check isolation only of clusters in fiducial region
3151 if(IsFiducialCutOn())
3153 Bool_t in = GetFiducialCut()->IsInFiducialCut(*aod->Momentum(),aod->GetDetector()) ;
3154 if(! in ) continue ;
3157 Float_t pt = aod->Pt();
3159 //If too small or too large pt, skip
3160 if(pt < GetMinPt() || pt > GetMaxPt() ) continue ;
3162 Int_t mcTag = aod->GetTag() ;
3163 Int_t mcIndex = GetMCIndex(mcTag);
3165 // --- In case of redoing isolation from delta AOD ----
3166 // Not standard case, not used since its implementation
3169 //Analysis of multiple IC at same time
3170 MakeSeveralICAnalysis(aod,mcIndex);
3174 // --- In case of redoing isolation multiple cuts ----
3178 //In case a more strict IC is needed in the produced AOD
3179 Bool_t isolated = kFALSE;
3180 Int_t n = 0, nfrac = 0;
3181 Float_t coneptsum = 0 ;
3183 //Recover reference arrays with clusters and tracks
3184 TObjArray * refclusters = aod->GetObjArray(GetAODObjArrayName()+"Clusters");
3185 TObjArray * reftracks = aod->GetObjArray(GetAODObjArrayName()+"Tracks");
3187 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
3188 GetReader(), GetCaloPID(),
3190 n,nfrac,coneptsum, isolated);
3193 Bool_t isolated = aod->IsIsolated();
3194 Bool_t decay = aod->IsTagged();
3195 Float_t energy = aod->E();
3196 Float_t phi = aod->Phi();
3197 Float_t eta = aod->Eta();
3199 if(GetDebug() > 0) printf(" AliAnaParticleIsolation::MakeAnalysisFillHistograms() - pt %1.1f, eta %1.1f, phi %1.1f, Isolated %d\n",
3200 pt, eta, phi, isolated);
3202 //---------------------------------------------------------------
3203 // Fill Shower shape and track matching histograms
3204 //---------------------------------------------------------------
3206 FillTrackMatchingShowerShapeControlHistograms(aod,mcIndex);
3208 //---------------------------------------------------------------
3209 // Fill pt/sum pT distribution of particles in cone or in UE band
3210 //---------------------------------------------------------------
3212 Float_t coneptsumCluster = 0;
3213 Float_t coneptsumTrack = 0;
3214 Float_t coneptsumCell = 0;
3215 Float_t etaBandptsumClusterNorm = 0;
3216 Float_t etaBandptsumTrackNorm = 0;
3218 CalculateTrackSignalInCone (aod,coneptsumTrack );
3219 CalculateCaloSignalInCone (aod,coneptsumCluster);
3220 if(fFillCellHistograms)
3221 CalculateCaloCellSignalInCone(aod,coneptsumCell );
3223 if(GetIsolationCut()->GetParticleTypeInCone()==AliIsolationCut::kNeutralAndCharged)
3225 fhConeSumPtClustervsTrack ->Fill(coneptsumCluster,coneptsumTrack);
3226 if(fFillCellHistograms)
3228 fhConeSumPtCellvsTrack ->Fill(coneptsumCell, coneptsumTrack);
3229 fhConeSumPtCellTrack ->Fill(pt, coneptsumTrack+coneptsumCell);
3230 fhConeSumPtCellTrackTrigEtaPhi->Fill(eta,phi,coneptsumTrack+coneptsumCell);
3234 fhConeSumPt ->Fill(pt, coneptsumTrack+coneptsumCluster);
3235 fhConeSumPtTrigEtaPhi ->Fill(eta,phi,coneptsumTrack+coneptsumCluster);
3238 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d Energy Sum in Isolation Cone %2.2f\n", iaod, coneptsumTrack+coneptsumCluster);
3240 //normalize phi/eta band per area unit
3241 if(fFillUEBandSubtractHistograms)
3242 CalculateNormalizeUEBandPerUnitArea(aod, coneptsumCluster, coneptsumCell, coneptsumTrack, etaBandptsumTrackNorm, etaBandptsumClusterNorm) ;
3244 // printf("Histograms analysis : cluster pt = %f, etaBandTrack = %f, etaBandCluster = %f, isolation = %d\n",aod->Pt(),etaBandptsumTrackNorm,etaBandptsumClusterNorm,aod->IsIsolated());
3247 //---------------------------------------------------------------
3248 // Isolated/ Non isolated histograms
3249 //---------------------------------------------------------------
3254 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d ISOLATED: fill histograms\n", iaod);
3256 fhEIso ->Fill(energy);
3258 fhPhiIso ->Fill(pt,phi);
3259 fhEtaIso ->Fill(pt,eta);
3260 fhEtaPhiIso ->Fill(eta,phi);
3264 // For histograms in arrays, index in the array, corresponding to any particle origin
3265 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
3267 fhPtIsoMC [kmcPhoton]->Fill(pt);
3268 fhPhiIsoMC[kmcPhoton]->Fill(pt,phi);
3269 fhEtaIsoMC[kmcPhoton]->Fill(pt,eta);
3272 fhPtIsoMC [mcIndex]->Fill(pt);
3273 fhPhiIsoMC[mcIndex]->Fill(pt,phi);
3274 fhEtaIsoMC[mcIndex]->Fill(pt,eta);
3275 }//Histograms with MC
3277 // Candidates tagged as decay in another analysis (AliAnaPi0EbE)
3278 if(decay && fFillTaggedDecayHistograms)
3280 fhPtDecayIso ->Fill(pt);
3281 fhEtaPhiDecayIso->Fill(eta,phi);
3285 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
3286 fhPtDecayIsoMC[kmcPhoton]->Fill(pt);
3288 fhPtDecayIsoMC[mcIndex]->Fill(pt);
3292 if(fFillNLMHistograms)
3293 fhPtNLocMaxIso ->Fill(pt,aod->GetFiducialArea()) ; // remember to change method name
3295 if(fFillHighMultHistograms)
3297 fhPtCentralityIso ->Fill(pt,GetEventCentrality()) ;
3298 fhPtEventPlaneIso ->Fill(pt,GetEventPlaneAngle() ) ;
3301 if(fFillPileUpHistograms)
3303 if(GetReader()->IsPileUpFromSPD()) { fhEIsoPileUp[0] ->Fill(energy) ; fhPtIsoPileUp[0]->Fill(pt) ; }
3304 if(GetReader()->IsPileUpFromEMCal()) { fhEIsoPileUp[1] ->Fill(energy) ; fhPtIsoPileUp[1]->Fill(pt) ; }
3305 if(GetReader()->IsPileUpFromSPDOrEMCal()) { fhEIsoPileUp[2] ->Fill(energy) ; fhPtIsoPileUp[2]->Fill(pt) ; }
3306 if(GetReader()->IsPileUpFromSPDAndEMCal()) { fhEIsoPileUp[3] ->Fill(energy) ; fhPtIsoPileUp[3]->Fill(pt) ; }
3307 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) { fhEIsoPileUp[4] ->Fill(energy) ; fhPtIsoPileUp[4]->Fill(pt) ; }
3308 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) { fhEIsoPileUp[5] ->Fill(energy) ; fhPtIsoPileUp[5]->Fill(pt) ; }
3309 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhEIsoPileUp[6] ->Fill(energy) ; fhPtIsoPileUp[6]->Fill(pt) ; }
3311 // Fill histograms to undertand pile-up before other cuts applied
3312 // Remember to relax time cuts in the reader
3313 FillPileUpHistograms(aod->GetCaloLabel(0));
3316 }//Isolated histograms
3317 else // NON isolated
3320 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d NOT ISOLATED, fill histograms\n", iaod);
3322 fhENoIso ->Fill(energy);
3323 fhPtNoIso ->Fill(pt);
3324 fhEtaPhiNoIso ->Fill(eta,phi);
3328 if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) )
3329 fhPtNoIsoMC[kmcPhoton]->Fill(pt);
3331 fhPtNoIsoMC[mcIndex]->Fill(pt);
3334 // Candidates tagged as decay in another analysis (AliAnaPi0EbE)
3335 if(decay && fFillTaggedDecayHistograms)
3337 fhPtDecayNoIso ->Fill(pt);
3338 fhEtaPhiDecayNoIso->Fill(eta,phi);
3342 if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) )
3343 fhPtDecayNoIsoMC[kmcPhoton]->Fill(pt);
3345 fhPtDecayNoIsoMC[mcIndex]->Fill(pt);
3349 if(fFillNLMHistograms)
3350 fhPtNLocMaxNoIso ->Fill(pt,aod->GetFiducialArea()); // remember to change method name
3352 if(fFillPileUpHistograms)
3354 if(GetReader()->IsPileUpFromSPD()) { fhENoIsoPileUp[0] ->Fill(energy) ; fhPtNoIsoPileUp[0]->Fill(pt) ; }
3355 if(GetReader()->IsPileUpFromEMCal()) { fhENoIsoPileUp[1] ->Fill(energy) ; fhPtNoIsoPileUp[1]->Fill(pt) ; }
3356 if(GetReader()->IsPileUpFromSPDOrEMCal()) { fhENoIsoPileUp[2] ->Fill(energy) ; fhPtNoIsoPileUp[2]->Fill(pt) ; }
3357 if(GetReader()->IsPileUpFromSPDAndEMCal()) { fhENoIsoPileUp[3] ->Fill(energy) ; fhPtNoIsoPileUp[3]->Fill(pt) ; }
3358 if(GetReader()->IsPileUpFromSPDAndNotEMCal()) { fhENoIsoPileUp[4] ->Fill(energy) ; fhPtNoIsoPileUp[4]->Fill(pt) ; }
3359 if(GetReader()->IsPileUpFromEMCalAndNotSPD()) { fhENoIsoPileUp[5] ->Fill(energy) ; fhPtNoIsoPileUp[5]->Fill(pt) ; }
3360 if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhENoIsoPileUp[6] ->Fill(energy) ; fhPtNoIsoPileUp[6]->Fill(pt) ; }
3367 //______________________________________________________________________
3368 void AliAnaParticleIsolation::FillAcceptanceHistograms()
3370 // Fill acceptance histograms if MC data is available
3371 // Only when particle is in the calorimeter. Rethink if CTS is used.
3373 if(GetDebug() > 0) printf("AliAnaParticleIsolation::FillAcceptanceHistograms() - Start \n");
3375 //Double_t photonY = -100 ;
3376 Double_t photonE = -1 ;
3377 Double_t photonPt = -1 ;
3378 Double_t photonPhi = 100 ;
3379 Double_t photonEta = -1 ;
3387 TParticle * primStack = 0;
3388 AliAODMCParticle * primAOD = 0;
3391 // Get the ESD MC particles container
3392 AliStack * stack = 0;
3393 if( GetReader()->ReadStack() )
3395 stack = GetMCStack();
3397 nprim = stack->GetNtrack();
3400 // Get the AOD MC particles container
3401 TClonesArray * mcparticles = 0;
3402 if( GetReader()->ReadAODMCParticles() )
3404 mcparticles = GetReader()->GetAODMCParticles();
3405 if( !mcparticles ) return;
3406 nprim = mcparticles->GetEntriesFast();
3409 for(Int_t i=0 ; i < nprim; i++)
3411 if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
3413 if(GetReader()->ReadStack())
3415 primStack = stack->Particle(i) ;
3418 printf("AliAnaParticleIsolation::FillAcceptanceHistograms() - ESD primaries pointer not available!!\n");
3422 pdg = primStack->GetPdgCode();
3423 status = primStack->GetStatusCode();
3425 if(primStack->Energy() == TMath::Abs(primStack->Pz())) continue ; //Protection against floating point exception
3427 //printf("i %d, %s %d %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
3428 // prim->GetName(), prim->GetPdgCode());
3430 //photonY = 0.5*TMath::Log((prim->Energy()-prim->Pz())/(prim->Energy()+prim->Pz())) ;
3433 primStack->Momentum(lv);
3438 primAOD = (AliAODMCParticle *) mcparticles->At(i);
3441 printf("AliAnaParticleIsolation::FillAcceptanceHistograms() - AOD primaries pointer not available!!\n");
3445 pdg = primAOD->GetPdgCode();
3446 status = primAOD->GetStatus();
3448 if(primAOD->E() == TMath::Abs(primAOD->Pz())) continue ; //Protection against floating point exception
3450 //photonY = 0.5*TMath::Log((prim->Energy()-prim->Pz())/(prim->Energy()+prim->Pz())) ;
3453 lv.SetPxPyPzE(primAOD->Px(),primAOD->Py(),primAOD->Pz(),primAOD->E());
3456 // Select only photons in the final state
3457 if(pdg != 22 ) continue ;
3459 // Consider only final state particles, but this depends on generator,
3460 // status 1 is the usual one, in case of not being ok, leave the possibility
3461 // to not consider this.
3462 if(status != 1 && GetMCAnalysisUtils()->GetMCGenerator()!="" ) continue ;
3464 // If too small or too large pt, skip, same cut as for data analysis
3465 photonPt = lv.Pt () ;
3467 if(photonPt < GetMinPt() || photonPt > GetMaxPt() ) continue ;
3470 photonEta = lv.Eta() ;
3471 photonPhi = lv.Phi() ;
3473 if(photonPhi < 0) photonPhi+=TMath::TwoPi();
3475 // Check if photons hit the Calorimeter acceptance
3476 if(IsRealCaloAcceptanceOn() && fIsoDetector!="CTS") // defined on base class
3478 if(GetReader()->ReadStack() &&
3479 !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(fIsoDetector, primStack)) continue ;
3480 if(GetReader()->ReadAODMCParticles() &&
3481 !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(fIsoDetector, primAOD )) continue ;
3484 // Check same fidutial borders as in data analysis on top of real acceptance if real was requested.
3485 if(!GetFiducialCut()->IsInFiducialCut(lv,fIsoDetector)) continue ;
3487 // Get tag of this particle photon from fragmentation, decay, prompt ...
3488 // Set the origin of the photon.
3489 tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader());
3491 if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
3493 // A conversion photon from a hadron, skip this kind of photon
3494 // printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
3495 // GetMCAnalysisUtils()->PrintMCTag(tag);
3501 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) )
3503 mcIndex = kmcPrimPrompt;
3505 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) )
3507 mcIndex = kmcPrimFrag ;
3509 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
3511 mcIndex = kmcPrimISR;
3513 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
3515 mcIndex = kmcPrimPi0Decay;
3517 else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
3518 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)))
3520 mcIndex = kmcPrimOtherDecay;
3524 // Other decay but from non final state particle
3525 mcIndex = kmcPrimOtherDecay;
3528 // ////////////////////ISO MC/////////////////////////
3529 Double_t sumPtInCone = 0; Double_t dR=0. ;
3530 TParticle * mcisopStack = 0;
3531 AliAODMCParticle * mcisopAOD = 0;
3532 Int_t partInConeStatus = -1, partInConeMother = -1;
3533 Double_t partInConePt = 0, partInConeEta = 0, partInConePhi = 0;
3535 for(Int_t ip = 0; ip < nprim ; ip++)
3539 if( GetReader()->ReadStack() )
3541 mcisopStack = static_cast<TParticle*>(stack->Particle(ip));
3542 if( !mcisopStack ) continue;
3543 partInConeStatus = mcisopStack->GetStatusCode();
3544 partInConeMother = mcisopStack->GetMother(0);
3545 partInConePt = mcisopStack->Pt();
3546 partInConeEta = mcisopStack->Eta();
3547 partInConePhi = mcisopStack->Phi();
3551 mcisopAOD = (AliAODMCParticle *) mcparticles->At(ip);
3552 if( !mcisopAOD ) continue;
3553 partInConeStatus = mcisopAOD->GetStatus();
3554 partInConeMother = mcisopAOD->GetMother();
3555 partInConePt = mcisopAOD->Pt();
3556 partInConeEta = mcisopAOD->Eta();
3557 partInConePhi = mcisopAOD->Phi();
3560 // Consider only final state particles, but this depends on generator,
3561 // status 1 is the usual one, in case of not being ok, leave the possibility
3562 // to not consider this.
3563 if( partInConeStatus != 1 && GetMCAnalysisUtils()->GetMCGenerator()!="" ) continue ;
3565 if( partInConeMother == i ) continue;
3567 if( partInConePt < GetReader()->GetCTSPtMin() ) continue;
3568 // Careful!!!, cut for TPC tracks and calorimeter clusters in cone can be different
3570 // TVector3 vmcv(mcisop->Px(),mcisop->Py(), mcisop->Pz());
3571 // if(vmcv.Perp()>1)
3575 // Add here Acceptance cut???, charged particles CTS fid cut, neutral Calo real acceptance.
3578 dR = GetIsolationCut()->Radius(photonEta, photonPhi, partInConeEta, partInConePhi);
3580 if(dR > GetIsolationCut()->GetConeSize())
3583 sumPtInCone += partInConePt;
3584 if(partInConePt > GetIsolationCut()->GetPtThreshold() &&
3585 partInConePt < GetIsolationCut()->GetPtThresholdMax()) npart++;
3588 ///////END ISO MC/////////////////////////
3590 // Fill the histograms, only those in the defined calorimeter acceptance
3592 fhEtaPrimMC[kmcPrimPhoton]->Fill(photonPt , photonEta) ;
3593 fhPhiPrimMC[kmcPrimPhoton]->Fill(photonPt , photonPhi) ;
3594 fhEPrimMC [kmcPrimPhoton]->Fill(photonE) ;
3596 fhEtaPrimMC[mcIndex]->Fill(photonPt , photonEta) ;
3597 fhPhiPrimMC[mcIndex]->Fill(photonPt , photonPhi) ;
3598 fhEPrimMC [mcIndex]->Fill(photonE ) ;
3601 Bool_t isolated = kFALSE;
3602 if(GetIsolationCut()->GetICMethod() == AliIsolationCut::kSumPtIC &&
3603 (sumPtInCone < GetIsolationCut()->GetSumPtThreshold() ||
3604 sumPtInCone > GetIsolationCut()->GetSumPtThresholdMax()))
3607 if(GetIsolationCut()->GetICMethod() == AliIsolationCut::kPtThresIC &&
3613 fhPtPrimMCiso [mcIndex] ->Fill(photonPt) ;
3614 fhPtPrimMCiso [kmcPrimPhoton]->Fill(photonPt) ;
3617 }//loop on primaries
3619 if(GetDebug() > 0) printf("AliAnaParticleIsolation::FillAcceptanceHistograms() - End \n");
3624 //_____________________________________________________________________________________
3625 void AliAnaParticleIsolation::MakeSeveralICAnalysis(AliAODPWG4ParticleCorrelation* ph,
3629 //Isolation Cut Analysis for both methods and different pt cuts and cones
3630 Float_t ptC = ph->Pt();
3631 Float_t etaC = ph->Eta();
3632 Float_t phiC = ph->Phi();
3633 Int_t tag = ph->GetTag();
3634 Bool_t decay = ph->IsTagged();
3637 printf("AliAnaParticleIsolation::MakeSeveralICAnalysis() - Isolate pT %2.2f, decay? %d\n",ptC, decay);
3639 //Keep original setting used when filling AODs, reset at end of analysis
3640 Float_t ptthresorg = GetIsolationCut()->GetPtThreshold();
3641 Float_t ptfracorg = GetIsolationCut()->GetPtFraction();
3642 Float_t ptsumcorg = GetIsolationCut()->GetSumPtThreshold();
3643 Float_t rorg = GetIsolationCut()->GetConeSize();
3645 Float_t coneptsum = 0 ;
3646 Int_t n [10][10];//[fNCones][fNPtThresFrac];
3647 Int_t nfrac[10][10];//[fNCones][fNPtThresFrac];
3648 Bool_t isolated = kFALSE;
3650 Int_t nFracCone = 0;
3652 // Fill hist with all particles before isolation criteria
3653 fhENoIso ->Fill(ph->E());
3654 fhPtNoIso ->Fill(ptC);
3655 fhEtaPhiNoIso->Fill(etaC,phiC);
3659 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
3660 fhPtNoIsoMC[kmcPhoton]->Fill(ptC);
3662 fhPtNoIsoMC[mcIndex]->Fill(ptC);
3665 // Candidates tagged as decay in another analysis (AliAnaPi0EbE)
3666 if(decay && fFillTaggedDecayHistograms)
3668 fhPtDecayNoIso ->Fill(ptC);
3669 fhEtaPhiDecayNoIso->Fill(etaC,phiC);
3673 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
3674 fhPtDecayNoIsoMC[kmcPhoton]->Fill(ptC);
3676 fhPtDecayNoIsoMC[mcIndex]->Fill(ptC);
3680 //Get vertex for photon momentum calculation
3681 Double_t vertex[] = {0,0,0} ; //vertex ;
3682 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
3683 GetReader()->GetVertex(vertex);
3685 //Loop on cone sizes
3686 for(Int_t icone = 0; icone<fNCones; icone++)
3688 //Recover reference arrays with clusters and tracks
3689 TObjArray * refclusters = ph->GetObjArray(GetAODObjArrayName()+"Clusters");
3690 TObjArray * reftracks = ph->GetObjArray(GetAODObjArrayName()+"Tracks");
3692 //If too small or too large pt, skip
3693 if(ptC < GetMinPt() || ptC > GetMaxPt() ) continue ;
3695 //In case a more strict IC is needed in the produced AOD
3697 nCone=0; nFracCone = 0; isolated = kFALSE; coneptsum = 0;
3699 GetIsolationCut()->SetSumPtThreshold(100);
3700 GetIsolationCut()->SetPtThreshold(100);
3701 GetIsolationCut()->SetPtFraction(100);
3702 GetIsolationCut()->SetConeSize(fConeSizes[icone]);
3704 // retreive pt tracks to fill histo vs. pt leading
3705 //Fill pt distribution of particles in cone
3706 //fhPtLeadingPt(),fhPerpSumPtLeadingPt(),fhPerpPtLeadingPt(),
3708 // Tracks in perpendicular cones
3709 Double_t sumptPerp = 0. ;
3710 TObjArray * trackList = GetCTSTracks() ;
3711 for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++)
3713 AliVTrack* track = (AliVTrack *) trackList->At(itrack);
3714 //fill the histograms at forward range
3717 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Track not available?");
3721 Double_t dPhi = phiC - track->Phi() + TMath::PiOver2();
3722 Double_t dEta = etaC - track->Eta();
3723 Double_t arg = dPhi*dPhi + dEta*dEta;
3724 if(TMath::Sqrt(arg) < fConeSizes[icone])
3726 fhPerpPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
3727 sumptPerp+=track->Pt();
3730 dPhi = phiC - track->Phi() - TMath::PiOver2();
3731 arg = dPhi*dPhi + dEta*dEta;
3732 if(TMath::Sqrt(arg) < fConeSizes[icone])
3734 fhPerpPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
3735 sumptPerp+=track->Pt();
3739 fhPerpSumPtLeadingPt[icone]->Fill(ptC,sumptPerp);
3741 // Tracks in isolation cone, pT distribution and sum
3742 if(reftracks && GetIsolationCut()->GetParticleTypeInCone()!= AliIsolationCut::kOnlyNeutral)
3744 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++)
3746 AliVTrack* track = (AliVTrack *) reftracks->At(itrack);
3748 Float_t rad = GetIsolationCut()->Radius(etaC, phiC, track->Eta(), track->Phi());
3750 if(rad > fConeSizes[icone]) continue ;
3752 fhPtLeadingPt[icone]->Fill(ptC, track->Pt());
3753 coneptsum += track->Pt();
3757 // Clusters in isolation cone, pT distribution and sum
3758 if(refclusters && GetIsolationCut()->GetParticleTypeInCone()!= AliIsolationCut::kOnlyCharged)
3760 TLorentzVector mom ;
3761 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++)
3763 AliVCluster* calo = (AliVCluster *) refclusters->At(icalo);
3765 calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
3767 Float_t rad = GetIsolationCut()->Radius(etaC, phiC, mom.Eta(), mom.Phi());
3769 if(rad > fConeSizes[icone]) continue ;
3771 fhPtLeadingPt[icone]->Fill(ptC, mom.Pt());
3772 coneptsum += mom.Pt();
3776 fhSumPtLeadingPt[icone]->Fill(ptC,coneptsum);
3780 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
3781 fhSumPtLeadingPtMC[kmcPhoton][icone]->Fill(ptC,coneptsum) ;
3783 fhSumPtLeadingPtMC[mcIndex][icone]->Fill(ptC,coneptsum) ;
3788 //Loop on pt thresholds
3789 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
3792 nfrac[icone][ipt]=0;
3793 GetIsolationCut()->SetPtThreshold(fPtThresholds[ipt]);
3794 GetIsolationCut()->SetPtFraction(fPtFractions[ipt]) ;
3795 GetIsolationCut()->SetSumPtThreshold(fSumPtThresholds[ipt]);
3797 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
3798 GetReader(), GetCaloPID(),
3800 n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);
3803 // Normal pT threshold cut
3807 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - cone size %1.1f, ptThres %1.1f, sumptThresh %1.1f\n",
3808 fConeSizes[icone],fPtThresholds[ipt],fSumPtThresholds[ipt]);
3809 printf("\t n %d, nfrac %d, coneptsum %2.2f\n",
3810 n[icone][ipt],nfrac[icone][ipt],coneptsum);
3812 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - pt %1.1f, eta %1.1f, phi %1.1f\n",ptC, etaC, phiC);
3815 if(n[icone][ipt] == 0)
3818 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling pt threshold loop\n");
3820 fhPtThresIsolated [icone][ipt]->Fill(ptC);
3821 fhEtaPhiPtThresIso[icone][ipt]->Fill(etaC,phiC);
3823 if( decay && fFillTaggedDecayHistograms )
3825 fhPtPtThresDecayIso [icone][ipt]->Fill(ptC);
3826 fhEtaPhiPtThresDecayIso[icone][ipt]->Fill(etaC,phiC);
3831 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) )
3832 fhPtThresIsolatedMC[kmcPhoton][icone][ipt]->Fill(ptC) ;
3834 fhPtThresIsolatedMC[mcIndex][icone][ipt]->Fill(ptC) ;
3839 // pt in cone fraction
3840 if(nfrac[icone][ipt] == 0)
3843 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling frac loop\n");
3845 fhPtFracIsolated [icone][ipt]->Fill(ptC);
3846 fhEtaPhiPtFracIso[icone][ipt]->Fill(etaC,phiC);
3848 if( decay && fFillTaggedDecayHistograms )
3850 fhPtPtFracDecayIso [icone][ipt]->Fill(ptC);
3851 fhEtaPhiPtFracDecayIso[icone][ipt]->Fill(etaC,phiC);
3856 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
3857 fhPtFracIsolatedMC[kmcPhoton][icone][ipt]->Fill(ptC) ;
3859 fhPtFracIsolatedMC[mcIndex][icone][ipt]->Fill(ptC) ;
3864 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - checking IC method : %i\n",GetIsolationCut()->GetICMethod());
3866 //Pt threshold on pt cand/ sum in cone histograms
3867 if(coneptsum < fSumPtThresholds[ipt])
3870 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling sum loop\n");
3872 fhSumPtIsolated [icone][ipt]->Fill(ptC) ;
3873 fhEtaPhiPtSumIso[icone][ipt]->Fill(etaC, phiC) ;
3875 if( decay && fFillTaggedDecayHistograms )
3877 fhPtPtSumDecayIso[icone][ipt]->Fill(ptC);
3878 fhEtaPhiPtSumDecayIso[icone][ipt]->Fill(etaC, phiC) ;
3883 if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
3884 fhSumPtIsolatedMC[kmcPhoton][icone][ipt]->Fill(ptC) ;
3886 fhSumPtIsolatedMC[mcIndex][icone][ipt]->Fill(ptC) ;
3890 // pt sum pt frac method
3891 // if( ((fPtFractions[ipt]*ptC < fSumPtThresholds[ipt]) && (coneptsum < fSumPtThresholds[ipt])) || ((fPtFractions[ipt]*ptC > fSumPtThresholds[ipt]) && (coneptsum < fPtFractions[ipt]*ptC)) )
3893 if(coneptsum < fPtFractions[ipt]*ptC)
3896 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling PtFrac PtSum loop\n");
3898 fhPtFracPtSumIso [icone][ipt]->Fill(ptC) ;
3899 fhEtaPhiFracPtSumIso[icone][ipt]->Fill(etaC,phiC) ;
3901 if( decay && fFillTaggedDecayHistograms )
3903 fhPtFracPtSumDecayIso [icone][ipt]->Fill(ptC);
3904 fhEtaPhiFracPtSumDecayIso[icone][ipt]->Fill(etaC,phiC);
3909 Float_t cellDensity = GetIsolationCut()->GetCellDensity( ph, GetReader());
3910 if(coneptsum < fSumPtThresholds[ipt]*cellDensity)
3913 printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling density loop\n");
3915 fhPtSumDensityIso [icone][ipt]->Fill(ptC) ;
3916 fhEtaPhiSumDensityIso[icone][ipt]->Fill(etaC,phiC) ;
3918 if( decay && fFillTaggedDecayHistograms )
3920 fhPtSumDensityDecayIso [icone][ipt]->Fill(ptC);
3921 fhEtaPhiSumDensityDecayIso[icone][ipt]->Fill(etaC, phiC);
3930 //Reset original parameters for AOD analysis
3931 GetIsolationCut()->SetPtThreshold(ptthresorg);
3932 GetIsolationCut()->SetPtFraction(ptfracorg);
3933 GetIsolationCut()->SetSumPtThreshold(ptsumcorg);
3934 GetIsolationCut()->SetConeSize(rorg);
3938 //_____________________________________________________________
3939 void AliAnaParticleIsolation::Print(const Option_t * opt) const
3942 //Print some relevant parameters set for the analysis
3946 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
3947 AliAnaCaloTrackCorrBaseClass::Print(" ");
3949 printf("ReMake Isolation = %d \n", fReMakeIC) ;
3950 printf("Make Several Isolation = %d \n", fMakeSeveralIC) ;
3951 printf("Calorimeter for isolation = %s \n", fCalorimeter.Data()) ;
3952 printf("Detector for candidate isolation = %s \n", fIsoDetector.Data()) ;
3956 printf("N Cone Sizes = %d\n", fNCones) ;
3957 printf("Cone Sizes = \n") ;
3958 for(Int_t i = 0; i < fNCones; i++)
3959 printf(" %1.2f;", fConeSizes[i]) ;
3962 printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ;
3963 printf(" pT thresholds = \n") ;
3964 for(Int_t i = 0; i < fNPtThresFrac; i++)
3965 printf(" %2.2f;", fPtThresholds[i]) ;
3969 printf(" pT fractions = \n") ;
3970 for(Int_t i = 0; i < fNPtThresFrac; i++)
3971 printf(" %2.2f;", fPtFractions[i]) ;
3975 printf("sum pT thresholds = \n") ;
3976 for(Int_t i = 0; i < fNPtThresFrac; i++)
3977 printf(" %2.2f;", fSumPtThresholds[i]) ;