1 - Rename AliAnaMaker and AliAnaBaseClass to AliAnaPartCorrMaker and AliAnaPartCorrB...
[u/mrichter/AliRoot.git] / PWG4 / AliNeutralMesonSelection.cxx
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 //_________________________________________________________________________
18 // Class that contains methods to select candidate pairs to neutral meson 
19 // 2 main selections, invariant mass around pi0 (also any other mass),
20 // apperture angle to distinguish from combinatorial.
21 //-- Author: Gustavo Conesa (INFN-LNF)
22
23 // --- ROOT system ---
24 #include <TLorentzVector.h>
25 #include <TH2.h>
26 #include <TList.h>
27 #include <TArrayD.h>
28  
29 //---- AliRoot system ----
30 #include "AliNeutralMesonSelection.h" 
31 #include "AliLog.h"
32
33 ClassImp(AliNeutralMesonSelection)
34   
35   
36 //____________________________________________________________________________
37   AliNeutralMesonSelection::AliNeutralMesonSelection() : 
38     TObject(), fM(0),
39     fInvMassMaxCut(0.), fInvMassMinCut(0.),
40     fAngleMaxParam(),  fKeepNeutralMesonHistos(0), 
41     fhAnglePairNoCut(0), fhAnglePairOpeningAngleCut(0), 
42     fhAnglePairAllCut(0), 
43     fhInvMassPairNoCut(0), fhInvMassPairOpeningAngleCut(0), 
44     fhInvMassPairAllCut(0) 
45 {
46   //Default Ctor
47   
48   //Initialize parameters
49   
50   // kGammaHadron and kGammaJet 
51   fAngleMaxParam.Set(4) ;
52   fAngleMaxParam.Reset(0.);
53   
54   //Initialize parameters
55   InitParameters();
56 }
57
58 //____________________________________________________________________________
59 AliNeutralMesonSelection::AliNeutralMesonSelection(const AliNeutralMesonSelection & g) :   
60   TObject(), fM(g.fM),
61   fInvMassMaxCut(g.fInvMassMaxCut), fInvMassMinCut(g.fInvMassMinCut),
62   fAngleMaxParam(g.fAngleMaxParam),
63   fKeepNeutralMesonHistos(g.fKeepNeutralMesonHistos),
64   fhAnglePairNoCut(g. fhAnglePairNoCut), 
65   fhAnglePairOpeningAngleCut(g. fhAnglePairOpeningAngleCut), 
66   fhAnglePairAllCut(g. fhAnglePairAllCut), 
67   fhInvMassPairNoCut(g.fhInvMassPairNoCut),  
68   fhInvMassPairOpeningAngleCut(g.fhInvMassPairOpeningAngleCut), 
69   fhInvMassPairAllCut(g.fhInvMassPairAllCut)
70 {
71   // cpy ctor
72 }
73
74 //_________________________________________________________________________
75 AliNeutralMesonSelection & AliNeutralMesonSelection::operator = (const AliNeutralMesonSelection & source)
76 {
77   // assignment operator
78   
79   if(this == &source)return *this;
80   ((TObject *)this)->operator=(source);
81
82   fM = source.fM ;
83   fInvMassMaxCut = source.fInvMassMaxCut ; 
84   fInvMassMinCut = source.fInvMassMinCut ;
85   fAngleMaxParam = source.fAngleMaxParam ;
86   fKeepNeutralMesonHistos = source.fKeepNeutralMesonHistos;
87  
88   fhAnglePairNoCut = source. fhAnglePairNoCut ; 
89   fhAnglePairOpeningAngleCut = source. fhAnglePairOpeningAngleCut ; 
90   fhAnglePairAllCut = source. fhAnglePairAllCut ; 
91   fhInvMassPairNoCut = source.fhInvMassPairNoCut ; 
92   fhInvMassPairOpeningAngleCut = source.fhInvMassPairOpeningAngleCut ; 
93   fhInvMassPairAllCut = source.fhInvMassPairAllCut ; 
94   
95   return *this;
96   
97 }
98
99 //____________________________________________________________________________
100 AliNeutralMesonSelection::~AliNeutralMesonSelection() 
101 {
102   //dtor
103
104   if(!fKeepNeutralMesonHistos){
105     //Histograms initialized and filled but not passed to output container
106     //delete here, I am not sure this is correct
107     
108     if(fhAnglePairNoCut) delete fhAnglePairNoCut;
109     if(fhAnglePairOpeningAngleCut) delete fhAnglePairOpeningAngleCut; 
110     if(fhAnglePairAllCut) delete fhAnglePairAllCut;
111     if(fhInvMassPairNoCut) delete fhInvMassPairNoCut;
112     if(fhInvMassPairOpeningAngleCut) delete fhInvMassPairOpeningAngleCut;
113     if(fhInvMassPairAllCut) delete fhInvMassPairAllCut; 
114
115   }
116   
117 }
118 //________________________________________________________________________
119 TList *  AliNeutralMesonSelection::GetCreateOutputObjects()
120 {  
121
122   // Create histograms to be saved in output file and 
123   // store them in outputContainer of the analysis class that calls this class.
124
125   TList * outputContainer = new TList() ; 
126   outputContainer->SetName("MesonDecayHistos") ; 
127   
128   fhAnglePairNoCut  = new TH2F
129     ("AnglePairNoCut",
130      "Angle between all #gamma pair vs E_{#pi^{0}}",200,0,50,200,0,0.2); 
131   fhAnglePairNoCut->SetYTitle("Angle (rad)");
132   fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
133   
134   fhAnglePairOpeningAngleCut  = new TH2F
135     ("AnglePairOpeningAngleCut",
136      "Angle between all #gamma pair (opening angle + azimuth cut) vs E_{#pi^{0}}"
137      ,200,0,50,200,0,0.2); 
138   fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)");
139   fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV)");
140   
141   fhAnglePairAllCut  = new TH2F
142     ("AnglePairAllCut",
143      "Angle between all #gamma pair (opening angle + inv mass cut+azimuth) vs E_{#pi^{0}}"
144      ,200,0,50,200,0,0.2); 
145   fhAnglePairAllCut->SetYTitle("Angle (rad)");
146   fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV)");    
147   
148   //
149   fhInvMassPairNoCut  = new TH2F
150     ("InvMassPairNoCut","Invariant Mass of all #gamma pair vs E_{#pi^{0}}",
151      120,0,120,360,0,0.5); 
152   fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})");
153   fhInvMassPairNoCut->SetXTitle("E_{ #pi^{0}} (GeV)");
154   
155   fhInvMassPairOpeningAngleCut  = new TH2F
156     ("InvMassPairOpeningAngleCut",
157      "Invariant Mass of #gamma pair (angle cut) vs E_{#pi^{0}}",
158      120,0,120,360,0,0.5); 
159   fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})");
160   fhInvMassPairOpeningAngleCut->SetXTitle(" E_{#pi^{0}}(GeV)");
161   
162   fhInvMassPairAllCut  = new TH2F
163     ("InvMassPairAllCut",
164      "Invariant Mass of #gamma pair (opening angle+invmass cut) vs E_{#pi^{0}}",
165      120,0,120,360,0,0.5); 
166   fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
167   fhInvMassPairAllCut->SetXTitle("E_{#pi^{0}}(GeV)");
168
169   outputContainer->Add(fhAnglePairNoCut) ; 
170   outputContainer->Add(fhAnglePairOpeningAngleCut) ;
171   outputContainer->Add(fhAnglePairAllCut) ; 
172   
173   outputContainer->Add(fhInvMassPairNoCut) ; 
174   outputContainer->Add(fhInvMassPairOpeningAngleCut) ; 
175   outputContainer->Add(fhInvMassPairAllCut) ; 
176   
177   return outputContainer;
178 }
179
180  //____________________________________________________________________________
181 void AliNeutralMesonSelection::InitParameters()
182 {
183  
184   //Initialize the parameters of the analysis.
185   fKeepNeutralMesonHistos = kFALSE ;
186
187   //-------------kHadron, kJetLeadCone-----------------
188   fAngleMaxParam.Set(4) ;
189   fAngleMaxParam.AddAt(0.4,0);//={0.4,-0.25,0.025,-2e-4};
190   fAngleMaxParam.AddAt(-0.25,1) ;
191   fAngleMaxParam.AddAt(0.025,2) ;
192   fAngleMaxParam.AddAt(-2e-4,3) ;
193
194   fInvMassMaxCut  = 0.16 ;
195   fInvMassMinCut  = 0.11 ;
196
197   fM = 0.1349766;//neutralMeson mass
198 }
199
200 //__________________________________________________________________________-
201 Bool_t AliNeutralMesonSelection::IsAngleInWindow(const Float_t angle,const Float_t e) const {
202   //Check if the opening angle of the candidate pairs is inside 
203   //our selection windowd
204
205   Bool_t result = kFALSE;
206   Double_t max =  fAngleMaxParam.At(0)*TMath::Exp(fAngleMaxParam.At(1)*e)
207     +fAngleMaxParam.At(2)+fAngleMaxParam.At(3)*e;
208   Double_t arg = (e*e-2*fM*fM)/(e*e);
209   Double_t min = 100. ;
210   if(arg>0.)
211     min = TMath::ACos(arg);
212
213   if((angle<max)&&(angle>=min))
214     result = kTRUE;
215  
216   return result;
217 }
218
219 //____________________________________________________________________________
220 Bool_t  AliNeutralMesonSelection::SelectPair(TLorentzVector gammai, TLorentzVector gammaj)  
221 {  
222   
223   //Search for the neutral pion within selection cuts
224   Bool_t goodpair = kFALSE ;
225   
226   Double_t pt  = (gammai+gammaj).Pt();
227   Double_t phi = (gammai+gammaj).Phi();
228   if(phi < 0)
229     phi+=TMath::TwoPi();
230   Double_t invmass = (gammai+gammaj).M();
231   Double_t angle   = gammaj.Angle(gammai.Vect());
232   Double_t e       = (gammai+gammaj).E();
233   
234   //Fill histograms with no cuts applied.
235   fhAnglePairNoCut->Fill(e,angle);
236   fhInvMassPairNoCut->Fill(e,invmass);
237     
238   //Cut on the aperture of the pair
239   if(IsAngleInWindow(angle,e)){
240     fhAnglePairOpeningAngleCut     ->Fill(e,angle);
241     fhInvMassPairOpeningAngleCut->Fill(e,invmass);
242     AliDebug(2,Form("Angle cut: pt %f, phi %f",pt,phi));
243     
244     //Cut on the invariant mass of the pair
245     if((invmass>fInvMassMinCut) && (invmass<fInvMassMaxCut)){ 
246       fhInvMassPairAllCut  ->Fill(e,invmass);
247       fhAnglePairAllCut       ->Fill(e,angle);
248       goodpair = kTRUE;
249       AliDebug(2,Form("IM cut: pt %f, phi %f",pt,phi));
250     }//(invmass>0.125) && (invmass<0.145)
251   }//Opening angle cut
252   
253   return goodpair; 
254   
255 }
256
257 //__________________________________________________________________
258 void AliNeutralMesonSelection::Print(const Option_t * opt) const
259 {
260
261   //Print some relevant parameters set for the analysis
262   if(! opt)
263     return;
264   
265   Info("Print", "%s %s", GetName(), GetTitle() ) ;
266
267   printf("mass : %f  \n", fM );
268   printf("Invariant mass limits : %f < m < %f \n", fInvMassMinCut , fInvMassMinCut );
269   printf("Angle selection param: \n");
270   printf("p0 :     %f", fAngleMaxParam.At(0));
271   printf("p1 :     %f", fAngleMaxParam.At(1));
272   printf("p2 :     %f", fAngleMaxParam.At(2));
273   printf("p3 :     %f", fAngleMaxParam.At(3));
274
275   printf("Keep Neutral Meson Histos = %d\n",fKeepNeutralMesonHistos);
276
277