]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/hfe/AliHFEV0cuts.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Warning fix (Laurent)
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEV0cuts.cxx
CommitLineData
c04c80e6 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// * 20/04/2010 *
16// Class for optimising and applying V0 cuts to obtain clean V0 samples
17// Compatible with ESDs only
18//
19// Authors:
20// Matus Kalisky <matus.kalisky@cern.ch>
21//
22
23#include "TDatabasePDG.h"
24
25#include "AliESDtrack.h"
26#include "AliMCEvent.h"
27#include "AliESDv0.h"
28#include "AliKFParticle.h"
29#include "AliKFVertex.h"
3a72645a 30#include "AliLog.h"
c04c80e6 31
32#include "AliHFEcollection.h"
33
34#include "AliHFEV0cuts.h"
35
36ClassImp(AliHFEV0cuts)
37
38//________________________________________________________________
39AliHFEV0cuts::AliHFEV0cuts():
40 fQA(NULL)
3a72645a 41 , fQAmc(NULL)
c04c80e6 42 , fMCEvent(NULL)
43 , fInputEvent(NULL)
44 , fPrimaryVertex(NULL)
3a72645a 45 , fCurrentV0id(0)
46 , fPdaughterPDG(0)
47 , fNdaughterPDG(0)
c04c80e6 48{
49
50 //
51 // Default constructor
52 //
53
54
55}
56//________________________________________________________________
57AliHFEV0cuts::~AliHFEV0cuts()
58{
59 //
60 // destructor
61 //
62 if (fQA) delete fQA;
3a72645a 63 if (fQAmc) delete fQAmc;
c04c80e6 64}
65
66//________________________________________________________________
67AliHFEV0cuts::AliHFEV0cuts(const AliHFEV0cuts &ref):
68 TObject(ref)
69 , fQA(NULL)
3a72645a 70 , fQAmc(NULL)
c04c80e6 71 , fMCEvent(NULL)
72 , fInputEvent(NULL)
73 , fPrimaryVertex(NULL)
3a72645a 74 , fCurrentV0id(0)
75 , fPdaughterPDG(0)
76 , fNdaughterPDG(0)
c04c80e6 77{
78 //
79 // Copy constructor
80 //
81 ref.Copy(*this);
82}
83//________________________________________________________________
84AliHFEV0cuts &AliHFEV0cuts::operator=(const AliHFEV0cuts &ref){
85 //
86 // Assignment operator
87 //
88 if(this != &ref)
89 ref.Copy(*this);
90 return *this;
91}
92//________________________________________________________________
93void AliHFEV0cuts::Copy(TObject &ref) const{
94 //
95 // Copy function
96 //
97 AliHFEV0cuts &target = dynamic_cast<AliHFEV0cuts &>(ref);
98
99 if(fQA) target.fQA = dynamic_cast<AliHFEcollection *>(fQA->Clone());
100
3a72645a 101 if(fQAmc) target.fQAmc = dynamic_cast<AliHFEcollection *>(fQAmc->Clone());
102
c04c80e6 103 if(target.fMCEvent) delete target.fMCEvent;
104 target.fMCEvent = new AliMCEvent;
105
106 if(target.fPrimaryVertex) delete target.fPrimaryVertex;
107 target.fPrimaryVertex = new AliKFVertex;
108
109 TObject::Copy(ref);
110
111}
112//___________________________________________________________________
113void AliHFEV0cuts::Init(const char* name){
114 //
115 // initialize the output objects and create histograms
116 //
117
118 //
119 // all the "h_cut_XXX" histograms hare cut value distributions:
120 // [0] for all candidates
121 // [1] jus before the cut on given variable was applied, but after all the previous cuts
122 //
123
124 fQA = new AliHFEcollection("fQA", name);
125
3a72645a 126 fQAmc = new AliHFEcollection("fQAmc", name);
c04c80e6 127
128 // common for all V0s
129 fQA->CreateTH2Fvector1(2, "h_all_AP", "armenteros plot for all V0 candidates", 200, -1, 1, 200, 0, 0.25);
130
131 // gammas
132 fQA->CreateTH1Fvector1(2, "h_cut_Gamma_CosPoint", "Gamma Cosine pointing angle; cos point. angle; counts", 100, 0, 0.1);
133 fQA->CreateTH1Fvector1(2, "h_cut_Gamma_DCA", "DCA between the gamma daughters; dca (cm); counts", 100, 0, 2);
3a72645a 134 fQA->CreateTH1Fvector1(2, "h_cut_Gamma_VtxR_old", "*old* Radius of the gamma conversion vertex; r (cm); counts", 1000, 0, 100);
c04c80e6 135 fQA->CreateTH1Fvector1(2, "h_cut_Gamma_VtxR", "Radius of the gamma conversion vertex; r (cm); counts", 1000, 0, 100);
136 fQA->CreateTH1Fvector1(2, "h_cut_Gamma_OA", "opening angle of the gamma products; opening angle (rad); counts", 100, 0, 1);
137 fQA->CreateTH1Fvector1(2, "h_cut_Gamma_PP", "gamma psi pair angle; psi pairangle (rad); counts", 100, 0, 2);
3a72645a 138 fQA->CreateTH1Fvector1(2, "h_cut_Gamma_Chi2", "gamma Chi2/NDF; Chi2/NDF; counts", 100, 0, 50);
139 fQA->CreateTH1Fvector1(7, "h_Gamma_Mass", "Invariant mass of gammas; mass (GeV/c^{2}); counts", 100, 0, 0.2);
c04c80e6 140
141
142 // kaons
143 fQA->CreateTH1Fvector1(2, "h_cut_K0_CosPoint", "K0 Cosine pointing angle; cos point. angle; counts", 100, 0, 0.1);
144 fQA->CreateTH1Fvector1(2, "h_cut_K0_DCA", "DCA between the K0 daughters; dca (cm); counts", 100, 0, 2);
145 fQA->CreateTH1Fvector1(2, "h_cut_K0_VtxR", "Radius of the K0 decay vertex; r (cm); counts", 1000, 0, 100);
3a72645a 146 fQA->CreateTH1Fvector1(2, "h_cut_K0_Chi2", "K0 Chi2/NDF; Chi2/NDF; counts", 100, 0, 50);
147 fQA->CreateTH1Fvector1(5, "h_K0_Mass", "Invariant mass of K0; mass (GeV/c^{2}); counts", 125, 0.45, 0.55);
c04c80e6 148
149 // lambda
150 fQA->CreateTH1Fvector1(2, "h_cut_L_CosPoint", "L Cosine pointing angle; cos point. angle; counts", 100, 0, 0.1);
151 fQA->CreateTH1Fvector1(2, "h_cut_L_DCA", "DCA between the L daughters; dca (cm); counts", 100, 0, 2);
152 fQA->CreateTH1Fvector1(2, "h_cut_L_VtxR", "Radius of the L decay vertex; r (cm); counts", 1000, 0, 100);
3a72645a 153 fQA->CreateTH1Fvector1(2, "h_cut_L_Chi2", "L Chi2/NDF; Chi2/NDF; counts", 100, 0, 50);
154 fQA->CreateTH1Fvector1(5, "h_L_Mass", "Invariant mass of L; mass (GeV/c^{2}); counts", 60, 1.1, 1.13);
155 fQA->CreateTH1Fvector1(5, "h_AL_Mass", "Invariant mass of anti L; mass (GeV/c^{2}); counts", 60, 1.1, 1.13);
c04c80e6 156
157 fQA->CreateTH2F("h_L_checks", "Lambda candidate check[0] -v- check[1]; check[0]; check[1]", 5, -0.75, 1.75, 6, -0.75, 1.75 );
158
159 // electrons
160 fQA->CreateTH1Fvector1(9, "h_Electron_P", "Momenta of conversion electrons -cuts-; P (GeV/c); counts", 50, 0.1, 20, 0);
161
162 // K0 pions
163 fQA->CreateTH1Fvector1(8, "h_PionK0_P", "Momenta of K0 pions -cuts-; P (GeV/c) counts;", 50, 0.1, 20, 0);
164
165 // L pions
3a72645a 166 fQA->CreateTH1Fvector1(9, "h_PionL_P", "Momenta of L pions -cuts-; P (GeV/c) counts;", 50, 0.1, 20, 0);
c04c80e6 167
168 // L protons
169 fQA->CreateTH1Fvector1(9, "h_ProtonL_P", "Momenta of L protons -cuts-; P (GeV/c) counts;", 50, 0.1, 20, 0);
170
171 // single track cuts
172 fQA->CreateTH1F("h_ST_NclsTPC", "Number of TPC clusters", 161, -1, 160);
173 fQA->CreateTH1F("h_ST_TPCrefit", "TPC refit", 2, -0.5, 1.5);
174 fQA->CreateTH1F("h_ST_chi2TPCcls", "chi2 per TPC cluster", 100, 0, 10);
175 fQA->CreateTH1F("h_ST_TPCclsR", "TPC cluster ratio", 120, -0.1, 1.1);
176 fQA->CreateTH1F("h_ST_kinks", "kinks", 2, -0.5, 1.5);
177 fQA->CreateTH1F("h_ST_pt", "track pt", 100, 0.1, 20, 0);
178 fQA->CreateTH1F("h_ST_eta", "track eta", 100, -1.5, 1.5);
179
180 //
181 // possibly new cuts
182 //
183
184 // Gamma
185 fQA->CreateTH2Fvector1(2, "h_cut_Gamma_OAvP", "open. ang. of the Gamma daughters versus Gamma mom; Gamma p (GeV/c); opening angle (pions) (rad)", 100, 0.1, 10, 200, 0., 0.2);
186 // K0
187 fQA->CreateTH2Fvector1(2, "h_cut_K0_OAvP", "open. ang. of the K0 daughters versus K0 momentum; K0 p (GeV/c); opening angle (pions) (rad)", 100, 0.1, 10, 100, 0, 3.5);
188 // Lambda
189 fQA->CreateTH2Fvector1(2, "h_cut_L_OAvP", "open. ang. of the L daughters versus L momentum; Lambda p (GeV/c); openeing angle pion-proton (rad)", 100, 0.1, 10, 100, 0, 3.5);
190 fQA->CreateTH2Fvector1(2, "h_cut_L_rdp_v_mp", "relative L daughter mom -v- mother mom; L mom (GeV/c); relative daughter mom p2/p1", 100, 0.1, 10, 100, 0, 1);
c04c80e6 191
192
193 // THnSparse histograms
194
195 // THnSparse for the K0 mass
196 // to be looked at after merging run by run
197 // axes: mass, pt, theta, phi
198 {
199 Int_t nBin[4] = {100, 10, 10, 18};
200 Double_t nMin[4] = {0.45, 0.1, 0., 0.};
201 Double_t nMax[4] = {0.55, 10., TMath::Pi(), 2*TMath::Pi()};
202 TString htitle = "K0 sparse; mass (GeV/c^{2}); p_{T} (GeV/c); theta (rad); phi(rad)";
203 fQA->CreateTHnSparse("hK0", htitle, 4, nBin, nMin, nMax);
204 fQA->BinLogAxis("hK0", 1);
205 }
206
3a72645a 207
208 //
209 // MC plots for checking and tuning the V0 cuts
210 //
211
212 const char *v0[4] = {"G", "K", "L"}; // to keep the names short
213 // number of V0s left after each cut step - for signal and background - within given mass window
214 for(Int_t i=0; i<3; ++i){
215 fQAmc->CreateTH1F(Form("h_%s_cuts_S", v0[i]), Form("h_%s_cuts_S", v0[i]), 10, -0.5, 9.5);
216 fQAmc->CreateTH1F(Form("h_%s_cuts_B", v0[i]), Form("h_%s_cuts_B", v0[i]), 10, -0.5, 9.5);
217 }
218
219 //
220 // cut distributions for signal and background
221 //
222
223 const Float_t pMin = 0.1;
224 const Float_t pMax = 10.;
225 const Int_t pN = 12;
226
227
228 // gamma signal
229 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_CosPoint_S", "S - Gamma Cosine pointing angle; mom (GeV/c); cos point. angle", pN, pMin, pMax, 50, 0, 0.1, 0);
230 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_DCA_S", "S - DCA between the gamma daughters; mom (GeV/c); dca (cm)", pN, pMin, pMax, 50, 0, 2, 0);
231 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_VtxR_S", "S - Radius of the gamma conversion vertex; mom (GeV/c); r (cm)", pN, pMin, pMax, 100, 0, 100, 0);
232 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_OA_S", "S - opening angle of the gamma products; mom (GeV/c); opening angle (rad)", pN, pMin, pMax, 50, 0, 0.3, 0);
233 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_PP_S", "S - gamma psi pair angle; mom (GeV/c); psi pairangle (rad)", pN, pMin, pMax, 50, 0, 0.5, 0);
234 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_Chi2_S", "S - gamma Chi2/NDF; mom (GeV/c); Chi2/NDF", pN, pMin, pMax, 50, 0, 100, 0);
235
236 fQAmc->CreateTH1Fvector1(8, "h_Gamma_Mass_S", "S - Invariant mass of gammas; mass (GeV/c^{2}); counts", 100, 0, 0.2);
237 // gamma background
238 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_CosPoint_B", "B - Gamma Cosine pointing angle; mom (GeV/c); cos point. angle", pN, pMin, pMax, 50, 0, 0.1, 0);
239 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_DCA_B", "B - DCA between the gamma daughters; mom (GeV/c); dca (cm)", pN, pMin, pMax, 50, 0, 2, 0);
240 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_VtxR_B", "B - Radius of the gamma conversion vertex; mom (GeV/c); r (cm)", pN, pMin, pMax, 100, 0, 100, 0);
241 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_OA_B", "B - opening angle of the gamma products; mom (GeV/c); opening angle (rad)", pN, pMin, pMax, 50, 0, 0.3, 0);
242 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_PP_B", "B - gamma psi pair angle; mom (GeV/c); psi pairangle (rad)", pN, pMin, pMax, 50, 0, 0.5, 0);
243 fQAmc->CreateTH2Fvector1(2, "h_cut_Gamma_Chi2_B", "B - gamma Chi2/NDF; mom (GeV/c); Chi2/NDF", pN, pMin, pMax, 50, 0, 100, 0);
244
245 fQAmc->CreateTH1Fvector1(8, "h_Gamma_Mass_B", "B - Invariant mass of gammas; mass (GeV/c^{2}); counts", 100, 0, 0.2);
246
247 // kaons signal
248 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_CosPoint_S", "S - K0 Cosine pointing angle; mom (GeV/c); cos point. angle", pN, pMin, pMax, 50, 0, 0.1, 0);
249 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_DCA_S", "S - DCA between the K0 daughters; mom (GeV/c); dca (cm)", pN, pMin, pMax, 50, 0, 2, 0);
250 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_VtxR_S", "S - Radius of the K0 decay vertex; mom (GeV/c); r (cm)", pN, pMin, pMax, 50, 0, 100, 0);
251 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_Chi2_S", "S - K0 Chi2/NDF; mom (GeV/c); Chi2/NDF", pN, pMin, pMax, 50, 0, 100, 0);
252 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_OA_S", "S - opening angle of the K0 pions; mom (GeV/c); opening angle (rad)", pN, pMin, pMax, 100, 0, 1, 0);
253
254 fQAmc->CreateTH1Fvector1(5, "h_K0_Mass_S", "S - Invariant mass of K0; mass (GeV/c^{2}); counts", 125, 0.45, 0.55);
255 // kaons background
256 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_CosPoint_B", "B - K0 Cosine pointing angle; mom (GeV/c); cos point. angle", pN, pMin, pMax, 50, 0, 0.1, 0);
257 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_DCA_B", "B - DCA between the K0 daughters; mom (GeV/c); dca (cm)", pN, pMin, pMax, 50, 0, 2, 0);
258 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_VtxR_B", "B - Radius of the K0 decay vertex; mom (GeV/c); r (cm)", pN, pMin, pMax, 50, 0, 100, 0);
259 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_Chi2_B", "B - K0 Chi2/NDF; mom (GeV/c); Chi2/NDF", pN, pMin, pMax, 50, 0, 100, 0);
260 fQAmc->CreateTH2Fvector1(2, "h_cut_K0_OA_B", "B - opening angle of the K0 pions; mom (GeV/c); opening angle (rad)", pN, pMin, pMax, 100, 0, 1, 0);
261
262 fQAmc->CreateTH1Fvector1(5, "h_K0_Mass_B", "B - Invariant mass of K0; mass (GeV/c^{2}); counts", 125, 0.45, 0.55);
263
264 // lambda signal
265 fQAmc->CreateTH2Fvector1(2, "h_cut_L_CosPoint_S", "S - L Cosine pointing angle; mom (GeV/c); cos point. angle", pN, pMin, pMax, 50, 0, 0.1, 0);
266 fQAmc->CreateTH2Fvector1(2, "h_cut_L_DCA_S", "S - DCA between the L daughters; mom (GeV/c); dca (cm)", pN, pMin, pMax, 50, 0, 2, 0);
267 fQAmc->CreateTH2Fvector1(2, "h_cut_L_VtxR_S", "S - Radius of the L decay vertex; mom (GeV/c); r (cm)", pN, pMin, pMax, 50, 0, 100, 0);
268 fQAmc->CreateTH2Fvector1(2, "h_cut_L_Chi2_S", "S - L Chi2/NDF; mom (GeV/c); Chi2/NDF", pN, pMin, pMax, 50, 0, 100, 0);
269 fQAmc->CreateTH2Fvector1(2, "h_cut_L_OA_S", "S - opening angle of the L p-p; mom (GeV/c); opening angle (rad)", pN, pMin, pMax, 100, 0, 1, 0);
270
271 fQAmc->CreateTH1Fvector1(5, "h_L_Mass_S", "S - Invariant mass of L; mass (GeV/c^{2}); counts", 60, 1.1, 1.13);
272 fQAmc->CreateTH1Fvector1(5, "h_AL_Mass_S", "S - Invariant mass of anti L; mass (GeV/c^{2}); counts", 60, 1.1, 1.13);
273 // lambda background
274 fQAmc->CreateTH2Fvector1(2, "h_cut_L_CosPoint_B", "B - L Cosine pointing angle; mom (GeV/c); cos point. angle", pN, pMin, pMax, 50, 0, 0.1, 0);
275 fQAmc->CreateTH2Fvector1(2, "h_cut_L_DCA_B", "B - DCA between the L daughters; mom (GeV/c); dca (cm)", pN, pMin, pMax, 50, 0, 2, 0);
276 fQAmc->CreateTH2Fvector1(2, "h_cut_L_VtxR_B", "B - Radius of the L decay vertex; mom (GeV/c); r (cm)", pN, pMin, pMax, 50, 0, 100, 0);
277 fQAmc->CreateTH2Fvector1(2, "h_cut_L_Chi2_B", "B - L Chi2/NDF; mom (GeV/c); Chi2/NDF", pN, pMin, pMax, 50, 0, 100, 0);
278 fQAmc->CreateTH2Fvector1(2, "h_cut_L_OA_B", "B - opening angle of the L p-p; mom (GeV/c); opening angle (rad)", pN, pMin, pMax, 100, 0, 1, 0);
279
280 fQAmc->CreateTH2Fvector1(2, "h_cut_L_rdp_v_mp_S", "S - relative L daughter mom -v- mother mom; L mom (GeV/c); relative daughter mom p2/p1", 100, 0.1, 10, 100, 0, 1);
281 fQAmc->CreateTH2Fvector1(2, "h_cut_L_rdp_v_mp_B", "B - relative L daughter mom -v- mother mom; L mom (GeV/c); relative daughter mom p2/p1", 100, 0.1, 10, 100, 0, 1);
282 fQAmc->CreateTH1Fvector1(5, "h_LAL_Mass_B", "B - Invariant mass of anti L; mass (GeV/c^{2}); counts", 60, 1.1, 1.13);
283
284
285 // MC tagged daughter track momentum distribution after each cut step
286// fQAmc->CreateTH1Fvector1(10, "h_electron_p_S", "h_electron_p_S", 20, 0.1, 20, 0);
287// fQAmc->CreateTH1Fvector1(10, "h_K0pion_p_S", "h_K0pion_p_S", 20, 0.1, 20, 0);
288// fQAmc->CreateTH1Fvector1(10, "h_Lpion_p_S", "h_Lpion_p_S", 20, 0.1, 20, 0);
289// fQAmc->CreateTH1Fvector1(10, "h_proton_p_S", "h_proton_p_S", 20, 0.1, 20, 0);
290
291 // V0 momnetum distribution of MC tagged signal and backglound after all cuts
292 fQAmc->CreateTH1F("h_gamma_p_S", "true gammas after all cuts", 20, 0.1, 10, 0);
293 fQAmc->CreateTH1F("h_gamma_p_B", "true gamma BG after all cuts", 20, 0.1, 10, 0);
294 fQAmc->CreateTH1F("h_K0_p_S", "true K0s after all cuts", 20, 0.1, 10, 0);
295 fQAmc->CreateTH1F("h_K0_p_B", "true K0 BG after all cuts", 20, 0.1, 10, 0);
296 fQAmc->CreateTH1F("h_lambda_p_S", "MC true lambdas after all cuts", 20, 0.1, 10, 0);
297 fQAmc->CreateTH1F("h_lambda_p_B", "MC true lambda BG after all cuts", 20, 0.1, 10, 0);
298 fQAmc->CreateTH1F("h_alambda_p_S", "MC true anti-lambdas after all cuts", 20, 0.1, 10, 0);
299 fQAmc->CreateTH1F("h_alambda_p_B", "MC true anti-lambda BG after all cuts", 20, 0.1, 10, 0);
300
301 // invariant mass ditributions for the V0 for different hypoteses (gamma, K0, L, AL)
302 fQAmc->CreateTH1F("h_Mass_gamma_as_K0","h_Mass_gamma_as_K0", 200, 0, 2);
303 fQAmc->CreateTH1F("h_Mass_gamma_as_L","h_Mass_gamma_as_L", 200, 0, 2);
304 fQAmc->CreateTH1F("h_Mass_K0_as_G", "h_Mass_K0_as_gamma", 200, 0, 2);
305 fQAmc->CreateTH1F("h_Mass_K0_as_L", "h_Mass_K0_as_Lambda", 200, 0, 2);
306 fQAmc->CreateTH1F("h_Mass_L_as_G", "h_Mass_L_as_gamma", 200, 0, 2);
307 fQAmc->CreateTH1F("h_Mass_L_as_K0", "h_Mass_L_as_K0", 200, 0, 2);
308
309 // Invariant mass distribution of MC tagged signal for diffrent momenta
310 fQAmc->CreateTH2F("h_gamma_MvP_S", "mc tagged gammas - signal; p (GeV/c); m (GeV/c^{2})", 12, 0.1, 20, 100, 0., 0.1, 0);
311 fQAmc->CreateTH2F("h_K0_MvP_S", "mc tagged K0s - signal; p (GeV/c); m (GeV/c^{2})", 12, 0.1, 20, 100, 0.45, 0.55, 0);
312 fQAmc->CreateTH2F("h_lambda_MvP_S", "mc tagged Lambdas - signal; p (GeV/c); m (GeV/c^{2})", 12, 0.1, 20, 100, 1.08, 1.14, 0);
313
314 // electrons
315 fQAmc->CreateTH1Fvector1(8, "h_Electron_P_S", "MC-S momenta of conversion electrons -cuts-; P (GeV/c); counts", 20, 0.1, 20, 0);
316 fQAmc->CreateTH1Fvector1(8, "h_Electron_P_B", "MC-B momenta of conversion electrons -cuts-; P (GeV/c); counts", 20, 0.1, 20, 0);
317
318 // K0 pions
319 fQAmc->CreateTH1Fvector1(7, "h_PionK0_P_S", "MC-S momenta of K0 pions -cuts-; P (GeV/c) counts;", 20, 0.1, 20, 0);
320 fQAmc->CreateTH1Fvector1(7, "h_PionK0_P_B", "MC-B momenta of K0 pions -cuts-; P (GeV/c) counts;", 20, 0.1, 20, 0);
321
322 // L pions
323 fQAmc->CreateTH1Fvector1(8, "h_PionL_P_S", "MC-S momenta of L pions -cuts-; P (GeV/c) counts;", 20, 0.1, 50, 0);
324 fQAmc->CreateTH1Fvector1(8, "h_PionL_P_B", "MC-B momenta of L pions -cuts-; P (GeV/c) counts;", 20, 0.1, 50, 0);
325
326 // L protons
327 fQAmc->CreateTH1Fvector1(8, "h_ProtonL_P_S", "MC-S momenta of L protons -cuts-; P (GeV/c) counts;", 20, 0.1, 20, 0);
328 fQAmc->CreateTH1Fvector1(8, "h_ProtonL_P_B", "MC-B momenta of L protons -cuts-; P (GeV/c) counts;", 20, 0.1, 20, 0);
329
330
331
332 // cut efficiencies
c04c80e6 333}
334//________________________________________________________________
335Bool_t AliHFEV0cuts::TrackCutsCommon(AliESDtrack* track){
336 //
337 // singe track cuts commom for all particle candidates
338 //
339
340 if(!track) return kFALSE;
3a72645a 341
c04c80e6 342
343 // status word
344 ULong_t status = track->GetStatus();
345
346
347 // No. of TPC clusters
348 fQA->Fill("h_ST_NclsTPC", track->GetTPCNcls());
3a72645a 349 if(track->GetTPCNcls() < 80) return kFALSE; //
c04c80e6 350
351 // TPC refit
352 if((status & AliESDtrack::kTPCrefit)){
353 fQA->Fill("h_ST_TPCrefit", 1);
354 }
355 if(!(status & AliESDtrack::kTPCrefit)){
356 fQA->Fill("h_ST_TPCrefit", 0);
357 return kFALSE;
358 }
359
360 // Chi2 per TPC cluster
361 Int_t nTPCclusters = track->GetTPCclusters(0);
362 Float_t chi2perTPCcluster = track->GetTPCchi2()/Float_t(nTPCclusters);
363 fQA->Fill("h_ST_chi2TPCcls", chi2perTPCcluster);
3a72645a 364 if(chi2perTPCcluster > 3.5) return kFALSE; // 4.0
c04c80e6 365
366 // TPC cluster ratio
367 Float_t cRatioTPC = track->GetTPCNclsF() > 0. ? static_cast<Float_t>(track->GetTPCNcls())/static_cast<Float_t> (track->GetTPCNclsF()) : 1.;
368 fQA->Fill("h_ST_TPCclsR", cRatioTPC);
369 if(cRatioTPC < 0.6) return kFALSE;
370
371 // kinks
372 fQA->Fill("h_ST_kinks", track->GetKinkIndex(0));
373 if(track->GetKinkIndex(0) != 0) return kFALSE;
374
375 // pt
376 fQA->Fill("h_ST_pt",track->Pt());
3a72645a 377 if(track->Pt() < 0.1 || track->Pt() > 100) return kFALSE; //
c04c80e6 378
379 // eta
380 fQA->Fill("h_ST_eta", track->Eta());
381 //if(TMath::Abs(track->Eta()) > 0.9) return kFALSE;
382
383 return kTRUE;
384}
385//________________________________________________________________
386Bool_t AliHFEV0cuts::V0CutsCommon(AliESDv0 *v0){
387 //
388 // V0 cuts common to all V0s
389 //
3a72645a 390
c04c80e6 391 AliESDtrack* dN, *dP;
392
393 dP = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(v0->GetPindex()));
394 dN = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(v0->GetNindex()));
395
396 if(!dN || !dP) return kFALSE;
397
398 Int_t qP = dP->Charge();
399 Int_t qN = dN->Charge();
400
401 if((qP*qN) != -1) return kFALSE;
402
403 return kTRUE;
404}
405//________________________________________________________________
406Bool_t AliHFEV0cuts::GammaCuts(AliESDv0 *v0){
407 //
408 // gamma cuts
409 //
410
411 if(!v0) return kFALSE;
412
3a72645a 413 if(fMCEvent){
414 if(1 == fCurrentV0id){
415 fQAmc->Fill("h_Mass_gamma_as_K0", v0->GetEffMass(2, 2));
416 fQAmc->Fill("h_Mass_gamma_as_L", v0->GetEffMass(2, 4));
417 fQAmc->Fill("h_Mass_gamma_as_L", v0->GetEffMass(4, 2));
418 }
419 }
420
c04c80e6 421 // loose cuts first
3a72645a 422 //if(LooseRejectK0(v0) || LooseRejectLambda(v0)) return kFALSE;
c04c80e6 423
424 AliVTrack* daughter[2];
425 Int_t pIndex = 0, nIndex = 0;
426 if(CheckSigns(v0)){
427 pIndex = v0->GetPindex();
428 nIndex = v0->GetNindex();
429 }
430 else{
431 pIndex = v0->GetNindex();
432 nIndex = v0->GetPindex();
433 }
434 daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(pIndex));
435 daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(nIndex));
436 if(!daughter[0] || !daughter[1]) return kFALSE;
437
438 AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kElectron), TMath::Abs(kElectron));
439 if(!kfMother) return kFALSE;
440
441 // production vertex is set in the 'CreateMotherParticle' function
3a72645a 442 //kfMother->SetMassConstraint(0, 0.001);
c04c80e6 443
444 AliESDtrack* d[2];
445 d[0] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(pIndex));
446 d[1] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(nIndex));
447
448 Float_t iMass = v0->GetEffMass(0, 0);
449 Float_t iP = v0->P();
450 Float_t p[2] = {d[0]->GetP(), d[1]->GetP()};
451
c04c80e6 452 // Cut values
3a72645a 453 const Double_t cutChi2NDF = 40.; // ORG [7.]
454 const Double_t cutCosPoint[2] = {0., 0.02}; // ORG [0., 0.03]
c04c80e6 455 const Double_t cutDCA[2] = {0., 0.25}; // ORG [0., 0.25]
3a72645a 456 const Double_t cutProdVtxR[2] = {8., 90.}; // ORG [6., 9999]
c04c80e6 457 const Double_t cutPsiPair[2] = {0., 0.05}; // ORG [0. 0.05]
3a72645a 458 const Double_t cutOAngle[2] = {0, 0.1}; // ORG [0., 0.1]
459 // mass cut
c04c80e6 460 const Double_t cutMass = 0.05; // ORG [0.05]
c04c80e6 461 // Values
462
463 // cos pointing angle
464 Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
465 cosPoint = TMath::ACos(cosPoint);
466
467 // DCA between daughters
468 Double_t dca = v0->GetDcaV0Daughters();
469
470 // Production vertex
471 Double_t x, y, z;
472 v0->GetXYZ(x,y,z);
473 Double_t r = TMath::Sqrt(x*x + y*y);
474
3a72645a 475 Double_t xy[2];
476 Double_t r2 = -1.;
477 if ( GetConvPosXY(d[0], d[1], xy) ){
478 r2 = TMath::Sqrt(xy[0]*xy[0] + xy[1]*xy[1]);
479 }
480
481
c04c80e6 482 // Opening angle
483 Double_t oAngle = OpenAngle(v0);
484
485 // psi pair
486 Double_t psiPair = PsiPair(v0);
487
488 // V0 chi2/ndf
489 Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
490
491 if(kfMother) delete kfMother;
492
493 //
494 // Apply the cuts, produce QA plots (with mass cut)
495 //
496 fQA->Fill("h_Gamma_Mass", 0, iMass);
3a72645a 497
498 // MC
499 if(fMCEvent){
500 if(1 == fCurrentV0id){
501 fQAmc->Fill("h_Gamma_Mass_S", 0, iMass);
502 fQAmc->Fill("h_gamma_MvP_S", iP, iMass);
503 }
504 else if(-2 != fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_B", 0, iMass);
505 }
506 // cut distributions
c04c80e6 507 if(iMass < cutMass){
c04c80e6 508 fQA->Fill("h_Electron_P", 0, p[0]);
509 fQA->Fill("h_Electron_P", 0, p[1]);
510 fQA->Fill("h_cut_Gamma_CosPoint", 0, cosPoint);
c04c80e6 511 fQA->Fill("h_cut_Gamma_DCA", 0, dca);
3a72645a 512 fQA->Fill("h_cut_Gamma_VtxR_old", 0, r);
513 fQA->Fill("h_cut_Gamma_VtxR", 0, r2);
c04c80e6 514 fQA->Fill("h_cut_Gamma_OA", 0, oAngle);
515 fQA->Fill("h_cut_Gamma_PP", 0, psiPair);
516 fQA->Fill("h_cut_Gamma_Chi2", 0, chi2ndf);
3a72645a 517 fQA->Fill("h_cut_Gamma_Chi2", 1, chi2ndf, iP);
c04c80e6 518 fQA->Fill("h_cut_Gamma_OAvP", 0, iP, oAngle);
3a72645a 519
520 if(fMCEvent){
521 // MC signal
522 if(1 == fCurrentV0id){
523 fQAmc->Fill("h_cut_Gamma_CosPoint_S", 0, iP, cosPoint);
524 fQAmc->Fill("h_cut_Gamma_DCA_S", 0, iP, dca);
525 fQAmc->Fill("h_cut_Gamma_VtxR_S", 0, iP, r2);
526 fQAmc->Fill("h_cut_Gamma_OA_S", 0, iP, oAngle);
527 fQAmc->Fill("h_cut_Gamma_PP_S", 0, iP, psiPair);
528 fQAmc->Fill("h_cut_Gamma_Chi2_S", 0, iP, chi2ndf);
529 fQAmc->Fill("h_cut_Gamma_Chi2_S", 1, iP, chi2ndf);
530 fQAmc->Fill("h_Electron_P_S", 0, p[0]);
531 fQAmc->Fill("h_Electron_P_S", 0, p[1]);
532 }
533 // MC background
534 else if(-2 != fCurrentV0id){
535 fQAmc->Fill("h_cut_Gamma_CosPoint_B", 0, iP, cosPoint);
536 fQAmc->Fill("h_cut_Gamma_DCA_B", 0, iP, dca);
537 fQAmc->Fill("h_cut_Gamma_VtxR_B", 0, iP, r2);
538 fQAmc->Fill("h_cut_Gamma_OA_B", 0, iP, oAngle);
539 fQAmc->Fill("h_cut_Gamma_PP_B", 0, iP, psiPair);
540 fQAmc->Fill("h_cut_Gamma_Chi2_B", 0, iP, chi2ndf);
541 fQAmc->Fill("h_cut_Gamma_Chi2_B", 1, iP, chi2ndf);
542 fQAmc->Fill("h_Electron_P_B", 0, p[0]);
543 fQAmc->Fill("h_Electron_P_B", 0, p[1]);
544 }
545 }
546 }
c04c80e6 547
3a72645a 548
549 //
550 // Chi2/NDF cut
551 //
552 if(chi2ndf > cutChi2NDF) return kFALSE;
c04c80e6 553 fQA->Fill("h_Gamma_Mass", 1, iMass);
554 if(iMass < cutMass){
3a72645a 555 fQA->Fill("h_cut_Gamma_CosPoint", 1, cosPoint);
c04c80e6 556 fQA->Fill("h_Electron_P", 1, p[0]);
557 fQA->Fill("h_Electron_P", 1, p[1]);
c04c80e6 558 }
3a72645a 559 if(fMCEvent){
560 if(1 == fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_S", 1, iMass);
561 else if(-2 != fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_B", 1, iMass);
562 if(iMass < cutMass){
563 if(1 == fCurrentV0id){
564 fQAmc->Fill("h_cut_Gamma_CosPoint_S", 1, iP, cosPoint);
565 fQAmc->Fill("h_Electron_P_S", 1, p[0]);
566 fQAmc->Fill("h_Electron_P_S", 1, p[1]);
567 }
568 else if(-2 != fCurrentV0id){
569 fQAmc->Fill("h_cut_Gamma_CosPoint_B", 1, iP, cosPoint);
570 fQAmc->Fill("h_Electron_P_B", 1, p[0]);
571 fQAmc->Fill("h_Electron_P_B", 1, p[1]);
572 }
573 }
574 }
575
576 //
577 // Cos point cut
578 //
579 if(cosPoint < cutCosPoint[0] || cosPoint > cutCosPoint[1]) return kFALSE;
c04c80e6 580 fQA->Fill("h_Gamma_Mass", 2, iMass);
581 if(iMass < cutMass){
582 fQA->Fill("h_Electron_P", 2, p[0]);
583 fQA->Fill("h_Electron_P", 2, p[1]);
3a72645a 584 fQA->Fill("h_cut_Gamma_DCA", 1, dca);
c04c80e6 585 }
3a72645a 586 if(fMCEvent){
587 if(1 == fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_S", 2, iMass);
588 else if(-2 != fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_B", 2, iMass);
589 if(iMass < cutMass){
590 if(1 == fCurrentV0id){
591 fQAmc->Fill("h_cut_Gamma_DCA_S", 1, iP, dca);
592 fQAmc->Fill("h_Electron_P_S", 2, p[0]);
593 fQAmc->Fill("h_Electron_P_S", 2, p[1]);
594
595 }
596 else if(-2 != fCurrentV0id){
597 fQAmc->Fill("h_cut_Gamma_DCA_B", 1, iP, dca);
598 fQAmc->Fill("h_Electron_P_B", 2, p[0]);
599 fQAmc->Fill("h_Electron_P_B", 2, p[1]);
600
601 }
602 }
603 }
604
605 //
606 // DCA cut
607 //
608 if(dca < cutDCA[0] || dca > cutDCA[1]) return kFALSE;
c04c80e6 609 fQA->Fill("h_Gamma_Mass", 3, iMass);
610 if(iMass < cutMass){
611 fQA->Fill("h_Electron_P", 3, p[0]);
612 fQA->Fill("h_Electron_P", 3, p[1]);
3a72645a 613 fQA->Fill("h_cut_Gamma_VtxR_old", 1, r);
614 fQA->Fill("h_cut_Gamma_VtxR", 1, r2);
615 }
616 if(fMCEvent){
617 if(1 == fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_S", 3, iMass);
618 else if(-2 != fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_B", 3, iMass);
619 if(iMass < cutMass){
620 if(1 == fCurrentV0id){
621 fQAmc->Fill("h_cut_Gamma_VtxR_S", 1, iP, r2);
622 fQAmc->Fill("h_Electron_P_S", 3, p[0]);
623 fQAmc->Fill("h_Electron_P_S", 3, p[1]);
624
625 }
626 else if(-2 != fCurrentV0id){
627 fQAmc->Fill("h_cut_Gamma_VtxR_B", 1, iP, r2);
628 fQAmc->Fill("h_Electron_P_B", 3, p[0]);
629 fQAmc->Fill("h_Electron_P_B", 3, p[1]);
630 }
631 }
c04c80e6 632 }
633
3a72645a 634 //
635 // Vertex radius cut
636 //
637 if(r < cutProdVtxR[0] || r > cutProdVtxR[1]) return kFALSE;
c04c80e6 638 fQA->Fill("h_Gamma_Mass", 4, iMass);
639 if(iMass < cutMass){
3a72645a 640 fQA->Fill("h_cut_Gamma_PP", 1, psiPair);
c04c80e6 641 fQA->Fill("h_Electron_P", 4, p[0]);
642 fQA->Fill("h_Electron_P", 4, p[1]);
3a72645a 643 }
644 if(fMCEvent){
645 if(1 == fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_S", 4, iMass);
646 else if(-2 != fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_B", 4, iMass);
647 if(iMass < cutMass){
648 if(1 == fCurrentV0id){
649 fQAmc->Fill("h_cut_Gamma_PP_S", 1, iP, psiPair);
650 fQAmc->Fill("h_Electron_P_S", 4, p[0]);
651 fQAmc->Fill("h_Electron_P_S", 4, p[1]);
652 }
653 else if(-2 != fCurrentV0id){
654 fQAmc->Fill("h_cut_Gamma_PP_B", 1, iP, psiPair);
655 fQAmc->Fill("h_Electron_P_B", 4, p[0]);
656 fQAmc->Fill("h_Electron_P_B", 4, p[1]);
657 }
658 }
c04c80e6 659 }
660
3a72645a 661
662 //
663 // PsiPair cut
664 //
c04c80e6 665 if(psiPair < cutPsiPair[0] || psiPair > cutPsiPair[1]) return kFALSE;
666 fQA->Fill("h_Gamma_Mass", 5, iMass);
667 if(iMass < cutMass){
3a72645a 668 fQA->Fill("h_cut_Gamma_OA", 1, oAngle);
669 fQA->Fill("h_cut_Gamma_OAvP", 1, iP, oAngle);
c04c80e6 670 fQA->Fill("h_Electron_P", 5, p[0]);
671 fQA->Fill("h_Electron_P", 5, p[1]);
3a72645a 672 }
673 if(fMCEvent){
674 if(1 == fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_S", 5, iMass);
675 else if(-2 != fCurrentV0id)fQAmc->Fill("h_Gamma_Mass_B", 5, iMass);
676
677 if(iMass < cutMass){
678 if(1 == fCurrentV0id){
679 fQAmc->Fill("h_cut_Gamma_OA_S", 1, iP, oAngle);
680 fQAmc->Fill("h_Electron_P_S", 5, p[0]);
681 fQAmc->Fill("h_Electron_P_S", 5, p[1]);
682 }
683 else if(-2 != fCurrentV0id){
684 fQAmc->Fill("h_cut_Gamma_OA_B", 1, iP, oAngle);
685 fQAmc->Fill("h_Electron_P_B", 5, p[0]);
686 fQAmc->Fill("h_Electron_P_B", 5, p[1]);
687 }
688 }
c04c80e6 689 }
690
3a72645a 691 //
692 // Opening angle cut (obsolete?)
693 //
694 if(oAngle < cutOAngle[0] || oAngle > cutOAngle[1]) return kFALSE;
c04c80e6 695 fQA->Fill("h_Gamma_Mass", 6, iMass);
696 if(iMass < cutMass){
697 fQA->Fill("h_Electron_P", 6, p[0]);
698 fQA->Fill("h_Electron_P", 6, p[1]);
c04c80e6 699 }
3a72645a 700 if(fMCEvent){
701 if(1 == fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_S", 6, iMass);
702 else if(-2 != fCurrentV0id) fQAmc->Fill("h_Gamma_Mass_B", 6, iMass);
703 if(iMass < cutMass){
704 if(1 == fCurrentV0id){
705 fQAmc->Fill("h_Electron_P_S", 6, p[0]);
706 fQAmc->Fill("h_Electron_P_S", 6, p[1]);
707 }
708 else if(-2 != fCurrentV0id){
709 fQAmc->Fill("h_Electron_P_B", 6, p[0]);
710 fQAmc->Fill("h_Electron_P_B", 6, p[1]);
711 }
712 }
c04c80e6 713 }
714
715
716 if(iMass > cutMass) return kFALSE;
717
718 // all cuts passed
3a72645a 719
720
721 // some MC stuff
722 //printf("**D: gamma V0id: %i, P: %i, N: %i \n", fCurrentV0id, fPdaughterPDG, fNdaughterPDG);
723 if(1 == fCurrentV0id){
724 fQAmc->Fill("h_gamma_p_S", iP);
725 fQAmc->Fill("h_Electron_P_S", 7, p[0]);
726 fQAmc->Fill("h_Electron_P_S", 7, p[1]);
727 }
728 else if (-2 != fCurrentV0id){
729 fQAmc->Fill("h_gamma_p_B", iP);
730 fQAmc->Fill("h_Electron_P_B", 7, p[0]);
731 fQAmc->Fill("h_Electron_P_B", 7, p[1]);
732 }
733
734
c04c80e6 735 return kTRUE;
736}
737//________________________________________________________________
738Bool_t AliHFEV0cuts::K0Cuts(AliESDv0 *v0){
739 //
740 // K0 cuts
741 //
742
743 if(!v0) return kFALSE;
744
3a72645a 745 if(fMCEvent){
746 if(2 == fCurrentV0id){
747 fQAmc->Fill("h_Mass_K0_as_G", v0->GetEffMass(0, 0));
748 fQAmc->Fill("h_Mass_K0_as_L", v0->GetEffMass(2, 4));
749 fQAmc->Fill("h_Mass_K0_as_L", v0->GetEffMass(4, 2));
750 }
751 }
752
753 //const Double_t cK0mass=TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(); // PDG K0s mass
c04c80e6 754 AliVTrack* daughter[2];
755 Int_t pIndex = 0, nIndex = 0;
756 if(CheckSigns(v0)){
757 pIndex = v0->GetPindex();
758 nIndex = v0->GetNindex();
759 }
760 else{
761 pIndex = v0->GetNindex();
762 nIndex = v0->GetPindex();
763 }
764
765 daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(pIndex));
766 daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(nIndex));
767 if(!daughter[0] || !daughter[1]) return kFALSE;
768
769 AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kPiPlus));
770 if(!kfMother) return kFALSE;
771 // production vertex is set in the 'CreateMotherParticle' function
3a72645a 772 //kfMother->SetMassConstraint(cK0mass, 0.);
c04c80e6 773
774 AliESDtrack* d[2];
775 d[0] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(pIndex));
776 d[1] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(nIndex));
777
778 Float_t iMass = v0->GetEffMass(2, 2);
779 Float_t iP = v0->P();
780 Float_t p[2] = {d[0]->GetP(), d[1]->GetP()};
781 Double_t theta = v0->Theta();
782 Double_t phi = v0->Phi();
783 Double_t pt = v0->Pt();
784 Double_t data[4] = {0., 0., 0., 0.};
3a72645a 785
c04c80e6 786 // Cut values
3a72645a 787 const Double_t cutChi2NDF = 40.; // ORG [7.]
788 const Double_t cutCosPoint[2] = {0., 0.02}; // ORG [0., 0.03]
c04c80e6 789 const Double_t cutDCA[2] = {0., 0.2}; // ORG [0., 0.1]
790 const Double_t cutProdVtxR[2] = {2.0, 30.}; // ORG [0., 8.1]
791 const Double_t cutMass[2] = {0.49, 0.51}; // ORG [0.485, 0.51]
3a72645a 792 //const Double_t cutOAngleP = (1.0/(iP + 0.3) - 0.1); // momentum dependent min. OAngle ~ 1/x
c04c80e6 793 // Values
794
795 // cos pointing angle
796 Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
797 cosPoint = TMath::ACos(cosPoint);
798
799 // DCA between daughters
800 Double_t dca = v0->GetDcaV0Daughters();
801
802 // Production vertex
803 Double_t x, y, z;
804 v0->GetXYZ(x,y,z);
805
806 Double_t r = TMath::Sqrt(x*x + y*y);
807
808 // V0 chi2/ndf
809 Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
810
811 if(kfMother) delete kfMother;
812
813 // Opening angle
814 Double_t oAngle = OpenAngle(v0);
815
c04c80e6 816 //
817 // Apply the cuts, produce QA plots (with mass cut)
818 //
819
820 fQA->Fill("h_K0_Mass", 0, iMass);
3a72645a 821 // MC
822 if(fMCEvent){
823 if(2 == fCurrentV0id){
824 fQAmc->Fill("h_K0_Mass_S", 0, iMass);
825 fQAmc->Fill("h_K0_MvP_S", iP, iMass);
826 }
827 else if(-2 != fCurrentV0id) fQAmc->Fill("h_K0_Mass_B", 0, iMass);
828 }
829
c04c80e6 830 if(iMass > cutMass[0] && iMass < cutMass[1]){
c04c80e6 831 fQA->Fill("h_PionK0_P", 0, p[0]);
832 fQA->Fill("h_PionK0_P", 0, p[1]);
833 fQA->Fill("h_cut_K0_CosPoint", 0, cosPoint);
c04c80e6 834 fQA->Fill("h_cut_K0_DCA", 0, dca);
835 fQA->Fill("h_cut_K0_VtxR", 0, r);
836 fQA->Fill("h_cut_K0_Chi2", 0, chi2ndf);
3a72645a 837 fQA->Fill("h_cut_K0_Chi2", 1, chi2ndf);
838 }
839
840 // MC
841 if(fMCEvent){
842 if(iMass > cutMass[0] && iMass < cutMass[1]){
843 if(2 == fCurrentV0id){
844 fQAmc->Fill("h_cut_K0_CosPoint_S", 0, iP, cosPoint);
845 fQAmc->Fill("h_cut_K0_DCA_S", 0, iP, dca);
846 fQAmc->Fill("h_cut_K0_VtxR_S", 0, iP, r);
847 fQAmc->Fill("h_cut_K0_Chi2_S", 0, iP, chi2ndf);
848 fQAmc->Fill("h_cut_K0_Chi2_S", 1, iP, chi2ndf);
849 fQAmc->Fill("h_cut_K0_OA_S", 0, iP, oAngle);
850 fQAmc->Fill("h_PionK0_P_S", 0, p[0]);
851 fQAmc->Fill("h_PionK0_P_S", 0, p[1]);
852 }
853 else if(-2 != fCurrentV0id){
854 fQAmc->Fill("h_cut_K0_CosPoint_B", 0, iP, cosPoint);
855 fQAmc->Fill("h_cut_K0_DCA_B", 0, iP, dca);
856 fQAmc->Fill("h_cut_K0_VtxR_B", 0, iP, r);
857 fQAmc->Fill("h_cut_K0_Chi2_B", 0, iP, chi2ndf);
858 fQAmc->Fill("h_cut_K0_Chi2_B", 1, iP, chi2ndf);
859 fQAmc->Fill("h_cut_K0_OA_B", 0, iP, oAngle);
860 fQAmc->Fill("h_PionK0_P_B", 0, p[0]);
861 fQAmc->Fill("h_PionK0_P_B", 0, p[1]);
862 }
863 }
c04c80e6 864 }
865
3a72645a 866 //
867 // Chi2/NDF cut
868 //
869 if(chi2ndf > cutChi2NDF) return kFALSE;
c04c80e6 870 fQA->Fill("h_K0_Mass", 1, iMass);
871 if(iMass > cutMass[0] && iMass < cutMass[1]){
3a72645a 872 fQA->Fill("h_cut_K0_CosPoint", 1, cosPoint);
c04c80e6 873 fQA->Fill("h_PionK0_P", 1, p[0]);
874 fQA->Fill("h_PionK0_P", 1, p[1]);
3a72645a 875 }
876 if(fMCEvent){
877 if(2 == fCurrentV0id) fQAmc->Fill("h_K0_Mass_S", 1, iMass);
878 else if(-2 != fCurrentV0id) fQAmc->Fill("h_K0_Mass_B", 1, iMass);
879 if(iMass > cutMass[0] && iMass < cutMass[1]){
880 if(2 == fCurrentV0id){
881 fQAmc->Fill("h_cut_K0_CosPoint_S", 1, iP, cosPoint);
882 fQAmc->Fill("h_PionK0_P_S", 1, p[0]);
883 fQAmc->Fill("h_PionK0_P_S", 1, p[1]);
884 }
885 else if(-2 != fCurrentV0id){
886 fQAmc->Fill("h_cut_K0_CosPoint_B", 1, iP, cosPoint);
887 fQAmc->Fill("h_PionK0_P_B", 1, p[0]);
888 fQAmc->Fill("h_PionK0_P_B", 1, p[1]);
889 }
890 }
c04c80e6 891 }
3a72645a 892
893 //
894 // Cos point cut
895 //
896 if(cosPoint < cutCosPoint[0] || cosPoint > cutCosPoint[1]) return kFALSE;
c04c80e6 897 fQA->Fill("h_K0_Mass", 2, iMass);
898 if(iMass > cutMass[0] && iMass < cutMass[1]){
899 fQA->Fill("h_PionK0_P", 2, p[0]);
900 fQA->Fill("h_PionK0_P", 2, p[1]);
3a72645a 901 fQA->Fill("h_cut_K0_DCA", 1, dca);
902 }
903 if(fMCEvent){
904 if(2 == fCurrentV0id) fQAmc->Fill("h_K0_Mass_S", 2, iMass);
905 else if(-2 != fCurrentV0id) fQAmc->Fill("h_K0_Mass_B", 2, iMass);
906 if(iMass > cutMass[0] && iMass < cutMass[1]){
907 if(2 == fCurrentV0id){
908 fQAmc->Fill("h_cut_K0_DCA_S", 1, iP, dca);
909 fQAmc->Fill("h_PionK0_P_S", 2, p[0]);
910 fQAmc->Fill("h_PionK0_P_S", 2, p[1]);
911 }
912 else if(-2 != fCurrentV0id){
913 fQAmc->Fill("h_cut_K0_DCA_B", 1, iP, dca);
914 fQAmc->Fill("h_PionK0_P_B", 2, p[0]);
915 fQAmc->Fill("h_PionK0_P_B", 2, p[1]);
916 }
917 }
c04c80e6 918 }
c04c80e6 919
3a72645a 920
921 //
922 // DCA cut
923 //
924 if(dca < cutDCA[0] || dca > cutDCA[1]) return kFALSE;
c04c80e6 925 fQA->Fill("h_K0_Mass", 3, iMass);
926 if(iMass > cutMass[0] && iMass < cutMass[1]){
927 fQA->Fill("h_PionK0_P", 3, p[0]);
928 fQA->Fill("h_PionK0_P", 3, p[1]);
3a72645a 929 fQA->Fill("h_cut_K0_VtxR", 1, r);
930 }
931 if(fMCEvent){
932 if(2 == fCurrentV0id) fQAmc->Fill("h_K0_Mass_S", 3, iMass);
933 else if(-2 != fCurrentV0id) fQAmc->Fill("h_K0_Mass_B", 3, iMass);
934 if(iMass > cutMass[0] && iMass < cutMass[1]){
935 if(2 == fCurrentV0id){
936 fQAmc->Fill("h_cut_K0_VtxR_S", 1, iP, r);
937 fQAmc->Fill("h_PionK0_P_S", 3, p[0]);
938 fQAmc->Fill("h_PionK0_P_S", 3, p[1]);
939 }
940 else if(-2 != fCurrentV0id){
941 fQAmc->Fill("h_cut_K0_VtxR_B", 1, iP, r);
942 fQAmc->Fill("h_PionK0_P_B", 3, p[0]);
943 fQAmc->Fill("h_PionK0_P_B", 3, p[1]);
944 }
945 }
c04c80e6 946 }
947
3a72645a 948
949 //
950 // Vertex R cut
951 //
952 if(r < cutProdVtxR[0] || r > cutProdVtxR[1]) return kFALSE;
c04c80e6 953 fQA->Fill("h_K0_Mass", 4, iMass);
954 if(iMass > cutMass[0] && iMass < cutMass[1]){
955 fQA->Fill("h_PionK0_P", 4, p[0]);
956 fQA->Fill("h_PionK0_P", 4, p[1]);
957 fQA->Fill("h_cut_K0_OAvP", 1, iP, oAngle);
c04c80e6 958 }
3a72645a 959 if(fMCEvent){
960 if(2 == fCurrentV0id) fQAmc->Fill("h_K0_Mass_S", 4, iMass);
961 else if(-2 != fCurrentV0id) fQAmc->Fill("h_K0_Mass_B", 4, iMass);
962 if(iMass > cutMass[0] && iMass < cutMass[1]){
963 if(2 == fCurrentV0id){
964 fQAmc->Fill("h_cut_K0_OA_S", 1, iP, oAngle);
965 fQAmc->Fill("h_PionK0_P_S", 4, p[0]);
966 fQAmc->Fill("h_PionK0_P_S", 4, p[1]);
967 }
968 else if(-2 != fCurrentV0id){
969 fQAmc->Fill("h_cut_K0_OA_B", 1, iP, oAngle);
970 fQAmc->Fill("h_PionK0_P_B", 4, p[0]);
971 fQAmc->Fill("h_PionK0_P_B", 4, p[1]);
972 }
973 }
c04c80e6 974 }
975
976 data[0] = iMass;
977 data[1] = pt;
978 data[2] = theta;
979 data[3] = phi;
980 //printf("-D: m: %f, pT: %f, theta: %f, phi: %f\n", invMass, mPt, theta, phi);
981 fQA->Fill("hK0", data);
982
983
984 if(iMass < cutMass[0] || iMass > cutMass[1]) return kFALSE;
985
986 // all cuts passed
987
3a72645a 988 // some MC stuff
989 if(2 == fCurrentV0id){
990 fQAmc->Fill("h_K0_p_S", iP);
991 fQAmc->Fill("h_PionK0_P_S", 5, p[0]);
992 fQAmc->Fill("h_PionK0_P_S", 5, p[1]);
993 }
994 else if (-2 != fCurrentV0id){
995 fQAmc->Fill("h_K0_p_B", iP);
996 fQAmc->Fill("h_PionK0_P_B", 5, p[0]);
997 fQAmc->Fill("h_PionK0_P_B", 5, p[1]);
998 }
999
c04c80e6 1000 return kTRUE;
1001}
1002//________________________________________________________________
1003Bool_t AliHFEV0cuts::LambdaCuts(AliESDv0 *v0, Bool_t &isLambda ){
1004 //
1005 // Lambda cuts - decision on Lambda - AntiLambda is taken too
1006 //
1007 // discrimination between lambda and antilambda - correlation of the following variables necessary:
1008 // - momentum of the proton AND momentum of the pion (proton momentum is allways larger)
1009 // - mass of the mother particle
1010
1011 if(!v0) return kFALSE;
1012
3a72645a 1013 if(fMCEvent){
1014 if(4 == fCurrentV0id){
1015 fQAmc->Fill("h_Mass_L_as_G", v0->GetEffMass(0, 0));
1016 fQAmc->Fill("h_Mass_L_as_K0", v0->GetEffMass(2, 0));
1017 }
1018 }
c04c80e6 1019 // loose cuts first
3a72645a 1020 //if(LooseRejectK0(v0) || LooseRejectGamma(v0)) return kFALSE;
c04c80e6 1021
1022 const Double_t cL0mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(); // PDG lambda mass
1023
1024 AliVTrack* daughter[2];
1025 Int_t pIndex = 0, nIndex = 0;
1026 Float_t mMass[2] = {-1., -1.};
1027 if(CheckSigns(v0)){
1028 pIndex = v0->GetPindex();
1029 nIndex = v0->GetNindex();
1030 mMass[0] = v0->GetEffMass(4, 2);
1031 mMass[1] = v0->GetEffMass(2, 4);
1032 }
1033 else{
1034 pIndex = v0->GetNindex();
1035 nIndex = v0->GetPindex();
1036 mMass[0] = v0->GetEffMass(2, 4);
1037 mMass[1] = v0->GetEffMass(4, 2);
1038 }
1039
1040 daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(pIndex));
1041 daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(nIndex));
1042 if(!daughter[0] || !daughter[1]) return kFALSE;
1043
1044 AliKFParticle *kfMother[2] = {0x0, 0x0};
1045 // Lambda
1046 kfMother[0] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kProton), TMath::Abs(kPiPlus));
1047 if(!kfMother[0]) return kFALSE;
1048
1049 // production vertex is set in the 'CreateMotherParticle' function
3a72645a 1050 //kfMother[0]->SetMassConstraint(cL0mass, 0.);
c04c80e6 1051
1052 // Anti Lambda
1053 kfMother[1] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kProton));
1054 if(!kfMother[1]) return kFALSE;
1055 // production vertex is set in the 'CreateMotherParticle' function
3a72645a 1056 //kfMother[1]->SetMassConstraint(cL0mass, 0.);
c04c80e6 1057
1058 Float_t dMass[2] = {TMath::Abs(mMass[0] - cL0mass), TMath::Abs(mMass[1] - cL0mass)};
1059
1060 AliESDtrack* d[2];
1061 d[0] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(pIndex));
1062 d[1] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(nIndex));
1063 if(!d[0] || !d[1]) return kFALSE;
1064
1065 Float_t p[2] = {d[0]->GetP(), d[1]->GetP()};
1066
1067 // check the 3 lambda - antilambda variables
1068 Int_t check[2] = {-1, -1}; // 0 : lambda, 1 : antilambda
1069 // 1) momentum of the daughter particles - proton is expected to have higher momentum than pion
1070 check[0] = (p[0] > p[1]) ? 0 : 1;
1071 // 2) mass of the mother particle
1072 check[1] = (dMass[0] < dMass[1]) ? 0 : 1;
1073 fQA->Fill("h_L_checks", check[0]*1.0, check[1]*1.0);
1074
1075 // if the two check do not agree
1076 if(check[0] != check[1]){
1077 if(kfMother[0]) delete kfMother[0];
1078 if(kfMother[1]) delete kfMother[1];
1079 return kFALSE;
1080 }
1081
1082 // now that the check[0] == check[1]
1083 const Int_t type = check[0];
1084
1085 Float_t iMass =0.;
1086 if(CheckSigns(v0)){
1087 iMass = (type == 0) ? v0->GetEffMass(4, 2) : v0->GetEffMass(2, 4);
1088 }
1089 else{
1090 iMass = (type == 0) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);
1091 }
1092 Float_t iP = v0->P();
1093
1094 // Cuts
3a72645a 1095 const Double_t cutChi2NDF = 40.; // ORG [5.]
1096 const Double_t cutCosPoint[2] = {0., 0.02}; // ORG [0., 0.03]
c04c80e6 1097 const Double_t cutDCA[2] = {0., 0.2}; // ORG [0., 0.2]
3a72645a 1098 const Double_t cutProdVtxR[2] = {2., 40.}; // ORG [0., 24.]
c04c80e6 1099 const Double_t cutMass[2] = {1.11, 1.12}; // ORG [1.11, 1.12]
c04c80e6 1100 // cundidate cuts
1101 // opening angle as a function of L momentum
3a72645a 1102 //const Double_t cutOAngleP = 0.3 - 0.2*iP; // momentum dependent min. OAngle linear cut
c04c80e6 1103 // relative daughter momentum versusu mother momentum
c04c80e6 1104
1105 // compute the cut values
1106
1107 // cos pointing angle
1108 Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
1109 cosPoint = TMath::ACos(cosPoint);
1110
1111 // DCA between daughters
1112 Double_t dca = v0->GetDcaV0Daughters();
1113
1114 // Production vertex
1115 Double_t x, y, z;
1116 v0->GetXYZ(x,y,z);
1117 Double_t r = TMath::Sqrt(x*x + y*y);
1118
3a72645a 1119 // proton - pion indices
c04c80e6 1120 Int_t ix[2] = {0, 1};
1121 if(1 == type){
1122 ix[0] = 1;
1123 ix[1] = 0;
1124 }
3a72645a 1125
1126 // proton - pion indices - based on MC truth
1127 // for background use the reconstructed indices
1128 Int_t ixMC[2] = {-1, -1}; // {proton, pion}
1129 if(fMCEvent){
1130 if(4 == fCurrentV0id){
1131 ixMC[0] = 0;
1132 ixMC[1] = 1;
1133 }
1134 else if(-4 == fCurrentV0id){
1135 ixMC[0] = 1;
1136 ixMC[1] = 0;
1137 }
1138 else{
1139 ixMC[0] = ix[0];
1140 ixMC[1] = ix[1];
1141 }
1142 }
1143
c04c80e6 1144 // V0 chi2/ndf
1145 Double_t chi2ndf = kfMother[type]->GetChi2()/kfMother[type]->GetNDF();
1146
1147 if(kfMother[0]) delete kfMother[0];
1148 if(kfMother[1]) delete kfMother[1];
1149
1150 // Opening angle
1151 Double_t oAngle = OpenAngle(v0);
1152
1153 // Relative daughter momentum
1154 Double_t rP = (0 == check[0]) ? p[1]/p[0] : p[0]/p[1];
1155
c04c80e6 1156
1157 //
1158 // Apply the cuts, produce QA plots (with mass cut)
1159 //
1160
1161 (type == 0) ? fQA->Fill("h_L_Mass", 0, iMass) : fQA->Fill("h_AL_Mass", 0, iMass);
1162
3a72645a 1163
1164
1165 // MC
1166 if(fMCEvent){
1167 if(4 == fCurrentV0id){
1168 fQAmc->Fill("h_L_Mass_S", 0, iMass);
1169 fQAmc->Fill("h_lambda_MvP_S", iP, iMass);
1170 }
1171 else if(-4 == fCurrentV0id){
1172 fQAmc->Fill("h_AL_Mass_S", 0, iMass);
1173 fQAmc->Fill("h_lambda_MvP_S", iP, iMass);
1174 }
1175 else if(-2 != fCurrentV0id) fQAmc->Fill("h_LAL_Mass_B", 0, iMass);
1176 }
1177
1178
c04c80e6 1179 if(iMass > cutMass[0] && iMass < cutMass[1]){
c04c80e6 1180 fQA->Fill("h_ProtonL_P", 0, p[ix[0]]);
1181 fQA->Fill("h_PionL_P", 0, p[ix[1]]);
1182 fQA->Fill("h_cut_L_Chi2", 0, chi2ndf);
3a72645a 1183 fQA->Fill("h_cut_L_Chi2", 1, chi2ndf);
c04c80e6 1184 fQA->Fill("h_cut_L_CosPoint", 0, cosPoint);
c04c80e6 1185 fQA->Fill("h_cut_L_DCA", 0, dca);
1186 fQA->Fill("h_cut_L_VtxR", 0, r);
1187 fQA->Fill("h_cut_L_OAvP", 0, iP, oAngle);
1188 fQA->Fill("h_cut_L_rdp_v_mp", 0, iP, rP);
c04c80e6 1189 }
3a72645a 1190 if(fMCEvent){
1191 if(iMass > cutMass[0] && iMass < cutMass[1]){
1192 if(4 == TMath::Abs(fCurrentV0id)){
1193 fQAmc->Fill("h_cut_L_Chi2_S", 0, iP, chi2ndf);
1194 fQAmc->Fill("h_cut_L_Chi2_S", 1, iP, chi2ndf);
1195 fQAmc->Fill("h_cut_L_CosPoint_S", 0, iP, cosPoint);
1196 fQAmc->Fill("h_cut_L_DCA_S", 0, iP, dca);
1197 fQAmc->Fill("h_cut_L_VtxR_S", 0, iP, r);
1198 fQAmc->Fill("h_cut_L_OA_S", 0, iP, oAngle);
1199 fQAmc->Fill("h_cut_L_rdp_v_mp_S", 0, iP, rP);
1200 fQAmc->Fill("h_ProtonL_P_S", 0, p[ixMC[0]]);
1201 fQAmc->Fill("h_PionL_P_S", 0, p[ixMC[1]]);
1202 }
1203 else if(-2 != fCurrentV0id){
1204 fQAmc->Fill("h_cut_L_Chi2_B", 0, iP, chi2ndf);
1205 fQAmc->Fill("h_cut_L_Chi2_B", 1, iP, chi2ndf);
1206 fQAmc->Fill("h_cut_L_CosPoint_B", 0, iP, cosPoint);
1207 fQAmc->Fill("h_cut_L_DCA_B", 0, iP, dca);
1208 fQAmc->Fill("h_cut_L_VtxR_B", 0, iP, r);
1209 fQAmc->Fill("h_cut_L_OA_B", 0, iP, oAngle);
1210 fQAmc->Fill("h_cut_L_rdp_v_mp_B", 0, iP, rP);
1211 fQAmc->Fill("h_ProtonL_P_B", 0, p[ixMC[0]]);
1212 fQAmc->Fill("h_PionL_P_B", 0, p[ixMC[1]]);
1213 }
1214 }
1215 }
1216 //
1217 // Chi2/NDF cut
1218 //
1219 if(chi2ndf > cutChi2NDF) return kFALSE;
c04c80e6 1220 (type == 0) ? fQA->Fill("h_L_Mass", 1, iMass) : fQA->Fill("h_AL_Mass", 1, iMass);
1221 if(iMass > cutMass[0] && iMass < cutMass[1]){
3a72645a 1222 fQA->Fill("h_cut_L_CosPoint", 1, cosPoint);
c04c80e6 1223 fQA->Fill("h_ProtonL_P", 1, p[ix[0]]);
1224 fQA->Fill("h_PionL_P", 1, p[ix[1]]);
c04c80e6 1225 }
3a72645a 1226 if(fMCEvent){
1227 if(4 == fCurrentV0id) fQAmc->Fill("h_L_Mass_S", 1, iMass);
1228 else if(-4 == fCurrentV0id) fQAmc->Fill("h_AL_Mass_S", 1, iMass);
1229 else if(-2 != fCurrentV0id) fQAmc->Fill("h_LAL_Mass_B", 1, iMass);
1230 if(iMass > cutMass[0] && iMass < cutMass[1]){
1231 if(4 == TMath::Abs(fCurrentV0id)){
1232 fQAmc->Fill("h_cut_L_CosPoint_S", 1, iP, cosPoint);
1233 fQAmc->Fill("h_ProtonL_P_S", 1, p[ixMC[0]]);
1234 fQAmc->Fill("h_PionL_P_S", 1, p[ixMC[1]]);
1235 }
1236 else if(-2 != fCurrentV0id){
1237 fQAmc->Fill("h_cut_L_CosPoint_B", 1, iP, cosPoint);
1238 fQAmc->Fill("h_ProtonL_P_B", 1, p[ixMC[0]]);
1239 fQAmc->Fill("h_PionL_P_B", 1, p[ixMC[1]]);
1240 }
1241 }
1242 }
1243
1244 //
1245 // Cos point cut
1246 //
1247 if(cosPoint < cutCosPoint[0] || cosPoint > cutCosPoint[1]) return kFALSE;
c04c80e6 1248 (type == 0) ? fQA->Fill("h_L_Mass", 2, iMass) : fQA->Fill("h_AL_Mass", 2, iMass);
1249 if(iMass > cutMass[0] && iMass < cutMass[1]){
1250 fQA->Fill("h_ProtonL_P", 2, p[ix[0]]);
1251 fQA->Fill("h_PionL_P", 2, p[ix[1]]);
3a72645a 1252 fQA->Fill("h_cut_L_DCA", 1, dca);
1253 }
1254 if(fMCEvent){
1255 if(4 == fCurrentV0id) fQAmc->Fill("h_L_Mass_S", 2, iMass);
1256 else if(-4 == fCurrentV0id) fQAmc->Fill("h_AL_Mass_S", 2, iMass);
1257 else if(-2 != fCurrentV0id) fQAmc->Fill("h_LAL_Mass_B", 2, iMass);
1258 if(iMass > cutMass[0] && iMass < cutMass[1]){
1259 if(4 == TMath::Abs(fCurrentV0id)){
1260 fQAmc->Fill("h_cut_L_DCA_S", 1, iP, dca);
1261 fQAmc->Fill("h_ProtonL_P_S", 2, p[ixMC[0]]);
1262 fQAmc->Fill("h_PionL_P_S", 2, p[ixMC[1]]);
1263 }
1264 else if(-2 != fCurrentV0id){
1265 fQAmc->Fill("h_cut_L_DCA_B", 1, iP, dca);
1266 fQAmc->Fill("h_ProtonL_P_B", 2, p[ixMC[0]]);
1267 fQAmc->Fill("h_PionL_P_B", 2, p[ixMC[1]]);
1268 }
1269 }
c04c80e6 1270 }
1271
3a72645a 1272 //
1273 // DCA cut
1274 //
1275 if(dca < cutDCA[0] || dca > cutDCA[1]) return kFALSE;
c04c80e6 1276 (type == 0) ? fQA->Fill("h_L_Mass", 3, iMass) : fQA->Fill("h_AL_Mass", 3, iMass);
1277 if(iMass > cutMass[0] && iMass < cutMass[1]){
1278 fQA->Fill("h_ProtonL_P", 3, p[ix[0]]);
1279 fQA->Fill("h_PionL_P", 3, p[ix[1]]);
3a72645a 1280 fQA->Fill("h_cut_L_VtxR", 1, r);
1281 }
1282 if(fMCEvent){
1283 if(4 == fCurrentV0id) fQAmc->Fill("h_L_Mass_S", 3, iMass);
1284 else if(-4 == fCurrentV0id) fQAmc->Fill("h_AL_Mass_S", 3, iMass);
1285 else if(-2 != fCurrentV0id) fQAmc->Fill("h_LAL_Mass_B", 3, iMass);
1286 if(iMass > cutMass[0] && iMass < cutMass[1]){
1287 if(4 == TMath::Abs(fCurrentV0id)){
1288 fQAmc->Fill("h_cut_L_VtxR_S", 1, iP, r);
1289 fQAmc->Fill("h_ProtonL_P_S", 3, p[ixMC[0]]);
1290 fQAmc->Fill("h_PionL_P_S", 3, p[ixMC[1]]);
1291 }
1292 else if(-2 != fCurrentV0id){
1293 fQAmc->Fill("h_cut_L_VtxR_B", 1, iP, r);
1294 fQAmc->Fill("h_ProtonL_P_B", 3, p[ixMC[0]]);
1295 fQAmc->Fill("h_PionL_P_B", 3, p[ixMC[1]]);
1296 }
1297 }
c04c80e6 1298 }
1299
3a72645a 1300 //
1301 // Vertex radius cut
1302 //
1303 if(r < cutProdVtxR[0] || r > cutProdVtxR[1]) return kFALSE;
c04c80e6 1304 (type == 0) ? fQA->Fill("h_L_Mass", 4, iMass) : fQA->Fill("h_AL_Mass", 4, iMass);
1305 if(iMass > cutMass[0] && iMass < cutMass[1]){
3a72645a 1306 fQA->Fill("h_cut_L_OAvP", 1, iP, oAngle);
c04c80e6 1307 fQA->Fill("h_ProtonL_P", 4, p[ix[0]]);
1308 fQA->Fill("h_PionL_P", 4, p[ix[1]]);
c04c80e6 1309 }
3a72645a 1310 if(fMCEvent){
1311 if(4 == fCurrentV0id) fQAmc->Fill("h_L_Mass_S", 4, iMass);
1312 else if(-4 == fCurrentV0id) fQAmc->Fill("h_AL_Mass_S", 4, iMass);
1313 else if(-2 != fCurrentV0id) fQAmc->Fill("h_LAL_Mass_B", 4, iMass);
1314 if(iMass > cutMass[0] && iMass < cutMass[1]){
1315 if(4 == TMath::Abs(fCurrentV0id)){
1316 fQAmc->Fill("h_cut_L_OA_S", 1, iP, oAngle);
1317 fQAmc->Fill("h_ProtonL_P_S", 4, p[ixMC[0]]);
1318 fQAmc->Fill("h_PionL_P_S", 4, p[ixMC[1]]);
1319 }
1320 else if(-2 != fCurrentV0id){
1321 fQAmc->Fill("h_cut_L_OA_B", 1, iP, oAngle);
1322 fQAmc->Fill("h_ProtonL_P_B", 4, p[ixMC[0]]);
1323 fQAmc->Fill("h_PionL_P_B", 4, p[ixMC[1]]);
1324 }
1325 }
c04c80e6 1326 }
1327
1328 if(iMass < cutMass[0] || iMass > cutMass[1]) {
1329 return kFALSE;
1330 }
1331
1332 // all cuts passed
1333
1334 // assign the lambda type value: Lambda: kTRUE, Anti-Lambda: kFALSE
1335 isLambda = (0 == type) ? kTRUE : kFALSE;
1336
3a72645a 1337 // some MC stuff
1338 if(4 == fCurrentV0id){
1339 fQAmc->Fill("h_lambda_p_S", iP);
1340 }
1341 else if(-4 == fCurrentV0id){
1342 fQAmc->Fill("h_alambda_p_S", iP);
1343 }
1344 else if (-2 != fCurrentV0id && 0 == type){
1345 fQAmc->Fill("h_lambda_p_B", iP);
1346 }
1347 else if(-2 != fCurrentV0id && 0 != type ){
1348 fQAmc->Fill("h_alambda_p_B", iP);
1349 }
1350 //
1351 if(4 == TMath::Abs(fCurrentV0id)){
1352 fQAmc->Fill("h_ProtonL_P_S", 5, p[ixMC[0]]);
1353 fQAmc->Fill("h_PionL_P_S", 5, p[ixMC[1]]);
1354 }
1355 else if(-2 != fCurrentV0id){
1356 fQAmc->Fill("h_ProtonL_P_B", 5, p[ixMC[0]]);
1357 fQAmc->Fill("h_PionL_P_B", 5, p[ixMC[1]]);
1358 }
1359
c04c80e6 1360 return kTRUE;
1361}
1362//________________________________________________________________
1363Double_t AliHFEV0cuts::OpenAngle(AliESDv0 *v0) const {
1364 //
1365 // Opening angle between two daughter tracks
1366 //
1367 Double_t mn[3] = {0,0,0};
1368 Double_t mp[3] = {0,0,0};
1369
1370
1371 v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
1372 v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
1373
1374
1375 Double_t openAngle = TMath::ACos((mp[0]*mn[0] + mp[1]*mn[1] + mp[2]*mn[2])/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1] + mp[2]*mp[2])*TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1] + mn[2]*mn[2])));
1376
1377 return TMath::Abs(openAngle);
1378}
1379//________________________________________________________________
1380Double_t AliHFEV0cuts::PsiPair(AliESDv0 *v0) {
1381 //
1382 // Angle between daughter momentum plane and plane
1383 //
1384
1385 if(!fInputEvent) return -1.;
1386
1387 Float_t magField = fInputEvent->GetMagneticField();
1388
1389 Int_t pIndex = -1;
1390 Int_t nIndex = -1;
1391 if(CheckSigns(v0)){
1392 pIndex = v0->GetPindex();
1393 nIndex = v0->GetNindex();
1394 }
1395 else{
1396 pIndex = v0->GetNindex();
1397 nIndex = v0->GetPindex();
1398 }
1399
1400
1401 AliESDtrack* daughter[2];
1402
1403 daughter[0] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(pIndex));
1404 daughter[1] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(nIndex));
1405
1406 Double_t x, y, z;
1407 v0->GetXYZ(x,y,z);//Reconstructed coordinates of V0; to be replaced by Markus Rammler's method in case of conversions!
1408
1409 Double_t mn[3] = {0,0,0};
1410 Double_t mp[3] = {0,0,0};
1411
1412
1413 v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
1414 v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
1415
1416
1417 Double_t deltat = 1.;
1418 deltat = TMath::ATan(mp[2]/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1])+1.e-13)) - TMath::ATan(mn[2]/(TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1])+1.e-13));//difference of angles of the two daughter tracks with z-axis
1419
1420 Double_t radiussum = TMath::Sqrt(x*x + y*y) + 50;//radius to which tracks shall be propagated
1421
1422 Double_t momPosProp[3];
1423 Double_t momNegProp[3];
1424
1425 AliExternalTrackParam pt(*daughter[0]), nt(*daughter[1]);
1426
1427 Double_t psiPair = 4.;
1428
1429 if(nt.PropagateTo(radiussum,magField) == 0)//propagate tracks to the outside
1430 psiPair = -5.;
1431 if(pt.PropagateTo(radiussum,magField) == 0)
1432 psiPair = -5.;
1433 pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation
1434 nt.GetPxPyPz(momNegProp);
1435
1436 Double_t pEle =
1437 TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
1438 Double_t pPos =
1439 TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
1440
1441 Double_t scalarproduct =
1442 momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta
1443
1444 Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks
1445
1446 psiPair = TMath::Abs(TMath::ASin(deltat/chipair));
1447
1448 return psiPair;
1449}
1450//________________________________________________________________
1451AliKFParticle *AliHFEV0cuts::CreateMotherParticle(AliVTrack* const pdaughter, AliVTrack* const ndaughter, Int_t pspec, Int_t nspec){
1452 //
1453 // Creates a mother particle
1454 //
1455 AliKFParticle pkfdaughter(*pdaughter, pspec);
1456 AliKFParticle nkfdaughter(*ndaughter, nspec);
1457
1458 // - check if the daughter particles are coming from the primary vertex
1459 // - check the number of tracks that take part in the creaton of primary vertex.
1460 // important: after removeal of candidate tracks there must be at least 2 tracks left
1461 // otherwise the primary vertex will be corrupted
1462
1463 // ESD Analyis
1464 //const AliESDVertex *esdvertex = dynamic_cast<const AliESDVertex *>(fInputEvent->GetPrimaryVertex());
1465 //if(!esdvertex) return NULL;
1466 //UShort_t *contrib = esdvertex->GetIndices();
1467
1468 //
1469 // not using the removal of the daughter track now
1470 //
1471 // Int_t nTracks = esdvertex->GetNIndices();
1472 // printf(" -D: N Vertex tracks: %i\n", nTracks);
1473 // printf(" -D: N Contributors: %i\n", fPrimaryVertex->GetNContributors());
1474 // Int_t nfound = 0;
1475 // for(Int_t id = 0; id < esdvertex->GetNIndices(); id++){
1476 // if(contrib[id] == pdaughter->GetID()){
1477 // if( (nTracks - nfound) <= 2 ) return NULL;
1478 // *fPrimaryVertex -= pkfdaughter;
1479 // removed[0] = kTRUE;
1480 // nfound++;
1481 // }
1482 // if(contrib[id] == ndaughter->GetID()){
1483 // if( (nTracks - nfound) <=2 ) return NULL;
1484 // *fPrimaryVertex -= nkfdaughter;
1485 // removed[1] = kTRUE;
1486 // nfound++;
1487 // }
1488 // if(nfound == 2) break;
1489 // }
1490
1491 // printf(" -D: n removed: %i\n", nfound);
1492
1493 // Create the mother particle
1494 AliKFParticle *m = new AliKFParticle(pkfdaughter, nkfdaughter);
3a72645a 1495 // DEBUG - testing
1496 if(TMath::Abs(kElectron) == pspec && TMath::Abs(kElectron) == nspec) m->SetMassConstraint(0, 0.001);
1497 else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(), 0.);
1498 else if(TMath::Abs(kProton) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
1499 else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kProton) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
1500 else{
1501 AliError("Wrong daughter ID - mass constraint can not be set");
1502 }
c04c80e6 1503
1504 AliKFVertex improvedVertex = *fPrimaryVertex;
1505 improvedVertex += *m;
1506 m->SetProductionVertex(improvedVertex);
1507
1508 // update 15/06/2010
1509 // mother particle will not be added to primary vertex but only to its copy
1510 // as this confilcts with calling
1511 // m->SetPrimaryVertex() function and
1512 // subsequently removing the mother particle afterwards
1513 // Sourse: Sergey Gorbunov
1514
1515 return m;
1516}
1517//_________________________________________________
1518Bool_t AliHFEV0cuts::LooseRejectK0(AliESDv0 * const v0) const {
1519 //
1520 // Reject K0 based on loose cuts
1521 //
1522 Double_t mass = v0->GetEffMass(AliPID::kPion, AliPID::kPion);
1523 if(mass > 0.494 && mass < 0.501) return kTRUE;
1524 return kFALSE;
1525}
1526
1527//_________________________________________________
1528Bool_t AliHFEV0cuts::LooseRejectLambda(AliESDv0 * const v0) const {
1529 //
1530 // Reject Lambda based on loose cuts
1531 //
1532 Double_t mass1 = v0->GetEffMass(AliPID::kPion, AliPID::kProton);
1533 Double_t mass2 = v0->GetEffMass(AliPID::kProton, AliPID::kPion);
1534
1535 if(mass1 > 1.1 && mass1 < 1.12) return kTRUE;
1536 if(mass2 > 1.1 && mass2 < 1.12) return kTRUE;
1537 return kFALSE;
1538}
1539
1540//_________________________________________________
1541Bool_t AliHFEV0cuts::LooseRejectGamma(AliESDv0 * const v0) const {
1542 //
1543 // Reject Lambda based on loose cuts
1544 //
1545
1546 Double_t mass = v0->GetEffMass(AliPID::kElectron, AliPID::kElectron);
1547
1548 if(mass < 0.02) return kTRUE;
1549 return kFALSE;
1550}
1551//___________________________________________________________________
1552void AliHFEV0cuts::Armenteros(AliESDv0 *v0, Float_t val[2]){
1553 //
1554 // computes the Armenteros variables for given V0
1555 // fills the histogram
1556 // returns the values via "val"
1557 //
1558
1559 Double_t mn[3] = {0,0,0};
1560 Double_t mp[3] = {0,0,0};
1561 Double_t mm[3] = {0,0,0};
1562
1563 if(CheckSigns(v0)){
1564 v0->GetNPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
1565 v0->GetPPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
1566 }
1567 else{
1568 v0->GetPPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
1569 v0->GetNPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
1570 }
1571 v0->GetPxPyPz(mm[0],mm[1],mm[2]); //reconstructed cartesian momentum components of mother
1572
1573 TVector3 vecN(mn[0],mn[1],mn[2]);
1574 TVector3 vecP(mp[0],mp[1],mp[2]);
1575 TVector3 vecM(mm[0],mm[1],mm[2]);
1576
1577 Double_t thetaP = acos((vecP * vecM)/(vecP.Mag() * vecM.Mag()));
1578 Double_t thetaN = acos((vecN * vecM)/(vecN.Mag() * vecM.Mag()));
1579
1580 Double_t alfa = ((vecP.Mag())*cos(thetaP)-(vecN.Mag())*cos(thetaN))/
1581 ((vecP.Mag())*cos(thetaP)+(vecN.Mag())*cos(thetaN)) ;
1582 Double_t qt = vecP.Mag()*sin(thetaP);
1583
1584 val[0] = alfa;
1585 val[1] = qt;
1586
1587}
1588//___________________________________________________________________
1589Bool_t AliHFEV0cuts::CheckSigns(AliESDv0* const v0){
1590 //
1591 // check wheter the sign was correctly applied to
1592 // V0 daughter tracks
1593 //
1594
1595 Bool_t correct = kFALSE;
1596
1597 Int_t pIndex = 0, nIndex = 0;
1598 pIndex = v0->GetPindex();
1599 nIndex = v0->GetNindex();
1600
1601 AliESDtrack* d[2];
1602 d[0] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(pIndex));
1603 d[1] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(nIndex));
1604
1605 Int_t sign[2];
1606 sign[0] = (int)d[0]->GetSign();
1607 sign[1] = (int)d[1]->GetSign();
1608
1609 if(-1 == sign[0] && 1 == sign[1]){
1610 correct = kFALSE;
1611 //v0->SetIndex(0, pIndex); // set the index of the negative v0 track
1612 //v0->SetIndex(1, nIndex); // set the index of the positive v0 track
1613 }
1614 else{
1615 correct = kTRUE;
1616 }
1617
1618 //pIndex = v0->GetPindex();
1619 //nIndex = v0->GetNindex();
1620 //printf("-D2: P: %i, N: %i\n", pIndex, nIndex);
1621
1622 return correct;
1623}
3a72645a 1624//___________________________________________________________________
1625Bool_t AliHFEV0cuts::GetConvPosXY(AliESDtrack * const ptrack, AliESDtrack * const ntrack, Double_t convpos[2]){
1626 //
1627 // recalculate the gamma conversion XY postition
1628 //
1629
1630 const Double_t b = fInputEvent->GetMagneticField();
1631
1632 Double_t helixcenterpos[2];
1633 GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos);
1634
1635 Double_t helixcenterneg[2];
1636 GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg);
1637
1638 Double_t poshelix[6];
1639 ptrack->GetHelixParameters(poshelix,b);
1640 Double_t posradius = TMath::Abs(1./poshelix[4]);
1641
1642 Double_t neghelix[6];
1643 ntrack->GetHelixParameters(neghelix,b);
1644 Double_t negradius = TMath::Abs(1./neghelix[4]);
1645
1646 Double_t xpos = helixcenterpos[0];
1647 Double_t ypos = helixcenterpos[1];
1648 Double_t xneg = helixcenterneg[0];
1649 Double_t yneg = helixcenterneg[1];
1650
1651 convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
1652 convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius);
1653
1654 return 1;
1655}
1656//___________________________________________________________________
1657Bool_t AliHFEV0cuts::GetHelixCenter(AliESDtrack * const track, Double_t b,Int_t charge, Double_t center[2]){
1658 // see header file for documentation
1659
1660 Double_t pi = TMath::Pi();
1661
1662 Double_t helix[6];
1663 track->GetHelixParameters(helix,b);
1664
1665 Double_t xpos = helix[5];
1666 Double_t ypos = helix[0];
1667 Double_t radius = TMath::Abs(1./helix[4]);
1668 Double_t phi = helix[2];
1669
1670 if(phi < 0){
1671 phi = phi + 2*pi;
1672 }
1673
1674 phi -= pi/2.;
1675 Double_t xpoint = radius * TMath::Cos(phi);
1676 Double_t ypoint = radius * TMath::Sin(phi);
1677
1678 if(b<0){
1679 if(charge > 0){
1680 xpoint = - xpoint;
1681 ypoint = - ypoint;
1682 }
1683
1684 if(charge < 0){
1685 xpoint = xpoint;
1686 ypoint = ypoint;
1687 }
1688 }
1689 if(b>0){
1690 if(charge > 0){
1691 xpoint = xpoint;
1692 ypoint = ypoint;
1693 }
1694
1695 if(charge < 0){
1696 xpoint = - xpoint;
1697 ypoint = - ypoint;
1698 }
1699 }
1700 center[0] = xpos + xpoint;
1701 center[1] = ypos + ypoint;
1702
1703 return 1;
1704}