1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //====================================================================================================================================================
18 // Description of an ALICE muon forward track, combining the information of the Muon Spectrometer and the Muon Forward Tracker
20 // Contact author: antonio.uras@cern.ch
22 //====================================================================================================================================================
25 #include "AliMUONTrack.h"
26 #include "AliMFTCluster.h"
27 #include "AliMUONVCluster.h"
28 #include "AliMUONTrackParam.h"
29 #include "AliMUONTrackExtrap.h"
30 #include "TClonesArray.h"
32 #include "TParticle.h"
33 #include "AliMuonForwardTrack.h"
34 #include "AliMFTConstants.h"
35 #include "TLorentzVector.h"
36 #include "TDatabasePDG.h"
37 #include "AliMUONConstants.h"
39 ClassImp(AliMuonForwardTrack)
41 //====================================================================================================================================================
43 AliMuonForwardTrack::AliMuonForwardTrack():
48 fNWrongClustersMC(-1),
55 // default constructor
57 for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = kFALSE;
58 for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
59 fParentMCLabel[iParent] = -1;
60 fParentPDGCode[iParent] = 0;
62 fMFTClusters = new TClonesArray("AliMFTCluster");
63 fMFTClusters -> SetOwner(kTRUE);
67 //====================================================================================================================================================
69 AliMuonForwardTrack::AliMuonForwardTrack(AliMUONTrack *MUONTrack):
74 fNWrongClustersMC(-1),
81 SetMUONTrack(MUONTrack);
82 for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = kFALSE;
83 for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
84 fParentMCLabel[iParent] = -1;
85 fParentPDGCode[iParent] = 0;
87 fMFTClusters = new TClonesArray("AliMFTCluster");
88 fMFTClusters -> SetOwner(kTRUE);
92 //====================================================================================================================================================
94 AliMuonForwardTrack::AliMuonForwardTrack(const AliMuonForwardTrack& track):
99 fNWrongClustersMC(track.fNWrongClustersMC),
100 fTrackMCId(track.fTrackMCId),
101 fKinem(track.fKinem),
102 fParamCovMatrix(track.fParamCovMatrix),
103 fRAtAbsorberEnd(track.fRAtAbsorberEnd)
107 fMUONTrack = new AliMUONTrack(*(track.fMUONTrack));
108 if (track.fMCTrackRef) fMCTrackRef = new TParticle(*(track.fMCTrackRef));
109 fMFTClusters = new TClonesArray(*(track.fMFTClusters));
110 fMFTClusters->SetOwner(kTRUE);
111 for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = (track.fPlaneExists)[iPlane];
112 for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
113 fParentMCLabel[iParent] = (track.fParentMCLabel)[iParent];
114 fParentPDGCode[iParent] = (track.fParentPDGCode)[iParent];
119 //====================================================================================================================================================
121 AliMuonForwardTrack& AliMuonForwardTrack::operator=(const AliMuonForwardTrack& track) {
123 // Asignment operator
125 // check assignement to self
126 if (this == &track) return *this;
128 // base class assignement
129 AliMUONTrack::operator=(track);
134 fMUONTrack = new AliMUONTrack(*(track.fMUONTrack));
135 if (track.fMCTrackRef) fMCTrackRef = new TParticle(*(track.fMCTrackRef));
136 fMFTClusters = new TClonesArray(*(track.fMFTClusters));
137 fMFTClusters->SetOwner(kTRUE);
138 fNWrongClustersMC = track.fNWrongClustersMC;
139 fTrackMCId = track.fTrackMCId;
140 fKinem = track.fKinem;
141 fParamCovMatrix = track.fParamCovMatrix;
142 fRAtAbsorberEnd = track.fRAtAbsorberEnd;
144 for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = (track.fPlaneExists)[iPlane];
145 for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
146 fParentMCLabel[iParent] = (track.fParentMCLabel)[iParent];
147 fParentPDGCode[iParent] = (track.fParentPDGCode)[iParent];
154 //====================================================================================================================================================
156 void AliMuonForwardTrack::Clear(const Option_t* /*opt*/) {
159 fMFTClusters -> Delete();
160 delete fMFTClusters; fMFTClusters = 0x0;
161 delete fMUONTrack; fMUONTrack = 0x0;
162 delete fMCTrackRef; fMCTrackRef = 0x0;
166 //====================================================================================================================================================
168 AliMuonForwardTrack::~AliMuonForwardTrack() {
172 fMFTClusters -> Delete();
177 //====================================================================================================================================================
179 void AliMuonForwardTrack::SetMUONTrack(AliMUONTrack *MUONTrack) {
182 AliInfo("fMUONTrack already exists, nothing will be done");
186 fMUONTrack = MUONTrack;
187 SetGlobalChi2(fMUONTrack->GetGlobalChi2());
191 //====================================================================================================================================================
193 void AliMuonForwardTrack::SetMCTrackRef(TParticle *MCTrackRef) {
196 AliInfo("fMCTrackRef already exists, nothing will be done");
200 fMCTrackRef = MCTrackRef;
204 //====================================================================================================================================================
206 void AliMuonForwardTrack::AddTrackParamAtMFTCluster(AliMUONTrackParam &trackParam, AliMFTCluster &mftCluster) {
208 AliDebug(1, Form("Before adding: this->fMFTClusters=%p has %d entries", this->fMFTClusters, this->fMFTClusters->GetEntries()));
209 Int_t iMFTCluster = this->fMFTClusters->GetEntries();
210 AliDebug(1, Form("mftCluster->GetX() = %f mftCluster->GetY() = %f mftCluster->GetErrX() = %f cmftCluster->GetErrY() = %f",
211 mftCluster.GetX(), mftCluster.GetY(), mftCluster.GetErrX(), mftCluster.GetErrY()));
212 AliMUONVCluster *muonCluster = (AliMUONVCluster*) mftCluster.CreateMUONCluster();
213 AliDebug(1, Form("Created MUON cluster %p", muonCluster));
214 trackParam.SetUniqueID(iMFTCluster); // we profit of this slot to store the reference to the corresponding MFTCluster
215 AliDebug(1, Form("Now adding muonCluster %p and trackParam %p",muonCluster, &trackParam));
216 AddTrackParamAtCluster(trackParam, *muonCluster, kTRUE);
217 // we pass the parameters this->GetTrackParamAtCluster()->First() to the Kalman Filter algorithm: they will be updated!!
218 Double_t chi2Kalman = RunKalmanFilter(*(AliMUONTrackParam*)(GetTrackParamAtCluster()->First()));
219 AliDebug(1, Form("Adding Kalman chi2 = %f to global chi2 = %f", chi2Kalman, GetGlobalChi2()));
220 Double_t newGlobalChi2 = GetGlobalChi2() + chi2Kalman;
221 mftCluster.SetLocalChi2(chi2Kalman);
222 mftCluster.SetTrackChi2(newGlobalChi2);
223 new ((*(this->fMFTClusters))[iMFTCluster]) AliMFTCluster(mftCluster);
224 AliDebug(1, Form("GetTrackParamAtCluster() = %p has %d entries while this->fMFTClusters=%p has %d entries",
225 GetTrackParamAtCluster(), GetTrackParamAtCluster()->GetEntries(), this->fMFTClusters, this->fMFTClusters->GetEntries()));
226 AliDebug(1, Form("muonCluster->GetZ() = %f, trackParam->GetZ() = %f",muonCluster->GetZ(), trackParam.GetZ()));
227 SetGlobalChi2(newGlobalChi2);
228 AliDebug(1, Form("New global chi2 = %f", GetGlobalChi2()));
229 ((AliMUONTrackParam*) GetTrackParamAtCluster()->First())->SetTrackChi2(newGlobalChi2);
230 for (Int_t iPar=0; iPar<GetTrackParamAtCluster()->GetEntries(); iPar++) {
231 AliDebug(1, Form("GetTrackParamAtCluster()->At(%d)->GetClusterPtr() = %p",
232 iPar, ((AliMUONTrackParam*)(GetTrackParamAtCluster()->At(iPar)))->GetClusterPtr()));
237 //====================================================================================================================================================
239 AliMUONTrackParam* AliMuonForwardTrack::GetTrackParamAtMUONCluster(Int_t iMUONCluster) {
241 if (iMUONCluster<0 || iMUONCluster>=GetNMUONClusters()) {
242 AliError("Invalid MUON cluster index. NULL pointer will be returned");
246 AliMUONTrackParam *trackParam = (AliMUONTrackParam*) fMUONTrack->GetTrackParamAtCluster()->At(iMUONCluster);
252 //====================================================================================================================================================
254 AliMUONTrackParam* AliMuonForwardTrack::GetTrackParamAtMFTCluster(Int_t iMFTCluster) {
256 if (iMFTCluster<0 || iMFTCluster>=GetNMFTClusters()) {
257 AliError("Invalid MFT cluster index. NULL pointer will be returned");
261 AliMUONTrackParam *trackParam = (AliMUONTrackParam*) GetTrackParamAtCluster()->At(iMFTCluster);
267 //====================================================================================================================================================
269 AliMUONVCluster* AliMuonForwardTrack::GetMUONCluster(Int_t iMUONCluster) {
271 if (iMUONCluster<0 || iMUONCluster>=GetNMUONClusters()) {
272 AliError("Invalid MUON cluster index. NULL pointer will be returned");
276 AliMUONTrackParam *trackParam = GetTrackParamAtMUONCluster(iMUONCluster);
277 AliMUONVCluster *muonCluster = trackParam->GetClusterPtr();
283 //====================================================================================================================================================
285 AliMFTCluster* AliMuonForwardTrack::GetMFTCluster(Int_t iMFTCluster) {
287 if (iMFTCluster<0 || iMFTCluster>=GetNMFTClusters()) {
288 AliError(Form("Invalid MFT cluster index (%d). GetNMFTClusters()=%d. NULL pointer will be returned", iMFTCluster, GetNMFTClusters()));
292 AliMUONTrackParam *trackParam = GetTrackParamAtMFTCluster(iMFTCluster);
293 AliMFTCluster *mftCluster = (AliMFTCluster*) this->fMFTClusters->At(trackParam->GetUniqueID());
299 //====================================================================================================================================================
301 Double_t AliMuonForwardTrack::RunKalmanFilter(AliMUONTrackParam &trackParamAtCluster) {
303 AliDebug(1, Form("Running Kalman filter for parameters %p (z = %f) and cluster %p (z = %f)",
304 &trackParamAtCluster, trackParamAtCluster.GetZ(), trackParamAtCluster.GetClusterPtr(), trackParamAtCluster.GetClusterPtr()->GetZ()));
306 // Compute new track parameters and their covariances including new cluster using kalman filter
307 // return the *additional* track chi2
309 // Get actual track parameters (p)
310 TMatrixD param(trackParamAtCluster.GetParameters());
312 // Get new cluster parameters (m)
313 AliMUONVCluster *cluster = trackParamAtCluster.GetClusterPtr();
314 AliDebug(1, Form("cluster->GetX() = %f cluster->GetY() = %f cluster->GetErrX() = %f cluster->GetErrY() = %f",
315 cluster->GetX(), cluster->GetY(), cluster->GetErrX(), cluster->GetErrY()));
316 TMatrixD clusterParam(5,1);
318 clusterParam(0,0) = cluster->GetX();
319 clusterParam(2,0) = cluster->GetY();
321 // Compute the current parameter weight (W)
322 TMatrixD paramWeight(trackParamAtCluster.GetCovariances());
323 if (paramWeight.Determinant() != 0) {
324 paramWeight.Invert();
326 Warning("RunKalmanFilter"," Determinant = 0");
330 // Compute the new cluster weight (U)
331 TMatrixD clusterWeight(5,5);
332 clusterWeight.Zero();
333 clusterWeight(0,0) = 1. / cluster->GetErrX2();
334 clusterWeight(2,2) = 1. / cluster->GetErrY2();
336 // Compute the new parameters covariance matrix ( (W+U)^-1 )
337 TMatrixD newParamCov(paramWeight,TMatrixD::kPlus,clusterWeight);
338 if (newParamCov.Determinant() != 0) {
339 newParamCov.Invert();
342 Warning("RunKalmanFilter"," Determinant = 0");
346 // Save the new parameters covariance matrix
347 trackParamAtCluster.SetCovariances(newParamCov);
349 // Compute the new parameters (p' = ((W+U)^-1)U(m-p) + p)
350 TMatrixD tmp(clusterParam,TMatrixD::kMinus,param);
351 TMatrixD tmp2(clusterWeight,TMatrixD::kMult,tmp); // U(m-p)
352 TMatrixD newParam(newParamCov,TMatrixD::kMult,tmp2); // ((W+U)^-1)U(m-p)
353 newParam += param; // ((W+U)^-1)U(m-p) + p
355 // Save the new parameters
356 trackParamAtCluster.SetParameters(newParam);
358 // Compute the additional chi2 (= ((p'-p)^-1)W(p'-p) + ((p'-m)^-1)U(p'-m))
359 tmp = newParam; // p'
360 tmp -= param; // (p'-p)
361 TMatrixD tmp3(paramWeight,TMatrixD::kMult,tmp); // W(p'-p)
362 TMatrixD addChi2Track(tmp,TMatrixD::kTransposeMult,tmp3); // ((p'-p)^-1)W(p'-p)
363 tmp = newParam; // p'
364 tmp -= clusterParam; // (p'-m)
365 TMatrixD tmp4(clusterWeight,TMatrixD::kMult,tmp); // U(p'-m)
366 addChi2Track += TMatrixD(tmp,TMatrixD::kTransposeMult,tmp4); // ((p'-p)^-1)W(p'-p) + ((p'-m)^-1)U(p'-m)
368 AliDebug(1,Form("Adding Kalman chi2 = %f",addChi2Track(0,0)));
370 return addChi2Track(0,0);
374 //====================================================================================================================================================
376 Double_t AliMuonForwardTrack::GetWeightedOffset(Double_t x, Double_t y, Double_t z) {
378 AliMUONTrackParam *param = GetTrackParamAtMFTCluster(0);
379 AliMUONTrackExtrap::ExtrapToZCov(param, z);
382 cov = param->GetCovariances();
384 TMatrixD covCoordinates(2,2);
385 covCoordinates(0,0) = cov(0,0);
386 covCoordinates(0,1) = cov(0,2);
387 covCoordinates(1,0) = cov(2,0);
388 covCoordinates(1,1) = cov(2,2);
390 TMatrixD covCoordinatesInverse = covCoordinates.Invert();
392 Double_t dX = param->GetNonBendingCoor() - x;
393 Double_t dY = param->GetBendingCoor() - y;
395 Double_t weightedOffset = TMath::Sqrt(0.5*(dX*dX*covCoordinatesInverse(0,0) +
396 dY*dY*covCoordinatesInverse(1,1) +
397 2.*dX*dY*covCoordinatesInverse(0,1)));
399 return weightedOffset;
403 //====================================================================================================================================================
405 Double_t AliMuonForwardTrack::GetOffsetX(Double_t x, Double_t z) {
407 AliMUONTrackParam *param = GetTrackParamAtMFTCluster(0);
408 AliMUONTrackExtrap::ExtrapToZCov(param, z);
409 Double_t dX = param->GetNonBendingCoor() - x;
414 //====================================================================================================================================================
416 Double_t AliMuonForwardTrack::GetOffsetY(Double_t y, Double_t z) {
418 AliMUONTrackParam *param = GetTrackParamAtMFTCluster(0);
419 AliMUONTrackExtrap::ExtrapToZCov(param, z);
420 Double_t dY = param->GetBendingCoor() - y;
425 //====================================================================================================================================================
427 Double_t AliMuonForwardTrack::GetOffset(Double_t x, Double_t y, Double_t z) {
429 AliMUONTrackParam *param = GetTrackParamAtMFTCluster(0);
430 AliMUONTrackExtrap::ExtrapToZCov(param, z);
431 Double_t dX = param->GetNonBendingCoor() - x;
432 Double_t dY = param->GetBendingCoor() - y;
433 return TMath::Sqrt(dX*dX + dY*dY);
437 //====================================================================================================================================================
439 Double_t AliMuonForwardTrack::GetDCA(Double_t x, Double_t y, Double_t z) {
441 // Distance of Closest Approach, according to the standard MUON terminology. Actually, the offset of the track w.r.t. the primary vertex,
442 // where the extrapolation of the track DOES NOT include the MFT information
444 AliMUONTrackParam param = (*((AliMUONTrackParam*)(GetTrackParamAtMUONCluster(0))));
446 AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¶m, z);
447 Double_t dX = param.GetNonBendingCoor() - x;
448 Double_t dY = param.GetBendingCoor() - y;
449 return TMath::Sqrt(dX*dX + dY*dY);
453 //====================================================================================================================================================
455 Double_t AliMuonForwardTrack::GetMomentumSpectrometer(Double_t z) {
457 // Momentum of the track at the primary vertex plane, where the extrapolation of the track DOES NOT include the MFT information
459 AliMUONTrackParam param = (*((AliMUONTrackParam*)(GetTrackParamAtMUONCluster(0))));
461 AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¶m, z);
466 //====================================================================================================================================================
468 Double_t AliMuonForwardTrack::GetThetaAbs() {
470 // it is the angle defined by the imaginary line goingo from the vertex to the exit point of the track at the end of the hadron absorber
472 Double_t z = AliMUONConstants::AbsZEnd();
474 AliMUONTrackParam param = (*((AliMUONTrackParam*)(GetTrackParamAtMUONCluster(0))));
476 AliMUONTrackExtrap::ExtrapToZ(¶m, z);
477 Double_t x = param.GetNonBendingCoor();
478 Double_t y = param.GetBendingCoor();
480 return 180. +TMath::ATan(TMath::Sqrt(x*x + y*y)/z)*TMath::RadToDeg();
484 //====================================================================================================================================================
486 Bool_t AliMuonForwardTrack::IsFromDirectResonance() {
488 if (!IsMuon()) return kFALSE;
490 for (Int_t i=GetFirstMotherID(); i>=0; i--) if (!IsPDGResonance(GetParentPDGCode(i))) return kFALSE; // the decay chain finds a non-resonance particle
496 //====================================================================================================================================================
498 Bool_t AliMuonForwardTrack::IsFromChainResonance() {
500 if (!IsMuon()) return kFALSE;
502 if (GetFirstMotherID() == 0) return kFALSE; // it is not a chain
504 if (!IsPDGResonance(GetParentPDGCode(GetFirstMotherID()))) return kFALSE; // primordial is not a resonance
506 for (Int_t i=GetFirstMotherID()-1; i>=0; i--) if (!IsPDGResonance(GetParentPDGCode(i))) return kTRUE; // the decay chain finds a non-resonance particle
512 //====================================================================================================================================================
514 Bool_t AliMuonForwardTrack::IsFromDirectCharm() {
516 if (!IsMuon()) return kFALSE;
518 for (Int_t i=GetFirstMotherID(); i>=0; i--) if (!IsPDGCharm(GetParentPDGCode(i))) return kFALSE; // the decay chain finds a non-charmed particle
524 //====================================================================================================================================================
526 Bool_t AliMuonForwardTrack::IsFromChainCharm() {
528 if (!IsMuon()) return kFALSE;
530 if (GetFirstMotherID() == 0) return kFALSE; // it is not a chain
532 if (!IsPDGCharm(GetParentPDGCode(GetFirstMotherID()))) return kFALSE; // primordial is not a charmed hadron
534 for (Int_t i=GetFirstMotherID()-1; i>=0; i--) if (!IsPDGCharm(GetParentPDGCode(i))) return kTRUE; // the decay chain finds a non-charmed particle
540 //====================================================================================================================================================
542 Bool_t AliMuonForwardTrack::IsFromDirectBeauty() {
544 if (!IsMuon()) return kFALSE;
546 for (Int_t i=GetFirstMotherID(); i>=0; i--) if (!IsPDGBeauty(GetParentPDGCode(i))) return kFALSE; // the decay chain finds a non-beauty particle
552 //====================================================================================================================================================
554 Bool_t AliMuonForwardTrack::IsFromChainBeauty() {
556 if (!IsMuon()) return kFALSE;
558 if (GetFirstMotherID() == 0) return kFALSE; // it is not a chain
560 if (!IsPDGBeauty(GetParentPDGCode(GetFirstMotherID()))) return kFALSE; // primordial is not a beauty hadron
562 for (Int_t i=GetFirstMotherID()-1; i>=0; i--) if (!IsPDGBeauty(GetParentPDGCode(i))) return kTRUE; // the decay chain finds a non-beauty particle
568 //====================================================================================================================================================
570 Bool_t AliMuonForwardTrack::IsMuon() {
572 if (IsFake()) return kFALSE;
574 if (TMath::Abs(fMCTrackRef->GetPdgCode())==13) return kTRUE;
580 //====================================================================================================================================================
582 Bool_t AliMuonForwardTrack::IsFake() {
584 if (!fMCTrackRef) return kTRUE;
590 //====================================================================================================================================================
592 Bool_t AliMuonForwardTrack::IsPDGResonance(Int_t pdg) {
594 // if (pdg<10) return kFALSE;
595 // Int_t id = pdg%100000;
596 // return (!((id-id%10)%110));
598 Bool_t result = kFALSE;
608 pdg == 100553 ) result = kTRUE;
614 //====================================================================================================================================================
616 Bool_t AliMuonForwardTrack::IsPDGCharm(Int_t pdg) {
618 Bool_t result = kFALSE;
620 if ( TMath::Abs(pdg) == 411 ||
621 TMath::Abs(pdg) == 421 ||
622 TMath::Abs(pdg) == 10411 ||
623 TMath::Abs(pdg) == 10421 ||
624 TMath::Abs(pdg) == 413 ||
625 TMath::Abs(pdg) == 423 ||
626 TMath::Abs(pdg) == 10413 ||
627 TMath::Abs(pdg) == 10423 ||
628 TMath::Abs(pdg) == 20413 ||
629 TMath::Abs(pdg) == 20423 ||
630 TMath::Abs(pdg) == 415 ||
631 TMath::Abs(pdg) == 425 ||
632 TMath::Abs(pdg) == 431 ||
633 TMath::Abs(pdg) == 10431 ||
634 TMath::Abs(pdg) == 433 ||
635 TMath::Abs(pdg) == 10433 ||
636 TMath::Abs(pdg) == 20433 ||
637 TMath::Abs(pdg) == 435 ) result = kTRUE;
643 //====================================================================================================================================================
645 Bool_t AliMuonForwardTrack::IsPDGBeauty(Int_t pdg) {
647 Bool_t result = kFALSE;
649 if ( TMath::Abs(pdg) == 511 ||
650 TMath::Abs(pdg) == 521 ||
651 TMath::Abs(pdg) == 10511 ||
652 TMath::Abs(pdg) == 10521 ||
653 TMath::Abs(pdg) == 513 ||
654 TMath::Abs(pdg) == 523 ||
655 TMath::Abs(pdg) == 10513 ||
656 TMath::Abs(pdg) == 10523 ||
657 TMath::Abs(pdg) == 20513 ||
658 TMath::Abs(pdg) == 20523 ||
659 TMath::Abs(pdg) == 515 ||
660 TMath::Abs(pdg) == 525 ||
661 TMath::Abs(pdg) == 531 ||
662 TMath::Abs(pdg) == 10531 ||
663 TMath::Abs(pdg) == 533 ||
664 TMath::Abs(pdg) == 10533 ||
665 TMath::Abs(pdg) == 20533 ||
666 TMath::Abs(pdg) == 535 ||
667 TMath::Abs(pdg) == 541 ||
668 TMath::Abs(pdg) == 10541 ||
669 TMath::Abs(pdg) == 543 ||
670 TMath::Abs(pdg) == 10543 ||
671 TMath::Abs(pdg) == 20543 ||
672 TMath::Abs(pdg) == 545 ) result = kTRUE;
678 //====================================================================================================================================================
680 Bool_t AliMuonForwardTrack::IsMuonFromBackground() {
682 Bool_t result = kFALSE;
684 if (!IsMuon()) return result;
686 if (!IsFromDirectResonance() && !IsFromChainResonance() && !IsFromDirectCharm() && !IsFromChainCharm() && !IsFromDirectBeauty() && !IsFromChainBeauty()) result = kTRUE;
688 if (result) AliDebug(1, Form("Muon comes from a background source %d", GetParentPDGCode(0)));
694 //====================================================================================================================================================
696 void AliMuonForwardTrack::EvalKinem(Double_t z) {
698 AliMUONTrackParam *param = GetTrackParamAtMFTCluster(0);
699 AliMUONTrackExtrap::ExtrapToZCov(param, z);
701 Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
702 Double_t energy = TMath::Sqrt(param->P()*param->P() + mMu*mMu);
703 fKinem.SetPxPyPzE(param->Px(), param->Py(), param->Pz(), energy);
706 fParamCovMatrix = param->GetCovariances();
710 //====================================================================================================================================================
712 Int_t AliMuonForwardTrack::GetFirstMotherID() {
714 Int_t motherLevel = 0;
715 Int_t motherMCLabel = GetParentMCLabel(motherLevel);
716 while (motherMCLabel >= 0) {
718 motherMCLabel = GetParentMCLabel(motherLevel);
720 return motherLevel-1;
724 //====================================================================================================================================================
726 void AliMuonForwardTrack::PrintHistory() {
728 if (IsFake()) printf("Track is a fake MUON\n");
730 TString history = "";
731 for (Int_t i=GetFirstMotherID(); i>=0; i--) {
732 history += TDatabasePDG::Instance()->GetParticle(GetParentPDGCode(i))->GetName();
735 history += TDatabasePDG::Instance()->GetParticle(fMCTrackRef->GetPdgCode())->GetName();
736 printf("%s\n",history.Data());
741 //====================================================================================================================================================