Enlarging the recoParam functionality and populating the track list (needed by HLT)
authorcblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 2 Jun 2008 15:36:18 +0000 (15:36 +0000)
committercblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 2 Jun 2008 15:36:18 +0000 (15:36 +0000)
TRD/AliTRDReconstructor.cxx
TRD/AliTRDReconstructor.h
TRD/AliTRDchamberTimeBin.cxx
TRD/AliTRDrecoParam.cxx
TRD/AliTRDrecoParam.h
TRD/AliTRDseedV1.cxx
TRD/AliTRDtracker.cxx
TRD/AliTRDtrackerV1.cxx
TRD/AliTRDtrackerV1.h
TRD/AliTRDtrackingChamber.cxx

index 4769387..8f39bb9 100644 (file)
 
 ClassImp(AliTRDReconstructor)
 
-Bool_t AliTRDReconstructor::fgkSeedingOn  = kFALSE;
-Int_t  AliTRDReconstructor::fgStreamLevel = 0;      // Stream (debug) level
 AliTRDrecoParam* AliTRDReconstructor::fgRecoParam = 0x0; 
 
-
 //_____________________________________________________________________________
 AliTRDReconstructor::~AliTRDReconstructor()                  { 
        if(fgRecoParam) delete fgRecoParam;
index 7103141..b887ecf 100644 (file)
@@ -40,17 +40,11 @@ class AliTRDReconstructor: public AliReconstructor
                                                                  , esd);                    }
   virtual void             FillESD(TTree *digitsTree, TTree *clusterTree, AliESDEvent *esd) const;
 
-  static  void             SetSeedingOn(Bool_t seeding)          { fgkSeedingOn  = seeding; }  
-  static  void             SetStreamLevel(Int_t level)           { fgStreamLevel = level;   }
   static  void             SetRecoParam(AliTRDrecoParam *reco)   { fgRecoParam   = reco;    }
 
-  static  Bool_t           SeedingOn()                           { return fgkSeedingOn;     }
-  static  Int_t            StreamLevel()                         { return fgStreamLevel;    }
 
  private:
 
-  static  Bool_t           fgkSeedingOn;  //  Set flag for seeding during reconstruction
-  static  Int_t            fgStreamLevel; //  Flag for streaming
   static  AliTRDrecoParam *fgRecoParam;   //  Reconstruction parameters
 
   ClassDef(AliTRDReconstructor,0)         //  Class for the TRD reconstruction
