X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=MUON%2FAliMUONTrack.cxx;h=066967b5c8163fe97faf7c0da0b32712d27775eb;hb=311d06910afa395419e973e632ba7527f365a3a7;hp=661271bcc9e3811a7bbc2bc42fba75d88d8cde0b;hpb=7771752ec862f31222329412e915e79197b8cea3;p=u%2Fmrichter%2FAliRoot.git diff --git a/MUON/AliMUONTrack.cxx b/MUON/AliMUONTrack.cxx index 661271bcc9e..066967b5c81 100644 --- a/MUON/AliMUONTrack.cxx +++ b/MUON/AliMUONTrack.cxx @@ -15,36 +15,31 @@ /* $Id$ */ -/////////////////////////////////////////////////// -// -// Reconstructed track -// in -// ALICE -// dimuon -// spectrometer -// -/////////////////////////////////////////////////// +//----------------------------------------------------------------------------- +// Class AliMUONTrack +//------------------- +// Reconstructed track in ALICE dimuon spectrometer +//----------------------------------------------------------------------------- #include "AliMUONTrack.h" #include "AliMUONTrackParam.h" #include "AliMUONHitForRec.h" +#include "AliMUONObjectPair.h" #include "AliMUONConstants.h" #include "AliMUONTrackExtrap.h" #include "AliLog.h" #include -#include #include +#include + /// \cond CLASSIMP ClassImp(AliMUONTrack) // Class implementation in ROOT context /// \endcond -const Double_t AliMUONTrack::fgkMaxTrackingDistanceBending = 2.; -const Double_t AliMUONTrack::fgkMaxTrackingDistanceNonBending = 2.; - //__________________________________________________________________________ AliMUONTrack::AliMUONTrack() : TObject(), @@ -52,176 +47,226 @@ AliMUONTrack::AliMUONTrack() fTrackParamAtHit(0x0), fHitForRecAtHit(0x0), fNTrackHits(0), - fExtrapTrackParam(), fFitWithVertex(kFALSE), fVertex(0x0), - fFitFMin(-1.), + fFitWithMCS(kFALSE), + fHitWeightsNonBending(0x0), + fHitWeightsBending(0x0), + fGlobalChi2(-1.), + fImproved(kFALSE), fMatchTrigger(-1), floTrgNum(-1), fChi2MatchTrigger(0.), fTrackID(0), - fHitsPatternInTrigCh(0) + fHitsPatternInTrigCh(0), + fLocalTrigger(234) { /// Default constructor } //__________________________________________________________________________ -AliMUONTrack::AliMUONTrack(AliMUONHitForRec* hitForRec1, AliMUONHitForRec* hitForRec2) +AliMUONTrack::AliMUONTrack(AliMUONObjectPair *segment) : TObject(), fTrackParamAtVertex(), fTrackParamAtHit(0x0), fHitForRecAtHit(0x0), fNTrackHits(0), - fExtrapTrackParam(), fFitWithVertex(kFALSE), fVertex(0x0), - fFitFMin(-1.), + fFitWithMCS(kFALSE), + fHitWeightsNonBending(0x0), + fHitWeightsBending(0x0), + fGlobalChi2(0.), + fImproved(kFALSE), fMatchTrigger(-1), floTrgNum(-1), fChi2MatchTrigger(0.), fTrackID(0), - fHitsPatternInTrigCh(0) + fHitsPatternInTrigCh(0), + fLocalTrigger(234) { /// Constructor from thw hitForRec's fTrackParamAtHit = new TClonesArray("AliMUONTrackParam",10); + fTrackParamAtHit->SetOwner(kTRUE); fHitForRecAtHit = new TClonesArray("AliMUONHitForRec",10); + fHitForRecAtHit->SetOwner(kTRUE); - if (!hitForRec1) return; //AZ + if (!segment) return; //AZ - // Add hits to the track - AddTrackParamAtHit(0,hitForRec1); - AddTrackParamAtHit(0,hitForRec2); + // Pointers to hits from the segment + AliMUONHitForRec* hit1 = (AliMUONHitForRec*) segment->First(); + AliMUONHitForRec* hit2 = (AliMUONHitForRec*) segment->Second(); - // sort TrackParamAtHit according to increasing -Z - fTrackParamAtHit->Sort(); + // check sorting in -Z (spectro z<0) + if (hit1->GetZ() < hit2->GetZ()) { + hit1 = hit2; + hit2 = (AliMUONHitForRec*) segment->First(); + } + + // order the hits into the track according to the station the segment belong to + //(the hit first attached is the one from which we will start the tracking procedure) + if (hit1->GetChamberNumber() == 8) { + AddTrackParamAtHit(0,hit1); + AddTrackParamAtHit(0,hit2); + } else { + AddTrackParamAtHit(0,hit2); + AddTrackParamAtHit(0,hit1); + } - // Set track parameters at first track hit AliMUONTrackParam* trackParamAtFirstHit = (AliMUONTrackParam*) fTrackParamAtHit->First(); AliMUONHitForRec* firstHit = trackParamAtFirstHit->GetHitForRecPtr(); - AliMUONHitForRec* lastHit = ((AliMUONTrackParam*) fTrackParamAtHit->Last())->GetHitForRecPtr(); + AliMUONTrackParam* trackParamAtLastHit = (AliMUONTrackParam*) fTrackParamAtHit->Last(); + AliMUONHitForRec* lastHit = trackParamAtLastHit->GetHitForRecPtr(); + + + // Compute track parameters Double_t dZ = firstHit->GetZ() - lastHit->GetZ(); // Non bending plane - Double_t nonBendingCoor = firstHit->GetNonBendingCoor(); - trackParamAtFirstHit->SetNonBendingCoor(nonBendingCoor); - trackParamAtFirstHit->SetNonBendingSlope((nonBendingCoor - lastHit->GetNonBendingCoor()) / dZ); + Double_t nonBendingCoor1 = firstHit->GetNonBendingCoor(); + Double_t nonBendingCoor2 = lastHit->GetNonBendingCoor(); + Double_t nonBendingSlope = (nonBendingCoor1 - nonBendingCoor2) / dZ; // Bending plane - Double_t bendingCoor = firstHit->GetBendingCoor(); - trackParamAtFirstHit->SetBendingCoor(bendingCoor); - Double_t bendingSlope = (bendingCoor - lastHit->GetBendingCoor()) / dZ; - trackParamAtFirstHit->SetBendingSlope(bendingSlope); + Double_t bendingCoor1 = firstHit->GetBendingCoor(); + Double_t bendingCoor2 = lastHit->GetBendingCoor(); + Double_t bendingSlope = (bendingCoor1 - bendingCoor2) / dZ; // Inverse bending momentum - Double_t bendingImpact = bendingCoor - firstHit->GetZ() * bendingSlope; + Double_t bendingImpact = bendingCoor1 - firstHit->GetZ() * bendingSlope; Double_t inverseBendingMomentum = 1. / AliMUONTrackExtrap::GetBendingMomentumFromImpactParam(bendingImpact); + + + // Set track parameters at first hit + trackParamAtFirstHit->SetNonBendingCoor(nonBendingCoor1); + trackParamAtFirstHit->SetNonBendingSlope(nonBendingSlope); + trackParamAtFirstHit->SetBendingCoor(bendingCoor1); + trackParamAtFirstHit->SetBendingSlope(bendingSlope); trackParamAtFirstHit->SetInverseBendingMomentum(inverseBendingMomentum); - // Evaluate covariances - TMatrixD *paramCov = trackParamAtFirstHit->GetCovariances(); - (*paramCov) = 0; + + // Set track parameters at last hit + trackParamAtLastHit->SetNonBendingCoor(nonBendingCoor2); + trackParamAtLastHit->SetNonBendingSlope(nonBendingSlope); + trackParamAtLastHit->SetBendingCoor(bendingCoor2); + trackParamAtLastHit->SetBendingSlope(bendingSlope); + trackParamAtLastHit->SetInverseBendingMomentum(inverseBendingMomentum); + + + // Compute and set track parameters covariances at first hit + TMatrixD paramCov1(5,5); + paramCov1.Zero(); // Non bending plane - (*paramCov)(0,0) = firstHit->GetNonBendingReso2(); - (*paramCov)(0,1) = firstHit->GetNonBendingReso2() / dZ; - (*paramCov)(1,0) = (*paramCov)(0,1); - (*paramCov)(1,1) = ( firstHit->GetNonBendingReso2() + lastHit->GetNonBendingReso2() ) / dZ / dZ; + paramCov1(0,0) = firstHit->GetNonBendingReso2(); + paramCov1(0,1) = firstHit->GetNonBendingReso2() / dZ; + paramCov1(1,0) = paramCov1(0,1); + paramCov1(1,1) = ( firstHit->GetNonBendingReso2() + lastHit->GetNonBendingReso2() ) / dZ / dZ; // Bending plane - (*paramCov)(2,2) = firstHit->GetBendingReso2(); - (*paramCov)(2,3) = firstHit->GetBendingReso2() / dZ; - (*paramCov)(3,2) = (*paramCov)(2,3); - (*paramCov)(3,3) = ( firstHit->GetBendingReso2() + lastHit->GetBendingReso2() ) / dZ / dZ; + paramCov1(2,2) = firstHit->GetBendingReso2(); + paramCov1(2,3) = firstHit->GetBendingReso2() / dZ; + paramCov1(3,2) = paramCov1(2,3); + paramCov1(3,3) = ( firstHit->GetBendingReso2() + lastHit->GetBendingReso2() ) / dZ / dZ; // Inverse bending momentum (50% error) - (*paramCov)(4,4) = 0.5*inverseBendingMomentum * 0.5*inverseBendingMomentum; + paramCov1(4,4) = 0.5*inverseBendingMomentum * 0.5*inverseBendingMomentum; + // Set covariances + trackParamAtFirstHit->SetCovariances(paramCov1); -} - - //__________________________________________________________________________ -AliMUONTrack::~AliMUONTrack() -{ - /// Destructor - if (fTrackParamAtHit) { - // delete the TClonesArray of pointers to TrackParam - delete fTrackParamAtHit; - fTrackParamAtHit = 0x0; - } - - if (fHitForRecAtHit) { - // delete the TClonesArray of pointers to HitForRec - delete fHitForRecAtHit; - fHitForRecAtHit = 0x0; - } - if (fVertex) { - // delete the vertex used during the tracking procedure - delete fVertex; - fVertex = 0x0; - } + // Compute and set track parameters covariances at last hit (as if the first hit did not exist) + TMatrixD paramCov2(5,5); + paramCov2.Zero(); + // Non bending plane + paramCov2(0,0) = paramCov1(0,0); + paramCov2(1,1) = 100.*paramCov1(1,1); + // Bending plane + paramCov2(2,2) = paramCov1(2,2); + paramCov2(3,3) = 100.*paramCov1(3,3); + // Inverse bending momentum + paramCov2(4,4) = paramCov1(4,4); + // Set covariances + trackParamAtLastHit->SetCovariances(paramCov2); + + + // Flag first hit as being removable + trackParamAtFirstHit->SetRemovable(kTRUE); + + // Flag last hit as being removable + trackParamAtLastHit->SetRemovable(kTRUE); + } //__________________________________________________________________________ -AliMUONTrack::AliMUONTrack (const AliMUONTrack& theMUONTrack) - : TObject(theMUONTrack), - fTrackParamAtVertex(theMUONTrack.fTrackParamAtVertex), +AliMUONTrack::AliMUONTrack (const AliMUONTrack& track) + : TObject(track), + fTrackParamAtVertex(track.fTrackParamAtVertex), fTrackParamAtHit(0x0), fHitForRecAtHit(0x0), - fNTrackHits(theMUONTrack.fNTrackHits), - fExtrapTrackParam(theMUONTrack.fExtrapTrackParam), - fFitWithVertex(theMUONTrack.fFitWithVertex), + fNTrackHits(track.fNTrackHits), + fFitWithVertex(track.fFitWithVertex), fVertex(0x0), - fFitFMin(theMUONTrack.fFitFMin), - fMatchTrigger(theMUONTrack.fMatchTrigger), - floTrgNum(theMUONTrack.floTrgNum), - fChi2MatchTrigger(theMUONTrack.fChi2MatchTrigger), - fTrackID(theMUONTrack.fTrackID), - fHitsPatternInTrigCh(theMUONTrack.fHitsPatternInTrigCh) + fFitWithMCS(track.fFitWithMCS), + fHitWeightsNonBending(0x0), + fHitWeightsBending(0x0), + fGlobalChi2(track.fGlobalChi2), + fImproved(track.fImproved), + fMatchTrigger(track.fMatchTrigger), + floTrgNum(track.floTrgNum), + fChi2MatchTrigger(track.fChi2MatchTrigger), + fTrackID(track.fTrackID), + fHitsPatternInTrigCh(track.fHitsPatternInTrigCh), + fLocalTrigger(track.fLocalTrigger) { ///copy constructor Int_t maxIndex = 0; // necessary to make a copy of the objects and not only the pointers in TClonesArray. - if (theMUONTrack.fTrackParamAtHit) { - maxIndex = (theMUONTrack.fTrackParamAtHit)->GetEntriesFast(); + if (track.fTrackParamAtHit) { + maxIndex = (track.fTrackParamAtHit)->GetEntriesFast(); fTrackParamAtHit = new TClonesArray("AliMUONTrackParam",maxIndex); for (Int_t index = 0; index < maxIndex; index++) { - new ((*fTrackParamAtHit)[index]) AliMUONTrackParam(*(AliMUONTrackParam*)theMUONTrack.fTrackParamAtHit->At(index)); + new ((*fTrackParamAtHit)[index]) AliMUONTrackParam(*(AliMUONTrackParam*)track.fTrackParamAtHit->At(index)); } } // necessary to make a copy of the objects and not only the pointers in TClonesArray. - if (theMUONTrack.fHitForRecAtHit) { - maxIndex = (theMUONTrack.fHitForRecAtHit)->GetEntriesFast(); + if (track.fHitForRecAtHit) { + maxIndex = (track.fHitForRecAtHit)->GetEntriesFast(); fHitForRecAtHit = new TClonesArray("AliMUONHitForRec",maxIndex); for (Int_t index = 0; index < maxIndex; index++) { - new ((*fHitForRecAtHit)[index]) AliMUONHitForRec(*(AliMUONHitForRec*)theMUONTrack.fHitForRecAtHit->At(index)); + new ((*fHitForRecAtHit)[index]) AliMUONHitForRec(*(AliMUONHitForRec*)track.fHitForRecAtHit->At(index)); } } // copy vertex used during the tracking procedure if any - if (theMUONTrack.fVertex) fVertex = new AliMUONHitForRec(*(theMUONTrack.fVertex)); + if (track.fVertex) fVertex = new AliMUONHitForRec(*(track.fVertex)); + + // copy hit weights matrices if any + if (track.fHitWeightsNonBending) fHitWeightsNonBending = new TMatrixD(*(track.fHitWeightsNonBending)); + if (track.fHitWeightsBending) fHitWeightsBending = new TMatrixD(*(track.fHitWeightsBending)); } //__________________________________________________________________________ -AliMUONTrack & AliMUONTrack::operator=(const AliMUONTrack& theMUONTrack) +AliMUONTrack & AliMUONTrack::operator=(const AliMUONTrack& track) { /// Asignment operator // check assignement to self - if (this == &theMUONTrack) + if (this == &track) return *this; // base class assignement - TObject::operator=(theMUONTrack); + TObject::operator=(track); - fTrackParamAtVertex = theMUONTrack.fTrackParamAtVertex; + fTrackParamAtVertex = track.fTrackParamAtVertex; Int_t maxIndex = 0; // necessary to make a copy of the objects and not only the pointers in TClonesArray. - if (theMUONTrack.fTrackParamAtHit) { + if (track.fTrackParamAtHit) { if (fTrackParamAtHit) fTrackParamAtHit->Clear(); else fTrackParamAtHit = new TClonesArray("AliMUONTrackParam",10); - maxIndex = (theMUONTrack.fTrackParamAtHit)->GetEntriesFast(); + maxIndex = (track.fTrackParamAtHit)->GetEntriesFast(); for (Int_t index = 0; index < maxIndex; index++) { new ((*fTrackParamAtHit)[fTrackParamAtHit->GetEntriesFast()]) - AliMUONTrackParam(*(AliMUONTrackParam*)(theMUONTrack.fTrackParamAtHit)->At(index)); + AliMUONTrackParam(*(AliMUONTrackParam*)(track.fTrackParamAtHit)->At(index)); } } else if (fTrackParamAtHit) { delete fTrackParamAtHit; @@ -229,13 +274,13 @@ AliMUONTrack & AliMUONTrack::operator=(const AliMUONTrack& theMUONTrack) } // necessary to make a copy of the objects and not only the pointers in TClonesArray. - if (theMUONTrack.fHitForRecAtHit) { + if (track.fHitForRecAtHit) { if (fHitForRecAtHit) fHitForRecAtHit->Clear(); else fHitForRecAtHit = new TClonesArray("AliMUONHitForRec",10); - maxIndex = (theMUONTrack.fHitForRecAtHit)->GetEntriesFast(); + maxIndex = (track.fHitForRecAtHit)->GetEntriesFast(); for (Int_t index = 0; index < maxIndex; index++) { new ((*fHitForRecAtHit)[fHitForRecAtHit->GetEntriesFast()]) - AliMUONHitForRec(*(AliMUONHitForRec*)(theMUONTrack.fHitForRecAtHit)->At(index)); + AliMUONHitForRec(*(AliMUONHitForRec*)(track.fHitForRecAtHit)->At(index)); } } else if (fHitForRecAtHit) { delete fHitForRecAtHit; @@ -243,30 +288,75 @@ AliMUONTrack & AliMUONTrack::operator=(const AliMUONTrack& theMUONTrack) } // copy vertex used during the tracking procedure if any. - if (theMUONTrack.fVertex) { - if (fVertex) *fVertex = *(theMUONTrack.fVertex); - else fVertex = new AliMUONHitForRec(*(theMUONTrack.fVertex)); + if (track.fVertex) { + if (fVertex) *fVertex = *(track.fVertex); + else fVertex = new AliMUONHitForRec(*(track.fVertex)); } else if (fVertex) { delete fVertex; fVertex = 0x0; } - fExtrapTrackParam = theMUONTrack.fExtrapTrackParam; + // copy hit weights matrix if any + if (track.fHitWeightsNonBending) { + if (fHitWeightsNonBending) { + fHitWeightsNonBending->ResizeTo(*(track.fHitWeightsNonBending)); + *fHitWeightsNonBending = *(track.fHitWeightsNonBending); + } else fHitWeightsNonBending = new TMatrixD(*(track.fHitWeightsNonBending)); + } else if (fHitWeightsNonBending) { + delete fHitWeightsNonBending; + fHitWeightsNonBending = 0x0; + } + + // copy hit weights matrix if any + if (track.fHitWeightsBending) { + if (fHitWeightsBending) { + fHitWeightsBending->ResizeTo(*(track.fHitWeightsBending)); + *fHitWeightsBending = *(track.fHitWeightsBending); + } else fHitWeightsBending = new TMatrixD(*(track.fHitWeightsBending)); + } else if (fHitWeightsBending) { + delete fHitWeightsBending; + fHitWeightsBending = 0x0; + } - fNTrackHits = theMUONTrack.fNTrackHits; - fFitWithVertex = theMUONTrack.fFitWithVertex; - fFitFMin = theMUONTrack.fFitFMin; - fMatchTrigger = theMUONTrack.fMatchTrigger; - floTrgNum = theMUONTrack.floTrgNum; - fChi2MatchTrigger = theMUONTrack.fChi2MatchTrigger; - fTrackID = theMUONTrack.fTrackID; - fHitsPatternInTrigCh = theMUONTrack.fHitsPatternInTrigCh; + fNTrackHits = track.fNTrackHits; + fFitWithVertex = track.fFitWithVertex; + fFitWithMCS = track.fFitWithMCS; + fGlobalChi2 = track.fGlobalChi2; + fImproved = track.fImproved; + fMatchTrigger = track.fMatchTrigger; + floTrgNum = track.floTrgNum; + fChi2MatchTrigger = track.fChi2MatchTrigger; + fTrackID = track.fTrackID; + fHitsPatternInTrigCh = track.fHitsPatternInTrigCh; + fLocalTrigger = track.fLocalTrigger; return *this; } //__________________________________________________________________________ -void AliMUONTrack::AddTrackParamAtHit(AliMUONTrackParam *trackParam, AliMUONHitForRec *hitForRec) +AliMUONTrack::~AliMUONTrack() +{ + /// Destructor + delete fTrackParamAtHit; + delete fHitForRecAtHit; + delete fVertex; + delete fHitWeightsNonBending; + delete fHitWeightsBending; +} + + //__________________________________________________________________________ +void AliMUONTrack::Clear(Option_t* opt) +{ + /// Clear arrays + if ( fTrackParamAtHit ) fTrackParamAtHit->Clear(opt); + if ( fHitForRecAtHit ) fHitForRecAtHit->Clear(opt); + delete fVertex; fVertex = 0x0; + delete fHitWeightsNonBending; fHitWeightsNonBending = 0x0; + delete fHitWeightsBending; fHitWeightsBending = 0x0; +} + + //__________________________________________________________________________ +void AliMUONTrack::AddTrackParamAtHit(const AliMUONTrackParam *trackParam, AliMUONHitForRec *hitForRec) { /// Add TrackParamAtHit if "trackParam" != NULL /// else create empty TrackParamAtHit and set the z position to the one of "hitForRec" if any @@ -284,12 +374,30 @@ void AliMUONTrack::AddTrackParamAtHit(AliMUONTrackParam *trackParam, AliMUONHitF } } else { trackParamAtHit = new ((*fTrackParamAtHit)[fNTrackHits]) AliMUONTrackParam(); - if (hitForRec) ((AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(fNTrackHits))->SetZ(hitForRec->GetZ()); + if (hitForRec) trackParamAtHit->SetZ(hitForRec->GetZ()); } if (hitForRec) trackParamAtHit->SetHitForRecPtr(hitForRec); fNTrackHits++; } + //__________________________________________________________________________ +void AliMUONTrack::RemoveTrackParamAtHit(AliMUONTrackParam *trackParam) +{ + /// Remove trackParam from the array of TrackParamAtHit + if (!fTrackParamAtHit) { + AliWarning("array fTrackParamAtHit does not exist"); + return; + } + + if (!fTrackParamAtHit->Remove(trackParam)) { + AliWarning("object to remove does not exist in array fTrackParamAtHit"); + return; + } + + fTrackParamAtHit->Compress(); + fNTrackHits--; +} + //__________________________________________________________________________ void AliMUONTrack::AddHitForRecAtHit(const AliMUONHitForRec *hitForRec) { @@ -304,13 +412,473 @@ void AliMUONTrack::AddHitForRecAtHit(const AliMUONHitForRec *hitForRec) } //__________________________________________________________________________ -void AliMUONTrack::SetVertex(AliMUONHitForRec* vertex) +void AliMUONTrack::UpdateTrackParamAtHit() +{ + /// Update track parameters at each attached hit + + if (fNTrackHits == 0) { + AliWarning("no hit attached to the track"); + return; + } + + Double_t z; + AliMUONTrackParam* startingTrackParam = (AliMUONTrackParam*) fTrackParamAtHit->First(); + AliMUONTrackParam* trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->After(startingTrackParam); + while (trackParamAtHit) { + + // save current z + z = trackParamAtHit->GetZ(); + + // reset track parameters and their covariances + trackParamAtHit->SetParameters(startingTrackParam->GetParameters()); + trackParamAtHit->SetZ(startingTrackParam->GetZ()); + + // extrapolation to the given z + AliMUONTrackExtrap::ExtrapToZ(trackParamAtHit, z); + + // prepare next step + startingTrackParam = trackParamAtHit; + trackParamAtHit = (AliMUONTrackParam*) (fTrackParamAtHit->After(trackParamAtHit)); + } + +} + + //__________________________________________________________________________ +void AliMUONTrack::UpdateCovTrackParamAtHit() +{ + /// Update track parameters and their covariances at each attached hit + + if (fNTrackHits == 0) { + AliWarning("no hit attached to the track"); + return; + } + + Double_t z; + AliMUONTrackParam* startingTrackParam = (AliMUONTrackParam*) fTrackParamAtHit->First(); + AliMUONTrackParam* trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->After(startingTrackParam); + while (trackParamAtHit) { + + // save current z + z = trackParamAtHit->GetZ(); + + // reset track parameters and their covariances + trackParamAtHit->SetParameters(startingTrackParam->GetParameters()); + trackParamAtHit->SetZ(startingTrackParam->GetZ()); + trackParamAtHit->SetCovariances(startingTrackParam->GetCovariances()); + + // extrapolation to the given z + AliMUONTrackExtrap::ExtrapToZCov(trackParamAtHit, z); + + // prepare next step + startingTrackParam = trackParamAtHit; + trackParamAtHit = (AliMUONTrackParam*) (fTrackParamAtHit->After(trackParamAtHit)); + } + +} + + //__________________________________________________________________________ +void AliMUONTrack::SetVertex(const AliMUONHitForRec* vertex) { /// Set the vertex used during the tracking procedure if (!fVertex) fVertex = new AliMUONHitForRec(*vertex); else *fVertex = *vertex; } + + //__________________________________________________________________________ +Bool_t AliMUONTrack::ComputeLocalChi2(Bool_t accountForMCS) +{ + /// Compute the removable hit contribution to the chi2 of the track + /// accounting for multiple scattering or not according to the flag + /// - Also recompute the weight matrices of the attached hits if accountForMCS=kTRUE + /// - Assume that track parameters at each hit are corrects + /// - Return kFALSE if computation failed + + // Check hits (if the first one exist, assume that the other ones exit too!) + AliMUONTrackParam* trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->First(); + if (!trackParamAtHit->GetHitForRecPtr()) { + AliWarning("hit is missing"); + return kFALSE; + } + + if (accountForMCS) { // Compute local chi2 taking into account multiple scattering effects + + // Compute MCS covariance matrix only once + TMatrixD mcsCovariances(fNTrackHits,fNTrackHits); + ComputeMCSCovariances(mcsCovariances); + + // Make sure hit weights are consistent with following calculations + if (!ComputeHitWeights(&mcsCovariances)) { + AliWarning("cannot take into account the multiple scattering effects"); + return ComputeLocalChi2(kFALSE); + } + + // Compute chi2 of the track + Double_t globalChi2 = ComputeGlobalChi2(kTRUE); + if (globalChi2 < 0.) return kFALSE; + + // Loop over removable hits and compute their local chi2 + AliMUONTrackParam* trackParamAtHit1; + AliMUONHitForRec *hitForRec, *discardedHit; + Int_t hitNumber1, hitNumber2, currentHitNumber1, currentHitNumber2; + TMatrixD hitWeightsNB(fNTrackHits-1,fNTrackHits-1); + TMatrixD hitWeightsB(fNTrackHits-1,fNTrackHits-1); + Double_t *dX = new Double_t[fNTrackHits-1]; + Double_t *dY = new Double_t[fNTrackHits-1]; + Double_t globalChi2b; + while (trackParamAtHit) { + + discardedHit = trackParamAtHit->GetHitForRecPtr(); + + // Recompute hit weights without the current hit + if (!ComputeHitWeights(hitWeightsNB, hitWeightsB, &mcsCovariances, discardedHit)) { + AliWarning("cannot take into account the multiple scattering effects"); + ComputeLocalChi2(kFALSE); + } + + // Compute track chi2 without the current hit + globalChi2b = 0.; + currentHitNumber1 = 0; + for (hitNumber1 = 0; hitNumber1 < fNTrackHits ; hitNumber1++) { + trackParamAtHit1 = (AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(hitNumber1); + hitForRec = trackParamAtHit1->GetHitForRecPtr(); + + if (hitForRec == discardedHit) continue; + + // Compute and save residuals + dX[currentHitNumber1] = hitForRec->GetNonBendingCoor() - trackParamAtHit1->GetNonBendingCoor(); + dY[currentHitNumber1] = hitForRec->GetBendingCoor() - trackParamAtHit1->GetBendingCoor(); + + currentHitNumber2 = 0; + for (hitNumber2 = 0; hitNumber2 < hitNumber1; hitNumber2++) { + hitForRec = ((AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(hitNumber2))->GetHitForRecPtr(); + + if (hitForRec == discardedHit) continue; + + // Add contribution from covariances + globalChi2b += (hitWeightsNB(currentHitNumber1, currentHitNumber2) + + hitWeightsNB(currentHitNumber2, currentHitNumber1)) * dX[currentHitNumber1] * dX[currentHitNumber2] + + (hitWeightsB(currentHitNumber1, currentHitNumber2) + + hitWeightsB(currentHitNumber2, currentHitNumber1)) * dY[currentHitNumber1] * dY[currentHitNumber2]; + + currentHitNumber2++; + } + + // Add contribution from variances + globalChi2b += hitWeightsNB(currentHitNumber1, currentHitNumber1) * dX[currentHitNumber1] * dX[currentHitNumber1] + + hitWeightsB(currentHitNumber1, currentHitNumber1) * dY[currentHitNumber1] * dY[currentHitNumber1]; + + currentHitNumber1++; + } + + // Set local chi2 + trackParamAtHit->SetLocalChi2(globalChi2 - globalChi2b); + + trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->After(trackParamAtHit); + } + + delete [] dX; + delete [] dY; + + } else { // without multiple scattering effects + + AliMUONHitForRec *discardedHit; + Double_t dX, dY; + while (trackParamAtHit) { + + discardedHit = trackParamAtHit->GetHitForRecPtr(); + + // Compute residuals + dX = discardedHit->GetNonBendingCoor() - trackParamAtHit->GetNonBendingCoor(); + dY = discardedHit->GetBendingCoor() - trackParamAtHit->GetBendingCoor(); + + // Set local chi2 + trackParamAtHit->SetLocalChi2(dX * dX / discardedHit->GetNonBendingReso2() + dY * dY / discardedHit->GetBendingReso2()); + + trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->After(trackParamAtHit); + } + + } + + return kTRUE; + +} + + //__________________________________________________________________________ +Double_t AliMUONTrack::ComputeGlobalChi2(Bool_t accountForMCS) +{ + /// Compute the chi2 of the track accounting for multiple scattering or not according to the flag + /// - Assume that track parameters at each hit are corrects + /// - Assume the hits weights matrices are corrects + /// - Return negative value if chi2 computation failed + + // Check hits (if the first one exist, assume that the other ones exit too!) + AliMUONTrackParam* trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->First(); + if (!trackParamAtHit->GetHitForRecPtr()) { + AliWarning("hit is missing"); + return -1.; + } + + Double_t chi2 = 0.; + + if (accountForMCS) { + + // Check the weight matrices. If weight matrices are not available compute chi2 without MCS + if (!fHitWeightsNonBending || !fHitWeightsBending) { + AliWarning("hit weights including multiple scattering effects are not available\n\t\t --> compute chi2 WITHOUT multiple scattering"); + return ComputeGlobalChi2(kFALSE); + } + if (fHitWeightsNonBending->GetNrows() != fNTrackHits || fHitWeightsBending->GetNcols() != fNTrackHits) { + AliWarning("hit weights including multiple scattering effects are not available\n\t\t --> compute chi2 WITHOUT multiple scattering"); + return ComputeGlobalChi2(kFALSE); + } + + // Compute chi2 + AliMUONHitForRec *hitForRec; + Double_t *dX = new Double_t[fNTrackHits]; + Double_t *dY = new Double_t[fNTrackHits]; + Int_t hitNumber1, hitNumber2; + for (hitNumber1 = 0; hitNumber1 < fNTrackHits ; hitNumber1++) { + trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(hitNumber1); + hitForRec = trackParamAtHit->GetHitForRecPtr(); + dX[hitNumber1] = hitForRec->GetNonBendingCoor() - trackParamAtHit->GetNonBendingCoor(); + dY[hitNumber1] = hitForRec->GetBendingCoor() - trackParamAtHit->GetBendingCoor(); + for (hitNumber2 = 0; hitNumber2 < hitNumber1; hitNumber2++) { + chi2 += ((*fHitWeightsNonBending)(hitNumber1, hitNumber2) + (*fHitWeightsNonBending)(hitNumber2, hitNumber1)) * dX[hitNumber1] * dX[hitNumber2] + + ((*fHitWeightsBending)(hitNumber1, hitNumber2) + (*fHitWeightsBending)(hitNumber2, hitNumber1)) * dY[hitNumber1] * dY[hitNumber2]; + } + chi2 += ((*fHitWeightsNonBending)(hitNumber1, hitNumber1) * dX[hitNumber1] * dX[hitNumber1]) + + ((*fHitWeightsBending)(hitNumber1, hitNumber1) * dY[hitNumber1] * dY[hitNumber1]); + } + delete [] dX; + delete [] dY; + + } else { + + AliMUONHitForRec *hitForRec; + Double_t dX, dY; + for (Int_t hitNumber = 0; hitNumber < fNTrackHits ; hitNumber++) { + trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(hitNumber); + hitForRec = trackParamAtHit->GetHitForRecPtr(); + dX = hitForRec->GetNonBendingCoor() - trackParamAtHit->GetNonBendingCoor(); + dY = hitForRec->GetBendingCoor() - trackParamAtHit->GetBendingCoor(); + chi2 += dX * dX / hitForRec->GetNonBendingReso2() + dY * dY / hitForRec->GetBendingReso2(); + } + + } + + return chi2; + +} + + //__________________________________________________________________________ +Bool_t AliMUONTrack::ComputeHitWeights(TMatrixD* mcsCovariances) +{ + /// Compute the weight matrices of the attached hits, in non bending and bending direction, + /// accounting for multiple scattering correlations and hits resolution + /// - Use the provided MCS covariance matrix if any (otherwise build it temporarily) + /// - Assume that track parameters at each hit are corrects + /// - Return kFALSE if computation failed + + // Alocate memory + if (!fHitWeightsNonBending) fHitWeightsNonBending = new TMatrixD(fNTrackHits,fNTrackHits); + if (!fHitWeightsBending) fHitWeightsBending = new TMatrixD(fNTrackHits,fNTrackHits); + + // Check hits (if the first one exist, assume that the other ones exit too!) + if (!((AliMUONTrackParam*) fTrackParamAtHit->First())->GetHitForRecPtr()) { + AliWarning("hit is missing"); + fHitWeightsNonBending->ResizeTo(0,0); + fHitWeightsBending->ResizeTo(0,0); + return kFALSE; + } + + // Compute weights matrices + if (!ComputeHitWeights(*fHitWeightsNonBending, *fHitWeightsBending, mcsCovariances)) return kFALSE; + + return kTRUE; + +} + + //__________________________________________________________________________ +Bool_t AliMUONTrack::ComputeHitWeights(TMatrixD& hitWeightsNB, TMatrixD& hitWeightsB, TMatrixD* mcsCovariances, AliMUONHitForRec* discardedHit) const +{ + /// Compute the weight matrices, in non bending and bending direction, + /// of the other attached hits assuming the discarded one does not exist + /// accounting for multiple scattering correlations and hits resolution + /// - Use the provided MCS covariance matrix if any (otherwise build it temporarily) + /// - Return kFALSE if computation failed + + // Check MCS covariance matrix and recompute it if need + Bool_t deleteMCSCov = kFALSE; + if (!mcsCovariances) { + mcsCovariances = new TMatrixD(fNTrackHits,fNTrackHits); + deleteMCSCov = kTRUE; + ComputeMCSCovariances(*mcsCovariances); + } + + // Resize the weights matrices; alocate memory + if (discardedHit) { + hitWeightsNB.ResizeTo(fNTrackHits-1,fNTrackHits-1); + hitWeightsB.ResizeTo(fNTrackHits-1,fNTrackHits-1); + } else { + hitWeightsNB.ResizeTo(fNTrackHits,fNTrackHits); + hitWeightsB.ResizeTo(fNTrackHits,fNTrackHits); + } + + // Define variables + AliMUONHitForRec *hitForRec1, *hitForRec2; + Int_t currentHitNumber1, currentHitNumber2; + + // Compute the covariance matrices + currentHitNumber1 = 0; + for (Int_t hitNumber1 = 0; hitNumber1 < fNTrackHits; hitNumber1++) { + hitForRec1 = ((AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(hitNumber1))->GetHitForRecPtr(); + + if (hitForRec1 == discardedHit) continue; + + // Loop over next hits + currentHitNumber2 = currentHitNumber1; + for (Int_t hitNumber2 = hitNumber1; hitNumber2 < fNTrackHits; hitNumber2++) { + hitForRec2 = ((AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(hitNumber2))->GetHitForRecPtr(); + + if (hitForRec2 == discardedHit) continue; + + // Fill with MCS covariances + hitWeightsNB(currentHitNumber1, currentHitNumber2) = (*mcsCovariances)(hitNumber1,hitNumber2); + + // Equal contribution from multiple scattering in non bending and bending directions + hitWeightsB(currentHitNumber1, currentHitNumber2) = hitWeightsNB(currentHitNumber1, currentHitNumber2); + + // Add contribution from hit resolution to diagonal element and symmetrize the matrix + if (currentHitNumber1 == currentHitNumber2) { + + // In non bending plane + hitWeightsNB(currentHitNumber1, currentHitNumber1) += hitForRec1->GetNonBendingReso2(); + // In bending plane + hitWeightsB(currentHitNumber1, currentHitNumber1) += hitForRec1->GetBendingReso2(); + + } else { + + // In non bending plane + hitWeightsNB(currentHitNumber2, currentHitNumber1) = hitWeightsNB(currentHitNumber1, currentHitNumber2); + // In bending plane + hitWeightsB(currentHitNumber2, currentHitNumber1) = hitWeightsB(currentHitNumber1, currentHitNumber2); + + } + + currentHitNumber2++; + } + + currentHitNumber1++; + } + + // Inversion of covariance matrices to get the weights + if (hitWeightsNB.Determinant() != 0 && hitWeightsB.Determinant() != 0) { + hitWeightsNB.Invert(); + hitWeightsB.Invert(); + } else { + AliWarning(" Determinant = 0"); + hitWeightsNB.ResizeTo(0,0); + hitWeightsB.ResizeTo(0,0); + if(deleteMCSCov) delete mcsCovariances; + return kFALSE; + } + + if(deleteMCSCov) delete mcsCovariances; + + return kTRUE; + +} + + //__________________________________________________________________________ +void AliMUONTrack::ComputeMCSCovariances(TMatrixD& mcsCovariances) const +{ + /// Compute the multiple scattering covariance matrix + /// (assume that track parameters at each hit are corrects) + + // Reset the size of the covariance matrix if needed + if (mcsCovariances.GetNrows() != fNTrackHits) mcsCovariances.ResizeTo(fNTrackHits,fNTrackHits); + + // Define variables + Int_t nChambers = AliMUONConstants::NTrackingCh(); + AliMUONTrackParam* trackParamAtHit; + AliMUONHitForRec *hitForRec; + AliMUONTrackParam extrapTrackParam; + Int_t currentChamber = 0, expectedChamber = 0, size = 0; + Double_t *mcsAngle2 = new Double_t[2*nChambers]; + Double_t *zMCS = new Double_t[2*nChambers]; + Int_t *indices = new Int_t[2*fNTrackHits]; + + // Compute multiple scattering dispersion angle at each chamber + // and save the z position where it is calculated + for (Int_t hitNumber = 0; hitNumber < fNTrackHits; hitNumber++) { + trackParamAtHit = (AliMUONTrackParam*) fTrackParamAtHit->UncheckedAt(hitNumber); + hitForRec = trackParamAtHit->GetHitForRecPtr(); + + // look for missing chambers if any + currentChamber = hitForRec->GetChamberNumber(); + while (currentChamber > expectedChamber) { + + // Save the z position where MCS dispersion is calculated + zMCS[size] = AliMUONConstants::DefaultChamberZ(expectedChamber); + + // Do not take into account MCS in chambers prior the first hit + if (hitNumber > 0) { + + // Get track parameters at missing chamber z + extrapTrackParam = *trackParamAtHit; + AliMUONTrackExtrap::ExtrapToZ(&extrapTrackParam, zMCS[expectedChamber]); + + // Save multiple scattering dispersion angle in missing chamber + mcsAngle2[size] = AliMUONTrackExtrap::GetMCSAngle2(extrapTrackParam,AliMUONConstants::ChamberThicknessInX0(),1.); + + } else mcsAngle2[size] = 0.; + + expectedChamber++; + size++; + } + + // Save z position where MCS dispersion is calculated + zMCS[size] = trackParamAtHit->GetZ(); + + // Save multiple scattering dispersion angle in current chamber + mcsAngle2[size] = AliMUONTrackExtrap::GetMCSAngle2(*trackParamAtHit,AliMUONConstants::ChamberThicknessInX0(),1.); + + // Save indice in zMCS array corresponding to the current cluster + indices[hitNumber] = size; + + expectedChamber = currentChamber + 1; + size++; + } + + // complete array of z if last hit is on the last but one chamber + if (currentChamber != nChambers-1) zMCS[size++] = AliMUONConstants::DefaultChamberZ(nChambers-1); + + // Compute the covariance matrix + for (Int_t hitNumber1 = 0; hitNumber1 < fNTrackHits; hitNumber1++) { + + for (Int_t hitNumber2 = hitNumber1; hitNumber2 < fNTrackHits; hitNumber2++) { + + // Initialization to 0 (diagonal plus upper triangular part) + mcsCovariances(hitNumber1,hitNumber2) = 0.; + + // Compute contribution from multiple scattering in upstream chambers + for (Int_t k = 0; k < indices[hitNumber1]; k++) { + mcsCovariances(hitNumber1,hitNumber2) += (zMCS[indices[hitNumber1]] - zMCS[k]) * (zMCS[indices[hitNumber2]] - zMCS[k]) * mcsAngle2[k]; + } + + // Symetrize the matrix + mcsCovariances(hitNumber2,hitNumber1) = mcsCovariances(hitNumber1,hitNumber2); + } + + } + + delete [] mcsAngle2; + delete [] zMCS; + delete [] indices; + +} + //__________________________________________________________________________ Int_t AliMUONTrack::HitsInCommon(AliMUONTrack* track) const { @@ -336,6 +904,16 @@ Int_t AliMUONTrack::HitsInCommon(AliMUONTrack* track) const return hitsInCommon; } + //__________________________________________________________________________ +Double_t AliMUONTrack::GetNormalizedChi2() const +{ + /// return the chi2 value divided by the number of degrees of freedom (or 1.e10 if ndf < 0) + + Double_t numberOfDegFree = (2. * fNTrackHits - 5.); + if (numberOfDegFree > 0.) return fGlobalChi2 / numberOfDegFree; + else return 1.e10; +} + //__________________________________________________________________________ Bool_t* AliMUONTrack::CompatibleTrack(AliMUONTrack * track, Double_t sigma2Cut) const { @@ -373,156 +951,6 @@ Bool_t* AliMUONTrack::CompatibleTrack(AliMUONTrack * track, Double_t sigma2Cut) } return nCompHit; -} - - //__________________________________________________________________________ -Double_t AliMUONTrack::TryOneHitForRec(AliMUONHitForRec* hitForRec) -{ -/// Test the compatibility between the track and the hitForRec: -/// return the corresponding Chi2 - - // Get track parameters and their covariances at the z position of hitForRec - AliMUONTrackParam extrapTrackParam(fExtrapTrackParam); - AliMUONTrackExtrap::ExtrapToZCov(&extrapTrackParam, hitForRec->GetZ()); - - // Set differences between trackParam and hitForRec in the bending and non bending directions - TMatrixD dPos(2,1); - dPos(0,0) = hitForRec->GetNonBendingCoor() - extrapTrackParam.GetNonBendingCoor(); - dPos(1,0) = hitForRec->GetBendingCoor() - extrapTrackParam.GetBendingCoor(); - - // quick test of hitForRec compatibility within a wide road of x*y = 10*1 cm2 to save computing time - if (TMath::Abs(dPos(0,0)) > fgkMaxTrackingDistanceNonBending || - TMath::Abs(dPos(1,0)) > fgkMaxTrackingDistanceBending) return 1.e10; - - // Set the error matrix from trackParam covariances and hitForRec resolution - TMatrixD* paramCov = extrapTrackParam.GetCovariances(); - TMatrixD error(2,2); - error(0,0) = (*paramCov)(0,0) + hitForRec->GetNonBendingReso2(); - error(0,1) = (*paramCov)(0,2); - error(1,0) = (*paramCov)(2,0); - error(1,1) = (*paramCov)(2,2) + hitForRec->GetBendingReso2(); - - // Invert the error matrix for Chi2 calculation - if (error.Determinant() != 0) { - error.Invert(); - } else { - AliWarning(" Determinant error=0"); - return 1.e10; - } - - // Compute the Chi2 value - TMatrixD tmp(error,TMatrixD::kMult,dPos); - TMatrixD result(dPos,TMatrixD::kTransposeMult,tmp); - - return result(0,0); - -} - - //__________________________________________________________________________ -Double_t AliMUONTrack::TryTwoHitForRec(AliMUONHitForRec* hitForRec1, AliMUONHitForRec* hitForRec2) -{ -/// Test the compatibility between the track and the 2 hitForRec together: -/// return the corresponding Chi2 accounting for covariances between the 2 hitForRec - - // Get track parameters and their covariances at the z position of the first hitForRec - AliMUONTrackParam extrapTrackParam1(fExtrapTrackParam); - AliMUONTrackExtrap::ExtrapToZCov(&extrapTrackParam1, hitForRec1->GetZ()); - - // Get track parameters at second hitForRec - AliMUONTrackParam extrapTrackParam2(extrapTrackParam1); - AliMUONTrackExtrap::ExtrapToZ(&extrapTrackParam2, hitForRec2->GetZ()); - - // Set differences between track and the 2 hitForRec in the bending and non bending directions - TMatrixD dPos(4,1); - dPos(0,0) = hitForRec1->GetNonBendingCoor() - extrapTrackParam1.GetNonBendingCoor(); - dPos(1,0) = hitForRec1->GetBendingCoor() - extrapTrackParam1.GetBendingCoor(); - dPos(2,0) = hitForRec2->GetNonBendingCoor() - extrapTrackParam2.GetNonBendingCoor(); - dPos(3,0) = hitForRec2->GetBendingCoor() - extrapTrackParam2.GetBendingCoor(); - - // quick tests of hitForRec compatibility within a wide road of x*y = 1*1 cm2 to save computing time - if (TMath::Abs(dPos(0,0)) > fgkMaxTrackingDistanceNonBending || - TMath::Abs(dPos(1,0)) > fgkMaxTrackingDistanceBending || - TMath::Abs(dPos(2,0)) > fgkMaxTrackingDistanceNonBending || - TMath::Abs(dPos(3,0)) > fgkMaxTrackingDistanceBending) return 1.e10; - - // Calculate the error matrix from the track parameter covariances at first hitForRec - TMatrixD error(4,4); - error = 0.; - if (extrapTrackParam1.CovariancesExist()) { - // Get the pointer to the parameter covariance matrix at first hitForRec - TMatrixD* paramCov = extrapTrackParam1.GetCovariances(); - - // Save track parameters at first hitForRec - AliMUONTrackParam extrapTrackParam1Save(extrapTrackParam1); - Double_t nonBendingCoor1 = extrapTrackParam1Save.GetNonBendingCoor(); - Double_t nonBendingSlope1 = extrapTrackParam1Save.GetNonBendingSlope(); - Double_t bendingCoor1 = extrapTrackParam1Save.GetBendingCoor(); - Double_t bendingSlope1 = extrapTrackParam1Save.GetBendingSlope(); - Double_t inverseBendingMomentum1 = extrapTrackParam1Save.GetInverseBendingMomentum(); - Double_t z1 = extrapTrackParam1Save.GetZ(); - - // Save track coordinates at second hitForRec - Double_t nonBendingCoor2 = extrapTrackParam2.GetNonBendingCoor(); - Double_t bendingCoor2 = extrapTrackParam2.GetBendingCoor(); - - // Calculate the jacobian related to the transformation between track parameters - // at first hitForRec and track coordinates at the 2 hitForRec z-position - TMatrixD jacob(4,5); - jacob = 0.; - // first derivative at the first hitForRec: - jacob(0,0) = 1.; // dx1/dx - jacob(1,2) = 1.; // dy1/dy - // first derivative at the second hitForRec: - Double_t dParam[5]; - for (Int_t i=0; i<5; i++) { - // Skip jacobian calculation for parameters with no associated error - if ((*paramCov)(i,i) == 0.) continue; - // Small variation of parameter i only - for (Int_t j=0; j<5; j++) { - if (j==i) { - dParam[j] = TMath::Sqrt((*paramCov)(i,i)); - if (j == 4) dParam[j] *= TMath::Sign(1.,-inverseBendingMomentum1); // variation always in the same direction - } else dParam[j] = 0.; - } - // Set new track parameters at first hitForRec - extrapTrackParam1Save.SetNonBendingCoor (nonBendingCoor1 + dParam[0]); - extrapTrackParam1Save.SetNonBendingSlope (nonBendingSlope1 + dParam[1]); - extrapTrackParam1Save.SetBendingCoor (bendingCoor1 + dParam[2]); - extrapTrackParam1Save.SetBendingSlope (bendingSlope1 + dParam[3]); - extrapTrackParam1Save.SetInverseBendingMomentum(inverseBendingMomentum1 + dParam[4]); - extrapTrackParam1Save.SetZ (z1); - // Extrapolate new track parameters to the z position of the second hitForRec - AliMUONTrackExtrap::ExtrapToZ(&extrapTrackParam1Save,hitForRec2->GetZ()); - // Calculate the jacobian - jacob(2,i) = (extrapTrackParam1Save.GetNonBendingCoor() - nonBendingCoor2) / dParam[i]; // dx2/dParami - jacob(3,i) = (extrapTrackParam1Save.GetBendingCoor() - bendingCoor2 ) / dParam[i]; // dy2/dParami - } - - // Calculate the error matrix - TMatrixD tmp((*paramCov),TMatrixD::kMultTranspose,jacob); - error = TMatrixD(jacob,TMatrixD::kMult,tmp); - } - - // Add hitForRec resolution to the error matrix - error(0,0) += hitForRec1->GetNonBendingReso2(); - error(1,1) += hitForRec1->GetBendingReso2(); - error(2,2) += hitForRec2->GetNonBendingReso2(); - error(3,3) += hitForRec2->GetBendingReso2(); - - // invert the error matrix for Chi2 calculation - if (error.Determinant() != 0) { - error.Invert(); - } else { - AliWarning(" Determinant error=0"); - return 1.e10; - } - - // Compute the Chi2 value - TMatrixD tmp2(error,TMatrixD::kMult,dPos); - TMatrixD result(dPos,TMatrixD::kTransposeMult,tmp2); - - return result(0,0); - } //__________________________________________________________________________ @@ -548,27 +976,32 @@ void AliMUONTrack::RecursiveDump(void) const } //_____________________________________________- -void AliMUONTrack::Print(Option_t* opt) const +void AliMUONTrack::Print(Option_t*) const { /// Printing Track information - /// "full" option for printing all the information about the track - TString sopt(opt); - sopt.ToUpper(); - - if ( sopt.Contains("FULL") ) { - cout << " No.Clusters=" << setw(2) << GetNTrackHits() << - // ", Bending P="<< setw(8) << setprecision(5) << 1./GetInverseBendingMomentum() << - //", NonBendSlope=" << setw(8) << setprecision(5) << GetNonBendingSlope()*180./TMath::Pi() << - //", BendSlope=" << setw(8) << setprecision(5) << GetBendingSlope()*180./TMath::Pi() << + + cout << " No.Clusters=" << setw(2) << GetNTrackHits() << ", Match2Trig=" << setw(1) << GetMatchTrigger() << ", LoTrgNum=" << setw(3) << GetLoTrgNum() << - ", Chi2-tracking-trigger=" << setw(8) << setprecision(5) << GetChi2MatchTrigger() << endl ; - GetTrackParamAtHit()->First()->Print("full"); - } - else { - cout << ""; - GetTrackParamAtHit()->First()->Print(""); + ", Chi2-tracking-trigger=" << setw(8) << setprecision(5) << GetChi2MatchTrigger(); + cout << Form(" HitTriggerPattern %x",fHitsPatternInTrigCh) << endl; + GetTrackParamAtHit()->First()->Print("FULL"); +} + +//__________________________________________________________________________ +void AliMUONTrack::SetLocalTrigger(Int_t loCirc, Int_t loStripX, Int_t loStripY, Int_t loDev, Int_t loLpt, Int_t loHpt) +{ + /// pack the local trigger information and store + + if (loCirc < 0 || loCirc > 233) return; + + fLocalTrigger = 0; + fLocalTrigger += loCirc; + fLocalTrigger += loStripX << 8; + fLocalTrigger += loStripY << 13; + fLocalTrigger += loDev << 17; + fLocalTrigger += loLpt << 22; + fLocalTrigger += loHpt << 24; - } - } +