PWG4/PartCorrDep/AliAnaParticleIsolation.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4  * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes hereby granted      *\r
   9  * without fee, provided that the above copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 /* $Id: AliAnaParticleIsolation.cxx 28688 2008-09-11 15:04:07Z gconesab $ */\r
  16 \r
  17 //_________________________________________________________________________\r
  18 // Class for analysis of particle isolation\r
  19 // Input is selected particles put in AOD branch (AliAODPWG4ParticleCorrelation)\r
  20 //\r
  21 //  Class created from old AliPHOSGammaJet \r
  22 //  (see AliRoot versions previous Release 4-09)\r
  23 //\r
  24 // -- Author: Gustavo Conesa (LNF-INFN) \r
  25 //////////////////////////////////////////////////////////////////////////////\r
  26   \r
  27   \r
  28 // --- ROOT system --- \r
  29 #include <TList.h>\r
  30 \r
  31 // --- Analysis system --- \r
  32 #include "AliAnaParticleIsolation.h" \r
  33 #include "AliLog.h"\r
  34 #include "AliCaloTrackReader.h"\r
  35 #include "AliIsolationCut.h"\r
  36 #include "AliNeutralMesonSelection.h"\r
  37 #include "AliAODPWG4ParticleCorrelation.h"\r
  38 #include "AliCaloPID.h"\r
  39 \r
  40 ClassImp(AliAnaParticleIsolation)\r
  41   \r
  42 //____________________________________________________________________________\r
  43   AliAnaParticleIsolation::AliAnaParticleIsolation() : \r
  44     AliAnaPartCorrBaseClass(), fCalorimeter(""), fVertex(),\r
  45         fReMakeIC(0), fMakeSeveralIC(0), fMakeInvMass(0),\r
  46     fhPtIso(0),fhPhiIso(0),fhEtaIso(0), fhConeSumPt(0), fhPtInCone(0),\r
  47     //Several IC\r
  48     fNCones(0),fNPtThresFrac(0), fConeSizes(),  fPtThresholds(),  fPtFractions(), \r
  49     //MC\r
  50         fhPtIsoPrompt(0),fhPhiIsoPrompt(0),fhEtaIsoPrompt(0), \r
  51     fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(),  fhPtSumIsolatedPrompt(),\r
  52         fhPtIsoFragmentation(0),fhPhiIsoFragmentation(0),fhEtaIsoFragmentation(0), \r
  53     fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(),  fhPtSumIsolatedFragmentation(),\r
  54         fhPtIsoPi0Decay(0),fhPhiIsoPi0Decay(0),fhEtaIsoPi0Decay(0),\r
  55     fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(),  fhPtSumIsolatedPi0Decay(),\r
  56         fhPtIsoOtherDecay(0),fhPhiIsoOtherDecay(0),fhEtaIsoOtherDecay(0), \r
  57     fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(),  fhPtSumIsolatedOtherDecay(),\r
  58         fhPtIsoConversion(0),fhPhiIsoConversion(0),fhEtaIsoConversion(0), \r
  59     fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(),  fhPtSumIsolatedConversion(),\r
  60         fhPtIsoUnknown(0),fhPhiIsoUnknown(0),fhEtaIsoUnknown(0), \r
  61     fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(),  fhPtSumIsolatedUnknown(),\r
  62         //Histograms settings\r
  63         fHistoNPtSumBins(0),    fHistoPtSumMax(0.),    fHistoPtSumMin(0.),\r
  64         fHistoNPtInConeBins(0), fHistoPtInConeMax(0.), fHistoPtInConeMin(0.)\r
  65 {\r
  66   //default ctor\r
  67   \r
  68   //Initialize parameters\r
  69   InitParameters();\r
  70   \r
  71   fVertex[0] = 0.;  fVertex[1] = 0.;  fVertex[2] = 0.;  \r
  72           \r
  73   for(Int_t i = 0; i < 5 ; i++){ \r
  74     fConeSizes[i] = 0 ; \r
  75     fhPtSumIsolated[i] = 0 ;  \r
  76     \r
  77     fhPtSumIsolatedPrompt[i] = 0 ;  \r
  78     fhPtSumIsolatedFragmentation[i] = 0 ;  \r
  79     fhPtSumIsolatedPi0Decay[i] = 0 ;  \r
  80     fhPtSumIsolatedOtherDecay[i] = 0 ;  \r
  81     fhPtSumIsolatedConversion[i] = 0 ;  \r
  82     fhPtSumIsolatedUnknown[i] = 0 ;  \r
  83     \r
  84     for(Int_t j = 0; j < 5 ; j++){ \r
  85       fhPtThresIsolated[i][j] = 0 ;  \r
  86       fhPtFracIsolated[i][j] = 0 ; \r
  87       \r
  88       fhPtThresIsolatedPrompt[i][j] = 0 ;  \r
  89       fhPtThresIsolatedFragmentation[i][j] = 0 ; \r
  90       fhPtThresIsolatedPi0Decay[i][j] = 0 ;  \r
  91       fhPtThresIsolatedOtherDecay[i][j] = 0 ;  \r
  92       fhPtThresIsolatedConversion[i][j] = 0 ;  \r
  93       fhPtThresIsolatedUnknown[i][j] = 0 ;  \r
  94   \r
  95       fhPtFracIsolatedPrompt[i][j] = 0 ;  \r
  96       fhPtFracIsolatedFragmentation[i][j] = 0 ;  \r
  97       fhPtFracIsolatedPi0Decay[i][j] = 0 ;  \r
  98       fhPtFracIsolatedOtherDecay[i][j] = 0 ;  \r
  99       fhPtFracIsolatedConversion[i][j] = 0 ;\r
 100       fhPtFracIsolatedUnknown[i][j] = 0 ;  \r
 101  \r
 102     }  \r
 103   } \r
 104   \r
 105   for(Int_t i = 0; i < 5 ; i++){ \r
 106     fPtFractions[i]=  0 ; \r
 107     fPtThresholds[i]= 0 ; \r
 108   } \r
 109 \r
 110 \r
 111 }\r
 112 \r
 113 //____________________________________________________________________________\r
 114 AliAnaParticleIsolation::AliAnaParticleIsolation(const AliAnaParticleIsolation & g) : \r
 115   AliAnaPartCorrBaseClass(g), fCalorimeter(g.fCalorimeter),fVertex(),\r
 116   fReMakeIC(g.fReMakeIC), fMakeSeveralIC(g.fMakeSeveralIC),  fMakeInvMass(g.fMakeInvMass),\r
 117   fhPtIso(g.fhPtIso),fhPhiIso(g.fhPhiIso),fhEtaIso(g.fhEtaIso), \r
 118   fhConeSumPt(g.fhConeSumPt), fhPtInCone(g.fhPtInCone),\r
 119   //Several IC\r
 120   fNCones(g.fNCones),fNPtThresFrac(g.fNPtThresFrac), fConeSizes(), fPtThresholds(),  fPtFractions(), \r
 121   fhPtThresIsolated(),  fhPtFracIsolated(), fhPtSumIsolated(),\r
 122   //MC\r
 123   fhPtIsoPrompt(g.fhPtIsoPrompt),fhPhiIsoPrompt(g.fhPhiIsoPrompt),fhEtaIsoPrompt(g.fhEtaIsoPrompt), \r
 124   fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(),  fhPtSumIsolatedPrompt(),\r
 125   fhPtIsoFragmentation(g.fhPtIsoFragmentation),fhPhiIsoFragmentation(g.fhPhiIsoFragmentation),fhEtaIsoFragmentation(g.fhEtaIsoFragmentation),\r
 126   fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(),  fhPtSumIsolatedFragmentation(),\r
 127   fhPtIsoPi0Decay(g.fhPtIsoPi0Decay),fhPhiIsoPi0Decay(g.fhPhiIsoPi0Decay),fhEtaIsoPi0Decay(g.fhEtaIsoPi0Decay), \r
 128   fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(),  fhPtSumIsolatedPi0Decay(),\r
 129   fhPtIsoOtherDecay(g.fhPtIsoOtherDecay),fhPhiIsoOtherDecay(g.fhPhiIsoOtherDecay),fhEtaIsoOtherDecay(g.fhEtaIsoOtherDecay), \r
 130   fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(),  fhPtSumIsolatedOtherDecay(),\r
 131   fhPtIsoConversion(g. fhPtIsoConversion),fhPhiIsoConversion(g.fhPhiIsoConversion),fhEtaIsoConversion(g.fhEtaIsoConversion), \r
 132   fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(),  fhPtSumIsolatedConversion(),\r
 133   fhPtIsoUnknown(g.fhPtIsoUnknown),fhPhiIsoUnknown(g.fhPhiIsoUnknown),fhEtaIsoUnknown(g.fhEtaIsoUnknown), \r
 134   fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(),  fhPtSumIsolatedUnknown(),\r
 135   //Histograms\r
 136   fHistoNPtSumBins(g.fHistoNPtSumBins), fHistoPtSumMax(g.fHistoPtSumMax), fHistoPtSumMin(g.fHistoPtSumMax),\r
 137   fHistoNPtInConeBins(g.fHistoNPtInConeBins), fHistoPtInConeMax(g.fHistoPtInConeMax), fHistoPtInConeMin(g.fHistoPtInConeMin)\r
 138 {  \r
 139   // cpy ctor\r
 140   \r
 141         fVertex[0] = g.fVertex[0];  \r
 142         fVertex[1] = g.fVertex[1];  \r
 143         fVertex[2] = g.fVertex[2];  \r
 144 \r
 145         //Several IC\r
 146         for(Int_t i = 0; i < fNCones ; i++){\r
 147                 fConeSizes[i] =  g.fConeSizes[i];\r
 148                 fhPtSumIsolated[i] = g.fhPtSumIsolated[i]; \r
 149                 \r
 150                 fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i]; \r
 151                 fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i]; \r
 152                 fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i]; \r
 153                 fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i]; \r
 154                 fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i]; \r
 155                 fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i]; \r
 156                 \r
 157                 for(Int_t j = 0; j < fNPtThresFrac ; j++){\r
 158                         fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j]; \r
 159                         fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j];\r
 160                         \r
 161                         fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j]; \r
 162                         fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j]; \r
 163                         fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j]; \r
 164                         fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j]; \r
 165                         fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j]; \r
 166                         fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j]; \r
 167                         \r
 168                         fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j]; \r
 169                         fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j]; \r
 170                         fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j]; \r
 171                         fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j]; \r
 172                         fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j]; \r
 173                         fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j]; \r
 174                         \r
 175                 } \r
 176         }\r
 177         \r
 178         for(Int_t i = 0; i < fNPtThresFrac ; i++){\r
 179                 fPtFractions[i]=   g.fPtFractions[i];\r
 180                 fPtThresholds[i]=   g.fPtThresholds[i];\r
 181         }\r
 182         \r
 183 }\r
 184 \r
 185 //_________________________________________________________________________\r
 186 AliAnaParticleIsolation & AliAnaParticleIsolation::operator = (const AliAnaParticleIsolation & g)\r
 187 {\r
 188         // assignment operator\r
 189   \r
 190         if(&g == this) return *this;\r
 191 \r
 192         fReMakeIC = g.fReMakeIC ;\r
 193         fMakeSeveralIC = g.fMakeSeveralIC ;\r
 194         fMakeInvMass = g.fMakeInvMass ;\r
 195         fCalorimeter = g.fCalorimeter ;\r
 196         fVertex[0] = g.fVertex[0];  \r
 197         fVertex[1] = g.fVertex[1];  \r
 198         fVertex[2] = g.fVertex[2]; \r
 199          \r
 200         fhConeSumPt = g.fhConeSumPt ;\r
 201         fhPtInCone = g.fhPtInCone;\r
 202         \r
 203         fhPtIso = g.fhPtIso ; \r
 204         fhPhiIso = g.fhPhiIso ;\r
 205         fhEtaIso = g.fhEtaIso ;\r
 206         \r
 207         fhPtIsoPrompt = g.fhPtIsoPrompt;\r
 208         fhPhiIsoPrompt = g.fhPhiIsoPrompt;\r
 209         fhEtaIsoPrompt = g.fhEtaIsoPrompt; \r
 210         fhPtIsoFragmentation = g.fhPtIsoFragmentation;\r
 211         fhPhiIsoFragmentation = g.fhPhiIsoFragmentation;\r
 212         fhEtaIsoFragmentation = g.fhEtaIsoFragmentation; \r
 213         fhPtIsoPi0Decay = g.fhPtIsoPi0Decay;\r
 214         fhPhiIsoPi0Decay = g.fhPhiIsoPi0Decay;\r
 215         fhEtaIsoPi0Decay = g.fhEtaIsoPi0Decay; \r
 216         fhPtIsoOtherDecay = g.fhPtIsoOtherDecay;\r
 217         fhPhiIsoOtherDecay = g.fhPhiIsoOtherDecay;\r
 218         fhEtaIsoOtherDecay = g.fhEtaIsoOtherDecay; \r
 219         fhPtIsoConversion = g. fhPtIsoConversion;\r
 220         fhPhiIsoConversion = g.fhPhiIsoConversion;\r
 221         fhEtaIsoConversion = g.fhEtaIsoConversion; \r
 222         fhPtIsoUnknown = g.fhPtIsoUnknown;\r
 223         fhPhiIsoUnknown = g.fhPhiIsoUnknown;\r
 224         fhEtaIsoUnknown = g.fhEtaIsoUnknown; \r
 225         \r
 226         //Several IC\r
 227         fNCones = g.fNCones ;\r
 228         fNPtThresFrac = g.fNPtThresFrac ; \r
 229         \r
 230         for(Int_t i = 0; i < fNCones ; i++){\r
 231                 fConeSizes[i] =  g.fConeSizes[i];\r
 232                 fhPtSumIsolated[i] = g.fhPtSumIsolated[i] ;\r
 233                 \r
 234                 fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i]; \r
 235                 fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i]; \r
 236                 fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i]; \r
 237                 fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i]; \r
 238                 fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i]; \r
 239                 fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i]; \r
 240                 \r
 241                 for(Int_t j = 0; j < fNPtThresFrac ; j++){\r
 242                         fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j] ;\r
 243                         fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j] ;\r
 244                         \r
 245                         fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j]; \r
 246                         fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j]; \r
 247                         fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j]; \r
 248                         fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j]; \r
 249                         fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j]; \r
 250                         fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j]; \r
 251                         \r
 252                         fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j]; \r
 253                         fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j]; \r
 254                         fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j]; \r
 255                         fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j]; \r
 256                         fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j]; \r
 257                         fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j]; \r
 258                         \r
 259                 }\r
 260         }\r
 261         \r
 262         for(Int_t i = 0; i < fNPtThresFrac ; i++){\r
 263                 fPtThresholds[i]=   g.fPtThresholds[i];\r
 264                 fPtFractions[i]=   g.fPtFractions[i];\r
 265         }\r
 266         \r
 267         \r
 268         fHistoNPtSumBins    = g.fHistoNPtSumBins;\r
 269         fHistoPtSumMax      = g.fHistoPtSumMax; \r
 270         fHistoPtSumMin      = g.fHistoPtSumMax;\r
 271         fHistoNPtInConeBins = g.fHistoNPtInConeBins;\r
 272         fHistoPtInConeMax   = g.fHistoPtInConeMax;\r
 273         fHistoPtInConeMin   = g.fHistoPtInConeMin;      \r
 274         \r
 275         return *this;\r
 276         \r
 277 }\r
 278 \r
 279 //____________________________________________________________________________\r
 280 AliAnaParticleIsolation::~AliAnaParticleIsolation() \r
 281 {\r
 282   //dtor\r
 283   //do not delete histograms\r
 284 \r
 285  delete [] fConeSizes ; \r
 286  delete [] fPtThresholds ; \r
 287  delete [] fPtFractions ; \r
 288  delete [] fVertex ;\r
 289 }\r
 290 \r
 291 //_________________________________________________________________________\r
 292 Bool_t AliAnaParticleIsolation::CheckInvMass(const Int_t iaod, const AliAODPWG4ParticleCorrelation * part1) const\r
 293 {\r
 294         // Search if there is a companion decay photon to the candidate \r
 295         // and discard it in such case\r
 296         // Use it only if isolation candidates are photons\r
 297         // Make sure that no selection on photon pt is done in the input aod photon list.\r
 298         TLorentzVector mom1 = *(part1->Momentum());\r
 299         TLorentzVector mom2 ;\r
 300         for(Int_t jaod = 0; jaod < GetInputAODBranch()->GetEntriesFast(); jaod++){\r
 301                 if(iaod == jaod) continue ;\r
 302                 AliAODPWG4ParticleCorrelation * part2 =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(jaod));\r
 303                 mom2 = *(part2->Momentum());\r
 304                 //Select good pair (good phi, pt cuts, aperture and invariant mass)\r
 305                 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2)){\r
 306                         if(GetDebug() > 1)printf("Selected gamma pair: pt %f, phi %f, eta%f",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());\r
 307                         return kTRUE ;\r
 308                 }\r
 309         }//loop\r
 310         \r
 311   return kFALSE;\r
 312 }\r
 313 \r
 314 //________________________________________________________________________\r
 315 TList *  AliAnaParticleIsolation::GetCreateOutputObjects()\r
 316 {  \r
 317         // Create histograms to be saved in output file and \r
 318         // store them in outputContainer\r
 319         TList * outputContainer = new TList() ; \r
 320         outputContainer->SetName("IsolatedParticleHistos") ; \r
 321         \r
 322         Int_t nptbins  = GetHistoNPtBins();\r
 323         Int_t nphibins = GetHistoNPhiBins();\r
 324         Int_t netabins = GetHistoNEtaBins();\r
 325         Float_t ptmax  = GetHistoPtMax();\r
 326         Float_t phimax = GetHistoPhiMax();\r
 327         Float_t etamax = GetHistoEtaMax();\r
 328         Float_t ptmin  = GetHistoPtMin();\r
 329         Float_t phimin = GetHistoPhiMin();\r
 330         Float_t etamin = GetHistoEtaMin();      \r
 331         \r
 332         Int_t nptsumbins    = fHistoNPtSumBins;\r
 333         Float_t ptsummax    = fHistoPtSumMax;\r
 334         Float_t ptsummin    = fHistoPtSumMin;   \r
 335         Int_t nptinconebins = fHistoNPtInConeBins;\r
 336         Float_t ptinconemax = fHistoPtInConeMax;\r
 337         Float_t ptinconemin = fHistoPtInConeMin;\r
 338         \r
 339         if(!fMakeSeveralIC){\r
 340                 \r
 341                 fhConeSumPt  = new TH2F\r
 342                 ("hConePtSum","#Sigma p_{T}  in cone ",nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 343                 fhConeSumPt->SetYTitle("#Sigma p_{T}");\r
 344                 fhConeSumPt->SetXTitle("p_{T #gamma} (GeV/c)");\r
 345                 outputContainer->Add(fhConeSumPt) ;\r
 346                 \r
 347                 fhPtInCone  = new TH2F\r
 348                 ("hPtInCone","p_{T} in cone ",nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);\r
 349                 fhPtInCone->SetYTitle("p_{T in cone} (GeV/c)");\r
 350                 fhPtInCone->SetXTitle("p_{T #gamma} (GeV/c)");\r
 351                 outputContainer->Add(fhPtInCone) ;\r
 352                 \r
 353                 fhPtIso  = new TH1F("hPtIso","Isolated Number of #gamma over calorimeter",nptbins,ptmin,ptmax); \r
 354                 fhPtIso->SetYTitle("N");\r
 355                 fhPtIso->SetXTitle("p_{T #gamma}(GeV/c)");\r
 356                 outputContainer->Add(fhPtIso) ; \r
 357                 \r
 358                 fhPhiIso  = new TH2F\r
 359                 ("hPhiIso","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
 360                 fhPhiIso->SetYTitle("#phi");\r
 361                 fhPhiIso->SetXTitle("p_{T #gamma} (GeV/c)");\r
 362                 outputContainer->Add(fhPhiIso) ; \r
 363                 \r
 364                 fhEtaIso  = new TH2F\r
 365                 ("hEtaIso","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
 366                 fhEtaIso->SetYTitle("#eta");\r
 367                 fhEtaIso->SetXTitle("p_{T #gamma} (GeV/c)");\r
 368                 outputContainer->Add(fhEtaIso) ;\r
 369                 \r
 370                 if(IsDataMC()){\r
 371                         \r
 372                         fhPtIsoPrompt  = new TH1F("hPtIsoPrompt","Isolated Number of #gamma over calorimeter",nptbins,ptmin,ptmax); \r
 373                         fhPtIsoPrompt->SetYTitle("N");\r
 374                         fhPtIsoPrompt->SetXTitle("p_{T #gamma}(GeV/c)");\r
 375                         outputContainer->Add(fhPtIsoPrompt) ; \r
 376                         \r
 377                         fhPhiIsoPrompt  = new TH2F\r
 378                         ("hPhiIsoPrompt","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
 379                         fhPhiIsoPrompt->SetYTitle("#phi");\r
 380                         fhPhiIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");\r
 381                         outputContainer->Add(fhPhiIsoPrompt) ; \r
 382                         \r
 383                         fhEtaIsoPrompt  = new TH2F\r
 384                         ("hEtaIsoPrompt","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
 385                         fhEtaIsoPrompt->SetYTitle("#eta");\r
 386                         fhEtaIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");\r
 387                         outputContainer->Add(fhEtaIsoPrompt) ;\r
 388                         \r
 389                         fhPtIsoFragmentation  = new TH1F("hPtIsoFragmentation","Isolated Number of #gamma over calorimeter",nptbins,ptmin,ptmax); \r
 390                         fhPtIsoFragmentation->SetYTitle("N");\r
 391                         fhPtIsoFragmentation->SetXTitle("p_{T #gamma}(GeV/c)");\r
 392                         outputContainer->Add(fhPtIsoFragmentation) ; \r
 393                         \r
 394                         fhPhiIsoFragmentation  = new TH2F\r
 395                         ("hPhiIsoFragmentation","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
 396                         fhPhiIsoFragmentation->SetYTitle("#phi");\r
 397                         fhPhiIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");\r
 398                         outputContainer->Add(fhPhiIsoFragmentation) ; \r
 399                         \r
 400                         fhEtaIsoFragmentation  = new TH2F\r
 401                         ("hEtaIsoFragmentation","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
 402                         fhEtaIsoFragmentation->SetYTitle("#eta");\r
 403                         fhEtaIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");\r
 404                         outputContainer->Add(fhEtaIsoFragmentation) ;\r
 405                         \r
 406                         fhPtIsoPi0Decay  = new TH1F("hPtIsoPi0Decay","Isolated Number of #gamma over calorimeter",nptbins,ptmin,ptmax); \r
 407                         fhPtIsoPi0Decay->SetYTitle("N");\r
 408                         fhPtIsoPi0Decay->SetXTitle("p_{T #gamma}(GeV/c)");\r
 409                         outputContainer->Add(fhPtIsoPi0Decay) ; \r
 410                         \r
 411                         fhPhiIsoPi0Decay  = new TH2F\r
 412                         ("hPhiIsoPi0Decay","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
 413                         fhPhiIsoPi0Decay->SetYTitle("#phi");\r
 414                         fhPhiIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");\r
 415                         outputContainer->Add(fhPhiIsoPi0Decay) ; \r
 416                         \r
 417                         fhEtaIsoPi0Decay  = new TH2F\r
 418                         ("hEtaIsoPi0Decay","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
 419                         fhEtaIsoPi0Decay->SetYTitle("#eta");\r
 420                         fhEtaIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");\r
 421                         outputContainer->Add(fhEtaIsoPi0Decay) ;\r
 422                         \r
 423                         fhPtIsoOtherDecay  = new TH1F("hPtIsoOtherDecay","Isolated Number of #gamma over calorimeter",nptbins,ptmin,ptmax); \r
 424                         fhPtIsoOtherDecay->SetYTitle("N");\r
 425                         fhPtIsoOtherDecay->SetXTitle("p_{T #gamma}(GeV/c)");\r
 426                         outputContainer->Add(fhPtIsoOtherDecay) ; \r
 427                         \r
 428                         fhPhiIsoOtherDecay  = new TH2F\r
 429                         ("hPhiIsoOtherDecay","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
 430                         fhPhiIsoOtherDecay->SetYTitle("#phi");\r
 431                         fhPhiIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");\r
 432                         outputContainer->Add(fhPhiIsoOtherDecay) ; \r
 433                         \r
 434                         fhEtaIsoOtherDecay  = new TH2F\r
 435                         ("hEtaIsoOtherDecay","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
 436                         fhEtaIsoOtherDecay->SetYTitle("#eta");\r
 437                         fhEtaIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");\r
 438                         outputContainer->Add(fhEtaIsoOtherDecay) ;\r
 439                         \r
 440                         fhPtIsoConversion  = new TH1F("hPtIsoConversion","Isolated Number of #gamma over calorimeter",nptbins,ptmin,ptmax); \r
 441                         fhPtIsoConversion->SetYTitle("N");\r
 442                         fhPtIsoConversion->SetXTitle("p_{T #gamma}(GeV/c)");\r
 443                         outputContainer->Add(fhPtIsoConversion) ; \r
 444                         \r
 445                         fhPhiIsoConversion  = new TH2F\r
 446                         ("hPhiIsoConversion","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
 447                         fhPhiIsoConversion->SetYTitle("#phi");\r
 448                         fhPhiIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");\r
 449                         outputContainer->Add(fhPhiIsoConversion) ; \r
 450                         \r
 451                         fhEtaIsoConversion  = new TH2F\r
 452                         ("hEtaIsoConversion","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
 453                         fhEtaIsoConversion->SetYTitle("#eta");\r
 454                         fhEtaIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");\r
 455                         outputContainer->Add(fhEtaIsoConversion) ;\r
 456                         \r
 457                         fhPtIsoUnknown  = new TH1F("hPtIsoUnknown","Isolated Number of #gamma over calorimeter",nptbins,ptmin,ptmax); \r
 458                         fhPtIsoUnknown->SetYTitle("N");\r
 459                         fhPtIsoUnknown->SetXTitle("p_{T #gamma}(GeV/c)");\r
 460                         outputContainer->Add(fhPtIsoUnknown) ; \r
 461                         \r
 462                         fhPhiIsoUnknown  = new TH2F\r
 463                         ("hPhiIsoUnknown","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); \r
 464                         fhPhiIsoUnknown->SetYTitle("#phi");\r
 465                         fhPhiIsoUnknown->SetXTitle("p_{T #gamma} (GeV/c)");\r
 466                         outputContainer->Add(fhPhiIsoUnknown) ; \r
 467                         \r
 468                         fhEtaIsoUnknown  = new TH2F\r
 469                         ("hEtaIsoUnknown","Isolated #phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); \r
 470                         fhEtaIsoUnknown->SetYTitle("#eta");\r
 471                         fhEtaIsoUnknown->SetXTitle("p_{T #gamma} (GeV/c)");\r
 472                         outputContainer->Add(fhEtaIsoUnknown) ;\r
 473                 }//Histos with MC\r
 474                 \r
 475         }\r
 476         \r
 477         if(fMakeSeveralIC){\r
 478                 char name[128];\r
 479                 char title[128];\r
 480                 for(Int_t icone = 0; icone<fNCones; icone++){\r
 481                         sprintf(name,"hPtSumIsolated_Cone_%d",icone);\r
 482                         sprintf(title,"Candidate cone sum p_{T} for cone size %d vs candidate p_{T}",icone);\r
 483                         fhPtSumIsolated[icone]  = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 484                         fhPtSumIsolated[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");\r
 485                         fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)");\r
 486                         outputContainer->Add(fhPtSumIsolated[icone]) ; \r
 487                         \r
 488                         if(IsDataMC()){\r
 489                                 sprintf(name,"hPtSumIsolatedPrompt_Cone_%d",icone);\r
 490                                 sprintf(title,"Candidate Prompt cone sum p_{T} for cone size %d vs candidate p_{T}",icone);\r
 491                                 fhPtSumIsolatedPrompt[icone]  = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 492                                 fhPtSumIsolatedPrompt[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");\r
 493                                 fhPtSumIsolatedPrompt[icone]->SetXTitle("p_{T} (GeV/c)");\r
 494                                 outputContainer->Add(fhPtSumIsolatedPrompt[icone]) ; \r
 495                                 \r
 496                                 sprintf(name,"hPtSumIsolatedFragmentation_Cone_%d",icone);\r
 497                                 sprintf(title,"Candidate Fragmentation cone sum p_{T} for cone size %d vs candidate p_{T}",icone);\r
 498                                 fhPtSumIsolatedFragmentation[icone]  = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 499                                 fhPtSumIsolatedFragmentation[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");\r
 500                                 fhPtSumIsolatedFragmentation[icone]->SetXTitle("p_{T} (GeV/c)");\r
 501                                 outputContainer->Add(fhPtSumIsolatedFragmentation[icone]) ; \r
 502                                 \r
 503                                 sprintf(name,"hPtSumIsolatedPi0Decay_Cone_%d",icone);\r
 504                                 sprintf(title,"Candidate Pi0Decay cone sum p_{T} for cone size %d vs candidate p_{T}",icone);\r
 505                                 fhPtSumIsolatedPi0Decay[icone]  = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 506                                 fhPtSumIsolatedPi0Decay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");\r
 507                                 fhPtSumIsolatedPi0Decay[icone]->SetXTitle("p_{T} (GeV/c)");\r
 508                                 outputContainer->Add(fhPtSumIsolatedPi0Decay[icone]) ; \r
 509                                 \r
 510                                 sprintf(name,"hPtSumIsolatedOtherDecay_Cone_%d",icone);\r
 511                                 sprintf(title,"Candidate OtherDecay cone sum p_{T} for cone size %d vs candidate p_{T}",icone);\r
 512                                 fhPtSumIsolatedOtherDecay[icone]  = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 513                                 fhPtSumIsolatedOtherDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");\r
 514                                 fhPtSumIsolatedOtherDecay[icone]->SetXTitle("p_{T} (GeV/c)");\r
 515                                 outputContainer->Add(fhPtSumIsolatedOtherDecay[icone]) ; \r
 516                                 \r
 517                                 sprintf(name,"hPtSumIsolatedConversion_Cone_%d",icone);\r
 518                                 sprintf(title,"Candidate Conversion cone sum p_{T} for cone size %d vs candidate p_{T}",icone);\r
 519                                 fhPtSumIsolatedConversion[icone]  = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 520                                 fhPtSumIsolatedConversion[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");\r
 521                                 fhPtSumIsolatedConversion[icone]->SetXTitle("p_{T} (GeV/c)");\r
 522                                 outputContainer->Add(fhPtSumIsolatedConversion[icone]) ; \r
 523                                 \r
 524                                 sprintf(name,"hPtSumIsolatedUnknown_Cone_%d",icone);\r
 525                                 sprintf(title,"Candidate Unknown cone sum p_{T} for cone size %d vs candidate p_{T}",icone);\r
 526                                 fhPtSumIsolatedUnknown[icone]  = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);\r
 527                                 fhPtSumIsolatedUnknown[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");\r
 528                                 fhPtSumIsolatedUnknown[icone]->SetXTitle("p_{T} (GeV/c)");\r
 529                                 outputContainer->Add(fhPtSumIsolatedUnknown[icone]) ; \r
 530                                 \r
 531                         }//Histos with MC\r
 532                         \r
 533                         for(Int_t ipt = 0; ipt<fNPtThresFrac;ipt++){ \r
 534                                 sprintf(name,"hPtThresIsol_Cone_%d_Pt%d",icone,ipt);\r
 535                                 sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 536                                 fhPtThresIsolated[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 537                                 fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 538                                 outputContainer->Add(fhPtThresIsolated[icone][ipt]) ; \r
 539                                 \r
 540                                 sprintf(name,"hPtFracIsol_Cone_%d_Pt%d",icone,ipt);\r
 541                                 sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 542                                 fhPtFracIsolated[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 543                                 fhPtFracIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 544                                 outputContainer->Add(fhPtFracIsolated[icone][ipt]) ; \r
 545                                 \r
 546                                 if(IsDataMC()){\r
 547                                         sprintf(name,"hPtThresIsolPrompt_Cone_%d_Pt%d",icone,ipt);\r
 548                                         sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 549                                         fhPtThresIsolatedPrompt[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 550                                         fhPtThresIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 551                                         outputContainer->Add(fhPtThresIsolatedPrompt[icone][ipt]) ; \r
 552                                         \r
 553                                         sprintf(name,"hPtFracIsolPrompt_Cone_%d_Pt%d",icone,ipt);\r
 554                                         sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 555                                         fhPtFracIsolatedPrompt[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 556                                         fhPtFracIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 557                                         outputContainer->Add(fhPtFracIsolatedPrompt[icone][ipt]) ; \r
 558                                         \r
 559                                         sprintf(name,"hPtThresIsolFragmentation_Cone_%d_Pt%d",icone,ipt);\r
 560                                         sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 561                                         fhPtThresIsolatedFragmentation[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 562                                         fhPtThresIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 563                                         outputContainer->Add(fhPtThresIsolatedFragmentation[icone][ipt]) ; \r
 564                                         \r
 565                                         sprintf(name,"hPtFracIsolFragmentation_Cone_%d_Pt%d",icone,ipt);\r
 566                                         sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 567                                         fhPtFracIsolatedFragmentation[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 568                                         fhPtFracIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 569                                         outputContainer->Add(fhPtFracIsolatedFragmentation[icone][ipt]) ; \r
 570                                         \r
 571                                         sprintf(name,"hPtThresIsolPi0Decay_Cone_%d_Pt%d",icone,ipt);\r
 572                                         sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 573                                         fhPtThresIsolatedPi0Decay[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 574                                         fhPtThresIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 575                                         outputContainer->Add(fhPtThresIsolatedPi0Decay[icone][ipt]) ; \r
 576                                         \r
 577                                         sprintf(name,"hPtFracIsolPi0Decay_Cone_%d_Pt%d",icone,ipt);\r
 578                                         sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 579                                         fhPtFracIsolatedPi0Decay[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 580                                         fhPtFracIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 581                                         outputContainer->Add(fhPtFracIsolatedPi0Decay[icone][ipt]) ; \r
 582                                         \r
 583                                         sprintf(name,"hPtThresIsolOtherDecay_Cone_%d_Pt%d",icone,ipt);\r
 584                                         sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 585                                         fhPtThresIsolatedOtherDecay[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 586                                         fhPtThresIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 587                                         outputContainer->Add(fhPtThresIsolatedOtherDecay[icone][ipt]) ; \r
 588                                         \r
 589                                         sprintf(name,"hPtFracIsolOtherDecay_Cone_%d_Pt%d",icone,ipt);\r
 590                                         sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 591                                         fhPtFracIsolatedOtherDecay[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 592                                         fhPtFracIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 593                                         outputContainer->Add(fhPtFracIsolatedOtherDecay[icone][ipt]) ;\r
 594                                         \r
 595                                         sprintf(name,"hPtThresIsolConversion_Cone_%d_Pt%d",icone,ipt);\r
 596                                         sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 597                                         fhPtThresIsolatedConversion[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 598                                         fhPtThresIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 599                                         outputContainer->Add(fhPtThresIsolatedConversion[icone][ipt]) ; \r
 600                                         \r
 601                                         sprintf(name,"hPtFracIsolConversion_Cone_%d_Pt%d",icone,ipt);\r
 602                                         sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 603                                         fhPtFracIsolatedConversion[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 604                                         fhPtFracIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 605                                         outputContainer->Add(fhPtFracIsolatedConversion[icone][ipt]) ;\r
 606                                         \r
 607                                         sprintf(name,"hPtThresIsolUnknown_Cone_%d_Pt%d",icone,ipt);\r
 608                                         sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 609                                         fhPtThresIsolatedUnknown[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 610                                         fhPtThresIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 611                                         outputContainer->Add(fhPtThresIsolatedUnknown[icone][ipt]) ; \r
 612                                         \r
 613                                         sprintf(name,"hPtFracIsolUnknown_Cone_%d_Pt%d",icone,ipt);\r
 614                                         sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);\r
 615                                         fhPtFracIsolatedUnknown[icone][ipt]  = new TH1F(name, title,nptbins,ptmin,ptmax);\r
 616                                         fhPtFracIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)");\r
 617                                         outputContainer->Add(fhPtFracIsolatedUnknown[icone][ipt]) ;  \r
 618                                         \r
 619                                 }//Histos with MC\r
 620                                 \r
 621                         }//icone loop\r
 622                 }//ipt loop\r
 623         }\r
 624         \r
 625         //Keep neutral meson selection histograms if requiered\r
 626         //Setting done in AliNeutralMesonSelection\r
 627         if(fMakeInvMass && GetNeutralMesonSelection()){\r
 628                 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;\r
 629                 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())\r
 630                         for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;\r
 631         }\r
 632         \r
 633         \r
 634         //Save parameters used for analysis\r
 635         TString parList ; //this will be list of parameters used for this analysis.\r
 636         char onePar[255] ;\r
 637         \r
 638         sprintf(onePar,"--- AliAnaParticleIsolation ---\n") ;\r
 639         parList+=onePar ;       \r
 640         sprintf(onePar,"Calorimeter: %s\n",fCalorimeter.Data()) ;\r
 641         parList+=onePar ;\r
 642         sprintf(onePar,"fReMakeIC =%d (Flag for reisolation during histogram filling) \n",fReMakeIC) ;\r
 643         parList+=onePar ;\r
 644         sprintf(onePar,"fMakeSeveralIC=%d (Flag for isolation with several cuts at the same time ) \n",fMakeSeveralIC) ;\r
 645         parList+=onePar ;\r
 646         sprintf(onePar,"fMakeInvMass=%d (Flag for rejection of candidates with a pi0 inv mass pair) \n",fMakeInvMass) ;\r
 647         parList+=onePar ;\r
 648     \r
 649         if(fMakeSeveralIC){\r
 650                 sprintf(onePar,"fNCones =%d (Number of cone sizes) \n",fNCones) ;\r
 651                 parList+=onePar ;\r
 652                 sprintf(onePar,"fNPtThresFrac=%d (Flag for isolation with several cuts at the same time ) \n",fNPtThresFrac) ;\r
 653                 parList+=onePar ;\r
 654                 \r
 655                 for(Int_t icone = 0; icone < fNCones ; icone++){\r
 656                         sprintf(onePar,"fConeSizes[%d]=%1.2f (isolation cone size) \n",icone, fConeSizes[icone]) ;\r
 657                         parList+=onePar ;       \r
 658                 }\r
 659                 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++){\r
 660                         sprintf(onePar,"fPtThresholds[%d]=%1.2f (isolation pt threshold) \n",ipt, fPtThresholds[ipt]) ;\r
 661                         parList+=onePar ;       \r
 662                 }\r
 663                 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++){\r
 664                         sprintf(onePar,"fPtFractions[%d]=%1.2f (isolation pt fraction threshold) \n",ipt, fPtFractions[ipt]) ;\r
 665                         parList+=onePar ;       \r
 666                 }               \r
 667         }\r
 668         \r
 669         //Get parameters set in base class.\r
 670         parList += GetBaseParametersList() ;\r
 671         \r
 672         //Get parameters set in IC class.\r
 673         if(!fMakeSeveralIC)parList += GetIsolationCut()->GetICParametersList() ;\r
 674         \r
 675         TObjString *oString= new TObjString(parList) ;\r
 676         outputContainer->Add(oString);\r
 677         \r
 678         \r
 679         return outputContainer ;\r
 680         \r
 681 }\r
 682 \r
 683 //__________________________________________________________________\r
 684 void  AliAnaParticleIsolation::MakeAnalysisFillAOD() \r
 685 {\r
 686         //Do analysis and fill aods\r
 687         //Search for the isolated photon in fCalorimeter with pt > GetMinPt()\r
 688         \r
 689         if(!GetInputAODBranch())\r
 690                 AliFatal(Form("ParticleIsolation::FillAOD: No input particles in AOD with name branch < %s > \n",GetInputAODName().Data()));\r
 691         \r
 692         Int_t n = 0, nfrac = 0;\r
 693         Bool_t isolated = kFALSE ; \r
 694         Float_t coneptsum = 0 ;\r
 695         TClonesArray * pl = new TClonesArray; \r
 696         \r
 697         //Select the calorimeter for candidate isolation with neutral particles\r
 698         if(fCalorimeter == "PHOS")\r
 699                 pl = GetAODPHOS();\r
 700         else if (fCalorimeter == "EMCAL")\r
 701                 pl = GetAODEMCAL();\r
 702         \r
 703         //Get vertex for photon momentum calculation\r
 704         if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(fVertex);\r
 705     \r
 706         //Loop on AOD branch, filled previously in AliAnaPhoton\r
 707         TLorentzVector mom ;\r
 708         \r
 709         for(Int_t iaod = 0; iaod < GetInputAODBranch()->GetEntriesFast(); iaod++){\r
 710                 AliAODPWG4ParticleCorrelation * aod =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));\r
 711                                 \r
 712                 //If too small or too large pt, skip\r
 713                 if(aod->Pt() < GetMinPt() || aod->Pt() > GetMaxPt() ) continue ; \r
 714                 \r
 715                 //Check invariant mass, if pi0, skip.\r
 716                 Bool_t decay = kFALSE ;\r
 717                 if(fMakeInvMass) decay = CheckInvMass(iaod,aod);\r
 718                 if(decay) continue ;\r
 719                 \r
 720                 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;\r
 721                 GetIsolationCut()->MakeIsolationCut((TSeqCollection*)GetAODCTS(), (TSeqCollection*)pl, \r
 722                                                                                                 fVertex, kTRUE, aod, n,nfrac,coneptsum, isolated);\r
 723                 aod->SetIsolated(isolated);\r
 724                 \r
 725         }//loop\r
 726         \r
 727         if(GetDebug() > 1) printf("End fill AODs ");  \r
 728         \r
 729 }\r
 730 \r
 731 //__________________________________________________________________\r
 732 void  AliAnaParticleIsolation::MakeAnalysisFillHistograms() \r
 733 {\r
 734         //Do analysis and fill histograms\r
 735         Int_t n = 0, nfrac = 0;\r
 736         Bool_t isolated = kFALSE ; \r
 737         Float_t coneptsum = 0 ;\r
 738         //cout<<"MakeAnalysisFillHistograms"<<endl;\r
 739         \r
 740         //Loop on stored AOD \r
 741         Int_t naod = GetInputAODBranch()->GetEntriesFast();\r
 742         if(GetDebug() > 0) printf("histo aod branch entries %d\n", naod);\r
 743         for(Int_t iaod = 0; iaod < naod ; iaod++){\r
 744                 AliAODPWG4ParticleCorrelation* aod =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));\r
 745                 Bool_t isolation = aod->IsIsolated();              \r
 746                 Float_t ptcluster  = aod->Pt();\r
 747                 Float_t phicluster = aod->Phi();\r
 748                 Float_t etacluster = aod->Eta();\r
 749                 \r
 750                 //is pt too small skip it\r
 751                 if(ptcluster < GetMinPt() || ptcluster > GetMaxPt() ) continue ; \r
 752                                 \r
 753                 if(fMakeSeveralIC) {\r
 754                         //Analysis of multiple IC at same time\r
 755                         MakeSeveralICAnalysis(aod);\r
 756                         continue;\r
 757                 }\r
 758                 else if(fReMakeIC){\r
 759                         //In case a more strict IC is needed in the produced AOD\r
 760                         n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;\r
 761                         GetIsolationCut()->MakeIsolationCut((TSeqCollection*)aod->GetRefIsolationConeTracks(), \r
 762                                                                                                 (TSeqCollection*)aod->GetRefIsolationConeClusters(), \r
 763                                                                                                 fVertex, kFALSE, aod,\r
 764                                                                                                 n,nfrac,coneptsum, isolated);\r
 765                         fhConeSumPt->Fill(ptcluster,coneptsum);       \r
 766                 }\r
 767                 \r
 768                 //Fill pt distribution of particles in cone\r
 769                 //Tracks\r
 770                 coneptsum=0;\r
 771                 for(Int_t itrack=0; itrack < (aod->GetRefIsolationConeTracks())->GetEntriesFast(); itrack++){\r
 772                         AliAODTrack* track = (AliAODTrack *)((aod->GetRefIsolationConeTracks())->At(itrack));\r
 773                         fhPtInCone->Fill(ptcluster,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));\r
 774                         coneptsum+=track->Pt();\r
 775                 }\r
 776                 //CaloClusters\r
 777                 TLorentzVector mom ;\r
 778                 for(Int_t icalo=0; icalo < (aod->GetRefIsolationConeClusters())->GetEntriesFast(); icalo++){\r
 779                         AliAODCaloCluster* calo = (AliAODCaloCluster *)((aod->GetRefIsolationConeClusters())->At(icalo));\r
 780                         calo->GetMomentum(mom,fVertex);//Assume that come from vertex in straight line\r
 781                         fhPtInCone->Fill(ptcluster, mom.Pt());\r
 782                         coneptsum+=mom.Pt();\r
 783                 }\r
 784                 \r
 785                 if(!fReMakeIC) fhConeSumPt->Fill(ptcluster,coneptsum);\r
 786                 \r
 787                 if(isolation){    \r
 788                         fhPtIso  ->Fill(ptcluster);\r
 789                         fhPhiIso ->Fill(ptcluster,phicluster);\r
 790                         fhEtaIso ->Fill(ptcluster,etacluster);\r
 791                         \r
 792                         if(IsDataMC()){\r
 793                                 Int_t tag =aod->GetTag();\r
 794                                 \r
 795                                 if(tag == AliCaloPID::kMCPrompt){\r
 796                                         fhPtIsoPrompt  ->Fill(ptcluster);\r
 797                                         fhPhiIsoPrompt ->Fill(ptcluster,phicluster);\r
 798                                         fhEtaIsoPrompt ->Fill(ptcluster,etacluster);\r
 799                                 }\r
 800                                 else if(tag == AliCaloPID::kMCFragmentation)\r
 801                                 {\r
 802                                         fhPtIsoFragmentation  ->Fill(ptcluster);\r
 803                                         fhPhiIsoFragmentation ->Fill(ptcluster,phicluster);\r
 804                                         fhEtaIsoFragmentation ->Fill(ptcluster,etacluster);\r
 805                                 }\r
 806                                 else if(tag == AliCaloPID::kMCPi0Decay)\r
 807                                 {\r
 808                                         fhPtIsoPi0Decay  ->Fill(ptcluster);\r
 809                                         fhPhiIsoPi0Decay ->Fill(ptcluster,phicluster);\r
 810                                         fhEtaIsoPi0Decay ->Fill(ptcluster,etacluster);\r
 811                                 }\r
 812                                 else if(tag == AliCaloPID::kMCEtaDecay || tag == AliCaloPID::kMCOtherDecay)\r
 813                                 {\r
 814                                         fhPtIsoOtherDecay  ->Fill(ptcluster);\r
 815                                         fhPhiIsoOtherDecay ->Fill(ptcluster,phicluster);\r
 816                                         fhEtaIsoOtherDecay ->Fill(ptcluster,etacluster);\r
 817                                 }\r
 818                                 else if(tag == AliCaloPID::kMCConversion)\r
 819                                 {\r
 820                                         fhPtIsoConversion  ->Fill(ptcluster);\r
 821                                         fhPhiIsoConversion ->Fill(ptcluster,phicluster);\r
 822                                         fhEtaIsoConversion ->Fill(ptcluster,etacluster);\r
 823                                 }\r
 824                                 else\r
 825                                 {\r
 826                                         fhPtIsoUnknown  ->Fill(ptcluster);\r
 827                                         fhPhiIsoUnknown ->Fill(ptcluster,phicluster);\r
 828                                         fhEtaIsoUnknown ->Fill(ptcluster,etacluster);\r
 829                                 }\r
 830                         }//Histograms with MC\r
 831                         \r
 832                 }//Isolated histograms\r
 833                 \r
 834         }// aod loop\r
 835         \r
 836 }\r
 837 \r
 838 //____________________________________________________________________________\r
 839 void AliAnaParticleIsolation::InitParameters()\r
 840 {\r
 841   \r
 842   //Initialize the parameters of the analysis.\r
 843   SetInputAODName("photons");\r
 844   \r
 845   fCalorimeter = "PHOS" ;\r
 846   fReMakeIC = kFALSE ;\r
 847   fMakeSeveralIC = kFALSE ;\r
 848   fMakeInvMass = kFALSE ;\r
 849   \r
 850  //----------- Several IC-----------------\r
 851   fNCones           = 5 ; \r
 852   fNPtThresFrac     = 5 ; \r
 853   fConeSizes[0]    = 0.1;  fConeSizes[1]     = 0.2;  fConeSizes[2]    = 0.3; fConeSizes[3]    = 0.4;  fConeSizes[4]    = 0.5;\r
 854   fPtThresholds[0] = 1.;   fPtThresholds[1] = 2.;    fPtThresholds[2] = 3.;  fPtThresholds[3] = 4.;   fPtThresholds[4] = 5.; \r
 855   fPtFractions[0]  = 0.05; fPtFractions[1]  = 0.075; fPtFractions[2]  = 0.1; fPtFractions[3]  = 1.25; fPtFractions[4]  = 1.5; \r
 856 \r
 857 //------------- Histograms settings -------\r
 858   fHistoNPtSumBins = 100 ;\r
 859   fHistoPtSumMax   = 50 ;\r
 860   fHistoPtSumMin   = 0.  ;\r
 861 \r
 862   fHistoNPtInConeBins = 100 ;\r
 863   fHistoPtInConeMax   = 50 ;\r
 864   fHistoPtInConeMin   = 0.  ;\r
 865 \r
 866 }\r
 867 \r
 868 //__________________________________________________________________\r
 869 void  AliAnaParticleIsolation::MakeSeveralICAnalysis(AliAODPWG4ParticleCorrelation* ph) \r
 870 {\r
 871         //Isolation Cut Analysis for both methods and different pt cuts and cones\r
 872         Float_t ptC = ph->Pt(); \r
 873         Int_t tag   = ph->GetTag();\r
 874         \r
 875         //Keep original setting used when filling AODs, reset at end of analysis  \r
 876         Float_t ptthresorg = GetIsolationCut()->GetPtThreshold();\r
 877         Float_t ptfracorg  = GetIsolationCut()->GetPtFraction();\r
 878         Float_t rorg       = GetIsolationCut()->GetConeSize();\r
 879         \r
 880         Float_t coneptsum = 0 ;  \r
 881         Int_t n[10][10];//[fNCones][fNPtThresFrac];\r
 882         Int_t nfrac[10][10];//[fNCones][fNPtThresFrac];\r
 883         Bool_t  isolated   = kFALSE;\r
 884         \r
 885         //Loop on cone sizes\r
 886         for(Int_t icone = 0; icone<fNCones; icone++){\r
 887                 GetIsolationCut()->SetConeSize(fConeSizes[icone]);\r
 888                 coneptsum = 0 ;\r
 889                 //Loop on ptthresholds\r
 890                 for(Int_t ipt = 0; ipt<fNPtThresFrac ;ipt++){\r
 891                         n[icone][ipt]=0;\r
 892                         nfrac[icone][ipt]=0;\r
 893                         GetIsolationCut()->SetPtThreshold(fPtThresholds[ipt]);\r
 894                         GetIsolationCut()->MakeIsolationCut((TSeqCollection*)ph->GetRefIsolationConeTracks(), \r
 895                                                                                                 (TSeqCollection*)ph->GetRefIsolationConeClusters(), \r
 896                                                                                                 fVertex, kFALSE, ph,\r
 897                                                                                                 n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);\r
 898                         //Normal ptThreshold cut\r
 899                         if(n[icone][ipt] == 0) {\r
 900                                 fhPtThresIsolated[icone][ipt]->Fill(ptC);\r
 901                                 if(IsDataMC()){\r
 902                                         if(tag == AliCaloPID::kMCPrompt) fhPtThresIsolatedPrompt[icone][ipt]->Fill(ptC) ;\r
 903                                         else if(tag == AliCaloPID::kMCConversion) fhPtThresIsolatedConversion[icone][ipt]->Fill(ptC) ;\r
 904                                         else if(tag == AliCaloPID::kMCFragmentation) fhPtThresIsolatedFragmentation[icone][ipt]->Fill(ptC) ;\r
 905                                         else if(tag == AliCaloPID::kMCPi0Decay) fhPtThresIsolatedPi0Decay[icone][ipt]->Fill(ptC) ;\r
 906                                         else if(tag == AliCaloPID::kMCOtherDecay || tag == AliCaloPID::kMCEtaDecay) fhPtThresIsolatedOtherDecay[icone][ipt]->Fill(ptC) ;\r
 907                                         else  fhPtThresIsolatedUnknown[icone][ipt]->Fill(ptC) ;\r
 908                                 }\r
 909                         }\r
 910                         \r
 911                         //Pt threshold on pt cand/ pt in cone fraction\r
 912                         if(nfrac[icone][ipt] == 0) {\r
 913                                 fhPtFracIsolated[icone][ipt]->Fill(ptC);\r
 914                                 if(IsDataMC()){\r
 915                                         if(tag == AliCaloPID::kMCPrompt) fhPtFracIsolatedPrompt[icone][ipt]->Fill(ptC) ;\r
 916                                         else if(tag == AliCaloPID::kMCConversion) fhPtFracIsolatedConversion[icone][ipt]->Fill(ptC) ;\r
 917                                         else if(tag == AliCaloPID::kMCFragmentation) fhPtFracIsolatedFragmentation[icone][ipt]->Fill(ptC) ;\r
 918                                         else if(tag == AliCaloPID::kMCPi0Decay) fhPtFracIsolatedPi0Decay[icone][ipt]->Fill(ptC) ;\r
 919                                         else if(tag == AliCaloPID::kMCOtherDecay || tag == AliCaloPID::kMCEtaDecay) fhPtFracIsolatedOtherDecay[icone][ipt]->Fill(ptC) ;\r
 920                                         else  fhPtFracIsolatedUnknown[icone][ipt]->Fill(ptC) ;\r
 921                                 }\r
 922                         }\r
 923                 }//pt thresh loop\r
 924                 \r
 925                 //Sum in cone histograms\r
 926                 fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ;\r
 927                 if(IsDataMC()){\r
 928                         if(tag == AliCaloPID::kMCPrompt) fhPtSumIsolatedPrompt[icone]->Fill(ptC,coneptsum) ;\r
 929                         else if(tag == AliCaloPID::kMCConversion) fhPtSumIsolatedConversion[icone]->Fill(ptC,coneptsum) ;\r
 930                         else if(tag == AliCaloPID::kMCFragmentation) fhPtSumIsolatedFragmentation[icone]->Fill(ptC,coneptsum) ;\r
 931                         else if(tag == AliCaloPID::kMCPi0Decay) fhPtSumIsolatedPi0Decay[icone]->Fill(ptC,coneptsum) ;\r
 932                         else if(tag == AliCaloPID::kMCOtherDecay || tag == AliCaloPID::kMCEtaDecay) fhPtSumIsolatedOtherDecay[icone]->Fill(ptC,coneptsum) ;\r
 933                         else  fhPtSumIsolatedUnknown[icone]->Fill(ptC,coneptsum) ;\r
 934                 }\r
 935                 \r
 936         }//cone size loop\r
 937         \r
 938         //Reset original parameters for AOD analysis\r
 939         GetIsolationCut()->SetPtThreshold(ptthresorg);\r
 940         GetIsolationCut()->SetPtFraction(ptfracorg);\r
 941         GetIsolationCut()->SetConeSize(rorg);\r
 942         \r
 943 }\r
 944 \r
 945 //__________________________________________________________________\r
 946 void AliAnaParticleIsolation::Print(const Option_t * opt) const\r
 947 {\r
 948         \r
 949         //Print some relevant parameters set for the analysis\r
 950         if(! opt)\r
 951                 return;\r
 952         \r
 953         printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
 954         \r
 955         printf("Min Gamma pT       =     %2.2f\n",  GetMinPt()) ;\r
 956         printf("Max Gamma pT       =     %3.2f\n",  GetMaxPt()) ;\r
 957         \r
 958         printf("ReMake Isolation          = %d \n",  fReMakeIC) ;\r
 959         printf("Make Several Isolation    = %d \n",  fMakeSeveralIC) ;\r
 960         printf("Calorimeter for isolation = %s \n",  fCalorimeter.Data()) ;\r
 961         \r
 962         if(fMakeSeveralIC){\r
 963                 printf("N Cone Sizes       =     %d\n", fNCones) ; \r
 964                 printf("Cone Sizes          =    \n") ;\r
 965                 for(Int_t i = 0; i < fNCones; i++)\r
 966                         printf("  %1.2f;",  fConeSizes[i]) ;\r
 967                 printf("    \n") ;\r
 968                 \r
 969                 printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ;\r
 970                 printf(" pT thresholds         =    \n") ;\r
 971                 for(Int_t i = 0; i < fNPtThresFrac; i++)\r
 972                         printf("   %2.2f;",  fPtThresholds[i]) ;\r
 973                 \r
 974                 printf("    \n") ;\r
 975                 \r
 976                 printf(" pT fractions         =    \n") ;\r
 977                 for(Int_t i = 0; i < fNPtThresFrac; i++)\r
 978                         printf("   %2.2f;",  fPtFractions[i]) ;\r
 979                 \r
 980         }  \r
 981         \r
 982         printf("Histograms: %3.1f < pT sum < %3.1f,  Nbin = %d\n", fHistoPtSumMin,  fHistoPtSumMax,  fHistoNPtSumBins);\r
 983         printf("Histograms: %3.1f < pT in cone < %3.1f, Nbin = %d\n", fHistoPtInConeMin, fHistoPtInConeMax, fHistoNPtInConeBins);\r
 984         \r
 985         printf("    \n") ;\r
 986         \r
 987 } \r
 988 \r