TRD/AliTRDchamberTimeBin.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4  * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes is hereby granted   *\r
   9  * without fee, provided that the above copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 \r
  16 /* $Id: AliTRDchamberTimeBin.cxx 23313 2008-01-11 14:56:43Z cblume $ */\r
  17 \r
  18 ///////////////////////////////////////////////////////////////////////////////\r
  19 //                                                                           //\r
  20 //  Organization of clusters at the level of 1 TRD chamber.                  //\r
  21 //  The data structure is used for tracking at the stack level.              //\r
  22 //                                                                           //\r
  23 //  Functionalities:                                                         //\r
  24 //   1. cluster organization and sorting                                     //\r
  25 //   2. fast data navigation                                                 //\r
  26 //                                                                           //\r
  27 //  Authors:                                                                 //\r
  28 //    Alex Bercuci <A.Bercuci@gsi.de>                                        //\r
  29 //    Markus Fasel <M.Fasel@gsi.de>                                          //\r
  30 //                                                                           //\r
  31 ///////////////////////////////////////////////////////////////////////////////\r
  32 \r
  33 #include <TObject.h>\r
  34 #include <TROOT.h>\r
  35 #include <TMath.h>\r
  36 #include <TStopwatch.h>\r
  37 #include <TTreeStream.h>\r
  38 \r
  39 #include "AliLog.h"\r
  40 \r
  41 #include "AliTRDcluster.h"\r
  42 #include "AliTRDchamberTimeBin.h"\r
  43 #include "AliTRDrecoParam.h"\r
  44 #include "AliTRDReconstructor.h"\r
  45 #include "AliTRDtrackerV1.h"\r
  46 \r
  47 \r
  48 ClassImp(AliTRDchamberTimeBin)\r
  49 \r
  50 //_____________________________________________________________________________\r
  51 AliTRDchamberTimeBin::AliTRDchamberTimeBin(Int_t plane, Int_t stack, Int_t sector, Double_t z0, Double_t zLength)\r
  52   :TObject()\r
  53   ,fOwner(kFALSE)\r
  54   ,fPlane(plane)\r
  55   ,fStack(stack)\r
  56   ,fSector(sector)\r
  57   ,fNRows(kMaxRows)\r
  58   ,fN(0)\r
  59   ,fX(0.)\r
  60   ,fZ0(z0)\r
  61   ,fZLength(zLength)\r
  62 {\r
  63   //\r
  64   // Default constructor (Only provided to use AliTRDchamberTimeBin with arrays)\r
  65   //\r
  66 \r
  67         for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r
  68         for(int ic=0; ic<kMaxClustersLayer; ic++){\r
  69                 fClusters[ic] = 0x0;\r
  70                 fIndex[ic]    = 0xffff;\r
  71         }\r
  72 }\r
  73 \r
  74 // //_____________________________________________________________________________\r
  75 // AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer, Double_t\r
  76 // z0, Double_t zLength, UChar_t stackNr):\r
  77 //      TObject()\r
  78 //      ,fOwner(kFALSE)\r
  79 //   ,fPlane(0xff)\r
  80 //   ,fStack(0xff)\r
  81 //   ,fSector(0xff)\r
  82 //   ,fNRows(kMaxRows)\r
  83 //   ,fN(0)\r
  84 //   ,fX(0.)\r
  85 //      ,fZ0(z0)\r
  86 //      ,fZLength(zLength)\r
  87 // {\r
  88 // // Standard constructor.\r
  89 // // Initialize also the underlying AliTRDpropagationLayer using the copy constructor.\r
  90 // \r
  91 //      SetT0(layer.IsT0());\r
  92 //      for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r
  93 //      for(int ic=0; ic<kMaxClustersLayer; ic++){\r
  94 //              fClusters[ic] = 0x0;\r
  95 //              fIndex[ic]    = 0xffff;\r
  96 //      }\r
  97 // }\r
  98 // \r
  99 // //_____________________________________________________________________________\r
 100 // AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer):\r
 101 //      TObject()\r
 102 //      ,fOwner(kFALSE)\r
 103 //   ,fPlane(0xff)\r
 104 //   ,fStack(0xff)\r
 105 //   ,fSector(0xff)\r
 106 //   ,fNRows(kMaxRows)\r
 107 //   ,fN(0)\r
 108 //   ,fX(0.)\r
 109 //      ,fZ0(0)\r
 110 //      ,fZLength(0)\r
 111 // {\r
 112 // // Standard constructor using only AliTRDpropagationLayer.\r
 113 //      \r
 114 //      SetT0(layer.IsT0());\r
 115 //      for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r
 116 //      for(int ic=0; ic<kMaxClustersLayer; ic++){\r
 117 //              fClusters[ic] = 0x0;\r
 118 //              fIndex[ic]    = 0xffff;\r
 119 //      }\r
 120 // }\r
 121 // //_____________________________________________________________________________\r
 122 // AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDpropagationLayer &layer)\r
 123 // {\r
 124 // // Assignment operator from an AliTRDpropagationLayer\r
 125 // \r
 126 //      if (this != &layer) layer.Copy(*this);\r
 127 //      return *this;\r
 128 // }\r
 129 // \r
 130 \r
 131 //_____________________________________________________________________________\r
 132 AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDchamberTimeBin &layer):\r
 133         TObject()\r
 134         ,fOwner(layer.fOwner)\r
 135   ,fPlane(layer.fPlane)\r
 136   ,fStack(layer.fStack)\r
 137   ,fSector(layer.fSector)\r
 138         ,fNRows(layer.fNRows)\r
 139   ,fN(layer.fN)\r
 140   ,fX(layer.fX)\r
 141         ,fZ0(layer.fZ0)\r
 142         ,fZLength(layer.fZLength)\r
 143 {\r
 144 // Copy Constructor (performs a deep copy)\r
 145         \r
 146         SetT0(layer.IsT0());\r
 147         for(int i=0; i<kMaxRows; i++) fPositions[i] = layer.fPositions[i];\r
 148         memcpy(&fClusters[0], &layer.fClusters[0], kMaxClustersLayer*sizeof(UChar_t));\r
 149         memcpy(&fIndex[0], &layer.fIndex[0], kMaxClustersLayer*sizeof(UInt_t));\r
 150 \r
 151 \r
 152 //      BuildIndices();\r
 153 }\r
 154 \r
 155 //_____________________________________________________________________________\r
 156 AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDchamberTimeBin &layer)\r
 157 {\r
 158 // Assignment operator\r
 159 \r
 160         if (this != &layer) layer.Copy(*this);\r
 161   return *this;\r
 162 }\r
 163 \r
 164 //_____________________________________________________________________________\r
 165 void AliTRDchamberTimeBin::Copy(TObject &o) const\r
 166 {\r
 167 // Copy method. Performs a deep copy of all data from this object to object o.\r
 168         \r
 169         AliTRDchamberTimeBin &layer = (AliTRDchamberTimeBin &)o;\r
 170         layer.fOwner       = kFALSE;\r
 171         layer.fPlane       = fPlane;\r
 172         layer.fStack       = fStack;\r
 173         layer.fSector      = fSector;\r
 174         layer.fNRows       = fNRows;\r
 175         layer.fN           = fN;\r
 176         layer.fX           = fX;\r
 177         layer.fZ0          = fZ0;\r
 178         layer.fZLength     = fZLength;\r
 179         layer.SetT0(IsT0());\r
 180         \r
 181         for(int i = 0; i < kMaxRows; i++) layer.fPositions[i] = 0;\r
 182         \r
 183         for(int i=0; i<kMaxRows; i++) layer.fPositions[i] = fPositions[i];\r
 184         memcpy(&layer.fClusters[0], &fClusters[0], kMaxClustersLayer*sizeof(UChar_t));\r
 185         memcpy(&layer.fIndex[0], &fIndex[0], kMaxClustersLayer*sizeof(UInt_t));\r
 186         \r
 187         TObject::Copy(layer); // copies everything into layer\r
 188         \r
 189 //      layer.BuildIndices();\r
 190 }\r
 191 \r
 192 //_____________________________________________________________________________\r
 193 AliTRDchamberTimeBin::~AliTRDchamberTimeBin()\r
 194 {\r
 195 // Destructor\r
 196         if(fOwner) for(int ic=0; ic<fN; ic++) delete fClusters[ic];\r
 197 }\r
 198 \r
 199 //_____________________________________________________________________________\r
 200 void AliTRDchamberTimeBin::SetRange(const Float_t z0, const Float_t zLength)\r
 201 {\r
 202 // Sets the range in z-direction\r
 203 //\r
 204 // Parameters:\r
 205 //   z0      : starting position of layer in the z direction\r
 206 //   zLength : length of layer in the z direction \r
 207 \r
 208         fZ0 = (z0 <= z0 + zLength) ? z0 : z0 + zLength;\r
 209         fZLength = TMath::Abs(zLength);\r
 210 }\r
 211 \r
 212 //_____________________________________________________________________________\r
 213 void AliTRDchamberTimeBin::InsertCluster(AliTRDcluster *c, UInt_t index) \r
 214 {\r
 215   //\r
 216   // Insert cluster in cluster array.\r
 217   // Clusters are sorted according to Y coordinate.  \r
 218   //\r
 219 \r
 220   //if (fTimeBinIndex < 0) { \r
 221     //AliWarning("Attempt to insert cluster into non-sensitive time bin!\n");\r
 222     //return;\r
 223   //}\r
 224 \r
 225   if (fN == (Int_t) kMaxClustersLayer) {\r
 226     //AliWarning("Too many clusters !\n"); \r
 227     return;\r
 228   }\r
 229 \r
 230   if (fN == 0) {\r
 231     fIndex[0]       = index; \r
 232     fClusters[fN++] = c; \r
 233     return;\r
 234   }\r
 235 \r
 236   Int_t i = Find(c->GetY());\r
 237   memmove(fClusters+i+1,fClusters+i,(fN-i)*sizeof(AliTRDcluster*));\r
 238   memmove(fIndex   +i+1,fIndex   +i,(fN-i)*sizeof(UInt_t)); \r
 239   fIndex[i]    = index; \r
 240   fClusters[i] = c; \r
 241   fN++;\r
 242 \r
 243 }  \r
 244 \r
 245 \r
 246 //_____________________________________________________________________________\r
 247 void AliTRDchamberTimeBin::BuildIndices(Int_t iter)\r
 248 {\r
 249 // Rearrangement of the clusters belonging to the propagation layer for the stack.\r
 250 //\r
 251 // Detailed description\r
 252 //\r
 253 // The array indices of all clusters in one PropagationLayer are stored in\r
 254 // array. The array is divided into several bins.\r
 255 // The clusters are sorted in increasing order of their y coordinate.\r
 256 //\r
 257 // Sorting algorithm: TreeSearch\r
 258 //\r
 259 \r
 260         if(!fN) return;\r
 261 \r
 262         // Select clusters that belong to the Stack\r
 263         Int_t nClStack = 0;                                     // Internal counter\r
 264         for(Int_t i = 0; i < fN; i++){\r
 265                 if(fClusters[i]->IsUsed()){\r
 266                         fClusters[i] = 0x0;\r
 267                         fIndex[i] = 0xffff;\r
 268                 } else nClStack++;\r
 269         }\r
 270         if(nClStack > kMaxClustersLayer) AliWarning(Form("Number of clusters in stack %d exceed buffer size %d. Truncating.", nClStack, kMaxClustersLayer));\r
 271                 \r
 272         // Nothing in this time bin. Reset indexes \r
 273         if(!nClStack){\r
 274                 fN = 0;\r
 275                 memset(&fPositions[0], 0xff, sizeof(UChar_t) * kMaxRows);\r
 276                 memset(&fClusters[0], 0x0, sizeof(AliTRDcluster*) * kMaxClustersLayer);\r
 277                 memset(&fIndex[0], 0xffff, sizeof(UInt_t) * kMaxClustersLayer);\r
 278                 return;\r
 279         }\r
 280         \r
 281         // Make a copy\r
 282         AliTRDcluster *helpCL[kMaxClustersLayer];\r
 283         Int_t helpInd[kMaxClustersLayer];\r
 284         nClStack = 0;\r
 285         for(Int_t i = 0; i < TMath::Min(fN, kMaxClustersLayer); i++){\r
 286                 if(!fClusters[i]) continue;\r
 287                 helpCL[nClStack]  = fClusters[i];\r
 288                 helpInd[nClStack] = fIndex[i];\r
 289                 fClusters[i]      = 0x0;\r
 290                 fIndex[i]         = 0xffff;\r
 291                 nClStack++;\r
 292         }\r
 293         \r
 294         // do clusters arrangement\r
 295         fX = 0.;\r
 296         fN =  nClStack;\r
 297         nClStack = 0;\r
 298         // Reset Positions array\r
 299         memset(fPositions, 0, sizeof(UChar_t)*kMaxRows);\r
 300         for(Int_t i = 0; i < fN; i++){\r
 301                 // boundary check\r
 302                 AliTRDcluster *cl = helpCL[i];\r
 303                 UChar_t rowIndex = cl->GetPadRow();\r
 304                 // Insert Leaf\r
 305                 Int_t pos = FindYPosition(cl->GetY(), rowIndex, i);\r
 306                 if(pos == -1){          // zbin is empty;\r
 307                         Int_t upper = (rowIndex == fNRows - 1) ? nClStack : fPositions[rowIndex + 1];\r
 308                         memmove(fClusters + upper + 1, fClusters + upper, (sizeof(AliTRDcluster *))*(nClStack-upper));\r
 309                         memmove(fIndex + upper + 1, fIndex + upper, (sizeof(UInt_t))*(nClStack-upper));\r
 310                         fClusters[upper] = cl;\r
 311                         fIndex[upper] = helpInd[i]; \r
 312                         // Move All pointer one position back\r
 313                         for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++;\r
 314                         nClStack++;\r
 315                 } else {                // zbin not empty\r
 316                         memmove(fClusters + pos + 2, fClusters + pos+1, (sizeof(AliTRDcluster *))*(nClStack-(pos+1)));\r
 317                         memmove(fIndex + pos + 2, fIndex + pos+1, (sizeof(UInt_t))*(nClStack-(pos+1)));\r
 318                         fClusters[pos + 1] = cl;        //fIndex[i];\r
 319                         fIndex[pos + 1] = helpInd[i];\r
 320                         // Move All pointer one position back\r
 321                         for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++; \r
 322                         nClStack++;\r
 323                 }\r
 324 \r
 325                 // calculate mean x\r
 326                 fX += cl->GetX(); \r
 327                 \r
 328                 // Debug Streaming\r
 329                 if(AliTRDtrackerV1::DebugStreamer() && AliTRDReconstructor::RecoParam()->GetStreamLevel() >= 3){\r
 330                         TTreeSRedirector &cstream = *AliTRDtrackerV1::DebugStreamer();\r
 331                         cstream << "BuildIndices"\r
 332                         << "Plane="    << fPlane\r
 333                         << "Stack="    << fStack\r
 334                         << "Sector="   << fSector\r
 335                         << "Iter="     << iter\r
 336                         << "C.="       << cl\r
 337                         << "rowIndex=" << rowIndex\r
 338                         << "\n";\r
 339                 }\r
 340         }\r
 341 \r
 342 //      AliInfo("Positions");\r
 343 //      for(int ir=0; ir<fNRows; ir++) printf("pos[%d] %d\n", ir, fPositions[ir]);\r
 344 \r
 345         fX /= fN;\r
 346 }\r
 347 \r
 348 //_____________________________________________________________________________\r
 349 Int_t AliTRDchamberTimeBin::Find(Float_t y) const\r
 350 {\r
 351   //\r
 352   // Returns index of the cluster nearest in Y    \r
 353   //\r
 354 \r
 355   if (fN <= 0) return 0;\r
 356   \r
 357   if (y <= fClusters[0]->GetY()) return 0;\r
 358   \r
 359   if (y >  fClusters[fN-1]->GetY()) return fN;\r
 360   \r
 361 \r
 362   Int_t b = 0;\r
 363   Int_t e = fN - 1;\r
 364   Int_t m = (b + e) / 2;\r
 365 \r
 366   for ( ; b < e; m = (b + e) / 2) {\r
 367     if (y > fClusters[m]->GetY()) b = m + 1;\r
 368     else e = m;\r
 369   }\r
 370 \r
 371   return m;\r
 372 }    \r
 373 \r
 374 //_____________________________________________________________________________\r
 375 Int_t AliTRDchamberTimeBin::FindYPosition(Double_t y, UChar_t z, Int_t nClusters) const\r
 376 {\r
 377 //\r
 378 // Tree search Algorithm to find the nearest left cluster for a given\r
 379 // y-position in a certain z-bin (in fact AVL-tree). \r
 380 // Making use of the fact that clusters are sorted in y-direction.\r
 381 //\r
 382 // Parameters:\r
 383 //   y : y position of the reference point in tracking coordinates\r
 384 //   z : z reference bin.\r
 385 //   nClusters : \r
 386 //\r
 387 // Output :\r
 388 // Index of the nearest left cluster in the StackLayer indexing (-1 if no clusters are found)\r
 389 //\r
 390 \r
 391         Int_t start = fPositions[z];            // starting Position of the bin\r
 392         Int_t upper = (Int_t)((z != fNRows - 1) ? fPositions[z+1] : nClusters); // ending Position of the bin \r
 393         Int_t end = upper - 1; // ending Position of the bin \r
 394         if(end < start) return -1; // Bin is empty\r
 395         Int_t middle = static_cast<Int_t>((start + end)/2);\r
 396         // 1st Part: climb down the tree: get the next cluster BEFORE ypos\r
 397         while(start + 1 < end){\r
 398                 if(y >= fClusters[middle]->GetY()) start = middle;\r
 399                 else end = middle;\r
 400                 middle = static_cast<Int_t>((start + end)/2);\r
 401         }\r
 402         if(y > fClusters[end]->GetY()) return end;\r
 403         return start;\r
 404 }\r
 405 \r
 406 //_____________________________________________________________________________\r
 407 Int_t AliTRDchamberTimeBin::FindNearestYCluster(Double_t y, UChar_t z) const\r
 408 {\r
 409 //\r
 410 // Tree search Algorithm to find the nearest cluster for a given\r
 411 // y-position in a certain z-bin (in fact AVL-tree). \r
 412 // Making use of the fact that clusters are sorted in y-direction.\r
 413 //\r
 414 // Parameters:\r
 415 //   y : y position of the reference point in tracking coordinates\r
 416 //   z : z reference bin.\r
 417 //\r
 418 // Output \r
 419 // Index of the nearest cluster in the StackLayer indexing (-1 if no clusters are found)\r
 420 //\r
 421 \r
 422         Int_t position = FindYPosition(y, z, fN);\r
 423         if(position == -1) return position; // bin empty\r
 424         // FindYPosition always returns the left Neighbor. We don't know if the left or the right Neighbor is nearest\r
 425         // to the Reference y-position, so test both\r
 426         Int_t upper = (Int_t)((z < fNRows-1) ? fPositions[z+1] : fN); // ending Position of the bin\r
 427         if((position + 1) < (upper)){\r
 428                 if(TMath::Abs(y - fClusters[position + 1]->GetY()) < TMath::Abs(y - fClusters[position]->GetY())) return position + 1;\r
 429                 else return position;\r
 430         }\r
 431         return position;\r
 432 }\r
 433 \r
 434 //_____________________________________________________________________________\r
 435 Int_t AliTRDchamberTimeBin::SearchNearestCluster(Double_t y, Double_t z, Double_t maxroady, Double_t maxroadz) const\r
 436 {\r
 437 //\r
 438 // Finds the nearest cluster from a given point in a defined range.\r
 439 // Distance is determined in a 2D space by the 2-Norm.\r
 440 //\r
 441 // Parameters :\r
 442 //   y : y position of the reference point in tracking coordinates\r
 443 //   z : z reference bin.\r
 444 //   maxroady : maximum searching distance in y direction\r
 445 //   maxroadz : maximum searching distance in z direction\r
 446 //\r
 447 // Output \r
 448 // Index of the nearest cluster in the StackLayer indexing (-1 if no cluster is found).\r
 449 // Cluster can be accessed with the operator[] or GetCluster(Int_t index)\r
 450 //\r
 451 // Detail description\r
 452 //\r
 453 // The following steps are perfomed:\r
 454 // 1. Get the expected z bins inside maxroadz.\r
 455 // 2. For each z bin find nearest y cluster.\r
 456 // 3. Select best candidate\r
 457 //\r
 458         Int_t   index   = -1;\r
 459         // initial minimal distance will be represented as ellipse: semi-major = z-direction\r
 460         // later 2-Norm will be used  \r
 461 //      Float_t nExcentricity = TMath::Sqrt(maxroadz*maxroadz - maxroad*maxroad)/maxroadz;\r
 462         Float_t mindist = maxroadz;\r
 463         \r
 464         // not very nice but unfortunately neccessarry: we have ho check the neighbors in both directions (+ and -) too. How \r
 465         // much neighbors depends on the Quotient maxroadz/fZLength   \r
 466         UChar_t maxRows = 3;\r
 467         UChar_t zpos[kMaxRows];\r
 468   // Float_t mindist = TMath::Sqrt(maxroad*maxroad + maxroadz*maxroadz);\r
 469 //      UChar_t myZbin = FindTreePosition(z, fZ0 + fZLength/2, fZLength/4, 8, 8, kFALSE);\r
 470         UChar_t myZbin = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - z)/fZLength * fNRows);\r
 471         if(z < fZ0) myZbin = fNRows - 1;\r
 472         if(z > fZ0 + fZLength) myZbin = 0;\r
 473         //printf("\n%f < %f < %f [%d]\n", fZ0, z, fZ0 + fZLength, myZbin);\r
 474         //for(int ic=0; ic<fN; ic++) printf("%d z = %f row %d\n", ic, fClusters[ic]->GetZ(), fClusters[ic]->GetPadRow());\r
 475 \r
 476         UChar_t nNeighbors = 0;\r
 477         for(UChar_t i = 0; i < maxRows; i++){\r
 478                 if((myZbin - 1 + i) < 0) continue;\r
 479                 if((myZbin - 1 + i) > fNRows - 1) break;\r
 480                 zpos[nNeighbors] = myZbin - 1 + i;\r
 481                 nNeighbors++;\r
 482         }\r
 483         Float_t ycl = 0, zcl = 0;\r
 484         for(UChar_t neighbor = 0; neighbor < nNeighbors; neighbor++){   // Always test the neighbors too\r
 485                 Int_t pos = FindNearestYCluster(y, zpos[neighbor]);\r
 486                 if(pos == -1) continue; // No cluster in bin\r
 487                 AliTRDcluster *c = (AliTRDcluster *) (fClusters[pos]);\r
 488                 if(c->IsUsed()) continue;               // we are only interested in unused clusters\r
 489                 ycl = c->GetY();\r
 490                 // Too far away in y-direction (Prearrangement)\r
 491                 if (TMath::Abs(ycl - y) > maxroady){ \r
 492                         //printf("y[%f] ycl[%f] roady[%f]\n", y, ycl, maxroady);\r
 493                         continue;\r
 494                 }\r
 495                 zcl = c->GetZ();\r
 496                 // Too far away in z-Direction\r
 497                 // (Prearrangement since we have not so many bins to test)\r
 498                 if (TMath::Abs(zcl - z) > maxroadz) continue;\r
 499                 \r
 500                 Float_t dist;           // distance defined as 2-Norm   \r
 501                 // if we havent found a Particle that is in the ellipse around (y,z) with maxroad as semi-minor and\r
 502                 // maxroadz as semi-major, we take the radius of the ellipse concerning the cluster as mindist, later we \r
 503                 // take the 2-Norm when we found a cluster inside the ellipse (The value 10000 is taken because it is surely\r
 504                 // large enough to be usable as an indicator whether we have found a nearer cluster or not)\r
 505 //              if(mindist > 10000.){\r
 506 //                      Float_t phi = ((zcl - z) == 0) ? TMath::Pi()/2 : TMath::ATan((ycl - y)/(zcl - z));\r
 507 //                      mindist = maxroad/TMath::Sqrt(1 - nExcentricity*nExcentricity * (TMath::Cos(phi))*(TMath::Cos(phi)));\r
 508 //              }\r
 509                 dist = TMath::Max(TMath::Abs(y-ycl),TMath::Abs(z-zcl)); // infinity Norm\r
 510 //              dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z));\r
 511                 if((Int_t)(dist * 100000) < (Int_t)(mindist * 100000)){\r
 512                 //if((dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z))) < mindist){\r
 513                         mindist = dist;\r
 514                         index   = pos;\r
 515                 }       \r
 516         }\r
 517         // This is the Array Position in fIndex2D of the Nearest cluster: if a\r
 518         // cluster is called, then the function has to retrieve the Information\r
 519         // which is Stored in the Array called, the function\r
 520         return index;\r
 521 }\r
 522 \r
 523 //_____________________________________________________________________________\r
 524 void AliTRDchamberTimeBin::BuildCond(AliTRDcluster *cl, Double_t *cond, UChar_t Layer, Double_t theta, Double_t phi)\r
 525 {\r
 526 // Helper function to calculate the area where to expect a cluster in THIS\r
 527 // layer. \r
 528 //\r
 529 // Parameters :\r
 530 //   cl    : \r
 531 //   cond  :\r
 532 //   Layer : \r
 533 //   theta : \r
 534 //   phi   : \r
 535 //\r
 536 // Detail description\r
 537 //\r
 538 // Helper function to calculate the area where to expect a cluster in THIS\r
 539 // layer. by using the information of a former cluster in another layer\r
 540 // and the angle in theta- and phi-direction between layer 0 and layer 3.\r
 541 // If the layer is zero, initial conditions are calculated. Otherwise a\r
 542 // linear interpolation is performed. \r
 543 //Begin_Html\r
 544 //<img src="gif/build_cond.gif">\r
 545 //End_Html\r
 546 //\r
 547 \r
 548         if(!AliTRDReconstructor::RecoParam()){\r
 549                 AliError("Reconstruction parameters not initialized.");\r
 550                 return;\r
 551         }\r
 552         \r
 553         if(Layer == 0){\r
 554                 cond[0] = cl->GetY();                   // center: y-Direction\r
 555                 cond[1] = cl->GetZ();                   // center: z-Direction\r
 556                 cond[2] = AliTRDReconstructor::RecoParam()->GetMaxPhi()   * (cl->GetX() - GetX()) + 1.0;                        // deviation: y-Direction\r
 557                 cond[3] = AliTRDReconstructor::RecoParam()->GetMaxTheta() * (cl->GetX() - GetX()) + 1.0;                        // deviation: z-Direction\r
 558         } else {\r
 559                 cond[0] = cl->GetY() + phi   * (GetX() - cl->GetX()); \r
 560                 cond[1] = cl->GetZ() + theta * (GetX() - cl->GetX());\r
 561                 cond[2] = AliTRDReconstructor::RecoParam()->GetRoad0y() + phi;\r
 562                 cond[3] = AliTRDReconstructor::RecoParam()->GetRoad0z();\r
 563         }\r
 564 }\r
 565 \r
 566 //_____________________________________________________________________________\r
 567 void AliTRDchamberTimeBin::GetClusters(Double_t *cond, Int_t *index, Int_t& ncl, Int_t BufferSize)\r
 568 {\r
 569 // Finds all clusters situated in this layer inside a rectangle  given by the center an ranges.\r
 570 //\r
 571 // Parameters :\r
 572 //   cond  :\r
 573 //   index : \r
 574 //   ncl :\r
 575 //   BufferSize   :\r
 576 //\r
 577 // Output :\r
 578 //\r
 579 // Detail description\r
 580 //\r
 581 // Function returs an array containing the indices in the stacklayer of\r
 582 // the clusters found an  the number of found clusters in the stacklayer\r
 583 \r
 584         ncl = 0;\r
 585         memset(index, 0, BufferSize*sizeof(Int_t));\r
 586         if(fN == 0) return;\r
 587                 \r
 588         //Boundary checks\r
 589         Double_t zvals[2];\r
 590         if(cond[1] - cond[3] > fZ0 + fZLength || cond[1] + cond[3] < fZ0) return; // We are outside of the chamvber\r
 591         zvals[0] = ((cond[1] - cond[3]) < fZ0) ? fZ0 : (cond[1] - cond[3]);\r
 592         zvals[1] = ((cond[1] + cond[3]) < fZ0 + fZLength) ? (cond[1] + cond[3]) : fZ0 + fZLength - 1.E-3;\r
 593 \r
 594         UChar_t zhi = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[0])/fZLength * fNRows);\r
 595         UChar_t zlo = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[1])/fZLength * fNRows);\r
 596         \r
 597         \r
 598 //      AliInfo(Form("yc[%f] zc[%f] dy[%f] dz[%f]", cond[0], cond[1], cond[2], cond[3]));\r
 599 //      PrintClusters();\r
 600 //      AliInfo(Form("zlo[%f] zhi[%f]", zvals[0], zvals[1]));\r
 601 //      AliInfo(Form("zlo[%d] zhi[%d]", zlo, zhi));\r
 602         \r
 603         //Preordering in Direction z saves a lot of loops (boundary checked)\r
 604         for(UChar_t z = zlo; z <= zhi; z++){\r
 605                 UInt_t upper = (z < fNRows-1) ? fPositions[z+1] : fN;\r
 606                 //AliInfo(Form("z[%d] y [%d %d]", z, fPositions[z], upper));\r
 607                 for(Int_t y = fPositions[z]; y < (Int_t)upper; y++){\r
 608                         if(ncl == BufferSize){\r
 609                                 AliWarning("Buffer size riched. Some clusters may be lost.");\r
 610                                 return; //Buffer filled\r
 611                         }\r
 612                         \r
 613                         if(fClusters[y]->GetY() > (cond[0] + cond[2])) break;                   // Abbortion conditions!!!\r
 614                         if(fClusters[y]->GetY() < (cond[0] - cond[2])) continue;        // Too small\r
 615                         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
 616                         index[ncl] = y;\r
 617                         ncl++;\r
 618                 }\r
 619         }\r
 620         if(ncl>fN) AliError(Form("Clusters found %d > %d (clusters in layer)", ncl, fN));\r
 621 }\r
 622 \r
 623 //_____________________________________________________________________________\r
 624 AliTRDcluster *AliTRDchamberTimeBin::GetNearestCluster(Double_t *cond)\r
 625 {\r
 626 // Function returning a pointer to the nearest cluster (nullpointer if not successfull).\r
 627 //\r
 628 // Parameters :\r
 629 //   cond  :\r
 630 //\r
 631 // Output :\r
 632 //   pointer to the nearest cluster (nullpointer if not successfull).\r
 633 // \r
 634 // Detail description\r
 635 //\r
 636 // returns a pointer to the nearest cluster (nullpointer if not\r
 637 // successfull) by the help of the method FindNearestCluster\r
 638         \r
 639         \r
 640         Double_t maxroad  = AliTRDReconstructor::RecoParam()->GetRoad2y();\r
 641         Double_t maxroadz = AliTRDReconstructor::RecoParam()->GetRoad2z();\r
 642         \r
 643         Int_t index = SearchNearestCluster(cond[0],cond[1],maxroad,maxroadz);\r
 644         AliTRDcluster *returnCluster = 0x0;\r
 645         if(index != -1) returnCluster = (AliTRDcluster *) fClusters[index];\r
 646         return returnCluster;\r
 647 }\r
 648 \r
 649 //_____________________________________________________________________________\r
 650 void AliTRDchamberTimeBin::PrintClusters() const\r
 651 {\r
 652 // Prints the position of each cluster in the stacklayer on the stdout\r
 653 //\r
 654         printf("\nnRows = %d\n", fNRows);\r
 655         printf("Z0 = %f\n", fZ0);\r
 656         printf("Z1 = %f\n", fZ0+fZLength);\r
 657         printf("clusters in AliTRDchamberTimeBin %d\n", fN);\r
 658         for(Int_t i = 0; i < fN; i++){\r
 659                 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
 660                 if(fClusters[i]->IsUsed()) printf("cluster allready used. rejected in search algorithm\n");\r
 661         }\r
 662 }\r