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