| f9cea31c |
1 | /************************************************************************** |
| 2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * |
| 3 | * * |
| 4 | * Author: The ALICE Off-line Project. * |
| 5 | * Contributors are mentioned in the code where appropriate. * |
| 6 | * * |
| 7 | * Permission to use, copy, modify and distribute this software and its * |
| 8 | * documentation strictly for non-commercial purposes is 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 | **************************************************************************/ |
| 15 | /* $Id$ */ |
| 16 | |
| 17 | /* History of cvs commits: |
| 18 | * |
| 19 | * $Log$ |
| 5e7c69e7 |
20 | * Revision 1.9 2007/11/17 16:39:49 gustavo |
| 21 | * removed deleting of not owned data and deleting of histograms which are exported to the output file (MG) |
| 22 | * |
| 3788def4 |
23 | * Revision 1.8 2007/10/29 13:48:42 gustavo |
| 24 | * Corrected coding violations |
| 25 | * |
| 4b707925 |
26 | * Revision 1.6 2007/08/17 12:40:04 schutz |
| 27 | * New analysis classes by Gustavo Conesa |
| 28 | * |
| bdcfac30 |
29 | * Revision 1.4.4.4 2007/07/26 10:32:09 schutz |
| 30 | * new analysis classes in the the new analysis framework |
| 2a1d8a29 |
31 | * |
| f9cea31c |
32 | * |
| 33 | */ |
| 34 | |
| 35 | //_________________________________________________________________________ |
| bdcfac30 |
36 | // Class for the prompt gamma analysis, isolation cut |
| f9cea31c |
37 | // |
| 38 | // Class created from old AliPHOSGammaJet |
| 39 | // (see AliRoot versions previous Release 4-09) |
| 40 | // |
| 2c8efea8 |
41 | // -- Author: Gustavo Conesa (LNF-INFN) |
| f9cea31c |
42 | ////////////////////////////////////////////////////////////////////////////// |
| bdcfac30 |
43 | |
| 44 | |
| 45 | // --- ROOT system --- |
| f9cea31c |
46 | #include <TParticle.h> |
| 47 | #include <TH2.h> |
| bdcfac30 |
48 | #include <TList.h> |
| f9cea31c |
49 | #include "Riostream.h" |
| 4b707925 |
50 | |
| 2c8efea8 |
51 | // --- Analysis system --- |
| 4b707925 |
52 | #include "AliAnaGammaDirect.h" |
| f9cea31c |
53 | #include "AliLog.h" |
| 2c8efea8 |
54 | #include "AliCaloTrackReader.h" |
| 55 | #include "AliIsolationCut.h" |
| 56 | #include "AliNeutralMesonSelection.h" |
| 4b707925 |
57 | |
| f9cea31c |
58 | ClassImp(AliAnaGammaDirect) |
| bdcfac30 |
59 | |
| f9cea31c |
60 | //____________________________________________________________________________ |
| bdcfac30 |
61 | AliAnaGammaDirect::AliAnaGammaDirect() : |
| 2c8efea8 |
62 | AliAnaBaseClass(), fDetector(""), fMakeIC(0), fReMakeIC(0), |
| 63 | fMakeSeveralIC(0), fMakeInvMass(0), |
| 64 | fhPtGamma(0),fhPhiGamma(0),fhEtaGamma(0), fhConeSumPt(0), |
| 65 | //Several IC |
| 66 | fNCones(0),fNPtThresFrac(0), fConeSizes(), fPtThresholds(), fPtFractions(), |
| 67 | fhPtThresIsolated(), fhPtFracIsolated(), fhPtSumIsolated(), |
| 68 | //MC |
| 69 | fhPtPrompt(0),fhPhiPrompt(0),fhEtaPrompt(0), |
| 70 | fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(), |
| 71 | fhPtFragmentation(0),fhPhiFragmentation(0),fhEtaFragmentation(0), |
| 72 | fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(), |
| 73 | fhPtPi0Decay(0),fhPhiPi0Decay(0),fhEtaPi0Decay(0), |
| 74 | fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(), |
| 75 | fhPtOtherDecay(0),fhPhiOtherDecay(0),fhEtaOtherDecay(0), |
| 76 | fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(), |
| 77 | fhPtConversion(0),fhPhiConversion(0),fhEtaConversion(0), |
| 78 | fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(), |
| 79 | fhPtUnknown(0),fhPhiUnknown(0),fhEtaUnknown(0), |
| 80 | fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(), fhPtSumIsolatedUnknown() |
| f9cea31c |
81 | { |
| bdcfac30 |
82 | //default ctor |
| 83 | |
| 463ee300 |
84 | //Initialize parameters |
| 85 | InitParameters(); |
| 86 | |
| f9cea31c |
87 | } |
| 88 | |
| f9cea31c |
89 | //____________________________________________________________________________ |
| 90 | AliAnaGammaDirect::AliAnaGammaDirect(const AliAnaGammaDirect & g) : |
| 2c8efea8 |
91 | AliAnaBaseClass(g), fDetector(g.fDetector), |
| 92 | fMakeIC(g.fMakeIC), fReMakeIC(g.fReMakeIC), |
| 93 | fMakeSeveralIC(g.fMakeSeveralIC), fMakeInvMass(g.fMakeInvMass), |
| 94 | fhPtGamma(g.fhPtGamma),fhPhiGamma(g.fhPhiGamma), |
| 95 | fhEtaGamma(g.fhEtaGamma), fhConeSumPt(g.fhConeSumPt), |
| 96 | //Several IC |
| 97 | fNCones(g.fNCones),fNPtThresFrac(g.fNPtThresFrac), fConeSizes(), fPtThresholds(), fPtFractions(), |
| 98 | fhPtThresIsolated(), fhPtFracIsolated(), fhPtSumIsolated(), |
| 99 | //MC |
| 100 | fhPtPrompt(g.fhPtPrompt),fhPhiPrompt(g.fhPhiPrompt),fhEtaPrompt(g.fhEtaPrompt), |
| 101 | fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(), |
| 102 | fhPtFragmentation(g.fhPtFragmentation),fhPhiFragmentation(g.fhPhiFragmentation),fhEtaFragmentation(g.fhEtaFragmentation), |
| 103 | fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(), |
| 104 | fhPtPi0Decay(g.fhPtPi0Decay),fhPhiPi0Decay(g.fhPhiPi0Decay),fhEtaPi0Decay(g.fhEtaPi0Decay), |
| 105 | fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(), |
| 106 | fhPtOtherDecay(g.fhPtOtherDecay),fhPhiOtherDecay(g.fhPhiOtherDecay),fhEtaOtherDecay(g.fhEtaOtherDecay), |
| 107 | fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(), |
| 108 | fhPtConversion(g. fhPtConversion),fhPhiConversion(g.fhPhiConversion),fhEtaConversion(g.fhEtaConversion), |
| 109 | fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(), |
| 110 | fhPtUnknown(g.fhPtUnknown),fhPhiUnknown(g.fhPhiUnknown),fhEtaUnknown(g.fhEtaUnknown), |
| 111 | fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(), fhPtSumIsolatedUnknown() |
| f9cea31c |
112 | { |
| 113 | // cpy ctor |
| bdcfac30 |
114 | |
| 2c8efea8 |
115 | //Several IC |
| bdcfac30 |
116 | for(Int_t i = 0; i < fNCones ; i++){ |
| 117 | fConeSizes[i] = g.fConeSizes[i]; |
| 118 | fhPtSumIsolated[i] = g.fhPtSumIsolated[i]; |
| 2c8efea8 |
119 | |
| 120 | fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i]; |
| 121 | fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i]; |
| 122 | fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i]; |
| 123 | fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i]; |
| 124 | fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i]; |
| 125 | fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i]; |
| 126 | |
| 127 | for(Int_t j = 0; j < fNPtThresFrac ; j++){ |
| bdcfac30 |
128 | fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j]; |
| 2c8efea8 |
129 | fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j]; |
| 130 | |
| 131 | fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j]; |
| 132 | fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j]; |
| 133 | fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j]; |
| 134 | fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j]; |
| 135 | fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j]; |
| 136 | fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j]; |
| 137 | |
| 138 | fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j]; |
| 139 | fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j]; |
| 140 | fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j]; |
| 141 | fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j]; |
| 142 | fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j]; |
| 143 | fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j]; |
| 144 | |
| 145 | } |
| bdcfac30 |
146 | } |
| 147 | |
| 2c8efea8 |
148 | for(Int_t i = 0; i < fNPtThresFrac ; i++){ |
| 149 | fPtFractions[i]= g.fPtFractions[i]; |
| bdcfac30 |
150 | fPtThresholds[i]= g.fPtThresholds[i]; |
| 2c8efea8 |
151 | } |
| 4b707925 |
152 | |
| f9cea31c |
153 | } |
| 154 | |
| 463ee300 |
155 | //_________________________________________________________________________ |
| 2c8efea8 |
156 | AliAnaGammaDirect & AliAnaGammaDirect::operator = (const AliAnaGammaDirect & g) |
| 463ee300 |
157 | { |
| 158 | // assignment operator |
| bdcfac30 |
159 | |
| 2c8efea8 |
160 | if(&g == this) return *this; |
| 161 | |
| 162 | fMakeIC = g.fMakeIC ; |
| 163 | fReMakeIC = g.fReMakeIC ; |
| 164 | fMakeSeveralIC = g.fMakeSeveralIC ; |
| 165 | fMakeInvMass = g.fMakeInvMass ; |
| 166 | fDetector = g.fDetector ; |
| 167 | |
| 168 | fhPtGamma = g.fhPtGamma ; |
| 169 | fhPhiGamma = g.fhPhiGamma ; |
| 170 | fhEtaGamma = g.fhEtaGamma ; |
| 171 | fhConeSumPt = g.fhConeSumPt ; |
| 172 | |
| 173 | fhPtPrompt = g.fhPtPrompt; |
| 174 | fhPhiPrompt = g.fhPhiPrompt; |
| 175 | fhEtaPrompt = g.fhEtaPrompt; |
| 176 | fhPtFragmentation = g.fhPtFragmentation; |
| 177 | fhPhiFragmentation = g.fhPhiFragmentation; |
| 178 | fhEtaFragmentation = g.fhEtaFragmentation; |
| 179 | fhPtPi0Decay = g.fhPtPi0Decay; |
| 180 | fhPhiPi0Decay = g.fhPhiPi0Decay; |
| 181 | fhEtaPi0Decay = g.fhEtaPi0Decay; |
| 182 | fhPtOtherDecay = g.fhPtOtherDecay; |
| 183 | fhPhiOtherDecay = g.fhPhiOtherDecay; |
| 184 | fhEtaOtherDecay = g.fhEtaOtherDecay; |
| 185 | fhPtConversion = g. fhPtConversion; |
| 186 | fhPhiConversion = g.fhPhiConversion; |
| 187 | fhEtaConversion = g.fhEtaConversion; |
| 188 | fhPtUnknown = g.fhPtUnknown; |
| 189 | fhPhiUnknown = g.fhPhiUnknown; |
| 190 | fhEtaUnknown = g.fhEtaUnknown; |
| 191 | |
| 192 | //Several IC |
| 193 | fNCones = g.fNCones ; |
| 194 | fNPtThresFrac = g.fNPtThresFrac ; |
| bdcfac30 |
195 | |
| 196 | for(Int_t i = 0; i < fNCones ; i++){ |
| 2c8efea8 |
197 | fConeSizes[i] = g.fConeSizes[i]; |
| 198 | fhPtSumIsolated[i] = g.fhPtSumIsolated[i] ; |
| 199 | |
| 200 | fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i]; |
| 201 | fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i]; |
| 202 | fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i]; |
| 203 | fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i]; |
| 204 | fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i]; |
| 205 | fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i]; |
| 206 | |
| 207 | for(Int_t j = 0; j < fNPtThresFrac ; j++){ |
| 208 | fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j] ; |
| 209 | fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j] ; |
| 210 | |
| 211 | fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j]; |
| 212 | fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j]; |
| 213 | fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j]; |
| 214 | fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j]; |
| 215 | fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j]; |
| 216 | fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j]; |
| 217 | |
| 218 | fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j]; |
| 219 | fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j]; |
| 220 | fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j]; |
| 221 | fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j]; |
| 222 | fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j]; |
| 223 | fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j]; |
| 224 | |
| 225 | } |
| bdcfac30 |
226 | } |
| 227 | |
| 2c8efea8 |
228 | for(Int_t i = 0; i < fNPtThresFrac ; i++){ |
| 229 | fPtThresholds[i]= g.fPtThresholds[i]; |
| 230 | fPtFractions[i]= g.fPtFractions[i]; |
| 231 | } |
| 4b707925 |
232 | |
| 463ee300 |
233 | return *this; |
| bdcfac30 |
234 | |
| 463ee300 |
235 | } |
| 236 | |
| f9cea31c |
237 | //____________________________________________________________________________ |
| 238 | AliAnaGammaDirect::~AliAnaGammaDirect() |
| 239 | { |
| 2c8efea8 |
240 | //dtor |
| 241 | //do not delete histograms |
| 242 | |
| 243 | delete [] fConeSizes ; |
| 244 | delete [] fPtThresholds ; |
| 245 | delete [] fPtFractions ; |
| 246 | |
| 247 | } |
| 248 | //_________________________________________________________________________ |
| 249 | Bool_t AliAnaGammaDirect::CheckInvMass(const Int_t icalo,const TLorentzVector mom, Double_t *vertex, TClonesArray * pl){ |
| 250 | //Search if there is a companion decay photon to the candidate |
| 251 | // and discard it in such case |
| 252 | TLorentzVector mom2 ; |
| 253 | for(Int_t jcalo = 0; jcalo < pl->GetEntries(); jcalo++){ |
| 254 | if(icalo==jcalo) continue ; |
| 255 | AliAODCaloCluster * calo = dynamic_cast<AliAODCaloCluster*> (pl->At(jcalo)); |
| 256 | |
| 257 | //Cluster selection, not charged, with photon id and in fidutial cut, fill TLorentz |
| 258 | if(!SelectCluster(calo, vertex, mom2)) continue ; |
| 259 | |
| 260 | //Select good pair (good phit, pt cuts, aperture and invariant mass) |
| 261 | if(GetNeutralMesonSelection()->SelectPair(mom, mom2)){ |
| 262 | if(GetDebug()>1)printf("Selected gamma pair: pt %f, phi %f, eta%f",(mom+mom2).Pt(), (mom+mom2).Phi(), (mom+mom2).Eta()); |
| 263 | return kTRUE ; |
| 264 | } |
| 265 | }//loop |
| 266 | |
| 267 | return kFALSE; |
| 4b707925 |
268 | |
| 2c8efea8 |
269 | } |
| 270 | //_________________________________________________________________________ |
| 271 | Int_t AliAnaGammaDirect::CheckOrigin(const Int_t label){ |
| 272 | //Play with the MC stack if available |
| 273 | //Check origin of the candidates, good for PYTHIA |
| 274 | |
| 275 | AliStack * stack = GetMCStack() ; |
| 276 | if(!stack) AliFatal("Stack is not available, check analysis settings in configuration file, STOP!!"); |
| 277 | |
| 278 | if(label >= 0 && label < stack->GetNtrack()){ |
| 279 | //Mother |
| 280 | TParticle * mom = stack->Particle(label); |
| 281 | Int_t mPdg = TMath::Abs(mom->GetPdgCode()); |
| 282 | Int_t mStatus = mom->GetStatusCode() ; |
| 283 | Int_t iParent = mom->GetFirstMother() ; |
| 284 | if(label < 8 ) printf("Mother is parton %d\n",iParent); |
| 285 | |
| 286 | //GrandParent |
| 287 | TParticle * parent = new TParticle ; |
| 288 | Int_t pPdg = -1; |
| 289 | Int_t pStatus =-1; |
| 290 | if(iParent > 0){ |
| 291 | parent = stack->Particle(iParent); |
| 292 | pPdg = TMath::Abs(parent->GetPdgCode()); |
| 293 | pStatus = parent->GetStatusCode(); |
| 294 | } |
| 295 | else |
| 296 | printf("Parent with label %d\n",iParent); |
| 297 | |
| 298 | //return tag |
| 299 | if(mPdg == 22){ |
| 300 | if(mStatus == 1){ |
| 301 | if(iParent < 8) { |
| 302 | if(pPdg == 22) return kPrompt; |
| 303 | else return kFragmentation; |
| 304 | } |
| 305 | else if(pStatus == 11){ |
| 306 | if(pPdg == 111) return kPi0Decay ; |
| 307 | else if (pPdg == 321) return kEtaDecay ; |
| 308 | else return kOtherDecay ; |
| 309 | } |
| 310 | }//Status 1 : Pythia generated |
| 311 | else if(mStatus == 0){ |
| 312 | if(pPdg ==22 || pPdg ==11) return kConversion ; |
| 313 | if(pPdg == 111) return kPi0Decay ; |
| 314 | else if (pPdg == 221) return kEtaDecay ; |
| 315 | else return kOtherDecay ; |
| 316 | }//status 0 : geant generated |
| 317 | }//Mother gamma |
| 318 | else if(mPdg == 111) return kPi0 ; |
| 319 | else if(mPdg == 221) return kEta ; |
| 320 | else if(mPdg ==11){ |
| 321 | if(mStatus == 0) return kConversion ; |
| 322 | else return kElectron ; |
| 323 | } |
| 324 | else return kUnknown; |
| 325 | }//Good label value |
| 326 | else{ |
| 327 | if(label < 0 ) printf("*** bad label or no stack ***: label %d \n", label); |
| 328 | if(label >= stack->GetNtrack()) printf("*** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack()); |
| 329 | return kUnknown; |
| 330 | }//Bad label |
| 331 | |
| 332 | return kUnknown; |
| 333 | |
| 2a1d8a29 |
334 | } |
| 335 | |
| 336 | //________________________________________________________________________ |
| bdcfac30 |
337 | TList * AliAnaGammaDirect::GetCreateOutputObjects() |
| 2a1d8a29 |
338 | { |
| 463ee300 |
339 | // Create histograms to be saved in output file and |
| bdcfac30 |
340 | // store them in outputContainer |
| 341 | TList * outputContainer = new TList() ; |
| 342 | outputContainer->SetName("DirectGammaHistos") ; |
| 2c8efea8 |
343 | |
| 2a1d8a29 |
344 | //Histograms of highest gamma identified in Event |
| 2c8efea8 |
345 | fhPtGamma = new TH1F("hPtGamma","Number of #gamma over calorimeter",240,0,120); |
| 346 | fhPtGamma->SetYTitle("N"); |
| 347 | fhPtGamma->SetXTitle("p_{T #gamma}(GeV/c)"); |
| 348 | outputContainer->Add(fhPtGamma) ; |
| 2a1d8a29 |
349 | |
| 350 | fhPhiGamma = new TH2F |
| 2c8efea8 |
351 | ("hPhiGamma","#phi_{#gamma}",200,0,120,200,0,7); |
| 2a1d8a29 |
352 | fhPhiGamma->SetYTitle("#phi"); |
| 353 | fhPhiGamma->SetXTitle("p_{T #gamma} (GeV/c)"); |
| bdcfac30 |
354 | outputContainer->Add(fhPhiGamma) ; |
| 2a1d8a29 |
355 | |
| 356 | fhEtaGamma = new TH2F |
| 2c8efea8 |
357 | ("hEtaGamma","#phi_{#gamma}",200,0,120,200,-0.8,0.8); |
| 2a1d8a29 |
358 | fhEtaGamma->SetYTitle("#eta"); |
| 359 | fhEtaGamma->SetXTitle("p_{T #gamma} (GeV/c)"); |
| bdcfac30 |
360 | outputContainer->Add(fhEtaGamma) ; |
| 361 | |
| 2c8efea8 |
362 | if(!fMakeSeveralIC){ |
| 363 | fhConeSumPt = new TH2F |
| 364 | ("hConePtSum","#Sigma p_{T} in cone ",200,0,120,100,0,100); |
| 365 | fhConeSumPt->SetYTitle("#Sigma p_{T}"); |
| 366 | fhConeSumPt->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 367 | outputContainer->Add(fhConeSumPt) ; |
| 368 | } |
| 4b707925 |
369 | |
| 2c8efea8 |
370 | if(IsDataMC()){ |
| 371 | |
| 372 | fhPtPrompt = new TH1F("hPtPrompt","Number of #gamma over calorimeter",240,0,120); |
| 373 | fhPtPrompt->SetYTitle("N"); |
| 374 | fhPtPrompt->SetXTitle("p_{T #gamma}(GeV/c)"); |
| 375 | outputContainer->Add(fhPtPrompt) ; |
| 376 | |
| 377 | fhPhiPrompt = new TH2F |
| 378 | ("hPhiPrompt","#phi_{#gamma}",200,0,120,200,0,7); |
| 379 | fhPhiPrompt->SetYTitle("#phi"); |
| 380 | fhPhiPrompt->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 381 | outputContainer->Add(fhPhiPrompt) ; |
| 382 | |
| 383 | fhEtaPrompt = new TH2F |
| 384 | ("hEtaPrompt","#phi_{#gamma}",200,0,120,200,-0.8,0.8); |
| 385 | fhEtaPrompt->SetYTitle("#eta"); |
| 386 | fhEtaPrompt->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 387 | outputContainer->Add(fhEtaPrompt) ; |
| 388 | |
| 389 | fhPtFragmentation = new TH1F("hPtFragmentation","Number of #gamma over calorimeter",240,0,120); |
| 390 | fhPtFragmentation->SetYTitle("N"); |
| 391 | fhPtFragmentation->SetXTitle("p_{T #gamma}(GeV/c)"); |
| 392 | outputContainer->Add(fhPtFragmentation) ; |
| 393 | |
| 394 | fhPhiFragmentation = new TH2F |
| 395 | ("hPhiFragmentation","#phi_{#gamma}",200,0,120,200,0,7); |
| 396 | fhPhiFragmentation->SetYTitle("#phi"); |
| 397 | fhPhiFragmentation->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 398 | outputContainer->Add(fhPhiFragmentation) ; |
| 399 | |
| 400 | fhEtaFragmentation = new TH2F |
| 401 | ("hEtaFragmentation","#phi_{#gamma}",200,0,120,200,-0.8,0.8); |
| 402 | fhEtaFragmentation->SetYTitle("#eta"); |
| 403 | fhEtaFragmentation->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 404 | outputContainer->Add(fhEtaFragmentation) ; |
| 405 | |
| 406 | fhPtPi0Decay = new TH1F("hPtPi0Decay","Number of #gamma over calorimeter",240,0,120); |
| 407 | fhPtPi0Decay->SetYTitle("N"); |
| 408 | fhPtPi0Decay->SetXTitle("p_{T #gamma}(GeV/c)"); |
| 409 | outputContainer->Add(fhPtPi0Decay) ; |
| 410 | |
| 411 | fhPhiPi0Decay = new TH2F |
| 412 | ("hPhiPi0Decay","#phi_{#gamma}",200,0,120,200,0,7); |
| 413 | fhPhiPi0Decay->SetYTitle("#phi"); |
| 414 | fhPhiPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 415 | outputContainer->Add(fhPhiPi0Decay) ; |
| 416 | |
| 417 | fhEtaPi0Decay = new TH2F |
| 418 | ("hEtaPi0Decay","#phi_{#gamma}",200,0,120,200,-0.8,0.8); |
| 419 | fhEtaPi0Decay->SetYTitle("#eta"); |
| 420 | fhEtaPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 421 | outputContainer->Add(fhEtaPi0Decay) ; |
| 422 | |
| 423 | fhPtOtherDecay = new TH1F("hPtOtherDecay","Number of #gamma over calorimeter",240,0,120); |
| 424 | fhPtOtherDecay->SetYTitle("N"); |
| 425 | fhPtOtherDecay->SetXTitle("p_{T #gamma}(GeV/c)"); |
| 426 | outputContainer->Add(fhPtOtherDecay) ; |
| 427 | |
| 428 | fhPhiOtherDecay = new TH2F |
| 429 | ("hPhiOtherDecay","#phi_{#gamma}",200,0,120,200,0,7); |
| 430 | fhPhiOtherDecay->SetYTitle("#phi"); |
| 431 | fhPhiOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 432 | outputContainer->Add(fhPhiOtherDecay) ; |
| 433 | |
| 434 | fhEtaOtherDecay = new TH2F |
| 435 | ("hEtaOtherDecay","#phi_{#gamma}",200,0,120,200,-0.8,0.8); |
| 436 | fhEtaOtherDecay->SetYTitle("#eta"); |
| 437 | fhEtaOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 438 | outputContainer->Add(fhEtaOtherDecay) ; |
| 439 | |
| 440 | fhPtConversion = new TH1F("hPtConversion","Number of #gamma over calorimeter",240,0,120); |
| 441 | fhPtConversion->SetYTitle("N"); |
| 442 | fhPtConversion->SetXTitle("p_{T #gamma}(GeV/c)"); |
| 443 | outputContainer->Add(fhPtConversion) ; |
| 444 | |
| 445 | fhPhiConversion = new TH2F |
| 446 | ("hPhiConversion","#phi_{#gamma}",200,0,120,200,0,7); |
| 447 | fhPhiConversion->SetYTitle("#phi"); |
| 448 | fhPhiConversion->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 449 | outputContainer->Add(fhPhiConversion) ; |
| 450 | |
| 451 | fhEtaConversion = new TH2F |
| 452 | ("hEtaConversion","#phi_{#gamma}",200,0,120,200,-0.8,0.8); |
| 453 | fhEtaConversion->SetYTitle("#eta"); |
| 454 | fhEtaConversion->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 455 | outputContainer->Add(fhEtaConversion) ; |
| 456 | |
| 457 | fhPtUnknown = new TH1F("hPtUnknown","Number of #gamma over calorimeter",240,0,120); |
| 458 | fhPtUnknown->SetYTitle("N"); |
| 459 | fhPtUnknown->SetXTitle("p_{T #gamma}(GeV/c)"); |
| 460 | outputContainer->Add(fhPtUnknown) ; |
| 461 | |
| 462 | fhPhiUnknown = new TH2F |
| 463 | ("hPhiUnknown","#phi_{#gamma}",200,0,120,200,0,7); |
| 464 | fhPhiUnknown->SetYTitle("#phi"); |
| 465 | fhPhiUnknown->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 466 | outputContainer->Add(fhPhiUnknown) ; |
| 467 | |
| 468 | fhEtaUnknown = new TH2F |
| 469 | ("hEtaUnknown","#phi_{#gamma}",200,0,120,200,-0.8,0.8); |
| 470 | fhEtaUnknown->SetYTitle("#eta"); |
| 471 | fhEtaUnknown->SetXTitle("p_{T #gamma} (GeV/c)"); |
| 472 | outputContainer->Add(fhEtaUnknown) ; |
| 473 | }//Histos with MC |
| 474 | |
| 475 | if(fMakeSeveralIC){ |
| bdcfac30 |
476 | char name[128]; |
| 477 | char title[128]; |
| 478 | for(Int_t icone = 0; icone<fNCones; icone++){ |
| 2c8efea8 |
479 | sprintf(name,"hPtSumIsolated_Cone_%d",icone); |
| bdcfac30 |
480 | sprintf(title,"Candidate cone sum p_{T} for cone size %d vs candidate p_{T}",icone); |
| 481 | fhPtSumIsolated[icone] = new TH2F(name, title,240,0,120,120,0,10); |
| 482 | fhPtSumIsolated[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); |
| 483 | fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)"); |
| 484 | outputContainer->Add(fhPtSumIsolated[icone]) ; |
| 2c8efea8 |
485 | |
| 486 | if(IsDataMC()){ |
| 487 | sprintf(name,"hPtSumIsolatedPrompt_Cone_%d",icone); |
| 488 | sprintf(title,"Candidate Prompt cone sum p_{T} for cone size %d vs candidate p_{T}",icone); |
| 489 | fhPtSumIsolatedPrompt[icone] = new TH2F(name, title,240,0,120,120,0,10); |
| 490 | fhPtSumIsolatedPrompt[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); |
| 491 | fhPtSumIsolatedPrompt[icone]->SetXTitle("p_{T} (GeV/c)"); |
| 492 | outputContainer->Add(fhPtSumIsolatedPrompt[icone]) ; |
| 493 | |
| 494 | sprintf(name,"hPtSumIsolatedFragmentation_Cone_%d",icone); |
| 495 | sprintf(title,"Candidate Fragmentation cone sum p_{T} for cone size %d vs candidate p_{T}",icone); |
| 496 | fhPtSumIsolatedFragmentation[icone] = new TH2F(name, title,240,0,120,120,0,10); |
| 497 | fhPtSumIsolatedFragmentation[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); |
| 498 | fhPtSumIsolatedFragmentation[icone]->SetXTitle("p_{T} (GeV/c)"); |
| 499 | outputContainer->Add(fhPtSumIsolatedFragmentation[icone]) ; |
| 500 | |
| 501 | sprintf(name,"hPtSumIsolatedPi0Decay_Cone_%d",icone); |
| 502 | sprintf(title,"Candidate Pi0Decay cone sum p_{T} for cone size %d vs candidate p_{T}",icone); |
| 503 | fhPtSumIsolatedPi0Decay[icone] = new TH2F(name, title,240,0,120,120,0,10); |
| 504 | fhPtSumIsolatedPi0Decay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); |
| 505 | fhPtSumIsolatedPi0Decay[icone]->SetXTitle("p_{T} (GeV/c)"); |
| 506 | outputContainer->Add(fhPtSumIsolatedPi0Decay[icone]) ; |
| 507 | |
| 508 | sprintf(name,"hPtSumIsolatedOtherDecay_Cone_%d",icone); |
| 509 | sprintf(title,"Candidate OtherDecay cone sum p_{T} for cone size %d vs candidate p_{T}",icone); |
| 510 | fhPtSumIsolatedOtherDecay[icone] = new TH2F(name, title,240,0,120,120,0,10); |
| 511 | fhPtSumIsolatedOtherDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); |
| 512 | fhPtSumIsolatedOtherDecay[icone]->SetXTitle("p_{T} (GeV/c)"); |
| 513 | outputContainer->Add(fhPtSumIsolatedOtherDecay[icone]) ; |
| 514 | |
| 515 | sprintf(name,"hPtSumIsolatedConversion_Cone_%d",icone); |
| 516 | sprintf(title,"Candidate Conversion cone sum p_{T} for cone size %d vs candidate p_{T}",icone); |
| 517 | fhPtSumIsolatedConversion[icone] = new TH2F(name, title,240,0,120,120,0,10); |
| 518 | fhPtSumIsolatedConversion[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); |
| 519 | fhPtSumIsolatedConversion[icone]->SetXTitle("p_{T} (GeV/c)"); |
| 520 | outputContainer->Add(fhPtSumIsolatedConversion[icone]) ; |
| 521 | |
| 522 | sprintf(name,"hPtSumIsolatedUnknown_Cone_%d",icone); |
| 523 | sprintf(title,"Candidate Unknown cone sum p_{T} for cone size %d vs candidate p_{T}",icone); |
| 524 | fhPtSumIsolatedUnknown[icone] = new TH2F(name, title,240,0,120,120,0,10); |
| 525 | fhPtSumIsolatedUnknown[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); |
| 526 | fhPtSumIsolatedUnknown[icone]->SetXTitle("p_{T} (GeV/c)"); |
| 527 | outputContainer->Add(fhPtSumIsolatedUnknown[icone]) ; |
| 528 | |
| 529 | }//Histos with MC |
| 530 | |
| 531 | for(Int_t ipt = 0; ipt<fNPtThresFrac;ipt++){ |
| 532 | sprintf(name,"hPtThresIsol_Cone_%d_Pt%d",icone,ipt); |
| bdcfac30 |
533 | sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 534 | fhPtThresIsolated[icone][ipt] = new TH1F(name, title,240,0,120); |
| 535 | fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 536 | outputContainer->Add(fhPtThresIsolated[icone][ipt]) ; |
| 2c8efea8 |
537 | |
| 538 | sprintf(name,"hPtFracIsol_Cone_%d_Pt%d",icone,ipt); |
| 539 | sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 540 | fhPtFracIsolated[icone][ipt] = new TH1F(name, title,240,0,120); |
| 541 | fhPtFracIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 542 | outputContainer->Add(fhPtFracIsolated[icone][ipt]) ; |
| 543 | |
| 544 | if(IsDataMC()){ |
| 545 | sprintf(name,"hPtThresIsolPrompt_Cone_%d_Pt%d",icone,ipt); |
| 546 | sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 547 | fhPtThresIsolatedPrompt[icone][ipt] = new TH1F(name, title,240,0,120); |
| 548 | fhPtThresIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 549 | outputContainer->Add(fhPtThresIsolatedPrompt[icone][ipt]) ; |
| 550 | |
| 551 | sprintf(name,"hPtFracIsolPrompt_Cone_%d_Pt%d",icone,ipt); |
| 552 | sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 553 | fhPtFracIsolatedPrompt[icone][ipt] = new TH1F(name, title,240,0,120); |
| 554 | fhPtFracIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 555 | outputContainer->Add(fhPtFracIsolatedPrompt[icone][ipt]) ; |
| 556 | |
| 557 | sprintf(name,"hPtThresIsolFragmentation_Cone_%d_Pt%d",icone,ipt); |
| 558 | sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 559 | fhPtThresIsolatedFragmentation[icone][ipt] = new TH1F(name, title,240,0,120); |
| 560 | fhPtThresIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 561 | outputContainer->Add(fhPtThresIsolatedFragmentation[icone][ipt]) ; |
| 562 | |
| 563 | sprintf(name,"hPtFracIsolFragmentation_Cone_%d_Pt%d",icone,ipt); |
| 564 | sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 565 | fhPtFracIsolatedFragmentation[icone][ipt] = new TH1F(name, title,240,0,120); |
| 566 | fhPtFracIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 567 | outputContainer->Add(fhPtFracIsolatedFragmentation[icone][ipt]) ; |
| 568 | |
| 569 | sprintf(name,"hPtThresIsolPi0Decay_Cone_%d_Pt%d",icone,ipt); |
| 570 | sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 571 | fhPtThresIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,240,0,120); |
| 572 | fhPtThresIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 573 | outputContainer->Add(fhPtThresIsolatedPi0Decay[icone][ipt]) ; |
| 574 | |
| 575 | sprintf(name,"hPtFracIsolPi0Decay_Cone_%d_Pt%d",icone,ipt); |
| 576 | sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 577 | fhPtFracIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,240,0,120); |
| 578 | fhPtFracIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 579 | outputContainer->Add(fhPtFracIsolatedPi0Decay[icone][ipt]) ; |
| 580 | |
| 581 | sprintf(name,"hPtThresIsolOtherDecay_Cone_%d_Pt%d",icone,ipt); |
| 582 | sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 583 | fhPtThresIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,240,0,120); |
| 584 | fhPtThresIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 585 | outputContainer->Add(fhPtThresIsolatedOtherDecay[icone][ipt]) ; |
| 586 | |
| 587 | sprintf(name,"hPtFracIsolOtherDecay_Cone_%d_Pt%d",icone,ipt); |
| 588 | sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 589 | fhPtFracIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,240,0,120); |
| 590 | fhPtFracIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 591 | outputContainer->Add(fhPtFracIsolatedOtherDecay[icone][ipt]) ; |
| 592 | |
| 593 | sprintf(name,"hPtThresIsolConversion_Cone_%d_Pt%d",icone,ipt); |
| 594 | sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 595 | fhPtThresIsolatedConversion[icone][ipt] = new TH1F(name, title,240,0,120); |
| 596 | fhPtThresIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 597 | outputContainer->Add(fhPtThresIsolatedConversion[icone][ipt]) ; |
| 598 | |
| 599 | sprintf(name,"hPtFracIsolConversion_Cone_%d_Pt%d",icone,ipt); |
| 600 | sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 601 | fhPtFracIsolatedConversion[icone][ipt] = new TH1F(name, title,240,0,120); |
| 602 | fhPtFracIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 603 | outputContainer->Add(fhPtFracIsolatedConversion[icone][ipt]) ; |
| 604 | |
| 605 | sprintf(name,"hPtThresIsolUnknown_Cone_%d_Pt%d",icone,ipt); |
| 606 | sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 607 | fhPtThresIsolatedUnknown[icone][ipt] = new TH1F(name, title,240,0,120); |
| 608 | fhPtThresIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 609 | outputContainer->Add(fhPtThresIsolatedUnknown[icone][ipt]) ; |
| 610 | |
| 611 | sprintf(name,"hPtFracIsolUnknown_Cone_%d_Pt%d",icone,ipt); |
| 612 | sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); |
| 613 | fhPtFracIsolatedUnknown[icone][ipt] = new TH1F(name, title,240,0,120); |
| 614 | fhPtFracIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); |
| 615 | outputContainer->Add(fhPtFracIsolatedUnknown[icone][ipt]) ; |
| 616 | |
| 617 | }//Histos with MC |
| 618 | |
| bdcfac30 |
619 | }//icone loop |
| 620 | }//ipt loop |
| f9cea31c |
621 | } |
| 622 | |
| 2c8efea8 |
623 | //Keep neutral meson selection histograms if requiered |
| 624 | //Setting done in AliNeutralMesonSelection |
| 625 | if(fMakeInvMass && GetNeutralMesonSelection()){ |
| 626 | TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ; |
| 627 | cout<<"NMSHistos "<< nmsHistos<<endl; |
| 628 | if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept()) |
| 629 | for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ; |
| 630 | } |
| 4b707925 |
631 | |
| bdcfac30 |
632 | return outputContainer ; |
| 463ee300 |
633 | |
| f9cea31c |
634 | } |
| 635 | |
| 2c8efea8 |
636 | //__________________________________________________________________ |
| 637 | void AliAnaGammaDirect::MakeAnalysisFillAOD() |
| f9cea31c |
638 | { |
| 2c8efea8 |
639 | //Do analysis and fill aods |
| 640 | //Search for the isolated photon in fDetector with pt > GetMinPt() |
| f9cea31c |
641 | |
| 2c8efea8 |
642 | //Fill CaloClusters if working with ESDs |
| 643 | //if(GetReader()->GetDataType() == AliCaloTrackReader::kESD) ConnectAODCaloClusters(); |
| 4b707925 |
644 | |
| 2c8efea8 |
645 | Int_t n = 0, nfrac = 0; |
| 646 | Bool_t isolated = kFALSE ; |
| 647 | Float_t coneptsum = 0 ; |
| 648 | TClonesArray * pl = new TClonesArray; |
| 649 | |
| 650 | //Get vertex for photon momentum calculation |
| 651 | Double_t vertex[]={0,0,0} ; //vertex ; |
| 652 | if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex); |
| 653 | |
| 654 | //Select the detector of the photon |
| 655 | if(fDetector == "PHOS") |
| 656 | pl = GetAODPHOS(); |
| 657 | else if (fDetector == "EMCAL") |
| 658 | pl = GetAODEMCAL(); |
| 659 | //cout<<"Number of entries "<<pl->GetEntries()<<endl; |
| f9cea31c |
660 | |
| 2c8efea8 |
661 | //Fill AODCaloClusters and AODParticleCorrelation with PHOS aods |
| 662 | TLorentzVector mom ; |
| 663 | for(Int_t icalo = 0; icalo < pl->GetEntries(); icalo++){ |
| 664 | AliAODCaloCluster * calo = dynamic_cast<AliAODCaloCluster*> (pl->At(icalo)); |
| 665 | |
| 666 | //Cluster selection, not charged, with photon id and in fidutial cut |
| 667 | if(!SelectCluster(calo,vertex,mom)) continue ; |
| 668 | |
| 669 | //If too small pt, skip |
| 670 | if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ; |
| 671 | |
| 672 | //Play with the MC stack if available |
| 673 | Int_t tag =-1; |
| 674 | if(IsDataMC()){ |
| 675 | //Check origin of the candidates |
| 676 | tag = CheckOrigin(calo->GetLabel(0)); |
| 677 | if(GetDebug() > 0) printf("Origin of candidate %d\n",tag); |
| 678 | }//Work with stack also |
| 679 | |
| 680 | //Check invariant mass |
| 681 | Bool_t decay = kFALSE ; |
| 682 | if(fMakeInvMass) decay = CheckInvMass(icalo,mom,vertex,pl); |
| 683 | if(decay) continue ; |
| 684 | |
| 685 | //Create AOD for analysis |
| 686 | AliAODParticleCorrelation ph = AliAODParticleCorrelation(mom); |
| 687 | ph.SetLabel(calo->GetLabel(0)); |
| 688 | ph.SetPdg(AliCaloPID::kPhoton); |
| 689 | ph.SetDetector(fDetector); |
| 690 | ph.SetTag(tag); |
| 691 | if(fMakeIC){ |
| 692 | n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0; |
| 693 | GetIsolationCut()->MakeIsolationCut((TSeqCollection*)GetAODCTS(), (TSeqCollection*)pl, |
| 694 | vertex, kTRUE, &ph,icalo,-1, |
| 695 | n,nfrac,coneptsum, isolated); |
| 696 | if(isolated){ |
| 697 | //cout<<"Isolated : E "<<mom.E()<<" ; Phi"<<mom.Phi()<< " ; Eta "<<mom.Eta()<<endl; |
| 698 | AddAODParticleCorrelation(ph); |
| 699 | } |
| 4b707925 |
700 | } |
| 2c8efea8 |
701 | else //Fill all if not isolation requested |
| 702 | AddAODParticleCorrelation(ph); |
| 4b707925 |
703 | |
| 2c8efea8 |
704 | }//loop |
| 705 | |
| 706 | if(GetDebug() > 1) printf("End fill AODs "); |
| 707 | |
| f9cea31c |
708 | } |
| 709 | |
| 2c8efea8 |
710 | //__________________________________________________________________ |
| 711 | void AliAnaGammaDirect::MakeAnalysisFillHistograms() |
| 712 | { |
| 713 | //Do analysis and fill histograms |
| 714 | Int_t n = 0, nfrac = 0; |
| 715 | Bool_t isolated = kFALSE ; |
| 716 | Float_t coneptsum = 0 ; |
| 717 | //cout<<"MakeAnalysisFillHistograms"<<endl; |
| 718 | |
| 719 | //Get vertex for photon momentum calculation |
| 720 | Double_t v[]={0,0,0} ; //vertex ; |
| 721 | if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(v); |
| 722 | |
| 723 | //Loop on stored AOD photons |
| 724 | Int_t naod = GetAODBranch()->GetEntriesFast(); |
| 725 | if(GetDebug() > 0) printf("histo aod branch entries %d\n", naod); |
| 726 | for(Int_t iaod = 0; iaod < naod ; iaod++){ |
| 727 | AliAODParticleCorrelation* ph = dynamic_cast<AliAODParticleCorrelation*> (GetAODBranch()->At(iaod)); |
| 728 | Int_t pdg = ph->GetPdg(); |
| 729 | |
| 730 | //Only isolate photons in detector fDetector |
| 731 | if(pdg != AliCaloPID::kPhoton || ph->GetDetector() != fDetector) continue; |
| 732 | |
| 733 | if(fMakeSeveralIC) { |
| 734 | //Analysis of multiple IC at same time |
| 735 | MakeSeveralICAnalysis(ph,v); |
| 736 | continue; |
| 737 | } |
| 738 | else if(fReMakeIC){ |
| 739 | //In case a more strict IC is needed in the produced AOD |
| 740 | n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0; |
| 741 | GetIsolationCut()->MakeIsolationCut((TSeqCollection*)ph->GetRefIsolationConeTracks(), |
| 742 | (TSeqCollection*)ph->GetRefIsolationConeClusters(), |
| 743 | v, kFALSE, ph,-1,-1, |
| 744 | n,nfrac,coneptsum, isolated); |
| 745 | } |
| 746 | |
| 747 | //Fill histograms (normal case) |
| 748 | if(!fReMakeIC || (fReMakeIC && isolated)){ |
| 749 | |
| 750 | //printf("Isolated Gamma: pt %f, phi %f, eta %f", ph->Pt(),ph->Phi(),ph->Eta()) ; |
| 751 | |
| 752 | //Fill prompt gamma histograms |
| 753 | Float_t ptcluster = ph->Pt(); |
| 754 | Float_t phicluster = ph->Phi(); |
| 755 | Float_t etacluster = ph->Eta(); |
| 756 | |
| 757 | fhPtGamma ->Fill(ptcluster); |
| 758 | fhPhiGamma ->Fill(ptcluster,phicluster); |
| 759 | fhEtaGamma ->Fill(ptcluster,etacluster); |
| 760 | fhConeSumPt->Fill(ptcluster,coneptsum); |
| 761 | |
| 762 | if(IsDataMC()){ |
| 763 | Int_t tag =ph->GetTag(); |
| 764 | |
| 765 | if(tag == kPrompt){ |
| 766 | fhPtPrompt ->Fill(ptcluster); |
| 767 | fhPhiPrompt ->Fill(ptcluster,phicluster); |
| 768 | fhEtaPrompt ->Fill(ptcluster,etacluster); |
| 769 | } |
| 770 | else if(tag==kFragmentation) |
| 771 | { |
| 772 | fhPtFragmentation ->Fill(ptcluster); |
| 773 | fhPhiFragmentation ->Fill(ptcluster,phicluster); |
| 774 | fhEtaFragmentation ->Fill(ptcluster,etacluster); |
| 775 | } |
| 776 | else if(tag==kPi0Decay) |
| 777 | { |
| 778 | fhPtPi0Decay ->Fill(ptcluster); |
| 779 | fhPhiPi0Decay ->Fill(ptcluster,phicluster); |
| 780 | fhEtaPi0Decay ->Fill(ptcluster,etacluster); |
| 781 | } |
| 782 | else if(tag==kEtaDecay || tag==kOtherDecay) |
| 783 | { |
| 784 | fhPtOtherDecay ->Fill(ptcluster); |
| 785 | fhPhiOtherDecay ->Fill(ptcluster,phicluster); |
| 786 | fhEtaOtherDecay ->Fill(ptcluster,etacluster); |
| 787 | } |
| 788 | else if(tag==kConversion) |
| 789 | { |
| 790 | fhPtConversion ->Fill(ptcluster); |
| 791 | fhPhiConversion ->Fill(ptcluster,phicluster); |
| 792 | fhEtaConversion ->Fill(ptcluster,etacluster); |
| 793 | } |
| 794 | else{ |
| 795 | |
| 796 | fhPtUnknown ->Fill(ptcluster); |
| 797 | fhPhiUnknown ->Fill(ptcluster,phicluster); |
| 798 | fhEtaUnknown ->Fill(ptcluster,etacluster); |
| 799 | } |
| 800 | }//Histograms with MC |
| 801 | |
| 802 | } |
| 803 | }// aod loop |
| 804 | |
| 805 | } |
| 806 | |
| 807 | //____________________________________________________________________________ |
| 2a1d8a29 |
808 | void AliAnaGammaDirect::InitParameters() |
| f9cea31c |
809 | { |
| 2c8efea8 |
810 | |
| 2a1d8a29 |
811 | //Initialize the parameters of the analysis. |
| f9cea31c |
812 | |
| 2c8efea8 |
813 | fDetector = "PHOS" ; |
| 814 | fMakeIC = kTRUE; |
| 815 | fReMakeIC = kFALSE ; |
| 816 | fMakeSeveralIC = kFALSE ; |
| 817 | fMakeInvMass = kFALSE ; |
| f9cea31c |
818 | |
| 2c8efea8 |
819 | //----------- Several IC----------------- |
| 4b707925 |
820 | fNCones = 5 ; |
| 2c8efea8 |
821 | fNPtThresFrac = 6 ; |
| 4b707925 |
822 | fConeSizes[0] = 0.1; fConeSizes[1] = 0.2; fConeSizes[2] = 0.3; fConeSizes[3] = 0.4; fConeSizes[4] = 0.5; |
| 3bb2c538 |
823 | fPtThresholds[0]=0.; fPtThresholds[1]=1.; fPtThresholds[2]=2.; fPtThresholds[3]=3.; fPtThresholds[4]=4.;fPtThresholds[5]=5.; |
| 2c8efea8 |
824 | fPtFractions[0]=0.05; fPtFractions[1]=0.075; fPtFractions[2]=0.1; fPtFractions[3]=1.25; fPtFractions[4]=1.5;fPtFractions[5]=2.; |
| bdcfac30 |
825 | |
| f9cea31c |
826 | } |
| 827 | |
| 828 | //__________________________________________________________________ |
| 2c8efea8 |
829 | void AliAnaGammaDirect::MakeSeveralICAnalysis(AliAODParticleCorrelation* ph, Double_t v[3]) |
| 830 | { |
| 831 | //Isolation Cut Analysis for both methods and different pt cuts and cones |
| 832 | Float_t ptC = ph->Pt(); |
| 833 | Float_t phiC = ph->Phi(); |
| 834 | Float_t etaC = ph->Eta(); |
| 835 | |
| 836 | fhPtGamma->Fill(ptC); |
| 837 | fhPhiGamma->Fill(ptC,ph->Phi()); |
| 838 | fhEtaGamma->Fill(ptC,ph->Eta()); |
| 839 | Int_t tag =ph->GetTag(); |
| 840 | |
| 841 | if(IsDataMC()){ |
| 842 | if(tag == kPrompt){ |
| 843 | fhPtPrompt ->Fill(ptC); |
| 844 | fhPhiPrompt ->Fill(ptC,phiC); |
| 845 | fhEtaPrompt ->Fill(ptC,etaC); |
| f9cea31c |
846 | } |
| 2c8efea8 |
847 | else if(tag==kFragmentation) |
| 848 | { |
| 849 | fhPtFragmentation ->Fill(ptC); |
| 850 | fhPhiFragmentation ->Fill(ptC,phiC); |
| 851 | fhEtaFragmentation ->Fill(ptC,etaC); |
| 852 | } |
| 853 | else if(tag==kPi0Decay) |
| 854 | { |
| 855 | fhPtPi0Decay ->Fill(ptC); |
| 856 | fhPhiPi0Decay ->Fill(ptC,phiC); |
| 857 | fhEtaPi0Decay ->Fill(ptC,etaC); |
| 858 | } |
| 859 | else if(tag==kEtaDecay || tag==kOtherDecay) |
| 860 | { |
| 861 | fhPtOtherDecay ->Fill(ptC); |
| 862 | fhPhiOtherDecay ->Fill(ptC,phiC); |
| 863 | fhEtaOtherDecay ->Fill(ptC,etaC); |
| 864 | } |
| 865 | else if(tag==kConversion) |
| 866 | { |
| 867 | fhPtConversion ->Fill(ptC); |
| 868 | fhPhiConversion ->Fill(ptC,phiC); |
| 869 | fhEtaConversion ->Fill(ptC,etaC); |
| 870 | } |
| 871 | else{ |
| f9cea31c |
872 | |
| 2c8efea8 |
873 | fhPtUnknown ->Fill(ptC); |
| 874 | fhPhiUnknown ->Fill(ptC,phiC); |
| 875 | fhEtaUnknown ->Fill(ptC,etaC); |
| 876 | } |
| 877 | }//Histograms with MC |
| 878 | //Keep original setting used when filling AODs, reset at end of analysis |
| 879 | Float_t ptthresorg = GetIsolationCut()->GetPtThreshold(); |
| 880 | Float_t ptfracorg = GetIsolationCut()->GetPtFraction(); |
| 881 | Float_t rorg = GetIsolationCut()->GetConeSize(); |
| 882 | |
| 883 | Float_t coneptsum = 0 ; |
| 884 | Int_t n[10][10];//[fNCones][fNPtThresFrac]; |
| 885 | Int_t nfrac[10][10];//[fNCones][fNPtThresFrac]; |
| 886 | Bool_t isolated = kFALSE; |
| 887 | |
| 888 | for(Int_t icone = 0; icone<fNCones; icone++){ |
| 889 | GetIsolationCut()->SetConeSize(fConeSizes[icone]); |
| 890 | coneptsum = 0 ; |
| 891 | for(Int_t ipt = 0; ipt<fNPtThresFrac ;ipt++){ |
| 892 | n[icone][ipt]=0; |
| 893 | nfrac[icone][ipt]=0; |
| 894 | GetIsolationCut()->SetPtThreshold(fPtThresholds[ipt]); |
| 895 | GetIsolationCut()->MakeIsolationCut((TSeqCollection*)ph->GetRefIsolationConeTracks(), |
| 896 | (TSeqCollection*)ph->GetRefIsolationConeClusters(), |
| 897 | v, kFALSE, ph,-1,-1, |
| 898 | n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated); |
| 899 | if(n[icone][ipt] == 0) { |
| 900 | fhPtThresIsolated[icone][ipt]->Fill(ptC); |
| 901 | if(IsDataMC()){ |
| 902 | if(tag==kPrompt) fhPtThresIsolatedPrompt[icone][ipt]->Fill(ptC,coneptsum) ; |
| 903 | else if(tag==kConversion) fhPtThresIsolatedConversion[icone][ipt]->Fill(ptC,coneptsum) ; |
| 904 | else if(tag==kFragmentation) fhPtThresIsolatedFragmentation[icone][ipt]->Fill(ptC,coneptsum) ; |
| 905 | else if(tag==kPi0Decay) fhPtThresIsolatedPi0Decay[icone][ipt]->Fill(ptC,coneptsum) ; |
| 906 | else if(tag==kOtherDecay || tag==kEtaDecay) fhPtThresIsolatedOtherDecay[icone][ipt]->Fill(ptC,coneptsum) ; |
| 907 | else fhPtThresIsolatedUnknown[icone][ipt]->Fill(ptC,coneptsum) ; |
| 908 | } |
| 909 | } |
| 910 | if(nfrac[icone][ipt] == 0) { |
| 911 | fhPtFracIsolated[icone][ipt]->Fill(ptC); |
| 912 | if(IsDataMC()){ |
| 913 | if(tag==kPrompt) fhPtFracIsolatedPrompt[icone][ipt]->Fill(ptC,coneptsum) ; |
| 914 | else if(tag==kConversion) fhPtFracIsolatedConversion[icone][ipt]->Fill(ptC,coneptsum) ; |
| 915 | else if(tag==kFragmentation) fhPtFracIsolatedFragmentation[icone][ipt]->Fill(ptC,coneptsum) ; |
| 916 | else if(tag==kPi0Decay) fhPtFracIsolatedPi0Decay[icone][ipt]->Fill(ptC,coneptsum) ; |
| 917 | else if(tag==kOtherDecay || tag==kEtaDecay) fhPtFracIsolatedOtherDecay[icone][ipt]->Fill(ptC,coneptsum) ; |
| 918 | else fhPtFracIsolatedUnknown[icone][ipt]->Fill(ptC,coneptsum) ; |
| 919 | } |
| f9cea31c |
920 | } |
| 2c8efea8 |
921 | }//pt thresh loop |
| 922 | fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ; |
| 923 | if(IsDataMC()){ |
| 924 | if(tag==kPrompt) fhPtSumIsolatedPrompt[icone]->Fill(ptC,coneptsum) ; |
| 925 | else if(tag==kConversion) fhPtSumIsolatedConversion[icone]->Fill(ptC,coneptsum) ; |
| 926 | else if(tag==kFragmentation) fhPtSumIsolatedFragmentation[icone]->Fill(ptC,coneptsum) ; |
| 927 | else if(tag==kPi0Decay) fhPtSumIsolatedPi0Decay[icone]->Fill(ptC,coneptsum) ; |
| 928 | else if(tag==kOtherDecay || tag==kEtaDecay) fhPtSumIsolatedOtherDecay[icone]->Fill(ptC,coneptsum) ; |
| 929 | else fhPtSumIsolatedUnknown[icone]->Fill(ptC,coneptsum) ; |
| f9cea31c |
930 | } |
| f9cea31c |
931 | |
| 2c8efea8 |
932 | }//cone size loop |
| 933 | |
| 934 | //Reset original parameters for AOD analysis |
| 935 | GetIsolationCut()->SetPtThreshold(ptthresorg); |
| 936 | GetIsolationCut()->SetPtFraction(ptfracorg); |
| 937 | GetIsolationCut()->SetConeSize(rorg); |
| f9cea31c |
938 | |
| bdcfac30 |
939 | } |
| 940 | |
| 4b707925 |
941 | //__________________________________________________________________ |
| bdcfac30 |
942 | void AliAnaGammaDirect::Print(const Option_t * opt) const |
| 943 | { |
| 944 | |
| f9cea31c |
945 | //Print some relevant parameters set for the analysis |
| 946 | if(! opt) |
| 947 | return; |
| bdcfac30 |
948 | |
| 2c8efea8 |
949 | printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ; |
| bdcfac30 |
950 | |
| 2c8efea8 |
951 | printf("Min Gamma pT = %2.2f\n", GetMinPt()) ; |
| 952 | printf("Max Gamma pT = %3.2f\n", GetMaxPt()) ; |
| 953 | |
| 954 | // if(IsCaloPIDOn())printf("Check PID \n") ; |
| 955 | // if(IsCaloPIDRecalculationOn())printf("Recalculate PID \n") ; |
| 956 | // if(IsFidutialCutOn())printf("Check Fidutial cut \n") ; |
| 957 | printf("Make Isolation = %d \n", fMakeIC) ; |
| 958 | printf("ReMake Isolation = %d \n", fReMakeIC) ; |
| 959 | printf("Make Several Isolation = %d \n", fMakeSeveralIC) ; |
| 960 | |
| 961 | if(fMakeSeveralIC){ |
| 962 | printf("N Cone Sizes = %d\n", fNCones) ; |
| 963 | printf("Cone Sizes = \n") ; |
| bdcfac30 |
964 | for(Int_t i = 0; i < fNCones; i++) |
| 2c8efea8 |
965 | printf(" %1.2f;", fConeSizes[i]) ; |
| bdcfac30 |
966 | printf(" \n") ; |
| 2c8efea8 |
967 | |
| 968 | printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ; |
| 969 | printf(" pT thresholds = \n") ; |
| 970 | for(Int_t i = 0; i < fNPtThresFrac; i++) |
| 971 | printf(" %2.2f;", fPtThresholds[i]) ; |
| 972 | |
| 973 | printf(" \n") ; |
| 974 | |
| 975 | printf(" pT fractions = \n") ; |
| 976 | for(Int_t i = 0; i < fNPtThresFrac; i++) |
| 977 | printf(" %2.2f;", fPtFractions[i]) ; |
| 978 | |
| 979 | } |
| 980 | |
| bdcfac30 |
981 | printf(" \n") ; |
| 982 | |
| 983 | } |
| 2c8efea8 |
984 | |
| 985 | //____________________________________________________________________________ |
| 986 | Bool_t AliAnaGammaDirect::SelectCluster(AliAODCaloCluster * calo, Double_t vertex[3], TLorentzVector & mom){ |
| 987 | //Select cluster depending on its pid and acceptance selections |
| 988 | |
| 989 | //Skip matched clusters with tracks |
| 990 | if(calo->GetNTracksMatched() > 0) return kFALSE ; |
| 991 | |
| 992 | //Check PID |
| 993 | calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line |
| 994 | Int_t pdg = AliCaloPID::kPhoton; |
| 995 | if(IsCaloPIDOn()){ |
| 996 | //Get most probable PID, 2 options check PID weights (in MC this option is mandatory) |
| 997 | //or redo PID, recommended option for EMCal. |
| 998 | if(!IsCaloPIDRecalculationOn() || GetReader()->GetDataType() == AliCaloTrackReader::kMC ) |
| 999 | pdg = GetCaloPID()->GetPdg(fDetector,calo->PID(),mom.E());//PID with weights |
| 1000 | else |
| 1001 | pdg = GetCaloPID()->GetPdg(fDetector,mom,calo->GetM02(),0,0,0,0);//PID recalculated |
| 1002 | |
| 1003 | if(GetDebug() > 1) printf("PDG of identified particle %d\n",pdg); |
| 1004 | |
| 1005 | //If it does not pass pid, skip |
| 1006 | if(pdg != AliCaloPID::kPhoton) return kFALSE ; |
| 1007 | } |
| 1008 | |
| 1009 | //Check acceptance selection |
| 1010 | if(IsFidutialCutOn()){ |
| 1011 | Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,fDetector) ; |
| 1012 | if(! in ) return kFALSE ; |
| 1013 | } |
| 1014 | |
| 1015 | if(GetDebug() > 1) printf("Correlation photon selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg); |
| 1016 | |
| 1017 | return kTRUE; |
| 1018 | |
| 1019 | } |
git://git.uio.no / u / mrichter / AliRoot.git / blame