1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //_________________________________________________________________________
17 // Class that contains methods to select candidate pairs to neutral meson
18 // 2 main selections, invariant mass around pi0 (also any other mass),
19 // apperture angle to distinguish from combinatorial.
20 //-- Author: Gustavo Conesa (INFN-LNF)
22 // --- ROOT system ---
23 #include <TLorentzVector.h>
28 //---- AliRoot system ----
29 #include "AliNeutralMesonSelection.h"
31 ClassImp(AliNeutralMesonSelection)
34 //______________________________________________________
35 AliNeutralMesonSelection::AliNeutralMesonSelection() :
36 TObject(), fAsymmetryCut(1),
37 fUseAsymmetryCut(0), fM(0),
38 fInvMassMaxCut(0.), fInvMassMinCut(0.),
39 fLeftBandMinCut(0.), fLeftBandMaxCut(0.),
40 fRightBandMinCut(0.), fRightBandMaxCut(0.),
41 fAngleMaxParam(), fUseAngleCut(0),
42 fKeepNeutralMesonHistos(0),
43 fParticle(""), fDecayBit(0),
45 fhAnglePairNoCut(0), fhAnglePairOpeningAngleCut(0),
46 fhAnglePairAsymmetryCut(0), fhAnglePairAllCut(0),
47 fhInvMassPairNoCut(0), fhInvMassPairOpeningAngleCut(0),
48 fhInvMassPairAsymmetryCut(0), fhInvMassPairAllCut(0),
49 fhAsymmetryNoCut(0), fhAsymmetryOpeningAngleCut(0),
51 // histogram ranges and bins
52 fHistoNEBins(0), fHistoEMax(0.), fHistoEMin(0.),
53 fHistoNAngleBins(0), fHistoAngleMax(0.), fHistoAngleMin(0.),
54 fHistoNIMBins(0), fHistoIMMax(0.), fHistoIMMin(0.)
58 //Initialize parameters
62 //_________________________________________________________
63 TList * AliNeutralMesonSelection::GetCreateOutputObjects()
65 // Create histograms to be saved in output file and
66 // store them in outputContainer of the analysis class that calls this class.
68 TList * outputContainer = new TList() ;
69 outputContainer->SetName("MesonDecayHistos") ;
71 if(fKeepNeutralMesonHistos){
73 outputContainer->SetOwner(kFALSE);
75 fhAnglePairNoCut = new TH2F
77 "Angle between all #gamma pair vs E_{#pi^{0}}",fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax);
78 fhAnglePairNoCut->SetYTitle("Angle (rad)");
79 fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
81 fhAsymmetryNoCut = new TH2F
82 ("AsymmetryNoCut","Asymmetry of all #gamma pair vs E_{#pi^{0}}",
83 fHistoNEBins,fHistoEMin,fHistoEMax,100,0,1);
84 fhAsymmetryNoCut->SetYTitle("Asymmetry");
85 fhAsymmetryNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
87 fhInvMassPairNoCut = new TH2F
88 ("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}",
89 fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax);
90 fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})");
91 fhInvMassPairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
93 outputContainer->Add(fhAnglePairNoCut) ;
94 outputContainer->Add(fhAsymmetryNoCut) ;
95 outputContainer->Add(fhInvMassPairNoCut) ;
98 fhAnglePairOpeningAngleCut = new TH2F
99 ("AnglePairOpeningAngleCut",
100 "Angle between all #gamma pair (opening angle) vs E_{#pi^{0}}"
101 ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax);
102 fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)");
103 fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV)");
105 fhAsymmetryOpeningAngleCut = new TH2F
106 ("AsymmetryOpeningAngleCut",
107 "Asymmetry of #gamma pair (angle cut) vs E_{#pi^{0}}",
108 fHistoNEBins,fHistoEMin,fHistoEMax,100,0,1);
109 fhAsymmetryOpeningAngleCut->SetYTitle("Asymmetry");
110 fhAsymmetryOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)");
112 fhInvMassPairOpeningAngleCut = new TH2F
113 ("InvMassPairOpeningAngleCut",
114 "Invariant Mass of #gamma pair (angle cut) vs E_{#pi^{0}}",
115 fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax);
116 fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})");
117 fhInvMassPairOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)");
119 outputContainer->Add(fhAnglePairOpeningAngleCut) ;
120 outputContainer->Add(fhAsymmetryOpeningAngleCut) ;
121 outputContainer->Add(fhInvMassPairOpeningAngleCut) ;
124 if(fUseAsymmetryCut) {
125 fhAnglePairAsymmetryCut = new TH2F
126 ("AnglePairAsymmetryCut",
127 "Angle between all #gamma pair (opening angle + asymetry cut) vs E_{#pi^{0}}"
128 ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax);
129 fhAnglePairAsymmetryCut->SetYTitle("Angle (rad)");
130 fhAnglePairAsymmetryCut->SetXTitle("E_{ #pi^{0}} (GeV)");
132 fhInvMassPairAsymmetryCut = new TH2F
133 ("InvMassPairAsymmetryCut",
134 "Invariant Mass of #gamma pair (opening angle + asymmetry) vs E_{#pi^{0}}",
135 fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax);
136 fhInvMassPairAsymmetryCut->SetYTitle("Invariant Mass (GeV/c^{2})");
137 fhInvMassPairAsymmetryCut->SetXTitle("E_{#pi^{0}}(GeV)");
139 outputContainer->Add(fhAnglePairAsymmetryCut) ;
140 outputContainer->Add(fhInvMassPairAsymmetryCut) ;
144 fhAnglePairAllCut = new TH2F
146 "Angle between all #gamma pair (opening angle + asymmetry + inv mass cut) vs E_{#pi^{0}}"
147 ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNAngleBins,fHistoAngleMin,fHistoAngleMax);
148 fhAnglePairAllCut->SetYTitle("Angle (rad)");
149 fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)");
151 fhInvMassPairAllCut = new TH2F
152 ("InvMassPairAllCut",
153 "Invariant Mass of #gamma pair (opening angle + asymmetry + invmass cut) vs E_{#pi^{0}}",
154 fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNIMBins,fHistoIMMin,fHistoIMMax);
155 fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
156 fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)");
158 fhAsymmetryAllCut = new TH2F
160 "Asymmetry of #gamma pair (opening angle+invmass cut) vs E_{#pi^{0}}",
161 fHistoNEBins,fHistoEMin,fHistoEMax,100,0,1);
162 fhAsymmetryAllCut->SetYTitle("Asymmetry");
163 fhAsymmetryAllCut->SetXTitle("E_{#pi^{0}}(GeV)");
165 outputContainer->Add(fhAnglePairAllCut) ;
166 outputContainer->Add(fhAsymmetryAllCut) ;
167 outputContainer->Add(fhInvMassPairAllCut) ;
171 return outputContainer;
175 //_____________________________________________
176 void AliNeutralMesonSelection::InitParameters()
179 //Initialize the parameters of the analysis.
180 fAngleMaxParam.Set(4) ;
181 fAngleMaxParam.Reset(0.);
185 //Histogrammes settings
190 fHistoNAngleBins = 200 ;
191 fHistoAngleMax = 0.5 ;
192 fHistoAngleMin = 0. ;
196 //______________________________________________________________________________
197 Bool_t AliNeutralMesonSelection::IsAngleInWindow(Float_t angle, Float_t e) const
200 // Check if the opening angle of the candidate pairs is inside
201 // our selection window
202 // Attention, only valid for Pi0, if needed for Eta need to revise max angle function or change parameters
204 Double_t max = fAngleMaxParam.At(0)*TMath::Exp(fAngleMaxParam.At(1)*e)
205 +fAngleMaxParam.At(2)+fAngleMaxParam.At(3)*e;
206 Double_t arg = (e*e-2*fM*fM)/(e*e);
207 Double_t min = 100. ;
209 min = TMath::ACos(arg)+fShiftMinAngle[0]+fShiftMinAngle[1]*e;
211 if((angle<max)&&(angle>=min)) return kTRUE ;
216 //_________________________________________________________________
217 Bool_t AliNeutralMesonSelection::SelectPair(TLorentzVector gammai,
218 TLorentzVector gammaj,
222 //Search for the neutral pion within selection cuts
224 // Double_t pt = (gammai+gammaj).Pt();
225 Double_t phi = (gammai+gammaj).Phi();
229 Double_t invmass = (gammai+gammaj).M();
230 Double_t angle = gammaj.Angle(gammai.Vect());
231 Double_t e = (gammai+gammaj).E();
232 Double_t asy = TMath::Abs((gammai-gammaj).E())/(gammai+gammaj).E();
234 //Fill histograms with no cuts applied.
235 if(fKeepNeutralMesonHistos)
237 fhAnglePairNoCut ->Fill(e,angle);
238 fhInvMassPairNoCut->Fill(e,invmass);
239 fhAsymmetryNoCut ->Fill(e,asy);
242 //Cut on the aperture of the pair
245 if(IsAngleInWindow(angle,e))
247 if(fKeepNeutralMesonHistos )
249 fhAnglePairOpeningAngleCut ->Fill(e,angle);
250 fhInvMassPairOpeningAngleCut->Fill(e,invmass);
251 fhAsymmetryOpeningAngleCut ->Fill(e,asy);
253 //AliDebug(2,Form("Angle cut: pt %f, phi %f",pt,phi));
254 } else return kFALSE;
260 if(fAsymmetryCut > asy)
262 if(fKeepNeutralMesonHistos)
264 fhInvMassPairAsymmetryCut->Fill(e,invmass);
265 fhAnglePairAsymmetryCut ->Fill(e,angle);
267 } else return kFALSE;
271 //Cut on the invariant mass of the pair
273 Float_t invmassmaxcut = fInvMassMaxCut;
274 Float_t invmassRightBandMinCut = fRightBandMinCut;
275 Float_t invmassRightBandMixCut = fRightBandMaxCut;
277 if(calo=="EMCAL" && e > 6.)
278 { // for EMCAL, pi0s, mass depends strongly with energy for e > 6, loose max cut
280 invmassmaxcut = (fInvMassMaxCutParam[0]+fInvMassMaxCut)+fInvMassMaxCutParam[1]*e+fInvMassMaxCutParam[2]*e*e;
281 invmassRightBandMinCut = (fInvMassMaxCutParam[0]+fRightBandMinCut)+fInvMassMaxCutParam[1]*e+fInvMassMaxCutParam[2]*e*e;
282 invmassRightBandMixCut = (fInvMassMaxCutParam[0]+fRightBandMaxCut)+fInvMassMaxCutParam[1]*e+fInvMassMaxCutParam[2]*e*e;
284 //printf("e %f, max cut %f, p00 %f,p0 %f,p1 %f,p2 %f\n",
285 // e,invmassmaxcut,fInvMassMaxCut,fInvMassMaxCutParam[0],fInvMassMaxCutParam[1],fInvMassMaxCutParam[2]);
288 // normal case, invariant mass selection around pi0/eta peak
289 if( !fParticle.Contains("SideBand") )
291 if( invmass > fInvMassMinCut && invmass < invmassmaxcut )
293 if(fKeepNeutralMesonHistos)
295 fhInvMassPairAllCut->Fill(e,invmass);
296 fhAnglePairAllCut ->Fill(e,angle);
297 fhAsymmetryAllCut ->Fill(e,asy);
300 //AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi));
303 }//(invmass>0.125) && (invmass<0.145)
310 else // select a band around pi0/eta
312 if((invmass > fLeftBandMinCut && invmass < fLeftBandMaxCut ) ||
313 (invmass > invmassRightBandMinCut && invmass < invmassRightBandMixCut))
315 if(fKeepNeutralMesonHistos)
317 fhInvMassPairAllCut->Fill(e,invmass);
318 fhAnglePairAllCut ->Fill(e,angle);
319 fhAsymmetryAllCut ->Fill(e,asy);
322 //AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi));
325 }//(invmass>0.125) && (invmass<0.145)
334 //_______________________________________________________________
335 void AliNeutralMesonSelection::SetParticle(TString particleName)
337 // Set some default parameters for selection of pi0 or eta
339 fParticle = particleName ;
342 if(particleName.Contains("Pi0"))
344 fHistoNIMBins = 150 ;
349 fInvMassMaxCut = 0.16 ; // GeV
350 fInvMassMinCut = 0.11 ; // GeV
352 fLeftBandMinCut = 0.05 ; // GeV
353 fLeftBandMaxCut = 0.09 ; // GeV
354 fRightBandMinCut = 0.160 ; // GeV
355 fRightBandMaxCut = 0.165 ; // GeV
357 fInvMassMaxCutParam[0] = 0.0 ;
358 fInvMassMaxCutParam[1] =-7.e-5 ;
359 fInvMassMaxCutParam[2] = 8.e-5 ;
361 fShiftMinAngle[0] =-0.03 ;
362 fShiftMinAngle[1] = 0.0025;
364 fAngleMaxParam.AddAt( 0.8, 0) ;
365 fAngleMaxParam.AddAt(-1, 1) ;
366 fAngleMaxParam.AddAt( 0.09, 2) ; //for pi0 shift, for eta maybe 0.09
367 fAngleMaxParam.AddAt(-2.e-3,3) ;
370 if(particleName.Contains("Side")) fDecayBit = kPi0Side;
373 else if(particleName.Contains("Eta"))
375 fHistoNIMBins = 200 ; // GeV
376 fHistoIMMax = 0.75 ; // GeV
377 fHistoIMMin = 0.35 ; // GeV
380 fInvMassMaxCut = 0.590 ; // GeV
381 fInvMassMinCut = 0.510 ; // GeV
383 fLeftBandMinCut = 0.450 ; // GeV
384 fLeftBandMaxCut = 0.500 ; // GeV
385 fRightBandMinCut = 0.600 ; // GeV
386 fRightBandMaxCut = 0.650 ; // GeV
388 fInvMassMaxCutParam[0] = 0.00 ;
389 fInvMassMaxCutParam[1] = 0.00 ;
390 fInvMassMaxCutParam[2] = 0.00 ;
392 fShiftMinAngle[0] =-0.03 ;
393 fShiftMinAngle[0] = 0.00 ;
395 fAngleMaxParam.AddAt( 0.80, 0) ; // Same as pi0
396 fAngleMaxParam.AddAt(-0.25, 1) ; // Same as pi0
397 fAngleMaxParam.AddAt( 0.12, 2) ; // Shifted with respect to pi0
398 fAngleMaxParam.AddAt(-5.e-4, 3) ; // Same as pi0
401 if(particleName.Contains("Side")) fDecayBit = kEtaSide;
404 printf("AliAnaNeutralMesonSelection::SetParticle(%s) *** Particle NOT defined (Pi0 or Eta), Pi0 settings by default *** \n",particleName.Data());
408 //______________________________________________________________
409 void AliNeutralMesonSelection::Print(const Option_t * opt) const
412 //Print some relevant parameters set for the analysis
416 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
418 printf("Particle %s, bit %d, mass : %f \n", fParticle.Data(), fDecayBit, fM );
419 printf("Invariant mass limits : %f < m < %f \n", fInvMassMinCut , fInvMassMinCut );
421 printf("Use asymmetry cut? : %d ; A < %f \n", fUseAngleCut, fAsymmetryCut );
423 printf("Use angle cut? : %d \n", fUseAngleCut );
426 printf("Angle selection param: \n");
427 printf("p0 : %f\n", fAngleMaxParam.At(0));
428 printf("p1 : %f\n", fAngleMaxParam.At(1));
429 printf("p2 : %f\n", fAngleMaxParam.At(2));
430 printf("p3 : %f\n", fAngleMaxParam.At(3));
431 printf("Min angle shift : p0 = %1.3f, p1 = %1.3f,\n", fShiftMinAngle[0],fShiftMinAngle[1]);
434 printf("Keep Neutral Meson Histos = %d\n",fKeepNeutralMesonHistos);
436 if(fKeepNeutralMesonHistos)
438 printf("Histograms: %3.1f < E < %3.1f, Nbin = %d\n", fHistoEMin, fHistoEMax, fHistoNEBins);
439 printf("Histograms: %3.1f < angle < %3.1f, Nbin = %d\n", fHistoAngleMin, fHistoAngleMax, fHistoNAngleBins);
440 printf("Histograms: %3.1f < IM < %3.1f, Nbin = %d\n", fHistoIMMin, fHistoIMMax, fHistoNIMBins);