git://git.uio.no / u/mrichter/AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
bug fix for new architecture
[u/mrichter/AliRoot.git] / STEER / AliAODRecoDecay.cxx
CommitLineData
7de7497b 1/**************************************************************************
2 * Copyright(c) 1998-2006, 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
16/////////////////////////////////////////////////////////////
17//
18// Base class for AOD reconstructed decay
19//
20// Author: A.Dainese, andrea.dainese@lnl.infn.it
21/////////////////////////////////////////////////////////////
22
23#include <TDatabasePDG.h>
24#include <TVector3.h>
bfd68286 25#include <TClonesArray.h>
e045ffda 26#include "AliLog.h"
27#include "AliVTrack.h"
bfd68286 28#include "AliAODMCParticle.h"
7de7497b 29#include "AliAODRecoDecay.h"
30
31ClassImp(AliAODRecoDecay)
32
33//--------------------------------------------------------------------------
34AliAODRecoDecay::AliAODRecoDecay() :
e045ffda 35 AliVTrack(),
7de7497b 36 fSecondaryVtx(0x0),
a11de4a0 37 fOwnSecondaryVtx(0x0),
7de7497b 38 fCharge(0),
39 fNProngs(0), fNDCA(0), fNPID(0),
40 fPx(0x0), fPy(0x0), fPz(0x0),
41 fd0(0x0),
42 fDCA(0x0),
43 fPID(0x0),
44 fEventNumber(-1),fRunNumber(-1)
45{
46 //
47 // Default Constructor
48 //
7de7497b 49}
50//--------------------------------------------------------------------------
51AliAODRecoDecay::AliAODRecoDecay(AliAODVertex *vtx2,Int_t nprongs,
52 Short_t charge,
53 Double_t *px,Double_t *py,Double_t *pz,
54 Double_t *d0) :
e045ffda 55 AliVTrack(),
7de7497b 56 fSecondaryVtx(vtx2),
a11de4a0 57 fOwnSecondaryVtx(0x0),
7de7497b 58 fCharge(charge),
59 fNProngs(nprongs), fNDCA(0), fNPID(0),
60 fPx(0x0), fPy(0x0), fPz(0x0),
61 fd0(0x0),
62 fDCA(0x0),
63 fPID(0x0),
64 fEventNumber(-1),fRunNumber(-1)
65{
66 //
67 // Constructor with AliAODVertex for decay vertex
68 //
69
70 fPx = new Double_t[GetNProngs()];
71 fPy = new Double_t[GetNProngs()];
72 fPz = new Double_t[GetNProngs()];
73 fd0 = new Double_t[GetNProngs()];
74 for(Int_t i=0; i<GetNProngs(); i++) {
75 fPx[i] = px[i];
76 fPy[i] = py[i];
77 fPz[i] = pz[i];
78 fd0[i] = d0[i];
79 }
80}
81//--------------------------------------------------------------------------
82AliAODRecoDecay::AliAODRecoDecay(AliAODVertex *vtx2,Int_t nprongs,
83 Short_t charge,
84 Double_t *d0) :
e045ffda 85 AliVTrack(),
7de7497b 86 fSecondaryVtx(vtx2),
a11de4a0 87 fOwnSecondaryVtx(0x0),
7de7497b 88 fCharge(charge),
89 fNProngs(nprongs), fNDCA(0), fNPID(0),
90 fPx(0x0), fPy(0x0), fPz(0x0),
91 fd0(0x0),
92 fDCA(0x0),
93 fPID(0x0),
94 fEventNumber(-1),fRunNumber(-1)
95{
96 //
97 // Constructor with AliAODVertex for decay vertex and without prongs momenta
98 //
99
100 fd0 = new Double_t[GetNProngs()];
101 for(Int_t i=0; i<GetNProngs(); i++) fd0[i] = d0[i];
102}
103//--------------------------------------------------------------------------
104AliAODRecoDecay::AliAODRecoDecay(const AliAODRecoDecay &source) :
e045ffda 105 AliVTrack(source),
7de7497b 106 fSecondaryVtx(source.fSecondaryVtx),
a11de4a0 107 fOwnSecondaryVtx(source.fOwnSecondaryVtx),
7de7497b 108 fCharge(source.fCharge),
109 fNProngs(source.fNProngs), fNDCA(source.fNDCA), fNPID(source.fNPID),
110 fPx(0x0), fPy(0x0), fPz(0x0),
111 fd0(0x0),
112 fDCA(0x0),
113 fPID(0x0),
114 fEventNumber(source.fEventNumber),fRunNumber(source.fRunNumber)
115{
116 //
117 // Copy constructor
118 //
119 if(source.GetNProngs()>0) {
58b0186f 120 fd0 = new Double32_t[GetNProngs()];
121 memcpy(fd0,source.fd0,GetNProngs()*sizeof(Double32_t));
7de7497b 122 if(source.fPx) {
58b0186f 123 fPx = new Double32_t[GetNProngs()];
124 fPy = new Double32_t[GetNProngs()];
125 fPz = new Double32_t[GetNProngs()];
126 memcpy(fPx,source.fPx,GetNProngs()*sizeof(Double32_t));
127 memcpy(fPy,source.fPy,GetNProngs()*sizeof(Double32_t));
128 memcpy(fPz,source.fPz,GetNProngs()*sizeof(Double32_t));
7de7497b 129 }
130 if(source.fPID) {
58b0186f 131 fPID = new Double32_t[fNPID];
132 memcpy(fPID,source.fPID,fNPID*sizeof(Double32_t));
7de7497b 133 }
134 if(source.fDCA) {
58b0186f 135 fDCA = new Double32_t[fNDCA];
136 memcpy(fDCA,source.fDCA,fNDCA*sizeof(Double32_t));
7de7497b 137 }
138 }
139}
140//--------------------------------------------------------------------------
141AliAODRecoDecay &AliAODRecoDecay::operator=(const AliAODRecoDecay &source)
142{
143 //
144 // assignment operator
145 //
146 if(&source == this) return *this;
147 fSecondaryVtx = source.fSecondaryVtx;
a11de4a0 148 fOwnSecondaryVtx = source.fOwnSecondaryVtx;
7de7497b 149 fCharge = source.fCharge;
150 fNProngs = source.fNProngs;
151 fNDCA = source.fNDCA;
152 fNPID = source.fNPID;
153 fEventNumber = source.fEventNumber;
154 fRunNumber = source.fRunNumber;
155 if(source.GetNProngs()>0) {
341952b8 156 if(fd0)delete [] fd0;
58b0186f 157 fd0 = new Double32_t[GetNProngs()];
158 memcpy(fd0,source.fd0,GetNProngs()*sizeof(Double32_t));
7de7497b 159 if(source.fPx) {
341952b8 160 if(fPx) delete [] fPx;
58b0186f 161 fPx = new Double32_t[GetNProngs()];
341952b8 162 if(fPy) delete [] fPy;
58b0186f 163 fPy = new Double32_t[GetNProngs()];
341952b8 164 if(fPz) delete [] fPz;
58b0186f 165 fPz = new Double32_t[GetNProngs()];
166 memcpy(fPx,source.fPx,GetNProngs()*sizeof(Double32_t));
167 memcpy(fPy,source.fPy,GetNProngs()*sizeof(Double32_t));
168 memcpy(fPz,source.fPz,GetNProngs()*sizeof(Double32_t));
7de7497b 169 }
170 if(source.fPID) {
341952b8 171 if(fPID) delete [] fPID;
58b0186f 172 fPID = new Double32_t[fNPID];
173 memcpy(fPID,source.fPID,fNPID*sizeof(Double32_t));
7de7497b 174 }
175 if(source.fDCA) {
341952b8 176 if(fDCA) delete [] fDCA;
58b0186f 177 fDCA = new Double32_t[fNDCA];
178 memcpy(fDCA,source.fDCA,fNDCA*sizeof(Double32_t));
7de7497b 179 }
180 }
181 return *this;
182}
183//--------------------------------------------------------------------------
184AliAODRecoDecay::~AliAODRecoDecay() {
185 //
186 // Default Destructor
187 //
a11de4a0 188 if(fPx) { delete [] fPx; fPx=NULL; }
189 if(fPy) { delete [] fPy; fPy=NULL; }
190 if(fPz) { delete [] fPz; fPz=NULL; }
191 if(fd0) { delete [] fd0; fd0=NULL; }
192 if(fPID) { delete [] fPID; fPID=NULL; }
193 if(fDCA) { delete [] fDCA; fDCA=NULL; }
194 if(fOwnSecondaryVtx) { delete fOwnSecondaryVtx; fOwnSecondaryVtx=NULL; }
7de7497b 195}
196//----------------------------------------------------------------------------
197Double_t AliAODRecoDecay::Alpha() const
198{
199 //
200 // Armenteros-Podolanski alpha for 2-prong decays
201 //
202 if(GetNProngs()!=2) {
203 printf("Can be called only for 2-prong decays");
204 return (Double_t)-99999.;
205 }
206 return 1.-2./(1.+QlProng(0)/QlProng(1));
207}
208//----------------------------------------------------------------------------
209Double_t AliAODRecoDecay::DecayLength(Double_t point[3]) const
210{
211 //
212 // Decay length assuming it is produced at "point" [cm]
213 //
214 return TMath::Sqrt((point[0]-GetSecVtxX())
215 *(point[0]-GetSecVtxX())
216 +(point[1]-GetSecVtxY())
217 *(point[1]-GetSecVtxY())
218 +(point[2]-GetSecVtxZ())
219 *(point[2]-GetSecVtxZ()));
220}
221//----------------------------------------------------------------------------
222Double_t AliAODRecoDecay::DecayLengthXY(Double_t point[3]) const
223{
224 //
225 // Decay length in XY assuming it is produced at "point" [cm]
226 //
227 return TMath::Sqrt((point[0]-GetSecVtxX())
228 *(point[0]-GetSecVtxX())
229 +(point[1]-GetSecVtxY())
230 *(point[1]-GetSecVtxY()));
231}
232//----------------------------------------------------------------------------
233Double_t AliAODRecoDecay::CosPointingAngle(Double_t point[3]) const
234{
235 //
236 // Cosine of pointing angle in space assuming it is produced at "point"
237 //
238 TVector3 mom(Px(),Py(),Pz());
239 TVector3 fline(GetSecVtxX()-point[0],
240 GetSecVtxY()-point[1],
241 GetSecVtxZ()-point[2]);
242
243 Double_t pta = mom.Angle(fline);
244
245 return TMath::Cos(pta);
246}
247//----------------------------------------------------------------------------
248Double_t AliAODRecoDecay::CosPointingAngleXY(Double_t point[3]) const
249{
250 //
251 // Cosine of pointing angle in transverse plane assuming it is produced
252 // at "point"
253 //
254 TVector3 momXY(Px(),Py(),0.);
255 TVector3 flineXY(GetSecVtxX()-point[0],
256 GetSecVtxY()-point[1],
257 0.);
258
259 Double_t ptaXY = momXY.Angle(flineXY);
260
261 return TMath::Cos(ptaXY);
262}
263//----------------------------------------------------------------------------
264Double_t AliAODRecoDecay::CosThetaStar(Int_t ip,UInt_t pdgvtx,UInt_t pdgprong0,UInt_t pdgprong1) const
265{
266 //
267 // Only for 2-prong decays:
268 // Cosine of decay angle (theta*) in the rest frame of the mother particle
269 // for prong ip (0 or 1) with mass hypotheses pdgvtx for mother particle,
270 // pdgprong0 for prong 0 and pdgprong1 for prong1
271 //
272 if(GetNProngs()!=2) {
273 printf("Can be called only for 2-prong decays");
274 return (Double_t)-99999.;
275 }
276 Double_t massvtx = TDatabasePDG::Instance()->GetParticle(pdgvtx)->Mass();
277 Double_t massp[2];
278 massp[0] = TDatabasePDG::Instance()->GetParticle(pdgprong0)->Mass();
279 massp[1] = TDatabasePDG::Instance()->GetParticle(pdgprong1)->Mass();
280
281 Double_t pStar = TMath::Sqrt(TMath::Power(massvtx*massvtx-massp[0]*massp[0]-massp[1]*massp[1],2.)-4.*massp[0]*massp[0]*massp[1]*massp[1])/(2.*massvtx);
282
283 Double_t beta = P()/E(pdgvtx);
284 Double_t gamma = E(pdgvtx)/massvtx;
285
286 Double_t cts = (QlProng(ip)/gamma-beta*TMath::Sqrt(pStar*pStar+massp[ip]*massp[ip]))/pStar;
287
288 return cts;
289}
290//---------------------------------------------------------------------------
291Double_t AliAODRecoDecay::Ct(UInt_t pdg,Double_t point[3]) const
292{
293 //
294 // Decay time * c assuming it is produced at "point" [cm]
295 //
296 Double_t mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
297 return DecayLength(point)*mass/P();
298}
299//---------------------------------------------------------------------------
300Double_t AliAODRecoDecay::E(UInt_t pdg) const
301{
302 //
303 // Energy
304 //
305 Double_t mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
306 return TMath::Sqrt(mass*mass+P()*P());
307}
308//---------------------------------------------------------------------------
309Double_t AliAODRecoDecay::EProng(Int_t ip,UInt_t pdg) const
310{
311 //
312 // Energy of ip-th prong
313 //
314 Double_t mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
315 return TMath::Sqrt(mass*mass+PProng(ip)*PProng(ip));
316}
e045ffda 317//--------------------------------------------------------------------------
318Bool_t AliAODRecoDecay::GetCovarianceXYZPxPyPz(Double_t cv[21]) const {
319 //
320 // This function returns the global covariance matrix of the track params
321 //
322 // Cov(x,x) ... : cv[0]
323 // Cov(y,x) ... : cv[1] cv[2]
324 // Cov(z,x) ... : cv[3] cv[4] cv[5]
325 // Cov(px,x)... : cv[6] cv[7] cv[8] cv[9]
326 // Cov(py,x)... : cv[10] cv[11] cv[12] cv[13] cv[14]
327 // Cov(pz,x)... : cv[15] cv[16] cv[17] cv[18] cv[19] cv[20]
7de7497b 328 //
e045ffda 329 // For XYZ we take the cov of the vertex, for PxPyPz we take the
330 // sum of the covs of PxPyPz from the daughters, for the moment
331 // we set the cov between position and momentum as the sum of
332 // the same cov from the daughters.
7de7497b 333 //
e045ffda 334
335 Int_t j;
336 for(j=0;j<21;j++) cv[j]=0.;
337
338 if(!GetNDaughters()) {
339 AliError("No daughters available");
340 return kFALSE;
341 }
342
343 Double_t v[6];
344 AliAODVertex *secv=GetSecondaryVtx();
345 if(!secv) {
346 AliError("Vertex covariance matrix not available");
347 return kFALSE;
348 }
349 if(!secv->GetCovMatrix(v)) {
350 AliError("Vertex covariance matrix not available");
351 return kFALSE;
352 }
353
354 Double_t p[21]; for(j=0;j<21;j++) p[j]=0.;
355 Bool_t error=kFALSE;
356 for(Int_t i=1; i<GetNDaughters(); i++) {
357 AliVTrack *daugh = (AliVTrack*)GetDaughter(i);
358 Double_t dcov[21];
359 if(!daugh->GetCovarianceXYZPxPyPz(dcov)) error=kTRUE;
360 for(j=0;j<21;j++) p[j] += dcov[j];
361 }
362 if(error) {
363 AliError("No covariance for at least one daughter")
364 return kFALSE;
365 }
366
367 for(j=0; j<21; j++) {
368 if(j<6) {
369 cv[j] = v[j];
370 } else {
371 cv[j] = p[j];
372 }
373 }
374
375 return kTRUE;
376}
7de7497b 377//----------------------------------------------------------------------------
e045ffda 378UChar_t AliAODRecoDecay::GetITSClusterMap() const {
379 //
380 // We take the logical AND of the daughters cluster maps
381 // (only if all daughters have the bit for given layer, we set the bit)
382 //
383 UChar_t map=0;
384
385 if(!GetNDaughters()) {
386 AliError("No daughters available");
387 return map;
388 }
389
390 for(Int_t l=0; l<12; l++) { // loop on ITS layers (from here we cannot know how many they are; let's put 12 to be conservative)
391 Int_t bit = 1;
392 for(Int_t i=0; i<GetNDaughters(); i++) {
393 AliVTrack *daugh = (AliVTrack*)GetDaughter(i);
394 if(!TESTBIT(daugh->GetITSClusterMap(),l)) bit=0;
395 }
396 if(bit) SETBIT(map,l);
397 }
398
399 return map;
400}
401//--------------------------------------------------------------------------
402ULong_t AliAODRecoDecay::GetStatus() const {
403 //
404 // Same as for ITSClusterMap
405 //
406 ULong_t status=0;
407
408 if(!GetNDaughters()) {
409 AliError("No daughters available");
410 return status;
411 }
412
413 AliVTrack *daugh0 = (AliVTrack*)GetDaughter(0);
414 status = status&(daugh0->GetStatus());
415
416 for(Int_t i=1; i<GetNDaughters(); i++) {
417 AliVTrack *daugh = (AliVTrack*)GetDaughter(i);
418 status = status&(daugh->GetStatus());
419 }
420
421 return status;
422}
423//--------------------------------------------------------------------------
7de7497b 424Double_t AliAODRecoDecay::ImpParXY(Double_t point[3]) const
425{
426 //
427 // Impact parameter in the bending plane of the particle
428 // w.r.t. to "point"
429 //
430 Double_t k = -(GetSecVtxX()-point[0])*Px()-(GetSecVtxY()-point[1])*Py();
431 k /= Pt()*Pt();
432 Double_t dx = GetSecVtxX()-point[0]+k*Px();
433 Double_t dy = GetSecVtxY()-point[1]+k*Py();
434 Double_t absImpPar = TMath::Sqrt(dx*dx+dy*dy);
435 TVector3 mom(Px(),Py(),Pz());
436 TVector3 fline(GetSecVtxX()-point[0],
437 GetSecVtxY()-point[1],
438 GetSecVtxZ()-point[2]);
439 TVector3 cross = mom.Cross(fline);
440 return (cross.Z()>0. ? absImpPar : -absImpPar);
441}
442//----------------------------------------------------------------------------
443Double_t AliAODRecoDecay::InvMass(Int_t npdg,UInt_t *pdg) const
444{
445 //
446 // Invariant mass for prongs mass hypotheses in pdg array
447 //
448 if(GetNProngs()!=npdg) {
449 printf("npdg != GetNProngs()");
450 return (Double_t)-99999.;
451 }
452 Double_t energysum = 0.;
453
454 for(Int_t i=0; i<GetNProngs(); i++) {
455 energysum += EProng(i,pdg[i]);
456 }
457
458 Double_t mass = TMath::Sqrt(energysum*energysum-P()*P());
459
460 return mass;
461}
462//----------------------------------------------------------------------------
463Double_t AliAODRecoDecay::InvMass2Prongs(Int_t ip1,Int_t ip2,
464 UInt_t pdg1,UInt_t pdg2) const
465{
466 //
467 // 2-prong(ip1,ip2) invariant mass for prongs mass hypotheses in pdg1,2
468 //
469 Double_t energysum = EProng(ip1,pdg1) + EProng(ip2,pdg2);
470 Double_t psum2 = (PxProng(ip1)+PxProng(ip2))*(PxProng(ip1)+PxProng(ip2))
471 +(PyProng(ip1)+PyProng(ip2))*(PyProng(ip1)+PyProng(ip2))
472 +(PzProng(ip1)+PzProng(ip2))*(PzProng(ip1)+PzProng(ip2));
473 Double_t mass = TMath::Sqrt(energysum*energysum-psum2);
474
475 return mass;
476}
bfd68286 477//----------------------------------------------------------------------------
478Int_t AliAODRecoDecay::MatchToMC(Int_t pdgabs,TClonesArray *mcArray) const
479{
480 //
481 // Check if this candidate is matched to a MC signal
482 // If no, return -1
483 // If yes, return label (>=0) of the AliAODMCParticle
484 //
485
486 if(!GetNDaughters()) {
487 AliError("No daughters available");
488 return -1;
489 }
490
491 Int_t *dgLabels = new Int_t[GetNDaughters()];
492
493 // loop on daughters and write the labels
494 for(Int_t i=0; i<GetNDaughters(); i++) {
495 AliVTrack *trk = (AliVTrack*)GetDaughter(i);
496 dgLabels[i] = trk->GetLabel();
497 }
498
499 Int_t labMother = MatchToMC(pdgabs,mcArray,dgLabels);
500
501 delete [] dgLabels; dgLabels=NULL;
502
503 return labMother;
504}
505//----------------------------------------------------------------------------
506Int_t AliAODRecoDecay::MatchToMC(Int_t pdgabs,TClonesArray *mcArray,
507 Int_t *dgLabels) const
508{
509 //
510 // Check if this candidate is matched to a MC signal
511 // If no, return -1
512 // If yes, return label (>=0) of the AliAODMCParticle
513 //
514
515 if(!GetNDaughters()) {
516 AliError("No daughters available");
517 return -1;
518 }
519
520 Int_t *labMom = new Int_t[GetNDaughters()];
521 Int_t i,lab,labMother,pdgMother;
522 AliAODMCParticle *part=0;
54ea2d63 523 Double_t pxSumDgs=0.,pySumDgs=0.,pzSumDgs=0.;
bfd68286 524
525 // loop on daughter labels
526 for(i=0; i<GetNDaughters(); i++) {
527 labMom[i]=-1;
528 lab = dgLabels[i];
529 if(lab<0) {
530 printf("daughter with negative label\n");
531 continue;
532 }
533 part = (AliAODMCParticle*)mcArray->At(lab);
534 if(!part) {
535 printf("no MC particle\n");
536 continue;
537 }
54ea2d63 538 // keep sum of daughters' momenta, to check for mom conservation
539 pxSumDgs += part->Px();
540 pySumDgs += part->Py();
541 pzSumDgs += part->Pz();
542 // for the J/psi, check that the daughters are electrons
543 if(pdgabs==443 && TMath::Abs(part->GetPdgCode())!=11) continue;
544
bfd68286 545 while(part->GetMother()>=0) {
546 labMother=part->GetMother();
547 part = (AliAODMCParticle*)mcArray->At(labMother);
548 if(!part) {
549 printf("no MC mother particle\n");
550 break;
551 }
552 pdgMother = TMath::Abs(part->GetPdgCode());
553 if(pdgMother==pdgabs) {
554 labMom[i]=labMother;
555 break;
556 }
557 }
558 }
559
560 // check if the candidate is signal
561 Bool_t isSignal=kTRUE;
562 labMother=labMom[0];
563 // all labels have to be the same and !=-1
564 for(i=0; i<GetNDaughters(); i++) {
565 if(labMom[i]==-1 || labMom[i]!=labMother) isSignal=kFALSE;
566 }
567
568 delete [] labMom; labMom=NULL;
569
570 if(!isSignal) return -1;
571
54ea2d63 572 /*
bfd68286 573 // check that the mother decayed in <GetNDaughters()> prongs
bfd68286 574 Int_t ndg = TMath::Abs(part->GetDaughter(1)-part->GetDaughter(0))+1;
bfd68286 575 if(ndg!=GetNDaughters()) return -1;
54ea2d63 576 AliAODMCParticle* p1=(AliAODMCParticle*)(mcArray->At(part->GetDaughter(1)));
577 AliAODMCParticle* p0=(AliAODMCParticle*)(mcArray->At(part->GetDaughter(0)));
578 printf("pdg %d %d %d %d %d %d\n",part->GetDaughter(1),part->GetDaughter(0),dgLabels[0],dgLabels[1],p0->GetPdgCode(),p1->GetPdgCode());
579 */
bfd68286 580
54ea2d63 581 // the above works only for non-resonant decays,
582 // it's better to check for mom conservation
583 part = (AliAODMCParticle*)mcArray->At(labMother);
584 Double_t pxMother = part->Px();
585 Double_t pyMother = part->Py();
586 Double_t pzMother = part->Pz();
a364c44a 587 // within 0.1%
588 if((TMath::Abs(pxMother-pxSumDgs)/(TMath::Abs(pxMother)+1.e-13)) > 0.001 &&
589 (TMath::Abs(pyMother-pySumDgs)/(TMath::Abs(pyMother)+1.e-13)) > 0.001 &&
590 (TMath::Abs(pzMother-pzSumDgs)/(TMath::Abs(pzMother)+1.e-13)) > 0.001)
591 return -1;
54ea2d63 592
bfd68286 593 return labMother;
594}
7de7497b 595//---------------------------------------------------------------------------
596void AliAODRecoDecay::Print(Option_t* /*option*/) const
597{
598 //
599 // Print some information
600 //
601 printf("AliAODRecoDecay with %d prongs\n",GetNProngs());
602 printf("Secondary Vertex: (%f, %f, %f)\n",GetSecVtxX(),GetSecVtxY(),GetSecVtxZ());
603
604 return;
605}
606//----------------------------------------------------------------------------
607Double_t AliAODRecoDecay::ProngsRelAngle(Int_t ip1,Int_t ip2) const
608{
609 //
610 // Relative angle between two prongs
611 //
612 TVector3 momA(PxProng(ip1),PyProng(ip1),PzProng(ip1));
613 TVector3 momB(PxProng(ip2),PyProng(ip2),PzProng(ip2));
614
615 Double_t angle = momA.Angle(momB);
616
617 return angle;
618}
619//----------------------------------------------------------------------------
620Double_t AliAODRecoDecay::QlProng(Int_t ip) const
621{
622 //
623 // Longitudinal momentum of prong w.r.t. to total momentum
624 //
625 TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
626 TVector3 momTot(Px(),Py(),Pz());
627
628 return mom.Dot(momTot)/momTot.Mag();
629}
630//----------------------------------------------------------------------------
631Double_t AliAODRecoDecay::QtProng(Int_t ip) const
632{
633 //
634 // Transverse momentum of prong w.r.t. to total momentum
635 //
636 TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
637 TVector3 momTot(Px(),Py(),Pz());
638
639 return mom.Perp(momTot);
640}
641//----------------------------------------------------------------------------
642Double_t AliAODRecoDecay::QlProngFlightLine(Int_t ip,Double_t point[3]) const
643{
644 //
645 // Longitudinal momentum of prong w.r.t. to flight line between "point"
646 // and fSecondaryVtx
647 //
648 TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
649 TVector3 fline(GetSecVtxX()-point[0],
650 GetSecVtxY()-point[1],
651 GetSecVtxZ()-point[2]);
652
653 return mom.Dot(fline)/fline.Mag();
654}
655//----------------------------------------------------------------------------
656Double_t AliAODRecoDecay::QtProngFlightLine(Int_t ip,Double_t point[3]) const
657{
658 //
659 // Transverse momentum of prong w.r.t. to flight line between "point" and
660 // fSecondaryVtx
661 //
662 TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
663 TVector3 fline(GetSecVtxX()-point[0],
664 GetSecVtxY()-point[1],
665 GetSecVtxZ()-point[2]);
666
667 return mom.Perp(fline);
668}
669//--------------------------------------------------------------------------