]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/AliMUONVTrackReconstructor.cxx
Default recoparam object
[u/mrichter/AliRoot.git] / MUON / AliMUONVTrackReconstructor.cxx
index 7c44f936b6975d403024135dafd0f07745befd2f..074c42da54911e6130c905880cac07209c6f44f0 100644 (file)
 /// * ValidateTracksWithTrigger to match tracker/trigger tracks
 ///
 /// Several options and adjustable parameters are available for both KALMAN and ORIGINAL
-/// tracking algorithms. They can be changed by using:
-/// AliMUONRecoParam *muonRecoParam = AliMUONRecoParam::GetLow(High)FluxParam();
-/// muonRecoParam->Set...(); // see methods in AliMUONRecoParam.h for details
-/// AliMUONReconstructor::SetRecoParam(muonRecoParam);
+/// tracking algorithms. They can be changed through the AliMUONRecoParam object
+/// set in the reconstruction macro or read from the CDB
+/// (see methods in AliMUONRecoParam.h file for details)
 ///
 /// Main parameters and options are:
 /// - *fgkSigmaToCutForTracking* : quality cut used to select new clusters to be
@@ -93,13 +92,16 @@ ClassImp(AliMUONVTrackReconstructor) // Class implementation in ROOT context
 /// \endcond
 
   //__________________________________________________________________________
