major TRD reconstruction update
authorabercuci <abercuci@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 8 Apr 2009 14:12:58 +0000 (14:12 +0000)
committerabercuci <abercuci@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 8 Apr 2009 14:12:58 +0000 (14:12 +0000)
- update clusterizer for flexible switch between methods to calculate
  r-phi cluster position. Currently the following methods are supported:
  - COG - center of gravity
  - LUT - look up table [default]
  - GAUS - gauss approximation [testing phase - to production]
- move position calculation related code from clusterizer to the cluster
  itself
- steering options for cluster position calculation implemented at Reconstructor level.
  - user help via static method AliTRDReconstructor::Options() function
- consistent pad row cross implementation in the tracklet. Reference
  radial position for pads crossing pad row is recalculated based on z
  cluster charge distribution.
- miscellaneous

17 files changed:
TRD/AliTRDReconstructor.cxx
TRD/AliTRDReconstructor.h
TRD/AliTRDarrayADC.h
TRD/AliTRDcluster.cxx
TRD/AliTRDcluster.h
TRD/AliTRDclusterizer.cxx
TRD/AliTRDclusterizer.h
TRD/AliTRDrecoParam.cxx
TRD/AliTRDrecoParam.h
TRD/AliTRDseedV1.cxx
TRD/AliTRDseedV1.h
TRD/AliTRDtrackV1.cxx
TRD/AliTRDtrackV1.h
TRD/AliTRDtrackerV1.cxx
TRD/AliTRDtransform.cxx
TRD/AliTRDtransform.h
TRD/qaRec/AliTRDresolution.cxx

index 738916c..c25b6a3 100644 (file)
 ClassImp(AliTRDReconstructor)
 
 TClonesArray *AliTRDReconstructor::fgClusters = 0x0;
+Char_t* AliTRDReconstructor::fgSteerNames[kNsteer] = {
+  "DigitsConversion       "
+ ,"Tail Cancellation      "
+ ,"Clusters LUT           "
+ ,"Clusters GAUSS         "
+ ,"Clusters Sharing       "
+ ,"NN PID                 "
+ ,"8 dEdx slices in ESD   "
+ ,"Write Clusters         "
+ ,"Write Online Tracklets "
+ ,"Drift Gas Argon        "
+ ,"Stand Alone Tracking   "
+ ,"Vertex Constrain       "
+ ,"Tracklet Improve       "
+ ,"HLT Mode              "
+ ,"Cosmic Reconstruction "
+};
+Char_t* AliTRDReconstructor::fgSteerFlags[kNsteer] = {
+  "dc"// digits conversion [false]
+ ,"tc"// apply tail cancellation [true]
+ ,"lt"// look-up-table for cluster shape in the r-phi direction
+ ,"gs"// gauss cluster shape in the r-phi direction
+ ,"sh"// cluster sharing between tracks
+ ,"nn"// PID method in reconstruction (NN) [true]
+ ,"8s"// 8 dEdx slices in ESD [true] 
+ ,"cw"// write clusters [true]
+ ,"tw"// write online tracklets [false]
+ ,"ar"// drift gas [false] - do not update the number of exponentials in the TC !
+ ,"sa"// track seeding (stand alone tracking) [true]
+ ,"vc"// vertex constrain on stand alone track finder [false]
+ ,"ti"// improve tracklets in stand alone track finder [true]
+ ,"hlt"// HLT reconstruction [false]
+ ,"cos"// Cosmic Reconstruction [false]
+};
+Char_t* AliTRDReconstructor::fgTaskNames[kNtasks] = {
+  "RawReader"
+ ,"Clusterizer"
+ ,"Tracker"
+ ,"PID"
+};
+Char_t* AliTRDReconstructor::fgTaskFlags[kNtasks] = {
+  "rr"
+ ,"cl"
+ ,"tr"
+ ,"pd"
+};
+
 //_____________________________________________________________________________
 AliTRDReconstructor::AliTRDReconstructor()
   :AliReconstructor()
@@ -65,7 +112,15 @@ AliTRDReconstructor::AliTRDReconstructor()
   // PID method in reconstruction (NN) [nn]
   SETFLG(fSteerParam, kSteerPID);
   // number of dEdx slices in the ESD track [8s]
-  //SETFLG(fSteerParam, kEightSlices);
+  SETFLG(fSteerParam, kEightSlices);
+  // vertex constrain for stand alone track finder
+  SETFLG(fSteerParam, kVertexConstrained);
+  // improve tracklets for stand alone track finder
+  SETFLG(fSteerParam, kImproveTracklet);
+  // use look up table for cluster r-phi position
+  SETFLG(fSteerParam, kLUT);
+  // use tail cancellation
+  SETFLG(fSteerParam, kTC);
 
   memset(fStreamLevel, 0, kNtasks*sizeof(UChar_t));
   memset(fDebugStream, 0, sizeof(TTreeSRedirector *) * kNtasks);
@@ -112,17 +167,7 @@ void AliTRDReconstructor::Init(){
   // Init Options
   //
   SetOption(GetOption());
-
-  AliInfo(Form("\tDigitsConversion       [dc] : %s", fSteerParam&kDigitsConversion?"yes":"no"));
-  AliInfo(Form("\tWrite Clusters         [cw] : %s", fSteerParam&kWriteClusters?"yes":"no"));
-  AliInfo(Form("\tWrite Online Tracklets [tw] : %s", fSteerParam&kWriteTracklets?"yes":"no"));
-  AliInfo(Form("\tDrift Gas Argon        [ar] : %s", fSteerParam&kDriftGas?"yes":"no"));
-  AliInfo(Form("\tStand Alone Tracking   [sa] : %s", fSteerParam&kSeeding?"yes":"no"));
-  AliInfo(Form("\tHLT         Tracking  [hlt] : %s", fSteerParam&kHLT?"yes":"no"));
-  AliInfo(Form("\tCosmic Reconstruction [cos] : %s", fSteerParam&kCosmic?"yes":"no"));
-  AliInfo(Form("\tNN PID                 [nn] : %s", fSteerParam&kSteerPID?"yes":"no"));
-  AliInfo(Form("\t8 dEdx slices in ESD   [8s] : %s", fSteerParam&kEightSlices?"yes":"no"));
-  AliInfo(Form("\tStreaming Levels            : Clusterizer[%d] Tracker[%d] PID[%d]", fStreamLevel[kClusterizer], fStreamLevel[kTracker], fStreamLevel[kPID]));
+  Options(fSteerParam, fStreamLevel);
 }
 
 //_____________________________________________________________________________
@@ -161,10 +206,10 @@ void AliTRDReconstructor::Reconstruct(AliRawReader *rawReader
   rawReader->Select("TRD");
 
   // New (fast) cluster finder
-  AliTRDclusterizer clusterer("clusterer","TRD clusterizer");
+  AliTRDclusterizer clusterer(fgTaskNames[kClusterizer], fgTaskNames[kClusterizer]);
   clusterer.SetReconstructor(this);
   clusterer.OpenOutput(clusterTree);
-  clusterer.SetAddLabels(kFALSE);
+  clusterer.SetUseLabels(kFALSE);
   clusterer.Raw2ClustersChamber(rawReader);
   
   if(IsWritingClusters()) return;
@@ -184,7 +229,7 @@ void AliTRDReconstructor::Reconstruct(TTree *digitsTree
 
   //AliInfo("Reconstruct TRD clusters from Digits [Digit TTree -> Cluster TTree]");
 
-  AliTRDclusterizer clusterer("clusterer","TRD clusterizer");
+  AliTRDclusterizer clusterer(fgTaskNames[kClusterizer], fgTaskNames[kClusterizer]);
   clusterer.SetReconstructor(this);
   clusterer.OpenOutput(clusterTree);
   clusterer.ReadDigits(digitsTree);
@@ -228,75 +273,41 @@ void AliTRDReconstructor::SetOption(Option_t *opt)
 {
 // Read option string into the steer param.
 //
-// Default steer param values
-//
-// digits conversion [dc] = false
-// drift gas [ar] = false - do not update the number of exponentials in the TC !
-// write clusters [cw] = true
-// write online tracklets [tw] = false
-// track seeding (stand alone tracking) [sa] = true
-// PID method in reconstruction (NN) [nn] = true
-// 8 dEdx slices in ESD [8s] = false 
-// HLT tracking [hlt] = false
-// Cosmic Reconstruction [cos] = false
-//
 
   AliReconstructor::SetOption(opt);
 
   TString s(opt);
   TObjArray *opar = s.Tokenize(",");
   for(Int_t ipar=0; ipar<opar->GetEntriesFast(); ipar++){
+    Bool_t PROCESSED = kFALSE;
     TString sopt(((TObjString*)(*opar)[ipar])->String());
-    if(sopt.Contains("dc")){
-      SETFLG(fSteerParam, kDigitsConversion);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kDigitsConversion);
-      continue;        
-    } else if(sopt.Contains("cw")){ 
-      SETFLG(fSteerParam, kWriteClusters);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kWriteClusters);
-      continue;
-    } else if(sopt.Contains("sa")){
-      SETFLG(fSteerParam, kSeeding);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kSeeding);
-      continue;
-    } else if(sopt.Contains("nn")){
-      SETFLG(fSteerParam, kSteerPID);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kSteerPID);
-      continue;
-    } else if(sopt.Contains("8s")){
-      SETFLG(fSteerParam, kEightSlices);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kEightSlices);
-      continue;
-    } else if(sopt.Contains("tw")){
-      SETFLG(fSteerParam, kWriteTracklets);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kWriteTracklets);
-      continue;        
-    } else if(sopt.Contains("ar")){
-      SETFLG(fSteerParam, kDriftGas);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kDriftGas);
-      continue;        
-    } else if(sopt.Contains("hlt")){
-      SETFLG(fSteerParam, kHLT);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kHLT);
-      continue;        
-    } else if(sopt.Contains("cos")){
-      SETFLG(fSteerParam, kCosmic);
-      if(sopt.Contains("!")) CLRFLG(fSteerParam, kCosmic);
-    } else if(sopt.Contains("sl")){
+    for(Int_t iopt=0; iopt<kNsteer; iopt++){
+      if(!sopt.Contains(fgSteerFlags[iopt])) continue;
+      SETFLG(fSteerParam, BIT(iopt));
+      if(sopt.Contains("!")) CLRFLG(fSteerParam, BIT(iopt));
+      PROCESSED = kTRUE;
+      break;   
+    }
+    if(PROCESSED) continue;
+
+    if(sopt.Contains("sl")){
       TObjArray *stl = sopt.Tokenize("_");
       if(stl->GetEntriesFast() < 3) continue;
       TString taskstr(((TObjString*)(*stl)[1])->String());
       TString levelstring(((TObjString*)(*stl)[2])->String());
-      // Set the stream Level
       Int_t level = levelstring.Atoi();
-      ETRDReconstructorTask task = kTracker;
-      if(taskstr.CompareTo("raw") == 0) task = kRawReader;     
-      else if(taskstr.CompareTo("cl") == 0) task = kClusterizer;       
-      else if(taskstr.CompareTo("tr") == 0) task = kTracker;
-      else if(taskstr.CompareTo("pi") == 0) task = kPID;
-      SetStreamLevel(level, task);
-      continue;
-    }
+
+      // Set the stream Level
+      PROCESSED = kFALSE;
+      for(Int_t it=0; it<kNtasks; it++){
+        if(taskstr.CompareTo(fgTaskFlags[it]) != 0) continue;
+        SetStreamLevel(level, ETRDReconstructorTask(it));
+        PROCESSED = kTRUE;
+      }
+    } 
+    if(PROCESSED) continue;
+
+    AliWarning(Form("Unknown option flag %s.", sopt.Data()));
   }
 }
 
@@ -305,15 +316,26 @@ void AliTRDReconstructor::SetStreamLevel(Int_t level, ETRDReconstructorTask task
   //
   // Set the Stream Level for one of the tasks Clusterizer, Tracker or PID
   //
-  TString taskname[4] = {"RawReader", "Clusterizer", "Tracker", "PID"};
   const Int_t minLevel[4] = {1, 1, 2, 1}; // the minimum debug level upon which a debug stream is created for different tasks
   //AliInfo(Form("Setting Stream Level for Task %s to %d", taskname.Data(),level));
   fStreamLevel[(Int_t)task] = level;
   // Initialize DebugStreamer if not yet done
   if(level >= minLevel[task] && !fDebugStream[task]){
     TDirectory *savedir = gDirectory;
-    fDebugStream[task] = new TTreeSRedirector(Form("TRD.%sDebug.root", taskname[task].Data()));
+    fDebugStream[task] = new TTreeSRedirector(Form("TRD.Debug%s.root", fgTaskNames[task]));
     savedir->cd();
     SETFLG(fSteerParam, kOwner);
   }
 }
+
+//_____________________________________________________________________________
+void AliTRDReconstructor::Options(UInt_t steer, UChar_t *stream)
+{
+  for(Int_t iopt=0; iopt<kNsteer; iopt++){
+    AliInfoGeneral("AliTRDReconstructor", Form(" %s[%s]%s", fgSteerNames[iopt], fgSteerFlags[iopt], steer ?(((steer>>iopt)&1)?" : ON":" : OFF"):""));
+  }
+  AliInfoGeneral("AliTRDReconstructor", " Debug Streaming"); 
+  for(Int_t it=0; it<kNtasks; it++) 
+    AliInfoGeneral("AliTRDReconstructor", Form(" %s [sl_%s] %d", fgTaskNames[it], fgTaskFlags[it], stream ? stream[it] : 0));
+}
+
index 61d6e35..69164dd 100644 (file)
@@ -24,15 +24,22 @@ class AliTRDReconstructor: public AliReconstructor
 public:
   enum ETRDReconstructorSteer {
     kDigitsConversion= BIT(0)
-    ,kWriteClusters  = BIT(1)
-    ,kSeeding        = BIT(2)
-    ,kSteerPID       = BIT(3)
-    ,kEightSlices    = BIT(4)
-    ,kWriteTracklets = BIT(5)
-    ,kDriftGas       = BIT(6)
-    ,kHLT            = BIT(7)
-    ,kCosmic         = BIT(8)
-    ,kOwner          = BIT(14)
+    ,kTC             = BIT(1) // tail cancelation
+    ,kLUT            = BIT(2) // look up table for cluster position determination 
+    ,kGAUS           = BIT(3) // look up table for cluster position determination 
+    ,kClusterSharing = BIT(4) // Toggle cluster sharing
+    ,kSteerPID       = BIT(5)
+    ,kEightSlices    = BIT(6)
+    ,kWriteClusters  = BIT(7)
+    ,kWriteTracklets = BIT(8)
+    ,kDriftGas       = BIT(9)
+    ,kSeeding        = BIT(10)
+    ,kVertexConstrained = BIT(11) // Perform vertex constrained fit
+    ,kImproveTracklet   = BIT(12) // Improve tracklet in the SA TRD track finder 
+    ,kHLT            = BIT(13)
+    ,kCosmic         = BIT(14)
+    ,kOwner          = BIT(15)
+    ,kNsteer         = 15       // number of tasks
   };
   enum ETRDReconstructorTask {
     kRawReader    = 0
@@ -66,14 +73,21 @@ public:
   static const AliTRDrecoParam* GetRecoParam() { return dynamic_cast<const AliTRDrecoParam*>(AliReconstructor::GetRecoParam(2)); }
   Int_t               GetStreamLevel(ETRDReconstructorTask task) const    { return fStreamLevel[task];} 
   inline void         GetTCParams(Double_t *par) const;
-  virtual Bool_t      HasDigitConversion() const                   { return fSteerParam&kDigitsConversion;  };
+  virtual Bool_t      HasDigitConversion() const { return fSteerParam&kDigitsConversion;  };
+  Bool_t              HasVertexConstrained() const { return fSteerParam&kVertexConstrained; }
+  Bool_t              HasImproveTracklets() const  { return fSteerParam&kImproveTracklet; }
   Bool_t              IsWritingClusters() const  { return fSteerParam&kWriteClusters;}
   Bool_t              IsWritingTracklets() const { return fSteerParam&kWriteTracklets;}
   Bool_t              IsHLT() const              { return fSteerParam&kHLT;}
   Bool_t              IsSeeding() const          { return fSteerParam&kSeeding;}
   Bool_t              IsCosmic() const           { return fSteerParam&kCosmic;}
   Bool_t              IsEightSlices() const      { return fSteerParam&kEightSlices;}
-
+  Bool_t              UseClusterSharing() const  { return fSteerParam&kClusterSharing;}
+  Bool_t              UseLUT() const             { return fSteerParam&kLUT;}
+  Bool_t              UseGAUS() const             { return fSteerParam&kGAUS;}
+  Bool_t              UseTailCancelation() const { return fSteerParam&kTC;}
+  
+  static void         Options(UInt_t steer=0, UChar_t *stream=0x0);
   virtual void        Reconstruct(AliRawReader *rawReader, TTree *clusterTree) const;
   virtual void        Reconstruct(TTree *digitsTree, TTree *clusterTree) const;
 
@@ -83,14 +97,18 @@ public:
   void                SetStreamLevel(Int_t level, ETRDReconstructorTask task= kTracker);
 
 private:
-  UChar_t           fStreamLevel[kNtasks];// stream level for each reconstruction task         
-  UInt_t            fSteerParam;          // steering flags
+  static Char_t    *fgSteerNames[kNsteer];//! steering names
+  static Char_t    *fgSteerFlags[kNsteer];//! steering flags
+  static Char_t    *fgTaskNames[kNtasks]; //! tasks names
+  static Char_t    *fgTaskFlags[kNtasks]; //! tasks flags
+  UChar_t           fStreamLevel[kNtasks];// stream level for each reconstruction task
+  UInt_t            fSteerParam;          // steering bits
   Double_t          fTCParams[8];         // Tail Cancellation parameters for drift gases 
   TTreeSRedirector *fDebugStream[kNtasks];// Debug Streamer container;
  
   static TClonesArray *fgClusters;    // list of clusters for local reconstructor
 
-  ClassDef(AliTRDReconstructor, 1)    //  Class for the TRD reconstruction
+  ClassDef(AliTRDReconstructor, 2)    //  Class for the TRD reconstruction
 
 };
 
index c31c6be..67e9d70 100644 (file)
@@ -31,6 +31,7 @@ class AliTRDarrayADC: public TObject
   Bool_t  HasData() const {return fNtime ? 1 : 0;};
   Short_t GetData(Int_t nrow, Int_t ncol, Int_t ntime) const
                        {return fADC[(nrow*fNcol+ncol)*fNtime+ntime];};
+  inline void GetData(Int_t r, Int_t c, Int_t t, Int_t n, Short_t *vals) const;
   Short_t GetDataB(Int_t nrow, Int_t ncol, Int_t ntime) const;
   UChar_t GetPadStatus(Int_t nrow, Int_t ncol, Int_t ntime) const;
   void    SetPadStatus(Int_t nrow, Int_t ncol, Int_t ntime, UChar_t status);
@@ -56,4 +57,10 @@ class AliTRDarrayADC: public TObject
   ClassDef(AliTRDarrayADC,1) //ADC container class
     
 };
+
+inline void AliTRDarrayADC::GetData(Int_t r, Int_t c, Int_t t, Int_t n, Short_t *vals) const
+{
+  for(Int_t ic=n, idx = (r*fNcol+c)*fNtime+t; ic--; idx+=fNtime) vals[ic] = fADC[idx];
+}
+
 #endif 
index a29a099..036f779 100644 (file)
@@ -31,6 +31,9 @@
 
 ClassImp(AliTRDcluster)
 
+const Int_t AliTRDcluster::fgkNlut = 128;
+Double_t *AliTRDcluster::fgLUT = 0x0;
+
 //___________________________________________________________________________
 AliTRDcluster::AliTRDcluster() 
   :AliCluster() 
@@ -55,6 +58,30 @@ AliTRDcluster::AliTRDcluster()
 }
 
 //___________________________________________________________________________
+AliTRDcluster::AliTRDcluster(Int_t det, UChar_t col, UChar_t row, UChar_t time, const Short_t *sig, UShort_t vid) 
+  :AliCluster() 
+  ,fPadCol(col)
+  ,fPadRow(row)
+  ,fPadTime(time)
+  ,fLocalTimeBin(0)
+  ,fNPads(0)
+  ,fClusterMasking(0)
+  ,fDetector(det)
+  ,fQ(0.)
+  ,fCenter(0.)
+{ 
+  //
+  // Constructor for self constructing cluster. In this approach the information is inserted gradualy into the 
+  // cluster and all dependencies are (re)calculated inside the cluster itself.
+  //
+  // A.Bercuci <A.Bercuci@gsi.de>
+
+  memcpy(&fSignals, sig, 7*sizeof(Short_t));
+  fQ = fSignals[2]+fSignals[3]+fSignals[4];
+  SetVolumeId(vid);
+}
+
+//___________________________________________________________________________
 AliTRDcluster::AliTRDcluster(Int_t det, Float_t q
                            , Float_t *pos, Float_t *sig
                            , Int_t *tracks, Char_t npads, Short_t *signals
@@ -239,26 +266,150 @@ Float_t AliTRDcluster::GetSumS() const
 }
 
 //___________________________________________________________________________
