Compilation on Windows/Cygwin
[u/mrichter/AliRoot.git] / PWG4 / AliNeutralMesonSelection.cxx
CommitLineData
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$
3788def4 20 * Revision 1.3 2007/10/29 13:48:42 gustavo
21 * Corrected coding violations
22 *
3bb2c538 23 * Revision 1.2 2007/08/17 12:40:04 schutz
24 * New analysis classes by Gustavo Conesa
25 *
bdcfac30 26 * Revision 1.1.2.1 2007/07/26 10:32:09 schutz
27 * new analysis classes in the the new analysis framework
28 *
29 *
30 */
31
32//_________________________________________________________________________
33// Class that contains methods to select candidate pairs to neutral meson
3bb2c538 34// 2 main selections, invariant mass around pi0 (also any other mass),
35// apperture angle to distinguish from combinatorial.
36// There is a 3rd cut based on the gamma correlation on phi or pt.
bdcfac30 37//-- Author: Gustavo Conesa (INFN-LNF)
38
39// --- ROOT system ---
40#include <TParticle.h>
41#include <TLorentzVector.h>
42#include <TH2.h>
43#include <TList.h>
3bb2c538 44#include <TArrayD.h>
45
bdcfac30 46//---- AliRoot system ----
47#include "AliNeutralMesonSelection.h"
bdcfac30 48#include "AliLog.h"
49
50ClassImp(AliNeutralMesonSelection)
51
52
53//____________________________________________________________________________
54 AliNeutralMesonSelection::AliNeutralMesonSelection() :
55 TObject(), fSelect(0), fM(0),
56 fInvMassMaxCut(0.), fInvMassMinCut(0.),
57 fAngleMaxParam(), fMinPt(0),
58 fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.),
59 fRatioMaxCut(0), fRatioMinCut(0), fKeepNeutralMesonHistos(0),
60 fhAnglePairNoCut(0), fhAnglePairCorrelationCut(0),
61 fhAnglePairOpeningAngleCut(0), fhAnglePairAllCut(0),
62 fhInvMassPairNoCut(0), fhInvMassPairCorrelationCut(0),
63 fhInvMassPairOpeningAngleCut(0), fhInvMassPairAllCut(0)
64{
65 //Default Ctor
66
67 //Initialize parameters
68
69 // kGammaHadron and kGammaJet
70 fAngleMaxParam.Set(4) ;
71 fAngleMaxParam.Reset(0.);
72
73 //Initialize parameters
74 InitParameters();
75}
76
77//____________________________________________________________________________
78AliNeutralMesonSelection::AliNeutralMesonSelection(const AliNeutralMesonSelection & g) :
79 TObject(),
80 fSelect(g.fSelect), fM(g.fM),
81 fInvMassMaxCut(g.fInvMassMaxCut), fInvMassMinCut(g.fInvMassMinCut),
82 fAngleMaxParam(g.fAngleMaxParam), fMinPt(g.fMinPt),
83 fDeltaPhiMaxCut(g.fDeltaPhiMaxCut), fDeltaPhiMinCut(g.fDeltaPhiMinCut),
84 fRatioMaxCut(g.fRatioMaxCut), fRatioMinCut(g.fRatioMinCut),
85 fKeepNeutralMesonHistos(g.fKeepNeutralMesonHistos),
86 fhAnglePairNoCut(g. fhAnglePairNoCut),
87 fhAnglePairCorrelationCut(g. fhAnglePairCorrelationCut),
88 fhAnglePairOpeningAngleCut(g. fhAnglePairOpeningAngleCut),
89 fhAnglePairAllCut(g. fhAnglePairAllCut),
90 fhInvMassPairNoCut(g.fhInvMassPairNoCut),
91 fhInvMassPairCorrelationCut(g.fhInvMassPairCorrelationCut),
92 fhInvMassPairOpeningAngleCut(g.fhInvMassPairOpeningAngleCut),
93 fhInvMassPairAllCut(g.fhInvMassPairAllCut)
94{
95 // cpy ctor
96}
97
98//_________________________________________________________________________
99AliNeutralMesonSelection & AliNeutralMesonSelection::operator = (const AliNeutralMesonSelection & source)
100{
101 // assignment operator
102
103 if(this == &source)return *this;
104 ((TObject *)this)->operator=(source);
105
106 fSelect = source.fSelect ;
107 fM = source.fM ;
108 fInvMassMaxCut = source.fInvMassMaxCut ;
109 fInvMassMinCut = source.fInvMassMinCut ;
110 fAngleMaxParam = source.fAngleMaxParam ;
111 fMinPt = source.fMinPt ;
112 fDeltaPhiMaxCut = source.fDeltaPhiMaxCut ;
113 fDeltaPhiMinCut = source.fDeltaPhiMinCut ;
114 fRatioMaxCut = source.fRatioMaxCut ;
115 fRatioMinCut = source.fRatioMinCut ;
116 fKeepNeutralMesonHistos = source.fKeepNeutralMesonHistos;
117
118 fhAnglePairNoCut = source. fhAnglePairNoCut ;
119 fhAnglePairCorrelationCut = source. fhAnglePairCorrelationCut ;
120 fhAnglePairOpeningAngleCut = source. fhAnglePairOpeningAngleCut ;
121 fhAnglePairAllCut = source. fhAnglePairAllCut ;
122 fhInvMassPairNoCut = source.fhInvMassPairNoCut ;
123 fhInvMassPairCorrelationCut = source.fhInvMassPairCorrelationCut ;
124 fhInvMassPairOpeningAngleCut = source.fhInvMassPairOpeningAngleCut ;
125 fhInvMassPairAllCut = source.fhInvMassPairAllCut ;
126
127 return *this;
128
129}
130
131//____________________________________________________________________________
132AliNeutralMesonSelection::~AliNeutralMesonSelection()
133{
3788def4 134 // Remove all pointers except analysis output pointers.
135
bdcfac30 136}
137
138
139
140//________________________________________________________________________
141TList * 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 //____________________________________________________________________________
217void 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//__________________________________________________________________________-
3bb2c538 243Bool_t AliNeutralMesonSelection::IsAngleInWindow(const Float_t angle,const Float_t e) const {
bdcfac30 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//____________________________________________________________________________
3bb2c538 262Bool_t AliNeutralMesonSelection::CutPtPhi(Double_t ptg, Double_t phig, Double_t pt, Double_t phi) const
bdcfac30 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//____________________________________________________________________________
280Bool_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//__________________________________________________________________
329void 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}