index 7e0fced..31b4f08 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-
-/* $Id: AliTRDchamberTimeBin.cxx 23313 2008-01-11 14:56:43Z cblume $ */
-
-///////////////////////////////////////////////////////////////////////////////
-//                                                                           //
-//  Organization of clusters at the level of 1 TRD chamber.                  //
-//  The data structure is used for tracking at the stack level.              //
-//                                                                           //
-//  Functionalities:                                                         //
-//   1. cluster organization and sorting                                     //
-//   2. fast data navigation                                                 //
-//                                                                           //
-//  Authors:                                                                 //
-//    Alex Bercuci <A.Bercuci@gsi.de>                                        //
-//    Markus Fasel <M.Fasel@gsi.de>                                          //
-//                                                                           //
-///////////////////////////////////////////////////////////////////////////////
-
-#include <TObject.h>
-#include <TROOT.h>
-#include <TMath.h>
-#include <TStopwatch.h>
-#include <TTreeStream.h>
-
-#include "AliLog.h"
-
-#include "AliTRDcluster.h"
-#include "AliTRDchamberTimeBin.h"
-#include "AliTRDrecoParam.h"
-#include "AliTRDReconstructor.h"
-#include "AliTRDtrackerV1.h"
-
-
-ClassImp(AliTRDchamberTimeBin)
-
-//_____________________________________________________________________________
-AliTRDchamberTimeBin::AliTRDchamberTimeBin(Int_t plane, Int_t stack, Int_t sector, Double_t z0, Double_t zLength)
-  :TObject()
-  ,fOwner(kFALSE)
-  ,fPlane(plane)
-  ,fStack(stack)
-  ,fSector(sector)
-  ,fNRows(kMaxRows)
-  ,fN(0)
-  ,fX(0.)
-  ,fZ0(z0)
-  ,fZLength(zLength)
-{
-  //
-  // Default constructor (Only provided to use AliTRDchamberTimeBin with arrays)
-  //
-
-       for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;
-       for(int ic=0; ic<kMaxClustersLayer; ic++){
-               fClusters[ic] = 0x0;
-               fIndex[ic]    = 0xffff;
-       }
-}
-
-// //_____________________________________________________________________________
-// AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer, Double_t
-// z0, Double_t zLength, UChar_t stackNr):
-//     TObject()
-//     ,fOwner(kFALSE)
-//   ,fPlane(0xff)
-//   ,fStack(0xff)
-//   ,fSector(0xff)
-//   ,fNRows(kMaxRows)
-//   ,fN(0)
-//   ,fX(0.)
-//     ,fZ0(z0)
-//     ,fZLength(zLength)
-// {
-// // Standard constructor.
-// // Initialize also the underlying AliTRDpropagationLayer using the copy constructor.
-// 
-//     SetT0(layer.IsT0());
-//     for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;
-//     for(int ic=0; ic<kMaxClustersLayer; ic++){
-//             fClusters[ic] = 0x0;
-//             fIndex[ic]    = 0xffff;
-//     }
-// }
-// 
-// //_____________________________________________________________________________
-// AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer):
-//     TObject()
-//     ,fOwner(kFALSE)
-//   ,fPlane(0xff)
-//   ,fStack(0xff)
-//   ,fSector(0xff)
-//   ,fNRows(kMaxRows)
-//   ,fN(0)
-//   ,fX(0.)
-//     ,fZ0(0)
-//     ,fZLength(0)
-// {
-// // Standard constructor using only AliTRDpropagationLayer.
-//     
-//     SetT0(layer.IsT0());
-//     for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;
-//     for(int ic=0; ic<kMaxClustersLayer; ic++){
-//             fClusters[ic] = 0x0;
-//             fIndex[ic]    = 0xffff;
-//     }
-// }
-// //_____________________________________________________________________________
-// AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDpropagationLayer &layer)
-// {
-// // Assignment operator from an AliTRDpropagationLayer
-// 
-//     if (this != &layer) layer.Copy(*this);
-//     return *this;
-// }
-// 
-
-//_____________________________________________________________________________
-AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDchamberTimeBin &layer):
-       TObject()
-       ,fOwner(layer.fOwner)
-  ,fPlane(layer.fPlane)
-  ,fStack(layer.fStack)
-  ,fSector(layer.fSector)
-       ,fNRows(layer.fNRows)
-  ,fN(layer.fN)
-  ,fX(layer.fX)
-       ,fZ0(layer.fZ0)
-       ,fZLength(layer.fZLength)
-{
-// Copy Constructor (performs a deep copy)
-       
-       SetT0(layer.IsT0());
-       for(int i=0; i<kMaxRows; i++) fPositions[i] = layer.fPositions[i];
-       memcpy(&fClusters[0], &layer.fClusters[0], kMaxClustersLayer*sizeof(UChar_t));
-       memcpy(&fIndex[0], &layer.fIndex[0], kMaxClustersLayer*sizeof(UInt_t));
-
-
-//     BuildIndices();
-}
-
-//_____________________________________________________________________________
-AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDchamberTimeBin &layer)
-{
-// Assignment operator
-
-       if (this != &layer) layer.Copy(*this);
-  return *this;
-}
-
-//_____________________________________________________________________________
-void AliTRDchamberTimeBin::Copy(TObject &o) const
-{
-// Copy method. Performs a deep copy of all data from this object to object o.
-       
-       AliTRDchamberTimeBin &layer = (AliTRDchamberTimeBin &)o;
-       layer.fOwner       = kFALSE;
-       layer.fPlane       = fPlane;
-       layer.fStack       = fStack;
-       layer.fSector      = fSector;
-       layer.fNRows       = fNRows;
-       layer.fN           = fN;
-       layer.fX           = fX;
-       layer.fZ0          = fZ0;
-       layer.fZLength     = fZLength;
-       layer.SetT0(IsT0());
-       
-       for(int i = 0; i < kMaxRows; i++) layer.fPositions[i] = 0;
-       
-       for(int i=0; i<kMaxRows; i++) layer.fPositions[i] = fPositions[i];
-       memcpy(&layer.fClusters[0], &fClusters[0], kMaxClustersLayer*sizeof(UChar_t));
-       memcpy(&layer.fIndex[0], &fIndex[0], kMaxClustersLayer*sizeof(UInt_t));
-       
-       TObject::Copy(layer); // copies everything into layer
-       
-//     layer.BuildIndices();
-}
-
-//_____________________________________________________________________________
-AliTRDchamberTimeBin::~AliTRDchamberTimeBin()
-{
-// Destructor
-       if(fOwner) for(int ic=0; ic<fN; ic++) delete fClusters[ic];
-}
-
-//_____________________________________________________________________________
-void AliTRDchamberTimeBin::SetRange(const Float_t z0, const Float_t zLength)
-{
-// Sets the range in z-direction
-//
-// Parameters:
-//   z0      : starting position of layer in the z direction
-//   zLength : length of layer in the z direction 
-
-       fZ0 = (z0 <= z0 + zLength) ? z0 : z0 + zLength;
-       fZLength = TMath::Abs(zLength);
-}
-
-//_____________________________________________________________________________
-void AliTRDchamberTimeBin::InsertCluster(AliTRDcluster *c, UInt_t index) 
-{
-  //
-  // Insert cluster in cluster array.
-  // Clusters are sorted according to Y coordinate.  
-  //
-
-  //if (fTimeBinIndex < 0) { 
-    //AliWarning("Attempt to insert cluster into non-sensitive time bin!\n");
-    //return;
-  //}
-
-  if (fN == (Int_t) kMaxClustersLayer) {
-    //AliWarning("Too many clusters !\n"); 
-    return;
-  }
-
-  if (fN == 0) {
-    fIndex[0]       = index; 
-    fClusters[fN++] = c; 
-    return;
-  }
-
-  Int_t i = Find(c->GetY());
-  memmove(fClusters+i+1,fClusters+i,(fN-i)*sizeof(AliTRDcluster*));
-  memmove(fIndex   +i+1,fIndex   +i,(fN-i)*sizeof(UInt_t)); 
-  fIndex[i]    = index; 
-  fClusters[i] = c; 
-  fN++;
-
-}  
-
-
-//_____________________________________________________________________________
-void AliTRDchamberTimeBin::BuildIndices(Int_t iter)
-{
-// Rearrangement of the clusters belonging to the propagation layer for the stack.
-//
-// Detailed description
-//
-// The array indices of all clusters in one PropagationLayer are stored in
-// array. The array is divided into several bins.
-// The clusters are sorted in increasing order of their y coordinate.
-//
-// Sorting algorithm: TreeSearch
-//
-
-       if(!fN) return;
-
-       // Select clusters that belong to the Stack
-       Int_t nClStack = 0;                                     // Internal counter
-       for(Int_t i = 0; i < fN; i++){
-               if(fClusters[i]->IsUsed()){
-                       fClusters[i] = 0x0;
-                       fIndex[i] = 0xffff;
-               } else nClStack++;
-       }
-       if(nClStack > kMaxClustersLayer) AliWarning(Form("Number of clusters in stack %d exceed buffer size %d. Truncating.", nClStack, kMaxClustersLayer));
-               
-       // Nothing in this time bin. Reset indexes 
-       if(!nClStack){
-               fN = 0;
-               memset(&fPositions[0], 0xff, sizeof(UChar_t) * kMaxRows);
-               memset(&fClusters[0], 0x0, sizeof(AliTRDcluster*) * kMaxClustersLayer);
-               memset(&fIndex[0], 0xffff, sizeof(UInt_t) * kMaxClustersLayer);
-               return;
-       }
-       
-       // Make a copy
-       AliTRDcluster *helpCL[kMaxClustersLayer];
-       Int_t helpInd[kMaxClustersLayer];
-       nClStack = 0;
-       for(Int_t i = 0; i < TMath::Min(fN, kMaxClustersLayer); i++){
-               if(!fClusters[i]) continue;
-               helpCL[nClStack]  = fClusters[i];
-               helpInd[nClStack] = fIndex[i];
-               fClusters[i]      = 0x0;
-               fIndex[i]         = 0xffff;
-               nClStack++;
-       }
-       
-       // do clusters arrangement
-       fX = 0.;
-       fN =  nClStack;
-       nClStack = 0;
-       // Reset Positions array
-       memset(fPositions, 0, sizeof(UChar_t)*kMaxRows);
-       for(Int_t i = 0; i < fN; i++){
-               // boundary check
-               AliTRDcluster *cl = helpCL[i];
-               UChar_t rowIndex = cl->GetPadRow();
-               // Insert Leaf
-               Int_t pos = FindYPosition(cl->GetY(), rowIndex, i);
-               if(pos == -1){          // zbin is empty;
-                       Int_t upper = (rowIndex == fNRows - 1) ? nClStack : fPositions[rowIndex + 1];
-                       memmove(fClusters + upper + 1, fClusters + upper, (sizeof(AliTRDcluster *))*(nClStack-upper));
-                       memmove(fIndex + upper + 1, fIndex + upper, (sizeof(UInt_t))*(nClStack-upper));
-                       fClusters[upper] = cl;
-                       fIndex[upper] = helpInd[i]; 
-                       // Move All pointer one position back
-                       for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++;
-                       nClStack++;
-               } else {                // zbin not empty
-                       memmove(fClusters + pos + 2, fClusters + pos+1, (sizeof(AliTRDcluster *))*(nClStack-(pos+1)));
-                       memmove(fIndex + pos + 2, fIndex + pos+1, (sizeof(UInt_t))*(nClStack-(pos+1)));
-                       fClusters[pos + 1] = cl;        //fIndex[i];
-                       fIndex[pos + 1] = helpInd[i];
-                       // Move All pointer one position back
-                       for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++; 
-                       nClStack++;
-               }
-
-               // calculate mean x
-               fX += cl->GetX(); 
-               
-               // Debug Streaming
-               if(AliTRDtrackerV1::DebugStreamer() && AliTRDReconstructor::StreamLevel() >= 3){
-                       TTreeSRedirector &cstream = *AliTRDtrackerV1::DebugStreamer();
-                       cstream << "BuildIndices"
-                       << "Plane="    << fPlane
-                       << "Stack="    << fStack
-                       << "Sector="   << fSector
-                       << "Iter="     << iter
-                       << "C.="       << cl
-                       << "rowIndex=" << rowIndex
-                       << "\n";
-               }
-       }
-
-//     AliInfo("Positions");
-//     for(int ir=0; ir<fNRows; ir++) printf("pos[%d] %d\n", ir, fPositions[ir]);
-
-       fX /= fN;
-}
-
-//_____________________________________________________________________________
-Int_t AliTRDchamberTimeBin::Find(Float_t y) const
-{
-  //
-  // Returns index of the cluster nearest in Y    
-  //
-
-  if (fN <= 0) return 0;
-  
-  if (y <= fClusters[0]->GetY()) return 0;
-  
-  if (y >  fClusters[fN-1]->GetY()) return fN;
-  
-
-  Int_t b = 0;
-  Int_t e = fN - 1;
-  Int_t m = (b + e) / 2;
-
-  for ( ; b < e; m = (b + e) / 2) {
-    if (y > fClusters[m]->GetY()) b = m + 1;
-    else e = m;
-  }
-
-  return m;
-}    
-
-//_____________________________________________________________________________
-Int_t AliTRDchamberTimeBin::FindYPosition(Double_t y, UChar_t z, Int_t nClusters) const
-{
-//
-// Tree search Algorithm to find the nearest left cluster for a given
-// y-position in a certain z-bin (in fact AVL-tree). 
-// Making use of the fact that clusters are sorted in y-direction.
-//
-// Parameters:
-//   y : y position of the reference point in tracking coordinates
-//   z : z reference bin.
-//   nClusters : 
-//
-// Output :
-// Index of the nearest left cluster in the StackLayer indexing (-1 if no clusters are found)
-//
-
-       Int_t start = fPositions[z];            // starting Position of the bin
-       Int_t upper = (Int_t)((z != fNRows - 1) ? fPositions[z+1] : nClusters); // ending Position of the bin 
-       Int_t end = upper - 1; // ending Position of the bin 
-       if(end < start) return -1; // Bin is empty
-       Int_t middle = static_cast<Int_t>((start + end)/2);
-       // 1st Part: climb down the tree: get the next cluster BEFORE ypos
-       while(start + 1 < end){
-               if(y >= fClusters[middle]->GetY()) start = middle;
-               else end = middle;
-               middle = static_cast<Int_t>((start + end)/2);
-       }
-       if(y > fClusters[end]->GetY()) return end;
-       return start;
-}
-
-//_____________________________________________________________________________
-Int_t AliTRDchamberTimeBin::FindNearestYCluster(Double_t y, UChar_t z) const
-{
-//
-// Tree search Algorithm to find the nearest cluster for a given
-// y-position in a certain z-bin (in fact AVL-tree). 
-// Making use of the fact that clusters are sorted in y-direction.
-//
-// Parameters:
-//   y : y position of the reference point in tracking coordinates
-//   z : z reference bin.
-//
-// Output 
-// Index of the nearest cluster in the StackLayer indexing (-1 if no clusters are found)
-//
-
-       Int_t position = FindYPosition(y, z, fN);
-       if(position == -1) return position; // bin empty
-       // FindYPosition always returns the left Neighbor. We don't know if the left or the right Neighbor is nearest
-       // to the Reference y-position, so test both
-       Int_t upper = (Int_t)((z < fNRows-1) ? fPositions[z+1] : fN); // ending Position of the bin
-       if((position + 1) < (upper)){
-               if(TMath::Abs(y - fClusters[position + 1]->GetY()) < TMath::Abs(y - fClusters[position]->GetY())) return position + 1;
-               else return position;
-       }
-       return position;
-}
-
-//_____________________________________________________________________________
-Int_t AliTRDchamberTimeBin::SearchNearestCluster(Double_t y, Double_t z, Double_t maxroady, Double_t maxroadz) const
-{
-//
-// Finds the nearest cluster from a given point in a defined range.
-// Distance is determined in a 2D space by the 2-Norm.
-//
-// Parameters :
-//   y : y position of the reference point in tracking coordinates
-//   z : z reference bin.
-//   maxroady : maximum searching distance in y direction
-//   maxroadz : maximum searching distance in z direction
-//
-// Output 
-// Index of the nearest cluster in the StackLayer indexing (-1 if no cluster is found).
-// Cluster can be accessed with the operator[] or GetCluster(Int_t index)
-//
-// Detail description
-//
-// The following steps are perfomed:
-// 1. Get the expected z bins inside maxroadz.
-// 2. For each z bin find nearest y cluster.
-// 3. Select best candidate
-//
-       Int_t   index   = -1;
-       // initial minimal distance will be represented as ellipse: semi-major = z-direction
-       // later 2-Norm will be used  
-//     Float_t nExcentricity = TMath::Sqrt(maxroadz*maxroadz - maxroad*maxroad)/maxroadz;
-       Float_t mindist = maxroadz;
-       
-       // not very nice but unfortunately neccessarry: we have ho check the neighbors in both directions (+ and -) too. How 
-       // much neighbors depends on the Quotient maxroadz/fZLength   
-       UChar_t maxRows = 3;
-       UChar_t zpos[kMaxRows];
-  // Float_t mindist = TMath::Sqrt(maxroad*maxroad + maxroadz*maxroadz);
-//     UChar_t myZbin = FindTreePosition(z, fZ0 + fZLength/2, fZLength/4, 8, 8, kFALSE);
-       UChar_t myZbin = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - z)/fZLength * fNRows);
-       if(z < fZ0) myZbin = fNRows - 1;
-       if(z > fZ0 + fZLength) myZbin = 0;
-       //printf("\n%f < %f < %f [%d]\n", fZ0, z, fZ0 + fZLength, myZbin);
-       //for(int ic=0; ic<fN; ic++) printf("%d z = %f row %d\n", ic, fClusters[ic]->GetZ(), fClusters[ic]->GetPadRow());
-
-       UChar_t nNeighbors = 0;
-       for(UChar_t i = 0; i < maxRows; i++){
-               if((myZbin - 1 + i) < 0) continue;
-               if((myZbin - 1 + i) > fNRows - 1) break;
-               zpos[nNeighbors] = myZbin - 1 + i;
-               nNeighbors++;
-       }
-       Float_t ycl = 0, zcl = 0;
-       for(UChar_t neighbor = 0; neighbor < nNeighbors; neighbor++){   // Always test the neighbors too
-               Int_t pos = FindNearestYCluster(y, zpos[neighbor]);
-               if(pos == -1) continue; // No cluster in bin
-               AliTRDcluster *c = (AliTRDcluster *) (fClusters[pos]);
-               if(c->IsUsed()) continue;               // we are only interested in unused clusters
-               ycl = c->GetY();
-               // Too far away in y-direction (Prearrangement)
-               if (TMath::Abs(ycl - y) > maxroady){ 
-                       //printf("y[%f] ycl[%f] roady[%f]\n", y, ycl, maxroady);
-                       continue;
-               }
-               zcl = c->GetZ();
-               // Too far away in z-Direction
-               // (Prearrangement since we have not so many bins to test)
-               if (TMath::Abs(zcl - z) > maxroadz) continue;
-               
-               Float_t dist;           // distance defined as 2-Norm   
-               // if we havent found a Particle that is in the ellipse around (y,z) with maxroad as semi-minor and
-               // maxroadz as semi-major, we take the radius of the ellipse concerning the cluster as mindist, later we 
-               // take the 2-Norm when we found a cluster inside the ellipse (The value 10000 is taken because it is surely
-               // large enough to be usable as an indicator whether we have found a nearer cluster or not)
-//             if(mindist > 10000.){
-//                     Float_t phi = ((zcl - z) == 0) ? TMath::Pi()/2 : TMath::ATan((ycl - y)/(zcl - z));
-//                     mindist = maxroad/TMath::Sqrt(1 - nExcentricity*nExcentricity * (TMath::Cos(phi))*(TMath::Cos(phi)));
-//             }
-               dist = TMath::Max(TMath::Abs(y-ycl),TMath::Abs(z-zcl)); // infinity Norm
-//             dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z));
-               if((Int_t)(dist * 100000) < (Int_t)(mindist * 100000)){
-               //if((dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z))) < mindist){
-                       mindist = dist;
-                       index   = pos;
-               }       
-       }
-       // This is the Array Position in fIndex2D of the Nearest cluster: if a
-       // cluster is called, then the function has to retrieve the Information
-       // which is Stored in the Array called, the function
-       return index;
-}
-
-//_____________________________________________________________________________
-void AliTRDchamberTimeBin::BuildCond(AliTRDcluster *cl, Double_t *cond, UChar_t Layer, Double_t theta, Double_t phi)
-{
-// Helper function to calculate the area where to expect a cluster in THIS
-// layer. 
-//
-// Parameters :
-//   cl    : 
-//   cond  :
-//   Layer : 
-//   theta : 
-//   phi   : 
-//
-// Detail description
-//
-// Helper function to calculate the area where to expect a cluster in THIS
-// layer. by using the information of a former cluster in another layer
-// and the angle in theta- and phi-direction between layer 0 and layer 3.
-// If the layer is zero, initial conditions are calculated. Otherwise a
-// linear interpolation is performed. 
-//Begin_Html
-//<img src="gif/build_cond.gif">
-//End_Html
-//
-
-       if(!AliTRDReconstructor::RecoParam()){
-               AliError("Reconstruction parameters not initialized.");
-               return;
-       }
-       
-       if(Layer == 0){
-               cond[0] = cl->GetY();                   // center: y-Direction
-               cond[1] = cl->GetZ();                   // center: z-Direction
-               cond[2] = AliTRDReconstructor::RecoParam()->GetMaxPhi()   * (cl->GetX() - GetX()) + 1.0;                        // deviation: y-Direction
-               cond[3] = AliTRDReconstructor::RecoParam()->GetMaxTheta() * (cl->GetX() - GetX()) + 1.0;                        // deviation: z-Direction
-       } else {
-               cond[0] = cl->GetY() + phi   * (GetX() - cl->GetX()); 
-               cond[1] = cl->GetZ() + theta * (GetX() - cl->GetX());
-               cond[2] = AliTRDReconstructor::RecoParam()->GetRoad0y() + phi;
-               cond[3] = AliTRDReconstructor::RecoParam()->GetRoad0z();
-       }
-}
-
-//_____________________________________________________________________________
-void AliTRDchamberTimeBin::GetClusters(Double_t *cond, Int_t *index, Int_t& ncl, Int_t BufferSize)
-{
-// Finds all clusters situated in this layer inside a rectangle  given by the center an ranges.
-//
-// Parameters :
-//   cond  :
-//   index : 
-//   ncl :
-//   BufferSize   :
-//
-// Output :
-//
-// Detail description
-//
-// Function returs an array containing the indices in the stacklayer of
-// the clusters found an  the number of found clusters in the stacklayer
-
-       ncl = 0;
-       memset(index, 0, BufferSize*sizeof(Int_t));
-       if(fN == 0) return;
-               
-       //Boundary checks
-       Double_t zvals[2];
-       if(cond[1] - cond[3] > fZ0 + fZLength || cond[1] + cond[3] < fZ0) return; // We are outside of the chamvber
-       zvals[0] = ((cond[1] - cond[3]) < fZ0) ? fZ0 : (cond[1] - cond[3]);
-       zvals[1] = ((cond[1] + cond[3]) < fZ0 + fZLength) ? (cond[1] + cond[3]) : fZ0 + fZLength - 1.E-3;
-
-       UChar_t zhi = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[0])/fZLength * fNRows);
-       UChar_t zlo = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[1])/fZLength * fNRows);
-       
-       
-//     AliInfo(Form("yc[%f] zc[%f] dy[%f] dz[%f]", cond[0], cond[1], cond[2], cond[3]));
-//     PrintClusters();
-//     AliInfo(Form("zlo[%f] zhi[%f]", zvals[0], zvals[1]));
-//     AliInfo(Form("zlo[%d] zhi[%d]", zlo, zhi));
-       
-       //Preordering in Direction z saves a lot of loops (boundary checked)
-       for(UChar_t z = zlo; z <= zhi; z++){
-               UInt_t upper = (z < fNRows-1) ? fPositions[z+1] : fN;
-               //AliInfo(Form("z[%d] y [%d %d]", z, fPositions[z], upper));
-               for(Int_t y = fPositions[z]; y < (Int_t)upper; y++){
-                       if(ncl == BufferSize){
-                               AliWarning("Buffer size riched. Some clusters may be lost.");
-                               return; //Buffer filled
-                       }
-                       
-                       if(fClusters[y]->GetY() > (cond[0] + cond[2])) break;                   // Abbortion conditions!!!
-                       if(fClusters[y]->GetY() < (cond[0] - cond[2])) continue;        // Too small
-                       if(((Int_t)((fClusters[y]->GetZ())*1000) < (Int_t)(zvals[0]*1000)) || ((Int_t)((fClusters[y]->GetZ())*1000) > (Int_t)(zvals[1]*1000))){/*printf("exit z\n"); TODO*/ continue;}
-                       index[ncl] = y;
-                       ncl++;
-               }
-       }
-       if(ncl>fN) AliError(Form("Clusters found %d > %d (clusters in layer)", ncl, fN));
-}
-
-//_____________________________________________________________________________
-AliTRDcluster *AliTRDchamberTimeBin::GetNearestCluster(Double_t *cond)
-{
-// Function returning a pointer to the nearest cluster (nullpointer if not successfull).
-//
-// Parameters :
-//   cond  :
-//
-// Output :
-//   pointer to the nearest cluster (nullpointer if not successfull).
-// 
-// Detail description
-//
-// returns a pointer to the nearest cluster (nullpointer if not
-// successfull) by the help of the method FindNearestCluster
-       
-       
-       Double_t maxroad  = AliTRDReconstructor::RecoParam()->GetRoad2y();
-       Double_t maxroadz = AliTRDReconstructor::RecoParam()->GetRoad2z();
-       
-       Int_t index = SearchNearestCluster(cond[0],cond[1],maxroad,maxroadz);
-       AliTRDcluster *returnCluster = 0x0;
-       if(index != -1) returnCluster = (AliTRDcluster *) fClusters[index];
-       return returnCluster;
-}
-
-//_____________________________________________________________________________
-void AliTRDchamberTimeBin::PrintClusters() const
-{
-// Prints the position of each cluster in the stacklayer on the stdout
-//
-       printf("\nnRows = %d\n", fNRows);
-       printf("Z0 = %f\n", fZ0);
-       printf("Z1 = %f\n", fZ0+fZLength);
-       printf("clusters in AliTRDchamberTimeBin %d\n", fN);
-       for(Int_t i = 0; i < fN; i++){
-               printf("AliTRDchamberTimeBin: index=%i, Cluster: X = %3.3f [%d] Y = %3.3f [%d] Z = %3.3f [%d]\n", i,  fClusters[i]->GetX(), fClusters[i]->GetLocalTimeBin(), fClusters[i]->GetY(), fClusters[i]->GetPadCol(), fClusters[i]->GetZ(), fClusters[i]->GetPadRow());
-               if(fClusters[i]->IsUsed()) printf("cluster allready used. rejected in search algorithm\n");
-       }
-}
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+/* $Id: AliTRDchamberTimeBin.cxx 23313 2008-01-11 14:56:43Z cblume $ */\r
+\r
+///////////////////////////////////////////////////////////////////////////////\r
+//                                                                           //\r
+//  Organization of clusters at the level of 1 TRD chamber.                  //\r
+//  The data structure is used for tracking at the stack level.              //\r
+//                                                                           //\r
+//  Functionalities:                                                         //\r
+//   1. cluster organization and sorting                                     //\r
+//   2. fast data navigation                                                 //\r
+//                                                                           //\r
+//  Authors:                                                                 //\r
+//    Alex Bercuci <A.Bercuci@gsi.de>                                        //\r
+//    Markus Fasel <M.Fasel@gsi.de>                                          //\r
+//                                                                           //\r
+///////////////////////////////////////////////////////////////////////////////\r
+\r
+#include <TObject.h>\r
+#include <TROOT.h>\r
+#include <TMath.h>\r
+#include <TStopwatch.h>\r
+#include <TTreeStream.h>\r
+\r
+#include "AliLog.h"\r
+\r
+#include "AliTRDcluster.h"\r
+#include "AliTRDchamberTimeBin.h"\r
+#include "AliTRDrecoParam.h"\r
+#include "AliTRDReconstructor.h"\r
+#include "AliTRDtrackerV1.h"\r
+\r
+\r
+ClassImp(AliTRDchamberTimeBin)\r
+\r
+//_____________________________________________________________________________\r
+AliTRDchamberTimeBin::AliTRDchamberTimeBin(Int_t plane, Int_t stack, Int_t sector, Double_t z0, Double_t zLength)\r
+  :TObject()\r
+  ,fOwner(kFALSE)\r
+  ,fPlane(plane)\r
+  ,fStack(stack)\r
+  ,fSector(sector)\r
+  ,fNRows(kMaxRows)\r
+  ,fN(0)\r
+  ,fX(0.)\r
+  ,fZ0(z0)\r
+  ,fZLength(zLength)\r
+{\r
+  //\r
+  // Default constructor (Only provided to use AliTRDchamberTimeBin with arrays)\r
+  //\r
+\r
+       for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r
+       for(int ic=0; ic<kMaxClustersLayer; ic++){\r
+               fClusters[ic] = 0x0;\r
+               fIndex[ic]    = 0xffff;\r
+       }\r
+}\r
+\r
+// //_____________________________________________________________________________\r
+// AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer, Double_t\r
+// z0, Double_t zLength, UChar_t stackNr):\r
+//     TObject()\r
+//     ,fOwner(kFALSE)\r
+//   ,fPlane(0xff)\r
+//   ,fStack(0xff)\r
+//   ,fSector(0xff)\r
+//   ,fNRows(kMaxRows)\r
+//   ,fN(0)\r
+//   ,fX(0.)\r
+//     ,fZ0(z0)\r
+//     ,fZLength(zLength)\r
+// {\r
+// // Standard constructor.\r
+// // Initialize also the underlying AliTRDpropagationLayer using the copy constructor.\r
+// \r
+//     SetT0(layer.IsT0());\r
+//     for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r
+//     for(int ic=0; ic<kMaxClustersLayer; ic++){\r
+//             fClusters[ic] = 0x0;\r
+//             fIndex[ic]    = 0xffff;\r
+//     }\r
+// }\r
+// \r
+// //_____________________________________________________________________________\r
+// AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer):\r
+//     TObject()\r
+//     ,fOwner(kFALSE)\r
+//   ,fPlane(0xff)\r
+//   ,fStack(0xff)\r
+//   ,fSector(0xff)\r
+//   ,fNRows(kMaxRows)\r
+//   ,fN(0)\r
+//   ,fX(0.)\r
+//     ,fZ0(0)\r
+//     ,fZLength(0)\r
+// {\r
+// // Standard constructor using only AliTRDpropagationLayer.\r
+//     \r
+//     SetT0(layer.IsT0());\r
+//     for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r
+//     for(int ic=0; ic<kMaxClustersLayer; ic++){\r
+//             fClusters[ic] = 0x0;\r
+//             fIndex[ic]    = 0xffff;\r
+//     }\r
+// }\r
+// //_____________________________________________________________________________\r
+// AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDpropagationLayer &layer)\r
+// {\r
+// // Assignment operator from an AliTRDpropagationLayer\r
+// \r
+//     if (this != &layer) layer.Copy(*this);\r
+//     return *this;\r
+// }\r
+// \r
+\r
+//_____________________________________________________________________________\r
+AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDchamberTimeBin &layer):\r
+       TObject()\r
+       ,fOwner(layer.fOwner)\r
+  ,fPlane(layer.fPlane)\r
+  ,fStack(layer.fStack)\r
+  ,fSector(layer.fSector)\r
+       ,fNRows(layer.fNRows)\r
+  ,fN(layer.fN)\r
+  ,fX(layer.fX)\r
+       ,fZ0(layer.fZ0)\r
+       ,fZLength(layer.fZLength)\r
+{\r
+// Copy Constructor (performs a deep copy)\r
+       \r
+       SetT0(layer.IsT0());\r
+       for(int i=0; i<kMaxRows; i++) fPositions[i] = layer.fPositions[i];\r
+       memcpy(&fClusters[0], &layer.fClusters[0], kMaxClustersLayer*sizeof(UChar_t));\r
+       memcpy(&fIndex[0], &layer.fIndex[0], kMaxClustersLayer*sizeof(UInt_t));\r
+\r
+\r
+//     BuildIndices();\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDchamberTimeBin &layer)\r
+{\r
+// Assignment operator\r
+\r
+       if (this != &layer) layer.Copy(*this);\r
+  return *this;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliTRDchamberTimeBin::Copy(TObject &o) const\r
+{\r
+// Copy method. Performs a deep copy of all data from this object to object o.\r
+       \r
+       AliTRDchamberTimeBin &layer = (AliTRDchamberTimeBin &)o;\r
+       layer.fOwner       = kFALSE;\r
+       layer.fPlane       = fPlane;\r
+       layer.fStack       = fStack;\r
+       layer.fSector      = fSector;\r
+       layer.fNRows       = fNRows;\r
+       layer.fN           = fN;\r
+       layer.fX           = fX;\r
+       layer.fZ0          = fZ0;\r
+       layer.fZLength     = fZLength;\r
+       layer.SetT0(IsT0());\r
+       \r
+       for(int i = 0; i < kMaxRows; i++) layer.fPositions[i] = 0;\r
+       \r
+       for(int i=0; i<kMaxRows; i++) layer.fPositions[i] = fPositions[i];\r
+       memcpy(&layer.fClusters[0], &fClusters[0], kMaxClustersLayer*sizeof(UChar_t));\r
+       memcpy(&layer.fIndex[0], &fIndex[0], kMaxClustersLayer*sizeof(UInt_t));\r
+       \r
+       TObject::Copy(layer); // copies everything into layer\r
+       \r
+//     layer.BuildIndices();\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AliTRDchamberTimeBin::~AliTRDchamberTimeBin()\r
+{\r
+// Destructor\r
+       if(fOwner) for(int ic=0; ic<fN; ic++) delete fClusters[ic];\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliTRDchamberTimeBin::SetRange(const Float_t z0, const Float_t zLength)\r
+{\r
+// Sets the range in z-direction\r
+//\r
+// Parameters:\r
+//   z0      : starting position of layer in the z direction\r
+//   zLength : length of layer in the z direction \r
+\r
+       fZ0 = (z0 <= z0 + zLength) ? z0 : z0 + zLength;\r
+       fZLength = TMath::Abs(zLength);\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliTRDchamberTimeBin::InsertCluster(AliTRDcluster *c, UInt_t index) \r
+{\r
+  //\r
+  // Insert cluster in cluster array.\r
+  // Clusters are sorted according to Y coordinate.  \r
+  //\r
+\r
+  //if (fTimeBinIndex < 0) { \r
+    //AliWarning("Attempt to insert cluster into non-sensitive time bin!\n");\r
+    //return;\r
+  //}\r
+\r
+  if (fN == (Int_t) kMaxClustersLayer) {\r
+    //AliWarning("Too many clusters !\n"); \r
+    return;\r
+  }\r
+\r
+  if (fN == 0) {\r
+    fIndex[0]       = index; \r
+    fClusters[fN++] = c; \r
+    return;\r
+  }\r
+\r
+  Int_t i = Find(c->GetY());\r
+  memmove(fClusters+i+1,fClusters+i,(fN-i)*sizeof(AliTRDcluster*));\r
+  memmove(fIndex   +i+1,fIndex   +i,(fN-i)*sizeof(UInt_t)); \r
+  fIndex[i]    = index; \r
+  fClusters[i] = c; \r
+  fN++;\r
+\r
+}  \r
+\r
+\r
+//_____________________________________________________________________________\r
+void AliTRDchamberTimeBin::BuildIndices(Int_t iter)\r
+{\r
+// Rearrangement of the clusters belonging to the propagation layer for the stack.\r
+//\r
+// Detailed description\r
+//\r
+// The array indices of all clusters in one PropagationLayer are stored in\r
+// array. The array is divided into several bins.\r
+// The clusters are sorted in increasing order of their y coordinate.\r
+//\r
+// Sorting algorithm: TreeSearch\r
+//\r
+\r
+       if(!fN) return;\r
+\r
+       // Select clusters that belong to the Stack\r
+       Int_t nClStack = 0;                                     // Internal counter\r
+       for(Int_t i = 0; i < fN; i++){\r
+               if(fClusters[i]->IsUsed()){\r
+                       fClusters[i] = 0x0;\r
+                       fIndex[i] = 0xffff;\r
+               } else nClStack++;\r
+       }\r
+       if(nClStack > kMaxClustersLayer) AliWarning(Form("Number of clusters in stack %d exceed buffer size %d. Truncating.", nClStack, kMaxClustersLayer));\r
+               \r
+       // Nothing in this time bin. Reset indexes \r
+       if(!nClStack){\r
+               fN = 0;\r
+               memset(&fPositions[0], 0xff, sizeof(UChar_t) * kMaxRows);\r
+               memset(&fClusters[0], 0x0, sizeof(AliTRDcluster*) * kMaxClustersLayer);\r
+               memset(&fIndex[0], 0xffff, sizeof(UInt_t) * kMaxClustersLayer);\r
+               return;\r
+       }\r
+       \r
+       // Make a copy\r
+       AliTRDcluster *helpCL[kMaxClustersLayer];\r
+       Int_t helpInd[kMaxClustersLayer];\r
+       nClStack = 0;\r
+       for(Int_t i = 0; i < TMath::Min(fN, kMaxClustersLayer); i++){\r
+               if(!fClusters[i]) continue;\r
+               helpCL[nClStack]  = fClusters[i];\r
+               helpInd[nClStack] = fIndex[i];\r
+               fClusters[i]      = 0x0;\r
+               fIndex[i]         = 0xffff;\r
+               nClStack++;\r
+       }\r
+       \r
+       // do clusters arrangement\r
+       fX = 0.;\r
+       fN =  nClStack;\r
+       nClStack = 0;\r
+       // Reset Positions array\r
+       memset(fPositions, 0, sizeof(UChar_t)*kMaxRows);\r
+       for(Int_t i = 0; i < fN; i++){\r
+               // boundary check\r
+               AliTRDcluster *cl = helpCL[i];\r
+               UChar_t rowIndex = cl->GetPadRow();\r
+               // Insert Leaf\r
+               Int_t pos = FindYPosition(cl->GetY(), rowIndex, i);\r
+               if(pos == -1){          // zbin is empty;\r
+                       Int_t upper = (rowIndex == fNRows - 1) ? nClStack : fPositions[rowIndex + 1];\r
+                       memmove(fClusters + upper + 1, fClusters + upper, (sizeof(AliTRDcluster *))*(nClStack-upper));\r
+                       memmove(fIndex + upper + 1, fIndex + upper, (sizeof(UInt_t))*(nClStack-upper));\r
+                       fClusters[upper] = cl;\r
+                       fIndex[upper] = helpInd[i]; \r
+                       // Move All pointer one position back\r
+                       for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++;\r
+                       nClStack++;\r
+               } else {                // zbin not empty\r
+                       memmove(fClusters + pos + 2, fClusters + pos+1, (sizeof(AliTRDcluster *))*(nClStack-(pos+1)));\r
+                       memmove(fIndex + pos + 2, fIndex + pos+1, (sizeof(UInt_t))*(nClStack-(pos+1)));\r
+                       fClusters[pos + 1] = cl;        //fIndex[i];\r
+                       fIndex[pos + 1] = helpInd[i];\r
+                       // Move All pointer one position back\r
+                       for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++; \r
+                       nClStack++;\r
+               }\r
+\r
+               // calculate mean x\r
+               fX += cl->GetX(); \r
+               \r
+               // Debug Streaming\r
+               if(AliTRDtrackerV1::DebugStreamer() && AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 3){\r
+                       TTreeSRedirector &cstream = *AliTRDtrackerV1::DebugStreamer();\r
+                       cstream << "BuildIndices"\r
+                       << "Plane="    << fPlane\r
+                       << "Stack="    << fStack\r
+                       << "Sector="   << fSector\r
+                       << "Iter="     << iter\r
+                       << "C.="       << cl\r
+                       << "rowIndex=" << rowIndex\r
+                       << "\n";\r
+               }\r
+       }\r
+\r
+//     AliInfo("Positions");\r
+//     for(int ir=0; ir<fNRows; ir++) printf("pos[%d] %d\n", ir, fPositions[ir]);\r
+\r
+       fX /= fN;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+Int_t AliTRDchamberTimeBin::Find(Float_t y) const\r
+{\r
+  //\r
+  // Returns index of the cluster nearest in Y    \r
+  //\r
+\r
+  if (fN <= 0) return 0;\r
+  \r
+  if (y <= fClusters[0]->GetY()) return 0;\r
+  \r
+  if (y >  fClusters[fN-1]->GetY()) return fN;\r
+  \r
+\r
+  Int_t b = 0;\r
+  Int_t e = fN - 1;\r
+  Int_t m = (b + e) / 2;\r
+\r
+  for ( ; b < e; m = (b + e) / 2) {\r
+    if (y > fClusters[m]->GetY()) b = m + 1;\r
+    else e = m;\r
+  }\r
+\r
+  return m;\r
+}    \r
+\r
+//_____________________________________________________________________________\r
+Int_t AliTRDchamberTimeBin::FindYPosition(Double_t y, UChar_t z, Int_t nClusters) const\r
+{\r
+//\r
+// Tree search Algorithm to find the nearest left cluster for a given\r
+// y-position in a certain z-bin (in fact AVL-tree). \r
+// Making use of the fact that clusters are sorted in y-direction.\r
+//\r
+// Parameters:\r
+//   y : y position of the reference point in tracking coordinates\r
+//   z : z reference bin.\r
+//   nClusters : \r
+//\r
+// Output :\r
+// Index of the nearest left cluster in the StackLayer indexing (-1 if no clusters are found)\r
+//\r
+\r
+       Int_t start = fPositions[z];            // starting Position of the bin\r
+       Int_t upper = (Int_t)((z != fNRows - 1) ? fPositions[z+1] : nClusters); // ending Position of the bin \r
+       Int_t end = upper - 1; // ending Position of the bin \r
+       if(end < start) return -1; // Bin is empty\r
+       Int_t middle = static_cast<Int_t>((start + end)/2);\r
+       // 1st Part: climb down the tree: get the next cluster BEFORE ypos\r
+       while(start + 1 < end){\r
+               if(y >= fClusters[middle]->GetY()) start = middle;\r
+               else end = middle;\r
+               middle = static_cast<Int_t>((start + end)/2);\r
+       }\r
+       if(y > fClusters[end]->GetY()) return end;\r
+       return start;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+Int_t AliTRDchamberTimeBin::FindNearestYCluster(Double_t y, UChar_t z) const\r
+{\r
+//\r
+// Tree search Algorithm to find the nearest cluster for a given\r
+// y-position in a certain z-bin (in fact AVL-tree). \r
+// Making use of the fact that clusters are sorted in y-direction.\r
+//\r
+// Parameters:\r
+//   y : y position of the reference point in tracking coordinates\r
+//   z : z reference bin.\r
+//\r
+// Output \r
+// Index of the nearest cluster in the StackLayer indexing (-1 if no clusters are found)\r
+//\r
+\r
+       Int_t position = FindYPosition(y, z, fN);\r
+       if(position == -1) return position; // bin empty\r
+       // FindYPosition always returns the left Neighbor. We don't know if the left or the right Neighbor is nearest\r
+       // to the Reference y-position, so test both\r
+       Int_t upper = (Int_t)((z < fNRows-1) ? fPositions[z+1] : fN); // ending Position of the bin\r
+       if((position + 1) < (upper)){\r
+               if(TMath::Abs(y - fClusters[position + 1]->GetY()) < TMath::Abs(y - fClusters[position]->GetY())) return position + 1;\r
+               else return position;\r
+       }\r
+       return position;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+Int_t AliTRDchamberTimeBin::SearchNearestCluster(Double_t y, Double_t z, Double_t maxroady, Double_t maxroadz) const\r
+{\r
+//\r
+// Finds the nearest cluster from a given point in a defined range.\r
+// Distance is determined in a 2D space by the 2-Norm.\r
+//\r
+// Parameters :\r
+//   y : y position of the reference point in tracking coordinates\r
+//   z : z reference bin.\r
+//   maxroady : maximum searching distance in y direction\r
+//   maxroadz : maximum searching distance in z direction\r
+//\r
+// Output \r
+// Index of the nearest cluster in the StackLayer indexing (-1 if no cluster is found).\r
+// Cluster can be accessed with the operator[] or GetCluster(Int_t index)\r
+//\r
+// Detail description\r
+//\r
+// The following steps are perfomed:\r
+// 1. Get the expected z bins inside maxroadz.\r
+// 2. For each z bin find nearest y cluster.\r
+// 3. Select best candidate\r
+//\r
+       Int_t   index   = -1;\r
+       // initial minimal distance will be represented as ellipse: semi-major = z-direction\r
+       // later 2-Norm will be used  \r
+//     Float_t nExcentricity = TMath::Sqrt(maxroadz*maxroadz - maxroad*maxroad)/maxroadz;\r
+       Float_t mindist = maxroadz;\r
+       \r
+       // not very nice but unfortunately neccessarry: we have ho check the neighbors in both directions (+ and -) too. How \r
+       // much neighbors depends on the Quotient maxroadz/fZLength   \r
+       UChar_t maxRows = 3;\r
+       UChar_t zpos[kMaxRows];\r
+  // Float_t mindist = TMath::Sqrt(maxroad*maxroad + maxroadz*maxroadz);\r
+//     UChar_t myZbin = FindTreePosition(z, fZ0 + fZLength/2, fZLength/4, 8, 8, kFALSE);\r
+       UChar_t myZbin = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - z)/fZLength * fNRows);\r
+       if(z < fZ0) myZbin = fNRows - 1;\r
+       if(z > fZ0 + fZLength) myZbin = 0;\r
+       //printf("\n%f < %f < %f [%d]\n", fZ0, z, fZ0 + fZLength, myZbin);\r
+       //for(int ic=0; ic<fN; ic++) printf("%d z = %f row %d\n", ic, fClusters[ic]->GetZ(), fClusters[ic]->GetPadRow());\r
+\r
+       UChar_t nNeighbors = 0;\r
+       for(UChar_t i = 0; i < maxRows; i++){\r
+               if((myZbin - 1 + i) < 0) continue;\r
+               if((myZbin - 1 + i) > fNRows - 1) break;\r
+               zpos[nNeighbors] = myZbin - 1 + i;\r
+               nNeighbors++;\r
+       }\r
+       Float_t ycl = 0, zcl = 0;\r
+       for(UChar_t neighbor = 0; neighbor < nNeighbors; neighbor++){   // Always test the neighbors too\r
+               Int_t pos = FindNearestYCluster(y, zpos[neighbor]);\r
+               if(pos == -1) continue; // No cluster in bin\r
+               AliTRDcluster *c = (AliTRDcluster *) (fClusters[pos]);\r
+               if(c->IsUsed()) continue;               // we are only interested in unused clusters\r
+               ycl = c->GetY();\r
+               // Too far away in y-direction (Prearrangement)\r
+               if (TMath::Abs(ycl - y) > maxroady){ \r
+                       //printf("y[%f] ycl[%f] roady[%f]\n", y, ycl, maxroady);\r
+                       continue;\r
+               }\r
+               zcl = c->GetZ();\r
+               // Too far away in z-Direction\r
+               // (Prearrangement since we have not so many bins to test)\r
+               if (TMath::Abs(zcl - z) > maxroadz) continue;\r
+               \r
+               Float_t dist;           // distance defined as 2-Norm   \r
+               // if we havent found a Particle that is in the ellipse around (y,z) with maxroad as semi-minor and\r
+               // maxroadz as semi-major, we take the radius of the ellipse concerning the cluster as mindist, later we \r
+               // take the 2-Norm when we found a cluster inside the ellipse (The value 10000 is taken because it is surely\r
+               // large enough to be usable as an indicator whether we have found a nearer cluster or not)\r
+//             if(mindist > 10000.){\r
+//                     Float_t phi = ((zcl - z) == 0) ? TMath::Pi()/2 : TMath::ATan((ycl - y)/(zcl - z));\r
+//                     mindist = maxroad/TMath::Sqrt(1 - nExcentricity*nExcentricity * (TMath::Cos(phi))*(TMath::Cos(phi)));\r
+//             }\r
+               dist = TMath::Max(TMath::Abs(y-ycl),TMath::Abs(z-zcl)); // infinity Norm\r
+//             dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z));\r
+               if((Int_t)(dist * 100000) < (Int_t)(mindist * 100000)){\r
+               //if((dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z))) < mindist){\r
+                       mindist = dist;\r
+                       index   = pos;\r
+               }       \r
+       }\r
+       // This is the Array Position in fIndex2D of the Nearest cluster: if a\r
+       // cluster is called, then the function has to retrieve the Information\r
+       // which is Stored in the Array called, the function\r
+       return index;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliTRDchamberTimeBin::BuildCond(AliTRDcluster *cl, Double_t *cond, UChar_t Layer, Double_t theta, Double_t phi)\r
+{\r
+// Helper function to calculate the area where to expect a cluster in THIS\r
+// layer. \r
+//\r
+// Parameters :\r
+//   cl    : \r
+//   cond  :\r
+//   Layer : \r
+//   theta : \r
+//   phi   : \r
+//\r
+// Detail description\r
+//\r
+// Helper function to calculate the area where to expect a cluster in THIS\r
+// layer. by using the information of a former cluster in another layer\r
+// and the angle in theta- and phi-direction between layer 0 and layer 3.\r
+// If the layer is zero, initial conditions are calculated. Otherwise a\r
+// linear interpolation is performed. \r
+//Begin_Html\r
+//<img src="gif/build_cond.gif">\r
+//End_Html\r
+//\r
+\r
+       if(!AliTRDReconstructor::RecoParam()){\r
+               AliError("Reconstruction parameters not initialized.");\r
+               return;\r
+       }\r
+       \r
+       if(Layer == 0){\r
+               cond[0] = cl->GetY();                   // center: y-Direction\r
+               cond[1] = cl->GetZ();                   // center: z-Direction\r
+               cond[2] = AliTRDReconstructor::RecoParam()->GetMaxPhi()   * (cl->GetX() - GetX()) + 1.0;                        // deviation: y-Direction\r
+               cond[3] = AliTRDReconstructor::RecoParam()->GetMaxTheta() * (cl->GetX() - GetX()) + 1.0;                        // deviation: z-Direction\r
+       } else {\r
+               cond[0] = cl->GetY() + phi   * (GetX() - cl->GetX()); \r
+               cond[1] = cl->GetZ() + theta * (GetX() - cl->GetX());\r
+               cond[2] = AliTRDReconstructor::RecoParam()->GetRoad0y() + phi;\r
+               cond[3] = AliTRDReconstructor::RecoParam()->GetRoad0z();\r
+       }\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliTRDchamberTimeBin::GetClusters(Double_t *cond, Int_t *index, Int_t& ncl, Int_t BufferSize)\r
+{\r
+// Finds all clusters situated in this layer inside a rectangle  given by the center an ranges.\r
+//\r
+// Parameters :\r
+//   cond  :\r
+//   index : \r
+//   ncl :\r
+//   BufferSize   :\r
+//\r
+// Output :\r
+//\r
+// Detail description\r
+//\r
+// Function returs an array containing the indices in the stacklayer of\r
+// the clusters found an  the number of found clusters in the stacklayer\r
+\r
+       ncl = 0;\r
+       memset(index, 0, BufferSize*sizeof(Int_t));\r
+       if(fN == 0) return;\r
+               \r
+       //Boundary checks\r
+       Double_t zvals[2];\r
+       if(cond[1] - cond[3] > fZ0 + fZLength || cond[1] + cond[3] < fZ0) return; // We are outside of the chamvber\r
+       zvals[0] = ((cond[1] - cond[3]) < fZ0) ? fZ0 : (cond[1] - cond[3]);\r
+       zvals[1] = ((cond[1] + cond[3]) < fZ0 + fZLength) ? (cond[1] + cond[3]) : fZ0 + fZLength - 1.E-3;\r
+\r
+       UChar_t zhi = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[0])/fZLength * fNRows);\r
+       UChar_t zlo = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[1])/fZLength * fNRows);\r
+       \r
+       \r
+//     AliInfo(Form("yc[%f] zc[%f] dy[%f] dz[%f]", cond[0], cond[1], cond[2], cond[3]));\r
+//     PrintClusters();\r
+//     AliInfo(Form("zlo[%f] zhi[%f]", zvals[0], zvals[1]));\r
+//     AliInfo(Form("zlo[%d] zhi[%d]", zlo, zhi));\r
+       \r
+       //Preordering in Direction z saves a lot of loops (boundary checked)\r
+       for(UChar_t z = zlo; z <= zhi; z++){\r
+               UInt_t upper = (z < fNRows-1) ? fPositions[z+1] : fN;\r
+               //AliInfo(Form("z[%d] y [%d %d]", z, fPositions[z], upper));\r
+               for(Int_t y = fPositions[z]; y < (Int_t)upper; y++){\r
+                       if(ncl == BufferSize){\r
+                               AliWarning("Buffer size riched. Some clusters may be lost.");\r
+                               return; //Buffer filled\r
+                       }\r
+                       \r
+                       if(fClusters[y]->GetY() > (cond[0] + cond[2])) break;                   // Abbortion conditions!!!\r
+                       if(fClusters[y]->GetY() < (cond[0] - cond[2])) continue;        // Too small\r
+                       if(((Int_t)((fClusters[y]->GetZ())*1000) < (Int_t)(zvals[0]*1000)) || ((Int_t)((fClusters[y]->GetZ())*1000) > (Int_t)(zvals[1]*1000))){/*printf("exit z\n"); TODO*/ continue;}\r
+                       index[ncl] = y;\r
+                       ncl++;\r
+               }\r
+       }\r
+       if(ncl>fN) AliError(Form("Clusters found %d > %d (clusters in layer)", ncl, fN));\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AliTRDcluster *AliTRDchamberTimeBin::GetNearestCluster(Double_t *cond)\r
+{\r
+// Function returning a pointer to the nearest cluster (nullpointer if not successfull).\r
+//\r
+// Parameters :\r
+//   cond  :\r
+//\r
+// Output :\r
+//   pointer to the nearest cluster (nullpointer if not successfull).\r
+// \r
+// Detail description\r
+//\r
+// returns a pointer to the nearest cluster (nullpointer if not\r
+// successfull) by the help of the method FindNearestCluster\r
+       \r
+       \r
+       Double_t maxroad  = AliTRDReconstructor::RecoParam()->GetRoad2y();\r
+       Double_t maxroadz = AliTRDReconstructor::RecoParam()->GetRoad2z();\r
+       \r
+       Int_t index = SearchNearestCluster(cond[0],cond[1],maxroad,maxroadz);\r
+       AliTRDcluster *returnCluster = 0x0;\r
+       if(index != -1) returnCluster = (AliTRDcluster *) fClusters[index];\r
+       return returnCluster;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliTRDchamberTimeBin::PrintClusters() const\r
+{\r
+// Prints the position of each cluster in the stacklayer on the stdout\r
+//\r
+       printf("\nnRows = %d\n", fNRows);\r
+       printf("Z0 = %f\n", fZ0);\r
+       printf("Z1 = %f\n", fZ0+fZLength);\r
+       printf("clusters in AliTRDchamberTimeBin %d\n", fN);\r
+       for(Int_t i = 0; i < fN; i++){\r
+               printf("AliTRDchamberTimeBin: index=%i, Cluster: X = %3.3f [%d] Y = %3.3f [%d] Z = %3.3f [%d]\n", i,  fClusters[i]->GetX(), fClusters[i]->GetLocalTimeBin(), fClusters[i]->GetY(), fClusters[i]->GetPadCol(), fClusters[i]->GetZ(), fClusters[i]->GetPadRow());\r
+               if(fClusters[i]->IsUsed()) printf("cluster allready used. rejected in search algorithm\n");\r
+       }\r
+}\r
index 87a9410..d24992e 100644 (file)
@@ -48,6 +48,9 @@ AliTRDrecoParam::AliTRDrecoParam()
   ,fkChi2Z(30./*14.*//*12.5*/)
   ,fkChi2Y(.25)
   ,fkTrackLikelihood(-15.)
+  ,fkStreamLevel(0)
+  ,fSeedingOn(kFALSE)
+  ,fVertexConstrained(kTRUE)
   ,fClusMaxThresh(4.5)
   ,fClusSigThresh(3.5)
   ,fLUTOn(kTRUE)
@@ -93,6 +96,8 @@ AliTRDrecoParam *AliTRDrecoParam::GetCosmicTestParam()
   AliTRDrawStreamBase::SetRawStreamVersion("TB");
   AliTRDrecoParam *par = new AliTRDrecoParam();
   par->SetADCbaseline(10);
+  par->SetSeedingOn(kTRUE);
+  par->SetVertexConstrained(kTRUE);
   return par;
 
 }
index fbc5c06..69d1575 100644 (file)
@@ -46,7 +46,11 @@ class AliTRDrecoParam : public AliDetectorRecoParam
        Double_t GetPlaneQualityThreshold() const { return fkPlaneQualityThreshold; }
 
        Double_t GetTrackLikelihood() const       { return fkTrackLikelihood;       }
+       Int_t    GetStreamLevel() const           { return fkStreamLevel;           }
 
+       Bool_t SeedingOn() const { return fSeedingOn; }
+       Bool_t IsVertexConstrained() const { return fVertexConstrained; }
+       
         Double_t GetClusMaxThresh() const         { return fClusMaxThresh;   };
         Double_t GetClusSigThresh() const         { return fClusSigThresh;   };
         Int_t    GetTCnexp() const                { return fTCnexp;          };
@@ -54,13 +58,16 @@ class AliTRDrecoParam : public AliDetectorRecoParam
         Bool_t   TCOn() const                     { return fTCOn;            };
 
         Int_t    GetADCbaseline() const           { return fADCbaseline;     };
-
+        
        static   AliTRDrecoParam *GetLowFluxParam();
         static   AliTRDrecoParam *GetHighFluxParam();
         static   AliTRDrecoParam *GetCosmicTestParam();
 
        void     SetClusterSharing(Bool_t share = kTRUE) { fkClusterSharing = share;  };
        void     SetPIDMethod(Int_t pid = 1)             { fkPIDMethod = pid ? 1 : 0; };
+       void     SetSeedingOn(Bool_t seedingOn = kTRUE)  { fSeedingOn = seedingOn; }
+       void     SetVertexConstrained(Bool_t vertexConstrained = kTRUE) { fVertexConstrained = vertexConstrained; }
+       void     SetStreamLevel(Int_t streamLevel= 1) { fkStreamLevel = streamLevel; }
         void     SetLUT(Int_t lutOn = 1)                 { fLUTOn           = lutOn;  };
         void     SetClusMaxThresh(Float_t thresh)        { fClusMaxThresh   = thresh; };
         void     SetClusSigThresh(Float_t thresh)        { fClusSigThresh   = thresh; };
@@ -90,6 +97,10 @@ class AliTRDrecoParam : public AliDetectorRecoParam
        Double_t  fkChi2Z;                 // Max chi2 on the z direction for seeding clusters fit
        Double_t  fkChi2Y;                 // Max chi2 on the y direction for seeding clusters Rieman fit
        Double_t  fkTrackLikelihood;       // Track likelihood for tracklets Rieman fit
+       Int_t     fkStreamLevel;                                         // Streaming Level in TRD Reconstruction
+       
+       Bool_t    fSeedingOn;                // Do stand alone tracking in the TRD
+       Bool_t    fVertexConstrained;      // Perform vertex constrained fit
 
         // Clusterization parameter
         Double_t  fClusMaxThresh;          // Threshold value for cluster maximum
index 04af0f3..3f45e8a 100644 (file)
@@ -355,7 +355,7 @@ Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
        }
 
        AliTRDchamberTimeBin *layer = 0x0;
-       if(AliTRDReconstructor::StreamLevel()>=7 && c){
+       if(AliTRDReconstructor::RecoParam()->GetStreamLevel()>=7 && c){
                TClonesArray clusters("AliTRDcluster", 24);
                clusters.SetOwner(kTRUE);
                AliTRDcluster *cc = 0x0;
@@ -423,7 +423,7 @@ Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
                        fZ[iTime]        = cl->GetZ();
                        ncl++;
                }
-       if(AliTRDReconstructor::StreamLevel()>=7) AliInfo(Form("iter = %d ncl [%d] = %d", iter, fPlane, ncl));
+       if(AliTRDReconstructor::RecoParam()->GetStreamLevel()>=7) AliInfo(Form("iter = %d ncl [%d] = %d", iter, fPlane, ncl));
                
                if(ncl>1){      
                        // calculate length of the time bin (calibration aware)
@@ -464,7 +464,7 @@ Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
                        
                        AliTRDseed::Update();
                }
-       if(AliTRDReconstructor::StreamLevel()>=7) AliInfo(Form("iter = %d nclFit [%d] = %d", iter, fPlane, fN2));
+       if(AliTRDReconstructor::RecoParam()->GetStreamLevel()>=7) AliInfo(Form("iter = %d nclFit [%d] = %d", iter, fPlane, fN2));
                
                if(IsOK()){
                        tquality = GetQuality(kZcorr);
index b092679..e0d182b 100644 (file)
@@ -50,6 +50,7 @@
 #include "AliTRDCommonParam.h"
 #include "AliTRDtracker.h"
 #include "AliTRDReconstructor.h"
+#include "AliTRDrecoParam.h"
 #include "AliTRDCalibraFillHisto.h"
 
 ClassImp(AliTRDtracker)
@@ -417,7 +418,7 @@ Int_t AliTRDtracker::PropagateBack(AliESDEvent *event)
        Int_t    nSeed   = event->GetNumberOfTracks();
        if(!nSeed){
                // run stand alone tracking
-               if (AliTRDReconstructor::SeedingOn()) Clusters2Tracks(event);
+               if (AliTRDReconstructor::RecoParam()->SeedingOn()) Clusters2Tracks(event);
                return 0;
        }
        
@@ -525,7 +526,7 @@ Int_t AliTRDtracker::PropagateBack(AliESDEvent *event)
                // Debug part of tracking
                TTreeSRedirector &cstream = *fDebugStreamer;
                Int_t eventNrInFile = event->GetEventNumberInFile(); // This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number.
-               if (AliTRDReconstructor::StreamLevel() > 0) {
+               if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0) {
                        if (track->GetBackupTrack()) {
                                cstream << "Tracks"
                                << "EventNrInFile="  << eventNrInFile
@@ -773,7 +774,7 @@ Int_t AliTRDtracker::RefitInward(AliESDEvent *event)
 
       // Add TRD track to ESDfriendTrack - maybe this tracks are
       // not useful for post-processing - TODO make decision
-      if (AliTRDReconstructor::StreamLevel() > 0)  {
+      if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0)  {
         seed->AddCalibObject(new AliTRDtrack(*pt2/*, kTRUE*/));
       }
       delete pt2;
@@ -781,7 +782,7 @@ Int_t AliTRDtracker::RefitInward(AliESDEvent *event)
     }
 
     // Add TRD track to ESDfriendTrack
-    if (AliTRDReconstructor::StreamLevel() > 0)  {
+    if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0)  {
       seed->AddCalibObject(new AliTRDtrack(*pt/*, kTRUE*/));
     }
     delete pt;
@@ -1559,7 +1560,7 @@ Int_t AliTRDtracker::Clusters2Tracks(AliESDEvent *esd)
               isFake = kTRUE;
            }
 
-           if (AliTRDReconstructor::StreamLevel() > 0) {
+           if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0) {
              if ((!isFake) || ((icl3%10) == 0)) {  // Debugging print
                TTreeSRedirector &cstream = *fDebugStreamer;
                cstream << "Seeds0"
@@ -2197,7 +2198,7 @@ Int_t AliTRDtracker::Clusters2Tracks(AliESDEvent *esd)
            if (1 || (!isFake)) {
              Float_t zvertex = GetZ();
              TTreeSRedirector &cstream = *fDebugStreamer;
-             if (AliTRDReconstructor::StreamLevel() > 0) {
+             if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0) {
                cstream << "Seeds1"
                        << "isFake="     << isFake
                        << "Vertex="     << zvertex
@@ -2418,7 +2419,7 @@ Int_t AliTRDtracker::Clusters2Tracks(AliESDEvent *esd)
          esdtrack.UpdateTrackParams(track,AliESDtrack::kTRDout);
          esdtrack.SetLabel(label);
          esd->AddTrack(&esdtrack);     
-         if (AliTRDReconstructor::StreamLevel() > 0) {
+         if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0) {
            cstream << "Tracks"
                    << "EventNrInFile="  << eventNrInFile
                    << "ESD.="     << &esdtrack
@@ -2428,7 +2429,7 @@ Int_t AliTRDtracker::Clusters2Tracks(AliESDEvent *esd)
          }
        }
 
-       if (AliTRDReconstructor::StreamLevel() > 0) {
+       if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0) {
          cstream << "Seeds2"
                  << "Iter="      << jter
                  << "Track.="    << track
@@ -3661,7 +3662,7 @@ Int_t AliTRDtracker::FindClusters(Int_t sector, Int_t t0, Int_t t1
   TGraph graphy(t1-t0,x,yt);
   TGraph graphz(t1-t0,x,zt);
 
-  if (AliTRDReconstructor::StreamLevel() > 0) {
+  if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0) {
     cstream << "tracklet"
            << "track.="      << track              // Track parameters
            << "tany="        << tany               // Tangent of the local track angle 
index 2768098..632aaec 100644 (file)
@@ -94,11 +94,17 @@ AliTRDtrackerV1::AliTRDtrackerV1()
   if (!AliTRDcalibDB::Instance()) {
     AliFatal("Could not get calibration object");
   }
-  fgNTimeBins = AliTRDcalibDB::Instance()->GetNumberOfTimeBins();
+  if(!fgNTimeBins) fgNTimeBins = AliTRDcalibDB::Instance()->GetNumberOfTimeBins();
 
   for (Int_t isector = 0; isector < AliTRDgeometry::kNsect; isector++) new(&fTrSec[isector]) AliTRDtrackingSector(fGeom, isector);
+
+  if (!AliTRDReconstructor::RecoParam()){
+       AliWarning("RecoParams not set in AliTRDReconstructor. Setting to default LowFluxParam.");
+    AliTRDReconstructor::SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
+  }
   
-  if(AliTRDReconstructor::StreamLevel() > 1){
+
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1){
     TDirectory *savedir = gDirectory; 
     fgDebugStreamer    = new TTreeSRedirector("TRD.TrackerDebug.root");
     savedir->cd();
@@ -252,7 +258,7 @@ Int_t AliTRDtrackerV1::PropagateBack(AliESDEvent *event)
   Int_t    nSeed   = event->GetNumberOfTracks();
   if(!nSeed){
     // run stand alone tracking
-    if (AliTRDReconstructor::SeedingOn()) Clusters2Tracks(event);
+    if (AliTRDReconstructor::RecoParam()->SeedingOn()) Clusters2Tracks(event);
     return 0;
   }
   
@@ -299,7 +305,7 @@ Int_t AliTRDtrackerV1::PropagateBack(AliESDEvent *event)
         track.UpdateESDtrack(seed);
         
         // Add TRD track to ESDfriendTrack
-        if (AliTRDReconstructor::StreamLevel() > 0 /*&& quality TODO*/){ 
+        if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0 /*&& quality TODO*/){ 
           AliTRDtrackV1 *calibTrack = new AliTRDtrackV1(track);
           calibTrack->SetOwner();
           seed->AddCalibObject(calibTrack);
@@ -388,7 +394,7 @@ Int_t AliTRDtrackerV1::PropagateBack(AliESDEvent *event)
         track.UpdateESDtrack(seed);
         
         // Add TRD track to ESDfriendTrack
-//                             if (AliTRDReconstructor::StreamLevel() > 0 /*&& quality TODO*/){ 
+//                             if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0 /*&& quality TODO*/){ 
 //                                     AliTRDtrackV1 *calibTrack = new AliTRDtrackV1(track);
 //                                     calibTrack->SetOwner();
 //                                     seed->AddCalibObject(calibTrack);
@@ -402,7 +408,7 @@ Int_t AliTRDtrackerV1::PropagateBack(AliESDEvent *event)
         track.UpdateESDtrack(seed);
         
         // Add TRD track to ESDfriendTrack
-//                             if (AliTRDReconstructor::StreamLevel() > 0 /*&& quality TODO*/){ 
+//                             if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0 /*&& quality TODO*/){ 
 //                                     AliTRDtrackV1 *calibTrack = new AliTRDtrackV1(track);
 //                                     calibTrack->SetOwner();
 //                                     seed->AddCalibObject(calibTrack);
@@ -556,7 +562,7 @@ Int_t AliTRDtrackerV1::FollowProlongation(AliTRDtrackV1 &t)
     }
   }
 
-  if(AliTRDReconstructor::StreamLevel() > 1){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1){
     Int_t index;
     for(int iplane=0; iplane<6; iplane++){
       AliTRDseedV1 *tracklet = GetTracklet(&t, iplane, index);
@@ -725,7 +731,7 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
     
   } // end planes loop
 
-  if(AliTRDReconstructor::StreamLevel() > 1){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1){
     TTreeSRedirector &cstreamer = *fgDebugStreamer;
     Int_t eventNumber = AliTRDtrackerDebug::GetEventNumber();
     //AliTRDtrackV1 *debugTrack = new AliTRDtrackV1(t);
@@ -879,7 +885,7 @@ Float_t AliTRDtrackerV1::FitTiltedRiemanConstraint(AliTRDseedV1 *tracklets, Doub
   for(Int_t ip = 0; ip < AliTRDtrackerV1::kNPlanes; ip++)
     tracklets[ip].SetCC(curvature);
 
-  if(AliTRDReconstructor::StreamLevel() >= 5){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 5){
     //Linear Model on z-direction
     Double_t xref = CalculateReferenceX(tracklets);            // Relative to the middle of the stack
     Double_t slope = fitter->GetParameter(2);
@@ -1051,7 +1057,7 @@ Float_t AliTRDtrackerV1::FitTiltedRieman(AliTRDseedV1 *tracklets, Bool_t sigErro
     tracklets[iLayer].SetChi2(chi2track);
   }
   
-  if(AliTRDReconstructor::StreamLevel() >=5){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >=5){
     TTreeSRedirector &cstreamer = *fgDebugStreamer;
     Int_t eventNumber                  = AliTRDtrackerDebug::GetEventNumber();
     Int_t candidateNumber      = AliTRDtrackerDebug::GetCandidateNumber();
@@ -1283,7 +1289,7 @@ Double_t AliTRDtrackerV1::FitRiemanTilt(AliTRDtrackV1 *track, AliTRDseedV1 *trac
     }
   }
   
-  if(AliTRDReconstructor::StreamLevel() >=5){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >=5){
     TTreeSRedirector &cstreamer = *fgDebugStreamer;
     Int_t eventNumber                  = AliTRDtrackerDebug::GetEventNumber();
     Int_t candidateNumber      = AliTRDtrackerDebug::GetCandidateNumber();
@@ -1568,6 +1574,30 @@ AliTRDseedV1* AliTRDtrackerV1::SetTracklet(AliTRDseedV1 *tracklet)
   return new ((*fTracklets)[nentries]) AliTRDseedV1(*tracklet);
 }
 
+//____________________________________________________________________
+AliTRDtrackV1* AliTRDtrackerV1::SetTrack(AliTRDtrackV1 *track)
+{
+  // Add this track to the list of tracks stored in the tracker
+  //
+  // Parameters
+  //   - track : pointer to the track to be added to the list
+  //
+  // Output
+  //   - the pointer added
+  //
+  // Detailed description
+  // Build the tracks list if it is not yet created (late initialization)
+  // and adds the new track to the list.
+  //
+  if(!fTracks){
+    fTracks = new TClonesArray("AliTRDtrackV1", AliTRDgeometry::Nsect()*kMaxTracksStack);
+    fTracks->SetOwner(kTRUE);
+  }
+  Int_t nentries = fTracks->GetEntriesFast();
+  return new ((*fTracks)[nentries]) AliTRDtrackV1(*track);
+}
+
+
 
 //____________________________________________________________________
 Int_t AliTRDtrackerV1::Clusters2TracksSM(Int_t sector, AliESDEvent *esd)
@@ -1661,7 +1691,7 @@ Int_t AliTRDtrackerV1::Clusters2TracksStack(AliTRDtrackingChamber **stack, TClon
   
   // Build initial seeding configurations
   Double_t quality = BuildSeedingConfigs(stack, configs);
-  if(AliTRDReconstructor::StreamLevel() > 1){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1){
     AliInfo(Form("Plane config %d %d %d Quality %f"
     , configs[0], configs[1], configs[2], quality));
   }
@@ -1680,7 +1710,7 @@ Int_t AliTRDtrackerV1::Clusters2TracksStack(AliTRDtrackingChamber **stack, TClon
       ntracks = MakeSeeds(stack, &sseed[6*ntracks], pars);
       if(ntracks == kMaxTracksStack) break;
     }
-    if(AliTRDReconstructor::StreamLevel() > 1) AliInfo(Form("Candidate TRD tracks %d in iteration %d.", ntracks, fSieveSeeding));
+    if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1) AliInfo(Form("Candidate TRD tracks %d in iteration %d.", ntracks, fSieveSeeding));
     
     if(!ntracks) break;
     
@@ -1844,7 +1874,7 @@ Int_t AliTRDtrackerV1::Clusters2TracksStack(AliTRDtrackingChamber **stack, TClon
   Int_t ich = 0; while(!(chamber = stack[ich])) ich++;
   trackParams[6] = fGeom->GetSector(chamber->GetDetector());/* *alpha+shift;   // Supermodule*/
 
-  if(AliTRDReconstructor::StreamLevel() > 1){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1){
     AliInfo(Form("Track %d [%d] nlayers %d trackQuality = %e nused %d, yref = %3.3f", itrack, trackIndex, nlayers, fTrackQuality[trackIndex], nused, trackParams[1]));
           
     Int_t nclusters = 0;
@@ -1903,7 +1933,7 @@ Int_t AliTRDtrackerV1::Clusters2TracksStack(AliTRDtrackingChamber **stack, TClon
   esdTrack.SetLabel(track->GetLabel());
   track->UpdateESDtrack(&esdTrack);
   // write ESD-friends if neccessary
-  if (AliTRDReconstructor::StreamLevel() > 0){
+  if (AliTRDReconstructor::RecoParam()->GetStreamLevel() > 0){
     //printf("Creating Calibrations Object\n");
     AliTRDtrackV1 *calibTrack = new AliTRDtrackV1(*track);
     calibTrack->SetOwner();
@@ -1931,7 +1961,7 @@ Int_t AliTRDtrackerV1::Clusters2TracksStack(AliTRDtrackingChamber **stack, TClon
       chamber->Build(fGeom);//Indices(fSieveSeeding);
     }
 
-    if(AliTRDReconstructor::StreamLevel() > 1){ 
+    if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1){ 
       AliInfo(Form("Sieve level %d Plane config %d %d %d Quality %f", fSieveSeeding, configs[0], configs[1], configs[2], quality));
     }
   } while(fSieveSeeding<10); // end stack clusters sieve
@@ -2076,7 +2106,7 @@ Int_t AliTRDtrackerV1::MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *ss
     padlength[iplane] = pp->GetLengthIPad();
   }
   
-  if(AliTRDReconstructor::StreamLevel() > 1){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() > 1){
     AliInfo(Form("Making seeds Stack[%d] Config[%d] Tracks[%d]...", istack, config, ntracks));
   }
 
@@ -2138,7 +2168,7 @@ Int_t AliTRDtrackerV1::MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *ss
   }
 
   Bool_t isFake = kFALSE;
-  if(AliTRDReconstructor::StreamLevel() >= 2){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 2){
     if (c[0]->GetLabel(0) != c[3]->GetLabel(0)) isFake = kTRUE;
     if (c[1]->GetLabel(0) != c[3]->GetLabel(0)) isFake = kTRUE;
     if (c[2]->GetLabel(0) != c[3]->GetLabel(0)) isFake = kTRUE;
@@ -2256,7 +2286,7 @@ Int_t AliTRDtrackerV1::MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *ss
 
   // AliInfo("Extrapolation done.");
   // Debug Stream containing all the 6 tracklets
-  if(AliTRDReconstructor::StreamLevel() >= 2){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 2){
     TTreeSRedirector &cstreamer = *fgDebugStreamer;
     TLinearFitter *tiltedRieman = GetTiltedRiemanFitter();
     Int_t eventNumber          = AliTRDtrackerDebug::GetEventNumber();
@@ -2288,7 +2318,10 @@ 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);
-  chi2Vals[1] = FitTiltedRiemanConstraint(&cseed[0], GetZ());
+  if(AliTRDReconstructor::RecoParam()->IsVertexConstrained())
+       chi2Vals[1] = FitTiltedRiemanConstraint(&cseed[0], GetZ()); // Do Vertex Constrained fit if desired
+  else
+       chi2Vals[1] = 1.;
   chi2Vals[2] = GetChi2Z(&cseed[0]) / TMath::Max((mlayers - 3.), 1.);
   // Chi2 definitions in testing stage
   //chi2Vals[2] = GetChi2ZTest(&cseed[0]);
@@ -2320,7 +2353,7 @@ Int_t AliTRDtrackerV1::MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *ss
     cseed[iLayer].SetChi2Z(chi2[1]);
   }
       
-  if(AliTRDReconstructor::StreamLevel() >= 2){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 2){
     TTreeSRedirector &cstreamer = *fgDebugStreamer;
     Int_t eventNumber = AliTRDtrackerDebug::GetEventNumber();
     Int_t candidateNumber = AliTRDtrackerDebug::GetCandidateNumber();
@@ -2395,24 +2428,20 @@ AliTRDtrackV1* AliTRDtrackerV1::MakeTrack(AliTRDseedV1 *seeds, Double_t *params)
   c[ 6] = 0.0; c[ 7] = 0.0; c[ 8] = 0.0;  c[ 9] = 0.1;
   c[10] = 0.0; c[11] = 0.0; c[12] = 0.0;  c[13] = 0.0; c[14] = params[5]*params[5]*0.01;
 
-  AliTRDtrackV1 *track = new AliTRDtrackV1(seeds, &params[1], c, params[0], params[6]*alpha+shift);
-  track->PropagateTo(params[0]-5.0);
-  track->ResetCovariance(1);
-  Int_t nc = FollowBackProlongation(*track);
-  if (nc < 30) {
-    delete track;
-    track = 0x0;
-  } else {
-    //track->CookdEdx();
-    //track->CookdEdxTimBin(-1);
-    track->CookLabel(.9);
-    // computes PID for track
-    track->CookPID();
-    // update calibration references using this track
-    if(calibra->GetHisto2d()) calibra->UpdateHistogramsV1(track);
-  }
-
-  return track;
+  AliTRDtrackV1 track(seeds, &params[1], c, params[0], params[6]*alpha+shift);
+  track.PropagateTo(params[0]-5.0);
+  track.ResetCovariance(1);
+  Int_t nc = FollowBackProlongation(track);
+  if (nc < 30) return 0x0;
+
+  AliTRDtrackV1 *ptrTrack = SetTrack(&track);
+  ptrTrack->CookLabel(.9);
+  // computes PID for track
+  ptrTrack->CookPID();
+  // update calibration references using this track
+  if(calibra->GetHisto2d()) calibra->UpdateHistogramsV1(ptrTrack);
+
+  return ptrTrack;
 }
 
 
@@ -2473,7 +2502,7 @@ Int_t AliTRDtrackerV1::ImproveSeedQuality(AliTRDtrackingChamber **stack, AliTRDs
     }
 
     chi2 = FitTiltedRieman(bseed, kTRUE);
-    if(AliTRDReconstructor::StreamLevel() >= 7){
+    if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 7){
       Int_t eventNumber                = AliTRDtrackerDebug::GetEventNumber();
       Int_t candidateNumber = AliTRDtrackerDebug::GetCandidateNumber();
       TLinearFitter *tiltedRieman = GetTiltedRiemanFitter();
@@ -2526,12 +2555,13 @@ Double_t AliTRDtrackerV1::CalculateTrackLikelihood(AliTRDseedV1 *tracklets, Doub
   sumdaf /= Float_t (nLayers - 2.0);
   
   Double_t likeChi2Z  = TMath::Exp(-chi2[2] * 0.14);                   // Chi2Z 
-  Double_t likeChi2TC = TMath::Exp(-chi2[1] * 0.677);                  // Constrained Tilted Riemann
+  Double_t likeChi2TC = (AliTRDReconstructor::RecoParam()->IsVertexConstrained()) ? 
+                                                                                       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 likeAF     = TMath::Exp(-sumdaf * 3.23);
   Double_t trackLikelihood     = likeChi2Z * likeChi2TR * likeAF;
 
-  if(AliTRDReconstructor::StreamLevel() >= 2){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 2){
     Int_t eventNumber = AliTRDtrackerDebug::GetEventNumber();
     Int_t candidateNumber = AliTRDtrackerDebug::GetCandidateNumber();
     TTreeSRedirector &cstreamer = *fgDebugStreamer;
@@ -2598,7 +2628,7 @@ Double_t AliTRDtrackerV1::CookLikelihood(AliTRDseedV1 *cseed, Int_t planes[4]
   Double_t like      = likea * likechi2y * likechi2z * likeN;
 
   //   AliInfo(Form("sumda(%f) chi2[0](%f) chi2[1](%f) likea(%f) likechi2y(%f) likechi2z(%f) nclusters(%d) likeN(%f)", sumda, chi2[0], chi2[1], likea, likechi2y, likechi2z, nclusters, likeN));
-  if(AliTRDReconstructor::StreamLevel() >= 2){
+  if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 2){
     Int_t eventNumber = AliTRDtrackerDebug::GetEventNumber();
     Int_t candidateNumber = AliTRDtrackerDebug::GetCandidateNumber();
     // The Debug Stream contains the seed 
index 042362d..7fed733 100644 (file)
@@ -130,6 +130,7 @@ protected:
   Int_t          MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *sseed, Int_t *ipar);
   AliTRDtrackV1* MakeTrack(AliTRDseedV1 *seeds, Double_t *params);
   static Int_t   PropagateToX(AliTRDtrackV1 &t, Double_t xToGo, Double_t maxStep);
+  AliTRDtrackV1* SetTrack(AliTRDtrackV1 *track);
   AliTRDseedV1*  SetTracklet(AliTRDseedV1 *tracklet);
 
 private:
index 2262930..5d3d1ae 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-
-/* $Id: AliTRDtrackingDebug.cxx 23810 2008-02-08 09:00:27Z hristov $ */
-
-////////////////////////////////////////////////////////////////////////////
-//                                                                        //
-//  Tracking in one chamber                                               //
-//                                                                        //
-//  Authors:                                                              //
-//    Alex Bercuci <A.Bercuci@gsi.de>                                     //
-//    Markus Fasel <M.Fasel@gsi.de>                                       //
-//                                                                        //
-////////////////////////////////////////////////////////////////////////////
-
-#include "AliTRDtrackingChamber.h"
-
-#include "TMath.h"
-#include "TMatrixTBase.h"
-#include <TTreeStream.h>
-
-#include "AliTRDReconstructor.h"
-#include "AliTRDrecoParam.h"
-#include "AliTRDtrackerV1.h"
-#include "AliTRDgeometry.h"
-#include "AliTRDpadPlane.h"
-#include "AliTRDcalibDB.h"
-
-ClassImp(AliTRDtrackingChamber)
-
-//_______________________________________________________
-AliTRDtrackingChamber::AliTRDtrackingChamber(Int_t det) :
-       fDetector(det)
-       ,fX0(0.)
-{}  
-
-//_______________________________________________________
-void AliTRDtrackingChamber::Clear(const Option_t *opt)
-{
-       for(Int_t itb=0; itb<kNTimeBins; itb++) fTB[itb].Clear(opt);
-}
-
-//_______________________________________________________
-void AliTRDtrackingChamber::InsertCluster(AliTRDcluster *c, Int_t index)
-{
-       fTB[c->GetLocalTimeBin()].InsertCluster(c, index);
-}
-
-//_______________________________________________________
-Bool_t AliTRDtrackingChamber::Build(AliTRDgeometry *geo)
-{
-// Init chamber and all time bins (AliTRDchamberTimeBin)
-// Calculates radial position of the chamber based on 
-// radial positions of the time bins (calibration/alignment aware)
-//
-       Int_t stack = geo->GetChamber(fDetector);
-       Int_t plane = geo->GetPlane(fDetector);
-       AliTRDpadPlane *pp = geo->GetPadPlane(plane, stack);
-       Double_t zl = pp->GetRow0ROC() - pp->GetRowEndROC();
-       Double_t z0 = geo->GetRow0(plane, stack, 0) - zl;
-       Int_t nrows = pp->GetNrows();
-       
-       Int_t index[50], jtb = 0;
-       for(Int_t itb=0; itb<kNTimeBins; itb++){ 
-               if(!fTB[itb]) continue;
-               fTB[itb].SetRange(z0, zl);
-               fTB[itb].SetNRows(nrows);
-               fTB[itb].BuildIndices();
-               index[jtb++] = itb;
-       }       
-       if(jtb<2) return kFALSE;
-       
-       
-       // ESTIMATE POSITION OF PAD PLANE FOR THIS CHAMBER
-       Double_t x0 = fTB[index[0]].GetX();
-       Double_t x1 = fTB[index[1]].GetX();
-       Double_t dx = (x0 - x1)/(index[1] - index[0]); 
-       fX0 = x0 + dx*(index[0] - (Int_t)AliTRDcalibDB::Instance()->GetT0Average(fDetector));   
-       return kTRUE;
-}
-       
-//_______________________________________________________      
-Int_t AliTRDtrackingChamber::GetNClusters() const
-{
-// Returns number of clusters in chamber
-//
-       Int_t n = 0;
-       for(Int_t itb=0; itb<kNTimeBins; itb++){ 
-               n += Int_t(fTB[itb]);
-       }
-       return n;       
-}      
-
-//_______________________________________________________
-Double_t AliTRDtrackingChamber::GetQuality()
-{
-  //
-  // Calculate chamber quality for seeding.
-  // 
-  //
-  // Parameters :
-  //   layers : Array of propagation layers for this plane.
-  //
-  // Output :
-  //   plane quality factor for seeding
-  // 
-  // Detailed description
-  //
-  // The quality of the plane for seeding is higher if:
-  //  1. the average timebin population is closer to an integer number
-  //  2. the distribution of clusters/timebin is closer to a uniform distribution.
-  //    - the slope of the first derivative of a parabolic fit is small or
-  //    - the slope of a linear fit is small
-  //
-
-       Int_t ncl   = 0;
-       Int_t nused = 0;
-       Int_t nClLayer;
-       for(int itb=0; itb<kNTimeBins; itb++){
-               if(!(nClLayer = fTB[itb].GetNClusters())) continue;
-               ncl += nClLayer;
-               for(Int_t incl = 0; incl < nClLayer; incl++){
-                       if((fTB[itb].GetCluster(incl))->IsUsed()) nused++;
-               }
-       }
-       
-       // calculate the deviation of the mean number of clusters from the
-       // closest integer values
-       Float_t nclMed = float(ncl-nused)/AliTRDtrackerV1::GetNTimeBins();
-       Int_t ncli = Int_t(nclMed);
-       Float_t nclDev = TMath::Abs(nclMed - TMath::Max(ncli, 1));
-       nclDev -= (nclDev>.5) && ncli ? 1. : 0.;
-       return TMath::Exp(-5.*TMath::Abs(nclDev));
-
-//     // get slope of the derivative
-//     if(!fitter.Eval()) return quality;
-//     fitter.PrintResults(3);
-//     Double_t a = fitter.GetParameter(1);
-// 
-//     printf("ncl_dev(%f)  a(%f)\n", ncl_dev, a);
-//     return quality*TMath::Exp(-a);
-
-}
-
-
-//_______________________________________________________
-AliTRDchamberTimeBin *AliTRDtrackingChamber::GetSeedingLayer(AliTRDgeometry *geo)
-{
-  //
-  // Creates a seeding layer
-  //
-       
-       // constants
-       const Int_t kMaxRows = 16;
-       const Int_t kMaxCols = 144;
-       const Int_t kMaxPads = 2304;
-               
-       // Get the geometrical data of the chamber
-       Int_t plane = geo->GetPlane(fDetector);
-       Int_t stack = geo->GetChamber(fDetector);
-       Int_t sector= geo->GetSector(fDetector);
-       AliTRDpadPlane *pp = geo->GetPadPlane(plane, stack);
-       Int_t nCols = pp->GetNcols();
-       Float_t ymin = TMath::Min(pp->GetCol0(), pp->GetColEnd());
-       Float_t ymax = TMath::Max(pp->GetCol0(), pp->GetColEnd());
-       Float_t zmin = TMath::Min(pp->GetRow0(), pp->GetRowEnd());
-       Float_t zmax = TMath::Max(pp->GetRow0(), pp->GetRowEnd());
-       Float_t z0 = -1., zl = -1.;
-       Int_t nRows = pp->GetNrows();
-       Float_t binlength = (ymax - ymin)/nCols; 
-       //AliInfo(Form("ymin(%f) ymax(%f) zmin(%f) zmax(%f) nRows(%d) binlength(%f)", ymin, ymax, zmin, zmax, nRows, binlength));
-       
-       // Fill the histogram
-       Int_t nClusters;        
-       Int_t *histogram[kMaxRows];                                                                                     // 2D-Histogram
-       Int_t hvals[kMaxPads];  memset(hvals, 0, sizeof(Int_t)*kMaxPads);       
-       Float_t *sigmas[kMaxRows];
-       Float_t svals[kMaxPads];        memset(svals, 0, sizeof(Float_t)*kMaxPads);     
-       AliTRDcluster *c = 0x0;
-       for(Int_t irs = 0; irs < kMaxRows; irs++){
-               histogram[irs] = &hvals[irs*kMaxCols];
-               sigmas[irs] = &svals[irs*kMaxCols];
-       }
-       for(Int_t iTime = 0; iTime < kNTimeBins; iTime++){
-               if(!(nClusters = fTB[iTime].GetNClusters())) continue;
-               z0 = fTB[iTime].GetZ0();
-               zl = fTB[iTime].GetDZ0();
-               for(Int_t incl = 0; incl < nClusters; incl++){
-                       c = fTB[iTime].GetCluster(incl);        
-                       histogram[c->GetPadRow()][c->GetPadCol()]++;
-                       sigmas[c->GetPadRow()][c->GetPadCol()] += c->GetSigmaZ2();
-               }
-       }
-       
-// Now I have everything in the histogram, do the selection
-       //Int_t nPads = nCols * nRows;
-       // This is what we are interested in: The center of gravity of the best candidates
-       Float_t cogyvals[kMaxPads]; memset(cogyvals, 0, sizeof(Float_t)*kMaxPads);
-       Float_t cogzvals[kMaxPads]; memset(cogzvals, 0, sizeof(Float_t)*kMaxPads);
-       Float_t *cogy[kMaxRows];
-       Float_t *cogz[kMaxRows];
-       
-       // Lookup-Table storing coordinates according to the bins
-       Float_t yLengths[kMaxCols];
-       Float_t zLengths[kMaxRows];
-       for(Int_t icnt = 0; icnt < nCols; icnt++){
-               yLengths[icnt] = pp->GetColPos(nCols - 1 - icnt) + binlength/2;
-       }
-       for(Int_t icnt = 0; icnt < nRows; icnt++){
-               zLengths[icnt] = pp->GetRowPos(icnt) - pp->GetRowSize(icnt)/2;
-       }
-
-       // A bitfield is used to mask the pads as usable
-       Short_t mask[kMaxCols]; memset(mask, 0 ,sizeof(Short_t) * kMaxCols);//bool mvals[kMaxPads];
-       for(UChar_t icount = 0; icount < nRows; icount++){
-               cogy[icount] = &cogyvals[icount*kMaxCols];
-               cogz[icount] = &cogzvals[icount*kMaxCols];
-       }
-       // In this array the array position of the best candidates will be stored
-       Int_t   cand[AliTRDtrackerV1::kMaxTracksStack];
-       Float_t sigcands[AliTRDtrackerV1::kMaxTracksStack];
-       
-       // helper variables
-       Int_t indices[kMaxPads]; memset(indices, -1, sizeof(Int_t)*kMaxPads);
-       Int_t nCandidates = 0;
-       Float_t norm, cogv;
-       // histogram filled -> Select best bins
-       Int_t nPads = nCols * nRows;
-       TMath::Sort(nPads, hvals, indices);                     // bins storing a 0 should not matter
-       // Set Threshold
-       Int_t maximum = hvals[indices[0]];      // best
-       Int_t threshold = Int_t(maximum * AliTRDReconstructor::RecoParam()->GetFindableClusters());
-       Int_t col, row, lower, lower1, upper, upper1;
-       for(Int_t ib = 0; ib < nPads; ib++){
-               if(nCandidates >= AliTRDtrackerV1::kMaxTracksStack){
-                       printf("Number of seed candidates %d exceeded maximum allowed per stack %d", nCandidates, AliTRDtrackerV1::kMaxTracksStack);
-                       break;
-               }
-               // Positions
-               row = indices[ib]/nCols;
-               col = indices[ib]%nCols;
-               // here will be the threshold condition:
-               if((mask[col] & (1 << row)) != 0) continue;             // Pad is masked: continue
-               if(histogram[row][col] < TMath::Max(threshold, 1)){     // of course at least one cluster is needed
-                       break;                  // number of clusters below threshold: break;
-               } 
-               // passing: Mark the neighbors
-               lower  = TMath::Max(col - 1, 0); upper  = TMath::Min(col + 2, nCols);
-               lower1 = TMath::Max(row - 1, 0); upper1 = TMath::Min(row + 2, nCols);
-               for(Int_t ic = lower; ic < upper; ++ic)
-                       for(Int_t ir = lower1; ir < upper1; ++ir){
-                               if(ic == col && ir == row) continue;
-                               mask[ic] |= (1 << ir);
-                       }
-               // Storing the position in an array
-               // testing for neigboring
-               cogv = 0;
-               norm = 0;
-               lower = TMath::Max(col - 1, 0);
-               upper = TMath::Min(col + 2, nCols);
-               for(Int_t inb = lower; inb < upper; ++inb){
-                       cogv += yLengths[inb] * histogram[row][inb];
-                       norm += histogram[row][inb];
-               }
-               cogy[row][col] = cogv / norm;
-               cogv = 0; norm = 0;
-               lower = TMath::Max(row - 1, 0);
-               upper = TMath::Min(row + 2, nRows);
-               for(Int_t inb = lower; inb < upper; ++inb){
-                       cogv += zLengths[inb] * histogram[inb][col];
-                       norm += histogram[inb][col];
-               }
-               cogz[row][col] = Float_t(cogv) /  norm;
-               // passed the filter
-               cand[nCandidates] = row*nCols + col;    // store the position of a passig candidate into an Array
-               sigcands[nCandidates] = sigmas[row][col] / histogram[row][col]; // never be a floating point exeption
-               // Analysis output
-               nCandidates++;
-       }
-       if(!nCandidates) return 0x0;
-       
-       Float_t pos[3], sig[2];
-       Short_t signal[7]; memset(&signal[0], 0, 7*sizeof(Short_t));
-       AliTRDchamberTimeBin *fakeLayer = new AliTRDchamberTimeBin(plane, stack, sector, z0, zl);
-       AliTRDcluster *cluster = 0x0;
-       if(nCandidates){
-               UInt_t fakeIndex = 0;
-               for(Int_t ican = 0; ican < nCandidates; ican++){
-                       row = cand[ican] / nCols;
-                       col = cand[ican] % nCols;
-                       //temporary
-                       Int_t n = 0; Double_t x = 0., y = 0., z = 0.;
-                       for(int itb=0; itb<kNTimeBins; itb++){
-                               if(!(nClusters = fTB[itb].GetNClusters())) continue;
-                               for(Int_t incl = 0; incl < nClusters; incl++){
-                                       c = fTB[itb].GetCluster(incl);  
-                                       if(c->GetPadRow() != row) continue;
-                                       if(TMath::Abs(c->GetPadCol() - col) > 2) continue;
-                                       x += c->GetX();
-                                       y += c->GetY();
-                                       z += c->GetZ();
-                                       n++;
-                               }
-                       }
-                       pos[0] = x/n;
-                       pos[1] = y/n;
-                       pos[2] = z/n;
-                       sig[0] = .02;
-                       sig[1] = sigcands[ican];
-                       cluster = new AliTRDcluster(fDetector, 0., pos, sig, 0x0, 3, signal, col, row, 0, 0, 0., 0);
-                       fakeLayer->InsertCluster(cluster, fakeIndex++);
-               }
-       }
-       fakeLayer->SetNRows(nRows);
-       fakeLayer->SetOwner();
-       fakeLayer->BuildIndices();
-       //fakeLayer->PrintClusters();
-       
-       if(AliTRDReconstructor::StreamLevel() >= 3){
-               //TMatrixD hist(nRows, nCols);
-               //for(Int_t i = 0; i < nRows; i++)
-               //      for(Int_t j = 0; j < nCols; j++)
-               //              hist(i,j) = histogram[i][j];
-               TTreeSRedirector &cstreamer = *AliTRDtrackerV1::DebugStreamer();
-               cstreamer << "GetSeedingLayer"
-               << "plane="      << plane
-               << "ymin="       << ymin
-               << "ymax="       << ymax
-               << "zmin="       << zmin
-               << "zmax="       << zmax
-               << "L.="         << fakeLayer
-               //<< "Histogram.=" << &hist
-               << "\n";
-       }
-       
-       return fakeLayer;
-}
-
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+/* $Id: AliTRDtrackingDebug.cxx 23810 2008-02-08 09:00:27Z hristov $ */\r
+\r
+////////////////////////////////////////////////////////////////////////////\r
+//                                                                        //\r
+//  Tracking in one chamber                                               //\r
+//                                                                        //\r
+//  Authors:                                                              //\r
+//    Alex Bercuci <A.Bercuci@gsi.de>                                     //\r
+//    Markus Fasel <M.Fasel@gsi.de>                                       //\r
+//                                                                        //\r
+////////////////////////////////////////////////////////////////////////////\r
+\r
+#include "AliTRDtrackingChamber.h"\r
+\r
+#include "TMath.h"\r
+#include "TMatrixTBase.h"\r
+#include <TTreeStream.h>\r
+\r
+#include "AliTRDReconstructor.h"\r
+#include "AliTRDrecoParam.h"\r
+#include "AliTRDtrackerV1.h"\r
+#include "AliTRDgeometry.h"\r
+#include "AliTRDpadPlane.h"\r
+#include "AliTRDcalibDB.h"\r
+\r
+ClassImp(AliTRDtrackingChamber)\r
+\r
+//_______________________________________________________\r
+AliTRDtrackingChamber::AliTRDtrackingChamber(Int_t det) :\r
+       fDetector(det)\r
+       ,fX0(0.)\r
+{}  \r
+\r
+//_______________________________________________________\r
+void AliTRDtrackingChamber::Clear(const Option_t *opt)\r
+{\r
+       for(Int_t itb=0; itb<kNTimeBins; itb++) fTB[itb].Clear(opt);\r
+}\r
+\r
+//_______________________________________________________\r
+void AliTRDtrackingChamber::InsertCluster(AliTRDcluster *c, Int_t index)\r
+{\r
+       fTB[c->GetLocalTimeBin()].InsertCluster(c, index);\r
+}\r
+\r
+//_______________________________________________________\r
+Bool_t AliTRDtrackingChamber::Build(AliTRDgeometry *geo)\r
+{\r
+// Init chamber and all time bins (AliTRDchamberTimeBin)\r
+// Calculates radial position of the chamber based on \r
+// radial positions of the time bins (calibration/alignment aware)\r
+//\r
+       Int_t stack = geo->GetChamber(fDetector);\r
+       Int_t plane = geo->GetPlane(fDetector);\r
+       AliTRDpadPlane *pp = geo->GetPadPlane(plane, stack);\r
+       Double_t zl = pp->GetRow0ROC() - pp->GetRowEndROC();\r
+       Double_t z0 = geo->GetRow0(plane, stack, 0) - zl;\r
+       Int_t nrows = pp->GetNrows();\r
+       \r
+       Int_t index[50], jtb = 0;\r
+       for(Int_t itb=0; itb<kNTimeBins; itb++){ \r
+               if(!fTB[itb]) continue;\r
+               fTB[itb].SetRange(z0, zl);\r
+               fTB[itb].SetNRows(nrows);\r
+               fTB[itb].BuildIndices();\r
+               index[jtb++] = itb;\r
+       }       \r
+       if(jtb<2) return kFALSE;\r
+       \r
+       \r
+       // ESTIMATE POSITION OF PAD PLANE FOR THIS CHAMBER\r
+       Double_t x0 = fTB[index[0]].GetX();\r
+       Double_t x1 = fTB[index[1]].GetX();\r
+       Double_t dx = (x0 - x1)/(index[1] - index[0]); \r
+       fX0 = x0 + dx*(index[0] - (Int_t)AliTRDcalibDB::Instance()->GetT0Average(fDetector));   \r
+       return kTRUE;\r
+}\r
+       \r
+//_______________________________________________________      \r
+Int_t AliTRDtrackingChamber::GetNClusters() const\r
+{\r
+// Returns number of clusters in chamber\r
+//\r
+       Int_t n = 0;\r
+       for(Int_t itb=0; itb<kNTimeBins; itb++){ \r
+               n += Int_t(fTB[itb]);\r
+       }\r
+       return n;       \r
+}      \r
+\r
+//_______________________________________________________\r
+Double_t AliTRDtrackingChamber::GetQuality()\r
+{\r
+  //\r
+  // Calculate chamber quality for seeding.\r
+  // \r
+  //\r
+  // Parameters :\r
+  //   layers : Array of propagation layers for this plane.\r
+  //\r
+  // Output :\r
+  //   plane quality factor for seeding\r
+  // \r
+  // Detailed description\r
+  //\r
+  // The quality of the plane for seeding is higher if:\r
+  //  1. the average timebin population is closer to an integer number\r
+  //  2. the distribution of clusters/timebin is closer to a uniform distribution.\r
+  //    - the slope of the first derivative of a parabolic fit is small or\r
+  //    - the slope of a linear fit is small\r
+  //\r
+\r
+       Int_t ncl   = 0;\r
+       Int_t nused = 0;\r
+       Int_t nClLayer;\r
+       for(int itb=0; itb<kNTimeBins; itb++){\r
+               if(!(nClLayer = fTB[itb].GetNClusters())) continue;\r
+               ncl += nClLayer;\r
+               for(Int_t incl = 0; incl < nClLayer; incl++){\r
+                       if((fTB[itb].GetCluster(incl))->IsUsed()) nused++;\r
+               }\r
+       }\r
+       \r
+       // calculate the deviation of the mean number of clusters from the\r
+       // closest integer values\r
+       Float_t nclMed = float(ncl-nused)/AliTRDtrackerV1::GetNTimeBins();\r
+       Int_t ncli = Int_t(nclMed);\r
+       Float_t nclDev = TMath::Abs(nclMed - TMath::Max(ncli, 1));\r
+       nclDev -= (nclDev>.5) && ncli ? 1. : 0.;\r
+       return TMath::Exp(-5.*TMath::Abs(nclDev));\r
+\r
+//     // get slope of the derivative\r
+//     if(!fitter.Eval()) return quality;\r
+//     fitter.PrintResults(3);\r
+//     Double_t a = fitter.GetParameter(1);\r
+// \r
+//     printf("ncl_dev(%f)  a(%f)\n", ncl_dev, a);\r
+//     return quality*TMath::Exp(-a);\r
+\r
+}\r
+\r
+\r
+//_______________________________________________________\r
+AliTRDchamberTimeBin *AliTRDtrackingChamber::GetSeedingLayer(AliTRDgeometry *geo)\r
+{\r
+  //\r
+  // Creates a seeding layer\r
+  //\r
+       \r
+       // constants\r
+       const Int_t kMaxRows = 16;\r
+       const Int_t kMaxCols = 144;\r
+       const Int_t kMaxPads = 2304;\r
+               \r
+       // Get the geometrical data of the chamber\r
+       Int_t plane = geo->GetPlane(fDetector);\r
+       Int_t stack = geo->GetChamber(fDetector);\r
+       Int_t sector= geo->GetSector(fDetector);\r
+       AliTRDpadPlane *pp = geo->GetPadPlane(plane, stack);\r
+       Int_t nCols = pp->GetNcols();\r
+       Float_t ymin = TMath::Min(pp->GetCol0(), pp->GetColEnd());\r
+       Float_t ymax = TMath::Max(pp->GetCol0(), pp->GetColEnd());\r
+       Float_t zmin = TMath::Min(pp->GetRow0(), pp->GetRowEnd());\r
+       Float_t zmax = TMath::Max(pp->GetRow0(), pp->GetRowEnd());\r
+       Float_t z0 = -1., zl = -1.;\r
+       Int_t nRows = pp->GetNrows();\r
+       Float_t binlength = (ymax - ymin)/nCols; \r
+       //AliInfo(Form("ymin(%f) ymax(%f) zmin(%f) zmax(%f) nRows(%d) binlength(%f)", ymin, ymax, zmin, zmax, nRows, binlength));\r
+       \r
+       // Fill the histogram\r
+       Int_t nClusters;        \r
+       Int_t *histogram[kMaxRows];                                                                                     // 2D-Histogram\r
+       Int_t hvals[kMaxPads];  memset(hvals, 0, sizeof(Int_t)*kMaxPads);       \r
+       Float_t *sigmas[kMaxRows];\r
+       Float_t svals[kMaxPads];        memset(svals, 0, sizeof(Float_t)*kMaxPads);     \r
+       AliTRDcluster *c = 0x0;\r
+       for(Int_t irs = 0; irs < kMaxRows; irs++){\r
+               histogram[irs] = &hvals[irs*kMaxCols];\r
+               sigmas[irs] = &svals[irs*kMaxCols];\r
+       }\r
+       for(Int_t iTime = 0; iTime < kNTimeBins; iTime++){\r
+               if(!(nClusters = fTB[iTime].GetNClusters())) continue;\r
+               z0 = fTB[iTime].GetZ0();\r
+               zl = fTB[iTime].GetDZ0();\r
+               for(Int_t incl = 0; incl < nClusters; incl++){\r
+                       c = fTB[iTime].GetCluster(incl);        \r
+                       histogram[c->GetPadRow()][c->GetPadCol()]++;\r
+                       sigmas[c->GetPadRow()][c->GetPadCol()] += c->GetSigmaZ2();\r
+               }\r
+       }\r
+       \r
+// Now I have everything in the histogram, do the selection\r
+       //Int_t nPads = nCols * nRows;\r
+       // This is what we are interested in: The center of gravity of the best candidates\r
+       Float_t cogyvals[kMaxPads]; memset(cogyvals, 0, sizeof(Float_t)*kMaxPads);\r
+       Float_t cogzvals[kMaxPads]; memset(cogzvals, 0, sizeof(Float_t)*kMaxPads);\r
+       Float_t *cogy[kMaxRows];\r
+       Float_t *cogz[kMaxRows];\r
+       \r
+       // Lookup-Table storing coordinates according to the bins\r
+       Float_t yLengths[kMaxCols];\r
+       Float_t zLengths[kMaxRows];\r
+       for(Int_t icnt = 0; icnt < nCols; icnt++){\r
+               yLengths[icnt] = pp->GetColPos(nCols - 1 - icnt) + binlength/2;\r
+       }\r
+       for(Int_t icnt = 0; icnt < nRows; icnt++){\r
+               zLengths[icnt] = pp->GetRowPos(icnt) - pp->GetRowSize(icnt)/2;\r
+       }\r
+\r
+       // A bitfield is used to mask the pads as usable\r
+       Short_t mask[kMaxCols]; memset(mask, 0 ,sizeof(Short_t) * kMaxCols);//bool mvals[kMaxPads];\r
+       for(UChar_t icount = 0; icount < nRows; icount++){\r
+               cogy[icount] = &cogyvals[icount*kMaxCols];\r
+               cogz[icount] = &cogzvals[icount*kMaxCols];\r
+       }\r
+       // In this array the array position of the best candidates will be stored\r
+       Int_t   cand[AliTRDtrackerV1::kMaxTracksStack];\r
+       Float_t sigcands[AliTRDtrackerV1::kMaxTracksStack];\r
+       \r
+       // helper variables\r
+       Int_t indices[kMaxPads]; memset(indices, -1, sizeof(Int_t)*kMaxPads);\r
+       Int_t nCandidates = 0;\r
+       Float_t norm, cogv;\r
+       // histogram filled -> Select best bins\r
+       Int_t nPads = nCols * nRows;\r
+       TMath::Sort(nPads, hvals, indices);                     // bins storing a 0 should not matter\r
+       // Set Threshold\r
+       Int_t maximum = hvals[indices[0]];      // best\r
+       Int_t threshold = Int_t(maximum * AliTRDReconstructor::RecoParam()->GetFindableClusters());\r
+       Int_t col, row, lower, lower1, upper, upper1;\r
+       for(Int_t ib = 0; ib < nPads; ib++){\r
+               if(nCandidates >= AliTRDtrackerV1::kMaxTracksStack){\r
+                       printf("Number of seed candidates %d exceeded maximum allowed per stack %d", nCandidates, AliTRDtrackerV1::kMaxTracksStack);\r
+                       break;\r
+               }\r
+               // Positions\r
+               row = indices[ib]/nCols;\r
+               col = indices[ib]%nCols;\r
+               // here will be the threshold condition:\r
+               if((mask[col] & (1 << row)) != 0) continue;             // Pad is masked: continue\r
+               if(histogram[row][col] < TMath::Max(threshold, 1)){     // of course at least one cluster is needed\r
+                       break;                  // number of clusters below threshold: break;\r
+               } \r
+               // passing: Mark the neighbors\r
+               lower  = TMath::Max(col - 1, 0); upper  = TMath::Min(col + 2, nCols);\r
+               lower1 = TMath::Max(row - 1, 0); upper1 = TMath::Min(row + 2, nCols);\r
+               for(Int_t ic = lower; ic < upper; ++ic)\r
+                       for(Int_t ir = lower1; ir < upper1; ++ir){\r
+                               if(ic == col && ir == row) continue;\r
+                               mask[ic] |= (1 << ir);\r
+                       }\r
+               // Storing the position in an array\r
+               // testing for neigboring\r
+               cogv = 0;\r
+               norm = 0;\r
+               lower = TMath::Max(col - 1, 0);\r
+               upper = TMath::Min(col + 2, nCols);\r
+               for(Int_t inb = lower; inb < upper; ++inb){\r
+                       cogv += yLengths[inb] * histogram[row][inb];\r
+                       norm += histogram[row][inb];\r
+               }\r
+               cogy[row][col] = cogv / norm;\r
+               cogv = 0; norm = 0;\r
+               lower = TMath::Max(row - 1, 0);\r
+               upper = TMath::Min(row + 2, nRows);\r
+               for(Int_t inb = lower; inb < upper; ++inb){\r
+                       cogv += zLengths[inb] * histogram[inb][col];\r
+                       norm += histogram[inb][col];\r
+               }\r
+               cogz[row][col] = Float_t(cogv) /  norm;\r
+               // passed the filter\r
+               cand[nCandidates] = row*nCols + col;    // store the position of a passig candidate into an Array\r
+               sigcands[nCandidates] = sigmas[row][col] / histogram[row][col]; // never be a floating point exeption\r
+               // Analysis output\r
+               nCandidates++;\r
+       }\r
+       if(!nCandidates) return 0x0;\r
+       \r
+       Float_t pos[3], sig[2];\r
+       Short_t signal[7]; memset(&signal[0], 0, 7*sizeof(Short_t));\r
+       AliTRDchamberTimeBin *fakeLayer = new AliTRDchamberTimeBin(plane, stack, sector, z0, zl);\r
+       AliTRDcluster *cluster = 0x0;\r
+       if(nCandidates){\r
+               UInt_t fakeIndex = 0;\r
+               for(Int_t ican = 0; ican < nCandidates; ican++){\r
+                       row = cand[ican] / nCols;\r
+                       col = cand[ican] % nCols;\r
+                       //temporary\r
+                       Int_t n = 0; Double_t x = 0., y = 0., z = 0.;\r
+                       for(int itb=0; itb<kNTimeBins; itb++){\r
+                               if(!(nClusters = fTB[itb].GetNClusters())) continue;\r
+                               for(Int_t incl = 0; incl < nClusters; incl++){\r
+                                       c = fTB[itb].GetCluster(incl);  \r
+                                       if(c->GetPadRow() != row) continue;\r
+                                       if(TMath::Abs(c->GetPadCol() - col) > 2) continue;\r
+                                       x += c->GetX();\r
+                                       y += c->GetY();\r
+                                       z += c->GetZ();\r
+                                       n++;\r
+                               }\r
+                       }\r
+                       pos[0] = x/n;\r
+                       pos[1] = y/n;\r
+                       pos[2] = z/n;\r
+                       sig[0] = .02;\r
+                       sig[1] = sigcands[ican];\r
+                       cluster = new AliTRDcluster(fDetector, 0., pos, sig, 0x0, 3, signal, col, row, 0, 0, 0., 0);\r
+                       fakeLayer->InsertCluster(cluster, fakeIndex++);\r
+               }\r
+       }\r
+       fakeLayer->SetNRows(nRows);\r
+       fakeLayer->SetOwner();\r
+       fakeLayer->BuildIndices();\r
+       //fakeLayer->PrintClusters();\r
+       \r
+       if(AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 3){\r
+               //TMatrixD hist(nRows, nCols);\r
+               //for(Int_t i = 0; i < nRows; i++)\r
+               //      for(Int_t j = 0; j < nCols; j++)\r
+               //              hist(i,j) = histogram[i][j];\r
+               TTreeSRedirector &cstreamer = *AliTRDtrackerV1::DebugStreamer();\r
+               cstreamer << "GetSeedingLayer"\r
+               << "plane="      << plane\r
+               << "ymin="       << ymin\r
+               << "ymax="       << ymax\r
+               << "zmin="       << zmin\r
+               << "zmax="       << zmax\r
+               << "L.="         << fakeLayer\r
+               //<< "Histogram.=" << &hist\r
+               << "\n";\r
+       }\r
+       \r
+       return fakeLayer;\r
+}\r
+\r