-Float_t AliTRDcluster::GetXcorr(Int_t tb)
+Double_t AliTRDcluster::GetSX(Int_t tb, Double_t z)
+{
+  if(tb<1 || tb>=24) return 10.; // return huge [10cm]
+  const Double_t sx[24][10]={
+    {0.000e+00, 9.352e-01, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 2.309e+00},
+    {8.387e-02, 8.718e-02, 8.816e-02, 9.444e-02, 9.993e-02, 1.083e-01, 1.161e-01, 1.280e-01, 1.417e-01, 1.406e-01},
+    {1.097e-01, 1.105e-01, 1.127e-01, 1.151e-01, 1.186e-01, 1.223e-01, 1.272e-01, 1.323e-01, 1.389e-01, 1.490e-01},
+    {1.407e-01, 1.404e-01, 1.414e-01, 1.430e-01, 1.429e-01, 1.449e-01, 1.476e-01, 1.494e-01, 1.515e-01, 1.589e-01},
+    {1.681e-01, 1.679e-01, 1.666e-01, 1.657e-01, 1.656e-01, 1.649e-01, 1.652e-01, 1.662e-01, 1.671e-01, 1.694e-01},
+    {1.745e-01, 1.737e-01, 1.707e-01, 1.690e-01, 1.643e-01, 1.610e-01, 1.612e-01, 1.628e-01, 1.638e-01, 1.659e-01},
+    {1.583e-01, 1.558e-01, 1.535e-01, 1.488e-01, 1.445e-01, 1.419e-01, 1.428e-01, 1.451e-01, 1.462e-01, 1.494e-01},
+    {1.414e-01, 1.391e-01, 1.368e-01, 1.300e-01, 1.256e-01, 1.259e-01, 1.285e-01, 1.326e-01, 1.358e-01, 1.406e-01},
+    {1.307e-01, 1.289e-01, 1.261e-01, 1.216e-01, 1.193e-01, 1.165e-01, 1.201e-01, 1.241e-01, 1.274e-01, 1.344e-01},
+    {1.251e-01, 1.227e-01, 1.208e-01, 1.155e-01, 1.110e-01, 1.116e-01, 1.133e-01, 1.187e-01, 1.229e-01, 1.308e-01},
+    {1.234e-01, 1.209e-01, 1.175e-01, 1.127e-01, 1.094e-01, 1.093e-01, 1.109e-01, 1.155e-01, 1.210e-01, 1.275e-01},
+    {1.215e-01, 1.187e-01, 1.156e-01, 1.108e-01, 1.070e-01, 1.065e-01, 1.090e-01, 1.134e-01, 1.196e-01, 1.251e-01},
+    {1.202e-01, 1.180e-01, 1.151e-01, 1.108e-01, 1.070e-01, 1.058e-01, 1.089e-01, 1.127e-01, 1.183e-01, 1.256e-01},
+    {1.207e-01, 1.176e-01, 1.142e-01, 1.109e-01, 1.072e-01, 1.069e-01, 1.088e-01, 1.122e-01, 1.182e-01, 1.252e-01},
+    {1.213e-01, 1.182e-01, 1.156e-01, 1.102e-01, 1.076e-01, 1.063e-01, 1.091e-01, 1.132e-01, 1.181e-01, 1.243e-01},
+    {1.205e-01, 1.180e-01, 1.150e-01, 1.104e-01, 1.072e-01, 1.063e-01, 1.083e-01, 1.132e-01, 1.183e-01, 1.243e-01},
+    {1.212e-01, 1.195e-01, 1.135e-01, 1.107e-01, 1.070e-01, 1.065e-01, 1.097e-01, 1.126e-01, 1.185e-01, 1.238e-01},
+    {1.201e-01, 1.184e-01, 1.155e-01, 1.111e-01, 1.088e-01, 1.075e-01, 1.089e-01, 1.131e-01, 1.189e-01, 1.237e-01},
+    {1.197e-01, 1.186e-01, 1.147e-01, 1.113e-01, 1.085e-01, 1.077e-01, 1.105e-01, 1.137e-01, 1.188e-01, 1.245e-01},
+    {1.213e-01, 1.194e-01, 1.154e-01, 1.114e-01, 1.091e-01, 1.082e-01, 1.098e-01, 1.140e-01, 1.194e-01, 1.247e-01},
+    {1.210e-01, 1.189e-01, 1.155e-01, 1.119e-01, 1.088e-01, 1.080e-01, 1.105e-01, 1.141e-01, 1.195e-01, 1.244e-01},
+    {1.196e-01, 1.189e-01, 1.145e-01, 1.105e-01, 1.095e-01, 1.083e-01, 1.087e-01, 1.121e-01, 1.173e-01, 1.208e-01},
+    {1.123e-01, 1.129e-01, 1.108e-01, 1.110e-01, 1.080e-01, 1.065e-01, 1.056e-01, 1.066e-01, 1.071e-01, 1.095e-01},
+    {1.136e-01, 1.135e-01, 1.130e-01, 1.122e-01, 1.113e-01, 1.071e-01, 1.041e-01, 1.025e-01, 1.014e-01, 9.973e-02}
+  };
+  if(z>=0. && z<.25) return sx[tb][Int_t(z/.025)];
+  
+  Double_t m = 1.e-8; for(Int_t id=10; id--;) if(sx[tb][id]>m) m=sx[tb][id];
+  return m;
+}
+
+//___________________________________________________________________________
+Double_t AliTRDcluster::GetSY(Int_t tb, Double_t z)
+{
+  if(tb<1 || tb>=24) return 10.; // return huge [10cm]
+  const Double_t sy[24][10]={
+    {0.000e+00, 2.610e-01, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 4.680e-01},
+    {3.019e-02, 3.036e-02, 3.131e-02, 3.203e-02, 3.294e-02, 3.407e-02, 3.555e-02, 3.682e-02, 3.766e-02, 3.824e-02},
+    {1.773e-02, 1.778e-02, 1.772e-02, 1.790e-02, 1.807e-02, 1.833e-02, 1.873e-02, 1.905e-02, 1.958e-02, 2.029e-02},
+    {1.774e-02, 1.772e-02, 1.746e-02, 1.738e-02, 1.756e-02, 1.756e-02, 1.739e-02, 1.720e-02, 1.743e-02, 1.769e-02},
+    {2.064e-02, 2.078e-02, 2.069e-02, 2.060e-02, 2.033e-02, 2.024e-02, 2.022e-02, 1.961e-02, 1.922e-02, 1.901e-02},
+    {2.382e-02, 2.379e-02, 2.371e-02, 2.333e-02, 2.318e-02, 2.285e-02, 2.255e-02, 2.244e-02, 2.174e-02, 2.132e-02},
+    {2.615e-02, 2.589e-02, 2.539e-02, 2.493e-02, 2.420e-02, 2.396e-02, 2.362e-02, 2.342e-02, 2.321e-02, 2.330e-02},
+    {2.640e-02, 2.638e-02, 2.577e-02, 2.548e-02, 2.477e-02, 2.436e-02, 2.416e-02, 2.401e-02, 2.399e-02, 2.402e-02},
+    {2.647e-02, 2.632e-02, 2.587e-02, 2.546e-02, 2.465e-02, 2.447e-02, 2.429e-02, 2.415e-02, 2.429e-02, 2.475e-02},
+    {2.657e-02, 2.637e-02, 2.580e-02, 2.525e-02, 2.492e-02, 2.441e-02, 2.446e-02, 2.441e-02, 2.478e-02, 2.491e-02},
+    {2.640e-02, 2.608e-02, 2.583e-02, 2.539e-02, 2.478e-02, 2.440e-02, 2.456e-02, 2.464e-02, 2.486e-02, 2.533e-02},
+    {2.636e-02, 2.630e-02, 2.584e-02, 2.542e-02, 2.483e-02, 2.451e-02, 2.449e-02, 2.467e-02, 2.496e-02, 2.554e-02},
+    {2.634e-02, 2.629e-02, 2.583e-02, 2.526e-02, 2.480e-02, 2.460e-02, 2.458e-02, 2.472e-02, 2.518e-02, 2.549e-02},
+    {2.629e-02, 2.621e-02, 2.581e-02, 2.527e-02, 2.480e-02, 2.458e-02, 2.451e-02, 2.485e-02, 2.516e-02, 2.547e-02},
+    {2.629e-02, 2.607e-02, 2.573e-02, 2.543e-02, 2.485e-02, 2.464e-02, 2.452e-02, 2.476e-02, 2.505e-02, 2.550e-02},
+    {2.635e-02, 2.613e-02, 2.578e-02, 2.523e-02, 2.491e-02, 2.465e-02, 2.470e-02, 2.467e-02, 2.515e-02, 2.564e-02},
+    {2.613e-02, 2.602e-02, 2.587e-02, 2.526e-02, 2.507e-02, 2.482e-02, 2.456e-02, 2.486e-02, 2.509e-02, 2.572e-02},
+    {2.620e-02, 2.599e-02, 2.563e-02, 2.528e-02, 2.484e-02, 2.462e-02, 2.464e-02, 2.476e-02, 2.513e-02, 2.571e-02},
+    {2.634e-02, 2.596e-02, 2.565e-02, 2.519e-02, 2.497e-02, 2.457e-02, 2.450e-02, 2.481e-02, 2.511e-02, 2.540e-02},
+    {2.593e-02, 2.589e-02, 2.563e-02, 2.511e-02, 2.472e-02, 2.453e-02, 2.452e-02, 2.474e-02, 2.501e-02, 2.543e-02},
+    {2.576e-02, 2.582e-02, 2.526e-02, 2.505e-02, 2.462e-02, 2.446e-02, 2.445e-02, 2.466e-02, 2.486e-02, 2.510e-02},
+    {2.571e-02, 2.549e-02, 2.533e-02, 2.501e-02, 2.453e-02, 2.443e-02, 2.445e-02, 2.450e-02, 2.448e-02, 2.469e-02},
+    {2.812e-02, 2.786e-02, 2.776e-02, 2.723e-02, 2.695e-02, 2.650e-02, 2.642e-02, 2.617e-02, 2.612e-02, 2.610e-02},
+    {3.251e-02, 3.267e-02, 3.223e-02, 3.183e-02, 3.125e-02, 3.106e-02, 3.067e-02, 3.010e-02, 2.936e-02, 2.927e-02}
+  };
+  if(z>=0. && z<.25) return sy[tb][Int_t(z/.025)];
+
+  Double_t m = 1.e-8; for(Int_t id=10; id--;) if(sy[tb][id]>m) m=sy[tb][id];
+
+  return m;
+}
+
+//___________________________________________________________________________
+Double_t AliTRDcluster::GetXcorr(Int_t tb, Double_t z)
 {
   // drift length correction [cm]
   // TODO to be parametrized in term of drift velocity
   // A.Bercuci (Mar 28 2009)
 
   if(tb<0 || tb>=24) return 0.;
-  const Float_t cx[] = {
-    1.6402e-01, 7.2917e-02,-6.7848e-02,-1.4529e-01,-1.6279e-01,-1.3116e-01,
-   -8.2712e-02,-4.9453e-02,-2.9501e-02,-1.4543e-02,-6.1749e-03, 3.9221e-04,
-    1.9711e-03, 2.7388e-03, 2.9070e-03, 3.4183e-03, 2.8014e-03, 1.9351e-03,
-    4.9252e-04, 4.5742e-04, 1.2263e-04,-1.2219e-02,-6.9658e-02,-1.6681e-01};
-  return cx[tb];
+  const Int_t nd = 5;
+  const Double_t dx[24][nd]={
+    {+1.747e-01,+3.195e-01,+1.641e-01,+1.607e-01,+6.002e-01},
+    {+5.468e-02,+5.760e-02,+6.365e-02,+8.003e-02,+1.067e-01},
+    {-6.327e-02,-6.339e-02,-6.423e-02,-6.900e-02,-7.949e-02},
+    {-1.417e-01,-1.424e-01,-1.450e-01,-1.465e-01,-1.514e-01},
+    {-1.637e-01,-1.619e-01,-1.622e-01,-1.613e-01,-1.648e-01},
+    {-1.386e-01,-1.334e-01,-1.261e-01,-1.276e-01,-1.314e-01},
+    {-8.799e-02,-8.299e-02,-7.861e-02,-8.038e-02,-8.436e-02},
+    {-5.139e-02,-4.849e-02,-4.641e-02,-4.965e-02,-5.286e-02},
+    {-2.927e-02,-2.773e-02,-2.807e-02,-3.021e-02,-3.378e-02},
+    {-1.380e-02,-1.229e-02,-1.335e-02,-1.547e-02,-1.984e-02},
+    {-4.168e-03,-4.601e-03,-5.462e-03,-8.164e-03,-1.035e-02},
+    {+2.044e-03,+1.889e-03,+9.603e-04,-1.342e-03,-3.736e-03},
+    {+3.568e-03,+3.581e-03,+2.391e-03,+2.942e-05,-1.585e-03},
+    {+4.403e-03,+4.571e-03,+3.509e-03,+8.703e-04,-1.425e-03},
+    {+4.941e-03,+4.808e-03,+3.284e-03,+1.105e-03,-1.208e-03},
+    {+5.124e-03,+5.022e-03,+4.305e-03,+2.023e-03,-1.145e-03},
+    {+4.882e-03,+4.008e-03,+3.408e-03,+7.886e-04,-1.356e-03},
+    {+3.852e-03,+3.539e-03,+2.057e-03,+1.670e-04,-1.993e-03},
+    {+2.154e-03,+2.111e-03,+5.723e-04,-1.254e-03,-3.256e-03},
+    {+1.755e-03,+2.101e-03,+9.516e-04,-1.649e-03,-3.394e-03},
+    {+1.617e-03,+1.662e-03,+4.169e-04,-9.843e-04,-4.309e-03},
+    {-9.204e-03,-9.069e-03,-1.182e-02,-1.458e-02,-1.880e-02},
+    {-6.727e-02,-6.820e-02,-6.804e-02,-7.134e-02,-7.615e-02},
+    {-1.802e-01,-1.733e-01,-1.633e-01,-1.601e-01,-1.632e-01}
+  };
+//   const Double_t dx[24][nd]={
+//     {+0.000e+00,+0.000e+00,+0.000e+00,+0.000e+00,+0.000e+00,+0.000e+00,+0.000e+00,+0.000e+00,+0.000e+00,+0.000e+00},
+//     {-2.763e-04,-2.380e-04,-6.286e-04,-9.424e-04,+1.046e-03,+1.932e-03,+1.620e-03,+1.951e-03,-1.321e-03,-1.115e-03},
+//     {-1.825e-03,-9.245e-04,-1.012e-03,-8.215e-04,+2.703e-05,+1.403e-03,+2.340e-03,+2.577e-03,+2.017e-03,+8.006e-04},
+//     {-3.070e-03,-8.563e-04,-1.257e-03,+8.491e-05,+4.503e-04,-2.467e-05,-1.793e-04,+5.085e-04,+1.321e-03,+4.056e-04},
+//     {-3.637e-03,-2.857e-03,-3.098e-03,-2.304e-03,-1.467e-03,-1.755e-03,+4.585e-04,+2.757e-03,+3.184e-03,+3.525e-03},
+//     {-9.884e-03,-7.695e-03,-7.290e-03,-3.990e-03,-9.982e-04,+2.226e-03,+3.375e-03,+6.419e-03,+7.209e-03,+6.891e-03},
+//     {-6.844e-03,-5.807e-03,-4.012e-03,-1.566e-03,+5.035e-04,+2.024e-03,+3.225e-03,+3.918e-03,+5.942e-03,+6.024e-03},
+//     {-2.628e-03,-2.201e-03,-4.562e-04,+9.832e-04,+3.411e-03,+2.062e-03,+1.526e-03,+9.350e-04,+8.842e-04,+1.007e-03},
+//     {+6.603e-04,+1.545e-03,+1.681e-03,+1.918e-03,+2.165e-03,+1.825e-03,+1.691e-03,-1.923e-04,+1.835e-04,-1.284e-03},
+//     {+1.895e-03,+1.586e-03,+2.000e-03,+3.537e-03,+2.526e-03,+1.316e-03,+8.229e-04,-7.671e-05,-2.175e-03,-3.528e-03},
+//     {+2.927e-03,+3.369e-03,+3.603e-03,+2.675e-03,+2.737e-03,+1.133e-03,+4.318e-04,-1.215e-03,-2.443e-03,-3.116e-03},
+//     {+3.072e-03,+3.564e-03,+3.612e-03,+3.149e-03,+2.768e-03,+1.186e-03,+3.083e-04,-1.447e-03,-2.480e-03,-3.263e-03},
+//     {+2.697e-03,+3.565e-03,+3.759e-03,+2.855e-03,+2.909e-03,+6.564e-04,-5.224e-04,-3.309e-04,-1.636e-03,-3.739e-03},
+//     {+3.633e-03,+3.232e-03,+3.727e-03,+3.024e-03,+3.365e-03,+1.598e-03,-6.903e-04,-1.039e-03,-3.176e-03,-4.472e-03},
+//     {+2.999e-03,+3.942e-03,+3.322e-03,+3.162e-03,+1.978e-03,+1.657e-03,-4.760e-04,-8.343e-04,-2.346e-03,-3.281e-03},
+//     {+3.734e-03,+3.098e-03,+3.435e-03,+2.512e-03,+2.651e-03,+1.745e-03,+9.424e-04,-1.404e-03,-3.177e-03,-4.444e-03},
+//     {+3.204e-03,+4.003e-03,+3.068e-03,+2.697e-03,+3.187e-03,+3.878e-04,-1.124e-04,-1.855e-03,-2.584e-03,-3.807e-03},
+//     {+2.653e-03,+3.631e-03,+2.327e-03,+3.460e-03,+1.810e-03,+1.244e-03,-3.651e-04,-2.664e-04,-2.307e-03,-3.642e-03},
+//     {+2.538e-03,+3.208e-03,+2.390e-03,+3.519e-03,+1.763e-03,+1.330e-04,+1.669e-04,-1.422e-03,-1.685e-03,-3.519e-03},
+//     {+2.605e-03,+2.465e-03,+2.771e-03,+2.966e-03,+2.361e-03,+6.029e-04,-4.435e-04,-1.876e-03,-1.694e-03,-3.757e-03},
+//     {+2.866e-03,+3.315e-03,+3.146e-03,+2.117e-03,+1.933e-03,+9.339e-04,+9.556e-04,-1.314e-03,-3.615e-03,-3.558e-03},
+//     {+4.002e-03,+3.543e-03,+3.631e-03,+4.127e-03,+1.919e-03,-2.852e-04,-9.484e-04,-2.060e-03,-4.477e-03,-5.491e-03},
+//     {+6.029e-03,+5.147e-03,+4.286e-03,+2.215e-03,+9.240e-04,-1.554e-03,-2.366e-03,-3.635e-03,-5.372e-03,-6.467e-03},
+//     {+3.941e-03,+3.995e-03,+5.638e-04,-3.332e-04,-2.539e-03,-3.764e-03,-3.647e-03,-4.900e-03,-5.414e-03,-5.202e-03}
+//   };
+  if(z>=0. && z<.25) return dx[tb][Int_t(z/.025)];
+
+  Double_t m = 0.; for(Int_t id=nd; id--;) m+=dx[tb][id];
+  return m/nd;
 }
 
 //___________________________________________________________________________
