Redesigning Original tracking classes (Philippe Pillot)
authormartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 25 Oct 2006 15:53:08 +0000 (15:53 +0000)
committermartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 25 Oct 2006 15:53:08 +0000 (15:53 +0000)
17 files changed:
MUON/AliMUONConstants.cxx
MUON/AliMUONConstants.h
MUON/AliMUONHitForRec.cxx
MUON/AliMUONHitForRec.h
MUON/AliMUONRecoCheck.cxx
MUON/AliMUONTrack.cxx
MUON/AliMUONTrack.h
MUON/AliMUONTrackHit.cxx [deleted file]
MUON/AliMUONTrackHit.h [deleted file]
MUON/AliMUONTrackK.cxx
MUON/AliMUONTrackParam.cxx
MUON/AliMUONTrackParam.h
MUON/AliMUONTrackReconstructor.cxx
MUON/AliMUONTrackReconstructor.h
MUON/AlirootRun_MUONtest.sh
MUON/MUONrecLinkDef.h
MUON/libMUONrec.pkg

index 84fcf96..6eb28db 100644 (file)
@@ -38,6 +38,7 @@ ClassImp(AliMUONConstants)
 
 Int_t   AliMUONConstants::fgNCh = 14;
 Int_t   AliMUONConstants::fgNTrackingCh = 10;
+Int_t   AliMUONConstants::fgNTrackingSt = 5;
 Int_t   AliMUONConstants::fgNTriggerCh = 4;
 Int_t   AliMUONConstants::fgNTriggerCircuit = 234;
 Int_t   AliMUONConstants::fgNDetElem = 228;
