X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;ds=sidebyside;f=MUON%2FAliMUONVTrackReconstructor.cxx;h=114b73d973455b988c66e310fdc2d95da5f50e8b;hb=a36106b26ac2ee368dce3bcdc773f5a01e12a0eb;hp=4367aa9ff13b1da0e0d41de6749fa9507f4965a0;hpb=3fe6651e174a9415758fafc2f9544e6aabaf4f32;p=u%2Fmrichter%2FAliRoot.git diff --git a/MUON/AliMUONVTrackReconstructor.cxx b/MUON/AliMUONVTrackReconstructor.cxx index 4367aa9ff13..114b73d9734 100644 --- a/MUON/AliMUONVTrackReconstructor.cxx +++ b/MUON/AliMUONVTrackReconstructor.cxx @@ -74,10 +74,19 @@ #include "AliMUONVTriggerStore.h" #include "AliMUONVTriggerTrackStore.h" #include "AliMUONRecoParam.h" +#include "AliMUONGeometryTransformer.h" +#include "AliMUONVDigit.h" #include "AliMpDEManager.h" #include "AliMpArea.h" +#include "AliMpDDLStore.h" +#include "AliMpVSegmentation.h" +#include "AliMpSegmentation.h" +#include "AliMpPad.h" +#include "AliMpDetElement.h" +#include "AliMpCathodType.h" + #include "AliLog.h" #include "AliCodeTimer.h" #include "AliTracker.h" @@ -89,19 +98,23 @@ #include +using std::cout; +using std::endl; /// \cond CLASSIMP ClassImp(AliMUONVTrackReconstructor) // Class implementation in ROOT context /// \endcond //__________________________________________________________________________ AliMUONVTrackReconstructor::AliMUONVTrackReconstructor(const AliMUONRecoParam* recoParam, - AliMUONVClusterServer* clusterServer) + AliMUONVClusterServer* clusterServer, + const AliMUONGeometryTransformer* transformer) : TObject(), fRecTracksPtr(0x0), fNRecTracks(0), fClusterServer(clusterServer), fkRecoParam(recoParam), -fMaxMCSAngle2(0.) +fkTransformer(transformer), +fMaxMCSAngle2(0x0) { /// Constructor for class AliMUONVTrackReconstructor /// WARNING: if clusterServer=0x0, no clusterization will be possible at this level @@ -116,7 +129,9 @@ fMaxMCSAngle2(0.) AliMUONTrackParam param; Double_t inverseBendingP = (GetRecoParam()->GetMinBendingMomentum() > 0.) ? 1./GetRecoParam()->GetMinBendingMomentum() : 1.; param.SetInverseBendingMomentum(inverseBendingP); - fMaxMCSAngle2 = AliMUONTrackExtrap::GetMCSAngle2(param, AliMUONConstants::ChamberThicknessInX0(), 1.); + fMaxMCSAngle2 = new Double_t [AliMUONConstants::NTrackingCh()]; + for (Int_t iCh=0; iChImproveTracks()) ImproveTracks(); // Remove connected tracks - if (GetRecoParam()->RemoveConnectedTracksInSt12()) RemoveConnectedTracks(kFALSE); - else RemoveConnectedTracks(kTRUE); + RemoveConnectedTracks(3, 4, kFALSE); + RemoveConnectedTracks(2, 2, kFALSE); + if (GetRecoParam()->RemoveConnectedTracksInSt12()) RemoveConnectedTracks(0, 1, kFALSE); // Fill AliMUONTrack data members Finalize(); + if (!GetRecoParam()->RemoveConnectedTracksInSt12()) TagConnectedTracks(0, 1, kTRUE); + + // Make sure there is no bad track left + RemoveBadTracks(); + + // Refit the reconstructed tracks with a different resolution for mono-cathod clusters + if (GetRecoParam()->DiscardMonoCathodClusters()) DiscardMonoCathodClusters(); // Add tracks to MUON data container - for (Int_t i=0; iAt(i); track->SetUniqueID(i+1); trackStore.Add(*track); } + } //__________________________________________________________________________ @@ -196,18 +221,31 @@ Bool_t AliMUONVTrackReconstructor::IsAcceptable(AliMUONTrackParam &trackParam) // MCS dipersion Double_t angleMCS2 = 0.; - if (AliMUONTrackExtrap::IsFieldON() && chamber < 6) - angleMCS2 = AliMUONTrackExtrap::GetMCSAngle2(trackParam, AliMUONConstants::ChamberThicknessInX0(), 1.); - else angleMCS2 = fMaxMCSAngle2; Double_t impactMCS2 = 0.; - if (!GetRecoParam()->SelectOnTrackSlope()) for (Int_t iCh=0; iCh<=chamber; iCh++) - impactMCS2 += AliMUONConstants::DefaultChamberZ(chamber) * AliMUONConstants::DefaultChamberZ(chamber) * angleMCS2; + if (AliMUONTrackExtrap::IsFieldON() && chamber < 6) { + + // track momentum is known + for (Int_t iCh=0; iCh<=chamber; iCh++) { + Double_t localMCS2 = AliMUONTrackExtrap::GetMCSAngle2(trackParam, AliMUONConstants::ChamberThicknessInX0(iCh), 1.); + angleMCS2 += localMCS2; + impactMCS2 += AliMUONConstants::DefaultChamberZ(chamber) * AliMUONConstants::DefaultChamberZ(chamber) * localMCS2; + } + + } else { + + // track momentum is unknown + for (Int_t iCh=0; iCh<=chamber; iCh++) { + angleMCS2 += fMaxMCSAngle2[iCh]; + impactMCS2 += AliMUONConstants::DefaultChamberZ(chamber) * AliMUONConstants::DefaultChamberZ(chamber) * fMaxMCSAngle2[iCh]; + } + + } // ------ track selection in non bending direction ------ if (GetRecoParam()->SelectOnTrackSlope()) { // check if non bending slope is within tolerances - Double_t nonBendingSlopeErr = TMath::Sqrt(kParamCov(1,1) + (chamber + 1.) * angleMCS2); + Double_t nonBendingSlopeErr = TMath::Sqrt(kParamCov(1,1) + angleMCS2); if ((TMath::Abs(trackParam.GetNonBendingSlope()) - sigmaCut * nonBendingSlopeErr) > GetRecoParam()->GetMaxNonBendingSlope()) return kFALSE; } else { @@ -233,7 +271,7 @@ Bool_t AliMUONVTrackReconstructor::IsAcceptable(AliMUONTrackParam &trackParam) if (GetRecoParam()->SelectOnTrackSlope()) { // check if bending slope is within tolerances - Double_t bendingSlopeErr = TMath::Sqrt(kParamCov(3,3) + (chamber + 1.) * angleMCS2); + Double_t bendingSlopeErr = TMath::Sqrt(kParamCov(3,3) + angleMCS2); if ((TMath::Abs(trackParam.GetBendingSlope()) - sigmaCut * bendingSlopeErr) > GetRecoParam()->GetMaxBendingSlope()) return kFALSE; } else { @@ -255,10 +293,9 @@ Bool_t AliMUONVTrackReconstructor::IsAcceptable(AliMUONTrackParam &trackParam) TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsBetweenChambers(const AliMUONVClusterStore& clusterStore, Int_t ch1, Int_t ch2) { /// To make the list of segments from the list of clusters in the 2 given chambers. - /// Return a new TClonesArray of segments. - /// It is the responsibility of the user to delete it afterward. + /// Return a TClonesArray of new segments (segments made in a previous call of this function are removed). AliDebug(1,Form("Enter MakeSegmentsBetweenChambers (1..) %d-%d", ch1+1, ch2+1)); - AliCodeTimerAuto(""); + AliCodeTimerAuto("",0); AliMUONVCluster *cluster1, *cluster2; AliMUONObjectPair *segment; @@ -267,9 +304,12 @@ TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsBetweenChambers(const AliM Double_t bendingSlope = 0., bendingSlopeErr = 0., bendingImpactParam = 0., bendingImpactParamErr = 0., bendingImpactParamErr2 = 0.; Double_t bendingMomentum = 0., bendingMomentumErr = 0.; Double_t bendingVertexDispersion2 = GetRecoParam()->GetBendingVertexDispersion() * GetRecoParam()->GetBendingVertexDispersion(); - Double_t impactMCS2 = 0; // maximum impact parameter dispersion**2 due to MCS in chamber - if (!GetRecoParam()->SelectOnTrackSlope() || AliMUONTrackExtrap::IsFieldON()) for (Int_t iCh=0; iCh<=ch1; iCh++) - impactMCS2 += AliMUONConstants::DefaultChamberZ(iCh) * AliMUONConstants::DefaultChamberZ(iCh) * fMaxMCSAngle2; + Double_t angleMCS2 = 0.; // maximum angular dispersion**2 due to MCS in chamber + Double_t impactMCS2 = 0.; // maximum impact parameter dispersion**2 due to MCS in chamber + for (Int_t iCh=0; iCh<=ch1; iCh++) { + angleMCS2 += fMaxMCSAngle2[iCh]; + impactMCS2 += AliMUONConstants::DefaultChamberZ(iCh) * AliMUONConstants::DefaultChamberZ(iCh) * fMaxMCSAngle2[iCh]; + } Double_t sigmaCut = GetRecoParam()->GetSigmaCutForTracking(); // Create iterators to loop over clusters in both chambers @@ -277,7 +317,8 @@ TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsBetweenChambers(const AliM TIter nextInCh2(clusterStore.CreateChamberIterator(ch2,ch2)); // list of segments - TClonesArray *segments = new TClonesArray("AliMUONObjectPair", 100); + static TClonesArray *segments = new TClonesArray("AliMUONObjectPair", 100); + segments->Clear("C"); // Loop over clusters in the first chamber of the station while ( ( cluster1 = static_cast(nextInCh1()) ) ) { @@ -296,7 +337,7 @@ TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsBetweenChambers(const AliM if (GetRecoParam()->SelectOnTrackSlope()) { // check if non bending slope is within tolerances - nonBendingSlopeErr = TMath::Sqrt((cluster1->GetErrX2() + cluster2->GetErrX2()) / dZ / dZ + (ch1 + 1.) * fMaxMCSAngle2); + nonBendingSlopeErr = TMath::Sqrt((cluster1->GetErrX2() + cluster2->GetErrX2()) / dZ / dZ + angleMCS2); if ((TMath::Abs(nonBendingSlope) - sigmaCut * nonBendingSlopeErr) > GetRecoParam()->GetMaxNonBendingSlope()) continue; } else { @@ -325,7 +366,7 @@ TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsBetweenChambers(const AliM if (GetRecoParam()->SelectOnTrackSlope()) { // check if bending slope is within tolerances - bendingSlopeErr = TMath::Sqrt((cluster1->GetErrY2() + cluster2->GetErrY2()) / dZ / dZ + (ch1 + 1.) * fMaxMCSAngle2); + bendingSlopeErr = TMath::Sqrt((cluster1->GetErrY2() + cluster2->GetErrY2()) / dZ / dZ + angleMCS2); if ((TMath::Abs(bendingSlope) - sigmaCut * bendingSlopeErr) > GetRecoParam()->GetMaxBendingSlope()) continue; } else { @@ -509,54 +550,109 @@ void AliMUONVTrackReconstructor::RemoveDoubleTracks() } //__________________________________________________________________________ -void AliMUONVTrackReconstructor::RemoveConnectedTracks(Bool_t inSt345) +void AliMUONVTrackReconstructor::RemoveBadTracks() { - /// To remove double tracks: - /// Tracks are considered identical if they share 1 cluster or more. - /// If inSt345=kTRUE only stations 3, 4 and 5 are considered. - /// Among two identical tracks, one keeps the track with the larger number of clusters - /// or, if these numbers are equal, the track with the minimum chi2. - AliMUONTrack *track1, *track2, *trackToRemove; - Int_t clustersInCommon, nClusters1, nClusters2; - Bool_t removedTrack1; + /// Remove tracks for which a problem occured somewhere during the tracking + + AliMUONTrack *track, *nextTrack; + Bool_t trackRemoved = kFALSE; + + track = (AliMUONTrack*) fRecTracksPtr->First(); + while (track) { + + nextTrack = (AliMUONTrack*) fRecTracksPtr->After(track); + + if (track->GetGlobalChi2() >= AliMUONTrack::MaxChi2()) { + AliWarning("problem occured somewhere during the tracking --> discard track"); + fRecTracksPtr->Remove(track); + fNRecTracks--; + trackRemoved = kTRUE; + } + + track = nextTrack; + + } + + // compress array of tracks if needed + if (trackRemoved) fRecTracksPtr->Compress(); + +} + + //__________________________________________________________________________ +void AliMUONVTrackReconstructor::RemoveConnectedTracks(Int_t stMin, Int_t stMax, Bool_t all) +{ + /// Find and remove tracks sharing 1 cluster or more in station(s) [stMin, stMax]. + /// For each couple of connected tracks, one removes the one with the smallest + /// number of clusters or with the highest chi2 value in case of equality. + /// If all=kTRUE: both tracks are removed. + + // tag the tracks to be removed + TagConnectedTracks(stMin, stMax, all); + + // remove them + Int_t nTracks = fRecTracksPtr->GetEntriesFast(); + for (Int_t i = 0; i < nTracks; i++) { + if (((AliMUONTrack*) fRecTracksPtr->UncheckedAt(i))->IsConnected()) { + fRecTracksPtr->RemoveAt(i); + fNRecTracks--; + } + } + + // remove holes in the array if any + fRecTracksPtr->Compress(); +} + + //__________________________________________________________________________ +void AliMUONVTrackReconstructor::TagConnectedTracks(Int_t stMin, Int_t stMax, Bool_t all) +{ + /// Find and tag tracks sharing 1 cluster or more in station(s) [stMin, stMax]. + /// For each couple of connected tracks, one tags the one with the smallest + /// number of clusters or with the highest chi2 value in case of equality. + /// If all=kTRUE: both tracks are tagged. + + AliMUONTrack *track1, *track2; + Int_t nClusters1, nClusters2; + Int_t nTracks = fRecTracksPtr->GetEntriesFast(); + + // reset the tags + for (Int_t i = 0; i < nTracks; i++) ((AliMUONTrack*) fRecTracksPtr->UncheckedAt(i))->Connected(kFALSE); + // Loop over first track of the pair - track1 = (AliMUONTrack*) fRecTracksPtr->First(); - while (track1) { - removedTrack1 = kFALSE; - nClusters1 = track1->GetNClusters(); + for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) { + track1 = (AliMUONTrack*) fRecTracksPtr->UncheckedAt(iTrack1); + // Loop over second track of the pair - track2 = (AliMUONTrack*) fRecTracksPtr->After(track1); - while (track2) { - nClusters2 = track2->GetNClusters(); - // number of clusters in common between two tracks - if (inSt345) clustersInCommon = track1->ClustersInCommonInSt345(track2); - else clustersInCommon = track1->ClustersInCommon(track2); - // check for identical tracks - if (clustersInCommon > 0) { - // decide which track to remove - if ((nClusters1 > nClusters2) || ((nClusters1 == nClusters2) && (track1->GetGlobalChi2() <= track2->GetGlobalChi2()))) { - // remove track2 and continue the second loop with the track next to track2 - trackToRemove = track2; - track2 = (AliMUONTrack*) fRecTracksPtr->After(track2); - fRecTracksPtr->Remove(trackToRemove); - fRecTracksPtr->Compress(); // this is essential to retrieve the TClonesArray afterwards - fNRecTracks--; - } else { - // else remove track1 and continue the first loop with the track next to track1 - trackToRemove = track1; - track1 = (AliMUONTrack*) fRecTracksPtr->After(track1); - fRecTracksPtr->Remove(trackToRemove); - fRecTracksPtr->Compress(); // this is essential to retrieve the TClonesArray afterwards - fNRecTracks--; - removedTrack1 = kTRUE; - break; - } - } else track2 = (AliMUONTrack*) fRecTracksPtr->After(track2); - } // track2 - if (removedTrack1) continue; - track1 = (AliMUONTrack*) fRecTracksPtr->After(track1); - } // track1 - return; + for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) { + track2 = (AliMUONTrack*) fRecTracksPtr->UncheckedAt(iTrack2); + + // check for connected tracks and tag them + if (track1->ClustersInCommon(track2, stMin, stMax) > 0) { + + if (all) { + + // tag both tracks + track1->Connected(); + track2->Connected(); + + } else { + + // tag only the worst track + nClusters1 = track1->GetNClusters(); + nClusters2 = track2->GetNClusters(); + if ((nClusters1 > nClusters2) || ((nClusters1 == nClusters2) && (track1->GetGlobalChi2() <= track2->GetGlobalChi2()))) { + track2->Connected(); + } else { + track1->Connected(); + } + + } + + } + + } + + } + } //__________________________________________________________________________ @@ -574,7 +670,7 @@ void AliMUONVTrackReconstructor::AskForNewClustersInChamber(const AliMUONTrackPa // extrapolate track parameters to the chamber AliMUONTrackParam extrapTrackParam(trackParam); - AliMUONTrackExtrap::ExtrapToZCov(&extrapTrackParam, AliMUONConstants::DefaultChamberZ(chamber)); + if (!AliMUONTrackExtrap::ExtrapToZCov(&extrapTrackParam, AliMUONConstants::DefaultChamberZ(chamber))) return; // build the searching area using the track and chamber resolutions and the maximum-distance-to-track value const TMatrixD& kParamCov = extrapTrackParam.GetCovariances(); @@ -616,11 +712,11 @@ Double_t AliMUONVTrackReconstructor::TryOneCluster(const AliMUONTrackParam &trac /// return trackParamAtCluster // extrapolate track parameters and covariances at the z position of the tested cluster + // and set pointer to cluster into trackParamAtCluster trackParamAtCluster = trackParam; - AliMUONTrackExtrap::ExtrapToZCov(&trackParamAtCluster, cluster->GetZ(), updatePropagator); - - // set pointer to cluster into trackParamAtCluster trackParamAtCluster.SetClusterPtr(cluster); + if (!AliMUONTrackExtrap::ExtrapToZCov(&trackParamAtCluster, cluster->GetZ(), updatePropagator)) + return 2.*AliMUONTrack::MaxChi2(); // Set differences between trackParam and cluster in the bending and non bending directions Double_t dX = cluster->GetX() - trackParamAtCluster.GetNonBendingCoor(); @@ -634,7 +730,7 @@ Double_t AliMUONVTrackReconstructor::TryOneCluster(const AliMUONTrackParam &trac Double_t det = sigmaX2 * sigmaY2 - covXY * covXY; // Compute chi2 - if (det == 0.) return 1.e10; + if (det == 0.) return 2.*AliMUONTrack::MaxChi2(); return (dX * dX * sigmaY2 + dY * dY * sigmaX2 - 2. * dX * dY * covXY) / det; } @@ -702,7 +798,7 @@ Double_t AliMUONVTrackReconstructor::TryTwoClustersFast(const AliMUONTrackParam // Compute chi2 Double_t detX = sigma2X1 * sigma2X2 - covX1X2 * covX1X2; Double_t detY = sigma2Y1 * sigma2Y2 - covY1Y2 * covY1Y2; - if (detX == 0. || detY == 0.) return 1.e10; + if (detX == 0. || detY == 0.) return 2.*AliMUONTrack::MaxChi2(); return (dX1 * dX1 * sigma2X2 + dX2 * dX2 * sigma2X1 - 2. * dX1 * dX2 * covX1X2) / detX + (dY1 * dY1 * sigma2Y2 + dY2 * dY2 * sigma2Y1 - 2. * dY1 * dY2 * covY1Y2) / detY; @@ -721,7 +817,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInChamber(AliMUONTrack &trac /// return kTRUE if new cluster(s) have been found (otherwise return kFALSE) AliDebug(1,Form("Enter FollowLinearTrackInChamber(1..) %d", nextChamber+1)); - Double_t chi2WithOneCluster = 1.e10; + Double_t chi2WithOneCluster = AliMUONTrack::MaxChi2(); Double_t maxChi2WithOneCluster = 2. * GetRecoParam()->GetSigmaCutForTracking() * GetRecoParam()->GetSigmaCutForTracking(); // 2 because 2 quantities in chi2 Double_t bestChi2WithOneCluster = maxChi2WithOneCluster; @@ -744,7 +840,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInChamber(AliMUONTrack &trac } // Add MCS effect - AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(),1.); + AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(trackParam.GetClusterPtr()->GetChamberId()),-1.); // Printout for debuging if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) { @@ -856,8 +952,8 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac ch2 = 2*nextStation+1; } - Double_t chi2WithOneCluster = 1.e10; - Double_t chi2WithTwoClusters = 1.e10; + Double_t chi2WithOneCluster = AliMUONTrack::MaxChi2(); + Double_t chi2WithTwoClusters = AliMUONTrack::MaxChi2(); Double_t maxChi2WithOneCluster = 2. * GetRecoParam()->GetSigmaCutForTracking() * GetRecoParam()->GetSigmaCutForTracking(); // 2 because 2 quantities in chi2 Double_t maxChi2WithTwoClusters = 4. * GetRecoParam()->GetSigmaCutForTracking() * @@ -891,7 +987,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac } // Add MCS effect - AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(),1.); + AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(trackParam.GetClusterPtr()->GetChamberId()),-1.); // Printout for debuging if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) { @@ -926,7 +1022,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac } // add MCS effect - AliMUONTrackExtrap::AddMCSEffect(&extrapTrackParamAtCluster2,AliMUONConstants::ChamberThicknessInX0(),1.); + AliMUONTrackExtrap::AddMCSEffect(&extrapTrackParamAtCluster2,AliMUONConstants::ChamberThicknessInX0(ch2),-1.); // reset cluster iterator of chamber 1 nextInCh1.Reset(); @@ -1029,7 +1125,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac AliMUONTrackExtrap::LinearExtrapToZCov(&trackParam, AliMUONConstants::DefaultChamberZ(ch2)); // add MCS effect for next step - AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(),1.); + AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(ch2),-1.); // reset cluster iterator of chamber 1 nextInCh1.Reset(); @@ -1172,19 +1268,134 @@ void AliMUONVTrackReconstructor::ImproveTracks() //__________________________________________________________________________ void AliMUONVTrackReconstructor::Finalize() { - /// Recompute track parameters and covariances at each attached cluster from those at the first one + /// Recompute track parameters and covariances at each attached cluster /// Set the label pointing to the corresponding MC track + /// Remove the track if finalization failed - AliMUONTrack *track; + AliMUONTrack *track, *nextTrack; + Bool_t trackRemoved = kFALSE; track = (AliMUONTrack*) fRecTracksPtr->First(); while (track) { - FinalizeTrack(*track); + nextTrack = (AliMUONTrack*) fRecTracksPtr->After(track); + + if (FinalizeTrack(*track)) track->FindMCLabel(); + else { + fRecTracksPtr->Remove(track); + fNRecTracks--; + trackRemoved = kTRUE; + } + + track = nextTrack; + + } + + // compress array of tracks if needed + if (trackRemoved) fRecTracksPtr->Compress(); + +} + +//__________________________________________________________________________ +void AliMUONVTrackReconstructor::DiscardMonoCathodClusters() +{ + /// Assign a different resolution to the mono-cathod clusters + /// in the direction of the missing plane and refit the track + /// Remove the track in case of failure + + if (!fkTransformer) AliFatal("missing geometry transformer"); + + AliMUONTrack *track, *nextTrack; + Bool_t trackRemoved = kFALSE; + + track = (AliMUONTrack*) fRecTracksPtr->First(); + while (track) { + + nextTrack = (AliMUONTrack*) fRecTracksPtr->After(track); + + ChangeMonoCathodClusterRes(*track); + + if (!RefitTrack(*track) || (GetRecoParam()->ImproveTracks() && !track->IsImproved())) { + fRecTracksPtr->Remove(track); + fNRecTracks--; + trackRemoved = kTRUE; + } + + track = nextTrack; + + } + + // compress array of tracks if needed + if (trackRemoved) fRecTracksPtr->Compress(); + +} + +//__________________________________________________________________________ +void AliMUONVTrackReconstructor::ChangeMonoCathodClusterRes(AliMUONTrack &track) +{ + /// Assign a different resolution to the mono-cathod clusters + /// in the direction of the missing plane and refit the track + + // Loop over clusters + AliMUONVCluster *cluster; + Int_t nClusters = track.GetNClusters(); + for (Int_t iCluster = 0; iCluster < nClusters; iCluster++) { + cluster = ((AliMUONTrackParam*) track.GetTrackParamAtCluster()->UncheckedAt(iCluster))->GetClusterPtr(); + + // do it only for stations 3, 4 & 5 + if (cluster->GetChamberId() < 4) continue; - track->FindMCLabel(); + // get the cathod corresponding to the bending/non-bending plane + Int_t deId = cluster->GetDetElemId(); + AliMpDetElement* de = AliMpDDLStore::Instance()->GetDetElement(deId, kFALSE); + if (!de) continue; + AliMp::CathodType cath1 = de->GetCathodType(AliMp::kBendingPlane); + AliMp::CathodType cath2 = de->GetCathodType(AliMp::kNonBendingPlane); - track = (AliMUONTrack*) fRecTracksPtr->After(track); + // get the corresponding segmentation + const AliMpVSegmentation* seg1 = AliMpSegmentation::Instance()->GetMpSegmentation(deId, cath1); + const AliMpVSegmentation* seg2 = AliMpSegmentation::Instance()->GetMpSegmentation(deId, cath2); + if (!seg1 || !seg2) continue; + + // get local coordinate of the cluster + Double_t lX,lY,lZ; + Double_t gX = cluster->GetX(); + Double_t gY = cluster->GetY(); + Double_t gZ = cluster->GetZ(); + fkTransformer->Global2Local(deId,gX,gY,gZ,lX,lY,lZ); + + // find pads below the cluster + AliMpPad pad1 = seg1->PadByPosition(lX, lY, kFALSE); + AliMpPad pad2 = seg2->PadByPosition(lX, lY, kFALSE); + + // build their ID if pads are valid + UInt_t padId1 = (pad1.IsValid()) ? AliMUONVDigit::BuildUniqueID(deId, pad1.GetManuId(), pad1.GetManuChannel(), cath1) : 0; + UInt_t padId2 = (pad2.IsValid()) ? AliMUONVDigit::BuildUniqueID(deId, pad2.GetManuId(), pad2.GetManuChannel(), cath2) : 0; + + // check if the cluster contains these pads + Bool_t hasNonBending = kFALSE; + Bool_t hasBending = kFALSE; + for (Int_t iDigit = 0; iDigit < cluster->GetNDigits(); iDigit++) { + + UInt_t digitId = cluster->GetDigitId(iDigit); + + if (digitId == padId1) { + + hasBending = kTRUE; + if (hasNonBending) break; + + } else if (digitId == padId2) { + + hasNonBending = kTRUE; + if (hasBending) break; + + } + + } + + // modify the cluster resolution if needed + if (!hasNonBending) cluster->SetErrXY(GetRecoParam()->GetMonoCathodClNonBendingRes(), cluster->GetErrY()); + if (!hasBending) cluster->SetErrXY(cluster->GetErrX(), GetRecoParam()->GetMonoCathodClBendingRes()); } @@ -1197,20 +1408,21 @@ void AliMUONVTrackReconstructor::ValidateTracksWithTrigger(AliMUONVTrackStore& t const AliMUONTrackHitPattern& trackHitPattern) { /// Try to match track from tracking system with trigger track - AliCodeTimerAuto(""); + AliCodeTimerAuto("",0); trackHitPattern.ExecuteValidation(trackStore, triggerTrackStore, triggerStore); } - //__________________________________________________________________________ + +//__________________________________________________________________________ void AliMUONVTrackReconstructor::EventReconstructTrigger(const AliMUONTriggerCircuit& circuit, const AliMUONVTriggerStore& triggerStore, AliMUONVTriggerTrackStore& triggerTrackStore) { - /// To make the trigger tracks from Local Trigger + /// Fill trigger track store from local trigger AliDebug(1, ""); - AliCodeTimerAuto(""); - + AliCodeTimerAuto("",0); + AliMUONGlobalTrigger* globalTrigger = triggerStore.Global(); UChar_t gloTrigPat = 0; @@ -1220,60 +1432,80 @@ void AliMUONVTrackReconstructor::EventReconstructTrigger(const AliMUONTriggerCir gloTrigPat = globalTrigger->GetGlobalResponse(); } + AliMUONTriggerTrack triggerTrack; + TIter next(triggerStore.CreateIterator()); AliMUONLocalTrigger* locTrg(0x0); - - Float_t z11 = AliMUONConstants::DefaultChamberZ(10); - Float_t z21 = AliMUONConstants::DefaultChamberZ(12); - - AliMUONTriggerTrack triggerTrack; while ( ( locTrg = static_cast(next()) ) ) { - Bool_t xTrig=locTrg->IsTrigX(); - Bool_t yTrig=locTrg->IsTrigY(); - - Int_t localBoardId = locTrg->LoCircuit(); - - if (xTrig && yTrig) + if ( locTrg->IsTrigX() && locTrg->IsTrigY() ) { // make Trigger Track if trigger in X and Y - Float_t y11 = circuit.GetY11Pos(localBoardId, locTrg->LoStripX()); - // need first to convert deviation to [0-30] - // (see AliMUONLocalTriggerBoard::LocalTrigger) - Int_t deviation = locTrg->GetDeviation(); - Int_t stripX21 = locTrg->LoStripX()+deviation+1; - Float_t y21 = circuit.GetY21Pos(localBoardId, stripX21); - Float_t x11 = circuit.GetX11Pos(localBoardId, locTrg->LoStripY()); - - AliDebug(1, Form(" MakeTriggerTrack %d %d %d %d %f %f %f \n",locTrg->LoCircuit(), - locTrg->LoStripX(),locTrg->LoStripX()+locTrg->LoDev()+1,locTrg->LoStripY(),y11, y21, x11)); - - Float_t thetax = TMath::ATan2( x11 , z11 ); - Float_t thetay = TMath::ATan2( (y21-y11) , (z21-z11) ); + if (TriggerToTrack(circuit, *locTrg, triggerTrack, gloTrigPat)) + triggerTrackStore.Add(triggerTrack); - CorrectThetaRange(thetax); - CorrectThetaRange(thetay); - - triggerTrack.SetX11(x11); - triggerTrack.SetY11(y11); - triggerTrack.SetThetax(thetax); - triggerTrack.SetThetay(thetay); - triggerTrack.SetGTPattern(gloTrigPat); - triggerTrack.SetLoTrgNum(localBoardId); - - triggerTrackStore.Add(triggerTrack); } // board is fired } // end of loop on Local Trigger } //__________________________________________________________________________ -void AliMUONVTrackReconstructor::CorrectThetaRange(Float_t& theta) +Bool_t AliMUONVTrackReconstructor::TriggerToTrack(const AliMUONTriggerCircuit& circuit, + const AliMUONLocalTrigger& locTrg, + AliMUONTriggerTrack& triggerTrack, + UChar_t globalTriggerPattern) { - /// The angles of the trigger tracks, obtained with ATan2, - /// have values around +pi and -pi. On the contrary, the angles - /// used in the tracker tracks have values around 0. - /// This function sets the same range for the trigger tracks angles. - if (theta < -TMath::PiOver2()) theta += TMath::Pi(); - else if(theta > TMath::PiOver2()) theta -= TMath::Pi(); + /// To make the trigger tracks from Local Trigger + const Double_t kTrigNonBendReso = AliMUONConstants::TriggerNonBendingReso(); + const Double_t kTrigBendReso = AliMUONConstants::TriggerBendingReso(); + const Double_t kSqrt12 = TMath::Sqrt(12.); + + TMatrixD trigCov(3,3); + + Int_t localBoardId = locTrg.LoCircuit(); + + Float_t y11 = circuit.GetY11Pos(localBoardId, locTrg.LoStripX()); + Float_t z11 = circuit.GetZ11Pos(localBoardId, locTrg.LoStripX()); + // need first to convert deviation to [0-30] + // (see AliMUONLocalTriggerBoard::LocalTrigger) + Int_t deviation = locTrg.GetDeviation(); + Int_t stripX21 = locTrg.LoStripX()+deviation+1; + Float_t y21 = circuit.GetY21Pos(localBoardId, stripX21); + Float_t z21 = circuit.GetZ21Pos(localBoardId, stripX21); + Float_t x11 = circuit.GetX11Pos(localBoardId, locTrg.LoStripY()); + + AliDebug(1, Form(" MakeTriggerTrack %3d %2d %2d %2d (%f %f %f) (%f %f)\n",locTrg.LoCircuit(), + locTrg.LoStripX(),locTrg.LoStripX()+deviation+1,locTrg.LoStripY(),x11, y11, z11, y21, z21)); + + if (TMath::Abs(z11) < 0.00001) return kFALSE; + + Double_t deltaZ = z11 - z21; + + Float_t slopeX = x11/z11; + Float_t slopeY = (y11-y21) / deltaZ; + + Float_t sigmaX = circuit.GetX11Width(localBoardId, locTrg.LoStripY()) / kSqrt12; + Float_t sigmaY = circuit.GetY11Width(localBoardId, locTrg.LoStripX()) / kSqrt12; + Float_t sigmaY21 = circuit.GetY21Width(localBoardId, locTrg.LoStripX()) / kSqrt12; + + trigCov.Zero(); + trigCov(0,0) = kTrigNonBendReso * kTrigNonBendReso + sigmaX * sigmaX; + trigCov(1,1) = kTrigBendReso * kTrigBendReso + sigmaY * sigmaY; + trigCov(2,2) = + (2. * kTrigBendReso * kTrigBendReso + sigmaY * sigmaY + sigmaY21 * sigmaY21 ) / deltaZ / deltaZ; + trigCov(1,2) = trigCov(2,1) = trigCov(1,1) / deltaZ; + + triggerTrack.SetX11(x11); + triggerTrack.SetY11(y11); + triggerTrack.SetZ11(z11); + triggerTrack.SetZ21(z21); + triggerTrack.SetSlopeX(slopeX); + triggerTrack.SetSlopeY(slopeY); + triggerTrack.SetGTPattern(globalTriggerPattern); + triggerTrack.SetLoTrgNum(localBoardId); + triggerTrack.SetCovariances(trigCov); + triggerTrack.SetUniqueID(locTrg.GetUniqueID()); + + return kTRUE; + }