-Float_t AliTRDcluster::GetYcorr(Int_t ly, Float_t y)
+Double_t AliTRDcluster::GetYcorr(Int_t ly, Float_t y)
 {
-  // PRF correction TODO to be replaced by the gaussian 
-  // approximation with full error parametrization and // moved to the clusterizer
+// PRF correction TODO to be replaced by the gaussian 
+// approximation with full error parametrization
   const Float_t cy[AliTRDgeometry::kNlayer][3] = {
     { 4.014e-04, 8.605e-03, -6.880e+00},
     {-3.061e-04, 9.663e-03, -6.789e+00},
@@ -272,7 +423,7 @@ Float_t AliTRDcluster::GetYcorr(Int_t ly, Float_t y)
 }
 
 //_____________________________________________________________________________
-Float_t AliTRDcluster::GetXloc(Double_t t0, Double_t vd, Double_t *const /*q*/, Double_t *const /*xq*/, Double_t z)
+Float_t AliTRDcluster::GetXloc(Double_t t0, Double_t vd, Double_t *const /*q*/, Double_t *const /*xq*/, Double_t /*z*/)
 {
 //
 // (Re)Calculate cluster position in the x direction in local chamber coordinates (with respect to the anode wire 
@@ -308,15 +459,28 @@ Float_t AliTRDcluster::GetXloc(Double_t t0, Double_t vd, Double_t *const /*q*/,
 
   AliTRDCommonParam *cp = AliTRDCommonParam::Instance(); 
   Double_t fFreq = cp->GetSamplingFrequency();
+
+  // calculate t0 corrected time bin
+  Double_t td = fPadTime - t0;
+  fLocalTimeBin = TMath::Nint(td);
   //drift time corresponding to the center of the time bin
-  Double_t td = (fPadTime + .5)/fFreq; // [us] 
+  td = (td + .5)/fFreq; // [us] 
   // correction for t0
   td -= t0;
   // calculate radial posion of clusters in the drift region
-  if(td < .2 || td > 2.4) return 0.; 
-  // correction for TRF rising time 0.2us
+  if(td < .2) return 0.;
+  // TRF rising time (fitted)
+  // It should be absorbed by the t0. For the moment t0 is 0 for simulations.
+  // A.Bercuci (Mar 26 2009)
   td -= 0.189;
 
+  // apply fitted correction 
+  Float_t x = td*vd + GetXcorr(fLocalTimeBin);
+  if(x>.5*AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght()) SetInChamber(kFALSE);
+
+  return x;
+
+/*
   // invert drift time function
   Double_t xM= AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght(),
            x = vd*td + .5*AliTRDgeometry::CamHght(), 
@@ -334,55 +498,22 @@ Float_t AliTRDcluster::GetXloc(Double_t t0, Double_t vd, Double_t *const /*q*/,
   }
 
   return x-.5*AliTRDgeometry::CamHght();
+*/
 }
 
 //_____________________________________________________________________________
-Float_t AliTRDcluster::GetYloc(Double_t s2, Double_t W, Double_t xd, Double_t wt, Double_t *const y1, Double_t *const y2)
+Float_t AliTRDcluster::GetYloc(Double_t y0, Double_t s2, Double_t W, Double_t *const y1, Double_t *const y2)
 {
-//
-// (Re)Calculate cluster position in the y direction in local chamber coordinates using all available information from tracking.
-//
-// Input parameters:
-//   s2 - sigma of gaussian parameterization (see bellow for the exact parameterization)
-//   W  - pad width
-//   xd - drift length (with respect to the anode wire) [cm]
-//   wt - omega*tau = tg(a_L)
-// Output values :
-//   y1 and y2 - partial positions based on 2 pads clusters
-//   return y position of the cluster from all information
-//
-// Estimation of y coordinate is based on the gaussian approximation of the PRF. Thus one may
-// calculate the y position knowing the signals q_i-1, q_i and q_i+1 in the 3 adiacent pads by:
-// BEGIN_LATEX
-// y = #frac{1}{w_{1}+w_{2}}#[]{w_{1}#(){y_{0}-#frac{W}{2}+#frac{s^{2}}{W}ln#frac{q_{i}}{q_{i-1}}}+w_{2}#(){y_{0}+ #frac{W}{2}+#frac{s^{2}}{W}ln#frac{q_{i+1}}{q_{i}}}}
-// END_LATEX
-// where W is the pad width, y_0 is the position of the center pad and s^2 is given by
-// BEGIN_LATEX
-// s^{2} = s^{2}_{0} + s^{2}_{diff} (x,B) + #frac{tg^{2}(#phi-#alpha_{L})*l^{2}}{12}
-// END_LATEX
-// with s_0 being the PRF for 0 drift and track incidence phi equal to the lorentz angle a_L and the diffusion term 
-// being described by:
-// BEGIN_LATEX
-// s_{diff} (x,B) = #frac{D_{L}#sqrt{x}}{1+#(){#omega#tau}^{2}}
-// END_LATEX
-// with x being the drift length. The weights w_1 and w_2 are taken to be q_i-1^2 and q_i+1^2 respectively
-// 
-// Authors
-// Alex Bercuci <A.Bercuci@gsi.de>
-// Theodor Rascanu <trascanu@stud.uni-frankfurt.de>
-//
-  Float_t y0 = GetY()-W*fCenter;
-  Double_t w1 = fSignals[2]*fSignals[2];
-  Double_t w2 = fSignals[4]*fSignals[4];
-  Float_t y1r = fSignals[2]>0 ? (y0 - .5*W + s2*TMath::Log(fSignals[3]/(Float_t)fSignals[2])/W) : 0.;
-  Float_t y2r = fSignals[4]>0 ? (y0 + .5*W + s2*TMath::Log(fSignals[4]/(Float_t)fSignals[3])/W) : 0.;
 
-  if(y1) (*y1) = y1r;
-  if(y2) (*y2) = y2r;
+  //printf("  s[%3d %3d %3d] w[%f %f] yr[%f %f]\n", fSignals[2], fSignals[3], fSignals[4], w1/(w1+w2), w2/(w1+w2), y1r*W, y2r*W);
+  if(IsRPhiMethod(kCOG)) GetDYcog();
+  else if(IsRPhiMethod(kLUT)) GetDYlut();
+  if(IsRPhiMethod(kGAUS)) GetDYgauss(s2/W/W, y1, y2);
 
-  Double_t ld  = TMath::Max(xd - 0.*AliTRDgeometry::CamHght(), 0.);
+  if(y1) (*y1)*=W;
+  if(y2) (*y2)*=W;
 
-  return (w1*y1r+w2*y2r)/(w1+w2) - ld*wt;
+  return y0-fCenter*W;
 }
 
 //_____________________________________________________________________________
@@ -414,7 +545,7 @@ Bool_t AliTRDcluster::IsEqual(const TObject *o) const
 //_____________________________________________________________________________
 void AliTRDcluster::Print(Option_t *o) const
 {
-  AliInfo(Form("Det[%3d] LTrC[%7.2f %7.2f %7.2f] Q[%f] Stat[in(%c) use(%c) sh(%c)]", 
+  AliInfo(Form("Det[%3d] LTrC[%+6.2f %+6.2f %+6.2f] Q[%5.1f] Stat[in(%c) use(%c) sh(%c)]", 
     fDetector, GetX(), GetY(), GetZ(), fQ, 
     IsInChamber() ? 'y' : 'n', IsUsed() ? 'y' : 'n', IsShared() ? 'y' : 'n'));
 
@@ -453,3 +584,258 @@ void AliTRDcluster::SetPadMaskedStatus(UChar_t status)
     if(TESTBIT(status, ipos))
       SETBIT(fClusterMasking, ipos + 3);
 }
+
+//___________________________________________________________________________
+Float_t AliTRDcluster::GetDYcog(Double_t *const, Double_t *const)
+{
+//
+// Get COG position
+// Used for clusters with more than 3 pads - where LUT not applicable
+//
+
+  Double_t sum = fSignals[1]
+                +fSignals[2]
+                +fSignals[3] 
+                +fSignals[4]
+                +fSignals[5];
+
+  // ???????????? CBL
+  // Go to 3 pad COG ????
+  // ???????????? CBL
+  fCenter = (0.0 * (-fSignals[1] + fSignals[5])
+                      + (-fSignals[2] + fSignals[4])) / sum;
+
+  return fCenter;
+}
+
+//___________________________________________________________________________
+Float_t AliTRDcluster::GetDYlut(Double_t *const, Double_t *const)
+{
+  //
+  // Calculates the cluster position using the lookup table.
+  // Method provided by Bogdan Vulpescu.
+  //
+
+  if(!fgLUT) FillLUT();
+
+  Double_t ampL = fSignals[2],
+           ampC = fSignals[3],
+           ampR = fSignals[4];
+  Int_t ilayer = AliTRDgeometry::GetLayer(fDetector);
+
+  Double_t x    = 0.0;
+  Double_t xmin, xmax, xwid;
+
+  Int_t    side = 0;
+  Int_t    ix;
+
+  Double_t xMin[AliTRDgeometry::kNlayer] = { 
+    0.006492, 0.006377, 0.006258, 0.006144, 0.006030, 0.005980 
+  };
+  Double_t xMax[AliTRDgeometry::kNlayer] = { 
+    0.960351, 0.965870, 0.970445, 0.974352, 0.977667, 0.996101 
+  };
+
+  if      (ampL > ampR) {
+    x    = (ampL - ampR) / ampC;
+    side = -1;
+  } 
+  else if (ampL < ampR) {
+    x    = (ampR - ampL) / ampC;
+    side = +1;
+  }
+
+  if (ampL != ampR) {
+
+    xmin = xMin[ilayer] + 0.000005;
+    xmax = xMax[ilayer] - 0.000005;
+    xwid = (xmax - xmin) / 127.0;
+
+    if      (x < xmin) fCenter = 0.0000;
+    else if (x > xmax) fCenter = side * 0.5000;
+    else {
+      ix      = (Int_t) ((x - xmin) / xwid);
+      fCenter = side * fgLUT[ilayer*fgkNlut+ix];
+    }
+  } else fCenter = 0.0;
+
+  return fCenter;
+}
+
+//___________________________________________________________________________
+Float_t AliTRDcluster::GetDYgauss(Double_t s2w, Double_t *const y1, Double_t *const y2)
+{
+//
+// (Re)Calculate cluster position in the y direction in local chamber coordinates using all available information from tracking.
+//
+// Input parameters:
+//   s2 - sigma of gaussian parameterization (see bellow for the exact parameterization)
+//   W  - pad width
+//   xd - drift length (with respect to the anode wire) [cm]
+//   wt - omega*tau = tg(a_L)
+// Output values :
+//   y1 and y2 - partial positions based on 2 pads clusters
+//   return y position of the cluster from all information
+//
+// Estimation of y coordinate is based on the gaussian approximation of the PRF. Thus one may
+// calculate the y position knowing the signals q_i-1, q_i and q_i+1 in the 3 adiacent pads by:
+// BEGIN_LATEX
+// y = #frac{1}{w_{1}+w_{2}}#[]{w_{1}#(){y_{0}-#frac{W}{2}+#frac{s^{2}}{W}ln#frac{q_{i}}{q_{i-1}}}+w_{2}#(){y_{0}+ #frac{W}{2}+#frac{s^{2}}{W}ln#frac{q_{i+1}}{q_{i}}}}
+// END_LATEX
+// where W is the pad width, y_0 is the position of the center pad and s^2 is given by
+// BEGIN_LATEX
+// s^{2} = s^{2}_{0} + s^{2}_{diff} (x,B) + #frac{tg^{2}(#phi-#alpha_{L})*l^{2}}{12}
+// END_LATEX
+// with s_0 being the PRF for 0 drift and track incidence phi equal to the lorentz angle a_L and the diffusion term 
+// being described by:
+// BEGIN_LATEX
+// s_{diff} (x,B) = #frac{D_{L}#sqrt{x}}{1+#(){#omega#tau}^{2}}
+// END_LATEX
+// with x being the drift length. The weights w_1 and w_2 are taken to be q_i-1^2 and q_i+1^2 respectively
+// 
+// Authors
+// Alex Bercuci <A.Bercuci@gsi.de>
+// Theodor Rascanu <trascanu@stud.uni-frankfurt.de>
+//
+  Double_t w1  = fSignals[2]*fSignals[2];
+  Double_t w2  = fSignals[4]*fSignals[4];
+  //Double_t s2w = s2/W/W;
+  Float_t y1r  = fSignals[2]>0 ? (-0.5 + s2w*TMath::Log(fSignals[3]/(Float_t)fSignals[2])) : 0.;
+  Float_t y2r  = fSignals[4]>0 ? (0.5 + s2w*TMath::Log(fSignals[4]/(Float_t)fSignals[3])) : 0.;
+
+  if(y1) (*y1) = y1r;
+  if(y2) (*y2) = y2r;
+
+  return fCenter      = (w1*y1r+w2*y2r)/(w1+w2);
+}
+
+
+
+//_____________________________________________________________________________
+void AliTRDcluster::FillLUT()
+{
+  //
+  // Create the LUT
+  //
+
+  // The lookup table from Bogdan
+  Float_t lut[AliTRDgeometry::kNlayer][fgkNlut] = {  
+    {
+      0.0070, 0.0150, 0.0224, 0.0298, 0.0374, 0.0454, 0.0533, 0.0611, 
+      0.0684, 0.0755, 0.0827, 0.0900, 0.0975, 0.1049, 0.1120, 0.1187, 
+      0.1253, 0.1318, 0.1385, 0.1453, 0.1519, 0.1584, 0.1646, 0.1704, 
+      0.1762, 0.1821, 0.1879, 0.1938, 0.1996, 0.2053, 0.2108, 0.2160, 
+      0.2210, 0.2260, 0.2310, 0.2361, 0.2411, 0.2461, 0.2509, 0.2557, 
+      0.2602, 0.2646, 0.2689, 0.2732, 0.2774, 0.2816, 0.2859, 0.2901, 
+      0.2942, 0.2983, 0.3022, 0.3061, 0.3099, 0.3136, 0.3172, 0.3207, 
+      0.3242, 0.3278, 0.3312, 0.3347, 0.3382, 0.3416, 0.3450, 0.3483, 
+      0.3515, 0.3547, 0.3579, 0.3609, 0.3639, 0.3669, 0.3698, 0.3727, 
+      0.3756, 0.3785, 0.3813, 0.3842, 0.3870, 0.3898, 0.3926, 0.3952, 
+      0.3979, 0.4005, 0.4032, 0.4057, 0.4082, 0.4108, 0.4132, 0.4157, 
+      0.4181, 0.4205, 0.4228, 0.4252, 0.4275, 0.4299, 0.4322, 0.4345, 
+      0.4367, 0.4390, 0.4412, 0.4434, 0.4456, 0.4478, 0.4499, 0.4520, 
+      0.4541, 0.4562, 0.4583, 0.4603, 0.4623, 0.4643, 0.4663, 0.4683, 
+      0.4702, 0.4722, 0.4741, 0.4758, 0.4774, 0.4790, 0.4805, 0.4824, 
+      0.4844, 0.4863, 0.4883, 0.4902, 0.4921, 0.4940, 0.4959, 0.4978 
+    },
+    {
+      0.0072, 0.0156, 0.0235, 0.0313, 0.0394, 0.0478, 0.0561, 0.0642, 
+      0.0718, 0.0792, 0.0868, 0.0947, 0.1025, 0.1101, 0.1172, 0.1241, 
+      0.1309, 0.1378, 0.1449, 0.1518, 0.1586, 0.1650, 0.1710, 0.1770, 
+      0.1830, 0.1891, 0.1952, 0.2011, 0.2070, 0.2125, 0.2177, 0.2229, 
+      0.2280, 0.2332, 0.2383, 0.2435, 0.2484, 0.2533, 0.2581, 0.2627, 
+      0.2670, 0.2714, 0.2757, 0.2799, 0.2842, 0.2884, 0.2927, 0.2968, 
+      0.3008, 0.3048, 0.3086, 0.3123, 0.3159, 0.3195, 0.3231, 0.3266, 
+      0.3301, 0.3335, 0.3370, 0.3404, 0.3438, 0.3471, 0.3504, 0.3536, 
+      0.3567, 0.3598, 0.3628, 0.3657, 0.3686, 0.3715, 0.3744, 0.3772, 
+      0.3800, 0.3828, 0.3856, 0.3884, 0.3911, 0.3938, 0.3965, 0.3991, 
+      0.4016, 0.4042, 0.4067, 0.4092, 0.4116, 0.4140, 0.4164, 0.4187, 
+      0.4211, 0.4234, 0.4257, 0.4280, 0.4302, 0.4325, 0.4347, 0.4369, 
+      0.4391, 0.4413, 0.4434, 0.4456, 0.4477, 0.4497, 0.4518, 0.4538, 
+      0.4558, 0.4578, 0.4598, 0.4618, 0.4637, 0.4656, 0.4675, 0.4694, 
+      0.4713, 0.4732, 0.4750, 0.4766, 0.4781, 0.4797, 0.4813, 0.4832, 
+      0.4851, 0.4870, 0.4888, 0.4906, 0.4925, 0.4942, 0.4960, 0.4978
+    },
+    {
+      0.0075, 0.0163, 0.0246, 0.0328, 0.0415, 0.0504, 0.0592, 0.0674, 
+      0.0753, 0.0832, 0.0914, 0.0996, 0.1077, 0.1154, 0.1225, 0.1296, 
+      0.1369, 0.1442, 0.1515, 0.1585, 0.1652, 0.1714, 0.1776, 0.1839, 
+      0.1902, 0.1965, 0.2025, 0.2085, 0.2141, 0.2194, 0.2247, 0.2299, 
+      0.2352, 0.2405, 0.2457, 0.2507, 0.2557, 0.2604, 0.2649, 0.2693, 
+      0.2737, 0.2780, 0.2823, 0.2867, 0.2909, 0.2951, 0.2992, 0.3033, 
+      0.3072, 0.3110, 0.3146, 0.3182, 0.3218, 0.3253, 0.3288, 0.3323, 
+      0.3357, 0.3392, 0.3426, 0.3459, 0.3492, 0.3524, 0.3555, 0.3586, 
+      0.3616, 0.3645, 0.3674, 0.3703, 0.3731, 0.3759, 0.3787, 0.3815, 
+      0.3843, 0.3870, 0.3897, 0.3925, 0.3950, 0.3976, 0.4002, 0.4027, 
+      0.4052, 0.4076, 0.4101, 0.4124, 0.4148, 0.4171, 0.4194, 0.4217, 
+      0.4239, 0.4262, 0.4284, 0.4306, 0.4328, 0.4350, 0.4371, 0.4393, 
+      0.4414, 0.4435, 0.4455, 0.4476, 0.4496, 0.4516, 0.4536, 0.4555, 
+      0.4575, 0.4594, 0.4613, 0.4632, 0.4650, 0.4669, 0.4687, 0.4705, 
+      0.4723, 0.4741, 0.4758, 0.4773, 0.4789, 0.4804, 0.4821, 0.4839, 
+      0.4857, 0.4875, 0.4893, 0.4910, 0.4928, 0.4945, 0.4961, 0.4978
+    },
+    {
+      0.0078, 0.0171, 0.0258, 0.0345, 0.0438, 0.0532, 0.0624, 0.0708, 
+      0.0791, 0.0875, 0.0962, 0.1048, 0.1130, 0.1206, 0.1281, 0.1356, 
+      0.1432, 0.1508, 0.1582, 0.1651, 0.1716, 0.1780, 0.1845, 0.1910, 
+      0.1975, 0.2038, 0.2099, 0.2155, 0.2210, 0.2263, 0.2317, 0.2371, 
+      0.2425, 0.2477, 0.2528, 0.2578, 0.2626, 0.2671, 0.2715, 0.2759, 
+      0.2803, 0.2846, 0.2890, 0.2933, 0.2975, 0.3016, 0.3056, 0.3095, 
+      0.3132, 0.3168, 0.3204, 0.3239, 0.3274, 0.3309, 0.3344, 0.3378, 
+      0.3412, 0.3446, 0.3479, 0.3511, 0.3543, 0.3574, 0.3603, 0.3633, 
+      0.3662, 0.3690, 0.3718, 0.3747, 0.3774, 0.3802, 0.3829, 0.3857, 
+      0.3883, 0.3910, 0.3936, 0.3962, 0.3987, 0.4012, 0.4037, 0.4061, 
+      0.4085, 0.4109, 0.4132, 0.4155, 0.4177, 0.4200, 0.4222, 0.4244, 
+      0.4266, 0.4288, 0.4309, 0.4331, 0.4352, 0.4373, 0.4394, 0.4414, 
+      0.4435, 0.4455, 0.4475, 0.4494, 0.4514, 0.4533, 0.4552, 0.4571, 
+      0.4590, 0.4608, 0.4626, 0.4645, 0.4662, 0.4680, 0.4698, 0.4715, 
+      0.4733, 0.4750, 0.4766, 0.4781, 0.4796, 0.4812, 0.4829, 0.4846, 
+      0.4863, 0.4880, 0.4897, 0.4914, 0.4930, 0.4946, 0.4963, 0.4979
+    },
+    {
+      0.0081, 0.0178, 0.0270, 0.0364, 0.0463, 0.0562, 0.0656, 0.0744, 
+      0.0831, 0.0921, 0.1013, 0.1102, 0.1183, 0.1261, 0.1339, 0.1419, 
+      0.1499, 0.1576, 0.1648, 0.1715, 0.1782, 0.1849, 0.1917, 0.1984, 
+      0.2048, 0.2110, 0.2167, 0.2223, 0.2278, 0.2333, 0.2389, 0.2444, 
+      0.2497, 0.2548, 0.2598, 0.2645, 0.2691, 0.2735, 0.2780, 0.2824, 
+      0.2868, 0.2912, 0.2955, 0.2997, 0.3038, 0.3078, 0.3116, 0.3152, 
+      0.3188, 0.3224, 0.3259, 0.3294, 0.3329, 0.3364, 0.3398, 0.3432, 
+      0.3465, 0.3497, 0.3529, 0.3561, 0.3591, 0.3620, 0.3649, 0.3677, 
+      0.3705, 0.3733, 0.3761, 0.3788, 0.3816, 0.3843, 0.3869, 0.3896, 
+      0.3922, 0.3948, 0.3973, 0.3998, 0.4022, 0.4047, 0.4070, 0.4094, 
+      0.4117, 0.4139, 0.4162, 0.4184, 0.4206, 0.4227, 0.4249, 0.4270, 
+      0.4291, 0.4313, 0.4334, 0.4354, 0.4375, 0.4395, 0.4415, 0.4435, 
+      0.4455, 0.4474, 0.4493, 0.4512, 0.4531, 0.4550, 0.4568, 0.4586, 
+      0.4604, 0.4622, 0.4639, 0.4657, 0.4674, 0.4691, 0.4708, 0.4725, 
+      0.4742, 0.4758, 0.4773, 0.4788, 0.4803, 0.4819, 0.4836, 0.4852, 
+      0.4869, 0.4885, 0.4901, 0.4917, 0.4933, 0.4948, 0.4964, 0.4979
+    },
+    {
+      0.0085, 0.0189, 0.0288, 0.0389, 0.0497, 0.0603, 0.0699, 0.0792, 
+      0.0887, 0.0985, 0.1082, 0.1170, 0.1253, 0.1336, 0.1421, 0.1505, 
+      0.1587, 0.1662, 0.1733, 0.1803, 0.1874, 0.1945, 0.2014, 0.2081, 
+      0.2143, 0.2201, 0.2259, 0.2316, 0.2374, 0.2431, 0.2487, 0.2541, 
+      0.2593, 0.2642, 0.2689, 0.2735, 0.2781, 0.2826, 0.2872, 0.2917, 
+      0.2961, 0.3003, 0.3045, 0.3086, 0.3125, 0.3162, 0.3198, 0.3235, 
+      0.3270, 0.3306, 0.3342, 0.3377, 0.3411, 0.3446, 0.3479, 0.3511, 
+      0.3543, 0.3575, 0.3605, 0.3634, 0.3663, 0.3691, 0.3720, 0.3748, 
+      0.3775, 0.3803, 0.3830, 0.3857, 0.3884, 0.3911, 0.3937, 0.3962, 
+      0.3987, 0.4012, 0.4036, 0.4060, 0.4084, 0.4107, 0.4129, 0.4152, 
+      0.4174, 0.4196, 0.4218, 0.4239, 0.4261, 0.4282, 0.4303, 0.4324, 
+      0.4344, 0.4365, 0.4385, 0.4405, 0.4425, 0.4445, 0.4464, 0.4483, 
+      0.4502, 0.4521, 0.4539, 0.4558, 0.4576, 0.4593, 0.4611, 0.4629, 
+      0.4646, 0.4663, 0.4680, 0.4697, 0.4714, 0.4730, 0.4747, 0.4759, 
+      0.4769, 0.4780, 0.4790, 0.4800, 0.4811, 0.4827, 0.4843, 0.4859, 
+      0.4874, 0.4889, 0.4905, 0.4920, 0.4935, 0.4950, 0.4965, 0.4979
+    }
+  }; 
+
+  if(!fgLUT) fgLUT = new Double_t[AliTRDgeometry::kNlayer*fgkNlut];
+
+  for (Int_t ilayer = 0; ilayer < AliTRDgeometry::kNlayer; ilayer++) {
+    for (Int_t ilut  = 0; ilut  < fgkNlut; ilut++  ) {
+      fgLUT[ilayer*fgkNlut+ilut] = lut[ilayer][ilut];
+    }
+  }
+}
+
index cfc732e..e4cbe44 100644 (file)
 
 class AliTRDcluster : public AliCluster {
 public:
-
   enum ETRDclusterStatus { 
     kInChamber = BIT(16) // Out of fiducial volume of chamber (signal tails)
    ,kFivePad   = BIT(17) // Deconvoluted clusters
+   ,kLUT       = BIT(18)
+   ,kGAUS      = BIT(19)
+   ,kCOG       = BIT(20)
   };
   enum ETRDclusterMask { 
     kMaskedLeft   = 0
@@ -27,6 +29,7 @@ public:
   };
 
   AliTRDcluster();
+  AliTRDcluster(Int_t det, UChar_t col, UChar_t row, UChar_t time, const Short_t *sig, UShort_t volid);
   AliTRDcluster(
     Int_t det, Float_t q, Float_t *pos, Float_t *sig,
     Int_t *tracks, Char_t npads, Short_t *signals,
@@ -44,6 +47,7 @@ public:
   Bool_t   IsShared() const                { return IsClusterShared();}
   Bool_t   IsUsed() const                  { return IsClusterUsed(); }
   Bool_t   IsFivePad() const               { return TestBit(kFivePad);}
+  inline Bool_t IsRPhiMethod(ETRDclusterStatus m) const;
 
   UChar_t  GetPadMaskedPosition() const    { return fClusterMasking & 7; }
   UChar_t  GetPadMaskedStatus() const      { return fClusterMasking >> 3; }
@@ -57,20 +61,24 @@ public:
   Int_t    GetPadTime() const              { return fPadTime;       }
   Short_t *GetSignals()                    { return fSignals;       }
   Float_t  GetSumS() const;
-  static Float_t  GetXcorr(Int_t tb);
-  static Float_t  GetYcorr(Int_t ly, Float_t y);
+  static Double_t  GetSX(Int_t tb, Double_t z=-1);
+  static Double_t  GetSY(Int_t tb, Double_t z=-1);
+  static Double_t  GetXcorr(Int_t tb, Double_t z=-1);
+  static Double_t  GetYcorr(Int_t ly, Float_t y);
   Float_t  GetXloc(Double_t t0, Double_t vd, Double_t *const q=0x0, Double_t *const xq=0x0, Double_t z = 0.2);
-  Float_t  GetYloc(Double_t s2, Double_t W, Double_t xd, Double_t wt, Double_t *const y1=0x0, Double_t *const y2=0x0);
+  Float_t  GetYloc(Double_t y0, Double_t s2, Double_t W, Double_t *const y1=0x0, Double_t *const y2=0x0);
 
   void     Print(Option_t* o="") const;
 
   void     SetLocalTimeBin(Char_t t)       { fLocalTimeBin = t;     }
   void     SetInChamber(Bool_t in = kTRUE) { SetBit(kInChamber,in); }
+  void     SetNPads(Int_t n)               { fNPads = n; }
   void     SetPadMaskedPosition(UChar_t position);
   void     SetPadMaskedStatus(UChar_t status);
   void     SetPadCol(UChar_t inPadCol){ fPadCol = inPadCol;}
   void     SetPadRow(UChar_t inPadRow){ fPadRow = inPadRow;}
   void     SetPadTime(UChar_t inPadTime){ fPadTime = inPadTime;}
+  inline void SetRPhiMethod(ETRDclusterStatus m);
   void     SetDetector(Short_t inDetector){ fDetector = inDetector;}
   void     SetQ(Float_t inQ){ fQ = inQ;}
   void     SetClusterMasking(UChar_t inClusterMasking){ fClusterMasking = inClusterMasking;}
@@ -90,6 +98,37 @@ protected:
   Float_t fQ;              //  Amplitude 
   Float_t fCenter;         //  Center of the cluster relative to the pad 
 
-  ClassDef(AliTRDcluster, 6)        //  Cluster for the TRD
+private:
+  Float_t  GetDYcog(Double_t *const y1=0x0, Double_t *const y2=0x0);
+  Float_t  GetDYlut(Double_t *const y1=0x0, Double_t *const y2=0x0);
+  Float_t  GetDYgauss(Double_t sw, Double_t *const y1=0x0, Double_t *const y2=0x0);
+  static void      FillLUT();
+
+  static const Int_t   fgkNlut;              //!  Number of bins of the LUT
+  static Double_t     *fgLUT;                //! The lookup table
+
+  ClassDef(AliTRDcluster, 7)        //  Cluster for the TRD
 };
+
+//________________________________________________
+inline Bool_t AliTRDcluster::IsRPhiMethod(ETRDclusterStatus m) const
+{
+  if(m==kLUT && TestBit(kLUT)) return kTRUE;
+  else if(m==kGAUS && TestBit(kGAUS)) return kTRUE;
+  else if(m==kCOG && (!TestBit(kLUT)&&!TestBit(kGAUS))) return kTRUE;
+
+  return kFALSE;
+}
+
+//________________________________________________
+inline void AliTRDcluster::SetRPhiMethod(ETRDclusterStatus m)
+{
+  switch(m){
+  case kCOG: SetBit(kCOG); break;
+  case kLUT: SetBit(kLUT); break;
+  case kGAUS: SetBit(kGAUS); break;
+  default: break;
+  }
+}
+
 #endif
index 6d8486e..203fd74 100644 (file)
@@ -70,8 +70,6 @@ AliTRDclusterizer::AliTRDclusterizer(const AliTRDReconstructor *const rec)
   ,fTrackletContainer(NULL)
   ,fRawVersion(2)
   ,fTransform(new AliTRDtransform(0))
-  ,fLUTbin(0)
-  ,fLUT(NULL)
   ,fDigits(NULL)
   ,fIndexes(NULL)
   ,fADCthresh(0)
@@ -97,7 +95,7 @@ AliTRDclusterizer::AliTRDclusterizer(const AliTRDReconstructor *const rec)
   // AliTRDclusterizer default constructor
   //
 
-  SetBit(kAddLabels, kTRUE);
+  SetBit(kLabels, kTRUE);
 
   AliTRDcalibDB *trd = 0x0;
   if (!(trd = AliTRDcalibDB::Instance())) {
@@ -129,8 +127,6 @@ AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title, co
   ,fTrackletContainer(NULL)
   ,fRawVersion(2)
   ,fTransform(new AliTRDtransform(0))
-  ,fLUTbin(0)
-  ,fLUT(NULL)
   ,fDigits(NULL)
   ,fIndexes(NULL)
   ,fADCthresh(0)
@@ -156,7 +152,7 @@ AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title, co
   // AliTRDclusterizer constructor
   //
 
-  SetBit(kAddLabels, kTRUE);
+  SetBit(kLabels, kTRUE);
 
   AliTRDcalibDB *trd = 0x0;
   if (!(trd = AliTRDcalibDB::Instance())) {
@@ -167,7 +163,7 @@ AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title, co
 
   fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
 
-  FillLUT();
+  //FillLUT();
 
 }
 
@@ -183,8 +179,6 @@ AliTRDclusterizer::AliTRDclusterizer(const AliTRDclusterizer &c)
   ,fTrackletContainer(NULL)
   ,fRawVersion(2)
   ,fTransform(NULL)
-  ,fLUTbin(0)
-  ,fLUT(0)
   ,fDigits(NULL)
   ,fIndexes(NULL)
   ,fADCthresh(0)
@@ -210,9 +204,9 @@ AliTRDclusterizer::AliTRDclusterizer(const AliTRDclusterizer &c)
   // AliTRDclusterizer copy constructor
   //
 
-  SetBit(kAddLabels, kTRUE);
+  SetBit(kLabels, kTRUE);
 
-  FillLUT();
+  //FillLUT();
 
 }
 
@@ -243,10 +237,10 @@ AliTRDclusterizer::~AliTRDclusterizer()
     fTransform     = NULL;
   }
 
-  if (fLUT) {
-    delete [] fLUT;
-    fLUT           = NULL;
-  }
+//   if (fLUT) {
+//     delete [] fLUT;
+//     fLUT           = NULL;
+//   }
 
 }
 
@@ -280,8 +274,6 @@ void AliTRDclusterizer::Copy(TObject &c) const
   ((AliTRDclusterizer &) c).fTrackletContainer = NULL;
   ((AliTRDclusterizer &) c).fRawVersion    = fRawVersion;
   ((AliTRDclusterizer &) c).fTransform     = NULL;
-  ((AliTRDclusterizer &) c).fLUTbin        = 0;
-  ((AliTRDclusterizer &) c).fLUT           = NULL;
   ((AliTRDclusterizer &) c).fDigits      = NULL;
   ((AliTRDclusterizer &) c).fIndexes       = NULL;
   ((AliTRDclusterizer &) c).fADCthresh     = 0;
@@ -561,8 +553,8 @@ Bool_t AliTRDclusterizer::MakeClusters()
   //
 
   // Propagate info from the digits manager
-  if (TestBit(kAddLabels)){
-    SetBit(kAddLabels, fDigitsManager->UsesDictionaries());
+  if (TestBit(kLabels)){
+    SetBit(kLabels, fDigitsManager->UsesDictionaries());
   }
   
   Bool_t fReturn = kTRUE;
@@ -580,7 +572,7 @@ Bool_t AliTRDclusterizer::MakeClusters()
   
     Bool_t fR = kFALSE;
     if (indexes->HasEntry()){
-      if (TestBit(kAddLabels)){
+      if (TestBit(kLabels)){
         for (Int_t iDict = 0; iDict < AliTRDdigitsManager::kNDict; iDict++){
           AliTRDarrayDictionary *tracksIn = 0; //mod
           tracksIn = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(i,iDict);  //mod
@@ -634,7 +626,7 @@ Bool_t AliTRDclusterizer::Raw2ClustersChamber(AliRawReader *rawReader)
     fDigitsManager->CreateArrays();
   }
 
-  fDigitsManager->SetUseDictionaries(TestBit(kAddLabels));
+  fDigitsManager->SetUseDictionaries(TestBit(kLabels));
 
   // tracklet container for raw tracklet writing
   if (!fTrackletContainer && fReconstructor->IsWritingTracklets()) {
@@ -741,24 +733,19 @@ Bool_t AliTRDclusterizer::MakeClusters(Int_t det)
   fDigits = (AliTRDarrayADC *) fDigitsManager->GetDigits(det); //mod     
   
   // This is to take care of switched off super modules
-  if (!fDigits->HasData()) 
-    {
-      return kFALSE;
-    }
+  if (!fDigits->HasData()) return kFALSE;
 
   fIndexes = fDigitsManager->GetIndexes(det);
-  if (fIndexes->IsAllocated() == kFALSE)
-    {
-      AliError("Indexes do not exist!");
-      return kFALSE;      
-    }
+  if (fIndexes->IsAllocated() == kFALSE) {
+    AliError("Indexes do not exist!");
+    return kFALSE;      
+  }
 
   AliTRDcalibDB  *calibration = AliTRDcalibDB::Instance();
-  if (!calibration) 
-    {
-      AliFatal("No AliTRDcalibDB instance available\n");
-      return kFALSE;  
-    }
+  if (!calibration) {
+    AliFatal("No AliTRDcalibDB instance available\n");
+    return kFALSE;  
+  }
 
   fADCthresh = 0; 
 
@@ -814,49 +801,46 @@ Bool_t AliTRDclusterizer::MakeClusters(Int_t det)
   // Calibration object with the pad status
   fCalPadStatusROC       = calibration->GetPadStatusROC(fDet);
   
-  SetBit(kIsLUT, fReconstructor->GetRecoParam()->IsLUT());
-  SetBit(kIsHLT, fReconstructor->IsHLT());
+  SetBit(kLUT, fReconstructor->UseLUT());
+  SetBit(kGAUS, fReconstructor->UseGAUS());
+  SetBit(kHLT, fReconstructor->IsHLT());
 
   firstClusterROC = -1;
   fClusterROC     =  0;
 
-  if(fReconstructor->GetRecoParam()->IsTailCancelation()){
-    // Apply the gain and the tail cancelation via digital filter
-    TailCancelation();
-  }
+  // Apply the gain and the tail cancelation via digital filter
+  if(fReconstructor->UseTailCancelation()) TailCancelation();
 
   MaxStruct curr, last;
   Int_t nMaximas = 0, nCorrupted = 0;
 
   // Here the clusterfining is happening
   
-  for(curr.Time = 0; curr.Time < fTimeTotal; curr.Time++)
-    while(fIndexes->NextRCIndex(curr.Row, curr.Col))
-      if(IsMaximum(curr, curr.padStatus, &curr.Signals[0]))
-       {
-         if(last.Row>-1)
-           {
-             if(curr.Time==last.Time && curr.Row==last.Row && curr.Col==last.Col+2)
-               FivePadCluster(last, curr);
-             CreateCluster(last);
-           }
-         last=curr; curr.FivePad=kFALSE;
-       }
-  if(last.Row>-1)
-    CreateCluster(last);
+  for(curr.Time = 0; curr.Time < fTimeTotal; curr.Time++){
+    while(fIndexes->NextRCIndex(curr.Row, curr.Col)){
+      //printf("\nCHECK r[%2d] c[%3d] t[%d]\n", curr.Row, curr.Col, curr.Time);
+      if(IsMaximum(curr, curr.padStatus, &curr.Signals[0])){
+        //printf("\tMAX s[%d %d %d %d %d %d %d]\n", curr.Signals[0], curr.Signals[1], curr.Signals[2], curr.Signals[3], curr.Signals[4], curr.Signals[5],  curr.Signals[6]);
+        if(last.Row>-1){
+          if(curr.Time==last.Time && curr.Row==last.Row && curr.Col==last.Col+2) FivePadCluster(last, curr);
+          CreateCluster(last);
+        }
+        last=curr; curr.FivePad=kFALSE;
+      }
+      //printf("\t--- s[%d %d %d %d %d %d %d]\n", curr.Signals[0], curr.Signals[1], curr.Signals[2], curr.Signals[3], curr.Signals[4], curr.Signals[5],  curr.Signals[6]);
+    }
+  }
+  if(last.Row>-1) CreateCluster(last);
 
   if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kClusterizer) > 2){
     (*fDebugStream) << "MakeClusters"
-                   << "Detector="   << det
-                   << "NMaxima="    << nMaximas
-                   << "NClusters="  << fClusterROC
-                   << "NCorrupted=" << nCorrupted
-                   << "\n";
-  }
-
-  if (TestBit(kAddLabels)) {
-    AddLabels();
+      << "Detector="   << det
+      << "NMaxima="    << nMaximas
+      << "NClusters="  << fClusterROC
+      << "NCorrupted=" << nCorrupted
+      << "\n";
   }
+  if (TestBit(kLabels)) AddLabels();
 
   return kTRUE;
 
@@ -869,47 +853,44 @@ Bool_t AliTRDclusterizer::IsMaximum(const MaxStruct &Max, UChar_t &padStatus, Sh
   // Returns true if this row,col,time combination is a maximum. 
   // Gives back the padStatus and the signals of the center pad and the two neighbouring pads.
   //
-
-  Signals[1] = fDigits->GetData(Max.Row, Max.Col, Max.Time);
-  if(Signals[1] < fMaxThresh) return kFALSE;
+  fDigits->GetData(Max.Row, Max.Col-3, Max.Time, 7, &Signals[0]);
+  if(Signals[3] < fMaxThresh) return kFALSE;
 
   Float_t  noiseMiddleThresh = fMinMaxCutSigma*fCalNoiseDetValue*fCalNoiseROC->GetValue(Max.Col, Max.Row);
-  if (Signals[1] < noiseMiddleThresh) return kFALSE;
+  if (Signals[3] < noiseMiddleThresh) return kFALSE;
 
   if (Max.Col + 1 >= fColMax || Max.Col < 1) return kFALSE;
 
-  UChar_t status[3]={fCalPadStatusROC->GetStatus(Max.Col-1, Max.Row), 
-                    fCalPadStatusROC->GetStatus(Max.Col,   Max.Row), 
-                    fCalPadStatusROC->GetStatus(Max.Col+1, Max.Row)};
-
-  Signals[0] = fDigits->GetData(Max.Row, Max.Col-1, Max.Time);
-  Signals[2] = fDigits->GetData(Max.Row, Max.Col+1, Max.Time);  
+  UChar_t status[3]={
+    fCalPadStatusROC->GetStatus(Max.Col-1, Max.Row)
+   ,fCalPadStatusROC->GetStatus(Max.Col,   Max.Row)
+   ,fCalPadStatusROC->GetStatus(Max.Col+1, Max.Row)
+  };
 
   if(!(status[0] | status[1] | status[2])) {//all pads are good
-    if ((Signals[2] <= Signals[1]) && (Signals[0] <  Signals[1])) {
-      if ((Signals[2] >= fSigThresh) || (Signals[0] >= fSigThresh)) {
-       Float_t  noiseSumThresh = fMinLeftRightCutSigma
-         * fCalNoiseDetValue
-         * fCalNoiseROC->GetValue(Max.Col, Max.Row);
-       if ((Signals[2]+Signals[0]+Signals[1]) < noiseSumThresh) return kFALSE;
-       padStatus = 0;
-       return kTRUE;
+    if ((Signals[4] <= Signals[3]) && (Signals[2] <  Signals[3])) {
+      if ((Signals[4] >= fSigThresh) || (Signals[2] >= fSigThresh)) {
+        Float_t  noiseSumThresh = fMinLeftRightCutSigma
+          * fCalNoiseDetValue
+          * fCalNoiseROC->GetValue(Max.Col, Max.Row);
+        if ((Signals[2]+Signals[3]+Signals[4]) < noiseSumThresh) return kFALSE;
+        padStatus = 0;
+        return kTRUE;
       }
     }
-  }
-  else { // at least one of the pads is bad, and reject candidates with more than 1 problematic pad
-    if (status[2] && (!(status[0] || status[1])) && Signals[1] > Signals[0] && Signals[0] >= fSigThresh) { 
-      Signals[2]=0;
+  } else { // at least one of the pads is bad, and reject candidates with more than 1 problematic pad
+    if (status[2] && (!(status[0] || status[1])) && Signals[3] > Signals[2] && Signals[2] >= fSigThresh) { 
+      Signals[4]=0;
       SetPadStatus(status[2], padStatus);
       return kTRUE;
     } 
-    else if (status[0] && (!(status[1] || status[2])) && Signals[1] >= Signals[2] && Signals[2] >= fSigThresh) {
-      Signals[0]=0;
+    else if (status[0] && (!(status[1] || status[2])) && Signals[3] >= Signals[4] && Signals[4] >= fSigThresh) {
+      Signals[2]=0;
       SetPadStatus(status[0], padStatus);
       return kTRUE;
     }
-    else if (status[1] && (!(status[0] || status[2])) && ((Signals[2] >= fSigThresh) || (Signals[0] >= fSigThresh))) {
-      Signals[1]=TMath::Nint(fMaxThresh);
+    else if (status[1] && (!(status[0] || status[2])) && ((Signals[4] >= fSigThresh) || (Signals[2] >= fSigThresh))) {
+      Signals[3]=TMath::Nint(fMaxThresh);
       SetPadStatus(status[1], padStatus);
       return kTRUE;
     }
@@ -924,7 +905,6 @@ Bool_t AliTRDclusterizer::FivePadCluster(MaxStruct &ThisMax, MaxStruct &Neighbou
   // Look for 5 pad cluster with minimum in the middle
   // Gives back the ratio
   //
-
   if (ThisMax.Col >= fColMax - 3) return kFALSE;
   if (ThisMax.Col < fColMax - 5){
     if (fDigits->GetData(ThisMax.Row, ThisMax.Col+4, ThisMax.Time) >= fSigThresh)
@@ -938,14 +918,16 @@ Bool_t AliTRDclusterizer::FivePadCluster(MaxStruct &ThisMax, MaxStruct &Neighbou
   //if (fSignalsThisMax[1] >= 0){ //TR: mod
   
   const Float_t kEpsilon = 0.01;
-  Double_t padSignal[5] = {ThisMax.Signals[0], ThisMax.Signals[1], ThisMax.Signals[2],
-                          NeighbourMax.Signals[1], NeighbourMax.Signals[2]};
+  Double_t padSignal[5] = {
+    ThisMax.Signals[2], ThisMax.Signals[3], ThisMax.Signals[4],
+    NeighbourMax.Signals[3], NeighbourMax.Signals[4]
+  };
   
   // Unfold the two maxima and set the signal on 
   // the overlapping pad to the ratio
   Float_t ratio = Unfold(kEpsilon,fLayer,padSignal);
-  ThisMax.Signals[2] = TMath::Nint(ThisMax.Signals[2]*ratio);
-  NeighbourMax.Signals[0] = TMath::Nint(NeighbourMax.Signals[0]*(1-ratio));
+  ThisMax.Signals[4] = TMath::Nint(ThisMax.Signals[4]*ratio);
+  NeighbourMax.Signals[2] = TMath::Nint(NeighbourMax.Signals[2]*(1-ratio));
   ThisMax.FivePad=kTRUE;
   NeighbourMax.FivePad=kTRUE;
   return kTRUE;
@@ -958,107 +940,41 @@ void AliTRDclusterizer::CreateCluster(const MaxStruct &Max)
   // Creates a cluster at the given position and saves it in fRecPoints
   //
 
-  // The position of the cluster in COL direction relative to the center pad (pad units)
-  Double_t clusterPosCol = 0.0;
-  if (TestBit(kIsLUT)) {
-    // Calculate the position of the cluster by using the
-    // lookup table method
-    clusterPosCol = LUTposition(fLayer,Max.Signals[0]
-                               ,Max.Signals[1]
-                               ,Max.Signals[2]);
-  } 
-  else {
-    // Calculate the position of the cluster by using the
-    // center of gravity method
-    const Int_t kNsig = 5;
-    Double_t padSignal[kNsig];
-    padSignal[1] = Max.Signals[0];
-    padSignal[2] = Max.Signals[1];
-    padSignal[3] = Max.Signals[2];
-    if(Max.Col > 2){
-      padSignal[0] = fDigits->GetData(Max.Row, Max.Col-2, Max.Time);
-      if(padSignal[0]>= padSignal[1])
-       padSignal[0] = 0;
-    }
-    if(Max.Col < fColMax - 3){
-      padSignal[4] = fDigits->GetData(Max.Row, Max.Col+2, Max.Time);
-      if(padSignal[4]>= padSignal[3])
-       padSignal[4] = 0;
-    }
-    clusterPosCol = GetCOG(padSignal);
-  }
+  AliTRDcluster cluster(fDet, ((UChar_t) Max.Col), ((UChar_t) Max.Row), ((UChar_t) Max.Time), Max.Signals, fVolid);
+  if(TestBit(kLUT)) cluster.SetRPhiMethod(AliTRDcluster::kLUT);
+  else if(TestBit(kGAUS)) cluster.SetRPhiMethod(AliTRDcluster::kGAUS);
+  else cluster.SetRPhiMethod(AliTRDcluster::kCOG);
 
-  // Count the number of pads in the cluster
+  cluster.SetFivePad(Max.FivePad);
+  // set pads status for the cluster
+  UChar_t maskPosition = GetCorruption(Max.padStatus);
+  if (maskPosition) { 
+    cluster.SetPadMaskedPosition(maskPosition);
+    cluster.SetPadMaskedStatus(GetPadStatus(Max.padStatus));
+  }
+  // set number of pads in the cluster
   Int_t nPadCount = 1;
   Short_t signals[7] = { 0, 0, 0, 0, 0, 0, 0 };
-
-  if(!TestBit(kIsHLT))CalcAdditionalInfo(Max, signals, nPadCount);
+  if(!TestBit(kHLT)) CalcAdditionalInfo(Max, signals, nPadCount);
+  cluster.SetNPads(nPadCount);
 
   // Transform the local cluster coordinates into calibrated 
   // space point positions defined in the local tracking system.
   // Here the calibration for T0, Vdrift and ExB is applied as well.
-  Double_t clusterXYZ[6];
-  clusterXYZ[0] = clusterPosCol;
-  clusterXYZ[1] = Max.Signals[2];
-  clusterXYZ[2] = Max.Signals[1];
-  clusterXYZ[3] = Max.Signals[0];
-  clusterXYZ[4] = 0.0;
-  clusterXYZ[5] = 0.0;
-  Int_t    clusterRCT[3];
-  clusterRCT[0] = Max.Row;
-  clusterRCT[1] = Max.Col;
-  clusterRCT[2] = 0;
-
-  Bool_t out = kTRUE;
-  if (fTransform->Transform(clusterXYZ,clusterRCT,((UInt_t) Max.Time),out,0)) {
-
-    Char_t  clusterTimeBin = ((Char_t) clusterRCT[2]);
-    Float_t clusterPos[3];
-    clusterPos[0] = clusterXYZ[0];
-    clusterPos[1] = clusterXYZ[1];
-    clusterPos[2] = clusterXYZ[2];
-    Float_t clusterSig[2];
-    clusterSig[0] = clusterXYZ[4];
-    clusterSig[1] = clusterXYZ[5];
-    Float_t clusterCharge  = clusterXYZ[3];
-    
-    AliTRDcluster cluster(
-                         fDet,
-                         clusterCharge, clusterPos, clusterSig,
-                         0x0,
-                         ((Char_t) nPadCount),
-                         signals,
-                         ((UChar_t) Max.Col), ((UChar_t) Max.Row), ((UChar_t) Max.Time),
-                         clusterTimeBin, clusterPosCol,
-                         fVolid);
-    
-    cluster.SetInChamber(!out);
-    cluster.SetFivePad(Max.FivePad);
-    
-    UChar_t maskPosition = GetCorruption(Max.padStatus);
-    if (maskPosition) { 
-      cluster.SetPadMaskedPosition(maskPosition);
-      cluster.SetPadMaskedStatus(GetPadStatus(Max.padStatus));
-    }
-    
-    // Temporarily store the Max.Row, column and time bin of the center pad
-    // Used to later on assign the track indices
-    cluster.SetLabel(Max.Row, 0);
-    cluster.SetLabel(Max.Col, 1);
-    cluster.SetLabel(Max.Time,2);
-    
-    AddClusterToArray(&cluster); //needs to be like that because HLT does things differently
-    
-    //AddCluster(Max,clusterXYZ,clusterTimeBin,signals,nPadCount,out,clusterPosCol);
-    // Store the index of the first cluster in the current ROC
-    if (firstClusterROC < 0) {
-      firstClusterROC = fNoOfClusters;
-    }
-    
-    fNoOfClusters++;
-    fClusterROC++;
-  }
-
+  if(!fTransform->Transform(&cluster)) return;
+  // Temporarily store the Max.Row, column and time bin of the center pad
+  // Used to later on assign the track indices
+  cluster.SetLabel(Max.Row, 0);
+  cluster.SetLabel(Max.Col, 1);
+  cluster.SetLabel(Max.Time,2);
+  //needed for HLT reconstruction
+  AddClusterToArray(&cluster); 
+
+  // Store the index of the first cluster in the current ROC
+  if (firstClusterROC < 0) firstClusterROC = fNoOfClusters;
+  
+  fNoOfClusters++;
+  fClusterROC++;
 }
 
 //_____________________________________________________________________________
@@ -1178,30 +1094,6 @@ Bool_t AliTRDclusterizer::AddLabels()
 }
 
 //_____________________________________________________________________________
-Double_t AliTRDclusterizer::GetCOG(Double_t signal[5]) const
-{
-  //
-  // Get COG position
-  // Used for clusters with more than 3 pads - where LUT not applicable
-  //
-
-  Double_t sum = signal[0]
-              + signal[1]
-              + signal[2] 
-              + signal[3]
-              + signal[4];
-
-  // ???????????? CBL
-  // Go to 3 pad COG ????
-  // ???????????? CBL
-  Double_t res = (0.0 * (-signal[0] + signal[4])
-                      + (-signal[1] + signal[3])) / sum;
-
-  return res;            
-
-}
-
-//_____________________________________________________________________________
 Float_t AliTRDclusterizer::Unfold(Double_t eps, Int_t layer, Double_t *padSignal) const
 {
   //
@@ -1431,202 +1323,3 @@ TClonesArray *AliTRDclusterizer::RecPoints()
 
 }
 
-//_____________________________________________________________________________
-void AliTRDclusterizer::FillLUT()
-{
-  //
-  // Create the LUT
-  //
-
-  const Int_t kNlut = 128;
-
-  fLUTbin = AliTRDgeometry::kNlayer * kNlut;
-
-  // The lookup table from Bogdan
-  Float_t lut[AliTRDgeometry::kNlayer][kNlut] = {  
-    {
-      0.0070, 0.0150, 0.0224, 0.0298, 0.0374, 0.0454, 0.0533, 0.0611, 
-      0.0684, 0.0755, 0.0827, 0.0900, 0.0975, 0.1049, 0.1120, 0.1187, 
-      0.1253, 0.1318, 0.1385, 0.1453, 0.1519, 0.1584, 0.1646, 0.1704, 
-      0.1762, 0.1821, 0.1879, 0.1938, 0.1996, 0.2053, 0.2108, 0.2160, 
-      0.2210, 0.2260, 0.2310, 0.2361, 0.2411, 0.2461, 0.2509, 0.2557, 
-      0.2602, 0.2646, 0.2689, 0.2732, 0.2774, 0.2816, 0.2859, 0.2901, 
-      0.2942, 0.2983, 0.3022, 0.3061, 0.3099, 0.3136, 0.3172, 0.3207, 
-      0.3242, 0.3278, 0.3312, 0.3347, 0.3382, 0.3416, 0.3450, 0.3483, 
-      0.3515, 0.3547, 0.3579, 0.3609, 0.3639, 0.3669, 0.3698, 0.3727, 
-      0.3756, 0.3785, 0.3813, 0.3842, 0.3870, 0.3898, 0.3926, 0.3952, 
-      0.3979, 0.4005, 0.4032, 0.4057, 0.4082, 0.4108, 0.4132, 0.4157, 
-      0.4181, 0.4205, 0.4228, 0.4252, 0.4275, 0.4299, 0.4322, 0.4345, 
-      0.4367, 0.4390, 0.4412, 0.4434, 0.4456, 0.4478, 0.4499, 0.4520, 
-      0.4541, 0.4562, 0.4583, 0.4603, 0.4623, 0.4643, 0.4663, 0.4683, 
-      0.4702, 0.4722, 0.4741, 0.4758, 0.4774, 0.4790, 0.4805, 0.4824, 
-      0.4844, 0.4863, 0.4883, 0.4902, 0.4921, 0.4940, 0.4959, 0.4978 
-    },
-    {
-      0.0072, 0.0156, 0.0235, 0.0313, 0.0394, 0.0478, 0.0561, 0.0642, 
-      0.0718, 0.0792, 0.0868, 0.0947, 0.1025, 0.1101, 0.1172, 0.1241, 
-      0.1309, 0.1378, 0.1449, 0.1518, 0.1586, 0.1650, 0.1710, 0.1770, 
-      0.1830, 0.1891, 0.1952, 0.2011, 0.2070, 0.2125, 0.2177, 0.2229, 
-      0.2280, 0.2332, 0.2383, 0.2435, 0.2484, 0.2533, 0.2581, 0.2627, 
-      0.2670, 0.2714, 0.2757, 0.2799, 0.2842, 0.2884, 0.2927, 0.2968, 
-      0.3008, 0.3048, 0.3086, 0.3123, 0.3159, 0.3195, 0.3231, 0.3266, 
-      0.3301, 0.3335, 0.3370, 0.3404, 0.3438, 0.3471, 0.3504, 0.3536, 
-      0.3567, 0.3598, 0.3628, 0.3657, 0.3686, 0.3715, 0.3744, 0.3772, 
-      0.3800, 0.3828, 0.3856, 0.3884, 0.3911, 0.3938, 0.3965, 0.3991, 
-      0.4016, 0.4042, 0.4067, 0.4092, 0.4116, 0.4140, 0.4164, 0.4187, 
-      0.4211, 0.4234, 0.4257, 0.4280, 0.4302, 0.4325, 0.4347, 0.4369, 
-      0.4391, 0.4413, 0.4434, 0.4456, 0.4477, 0.4497, 0.4518, 0.4538, 
-      0.4558, 0.4578, 0.4598, 0.4618, 0.4637, 0.4656, 0.4675, 0.4694, 
-      0.4713, 0.4732, 0.4750, 0.4766, 0.4781, 0.4797, 0.4813, 0.4832, 
-      0.4851, 0.4870, 0.4888, 0.4906, 0.4925, 0.4942, 0.4960, 0.4978
-    },
-    {
-      0.0075, 0.0163, 0.0246, 0.0328, 0.0415, 0.0504, 0.0592, 0.0674, 
-      0.0753, 0.0832, 0.0914, 0.0996, 0.1077, 0.1154, 0.1225, 0.1296, 
-      0.1369, 0.1442, 0.1515, 0.1585, 0.1652, 0.1714, 0.1776, 0.1839, 
-      0.1902, 0.1965, 0.2025, 0.2085, 0.2141, 0.2194, 0.2247, 0.2299, 
-      0.2352, 0.2405, 0.2457, 0.2507, 0.2557, 0.2604, 0.2649, 0.2693, 
-      0.2737, 0.2780, 0.2823, 0.2867, 0.2909, 0.2951, 0.2992, 0.3033, 
-      0.3072, 0.3110, 0.3146, 0.3182, 0.3218, 0.3253, 0.3288, 0.3323, 
-      0.3357, 0.3392, 0.3426, 0.3459, 0.3492, 0.3524, 0.3555, 0.3586, 
-      0.3616, 0.3645, 0.3674, 0.3703, 0.3731, 0.3759, 0.3787, 0.3815, 
-      0.3843, 0.3870, 0.3897, 0.3925, 0.3950, 0.3976, 0.4002, 0.4027, 
-      0.4052, 0.4076, 0.4101, 0.4124, 0.4148, 0.4171, 0.4194, 0.4217, 
-      0.4239, 0.4262, 0.4284, 0.4306, 0.4328, 0.4350, 0.4371, 0.4393, 
-      0.4414, 0.4435, 0.4455, 0.4476, 0.4496, 0.4516, 0.4536, 0.4555, 
-      0.4575, 0.4594, 0.4613, 0.4632, 0.4650, 0.4669, 0.4687, 0.4705, 
-      0.4723, 0.4741, 0.4758, 0.4773, 0.4789, 0.4804, 0.4821, 0.4839, 
-      0.4857, 0.4875, 0.4893, 0.4910, 0.4928, 0.4945, 0.4961, 0.4978
-    },
-    {
-      0.0078, 0.0171, 0.0258, 0.0345, 0.0438, 0.0532, 0.0624, 0.0708, 
-      0.0791, 0.0875, 0.0962, 0.1048, 0.1130, 0.1206, 0.1281, 0.1356, 
-      0.1432, 0.1508, 0.1582, 0.1651, 0.1716, 0.1780, 0.1845, 0.1910, 
-      0.1975, 0.2038, 0.2099, 0.2155, 0.2210, 0.2263, 0.2317, 0.2371, 
-      0.2425, 0.2477, 0.2528, 0.2578, 0.2626, 0.2671, 0.2715, 0.2759, 
-      0.2803, 0.2846, 0.2890, 0.2933, 0.2975, 0.3016, 0.3056, 0.3095, 
-      0.3132, 0.3168, 0.3204, 0.3239, 0.3274, 0.3309, 0.3344, 0.3378, 
-      0.3412, 0.3446, 0.3479, 0.3511, 0.3543, 0.3574, 0.3603, 0.3633, 
-      0.3662, 0.3690, 0.3718, 0.3747, 0.3774, 0.3802, 0.3829, 0.3857, 
-      0.3883, 0.3910, 0.3936, 0.3962, 0.3987, 0.4012, 0.4037, 0.4061, 
-      0.4085, 0.4109, 0.4132, 0.4155, 0.4177, 0.4200, 0.4222, 0.4244, 
-      0.4266, 0.4288, 0.4309, 0.4331, 0.4352, 0.4373, 0.4394, 0.4414, 
-      0.4435, 0.4455, 0.4475, 0.4494, 0.4514, 0.4533, 0.4552, 0.4571, 
-      0.4590, 0.4608, 0.4626, 0.4645, 0.4662, 0.4680, 0.4698, 0.4715, 
-      0.4733, 0.4750, 0.4766, 0.4781, 0.4796, 0.4812, 0.4829, 0.4846, 
-      0.4863, 0.4880, 0.4897, 0.4914, 0.4930, 0.4946, 0.4963, 0.4979
-    },
-    {
-      0.0081, 0.0178, 0.0270, 0.0364, 0.0463, 0.0562, 0.0656, 0.0744, 
-      0.0831, 0.0921, 0.1013, 0.1102, 0.1183, 0.1261, 0.1339, 0.1419, 
-      0.1499, 0.1576, 0.1648, 0.1715, 0.1782, 0.1849, 0.1917, 0.1984, 
-      0.2048, 0.2110, 0.2167, 0.2223, 0.2278, 0.2333, 0.2389, 0.2444, 
-      0.2497, 0.2548, 0.2598, 0.2645, 0.2691, 0.2735, 0.2780, 0.2824, 
-      0.2868, 0.2912, 0.2955, 0.2997, 0.3038, 0.3078, 0.3116, 0.3152, 
-      0.3188, 0.3224, 0.3259, 0.3294, 0.3329, 0.3364, 0.3398, 0.3432, 
-      0.3465, 0.3497, 0.3529, 0.3561, 0.3591, 0.3620, 0.3649, 0.3677, 
-      0.3705, 0.3733, 0.3761, 0.3788, 0.3816, 0.3843, 0.3869, 0.3896, 
-      0.3922, 0.3948, 0.3973, 0.3998, 0.4022, 0.4047, 0.4070, 0.4094, 
-      0.4117, 0.4139, 0.4162, 0.4184, 0.4206, 0.4227, 0.4249, 0.4270, 
-      0.4291, 0.4313, 0.4334, 0.4354, 0.4375, 0.4395, 0.4415, 0.4435, 
-      0.4455, 0.4474, 0.4493, 0.4512, 0.4531, 0.4550, 0.4568, 0.4586, 
-      0.4604, 0.4622, 0.4639, 0.4657, 0.4674, 0.4691, 0.4708, 0.4725, 
-      0.4742, 0.4758, 0.4773, 0.4788, 0.4803, 0.4819, 0.4836, 0.4852, 
-      0.4869, 0.4885, 0.4901, 0.4917, 0.4933, 0.4948, 0.4964, 0.4979
-    },
-    {
-      0.0085, 0.0189, 0.0288, 0.0389, 0.0497, 0.0603, 0.0699, 0.0792, 
-      0.0887, 0.0985, 0.1082, 0.1170, 0.1253, 0.1336, 0.1421, 0.1505, 
-      0.1587, 0.1662, 0.1733, 0.1803, 0.1874, 0.1945, 0.2014, 0.2081, 
-      0.2143, 0.2201, 0.2259, 0.2316, 0.2374, 0.2431, 0.2487, 0.2541, 
-      0.2593, 0.2642, 0.2689, 0.2735, 0.2781, 0.2826, 0.2872, 0.2917, 
-      0.2961, 0.3003, 0.3045, 0.3086, 0.3125, 0.3162, 0.3198, 0.3235, 
-      0.3270, 0.3306, 0.3342, 0.3377, 0.3411, 0.3446, 0.3479, 0.3511, 
-      0.3543, 0.3575, 0.3605, 0.3634, 0.3663, 0.3691, 0.3720, 0.3748, 
-      0.3775, 0.3803, 0.3830, 0.3857, 0.3884, 0.3911, 0.3937, 0.3962, 
-      0.3987, 0.4012, 0.4036, 0.4060, 0.4084, 0.4107, 0.4129, 0.4152, 
-      0.4174, 0.4196, 0.4218, 0.4239, 0.4261, 0.4282, 0.4303, 0.4324, 
-      0.4344, 0.4365, 0.4385, 0.4405, 0.4425, 0.4445, 0.4464, 0.4483, 
-      0.4502, 0.4521, 0.4539, 0.4558, 0.4576, 0.4593, 0.4611, 0.4629, 
-      0.4646, 0.4663, 0.4680, 0.4697, 0.4714, 0.4730, 0.4747, 0.4759, 
-      0.4769, 0.4780, 0.4790, 0.4800, 0.4811, 0.4827, 0.4843, 0.4859, 
-      0.4874, 0.4889, 0.4905, 0.4920, 0.4935, 0.4950, 0.4965, 0.4979
-    }
-  }; 
-
-  if (fLUT) {
-    delete [] fLUT;
-  }
-  fLUT = new Double_t[fLUTbin];
-
-  for (Int_t ilayer = 0; ilayer < AliTRDgeometry::kNlayer; ilayer++) {
-    for (Int_t ilut  = 0; ilut  <  kNlut; ilut++  ) {
-      fLUT[ilayer*kNlut+ilut] = lut[ilayer][ilut];
-    }
-  }
-
-}
-
-//_____________________________________________________________________________
-Double_t AliTRDclusterizer::LUTposition(Int_t ilayer
-                                      , Double_t ampL
-                                      , Double_t ampC
-                                      , Double_t ampR) const
-{
-  //
-  // Calculates the cluster position using the lookup table.
-  // Method provided by Bogdan Vulpescu.
-  //
-
-  const Int_t kNlut = 128;
-
-  Double_t pos;
-  Double_t x    = 0.0;
-  Double_t xmin;
-  Double_t xmax;
-  Double_t xwid;
-
-  Int_t    side = 0;
-  Int_t    ix;
-
-  Double_t xMin[AliTRDgeometry::kNlayer] = { 0.006492, 0.006377, 0.006258
-                                          , 0.006144, 0.006030, 0.005980 };
-  Double_t xMax[AliTRDgeometry::kNlayer] = { 0.960351, 0.965870, 0.970445
-                                          , 0.974352, 0.977667, 0.996101 };
-
-  if      (ampL > ampR) {
-    x    = (ampL - ampR) / ampC;
-    side = -1;
-  } 
-  else if (ampL < ampR) {
-    x    = (ampR - ampL) / ampC;
-    side = +1;
-  }
-
-  if (ampL != ampR) {
-
-    xmin = xMin[ilayer] + 0.000005;
-    xmax = xMax[ilayer] - 0.000005;
-    xwid = (xmax - xmin) / 127.0;
-
-    if      (x < xmin) {
-      pos = 0.0000;
-    } 
-    else if (x > xmax) {
-      pos = side * 0.5000;
-    } 
-    else {
-      ix  = (Int_t) ((x - xmin) / xwid);
-      pos = side * fLUT[ilayer*kNlut+ix];
-    }
-      
-  } 
-  else {
-
-    pos = 0.0;
-
-  }
-
-  return pos;
-
-}
index ee10bea..0814879 100644 (file)
@@ -39,10 +39,11 @@ class AliTRDclusterizer : public TNamed
 
   // steering flags
   enum{
-    kIsHLT = BIT(12),    //  are we just in HLT?
-    kIsLUT = BIT(13),    //  are we using look up table for center of the cluster?
-    kOwner = BIT(14),    //  is the clusterizer owner of the clusters?
-    kAddLabels  = BIT(15)   //  Should clusters have MC labels?
+    kOwner  = BIT(14)  //  toggle cluster ownership
+   ,kLabels = BIT(15)  //  toggle MC labels for clusters
+   ,kHLT    = BIT(16)  //  HLT mode
+   ,kLUT    = BIT(17)  //  using look up table for cluster's r-phi position
+   ,kGAUS   = BIT(18)  //  using gauss approx. for cluster's r-phi position
   };
 
   struct MaxStruct
@@ -51,13 +52,13 @@ class AliTRDclusterizer : public TNamed
     Int_t       Col;
     Int_t       Time;
     UChar_t     padStatus;
-    Short_t     Signals[3];
+    Short_t     Signals[7];
     Bool_t      FivePad;
   MaxStruct():Row(-1),Col(0),Time(0),padStatus(0),FivePad(kFALSE)
-      {}
+      {memset(Signals, 0, 7*sizeof(Short_t));}
     MaxStruct &operator=(const MaxStruct &a)
     {Row=a.Row; Col=a.Col; Time=a.Time; padStatus=a.padStatus; FivePad=a.FivePad;
-     memcpy(Signals, a.Signals, 3*sizeof(Signals[0])); return *this;}
+     memcpy(Signals, a.Signals, 7*sizeof(Short_t)); return *this;}
   };
   
   AliTRDclusterizer(const AliTRDReconstructor *const rec = 0x0);
@@ -89,7 +90,7 @@ class AliTRDclusterizer : public TNamed
   Bool_t   MakeClusters(Int_t det);
 
   Bool_t   AddLabels();
-  Bool_t   SetAddLabels(const Bool_t kset) { SetBit(kAddLabels, kset); return TestBit(kAddLabels);  } // should we assign labels to clusters
+  Bool_t   SetUseLabels(const Bool_t kset) { SetBit(kLabels, kset); return TestBit(kLabels);  } // should we assign labels to clusters
   void     SetRawVersion(const Int_t iver) { fRawVersion = iver; } // set the expected raw data version
   void             SetReconstructor(const AliTRDReconstructor *rec) {fReconstructor = rec;}
   static UChar_t   GetStatus(Short_t &signal);
@@ -98,16 +99,15 @@ class AliTRDclusterizer : public TNamed
   Bool_t   IsClustersOwner() const {return TestBit(kOwner);}
   virtual void     SetClustersOwner(Bool_t own=kTRUE) {SetBit(kOwner, own); if(!own) {fRecPoints = 0x0; fNoOfClusters=0;} }
 
- protected:
+protected:
 
   void             DeConvExp (const Double_t *const source, Double_t *const target
                             ,const Int_t nTimeTotal, const Int_t nexp);
   void             TailCancelation();
 
   Float_t  Unfold(Double_t eps, Int_t layer, Double_t *padSignal) const;
-  Double_t GetCOG(Double_t signal[5]) const; 
-  void             FillLUT();
-  Double_t LUTposition(Int_t ilayer, Double_t ampL, Double_t ampC, Double_t ampR) const;
+  //Double_t GetCOG(Double_t signal[5]) const; 
+  //Double_t LUTposition(Int_t ilayer, Double_t ampL, Double_t ampC, Double_t ampR) const;
   
   void             SetPadStatus(const UChar_t status, UChar_t &encoding);
   UChar_t          GetPadStatus(UChar_t encoding) const;
@@ -134,9 +134,6 @@ class AliTRDclusterizer : public TNamed
 
   AliTRDtransform     *fTransform;           //! Transforms the reconstructed space points
 
-  Int_t                fLUTbin;              //  Number of bins of the LUT
-  Double_t            *fLUT;                 //! The lookup table
-
   AliTRDarrayADC      *fDigits;               // Array holding the digits
   AliTRDSignalIndex   *fIndexes;              // Array holding the indexes to the digits
   Float_t              fADCthresh;            // ADC thresholds: There is no ADC threshold anymore, and simParam should not be used in clusterizer. KO
index 41f87c8..63b6142 100644 (file)
@@ -66,12 +66,6 @@ AliTRDrecoParam::AliTRDrecoParam()
   //
   // Default constructor
   //
-  SetTailCancelation();
-  SetLUT();
-  SetClusterSharing(kFALSE);
-  SetVertexConstrained();
-  SetImproveTracklets(kFALSE);
-
   fSysCovMatrix[0] = 0.; // y direction (1 cm)
   fSysCovMatrix[1] = 0.; // z direction (1 cm)
   fSysCovMatrix[2] = 0.; // snp
@@ -117,12 +111,6 @@ AliTRDrecoParam::AliTRDrecoParam(const AliTRDrecoParam &ref)
   //
   // Copy constructor
   //
-  SetClusterSharing(ref.IsClusterSharing());
-  SetVertexConstrained(ref.IsVertexConstrained());
-  SetLUT(ref.IsLUT());
-  SetTailCancelation(ref.IsTailCancelation());
-  SetImproveTracklets(ref.HasImproveTracklets());
-
   memcpy(fSysCovMatrix, ref.fSysCovMatrix, 5*sizeof(Double_t));
   memcpy(fPIDThreshold, ref.fPIDThreshold, AliTRDCalPID::kNMom*sizeof(Double_t));
 }
@@ -149,7 +137,6 @@ AliTRDrecoParam *AliTRDrecoParam::GetHighFluxParam()
 
   AliTRDrecoParam *rec = new AliTRDrecoParam();
   rec->fkdNchdy = 4000.; // PbPb in TRD
-  rec->SetImproveTracklets(kTRUE);
 
   return rec;
 
@@ -166,11 +153,10 @@ AliTRDrecoParam *AliTRDrecoParam::GetCosmicTestParam()
   AliTRDrecoParam *par = new AliTRDrecoParam();
   par->fSysCovMatrix[0] = 1.; // y direction (1 cm)
   par->fSysCovMatrix[1] = 1.; // z direction (1 cm)
-  par->SetVertexConstrained(kFALSE);
-  par->fkChi2YSlope       = 0.11853;
-  par->fkChi2ZSlope       = 0.04527;
-       par->fkChi2YCut                                 = 1.;
-       par->fkPhiSlope                                 = 3.17954;
+  par->fkChi2YSlope     = 0.11853;
+  par->fkChi2ZSlope     = 0.04527;
+  par->fkChi2YCut       = 1.;
+  par->fkPhiSlope       = 3.17954;
   par->fkMaxTheta       = 2.1445;
   par->fkMaxPhi         = 2.7475;
   par->fkNMeanClusters  = 12.89;
index 36b1d60..b3cc9e5 100644 (file)
@@ -60,12 +60,6 @@ public:
   static   AliTRDrecoParam *GetHighFluxParam();
   static   AliTRDrecoParam *GetCosmicTestParam();
 
-  Bool_t   IsClusterSharing() const         { return TestBit(kClusterSharing);}
-  Bool_t   IsLUT() const                    { return TestBit(kLUT);}
-  Bool_t   IsTailCancelation() const        { return TestBit(kTC);}
-  Bool_t   IsVertexConstrained() const      { return TestBit(kVertexConstrained); }
-  Bool_t   HasImproveTracklets() const       { return TestBit(kImproveTracklet); }
-
   void     SetMaxTheta(Double_t maxTheta) {fkMaxTheta = maxTheta;}
   void     SetMaxPhi(Double_t maxPhi) {fkMaxPhi = maxPhi;}
   void     SetFindableClusters(Double_t r) {fkFindable = r;}
@@ -77,15 +71,10 @@ public:
   void     SetPhiSlope(Double_t phiSlope) {fkPhiSlope = phiSlope;}
   void     SetNMeanClusters(Double_t meanNclusters) {fkNMeanClusters = meanNclusters;}
   void     SetNSigmaClusters(Double_t sigmaNclusters) {fkNSigmaClusters = sigmaNclusters;} 
-  void     SetClusterSharing(Bool_t share = kTRUE)            { SetBit(kClusterSharing, share);}
-  void     SetImproveTracklets(Bool_t improve = kTRUE)         { SetBit(kImproveTracklet, improve);}
-  void     SetVertexConstrained(Bool_t vc = kTRUE)            { SetBit(kVertexConstrained, vc); }
-  void     SetLUT(Bool_t lut = kTRUE)                         { SetBit(kLUT, lut);};
   void     SetMinMaxCutSigma(Float_t minMaxCutSigma)          { fMinMaxCutSigma   = minMaxCutSigma; }
   void     SetMinLeftRightCutSigma(Float_t minLeftRightCutSigma) { fMinLeftRightCutSigma   = minLeftRightCutSigma; };
   void     SetClusMaxThresh(Float_t thresh)                   { fClusMaxThresh   = thresh; };
   void     SetClusSigThresh(Float_t thresh)                   { fClusSigThresh   = thresh; };
-  void     SetTailCancelation(Bool_t tc = kTRUE)              { SetBit(kTC, tc);  };
   inline void SetPIDThreshold(Double_t *pid);
   void     SetNexponential(Int_t nexp)                        { fTCnexp          = nexp;   };
   inline void SetSysCovMatrix(Double_t *sys);
@@ -93,13 +82,6 @@ public:
   void     SetNumberOfPostsamples(Int_t n)                    { fNumberOfPostsamples = n;}
 
 private:
-  enum{
-    kTC                = BIT(14) // tail cancelation
-   ,kLUT               = BIT(15) // look up table for cluster position determination 
-   ,kClusterSharing    = BIT(16) // Toggle cluster sharing
-   ,kVertexConstrained = BIT(17) // Perform vertex constrained fit
-   ,kImproveTracklet   = BIT(18) // Improve tracklet in the SA TRD track finder 
-  };
   // Physics reference values for TRD
   Double_t  fkdNchdy;                // dNch/dy
   Double_t  fkMaxTheta;              // Maximum theta
index fac4bd3..e6b00a8 100644 (file)
@@ -594,7 +594,8 @@ void AliTRDseedV1::GetCovAt(Double_t x, Double_t *cov) const
 
   Double_t xr     = fX0-x; 
   Double_t sy2    = fCov[0] +2.*xr*fCov[1] + xr*xr*fCov[2];
-  Double_t sz2    = GetPadLength()*GetPadLength()/12.;
+  Double_t sz2    = fS2Z;
+  //GetPadLength()*GetPadLength()/12.;
 
   // insert systematic uncertainties
   if(fReconstructor){
@@ -609,6 +610,8 @@ void AliTRDseedV1::GetCovAt(Double_t x, Double_t *cov) const
   cov[0] = (sy2+t2*sz2)*correction;
   cov[1] = GetTilt()*(sz2 - sy2)*correction;
   cov[2] = (t2*sy2+sz2)*correction;
+
+  //printf("C(%6.1f %+6.3f %6.1f)  [%s]\n", 1.e4*TMath::Sqrt(cov[0]), cov[1], 1.e4*TMath::Sqrt(cov[2]), IsRowCross()?" RC ":"-");
 }
 
 //____________________________________________________________
@@ -682,6 +685,15 @@ Double_t AliTRDseedV1::GetCovInv(Double_t *c, Double_t *d)
 }
 
 //____________________________________________________________________
+UShort_t AliTRDseedV1::GetVolumeId() const
+{
+  Int_t ic=0;
+  while(ic<kNclusters && !fClusters[ic]) ic++;
+  return fClusters[ic] ? fClusters[ic]->GetVolumeId() : 0;
+}
+
+
+//____________________________________________________________________
 void AliTRDseedV1::Calibrate()
 {
 // Retrieve calibration and position parameters from OCDB. 
@@ -758,135 +770,6 @@ void AliTRDseedV1::SetPadPlane(AliTRDpadPlane *p)
 }
 
 
-// //____________________________________________________________________
-// Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t quality, Bool_t kZcorr, AliTRDcluster *c)
-// {
-//   //
-//   // Iterative process to register clusters to the seed.
-//   // In iteration 0 we try only one pad-row and if quality not
-//   // sufficient we try 2 pad-rows (about 5% of tracks cross 2 pad-rows)
-//   //
-//   // debug level 7
-//   //
-//   
-//   if(!fReconstructor->GetRecoParam() ){
-//     AliError("Seed can not be used without a valid RecoParam.");
-//     return kFALSE;
-//   }
-// 
-//   AliTRDchamberTimeBin *layer = 0x0;
-//   if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7){
-//     AliTRDtrackingChamber ch(*chamber);
-//     ch.SetOwner(); 
-//     TTreeSRedirector &cstreamer = *fReconstructor->GetDebugStream(AliTRDReconstructor::kTracker);
-//     cstreamer << "AttachClustersIter"
-//       << "chamber.="   << &ch
-//       << "tracklet.="  << this
-//       << "\n";      
-//   }
-// 
-//   Float_t  tquality;
-//   Double_t kroady = fReconstructor->GetRecoParam() ->GetRoad1y();
-//   Double_t kroadz = GetPadLength() * .5 + 1.;
-//   
-//   // initialize configuration parameters
-//   Float_t zcorr = kZcorr ? GetTilt() * (fZfit[0] - fZref[0]) : 0.;
-//   Int_t   niter = kZcorr ? 1 : 2;
-//   
-//   Double_t yexp, zexp;
-//   Int_t ncl = 0;
-//   // start seed update
-//   for (Int_t iter = 0; iter < niter; iter++) {
-//     ncl = 0;
-//     for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
-//       if(!(layer = chamber->GetTB(iTime))) continue;
-//       if(!Int_t(*layer)) continue;
-//       
-//       // define searching configuration
-//       Double_t dxlayer = layer->GetX() - fX0;
-//       if(c){
-//         zexp = c->GetZ();
-//         //Try 2 pad-rows in second iteration
-//         if (iter > 0) {
-//           zexp = fZref[0] + fZref[1] * dxlayer - zcorr;
-//           if (zexp > c->GetZ()) zexp = c->GetZ() + GetPadLength()*0.5;
-//           if (zexp < c->GetZ()) zexp = c->GetZ() - GetPadLength()*0.5;
-//         }
-//       } else zexp = fZref[0] + (kZcorr ? fZref[1] * dxlayer : 0.);
-//       yexp  = fYref[0] + fYref[1] * dxlayer - zcorr;
-//       
-//       // Get and register cluster
-//       Int_t    index = layer->SearchNearestCluster(yexp, zexp, kroady, kroadz);
-//       if (index < 0) continue;
-//       AliTRDcluster *cl = (*layer)[index];
-//       
-//       fIndexes[iTime]  = layer->GetGlobalIndex(index);
-//       fClusters[iTime] = cl;
-// //       fY[iTime]        = cl->GetY();
-// //       fZ[iTime]        = cl->GetZ();
-//       ncl++;
-//     }
-//     if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7) AliInfo(Form("iter = %d ncl [%d] = %d", iter, fDet, ncl));
-//     
-//     if(ncl>1){      
-//       // calculate length of the time bin (calibration aware)
-//       Int_t irp = 0; Float_t x[2]={0., 0.}; Int_t tb[2] = {0,0};
-//       for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
-//         if(!fClusters[iTime]) continue;
-//         x[irp]  = fClusters[iTime]->GetX();
-//         tb[irp] = iTime;
-//         irp++;
-//         if(irp==2) break;
-//       } 
-//       Int_t dtb = tb[1] - tb[0];
-//       fdX = dtb ? (x[0] - x[1]) / dtb : 0.15;
-// 
-//       // update X0 from the clusters (calibration/alignment aware)
-//       for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
-//         if(!(layer = chamber->GetTB(iTime))) continue;
-//         if(!layer->IsT0()) continue;
-//         if(fClusters[iTime]){ 
-//           fX0 = fClusters[iTime]->GetX();
-//           break;
-//         } else { // we have to infere the position of the anode wire from the other clusters
-//           for (Int_t jTime = iTime+1; jTime < AliTRDtrackerV1::GetNTimeBins(); jTime++) {
-//             if(!fClusters[jTime]) continue;
-//             fX0 = fClusters[jTime]->GetX() + fdX * (jTime - iTime);
-//             break;
-//           }
-//         }
-//       }     
-//       
-//       // update YZ reference point
-//       // TODO
-//       
-//       // update x reference positions (calibration/alignment aware)
-// //       for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
-// //         if(!fClusters[iTime]) continue;
-// //         fX[iTime] = fX0 - fClusters[iTime]->GetX();
-// //       } 
-//       
-//       FitMI();
-//     }
-//     if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7) AliInfo(Form("iter = %d nclFit [%d] = %d", iter, fDet, fN2));
-//     
-//     if(IsOK()){
-//       tquality = GetQuality(kZcorr);
-//       if(tquality < quality) break;
-//       else quality = tquality;
-//     }
-//     kroadz *= 2.;
-//   } // Loop: iter
-//   if (!IsOK()) return kFALSE;
-// 
-//   if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=1) CookLabels();
-// 
-//   // load calibration params
-//   Calibrate();
-//   UpdateUsed();
-//   return kTRUE;     
-// }
-
 //____________________________________________________________________
 Bool_t AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *chamber, Bool_t tilt)
 {
@@ -1114,7 +997,7 @@ void AliTRDseedV1::Bootstrap(const AliTRDReconstructor *rec)
 
 
 //____________________________________________________________________
-Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
+Bool_t AliTRDseedV1::Fit(Bool_t tilt, Bool_t zcorr)
 {
   //
   // Linear fit of the tracklet
@@ -1130,10 +1013,7 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
   // 3. Do a Least Square Fit to the data
   //
 
-  if(!IsCalibrated()){
-    AliWarning("Tracklet fit failed. Call Calibrate().");
-    return kFALSE;
-  }
+  if(!IsCalibrated()) Calibrate();
 
   const Int_t kClmin = 8;
 
@@ -1151,22 +1031,6 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
     {3.439e-02,-3.601e-04, 3.883e-01, 2.623e-02},
     {3.510e-02, 2.066e-03, 3.651e-01, 2.588e-02},
   };
-  // 3. sy parallel to the track
-  const Float_t sy0 =  2.649e-02; // [cm]
-  const Float_t sya = -8.864e-04; // [cm]
-  const Float_t syb = -2.435e-01; // [cm]
-
-  // 4. sx parallel to the track
-  const Float_t sxgc = 5.427e-02;
-  const Float_t sxgm = 7.783e-01;
-  const Float_t sxgs = 2.743e-01;
-  const Float_t sxe0 =-2.065e+00;
-  const Float_t sxe1 =-2.978e-02;
-
-  // 5. sx perpendicular to the track
-//   const Float_t sxd0 = 1.881e-02;
-//   const Float_t sxd1 =-4.101e-01;
-//   const Float_t sxd2 = 1.572e+00;
 
   // get track direction
   Double_t y0   = fYref[0];
@@ -1181,8 +1045,7 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
 
   //AliTRDtrackerV1::AliTRDLeastSquare fitterZ;
   TLinearFitter  fitterY(1, "pol1");
-  // convertion factor from square to gauss distribution for sigma
-  //Double_t convert = 1./TMath::Sqrt(12.);
+  TLinearFitter  fitterZ(1, "pol1");
   
   // book cluster information
   Double_t qc[kNclusters], xc[kNclusters], yc[kNclusters], zc[kNclusters], sy[kNclusters];
@@ -1196,29 +1059,23 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
     yc[ic]  = 999.;
     zc[ic]  = 999.;
     sy[ic]  = 0.;
-    //sz[ic]  = 0.;
     if(!(c = (*jc))) continue;
     if(!c->IsInChamber()) continue;
 
     Float_t w = 1.;
     if(c->GetNPads()>4) w = .5;
     if(c->GetNPads()>5) w = .2;
+    Int_t tb = c->GetLocalTimeBin();
 
-    //zRow[fN] = c->GetPadRow();
     qc[n]   = TMath::Abs(c->GetQ());
-    // correct cluster position for PRF and v drift
+    // Radial cluster position
     //Int_t jc = TMath::Max(fN-3, 0);
     //xc[fN]   = c->GetXloc(fT0, fVD, &qc[jc], &xc[jc]/*, z0 - c->GetX()*dzdx*/);
-    //Double_t s2 = fS2PRF + fDiffL*fDiffL*xc[fN]/(1.+2.*exb2)+tgg*xc[fN]*xc[fN]*exb2/12.;
+    xc[n]   = fX0 - c->GetX();
+
+    //Double_t s2 = fS2PRF + fDiffL*fDiffL*xc[n]/(1.+2.*exb2)+tgg*xc[n]*xc[n]*exb2/12.;
     //yc[fN]   = c->GetYloc(s2, GetPadWidth(), xc[fN], fExB);
-    
-    // uncalibrated cluster correction 
-    // TODO remove
-    Double_t x, y; 
-    //GetClusterXY(c, x, y);
-    x = c->GetX(); y = c->GetY();
-    xc[n]   = fX0 - x;
-    yc[n]   = y;
+    yc[n]   = c->GetY()-AliTRDcluster::GetYcorr(ily, c->GetCenter());
     zc[n]   = c->GetZ();
 
     // extrapolated y value for the track
@@ -1229,9 +1086,8 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
     if(tilt) yc[n] -= GetTilt()*(zc[n] - zt); 
 
     // ELABORATE CLUSTER ERROR
-    // TODO to be moved to AliTRDcluster
     // basic y error (|| to track).
-    sy[n]  = xc[n] < AliTRDgeometry::CamHght() ? 2. : sy0 + sya*TMath::Exp(1./(xc[n]+syb));
+    sy[n]  = AliTRDcluster::GetSY(tb, zcorr?zt:-1.);
     //printf("cluster[%d]\n\tsy[0] = %5.3e [um]\n", fN,  sy[fN]*1.e4);
     // y error due to total charge
     sy[n] += sqb*(1./qc[n] - sq0inv);
@@ -1244,26 +1100,25 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
 
     // ADD ERROR ON x
     // error of drift length parallel to the track
-    Double_t sx = sxgc*TMath::Gaus(xc[n], sxgm, sxgs) + TMath::Exp(sxe0+sxe1*xc[n]); // [cm]
+    Double_t sx = AliTRDcluster::GetSX(tb, zcorr?zt:-1.); // [cm]
     //printf("\tsx[0] = %5.3e [um]\n", sx*1.e4);
-    // error of drift length perpendicular to the track
-    //sx += sxd0 + sxd1*d + sxd2*d*d;
     sx *= sx; // square sx
 
     // add error from ExB 
-    if(errors>0) sy[n] += fExB*fExB*sx;
+    sy[n] += fExB*fExB*sx;
     //printf("\tsy[3] = %5.3e [um^2]\n", sy[fN]*1.e8);
 
     // global radial error due to misalignment/miscalibration
     Double_t sx0  = 0.; sx0 *= sx0;
     // add sx contribution to sy due to track angle
-    if(errors>1) sy[n] += tgg*(sx+sx0);
+    sy[n] += tgg*(sx+sx0);
     // TODO we should add tilt pad correction here
     //printf("\tsy[4] = %5.3e [um^2]\n", sy[fN]*1.e8);
     c->SetSigmaY2(sy[n]);
 
     sy[n]  = TMath::Sqrt(sy[n]);
     fitterY.AddPoint(&xc[n], yc[n], sy[n]);
+    fitterZ.AddPoint(&xc[n], qc[n], 1.);
     n++;
   }
   // to few clusters
@@ -1280,61 +1135,36 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
   fCov[2] = p[3]; // variance of dydx
   // the ref radial position is set at the minimum of 
   // the y variance of the tracklet
-  fX   = -fCov[1]/fCov[2]; //fXref = fX0 - fXref;
-  fS2Y = fCov[0] +2.*fX*fCov[1] + fX*fX*fCov[2];
+  fX   = -fCov[1]/fCov[2];
 
   // fit XZ
-  if(IsRowCross()){ 
-    // TODO pad row cross position estimation !!!
-    //AliInfo(Form("Padrow cross in detector %d", fDet));
-    fZfit[0] = .5*(zc[0]+zc[n-1]); fZfit[1] = 0.;
+  if(IsRowCross()){
+    Int_t ic=n=kNclusters-1; jc = &fClusters[ic];
+    for(; ic>kNtb; ic--, --jc){
+      if(!(c = (*jc))) continue;
+      if(!c->IsInChamber()) continue;
+      qc[n]   = TMath::Abs(c->GetQ());
+      xc[n]   = fX0 - c->GetX();
+      zc[n]   = c->GetZ();
+      fitterZ.AddPoint(&xc[n], -qc[n], 1.);
+      n--;
+    }
+    // fit XZ
+    fitterZ.Eval();
+    if(fitterZ.GetParameter(1)!=0.){ 
+      fX = -fitterZ.GetParameter(0)/fitterZ.GetParameter(1);
+      fX=(fX<0.)?0.:fX;
+      Float_t dl = .5*AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght();
+      fX=(fX> dl)?dl:fX;
+    }
+
+    fZfit[0] = .5*(zc[0]+zc[kNclusters-1]); fZfit[1] = 0.;
     fS2Z     = 0.02+1.55*fZref[1]; fS2Z *= fS2Z;
   } else {
     fZfit[0] = zc[0]; fZfit[1] = 0.;
     fS2Z     = GetPadLength()*GetPadLength()/12.;
   }
-
-
-//   // determine z offset of the fit
-//   Float_t zslope = 0.;
-//   Int_t nchanges = 0, nCross = 0;
-//   if(nz==2){ // tracklet is crossing pad row
-//     // Find the break time allowing one chage on pad-rows
-//     // with maximal number of accepted clusters
-//     Int_t padRef = zRow[0];
-//     for (Int_t ic=1; ic<fN; ic++) {
-//       if(zRow[ic] == padRef) continue;
-//       
-//       // debug
-//       if(zRow[ic-1] == zRow[ic]){
-//         printf("ERROR in pad row change!!!\n");
-//       }
-//     
-//       // evaluate parameters of the crossing point
-//       Float_t sx = (xc[ic-1] - xc[ic])*convert;
-//       fCross[0] = .5 * (xc[ic-1] + xc[ic]);
-//       fCross[2] = .5 * (zc[ic-1] + zc[ic]);
-//       fCross[3] = TMath::Max(dzdx * sx, .01);
-//       zslope    = zc[ic-1] > zc[ic] ? 1. : -1.;
-//       padRef    = zRow[ic];
-//       nCross    = ic;
-//       nchanges++;
-//     }
-//   }
-// 
-//   // condition on nCross and reset nchanges TODO
-// 
-//   if(nchanges==1){
-//     if(dzdx * zslope < 0.){
-//       AliInfo("Tracklet-Track mismatch in dzdx. TODO.");
-//     }
-// 
-// 
-//     //zc[nc] = fitterZ.GetFunctionParameter(0); 
-//     fCross[1] = fYfit[0] - fCross[0] * fYfit[1];
-//     fCross[0] = fX0 - fCross[0];
-//   }
-
+  fS2Y = fCov[0] +2.*fX*fCov[1] + fX*fX*fCov[2];
   return kTRUE;
 }
 
@@ -1598,8 +1428,9 @@ void AliTRDseedV1::Print(Option_t *o) const
   // Printing the seedstatus
   //
 
-  AliInfo(Form("Det[%3d] X0[%7.2f] Pad[L[%5.2f] W[%5.2f] Tilt[%+6.2f]]", fDet, fX0, GetPadLength(), GetPadWidth(), GetTilt()));
+  AliInfo(Form("Det[%3d] X0[%7.2f] Pad{L[%5.2f] W[%5.2f] Tilt[%+6.2f]}", fDet, fX0, GetPadLength(), GetPadWidth(), GetTilt()));
   AliInfo(Form("N[%2d] Nused[%2d] Nshared[%2d] [%d]", GetN(), GetNUsed(), GetNShared(), fN));
+  AliInfo(Form("FLAGS : RC[%c] Kink[%c] SA[%c]", IsRowCross()?'y':'n', IsKink()?'y':'n', IsStandAlone()?'y':'n'));
 
   Double_t cov[3], x=GetX();
   GetCovAt(x, cov);
index 5645a35..c654d03 100644 (file)
@@ -76,8 +76,7 @@ public:
   void      CookdEdx(Int_t nslices);
   void      CookLabels();
   Bool_t    CookPID();
-  Bool_t    Fit(Bool_t tilt=kTRUE, Int_t errors = 2);
-//   void      FitMI();
+  Bool_t    Fit(Bool_t tilt=kFALSE, Bool_t zcorr=kFALSE);
   Bool_t    Init(AliTRDtrackV1 *track);
   inline void      Init(const AliRieman *fit);
   Bool_t    IsEqual(const TObject *inTracklet) const;
@@ -131,6 +130,7 @@ public:
   Float_t   GetTgl() const           { return fZref[1];}
   Float_t   GetTilt() const          { return fPad[2];}
   UInt_t    GetTrackletWord() const  { return 0;}
+  UShort_t  GetVolumeId() const;
   Float_t   GetX0() const            { return fX0;}
   Float_t   GetX() const             { return fX0 - fX;}
   Float_t   GetY() const             { return fYfit[0] - fYfit[1] * fX;}
index cd86f22..32138cb 100644 (file)
@@ -400,7 +400,7 @@ Int_t  AliTRDtrackV1::GetClusterIndex(Int_t id) const
 }
 
 //_______________________________________________________________
-Double_t AliTRDtrackV1::GetPredictedChi2(const AliTRDseedV1 *trklt) const
+Double_t AliTRDtrackV1::GetPredictedChi2(const AliTRDseedV1 *trklt, Double_t *cov) const
 {
 // Compute chi2 between tracklet and track. The value is calculated at the radial position of the track
 // equal to the reference radial position of the tracklet (see AliTRDseedV1)
@@ -412,9 +412,8 @@ Double_t AliTRDtrackV1::GetPredictedChi2(const AliTRDseedV1 *trklt) const
 // where X=(y z), the position of the track/tracklet in the yz plane
 // 
 
-  Double_t x = GetX();
-  Double_t p[2]   = { trklt->GetYat(x), trklt->GetZat(x)};
-  Double_t cov[3]; trklt->GetCovAt(x, cov);
+  Double_t p[2]   = { trklt->GetY(), trklt->GetZ()};
+  trklt->GetCovAt(trklt->GetX(), cov);
   return AliExternalTrackParam::GetPredictedChi2(p, cov);
 }
 
@@ -774,45 +773,16 @@ void AliTRDtrackV1::UnsetTracklet(Int_t plane)
 
 
 //_______________________________________________________________
-Bool_t  AliTRDtrackV1::Update(AliTRDseedV1 *trklt, Double_t chisq)
+Bool_t  AliTRDtrackV1::Update(Double_t *p, Double_t *cov, Double_t chi2)
 {
   //
-  // Update track and tracklet parameters 
+  // Update track 
   //
-  
-  Double_t x      = GetX();
-  Double_t p[2]   = { trklt->GetYat(x)
-                    , trklt->GetZat(x) };
-  Double_t cov[3]/*, covR[3], cov0[3]*/;
-
-//   printf("\tD[%3d] Ly[%d] Trk: x[%f] y[%f] z[%f]\n", trklt->GetDetector(), trklt->GetPlane(), GetX(), GetY(), GetZ());
-// //   
-//   Double_t xref = trklt->GetXref();
-//   trklt->GetCovAt(xref, covR);
-//   printf("xr=%5.3f y=%f+-%f z=%f+-%f (covYZ=%f)\n", xref, trklt->GetYat(xref), TMath::Sqrt(covR[0]), trklt->GetZat(xref), TMath::Sqrt(covR[2]), covR[1]);
-// 
-//   Double_t x0 = trklt->GetX0();
-//   trklt->GetCovAt(x0, cov0);
-//   printf("x0=%5.3f y=%f+-%f z=%f+-%f (covYZ=%f)\n", x0, trklt->GetYat(x0), TMath::Sqrt(cov0[0]), trklt->GetZat(x0), TMath::Sqrt(cov0[2]), cov0[1]);
-// 
-//   trklt->GetCovAt(x, cov);
-//   printf("x =%5.3f y=%f+-%f z=%f+-%f (covYZ=%f)\n", x, p[0], TMath::Sqrt(cov[0]), p[1], TMath::Sqrt(cov[2]), cov[1]);
-// 
-//   const Double_t *cc = GetCovariance();
-//   printf("yklm[0] = %f +- %f\n", GetY(), TMath::Sqrt(cc[0]));
-
-  trklt->GetCovAt(x, cov); 
   if(!AliExternalTrackParam::Update(p, cov)) return kFALSE;
-//   cc = GetCovariance();
-//   printf("yklm[1] = %f +- %f\n", GetY(), TMath::Sqrt(cc[0]));
 
-  AliTRDcluster *c = 0x0;
-  Int_t ic = 0; while(!(c = trklt->GetClusters(ic))) ic++;
-  AliTracker::FillResiduals(this, p, cov, c->GetVolumeId());
-  
   // Register info to track
   SetNumberOfClusters();
-  SetChi2(GetChi2() + chisq);
+  SetChi2(GetChi2() + chi2);
   return kTRUE;
 }
 
index 12879ab..eafe6fd 100644 (file)
@@ -89,7 +89,7 @@ public:
   Double_t       GetPIDsignal() const   { return 0.;}
   Double_t       GetPID(Int_t is) const { return (is >=0 && is < AliPID::kSPECIES) ? fPID[is] : -1.;}
   UChar_t        GetNumberOfTrackletsPID() const;
-  Double_t       GetPredictedChi2(const AliTRDseedV1 *tracklet) const;
+  Double_t       GetPredictedChi2(const AliTRDseedV1 *tracklet, Double_t *cov) const;
   Double_t       GetPredictedChi2(const AliCluster* /*c*/) const                   { return 0.0; }
   Int_t          GetProlongation(Double_t xk, Double_t &y, Double_t &z);
   inline UChar_t GetStatusTRD(Int_t ly=-1) const;
@@ -132,8 +132,7 @@ public:
   inline void    SetReconstructor(const AliTRDReconstructor *rec);
   inline Float_t StatusForTOF();
   void           UnsetTracklet(Int_t plane);
-  Bool_t         Update(AliTRDseedV1 *tracklet, Double_t chi2);
-  //Bool_t         Update(const AliTRDcluster *c, Double_t chi2, Int_t index, Double_t h01){ return AliTRDtrack::Update(c,chi2,index,h01); };
+  Bool_t         Update(Double_t *p, Double_t *cov, Double_t chi2);
   Bool_t         Update(const AliCluster *, Double_t, Int_t)                        { return kFALSE; };
   void           UpdateESDtrack(AliESDtrack *t);
 
index 2236e1b..73547b8 100644 (file)
@@ -598,8 +598,10 @@ Int_t AliTRDtrackerV1::FollowProlongation(AliTRDtrackV1 &t)
       kStoreIn = kFALSE;
     }
 
-    Double_t maxChi2 = t.GetPredictedChi2(tracklet);
-    if (maxChi2 < 1e+10 && t.Update(tracklet, maxChi2)){ 
+    Double_t cov[3]; tracklet->GetCovAt(x, cov);
+    Double_t p[2] = { tracklet->GetY(), tracklet->GetZ()};
+    Double_t chi2 = ((AliExternalTrackParam)t).GetPredictedChi2(p, cov);
+    if (chi2 < 1e+10 && t.Update(p, cov, chi2)){ 
       nClustersExpected += tracklet->GetN();
     }
   }
@@ -841,16 +843,21 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
     }
 
     // update Kalman with the TRD measurement
-    Double_t chi2 = t.GetPredictedChi2(ptrTracklet);
-    if(chi2>1e+10){
+    Double_t cov[3]; ptrTracklet->GetCovAt(x, cov);
+    Double_t p[2] = { ptrTracklet->GetY(), ptrTracklet->GetZ()};
+    Double_t chi2 = ((AliExternalTrackParam)t).GetPredictedChi2(p, cov);
+    if(chi2>1e+10){ // TODO
       t.SetStatus(AliTRDtrackV1::kChi2, ily);
       continue; 
     }
-    if(!t.Update(ptrTracklet, chi2)) {
+    if(!t.Update(p, cov, chi2)) {
       n=-1; 
       t.SetStatus(AliTRDtrackV1::kUpdate);
       break;
     }
+    // fill residuals ?!
+    AliTracker::FillResiduals(&t, p, cov, ptrTracklet->GetVolumeId());
+  
 
     // load tracklet to the tracker
     ptrTracklet->UpDate(&t);
@@ -1542,8 +1549,10 @@ Double_t AliTRDtrackerV1::FitKalman(AliTRDtrackV1 *track, AliTRDseedV1 *tracklet
     if (!AdjustSector(track)) break;
   
     //Update track
-    Double_t chi2 = track->GetPredictedChi2(ptrTracklet);
-    if(chi2<1e+10) track->Update(ptrTracklet, chi2);
+    Double_t cov[3]; ptrTracklet->GetCovAt(x, cov);
+    Double_t p[2] = { ptrTracklet->GetY(), ptrTracklet->GetZ()};
+    Double_t chi2 = ((AliExternalTrackParam*)track)->GetPredictedChi2(p, cov);
+    if(chi2<1e+10) track->Update(p, cov, chi2);
     if(!up) continue;
 
                //Reset material budget if 2 consecutive gold
@@ -2615,7 +2624,7 @@ Int_t AliTRDtrackerV1::MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *ss
               << "\n";
         }
               
-        if(fReconstructor->GetRecoParam()->HasImproveTracklets() && ImproveSeedQuality(stack, cseed) < 4){
+        if(fReconstructor->HasImproveTracklets() && ImproveSeedQuality(stack, cseed) < 4){
           AliTRDtrackerDebug::SetCandidateNumber(AliTRDtrackerDebug::GetCandidateNumber() + 1);
           continue;
         }
@@ -2629,7 +2638,7 @@ Int_t AliTRDtrackerV1::MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *ss
         // do the final track fitting (Once with vertex constraint and once without vertex constraint)
         Double_t chi2Vals[3];
         chi2Vals[0] = FitTiltedRieman(&cseed[0], kFALSE);
-        if(fReconstructor->GetRecoParam()->IsVertexConstrained())
+        if(fReconstructor->HasVertexConstrained())
           chi2Vals[1] = FitTiltedRiemanConstraint(&cseed[0], GetZ()); // Do Vertex Constrained fit if desired
         else
           chi2Vals[1] = 1.;
@@ -2887,7 +2896,7 @@ Double_t AliTRDtrackerV1::CalculateTrackLikelihood(AliTRDseedV1 *tracklets, Doub
   chi2phi /= Float_t (nLayers - 2.0);
   
   Double_t likeChi2Z  = TMath::Exp(-chi2[2] * 0.14);                   // Chi2Z 
-  Double_t likeChi2TC = (fReconstructor->GetRecoParam() ->IsVertexConstrained()) ? 
+  Double_t likeChi2TC = (fReconstructor->HasVertexConstrained()) ? 
                                                                                        TMath::Exp(-chi2[1] * 0.677) : 1;                       // Constrained Tilted Riemann
   Double_t likeChi2TR = TMath::Exp(-chi2[0] * 0.78);                   // Non-constrained Tilted Riemann
   Double_t likeChi2Phi= TMath::Exp(-chi2phi * 3.23);
index f0cfd9f..416da67 100644 (file)
@@ -201,7 +201,7 @@ void AliTRDtransform::SetDetector(Int_t det)
 }
 
 //_____________________________________________________________________________
-Bool_t AliTRDtransform::Transform(Double_t *x, Int_t *i, UInt_t time, Bool_t &out, Int_t  /*coordinateType*/)
+Bool_t AliTRDtransform::Transform(AliTRDcluster *c)
 {
   //
   // Transforms the local cluster coordinates into calibrated 
@@ -209,110 +209,54 @@ Bool_t AliTRDtransform::Transform(Double_t *x, Int_t *i, UInt_t time, Bool_t &ou
   //
   // Here the calibration for T0, Vdrift and ExB is applied as well.
   //
-  // Input:
-  //   x[0] = COL-position relative to the center pad (pad units)
-  //   x[1] = cluster signal in left pad
-  //   x[2] = cluster signal in middle pad
-  //   x[3] = cluster signal in right pad
-  //   i[0] = ROW pad number
-  //   i[1] = COL pad number
-  //   time = time bin number (uncalibrated for t0)
-  //
-  // Output:
-  //   x[0] = X-positions in tracking CS
-  //   x[1] = Y-positions in tracking CS
-  //   x[2] = Z-positions in tracking CS
-  //   x[3] = total cluster charge
-  //   x[4] = error in Y-direction
-  //   x[5] = error in Z-direction
-  //   i[2] = time bin number (calibrated for t0)
-  //
+  // Input: Cluster in the local chamber coordinates
+  // Output: Tracking cluster
 
-  Double_t posLocal[3];
-  Double_t posTracking[3];
+  if (!fMatrix) return kFALSE;
 
-  Int_t row = i[0];
-  Int_t col = i[1];
 
-  // Parameters to adjust the X position of clusters
+  // Parameters to adjust the X position of clusters in the alignable volume
   const Double_t kX0shift = AliTRDgeometry::AnodePos(); //[cm]
-  // TRF rising time (fitted)
-  // It should be absorbed by the t0. For the moment t0 is 0 for simulations.
-  // A.Bercuci (Mar 26 2009)
-  const Double_t kT0shift = 0.189;        //[us]
-
-  if (!fMatrix) {
-
-    x[0] = 0.0;
-    x[1] = 0.0;
-    x[2] = 0.0;
-    x[3] = 0.0;
-    x[4] = 0.0;
-    x[5] = 0.0;
-    i[2] = 0;
 
-    return kFALSE;
-
-  }
-  else {
  
-    // Calibration values
-    Double_t vdrift      = fCalVdriftDetValue * fCalVdriftROC->GetValue(col,row);
-    Double_t t0          = fCalT0DetValue     + fCalT0ROC->GetValue(col,row);
-
-    // T0 correction
-    Double_t timeT0Cal   = time - t0;
-    // Calculate the X-position,
-    Double_t xLocal      = ((timeT0Cal + 0.5) / fSamplingFrequency - kT0shift) * vdrift + AliTRDcluster::GetXcorr(TMath::Nint(timeT0Cal)); 
-
-    // Length of the amplification region
-    //Double_t ampLength   = (Double_t) AliTRDgeometry::CamHght();
-    // The drift distance
-    Double_t driftLength = TMath::Max(xLocal,0.0);
-    // ExB correction
-    Double_t exbCorr     = AliTRDCommonParam::Instance()->GetOmegaTau(vdrift);
-
-    // Pad dimensions
-    Double_t rowSize     = fPadPlane->GetRowSize(row);
-    Double_t colSize     = fPadPlane->GetColSize(col);
-
-    // Invert the X-position,
-    // apply ExB correction to the Y-position
-    // and move to the Z-position relative to the middle of the chamber
-    posLocal[0] = kX0shift-xLocal;
-    posLocal[1] = (fPadPlane->GetColPos(col) + (0.5 + x[0]) * colSize) - driftLength * exbCorr;
-    posLocal[2] = (fPadPlane->GetRowPos(row) -  0.5  * rowSize) - fZShiftIdeal;
-
-    // Go to tracking coordinates
-    fMatrix->LocalToMaster(posLocal,posTracking);
-
-    // The total charge of the cluster
-    Double_t q0             = x[1];
-    Double_t q1             = x[2];
-    Double_t q2             = x[3];
-    Double_t clusterCharge  = q0 + q1 + q2; 
-    Double_t clusterSigmaY2 = 0.0;
-    if (clusterCharge > 0.0) {
-      clusterSigmaY2 = (q1 * (q0 + q2) + 4.0 * q0 * q2)
-                    / (clusterCharge*clusterCharge);
-    }
-
-    // Output values
-    x[0] = posTracking[0];
-    x[1] = posTracking[1];
-    x[2] = posTracking[2];
-    x[3] = clusterCharge;
-    x[4] = colSize*colSize * (clusterSigmaY2 + 1.0/12.0);
-    x[5] = rowSize*rowSize / 12.0;                                       
-    i[2] = TMath::Nint(timeT0Cal);
-               
-    // For TRD tracking calibration awareness
-    out  = ((i[2] < 0) || (i[2] > Int_t(3.5 * fSamplingFrequency/vdrift))) ? kTRUE : kFALSE; 
-
-    return kTRUE;
-
-  }
-
+  // Retrieve calibration values
+  Int_t col = c->GetPadCol(), row = c->GetPadRow();
+  // drift velocity
+  Double_t vd  = fCalVdriftDetValue * fCalVdriftROC->GetValue(col,row);
+  // t0
+  Double_t t0  = fCalT0DetValue     + fCalT0ROC->GetValue(col,row);
+  // ExB correction
+  Double_t exb = AliTRDCommonParam::Instance()->GetOmegaTau(vd);
+
+  Float_t x = c->GetXloc(t0, vd);
+
+  // pad response width with diffusion corrections
+  Double_t s2  = AliTRDcalibDB::Instance()->GetPRFWidth(fDetector, col, row); s2 *= s2; 
+  // Float_t dl, dt;
+  // AliTRDCommonParam::Instance()->GetDiffCoeff(dl, dt, vd);
+  // s2 += dl*dl*x/(1.+2.*exb*exb);
+  s2 -= - 1.5e-1;
+
+  // Pad dimensions
+  Double_t rs = fPadPlane->GetRowSize(row);
+  Double_t cs = fPadPlane->GetColSize(col);
+  Double_t y0 = fPadPlane->GetColPos(col) + .5*cs;
+  Double_t loc[] = {
+    kX0shift-x,                    // Invert the X-position,
+    c->GetYloc(y0, s2, cs) - x*exb,// apply ExB correction
+    fPadPlane->GetRowPos(row) - .5*rs - fZShiftIdeal // move the Z-position relative to the middle of the chamber
+  };
+
+  // Go to tracking coordinates
+  Double_t trk[3];
+  fMatrix->LocalToMaster(loc, trk);
+
+  // store tracking values
+  c->SetX(trk[0]);c->SetY(trk[1]);c->SetZ(trk[2]);
+  c->SetSigmaY2(s2);
+  c->SetSigmaZ2(fPadPlane->GetRowSize(row)*fPadPlane->GetRowSize(row)/12.);
+  
+  return kTRUE;
 }
 
 //_____________________________________________________________________________
@@ -325,33 +269,6 @@ void AliTRDtransform::Recalibrate(AliTRDcluster *c, Bool_t setDet)
   // be used.
   //
 
-  if (setDet) {
-    SetDetector(c->GetDetector());
-  }
-
-  // Transform the local cluster coordinates into recalibrated 
-  // space point positions defined in the local tracking system.
-  // Here the calibration for T0, Vdrift and ExB is applied as well.
-  Double_t clusterXYZ[6];
-  clusterXYZ[0] = c->GetCenter();
-  clusterXYZ[1] = 0.0;
-  clusterXYZ[2] = 0.0;
-  clusterXYZ[3] = 0.0;
-  clusterXYZ[4] = 0.0;
-  clusterXYZ[5] = 0.0;
-  Int_t    clusterRCT[3];
-  clusterRCT[0] = c->GetPadRow();
-  clusterRCT[1] = c->GetPadCol();
-  clusterRCT[2] = 0;
-  Int_t time    = c->GetPadTime();
-  Bool_t out;
-  Transform(clusterXYZ,clusterRCT,((UInt_t) time), out, 0);
-
-  // Set the recalibrated coordinates
-  c->SetX(clusterXYZ[0]);
-  c->SetY(clusterXYZ[1]);
-  c->SetZ(clusterXYZ[2]);
-  c->SetLocalTimeBin(((Char_t) clusterRCT[2]));
-  c->SetInChamber(!out);
-
+  if (setDet) SetDetector(c->GetDetector());
+  Transform(c);
 }
index 6cac4de..65a3d12 100644 (file)
@@ -37,11 +37,7 @@ class AliTRDtransform : public TObject {
   AliTRDtransform &operator=(const AliTRDtransform &t) { *(new(this) AliTRDtransform(t));
                                                           return *this; }
 
-  virtual Bool_t   Transform(Double_t *x
-                           , Int_t    *i
-                           , UInt_t    time
-                           , Bool_t   &out
-                           , Int_t     coordinateType);
+  virtual Bool_t   Transform(AliTRDcluster *c);
   virtual void     Recalibrate(AliTRDcluster *c, Bool_t setDet = kTRUE);
 
           void     SetDetector(Int_t det);
index 92d410e..a3d0ed6 100644 (file)
@@ -71,6 +71,8 @@
 #include "AliTrackPointArray.h"
 #include "AliCDBManager.h"
 
+#include "AliTRDcalibDB.h"
+#include "AliTRDCommonParam.h"
 #include "AliTRDSimParam.h"
 #include "AliTRDgeometry.h"
 #include "AliTRDpadPlane.h"
@@ -384,6 +386,12 @@ TH1* AliTRDresolution::PlotTrackIn(const AliTRDtrackV1 *track)
       COV(ir,ic) = (*pc); COV(ic,ir) = (*pc);
     }
   }
+//   TMatrixDSymEigen eigen(COV);
+//   TVectorD evals = eigen.GetEigenValues();
+//   TMatrixDSym evalsm(kNPAR);
+//   for(Int_t ir=0; ir<kNPAR; ir++) for(Int_t ic=0; ic<kNPAR; ic++) evalsm(ir,ic) = (ir==ic ? evals(ir): 0.);
+//   TMatrixD evecs = eigen.GetEigenVectors();
+//   TMatrixD sqrcov(evecs, TMatrixD::kMult, TMatrixD(evalsm, TMatrixD::kMult, evecs.T()));
 
   // calculate delta
   parR[0]-=y0;parR[1]-=z0;
@@ -513,7 +521,7 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
     AliTRDseedV1 tt(*fTracklet);
     tt.SetZref(0, z0 - (x0-xAnode)*dzdx0);
     tt.SetZref(1, dzdx0); 
-    tt.Fit(kTRUE);
+    tt.Fit(kTRUE, kTRUE);
     x= tt.GetX();y= tt.GetY();z= tt.GetZ();
     dydx = tt.GetYfit(1);
     dx = x0 - x;
@@ -557,13 +565,20 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
     AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
     Float_t zr0 = pp->GetRow0() + pp->GetAnodeWireOffset();
     Float_t tilt = fTracklet->GetTilt();
+    Double_t exb = AliTRDCommonParam::Instance()->GetOmegaTau(1.5);
 
     AliTRDcluster *c = 0x0;
     fTracklet->ResetClusterIter(kFALSE);
     while((c = fTracklet->PrevCluster())){
       Float_t  q = TMath::Abs(c->GetQ());
-      //AliTRDseedV1::GetClusterXY(c,x,y);
       x = c->GetX(); y = c->GetY();
+//       Int_t col = c->GetPadCol();
+//       Int_t row = c->GetPadRow();
+//       Double_t pw = pp->GetColSize(col);
+//       Double_t y0 = (pp->GetColPos(col) + 0.5) * pw;
+//       Double_t s2 = AliTRDcalibDB::Instance()->GetPRFWidth(det, col, row); s2 *= s2; s2 -= - 1.5e-1;
+//       y = c->GetYloc(y0, s2, pw); y-=(xAnode-x)*exb;
+
       z = c->GetZ();
       dx = x0 - x; 
       yt = y0 - dx*dydx0;
@@ -995,7 +1010,7 @@ TObjArray* AliTRDresolution::Histos()
 
   // tracklet y resolution [0]
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltZ"))){
-    h = new TH2I("hMCtrkltZ", "Tracklet Resolution (Z)", 31, -.48, .48, 100, -.5, .5);
+    h = new TH2I("hMCtrkltZ", "Tracklet Resolution (Z)", 50, -1., 1., 100, -.5, .5);
     h->GetXaxis()->SetTitle("tg(#theta)");
     h->GetYaxis()->SetTitle("#Delta z [cm]");
     h->GetZaxis()->SetTitle("entries");
@@ -1004,7 +1019,7 @@ TObjArray* AliTRDresolution::Histos()
 
   // tracklet y resolution [0]
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltZPull"))){
-    h = new TH2I("hMCtrkltZPull", "Tracklet Pulls (Z)", 31, -.48, .48, 100, -3.5, 3.5);
+    h = new TH2I("hMCtrkltZPull", "Tracklet Pulls (Z)", 31, -1., 1., 100, -3.5, 3.5);
     h->GetXaxis()->SetTitle("tg(#theta)");
     h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
     h->GetZaxis()->SetTitle("entries");
@@ -1029,7 +1044,7 @@ TObjArray* AliTRDresolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackY);
 
-  // Kalman track y resolution
+  // Kalman track y pulls
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkYPull"))){
     h = new TH2I("hMCtrkYPull", "Kalman Track Pulls (Y)", 31, -.48, .48, 100, -3.5, 3.5);
     h->GetXaxis()->SetTitle("tg(#phi)");
@@ -1038,7 +1053,7 @@ TObjArray* AliTRDresolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackYPull);
 
-  // Kalman track Z resolution
+  // Kalman track Z resolution [IN]
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZIn"))){
     h = new TH2I("hMCtrkZIn", "Kalman Track Resolution (Zin)", 20, -1., 1., 100, -1., 1.);
     h->GetXaxis()->SetTitle("tg(#theta)");
@@ -1047,7 +1062,7 @@ TObjArray* AliTRDresolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackZIn);
 
-  // Kalman track Z resolution
+  // Kalman track Z resolution [OUT]
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZOut"))){
     h = new TH2I("hMCtrkZOut", "Kalman Track Resolution (Zout)", 20, -1., 1., 100, -1., 1.);
     h->GetXaxis()->SetTitle("tg(#theta)");
@@ -1056,7 +1071,7 @@ TObjArray* AliTRDresolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackZOut);
 
-  // Kalman track Z resolution
+  // Kalman track Z pulls
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZInPull"))){
     h = new TH2I("hMCtrkZInPull", "Kalman Track Pulls (Zin)", 20, -1., 1., 100, -4.5, 4.5);
     h->GetXaxis()->SetTitle("tg(#theta)");
@@ -1065,7 +1080,7 @@ TObjArray* AliTRDresolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackZInPull);
 
-  // Kalman track Z resolution
+  // Kalman track Z pulls
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZOutPull"))){
     h = new TH2I("hMCtrkZOutPull", "Kalman Track Pulls (Zout)", 20, -1., 1., 100, -4.5, 4.5);
     h->GetXaxis()->SetTitle("tg(#theta)");
@@ -1076,7 +1091,7 @@ TObjArray* AliTRDresolution::Histos()
 
   // Kalman track Pt resolution
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkPt"))){
-    h = new TH2I("hMCtrkPt", "Kalman Track Resolution (Pt)", 100, 0., 2., 150, -5., 20.);
+    h = new TH2I("hMCtrkPt", "Kalman Track Resolution (Pt)", 100, 0., 2., 150, -1., 10.);
     h->GetXaxis()->SetTitle("1/p_{t}");
     h->GetYaxis()->SetTitle("#Delta p_{t} [GeV/c]");
     h->GetZaxis()->SetTitle("entries");