]>
Commit | Line | Data |
---|---|---|
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 | } |