@@ -79,6 +80,9 @@ Float_t  AliMUONConstants::fgDmax[7]  = {  176.6, 229.0, 308.84, 418.2,  522.0,
  
 Int_t    AliMUONConstants::fgMaxZoom = 20;
 
+// Defaults parameters for track reconstruction
+Double_t  AliMUONConstants::fgDefaultChamberThicknessInX0 = 0.03;
+
 //______________________________________________________________________________
 Int_t AliMUONConstants::ChamberNumber(Float_t z) 
 {
index 2f66c35..2710dc2 100644 (file)
@@ -18,6 +18,8 @@ class AliMUONConstants : public TObject {
     static Int_t    NCh() {return fgNCh;}
     /// Return number of tracking chambers
     static Int_t    NTrackingCh() {return fgNTrackingCh;}
+    /// Return number of tracking stations
+    static Int_t    NTrackingSt() {return fgNTrackingSt;}
     /// Return number of trigger chambers
     static Int_t    NTriggerCh() {return fgNTriggerCh;}
     /// Return number of trigger circuits
@@ -66,6 +68,8 @@ class AliMUONConstants : public TObject {
     static Float_t Pitch()    {return fgPitch;}
     /// Return wire pitch for Station 1 & 2
     static Float_t PitchSt1() {return fgPitchSt1;}
+    /// Return chamber thickness in X0
+    static Double_t DefaultChamberThicknessInX0() {return fgDefaultChamberThicknessInX0;}
     /// Return Trigger ToF Limit (75 ns)
     static Float_t TriggerTofLimit() {return fgkTriggerTofLimit;}
  
@@ -78,6 +82,7 @@ class AliMUONConstants : public TObject {
  private:
     static Int_t  fgNCh;                ///<  Number of Chambers    
     static Int_t  fgNTrackingCh;        ///<  Number of Tracking Chambers
+    static Int_t  fgNTrackingSt;        ///<  Number of Tracking Stations
     static Int_t  fgNTriggerCh;         ///<  Number of Trigger Chambers
     static Int_t  fgNTriggerCircuit;    ///<  Number of Trigger Circuits
     static Int_t  fgNDetElem;           ///<  Number of Detection Elements.
@@ -102,6 +107,8 @@ class AliMUONConstants : public TObject {
     static Float_t  fgPitch;             ///< Wire pitch for St2 & Slats
     static Float_t  fgPitchSt1;          ///< Wire pitch for Station 1
 
+    static Double_t  fgDefaultChamberThicknessInX0; ///< default chamber thickness in X0 for reconstruction
+    
     static Int_t    fgMaxZoom;           ///< Maximum Zoom for event display
     static Float_t fgkTriggerTofLimit;   ///< Particle above this threshold are discarded in trigger algorithm
     
index 6f13bdd..90542b4 100644 (file)
@@ -46,8 +46,6 @@ AliMUONHitForRec::AliMUONHitForRec()
     fTrackRefSignal(0),
     fIndexOfFirstSegment(0),
     fNSegments(0),
-    fFirstTrackHitPtr(0x0),
-    fLastTrackHitPtr(0x0),
     fNTrackHits(0)
 {
 /// Default Constructor
@@ -69,8 +67,6 @@ AliMUONHitForRec::AliMUONHitForRec(AliTrackReference* theGhit)
     fTrackRefSignal(0),
     fIndexOfFirstSegment(-1),
     fNSegments(0),
-    fFirstTrackHitPtr(0x0),
-    fLastTrackHitPtr(0x0),
     fNTrackHits(0)
 {
 /// Constructor for AliMUONHitForRec from a track ref. hit.
@@ -119,8 +115,6 @@ AliMUONHitForRec::AliMUONHitForRec(AliMUONRawCluster* theRawCluster)
     fTrackRefSignal(-1),
     fIndexOfFirstSegment(-1),
     fNSegments(0),
-    fFirstTrackHitPtr(0x0),
-    fLastTrackHitPtr(0x0),
     fNTrackHits(0)
 {
 /// Constructor for AliMUONHitForRec from a raw cluster.
@@ -146,8 +140,6 @@ AliMUONHitForRec::AliMUONHitForRec (const AliMUONHitForRec& theMUONHitForRec)
     fTrackRefSignal(theMUONHitForRec.fTrackRefSignal),
     fIndexOfFirstSegment(theMUONHitForRec.fIndexOfFirstSegment),
     fNSegments(theMUONHitForRec.fNSegments),
-    fFirstTrackHitPtr(theMUONHitForRec.fFirstTrackHitPtr),
-    fLastTrackHitPtr(theMUONHitForRec.fLastTrackHitPtr),
     fNTrackHits(theMUONHitForRec.fNTrackHits)
 {
 /// Copy constructor
@@ -171,8 +163,6 @@ AliMUONHitForRec & AliMUONHitForRec::operator=(const AliMUONHitForRec& theMUONHi
   fTrackRefSignal = theMUONHitForRec.fTrackRefSignal;
   fIndexOfFirstSegment = theMUONHitForRec.fIndexOfFirstSegment;
   fNSegments = theMUONHitForRec.fNSegments;
-  fFirstTrackHitPtr = theMUONHitForRec.fFirstTrackHitPtr;
-  fLastTrackHitPtr = theMUONHitForRec.fLastTrackHitPtr;
   fNTrackHits = theMUONHitForRec.fNTrackHits;
   return *this;
 }
index 395aabf..02c034b 100644 (file)
@@ -51,10 +51,6 @@ class AliMUONHitForRec : public TObject {
   void SetIndexOfFirstSegment(Int_t IndexOfFirstSegment) { fIndexOfFirstSegment = IndexOfFirstSegment;}
   Int_t GetNSegments(void) const { return fNSegments;}
   void SetNSegments(Int_t NSegments) { fNSegments = NSegments;}
-  AliMUONTrackHit* GetFirstTrackHitPtr(void) const { return fFirstTrackHitPtr;}
-  void SetFirstTrackHitPtr(AliMUONTrackHit* FirstTrackHitPtr) { fFirstTrackHitPtr = FirstTrackHitPtr;}
-  AliMUONTrackHit* GetLastTrackHitPtr(void) const { return fLastTrackHitPtr;}
-  void SetLastTrackHitPtr(AliMUONTrackHit* LastTrackHitPtr) { fLastTrackHitPtr = LastTrackHitPtr;}
   Int_t GetNTrackHits(void) const { return fNTrackHits;}
   void SetNTrackHits(Int_t NTrackHits) { fNTrackHits = NTrackHits;}
 
@@ -88,9 +84,6 @@ class AliMUONHitForRec : public TObject {
   Int_t fIndexOfFirstSegment; //!<  index of first Segment
   Int_t fNSegments; //!<  number of Segments
 
-  // links forward to reconstructed track hits
-  AliMUONTrackHit *fFirstTrackHitPtr ; //!<  pointer to first TrackHit made with HitForRec
-  AliMUONTrackHit *fLastTrackHitPtr ; //!<  pointer to last TrackHit made with HitForRec
   Int_t fNTrackHits; //!<  number of TrackHit's made with HitForRec
   
   ClassDef(AliMUONHitForRec, 2) // Hit for reconstruction in ALICE dimuon spectrometer
index a6a3ca6..5fd1199 100644 (file)
@@ -179,7 +179,7 @@ void AliMUONRecoCheck::MakeTrackRef()
        else iChamber = AliMUONConstants::ChamberNumber(z);
        hitForRec->SetChamberNumber(iChamber);
 
-       muonTrack->AddTrackParamAtHit(trackParam);
+       muonTrack->AddTrackParamAtHit(trackParam,hitForRec);
        muonTrack->AddHitForRecAtHit(hitForRec);
        muonTrack->SetTrackID(track);
 
@@ -386,7 +386,7 @@ void AliMUONRecoCheck::CleanMuonTrackRef()
        trackParam->SetInverseBendingMomentum(1./bendingMomentum);
 
       trackNew->AddHitForRecAtHit(hitForRec);
-      trackNew->AddTrackParamAtHit(trackParam);
+      trackNew->AddTrackParamAtHit(trackParam,hitForRec);
       
       iHit1++;
     } // end iHit1
index 1d552e6..510c15e 100644 (file)
 
 #include "AliMUONTrack.h"
 
-#include "AliMUONTrackReconstructor.h" 
+#include "AliMUONTrackParam.h" 
 #include "AliMUONHitForRec.h" 
 #include "AliMUONSegment.h" 
-#include "AliMUONTrackHit.h"
 #include "AliMUONTriggerTrack.h"
 #include "AliMUONConstants.h"
 
 #include <TMath.h>
 #include <TMatrixD.h>
 #include <TObjArray.h>
-#include <TVirtualFitter.h>
 
 #include <stdlib.h> // for exit()
 
-// Functions to be minimized with Minuit
-void TrackChi2(Int_t &NParam, Double_t *Gradient, Double_t &Chi2, Double_t *Param, Int_t Flag);
-void TrackChi2MCS(Int_t &NParam, Double_t *Gradient, Double_t &Chi2, Double_t *Param, Int_t Flag);
-
-void mnvertLocal(Double_t* a, Int_t l, Int_t m, Int_t n, Int_t& ifail);
-
-Double_t MultipleScatteringAngle2(AliMUONTrackHit *TrackHit);
-
-TVirtualFitter* AliMUONTrack::fgFitter = NULL; 
-
 ClassImp(AliMUONTrack) // Class implementation in ROOT context
 
 //__________________________________________________________________________
 AliMUONTrack::AliMUONTrack()
   : TObject(),
-    fTrackReconstructor(0x0),
     fTrackParamAtVertex(),
     fTrackParamAtHit(0x0),
     fHitForRecAtHit(0x0),
-    fTrackHitsPtr(0x0),
     fNTrackHits(0),
-    fFitMCS(0),
-    fFitNParam(0),
-    fFitStart(0),
     fFitFMin(-1.),
     fMatchTrigger(kFALSE),
     fChi2MatchTrigger(0.),
     fTrackID(0)
 {
   // Default constructor
-  fgFitter = 0;
 }
 
   //__________________________________________________________________________
-AliMUONTrack::AliMUONTrack(AliMUONSegment* BegSegment, AliMUONSegment* EndSegment, AliMUONTrackReconstructor* TrackReconstructor)
+AliMUONTrack::AliMUONTrack(AliMUONSegment* BegSegment, AliMUONSegment* EndSegment)
   : TObject(),
-    fTrackReconstructor(TrackReconstructor),
     fTrackParamAtVertex(),
     fTrackParamAtHit(0x0),
     fHitForRecAtHit(0x0),
-    fTrackHitsPtr(new TObjArray(10)),
     fNTrackHits(0),
-    fFitMCS(0),
-    fFitNParam(3),
-    fFitStart(1),
     fFitFMin(-1.),
     fMatchTrigger(kFALSE),
     fChi2MatchTrigger(0.),
@@ -96,29 +73,25 @@ AliMUONTrack::AliMUONTrack(AliMUONSegment* BegSegment, AliMUONSegment* EndSegmen
 {
   // Constructor from two Segment's
 
-  if (BegSegment) { //AZ
-    AddSegment(BegSegment); // add hits from BegSegment
-    AddSegment(EndSegment); // add hits from EndSegment
-    fTrackHitsPtr->Sort(); // sort TrackHits according to increasing Z
-    SetTrackParamAtVertex(); // set track parameters at vertex
-  }
   fTrackParamAtHit = new TClonesArray("AliMUONTrackParam",10);
   fHitForRecAtHit = new TClonesArray("AliMUONHitForRec",10);
-  return;
+  
+  if (BegSegment) { //AZ
+    AddTrackParamAtHit(0,BegSegment->GetHitForRec1());
+    AddTrackParamAtHit(0,BegSegment->GetHitForRec2());
+    AddTrackParamAtHit(0,EndSegment->GetHitForRec1());
+    AddTrackParamAtHit(0,EndSegment->GetHitForRec2());
+    fTrackParamAtHit->Sort(); // sort TrackParamAtHit according to increasing Z
+  }
 }
 
   //__________________________________________________________________________
-AliMUONTrack::AliMUONTrack(AliMUONSegment* Segment, AliMUONHitForRec* HitForRec, AliMUONTrackReconstructor* TrackReconstructor)
+AliMUONTrack::AliMUONTrack(AliMUONSegment* Segment, AliMUONHitForRec* HitForRec)
   : TObject(),
-    fTrackReconstructor(TrackReconstructor),
     fTrackParamAtVertex(),
     fTrackParamAtHit(0x0),
     fHitForRecAtHit(0x0),
-    fTrackHitsPtr(new TObjArray(10)),
     fNTrackHits(0),
-    fFitMCS(0),
-    fFitNParam(3),
-    fFitStart(1),
     fFitFMin(-1.),
     fMatchTrigger(kFALSE),
     fChi2MatchTrigger(0.),
@@ -126,25 +99,19 @@ AliMUONTrack::AliMUONTrack(AliMUONSegment* Segment, AliMUONHitForRec* HitForRec,
 {
   // Constructor from one Segment and one HitForRec
 
-  AddSegment(Segment); // add hits from Segment
-  AddHitForRec(HitForRec); // add HitForRec
-  fTrackHitsPtr->Sort(); // sort TrackHits according to increasing Z
-  SetTrackParamAtVertex(); // set track parameters at vertex
   fTrackParamAtHit = new TClonesArray("AliMUONTrackParam",10);
   fHitForRecAtHit = new TClonesArray("AliMUONHitForRec",10);
-  return;
+  
+  AddTrackParamAtHit(0,Segment->GetHitForRec1());
+  AddTrackParamAtHit(0,Segment->GetHitForRec2());
+  AddTrackParamAtHit(0,HitForRec);
+  fTrackParamAtHit->Sort(); // sort TrackParamAtHit according to increasing Z
 }
 
   //__________________________________________________________________________
 AliMUONTrack::~AliMUONTrack()
 {
   // Destructor
-  if (fTrackHitsPtr) {
-    // delete the TObjArray of pointers to TrackHit's   
-    delete fTrackHitsPtr; 
-    fTrackHitsPtr = NULL;
-  }
-  
   if (fTrackParamAtHit) {
     // delete the TClonesArray of pointers to TrackParam
     delete fTrackParamAtHit;
@@ -161,51 +128,33 @@ AliMUONTrack::~AliMUONTrack()
   //__________________________________________________________________________
 AliMUONTrack::AliMUONTrack (const AliMUONTrack& theMUONTrack)
   : TObject(theMUONTrack),
-    fTrackReconstructor(theMUONTrack.fTrackReconstructor),
     fTrackParamAtVertex(theMUONTrack.fTrackParamAtVertex),
     fTrackParamAtHit(0x0),
     fHitForRecAtHit(0x0),
-    fTrackHitsPtr(0x0),
     fNTrackHits(theMUONTrack.fNTrackHits),
-    fFitMCS(theMUONTrack.fFitMCS),
-    fFitNParam(theMUONTrack.fFitNParam),
-    fFitStart(theMUONTrack.fFitStart),
     fFitFMin(theMUONTrack.fFitFMin),
     fMatchTrigger(theMUONTrack.fMatchTrigger),
     fChi2MatchTrigger(theMUONTrack.fChi2MatchTrigger),
     fTrackID(theMUONTrack.fTrackID)
 {
-  //fTrackReconstructor = new AliMUONTrackReconstructor(*MUONTrack.fTrackReconstructor);
-                               // is it right ?
-                              // NO, because it would use dummy copy constructor
-                              // and AliMUONTrack is not the owner of its TrackReconstructor 
-  //fTrackParamAtVertex = theMUONTrack.fTrackParamAtVertex;
-
- // necessary to make a copy of the objects and not only the pointers in TObjArray.
-  if (theMUONTrack.fTrackHitsPtr) {
-    fTrackHitsPtr  =  new TObjArray(10);
-    for (Int_t index = 0; index < (theMUONTrack.fTrackHitsPtr)->GetEntriesFast(); index++) {
-      AliMUONTrackHit *trackHit = new AliMUONTrackHit(*(AliMUONTrackHit*)(theMUONTrack.fTrackHitsPtr)->At(index));
-      fTrackHitsPtr->Add(trackHit);
-    }
-    fTrackHitsPtr->SetOwner(); // nedeed for deleting TClonesArray
-  }
-
+  
+  Int_t maxIndex = 0;
+  
   // necessary to make a copy of the objects and not only the pointers in TClonesArray.
   if (theMUONTrack.fTrackParamAtHit) {
-    fTrackParamAtHit  =  new TClonesArray("AliMUONTrackParam",10);
-    for (Int_t index = 0; index < (theMUONTrack.fTrackParamAtHit)->GetEntriesFast(); index++) {
-      {new ((*fTrackParamAtHit)[fTrackParamAtHit->GetEntriesFast()]) 
-         AliMUONTrackParam(*(AliMUONTrackParam*)(theMUONTrack.fTrackParamAtHit)->At(index));}
+    maxIndex = (theMUONTrack.fTrackParamAtHit)->GetEntriesFast();
+    fTrackParamAtHit  =  new TClonesArray("AliMUONTrackParam",maxIndex);
+    for (Int_t index = 0; index < maxIndex; index++) {
+      new ((*fTrackParamAtHit)[index]) AliMUONTrackParam(*(AliMUONTrackParam*)theMUONTrack.fTrackParamAtHit->At(index));
     }
   }  
 
   // necessary to make a copy of the objects and not only the pointers in TClonesArray.
   if (theMUONTrack.fHitForRecAtHit) {
-    fHitForRecAtHit  =  new TClonesArray("AliMUONHitForRec",10);
-    for (Int_t index = 0; index < (theMUONTrack.fHitForRecAtHit)->GetEntriesFast(); index++) {
-      {new ((*fHitForRecAtHit)[fHitForRecAtHit->GetEntriesFast()]) 
-         AliMUONHitForRec(*(AliMUONHitForRec*)(theMUONTrack.fHitForRecAtHit)->At(index));}
+    maxIndex = (theMUONTrack.fHitForRecAtHit)->GetEntriesFast();
+    fHitForRecAtHit  =  new TClonesArray("AliMUONHitForRec",maxIndex);
+    for (Int_t index = 0; index < maxIndex; index++) {
+      new ((*fHitForRecAtHit)[index]) AliMUONHitForRec(*(AliMUONHitForRec*)theMUONTrack.fHitForRecAtHit->At(index));
     }
   }  
 
@@ -222,48 +171,34 @@ AliMUONTrack & AliMUONTrack::operator=(const AliMUONTrack& theMUONTrack)
   // base class assignement
   TObject::operator=(theMUONTrack);
 
-  // fTrackReconstructor =  new AliMUONTrackReconstructor(*MUONTrack.fTrackReconstructor); 
-  // is it right ?
-  // is it right ? NO because it would use dummy copy constructor
-  fTrackReconstructor =  theMUONTrack.fTrackReconstructor;
-  fTrackParamAtVertex =  theMUONTrack.fTrackParamAtVertex;
-
- // necessary to make a copy of the objects and not only the pointers in TObjArray.
-  fTrackHitsPtr = 0;
-  if (theMUONTrack.fTrackHitsPtr) {
-    fTrackHitsPtr  =  new TObjArray(10);
-    for (Int_t index = 0; index < (theMUONTrack.fTrackHitsPtr)->GetEntriesFast(); index++) {
-      AliMUONTrackHit *trackHit = new AliMUONTrackHit(*(AliMUONTrackHit*)(theMUONTrack.fTrackHitsPtr)->At(index));
-      fTrackHitsPtr->Add(trackHit);
-    }
-    fTrackHitsPtr->SetOwner();  // nedeed for deleting TClonesArray
-  }  
+  fTrackParamAtVertex = theMUONTrack.fTrackParamAtVertex;
 
+  Int_t maxIndex = 0;
+  
   // necessary to make a copy of the objects and not only the pointers in TClonesArray.
   fTrackParamAtHit = 0;
   if (theMUONTrack.fTrackParamAtHit) {
     fTrackParamAtHit  =  new TClonesArray("AliMUONTrackParam",10);
-    for (Int_t index = 0; index < (theMUONTrack.fTrackParamAtHit)->GetEntriesFast(); index++) {
-      {new ((*fTrackParamAtHit)[fTrackParamAtHit->GetEntriesFast()]) 
-         AliMUONTrackParam(*(AliMUONTrackParam*)(theMUONTrack.fTrackParamAtHit)->At(index));}
+    maxIndex = (theMUONTrack.fTrackParamAtHit)->GetEntriesFast();
+    for (Int_t index = 0; index < maxIndex; index++) {
+      new ((*fTrackParamAtHit)[fTrackParamAtHit->GetEntriesFast()])
+       AliMUONTrackParam(*(AliMUONTrackParam*)(theMUONTrack.fTrackParamAtHit)->At(index));
     }
   }  
 
   // necessary to make a copy of the objects and not only the pointers in TClonesArray.
   fHitForRecAtHit = 0;
-  if (theMUONTrack.fHitForRecAtHit) {  
+  if (theMUONTrack.fHitForRecAtHit) {
     fHitForRecAtHit  =  new TClonesArray("AliMUONHitForRec",10);
-    for (Int_t index = 0; index < (theMUONTrack.fHitForRecAtHit)->GetEntriesFast(); index++) {
-      {new ((*fHitForRecAtHit)[fHitForRecAtHit->GetEntriesFast()]) 
-         AliMUONHitForRec(*(AliMUONHitForRec*)(theMUONTrack.fHitForRecAtHit)->At(index));}
+    maxIndex = (theMUONTrack.fHitForRecAtHit)->GetEntriesFast();
+    for (Int_t index = 0; index < maxIndex; index++) {
+      new ((*fHitForRecAtHit)[fHitForRecAtHit->GetEntriesFast()])
+       AliMUONHitForRec(*(AliMUONHitForRec*)(theMUONTrack.fHitForRecAtHit)->At(index));
     }
   }  
-
+  
   fNTrackHits         =  theMUONTrack.fNTrackHits;
-  fFitMCS             =  theMUONTrack.fFitMCS;
-  fFitNParam          =  theMUONTrack.fFitNParam;
   fFitFMin            =  theMUONTrack.fFitFMin;
-  fFitStart           =  theMUONTrack.fFitStart;
   fMatchTrigger       =  theMUONTrack.fMatchTrigger;
   fChi2MatchTrigger   =  theMUONTrack.fChi2MatchTrigger;
   fTrackID            =  theMUONTrack.fTrackID;
@@ -272,110 +207,34 @@ AliMUONTrack & AliMUONTrack::operator=(const AliMUONTrack& theMUONTrack)
 }
 
   //__________________________________________________________________________
-void AliMUONTrack::Remove()
+void AliMUONTrack::AddTrackParamAtHit(AliMUONTrackParam *trackParam, AliMUONHitForRec *hitForRec) 
 {
-  // Remove current track from array of tracks,
-  // and corresponding track hits from array of track hits.
-  // Compress the TClonesArray it belongs to.
-  AliMUONTrackHit *nextTrackHit;
-  AliMUONTrackReconstructor *eventRec = this->fTrackReconstructor;
-  TClonesArray *trackHitsPtr = eventRec->GetRecTrackHitsPtr();
-  // Loop over all track hits of track
-  AliMUONTrackHit *trackHit = (AliMUONTrackHit*) fTrackHitsPtr->First();
-  while (trackHit) {
-    nextTrackHit = (AliMUONTrackHit*) fTrackHitsPtr->After(trackHit);
-    // Remove TrackHit from event TClonesArray.
-    // Destructor is called,
-    // hence links between HitForRec's and TrackHit's are updated
-    trackHitsPtr->Remove(trackHit);
-    trackHit = nextTrackHit;
-  }
-  // Remove the track from event TClonesArray
-  // Destructor is called,
-  // hence space for TObjArray of pointers to TrackHit's is freed
-  eventRec->GetRecTracksPtr()->Remove(this);
-  // Number of tracks decreased by 1
-  eventRec->SetNRecTracks(eventRec->GetNRecTracks() - 1);
-  // Compress event TClonesArray of Track's:
-  // this is essential to retrieve the TClonesArray afterwards
-  eventRec->GetRecTracksPtr()->Compress();
-  // Compress event TClonesArray of TrackHit's:
-  // this is probably also essential to retrieve the TClonesArray afterwards
-  trackHitsPtr->Compress();
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::SetFitMCS(Int_t FitMCS)
-{
-  // Set multiple Coulomb scattering option for track fit "fFitMCS"
-  // from "FitMCS" argument: 0 without, 1 with
-  if ((FitMCS == 0) || (FitMCS == 1)) fFitMCS = FitMCS;
-  // better implementation with enum(with, without) ????
-  else {
-    cout << "ERROR in AliMUONTrack::SetFitMCS(FitMCS)" << endl;
-    cout << "FitMCS = " << FitMCS << " is neither 0 nor 1" << endl;
-    exit(0);
-  }
-  return;
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::SetFitNParam(Int_t FitNParam)
-{
-  // Set number of parameters for track fit "fFitNParam" from "FitNParam":
-  // 3 for momentum, 5 for momentum and position
-  if ((FitNParam == 3) || (FitNParam == 5)) fFitNParam = FitNParam;
-  else {
-    cout << "ERROR in AliMUONTrack::SetFitNParam(FitNParam)" << endl;
-    cout << "FitNParam = " << FitNParam << " is neither 3 nor 5" << endl;
-    exit(0);
+  // Add TrackParamAtHit if "trackParam" != NULL else create empty TrackParamAtHit
+  // Update link to HitForRec if "hitForRec" != NULL
+  if (!fTrackParamAtHit) {
+    fTrackParamAtHit = new TClonesArray("AliMUONTrackParam",10);  
+    fNTrackHits = 0;
   }
-  return;
+  AliMUONTrackParam* TrackParamAtHit;
+  if (trackParam) TrackParamAtHit = new ((*fTrackParamAtHit)[fNTrackHits]) AliMUONTrackParam(*trackParam);
+  else TrackParamAtHit = new ((*fTrackParamAtHit)[fNTrackHits]) AliMUONTrackParam();
+  if (hitForRec) TrackParamAtHit->SetHitForRecPtr(hitForRec);
+  fNTrackHits++;
 }
 
   //__________________________________________________________________________
-void AliMUONTrack::SetFitStart(Int_t FitStart)
+void AliMUONTrack::AddHitForRecAtHit(const AliMUONHitForRec *hitForRec) 
 {
-  // Set multiple Coulomb scattering option for track fit "fFitStart"
-  // from "FitStart" argument: 0 without, 1 with
-  if ((FitStart == 0) || (FitStart == 1)) fFitStart = FitStart;
-  // better implementation with enum(vertex, firstHit) ????
-  else {
-    cout << "ERROR in AliMUONTrack::SetFitStart(FitStart)" << endl;
-    cout << "FitStart = " << FitStart << " is neither 0 nor 1" << endl;
-    exit(0);
-  }
-  return;
+  if (!fHitForRecAtHit)
+    fHitForRecAtHit = new TClonesArray("AliMUONHitForRec",10); 
+  
+  if (!hitForRec)
+    AliFatal("AliMUONTrack::AddHitForRecAtHit: hitForRec == NULL");
+  
+  new ((*fHitForRecAtHit)[fHitForRecAtHit->GetEntriesFast()]) AliMUONHitForRec(*hitForRec);
 }
 
   //__________________________________________________________________________
-AliMUONTrackParam* AliMUONTrack::GetTrackParamAtFirstHit(void) const {
-  // Get pointer to TrackParamAtFirstHit
-  return ((AliMUONTrackHit*) (fTrackHitsPtr->First()))->GetTrackParam();}
-
-  //__________________________________________________________________________
-void AliMUONTrack::RecursiveDump(void) const
-{
-  // Recursive dump of AliMUONTrack, i.e. with dump of TrackHit's and HitForRec's
-  AliMUONTrackHit *trackHit;
-  AliMUONHitForRec *hitForRec;
-  cout << "Recursive dump of Track: " << this << endl;
-  // Track
-  this->Dump();
-  for (Int_t trackHitIndex = 0; trackHitIndex < fNTrackHits; trackHitIndex++) {
-    trackHit = (AliMUONTrackHit*) ((*fTrackHitsPtr)[trackHitIndex]);
-    // TrackHit
-    cout << "TrackHit: " << trackHit << " (index: " << trackHitIndex << ")" << endl;
-    trackHit->Dump();
-    hitForRec = trackHit->GetHitForRecPtr();
-    // HitForRec
-    cout << "HitForRec: " << hitForRec << endl;
-    hitForRec->Dump();
-  }
-  return;
-}
-  
-  //__________________________________________________________________________
 Bool_t* AliMUONTrack::CompatibleTrack(AliMUONTrack * Track, Double_t Sigma2Cut) const
 {
   // Return kTRUE/kFALSE for each chamber if hit is compatible or not 
@@ -421,585 +280,44 @@ Int_t AliMUONTrack::HitsInCommon(AliMUONTrack* Track) const
   // between the current track ("this")
   // and the track pointed to by "Track".
   Int_t hitsInCommon = 0;
-  AliMUONTrackHit *trackHit1, *trackHit2;
+  AliMUONTrackParam *trackParamAtHit1, *trackParamAtHit2;
   // Loop over hits of first track
-  trackHit1 = (AliMUONTrackHit*) this->GetTrackHitsPtr()->First();
-  while (trackHit1) {
+  trackParamAtHit1 = (AliMUONTrackParam*) this->fTrackParamAtHit->First();
+  while (trackParamAtHit1) {
     // Loop over hits of second track
-    trackHit2 = (AliMUONTrackHit*) Track->GetTrackHitsPtr()->First();
-    while (trackHit2) {
-      // Increment "hitsInCommon" if both TrackHits point to the same HitForRec
-      if ( (trackHit1->GetHitForRecPtr()) ==
-          (trackHit2->GetHitForRecPtr())    ) hitsInCommon++;
-      trackHit2 = (AliMUONTrackHit*) Track->GetTrackHitsPtr()->After(trackHit2);
-    } // trackHit2
-    trackHit1 = (AliMUONTrackHit*) this->GetTrackHitsPtr()->After(trackHit1);
-  } // trackHit1
+    trackParamAtHit2 = (AliMUONTrackParam*) Track->fTrackParamAtHit->First();
+    while (trackParamAtHit2) {
+      // Increment "hitsInCommon" if both TrackParamAtHits point to the same HitForRec
+      if ((trackParamAtHit1->GetHitForRecPtr()) == (trackParamAtHit2->GetHitForRecPtr())) hitsInCommon++;
+      trackParamAtHit2 = (AliMUONTrackParam*) Track->fTrackParamAtHit->After(trackParamAtHit2);
+    } // trackParamAtHit2
+    trackParamAtHit1 = (AliMUONTrackParam*) this->fTrackParamAtHit->After(trackParamAtHit1);
+  } // trackParamAtHit1
   return hitsInCommon;
 }
 
   //__________________________________________________________________________
-void AliMUONTrack::MatchTriggerTrack(TClonesArray *triggerTrackArray)
-{
-  // Match this track with one trigger track if possible
-  AliMUONTrackParam trackParam; 
-  AliMUONTriggerTrack *triggerTrack;
-  Double_t xTrack, yTrack, ySlopeTrack, dTrigTrackMin2, dTrigTrack2;
-  Double_t nSigmaCut2;
-
-  Double_t distSigma[3]={1,1,0.02}; // sigma of distributions (trigger-track) X,Y,slopeY
-  Double_t distTriggerTrack[3] = {0,0,0};
-
-  fMatchTrigger = kFALSE;
-  fChi2MatchTrigger = 0;
-
-  trackParam = *((AliMUONTrackParam*) fTrackParamAtHit->Last()); 
-  trackParam.ExtrapToZ(AliMUONConstants::DefaultChamberZ(10)); // extrap to 1st trigger chamber
-
-  nSigmaCut2 =  fTrackReconstructor->GetMaxSigma2Distance(); // nb of sigma**2 for cut
-  xTrack = trackParam.GetNonBendingCoor();
-  yTrack = trackParam.GetBendingCoor();
-  ySlopeTrack = trackParam.GetBendingSlope();
-  dTrigTrackMin2 = 999;
-  
-  triggerTrack = (AliMUONTriggerTrack*) triggerTrackArray->First();
-  while(triggerTrack){
-    distTriggerTrack[0] = (triggerTrack->GetX11()-xTrack)/distSigma[0];
-    distTriggerTrack[1] = (triggerTrack->GetY11()-yTrack)/distSigma[1];
-    distTriggerTrack[2] = (TMath::Tan(triggerTrack->GetThetay())-ySlopeTrack)/distSigma[2];
-    dTrigTrack2 = 0;
-    for (Int_t iVar = 0; iVar < 3; iVar++)
-      dTrigTrack2 += distTriggerTrack[iVar]*distTriggerTrack[iVar];
-    if (dTrigTrack2 < dTrigTrackMin2 && dTrigTrack2 < nSigmaCut2) {
-      dTrigTrackMin2 = dTrigTrack2;
-      fMatchTrigger = kTRUE;
-      fChi2MatchTrigger =  dTrigTrack2/3.; // Normalized Chi2, 3 variables (X,Y,slopeY)
-    }
-    triggerTrack = (AliMUONTriggerTrack*) triggerTrackArray->After(triggerTrack);
-  }
-
-}
-  //__________________________________________________________________________
-void AliMUONTrack::Fit()
-{
-  // Fit the current track ("this"),
-  // with or without multiple Coulomb scattering according to "fFitMCS",
-  // with the number of parameters given by "fFitNParam"
-  // (3 if one keeps X and Y fixed in "TrackParam", 5 if one lets them vary),
-  // starting, according to "fFitStart",
-  // with track parameters at vertex or at the first TrackHit.
-  // "fFitMCS", "fFitNParam" and "fFitStart" have to be set before
-  // by calling the corresponding Set methods.
-  Double_t arg[1], benC, errorParam, invBenP, lower, nonBenC, upper, x, y;
-  char parName[50];
-  AliMUONTrackParam *trackParam;
-  // Check if Minuit is initialized...
-  fgFitter = TVirtualFitter::Fitter(this); // add 3 or 5 for the maximum number of parameters ???
-  fgFitter->Clear(); // necessary ???? probably yes
-  // how to reset the printout number at every fit ????
-  // is there any risk to leave it like that ????
-  // how to go faster ???? choice of Minuit parameters like EDM ????
-  // choice of function to be minimized according to fFitMCS
-  if (fFitMCS == 0) fgFitter->SetFCN(TrackChi2);
-  else fgFitter->SetFCN(TrackChi2MCS);
-  // Switch off printout
-  arg[0] = -1;
-  fgFitter->ExecuteCommand("SET PRINT", arg, 1); // More printing !!!!
-  // No warnings
-  fgFitter->ExecuteCommand("SET NOW", arg, 0);
-  // Parameters according to "fFitStart"
-  // (should be a function to be used at every place where needed ????)
-  if (fFitStart == 0) trackParam = &fTrackParamAtVertex;
-  else trackParam = this->GetTrackParamAtFirstHit();
-  // set first 3 Minuit parameters
-  // could be tried with no limits for the search (min=max=0) ????
-  fgFitter->SetParameter(0, "InvBenP",
-                        trackParam->GetInverseBendingMomentum(),
-                        0.003, -0.4, 0.4);
-  fgFitter->SetParameter(1, "BenS",
-                        trackParam->GetBendingSlope(),
-                        0.001, -0.5, 0.5);
-  fgFitter->SetParameter(2, "NonBenS",
-                        trackParam->GetNonBendingSlope(),
-                        0.001, -0.5, 0.5);
-  if (fFitNParam == 5) {
-    // set last 2 Minuit parameters
-    // mandatory limits in Bending to avoid NaN values of parameters
-    fgFitter->SetParameter(3, "X",
-                          trackParam->GetNonBendingCoor(),
-                          0.03, -500.0, 500.0);
-    // mandatory limits in non Bending to avoid NaN values of parameters
-    fgFitter->SetParameter(4, "Y",
-                          trackParam->GetBendingCoor(),
-                          0.10, -500.0, 500.0);
-  }
-  // search without gradient calculation in the function
-  fgFitter->ExecuteCommand("SET NOGRADIENT", arg, 0);
-  // minimization
-  fgFitter->ExecuteCommand("MINIMIZE", arg, 0);
-  // exit from Minuit
-  //  fgFitter->ExecuteCommand("EXIT", arg, 0); // necessary ????
-  // get results into "invBenP", "benC", "nonBenC" ("x", "y")
-  fgFitter->GetParameter(0, parName, invBenP, errorParam, lower, upper);
-  fgFitter->GetParameter(1, parName, benC, errorParam, lower, upper);
-  fgFitter->GetParameter(2, parName, nonBenC, errorParam, lower, upper);
-  if (fFitNParam == 5) {
-    fgFitter->GetParameter(3, parName, x, errorParam, lower, upper);
-    fgFitter->GetParameter(4, parName, y, errorParam, lower, upper);
-  }
-  // result of the fit into track parameters
-  trackParam->SetInverseBendingMomentum(invBenP);
-  trackParam->SetBendingSlope(benC);
-  trackParam->SetNonBendingSlope(nonBenC);
-  if (fFitNParam == 5) {
-    trackParam->SetNonBendingCoor(x);
-    trackParam->SetBendingCoor(y);
-  }
-  // global result of the fit
-  Double_t fedm, errdef;
-  Int_t npari, nparx;
-  fgFitter->GetStats(fFitFMin, fedm, errdef, npari, nparx);
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::AddSegment(AliMUONSegment* Segment)
-{
-  // Add Segment to the track
-  AddHitForRec(Segment->GetHitForRec1()); // 1st hit
-  AddHitForRec(Segment->GetHitForRec2()); // 2nd hit
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::AddHitForRec(AliMUONHitForRec* HitForRec)
-{
-  // Add HitForRec to the track:
-  // actual TrackHit into TClonesArray of TrackHit's for the event;
-  // pointer to actual TrackHit in TObjArray of pointers to TrackHit's for the track
-  TClonesArray *recTrackHitsPtr = this->fTrackReconstructor->GetRecTrackHitsPtr();
-  Int_t eventTrackHits = this->fTrackReconstructor->GetNRecTrackHits();
-  // event
-  AliMUONTrackHit* trackHit =
-    new ((*recTrackHitsPtr)[eventTrackHits]) AliMUONTrackHit(HitForRec);
-  this->fTrackReconstructor->SetNRecTrackHits(eventTrackHits + 1);
-  // track
-  if (fTrackHitsPtr->IsOwner()) AliFatal("fTrackHitsPtr is owner");
-  fTrackHitsPtr->Add(trackHit);
-  fNTrackHits++;
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::SetTrackParamAtHit(Int_t indexHit, AliMUONTrackParam *TrackParam) const
-{
-  // Set track parameters at TrackHit with index "indexHit"
-  // from the track parameters pointed to by "TrackParam".
-  //PH  AliMUONTrackHit* trackHit = (AliMUONTrackHit*) ((*fTrackHitsPtr)[indexHit]);
-  AliMUONTrackHit* trackHit = (AliMUONTrackHit*) (fTrackHitsPtr->At(indexHit));
-  trackHit->SetTrackParam(TrackParam);
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::SetTrackParamAtVertex()
-{
-  // Set track parameters at vertex.
-  // TrackHit's are assumed to be only in stations(1..) 4 and 5,
-  // and sorted according to increasing Z..
-  // Parameters are calculated from information in HitForRec's
-  // of first and last TrackHit's.
-  AliMUONTrackParam *trackParam =
-    &fTrackParamAtVertex; // pointer to track parameters
-  // Pointer to HitForRec of first TrackHit
-  AliMUONHitForRec *firstHit =
-    ((AliMUONTrackHit*) (fTrackHitsPtr->First()))->GetHitForRecPtr();
-  // Pointer to HitForRec of last TrackHit
-  AliMUONHitForRec *lastHit =
-    ((AliMUONTrackHit*) (fTrackHitsPtr->Last()))->GetHitForRecPtr();
-  // Z difference between first and last hits
-  Double_t deltaZ = firstHit->GetZ() - lastHit->GetZ();
-  // bending slope in stations(1..) 4 and 5
-  Double_t bendingSlope =
-    (firstHit->GetBendingCoor() - lastHit->GetBendingCoor()) / deltaZ;
-  trackParam->SetBendingSlope(bendingSlope);
-  // impact parameter
-  Double_t impactParam =
-    firstHit->GetBendingCoor() - bendingSlope * firstHit->GetZ(); // same if from firstHit and  lastHit ????
-  // signed bending momentum
-  Double_t signedBendingMomentum =
-    fTrackReconstructor->GetBendingMomentumFromImpactParam(impactParam);
-  trackParam->SetInverseBendingMomentum(1.0 / signedBendingMomentum);
-  // bending slope at vertex
-  trackParam->
-    SetBendingSlope(bendingSlope +
-                   impactParam / fTrackReconstructor->GetSimpleBPosition());
-  // non bending slope
-  Double_t nonBendingSlope =
-    (firstHit->GetNonBendingCoor() - lastHit->GetNonBendingCoor()) / deltaZ;
-  trackParam->SetNonBendingSlope(nonBendingSlope);
-  // vertex coordinates at (0,0,0)
-  trackParam->SetZ(0.0);
-  trackParam->SetBendingCoor(0.0);
-  trackParam->SetNonBendingCoor(0.0);
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::AddTrackParamAtHit(const AliMUONTrackParam *trackParam) 
-{
-  // Add track paremeters
-
-  if (!fTrackParamAtHit)
-    fTrackParamAtHit = new TClonesArray("AliMUONTrackParam",10);  
-
-  new ((*fTrackParamAtHit)[fTrackParamAtHit->GetEntriesFast()]) AliMUONTrackParam(*trackParam);
-}
-
-  //__________________________________________________________________________
-void AliMUONTrack::AddHitForRecAtHit(const AliMUONHitForRec *hitForRec) 
-{
-  if (!fHitForRecAtHit)
-    fHitForRecAtHit = new TClonesArray("AliMUONHitForRec",10); 
-
-  new ((*fHitForRecAtHit)[fHitForRecAtHit->GetEntriesFast()]) AliMUONHitForRec(*hitForRec);
-}
-
-  //__________________________________________________________________________
-void TrackChi2(Int_t &NParam, Double_t * /*Gradient*/, Double_t &Chi2, Double_t *Param, Int_t /*Flag*/)
+void AliMUONTrack::RecursiveDump(void) const
 {
-  // Return the "Chi2" to be minimized with Minuit for track fitting,
-  // with "NParam" parameters
-  // and their current values in array pointed to by "Param".
-  // Assumes that the track hits are sorted according to increasing Z.
-  // Track parameters at each TrackHit are updated accordingly.
-  // Multiple Coulomb scattering is not taken into account
-  AliMUONTrack *trackBeingFitted;
-  AliMUONTrackHit* hit;
-  AliMUONTrackParam param1;
-  Int_t hitNumber;
-  Double_t zHit;
-  Chi2 = 0.0; // initialize Chi2
-  // copy of track parameters to be fitted
-  trackBeingFitted = (AliMUONTrack*) AliMUONTrack::Fitter()->GetObjectFit();
-  if (trackBeingFitted->GetFitStart() == 0)
-    param1 = *(trackBeingFitted->GetTrackParamAtVertex());
-  else param1 = *(trackBeingFitted->GetTrackParamAtFirstHit());
-  // Minuit parameters to be fitted into this copy
-  param1.SetInverseBendingMomentum(Param[0]);
-  param1.SetBendingSlope(Param[1]);
-  param1.SetNonBendingSlope(Param[2]);
-  if (NParam == 5) {
-    param1.SetNonBendingCoor(Param[3]);
-    param1.SetBendingCoor(Param[4]);
-  }
-  // Follow track through all planes of track hits
-  for (hitNumber = 0; hitNumber < trackBeingFitted->GetNTrackHits(); hitNumber++) {
-    hit = (AliMUONTrackHit*) (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber];
-    zHit = hit->GetHitForRecPtr()->GetZ();
-    // do something special if 2 hits with same Z ????
-    // security against infinite loop ????
-    (&param1)->ExtrapToZ(zHit); // extrapolation
-    hit->SetTrackParam(&param1);
-    // Increment Chi2
-    // done hit per hit, with hit resolution,
-    // and not with point and angle like in "reco_muon.F" !!!!
-    // Needs to add multiple scattering contribution ????
-    Double_t dX =
-      hit->GetHitForRecPtr()->GetNonBendingCoor() - (&param1)->GetNonBendingCoor();
-    Double_t dY =
-      hit->GetHitForRecPtr()->GetBendingCoor() - (&param1)->GetBendingCoor();
-    Chi2 =
-      Chi2 +
-      dX * dX / hit->GetHitForRecPtr()->GetNonBendingReso2() +
-      dY * dY / hit->GetHitForRecPtr()->GetBendingReso2();
+  // Recursive dump of AliMUONTrack, i.e. with dump of TrackParamAtHit's and attached HitForRec's
+  AliMUONTrackParam *trackParamAtHit;
+  AliMUONHitForRec *hitForRec;
+  cout << "Recursive dump of Track: " << this << endl;
+  // Track
+  this->Dump();
+  for (Int_t trackHitIndex = 0; trackHitIndex < fNTrackHits; trackHitIndex++) {
+    trackParamAtHit = (AliMUONTrackParam*) ((*fTrackParamAtHit)[trackHitIndex]);
+    // TrackHit
+    cout << "TrackParamAtHit: " << trackParamAtHit << " (index: " << trackHitIndex << ")" << endl;
+    trackParamAtHit->Dump();
+    hitForRec = trackParamAtHit->GetHitForRecPtr();
+    // HitForRec
+    cout << "HitForRec: " << hitForRec << endl;
+    hitForRec->Dump();
   }
+  return;
 }
-
-  //__________________________________________________________________________
-void TrackChi2MCS(Int_t &NParam, Double_t * /*Gradient*/, Double_t &Chi2, Double_t *Param, Int_t /*Flag*/)
-{
-  // Return the "Chi2" to be minimized with Minuit for track fitting,
-  // with "NParam" parameters
-  // and their current values in array pointed to by "Param".
-  // Assumes that the track hits are sorted according to increasing Z.
-  // Track parameters at each TrackHit are updated accordingly.
-  // Multiple Coulomb scattering is taken into account with covariance matrix.
-  AliMUONTrack *trackBeingFitted;
-  AliMUONTrackParam param1;
-  Chi2 = 0.0; // initialize Chi2
-  // copy of track parameters to be fitted
-  trackBeingFitted = (AliMUONTrack*) AliMUONTrack::Fitter()->GetObjectFit();
-  if (trackBeingFitted->GetFitStart() == 0)
-    param1 = *(trackBeingFitted->GetTrackParamAtVertex());
-  else param1 = *(trackBeingFitted->GetTrackParamAtFirstHit());
-  // Minuit parameters to be fitted into this copy
-  param1.SetInverseBendingMomentum(Param[0]);
-  param1.SetBendingSlope(Param[1]);
-  param1.SetNonBendingSlope(Param[2]);
-  if (NParam == 5) {
-    param1.SetNonBendingCoor(Param[3]);
-    param1.SetBendingCoor(Param[4]);
-  }
-
-  AliMUONTrackHit *hit;
-  Int_t chCurrent, chPrev = 0, hitNumber, hitNumber1, hitNumber2, hitNumber3;
-  Double_t z, z1, z2, z3;
-  AliMUONTrackHit *hit1, *hit2, *hit3;
-  Double_t hbc1, hbc2, pbc1, pbc2;
-  Double_t hnbc1, hnbc2, pnbc1, pnbc2;
-  Int_t numberOfHit = trackBeingFitted->GetNTrackHits();
-  TMatrixD *covBending = new TMatrixD(numberOfHit, numberOfHit);
-  TMatrixD *covNonBending = new TMatrixD(numberOfHit, numberOfHit);
-  Double_t *msa2 = new Double_t[numberOfHit];
-
-  // Predicted coordinates and  multiple scattering angles are first calculated
-  for (hitNumber = 0; hitNumber < numberOfHit; hitNumber++) {
-    hit = (AliMUONTrackHit*) (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber];
-    z = hit->GetHitForRecPtr()->GetZ();
-    // do something special if 2 hits with same Z ????
-    // security against infinite loop ????
-    (&param1)->ExtrapToZ(z); // extrapolation
-    hit->SetTrackParam(&param1);
-    // square of multiple scattering angle at current hit, with one chamber
-    msa2[hitNumber] = MultipleScatteringAngle2(hit);
-    // correction for eventual missing hits or multiple hits in a chamber,
-    // according to the number of chambers
-    // between the current hit and the previous one
-    chCurrent = hit->GetHitForRecPtr()->GetChamberNumber();
-    if (hitNumber > 0) msa2[hitNumber] = msa2[hitNumber] * (chCurrent - chPrev);
-    chPrev = chCurrent;
-  }
-
-  // Calculates the covariance matrix
-  for (hitNumber1 = 0; hitNumber1 < numberOfHit; hitNumber1++) { 
-    hit1 = (AliMUONTrackHit*) (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber1];
-    z1 = hit1->GetHitForRecPtr()->GetZ();
-    for (hitNumber2 = hitNumber1; hitNumber2 < numberOfHit; hitNumber2++) {
-      hit2 = (AliMUONTrackHit*) (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber2];
-      z2 = hit2->GetHitForRecPtr()->GetZ();
-      // initialization to 0 (diagonal plus upper triangular part)
-      (*covBending)(hitNumber2, hitNumber1) = 0.0;
-      // contribution from multiple scattering in bending plane:
-      // loop over upstream hits
-      for (hitNumber3 = 0; hitNumber3 < hitNumber1; hitNumber3++) {    
-       hit3 = (AliMUONTrackHit*)
-         (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber3];
-       z3 = hit3->GetHitForRecPtr()->GetZ();
-       (*covBending)(hitNumber2, hitNumber1) =
-         (*covBending)(hitNumber2, hitNumber1) +
-         ((z1 - z3) * (z2 - z3) * msa2[hitNumber3]); 
-      }
-      // equal contribution from multiple scattering in non bending plane
-      (*covNonBending)(hitNumber2, hitNumber1) =
-       (*covBending)(hitNumber2, hitNumber1);
-      if (hitNumber1 == hitNumber2) {
-       // Diagonal elements: add contribution from position measurements
-       // in bending plane
-       (*covBending)(hitNumber2, hitNumber1) =
-         (*covBending)(hitNumber2, hitNumber1) +
-         hit1->GetHitForRecPtr()->GetBendingReso2();
-       // and in non bending plane
-       (*covNonBending)(hitNumber2, hitNumber1) =
-         (*covNonBending)(hitNumber2, hitNumber1) +
-         hit1->GetHitForRecPtr()->GetNonBendingReso2();
-      }
-      else {
-       // Non diagonal elements: symmetrization
-       // for bending plane
-       (*covBending)(hitNumber1, hitNumber2) =
-         (*covBending)(hitNumber2, hitNumber1);
-       // and non bending plane
-       (*covNonBending)(hitNumber1, hitNumber2) =
-         (*covNonBending)(hitNumber2, hitNumber1);
-      }
-    } // for (hitNumber2 = hitNumber1;...
-  } // for (hitNumber1 = 0;...
-    
-  // Inversion of covariance matrices
-  // with "mnvertLocal", local "mnvert" function of Minuit.
-  // One cannot use directly "mnvert" since "TVirtualFitter" does not know it.
-  // One will have to replace this local function by the right inversion function
-  // from a specialized Root package for symmetric positive definite matrices,
-  // when available!!!!
-  Int_t ifailBending;
-  mnvertLocal(&((*covBending)(0,0)), numberOfHit, numberOfHit, numberOfHit,
-             ifailBending);
-  Int_t ifailNonBending;
-  mnvertLocal(&((*covNonBending)(0,0)), numberOfHit, numberOfHit, numberOfHit,
-             ifailNonBending);
-
-  // It would be worth trying to calculate the inverse of the covariance matrix
-  // only once per fit, since it cannot change much in principle,
-  // and it would save a lot of computing time !!!!
   
-  // Calculates Chi2
-  if ((ifailBending == 0) && (ifailNonBending == 0)) {
-    // with Multiple Scattering if inversion correct
-    for (hitNumber1 = 0; hitNumber1 < numberOfHit ; hitNumber1++) { 
-      hit1 = (AliMUONTrackHit*) (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber1];
-      hbc1 = hit1->GetHitForRecPtr()->GetBendingCoor();
-      pbc1 = hit1->GetTrackParam()->GetBendingCoor();
-      hnbc1 = hit1->GetHitForRecPtr()->GetNonBendingCoor();
-      pnbc1 = hit1->GetTrackParam()->GetNonBendingCoor();
-      for (hitNumber2 = 0; hitNumber2 < numberOfHit; hitNumber2++) {
-       hit2 = (AliMUONTrackHit*)
-         (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber2];
-       hbc2 = hit2->GetHitForRecPtr()->GetBendingCoor();
-       pbc2 = hit2->GetTrackParam()->GetBendingCoor();
-       hnbc2 = hit2->GetHitForRecPtr()->GetNonBendingCoor();
-       pnbc2 = hit2->GetTrackParam()->GetNonBendingCoor();
-       Chi2 = Chi2 +
-         ((*covBending)(hitNumber2, hitNumber1) *
-          (hbc1 - pbc1) * (hbc2 - pbc2)) +
-         ((*covNonBending)(hitNumber2, hitNumber1) *
-          (hnbc1 - pnbc1) * (hnbc2 - pnbc2));
-      }
-    }
-  } else {
-    // without Multiple Scattering if inversion impossible
-    for (hitNumber1 = 0; hitNumber1 < numberOfHit ; hitNumber1++) { 
-      hit1 = (AliMUONTrackHit*) (*(trackBeingFitted->GetTrackHitsPtr()))[hitNumber1];
-      hbc1 = hit1->GetHitForRecPtr()->GetBendingCoor();
-      pbc1 = hit1->GetTrackParam()->GetBendingCoor();
-      hnbc1 = hit1->GetHitForRecPtr()->GetNonBendingCoor();
-      pnbc1 = hit1->GetTrackParam()->GetNonBendingCoor();
-      Chi2 = Chi2 + 
-       ((hbc1 - pbc1) * (hbc1 - pbc1) /
-        hit1->GetHitForRecPtr()->GetBendingReso2()) +
-       ((hnbc1 - pnbc1) * (hnbc1 - pnbc1) /
-        hit1->GetHitForRecPtr()->GetNonBendingReso2());
-    }
-  }
-  
-  delete covBending;
-  delete covNonBending;
-  delete [] msa2;
-}
-
-Double_t MultipleScatteringAngle2(AliMUONTrackHit *TrackHit)
-{
-  // Returns square of multiple Coulomb scattering angle
-  // at TrackHit pointed to by "TrackHit"
-  Double_t slopeBending, slopeNonBending, radiationLength, inverseBendingMomentum2, inverseTotalMomentum2;
-  Double_t varMultipleScatteringAngle;
-  AliMUONTrack *trackBeingFitted = (AliMUONTrack*) AliMUONTrack::Fitter()->GetObjectFit();
-  AliMUONTrackParam *param = TrackHit->GetTrackParam();
-  // Better implementation in AliMUONTrack class ????
-  slopeBending = param->GetBendingSlope();
-  slopeNonBending = param->GetNonBendingSlope();
-  // thickness in radiation length for the current track,
-  // taking local angle into account
-  radiationLength =
-    trackBeingFitted->GetTrackReconstructor()->GetChamberThicknessInX0() *
-    TMath::Sqrt(1.0 +
-               slopeBending * slopeBending + slopeNonBending * slopeNonBending);
-  inverseBendingMomentum2 = 
-    param->GetInverseBendingMomentum() * param->GetInverseBendingMomentum();
-  inverseTotalMomentum2 =
-    inverseBendingMomentum2 * (1.0 + slopeBending * slopeBending) /
-    (1.0 + slopeBending *slopeBending + slopeNonBending * slopeNonBending); 
-  varMultipleScatteringAngle = 0.0136 * (1.0 + 0.038 * TMath::Log(radiationLength));
-  // The velocity is assumed to be 1 !!!!
-  varMultipleScatteringAngle = inverseTotalMomentum2 * radiationLength *
-    varMultipleScatteringAngle * varMultipleScatteringAngle;
-  return varMultipleScatteringAngle;
-}
-
-//______________________________________________________________________________
- void mnvertLocal(Double_t *a, Int_t l, Int_t, Int_t n, Int_t &ifail)
-{
-//*-*-*-*-*-*-*-*-*-*-*-*Inverts a symmetric matrix*-*-*-*-*-*-*-*-*-*-*-*-*
-//*-*                    ==========================
-//*-*        inverts a symmetric matrix.   matrix is first scaled to
-//*-*        have all ones on the diagonal (equivalent to change of units)
-//*-*        but no pivoting is done since matrix is positive-definite.
-//*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
-
-  // taken from TMinuit package of Root (l>=n)
-  // fVERTs, fVERTq and fVERTpp changed to localVERTs, localVERTq and localVERTpp
-  //  Double_t localVERTs[n], localVERTq[n], localVERTpp[n];
-  Double_t * localVERTs = new Double_t[n];
-  Double_t * localVERTq = new Double_t[n];
-  Double_t * localVERTpp = new Double_t[n];
-  // fMaxint changed to localMaxint
-  Int_t localMaxint = n;
-
-    /* System generated locals */
-    Int_t aOffset;
-
-    /* Local variables */
-    Double_t si;
-    Int_t i, j, k, kp1, km1;
-
-    /* Parameter adjustments */
-    aOffset = l + 1;
-    a -= aOffset;
-
-    /* Function Body */
-    ifail = 0;
-    if (n < 1) goto L100;
-    if (n > localMaxint) goto L100;
-//*-*-                  scale matrix by sqrt of diag elements
-    for (i = 1; i <= n; ++i) {
-        si = a[i + i*l];
-        if (si <= 0) goto L100;
-        localVERTs[i-1] = 1 / TMath::Sqrt(si);
-    }
-    for (i = 1; i <= n; ++i) {
-        for (j = 1; j <= n; ++j) {
-            a[i + j*l] = a[i + j*l]*localVERTs[i-1]*localVERTs[j-1];
-        }
-    }
-//*-*-                                       . . . start main loop . . . .
-    for (i = 1; i <= n; ++i) {
-        k = i;
-//*-*-                  preparation for elimination step1
-        if (a[k + k*l] != 0) localVERTq[k-1] = 1 / a[k + k*l];
-        else goto L100;
-        localVERTpp[k-1] = 1;
-        a[k + k*l] = 0;
-        kp1 = k + 1;
-        km1 = k - 1;
-        if (km1 < 0) goto L100;
-        else if (km1 == 0) goto L50;
-        else               goto L40;
-L40:
-        for (j = 1; j <= km1; ++j) {
-            localVERTpp[j-1] = a[j + k*l];
-            localVERTq[j-1]  = a[j + k*l]*localVERTq[k-1];
-            a[j + k*l]   = 0;
-        }
-L50:
-        if (k - n < 0) goto L51;
-        else if (k - n == 0) goto L60;
-        else                goto L100;
-L51:
-        for (j = kp1; j <= n; ++j) {
-            localVERTpp[j-1] = a[k + j*l];
-            localVERTq[j-1]  = -a[k + j*l]*localVERTq[k-1];
-            a[k + j*l]   = 0;
-        }
-//*-*-                  elimination proper
-L60:
-        for (j = 1; j <= n; ++j) {
-            for (k = j; k <= n; ++k) { a[j + k*l] += localVERTpp[j-1]*localVERTq[k-1]; }
-        }
-    }
-//*-*-                  elements of left diagonal and unscaling
-    for (j = 1; j <= n; ++j) {
-        for (k = 1; k <= j; ++k) {
-            a[k + j*l] = a[k + j*l]*localVERTs[k-1]*localVERTs[j-1];
-            a[j + k*l] = a[k + j*l];
-        }
-    }
-    delete [] localVERTs;
-    delete [] localVERTq;
-    delete [] localVERTpp;
-    return;
-//*-*-                  failure return
-L100:
-    delete [] localVERTs;
-    delete [] localVERTq;
-    delete [] localVERTpp;
-    ifail = 1;
-} /* mnvertLocal */
-
 //_____________________________________________-
 void AliMUONTrack::Print(Option_t* opt) const
 {
index 50d0966..007f7f4 100644 (file)
 #include <TClonesArray.h>
 
 #include "AliMUONTrackParam.h" // object belongs to the class
-#include "AliMUONHitForRec.h"  // object belongs to the class
 
-//const Int_t kMaxTrackingChamber=10;
-        // not used
-
-class TObjArray;
-class TVirtualFitter;
-class AliMUONTrackReconstructor;
 class AliMUONHitForRec;
 class AliMUONSegment;
 
@@ -36,74 +29,55 @@ class AliMUONTrack : public TObject
   AliMUONTrack (const AliMUONTrack& AliMUONTrack); // copy constructor
   AliMUONTrack& operator=(const AliMUONTrack& AliMUONTrack); // assignment operator
 
-  AliMUONTrack(AliMUONSegment* BegSegment, AliMUONSegment* EndSegment, AliMUONTrackReconstructor* TrackReconstructor); // Constructor from two Segment's
-  AliMUONTrack(AliMUONSegment* Segment, AliMUONHitForRec* HitForRec, AliMUONTrackReconstructor* TrackReconstructor); // Constructor from one Segment and one HitForRec
-  void Remove(void);
+  AliMUONTrack(AliMUONSegment* BegSegment, AliMUONSegment* EndSegment); // Constructor from two Segment's
+  AliMUONTrack(AliMUONSegment* Segment, AliMUONHitForRec* HitForRec); // Constructor from one Segment and one HitForRec
 
-  AliMUONTrackReconstructor* GetTrackReconstructor(void) const {return fTrackReconstructor;}
   AliMUONTrackParam*         GetTrackParamAtVertex(void) {return &fTrackParamAtVertex;}
-  void                       SetTrackParamAtVertex(void); // Set track parameters at vertex from last stations 4 & 5
   void                       SetTrackParamAtVertex(AliMUONTrackParam* TrackParam) {fTrackParamAtVertex = *TrackParam;}
+
   TClonesArray*              GetTrackParamAtHit(void) const {return fTrackParamAtHit;}
-  TClonesArray*              GetHitForRecAtHit(void) const {return fHitForRecAtHit;}
   void                       ResetTrackParamAtHit(void) { fTrackParamAtHit->Delete(); }
+  void                       AddTrackParamAtHit(AliMUONTrackParam *trackParam, AliMUONHitForRec *hitForRec); 
+  
+  TClonesArray*              GetHitForRecAtHit(void) const {return fHitForRecAtHit;}
   void                       ResetHitForRecAtHit(void) { fHitForRecAtHit->Delete(); }
-  void                       AddTrackParamAtHit(const AliMUONTrackParam *trackParam); 
   void                       AddHitForRecAtHit(const AliMUONHitForRec *hitForRec); 
 
-  TObjArray*                 GetTrackHitsPtr(void) const {return fTrackHitsPtr;}
   Int_t                      GetNTrackHits(void) const {return fNTrackHits;}
   void                       SetNTrackHits(Int_t nTrackHits) {fNTrackHits = nTrackHits;}
-  Int_t                      GetFitMCS(void) const {return fFitMCS;}
-  Int_t                      GetFitNParam(void) const {return fFitNParam;}
-  Int_t                      GetFitStart(void) const {return fFitStart;}
+
   Double_t                   GetFitFMin(void) const {return fFitFMin;}
+  void                       SetFitFMin(Double_t chi2) { fFitFMin = chi2; } // set Chi2
   Bool_t                     GetMatchTrigger(void) const {return fMatchTrigger;}
+  void                      SetMatchTrigger(Bool_t MatchTrigger) {fMatchTrigger = MatchTrigger;}
   Double_t                   GetChi2MatchTrigger(void) const {return fChi2MatchTrigger;}
-  void                       SetFitMCS(Int_t FitMCS);
-  void                       SetFitNParam(Int_t FitNParam);
-  void                       SetFitStart(Int_t FitStart);
-  void                       SetFitFMin(Double_t chi2) { fFitFMin = chi2; } // set Chi2
-
-  AliMUONTrackParam*         GetTrackParamAtFirstHit(void) const;
-
-  void                       RecursiveDump(void) const; // Recursive dump (with track hits)
-  void                       Fit(); // Fit
-  void                       AddSegment(AliMUONSegment* Segment); // Add Segment
-  void                       AddHitForRec(AliMUONHitForRec* HitForRec); // Add HitForRec
-  void                       SetTrackParamAtHit(Int_t indexHit, AliMUONTrackParam *TrackParam) const;
+  void                       SetChi2MatchTrigger(Double_t Chi2MatchTrigger) {fChi2MatchTrigger = Chi2MatchTrigger;}
+  
   Int_t                      HitsInCommon(AliMUONTrack* Track) const;
-  void                       MatchTriggerTrack(TClonesArray* TriggerTrackArray);
   Bool_t*                    CompatibleTrack(AliMUONTrack* Track, Double_t Sigma2Cut) const; // return array of compatible chamber
   
   Int_t                      GetTrackID() const {return fTrackID;}
   void                       SetTrackID(Int_t trackID) {fTrackID = trackID;}
 
-  virtual void               Print(Option_t* opt="") const;
-
-  static TVirtualFitter*     Fitter(void) {return fgFitter;}
+  void                       RecursiveDump(void) const; // Recursive dump (with track hits)
 
+  virtual void               Print(Option_t* opt="") const;
 
 
- protected:
  private:
-  static TVirtualFitter* fgFitter; //!< Pointer to track fitter
-  AliMUONTrackReconstructor* fTrackReconstructor; //!< Pointer to TrackReconstructor
   AliMUONTrackParam fTrackParamAtVertex; ///< Track parameters at vertex
   TClonesArray *fTrackParamAtHit; ///< Track parameters at hit
   TClonesArray *fHitForRecAtHit; ///< Cluster parameters at hit
-  TObjArray *fTrackHitsPtr; //!<  Pointer to array of pointers to TrackHit's
   Int_t fNTrackHits; ///< Number of TrackHit's
-  Int_t fFitMCS; ///< 0(1) for fit without(with) multiple Coulomb scattering
-  Int_t fFitNParam; ///< 3(5) for fit with 3(5) parameters
-  Int_t fFitStart; ///< 0 or 1 for fit starting from parameters at vertex (0) or at first TrackHit(1)
+  
   Double_t fFitFMin; ///< minimum value of the function minimized by the fit
   Bool_t fMatchTrigger; ///< 1 if track matches with trigger track, 0 if not
   Double_t fChi2MatchTrigger; ///< chi2 of trigger/track matching 
+  
   Int_t fTrackID; ///< track ID = track number in TrackRefs
-
-  ClassDef(AliMUONTrack, 2) // Reconstructed track in ALICE dimuon spectrometer
-    };
+  
+  
+  ClassDef(AliMUONTrack, 3) // Reconstructed track in ALICE dimuon spectrometer
+};
        
 #endif
diff --git a/MUON/AliMUONTrackHit.cxx b/MUON/AliMUONTrackHit.cxx
deleted file mode 100644 (file)
index 06e7439..0000000
+++ /dev/null
@@ -1,138 +0,0 @@
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-
-/* $Id$ */
-
-///////////////////////////////////////////////////////
-//
-// Reconstructed track hit
-// in
-// ALICE
-// dimuon
-// spectrometer
-//
-///////////////////////////////////////////////////////
-
-#include "AliMUONTrackHit.h" 
-#include "AliMUONHitForRec.h" 
-#include "AliLog.h" 
-
-ClassImp(AliMUONTrackHit) // Class implementation in ROOT context
-
-  //__________________________________________________________________________
-AliMUONTrackHit::AliMUONTrackHit()
-  : TObject(),
-    fTrackParam(),
-    fHitForRecPtr(0x0),
-    fNextTrackHitWithSameHitForRec(0x0),
-    fPrevTrackHitWithSameHitForRec(0x0)
-{
-/// Default constructor
-
-}
-  //__________________________________________________________________________
-AliMUONTrackHit::AliMUONTrackHit (const AliMUONTrackHit& theMUONTrackHit)
-  : TObject(theMUONTrackHit),
-    fTrackParam(theMUONTrackHit.fTrackParam),
-    fHitForRecPtr(theMUONTrackHit.fHitForRecPtr),
-    fNextTrackHitWithSameHitForRec(theMUONTrackHit.fNextTrackHitWithSameHitForRec),
-    fPrevTrackHitWithSameHitForRec(theMUONTrackHit.fPrevTrackHitWithSameHitForRec)
-{
-/// Copy constructor
-
-}
-  //__________________________________________________________________________
-AliMUONTrackHit & AliMUONTrackHit::operator=(const AliMUONTrackHit& theMUONTrackHit)
-{
-/// Assignment operator
-
-  // check assignement to self
-  if (this == &theMUONTrackHit)
-    return *this;
-
-  // base class assignement
-  TObject::operator=(theMUONTrackHit);
-
-  fTrackParam                    =  theMUONTrackHit.fTrackParam;
-  fHitForRecPtr                  =  theMUONTrackHit.fHitForRecPtr;
-  fNextTrackHitWithSameHitForRec = theMUONTrackHit.fNextTrackHitWithSameHitForRec;
-  fPrevTrackHitWithSameHitForRec = theMUONTrackHit.fPrevTrackHitWithSameHitForRec;
-
-  return *this;
-
-}
-  //__________________________________________________________________________
-AliMUONTrackHit::AliMUONTrackHit(AliMUONHitForRec* Hit)
-  : TObject(),
-    fTrackParam(),
-    fHitForRecPtr(Hit),
-    fNextTrackHitWithSameHitForRec(0x0),
-    fPrevTrackHitWithSameHitForRec(0x0)
-{
-/// Constructor from the HitForRec pointed to by "Hit"
-
-  // links from/to HitForRec
-  if (Hit->GetNTrackHits() == 0) {
-    fPrevTrackHitWithSameHitForRec = NULL;
-    Hit->SetFirstTrackHitPtr(this);
-  }
-  else {
-    fPrevTrackHitWithSameHitForRec = Hit->GetLastTrackHitPtr();
-    fNextTrackHitWithSameHitForRec = NULL;
-  }
-  Hit->SetLastTrackHitPtr(this);
-  fNextTrackHitWithSameHitForRec = NULL;
-  Hit->SetNTrackHits(Hit->GetNTrackHits() + 1);
-}
-  //__________________________________________________________________________
-AliMUONTrackHit::~AliMUONTrackHit()
-{
-/// Destructor
-/// Update links between HitForRec's and TrackHit's
-/// connected to the current TrackHit being removed.
-
-  AliMUONHitForRec *hit = fHitForRecPtr; // pointer to HitForRec
-  // remove current TrackHit in HitForRec links
-  if (this == hit->GetFirstTrackHitPtr())
-    hit->SetFirstTrackHitPtr(fNextTrackHitWithSameHitForRec); // if first
-  if (this == hit->GetLastTrackHitPtr())
-    hit->SetLastTrackHitPtr(fPrevTrackHitWithSameHitForRec); // if last
-  hit->SetNTrackHits(hit->GetNTrackHits() - 1); // decrement NTrackHits of hit
-  // update link to next TrackHit of previous TrackHit
-  if (fPrevTrackHitWithSameHitForRec != NULL)
-    fPrevTrackHitWithSameHitForRec->
-      SetNextTrackHitWithSameHitForRec(fNextTrackHitWithSameHitForRec);
-  // update link to previous TrackHit of next TrackHit
-  if (fNextTrackHitWithSameHitForRec)
-    fNextTrackHitWithSameHitForRec->
-      SetPrevTrackHitWithSameHitForRec(fPrevTrackHitWithSameHitForRec);
-  // to be checked thoroughly !!!!
-  // with Root counter of AliMUONTrackHit objects,
-  // with loop over all these links after the update
-}
-
-  //__________________________________________________________________________
-Int_t AliMUONTrackHit::Compare(const TObject* TrackHit) const
-{
-/// "Compare" function to sort with decreasing Z (spectro. muon Z <0).
-/// Returns 1 (0, -1) if Z of current TrackHit
-/// is smaller than (equal to, larger than) Z of TrackHit
-
-  if (fHitForRecPtr->GetZ() <
-      ((AliMUONTrackHit*)TrackHit)->fHitForRecPtr->GetZ()) return(1);
-  else if (fHitForRecPtr->GetZ() ==
-          ((AliMUONTrackHit*)TrackHit)->fHitForRecPtr->GetZ()) return( 0);
-  else return(-1);
-}
diff --git a/MUON/AliMUONTrackHit.h b/MUON/AliMUONTrackHit.h
deleted file mode 100644 (file)
index c646716..0000000
+++ /dev/null
@@ -1,55 +0,0 @@
-#ifndef ALIMUONTRACKHIT_H
-#define ALIMUONTRACKHIT_H
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice                               */
-
-/*$Id$*/
-// Revision of includes 07/05/2004
-
-/// \ingroup rec
-/// \class AliMUONTrackHit
-/// \brief Reconstructed track hit in ALICE dimuon spectrometer
-///
-////////////////////////////////////////////////////////
-/// Reconstructed track hit in ALICE dimuon spectrometer
-////////////////////////////////////////////////////////
-
-#include <TObject.h>
-#include "AliMUONTrackParam.h" // object belongs to the class
-
-class AliMUONHitForRec;
-
-class AliMUONTrackHit : public TObject 
-{
- public:
-  AliMUONTrackHit(); // Constructor
-  virtual ~AliMUONTrackHit(); // Destructor
-  AliMUONTrackHit (const AliMUONTrackHit& AliMUONTrackHit); // copy constructor
-  AliMUONTrackHit& operator=(const AliMUONTrackHit& AliMUONTrackHit); // assignment operator
-  AliMUONTrackHit(AliMUONHitForRec* Hit); // Constructor from one HitForRec
-
-  // Inline functions for Get and Set
-  AliMUONHitForRec* GetHitForRecPtr(void) const {return fHitForRecPtr;}
-  AliMUONTrackParam* GetTrackParam(void) {return &(fTrackParam);}
-  void SetTrackParam(AliMUONTrackParam* TrackParam) {fTrackParam = *TrackParam;}
-
-  // What is necessary for sorting TClonesArray's; sufficient too ????
-  Bool_t IsSortable () const {
-    // necessary for sorting TClonesArray of TrackHit's
-    return kTRUE; }
-  Int_t Compare(const TObject* TrackHit) const; // "Compare" function for sorting
-
-
- private:
-  void SetNextTrackHitWithSameHitForRec(AliMUONTrackHit *Next) {fNextTrackHitWithSameHitForRec = Next;}
-  void SetPrevTrackHitWithSameHitForRec(AliMUONTrackHit *Prev) {fPrevTrackHitWithSameHitForRec = Prev;}
-
-  AliMUONTrackParam fTrackParam; ///< Track parameters
-  AliMUONHitForRec *fHitForRecPtr; ///< Pointer to HitForRec
-  AliMUONTrackHit *fNextTrackHitWithSameHitForRec; ///< Pointer to next track hit with same HitForRec
-  AliMUONTrackHit *fPrevTrackHitWithSameHitForRec; ///< Pointer to previous track hit with same HitForRec
-
-  ClassDef(AliMUONTrackHit, 1) // Reconstructed track hit in ALICE dimuon spectrometer
-    };
-       
-#endif
index fc1afb1..305a0da 100644 (file)
@@ -135,7 +135,7 @@ AliMUONTrackK::AliMUONTrackK(AliMUONTrackReconstructor *TrackReconstructor, TClo
   //__________________________________________________________________________
 AliMUONTrackK::AliMUONTrackK(AliMUONSegment *segment)
   //: AliMUONTrack(segment, segment, fgTrackReconstructor)
-  : AliMUONTrack(NULL, segment, fgTrackReconstructor),
+  : AliMUONTrack(NULL, segment),
     fStartSegment(segment),
     fPosition(0.),
     fPositionNew(0.),
@@ -1276,8 +1276,6 @@ void AliMUONTrackK::FillMUONTrack(void)
 {
 /// Compute track parameters at hit positions (as for AliMUONTrack)
 
-  // Set number of hits per track
-  SetNTrackHits(fNmbTrackHits);
   // Set Chi2
   SetFitFMin(fChi2);
 
@@ -1295,7 +1293,7 @@ void AliMUONTrackK::FillMUONTrack(void)
       // Take saved info
       SetTrackParam(&trackParam, (TMatrixD*)fParSmooth->UncheckedAt(i), ((AliMUONHitForRec*)((*fTrackHits)[i]))->GetZ());
     }
-    AddTrackParamAtHit(&trackParam); 
+    AddTrackParamAtHit(&trackParam,(AliMUONHitForRec*)fTrackHits->UncheckedAt(i));
     // Fill array of HitForRec's
     AddHitForRecAtHit((AliMUONHitForRec*)fTrackHits->UncheckedAt(i)); 
   }
index e3dba86..156fd69 100644 (file)
@@ -33,6 +33,7 @@
 #include "AliMagF.h" 
 #include "AliLog.h" 
 #include "AliTracker.h"
+#include "AliMUONHitForRec.h"
 
 ClassImp(AliMUONTrackParam) // Class implementation in ROOT context
 
@@ -44,7 +45,9 @@ AliMUONTrackParam::AliMUONTrackParam()
     fNonBendingSlope(0.),
     fZ(0.),
     fBendingCoor(0.),
-    fNonBendingCoor(0.)
+    fNonBendingCoor(0.),
+    fkField(0x0),
+    fHitForRecPtr(0x0)
 {
 // Constructor
   // get field from outside
@@ -53,8 +56,23 @@ AliMUONTrackParam::AliMUONTrackParam()
 }
 
   //_________________________________________________________________________
-AliMUONTrackParam& 
-AliMUONTrackParam::operator=(const AliMUONTrackParam& theMUONTrackParam)
+AliMUONTrackParam::AliMUONTrackParam(const AliMUONTrackParam& theMUONTrackParam)
+  : TObject(theMUONTrackParam),
+    fInverseBendingMomentum(theMUONTrackParam.fInverseBendingMomentum), 
+    fBendingSlope(theMUONTrackParam.fBendingSlope),
+    fNonBendingSlope(theMUONTrackParam.fNonBendingSlope),
+    fZ(theMUONTrackParam.fZ),
+    fBendingCoor(theMUONTrackParam.fBendingCoor),
+    fNonBendingCoor(theMUONTrackParam.fNonBendingCoor),
+    fkField(theMUONTrackParam.fkField),
+    fHitForRecPtr(theMUONTrackParam.fHitForRecPtr)
+{
+  // Copy constructor
+
+}
+
+  //_________________________________________________________________________
+AliMUONTrackParam& AliMUONTrackParam::operator=(const AliMUONTrackParam& theMUONTrackParam)
 {
   // Asignment operator
   if (this == &theMUONTrackParam)
@@ -69,21 +87,51 @@ AliMUONTrackParam::operator=(const AliMUONTrackParam& theMUONTrackParam)
   fZ                      =  theMUONTrackParam.fZ; 
   fBendingCoor            =  theMUONTrackParam.fBendingCoor; 
   fNonBendingCoor         =  theMUONTrackParam.fNonBendingCoor;
+  fkField                 =  theMUONTrackParam.fkField;
+  fHitForRecPtr           =  theMUONTrackParam.fHitForRecPtr;
 
   return *this;
 }
-  //_________________________________________________________________________
-AliMUONTrackParam::AliMUONTrackParam(const AliMUONTrackParam& theMUONTrackParam)
-  : TObject(theMUONTrackParam),
-    fInverseBendingMomentum(theMUONTrackParam.fInverseBendingMomentum), 
-    fBendingSlope(theMUONTrackParam.fBendingSlope),
-    fNonBendingSlope(theMUONTrackParam.fNonBendingSlope),
-    fZ(theMUONTrackParam.fZ),
-    fBendingCoor(theMUONTrackParam.fBendingCoor),
-    fNonBendingCoor(theMUONTrackParam.fNonBendingCoor)
+
+  //__________________________________________________________________________
+AliMUONTrackParam::~AliMUONTrackParam()
 {
-  // Copy constructor
+/// Destructor
+/// Update the number of TrackHit's connected to the attached HitForRec if any
+  if (fHitForRecPtr) fHitForRecPtr->SetNTrackHits(fHitForRecPtr->GetNTrackHits() - 1); // decrement NTrackHits of hit
+}
 
+  //__________________________________________________________________________
+void AliMUONTrackParam::SetTrackParam(AliMUONTrackParam& theMUONTrackParam)
+{
+  /// Set track parameters from "TrackParam" leaving pointer to fHitForRecPtr unchanged
+  fInverseBendingMomentum =  theMUONTrackParam.fInverseBendingMomentum; 
+  fBendingSlope           =  theMUONTrackParam.fBendingSlope; 
+  fNonBendingSlope        =  theMUONTrackParam.fNonBendingSlope; 
+  fZ                      =  theMUONTrackParam.fZ; 
+  fBendingCoor            =  theMUONTrackParam.fBendingCoor; 
+  fNonBendingCoor         =  theMUONTrackParam.fNonBendingCoor;
+  
+}
+
+  //__________________________________________________________________________
+AliMUONHitForRec* AliMUONTrackParam::GetHitForRecPtr(void) const
+{
+/// return pointer to HitForRec attached to the current TrackParam
+/// this method should not be called when fHitForRecPtr == NULL
+  if (!fHitForRecPtr) AliWarning("AliMUONTrackParam::GetHitForRecPtr: fHitForRecPtr == NULL");
+  return fHitForRecPtr;
+}
+
+  //__________________________________________________________________________
+Int_t AliMUONTrackParam::Compare(const TObject* TrackParam) const
+{
+/// "Compare" function to sort with decreasing Z (spectro. muon Z <0).
+/// Returns 1 (0, -1) if Z of current TrackHit
+/// is smaller than (equal to, larger than) Z of TrackHit
+  if (fHitForRecPtr->GetZ() < ((AliMUONTrackParam*)TrackParam)->fHitForRecPtr->GetZ()) return(1);
+  else if (fHitForRecPtr->GetZ() == ((AliMUONTrackParam*)TrackParam)->fHitForRecPtr->GetZ()) return(0);
+  else return(-1);
 }
 
   //_________________________________________________________________________
index 694ebe5..4210654 100644 (file)
@@ -15,6 +15,7 @@
 ////////////////////////////////////////////////////
 
 #include <TObject.h>
+#include "AliMUONHitForRec.h"
 
 class AliESDMuonTrack;
 class AliMagF;
@@ -23,10 +24,10 @@ class AliMUONTrackParam : public TObject
 {
  public:
   AliMUONTrackParam(); // Constructor
-  virtual ~AliMUONTrackParam(){} // Destructor
+  virtual ~AliMUONTrackParam(); // Destructor
   
-  AliMUONTrackParam(const AliMUONTrackParam& rhs);// copy constructor (should be added per default !)
-  AliMUONTrackParam& operator=(const  AliMUONTrackParam& rhs);// (should be added per default !)
+  AliMUONTrackParam(const AliMUONTrackParam& theMUONTrackParam);
+  AliMUONTrackParam& operator=(const  AliMUONTrackParam& theMUONTrackParam);
 
   void GetParamFrom(const AliESDMuonTrack& esdMuonTrack);
   void SetParamFor(AliESDMuonTrack& esdMuonTrack);
@@ -45,11 +46,18 @@ class AliMUONTrackParam : public TObject
   void     SetBendingCoor(Double_t BendingCoor) {fBendingCoor = BendingCoor;}
   Double_t GetNonBendingCoor(void) const {return fNonBendingCoor;}
   void     SetNonBendingCoor(Double_t NonBendingCoor) {fNonBendingCoor = NonBendingCoor;}
+  void              SetTrackParam(AliMUONTrackParam& TrackParam);
+  AliMUONHitForRec* GetHitForRecPtr(void) const;
+  void              SetHitForRecPtr(AliMUONHitForRec* HitForRec) {fHitForRecPtr = HitForRec;}
+  
   Double_t Px() const;  // return px
   Double_t Py() const;  // return py
   Double_t Pz() const;  // return pz
   Double_t P()  const;  // return total momentum
 
+  Bool_t IsSortable () const {return kTRUE;} // necessary for sorting TClonesArray of TrackHit's
+  Int_t Compare(const TObject* TrackParam) const; // "Compare" function for sorting
+
   void ExtrapToZ(Double_t Z);
   void ExtrapToStation(Int_t Station, AliMUONTrackParam *TrackParam);
   void ExtrapToVertex(Double_t xVtx, Double_t yVtx, Double_t zVtx);  // extrapolation to vertex through the absorber (with true vertex) 
@@ -86,8 +94,10 @@ class AliMUONTrackParam : public TObject
   void GetFromGeant3Parameters(Double_t *VGeant3, Double_t Charge);
 
   void GetField(Double_t *Position, Double_t *Field) const;
-
-  ClassDef(AliMUONTrackParam, 1) // Track parameters in ALICE dimuon spectrometer
+  
+  AliMUONHitForRec *fHitForRecPtr; //!< Pointer to associated HitForRec if any
+  
+  ClassDef(AliMUONTrackParam, 2) // Track parameters in ALICE dimuon spectrometer
     };
        
 #endif
index 449b3ba..701a597 100644 (file)
 #include "AliMUONGlobalTrigger.h"
 #include "AliMUONSegment.h"
 #include "AliMUONTrack.h"
-#include "AliMUONTrackHit.h"
 #include "AliMagF.h"
 #include "AliMUONTrackK.h" 
 #include "AliLog.h"
 #include "AliTracker.h"
+#include <TVirtualFitter.h>
 
 //************* Defaults parameters for reconstruction
 const Double_t AliMUONTrackReconstructor::fgkDefaultMinBendingMomentum = 3.0;
@@ -62,7 +62,6 @@ const Double_t AliMUONTrackReconstructor::fgkDefaultMaxChi2 = 100.0;
 const Double_t AliMUONTrackReconstructor::fgkDefaultMaxSigma2Distance = 16.0;
 const Double_t AliMUONTrackReconstructor::fgkDefaultBendingResolution = 0.01;
 const Double_t AliMUONTrackReconstructor::fgkDefaultNonBendingResolution = 0.144;
-const Double_t AliMUONTrackReconstructor::fgkDefaultChamberThicknessInX0 = 0.03;
 // Simple magnetic field:
 // Value taken from macro MUONtracking.C: 0.7 T, hence 7 kG
 // Length and Position from reco_muon.F, with opposite sign:
@@ -74,6 +73,16 @@ const Double_t AliMUONTrackReconstructor::fgkDefaultSimpleBLength = 428.0;
 const Double_t AliMUONTrackReconstructor::fgkDefaultSimpleBPosition = 1019.0;
 const Double_t AliMUONTrackReconstructor::fgkDefaultEfficiency = 0.95;
 
+TVirtualFitter* AliMUONTrackReconstructor::fgFitter = NULL; 
+
+// Functions to be minimized with Minuit
+void TrackChi2(Int_t &NParam, Double_t *Gradient, Double_t &Chi2, Double_t *Param, Int_t Flag);
+void TrackChi2MCS(Int_t &NParam, Double_t *Gradient, Double_t &Chi2, Double_t *Param, Int_t Flag);
+
+void mnvertLocal(Double_t* a, Int_t l, Int_t m, Int_t n, Int_t& ifail);
+
+Double_t MultipleScatteringAngle2(AliMUONTrackParam *param);
+
 ClassImp(AliMUONTrackReconstructor) // Class implementation in ROOT context
 
 //__________________________________________________________________________
@@ -84,24 +93,39 @@ AliMUONTrackReconstructor::AliMUONTrackReconstructor(AliMUONData* data)
     fMaxBendingMomentum(fgkDefaultMaxBendingMomentum),
     fMaxChi2(fgkDefaultMaxChi2),
     fMaxSigma2Distance(fgkDefaultMaxSigma2Distance),
+    fRMin(0x0),
+    fRMax(0x0),
+    fSegmentMaxDistBending(0x0),
+    fSegmentMaxDistNonBending(0x0),
     fBendingResolution(fgkDefaultBendingResolution),
     fNonBendingResolution(fgkDefaultNonBendingResolution),
-    fChamberThicknessInX0(fgkDefaultChamberThicknessInX0),
+    fChamberThicknessInX0(AliMUONConstants::DefaultChamberThicknessInX0()),
     fSimpleBValue(fgkDefaultSimpleBValue),
     fSimpleBLength(fgkDefaultSimpleBLength),
     fSimpleBPosition(fgkDefaultSimpleBPosition),
     fEfficiency(fgkDefaultEfficiency),
     fHitsForRecPtr(0x0),
     fNHitsForRec(0),
+    fNHitsForRecPerChamber(0x0),
+    fIndexOfFirstHitForRecPerChamber(0x0),
+    fSegmentsPtr(0x0),
+    fNSegments(0x0),
     fRecTracksPtr(0x0),
     fNRecTracks(0),
-    fRecTrackHitsPtr(0x0),
-    fNRecTrackHits(0),
     fMUONData(data),
     fMuons(0),
     fTriggerTrack(new AliMUONTriggerTrack()),
     fTriggerCircuit(0x0)
 {
+  
+  // Memory allocation
+  fRMin = new Double_t[AliMUONConstants::NTrackingCh()];
+  fRMax = new Double_t[AliMUONConstants::NTrackingCh()];
+  fSegmentMaxDistBending = new Double_t[AliMUONConstants::NTrackingSt()];
+  fSegmentMaxDistNonBending = new Double_t[AliMUONConstants::NTrackingSt()];
+  fNHitsForRecPerChamber = new Int_t[AliMUONConstants::NTrackingCh()];
+  fIndexOfFirstHitForRecPerChamber = new Int_t[AliMUONConstants::NTrackingCh()];
+  
   // Constructor for class AliMUONTrackReconstructor
   SetReconstructionParametersToDefaults();
 
@@ -111,7 +135,9 @@ AliMUONTrackReconstructor::AliMUONTrackReconstructor(AliMUONData* data)
 
   // Memory allocation for the TClonesArray's of segments in stations
   // Is 2000 the right size ????
-  for (Int_t st = 0; st < AliMUONConstants::NTrackingCh()/2; st++) {
+  fSegmentsPtr = new TClonesArray*[AliMUONConstants::NTrackingSt()];
+  fNSegments = new Int_t[AliMUONConstants::NTrackingSt()];
+  for (Int_t st = 0; st < AliMUONConstants::NTrackingSt(); st++) {
     fSegmentsPtr[st] = new TClonesArray("AliMUONSegment", 2000);
     fNSegments[st] = 0; // really needed or GetEntriesFast sufficient ????
   }
@@ -119,10 +145,6 @@ AliMUONTrackReconstructor::AliMUONTrackReconstructor(AliMUONData* data)
   // Is 10 the right size ????
   fRecTracksPtr = new TClonesArray("AliMUONTrack", 10);
 
-  // Memory allocation for the TClonesArray of hits on reconstructed tracks
-  // Is 100 the right size ????
-  fRecTrackHitsPtr = new TClonesArray("AliMUONTrackHit", 100);
-
   const AliMagF* kField = AliTracker::GetFieldMap();
   if (!kField) AliFatal("No field available");
  // Sign of fSimpleBValue according to sign of Bx value at (50,50,-950).
@@ -142,11 +164,19 @@ AliMUONTrackReconstructor::AliMUONTrackReconstructor(AliMUONData* data)
 AliMUONTrackReconstructor::~AliMUONTrackReconstructor(void)
 {
   // Destructor for class AliMUONTrackReconstructor
-  delete fHitsForRecPtr; // Correct destruction of everything ???? or delete [] ????
-  for (Int_t st = 0; st < AliMUONConstants::NTrackingCh()/2; st++)
-    delete fSegmentsPtr[st]; // Correct destruction of everything ????
-
+  delete [] fRMin;
+  delete [] fRMax;
+  delete [] fSegmentMaxDistBending;
+  delete [] fSegmentMaxDistNonBending;
+  delete [] fNHitsForRecPerChamber;
+  delete [] fIndexOfFirstHitForRecPerChamber;
   delete fTriggerTrack;
+  delete fHitsForRecPtr;
+  // Correct destruction of everything ????
+  for (Int_t st = 0; st < AliMUONConstants::NTrackingSt(); st++) delete fSegmentsPtr[st];
+  delete [] fSegmentsPtr;
+  delete [] fNSegments;
+  delete fRecTracksPtr;
 }
   //__________________________________________________________________________
 void AliMUONTrackReconstructor::SetReconstructionParametersToDefaults(void)
@@ -221,7 +251,6 @@ void AliMUONTrackReconstructor::EventReconstruct(void)
   // To reconstruct one event
   AliDebug(1,"Enter EventReconstruct");
   ResetTracks(); //AZ
-  ResetTrackHits(); //AZ 
   ResetSegments(); //AZ
   ResetHitsForRec(); //AZ
   MakeEventToBeReconstructed();
@@ -229,10 +258,10 @@ void AliMUONTrackReconstructor::EventReconstruct(void)
   MakeTracks();
   if (fMUONData->IsTriggerTrackBranchesInTree()) 
     ValidateTracksWithTrigger(); 
-
+  
   // Add tracks to MUON data container 
-  for(Int_t i=0; i<GetNRecTracks(); i++) {
-    AliMUONTrack * track = (AliMUONTrack*) GetRecTracksPtr()->At(i);
+  for(Int_t i=0; i<fNRecTracks; i++) {
+    AliMUONTrack * track = (AliMUONTrack*) fRecTracksPtr->At(i);
     fMUONData->AddRecTrack(*track);
   }
 }
@@ -281,16 +310,6 @@ void AliMUONTrackReconstructor::ResetTracks(void)
   return;
 }
 
-
-  //__________________________________________________________________________
-void AliMUONTrackReconstructor::ResetTrackHits(void)
-{
-  // To reset the TClonesArray of hits on reconstructed tracks
-  if (fRecTrackHitsPtr) fRecTrackHitsPtr->Clear();
-  fNRecTrackHits = 0;
-  return;
-}
-
   //__________________________________________________________________________
 void AliMUONTrackReconstructor::MakeEventToBeReconstructed(void)
 {
@@ -575,7 +594,6 @@ void AliMUONTrackReconstructor::MakeTracks(void)
   AliDebug(1,"Enter MakeTracks");
   // The order may be important for the following Reset's
   //AZ ResetTracks();
-  //AZ ResetTrackHits();
   if (fTrackMethod != 1) { //AZ - Kalman filter
     MakeTrackCandidatesK();
     if (fRecTracksPtr->GetEntriesFast() == 0) return;
@@ -594,7 +612,6 @@ void AliMUONTrackReconstructor::MakeTracks(void)
     FollowTracks();
     // Remove double tracks
     RemoveDoubleTracks();
-    UpdateTrackParamAtHit();
     UpdateHitForRecAtHit();
   }
   return;
@@ -605,15 +622,49 @@ void AliMUONTrackReconstructor::ValidateTracksWithTrigger(void)
 {
   // Try to match track from tracking system with trigger track
   AliMUONTrack *track;
-  TClonesArray *recTriggerTracks;
+  AliMUONTrackParam trackParam; 
+  AliMUONTriggerTrack *triggerTrack;
   
   fMUONData->SetTreeAddress("RL");
   fMUONData->GetRecTriggerTracks();
-  recTriggerTracks = fMUONData->RecTriggerTracks();
-
+  TClonesArray *recTriggerTracks = fMUONData->RecTriggerTracks();
+  
+  Bool_t MatchTrigger;
+  Double_t distSigma[3]={1,1,0.02}; // sigma of distributions (trigger-track) X,Y,slopeY
+  Double_t distTriggerTrack[3];
+  Double_t Chi2MatchTrigger, xTrack, yTrack, ySlopeTrack, dTrigTrackMin2, dTrigTrack2 = 0.;
+  
   track = (AliMUONTrack*) fRecTracksPtr->First();
   while (track) {
-    track->MatchTriggerTrack(recTriggerTracks);
+    MatchTrigger = kFALSE;
+    Chi2MatchTrigger = 0.;
+    
+    trackParam = *((AliMUONTrackParam*) (track->GetTrackParamAtHit()->Last()));
+    trackParam.ExtrapToZ(AliMUONConstants::DefaultChamberZ(10)); // extrap to 1st trigger chamber
+    
+    xTrack = trackParam.GetNonBendingCoor();
+    yTrack = trackParam.GetBendingCoor();
+    ySlopeTrack = trackParam.GetBendingSlope();
+    dTrigTrackMin2 = 999.;
+  
+    triggerTrack = (AliMUONTriggerTrack*) recTriggerTracks->First();
+    while(triggerTrack){
+      distTriggerTrack[0] = (triggerTrack->GetX11()-xTrack)/distSigma[0];
+      distTriggerTrack[1] = (triggerTrack->GetY11()-yTrack)/distSigma[1];
+      distTriggerTrack[2] = (TMath::Tan(triggerTrack->GetThetay())-ySlopeTrack)/distSigma[2];
+      dTrigTrack2 = 0.;
+      for (Int_t iVar = 0; iVar < 3; iVar++) dTrigTrack2 += distTriggerTrack[iVar]*distTriggerTrack[iVar];
+      if (dTrigTrack2 < dTrigTrackMin2 && dTrigTrack2 < GetMaxSigma2Distance()) {
+        dTrigTrackMin2 = dTrigTrack2;
+        MatchTrigger = kTRUE;
+        Chi2MatchTrigger =  dTrigTrack2/3.; // Normalized Chi2, 3 variables (X,Y,slopeY)
+      }
+      triggerTrack = (AliMUONTriggerTrack*) recTriggerTracks->After(triggerTrack);
+    }
+    
+    track->SetMatchTrigger(MatchTrigger);
+    track->SetChi2MatchTrigger(Chi2MatchTrigger);
+    
     track = (AliMUONTrack*) fRecTracksPtr->After(track);
   }
 
@@ -697,6 +748,39 @@ Bool_t AliMUONTrackReconstructor::MakeTriggerTracks(void)
 }
 
   //__________________________________________________________________________
+void AliMUONTrackReconstructor::MakeTrackCandidates(void)
+{
+  // To make track candidates
+  // with at least 3 aligned points in stations(1..) 4 and 5
+  // (two Segment's or one Segment and one HitForRec)
+  Int_t begStation, iBegSegment, nbCan1Seg1Hit, nbCan2Seg;
+  AliMUONSegment *begSegment;
+  AliDebug(1,"Enter MakeTrackCandidates");
+  // Loop over stations(1..) 5 and 4 for the beginning segment
+  for (begStation = 4; begStation > 2; begStation--) {
+    // Loop over segments in the beginning station
+    for (iBegSegment = 0; iBegSegment < fNSegments[begStation]; iBegSegment++) {
+      // pointer to segment
+      begSegment = (AliMUONSegment*) ((*fSegmentsPtr[begStation])[iBegSegment]);
+      AliDebug(2,Form("Look for TrackCandidate's with Segment %d  in Station(0..) %d", iBegSegment, begStation));
+      // Look for track candidates with two segments,
+      // "begSegment" and all compatible segments in other station.
+      // Only for beginning station(1..) 5
+      // because candidates with 2 segments have to looked for only once.
+      if (begStation == 4)
+       nbCan2Seg = MakeTrackCandidatesWithTwoSegments(begSegment);
+      // Look for track candidates with one segment and one point,
+      // "begSegment" and all compatible HitForRec's in other station.
+      // Only if "begSegment" does not belong already to a track candidate.
+      // Is that a too strong condition ????
+      if (!(begSegment->GetInTrack()))
+       nbCan1Seg1Hit = MakeTrackCandidatesWithOneSegmentAndOnePoint(begSegment);
+    } // for (iBegSegment = 0;...
+  } // for (begStation = 4;...
+  return;
+}
+
+  //__________________________________________________________________________
 Int_t AliMUONTrackReconstructor::MakeTrackCandidatesWithTwoSegments(AliMUONSegment *BegSegment)
 {
   // To make track candidates with two segments in stations(1..) 4 and 5,
@@ -735,8 +819,9 @@ Int_t AliMUONTrackReconstructor::MakeTrackCandidatesWithTwoSegments(AliMUONSegme
       // to be done in track constructor ????
       BegSegment->SetInTrack(kTRUE);
       endSegment->SetInTrack(kTRUE);
-      recTrack = new ((*fRecTracksPtr)[fNRecTracks])
-       AliMUONTrack(BegSegment, endSegment, this);
+      recTrack = new ((*fRecTracksPtr)[fNRecTracks]) AliMUONTrack(BegSegment, endSegment);
+      // Set track parameters at vertex from last stations 4 & 5
+      CalcTrackParamAtVertex(recTrack);
       fNRecTracks++;
       if (AliLog::GetGlobalDebugLevel() > 1) recTrack->RecursiveDump();
       // increment number of track candidates with 2 segments
@@ -793,8 +878,9 @@ Int_t AliMUONTrackReconstructor::MakeTrackCandidatesWithOneSegmentAndOnePoint(Al
        // flag for beginning segments in one track:
        // to be done in track constructor ????
        BegSegment->SetInTrack(kTRUE);
-       recTrack = new ((*fRecTracksPtr)[fNRecTracks])
-         AliMUONTrack(BegSegment, hit, this);
+       recTrack = new ((*fRecTracksPtr)[fNRecTracks]) AliMUONTrack(BegSegment, hit);
+        // Set track parameters at vertex from last stations 4 & 5
+        CalcTrackParamAtVertex(recTrack);
        // the right place to eliminate "double counting" ???? how ????
        fNRecTracks++;
        if (AliLog::GetGlobalDebugLevel() > 1) recTrack->RecursiveDump();
@@ -808,36 +894,35 @@ Int_t AliMUONTrackReconstructor::MakeTrackCandidatesWithOneSegmentAndOnePoint(Al
 }
 
   //__________________________________________________________________________
-void AliMUONTrackReconstructor::MakeTrackCandidates(void)
+void AliMUONTrackReconstructor::CalcTrackParamAtVertex(AliMUONTrack *Track)
 {
-  // To make track candidates
-  // with at least 3 aligned points in stations(1..) 4 and 5
-  // (two Segment's or one Segment and one HitForRec)
-  Int_t begStation, iBegSegment, nbCan1Seg1Hit, nbCan2Seg;
-  AliMUONSegment *begSegment;
-  AliDebug(1,"Enter MakeTrackCandidates");
-  // Loop over stations(1..) 5 and 4 for the beginning segment
-  for (begStation = 4; begStation > 2; begStation--) {
-    // Loop over segments in the beginning station
-    for (iBegSegment = 0; iBegSegment < fNSegments[begStation]; iBegSegment++) {
-      // pointer to segment
-      begSegment = (AliMUONSegment*) ((*fSegmentsPtr[begStation])[iBegSegment]);
-      AliDebug(2,Form("Look for TrackCandidate's with Segment %d  in Station(0..) %d", iBegSegment, begStation));
-      // Look for track candidates with two segments,
-      // "begSegment" and all compatible segments in other station.
-      // Only for beginning station(1..) 5
-      // because candidates with 2 segments have to looked for only once.
-      if (begStation == 4)
-       nbCan2Seg = MakeTrackCandidatesWithTwoSegments(begSegment);
-      // Look for track candidates with one segment and one point,
-      // "begSegment" and all compatible HitForRec's in other station.
-      // Only if "begSegment" does not belong already to a track candidate.
-      // Is that a too strong condition ????
-      if (!(begSegment->GetInTrack()))
-       nbCan1Seg1Hit = MakeTrackCandidatesWithOneSegmentAndOnePoint(begSegment);
-    } // for (iBegSegment = 0;...
-  } // for (begStation = 4;...
-  return;
+  // Set track parameters at vertex.
+  // TrackHit's are assumed to be only in stations(1..) 4 and 5,
+  // and sorted according to increasing Z..
+  // Parameters are calculated from information in HitForRec's
+  // of first and last TrackHit's.
+  AliMUONTrackParam *trackParamVertex = Track->GetTrackParamAtVertex(); // pointer to track parameters at vertex
+  // Pointer to HitForRec attached to first TrackParamAtHit
+  AliMUONHitForRec *firstHit = ((AliMUONTrackParam*) (Track->GetTrackParamAtHit()->First()))->GetHitForRecPtr();
+  // Pointer to HitForRec attached to last TrackParamAtHit
+  AliMUONHitForRec *lastHit = ((AliMUONTrackParam*) (Track->GetTrackParamAtHit()->Last()))->GetHitForRecPtr();
+  // Z difference between first and last hits
+  Double_t deltaZ = firstHit->GetZ() - lastHit->GetZ();
+  // bending slope in stations(1..) 4 and 5
+  Double_t bendingSlope = (firstHit->GetBendingCoor() - lastHit->GetBendingCoor()) / deltaZ;
+  trackParamVertex->SetBendingSlope(bendingSlope);
+  // impact parameter
+  Double_t impactParam = firstHit->GetBendingCoor() - bendingSlope * firstHit->GetZ();
+  // signed bending momentum
+  trackParamVertex->SetInverseBendingMomentum(1.0 / GetBendingMomentumFromImpactParam(impactParam));
+  // bending slope at vertex
+  trackParamVertex->SetBendingSlope(bendingSlope + impactParam / GetSimpleBPosition());
+  // non bending slope
+  trackParamVertex->SetNonBendingSlope((firstHit->GetNonBendingCoor() - lastHit->GetNonBendingCoor()) / deltaZ);
+  // vertex coordinates at (0,0,0)
+  trackParamVertex->SetZ(0.0);
+  trackParamVertex->SetBendingCoor(0.0);
+  trackParamVertex->SetNonBendingCoor(0.0);
 }
 
   //__________________________________________________________________________
@@ -862,15 +947,13 @@ void AliMUONTrackReconstructor::FollowTracks(void)
   track = (AliMUONTrack*) fRecTracksPtr->First();
   trackIndex = -1;
   while (track) {
-    // Follow function for each track candidate ????
     trackIndex++;
     nextTrack = (AliMUONTrack*) fRecTracksPtr->After(track); // prepare next track
     AliDebug(2,Form("FollowTracks: track candidate(0..): %d", trackIndex));
-    // Fit track candidate
-    track->SetFitMCS(0); // without multiple Coulomb scattering
-    track->SetFitNParam(3); // with 3 parameters (X = Y = 0)
-    track->SetFitStart(0); // from parameters at vertex
-    track->Fit();
+    // Fit track candidate from parameters at vertex
+    // -> with 3 parameters (X_vertex and Y_vertex are fixed)
+    // without multiple Coulomb scattering
+    Fit(track,0,0);
     if (AliLog::GetGlobalDebugLevel()> 2) {
       cout << "FollowTracks: track candidate(0..): " << trackIndex
           << " after fit in stations(0..) 3 and 4" << endl;
@@ -881,8 +964,7 @@ void AliMUONTrackReconstructor::FollowTracks(void)
     // otherwise: majority coincidence 2 !!!!
     for (station = 2; station >= 0; station--) {
       // Track parameters at first track hit (smallest Z)
-      trackParam1 = ((AliMUONTrackHit*)
-                    (track->GetTrackHitsPtr()->First()))->GetTrackParam();
+      trackParam1 = (AliMUONTrackParam*) (track->GetTrackParamAtHit()->First());
       // extrapolation to station
       trackParam1->ExtrapToStation(station, trackParam);
       extrapSegment = new AliMUONSegment(); //  empty segment
@@ -900,11 +982,9 @@ void AliMUONTrackReconstructor::FollowTracks(void)
       dZ3 = AliMUONConstants::DefaultChamberZ(2 * station) -
            AliMUONConstants::DefaultChamberZ(2 * station + 1);
       extrapSegment->SetBendingCoorReso2(fBendingResolution * fBendingResolution);
-      extrapSegment->
-       SetNonBendingCoorReso2(fNonBendingResolution * fNonBendingResolution);
-      extrapSegment->UpdateFromStationTrackParam
-       (trackParam, mcsFactor, dZ1, dZ2, dZ3, station,
-        trackParam1->GetInverseBendingMomentum());
+      extrapSegment->SetNonBendingCoorReso2(fNonBendingResolution * fNonBendingResolution);
+      extrapSegment->UpdateFromStationTrackParam(trackParam, mcsFactor, dZ1, dZ2, dZ3, station,
+                                                trackParam1->GetInverseBendingMomentum());
       bestChi2 = 5.0;
       bestSegment = NULL;
       if (AliLog::GetGlobalDebugLevel() > 2) {
@@ -921,14 +1001,13 @@ void AliMUONTrackReconstructor::FollowTracks(void)
        segment = (AliMUONSegment*) ((*fSegmentsPtr[station])[iSegment]);
        // correction of corrected segment (fBendingCoor and fNonBendingCoor)
        // according to real Z value of "segment" and slopes of "extrapSegment"
-       (&(trackParam[0]))->ExtrapToZ(segment->GetZ());
-       (&(trackParam[1]))->ExtrapToZ(segment->GetZ());
+       trackParam[0].ExtrapToZ(segment->GetZ());
+       trackParam[1].ExtrapToZ(segment->GetZ()); // now same as trackParam[0] !?!?!?!?!?!
        extrapSegment->SetBendingCoor((&(trackParam[0]))->GetBendingCoor());
        extrapSegment->SetNonBendingCoor((&(trackParam[0]))->GetNonBendingCoor());
        extrapSegment->SetBendingSlope((&(trackParam[0]))->GetBendingSlope());
        extrapSegment->SetNonBendingSlope((&(trackParam[0]))->GetNonBendingSlope());
-       chi2 = segment->
-         NormalizedChi2WithSegment(extrapSegment, maxSigma2Distance);
+       chi2 = segment->NormalizedChi2WithSegment(extrapSegment, maxSigma2Distance);
        if (chi2 < bestChi2) {
          // update best Chi2 and Segment if better found
          bestSegment = segment;
@@ -937,16 +1016,13 @@ void AliMUONTrackReconstructor::FollowTracks(void)
       }
       if (bestSegment) {
        // best segment found: add it to track candidate
-       (&(trackParam[0]))->ExtrapToZ(bestSegment->GetZ());
-       (&(trackParam[1]))->ExtrapToZ(bestSegment->GetZ());
-       track->AddSegment(bestSegment);
-       // set track parameters at these two TrakHit's
-       track->SetTrackParamAtHit(track->GetNTrackHits() - 2, &(trackParam[0]));
-       track->SetTrackParamAtHit(track->GetNTrackHits() - 1, &(trackParam[1]));
+       trackParam[0].ExtrapToZ(bestSegment->GetZ());
+       track->AddTrackParamAtHit(&(trackParam[0]),bestSegment->GetHitForRec1());
+       trackParam[1].ExtrapToZ(bestSegment->GetZ()); // now same as trackParam[0] !?!?!?!?!?!
+       track->AddTrackParamAtHit(&(trackParam[1]),bestSegment->GetHitForRec2());
        AliDebug(3, Form("FollowTracks: track candidate(0..): %d  Added segment in station(0..): %d", trackIndex, station));
        if (AliLog::GetGlobalDebugLevel()>2) track->RecursiveDump();
-      }
-      else {
+      } else {
        // No best segment found:
        // Look for best compatible HitForRec in station:
        // should consider all possibilities ????
@@ -960,17 +1036,12 @@ void AliMUONTrackReconstructor::FollowTracks(void)
        // Loop over chambers of the station
        for (chInStation = 0; chInStation < 2; chInStation++) {
          ch = 2 * station + chInStation;
-         for (iHit = fIndexOfFirstHitForRecPerChamber[ch];
-              iHit < fIndexOfFirstHitForRecPerChamber[ch] +
-                fNHitsForRecPerChamber[ch];
-              iHit++) {
+         for (iHit = fIndexOfFirstHitForRecPerChamber[ch]; iHit < fIndexOfFirstHitForRecPerChamber[ch]+fNHitsForRecPerChamber[ch]; iHit++) {
            hit = (AliMUONHitForRec*) ((*fHitsForRecPtr)[iHit]);
            // coordinates of extrapolated hit
-           (&(trackParam[chInStation]))->ExtrapToZ(hit->GetZ());
-           extrapHit->
-             SetBendingCoor((&(trackParam[chInStation]))->GetBendingCoor());
-           extrapHit->
-             SetNonBendingCoor((&(trackParam[chInStation]))->GetNonBendingCoor());
+           trackParam[chInStation].ExtrapToZ(hit->GetZ());
+           extrapHit->SetBendingCoor((&(trackParam[chInStation]))->GetBendingCoor());
+           extrapHit->SetNonBendingCoor((&(trackParam[chInStation]))->GetNonBendingCoor());
            // resolutions from "extrapSegment"
            extrapHit->SetBendingReso2(extrapSegment->GetBendingCoorReso2());
            extrapHit->SetNonBendingReso2(extrapSegment->GetNonBendingCoorReso2());
@@ -987,21 +1058,18 @@ void AliMUONTrackReconstructor::FollowTracks(void)
        }
        if (bestHit) {
          // best hit found: add it to track candidate
-         (&(trackParam[chBestHit]))->ExtrapToZ(bestHit->GetZ());
-         track->AddHitForRec(bestHit);
-         // set track parameters at this TrackHit
-         track->SetTrackParamAtHit(track->GetNTrackHits() - 1,
-                                   &(trackParam[chBestHit]));
+         trackParam[chBestHit].ExtrapToZ(bestHit->GetZ());
+         track->AddTrackParamAtHit(&(trackParam[chBestHit]),bestHit);
          if (AliLog::GetGlobalDebugLevel() > 2) {
            cout << "FollowTracks: track candidate(0..): " << trackIndex
                 << " Added hit in station(0..): " << station << endl;
            track->RecursiveDump();
          }
-       }
-       else {
+       } else {
          // Remove current track candidate
          // and corresponding TrackHit's, ...
-         track->Remove();
+         fRecTracksPtr->Remove(track);
+         fNRecTracks--;
          delete extrapSegment;
          delete extrapHit;
          break; // stop the search for this candidate:
@@ -1010,28 +1078,26 @@ void AliMUONTrackReconstructor::FollowTracks(void)
        delete extrapHit;
       }
       delete extrapSegment;
-      // Sort track hits according to increasing Z
-      track->GetTrackHitsPtr()->Sort();
+      // Sort TrackParamAtHit according to increasing Z
+      track->GetTrackParamAtHit()->Sort();
       // Update track parameters at first track hit (smallest Z)
-      trackParam1 = ((AliMUONTrackHit*)
-                    (track->GetTrackHitsPtr()->First()))->GetTrackParam();
+      trackParam1 = (AliMUONTrackParam*) (track->GetTrackParamAtHit()->First());
       bendingMomentum = 0.;
       if (TMath::Abs(trackParam1->GetInverseBendingMomentum()) > 0.)
        bendingMomentum = TMath::Abs(1/(trackParam1->GetInverseBendingMomentum()));
       // Track removed if bendingMomentum not in window [min, max]
       if ((bendingMomentum < fMinBendingMomentum) || (bendingMomentum > fMaxBendingMomentum)) {
-       track->Remove();
+       fRecTracksPtr->Remove(track);
+       fNRecTracks--;
        break; // stop the search for this candidate:
        // exit from the loop over station 
       }
-      // Track fit
+      // Track fit from parameters at first hit
+      // -> with 5 parameters (momentum and position)
       // with multiple Coulomb scattering if all stations
-      if (station == 0) track->SetFitMCS(1);
+      if (station == 0) Fit(track,1,1);
       // without multiple Coulomb scattering if not all stations
-      else track->SetFitMCS(0);
-      track->SetFitNParam(5);  // with 5 parameters (momentum and position)
-      track->SetFitStart(1);  // from parameters at first hit
-      track->Fit();
+      else Fit(track,1,0);
       Double_t numberOfDegFree = (2.0 * track->GetNTrackHits() - 5);
       if (numberOfDegFree > 0) {
         chi2Norm =  track->GetFitFMin() / numberOfDegFree;
@@ -1040,7 +1106,8 @@ void AliMUONTrackReconstructor::FollowTracks(void)
       }
       // Track removed if normalized chi2 too high
       if (chi2Norm > fMaxChi2) {
-       track->Remove();
+       fRecTracksPtr->Remove(track);
+       fNRecTracks--;
        break; // stop the search for this candidate:
        // exit from the loop over station 
       }
@@ -1052,7 +1119,7 @@ void AliMUONTrackReconstructor::FollowTracks(void)
       // Track extrapolation to the vertex through the absorber (Branson)
       // after going through the first station
       if (station == 0) {
-       trackParamVertex = *trackParam1;
+       trackParamVertex = *((AliMUONTrackParam*) (track->GetTrackParamAtHit()->First()));
        (&trackParamVertex)->ExtrapToVertex(0.,0.,0.);
        track->SetTrackParamAtVertex(&trackParamVertex);
        if (AliLog::GetGlobalDebugLevel() > 0) {
@@ -1065,12 +1132,414 @@ void AliMUONTrackReconstructor::FollowTracks(void)
     // go really to next track
     track = nextTrack;
   } // while (track)
-  // Compression of track array (necessary after Remove ????)
+  // Compression of track array (necessary after Remove)
   fRecTracksPtr->Compress();
   return;
 }
 
   //__________________________________________________________________________
+void AliMUONTrackReconstructor::Fit(AliMUONTrack *Track, Int_t FitStart, Int_t FitMCS)
+{
+  // Fit the track "Track",
+  // with or without multiple Coulomb scattering according to "FitMCS",
+  // starting, according to "FitStart",
+  // with track parameters at vertex or at the first TrackHit.
+  
+  if ((FitStart != 0) && (FitStart != 1)) {
+    cout << "ERROR in AliMUONTrackReconstructor::Fit(...)" << endl;
+    cout << "FitStart = " << FitStart << " is neither 0 nor 1" << endl;
+    exit(0);
+  }
+  if ((FitMCS != 0) && (FitMCS != 1)) {
+    cout << "ERROR in AliMUONTrackReconstructor::Fit(...)" << endl;
+    cout << "FitMCS = " << FitMCS << " is neither 0 nor 1" << endl;
+    exit(0);
+  }
+  
+  Double_t arg[1], benC, errorParam, invBenP, lower, nonBenC, upper, x, y;
+  char parName[50];
+  AliMUONTrackParam *trackParam;
+  // Check if Minuit is initialized...
+  fgFitter = TVirtualFitter::Fitter(Track,5);
+  fgFitter->Clear(); // necessary ???? probably yes
+  // how to reset the printout number at every fit ????
+  // is there any risk to leave it like that ????
+  // how to go faster ???? choice of Minuit parameters like EDM ????
+  // choice of function to be minimized according to fFitMCS
+  if (FitMCS == 0) fgFitter->SetFCN(TrackChi2);
+  else fgFitter->SetFCN(TrackChi2MCS);
+  // Switch off printout
+  arg[0] = -1;
+  fgFitter->ExecuteCommand("SET PRINT", arg, 1); // More printing !!!!
+  // No warnings
+  fgFitter->ExecuteCommand("SET NOW", arg, 0);
+  // Parameters according to "fFitStart"
+  // (should be a function to be used at every place where needed ????)
+  if (FitStart == 0) trackParam = Track->GetTrackParamAtVertex();
+  else trackParam = (AliMUONTrackParam*) (Track->GetTrackParamAtHit()->First());
+  // set first 3 Minuit parameters
+  // could be tried with no limits for the search (min=max=0) ????
+  fgFitter->SetParameter(0, "InvBenP",
+                        trackParam->GetInverseBendingMomentum(),
+                        0.003, -0.4, 0.4);
+  fgFitter->SetParameter(1, "BenS",
+                        trackParam->GetBendingSlope(),
+                        0.001, -0.5, 0.5);
+  fgFitter->SetParameter(2, "NonBenS",
+                        trackParam->GetNonBendingSlope(),
+                        0.001, -0.5, 0.5);
+  if (FitStart == 1) {
+    // set last 2 Minuit parameters when we start from first track hit
+    // mandatory limits in Bending to avoid NaN values of parameters
+    fgFitter->SetParameter(3, "X",
+                          trackParam->GetNonBendingCoor(),
+                          0.03, -500.0, 500.0);
+    // mandatory limits in non Bending to avoid NaN values of parameters
+    fgFitter->SetParameter(4, "Y",
+                          trackParam->GetBendingCoor(),
+                          0.10, -500.0, 500.0);
+  }
+  // search without gradient calculation in the function
+  fgFitter->ExecuteCommand("SET NOGRADIENT", arg, 0);
+  // minimization
+  fgFitter->ExecuteCommand("MINIMIZE", arg, 0);
+  // exit from Minuit
+  //  fgFitter->ExecuteCommand("EXIT", arg, 0); // necessary ????
+  // get results into "invBenP", "benC", "nonBenC" ("x", "y")
+  fgFitter->GetParameter(0, parName, invBenP, errorParam, lower, upper);
+  trackParam->SetInverseBendingMomentum(invBenP);
+  fgFitter->GetParameter(1, parName, benC, errorParam, lower, upper);
+  trackParam->SetBendingSlope(benC);
+  fgFitter->GetParameter(2, parName, nonBenC, errorParam, lower, upper);
+  trackParam->SetNonBendingSlope(nonBenC);
+  if (FitStart == 1) {
+    fgFitter->GetParameter(3, parName, x, errorParam, lower, upper);
+    trackParam->SetNonBendingCoor(x);
+    fgFitter->GetParameter(4, parName, y, errorParam, lower, upper);
+    trackParam->SetBendingCoor(y);
+  }
+  // global result of the fit
+  Double_t fedm, errdef, FitFMin;
+  Int_t npari, nparx;
+  fgFitter->GetStats(FitFMin, fedm, errdef, npari, nparx);
+  Track->SetFitFMin(FitFMin);
+}
+
+  //__________________________________________________________________________
+void TrackChi2(Int_t &NParam, Double_t * /*Gradient*/, Double_t &Chi2, Double_t *Param, Int_t /*Flag*/)
+{
+  // Return the "Chi2" to be minimized with Minuit for track fitting,
+  // with "NParam" parameters
+  // and their current values in array pointed to by "Param".
+  // Assumes that the track hits are sorted according to increasing Z.
+  // Track parameters at each TrackHit are updated accordingly.
+  // Multiple Coulomb scattering is not taken into account
+  AliMUONTrack *trackBeingFitted;
+  AliMUONTrackParam param1;
+  AliMUONTrackParam* TrackParamAtHit;
+  AliMUONHitForRec* HitForRec;
+  Chi2 = 0.0; // initialize Chi2
+  // copy of track parameters to be fitted
+  trackBeingFitted = (AliMUONTrack*) AliMUONTrackReconstructor::Fitter()->GetObjectFit();
+  // 3 parameters means fit track candidate from parameters at vertex (X_vertex and Y_vertex are fixed)
+  if (NParam == 3) param1 = *(trackBeingFitted->GetTrackParamAtVertex());
+  else param1 = *((AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->First()));
+  // Minuit parameters to be fitted into this copy
+  param1.SetInverseBendingMomentum(Param[0]);
+  param1.SetBendingSlope(Param[1]);
+  param1.SetNonBendingSlope(Param[2]);
+  if (NParam == 5) {
+    param1.SetNonBendingCoor(Param[3]);
+    param1.SetBendingCoor(Param[4]);
+  }
+  // Follow track through all planes of track hits
+  TrackParamAtHit = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->First());
+  while (TrackParamAtHit) {
+    HitForRec = TrackParamAtHit->GetHitForRecPtr();
+    // extrapolation to the plane of the HitForRec attached to the current TrackParamAtHit
+    param1.ExtrapToZ(HitForRec->GetZ());
+    // update track parameters of the current hit
+    TrackParamAtHit->SetTrackParam(param1);
+    // Increment Chi2
+    // done hit per hit, with hit resolution,
+    // and not with point and angle like in "reco_muon.F" !!!!
+    // Needs to add multiple scattering contribution ????
+    Double_t dX = HitForRec->GetNonBendingCoor() - param1.GetNonBendingCoor();
+    Double_t dY = HitForRec->GetBendingCoor() - param1.GetBendingCoor();
+    Chi2 = Chi2 + dX * dX / HitForRec->GetNonBendingReso2() + dY * dY / HitForRec->GetBendingReso2();
+    TrackParamAtHit = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->After(TrackParamAtHit));
+  }
+}
+
+  //__________________________________________________________________________
+void TrackChi2MCS(Int_t &NParam, Double_t * /*Gradient*/, Double_t &Chi2, Double_t *Param, Int_t /*Flag*/)
+{
+  // Return the "Chi2" to be minimized with Minuit for track fitting,
+  // with "NParam" parameters
+  // and their current values in array pointed to by "Param".
+  // Assumes that the track hits are sorted according to increasing Z.
+  // Track parameters at each TrackHit are updated accordingly.
+  // Multiple Coulomb scattering is taken into account with covariance matrix.
+  AliMUONTrack *trackBeingFitted;
+  AliMUONTrackParam param1;
+  AliMUONTrackParam* TrackParamAtHit;
+  AliMUONHitForRec* HitForRec;
+  Chi2 = 0.0; // initialize Chi2
+  // copy of track parameters to be fitted
+  trackBeingFitted = (AliMUONTrack*) AliMUONTrackReconstructor::Fitter()->GetObjectFit();
+  // 3 parameters means fit track candidate from parameters at vertex (X_vertex and Y_vertex are fixed)
+  if (NParam == 3) param1 = *(trackBeingFitted->GetTrackParamAtVertex());
+  else param1 = *((AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->First()));
+  // Minuit parameters to be fitted into this copy
+  param1.SetInverseBendingMomentum(Param[0]);
+  param1.SetBendingSlope(Param[1]);
+  param1.SetNonBendingSlope(Param[2]);
+  if (NParam == 5) {
+    param1.SetNonBendingCoor(Param[3]);
+    param1.SetBendingCoor(Param[4]);
+  }
+
+  Int_t chCurrent, chPrev = 0, hitNumber, hitNumber1, hitNumber2, hitNumber3;
+  Double_t z1, z2, z3;
+  AliMUONTrackParam *TrackParamAtHit1, *TrackParamAtHit2, *TrackParamAtHit3;
+  AliMUONHitForRec *HitForRec1, *HitForRec2;
+  Double_t hbc1, hbc2, pbc1, pbc2;
+  Double_t hnbc1, hnbc2, pnbc1, pnbc2;
+  Int_t numberOfHit = trackBeingFitted->GetNTrackHits();
+  TMatrixD *covBending = new TMatrixD(numberOfHit, numberOfHit);
+  TMatrixD *covNonBending = new TMatrixD(numberOfHit, numberOfHit);
+  Double_t *msa2 = new Double_t[numberOfHit];
+
+  // Predicted coordinates and  multiple scattering angles are first calculated
+  for (hitNumber = 0; hitNumber < numberOfHit; hitNumber++) {
+    TrackParamAtHit = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->UncheckedAt(hitNumber));
+    HitForRec = TrackParamAtHit->GetHitForRecPtr();
+    // extrapolation to the plane of the HitForRec attached to the current TrackParamAtHit
+    param1.ExtrapToZ(HitForRec->GetZ());
+    // update track parameters of the current hit
+    TrackParamAtHit->SetTrackParam(param1);
+    // square of multiple scattering angle at current hit, with one chamber
+    msa2[hitNumber] = MultipleScatteringAngle2(&param1);
+    // correction for eventual missing hits or multiple hits in a chamber,
+    // according to the number of chambers
+    // between the current hit and the previous one
+    chCurrent = HitForRec->GetChamberNumber();
+    if (hitNumber > 0) msa2[hitNumber] = msa2[hitNumber] * (chCurrent - chPrev);
+    chPrev = chCurrent;
+  }
+
+  // Calculates the covariance matrix
+  for (hitNumber1 = 0; hitNumber1 < numberOfHit; hitNumber1++) { 
+    TrackParamAtHit1 = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->UncheckedAt(hitNumber1));
+    HitForRec1 = TrackParamAtHit1->GetHitForRecPtr();
+    z1 = HitForRec1->GetZ();
+    for (hitNumber2 = hitNumber1; hitNumber2 < numberOfHit; hitNumber2++) {
+      TrackParamAtHit2 = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->UncheckedAt(hitNumber2));
+      z2 = TrackParamAtHit2->GetHitForRecPtr()->GetZ();
+      // initialization to 0 (diagonal plus upper triangular part)
+      (*covBending)(hitNumber2, hitNumber1) = 0.0;
+      // contribution from multiple scattering in bending plane:
+      // loop over upstream hits
+      for (hitNumber3 = 0; hitNumber3 < hitNumber1; hitNumber3++) {    
+        TrackParamAtHit3 = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->UncheckedAt(hitNumber3));
+       z3 = TrackParamAtHit3->GetHitForRecPtr()->GetZ();
+       (*covBending)(hitNumber2, hitNumber1) = (*covBending)(hitNumber2, hitNumber1) + ((z1 - z3) * (z2 - z3) * msa2[hitNumber3]); 
+      }
+      // equal contribution from multiple scattering in non bending plane
+      (*covNonBending)(hitNumber2, hitNumber1) = (*covBending)(hitNumber2, hitNumber1);
+      if (hitNumber1 == hitNumber2) {
+       // Diagonal elements: add contribution from position measurements
+       // in bending plane
+       (*covBending)(hitNumber2, hitNumber1) = (*covBending)(hitNumber2, hitNumber1) + HitForRec1->GetBendingReso2();
+       // and in non bending plane
+       (*covNonBending)(hitNumber2, hitNumber1) = (*covNonBending)(hitNumber2, hitNumber1) + HitForRec1->GetNonBendingReso2();
+      } else {
+       // Non diagonal elements: symmetrization
+       // for bending plane
+       (*covBending)(hitNumber1, hitNumber2) = (*covBending)(hitNumber2, hitNumber1);
+       // and non bending plane
+       (*covNonBending)(hitNumber1, hitNumber2) = (*covNonBending)(hitNumber2, hitNumber1);
+      }
+    } // for (hitNumber2 = hitNumber1;...
+  } // for (hitNumber1 = 0;...
+    
+  // Inversion of covariance matrices
+  // with "mnvertLocal", local "mnvert" function of Minuit.
+  // One cannot use directly "mnvert" since "TVirtualFitter" does not know it.
+  // One will have to replace this local function by the right inversion function
+  // from a specialized Root package for symmetric positive definite matrices,
+  // when available!!!!
+  Int_t ifailBending;
+  mnvertLocal(&((*covBending)(0,0)), numberOfHit, numberOfHit, numberOfHit, ifailBending);
+  Int_t ifailNonBending;
+  mnvertLocal(&((*covNonBending)(0,0)), numberOfHit, numberOfHit, numberOfHit, ifailNonBending);
+
+  // It would be worth trying to calculate the inverse of the covariance matrix
+  // only once per fit, since it cannot change much in principle,
+  // and it would save a lot of computing time !!!!
+  
+  // Calculates Chi2
+  if ((ifailBending == 0) && (ifailNonBending == 0)) {
+    // with Multiple Scattering if inversion correct
+    for (hitNumber1 = 0; hitNumber1 < numberOfHit ; hitNumber1++) { 
+      TrackParamAtHit1 = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->UncheckedAt(hitNumber1));
+      HitForRec1 = TrackParamAtHit1->GetHitForRecPtr();
+      hbc1 = HitForRec1->GetBendingCoor();
+      pbc1 = TrackParamAtHit1->GetBendingCoor();
+      hnbc1 = HitForRec1->GetNonBendingCoor();
+      pnbc1 = TrackParamAtHit1->GetNonBendingCoor();
+      for (hitNumber2 = 0; hitNumber2 < numberOfHit; hitNumber2++) {
+       TrackParamAtHit2 = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->UncheckedAt(hitNumber2));
+        HitForRec2 = TrackParamAtHit2->GetHitForRecPtr();
+       hbc2 = HitForRec2->GetBendingCoor();
+       pbc2 = TrackParamAtHit2->GetBendingCoor();
+       hnbc2 = HitForRec2->GetNonBendingCoor();
+       pnbc2 = TrackParamAtHit2->GetNonBendingCoor();
+       Chi2 += ((*covBending)(hitNumber2, hitNumber1) * (hbc1 - pbc1) * (hbc2 - pbc2)) +
+               ((*covNonBending)(hitNumber2, hitNumber1) * (hnbc1 - pnbc1) * (hnbc2 - pnbc2));
+      }
+    }
+  } else {
+    // without Multiple Scattering if inversion impossible
+    for (hitNumber1 = 0; hitNumber1 < numberOfHit ; hitNumber1++) { 
+      TrackParamAtHit1 = (AliMUONTrackParam*) (trackBeingFitted->GetTrackParamAtHit()->UncheckedAt(hitNumber1));
+      HitForRec1 = TrackParamAtHit1->GetHitForRecPtr();
+      hbc1 = HitForRec1->GetBendingCoor();
+      pbc1 = TrackParamAtHit1->GetBendingCoor();
+      hnbc1 = HitForRec1->GetNonBendingCoor();
+      pnbc1 = TrackParamAtHit1->GetNonBendingCoor();
+      Chi2 += ((hbc1 - pbc1) * (hbc1 - pbc1) / HitForRec1->GetBendingReso2()) +
+             ((hnbc1 - pnbc1) * (hnbc1 - pnbc1) / HitForRec1->GetNonBendingReso2());
+    }
+  }
+  
+  delete covBending;
+  delete covNonBending;
+  delete [] msa2;
+}
+
+Double_t MultipleScatteringAngle2(AliMUONTrackParam *param)
+{
+  // Returns square of multiple Coulomb scattering angle
+  // from TrackParamAtHit pointed to by "param"
+  Double_t slopeBending, slopeNonBending, radiationLength, inverseBendingMomentum2, inverseTotalMomentum2;
+  Double_t varMultipleScatteringAngle;
+  // Better implementation in AliMUONTrack class ????
+  slopeBending = param->GetBendingSlope();
+  slopeNonBending = param->GetNonBendingSlope();
+  // thickness in radiation length for the current track,
+  // taking local angle into account
+  radiationLength = AliMUONConstants::DefaultChamberThicknessInX0() *
+                   TMath::Sqrt(1.0 + slopeBending*slopeBending + slopeNonBending*slopeNonBending);
+  inverseBendingMomentum2 =  param->GetInverseBendingMomentum() * param->GetInverseBendingMomentum();
+  inverseTotalMomentum2 = inverseBendingMomentum2 * (1.0 + slopeBending * slopeBending) /
+                         (1.0 + slopeBending *slopeBending + slopeNonBending * slopeNonBending); 
+  varMultipleScatteringAngle = 0.0136 * (1.0 + 0.038 * TMath::Log(radiationLength));
+  // The velocity is assumed to be 1 !!!!
+  varMultipleScatteringAngle = inverseTotalMomentum2 * radiationLength * varMultipleScatteringAngle * varMultipleScatteringAngle;
+  return varMultipleScatteringAngle;
+}
+
+//______________________________________________________________________________
+ void mnvertLocal(Double_t *a, Int_t l, Int_t, Int_t n, Int_t &ifail)
+{
+//*-*-*-*-*-*-*-*-*-*-*-*Inverts a symmetric matrix*-*-*-*-*-*-*-*-*-*-*-*-*
+//*-*                    ==========================
+//*-*        inverts a symmetric matrix.   matrix is first scaled to
+//*-*        have all ones on the diagonal (equivalent to change of units)
+//*-*        but no pivoting is done since matrix is positive-definite.
+//*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
+
+  // taken from TMinuit package of Root (l>=n)
+  // fVERTs, fVERTq and fVERTpp changed to localVERTs, localVERTq and localVERTpp
+  //  Double_t localVERTs[n], localVERTq[n], localVERTpp[n];
+  Double_t * localVERTs = new Double_t[n];
+  Double_t * localVERTq = new Double_t[n];
+  Double_t * localVERTpp = new Double_t[n];
+  // fMaxint changed to localMaxint
+  Int_t localMaxint = n;
+
+    /* System generated locals */
+    Int_t aOffset;
+
+    /* Local variables */
+    Double_t si;
+    Int_t i, j, k, kp1, km1;
+
+    /* Parameter adjustments */
+    aOffset = l + 1;
+    a -= aOffset;
+
+    /* Function Body */
+    ifail = 0;
+    if (n < 1) goto L100;
+    if (n > localMaxint) goto L100;
+//*-*-                  scale matrix by sqrt of diag elements
+    for (i = 1; i <= n; ++i) {
+        si = a[i + i*l];
+        if (si <= 0) goto L100;
+        localVERTs[i-1] = 1 / TMath::Sqrt(si);
+    }
+    for (i = 1; i <= n; ++i) {
+        for (j = 1; j <= n; ++j) {
+            a[i + j*l] = a[i + j*l]*localVERTs[i-1]*localVERTs[j-1];
+        }
+    }
+//*-*-                                       . . . start main loop . . . .
+    for (i = 1; i <= n; ++i) {
+        k = i;
+//*-*-                  preparation for elimination step1
+        if (a[k + k*l] != 0) localVERTq[k-1] = 1 / a[k + k*l];
+        else goto L100;
+        localVERTpp[k-1] = 1;
+        a[k + k*l] = 0;
+        kp1 = k + 1;
+        km1 = k - 1;
+        if (km1 < 0) goto L100;
+        else if (km1 == 0) goto L50;
+        else               goto L40;
+L40:
+        for (j = 1; j <= km1; ++j) {
+            localVERTpp[j-1] = a[j + k*l];
+            localVERTq[j-1]  = a[j + k*l]*localVERTq[k-1];
+            a[j + k*l]   = 0;
+        }
+L50:
+        if (k - n < 0) goto L51;
+        else if (k - n == 0) goto L60;
+        else                goto L100;
+L51:
+        for (j = kp1; j <= n; ++j) {
+            localVERTpp[j-1] = a[k + j*l];
+            localVERTq[j-1]  = -a[k + j*l]*localVERTq[k-1];
+            a[k + j*l]   = 0;
+        }
+//*-*-                  elimination proper
+L60:
+        for (j = 1; j <= n; ++j) {
+            for (k = j; k <= n; ++k) { a[j + k*l] += localVERTpp[j-1]*localVERTq[k-1]; }
+        }
+    }
+//*-*-                  elements of left diagonal and unscaling
+    for (j = 1; j <= n; ++j) {
+        for (k = 1; k <= j; ++k) {
+            a[k + j*l] = a[k + j*l]*localVERTs[k-1]*localVERTs[j-1];
+            a[j + k*l] = a[k + j*l];
+        }
+    }
+    delete [] localVERTs;
+    delete [] localVERTq;
+    delete [] localVERTpp;
+    return;
+//*-*-                  failure return
+L100:
+    delete [] localVERTs;
+    delete [] localVERTq;
+    delete [] localVERTpp;
+    ifail = 1;
+} /* mnvertLocal */
+
+  //__________________________________________________________________________
 void AliMUONTrackReconstructor::RemoveDoubleTracks(void)
 {
   // To remove double tracks.
@@ -1079,78 +1548,62 @@ void AliMUONTrackReconstructor::RemoveDoubleTracks(void)
   // Among two identical tracks, one keeps the track with the larger number of hits
   // or, if these numbers are equal, the track with the minimum Chi2.
   AliMUONTrack *track1, *track2, *trackToRemove;
-  Bool_t identicalTracks;
   Int_t hitsInCommon, nHits1, nHits2;
-  identicalTracks = kTRUE;
-  while (identicalTracks) {
-    identicalTracks = kFALSE;
-    // Loop over first track of the pair
-    track1 = (AliMUONTrack*) fRecTracksPtr->First();
-    while (track1 && (!identicalTracks)) {
-      nHits1 = track1->GetNTrackHits();
-      // Loop over second track of the pair
-      track2 = (AliMUONTrack*) fRecTracksPtr->After(track1);
-      while (track2 && (!identicalTracks)) {
-       nHits2 = track2->GetNTrackHits();
-       // number of hits in common between two tracks
-       hitsInCommon = track1->HitsInCommon(track2);
-       // check for identical tracks
-       if ((4 * hitsInCommon) >= (nHits1 + nHits2)) {
-         identicalTracks = kTRUE;
-         // decide which track to remove
-         if (nHits1 > nHits2) trackToRemove = track2;
-         else if (nHits1 < nHits2) trackToRemove = track1;
-         else if ((track1->GetFitFMin()) < (track2->GetFitFMin()))
-           trackToRemove = track2;
-         else trackToRemove = track1;
-         // remove it
-         trackToRemove->Remove();
-       }
-       track2 = (AliMUONTrack*) fRecTracksPtr->After(track2);
-      } // track2
-      track1 = (AliMUONTrack*) fRecTracksPtr->After(track1);
-    } // track1
-  }
-  return;
-}
-
-  //__________________________________________________________________________
-void AliMUONTrackReconstructor::UpdateTrackParamAtHit()
-{
-  // Set track parameters after track fitting. Fill fTrackParamAtHit of AliMUONTrack's
-  AliMUONTrack *track;
-  AliMUONTrackHit *trackHit;
-  AliMUONTrackParam *trackParam;
-  track = (AliMUONTrack*) fRecTracksPtr->First();
-  while (track) {
-    trackHit = (AliMUONTrackHit*) (track->GetTrackHitsPtr())->First();
-    while (trackHit) {
-      trackParam = trackHit->GetTrackParam();
-      track->AddTrackParamAtHit(trackParam);
-      trackHit = (AliMUONTrackHit*) (track->GetTrackHitsPtr())->After(trackHit); 
-    } // trackHit    
-    track = (AliMUONTrack*) fRecTracksPtr->After(track);
-  } // track
+  Bool_t removedTrack1;
+  // Loop over first track of the pair
+  track1 = (AliMUONTrack*) fRecTracksPtr->First();
+  while (track1) {
+    removedTrack1 = kFALSE;
+    nHits1 = track1->GetNTrackHits();
+    // Loop over second track of the pair
+    track2 = (AliMUONTrack*) fRecTracksPtr->After(track1);
+    while (track2) {
+      nHits2 = track2->GetNTrackHits();
+      // number of hits in common between two tracks
+      hitsInCommon = track1->HitsInCommon(track2);
+      // check for identical tracks
+      if ((4 * hitsInCommon) >= (nHits1 + nHits2)) {
+        // decide which track to remove
+        if ((nHits1 > nHits2) || ((nHits1 == nHits2) && (track1->GetFitFMin() < track2->GetFitFMin()))) {
+         // remove track2 and continue the second loop with the track next to track2
+         trackToRemove = track2;
+         track2 = (AliMUONTrack*) fRecTracksPtr->After(track2);
+         fRecTracksPtr->Remove(trackToRemove);
+         fNRecTracks--;
+         fRecTracksPtr->Compress(); // this is essential to retrieve the TClonesArray afterwards
+        } 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);
+         fNRecTracks--;
+         fRecTracksPtr->Compress(); // this is essential to retrieve the TClonesArray afterwards
+         removedTrack1 = kTRUE;
+         break;
+        }
+      } else track2 = (AliMUONTrack*) fRecTracksPtr->After(track2);
+    } // track2
+    if (removedTrack1) continue;
+    track1 = (AliMUONTrack*) fRecTracksPtr->After(track1);
+  } // track1
   return;
 }
 
   //__________________________________________________________________________
 void AliMUONTrackReconstructor::UpdateHitForRecAtHit()
 {
-  // Set cluster parameterss after track fitting. Fill fHitForRecAtHit of AliMUONTrack's
+  // Set cluster parameters after track fitting. Fill fHitForRecAtHit of AliMUONTrack's
   AliMUONTrack *track;
-  AliMUONTrackHit *trackHit;
-  AliMUONHitForRec *hitForRec;
+  AliMUONTrackParam *trackParamAtHit;
   track = (AliMUONTrack*) fRecTracksPtr->First();
   while (track) {
-    trackHit = (AliMUONTrackHit*) (track->GetTrackHitsPtr())->First();
-    while (trackHit) {
-      hitForRec = trackHit->GetHitForRecPtr();
-      track->AddHitForRecAtHit(hitForRec);
-      trackHit = (AliMUONTrackHit*) (track->GetTrackHitsPtr())->After(trackHit); 
-    } // trackHit    
+    trackParamAtHit = (AliMUONTrackParam*) (track->GetTrackParamAtHit()->First());
+    while (trackParamAtHit) {
+      track->AddHitForRecAtHit(trackParamAtHit->GetHitForRecPtr());
+      trackParamAtHit = (AliMUONTrackParam*) (track->GetTrackParamAtHit()->After(trackParamAtHit)); 
+    }
     track = (AliMUONTrack*) fRecTracksPtr->After(track);
-  } // track
+  }
   return;
 }
 
@@ -1207,8 +1660,7 @@ void AliMUONTrackReconstructor::EventDump(void)
                     z, x, y, pX, pY, pZ, c));
 
     // track parameters at first hit
-    trackParam1 = ((AliMUONTrackHit*)
-                  (track->GetTrackHitsPtr()->First()))->GetTrackParam();
+    trackParam1 = (AliMUONTrackParam*) track->GetTrackParamAtHit()->First();
     x = trackParam1->GetNonBendingCoor();
     y = trackParam1->GetBendingCoor();
     z = trackParam1->GetZ();
index b085328..774abca 100644 (file)
@@ -24,6 +24,8 @@ class TTree;
 class AliMUONData;
 class AliRunLoader;
 class AliMUONTriggerTrack;
+class AliMUONTrack;
+class TVirtualFitter;
 
 class AliMUONTrackReconstructor : public TObject {
 
@@ -56,6 +58,7 @@ class AliMUONTrackReconstructor : public TObject {
   Double_t GetEfficiency(void) const {return fEfficiency;}
   void SetEfficiency(Double_t Efficiency) {fEfficiency = Efficiency;}
   void SetReconstructionParametersToDefaults(void);
+  static TVirtualFitter* Fitter(void) {return fgFitter;}
 
   // Hits for reconstruction
   Int_t GetNHitsForRec(void) const {return fNHitsForRec;} // Number
@@ -65,11 +68,6 @@ class AliMUONTrackReconstructor : public TObject {
   void SetNRecTracks(Int_t NRecTracks) {fNRecTracks = NRecTracks;}
   TClonesArray* GetRecTracksPtr(void) const {return fRecTracksPtr;} // Array
  
-  // Hits on reconstructed tracks
-  Int_t GetNRecTrackHits() const {return fNRecTrackHits;} // Number
-  void SetNRecTrackHits(Int_t NRecTrackHits) {fNRecTrackHits = NRecTrackHits;}
-  TClonesArray* GetRecTrackHitsPtr(void) const {return fRecTrackHitsPtr;} // Array
-
   // Functions
   Double_t GetImpactParamFromBendingMomentum(Double_t BendingMomentum) const;
   Double_t GetBendingMomentumFromImpactParam(Double_t ImpactParam) const;
@@ -90,10 +88,6 @@ class AliMUONTrackReconstructor : public TObject {
 
  private:
 
-  // Constants which should be elsewhere ????
-  static const Int_t fgkMaxMuonTrackingChambers = 10; ///< Max number of Muon tracking chambers
-  static const Int_t fgkMaxMuonTrackingStations = 5; ///< Max number of Muon tracking stations
-
   // Defaults parameters for reconstruction
   static const Double_t fgkDefaultMinBendingMomentum; ///< default min. bending momentum for reconstruction
   static const Double_t fgkDefaultMaxBendingMomentum; ///< default max. bending momentum for reconstruction
@@ -101,7 +95,6 @@ class AliMUONTrackReconstructor : public TObject {
   static const Double_t fgkDefaultMaxSigma2Distance; ///< default square of max. distance for window size 
   static const Double_t fgkDefaultBendingResolution; ///< default bending coordinate resolution for reconstruction 
   static const Double_t fgkDefaultNonBendingResolution; ///< default non bending coordinate resolution for reconstruction
-  static const Double_t fgkDefaultChamberThicknessInX0; ///< default chamber thickness in X0 for reconstruction
   // Simple magnetic field:
   // Value taken from macro MUONtracking.C: 0.7 T, hence 7 kG
   // Length and Position from reco_muon.F, with opposite sign:
@@ -113,6 +106,7 @@ class AliMUONTrackReconstructor : public TObject {
   static const Double_t fgkDefaultSimpleBPosition; ///< default position of magnetic field (dipole)
   static const Double_t fgkDefaultEfficiency; ///< default chamber efficiency for track ref. hits recontruction
 
+  static TVirtualFitter* fgFitter; //!< Pointer to track fitter
 
   Int_t fTrackMethod; ///< AZ - tracking method
 
@@ -122,10 +116,10 @@ class AliMUONTrackReconstructor : public TObject {
   Double_t fMaxBendingMomentum; ///< maximum value (GeV/c) of momentum in bending plane
   Double_t fMaxChi2; ///< maximum Chi2 per degree of Freedom
   Double_t fMaxSigma2Distance; ///< maximum square distance in units of the variance (maximum chi2)
-  Double_t fRMin[fgkMaxMuonTrackingChambers]; ///< minimum radius (cm)
-  Double_t fRMax[fgkMaxMuonTrackingChambers]; ///< maximum radius (cm)
-  Double_t fSegmentMaxDistBending[fgkMaxMuonTrackingStations]; ///< maximum distance (cm) for segments in bending plane
-  Double_t fSegmentMaxDistNonBending[fgkMaxMuonTrackingStations]; ///< maximum distance (cm) for segments in non bending plane
+  Double_t* fRMin; ///< minimum radius (cm)
+  Double_t* fRMax; ///< maximum radius (cm)
+  Double_t* fSegmentMaxDistBending; ///< maximum distance (cm) for segments in bending plane
+  Double_t* fSegmentMaxDistNonBending; ///< maximum distance (cm) for segments in non bending plane
   Double_t fBendingResolution; ///< chamber resolution (cm) in bending plane
   Double_t fNonBendingResolution; ///< chamber resolution (cm) in non bending plane
   Double_t fChamberThicknessInX0; ///< chamber thickness in number of radiation lengths
@@ -136,24 +130,20 @@ class AliMUONTrackReconstructor : public TObject {
   Double_t fEfficiency; ///< chamber efficiency (used for track ref. hits only)
   
   // Hits for reconstruction (should be in AliMUON ????)
-  TClonesArray *fHitsForRecPtr; ///< pointer to the array of hits for reconstruction
+  TClonesArray* fHitsForRecPtr; ///< pointer to the array of hits for reconstruction
   Int_t fNHitsForRec; ///< number of hits for reconstruction
   // Information per chamber (should be in AliMUONChamber ????)
-  Int_t fNHitsForRecPerChamber[fgkMaxMuonTrackingChambers]; ///< number of HitsForRec
-  Int_t fIndexOfFirstHitForRecPerChamber[fgkMaxMuonTrackingChambers]; ///< index (0...) of first HitForRec
+  Int_t* fNHitsForRecPerChamber; ///< number of HitsForRec
+  Int_t* fIndexOfFirstHitForRecPerChamber; ///< index (0...) of first HitForRec
 
   // Segments inside a station
-  TClonesArray *fSegmentsPtr[fgkMaxMuonTrackingStations]; ///< array of pointers to the segments for each station
-  Int_t fNSegments[fgkMaxMuonTrackingStations]; ///< number of segments for each station
+  TClonesArray** fSegmentsPtr; ///< array of pointers to the segments for each station
+  Int_t* fNSegments; ///< number of segments for each station
 
   // Reconstructed tracks
   TClonesArray *fRecTracksPtr; ///< pointer to array of reconstructed tracks
   Int_t fNRecTracks; ///< number of reconstructed tracks
 
-  // Track hits on reconstructed tracks
-  TClonesArray *fRecTrackHitsPtr; ///< pointer to array of hits on reconstructed tracks
-  Int_t fNRecTrackHits; ///< number of hits on reconstructed tracks
-
   // data container
   AliMUONData* fMUONData; ///< Data container for MUON subsystem 
 
@@ -175,14 +165,14 @@ class AliMUONTrackReconstructor : public TObject {
   void MakeSegmentsPerStation(Int_t Station);
   void MakeTracks(void);
   Bool_t MakeTriggerTracks(void);
-  void ResetTrackHits(void);
   void ResetTracks(void);
+  void MakeTrackCandidates(void);
   Int_t MakeTrackCandidatesWithTwoSegments(AliMUONSegment *BegSegment);
   Int_t MakeTrackCandidatesWithOneSegmentAndOnePoint(AliMUONSegment *BegSegment);
-  void MakeTrackCandidates(void);
+  void CalcTrackParamAtVertex(AliMUONTrack *Track);
   void FollowTracks(void);
+  void Fit(AliMUONTrack *Track, Int_t FitStart, Int_t FitMCS);
   void RemoveDoubleTracks(void);
-  void UpdateTrackParamAtHit(void);
   void UpdateHitForRecAtHit(void);
   void ValidateTracksWithTrigger(void);
 
index ef8b70f..36679d9 100755 (executable)
@@ -41,7 +41,10 @@ MuonRec.SetRunTracking("");
 MuonRec.SetRunVertexFinder(kFALSE);
 MuonRec.SetRunLocalReconstruction("MUON");
 MuonRec.SetFillESD("MUON");
-MuonRec.SetLoadAlignData("MUON")
+MuonRec.SetLoadAlignData("MUON");
+// Uncoment following line to run reconstruction with the orginal tracking method
+// instead of the kalman one (default)
+// MuonRec.SetOption("MUON","Original");
 MuonRec.Run();
 .q
 EOF
index 5954aa8..128daca 100644 (file)
@@ -18,7 +18,6 @@
 #pragma link C++ class AliMUONTriggerTrack+; 
 #pragma link C++ class AliMUONRecoTrack+; 
 #pragma link C++ class AliMUONHitForRec+; 
-#pragma link C++ class AliMUONTrackHit+; 
 #pragma link C++ class AliMUONSegment+;
 #pragma link C++ class AliMUONClusterDrawAZ+;
 #pragma link C++ class AliMUONDetElement+; 
index 8c3c965..2520ff2 100644 (file)
@@ -14,7 +14,6 @@ SRCS:= AliMUONClusterReconstructor.cxx \
        AliMUONTriggerTrack.cxx \
        AliMUONRecoTrack.cxx \
        AliMUONHitForRec.cxx \
-       AliMUONTrackHit.cxx \
        AliMUONSegment.cxx \
        AliMUONClusterDrawAZ.cxx \
        AliMUONDetElement.cxx \