#include "TMatrixD.h"
#include "TParticle.h"
#include "AliMuonForwardTrack.h"
+#include "AliMFTConstants.h"
+#include "TLorentzVector.h"
+#include "TDatabasePDG.h"
+#include "AliMUONConstants.h"
ClassImp(AliMuonForwardTrack)
AliMUONTrack(),
fMUONTrack(0),
fMCTrackRef(0),
- fMFTClusters(0)
+ fMFTClusters(0),
+ fNWrongClustersMC(-1),
+ fTrackMCId(-1),
+ fKinem(0,0,0,0),
+ fParamCovMatrix(5,5),
+ fRAtAbsorberEnd(0)
{
// default constructor
- for (Int_t iPlane=0; iPlane<fMaxNPlanesMFT; iPlane++) fPlaneExists[iPlane] = kFALSE;
+
+ for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = kFALSE;
+ for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
+ fParentMCLabel[iParent] = -1;
+ fParentPDGCode[iParent] = 0;
+ }
fMFTClusters = new TClonesArray("AliMFTCluster");
+ fMFTClusters -> SetOwner(kTRUE);
}
AliMUONTrack(),
fMUONTrack(0),
fMCTrackRef(0),
- fMFTClusters(0)
+ fMFTClusters(0),
+ fNWrongClustersMC(-1),
+ fTrackMCId(-1),
+ fKinem(0,0,0,0),
+ fParamCovMatrix(5,5),
+ fRAtAbsorberEnd(0)
{
SetMUONTrack(MUONTrack);
- for (Int_t iPlane=0; iPlane<fMaxNPlanesMFT; iPlane++) fPlaneExists[iPlane] = kFALSE;
+ for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = kFALSE;
+ for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
+ fParentMCLabel[iParent] = -1;
+ fParentPDGCode[iParent] = 0;
+ }
fMFTClusters = new TClonesArray("AliMFTCluster");
+ fMFTClusters -> SetOwner(kTRUE);
}
AliMuonForwardTrack::AliMuonForwardTrack(const AliMuonForwardTrack& track):
AliMUONTrack(track),
- fMUONTrack(track.fMUONTrack),
- fMCTrackRef(track.fMCTrackRef),
- fMFTClusters(track.fMFTClusters)
+ fMUONTrack(0x0),
+ fMCTrackRef(0x0),
+ fMFTClusters(0x0),
+ fNWrongClustersMC(track.fNWrongClustersMC),
+ fTrackMCId(track.fTrackMCId),
+ fKinem(track.fKinem),
+ fParamCovMatrix(track.fParamCovMatrix),
+ fRAtAbsorberEnd(track.fRAtAbsorberEnd)
{
// copy constructor
- for (Int_t iPlane=0; iPlane<fMaxNPlanesMFT; iPlane++) fPlaneExists[iPlane] = (track.fPlaneExists)[iPlane];
+ fMUONTrack = new AliMUONTrack(*(track.fMUONTrack));
+ if (track.fMCTrackRef) fMCTrackRef = new TParticle(*(track.fMCTrackRef));
+ fMFTClusters = new TClonesArray(*(track.fMFTClusters));
+ fMFTClusters->SetOwner(kTRUE);
+ for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = (track.fPlaneExists)[iPlane];
+ for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
+ fParentMCLabel[iParent] = (track.fParentMCLabel)[iParent];
+ fParentPDGCode[iParent] = (track.fParentPDGCode)[iParent];
+ }
}
AliMUONTrack::operator=(track);
// clear memory
- Clear();
+ Clear("");
- fMUONTrack = track.fMUONTrack;
- fMCTrackRef = track.fMCTrackRef;
- fMFTClusters = track.fMFTClusters;
-
- for (Int_t iPlane=0; iPlane<fMaxNPlanesMFT; iPlane++) fPlaneExists[iPlane] = (track.fPlaneExists)[iPlane];
+ fMUONTrack = new AliMUONTrack(*(track.fMUONTrack));
+ if (track.fMCTrackRef) fMCTrackRef = new TParticle(*(track.fMCTrackRef));
+ fMFTClusters = new TClonesArray(*(track.fMFTClusters));
+ fMFTClusters->SetOwner(kTRUE);
+ fNWrongClustersMC = track.fNWrongClustersMC;
+ fTrackMCId = track.fTrackMCId;
+ fKinem = track.fKinem;
+ fParamCovMatrix = track.fParamCovMatrix;
+ fRAtAbsorberEnd = track.fRAtAbsorberEnd;
+
+ for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) fPlaneExists[iPlane] = (track.fPlaneExists)[iPlane];
+ for (Int_t iParent=0; iParent<fgkNParentsMax; iParent++) {
+ fParentMCLabel[iParent] = (track.fParentMCLabel)[iParent];
+ fParentPDGCode[iParent] = (track.fParentPDGCode)[iParent];
+ }
return *this;
//====================================================================================================================================================
+void AliMuonForwardTrack::Clear(const Option_t* /*opt*/) {
+
+ // Clear arrays
+ fMFTClusters -> Delete();
+ delete fMFTClusters; fMFTClusters = 0x0;
+ delete fMUONTrack; fMUONTrack = 0x0;
+ delete fMCTrackRef; fMCTrackRef = 0x0;
+
+}
+
+//====================================================================================================================================================
+
+AliMuonForwardTrack::~AliMuonForwardTrack() {
+
+ delete fMUONTrack;
+ delete fMCTrackRef;
+ fMFTClusters -> Delete();
+ delete fMFTClusters;
+
+}
+
+//====================================================================================================================================================
+
void AliMuonForwardTrack::SetMUONTrack(AliMUONTrack *MUONTrack) {
if (fMUONTrack) {
}
fMUONTrack = MUONTrack;
+ SetGlobalChi2(fMUONTrack->GetGlobalChi2());
}
void AliMuonForwardTrack::AddTrackParamAtMFTCluster(AliMUONTrackParam &trackParam, AliMFTCluster &mftCluster) {
- AliDebug(1, Form("fMFTClusters=%p has %d entries", fMFTClusters, fMFTClusters->GetEntries()));
- Int_t iMFTCluster = fMFTClusters->GetEntries();
+ AliDebug(1, Form("Before adding: this->fMFTClusters=%p has %d entries", this->fMFTClusters, this->fMFTClusters->GetEntries()));
+ Int_t iMFTCluster = this->fMFTClusters->GetEntries();
AliDebug(1, Form("mftCluster->GetX() = %f mftCluster->GetY() = %f mftCluster->GetErrX() = %f cmftCluster->GetErrY() = %f",
mftCluster.GetX(), mftCluster.GetY(), mftCluster.GetErrX(), mftCluster.GetErrY()));
AliMUONVCluster *muonCluster = (AliMUONVCluster*) mftCluster.CreateMUONCluster();
trackParam.SetUniqueID(iMFTCluster); // we profit of this slot to store the reference to the corresponding MFTCluster
AliDebug(1, Form("Now adding muonCluster %p and trackParam %p",muonCluster, &trackParam));
AddTrackParamAtCluster(trackParam, *muonCluster, kTRUE);
- AliDebug(1, Form("GetTrackParamAtCluster() = %p has %d entries",GetTrackParamAtCluster(), GetTrackParamAtCluster()->GetEntries()));
// we pass the parameters this->GetTrackParamAtCluster()->First() to the Kalman Filter algorithm: they will be updated!!
Double_t chi2Kalman = RunKalmanFilter(*(AliMUONTrackParam*)(GetTrackParamAtCluster()->First()));
+ AliDebug(1, Form("Adding Kalman chi2 = %f to global chi2 = %f", chi2Kalman, GetGlobalChi2()));
Double_t newGlobalChi2 = GetGlobalChi2() + chi2Kalman;
mftCluster.SetLocalChi2(chi2Kalman);
mftCluster.SetTrackChi2(newGlobalChi2);
- new ((*fMFTClusters)[iMFTCluster]) AliMFTCluster(mftCluster);
+ new ((*(this->fMFTClusters))[iMFTCluster]) AliMFTCluster(mftCluster);
+ AliDebug(1, Form("GetTrackParamAtCluster() = %p has %d entries while this->fMFTClusters=%p has %d entries",
+ GetTrackParamAtCluster(), GetTrackParamAtCluster()->GetEntries(), this->fMFTClusters, this->fMFTClusters->GetEntries()));
AliDebug(1, Form("muonCluster->GetZ() = %f, trackParam->GetZ() = %f",muonCluster->GetZ(), trackParam.GetZ()));
SetGlobalChi2(newGlobalChi2);
+ AliDebug(1, Form("New global chi2 = %f", GetGlobalChi2()));
((AliMUONTrackParam*) GetTrackParamAtCluster()->First())->SetTrackChi2(newGlobalChi2);
for (Int_t iPar=0; iPar<GetTrackParamAtCluster()->GetEntries(); iPar++) {
AliDebug(1, Form("GetTrackParamAtCluster()->At(%d)->GetClusterPtr() = %p",
AliMFTCluster* AliMuonForwardTrack::GetMFTCluster(Int_t iMFTCluster) {
if (iMFTCluster<0 || iMFTCluster>=GetNMFTClusters()) {
- AliError("Invalid MFT cluster index. NULL pointer will be returned");
+ AliError(Form("Invalid MFT cluster index (%d). GetNMFTClusters()=%d. NULL pointer will be returned", iMFTCluster, GetNMFTClusters()));
return NULL;
}
AliMUONTrackParam *trackParam = GetTrackParamAtMFTCluster(iMFTCluster);
- AliMFTCluster *mftCluster = (AliMFTCluster*) fMFTClusters->At(trackParam->GetUniqueID());
+ AliMFTCluster *mftCluster = (AliMFTCluster*) this->fMFTClusters->At(trackParam->GetUniqueID());
return mftCluster;
//====================================================================================================================================================
+Double_t AliMuonForwardTrack::GetDCA(Double_t x, Double_t y, Double_t z) {
+
+ // Distance of Closest Approach, according to the standard MUON terminology. Actually, the offset of the track w.r.t. the primary vertex,
+ // where the extrapolation of the track DOES NOT include the MFT information
+
+ AliMUONTrackParam param = (*((AliMUONTrackParam*)(GetTrackParamAtMUONCluster(0))));
+
+ AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¶m, z);
+ Double_t dX = param.GetNonBendingCoor() - x;
+ Double_t dY = param.GetBendingCoor() - y;
+ return TMath::Sqrt(dX*dX + dY*dY);
+
+}
+
+//====================================================================================================================================================
+
+Double_t AliMuonForwardTrack::GetMomentumSpectrometer(Double_t z) {
+
+ // Momentum of the track at the primary vertex plane, where the extrapolation of the track DOES NOT include the MFT information
+
+ AliMUONTrackParam param = (*((AliMUONTrackParam*)(GetTrackParamAtMUONCluster(0))));
+
+ AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¶m, z);
+ return param.P();
+
+}
+
+//====================================================================================================================================================
+
+Double_t AliMuonForwardTrack::GetThetaAbs() {
+
+ // 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
+
+ Double_t z = AliMUONConstants::AbsZEnd();
+
+ AliMUONTrackParam param = (*((AliMUONTrackParam*)(GetTrackParamAtMUONCluster(0))));
+
+ AliMUONTrackExtrap::ExtrapToZ(¶m, z);
+ Double_t x = param.GetNonBendingCoor();
+ Double_t y = param.GetBendingCoor();
+
+ return 180. +TMath::ATan(TMath::Sqrt(x*x + y*y)/z)*TMath::RadToDeg();
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsFromDirectResonance() {
+
+ if (!IsMuon()) return kFALSE;
+
+ for (Int_t i=GetFirstMotherID(); i>=0; i--) if (!IsPDGResonance(GetParentPDGCode(i))) return kFALSE; // the decay chain finds a non-resonance particle
+
+ return kTRUE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsFromChainResonance() {
+
+ if (!IsMuon()) return kFALSE;
+
+ if (GetFirstMotherID() == 0) return kFALSE; // it is not a chain
+
+ if (!IsPDGResonance(GetParentPDGCode(GetFirstMotherID()))) return kFALSE; // primordial is not a resonance
+
+ for (Int_t i=GetFirstMotherID()-1; i>=0; i--) if (!IsPDGResonance(GetParentPDGCode(i))) return kTRUE; // the decay chain finds a non-resonance particle
+
+ return kFALSE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsFromDirectCharm() {
+
+ if (!IsMuon()) return kFALSE;
+
+ for (Int_t i=GetFirstMotherID(); i>=0; i--) if (!IsPDGCharm(GetParentPDGCode(i))) return kFALSE; // the decay chain finds a non-charmed particle
+
+ return kTRUE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsFromChainCharm() {
+
+ if (!IsMuon()) return kFALSE;
+
+ if (GetFirstMotherID() == 0) return kFALSE; // it is not a chain
+
+ if (!IsPDGCharm(GetParentPDGCode(GetFirstMotherID()))) return kFALSE; // primordial is not a charmed hadron
+
+ for (Int_t i=GetFirstMotherID()-1; i>=0; i--) if (!IsPDGCharm(GetParentPDGCode(i))) return kTRUE; // the decay chain finds a non-charmed particle
+
+ return kFALSE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsFromDirectBeauty() {
+
+ if (!IsMuon()) return kFALSE;
+
+ for (Int_t i=GetFirstMotherID(); i>=0; i--) if (!IsPDGBeauty(GetParentPDGCode(i))) return kFALSE; // the decay chain finds a non-beauty particle
+
+ return kTRUE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsFromChainBeauty() {
+
+ if (!IsMuon()) return kFALSE;
+
+ if (GetFirstMotherID() == 0) return kFALSE; // it is not a chain
+
+ if (!IsPDGBeauty(GetParentPDGCode(GetFirstMotherID()))) return kFALSE; // primordial is not a beauty hadron
+
+ for (Int_t i=GetFirstMotherID()-1; i>=0; i--) if (!IsPDGBeauty(GetParentPDGCode(i))) return kTRUE; // the decay chain finds a non-beauty particle
+
+ return kFALSE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsMuon() {
+
+ if (IsFake()) return kFALSE;
+
+ if (TMath::Abs(fMCTrackRef->GetPdgCode())==13) return kTRUE;
+
+ return kFALSE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsFake() {
+
+ if (!fMCTrackRef) return kTRUE;
+
+ return kFALSE;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsPDGResonance(Int_t pdg) {
+
+ // if (pdg<10) return kFALSE;
+ // Int_t id = pdg%100000;
+ // return (!((id-id%10)%110));
+
+ Bool_t result = kFALSE;
+
+ if ( pdg == 113 ||
+ pdg == 221 ||
+ pdg == 223 ||
+ pdg == 331 ||
+ pdg == 333 ||
+ pdg == 443 ||
+ pdg == 100443 ||
+ pdg == 553 ||
+ pdg == 100553 ) result = kTRUE;
+
+ return result;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsPDGCharm(Int_t pdg) {
+
+ Bool_t result = kFALSE;
+
+ if ( TMath::Abs(pdg) == 411 ||
+ TMath::Abs(pdg) == 421 ||
+ TMath::Abs(pdg) == 10411 ||
+ TMath::Abs(pdg) == 10421 ||
+ TMath::Abs(pdg) == 413 ||
+ TMath::Abs(pdg) == 423 ||
+ TMath::Abs(pdg) == 10413 ||
+ TMath::Abs(pdg) == 10423 ||
+ TMath::Abs(pdg) == 20413 ||
+ TMath::Abs(pdg) == 20423 ||
+ TMath::Abs(pdg) == 415 ||
+ TMath::Abs(pdg) == 425 ||
+ TMath::Abs(pdg) == 431 ||
+ TMath::Abs(pdg) == 10431 ||
+ TMath::Abs(pdg) == 433 ||
+ TMath::Abs(pdg) == 10433 ||
+ TMath::Abs(pdg) == 20433 ||
+ TMath::Abs(pdg) == 435 ) result = kTRUE;
+
+ return result;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsPDGBeauty(Int_t pdg) {
+
+ Bool_t result = kFALSE;
+
+ if ( TMath::Abs(pdg) == 511 ||
+ TMath::Abs(pdg) == 521 ||
+ TMath::Abs(pdg) == 10511 ||
+ TMath::Abs(pdg) == 10521 ||
+ TMath::Abs(pdg) == 513 ||
+ TMath::Abs(pdg) == 523 ||
+ TMath::Abs(pdg) == 10513 ||
+ TMath::Abs(pdg) == 10523 ||
+ TMath::Abs(pdg) == 20513 ||
+ TMath::Abs(pdg) == 20523 ||
+ TMath::Abs(pdg) == 515 ||
+ TMath::Abs(pdg) == 525 ||
+ TMath::Abs(pdg) == 531 ||
+ TMath::Abs(pdg) == 10531 ||
+ TMath::Abs(pdg) == 533 ||
+ TMath::Abs(pdg) == 10533 ||
+ TMath::Abs(pdg) == 20533 ||
+ TMath::Abs(pdg) == 535 ||
+ TMath::Abs(pdg) == 541 ||
+ TMath::Abs(pdg) == 10541 ||
+ TMath::Abs(pdg) == 543 ||
+ TMath::Abs(pdg) == 10543 ||
+ TMath::Abs(pdg) == 20543 ||
+ TMath::Abs(pdg) == 545 ) result = kTRUE;
+
+ return result;
+
+}
+
+//====================================================================================================================================================
+
+Bool_t AliMuonForwardTrack::IsMuonFromBackground() {
+
+ Bool_t result = kFALSE;
+
+ if (!IsMuon()) return result;
+
+ if (!IsFromDirectResonance() && !IsFromChainResonance() && !IsFromDirectCharm() && !IsFromChainCharm() && !IsFromDirectBeauty() && !IsFromChainBeauty()) result = kTRUE;
+
+ if (result) AliDebug(1, Form("Muon comes from a background source %d", GetParentPDGCode(0)));
+
+ return result;
+
+}
+
+//====================================================================================================================================================
+
+void AliMuonForwardTrack::EvalKinem(Double_t z) {
+
+ AliMUONTrackParam *param = GetTrackParamAtMFTCluster(0);
+ AliMUONTrackExtrap::ExtrapToZCov(param, z);
+
+ Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
+ Double_t energy = TMath::Sqrt(param->P()*param->P() + mMu*mMu);
+ fKinem.SetPxPyPzE(param->Px(), param->Py(), param->Pz(), energy);
+
+ TMatrixD cov(5,5);
+ fParamCovMatrix = param->GetCovariances();
+
+}
+
+//====================================================================================================================================================
+
+Int_t AliMuonForwardTrack::GetFirstMotherID() {
+
+ Int_t motherLevel = 0;
+ Int_t motherMCLabel = GetParentMCLabel(motherLevel);
+ while (motherMCLabel >= 0) {
+ motherLevel++;
+ motherMCLabel = GetParentMCLabel(motherLevel);
+ }
+ return motherLevel-1;
+
+}
+
+//====================================================================================================================================================
+
+void AliMuonForwardTrack::PrintHistory() {
+
+ if (IsFake()) printf("Track is a fake MUON\n");
+ else {
+ TString history = "";
+ for (Int_t i=GetFirstMotherID(); i>=0; i--) {
+ history += TDatabasePDG::Instance()->GetParticle(GetParentPDGCode(i))->GetName();
+ history += " -> ";
+ }
+ history += TDatabasePDG::Instance()->GetParticle(fMCTrackRef->GetPdgCode())->GetName();
+ printf("%s\n",history.Data());
+ }
+
+}
+
+//====================================================================================================================================================