]>
Commit | Line | Data |
---|---|---|
bdcfac30 | 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$ | |
20 | * Revision 1.1.2.1 2007/07/26 10:32:09 schutz | |
21 | * new analysis classes in the the new analysis framework | |
22 | * | |
23 | * | |
24 | */ | |
25 | ||
26 | //_________________________________________________________________________ | |
27 | // Class that contains methods to select candidate pairs to neutral meson | |
28 | //-- Author: Gustavo Conesa (INFN-LNF) | |
29 | ||
30 | // --- ROOT system --- | |
31 | #include <TParticle.h> | |
32 | #include <TLorentzVector.h> | |
33 | #include <TH2.h> | |
34 | #include <TList.h> | |
35 | ||
36 | //---- AliRoot system ---- | |
37 | #include "AliNeutralMesonSelection.h" | |
38 | #include "Riostream.h" | |
39 | #include "AliLog.h" | |
40 | ||
41 | ClassImp(AliNeutralMesonSelection) | |
42 | ||
43 | ||
44 | //____________________________________________________________________________ | |
45 | AliNeutralMesonSelection::AliNeutralMesonSelection() : | |
46 | TObject(), fSelect(0), fM(0), | |
47 | fInvMassMaxCut(0.), fInvMassMinCut(0.), | |
48 | fAngleMaxParam(), fMinPt(0), | |
49 | fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), | |
50 | fRatioMaxCut(0), fRatioMinCut(0), fKeepNeutralMesonHistos(0), | |
51 | fhAnglePairNoCut(0), fhAnglePairCorrelationCut(0), | |
52 | fhAnglePairOpeningAngleCut(0), fhAnglePairAllCut(0), | |
53 | fhInvMassPairNoCut(0), fhInvMassPairCorrelationCut(0), | |
54 | fhInvMassPairOpeningAngleCut(0), fhInvMassPairAllCut(0) | |
55 | { | |
56 | //Default Ctor | |
57 | ||
58 | //Initialize parameters | |
59 | ||
60 | // kGammaHadron and kGammaJet | |
61 | fAngleMaxParam.Set(4) ; | |
62 | fAngleMaxParam.Reset(0.); | |
63 | ||
64 | //Initialize parameters | |
65 | InitParameters(); | |
66 | } | |
67 | ||
68 | //____________________________________________________________________________ | |
69 | AliNeutralMesonSelection::AliNeutralMesonSelection(const AliNeutralMesonSelection & g) : | |
70 | TObject(), | |
71 | fSelect(g.fSelect), fM(g.fM), | |
72 | fInvMassMaxCut(g.fInvMassMaxCut), fInvMassMinCut(g.fInvMassMinCut), | |
73 | fAngleMaxParam(g.fAngleMaxParam), fMinPt(g.fMinPt), | |
74 | fDeltaPhiMaxCut(g.fDeltaPhiMaxCut), fDeltaPhiMinCut(g.fDeltaPhiMinCut), | |
75 | fRatioMaxCut(g.fRatioMaxCut), fRatioMinCut(g.fRatioMinCut), | |
76 | fKeepNeutralMesonHistos(g.fKeepNeutralMesonHistos), | |
77 | fhAnglePairNoCut(g. fhAnglePairNoCut), | |
78 | fhAnglePairCorrelationCut(g. fhAnglePairCorrelationCut), | |
79 | fhAnglePairOpeningAngleCut(g. fhAnglePairOpeningAngleCut), | |
80 | fhAnglePairAllCut(g. fhAnglePairAllCut), | |
81 | fhInvMassPairNoCut(g.fhInvMassPairNoCut), | |
82 | fhInvMassPairCorrelationCut(g.fhInvMassPairCorrelationCut), | |
83 | fhInvMassPairOpeningAngleCut(g.fhInvMassPairOpeningAngleCut), | |
84 | fhInvMassPairAllCut(g.fhInvMassPairAllCut) | |
85 | { | |
86 | // cpy ctor | |
87 | } | |
88 | ||
89 | //_________________________________________________________________________ | |
90 | AliNeutralMesonSelection & AliNeutralMesonSelection::operator = (const AliNeutralMesonSelection & source) | |
91 | { | |
92 | // assignment operator | |
93 | ||
94 | if(this == &source)return *this; | |
95 | ((TObject *)this)->operator=(source); | |
96 | ||
97 | fSelect = source.fSelect ; | |
98 | fM = source.fM ; | |
99 | fInvMassMaxCut = source.fInvMassMaxCut ; | |
100 | fInvMassMinCut = source.fInvMassMinCut ; | |
101 | fAngleMaxParam = source.fAngleMaxParam ; | |
102 | fMinPt = source.fMinPt ; | |
103 | fDeltaPhiMaxCut = source.fDeltaPhiMaxCut ; | |
104 | fDeltaPhiMinCut = source.fDeltaPhiMinCut ; | |
105 | fRatioMaxCut = source.fRatioMaxCut ; | |
106 | fRatioMinCut = source.fRatioMinCut ; | |
107 | fKeepNeutralMesonHistos = source.fKeepNeutralMesonHistos; | |
108 | ||
109 | fhAnglePairNoCut = source. fhAnglePairNoCut ; | |
110 | fhAnglePairCorrelationCut = source. fhAnglePairCorrelationCut ; | |
111 | fhAnglePairOpeningAngleCut = source. fhAnglePairOpeningAngleCut ; | |
112 | fhAnglePairAllCut = source. fhAnglePairAllCut ; | |
113 | fhInvMassPairNoCut = source.fhInvMassPairNoCut ; | |
114 | fhInvMassPairCorrelationCut = source.fhInvMassPairCorrelationCut ; | |
115 | fhInvMassPairOpeningAngleCut = source.fhInvMassPairOpeningAngleCut ; | |
116 | fhInvMassPairAllCut = source.fhInvMassPairAllCut ; | |
117 | ||
118 | return *this; | |
119 | ||
120 | } | |
121 | ||
122 | //____________________________________________________________________________ | |
123 | AliNeutralMesonSelection::~AliNeutralMesonSelection() | |
124 | { | |
125 | // Remove all pointers | |
126 | ||
127 | delete fhAnglePairNoCut ; | |
128 | delete fhAnglePairCorrelationCut ; | |
129 | delete fhAnglePairOpeningAngleCut ; | |
130 | delete fhAnglePairAllCut ; | |
131 | delete fhInvMassPairNoCut ; | |
132 | delete fhInvMassPairCorrelationCut ; | |
133 | delete fhInvMassPairOpeningAngleCut ; | |
134 | delete fhInvMassPairAllCut ; | |
135 | ||
136 | } | |
137 | ||
138 | ||
139 | ||
140 | //________________________________________________________________________ | |
141 | TList * AliNeutralMesonSelection::GetCreateOutputObjects() | |
142 | { | |
143 | ||
144 | // Create histograms to be saved in output file and | |
145 | // store them in outputContainer | |
146 | TList * outputContainer = new TList() ; | |
147 | outputContainer->SetName("MesonDecayHistos") ; | |
148 | ||
149 | fhAnglePairNoCut = new TH2F | |
150 | ("AnglePairNoCut", | |
151 | "Angle between all #gamma pair vs E_{#pi^{0}}",200,0,50,200,0,0.2); | |
152 | fhAnglePairNoCut->SetYTitle("Angle (rad)"); | |
153 | fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)"); | |
154 | ||
155 | fhAnglePairOpeningAngleCut = new TH2F | |
156 | ("AnglePairOpeningAngleCut", | |
157 | "Angle between all #gamma pair (opening angle + azimuth cut) vs E_{#pi^{0}}" | |
158 | ,200,0,50,200,0,0.2); | |
159 | fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)"); | |
160 | fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV)"); | |
161 | ||
162 | fhAnglePairAllCut = new TH2F | |
163 | ("AnglePairAllCut", | |
164 | "Angle between all #gamma pair (opening angle + inv mass cut+azimuth) vs E_{#pi^{0}}" | |
165 | ,200,0,50,200,0,0.2); | |
166 | fhAnglePairAllCut->SetYTitle("Angle (rad)"); | |
167 | fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)"); | |
168 | ||
169 | // | |
170 | fhInvMassPairNoCut = new TH2F | |
171 | ("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}", | |
172 | 120,0,120,360,0,0.5); | |
173 | fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})"); | |
174 | fhInvMassPairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)"); | |
175 | ||
176 | fhInvMassPairOpeningAngleCut = new TH2F | |
177 | ("InvMassPairOpeningAngleCut", | |
178 | "Invariant Mass of #gamma pair (angle cut) vs E_{#pi^{0}}", | |
179 | 120,0,120,360,0,0.5); | |
180 | fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})"); | |
181 | fhInvMassPairOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)"); | |
182 | ||
183 | fhInvMassPairAllCut = new TH2F | |
184 | ("InvMassPairAllCut", | |
185 | "Invariant Mass of #gamma pair (opening angle+invmass cut) vs E_{#pi^{0}}", | |
186 | 120,0,120,360,0,0.5); | |
187 | fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})"); | |
188 | fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)"); | |
189 | ||
190 | fhAnglePairCorrelationCut = new TH2F | |
191 | ("AnglePairCorrelationCut", | |
192 | "Angle between correlated #gamma pair vs E_{#pi^{0}}",200,0,50,200,0,0.2); | |
193 | fhAnglePairCorrelationCut->SetYTitle("Angle (rad)"); | |
194 | fhAnglePairCorrelationCut->SetXTitle("E_{ #pi^{0}} (GeV)"); | |
195 | ||
196 | fhInvMassPairCorrelationCut = new TH2F | |
197 | ("InvMassPairCorrelationCut","Invariant Mass of correlated #gamma pair vs E_{#pi^{0}}", | |
198 | 120,0,120,360,0,0.5); | |
199 | fhInvMassPairCorrelationCut->SetYTitle("Invariant Mass (GeV/c^{2})"); | |
200 | fhInvMassPairCorrelationCut->SetXTitle("E_{ #pi^{0}} (GeV)"); | |
201 | ||
202 | outputContainer->Add(fhAnglePairNoCut) ; | |
203 | outputContainer->Add(fhAnglePairOpeningAngleCut) ; | |
204 | outputContainer->Add(fhAnglePairAllCut) ; | |
205 | ||
206 | outputContainer->Add(fhInvMassPairNoCut) ; | |
207 | outputContainer->Add(fhInvMassPairOpeningAngleCut) ; | |
208 | outputContainer->Add(fhInvMassPairAllCut) ; | |
209 | ||
210 | outputContainer->Add(fhAnglePairCorrelationCut) ; | |
211 | outputContainer->Add(fhInvMassPairCorrelationCut) ; | |
212 | ||
213 | return outputContainer; | |
214 | } | |
215 | ||
216 | //____________________________________________________________________________ | |
217 | void AliNeutralMesonSelection::InitParameters() | |
218 | { | |
219 | ||
220 | //Initialize the parameters of the analysis. | |
221 | fKeepNeutralMesonHistos = kTRUE ; | |
222 | ||
223 | //-------------kHadron, kJetLeadCone----------------- | |
224 | fAngleMaxParam.Set(4) ; | |
225 | fAngleMaxParam.AddAt(0.4,0);//={0.4,-0.25,0.025,-2e-4}; | |
226 | fAngleMaxParam.AddAt(-0.25,1) ; | |
227 | fAngleMaxParam.AddAt(0.025,2) ; | |
228 | fAngleMaxParam.AddAt(-2e-4,3) ; | |
229 | ||
230 | fInvMassMaxCut = 0.16 ; | |
231 | fInvMassMinCut = 0.11 ; | |
232 | ||
233 | fM = 0.1349766;//neutralMeson mass | |
234 | ||
235 | fMinPt = 0. ; | |
236 | fDeltaPhiMaxCut = 4.5; | |
237 | fDeltaPhiMinCut = 1.5 ; | |
238 | fRatioMaxCut = 1.0 ; | |
239 | fRatioMinCut = 0.1 ; | |
240 | } | |
241 | ||
242 | //__________________________________________________________________________- | |
243 | Bool_t AliNeutralMesonSelection::IsAngleInWindow(const Float_t angle,const Float_t e) { | |
244 | //Check if the opening angle of the candidate pairs is inside | |
245 | //our selection windowd | |
246 | ||
247 | Bool_t result = kFALSE; | |
248 | Double_t max = fAngleMaxParam.At(0)*TMath::Exp(fAngleMaxParam.At(1)*e) | |
249 | +fAngleMaxParam.At(2)+fAngleMaxParam.At(3)*e; | |
250 | Double_t arg = (e*e-2*fM*fM)/(e*e); | |
251 | Double_t min = 100. ; | |
252 | if(arg>0.) | |
253 | min = TMath::ACos(arg); | |
254 | ||
255 | if((angle<max)&&(angle>=min)) | |
256 | result = kTRUE; | |
257 | ||
258 | return result; | |
259 | } | |
260 | ||
261 | //____________________________________________________________________________ | |
262 | Bool_t AliNeutralMesonSelection::CutPtPhi(Double_t ptg, Double_t phig, Double_t pt, Double_t phi) | |
263 | { | |
264 | //Select pair if delta | |
265 | Bool_t cut = kFALSE ; | |
266 | ||
267 | if(fSelect == kNoSelectPhiPt) cut = kTRUE ; | |
268 | else if((phig-phi) > fDeltaPhiMinCut && ((phig-phi) < fDeltaPhiMaxCut)){ | |
269 | //Cut on pt | |
270 | if((fSelect == kSelectPhiPtRatio && ptg > 0. && pt/ptg > fRatioMinCut && pt/ptg < fRatioMaxCut) || | |
271 | (fSelect == kSelectPhiMinPt && pt > fMinPt) ) cut = kTRUE ; | |
272 | } | |
273 | else cut = kFALSE ; | |
274 | ||
275 | return cut ; | |
276 | ||
277 | } | |
278 | ||
279 | //____________________________________________________________________________ | |
280 | Bool_t AliNeutralMesonSelection::SelectPair(TParticle * pGamma, TLorentzVector gammai, TLorentzVector gammaj) | |
281 | { | |
282 | ||
283 | //Search for the neutral pion within selection cuts | |
284 | ||
285 | Double_t ptg = pGamma->Pt(); | |
286 | Double_t phig = pGamma->Phi() ; | |
287 | Bool_t goodpair = kFALSE ; | |
288 | ||
289 | Double_t pt = (gammai+gammaj).Pt(); | |
290 | Double_t phi = (gammai+gammaj).Phi(); | |
291 | if(phi < 0) | |
292 | phi+=TMath::TwoPi(); | |
293 | Double_t invmass = (gammai+gammaj).M(); | |
294 | Double_t angle = gammaj.Angle(gammai.Vect()); | |
295 | Double_t e = (gammai+gammaj).E(); | |
296 | ||
297 | //Fill histograms with no cuts applied. | |
298 | fhAnglePairNoCut->Fill(e,angle); | |
299 | fhInvMassPairNoCut->Fill(e,invmass); | |
300 | ||
301 | //Cut on phig-phi meson | |
302 | if(CutPtPhi(ptg, phig, pt, phi)){ | |
303 | ||
304 | fhAnglePairCorrelationCut ->Fill(e,angle); | |
305 | fhInvMassPairCorrelationCut->Fill(e,invmass); | |
306 | ||
307 | //Cut on the aperture of the pair | |
308 | if(IsAngleInWindow(angle,e)){ | |
309 | fhAnglePairOpeningAngleCut ->Fill(e,angle); | |
310 | fhInvMassPairOpeningAngleCut->Fill(e,invmass); | |
311 | AliDebug(2,Form("Angle cut: pt %f, phi %f",pt,phi)); | |
312 | ||
313 | //Cut on the invariant mass of the pair | |
314 | if((invmass>fInvMassMinCut) && (invmass<fInvMassMaxCut)){ | |
315 | fhInvMassPairAllCut ->Fill(e,invmass); | |
316 | fhAnglePairAllCut ->Fill(e,angle); | |
317 | goodpair = kTRUE; | |
318 | AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi)); | |
319 | }//(invmass>0.125) && (invmass<0.145) | |
320 | }//Opening angle cut | |
321 | } // cut on pt and phi | |
322 | ||
323 | ||
324 | return goodpair; | |
325 | ||
326 | } | |
327 | ||
328 | //__________________________________________________________________ | |
329 | void AliNeutralMesonSelection::Print(const Option_t * opt) const | |
330 | { | |
331 | ||
332 | //Print some relevant parameters set for the analysis | |
333 | if(! opt) | |
334 | return; | |
335 | ||
336 | Info("Print", "%s %s", GetName(), GetTitle() ) ; | |
337 | ||
338 | printf("mass : %f \n", fM ); | |
339 | printf("Invariant mass limits : %f < m < %f \n", fInvMassMinCut , fInvMassMinCut ); | |
340 | printf("Angle selection param: \n"); | |
341 | printf("p0 : %f", fAngleMaxParam.At(0)); | |
342 | printf("p1 : %f", fAngleMaxParam.At(1)); | |
343 | printf("p2 : %f", fAngleMaxParam.At(2)); | |
344 | printf("p3 : %f", fAngleMaxParam.At(3)); | |
345 | ||
346 | printf("pT meson > %f\n", fMinPt) ; | |
347 | printf("Phi gamma-meson < %f\n", fDeltaPhiMaxCut) ; | |
348 | printf("Phi gamma-meson > %f\n", fDeltaPhiMinCut) ; | |
349 | printf("pT meson / pT Gamma < %f\n", fRatioMaxCut) ; | |
350 | printf("pT meson / pT Gamma > %f\n", fRatioMinCut) ; | |
351 | printf("Keep Neutral Meson Histos = %d\n",fKeepNeutralMesonHistos); | |
352 | ||
353 | } |