-AliMUONVTrackReconstructor::AliMUONVTrackReconstructor(AliMUONVClusterServer& clusterServer)
-  : TObject(),
-    fRecTracksPtr(0x0),
-    fNRecTracks(0),
-  fClusterServer(clusterServer)
+AliMUONVTrackReconstructor::AliMUONVTrackReconstructor(const AliMUONRecoParam* recoParam,
+                                                       AliMUONVClusterServer* clusterServer)
+: TObject(),
+fRecTracksPtr(0x0),
+fNRecTracks(0),
+fClusterServer(clusterServer),
+fRecoParam(recoParam)
 {
   /// Constructor for class AliMUONVTrackReconstructor
+  /// WARNING: if clusterServer=0x0, no clusterization will be possible at this level
   
   // Memory allocation for the TClonesArray of reconstructed tracks
   fRecTracksPtr = new TClonesArray("AliMUONTrack", 100);
@@ -136,9 +138,12 @@ void AliMUONVTrackReconstructor::EventReconstruct(AliMUONVClusterStore& clusterS
   // Reset array of tracks
   ResetTracks();
   
-  // Look for candidates from at least 3 aligned points in stations(1..) 4 and 5
+  // Look for candidates from clusters in stations(1..) 4 and 5
   MakeTrackCandidates(clusterStore);
   
+  // Look for extra candidates from clusters in stations(1..) 4 and 5
+  if (GetRecoParam()->MakeMoreTrackCandidates()) MakeMoreTrackCandidates(clusterStore);
+  
   // Stop tracking if no candidate found
   if (fRecTracksPtr->GetEntriesFast() == 0) return;
   
@@ -146,10 +151,10 @@ void AliMUONVTrackReconstructor::EventReconstruct(AliMUONVClusterStore& clusterS
   FollowTracks(clusterStore);
   
   // Complement the reconstructed tracks
-  if (AliMUONReconstructor::GetRecoParam()->ComplementTracks()) ComplementTracks(clusterStore);
+  if (GetRecoParam()->ComplementTracks()) ComplementTracks(clusterStore);
   
   // Improve the reconstructed tracks
-  if (AliMUONReconstructor::GetRecoParam()->ImproveTracks()) ImproveTracks();
+  if (GetRecoParam()->ImproveTracks()) ImproveTracks();
   
   // Remove double tracks
   RemoveDoubleTracks();
@@ -161,23 +166,22 @@ void AliMUONVTrackReconstructor::EventReconstruct(AliMUONVClusterStore& clusterS
   for (Int_t i=0; i<fNRecTracks; ++i) 
   {
     AliMUONTrack * track = (AliMUONTrack*) fRecTracksPtr->At(i);
+    track->SetUniqueID(i+1);
     trackStore.Add(*track);
   }
 }
 
   //__________________________________________________________________________
-TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsInStation(const AliMUONVClusterStore& clusterStore, Int_t station)
+TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsBetweenChambers(const AliMUONVClusterStore& clusterStore, Int_t ch1, Int_t ch2)
 {
-  /// To make the list of segments in station(0..) "Station" from the list of clusters to be reconstructed.
+  /// 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.
-  AliDebug(1,Form("Enter MakeSegmentsPerStation (1..) %d",station+1));
+  AliDebug(1,Form("Enter MakeSegmentsBetweenChambers (1..) %d-%d", ch1+1, ch2+1));
   
   AliMUONVCluster *cluster1, *cluster2;
   AliMUONObjectPair *segment;
-  Double_t bendingSlope = 0, impactParam = 0., bendingMomentum = 0.; // to avoid compilation warning
-  Int_t ch1 = 2 * station;
-  Int_t ch2 = ch1 + 1;
+  Double_t nonBendingSlope = 0, bendingSlope = 0, impactParam = 0., bendingMomentum = 0.; // to avoid compilation warning
   
   // Create iterators to loop over clusters in both chambers
   TIter nextInCh1(clusterStore.CreateChamberIterator(ch1,ch1));
@@ -185,7 +189,6 @@ TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsInStation(const AliMUONVCl
   
   // list of segments
   TClonesArray *segments = new TClonesArray("AliMUONObjectPair", 100);
-  segments->SetOwner(kTRUE);
   
   // Loop over clusters in the first chamber of the station
   while ( ( cluster1 = static_cast<AliMUONVCluster*>(nextInCh1()) ) ) {
@@ -196,32 +199,45 @@ TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsInStation(const AliMUONVCl
     // Loop over clusters in the second chamber of the station
     while ( ( cluster2 = static_cast<AliMUONVCluster*>(nextInCh2()) ) ) {
       
+      // non bending slope
+      nonBendingSlope = (cluster1->GetX() - cluster2->GetX()) / (cluster1->GetZ() - cluster2->GetZ());
+      
+      // check if non bending slope is within tolerances
+      if (TMath::Abs(nonBendingSlope) > GetRecoParam()->GetMaxNonBendingSlope()) continue;
+      
       // bending slope
       bendingSlope = (cluster1->GetY() - cluster2->GetY()) / (cluster1->GetZ() - cluster2->GetZ());
       
-      // impact parameter
-      impactParam = cluster1->GetY() - cluster1->GetZ() * bendingSlope;
-     
-      // absolute value of bending momentum
-      bendingMomentum = TMath::Abs(AliMUONTrackExtrap::GetBendingMomentumFromImpactParam(impactParam));
-      
-      // check for bending momentum within tolerances
-      if ((bendingMomentum < AliMUONReconstructor::GetRecoParam()->GetMaxBendingMomentum()) &&
-         (bendingMomentum > AliMUONReconstructor::GetRecoParam()->GetMinBendingMomentum())) {
-        
-       // make new segment
-        segment = new ((*segments)[segments->GetLast()+1]) AliMUONObjectPair(cluster1, cluster2, kFALSE, kFALSE);
-        
-       // Printout for debug
-       if (AliLog::GetGlobalDebugLevel() > 1) {
-          cout << "segmentIndex(0...): " << segments->GetLast() << endl;
-          segment->Dump();
-          cout << "Cluster in first chamber" << endl;
-          cluster1->Print();
-          cout << "Cluster in second chamber" << endl;
-          cluster2->Print();
-        }
+      // check the bending momentum of the bending slope depending if the field is ON or OFF
+      if (AliMUONTrackExtrap::IsFieldON()) {
+       
+       // impact parameter
+       impactParam = cluster1->GetY() - cluster1->GetZ() * bendingSlope;
+       
+       // absolute value of bending momentum
+       bendingMomentum = TMath::Abs(AliMUONTrackExtrap::GetBendingMomentumFromImpactParam(impactParam));
+       
+       // check if bending momentum is within tolerances
+       if (bendingMomentum < GetRecoParam()->GetMinBendingMomentum() ||
+           bendingMomentum > GetRecoParam()->GetMaxBendingMomentum()) continue;
        
+      } else {
+       
+       // check if non bending slope is within tolerances
+       if (TMath::Abs(bendingSlope) > GetRecoParam()->GetMaxBendingSlope()) continue;
+      
+      }
+      // make new segment
+      segment = new ((*segments)[segments->GetLast()+1]) AliMUONObjectPair(cluster1, cluster2, kFALSE, kFALSE);
+      
+      // Printout for debug
+      if (AliLog::GetGlobalDebugLevel() > 1) {
+       cout << "segmentIndex(0...): " << segments->GetLast() << endl;
+       segment->Dump();
+       cout << "Cluster in first chamber" << endl;
+       cluster1->Print();
+       cout << "Cluster in second chamber" << endl;
+       cluster2->Print();
       }
       
     }
@@ -229,7 +245,7 @@ TClonesArray* AliMUONVTrackReconstructor::MakeSegmentsInStation(const AliMUONVCl
   }
   
   // Printout for debug
-  AliDebug(1,Form("Station: %d  NSegments:  %d ", station+1, segments->GetEntriesFast()));
+  AliDebug(1,Form("chambers%d-%d: NSegments =  %d ", ch1+1, ch2+1, segments->GetEntriesFast()));
   
   return segments;
 }
@@ -333,39 +349,49 @@ void AliMUONVTrackReconstructor::RemoveDoubleTracks()
 }
 
   //__________________________________________________________________________
-void AliMUONVTrackReconstructor::AskForNewClustersInStation(const AliMUONTrackParam &trackParam,
-                                                           AliMUONVClusterStore& clusterStore, Int_t station)
+void AliMUONVTrackReconstructor::AskForNewClustersInChamber(const AliMUONTrackParam &trackParam,
+                                                           AliMUONVClusterStore& clusterStore, Int_t chamber)
 {
   /// Ask the clustering to reconstruct new clusters around the track candidate position
-  /// in the 2 chambers of the given station
   
-  // maximum shift of the searching area due to distance between detection elements and the track slope
-  static const Double_t kgMaxShift = 2.; // 2 cm
+  // check if the current chamber is useable
+  if (!fClusterServer || !GetRecoParam()->UseChamber(chamber)) return;
+  
+  // maximum distance between the center of the chamber and a detection element
+  // (accounting for the inclination of the chamber)
+  static const Double_t kMaxDZ = 15.; // 15 cm
   
-  // extrapolate track parameters to the second chamber of the station
+  // extrapolate track parameters to the chamber
   AliMUONTrackParam extrapTrackParam(trackParam);
-  AliMUONTrackExtrap::ExtrapToZ(&extrapTrackParam, AliMUONConstants::DefaultChamberZ(2*station+1));
+  AliMUONTrackExtrap::ExtrapToZCov(&extrapTrackParam, AliMUONConstants::DefaultChamberZ(chamber));
   
-  // build the searching area
+  // build the searching area using the track resolution and the maximum-distance-to-track value
+  const TMatrixD& kParamCov = extrapTrackParam.GetCovariances();
+  Double_t errX2 = kParamCov(0,0) + kMaxDZ * kMaxDZ * kParamCov(1,1) + 2. * kMaxDZ * TMath::Abs(kParamCov(0,1));
+  Double_t errY2 = kParamCov(2,2) + kMaxDZ * kMaxDZ * kParamCov(3,3) + 2. * kMaxDZ * TMath::Abs(kParamCov(2,3));
+  Double_t dX = TMath::Abs(trackParam.GetNonBendingSlope()) * kMaxDZ +
+               GetRecoParam()->GetMaxNonBendingDistanceToTrack() +
+               GetRecoParam()->GetSigmaCutForTracking() * TMath::Sqrt(errX2);
+  Double_t dY = TMath::Abs(trackParam.GetBendingSlope()) * kMaxDZ +
+               GetRecoParam()->GetMaxBendingDistanceToTrack() +
+               GetRecoParam()->GetSigmaCutForTracking() * TMath::Sqrt(errY2);
+  TVector2 dimensions(dX, dY);
   TVector2 position(extrapTrackParam.GetNonBendingCoor(), extrapTrackParam.GetBendingCoor());
-  TVector2 dimensions(AliMUONReconstructor::GetRecoParam()->GetMaxNonBendingDistanceToTrack() + kgMaxShift,
-                     AliMUONReconstructor::GetRecoParam()->GetMaxBendingDistanceToTrack() + kgMaxShift);
-  AliMpArea area2(position, dimensions);
+  AliMpArea area(position, dimensions);
   
   // ask to cluterize in the given area of the given chamber
-  fClusterServer.Clusterize(2*station+1, clusterStore, area2);
-  
-  // extrapolate track parameters to the first chamber of the station
-  AliMUONTrackExtrap::ExtrapToZ(&extrapTrackParam, AliMUONConstants::DefaultChamberZ(2*station));
-  
-  // build the searching area
-  position.Set(extrapTrackParam.GetNonBendingCoor(), extrapTrackParam.GetBendingCoor());
-//  dimensions.Set(AliMUONReconstructor::GetRecoParam()->GetMaxNonBendingDistanceToTrack() + kgMaxShift,
-//              AliMUONReconstructor::GetRecoParam()->GetMaxBendingDistanceToTrack() + kgMaxShift);
-  AliMpArea area1(position, dimensions);
+  fClusterServer->Clusterize(chamber, clusterStore, area);
   
-  // ask to cluterize in the given area of the given chamber
-  fClusterServer.Clusterize(2*station, clusterStore, area1);
+}
+
+  //__________________________________________________________________________
+void AliMUONVTrackReconstructor::AskForNewClustersInStation(const AliMUONTrackParam &trackParam,
+                                                           AliMUONVClusterStore& clusterStore, Int_t station)
+{
+  /// Ask the clustering to reconstruct new clusters around the track candidate position
+  /// in the 2 chambers of the given station
+  AskForNewClustersInChamber(trackParam, clusterStore, 2*station+1);
+  AskForNewClustersInChamber(trackParam, clusterStore, 2*station);
 }
 
   //__________________________________________________________________________
@@ -400,16 +426,24 @@ Double_t AliMUONVTrackReconstructor::TryOneCluster(const AliMUONTrackParam &trac
   //__________________________________________________________________________
 Bool_t AliMUONVTrackReconstructor::TryOneClusterFast(const AliMUONTrackParam &trackParam, AliMUONVCluster* cluster)
 {
-/// Test the compatibility between the track and the cluster within a wide fix window
-/// assuming linear propagation of the track:
+/// Test the compatibility between the track and the cluster
+/// given the track resolution + the maximum-distance-to-track value
+/// and assuming linear propagation of the track:
 /// return kTRUE if they are compatibles
   
   Double_t dZ = cluster->GetZ() - trackParam.GetZ();
   Double_t dX = cluster->GetX() - (trackParam.GetNonBendingCoor() + trackParam.GetNonBendingSlope() * dZ);
   Double_t dY = cluster->GetY() - (trackParam.GetBendingCoor() + trackParam.GetBendingSlope() * dZ);
+  const TMatrixD& kParamCov = trackParam.GetCovariances();
+  Double_t errX2 = kParamCov(0,0) + dZ * dZ * kParamCov(1,1) + 2. * dZ * kParamCov(0,1);
+  Double_t errY2 = kParamCov(2,2) + dZ * dZ * kParamCov(3,3) + 2. * dZ * kParamCov(2,3);
+
+  Double_t dXmax = GetRecoParam()->GetSigmaCutForTracking() * TMath::Sqrt(errX2) +
+                   GetRecoParam()->GetMaxNonBendingDistanceToTrack();
+  Double_t dYmax = GetRecoParam()->GetSigmaCutForTracking() * TMath::Sqrt(errY2) +
+                  GetRecoParam()->GetMaxBendingDistanceToTrack();
   
-  if (TMath::Abs(dX) > AliMUONReconstructor::GetRecoParam()->GetMaxNonBendingDistanceToTrack() ||
-      TMath::Abs(dY) > AliMUONReconstructor::GetRecoParam()->GetMaxBendingDistanceToTrack()) return kFALSE;
+  if (TMath::Abs(dX) > dXmax || TMath::Abs(dY) > dYmax) return kFALSE;
   
   return kTRUE;
   
@@ -459,6 +493,127 @@ Double_t AliMUONVTrackReconstructor::TryTwoClustersFast(const AliMUONTrackParam
 }
 
   //__________________________________________________________________________
+Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInChamber(AliMUONTrack &trackCandidate, const AliMUONVClusterStore& clusterStore,
+                                                             Int_t nextChamber)
+{
+  /// Follow trackCandidate in chamber(0..) nextChamber assuming linear propagation, and search for compatible cluster(s)
+  /// Keep all possibilities or only the best one(s) according to the flag fgkTrackAllTracks:
+  /// kTRUE:  duplicate "trackCandidate" if there are several possibilities and add the new tracks at the end of
+  ///         fRecTracksPtr to avoid conficts with other track candidates at this current stage of the tracking procedure.
+  ///         Remove the obsolete "trackCandidate" at the end.
+  /// kFALSE: add only the best cluster(s) to the "trackCandidate". Try to add a couple of clusters in priority.
+  /// 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 maxChi2WithOneCluster = 2. * GetRecoParam()->GetSigmaCutForTracking() *
+                                       GetRecoParam()->GetSigmaCutForTracking(); // 2 because 2 quantities in chi2
+  Double_t bestChi2WithOneCluster = maxChi2WithOneCluster;
+  Bool_t foundOneCluster = kFALSE;
+  AliMUONTrack *newTrack = 0x0;
+  AliMUONVCluster *cluster;
+  AliMUONTrackParam trackParam;
+  AliMUONTrackParam extrapTrackParamAtCluster;
+  AliMUONTrackParam bestTrackParamAtCluster;
+  
+  // Get track parameters according to the propagation direction
+  if (nextChamber > 7) trackParam = *(AliMUONTrackParam*)trackCandidate.GetTrackParamAtCluster()->Last();
+  else trackParam = *(AliMUONTrackParam*)trackCandidate.GetTrackParamAtCluster()->First();
+  
+  // Printout for debuging
+  if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 2) || (AliLog::GetGlobalDebugLevel() >= 2)) {
+    cout<<endl<<"Track parameters and covariances at first cluster:"<<endl;
+    trackParam.GetParameters().Print();
+    trackParam.GetCovariances().Print();
+  }
+  
+  // Add MCS effect
+  AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(),1.);
+  
+  // Printout for debuging
+  if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
+    cout << "FollowLinearTrackInChamber: look for cluster in chamber(1..): " << nextChamber+1 << endl;
+  }
+  
+  // Create iterators to loop over clusters in chamber
+  TIter next(clusterStore.CreateChamberIterator(nextChamber,nextChamber));
+  
+  // look for candidates in chamber
+  while ( ( cluster = static_cast<AliMUONVCluster*>(next()) ) ) {
+    
+    // try to add the current cluster fast
+    if (!TryOneClusterFast(trackParam, cluster)) continue;
+    
+    // try to add the current cluster accuratly
+    extrapTrackParamAtCluster = trackParam;
+    AliMUONTrackExtrap::LinearExtrapToZ(&extrapTrackParamAtCluster, cluster->GetZ());
+    chi2WithOneCluster = TryOneCluster(extrapTrackParamAtCluster, cluster, extrapTrackParamAtCluster);
+    
+    // if good chi2 then consider to add cluster
+    if (chi2WithOneCluster < maxChi2WithOneCluster) {
+      foundOneCluster = kTRUE;
+      
+      // Printout for debuging
+      if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
+       cout << "FollowLinearTrackInChamber: found one cluster in chamber(1..): " << nextChamber+1
+       << " (Chi2 = " << chi2WithOneCluster << ")" << endl;
+       cluster->Print();
+      }
+      
+      if (GetRecoParam()->TrackAllTracks()) {
+       // copy trackCandidate into a new track put at the end of fRecTracksPtr and add the new cluster
+       newTrack = new ((*fRecTracksPtr)[fRecTracksPtr->GetLast()+1]) AliMUONTrack(trackCandidate);
+       if (GetRecoParam()->RequestStation(nextChamber/2))
+         extrapTrackParamAtCluster.SetRemovable(kFALSE);
+       else extrapTrackParamAtCluster.SetRemovable(kTRUE);
+       newTrack->AddTrackParamAtCluster(extrapTrackParamAtCluster,*cluster);
+       fNRecTracks++;
+       
+       // Printout for debuging
+       if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
+         cout << "FollowLinearTrackInChamber: added one cluster in chamber(1..): " << nextChamber+1 << endl;
+         if (AliLog::GetGlobalDebugLevel() >= 3) newTrack->RecursiveDump();
+       }
+       
+      } else if (chi2WithOneCluster < bestChi2WithOneCluster) {
+       // keep track of the best cluster
+       bestChi2WithOneCluster = chi2WithOneCluster;
+       bestTrackParamAtCluster = extrapTrackParamAtCluster;
+      }
+      
+    }
+    
+  }
+  
+  // fill out the best track if required else clean up the fRecTracksPtr array
+  if (!GetRecoParam()->TrackAllTracks()) {
+    if (foundOneCluster) {
+      if (GetRecoParam()->RequestStation(nextChamber/2))
+       bestTrackParamAtCluster.SetRemovable(kFALSE);
+      else bestTrackParamAtCluster.SetRemovable(kTRUE);
+      trackCandidate.AddTrackParamAtCluster(bestTrackParamAtCluster,*(bestTrackParamAtCluster.GetClusterPtr()));
+      
+      // Printout for debuging
+      if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
+        cout << "FollowLinearTrackInChamber: added the best cluster in chamber(1..): " << bestTrackParamAtCluster.GetClusterPtr()->GetChamberId()+1 << endl;
+        if (AliLog::GetGlobalDebugLevel() >= 3) newTrack->RecursiveDump();
+      }
+      
+    } else return kFALSE;
+    
+  } else if (foundOneCluster) {
+    
+    // remove obsolete track
+    fRecTracksPtr->Remove(&trackCandidate);
+    fNRecTracks--;
+    
+  } else return kFALSE;
+  
+  return kTRUE;
+  
+}
+
+//__________________________________________________________________________
 Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trackCandidate, const AliMUONVClusterStore& clusterStore,
                                                              Int_t nextStation)
 {
@@ -485,16 +640,17 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
   
   Double_t chi2WithOneCluster = 1.e10;
   Double_t chi2WithTwoClusters = 1.e10;
-  Double_t maxChi2WithOneCluster = 2. * AliMUONReconstructor::GetRecoParam()->GetSigmaCutForTracking() *
-                                       AliMUONReconstructor::GetRecoParam()->GetSigmaCutForTracking(); // 2 because 2 quantities in chi2
-  Double_t maxChi2WithTwoClusters = 4. * AliMUONReconstructor::GetRecoParam()->GetSigmaCutForTracking() *
-                                        AliMUONReconstructor::GetRecoParam()->GetSigmaCutForTracking(); // 4 because 4 quantities in chi2
+  Double_t maxChi2WithOneCluster = 2. * GetRecoParam()->GetSigmaCutForTracking() *
+                                       GetRecoParam()->GetSigmaCutForTracking(); // 2 because 2 quantities in chi2
+  Double_t maxChi2WithTwoClusters = 4. * GetRecoParam()->GetSigmaCutForTracking() *
+                                        GetRecoParam()->GetSigmaCutForTracking(); // 4 because 4 quantities in chi2
   Double_t bestChi2WithOneCluster = maxChi2WithOneCluster;
   Double_t bestChi2WithTwoClusters = maxChi2WithTwoClusters;
   Bool_t foundOneCluster = kFALSE;
   Bool_t foundTwoClusters = kFALSE;
   AliMUONTrack *newTrack = 0x0;
   AliMUONVCluster *clusterCh1, *clusterCh2;
+  AliMUONTrackParam trackParam;
   AliMUONTrackParam extrapTrackParamAtCluster1;
   AliMUONTrackParam extrapTrackParamAtCluster2;
   AliMUONTrackParam bestTrackParamAtCluster1;
@@ -505,8 +661,16 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
   for (Int_t i = 0; i < nClusters; i++) clusterCh1Used[i] = kFALSE;
   Int_t iCluster1;
   
-  // Get track parameters
-  AliMUONTrackParam trackParam(*(AliMUONTrackParam*)trackCandidate.GetTrackParamAtCluster()->First());
+  // Get track parameters according to the propagation direction
+  if (nextStation==4) trackParam = *(AliMUONTrackParam*)trackCandidate.GetTrackParamAtCluster()->Last();
+  else trackParam = *(AliMUONTrackParam*)trackCandidate.GetTrackParamAtCluster()->First();
+  
+  // Printout for debuging
+  if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 2) || (AliLog::GetGlobalDebugLevel() >= 2)) {
+    cout<<endl<<"Track parameters and covariances at first cluster:"<<endl;
+    trackParam.GetParameters().Print();
+    trackParam.GetCovariances().Print();
+  }
   
   // Add MCS effect
   AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(),1.);
