]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/vertexingHF/AliAODRecoCascadeHF.cxx
Updates for proper treatment of fQuality flag in centrality selection and for usage...
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAODRecoCascadeHF.cxx
CommitLineData
05d80dd6 1/**************************************************************************
2 * Copyright(c) 1998-2008, 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
27de2dfb 16/* $Id$ */
17
05d80dd6 18/////////////////////////////////////////////////////////////
19//
20// Class for AOD reconstructed heavy-flavour cascades
21//
c3fe1eed 22// Author: X-M. Zhang, zhangxm@iopp.ccnu.edu.cn
05d80dd6 23/////////////////////////////////////////////////////////////
24
25#include <TVector3.h>
26#include <TDatabasePDG.h>
b3999999 27#include <TClonesArray.h>
28#include "AliAODMCParticle.h"
05d80dd6 29#include "AliAODRecoDecay.h"
30#include "AliAODVertex.h"
31#include "AliAODRecoDecayHF2Prong.h"
32#include "AliAODRecoCascadeHF.h"
33
34ClassImp(AliAODRecoCascadeHF)
35//-----------------------------------------------------------------------------
36
37AliAODRecoCascadeHF::AliAODRecoCascadeHF() :
dcfa35b3 38 AliAODRecoDecayHF2Prong()
05d80dd6 39{
40 //
41 // Default Constructor
42 //
43}
44//-----------------------------------------------------------------------------
45AliAODRecoCascadeHF::AliAODRecoCascadeHF(AliAODVertex *vtx2, Short_t charge,
46 Double_t *px, Double_t *py, Double_t *pz,
47 Double_t *d0, Double_t *d0err, Double_t dca) :
dcfa35b3 48 AliAODRecoDecayHF2Prong(vtx2, px, py, pz, d0, d0err, dca)
05d80dd6 49{
50 //
51 // Constructor with AliAODVertex for decay vertex
52 //
53 SetCharge(charge);
54}
55//-----------------------------------------------------------------------------
56AliAODRecoCascadeHF::AliAODRecoCascadeHF(AliAODVertex *vtx2, Short_t charge,
57 Double_t *d0, Double_t *d0err, Double_t dca) :
dcfa35b3 58 AliAODRecoDecayHF2Prong(vtx2, d0, d0err, dca)
05d80dd6 59{
60 //
61 // Constructor with decay vertex and without prongs momenta
62 //
63 SetCharge(charge);
64}
65//-----------------------------------------------------------------------------
66AliAODRecoCascadeHF::AliAODRecoCascadeHF(const AliAODRecoCascadeHF &source) :
dcfa35b3 67 AliAODRecoDecayHF2Prong(source)
05d80dd6 68{
69 //
70 // Copy constructor
71 //
72}
73//-----------------------------------------------------------------------------
74AliAODRecoCascadeHF &AliAODRecoCascadeHF::operator=(const AliAODRecoCascadeHF &source)
75{
76 //
77 // assignment operator
78 //
79 if(&source == this) return *this;
80
81 AliAODRecoDecayHF2Prong::operator=(source);
82
05d80dd6 83 return *this;
84}
85//-----------------------------------------------------------------------------
86AliAODRecoCascadeHF::~AliAODRecoCascadeHF()
87{
88 //
89 // Default Destructor
90 //
91}
92//-----------------------------------------------------------------------------
93Double_t AliAODRecoCascadeHF::InvMassDstarKpipi() const
94{
95 //
96 // 3 prong invariant mass of the D0 daughters and the soft pion
97 //
78489ccc 98 Double_t e[3];
99 if (Charge()>0){
100 e[0]=Get2Prong()->EProng(0,211);
101 e[1]=Get2Prong()->EProng(1,321);
102 }else{
103 e[0]=Get2Prong()->EProng(0,321);
104 e[1]=Get2Prong()->EProng(1,211);
105 }
106 e[2]=EProng(0,211);
05d80dd6 107
78489ccc 108 Double_t esum = e[0]+e[1]+e[2];
109 Double_t minv = TMath::Sqrt(esum*esum-P()*P());
05d80dd6 110
78489ccc 111 return minv;
05d80dd6 112}
b3999999 113//----------------------------------------------------------------------------
114Int_t AliAODRecoCascadeHF::MatchToMC(Int_t pdgabs,Int_t pdgabs2prong,
383937f0 115 Int_t *pdgDg,Int_t *pdgDg2prong,
b3999999 116 TClonesArray *mcArray) const
117{
118 //
119 // Check if this candidate is matched to a MC signal
120 // If no, return -1
121 // If yes, return label (>=0) of the AliAODMCParticle
122 //
123
8862e6bf 124 Int_t ndg=GetNDaughters();
125 if(!ndg) {
b3999999 126 AliError("No daughters available");
127 return -1;
128 }
129
130 AliAODRecoDecayHF2Prong *the2Prong = Get2Prong();
383937f0 131 Int_t lab2Prong = the2Prong->MatchToMC(pdgabs2prong,mcArray,2,pdgDg2prong);
b3999999 132
133 if(lab2Prong<0) return -1;
134
e11ae259 135 Int_t dgLabels[10]={0,0,0,0,0,0,0,0,0,0};
b3999999 136
137 // loop on daughters and write labels
8862e6bf 138 for(Int_t i=0; i<ndg; i++) {
b3999999 139 AliVTrack *trk = (AliVTrack*)GetDaughter(i);
140 Int_t lab = trk->GetLabel();
141 if(lab==-1) { // this daughter is the 2prong
142 lab=lab2Prong;
143 } else if(lab<-1) {
144 printf("daughter with negative label\n");
145 continue;
146 }
147 dgLabels[i] = lab;
148 }
149
383937f0 150 return AliAODRecoDecay::MatchToMC(pdgabs,mcArray,dgLabels,2,2,pdgDg);
b3999999 151}
05d80dd6 152//-----------------------------------------------------------------------------
153Bool_t AliAODRecoCascadeHF::SelectDstar(const Double_t *cutsDstar,
154 const Double_t *cutsD0,
155 Bool_t testD0) const
156{
157 //
158 // cutsDstar[0] = inv. mass half width of D* [GeV]
159 // cutsDstar[1] = half width of (M_Kpipi-M_D0) [GeV]
160 // cutsDstar[2] = PtMin of pi_s [GeV/c]
161 // cutsDstar[3] = PtMax of pi_s [GeV/c]
162 // cutsDstar[4] = theta, angle between the pi_s and decay plane of the D0 [rad]
163 //
164 // cutsD0[0] = inv. mass half width [GeV]
165 // cutsD0[1] = dca [cm]
166 // cutsD0[2] = cosThetaStar
167 // cutsD0[3] = pTK [GeV/c]
168 // cutsD0[4] = pTPi [GeV/c]
169 // cutsD0[5] = d0K [cm] upper limit!
170 // cutsD0[6] = d0Pi [cm] upper limit!
171 // cutsD0[7] = d0d0 [cm^2]
172 // cutsD0[8] = cosThetaPoint
173
174
175 // check that the D0 passes the cuts
176 // (if we have a D*+, it has to pass as D0,
177 // if we have a D*-, it has to pass as D0bar)
178
179 if(testD0) {
180 Int_t okD0=0,okD0bar=0;
181 Get2Prong()->SelectD0(cutsD0,okD0,okD0bar);
182 if((Charge()==+1 && !okD0) || (Charge()==-1 && !okD0bar)) return kFALSE;
183 }
184
185 if( (PtProng(0)<cutsDstar[2]) || (PtProng(0)>cutsDstar[3]) ) return kFALSE;
186
187 Double_t mDstar = TDatabasePDG::Instance()->GetParticle(413)->Mass();
188 Double_t invmDstar = InvMassDstarKpipi();
189 if(TMath::Abs(mDstar-invmDstar)>cutsDstar[0]) return kFALSE;
190
191 Double_t mD0 = TDatabasePDG::Instance()->GetParticle(421)->Mass();
192 if(TMath::Abs((mDstar-mD0)-DeltaInvMass())>cutsDstar[1]) return kFALSE;
193
3f6adfc9 194 Double_t theta = AngleD0dkpPisoft();
05d80dd6 195 if(theta>cutsDstar[4]) return kFALSE;
05d80dd6 196
197 return kTRUE;
198}
a07ad8e0 199//-----------------------------------------------------------------------------
200Bool_t AliAODRecoCascadeHF::SelectLctoV0(const Double_t *cutsLctoV0,
201 Bool_t okLck0sp, Bool_t okLcLpi) const
202{
203 // cuts on Lambdac candidates to V0+bachelor
204 // (to be passed to AliAODRecoDecayHF3Prong::SelectLctoV0())
205 // 0 = inv. mass half width in K0s hypothesis [GeV]
206 // 1 = inv. mass half width in Lambda hypothesis [GeV]
207 // 2 = inv. mass V0 in K0s hypothesis half width [GeV]
208 // 3 = inv. mass V0 in Lambda hypothesis half width [GeV]
209 // 4 = pT min Bachelor track [GeV/c]
210 // 5 = pT min V0-Positive track [GeV/c]
211 // 6 = pT min V0-Negative track [GeV/c]
212 // 7 = dca cut on the V0 (cm)
213 // 8 = dca cut on the cascade (cm)
214
215// if ( !Getv0() || !Getv0PositiveTrack() || !Getv0NegativeTrack() )
216// { AliInfo(Form("Not adapted for ESDv0s, return true...")); return false; }
217
218 Double_t mLck0sp,mLcLpi;
219 okLck0sp=1; okLcLpi=1;
220
221 Double_t mLcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
222 Double_t mk0sPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
223 Double_t mLPDG = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
224
225 // k0s + p
226 double mk0s = Getv0()->MassK0Short();
227 mLck0sp = InvMassLctoK0sP();
228
229 // lambda + pi
230 double mlambda = Getv0()->MassLambda();
231 double malambda = Getv0()->MassAntiLambda();
232 mLcLpi = InvMassLctoLambdaPi();
233
234 // cut on Lc mass
235 // with k0s p hypothesis
236 if(TMath::Abs(mLck0sp-mLcPDG)>cutsLctoV0[0]) okLck0sp = 0;
237 // with Lambda pi hypothesis
238 if(TMath::Abs(mLcLpi-mLcPDG)>cutsLctoV0[1]) okLcLpi = 0;
239
240 // cuts on the v0 mass
241 if(TMath::Abs(mk0s-mk0sPDG)>cutsLctoV0[2]) okLck0sp = 0;
242 if( TMath::Abs(mlambda-mLPDG)>cutsLctoV0[3] &&
243 TMath::Abs(malambda-mLPDG)>cutsLctoV0[3] ) okLcLpi = 0;
244
245 if(!okLck0sp && !okLcLpi) return 0;
246
247 // cuts on the minimum pt of the tracks
248 if(TMath::Abs(GetBachelor()->Pt()) < cutsLctoV0[4]) return 0;
249 if(TMath::Abs(Getv0PositiveTrack()->Pt()) < cutsLctoV0[5]) return 0;
250 if(TMath::Abs(Getv0NegativeTrack()->Pt()) < cutsLctoV0[6]) return 0;
251
252 // cut on the v0 dca
253 if(TMath::Abs(Getv0()->DcaV0Daughters()) > cutsLctoV0[7]) return 0;
254
255 // cut on the cascade dca
256 if( TMath::Abs(GetDCA(0))>cutsLctoV0[8] ||
257 TMath::Abs(Getv0()->DcaPosToPrimVertex())>cutsLctoV0[8] ||
258 TMath::Abs(Getv0()->DcaNegToPrimVertex())>cutsLctoV0[8] ) return 0;
259
260 return true;
261
262}
263//-----------------------------------------------------------------------------
3f6adfc9 264Double_t AliAODRecoCascadeHF::AngleD0dkpPisoft() const {
265 //
266 // Angle of soft pion to D0 decay plane
267 //
268
269 TVector3 p3Trk0(Get2Prong()->PxProng(0),Get2Prong()->PyProng(0),Get2Prong()->PzProng(0)); // from D0
270 TVector3 p3Trk1(Get2Prong()->PxProng(1),Get2Prong()->PyProng(1),Get2Prong()->PzProng(1)); // from D0
271 TVector3 p3Trk2(PxProng(0),PyProng(0),PzProng(0)); // pi_s
272
273 TVector3 perp = p3Trk0.Cross(p3Trk1);
274 Double_t theta = p3Trk2.Angle(perp);
275 if(theta>(TMath::Pi()-theta)) theta = TMath::Pi() - theta;
276 theta = TMath::Pi()/2. - theta;
277
278 return theta;
279}
ccfae9f3 280//-----------------------------------------------------------------------------
281Bool_t AliAODRecoCascadeHF::TrigonometricalCut() const {
282 //
283 // Trigonometrical constraint
284 //
285 TVector3 p3Trk0(Get2Prong()->PxProng(0),Get2Prong()->PyProng(0),Get2Prong()->PzProng(0)); // from D0
286 TVector3 p3Trk1(Get2Prong()->PxProng(1),Get2Prong()->PyProng(1),Get2Prong()->PzProng(1)); // from D0
287 TVector3 p3Trk2(PxProng(0),PyProng(0),PzProng(0)); // pi_s
288
289 Double_t alpha = p3Trk0.Angle(p3Trk2);
290 Double_t beta = p3Trk1.Angle(p3Trk2);
291
292 Double_t cosphi01 = TMath::Cos(alpha) / TMath::Cos(AngleD0dkpPisoft());
293 Double_t cosphi02 = TMath::Cos(beta) / TMath::Cos(AngleD0dkpPisoft());
294
295 Double_t phi01 = TMath::ACos(cosphi01);
296 Double_t phi02 = TMath::ACos(cosphi02);
297 Double_t phi00 = p3Trk0.Angle(p3Trk1);
298
299 if((phi01>phi00) || (phi02>phi00)) return kFALSE;
300 return kTRUE;
301}