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 // --- Analysis system ---
37 #include "AliAnaParticleIsolation.h"
38 #include "AliCaloTrackReader.h"
39 #include "AliIsolationCut.h"
40 #include "AliNeutralMesonSelection.h"
41 #include "AliAODPWG4ParticleCorrelation.h"
42 #include "AliMCAnalysisUtils.h"
43 #include "AliVTrack.h"
44 #include "AliVCluster.h"
45 #include "AliESDEvent.h"
46 #include "AliAODEvent.h"
48 ClassImp(AliAnaParticleIsolation)
50 //______________________________________________________________________________
51 AliAnaParticleIsolation::AliAnaParticleIsolation() :
52 AliAnaCaloTrackCorrBaseClass(), fCalorimeter(""),
53 fReMakeIC(0), fMakeSeveralIC(0),
54 fFillPileUpHistograms(0),
55 fFillTMHisto(0), fFillSSHisto(0),
57 fNCones(0), fNPtThresFrac(0),
58 fConeSizes(), fPtThresholds(),
59 fPtFractions(), fSumPtThresholds(),
61 fhEIso(0), fhPtIso(0),
62 fhPhiIso(0), fhEtaIso(0), fhEtaPhiIso(0),
64 fhPtNoIso(0), fhPtDecayIso(0), fhPtDecayNoIso(0),
65 fhEtaPhiDecayIso(0), fhEtaPhiDecayNoIso(0),
66 fhConeSumPt(0), fhPtInCone(0), fhPtInConeCent(0),
67 fhFRConeSumPt(0), fhPtInFRCone(0), fhPhiUEConeSumPt(0),
68 fhEtaUEConeSumPt(0), fhEtaBand(0), fhPhiBand(0),
69 fhConeSumPtEtaUESub(0), fhConeSumPtPhiUESub(0),
71 fhPtIsoPrompt(0), fhPhiIsoPrompt(0), fhEtaIsoPrompt(0),
72 fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(),
73 fhPtIsoFragmentation(0), fhPhiIsoFragmentation(0), fhEtaIsoFragmentation(0),
74 fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(),
75 fhPtIsoPi0(0), fhPhiIsoPi0(0), fhEtaIsoPi0(0),
76 fhPtThresIsolatedPi0(), fhPtFracIsolatedPi0(), fhPtSumIsolatedPi0(),
77 fhPtIsoPi0Decay(0), fhPhiIsoPi0Decay(0), fhEtaIsoPi0Decay(0),
78 fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(),
79 fhPtIsoEtaDecay(0), fhPhiIsoEtaDecay(0), fhEtaIsoEtaDecay(0),
80 fhPtThresIsolatedEtaDecay(), fhPtFracIsolatedEtaDecay(), fhPtSumIsolatedEtaDecay(),
81 fhPtIsoOtherDecay(0), fhPhiIsoOtherDecay(0), fhEtaIsoOtherDecay(0),
82 fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(),
83 //fhPtIsoConversion(0), fhPhiIsoConversion(0), fhEtaIsoConversion(0),
84 //fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(),
85 fhPtIsoHadron(0), fhPhiIsoHadron(0), fhEtaIsoHadron(0),
86 fhPtThresIsolatedHadron(), fhPtFracIsolatedHadron(), fhPtSumIsolatedHadron(),
87 fhPtNoIsoPi0(0), fhPtNoIsoPi0Decay(0),
88 fhPtNoIsoEtaDecay(0), fhPtNoIsoOtherDecay(0),
89 fhPtNoIsoPrompt(0), fhPtIsoMCPhoton(0), fhPtNoIsoMCPhoton(0),
90 //fhPtNoIsoConversion(0),
91 fhPtNoIsoFragmentation(0), fhPtNoIsoHadron(0),
93 fhSumPtLeadingPt(), fhPtLeadingPt(),
94 fhFRSumPtLeadingPt(), fhFRPtLeadingPt(),
95 fhPtThresIsolated(), fhPtFracIsolated(), fhPtSumIsolated(),
96 fhEtaPhiPtThresIso(), fhEtaPhiPtThresDecayIso(), fhPtPtThresDecayIso(),
97 fhEtaPhiPtFracIso(), fhEtaPhiPtFracDecayIso(), fhPtPtFracDecayIso(),
98 fhPtPtSumDecayIso(), fhEtaPhiSumDensityIso(), fhEtaPhiSumDensityDecayIso(),
99 fhPtSumDensityIso(), fhPtSumDensityDecayIso(),
100 fhPtFracPtSumIso(), fhPtFracPtSumDecayIso(),
101 fhEtaPhiFracPtSumIso(), fhEtaPhiFracPtSumDecayIso(),
102 // Cluster control histograms
103 fhTrackMatchedDEta(), fhTrackMatchedDPhi(), fhTrackMatchedDEtaDPhi(),
104 fhdEdx(), fhEOverP(), fhTrackMatchedMCParticle(),
105 fhELambda0() , fhELambda1(), fhELambda0SSBkg(),
106 fhELambda0TRD(), fhELambda1TRD(),
107 fhELambda0MCPhoton(), fhELambda0MCPi0(), fhELambda0MCPi0Decay(),
108 fhELambda0MCEtaDecay(), fhELambda0MCOtherDecay(), fhELambda0MCHadron(),
109 // Number of local maxima in cluster
111 fhELambda0LocMax1(), fhELambda1LocMax1(),
112 fhELambda0LocMax2(), fhELambda1LocMax2(),
113 fhELambda0LocMaxN(), fhELambda1LocMaxN(),
115 fhTimeENoCut(0), fhTimeESPD(0), fhTimeESPDMulti(0),
116 fhTimeNPileUpVertSPD(0), fhTimeNPileUpVertTrack(0),
117 fhTimeNPileUpVertContributors(0),
118 fhTimePileUpMainVertexZDistance(0), fhTimePileUpMainVertexZDiamond(0),
119 // Histograms settings
120 fHistoNPtSumBins(0), fHistoPtSumMax(0.), fHistoPtSumMin(0.),
121 fHistoNPtInConeBins(0), fHistoPtInConeMax(0.), fHistoPtInConeMin(0.)
125 //Initialize parameters
128 for(Int_t i = 0; i < 5 ; i++)
132 fhPtSumIsolatedPrompt [i] = 0 ;
133 fhPtSumIsolatedFragmentation[i] = 0 ;
134 fhPtSumIsolatedPi0Decay [i] = 0 ;
135 fhPtSumIsolatedPi0 [i] = 0 ;
136 fhPtSumIsolatedEtaDecay [i] = 0 ;
137 fhPtSumIsolatedOtherDecay [i] = 0 ;
138 // fhPtSumIsolatedConversion [i] = 0 ;
139 fhPtSumIsolatedHadron [i] = 0 ;
141 for(Int_t j = 0; j < 5 ; j++)
143 fhPtThresIsolated [i][j] = 0 ;
144 fhPtFracIsolated [i][j] = 0 ;
145 fhPtSumIsolated [i][j] = 0 ;
147 fhEtaPhiPtThresIso [i][j] = 0 ;
148 fhEtaPhiPtThresDecayIso[i][j] = 0 ;
149 fhPtPtThresDecayIso [i][j] = 0 ;
151 fhEtaPhiPtFracIso [i][j] = 0 ;
152 fhEtaPhiPtFracDecayIso [i][j] = 0 ;
153 fhPtPtFracDecayIso [i][j] = 0 ;
154 fhPtPtSumDecayIso [i][j] = 0 ;
155 fhPtSumDensityIso [i][j] = 0 ;
156 fhPtSumDensityDecayIso [i][j] = 0 ;
157 fhEtaPhiSumDensityIso [i][j] = 0 ;
158 fhEtaPhiSumDensityDecayIso [i][j] = 0 ;
159 fhPtFracPtSumIso [i][j] = 0 ;
160 fhPtFracPtSumDecayIso [i][j] = 0 ;
161 fhEtaPhiFracPtSumIso [i][j] = 0 ;
162 fhEtaPhiFracPtSumDecayIso [i][j] = 0 ;
164 fhPtThresIsolatedPrompt [i][j] = 0 ;
165 fhPtThresIsolatedFragmentation[i][j] = 0 ;
166 fhPtThresIsolatedPi0Decay [i][j] = 0 ;
167 fhPtThresIsolatedPi0 [i][j] = 0 ;
168 fhPtThresIsolatedEtaDecay [i][j] = 0 ;
169 fhPtThresIsolatedOtherDecay [i][j] = 0 ;
170 // fhPtThresIsolatedConversion [i][j] = 0 ;
171 fhPtThresIsolatedHadron [i][j] = 0 ;
173 fhPtFracIsolatedPrompt [i][j] = 0 ;
174 fhPtFracIsolatedFragmentation [i][j] = 0 ;
175 fhPtFracIsolatedPi0 [i][j] = 0 ;
176 fhPtFracIsolatedPi0Decay [i][j] = 0 ;
177 fhPtFracIsolatedEtaDecay [i][j] = 0 ;
178 fhPtFracIsolatedOtherDecay [i][j] = 0 ;
179 // fhPtFracIsolatedConversion [i][j] = 0 ;
180 fhPtFracIsolatedHadron [i][j] = 0 ;
185 for(Int_t i = 0; i < 5 ; i++)
187 fPtFractions [i] = 0 ;
188 fPtThresholds [i] = 0 ;
189 fSumPtThresholds[i] = 0 ;
193 for(Int_t i = 0; i < 2 ; i++)
195 fhTrackMatchedDEta[i] = 0 ; fhTrackMatchedDPhi[i] = 0 ; fhTrackMatchedDEtaDPhi [i] = 0 ;
196 fhdEdx [i] = 0 ; fhEOverP [i] = 0 ; fhTrackMatchedMCParticle[i] = 0 ;
197 fhELambda0 [i] = 0 ; fhELambda1 [i] = 0 ;
198 fhELambda0TRD [i] = 0 ; fhELambda1TRD [i] = 0 ;
200 fhELambda0MCPhoton [i] = 0 ; fhELambda0MCPi0 [i] = 0 ; fhELambda0MCPi0Decay[i] = 0 ;
201 fhELambda0MCEtaDecay[i] = 0 ; fhELambda0MCOtherDecay[i] = 0 ; fhELambda0MCHadron [i] = 0 ;
204 // Number of local maxima in cluster
206 fhELambda0LocMax1[i] = 0 ; fhELambda1LocMax1[i] = 0 ;
207 fhELambda0LocMax2[i] = 0 ; fhELambda1LocMax2[i] = 0 ;
208 fhELambda0LocMaxN[i] = 0 ; fhELambda1LocMaxN[i] = 0 ;
214 //_________________________________________________________________
215 void AliAnaParticleIsolation::FillPileUpHistograms(Int_t clusterID)
217 // Fill some histograms to understand pile-up
218 if(!fFillPileUpHistograms) return;
222 printf("AliAnaParticleIsolation::FillPileUpHistograms(), ID of cluster = %d, not possible! ", clusterID);
227 TObjArray* clusters = 0x0;
228 if (fCalorimeter == "EMCAL") clusters = GetEMCALClusters();
229 else if(fCalorimeter == "PHOS" ) clusters = GetPHOSClusters();
232 Float_t time = -1000;
236 AliVCluster *cluster = FindCluster(clusters,clusterID,iclus);
237 energy = cluster->E();
238 time = cluster->GetTOF()*1e9;
241 //printf("E %f, time %f\n",energy,time);
242 AliVEvent * event = GetReader()->GetInputEvent();
244 fhTimeENoCut->Fill(energy,time);
245 if(GetReader()->IsPileUpFromSPD()) fhTimeESPD ->Fill(energy,time);
246 if(event->IsPileupFromSPDInMultBins()) fhTimeESPDMulti->Fill(energy,time);
248 if(energy < 8) return; // Fill time figures for high energy clusters not too close to trigger threshold
250 AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
251 AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
253 // N pile up vertices
254 Int_t nVerticesSPD = -1;
255 Int_t nVerticesTracks = -1;
259 nVerticesSPD = esdEv->GetNumberOfPileupVerticesSPD();
260 nVerticesTracks = esdEv->GetNumberOfPileupVerticesTracks();
265 nVerticesSPD = aodEv->GetNumberOfPileupVerticesSPD();
266 nVerticesTracks = aodEv->GetNumberOfPileupVerticesTracks();
269 fhTimeNPileUpVertSPD ->Fill(time,nVerticesSPD);
270 fhTimeNPileUpVertTrack->Fill(time,nVerticesTracks);
272 //printf("Is SPD %d, Is SPD Multi %d, n spd %d, n track %d\n",
273 // GetReader()->IsPileUpFromSPD(),event->IsPileupFromSPDInMultBins(),nVerticesSPD,nVerticesTracks);
276 Float_t z1 = -1, z2 = -1;
278 for(Int_t iVert=0; iVert<nVerticesSPD;iVert++)
282 const AliESDVertex* pv=esdEv->GetPileupVertexSPD(iVert);
283 ncont=pv->GetNContributors();
284 z1 = esdEv->GetPrimaryVertexSPD()->GetZ();
286 diamZ = esdEv->GetDiamondZ();
290 AliAODVertex *pv=aodEv->GetVertex(iVert);
291 if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
292 ncont=pv->GetNContributors();
293 z1=aodEv->GetPrimaryVertexSPD()->GetZ();
295 diamZ = aodEv->GetDiamondZ();
298 Double_t distZ = TMath::Abs(z2-z1);
299 diamZ = TMath::Abs(z2-diamZ);
301 fhTimeNPileUpVertContributors ->Fill(time,ncont);
302 fhTimePileUpMainVertexZDistance->Fill(time,distZ);
303 fhTimePileUpMainVertexZDiamond ->Fill(time,diamZ);
308 //________________________________________________________________________________________________
309 void AliAnaParticleIsolation::FillTrackMatchingShowerShapeControlHistograms(const Bool_t isolated,
310 const Int_t clusterID,
313 const TObjArray * plCTS,
314 const TObjArray * plNe,
315 AliAODPWG4ParticleCorrelation *pCandidate,
316 const AliCaloTrackReader * reader,
317 const AliCaloPID * pid)
319 // Fill Track matching and Shower Shape control histograms
320 if(!fFillTMHisto && !fFillSSHisto) return;
324 printf("AliAnaParticleIsolation::FillTrackMatchingShowerShapeControlHistograms(), ID of cluster = %d, not possible! ", clusterID);
329 TObjArray* clusters = 0x0;
330 if (fCalorimeter == "EMCAL") clusters = GetEMCALClusters();
331 else if(fCalorimeter == "PHOS" ) clusters = GetPHOSClusters();
336 AliVCluster *cluster = FindCluster(clusters,clusterID,iclus);
337 Float_t energy = cluster->E();
341 fhELambda0[isolated]->Fill(energy, cluster->GetM02() );
342 fhELambda1[isolated]->Fill(energy, cluster->GetM20() );
346 if (GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt) ||
347 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation)) fhELambda0MCPhoton [isolated]->Fill(energy, cluster->GetM02());
348 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0)) fhELambda0MCPi0 [isolated]->Fill(energy, cluster->GetM02());
349 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay)) fhELambda0MCPi0Decay [isolated]->Fill(energy, cluster->GetM02());
350 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay)) fhELambda0MCEtaDecay [isolated]->Fill(energy, cluster->GetM02());
351 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay)) fhELambda0MCOtherDecay[isolated]->Fill(energy, cluster->GetM02());
353 // else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion)) fhPtNoIsoConversion ->Fill(energy, cluster->GetM02());
354 else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron)) fhELambda0MCHadron [isolated]->Fill(energy, cluster->GetM02());
358 if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5) // TO DO: CHANGE FOR 2012
360 fhELambda0TRD[isolated]->Fill(energy, cluster->GetM02() );
361 fhELambda1TRD[isolated]->Fill(energy, cluster->GetM20() );
364 fhNLocMax[isolated]->Fill(energy,nMaxima);
365 if (nMaxima==1) { fhELambda0LocMax1[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMax1[isolated]->Fill(energy,cluster->GetM20()); }
366 else if(nMaxima==2) { fhELambda0LocMax2[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMax2[isolated]->Fill(energy,cluster->GetM20()); }
367 else { fhELambda0LocMaxN[isolated]->Fill(energy,cluster->GetM02()); fhELambda1LocMaxN[isolated]->Fill(energy,cluster->GetM20()); }
371 //Analyse non-isolated events
374 Bool_t iso = kFALSE ;
375 Float_t coneptsum = 0 ;
376 GetIsolationCut()->SetPtThresholdMax(1.);
377 GetIsolationCut()->MakeIsolationCut(plCTS, plNe,
379 kFALSE, pCandidate, "",
380 n,nfrac,coneptsum, iso);
381 if (!iso) fhELambda0SSBkg->Fill(energy, cluster->GetM02());
384 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Energy Sum in Isolation Cone %2.2f\n", coneptsum);
386 GetIsolationCut()->SetPtThresholdMax(10000.);
393 Float_t dZ = cluster->GetTrackDz();
394 Float_t dR = cluster->GetTrackDx();
396 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
398 dR = 2000., dZ = 2000.;
399 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
402 //printf("ParticleIsolation: dPhi %f, dEta %f\n",dR,dZ);
403 if(fhTrackMatchedDEta && TMath::Abs(dR) < 999)
405 fhTrackMatchedDEta[isolated]->Fill(energy,dZ);
406 fhTrackMatchedDPhi[isolated]->Fill(energy,dR);
407 if(energy > 0.5) fhTrackMatchedDEtaDPhi[isolated]->Fill(dZ,dR);
410 // Check dEdx and E/p of matched clusters
412 if(TMath::Abs(dZ) < 0.05 && TMath::Abs(dR) < 0.05)
415 AliVTrack *track = GetCaloUtils()->GetMatchedTrack(cluster, GetReader()->GetInputEvent());
419 Float_t dEdx = track->GetTPCsignal();
420 fhdEdx[isolated]->Fill(cluster->E(), dEdx);
422 Float_t eOverp = cluster->E()/track->P();
423 fhEOverP[isolated]->Fill(cluster->E(), eOverp);
426 // printf("AliAnaParticleIsolation::FillTrackMatchingShowerShapeHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
431 if ( !GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion) )
433 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
434 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 2.5 );
435 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 0.5 );
436 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 1.5 );
437 else fhTrackMatchedMCParticle[isolated]->Fill(energy, 3.5 );
442 if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ||
443 GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 6.5 );
444 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 4.5 );
445 else if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) fhTrackMatchedMCParticle[isolated]->Fill(energy, 5.5 );
446 else fhTrackMatchedMCParticle[isolated]->Fill(energy, 7.5 );
455 } // clusters array available
459 //______________________________________________________
460 TObjString * AliAnaParticleIsolation::GetAnalysisCuts()
462 //Save parameters used for analysis
463 TString parList ; //this will be list of parameters used for this analysis.
464 const Int_t buffersize = 255;
465 char onePar[buffersize] ;
467 snprintf(onePar, buffersize,"--- AliAnaParticleIsolation ---\n") ;
469 snprintf(onePar, buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
471 snprintf(onePar, buffersize,"fReMakeIC =%d (Flag for reisolation during histogram filling) \n",fReMakeIC) ;
473 snprintf(onePar, buffersize,"fMakeSeveralIC=%d (Flag for isolation with several cuts at the same time ) \n",fMakeSeveralIC) ;
475 snprintf(onePar, buffersize,"fFillTMHisto=%d (Flag for track matching histograms) \n",fFillTMHisto) ;
477 snprintf(onePar, buffersize,"fFillSSHisto=%d (Flag for shower shape histograms) \n",fFillSSHisto) ;
482 snprintf(onePar, buffersize,"fNCones =%d (Number of cone sizes) \n",fNCones) ;
484 snprintf(onePar, buffersize,"fNPtThresFrac=%d (Flag for isolation with several cuts at the same time ) \n",fNPtThresFrac) ;
487 for(Int_t icone = 0; icone < fNCones ; icone++)
489 snprintf(onePar, buffersize,"fConeSizes[%d]=%1.2f (isolation cone size) \n",icone, fConeSizes[icone]) ;
492 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
494 snprintf(onePar, buffersize,"fPtThresholds[%d]=%1.2f (isolation pt threshold) \n",ipt, fPtThresholds[ipt]) ;
497 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
499 snprintf(onePar, buffersize,"fPtFractions[%d]=%1.2f (isolation pt fraction threshold) \n",ipt, fPtFractions[ipt]) ;
502 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++)
504 snprintf(onePar, buffersize,"fSumPtThresholds[%d]=%1.2f (isolation sum pt threshold) \n",ipt, fSumPtThresholds[ipt]) ;
509 //Get parameters set in base class.
510 parList += GetBaseParametersList() ;
512 //Get parameters set in IC class.
513 if(!fMakeSeveralIC)parList += GetIsolationCut()->GetICParametersList() ;
515 return new TObjString(parList) ;
519 //________________________________________________________
520 TList * AliAnaParticleIsolation::GetCreateOutputObjects()
522 // Create histograms to be saved in output file and
523 // store them in outputContainer
524 TList * outputContainer = new TList() ;
525 outputContainer->SetName("IsolatedParticleHistos") ;
527 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins();
528 Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins();
529 Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
530 Float_t ptmax = GetHistogramRanges()->GetHistoPtMax();
531 Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();
532 Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
533 Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
534 Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
535 Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
536 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins();
537 Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax();
538 Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
539 Int_t ntimebins= GetHistogramRanges()->GetHistoTimeBins();
540 Float_t timemax = GetHistogramRanges()->GetHistoTimeMax();
541 Float_t timemin = GetHistogramRanges()->GetHistoTimeMin();
543 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
544 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
545 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
546 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
547 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
548 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
550 Int_t ndedxbins = GetHistogramRanges()->GetHistodEdxBins();
551 Float_t dedxmax = GetHistogramRanges()->GetHistodEdxMax();
552 Float_t dedxmin = GetHistogramRanges()->GetHistodEdxMin();
553 Int_t nPoverEbins = GetHistogramRanges()->GetHistoPOverEBins();
554 Float_t pOverEmax = GetHistogramRanges()->GetHistoPOverEMax();
555 Float_t pOverEmin = GetHistogramRanges()->GetHistoPOverEMin();
557 Int_t nptsumbins = fHistoNPtSumBins;
558 Float_t ptsummax = fHistoPtSumMax;
559 Float_t ptsummin = fHistoPtSumMin;
560 Int_t nptinconebins = fHistoNPtInConeBins;
561 Float_t ptinconemax = fHistoPtInConeMax;
562 Float_t ptinconemin = fHistoPtInConeMin;
564 Float_t ptthre = GetIsolationCut()->GetPtThreshold();
565 Float_t ptfrac = GetIsolationCut()->GetPtFraction();
566 Float_t r = GetIsolationCut()->GetConeSize();
570 TString hName [] = {"NoIso",""};
571 TString hTitle[] = {"Not isolated" ,"isolated"};
574 fhELambda0SSBkg = new TH2F
575 ("hELambda0SSBkg","Non isolated clusters : E vs #lambda_{0}",nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
576 fhELambda0SSBkg->SetYTitle("#lambda_{0}^{2}");
577 fhELambda0SSBkg->SetXTitle("E (GeV)");
578 outputContainer->Add(fhELambda0SSBkg) ;
581 for(Int_t iso = 0; iso < 2; iso++)
585 fhTrackMatchedDEta[iso] = new TH2F
586 (Form("hTrackMatchedDEta%s",hName[iso].Data()),
587 Form("%s - d#eta of cluster-track vs cluster energy for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
588 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
589 fhTrackMatchedDEta[iso]->SetYTitle("d#eta");
590 fhTrackMatchedDEta[iso]->SetXTitle("E_{cluster} (GeV)");
592 fhTrackMatchedDPhi[iso] = new TH2F
593 (Form("hTrackMatchedDPhi%s",hName[iso].Data()),
594 Form("%s - d#phi of cluster-track vs cluster energy for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
595 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
596 fhTrackMatchedDPhi[iso]->SetYTitle("d#phi (rad)");
597 fhTrackMatchedDPhi[iso]->SetXTitle("E_{cluster} (GeV)");
599 fhTrackMatchedDEtaDPhi[iso] = new TH2F
600 (Form("hTrackMatchedDEtaDPhi%s",hName[iso].Data()),
601 Form("%s - d#eta vs d#phi of cluster-track for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
602 nresetabins,resetamin,resetamax,nresphibins,resphimin,resphimax);
603 fhTrackMatchedDEtaDPhi[iso]->SetYTitle("d#phi (rad)");
604 fhTrackMatchedDEtaDPhi[iso]->SetXTitle("d#eta");
606 outputContainer->Add(fhTrackMatchedDEta[iso]) ;
607 outputContainer->Add(fhTrackMatchedDPhi[iso]) ;
608 outputContainer->Add(fhTrackMatchedDEtaDPhi[iso]) ;
610 fhdEdx[iso] = new TH2F
611 (Form("hdEdx%s",hName[iso].Data()),
612 Form("%s - Matched track <dE/dx> vs cluster E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
613 nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
614 fhdEdx[iso]->SetXTitle("E (GeV)");
615 fhdEdx[iso]->SetYTitle("<dE/dx>");
616 outputContainer->Add(fhdEdx[iso]);
618 fhEOverP[iso] = new TH2F
619 (Form("hEOverP%s",hName[iso].Data()),
620 Form("%s - Matched track E/p vs cluster E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
621 nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
622 fhEOverP[iso]->SetXTitle("E (GeV)");
623 fhEOverP[iso]->SetYTitle("E/p");
624 outputContainer->Add(fhEOverP[iso]);
628 fhTrackMatchedMCParticle[iso] = new TH2F
629 (Form("hTrackMatchedMCParticle%s",hName[iso].Data()),
630 Form("%s - Origin of particle vs energy vs cluster E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",hTitle[iso].Data(),r,ptthre,ptfrac),
631 nptbins,ptmin,ptmax,8,0,8);
632 fhTrackMatchedMCParticle[iso]->SetXTitle("E (GeV)");
633 //fhTrackMatchedMCParticle[iso]->SetYTitle("Particle type");
635 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(1 ,"Photon");
636 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(2 ,"Electron");
637 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(3 ,"Meson Merged");
638 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(4 ,"Rest");
639 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(5 ,"Conv. Photon");
640 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(6 ,"Conv. Electron");
641 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(7 ,"Conv. Merged");
642 fhTrackMatchedMCParticle[iso]->GetYaxis()->SetBinLabel(8 ,"Conv. Rest");
644 outputContainer->Add(fhTrackMatchedMCParticle[iso]);
650 fhELambda0[iso] = new TH2F
651 (Form("hELambda0%s",hName[iso].Data()),
652 Form("%s cluster : E vs #lambda_{0}",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
653 fhELambda0[iso]->SetYTitle("#lambda_{0}^{2}");
654 fhELambda0[iso]->SetXTitle("E (GeV)");
655 outputContainer->Add(fhELambda0[iso]) ;
659 fhELambda0MCPhoton[iso] = new TH2F
660 (Form("hELambda0%s_MCPhoton",hName[iso].Data()),
661 Form("%s cluster : E vs #lambda_{0}: Origin is final state photon",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
662 fhELambda0MCPhoton[iso]->SetYTitle("#lambda_{0}^{2}");
663 fhELambda0MCPhoton[iso]->SetXTitle("E (GeV)");
664 outputContainer->Add(fhELambda0MCPhoton[iso]) ;
666 fhELambda0MCPi0[iso] = new TH2F
667 (Form("hELambda0%s_MCPi0",hName[iso].Data()),
668 Form("%s cluster : E vs #lambda_{0}: Origin is pi0 (2 #gamma)",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
669 fhELambda0MCPi0[iso]->SetYTitle("#lambda_{0}^{2}");
670 fhELambda0MCPi0[iso]->SetXTitle("E (GeV)");
671 outputContainer->Add(fhELambda0MCPi0[iso]) ;
673 fhELambda0MCPi0Decay[iso] = new TH2F
674 (Form("hELambda0%s_MCPi0Decay",hName[iso].Data()),
675 Form("%s cluster : E vs #lambda_{0}: Origin is pi0 decay",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
676 fhELambda0MCPi0Decay[iso]->SetYTitle("#lambda_{0}^{2}");
677 fhELambda0MCPi0Decay[iso]->SetXTitle("E (GeV)");
678 outputContainer->Add(fhELambda0MCPi0Decay[iso]) ;
680 fhELambda0MCEtaDecay[iso] = new TH2F
681 (Form("hELambda0%s_MCEtaDecay",hName[iso].Data()),
682 Form("%s cluster : E vs #lambda_{0}: Origin is eta decay",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
683 fhELambda0MCEtaDecay[iso]->SetYTitle("#lambda_{0}^{2}");
684 fhELambda0MCEtaDecay[iso]->SetXTitle("E (GeV)");
685 outputContainer->Add(fhELambda0MCEtaDecay[iso]) ;
687 fhELambda0MCOtherDecay[iso] = new TH2F
688 (Form("hELambda0%s_MCOtherDecay",hName[iso].Data()),
689 Form("%s cluster : E vs #lambda_{0}: Origin is other decay",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
690 fhELambda0MCOtherDecay[iso]->SetYTitle("#lambda_{0}^{2}");
691 fhELambda0MCOtherDecay[iso]->SetXTitle("E (GeV)");
692 outputContainer->Add(fhELambda0MCOtherDecay[iso]) ;
694 fhELambda0MCHadron[iso] = new TH2F
695 (Form("hELambda0%s_MCHadron",hName[iso].Data()),
696 Form("%s cluster : E vs #lambda_{0}: Origin is hadron",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
697 fhELambda0MCHadron[iso]->SetYTitle("#lambda_{0}^{2}");
698 fhELambda0MCHadron[iso]->SetXTitle("E (GeV)");
699 outputContainer->Add(fhELambda0MCHadron[iso]) ;
702 fhELambda1[iso] = new TH2F
703 (Form("hELambda1%s",hName[iso].Data()),
704 Form("%s cluster: E vs #lambda_{1}",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
705 fhELambda1[iso]->SetYTitle("#lambda_{1}^{2}");
706 fhELambda1[iso]->SetXTitle("E (GeV)");
707 outputContainer->Add(fhELambda1[iso]) ;
709 if(fCalorimeter=="EMCAL")
711 fhELambda0TRD[iso] = new TH2F
712 (Form("hELambda0TRD%s",hName[iso].Data()),
713 Form("%s cluster: E vs #lambda_{0}, SM behind TRD",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
714 fhELambda0TRD[iso]->SetYTitle("#lambda_{0}^{2}");
715 fhELambda0TRD[iso]->SetXTitle("E (GeV)");
716 outputContainer->Add(fhELambda0TRD[iso]) ;
718 fhELambda1TRD[iso] = new TH2F
719 (Form("hELambda1TRD%s",hName[iso].Data()),
720 Form("%s cluster: E vs #lambda_{1}, SM behind TRD",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
721 fhELambda1TRD[iso]->SetYTitle("#lambda_{1}^{2}");
722 fhELambda1TRD[iso]->SetXTitle("E (GeV)");
723 outputContainer->Add(fhELambda1TRD[iso]) ;
726 fhNLocMax[iso] = new TH2F
727 (Form("hNLocMax%s",hName[iso].Data()),
728 Form("%s - Number of local maxima in cluster",hTitle[iso].Data()),
729 nptbins,ptmin,ptmax,10,0,10);
730 fhNLocMax[iso]->SetYTitle("N maxima");
731 fhNLocMax[iso]->SetXTitle("E (GeV)");
732 outputContainer->Add(fhNLocMax[iso]) ;
734 fhELambda0LocMax1[iso] = new TH2F
735 (Form("hELambda0LocMax1%s",hName[iso].Data()),
736 Form("%s cluster (#eta) pairs: E vs #lambda_{0}, 1 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
737 fhELambda0LocMax1[iso]->SetYTitle("#lambda_{0}^{2}");
738 fhELambda0LocMax1[iso]->SetXTitle("E (GeV)");
739 outputContainer->Add(fhELambda0LocMax1[iso]) ;
741 fhELambda1LocMax1[iso] = new TH2F
742 (Form("hELambda1LocMax1%s",hName[iso].Data()),
743 Form("%s cluster (#eta) pairs: E vs #lambda_{1}, 1 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
744 fhELambda1LocMax1[iso]->SetYTitle("#lambda_{1}^{2}");
745 fhELambda1LocMax1[iso]->SetXTitle("E (GeV)");
746 outputContainer->Add(fhELambda1LocMax1[iso]) ;
748 fhELambda0LocMax2[iso] = new TH2F
749 (Form("hELambda0LocMax2%s",hName[iso].Data()),
750 Form("%s cluster (#eta) pairs: E vs #lambda_{0}, 2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
751 fhELambda0LocMax2[iso]->SetYTitle("#lambda_{0}^{2}");
752 fhELambda0LocMax2[iso]->SetXTitle("E (GeV)");
753 outputContainer->Add(fhELambda0LocMax2[iso]) ;
755 fhELambda1LocMax2[iso] = new TH2F
756 (Form("hELambda1LocMax2%s",hName[iso].Data()),
757 Form("%s cluster (#eta) pairs: E vs #lambda_{1}, 2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
758 fhELambda1LocMax2[iso]->SetYTitle("#lambda_{1}^{2}");
759 fhELambda1LocMax2[iso]->SetXTitle("E (GeV)");
760 outputContainer->Add(fhELambda1LocMax2[iso]) ;
762 fhELambda0LocMaxN[iso] = new TH2F
763 ( Form("hELambda0LocMaxN%s",hName[iso].Data()),
764 Form("%s cluster (#eta) pairs: E vs #lambda_{0}, N>2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
765 fhELambda0LocMaxN[iso]->SetYTitle("#lambda_{0}^{2}");
766 fhELambda0LocMaxN[iso]->SetXTitle("E (GeV)");
767 outputContainer->Add(fhELambda0LocMaxN[iso]) ;
769 fhELambda1LocMaxN[iso] = new TH2F
770 (Form("hELambda1LocMaxN%s",hName[iso].Data()),
771 Form("%s cluster (#eta) pairs: E vs #lambda_{1}, N>2 Local maxima",hTitle[iso].Data()),nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
772 fhELambda1LocMaxN[iso]->SetYTitle("#lambda_{1}^{2}");
773 fhELambda1LocMaxN[iso]->SetXTitle("E (GeV)");
774 outputContainer->Add(fhELambda1LocMaxN[iso]) ;
777 } // control histograms for isolated and non isolated objects
779 fhConeSumPt = new TH2F("hConePtSum",
780 Form("#Sigma p_{T} in isolation cone for R = %2.2f",r),
781 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
782 fhConeSumPt->SetYTitle("#Sigma p_{T}");
783 fhConeSumPt->SetXTitle("p_{T} (GeV/c)");
784 outputContainer->Add(fhConeSumPt) ;
786 fhPtInCone = new TH2F("hPtInCone",
787 Form("p_{T} in isolation cone for R = %2.2f",r),
788 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
789 fhPtInCone->SetYTitle("p_{T in cone} (GeV/c)");
790 fhPtInCone->SetXTitle("p_{T} (GeV/c)");
791 outputContainer->Add(fhPtInCone) ;
793 fhPtInConeCent = new TH2F("hPtInConeCent",
794 Form("p_{T} in isolation cone for R = %2.2f",r),
795 100,0,100,nptinconebins,ptinconemin,ptinconemax);
796 fhPtInConeCent->SetYTitle("p_{T in cone} (GeV/c)");
797 fhPtInConeCent->SetXTitle("centrality");
798 outputContainer->Add(fhPtInConeCent) ;
800 fhFRConeSumPt = new TH2F("hFRConePtSum",
801 Form("#Sigma p_{T} in the forward region isolation cone for R = %2.2f",r),
802 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
803 fhFRConeSumPt->SetYTitle("#Sigma p_{T}");
804 fhFRConeSumPt->SetXTitle("p_{T} (GeV/c)");
805 outputContainer->Add(fhFRConeSumPt) ;
807 fhPtInFRCone = new TH2F("hPtInFRCone",
808 Form("p_{T} in forward region isolation cone for R = %2.2f",r),
809 nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
810 fhPtInFRCone->SetYTitle("p_{T in cone} (GeV/c)");
811 fhPtInFRCone->SetXTitle("p_{T} (GeV/c)");
812 outputContainer->Add(fhPtInFRCone) ;
814 fhPhiUEConeSumPt = new TH2F("hPhiUEConeSumPt",
815 Form("p_{T} in phi band around isolation cone for R = %2.2f",r),
816 nptbins,ptmin,ptmax,4*nptinconebins,ptinconemin,4*ptinconemax);
817 fhPhiUEConeSumPt->SetYTitle("p_{T in band} (GeV/c)");
818 fhPhiUEConeSumPt->SetXTitle("p_{T} (GeV/c)");
819 outputContainer->Add(fhPhiUEConeSumPt) ;
821 fhEtaUEConeSumPt = new TH2F("hEtaUEConeSumPt",
822 Form("p_{T} in eta band around isolation cone for R = %2.2f",r),
823 nptbins,ptmin,ptmax,4*nptinconebins,ptinconemin,4*ptinconemax);
824 fhEtaUEConeSumPt->SetYTitle("p_{T in band} (GeV/c)");
825 fhEtaUEConeSumPt->SetXTitle("p_{T} (GeV/c)");
826 outputContainer->Add(fhEtaUEConeSumPt) ;
828 fhEtaBand = new TH2F("fhEtaBand",
829 Form("Eta/Phi of particle in Eta band isolation cone for R = %2.2f",r),
830 netabins,etamin,etamax,nphibins,phimin,phimax);
831 fhEtaBand->SetXTitle("#eta");
832 fhEtaBand->SetYTitle("#phi");
833 outputContainer->Add(fhEtaBand) ;
835 fhPhiBand = new TH2F("fhPhiBand",
836 Form("Eta/Phi of particle in Phi band isolation cone for R = %2.2f",r),
837 netabins,etamin,etamax,nphibins,phimin,phimax);
838 fhPhiBand->SetXTitle("#eta");
839 fhPhiBand->SetYTitle("#phi");
840 outputContainer->Add(fhPhiBand) ;
842 fhConeSumPtEtaUESub = new TH2F("hConeSumPtEtaUESub",
843 Form("#Sigma p_{T} after bkg subtraction from eta band in the isolation cone for R = %2.2f",r),
844 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
845 fhConeSumPtEtaUESub->SetYTitle("#Sigma p_{T}");
846 fhConeSumPtEtaUESub->SetXTitle("p_{T} (GeV/c)");
847 outputContainer->Add(fhConeSumPtEtaUESub) ;
849 fhConeSumPtPhiUESub = new TH2F("hConeSumPtPhiUESub",
850 Form("#Sigma p_{T} after bkg subtraction from phi band in the isolation cone for R = %2.2f",r),
851 nptbins,ptmin,ptmax,2*nptsumbins,-ptsummax,ptsummax);
852 fhConeSumPtPhiUESub->SetYTitle("#Sigma p_{T}");
853 fhConeSumPtPhiUESub->SetXTitle("p_{T} (GeV/c)");
854 outputContainer->Add(fhConeSumPtPhiUESub) ;
856 fhEIso = new TH1F("hE",
857 Form("Number of isolated particles vs E for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
858 nptbins,ptmin,ptmax);
859 fhEIso->SetYTitle("dN / dE");
860 fhEIso->SetXTitle("E (GeV/c)");
861 outputContainer->Add(fhEIso) ;
863 fhPtIso = new TH1F("hPt",
864 Form("Number of isolated particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
865 nptbins,ptmin,ptmax);
866 fhPtIso->SetYTitle("dN / p_{T}");
867 fhPtIso->SetXTitle("p_{T} (GeV/c)");
868 outputContainer->Add(fhPtIso) ;
870 fhPhiIso = new TH2F("hPhi",
871 Form("Number of isolated particles vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
872 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
873 fhPhiIso->SetYTitle("#phi");
874 fhPhiIso->SetXTitle("p_{T} (GeV/c)");
875 outputContainer->Add(fhPhiIso) ;
877 fhEtaIso = new TH2F("hEta",
878 Form("Number of isolated particles vs #eta for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
879 nptbins,ptmin,ptmax,netabins,etamin,etamax);
880 fhEtaIso->SetYTitle("#eta");
881 fhEtaIso->SetXTitle("p_{T} (GeV/c)");
882 outputContainer->Add(fhEtaIso) ;
884 fhEtaPhiIso = new TH2F("hEtaPhiIso",
885 Form("Number of isolated particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
886 netabins,etamin,etamax,nphibins,phimin,phimax);
887 fhEtaPhiIso->SetXTitle("#eta");
888 fhEtaPhiIso->SetYTitle("#phi");
889 outputContainer->Add(fhEtaPhiIso) ;
891 fhPtDecayIso = new TH1F("hPtDecayIso",
892 Form("Number of isolated #pi^{0} decay particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
893 nptbins,ptmin,ptmax);
894 fhPtDecayIso->SetYTitle("N");
895 fhPtDecayIso->SetXTitle("p_{T}(GeV/c)");
896 outputContainer->Add(fhPtDecayIso) ;
898 fhEtaPhiDecayIso = new TH2F("hEtaPhiDecayIso",
899 Form("Number of isolated Pi0 decay particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
900 netabins,etamin,etamax,nphibins,phimin,phimax);
901 fhEtaPhiDecayIso->SetXTitle("#eta");
902 fhEtaPhiDecayIso->SetYTitle("#phi");
903 outputContainer->Add(fhEtaPhiDecayIso) ;
907 fhPtIsoPrompt = new TH1F("hPtMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax);
908 fhPtIsoPrompt->SetYTitle("N");
909 fhPtIsoPrompt->SetXTitle("p_{T #gamma}(GeV/c)");
910 outputContainer->Add(fhPtIsoPrompt) ;
912 fhPhiIsoPrompt = new TH2F
913 ("hPhiMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
914 fhPhiIsoPrompt->SetYTitle("#phi");
915 fhPhiIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
916 outputContainer->Add(fhPhiIsoPrompt) ;
918 fhEtaIsoPrompt = new TH2F
919 ("hEtaMCPrompt","Number of isolated prompt #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
920 fhEtaIsoPrompt->SetYTitle("#eta");
921 fhEtaIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
922 outputContainer->Add(fhEtaIsoPrompt) ;
924 fhPtIsoFragmentation = new TH1F("hPtMCFragmentation","Number of isolated #gamma",nptbins,ptmin,ptmax);
925 fhPtIsoFragmentation->SetYTitle("N");
926 fhPtIsoFragmentation->SetXTitle("p_{T #gamma}(GeV/c)");
927 outputContainer->Add(fhPtIsoFragmentation) ;
929 fhPhiIsoFragmentation = new TH2F
930 ("hPhiMCFragmentation","Number of isolated fragmentation #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
931 fhPhiIsoFragmentation->SetYTitle("#phi");
932 fhPhiIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
933 outputContainer->Add(fhPhiIsoFragmentation) ;
935 fhEtaIsoFragmentation = new TH2F
936 ("hEtaMCFragmentation","Number of isolated fragmentation #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
937 fhEtaIsoFragmentation->SetYTitle("#eta");
938 fhEtaIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
939 outputContainer->Add(fhEtaIsoFragmentation) ;
941 fhPtIsoPi0 = new TH1F("hPtMCPi0","Number of isolated #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax);
942 fhPtIsoPi0->SetYTitle("N");
943 fhPtIsoPi0->SetXTitle("p_{T #gamma}(GeV/c)");
944 outputContainer->Add(fhPtIsoPi0) ;
946 fhPhiIsoPi0 = new TH2F
947 ("hPhiMCPi0","Number of isolated #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
948 fhPhiIsoPi0->SetYTitle("#phi");
949 fhPhiIsoPi0->SetXTitle("p_{T #gamma} (GeV/c)");
950 outputContainer->Add(fhPhiIsoPi0) ;
952 fhEtaIsoPi0 = new TH2F
953 ("hEtaMCPi0","Number of isolated #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax,netabins,etamin,etamax);
954 fhEtaIsoPi0->SetYTitle("#eta");
955 fhEtaIsoPi0->SetXTitle("p_{T #gamma} (GeV/c)");
956 outputContainer->Add(fhEtaIsoPi0) ;
958 fhPtIsoPi0Decay = new TH1F("hPtMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax);
959 fhPtIsoPi0Decay->SetYTitle("N");
960 fhPtIsoPi0Decay->SetXTitle("p_{T #gamma}(GeV/c)");
961 outputContainer->Add(fhPtIsoPi0Decay) ;
963 fhPhiIsoPi0Decay = new TH2F
964 ("hPhiMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
965 fhPhiIsoPi0Decay->SetYTitle("#phi");
966 fhPhiIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
967 outputContainer->Add(fhPhiIsoPi0Decay) ;
969 fhEtaIsoPi0Decay = new TH2F
970 ("hEtaMCPi0Decay","Number of isolated #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
971 fhEtaIsoPi0Decay->SetYTitle("#eta");
972 fhEtaIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
973 outputContainer->Add(fhEtaIsoPi0Decay) ;
975 fhPtIsoEtaDecay = new TH1F("hPtMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax);
976 fhPtIsoEtaDecay->SetYTitle("N");
977 fhPtIsoEtaDecay->SetXTitle("p_{T #gamma}(GeV/c)");
978 outputContainer->Add(fhPtIsoEtaDecay) ;
980 fhPhiIsoEtaDecay = new TH2F
981 ("hPhiMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
982 fhPhiIsoEtaDecay->SetYTitle("#phi");
983 fhPhiIsoEtaDecay->SetXTitle("p_{T #gamma} (GeV/c)");
984 outputContainer->Add(fhPhiIsoEtaDecay) ;
986 fhEtaIsoEtaDecay = new TH2F
987 ("hEtaMCEtaDecay","Number of isolated #gamma from #eta decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
988 fhEtaIsoEtaDecay->SetYTitle("#eta");
989 fhEtaIsoEtaDecay->SetXTitle("p_{T #gamma} (GeV/c)");
990 outputContainer->Add(fhEtaIsoEtaDecay) ;
992 fhPtIsoOtherDecay = new TH1F("hPtMCOtherDecay","Number of isolated #gamma from non-#pi^{0} decay",nptbins,ptmin,ptmax);
993 fhPtIsoOtherDecay->SetYTitle("N");
994 fhPtIsoOtherDecay->SetXTitle("p_{T #gamma}(GeV/c)");
995 outputContainer->Add(fhPtIsoOtherDecay) ;
997 fhPhiIsoOtherDecay = new TH2F
998 ("hPhiMCOtherDecay","Number of isolated #gamma from non-#pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
999 fhPhiIsoOtherDecay->SetYTitle("#phi");
1000 fhPhiIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
1001 outputContainer->Add(fhPhiIsoOtherDecay) ;
1003 fhEtaIsoOtherDecay = new TH2F
1004 ("hEtaMCOtherDecay","Number of isolated #gamma non-#pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1005 fhEtaIsoOtherDecay->SetYTitle("#eta");
1006 fhEtaIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
1007 outputContainer->Add(fhEtaIsoOtherDecay) ;
1009 // fhPtIsoConversion = new TH1F("hPtMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax);
1010 // fhPtIsoConversion->SetYTitle("N");
1011 // fhPtIsoConversion->SetXTitle("p_{T #gamma}(GeV/c)");
1012 // outputContainer->Add(fhPtIsoConversion) ;
1014 // fhPhiIsoConversion = new TH2F
1015 // ("hPhiMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1016 // fhPhiIsoConversion->SetYTitle("#phi");
1017 // fhPhiIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
1018 // outputContainer->Add(fhPhiIsoConversion) ;
1020 // fhEtaIsoConversion = new TH2F
1021 // ("hEtaMCConversion","Number of isolated converted #gamma",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1022 // fhEtaIsoConversion->SetYTitle("#eta");
1023 // fhEtaIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
1024 // outputContainer->Add(fhEtaIsoConversion) ;
1026 fhPtIsoHadron = new TH1F("hPtMCHadron","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax);
1027 fhPtIsoHadron->SetYTitle("N");
1028 fhPtIsoHadron->SetXTitle("p_{T}(GeV/c)");
1029 outputContainer->Add(fhPtIsoHadron) ;
1031 fhPhiIsoHadron = new TH2F
1032 ("hPhiMCHadron","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1033 fhPhiIsoHadron->SetYTitle("#phi");
1034 fhPhiIsoHadron->SetXTitle("p_{T} (GeV/c)");
1035 outputContainer->Add(fhPhiIsoHadron) ;
1037 fhEtaIsoHadron = new TH2F
1038 ("hEtaMCHadron","Number of isolated non-#gamma particles",nptbins,ptmin,ptmax,netabins,etamin,etamax);
1039 fhEtaIsoHadron->SetYTitle("#eta");
1040 fhEtaIsoHadron->SetXTitle("p_{T} (GeV/c)");
1041 outputContainer->Add(fhEtaIsoHadron) ;
1047 // Not Isolated histograms, reference histograms
1049 fhPtNoIso = new TH1F("hPtNoIso",
1050 Form("Number of not isolated leading particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1051 nptbins,ptmin,ptmax);
1052 fhPtNoIso->SetYTitle("N");
1053 fhPtNoIso->SetXTitle("p_{T}(GeV/c)");
1054 outputContainer->Add(fhPtNoIso) ;
1057 fhEtaPhiNoIso = new TH2F("hEtaPhiNoIso",
1058 Form("Number of not isolated leading particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1059 netabins,etamin,etamax,nphibins,phimin,phimax);
1060 fhEtaPhiNoIso->SetXTitle("#eta");
1061 fhEtaPhiNoIso->SetYTitle("#phi");
1062 outputContainer->Add(fhEtaPhiNoIso) ;
1064 fhPtDecayNoIso = new TH1F("hPtDecayNoIso",
1065 Form("Number of not isolated leading pi0 decay particles vs p_{T} for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1066 nptbins,ptmin,ptmax);
1067 fhPtDecayNoIso->SetYTitle("N");
1068 fhPtDecayNoIso->SetXTitle("p_{T}(GeV/c)");
1069 outputContainer->Add(fhPtDecayNoIso) ;
1071 fhEtaPhiDecayNoIso = new TH2F("hEtaPhiDecayNoIso",
1072 Form("Number of not isolated leading Pi0 decay particles #eta vs #phi for R = %2.2f, p_{T}^{th} = %2.2f, p_{T}^{fr} = %2.2f",r,ptthre,ptfrac),
1073 netabins,etamin,etamax,nphibins,phimin,phimax);
1074 fhEtaPhiDecayNoIso->SetXTitle("#eta");
1075 fhEtaPhiDecayNoIso->SetYTitle("#phi");
1076 outputContainer->Add(fhEtaPhiDecayNoIso) ;
1082 fhPtNoIsoPi0 = new TH1F
1083 ("hPtNoIsoPi0","Number of not isolated leading #gamma from #pi^{0} (2 #gamma)",nptbins,ptmin,ptmax);
1084 fhPtNoIsoPi0->SetYTitle("N");
1085 fhPtNoIsoPi0->SetXTitle("p_{T} (GeV/c)");
1086 outputContainer->Add(fhPtNoIsoPi0) ;
1088 fhPtNoIsoPi0Decay = new TH1F
1089 ("hPtNoIsoPi0Decay","Number of not isolated leading #gamma from #pi^{0} decay",nptbins,ptmin,ptmax);
1090 fhPtNoIsoPi0Decay->SetYTitle("N");
1091 fhPtNoIsoPi0Decay->SetXTitle("p_{T} (GeV/c)");
1092 outputContainer->Add(fhPtNoIsoPi0Decay) ;
1094 fhPtNoIsoEtaDecay = new TH1F
1095 ("hPtNoIsoEtaDecay","Number of not isolated leading #gamma from eta decay",nptbins,ptmin,ptmax);
1096 fhPtNoIsoEtaDecay->SetYTitle("N");
1097 fhPtNoIsoEtaDecay->SetXTitle("p_{T} (GeV/c)");
1098 outputContainer->Add(fhPtNoIsoEtaDecay) ;
1100 fhPtNoIsoOtherDecay = new TH1F
1101 ("hPtNoIsoOtherDecay","Number of not isolated leading #gamma from other decay",nptbins,ptmin,ptmax);
1102 fhPtNoIsoOtherDecay->SetYTitle("N");
1103 fhPtNoIsoOtherDecay->SetXTitle("p_{T} (GeV/c)");
1104 outputContainer->Add(fhPtNoIsoOtherDecay) ;
1106 fhPtNoIsoPrompt = new TH1F
1107 ("hPtNoIsoPrompt","Number of not isolated leading prompt #gamma",nptbins,ptmin,ptmax);
1108 fhPtNoIsoPrompt->SetYTitle("N");
1109 fhPtNoIsoPrompt->SetXTitle("p_{T} (GeV/c)");
1110 outputContainer->Add(fhPtNoIsoPrompt) ;
1112 fhPtIsoMCPhoton = new TH1F
1113 ("hPtIsoMCPhoton","Number of isolated leading #gamma",nptbins,ptmin,ptmax);
1114 fhPtIsoMCPhoton->SetYTitle("N");
1115 fhPtIsoMCPhoton->SetXTitle("p_{T} (GeV/c)");
1116 outputContainer->Add(fhPtIsoMCPhoton) ;
1118 fhPtNoIsoMCPhoton = new TH1F
1119 ("hPtNoIsoMCPhoton","Number of not isolated leading #gamma",nptbins,ptmin,ptmax);
1120 fhPtNoIsoMCPhoton->SetYTitle("N");
1121 fhPtNoIsoMCPhoton->SetXTitle("p_{T} (GeV/c)");
1122 outputContainer->Add(fhPtNoIsoMCPhoton) ;
1124 // fhPtNoIsoConversion = new TH1F
1125 // ("hPtNoIsoConversion","Number of not isolated leading conversion #gamma",nptbins,ptmin,ptmax);
1126 // fhPtNoIsoConversion->SetYTitle("N");
1127 // fhPtNoIsoConversion->SetXTitle("p_{T} (GeV/c)");
1128 // outputContainer->Add(fhPtNoIsoConversion) ;
1130 fhPtNoIsoFragmentation = new TH1F
1131 ("hPtNoIsoFragmentation","Number of not isolated leading fragmentation #gamma",nptbins,ptmin,ptmax);
1132 fhPtNoIsoFragmentation->SetYTitle("N");
1133 fhPtNoIsoFragmentation->SetXTitle("p_{T} (GeV/c)");
1134 outputContainer->Add(fhPtNoIsoFragmentation) ;
1136 fhPtNoIsoHadron = new TH1F
1137 ("hPtNoIsoHadron","Number of not isolated leading hadrons",nptbins,ptmin,ptmax);
1138 fhPtNoIsoHadron->SetYTitle("N");
1139 fhPtNoIsoHadron->SetXTitle("p_{T} (GeV/c)");
1140 outputContainer->Add(fhPtNoIsoHadron) ;
1147 const Int_t buffersize = 255;
1148 char name[buffersize];
1149 char title[buffersize];
1150 for(Int_t icone = 0; icone<fNCones; icone++)
1152 // sum pt in cone vs. pt leading
1153 snprintf(name, buffersize,"hSumPtLeadingPt_Cone_%d",icone);
1154 snprintf(title, buffersize,"#Sigma p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1155 fhSumPtLeadingPt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1156 fhSumPtLeadingPt[icone] ->SetYTitle("#sum_{cone}p_{T} (GeV/c)");//#Sigma p_{T}
1157 fhSumPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1158 outputContainer->Add(fhSumPtLeadingPt[icone]) ;
1160 // pt in cone vs. pt leading
1161 snprintf(name, buffersize,"hPtLeadingPt_Cone_%d",icone);
1162 snprintf(title, buffersize,"p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1163 fhPtLeadingPt[icone] = new TH2F(name, title, nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1164 fhPtLeadingPt[icone] ->SetYTitle("p_{T}^{cone} (GeV/c)");
1165 fhPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1166 outputContainer->Add(fhPtLeadingPt[icone]) ;
1168 // sum pt in cone vs. pt leading in the forward region (for background subtraction studies)
1169 snprintf(name, buffersize,"hFRSumPtLeadingPt_Cone_%d",icone);
1170 snprintf(title, buffersize,"#Sigma p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1171 fhFRSumPtLeadingPt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1172 fhFRSumPtLeadingPt[icone] ->SetYTitle("#sum_{cone}p_{T} (GeV/c)");//#Sigma p_{T}
1173 fhFRSumPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1174 outputContainer->Add(fhFRSumPtLeadingPt[icone]) ;
1176 // pt in cone vs. pt leading in the forward region (for background subtraction studies)
1177 snprintf(name, buffersize,"hFRPtLeadingPt_Cone_%d",icone);
1178 snprintf(title, buffersize,"p_{T} in isolation cone for R = %2.2f",fConeSizes[icone]);
1179 fhFRPtLeadingPt[icone] = new TH2F(name, title, nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
1180 fhFRPtLeadingPt[icone] ->SetYTitle("p_{T}^{cone} (GeV/c)");
1181 fhFRPtLeadingPt[icone] ->SetXTitle("p_{T}^{leading} (GeV/c)");
1182 outputContainer->Add(fhFRPtLeadingPt[icone]) ;
1187 snprintf(name, buffersize,"hPtSumPrompt_Cone_%d",icone);
1188 snprintf(title, buffersize,"Candidate Prompt cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1189 fhPtSumIsolatedPrompt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1190 fhPtSumIsolatedPrompt[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1191 fhPtSumIsolatedPrompt[icone]->SetXTitle("p_{T} (GeV/c)");
1192 outputContainer->Add(fhPtSumIsolatedPrompt[icone]) ;
1194 snprintf(name, buffersize,"hPtSumFragmentation_Cone_%d",icone);
1195 snprintf(title, buffersize,"Candidate Fragmentation cone sum p_{T} for R = %2.2fvs candidate p_{T}",fConeSizes[icone]);
1196 fhPtSumIsolatedFragmentation[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1197 fhPtSumIsolatedFragmentation[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1198 fhPtSumIsolatedFragmentation[icone]->SetXTitle("p_{T} (GeV/c)");
1199 outputContainer->Add(fhPtSumIsolatedFragmentation[icone]) ;
1201 snprintf(name, buffersize,"hPtSumPi0_Cone_%d",icone);
1202 snprintf(title, buffersize,"Candidate Pi0 cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1203 fhPtSumIsolatedPi0[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1204 fhPtSumIsolatedPi0[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1205 fhPtSumIsolatedPi0[icone]->SetXTitle("p_{T} (GeV/c)");
1206 outputContainer->Add(fhPtSumIsolatedPi0[icone]) ;
1208 snprintf(name, buffersize,"hPtSumPi0Decay_Cone_%d",icone);
1209 snprintf(title, buffersize,"Candidate Pi0Decay cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1210 fhPtSumIsolatedPi0Decay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1211 fhPtSumIsolatedPi0Decay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1212 fhPtSumIsolatedPi0Decay[icone]->SetXTitle("p_{T} (GeV/c)");
1213 outputContainer->Add(fhPtSumIsolatedPi0Decay[icone]) ;
1215 snprintf(name, buffersize,"hPtSumEtaDecay_Cone_%d",icone);
1216 snprintf(title, buffersize,"Candidate EtaDecay cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1217 fhPtSumIsolatedEtaDecay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1218 fhPtSumIsolatedEtaDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1219 fhPtSumIsolatedEtaDecay[icone]->SetXTitle("p_{T} (GeV/c)");
1220 outputContainer->Add(fhPtSumIsolatedEtaDecay[icone]) ;
1222 snprintf(name, buffersize,"hPtSumOtherDecay_Cone_%d",icone);
1223 snprintf(title, buffersize,"Candidate OtherDecay cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1224 fhPtSumIsolatedOtherDecay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1225 fhPtSumIsolatedOtherDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1226 fhPtSumIsolatedOtherDecay[icone]->SetXTitle("p_{T} (GeV/c)");
1227 outputContainer->Add(fhPtSumIsolatedOtherDecay[icone]) ;
1229 // snprintf(name, buffersize,"hPtSumConversion_Cone_%d",icone);
1230 // snprintf(title, buffersize,"Candidate Conversion cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1231 // fhPtSumIsolatedConversion[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1232 // fhPtSumIsolatedConversion[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1233 // fhPtSumIsolatedConversion[icone]->SetXTitle("p_{T} (GeV/c)");
1234 // outputContainer->Add(fhPtSumIsolatedConversion[icone]) ;
1236 snprintf(name, buffersize,"hPtSumHadron_Cone_%d",icone);
1237 snprintf(title, buffersize,"Candidate Hadron cone sum p_{T} for R = %2.2f vs candidate p_{T}",fConeSizes[icone]);
1238 fhPtSumIsolatedHadron[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
1239 fhPtSumIsolatedHadron[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1240 fhPtSumIsolatedHadron[icone]->SetXTitle("p_{T} (GeV/c)");
1241 outputContainer->Add(fhPtSumIsolatedHadron[icone]) ;
1245 for(Int_t ipt = 0; ipt<fNPtThresFrac;ipt++)
1248 snprintf(name, buffersize,"hPtThres_Cone_%d_Pt%d",icone,ipt);
1249 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1250 fhPtThresIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1251 fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1252 outputContainer->Add(fhPtThresIsolated[icone][ipt]) ;
1254 snprintf(name, buffersize,"hPtFrac_Cone_%d_Pt%d",icone,ipt);
1255 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1256 fhPtFracIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1257 fhPtFracIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1258 outputContainer->Add(fhPtFracIsolated[icone][ipt]) ;
1261 snprintf(name, buffersize,"hPtSum_Cone_%d_Pt%d",icone,ipt);
1262 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1263 fhPtSumIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1264 // fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma p_{T} (GeV/c)");
1265 fhPtSumIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1266 outputContainer->Add(fhPtSumIsolated[icone][ipt]) ;
1268 snprintf(name, buffersize,"hPtSumDensity_Cone_%d_Pt%d",icone,ipt);
1269 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for density in R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1270 fhPtSumDensityIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1271 //fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma p_{T} (GeV/c)");
1272 fhPtSumDensityIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1273 outputContainer->Add(fhPtSumDensityIso[icone][ipt]) ;
1275 snprintf(name, buffersize,"hPtFracPtSum_Cone_%d_Pt%d",icone,ipt);
1276 snprintf(title, buffersize,"Isolated candidate p_{T} distribution for PtFracPtSum in R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1277 fhPtFracPtSumIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1278 //fhPtSumIsolated[icone][ipt]->SetYTitle("#Sigma p_{T} (GeV/c)");
1279 fhPtFracPtSumIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1280 outputContainer->Add(fhPtFracPtSumIso[icone][ipt]) ;
1282 // pt decays isolated
1283 snprintf(name, buffersize,"hPtThres_Decay_Cone_%d_Pt%d",icone,ipt);
1284 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1285 fhPtPtThresDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1286 fhPtPtThresDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1287 outputContainer->Add(fhPtPtThresDecayIso[icone][ipt]) ;
1289 snprintf(name, buffersize,"hPtFrac_Decay_Cone_%d_Pt%d",icone,ipt);
1290 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1291 fhPtPtFracDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1292 fhPtPtFracDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1293 outputContainer->Add(fhPtPtFracDecayIso[icone][ipt]) ;
1295 snprintf(name, buffersize,"hPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1296 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1297 fhPtPtSumDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1298 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1299 fhPtPtSumDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1300 outputContainer->Add(fhPtPtSumDecayIso[icone][ipt]) ;
1302 snprintf(name, buffersize,"hPtSumDensity_Decay_Cone_%d_Pt%d",icone,ipt);
1303 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for density in R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1304 fhPtSumDensityDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1305 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1306 fhPtSumDensityDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1307 outputContainer->Add(fhPtSumDensityDecayIso[icone][ipt]) ;
1309 snprintf(name, buffersize,"hPtFracPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1310 snprintf(title, buffersize,"Isolated decay candidate p_{T} distribution for PtFracPtSum in R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1311 fhPtFracPtSumDecayIso[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);//,nptsumbins,ptsummin,ptsummax);
1312 // fhPtPtSumDecayIso[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
1313 fhPtFracPtSumDecayIso[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1314 outputContainer->Add(fhPtFracPtSumDecayIso[icone][ipt]) ;
1318 snprintf(name, buffersize,"hEtaPhiPtThres_Cone_%d_Pt%d",icone,ipt);
1319 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1320 fhEtaPhiPtThresIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1321 fhEtaPhiPtThresIso[icone][ipt]->SetXTitle("#eta");
1322 fhEtaPhiPtThresIso[icone][ipt]->SetYTitle("#phi");
1323 outputContainer->Add(fhEtaPhiPtThresIso[icone][ipt]) ;
1325 snprintf(name, buffersize,"hEtaPhiPtFrac_Cone_%d_Pt%d",icone,ipt);
1326 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1327 fhEtaPhiPtFracIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1328 fhEtaPhiPtFracIso[icone][ipt]->SetXTitle("#eta");
1329 fhEtaPhiPtFracIso[icone][ipt]->SetYTitle("#phi");
1330 outputContainer->Add(fhEtaPhiPtFracIso[icone][ipt]) ;
1332 snprintf(name, buffersize,"hEtaPhiPtSum_Cone_%d_Pt%d",icone,ipt);
1333 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1334 fhEtaPhiPtSumIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1335 fhEtaPhiPtSumIso[icone][ipt]->SetXTitle("#eta");
1336 fhEtaPhiPtSumIso[icone][ipt]->SetYTitle("#phi");
1337 outputContainer->Add(fhEtaPhiPtSumIso[icone][ipt]) ;
1339 snprintf(name, buffersize,"hEtaPhiSumDensity_Cone_%d_Pt%d",icone,ipt);
1340 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for density R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1341 fhEtaPhiSumDensityIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1342 fhEtaPhiSumDensityIso[icone][ipt]->SetXTitle("#eta");
1343 fhEtaPhiSumDensityIso[icone][ipt]->SetYTitle("#phi");
1344 outputContainer->Add(fhEtaPhiSumDensityIso[icone][ipt]) ;
1346 snprintf(name, buffersize,"hEtaPhiFracPtSum_Cone_%d_Pt%d",icone,ipt);
1347 snprintf(title, buffersize,"Isolated candidate #eta:#phi distribution for FracPtSum R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1348 fhEtaPhiFracPtSumIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1349 fhEtaPhiFracPtSumIso[icone][ipt]->SetXTitle("#eta");
1350 fhEtaPhiFracPtSumIso[icone][ipt]->SetYTitle("#phi");
1351 outputContainer->Add(fhEtaPhiFracPtSumIso[icone][ipt]) ;
1354 snprintf(name, buffersize,"hEtaPhiPtThres_Decay_Cone_%d_Pt%d",icone,ipt);
1355 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for R = %2.2f and p_{T}^{th} = %2.2f GeV/c",fConeSizes[icone],fPtThresholds[ipt]);
1356 fhEtaPhiPtThresDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1357 fhEtaPhiPtThresDecayIso[icone][ipt]->SetXTitle("#eta");
1358 fhEtaPhiPtThresDecayIso[icone][ipt]->SetYTitle("#phi");
1359 outputContainer->Add(fhEtaPhiPtThresDecayIso[icone][ipt]) ;
1361 snprintf(name, buffersize,"hEtaPhiPtFrac_Decay_Cone_%d_Pt%d",icone,ipt);
1362 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1363 fhEtaPhiPtFracDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1364 fhEtaPhiPtFracDecayIso[icone][ipt]->SetXTitle("#eta");
1365 fhEtaPhiPtFracDecayIso[icone][ipt]->SetYTitle("#phi");
1366 outputContainer->Add(fhEtaPhiPtFracDecayIso[icone][ipt]) ;
1369 snprintf(name, buffersize,"hEtaPhiPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1370 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1371 fhEtaPhiPtSumDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1372 fhEtaPhiPtSumDecayIso[icone][ipt]->SetXTitle("#eta");
1373 fhEtaPhiPtSumDecayIso[icone][ipt]->SetYTitle("#phi");
1374 outputContainer->Add(fhEtaPhiPtSumDecayIso[icone][ipt]) ;
1376 snprintf(name, buffersize,"hEtaPhiSumDensity_Decay_Cone_%d_Pt%d",icone,ipt);
1377 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for density R = %2.2f and p_{T}^{sum} = %2.2f GeV/c",fConeSizes[icone],fSumPtThresholds[ipt]);
1378 fhEtaPhiSumDensityDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1379 fhEtaPhiSumDensityDecayIso[icone][ipt]->SetXTitle("#eta");
1380 fhEtaPhiSumDensityDecayIso[icone][ipt]->SetYTitle("#phi");
1381 outputContainer->Add(fhEtaPhiSumDensityDecayIso[icone][ipt]) ;
1383 snprintf(name, buffersize,"hEtaPhiFracPtSum_Decay_Cone_%d_Pt%d",icone,ipt);
1384 snprintf(title, buffersize,"Isolated decay candidate #eta:#phi distribution for FracPtSum R = %2.2f and p_{T}^{fr} = %2.2f GeV/c",fConeSizes[icone],fPtFractions[ipt]);
1385 fhEtaPhiFracPtSumDecayIso[icone][ipt] = new TH2F(name, title,netabins,etamin,etamax,nphibins,phimin,phimax);
1386 fhEtaPhiFracPtSumDecayIso[icone][ipt]->SetXTitle("#eta");
1387 fhEtaPhiFracPtSumDecayIso[icone][ipt]->SetYTitle("#phi");
1388 outputContainer->Add(fhEtaPhiFracPtSumDecayIso[icone][ipt]) ;
1393 snprintf(name, buffersize,"hPtThresMCPrompt_Cone_%d_Pt%d",icone,ipt);
1394 snprintf(title, buffersize,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1395 fhPtThresIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1396 fhPtThresIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1397 outputContainer->Add(fhPtThresIsolatedPrompt[icone][ipt]) ;
1399 snprintf(name, buffersize,"hPtFracMCPrompt_Cone_%d_Pt%d",icone,ipt);
1400 snprintf(title, buffersize,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1401 fhPtFracIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1402 fhPtFracIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1403 outputContainer->Add(fhPtFracIsolatedPrompt[icone][ipt]) ;
1405 snprintf(name, buffersize,"hPtThresMCFragmentation_Cone_%d_Pt%d",icone,ipt);
1406 snprintf(title, buffersize,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1407 fhPtThresIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1408 fhPtThresIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1409 outputContainer->Add(fhPtThresIsolatedFragmentation[icone][ipt]) ;
1411 snprintf(name, buffersize,"hPtFracMCFragmentation_Cone_%d_Pt%d",icone,ipt);
1412 snprintf(title, buffersize,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1413 fhPtFracIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1414 fhPtFracIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1415 outputContainer->Add(fhPtFracIsolatedFragmentation[icone][ipt]) ;
1417 snprintf(name, buffersize,"hPtThresMCPi0_Cone_%d_Pt%d",icone,ipt);
1418 snprintf(title, buffersize,"Isolated candidate Pi0 p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1419 fhPtThresIsolatedPi0[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1420 fhPtThresIsolatedPi0[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1421 outputContainer->Add(fhPtThresIsolatedPi0[icone][ipt]) ;
1423 snprintf(name, buffersize,"hPtFracMCPi0_Cone_%d_Pt%d",icone,ipt);
1424 snprintf(title, buffersize,"Isolated candidate Pi0 p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1425 fhPtFracIsolatedPi0[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1426 fhPtFracIsolatedPi0[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1427 outputContainer->Add(fhPtFracIsolatedPi0[icone][ipt]) ;
1429 snprintf(name, buffersize,"hPtThresMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
1430 snprintf(title, buffersize,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1431 fhPtThresIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1432 fhPtThresIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1433 outputContainer->Add(fhPtThresIsolatedPi0Decay[icone][ipt]) ;
1435 snprintf(name, buffersize,"hPtFracMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
1436 snprintf(title, buffersize,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1437 fhPtFracIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1438 fhPtFracIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1439 outputContainer->Add(fhPtFracIsolatedPi0Decay[icone][ipt]) ;
1441 snprintf(name, buffersize,"hPtThresMCEtaDecay_Cone_%d_Pt%d",icone,ipt);
1442 snprintf(title, buffersize,"Isolated candidate EtaDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1443 fhPtThresIsolatedEtaDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1444 fhPtThresIsolatedEtaDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1445 outputContainer->Add(fhPtThresIsolatedEtaDecay[icone][ipt]) ;
1447 snprintf(name, buffersize,"hPtFracMCEtaDecay_Cone_%d_Pt%d",icone,ipt);
1448 snprintf(title, buffersize,"Isolated candidate EtaDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1449 fhPtFracIsolatedEtaDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1450 fhPtFracIsolatedEtaDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1451 outputContainer->Add(fhPtFracIsolatedEtaDecay[icone][ipt]) ;
1454 snprintf(name, buffersize,"hPtThresMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
1455 snprintf(title, buffersize,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1456 fhPtThresIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1457 fhPtThresIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1458 outputContainer->Add(fhPtThresIsolatedOtherDecay[icone][ipt]) ;
1460 snprintf(name, buffersize,"hPtFracMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
1461 snprintf(title, buffersize,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1462 fhPtFracIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1463 fhPtFracIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1464 outputContainer->Add(fhPtFracIsolatedOtherDecay[icone][ipt]) ;
1466 // snprintf(name, buffersize,"hPtThresMCConversion_Cone_%d_Pt%d",icone,ipt);
1467 // snprintf(title, buffersize,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1468 // fhPtThresIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1469 // fhPtThresIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1470 // outputContainer->Add(fhPtThresIsolatedConversion[icone][ipt]) ;
1472 // snprintf(name, buffersize,"hPtFracMCConversion_Cone_%d_Pt%d",icone,ipt);
1473 // snprintf(title, buffersize,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1474 // fhPtFracIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1475 // fhPtFracIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1476 // outputContainer->Add(fhPtFracIsolatedConversion[icone][ipt]) ;
1478 snprintf(name, buffersize,"hPtThresMCHadron_Cone_%d_Pt%d",icone,ipt);
1479 snprintf(title, buffersize,"Isolated candidate Hadron p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1480 fhPtThresIsolatedHadron[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1481 fhPtThresIsolatedHadron[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1482 outputContainer->Add(fhPtThresIsolatedHadron[icone][ipt]) ;
1484 snprintf(name, buffersize,"hPtFracMCHadron_Cone_%d_Pt%d",icone,ipt);
1485 snprintf(title, buffersize,"Isolated candidate Hadron p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
1486 fhPtFracIsolatedHadron[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
1487 fhPtFracIsolatedHadron[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
1488 outputContainer->Add(fhPtFracIsolatedHadron[icone][ipt]) ;
1495 if(fFillPileUpHistograms)
1497 fhTimeENoCut = new TH2F ("hTimeE_NoCut","time of cluster vs E of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1498 fhTimeENoCut->SetXTitle("E (GeV)");
1499 fhTimeENoCut->SetYTitle("time (ns)");
1500 outputContainer->Add(fhTimeENoCut);
1502 fhTimeESPD = new TH2F ("hTimeE_SPD","time of cluster vs E of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1503 fhTimeESPD->SetXTitle("E (GeV)");
1504 fhTimeESPD->SetYTitle("time (ns)");
1505 outputContainer->Add(fhTimeESPD);
1507 fhTimeESPDMulti = new TH2F ("hTimeE_SPDMulti","time of cluster vs E of clusters, SPD multi cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
1508 fhTimeESPDMulti->SetXTitle("E (GeV)");
1509 fhTimeESPDMulti->SetYTitle("time (ns)");
1510 outputContainer->Add(fhTimeESPDMulti);
1512 fhTimeNPileUpVertSPD = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
1513 fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
1514 fhTimeNPileUpVertSPD->SetXTitle("time (ns)");
1515 outputContainer->Add(fhTimeNPileUpVertSPD);
1517 fhTimeNPileUpVertTrack = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 50,0,50 );
1518 fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
1519 fhTimeNPileUpVertTrack->SetXTitle("time (ns)");
1520 outputContainer->Add(fhTimeNPileUpVertTrack);
1522 fhTimeNPileUpVertContributors = new TH2F ("hTime_NPileUpVertContributors","time of cluster vs N constributors to pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
1523 fhTimeNPileUpVertContributors->SetYTitle("# vertex ");
1524 fhTimeNPileUpVertContributors->SetXTitle("time (ns)");
1525 outputContainer->Add(fhTimeNPileUpVertContributors);
1527 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);
1528 fhTimePileUpMainVertexZDistance->SetYTitle("distance Z (cm) ");
1529 fhTimePileUpMainVertexZDistance->SetXTitle("time (ns)");
1530 outputContainer->Add(fhTimePileUpMainVertexZDistance);
1532 fhTimePileUpMainVertexZDiamond = new TH2F ("hTime_PileUpMainVertexZDiamond","time of cluster vs distance in Z pile-up SPD vertex - z diamond",ntimebins,timemin,timemax,100,0,50);
1533 fhTimePileUpMainVertexZDiamond->SetYTitle("diamond distance Z (cm) ");
1534 fhTimePileUpMainVertexZDiamond->SetXTitle("time (ns)");
1535 outputContainer->Add(fhTimePileUpMainVertexZDiamond);
1539 return outputContainer ;
1543 //__________________________________
1544 void AliAnaParticleIsolation::Init()
1546 // Do some checks and init stuff
1548 // In case of several cone and thresholds analysis, open the cuts for the filling of the
1549 // track and cluster reference arrays in cone when done in the MakeAnalysisFillAOD().
1550 // The different cones, thresholds are tested for this list of tracks, clusters.
1553 printf("AliAnaParticleIsolation::Init() - Open default isolation cuts for multiple Isolation analysis\n");
1554 GetIsolationCut()->SetPtThreshold(100);
1555 GetIsolationCut()->SetPtFraction(100);
1556 GetIsolationCut()->SetConeSize(1);
1560 //____________________________________________
1561 void AliAnaParticleIsolation::InitParameters()
1564 //Initialize the parameters of the analysis.
1565 SetInputAODName("PWG4Particle");
1566 SetAODObjArrayName("IsolationCone");
1567 AddToHistogramsName("AnaIsolation_");
1569 fCalorimeter = "PHOS" ;
1570 fReMakeIC = kFALSE ;
1571 fMakeSeveralIC = kFALSE ;
1573 //----------- Several IC-----------------
1576 fConeSizes [0] = 0.1; fConeSizes [1] = 0.2; fConeSizes [2] = 0.3; fConeSizes [3] = 0.4; fConeSizes [4] = 0.5;
1577 fPtThresholds [0] = 1.; fPtThresholds [1] = 2.; fPtThresholds [2] = 3.; fPtThresholds [3] = 4.; fPtThresholds [4] = 5.;
1578 fPtFractions [0] = 0.05; fPtFractions [1] = 0.075; fPtFractions [2] = 0.1; fPtFractions [3] = 1.25; fPtFractions [4] = 1.5;
1579 fSumPtThresholds[0] = 1.; fSumPtThresholds[1] = 2.; fSumPtThresholds[2] = 3.; fSumPtThresholds[3] = 4.; fSumPtThresholds[4] = 5.;
1581 //------------- Histograms settings -------
1582 fHistoNPtSumBins = 100 ;
1583 fHistoPtSumMax = 50 ;
1584 fHistoPtSumMin = 0. ;
1586 fHistoNPtInConeBins = 100 ;
1587 fHistoPtInConeMax = 50 ;
1588 fHistoPtInConeMin = 0. ;
1592 //__________________________________________________
1593 void AliAnaParticleIsolation::MakeAnalysisFillAOD()
1595 //Do analysis and fill aods
1596 //Search for the isolated photon in fCalorimeter with pt > GetMinPt()
1598 if(!GetInputAODBranch())
1600 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
1604 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation"))
1606 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
1610 Int_t n = 0, nfrac = 0;
1611 Bool_t isolated = kFALSE ;
1612 Float_t coneptsum = 0 ;
1613 TObjArray * pl = 0x0; ;
1615 //Select the calorimeter for candidate isolation with neutral particles
1616 if (fCalorimeter == "PHOS" )
1617 pl = GetPHOSClusters();
1618 else if (fCalorimeter == "EMCAL")
1619 pl = GetEMCALClusters();
1621 //Loop on AOD branch, filled previously in AliAnaPhoton, find leading particle to do isolation only with it
1622 Double_t ptLeading = 0. ;
1623 Int_t idLeading = -1 ;
1624 TLorentzVector mom ;
1625 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1626 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Input aod branch entries %d\n", naod);
1628 for(Int_t iaod = 0; iaod < naod; iaod++)
1630 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1632 //If too small or too large pt, skip
1633 if(aodinput->Pt() < GetMinPt() || aodinput->Pt() > GetMaxPt() ) continue ;
1635 //check if it is low pt trigger particle
1636 if((aodinput->Pt() < GetIsolationCut()->GetPtThreshold() ||
1637 aodinput->Pt() < GetIsolationCut()->GetSumPtThreshold()) &&
1640 continue ; //trigger should not come from underlying event
1643 //vertex cut in case of mixing
1644 Int_t check = CheckMixedEventVertex(aodinput->GetCaloLabel(0), aodinput->GetTrackLabel(0));
1645 if(check == 0) continue;
1646 if(check == -1) return;
1648 //find the leading particles with highest momentum
1649 if ( aodinput->Pt() > ptLeading )
1651 ptLeading = aodinput->Pt() ;
1655 aodinput->SetLeadingParticle(kFALSE);
1657 }//finish searching for leading trigger particle
1659 // Check isolation of leading particle
1660 if(idLeading < 0) return;
1662 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(idLeading));
1663 aodinput->SetLeadingParticle(kTRUE);
1665 // Check isolation only of clusters in fiducial region
1666 if(IsFiducialCutOn())
1668 Bool_t in = GetFiducialCut()->IsInFiducialCut(*aodinput->Momentum(),aodinput->GetDetector()) ;
1672 //After cuts, study isolation
1673 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
1674 GetIsolationCut()->MakeIsolationCut(GetCTSTracks(),pl,
1675 GetReader(), GetCaloPID(),
1676 kTRUE, aodinput, GetAODObjArrayName(),
1677 n,nfrac,coneptsum, isolated);
1679 if(!fMakeSeveralIC) aodinput->SetIsolated(isolated);
1683 if(isolated)printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() : Particle %d IS ISOLATED \n",idLeading);
1684 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - End fill AODs \n");
1689 //_________________________________________________________
1690 void AliAnaParticleIsolation::MakeAnalysisFillHistograms()
1692 //Do analysis and fill histograms
1694 Int_t n = 0, nfrac = 0;
1695 Bool_t isolated = kFALSE ;
1696 Float_t coneptsum = 0 ;
1697 Float_t etaUEptsum = 0 ;
1698 Float_t phiUEptsum = 0 ;
1700 //Loop on stored AOD
1701 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1702 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Histo aod branch entries %d\n", naod);
1704 //Get vertex for photon momentum calculation
1705 Double_t vertex[] = {0,0,0} ; //vertex ;
1706 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1708 GetReader()->GetVertex(vertex);
1711 for(Int_t iaod = 0; iaod < naod ; iaod++)
1713 AliAODPWG4ParticleCorrelation* aod = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1715 if(!aod->IsLeadingParticle()) continue; // Try to isolate only leading cluster or track
1717 // Check isolation only of clusters in fiducial region
1718 if(IsFiducialCutOn())
1720 Bool_t in = GetFiducialCut()->IsInFiducialCut(*aod->Momentum(),aod->GetDetector()) ;
1721 if(! in ) continue ;
1724 Bool_t isolation = aod->IsIsolated();
1725 Bool_t decay = aod->IsTagged();
1726 Float_t energy = aod->E();
1727 Float_t pt = aod->Pt();
1728 Float_t phi = aod->Phi();
1729 Float_t eta = aod->Eta();
1730 Float_t conesize = GetIsolationCut()->GetConeSize();
1732 //Recover reference arrays with clusters and tracks
1733 TObjArray * refclusters = aod->GetObjArray(GetAODObjArrayName()+"Clusters");
1734 TObjArray * reftracks = aod->GetObjArray(GetAODObjArrayName()+"Tracks");
1736 //If too small or too large pt, skip
1737 if(pt < GetMinPt() || pt > GetMaxPt() ) continue ;
1739 // --- In case of redoing isolation from delta AOD ----
1743 //Analysis of multiple IC at same time
1744 MakeSeveralICAnalysis(aod);
1749 //In case a more strict IC is needed in the produced AOD
1750 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
1751 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
1752 GetReader(), GetCaloPID(),
1754 n,nfrac,coneptsum, isolated);
1755 fhConeSumPt->Fill(pt,coneptsum);
1756 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Energy Sum in Isolation Cone %2.2f\n", coneptsum);
1759 // ---------------------------------------------------
1761 //Fill pt distribution of particles in cone
1764 Double_t sumptFR = 0. ;
1765 TObjArray * trackList = GetCTSTracks() ;
1766 for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++)
1768 AliVTrack* track = (AliVTrack *) trackList->At(itrack);
1772 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Track not available?");
1776 //fill histogram for UE in phi band
1777 if(track->Eta() > (eta-conesize) && track->Eta() < (eta+conesize))
1779 phiUEptsum+=track->Pt();
1780 fhPhiBand->Fill(track->Eta(),track->Phi());
1782 //fill histogram for UE in eta band in EMCal acceptance
1783 if(track->Phi() > (phi-conesize) && track->Phi() < (phi+conesize) && track->Eta() > -0.6 && track->Eta() < 0.6)
1785 etaUEptsum+=track->Pt();
1786 fhEtaBand->Fill(track->Eta(),track->Phi());
1790 //fill the histograms at forward range
1791 Double_t dPhi = phi - track->Phi() + TMath::PiOver2();
1792 Double_t dEta = eta - track->Eta();
1793 Double_t arg = dPhi*dPhi + dEta*dEta;
1794 if(TMath::Sqrt(arg) < conesize)
1796 fhPtInFRCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1797 sumptFR+=track->Pt();
1800 dPhi = phi - track->Phi() - TMath::PiOver2();
1801 arg = dPhi*dPhi + dEta*dEta;
1802 if(TMath::Sqrt(arg) < conesize)
1804 fhPtInFRCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1805 sumptFR+=track->Pt();
1809 fhFRConeSumPt->Fill(pt,sumptFR);
1812 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++)
1814 AliVTrack* track = (AliVTrack *) reftracks->At(itrack);
1815 fhPtInCone->Fill(pt,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1816 if (GetEventCentrality()) fhPtInConeCent->Fill(GetEventCentrality(),TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
1817 coneptsum+=track->Pt();
1824 TLorentzVector mom ;
1825 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++)
1827 AliVCluster* calo = (AliVCluster *) refclusters->At(icalo);
1828 calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
1830 //fill histogram for UE in phi band
1831 if(mom.Eta() > (eta-conesize) && mom.Eta() < (eta+conesize))
1833 phiUEptsum+=mom.Pt();
1834 fhPhiBand->Fill(mom.Eta(),mom.Phi());
1836 //fill histogram for UE in eta band in EMCal acceptance
1837 if(mom.Phi() > (phi-conesize) && mom.Phi() < (phi+conesize))
1839 etaUEptsum+=mom.Pt();
1840 fhEtaBand->Fill(mom.Eta(),mom.Phi());
1843 fhPtInCone->Fill(pt, mom.Pt());
1844 if (GetEventCentrality()) fhPtInConeCent->Fill(GetEventCentrality(),mom.Pt());
1845 coneptsum+=mom.Pt();
1849 //normalize phi/eta band per area unit
1850 fhPhiUEConeSumPt->Fill(pt, phiUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*TMath::Pi()));
1851 fhEtaUEConeSumPt->Fill(pt, etaUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*0.6));
1853 Double_t SumPhiUESub = coneptsum-(phiUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*TMath::Pi()));
1854 Double_t SumEtaUESub = coneptsum-(etaUEptsum*(TMath::Pi()*conesize*conesize)/(2*conesize*2*0.6));
1856 fhConeSumPtPhiUESub->Fill(pt,SumPhiUESub);
1857 fhConeSumPtEtaUESub->Fill(pt,SumEtaUESub);
1860 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d Energy Sum in Isolation Cone %2.2f\n", iaod, coneptsum);
1862 if(!fReMakeIC) fhConeSumPt->Fill(pt,coneptsum);
1864 Int_t mcTag = aod->GetTag() ;
1865 Int_t clID = aod->GetCaloLabel(0) ;
1867 if(GetDebug() > 0) printf(" AliAnaParticleIsolation::MakeAnalysisFillHistograms() - pt %1.1f, eta %1.1f, phi %1.1f\n",pt, eta, phi);
1869 FillTrackMatchingShowerShapeControlHistograms(isolation, clID,aod->GetFiducialArea(),mcTag,reftracks,refclusters,aod,GetReader(), GetCaloPID());
1873 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d ISOLATED: fill histograms\n", iaod);
1875 // Fill histograms to undertand pile-up before other cuts applied
1876 // Remember to relax time cuts in the reader
1877 FillPileUpHistograms(clID);
1879 fhEIso ->Fill(energy);
1881 fhPhiIso ->Fill(pt,phi);
1882 fhEtaIso ->Fill(pt,eta);
1883 fhEtaPhiIso ->Fill(eta,phi);
1887 fhPtDecayIso->Fill(pt);
1888 fhEtaPhiDecayIso->Fill(eta,phi);
1894 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
1896 fhPtIsoMCPhoton ->Fill(pt);
1899 if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt))
1901 fhPtIsoPrompt ->Fill(pt);
1902 fhPhiIsoPrompt ->Fill(pt,phi);
1903 fhEtaIsoPrompt ->Fill(pt,eta);
1905 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation))
1907 fhPtIsoFragmentation ->Fill(pt);
1908 fhPhiIsoFragmentation ->Fill(pt,phi);
1909 fhEtaIsoFragmentation ->Fill(pt,eta);
1911 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0))
1913 fhPtIsoPi0 ->Fill(pt);
1914 fhPhiIsoPi0 ->Fill(pt,phi);
1915 fhEtaIsoPi0 ->Fill(pt,eta);
1917 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay))
1919 fhPtIsoPi0Decay ->Fill(pt);
1920 fhPhiIsoPi0Decay ->Fill(pt,phi);
1921 fhEtaIsoPi0Decay ->Fill(pt,eta);
1923 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay))
1925 fhPtIsoEtaDecay ->Fill(pt);
1926 fhPhiIsoEtaDecay ->Fill(pt,phi);
1927 fhEtaIsoEtaDecay ->Fill(pt,eta);
1929 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay))
1931 fhPtIsoOtherDecay ->Fill(pt);
1932 fhPhiIsoOtherDecay ->Fill(pt,phi);
1933 fhEtaIsoOtherDecay ->Fill(pt,eta);
1935 // else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion))
1937 // fhPtIsoConversion ->Fill(pt);
1938 // fhPhiIsoConversion ->Fill(pt,phi);
1939 // fhEtaIsoConversion ->Fill(pt,eta);
1941 else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron))// anything else but electrons
1943 fhPtIsoHadron ->Fill(pt);
1944 fhPhiIsoHadron ->Fill(pt,phi);
1945 fhEtaIsoHadron ->Fill(pt,eta);
1947 }//Histograms with MC
1949 }//Isolated histograms
1950 else // NON isolated
1952 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d NOT ISOLATED, fill histograms\n", iaod);
1954 fhPtNoIso ->Fill(pt);
1955 fhEtaPhiNoIso->Fill(eta,phi);
1959 fhPtDecayNoIso ->Fill(pt);
1960 fhEtaPhiDecayNoIso->Fill(eta,phi);
1965 if (GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton)) fhPtNoIsoMCPhoton ->Fill(pt);
1966 if (GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0)) fhPtNoIsoPi0 ->Fill(pt);
1967 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtNoIsoPi0Decay ->Fill(pt);
1968 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtNoIsoEtaDecay ->Fill(pt);
1969 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtNoIsoOtherDecay ->Fill(pt);
1970 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt)) fhPtNoIsoPrompt ->Fill(pt);
1971 else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation)) fhPtNoIsoFragmentation->Fill(pt);
1972 // else if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCConversion)) fhPtNoIsoConversion ->Fill(pt);
1973 else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron)) fhPtNoIsoHadron ->Fill(pt);
1981 //_____________________________________________________________________________________
1982 void AliAnaParticleIsolation::MakeSeveralICAnalysis(AliAODPWG4ParticleCorrelation* ph)
1985 //Isolation Cut Analysis for both methods and different pt cuts and cones
1986 Float_t ptC = ph->Pt();
1987 Float_t etaC = ph->Eta();
1988 Float_t phiC = ph->Phi();
1989 Int_t tag = ph->GetTag();
1990 Bool_t decay = ph->IsTagged();
1992 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeSeveralICAnalysis() - Isolate pT %2.2f\n",ptC);
1994 //Keep original setting used when filling AODs, reset at end of analysis
1995 Float_t ptthresorg = GetIsolationCut()->GetPtThreshold();
1996 Float_t ptfracorg = GetIsolationCut()->GetPtFraction();
1997 Float_t rorg = GetIsolationCut()->GetConeSize();
1999 Float_t coneptsum = 0 ;
2000 Int_t n [10][10];//[fNCones][fNPtThresFrac];
2001 Int_t nfrac[10][10];//[fNCones][fNPtThresFrac];
2002 Bool_t isolated = kFALSE;
2004 Int_t nFracCone = 0;
2006 // fill hist with all particles before isolation criteria
2007 fhPtNoIso ->Fill(ptC);
2008 fhEtaPhiNoIso->Fill(etaC,phiC);
2012 if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)) fhPtNoIsoMCPhoton ->Fill(ptC);
2013 if (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtNoIsoPi0 ->Fill(ptC);
2014 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtNoIsoPi0Decay ->Fill(ptC);
2015 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtNoIsoEtaDecay ->Fill(ptC);
2016 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtNoIsoOtherDecay ->Fill(ptC);
2017 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtNoIsoPrompt ->Fill(ptC);
2018 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtNoIsoFragmentation->Fill(ptC);
2019 // else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtNoIsoConversion ->Fill(ptC);
2020 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtNoIsoHadron ->Fill(ptC);
2025 fhPtDecayNoIso ->Fill(ptC);
2026 fhEtaPhiDecayNoIso->Fill(etaC,phiC);
2028 //Get vertex for photon momentum calculation
2029 Double_t vertex[] = {0,0,0} ; //vertex ;
2030 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
2032 GetReader()->GetVertex(vertex);
2035 //Loop on cone sizes
2036 for(Int_t icone = 0; icone<fNCones; icone++)
2038 //Recover reference arrays with clusters and tracks
2039 TObjArray * refclusters = ph->GetObjArray(GetAODObjArrayName()+"Clusters");
2040 TObjArray * reftracks = ph->GetObjArray(GetAODObjArrayName()+"Tracks");
2042 //If too small or too large pt, skip
2043 if(ptC < GetMinPt() || ptC > GetMaxPt() ) continue ;
2045 //In case a more strict IC is needed in the produced AOD
2047 nCone=0; nFracCone = 0; isolated = kFALSE; coneptsum = 0;
2049 GetIsolationCut()->SetSumPtThreshold(100);
2050 GetIsolationCut()->SetPtThreshold(100);
2051 GetIsolationCut()->SetPtFraction(100);
2052 GetIsolationCut()->SetConeSize(fConeSizes[icone]);
2053 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
2054 GetReader(), GetCaloPID(),
2056 nCone,nFracCone,coneptsum, isolated);
2059 fhSumPtLeadingPt[icone]->Fill(ptC,coneptsum);
2061 // retreive pt tracks to fill histo vs. pt leading
2062 //Fill pt distribution of particles in cone
2063 //fhPtLeadingPt(),fhFRSumPtLeadingPt(),fhFRPtLeadingPt(),
2067 Double_t sumptFR = 0. ;
2068 TObjArray * trackList = GetCTSTracks() ;
2069 for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++)
2071 AliVTrack* track = (AliVTrack *) trackList->At(itrack);
2072 //fill the histograms at forward range
2075 printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Track not available?");
2079 Double_t dPhi = phiC - track->Phi() + TMath::PiOver2();
2080 Double_t dEta = etaC - track->Eta();
2081 Double_t arg = dPhi*dPhi + dEta*dEta;
2082 if(TMath::Sqrt(arg) < fConeSizes[icone])
2084 fhFRPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
2085 sumptFR+=track->Pt();
2088 dPhi = phiC - track->Phi() - TMath::PiOver2();
2089 arg = dPhi*dPhi + dEta*dEta;
2090 if(TMath::Sqrt(arg) < fConeSizes[icone])
2092 fhFRPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
2093 sumptFR+=track->Pt();
2096 fhFRSumPtLeadingPt[icone]->Fill(ptC,sumptFR);
2099 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++)
2101 AliVTrack* track = (AliVTrack *) reftracks->At(itrack);
2102 fhPtLeadingPt[icone]->Fill(ptC,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
2103 coneptsum+=track->Pt();
2109 TLorentzVector mom ;
2110 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++)
2112 AliVCluster* calo = (AliVCluster *) refclusters->At(icalo);
2113 calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
2115 fhPtLeadingPt[icone]->Fill(ptC, mom.Pt());
2116 coneptsum+=mom.Pt();
2122 //Loop on ptthresholds
2123 for(Int_t ipt = 0; ipt<fNPtThresFrac ;ipt++)
2126 nfrac[icone][ipt]=0;
2127 GetIsolationCut()->SetPtThreshold(fPtThresholds[ipt]);
2128 GetIsolationCut()->SetPtFraction(fPtFractions[ipt]) ;
2129 GetIsolationCut()->SetSumPtThreshold(fSumPtThresholds[ipt]);
2131 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters,
2132 GetReader(), GetCaloPID(),
2134 n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);
2136 if(!isolated) continue;
2137 //Normal ptThreshold cut
2139 if(GetDebug() > 0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - cone size %1.1f, ptThres %1.1f, sumptThresh %1.1f, n %d, nfrac %d, coneptsum %2.2f, isolated %d\n",
2140 fConeSizes[icone],fPtThresholds[ipt],fSumPtThresholds[ipt],n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);
2141 if(GetDebug() > 0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - pt %1.1f, eta %1.1f, phi %1.1f\n",ptC, etaC, phiC);
2143 if(n[icone][ipt] == 0)
2145 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling pt threshold loop\n");
2146 fhPtThresIsolated[icone][ipt]->Fill(ptC);
2147 fhEtaPhiPtThresIso[icone][ipt]->Fill(etaC,phiC);
2151 fhPtPtThresDecayIso[icone][ipt]->Fill(ptC);
2152 // fhEtaPhiPtThresDecayIso[icone][ipt]->Fill(etaC,phiC);
2157 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtThresIsolatedPrompt[icone][ipt] ->Fill(ptC) ;
2158 // else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtThresIsolatedConversion[icone][ipt] ->Fill(ptC) ;
2159 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtThresIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
2160 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtThresIsolatedPi0[icone][ipt] ->Fill(ptC) ;
2161 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtThresIsolatedPi0Decay[icone][ipt] ->Fill(ptC) ;
2162 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtThresIsolatedEtaDecay[icone][ipt] ->Fill(ptC) ;
2163 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtThresIsolatedOtherDecay[icone][ipt] ->Fill(ptC) ;
2164 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtThresIsolatedHadron[icone][ipt] ->Fill(ptC) ;
2168 // pt in cone fraction
2169 if(nfrac[icone][ipt] == 0)
2171 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling frac loop\n");
2172 fhPtFracIsolated[icone][ipt]->Fill(ptC);
2173 fhEtaPhiPtFracIso[icone][ipt]->Fill(etaC,phiC);
2177 fhPtPtFracDecayIso[icone][ipt]->Fill(ptC);
2178 fhEtaPhiPtFracDecayIso[icone][ipt]->Fill(etaC,phiC);
2183 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtFracIsolatedPrompt[icone][ipt] ->Fill(ptC) ;
2184 // else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtFracIsolatedConversion[icone][ipt] ->Fill(ptC) ;
2185 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtFracIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
2186 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtFracIsolatedPi0[icone][ipt] ->Fill(ptC) ;
2187 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtFracIsolatedPi0Decay[icone][ipt] ->Fill(ptC) ;
2188 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtFracIsolatedEtaDecay[icone][ipt] ->Fill(ptC) ;
2189 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtFracIsolatedOtherDecay[icone][ipt] ->Fill(ptC) ;
2190 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtFracIsolatedHadron[icone][ipt]->Fill(ptC) ;
2194 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - checking IC method : %i\n",GetIsolationCut()->GetICMethod());
2196 //Pt threshold on pt cand/ sum in cone histograms
2197 if(coneptsum<fSumPtThresholds[ipt])
2198 {// if((GetIsolationCut()->GetICMethod())==1){//kSumPtIC){
2199 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling sum loop\n");
2200 fhPtSumIsolated[icone][ipt]->Fill(ptC) ;
2201 fhEtaPhiPtSumIso[icone][ipt]->Fill(etaC, phiC) ;
2204 fhPtPtSumDecayIso[icone][ipt]->Fill(ptC);
2205 fhEtaPhiPtSumDecayIso[icone][ipt]->Fill(etaC, phiC) ;
2209 // pt sum pt frac method
2210 // if( ((fPtFractions[ipt]*ptC < fSumPtThresholds[ipt]) && (coneptsum < fSumPtThresholds[ipt])) || ((fPtFractions[ipt]*ptC > fSumPtThresholds[ipt]) && (coneptsum < fPtFractions[ipt]*ptC)) )
2212 if(coneptsum < fPtFractions[ipt]*ptC)
2214 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling PtFrac PtSum loop\n");
2215 fhPtFracPtSumIso[icone][ipt]->Fill(ptC) ;
2216 fhEtaPhiFracPtSumIso[icone][ipt]->Fill(etaC,phiC) ;
2220 fhPtFracPtSumDecayIso[icone][ipt]->Fill(ptC);
2221 fhEtaPhiFracPtSumDecayIso[icone][ipt]->Fill(etaC,phiC);
2226 Float_t cellDensity = GetIsolationCut()->GetCellDensity( ph, GetReader());
2227 if(coneptsum<fSumPtThresholds[ipt]*cellDensity)
2228 {//(GetIsolationCut()->GetICMethod())==4){//kSumDensityIC) {
2229 if(GetDebug()>0) printf(" AliAnaParticleIsolation::MakeSeveralICAnalysis() - filling density loop\n");
2230 fhPtSumDensityIso[icone][ipt]->Fill(ptC) ;
2231 fhEtaPhiSumDensityIso[icone][ipt]->Fill(etaC,phiC) ;
2235 fhPtSumDensityDecayIso[icone][ipt]->Fill(ptC);
2236 fhEtaPhiSumDensityDecayIso[icone][ipt]->Fill(etaC, phiC);
2244 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtSumIsolatedPrompt[icone] ->Fill(ptC,coneptsum) ;
2245 // else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtSumIsolatedConversion[icone] ->Fill(ptC,coneptsum) ;
2246 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtSumIsolatedFragmentation[icone]->Fill(ptC,coneptsum) ;
2247 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) fhPtSumIsolatedPi0[icone] ->Fill(ptC,coneptsum) ;
2248 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtSumIsolatedPi0Decay[icone] ->Fill(ptC,coneptsum) ;
2249 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtSumIsolatedEtaDecay[icone] ->Fill(ptC,coneptsum) ;
2250 else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) fhPtSumIsolatedOtherDecay[icone] ->Fill(ptC,coneptsum) ;
2251 else if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) fhPtSumIsolatedHadron[icone]->Fill(ptC,coneptsum) ;
2256 //Reset original parameters for AOD analysis
2257 GetIsolationCut()->SetPtThreshold(ptthresorg);
2258 GetIsolationCut()->SetPtFraction(ptfracorg);
2259 GetIsolationCut()->SetConeSize(rorg);
2263 //_____________________________________________________________
2264 void AliAnaParticleIsolation::Print(const Option_t * opt) const
2267 //Print some relevant parameters set for the analysis
2271 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
2272 AliAnaCaloTrackCorrBaseClass::Print(" ");
2274 printf("ReMake Isolation = %d \n", fReMakeIC) ;
2275 printf("Make Several Isolation = %d \n", fMakeSeveralIC) ;
2276 printf("Calorimeter for isolation = %s \n", fCalorimeter.Data()) ;
2280 printf("N Cone Sizes = %d\n", fNCones) ;
2281 printf("Cone Sizes = \n") ;
2282 for(Int_t i = 0; i < fNCones; i++)
2283 printf(" %1.2f;", fConeSizes[i]) ;
2286 printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ;
2287 printf(" pT thresholds = \n") ;
2288 for(Int_t i = 0; i < fNPtThresFrac; i++)
2289 printf(" %2.2f;", fPtThresholds[i]) ;
2293 printf(" pT fractions = \n") ;
2294 for(Int_t i = 0; i < fNPtThresFrac; i++)
2295 printf(" %2.2f;", fPtFractions[i]) ;
2299 printf("sum pT thresholds = \n") ;
2300 for(Int_t i = 0; i < fNPtThresFrac; i++)
2301 printf(" %2.2f;", fSumPtThresholds[i]) ;
2306 printf("Histograms: %3.1f < pT sum < %3.1f, Nbin = %d\n", fHistoPtSumMin, fHistoPtSumMax, fHistoNPtSumBins );
2307 printf("Histograms: %3.1f < pT in cone < %3.1f, Nbin = %d\n", fHistoPtInConeMin, fHistoPtInConeMax, fHistoNPtInConeBins);