@@ -539,6 +703,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
       if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
         cout << "FollowLinearTrackInStation: found one cluster in chamber(1..): " << ch2+1
             << " (Chi2 = " << chi2WithOneCluster << ")" << endl;
+       clusterCh2->Print();
         cout << "                      look for second clusters in chamber(1..): " << ch1+1 << " ..." << endl;
       }
       
@@ -567,17 +732,17 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
          // Printout for debuging
          if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
            cout << "FollowLinearTrackInStation: found one cluster in chamber(1..): " << ch1+1
-                << " (Chi2 = " << chi2WithTwoClusters << ")" << endl;
+                << " (Chi2 = " << chi2WithTwoClusters << ")" << endl;
+           clusterCh1->Print();
          }
          
-         if (AliMUONReconstructor::GetRecoParam()->TrackAllTracks()) {
+         if (GetRecoParam()->TrackAllTracks()) {
            // copy trackCandidate into a new track put at the end of fRecTracksPtr and add the new clusters
             newTrack = new ((*fRecTracksPtr)[fRecTracksPtr->GetLast()+1]) AliMUONTrack(trackCandidate);
            extrapTrackParamAtCluster1.SetRemovable(kTRUE);
            newTrack->AddTrackParamAtCluster(extrapTrackParamAtCluster1,*clusterCh1);
            extrapTrackParamAtCluster2.SetRemovable(kTRUE);
            newTrack->AddTrackParamAtCluster(extrapTrackParamAtCluster2,*clusterCh2);
-           newTrack->GetTrackParamAtCluster()->Sort();
            fNRecTracks++;
            
            // Tag clusterCh1 as used
@@ -604,12 +769,13 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
       if (!foundSecondCluster) {
        foundOneCluster = kTRUE;
         
-       if (AliMUONReconstructor::GetRecoParam()->TrackAllTracks()) {
+       if (GetRecoParam()->TrackAllTracks()) {
          // copy trackCandidate into a new track put at the end of fRecTracksPtr and add the new cluster
           newTrack = new ((*fRecTracksPtr)[fRecTracksPtr->GetLast()+1]) AliMUONTrack(trackCandidate);
-         extrapTrackParamAtCluster2.SetRemovable(kFALSE);
+         if (GetRecoParam()->RequestStation(nextStation))
+           extrapTrackParamAtCluster2.SetRemovable(kFALSE);
+         else extrapTrackParamAtCluster2.SetRemovable(kTRUE);
          newTrack->AddTrackParamAtCluster(extrapTrackParamAtCluster2,*clusterCh2);
-         newTrack->GetTrackParamAtCluster()->Sort();
          fNRecTracks++;
          
          // Printout for debuging
@@ -632,7 +798,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
   
   // look for candidates in chamber 1 not already attached to a track
   // if we want to keep all possible tracks or if no good couple of clusters has been found
-  if (AliMUONReconstructor::GetRecoParam()->TrackAllTracks() || !foundTwoClusters) {
+  if (GetRecoParam()->TrackAllTracks() || !foundTwoClusters) {
     
     // Printout for debuging
     if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
@@ -644,7 +810,7 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
     
     // add MCS effect for next step
     AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(),1.);
-      
+    
     // reset cluster iterator of chamber 1
     nextInCh1.Reset();
     iCluster1 = -1;
@@ -657,12 +823,12 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
       
       // try to add the current cluster fast
       if (!TryOneClusterFast(trackParam, clusterCh1)) continue;
-       
+      
       // try to add the current cluster accuratly
       extrapTrackParamAtCluster1 = trackParam;
       AliMUONTrackExtrap::LinearExtrapToZ(&extrapTrackParamAtCluster1, clusterCh1->GetZ());
       chi2WithOneCluster = TryOneCluster(extrapTrackParamAtCluster1, clusterCh1, extrapTrackParamAtCluster1);
-    
+      
       // if good chi2 then consider to add clusterCh1
       // We do not try to attach a cluster in the other chamber too since it has already been done above
       if (chi2WithOneCluster < maxChi2WithOneCluster) {
@@ -672,14 +838,16 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
        if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
          cout << "FollowLinearTrackInStation: found one cluster in chamber(1..): " << ch1+1
               << " (Chi2 = " << chi2WithOneCluster << ")" << endl;
+         clusterCh1->Print();
        }
        
-       if (AliMUONReconstructor::GetRecoParam()->TrackAllTracks()) {
+       if (GetRecoParam()->TrackAllTracks()) {
          // copy trackCandidate into a new track put at the end of fRecTracksPtr and add the new cluster
          newTrack = new ((*fRecTracksPtr)[fRecTracksPtr->GetLast()+1]) AliMUONTrack(trackCandidate);
-         extrapTrackParamAtCluster1.SetRemovable(kFALSE);
+         if (GetRecoParam()->RequestStation(nextStation))
+           extrapTrackParamAtCluster1.SetRemovable(kFALSE);
+         else extrapTrackParamAtCluster1.SetRemovable(kTRUE);
          newTrack->AddTrackParamAtCluster(extrapTrackParamAtCluster1,*clusterCh1);
-         newTrack->GetTrackParamAtCluster()->Sort();
          fNRecTracks++;
          
          // Printout for debuging
@@ -701,13 +869,12 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
   }
   
   // fill out the best track if required else clean up the fRecTracksPtr array
-  if (!AliMUONReconstructor::GetRecoParam()->TrackAllTracks()) {
+  if (!GetRecoParam()->TrackAllTracks()) {
     if (foundTwoClusters) {
       bestTrackParamAtCluster1.SetRemovable(kTRUE);
       trackCandidate.AddTrackParamAtCluster(bestTrackParamAtCluster1,*(bestTrackParamAtCluster1.GetClusterPtr()));
       bestTrackParamAtCluster2.SetRemovable(kTRUE);
       trackCandidate.AddTrackParamAtCluster(bestTrackParamAtCluster2,*(bestTrackParamAtCluster2.GetClusterPtr()));
-      trackCandidate.GetTrackParamAtCluster()->Sort();
       
       // Printout for debuging
       if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
@@ -716,9 +883,10 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
       }
       
     } else if (foundOneCluster) {
-      bestTrackParamAtCluster1.SetRemovable(kFALSE);
+      if (GetRecoParam()->RequestStation(nextStation))
+       bestTrackParamAtCluster1.SetRemovable(kFALSE);
+      else bestTrackParamAtCluster1.SetRemovable(kTRUE);
       trackCandidate.AddTrackParamAtCluster(bestTrackParamAtCluster1,*(bestTrackParamAtCluster1.GetClusterPtr()));
-      trackCandidate.GetTrackParamAtCluster()->Sort();
       
       // Printout for debuging
       if ((AliLog::GetDebugLevel("MUON","AliMUONVTrackReconstructor") >= 1) || (AliLog::GetGlobalDebugLevel() >= 1)) {
@@ -747,173 +915,65 @@ Bool_t AliMUONVTrackReconstructor::FollowLinearTrackInStation(AliMUONTrack &trac
   
 }
 
-  //__________________________________________________________________________
-void AliMUONVTrackReconstructor::ValidateTracksWithTrigger(AliMUONVTrackStore& trackStore,
-                                                           const AliMUONVTriggerTrackStore& triggerTrackStore,
-                                                           const AliMUONVTriggerStore& triggerStore,
-                                                           const AliMUONTrackHitPattern& trackHitPattern)
+//__________________________________________________________________________
+void AliMUONVTrackReconstructor::ImproveTracks()
 {
-  /// Try to match track from tracking system with trigger track
-  AliCodeTimerAuto("");
-
-  const Double_t kDeltaZ = TMath::Abs(AliMUONConstants::DefaultChamberZ(12) - AliMUONConstants::DefaultChamberZ(10));
-
-  //static const Double_t kDistSigma[3]={1,1,0.02}; // sigma of distributions (trigger-track) X,Y,slopeY
-  // sigma of distributions (trigger-track) X,Y,slopeY
-  const Double_t kDistSigma[3]={AliMUONConstants::TriggerNonBendingReso(),
-                               AliMUONConstants::TriggerBendingReso(),
-                               1.414 * AliMUONConstants::TriggerBendingReso()/kDeltaZ};
-
-  const Double_t kTrigNonBendReso = AliMUONConstants::TriggerNonBendingReso();
-  const Double_t kTrigBendReso = AliMUONConstants::TriggerBendingReso();
-  const Double_t kTrigSlopeBendReso = 1.414 * AliMUONConstants::TriggerBendingReso()/kDeltaZ;
-  const Double_t kTrigCovSlopeBend = - kTrigBendReso * kTrigBendReso / kDeltaZ;
-
-  // Covariance matrix 3x3 (X,Y,slopeY) for trigger tracks
-  TMatrixD trigCov(3,3);
-  trigCov.Zero();
-  trigCov(0,0) = kTrigNonBendReso * kTrigNonBendReso;
-  trigCov(1,1) = kTrigBendReso * kTrigBendReso;
-  trigCov(2,2) = kTrigSlopeBendReso * kTrigSlopeBendReso;
-  trigCov(1,2) = trigCov(2,1) = kTrigCovSlopeBend;
-
-  Int_t matchTrigger;
-  Int_t loTrgNum(-1);
-  Double_t distTriggerTrack[3], sigma2[3];
-  Double_t xTrack, yTrack, ySlopeTrack, chi2MatchTrigger, minChi2MatchTrigger, chi2;
-
-  TIter itTrack(trackStore.CreateIterator());
-  AliMUONTrack* track;
-
-  const Float_t kZFilterOut = AliMUONConstants::MuonFilterZEnd();
-  const Float_t kFilterThickness = TMath::Abs(kZFilterOut-AliMUONConstants::MuonFilterZBeg()); // cm
-  const Int_t kFirstTrigCh = AliMUONConstants::NTrackingCh();
+  /// Improve tracks by removing clusters with local chi2 highter than the defined cut
+  /// Recompute track parameters and covariances at the remaining clusters
+  AliDebug(1,"Enter ImproveTracks");
   
-  while ( ( track = static_cast<AliMUONTrack*>(itTrack()) ) )
-  {
-    matchTrigger = 0;
-    chi2MatchTrigger = 0.;
-    loTrgNum = -1;
-    Int_t doubleMatch=-1; // Check if track matches 2 trigger tracks
-    Double_t doubleChi2 = -1.;
+  AliMUONTrack *track, *nextTrack;
+  
+  track = (AliMUONTrack*) fRecTracksPtr->First();
+  while (track) {
     
-    AliMUONTrackParam trackParam(*((AliMUONTrackParam*) (track->GetTrackParamAtCluster()->Last())));
-
-    AliMUONTrackExtrap::ExtrapToZCov(&trackParam, kZFilterOut); // Extrap to muon filter end
-    AliMUONTrackExtrap::AddMCSEffect(&trackParam, kFilterThickness, AliMUONConstants::MuonFilterX0()); // Add MCS effects
-    AliMUONTrackExtrap::ExtrapToZCov(&trackParam, AliMUONConstants::DefaultChamberZ(kFirstTrigCh)); // extrap to 1st trigger chamber
-
-    const TMatrixD& kParamCov = trackParam.GetCovariances();
+    // prepare next track in case the actual track is suppressed
+    nextTrack = (AliMUONTrack*) fRecTracksPtr->After(track);
     
-    xTrack = trackParam.GetNonBendingCoor();
-    yTrack = trackParam.GetBendingCoor();
-    ySlopeTrack = trackParam.GetBendingSlope();
-
-    // Covariance matrix 3x3 (X,Y,slopeY) for tracker tracks
-    TMatrixD trackCov(3,3);
-    trackCov.Zero();
-    trackCov(0,0) = kParamCov(0,0);
-    trackCov(1,1) = kParamCov(2,2);
-    trackCov(2,2) = kParamCov(3,3);
-    trackCov(1,2) = kParamCov(2,3);
-    trackCov(2,1) = kParamCov(3,2);
-
-    TMatrixD sumCov(trackCov,TMatrixD::kPlus,trigCov);
-
-    Bool_t isCovOK = kTRUE;
-
-    if (sumCov.Determinant() != 0) {
-      sumCov.Invert();
-    } else {
-      AliWarning(" Determinant = 0");
-      isCovOK = kFALSE;
-      sigma2[0] = kParamCov(0,0);
-      sigma2[1] = kParamCov(2,2);
-      sigma2[2] = kParamCov(3,3);
-      for (Int_t iVar = 0; iVar < 3; iVar++) sigma2[iVar] += kDistSigma[iVar] * kDistSigma[iVar];
+    ImproveTrack(*track);
+    
+    // remove track if improvement failed
+    if (!track->IsImproved()) {
+      fRecTracksPtr->Remove(track);
+      fNRecTracks--;
     }
+    
+    track = nextTrack;
+  }
+  
+  // compress the array in case of some tracks have been removed
+  fRecTracksPtr->Compress();
+  
+}
 
-    minChi2MatchTrigger = 999.;
-
-    AliMUONTriggerTrack *triggerTrack;
-    TIter itTriggerTrack(triggerTrackStore.CreateIterator());
-    while ( ( triggerTrack = static_cast<AliMUONTriggerTrack*>(itTriggerTrack() ) ) )
-    {
-      distTriggerTrack[0] = triggerTrack->GetX11()-xTrack;
-      distTriggerTrack[1] = triggerTrack->GetY11()-yTrack;
-      distTriggerTrack[2] = TMath::Tan(triggerTrack->GetThetay())-ySlopeTrack;
-
-      if(isCovOK){
-       TMatrixD paramDiff(3,1);
-       for(Int_t iVar = 0; iVar < 3; iVar++)
-         paramDiff(iVar,0) = distTriggerTrack[iVar];
-       
-       TMatrixD tmp(sumCov,TMatrixD::kMult,paramDiff);
-       TMatrixD chi2M(paramDiff,TMatrixD::kTransposeMult,tmp);
-       chi2 = chi2M(0,0);
-      }
-      else {
-       chi2 = 0.;
-       for (Int_t iVar = 0; iVar < 3; iVar++) chi2 += distTriggerTrack[iVar]*distTriggerTrack[iVar]/sigma2[iVar];
-      }
-
-      chi2 /= 3.; // Normalized Chi2: 3 degrees of freedom (X,Y,slopeY)
-      if (chi2 < AliMUONReconstructor::GetRecoParam()->GetMaxNormChi2MatchTrigger()) 
-      {
-        Bool_t isDoubleTrack = (TMath::Abs(chi2 - minChi2MatchTrigger)<1.);
-        if (chi2 < minChi2MatchTrigger && chi2 < AliMUONReconstructor::GetRecoParam()->GetMaxNormChi2MatchTrigger()) 
-        {
-          if(isDoubleTrack)
-          {
-            doubleMatch = loTrgNum;
-            doubleChi2 = chi2MatchTrigger;
-          }
-          minChi2MatchTrigger = chi2;
-          chi2MatchTrigger = chi2;
-          loTrgNum = triggerTrack->GetLoTrgNum();
-          AliMUONLocalTrigger* locTrg = triggerStore.FindLocal(loTrgNum);
-          matchTrigger=1;
-          if(locTrg->LoLpt()>0)matchTrigger=2;
-          if(locTrg->LoHpt()>0)matchTrigger=3;
-        }
-        else if(isDoubleTrack) 
-        {
-          doubleMatch = triggerTrack->GetLoTrgNum();
-          doubleChi2 = chi2;
-        }
-      }
-    }
-    if(doubleMatch>=0)
-    { // If two trigger tracks match, select the one passing more trigger cuts
-      AliDebug(1, Form("Two candidates found: %i and %i",loTrgNum,doubleMatch));
-      AliMUONLocalTrigger* locTrg1 = triggerStore.FindLocal(doubleMatch);
-      if((locTrg1->LoLpt()>0 && matchTrigger<2) || (locTrg1->LoHpt() && matchTrigger<3))
-      {
-        if(locTrg1->LoHpt()>0)matchTrigger=3;
-        else matchTrigger=2;
-        loTrgNum = doubleMatch;
-        chi2MatchTrigger=doubleChi2;
-      }
-    }
+//__________________________________________________________________________
+void AliMUONVTrackReconstructor::Finalize()
+{
+  /// Recompute track parameters and covariances at each attached cluster from those at the first one
+  
+  AliMUONTrack *track;
+  
+  track = (AliMUONTrack*) fRecTracksPtr->First();
+  while (track) {
     
-    track->SetMatchTrigger(matchTrigger);
-    track->SetLoTrgNum(loTrgNum);
-    track->SetChi2MatchTrigger(chi2MatchTrigger);
+    FinalizeTrack(*track);
     
-    AliMUONLocalTrigger* locTrg = static_cast<AliMUONLocalTrigger*>(triggerStore.FindLocal(loTrgNum));
+    track = (AliMUONTrack*) fRecTracksPtr->After(track);
     
-    if (locTrg)
-    {    
-      track->SetLocalTrigger(locTrg->LoCircuit(),
-                             locTrg->LoStripX(),
-                             locTrg->LoStripY(),
-                             locTrg->LoDev(),
-                             locTrg->LoLpt(),
-                             locTrg->LoHpt());
-    }    
-  }  
+  }
   
-  trackHitPattern.GetHitPattern(trackStore,triggerStore);
+}
+
+//__________________________________________________________________________
+void AliMUONVTrackReconstructor::ValidateTracksWithTrigger(AliMUONVTrackStore& trackStore,
+                                                           const AliMUONVTriggerTrackStore& triggerTrackStore,
+                                                           const AliMUONVTriggerStore& triggerStore,
+                                                           const AliMUONTrackHitPattern& trackHitPattern)
+{
+  /// Try to match track from tracking system with trigger track
+  AliCodeTimerAuto("");
+
+  trackHitPattern.ExecuteValidation(trackStore, triggerTrackStore, triggerStore);
 }
 
   //__________________________________________________________________________
@@ -944,16 +1004,10 @@ void AliMUONVTrackReconstructor::EventReconstructTrigger(const AliMUONTriggerCir
   
   while ( ( locTrg = static_cast<AliMUONLocalTrigger*>(next()) ) )
   {
-    Bool_t xTrig=kFALSE;
-    Bool_t yTrig=kFALSE;
+    Bool_t xTrig=locTrg->IsTrigX();
+    Bool_t yTrig=locTrg->IsTrigY();
     
     Int_t localBoardId = locTrg->LoCircuit();
-    if ( locTrg->LoSdev()==1 && locTrg->LoDev()==0 && 
-         locTrg->LoStripX()==0) xTrig=kFALSE; // no trigger in X
-    else xTrig=kTRUE;                         // trigger in X
-    if (locTrg->LoTrigY()==1 && 
-        locTrg->LoStripY()==15 ) yTrig = kFALSE; // no trigger in Y
-    else yTrig = kTRUE;                          // trigger in Y
     
     if (xTrig && yTrig) 
     { // make Trigger Track if trigger in X and Y
@@ -961,13 +1015,7 @@ void AliMUONVTrackReconstructor::EventReconstructTrigger(const AliMUONTriggerCir
       Float_t y11 = circuit.GetY11Pos(localBoardId, locTrg->LoStripX()); 
       // need first to convert deviation to [0-30] 
       // (see AliMUONLocalTriggerBoard::LocalTrigger)
-      Int_t deviation = locTrg->LoDev(); 
-      Int_t sign = 0;
-      if ( !locTrg->LoSdev() &&  deviation ) sign=-1;
-      if ( !locTrg->LoSdev() && !deviation ) sign= 0;
-      if (  locTrg->LoSdev() == 1 )          sign=+1;
-      deviation *= sign;
-      deviation += 15;
+      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());