HLT/hough/AliL3Hough.cxx

   1 // @(#) $Id$
   2
   3 // Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
   4 //*-- Copyright &copy ALICE HLT Group
   5
   6 #include "AliL3StandardIncludes.h"
   7 #include <sys/time.h>
   8
   9 #include "AliL3Logging.h"
  10 #include "AliL3HoughMerger.h"
  11 #include "AliL3HoughIntMerger.h"
  12 #include "AliL3HoughGlobalMerger.h"
  13 #include "AliL3Histogram.h"
  14 #include "AliL3Hough.h"
  15 #include "AliL3HoughTransformer.h"
  16 #include "AliL3HoughClusterTransformer.h"
  17 #include "AliL3HoughTransformerLUT.h"
  18 #include "AliL3HoughTransformerVhdl.h"
  19 #include "AliL3HoughTransformerRow.h"
  20 #include "AliL3HoughMaxFinder.h"
  21 #include "AliL3Benchmark.h"
  22 #ifdef use_aliroot
  23 #include "AliL3FileHandler.h"
  24 #else
  25 #include "AliL3MemHandler.h"
  26 #endif
  27 #include "AliL3DataHandler.h"
  28 #include "AliL3DigitData.h"
  29 #include "AliL3HoughEval.h"
  30 #include "AliL3Transform.h"
  31 #include "AliL3TrackArray.h"
  32 #include "AliL3HoughTrack.h"
  33 #include "AliL3DDLDataFileHandler.h"
  34
  35 #include "TThread.h"
  36
  37 #if __GNUC__ == 3
  38 using namespace std;
  39 #endif
  40
  41 /** /class AliL3Hough
  42 //<pre>
  43 //_____________________________________________________________
  44 // AliL3Hough
  45 //
  46 // Interface class for the Hough transform
  47 //
  48 // Example how to use:
  49 //
  50 // AliL3Hough *hough = new AliL3Hough(path,kTRUE,NumberOfEtaSegments);
  51 // hough->ReadData(slice);
  52 // hough->Transform();
  53 // hough->FindTrackCandidates();
  54 //
  55 // AliL3TrackArray *tracks = hough->GetTracks(patch);
  56 //
  57 //</pre>
  58 */
  59
  60 ClassImp(AliL3Hough)
  61
  62 AliL3Hough::AliL3Hough()
  63 {
  64   //Constructor
  65
  66   fBinary        = kFALSE;
  67   fAddHistograms = kFALSE;
  68   fDoIterative   = kFALSE;
  69   fWriteDigits   = kFALSE;
  70   fUse8bits      = kFALSE;
  71
  72   fMemHandler       = 0;
  73   fHoughTransformer = 0;
  74   fEval             = 0;
  75   fPeakFinder       = 0;
  76   fTracks           = 0;
  77   fGlobalTracks     = 0;
  78   fMerger           = 0;
  79   fInterMerger      = 0;
  80   fGlobalMerger     = 0;
  81   fBenchmark        = 0;
  82
  83   fNEtaSegments     = 0;
  84   fNPatches         = 0;
  85   fLastPatch        =-1;
  86   fVersion          = 0;
  87   fCurrentSlice     = 0;
  88   fEvent            = 0;
  89
  90   fKappaSpread = 6;
  91   fPeakRatio   = 0.5;
  92   fInputFile   = 0;
  93   fInputPtr    = 0;
  94
  95   SetTransformerParams();
  96   SetThreshold();
  97   SetNSaveIterations();
  98   SetPeakThreshold();
  99 #ifdef use_aliroot
 100   //just be sure that index is empty for new event
 101     AliL3FileHandler::CleanStaticIndex();
 102 #ifdef use_newio
 103     fRunLoader = 0;
 104 #endif
 105 #endif
 106   fThread = 0;
 107 }
 108
 109 AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t netasegments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr)
 110 {
 111   //Normal constructor
 112   fBinary = binary;
 113   strcpy(fPath,path);
 114   fNEtaSegments  = netasegments;
 115   fAddHistograms = kFALSE;
 116   fDoIterative   = kFALSE;
 117   fWriteDigits   = kFALSE;
 118   fUse8bits      = bit8;
 119   fVersion       = tv;
 120   fKappaSpread=6;
 121   fPeakRatio=0.5;
 122   if(!fBinary) {
 123     if(infile) {
 124       fInputFile = infile;
 125       fInputPtr = 0;
 126     }
 127     else {
 128       fInputFile = 0;
 129       fInputPtr = ptr;
 130     }
 131   }
 132   else {
 133     fInputFile = 0;
 134     fInputPtr = 0;
 135   }
 136 #ifdef use_aliroot
 137   //just be sure that index is empty for new event
 138     AliL3FileHandler::CleanStaticIndex();
 139 #ifdef use_newio
 140     fRunLoader = 0;
 141 #endif
 142 #endif
 143   fThread = 0;
 144 }
 145
 146 AliL3Hough::~AliL3Hough()
 147 {
 148   //dtor
 149
 150   CleanUp();
 151   if(fMerger)
 152     delete fMerger;
 153   //cout << "Cleaned class merger " << endl;
 154   if(fInterMerger)
 155     delete fInterMerger;
 156   //cout << "Cleaned class inter " << endl;
 157   if(fPeakFinder)
 158     delete fPeakFinder;
 159   //cout << "Cleaned class peak " << endl;
 160   if(fGlobalMerger)
 161     delete fGlobalMerger;
 162   //cout << "Cleaned class global " << endl;
 163   if(fBenchmark)
 164     delete fBenchmark;
 165   //cout << "Cleaned class bench " << endl;
 166   if(fGlobalTracks)
 167     delete fGlobalTracks;
 168   //cout << "Cleaned class globaltracks " << endl;
 169   if(fThread) {
 170     //    fThread->Delete();
 171     delete fThread;
 172     fThread = 0;
 173   }
 174 }
 175
 176 void AliL3Hough::CleanUp()
 177 {
 178   //Cleanup memory
 179
 180   for(Int_t i=0; i<fNPatches; i++)
 181     {
 182       if(fTracks[i]) delete fTracks[i];
 183       //cout << "Cleaned tracks " << i << endl;
 184       if(fEval[i]) delete fEval[i];
 185       //cout << "Cleaned eval " << i << endl;
 186       if(fHoughTransformer[i]) delete fHoughTransformer[i];
 187       //cout << "Cleaned traf " << i << endl;
 188       if(fMemHandler[i]) delete fMemHandler[i];
 189       //cout << "Cleaned mem " << i << endl;
 190     }
 191
 192   if(fTracks) delete [] fTracks;
 193   //cout << "Cleaned class tracks " << endl;
 194   if(fEval) delete [] fEval;
 195   //cout << "Cleaned class eval " << endl;
 196   if(fHoughTransformer) delete [] fHoughTransformer;
 197   //cout << "Cleaned cleass trafo " << endl;
 198   if(fMemHandler) delete [] fMemHandler;
 199   //cout << "Cleaned class mem " << endl;
 200 }
 201
 202 void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t netasegments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr,Float_t zvertex)
 203 {
 204   //Normal init of the AliL3Hough
 205   fBinary = binary;
 206   strcpy(fPath,path);
 207   fNEtaSegments = netasegments;
 208   fWriteDigits  = kFALSE;
 209   fUse8bits     = bit8;
 210   fVersion      = tv;
 211   if(!fBinary) {
 212     if(infile) {
 213       fInputFile = infile;
 214       fInputPtr = 0;
 215     }
 216     else {
 217       fInputFile = 0;
 218       fInputPtr = ptr;
 219     }
 220   }
 221   else {
 222     fInputFile = 0;
 223     fInputPtr = 0;
 224   }
 225   fZVertex = zvertex;
 226
 227   Init(); //do the rest
 228 }
 229
 230 void AliL3Hough::Init(Bool_t doit, Bool_t addhists)
 231 {
 232   // Init
 233   fDoIterative   = doit;
 234   fAddHistograms = addhists;
 235
 236   fNPatches = AliL3Transform::GetNPatches();
 237   fHoughTransformer = new AliL3HoughBaseTransformer*[fNPatches];
 238   fMemHandler = new AliL3MemHandler*[fNPatches];
 239
 240   fTracks = new AliL3TrackArray*[fNPatches];
 241   fEval = new AliL3HoughEval*[fNPatches];
 242
 243   fGlobalTracks = new AliL3TrackArray("AliL3HoughTrack");
 244
 245   AliL3HoughBaseTransformer *lasttransformer = 0;
 246
 247   for(Int_t i=0; i<fNPatches; i++)
 248     {
 249       switch (fVersion){ //choose Transformer
 250       case 1:
 251         fHoughTransformer[i] = new AliL3HoughTransformerLUT(0,i,fNEtaSegments);
 252         break;
 253       case 2:
 254         fHoughTransformer[i] = new AliL3HoughClusterTransformer(0,i,fNEtaSegments);
 255         break;
 256       case 3:
 257         fHoughTransformer[i] = new AliL3HoughTransformerVhdl(0,i,fNEtaSegments,fNSaveIterations);
 258         break;
 259       case 4:
 260         fHoughTransformer[i] = new AliL3HoughTransformerRow(0,i,fNEtaSegments,kFALSE,fZVertex);
 261         break;
 262       default:
 263         fHoughTransformer[i] = new AliL3HoughTransformer(0,i,fNEtaSegments,kFALSE,kFALSE);
 264       }
 265
 266       fHoughTransformer[i]->SetLastTransformer(lasttransformer);
 267       lasttransformer = fHoughTransformer[i];
 268       //      fHoughTransformer[i]->CreateHistograms(fNBinX[i],fLowPt[i],fNBinY[i],-fPhi[i],fPhi[i]);
 269       fHoughTransformer[i]->CreateHistograms(fNBinX[i],-fLowPt[i],fLowPt[i],fNBinY[i],-fPhi[i],fPhi[i]);
 270       //fHoughTransformer[i]->CreateHistograms(fLowPt[i],fUpperPt[i],fPtRes[i],fNBinY[i],fPhi[i]);
 271
 272       fHoughTransformer[i]->SetLowerThreshold(fThreshold[i]);
 273       fHoughTransformer[i]->SetUpperThreshold(100);
 274
 275       LOG(AliL3Log::kInformational,"AliL3Hough::Init","Version")
 276         <<"Initializing Hough transformer version "<<fVersion<<ENDLOG;
 277
 278       fEval[i] = new AliL3HoughEval();
 279       fTracks[i] = new AliL3TrackArray("AliL3HoughTrack");
 280       if(fUse8bits)
 281         fMemHandler[i] = new AliL3DataHandler();
 282       else
 283 #ifdef use_aliroot
 284         {
 285           if(!fInputFile) {
 286             if(!fInputPtr) {
 287               /* In case of reading digits file */
 288               fMemHandler[i] = new AliL3FileHandler(kTRUE); //use static index
 289               if(!fBinary) {
 290 #if use_newio
 291                 if(!fRunLoader) {
 292 #endif
 293                   Char_t filename[1024];
 294                   sprintf(filename,"%s/digitfile.root",fPath);
 295                   fMemHandler[i]->SetAliInput(filename);
 296 #if use_newio
 297                 }
 298                 else {
 299                   fMemHandler[i]->SetAliInput(fRunLoader);
 300                 }
 301 #endif
 302               }
 303             }
 304             else {
 305               /* In case of reading from DATE */
 306               fMemHandler[i] = new AliL3DDLDataFileHandler();
 307               fMemHandler[i]->SetReaderInput(fInputPtr,-1);
 308             }
 309           }
 310           else {
 311             /* In case of reading rawdata from ROOT file */
 312             fMemHandler[i] = new AliL3DDLDataFileHandler();
 313             fMemHandler[i]->SetReaderInput(fInputFile);
 314           }
 315         }
 316 #else
 317       fMemHandler[i] = new AliL3MemHandler();
 318 #endif
 319     }
 320
 321   fPeakFinder = new AliL3HoughMaxFinder("KappaPhi",50000);
 322   fMerger = new AliL3HoughMerger(fNPatches);
 323   fInterMerger = new AliL3HoughIntMerger();
 324   fGlobalMerger = 0;
 325   fBenchmark = new AliL3Benchmark();
 326 }
 327
 328 void AliL3Hough::SetTransformerParams(Float_t ptres,Float_t ptmin,Float_t ptmax,Int_t ny,Int_t patch)
 329 {
 330   // Setup the parameters for the Hough Transformer
 331   Int_t mrow;
 332   Float_t psi=0;
 333   if(patch==-1)
 334     mrow = 80;
 335   else
 336     mrow = AliL3Transform::GetLastRow(patch);
 337   if(ptmin)
 338     {
 339       Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2));
 340       Double_t kappa = -1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/ptmin;
 341       psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2);
 342       cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl;
 343     }
 344
 345   if(patch==-1)
 346     {
 347       Int_t i=0;
 348       while(i < 6)
 349         {
 350           fPtRes[i] = ptres;
 351           fLowPt[i] = ptmin;
 352           fUpperPt[i] = ptmax;
 353           fNBinY[i] = ny;
 354           fPhi[i] = psi;
 355           fNBinX[i]=0;
 356           i++;
 357         }
 358       return;
 359     }
 360
 361   fPtRes[patch] = ptres;
 362   fLowPt[patch] = ptmin;
 363   fUpperPt[patch] = ptmax;
 364   fNBinY[patch] = ny;
 365   fPhi[patch] = psi;
 366 }
 367 /*
 368 void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t patch)
 369 {
 370   // Setup the parameters for the Hough Transformer
 371
 372   Int_t mrow=80;
 373   Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2));
 374   Double_t kappa = -1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/ptmin;
 375   Double_t psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2);
 376   cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl;
 377
 378   Int_t i=0;
 379   while(i < 6)
 380     {
 381       fLowPt[i] = ptmin;
 382       fNBinY[i] = ny;
 383       fNBinX[i] = nx;
 384       fPhi[i] = psi;
 385       i++;
 386     }
 387 }
 388 */
 389 void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t /*patch*/)
 390 {
 391   // Setup the parameters for the Hough Transformer
 392
 393
 394   Int_t mrow=79;
 395   Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2));
 396   Double_t alpha1 = AliL3Transform::GetMaxY(mrow)/pow(lineradius,2);
 397   Double_t kappa = 1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/ptmin;
 398   Double_t psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2);
 399   //  cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl;
 400   AliL3HoughTrack track;
 401   track.SetTrackParameters(kappa,psi,1);
 402   Float_t hit[3];
 403   Int_t mrow2 = 158;
 404   track.GetCrossingPoint(mrow2,hit);
 405   Double_t lineradius2 = sqrt(pow(AliL3Transform::Row2X(mrow2),2) + pow(AliL3Transform::GetMaxY(mrow2),2));
 406   Double_t alpha2 = hit[1]/pow(lineradius2,2);
 407   //  cout<<"Calculated alphas range "<<alpha1<<" "<<alpha2<<" in patch "<<patch<<endl;
 408
 409   Int_t i=0;
 410   while(i < 6)
 411     {
 412       fLowPt[i] = 1.15*alpha1;
 413       fNBinY[i] = ny;
 414       fNBinX[i] = nx;
 415       fPhi[i] = 1.15*alpha2;
 416       i++;
 417     }
 418 }
 419
 420 void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t lpt,Float_t phi)
 421 {
 422   Int_t i=0;
 423   while(i < 6)
 424     {
 425       fLowPt[i] = lpt;
 426       fNBinY[i] = ny;
 427       fNBinX[i] = nx;
 428       fPhi[i] = phi;
 429       i++;
 430     }
 431 }
 432
 433 void AliL3Hough::SetThreshold(Int_t t3,Int_t patch)
 434 {
 435   // Set digits threshold
 436   if(patch==-1)
 437     {
 438       Int_t i=0;
 439       while(i < 6)
 440         fThreshold[i++]=t3;
 441       return;
 442     }
 443   fThreshold[patch]=t3;
 444 }
 445
 446 void AliL3Hough::SetPeakThreshold(Int_t threshold,Int_t patch)
 447 {
 448   // Set Peak Finder threshold
 449   if(patch==-1)
 450     {
 451       Int_t i=0;
 452       while(i < 6)
 453         fPeakThreshold[i++]=threshold;
 454       return;
 455     }
 456   fPeakThreshold[patch]=threshold;
 457 }
 458
 459 void AliL3Hough::DoBench(Char_t *name)
 460 {
 461   fBenchmark->Analyze(name);
 462 }
 463
 464 void AliL3Hough::Process(Int_t minslice,Int_t maxslice)
 465 {
 466   //Process all slices [minslice,maxslice].
 467   fGlobalMerger = new AliL3HoughGlobalMerger(minslice,maxslice);
 468
 469   for(Int_t i=minslice; i<=maxslice; i++)
 470     {
 471       ReadData(i);
 472       Transform();
 473       if(fAddHistograms) {
 474         if(fVersion != 4)
 475           AddAllHistograms();
 476         else
 477           AddAllHistogramsRows();
 478       }
 479       FindTrackCandidates();
 480       //Evaluate();
 481       //fGlobalMerger->FillTracks(fTracks[0],i);
 482     }
 483 }
 484
 485 void AliL3Hough::ReadData(Int_t slice,Int_t eventnr)
 486 {
 487   //Read data from files, binary or root.
 488
 489 #ifdef use_aliroot
 490   if(fEvent!=eventnr) //just be sure that index is empty for new event
 491     AliL3FileHandler::CleanStaticIndex();
 492 #endif
 493   fCurrentSlice = slice;
 494
 495   for(Int_t i=0; i<fNPatches; i++)
 496     {
 497       fMemHandler[i]->Free();
 498       UInt_t ndigits=0;
 499       AliL3DigitRowData *digits =0;
 500       Char_t name[256];
 501       fMemHandler[i]->Init(slice,i);
 502       if(fBinary)//take input data from binary files
 503         {
 504           if(fUse8bits)
 505             sprintf(name,"%s/binaries/digits_c8_%d_%d_%d.raw",fPath,eventnr,slice,i);
 506           else
 507             sprintf(name,"%s/binaries/digits_%d_%d_%d.raw",fPath,eventnr,slice,i);
 508
 509           fMemHandler[i]->SetBinaryInput(name);
 510           digits = (AliL3DigitRowData *)fMemHandler[i]->CompBinary2Memory(ndigits);
 511           fMemHandler[i]->CloseBinaryInput();
 512         }
 513       else //read data from root file
 514         {
 515 #ifdef use_aliroot
 516           if(fEvent!=eventnr)
 517             fMemHandler[i]->FreeDigitsTree();//or else the new event is not loaded
 518           digits=(AliL3DigitRowData *)fMemHandler[i]->AliAltroDigits2Memory(ndigits,eventnr);
 519 #else
 520           cerr<<"You cannot read from rootfile now"<<endl;
 521 #endif
 522         }
 523
 524       //Set the pointer to the TPCRawStream in case of fast raw data reading
 525       fHoughTransformer[i]->SetTPCRawStream(fMemHandler[i]->GetTPCRawStream());
 526
 527       //set input data and init transformer
 528       fHoughTransformer[i]->SetInputData(ndigits,digits);
 529       fHoughTransformer[i]->Init(slice,i,fNEtaSegments);
 530     }
 531
 532   fEvent=eventnr;
 533 }
 534
 535 void AliL3Hough::Transform(Int_t *rowrange)
 536 {
 537   //Transform all data given to the transformer within the given slice
 538   //(after ReadData(slice))
 539
 540   Double_t initTime,cpuTime;
 541   initTime = GetCpuTime();
 542   Int_t patchorder[6] = {5,2,0,1,3,4}; //The order in which patches are processed
 543   //  Int_t patchorder[6] = {0,1,2,3,4,5}; //The order in which patches are processed
 544   //  Int_t patchorder[6] = {5,4,3,2,1,0}; //The order in which patches are processed
 545   //  Int_t patchorder[6] = {5,2,4,3,1,0}; //The order in which patches are processed
 546   fLastPatch=-1;
 547   for(Int_t i=0; i<fNPatches; i++)
 548     {
 549       // In case of Row transformer reset the arrays only once
 550       if((fVersion != 4) || (i == 0)) {
 551         fBenchmark->Start("Hough Reset");
 552         fHoughTransformer[0]->Reset();//Reset the histograms
 553         fBenchmark->Stop("Hough Reset");
 554       }
 555       fBenchmark->Start("Hough Transform");
 556       PrepareForNextPatch(patchorder[i]);
 557       if(!rowrange) {
 558         char buf[256];
 559         sprintf(buf,"Patch %d",patchorder[i]);
 560         fBenchmark->Start(buf);
 561         fHoughTransformer[patchorder[i]]->SetLastPatch(fLastPatch);
 562         fHoughTransformer[patchorder[i]]->TransformCircle();
 563         fBenchmark->Stop(buf);
 564       }
 565       else
 566         fHoughTransformer[i]->TransformCircleC(rowrange,1);
 567       fBenchmark->Stop("Hough Transform");
 568       fLastPatch=patchorder[i];
 569     }
 570   cpuTime = GetCpuTime() - initTime;
 571   LOG(AliL3Log::kInformational,"AliL3Hough::Transform()","Timing")
 572     <<"Transform done in average per patch of "<<cpuTime*1000/fNPatches<<" ms"<<ENDLOG;
 573 }
 574
 575 void AliL3Hough::MergePatches()
 576 {
 577   // Merge patches if they are not summed
 578   if(fAddHistograms) //Nothing to merge here
 579     return;
 580   fMerger->MergePatches(kTRUE);
 581 }
 582
 583 void AliL3Hough::MergeInternally()
 584 {
 585   // Merge patches internally
 586   if(fAddHistograms)
 587     fInterMerger->FillTracks(fTracks[0]);
 588   else
 589     fInterMerger->FillTracks(fMerger->GetOutTracks());
 590
 591   fInterMerger->MMerge();
 592 }
 593
 594 void AliL3Hough::ProcessSliceIter()
 595 {
 596   //Process current slice (after ReadData(slice)) iteratively.
 597
 598   if(!fAddHistograms)
 599     {
 600       for(Int_t i=0; i<fNPatches; i++)
 601         {
 602           ProcessPatchIter(i);
 603           fMerger->FillTracks(fTracks[i],i); //Copy tracks to merger
 604         }
 605     }
 606   else
 607     {
 608       for(Int_t i=0; i<10; i++)
 609         {
 610           Transform();
 611           AddAllHistograms();
 612           InitEvaluate();
 613           AliL3HoughBaseTransformer *tr = fHoughTransformer[0];
 614           for(Int_t j=0; j<fNEtaSegments; j++)
 615             {
 616               AliL3Histogram *hist = tr->GetHistogram(j);
 617               if(hist->GetNEntries()==0) continue;
 618               fPeakFinder->Reset();
 619               fPeakFinder->SetHistogram(hist);
 620               fPeakFinder->FindAbsMaxima();
 621               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[0]->NextTrack();
 622               track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0));
 623               track->SetEtaIndex(j);
 624               track->SetEta(tr->GetEta(j,fCurrentSlice));
 625               for(Int_t k=0; k<fNPatches; k++)
 626                 {
 627                   fEval[i]->SetNumOfPadsToLook(2);
 628                   fEval[i]->SetNumOfRowsToMiss(2);
 629                   fEval[i]->RemoveFoundTracks();
 630                   /*
 631                   Int_t nrows=0;
 632                   if(!fEval[i]->LookInsideRoad(track,nrows))
 633                     {
 634                       fTracks[0]->Remove(fTracks[0]->GetNTracks()-1);
 635                       fTracks[0]->Compress();
 636                     }
 637                   */
 638                 }
 639             }
 640
 641         }
 642
 643     }
 644 }
 645
 646 void AliL3Hough::ProcessPatchIter(Int_t patch)
 647 {
 648   //Process patch in a iterative way.
 649   //transform + peakfinding + evaluation + transform +...
 650
 651   Int_t numoftries = 5;
 652   AliL3HoughBaseTransformer *tr = fHoughTransformer[patch];
 653   AliL3TrackArray *tracks = fTracks[patch];
 654   tracks->Reset();
 655   AliL3HoughEval *ev = fEval[patch];
 656   ev->InitTransformer(tr);
 657   //ev->RemoveFoundTracks();
 658   ev->SetNumOfRowsToMiss(3);
 659   ev->SetNumOfPadsToLook(2);
 660   AliL3Histogram *hist;
 661   for(Int_t t=0; t<numoftries; t++)
 662     {
 663       tr->Reset();
 664       tr->TransformCircle();
 665       for(Int_t i=0; i<fNEtaSegments; i++)
 666         {
 667           hist = tr->GetHistogram(i);
 668           if(hist->GetNEntries()==0) continue;
 669           fPeakFinder->Reset();
 670           fPeakFinder->SetHistogram(hist);
 671           fPeakFinder->FindAbsMaxima();
 672           //fPeakFinder->FindPeak1();
 673           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
 674           track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0));
 675           track->SetEtaIndex(i);
 676           track->SetEta(tr->GetEta(i,fCurrentSlice));
 677           /*
 678           Int_t nrows=0;
 679           if(!ev->LookInsideRoad(track,nrows))
 680             {
 681               tracks->Remove(tracks->GetNTracks()-1);
 682               tracks->Compress();
 683             }
 684           */
 685         }
 686     }
 687   fTracks[0]->QSort();
 688   LOG(AliL3Log::kInformational,"AliL3Hough::ProcessPatch","NTracks")
 689     <<AliL3Log::kDec<<"Found "<<tracks->GetNTracks()<<" tracks in patch "<<patch<<ENDLOG;
 690 }
 691
 692 void AliL3Hough::AddAllHistograms()
 693 {
 694   //Add the histograms within one etaslice.
 695   //Resulting histogram are in patch=0.
 696
 697   Double_t initTime,cpuTime;
 698   initTime = GetCpuTime();
 699   fBenchmark->Start("Add Histograms");
 700   for(Int_t i=0; i<fNEtaSegments; i++)
 701     {
 702       AliL3Histogram *hist0 = fHoughTransformer[0]->GetHistogram(i);
 703       for(Int_t j=1; j<fNPatches; j++)
 704         {
 705           AliL3Histogram *hist = fHoughTransformer[j]->GetHistogram(i);
 706           hist0->Add(hist);
 707         }
 708     }
 709   fBenchmark->Stop("Add Histograms");
 710   fAddHistograms = kTRUE;
 711   cpuTime = GetCpuTime() - initTime;
 712   LOG(AliL3Log::kInformational,"AliL3Hough::AddAllHistograms()","Timing")
 713     <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG;
 714 }
 715
 716 void AliL3Hough::AddAllHistogramsRows()
 717 {
 718   //Add the histograms within one etaslice.
 719   //Resulting histogram are in patch=0.
 720
 721   Double_t initTime,cpuTime;
 722   initTime = GetCpuTime();
 723   fBenchmark->Start("Add HistogramsRows");
 724
 725   UChar_t lastpatchlastrow = AliL3Transform::GetLastRowOnDDL(fLastPatch)+1;
 726
 727   UChar_t *tracklastrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetTrackLastRow();
 728
 729   for(Int_t i=0; i<fNEtaSegments; i++)
 730     {
 731       UChar_t *gapcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetGapCount(i);
 732       UChar_t *currentrowcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetCurrentRowCount(i);
 733
 734       AliL3Histogram *hist = fHoughTransformer[0]->GetHistogram(i);
 735       Int_t xmin = hist->GetFirstXbin();
 736       Int_t xmax = hist->GetLastXbin();
 737       Int_t ymin = hist->GetFirstYbin();
 738       Int_t ymax = hist->GetLastYbin();
 739       Int_t nxbins = hist->GetNbinsX()+2;
 740
 741       for(Int_t ybin=ymin; ybin<=ymax; ybin++)
 742         {
 743           for(Int_t xbin=xmin; xbin<=xmax; xbin++)
 744             {
 745               Int_t bin = xbin + ybin*nxbins; //Int_t bin = hist->GetBin(xbin,ybin);
 746               if(gapcount[bin] < MAX_N_GAPS) {
 747                 if(tracklastrow[bin] > lastpatchlastrow) {
 748                   if(lastpatchlastrow > currentrowcount[bin])
 749                     gapcount[bin] += (lastpatchlastrow-currentrowcount[bin]-1);
 750                 }
 751                 else {
 752                   if(tracklastrow[bin] > currentrowcount[bin])
 753                     gapcount[bin] += (tracklastrow[bin]-currentrowcount[bin]-1);
 754                 }
 755                 if(gapcount[bin] < MAX_N_GAPS)
 756                   hist->AddBinContent(bin,(159-gapcount[bin]));
 757               }
 758             }
 759         }
 760     }
 761
 762   fBenchmark->Stop("Add HistogramsRows");
 763   fAddHistograms = kTRUE;
 764   cpuTime = GetCpuTime() - initTime;
 765   LOG(AliL3Log::kInformational,"AliL3Hough::AddAllHistogramsRows()","Timing")
 766     <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG;
 767 }
 768
 769 void AliL3Hough::PrepareForNextPatch(Int_t nextpatch)
 770 {
 771   char buf[256];
 772   sprintf(buf,"Prepare For Patch %d",nextpatch);
 773   fBenchmark->Start(buf);
 774
 775   UChar_t lastpatchlastrow;
 776   if(fLastPatch == -1)
 777     lastpatchlastrow = 0;
 778   else
 779     lastpatchlastrow = AliL3Transform::GetLastRowOnDDL(fLastPatch)+1;
 780   UChar_t nextpatchfirstrow;
 781   if(nextpatch==0)
 782     nextpatchfirstrow = 0;
 783   else
 784     nextpatchfirstrow = AliL3Transform::GetFirstRowOnDDL(nextpatch)-1;
 785
 786   UChar_t *trackfirstrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetTrackFirstRow();
 787   UChar_t *tracklastrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetTrackLastRow();
 788
 789   for(Int_t i=0; i<fNEtaSegments; i++)
 790     {
 791       UChar_t *gapcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetGapCount(i);
 792       UChar_t *currentrowcount = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetCurrentRowCount(i);
 793       UChar_t *prevbin = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetPrevBin(i);
 794       UChar_t *nextbin = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetNextBin(i);
 795       UChar_t *nextrow = ((AliL3HoughTransformerRow *)fHoughTransformer[0])->GetNextRow(i);
 796
 797       AliL3Histogram *hist = fHoughTransformer[0]->GetHistogram(i);
 798       Int_t xmin = hist->GetFirstXbin();
 799       Int_t xmax = hist->GetLastXbin();
 800       Int_t ymin = hist->GetFirstYbin();
 801       Int_t ymax = hist->GetLastYbin();
 802       Int_t nxbins = hist->GetNbinsX()+2;
 803
 804       UChar_t lastyvalue = 0;
 805       Int_t endybin = ymin - 1;
 806       for(Int_t ybin=ymin; ybin<=ymax; ybin++)
 807         {
 808           UChar_t lastxvalue = 0;
 809           UChar_t maxvalue = 0;
 810           Int_t endxbin = xmin - 1;
 811           for(Int_t xbin=xmin; xbin<=xmax; xbin++)
 812             {
 813               Int_t bin = xbin + ybin*nxbins;
 814               UChar_t value = 0;
 815               if(gapcount[bin] < MAX_N_GAPS) {
 816                 value = 1;
 817                 maxvalue = 1;
 818                 if(tracklastrow[bin] > lastpatchlastrow) {
 819                   if(lastpatchlastrow > currentrowcount[bin])
 820                     gapcount[bin] += (lastpatchlastrow-currentrowcount[bin]-1);
 821                 }
 822                 else {
 823                   if(tracklastrow[bin] > currentrowcount[bin])
 824                     gapcount[bin] += (tracklastrow[bin]-currentrowcount[bin]-1);
 825                 }
 826                 if(trackfirstrow[bin] < nextpatchfirstrow)
 827                   currentrowcount[bin] = nextpatchfirstrow;
 828                 else
 829                   currentrowcount[bin] = trackfirstrow[bin];
 830               }
 831               if(fLastPatch != -1) {
 832                 if(value > 0)
 833                   {
 834                     nextbin[xbin + ybin*nxbins] = (UChar_t)xbin;
 835                     prevbin[xbin + ybin*nxbins] = (UChar_t)xbin;
 836                     if(value > lastxvalue)
 837                       {
 838                         UChar_t *tempnextbin = nextbin + endxbin + 1 + ybin*nxbins;
 839                         memset(tempnextbin,(UChar_t)xbin,xbin-endxbin-1);
 840                       }
 841                     endxbin = xbin;
 842                   }
 843                 else
 844                   {
 845                     prevbin[xbin + ybin*nxbins] = (UChar_t)endxbin;
 846                   }
 847                 lastxvalue = value;
 848               }
 849             }
 850           if(fLastPatch != -1) {
 851             UChar_t *tempnextbin = nextbin + endxbin + 1 + ybin*nxbins;
 852             memset(tempnextbin,(UChar_t)(xmax+1),xmax-endxbin);
 853             if(maxvalue > 0)
 854               {
 855                 nextrow[ybin] = (UChar_t)ybin;
 856                 if(maxvalue > lastyvalue)
 857                   {
 858                     UChar_t *tempnextrow = nextrow + endybin + 1;
 859                     memset(tempnextrow,(UChar_t)ybin,ybin-endybin-1);
 860                   }
 861                 endybin = ybin;
 862               }
 863             lastyvalue = maxvalue;
 864           }
 865         }
 866       if(fLastPatch != -1) {
 867         UChar_t *tempnextrow = nextrow + endybin + 1;
 868         memset(tempnextrow,(UChar_t)(ymax+1),ymax-endybin);
 869       }
 870     }
 871
 872   fBenchmark->Stop(buf);
 873 }
 874
 875 void AliL3Hough::AddTracks()
 876 {
 877   // Add current slice slice tracks to the global list of found tracks
 878   if(!fTracks[0])
 879     {
 880       cerr<<"AliL3Hough::AddTracks : No tracks"<<endl;
 881       return;
 882     }
 883   AliL3TrackArray *tracks = fTracks[0];
 884   for(Int_t i=0; i<tracks->GetNTracks(); i++)
 885     {
 886       AliL3Track *track = tracks->GetCheckedTrack(i);
 887       if(!track) continue;
 888       if(track->GetNHits()!=1) cerr<<"NHITS "<<track->GetNHits()<<endl;
 889       UInt_t *ids = track->GetHitNumbers();
 890       ids[0] = (fCurrentSlice&0x7f)<<25;
 891     }
 892
 893   fGlobalTracks->AddTracks(fTracks[0],0,fCurrentSlice);
 894 }
 895
 896 void AliL3Hough::FindTrackCandidatesRow()
 897 {
 898   // Find AliL3HoughTransformerRow track candidates
 899   if(fVersion != 4) {
 900     LOG(AliL3Log::kError,"AliL3Hough::FindTrackCandidatesRow()","")
 901       <<"Incompatible Peak Finder version!"<<ENDLOG;
 902     return;
 903   }
 904
 905   //Look for peaks in histograms, and find the track candidates
 906   Int_t npatches;
 907   if(fAddHistograms)
 908     npatches = 1; //Histograms have been added.
 909   else
 910     npatches = fNPatches;
 911
 912   Double_t initTime,cpuTime;
 913   initTime = GetCpuTime();
 914   fBenchmark->Start("Find Maxima");
 915   for(Int_t i=0; i<npatches; i++)
 916     {
 917       AliL3HoughBaseTransformer *tr = fHoughTransformer[i];
 918       fTracks[i]->Reset();
 919       fPeakFinder->Reset();
 920
 921       for(Int_t j=0; j<fNEtaSegments; j++)
 922         {
 923           AliL3Histogram *hist = tr->GetHistogram(j);
 924           if(hist->GetNEntries()==0) continue;
 925           fPeakFinder->SetHistogram(hist);
 926           fPeakFinder->SetEtaSlice(j);
 927           fPeakFinder->SetTrackLUTs(((AliL3HoughTransformerRow *)tr)->fTrackNRows,((AliL3HoughTransformerRow *)tr)->fTrackFirstRow,((AliL3HoughTransformerRow *)tr)->fTrackLastRow);
 928 #ifdef do_mc
 929           LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","")
 930             <<"Starting "<<j<<" etaslice"<<ENDLOG;
 931 #endif
 932           fPeakFinder->SetThreshold(fPeakThreshold[i]);
 933           fPeakFinder->FindAdaptedRowPeaks(1,0,0);//Maxima finder for HoughTransformerRow
 934
 935           //fPeakFinder->FindMaxima(fPeakThreshold[i]); //Simple maxima finder
 936         }
 937
 938       for(Int_t k=0; k<fPeakFinder->GetEntries(); k++)
 939         {
 940           if(fPeakFinder->GetWeight(k) < 0) continue;
 941           AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack();
 942           Float_t psi = atan((fPeakFinder->GetXPeak(k)-fPeakFinder->GetYPeak(k))/(AliL3HoughTransformerRow::GetBeta1()-AliL3HoughTransformerRow::GetBeta2()));
 943           Float_t kappa = 2.0*(fPeakFinder->GetXPeak(k)*cos(psi)-AliL3HoughTransformerRow::GetBeta1()*sin(psi));
 944           //          track->SetTrackParameters(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetWeight(k));
 945           track->SetTrackParameters(kappa,psi,fPeakFinder->GetWeight(k));
 946           track->SetBinXY(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetXPeakSize(k),fPeakFinder->GetYPeakSize(k));
 947           Int_t etaindex = (fPeakFinder->GetStartEta(k)+fPeakFinder->GetEndEta(k))/2;
 948           track->SetEtaIndex(etaindex);
 949           Float_t starteta = tr->GetEta(fPeakFinder->GetStartEta(k),fCurrentSlice);
 950           Float_t endeta = tr->GetEta(fPeakFinder->GetEndEta(k),fCurrentSlice);
 951           track->SetEta((starteta+endeta)/2.0);
 952           track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));
 953           track->SetSector(fCurrentSlice);
 954           track->SetSlice(fCurrentSlice);
 955 #ifdef do_mc
 956           Int_t label = tr->GetTrackID(etaindex,fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k));
 957           track->SetMCid(label);
 958           //      cout<<"Track found with label "<<label<<" at "<<fPeakFinder->GetXPeak(k)<<" "<<fPeakFinder->GetYPeak(k)<<" with weight "<<fPeakFinder->GetWeight(k)<<endl;
 959 #endif
 960         }
 961       LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","")
 962         <<"Found "<<fTracks[i]->GetNTracks()<<" tracks in slice "<<fCurrentSlice<<ENDLOG;
 963       fTracks[i]->QSort();
 964     }
 965   fBenchmark->Stop("Find Maxima");
 966   cpuTime = GetCpuTime() - initTime;
 967   LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","Timing")
 968     <<"Maxima finding done in "<<cpuTime*1000<<" ms"<<ENDLOG;
 969 }
 970
 971 void AliL3Hough::FindTrackCandidates()
 972 {
 973   // Find AliL3HoughTransformer track candidates
 974   if(fVersion == 4) {
 975     LOG(AliL3Log::kError,"AliL3Hough::FindTrackCandidatesRow()","")
 976       <<"Incompatible Peak Finder version!"<<ENDLOG;
 977     return;
 978   }
 979
 980   Int_t npatches;
 981   if(fAddHistograms)
 982     npatches = 1; //Histograms have been added.
 983   else
 984     npatches = fNPatches;
 985
 986   Double_t initTime,cpuTime;
 987   initTime = GetCpuTime();
 988   fBenchmark->Start("Find Maxima");
 989   for(Int_t i=0; i<npatches; i++)
 990     {
 991       AliL3HoughBaseTransformer *tr = fHoughTransformer[i];
 992       fTracks[i]->Reset();
 993
 994       for(Int_t j=0; j<fNEtaSegments; j++)
 995         {
 996           AliL3Histogram *hist = tr->GetHistogram(j);
 997           if(hist->GetNEntries()==0) continue;
 998           fPeakFinder->Reset();
 999           fPeakFinder->SetHistogram(hist);
1000 #ifdef do_mc
1001           cout<<"Starting "<<j<<" etaslice"<<endl;
1002 #endif
1003           fPeakFinder->SetThreshold(fPeakThreshold[i]);
1004           fPeakFinder->FindAdaptedPeaks(fKappaSpread,fPeakRatio);
1005
1006           for(Int_t k=0; k<fPeakFinder->GetEntries(); k++)
1007             {
1008               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack();
1009               track->SetTrackParameters(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetWeight(k));
1010               track->SetEtaIndex(j);
1011               track->SetEta(tr->GetEta(j,fCurrentSlice));
1012               track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));
1013             }
1014         }
1015       cout<<"Found "<<fTracks[i]->GetNTracks()<<" tracks in patch "<<i<<endl;
1016       fTracks[i]->QSort();
1017     }
1018   fBenchmark->Stop("Find Maxima");
1019   cpuTime = GetCpuTime() - initTime;
1020   LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","Timing")
1021     <<"Maxima finding done in "<<cpuTime*1000<<" ms"<<ENDLOG;
1022 }
1023
1024 void AliL3Hough::InitEvaluate()
1025 {
1026   //Pass the transformer objects to the AliL3HoughEval objects:
1027   //This will provide the evaluation objects with all the necessary
1028   //data and parameters it needs.
1029
1030   for(Int_t i=0; i<fNPatches; i++)
1031     fEval[i]->InitTransformer(fHoughTransformer[i]);
1032 }
1033
1034 Int_t AliL3Hough::Evaluate(Int_t roadwidth,Int_t nrowstomiss)
1035 {
1036   //Evaluate the tracks, by looking along the road in the raw data.
1037   //If track does not cross all padrows - rows2miss, it is removed from the arrray.
1038   //If histograms were not added, the check is done locally in patch,
1039   //meaning that nrowstomiss is the number of padrows the road can miss with respect
1040   //to the number of rows in the patch.
1041   //If the histograms were added, the comparison is done globally in the _slice_,
1042   //meaing that nrowstomiss is the number of padrows the road can miss with
1043   //respect to the total number of padrows in the slice.
1044   //
1045   //Return value = number of tracks which were removed (only in case of fAddHistograms)
1046
1047   if(!fTracks[0])
1048     {
1049       LOG(AliL3Log::kError,"AliL3Hough::Evaluate","Track Array")
1050         <<"No tracks to work with..."<<ENDLOG;
1051       return 0;
1052     }
1053
1054   Int_t removedtracks=0;
1055   AliL3TrackArray *tracks=0;
1056
1057   if(fAddHistograms)
1058     {
1059       tracks = fTracks[0];
1060       for(Int_t i=0; i<tracks->GetNTracks(); i++)
1061         {
1062           AliL3Track *track = tracks->GetCheckedTrack(i);
1063           if(!track) continue;
1064           track->SetNHits(0);
1065         }
1066     }
1067
1068   for(Int_t i=0; i<fNPatches; i++)
1069     EvaluatePatch(i,roadwidth,nrowstomiss);
1070
1071   //Here we check the tracks globally;
1072   //how many good rows (padrows with signal)
1073   //did it cross in the slice
1074   if(fAddHistograms)
1075     {
1076       for(Int_t j=0; j<tracks->GetNTracks(); j++)
1077         {
1078           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
1079
1080           if(track->GetNHits() < AliL3Transform::GetNRows() - nrowstomiss)
1081             {
1082               tracks->Remove(j);
1083               removedtracks++;
1084             }
1085         }
1086       tracks->Compress();
1087       tracks->QSort();
1088     }
1089
1090   return removedtracks;
1091 }
1092
1093 void AliL3Hough::EvaluatePatch(Int_t i,Int_t roadwidth,Int_t nrowstomiss)
1094 {
1095   //Evaluate patch i.
1096
1097   fEval[i]->InitTransformer(fHoughTransformer[i]);
1098   fEval[i]->SetNumOfPadsToLook(roadwidth);
1099   fEval[i]->SetNumOfRowsToMiss(nrowstomiss);
1100   //fEval[i]->RemoveFoundTracks();
1101
1102   AliL3TrackArray *tracks=0;
1103
1104   if(!fAddHistograms)
1105     tracks = fTracks[i];
1106   else
1107     tracks = fTracks[0];
1108
1109   Int_t nrows=0;
1110   for(Int_t j=0; j<tracks->GetNTracks(); j++)
1111     {
1112       AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
1113       if(!track)
1114         {
1115           LOG(AliL3Log::kWarning,"AliL3Hough::EvaluatePatch","Track array")
1116             <<"Track object missing!"<<ENDLOG;
1117           continue;
1118         }
1119       nrows=0;
1120       Int_t rowrange[2] = {AliL3Transform::GetFirstRow(i),AliL3Transform::GetLastRow(i)};
1121       Bool_t result = fEval[i]->LookInsideRoad(track,nrows,rowrange);
1122       if(fAddHistograms)
1123         {
1124           Int_t pre=track->GetNHits();
1125           track->SetNHits(pre+nrows);
1126         }
1127       else//the track crossed too few good padrows (padrows with signal) in the patch, so remove it
1128         {
1129           if(result == kFALSE)
1130             tracks->Remove(j);
1131         }
1132     }
1133
1134   tracks->Compress();
1135
1136 }
1137
1138 void AliL3Hough::MergeEtaSlices()
1139 {
1140   //Merge tracks found in neighbouring eta slices.
1141   //Removes the track with the lower weight.
1142
1143   fBenchmark->Start("Merge Eta-slices");
1144   AliL3TrackArray *tracks = fTracks[0];
1145   if(!tracks)
1146     {
1147       cerr<<"AliL3Hough::MergeEtaSlices : No tracks "<<endl;
1148       return;
1149     }
1150   for(Int_t j=0; j<tracks->GetNTracks(); j++)
1151     {
1152       AliL3HoughTrack *track1 = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
1153       if(!track1) continue;
1154       for(Int_t k=j+1; k<tracks->GetNTracks(); k++)
1155         {
1156           AliL3HoughTrack *track2 = (AliL3HoughTrack*)tracks->GetCheckedTrack(k);
1157           if(!track2) continue;
1158           if(abs(track1->GetEtaIndex() - track2->GetEtaIndex()) != 1) continue;
1159           if(fabs(track1->GetKappa()-track2->GetKappa()) < 0.006 &&
1160              fabs(track1->GetPsi()- track2->GetPsi()) < 0.1)
1161             {
1162               //cout<<"Merging track in slices "<<track1->GetEtaIndex()<<" "<<track2->GetEtaIndex()<<endl;
1163               if(track1->GetWeight() > track2->GetWeight())
1164                 tracks->Remove(k);
1165               else
1166                 tracks->Remove(j);
1167             }
1168         }
1169     }
1170   fBenchmark->Stop("Merge Eta-slices");
1171   tracks->Compress();
1172 }
1173
1174 void AliL3Hough::WriteTracks(Char_t *path)
1175 {
1176   // Write found tracks into file
1177   //cout<<"AliL3Hough::WriteTracks : Sorting the tracsk"<<endl;
1178   //fGlobalTracks->QSort();
1179
1180   Char_t filename[1024];
1181   sprintf(filename,"%s/tracks_%d.raw",path,fEvent);
1182   AliL3MemHandler mem;
1183   mem.SetBinaryOutput(filename);
1184   mem.TrackArray2Binary(fGlobalTracks);
1185   mem.CloseBinaryOutput();
1186   fGlobalTracks->Reset();
1187 }
1188
1189 void AliL3Hough::WriteTracks(Int_t slice,Char_t *path)
1190 {
1191   // Write found tracks slice by slice into file
1192
1193   AliL3MemHandler mem;
1194   Char_t fname[100];
1195   if(fAddHistograms)
1196     {
1197       sprintf(fname,"%s/tracks_ho_%d_%d.raw",path,fEvent,slice);
1198       mem.SetBinaryOutput(fname);
1199       mem.TrackArray2Binary(fTracks[0]);
1200       mem.CloseBinaryOutput();
1201     }
1202   else
1203     {
1204       for(Int_t i=0; i<fNPatches; i++)
1205         {
1206           sprintf(fname,"%s/tracks_ho_%d_%d_%d.raw",path,fEvent,slice,i);
1207           mem.SetBinaryOutput(fname);
1208           mem.TrackArray2Binary(fTracks[i]);
1209           mem.CloseBinaryOutput();
1210         }
1211     }
1212 }
1213
1214 void AliL3Hough::WriteDigits(Char_t *outfile)
1215 {
1216   //Write the current data to a new rootfile.
1217 #ifdef use_aliroot
1218
1219   for(Int_t i=0; i<fNPatches; i++)
1220     {
1221       AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer[i]->GetDataPointer();
1222       fMemHandler[i]->AliDigits2RootFile(tempPt,outfile);
1223     }
1224 #else
1225   cerr<<"AliL3Hough::WriteDigits : You need to compile with AliROOT!"<<endl;
1226   return;
1227 #endif
1228 }
1229
1230 Double_t AliL3Hough::GetCpuTime()
1231 {
1232   //Return the Cputime in seconds.
1233  struct timeval tv;
1234  gettimeofday( &tv, NULL );
1235  return tv.tv_sec+(((Double_t)tv.tv_usec)/1000000.);
1236 }
1237
1238 void *AliL3Hough::ProcessInThread(void *args)
1239 {
1240   AliL3Hough *instance = (AliL3Hough *)args;
1241   Int_t minslice = instance->GetMinSlice();
1242   Int_t maxslice = instance->GetMaxSlice();
1243   for(Int_t i=minslice; i<=maxslice; i++)
1244     {
1245       instance->ReadData(i,0);
1246       instance->Transform();
1247       instance->AddAllHistogramsRows();
1248       instance->FindTrackCandidatesRow();
1249       instance->AddTracks();
1250     }
1251   return (void *)0;
1252 }
1253
1254 void AliL3Hough::StartProcessInThread(Int_t minslice,Int_t maxslice)
1255 {
1256   if(!fThread) {
1257     char buf[255];
1258     sprintf(buf,"houghtrans_%d_%d",minslice,maxslice);
1259     SetMinMaxSlices(minslice,maxslice);
1260     //    fThread = new TThread(buf,(void (*) (void *))&ProcessInThread,(void *)this);
1261     fThread = new TThread(buf,&ProcessInThread,(void *)this);
1262     fThread->Run();
1263   }
1264   return;
1265 }
1266
1267 Int_t AliL3Hough::WaitForThreadFinish()
1268 {
1269   return TThread::Join(fThread->GetId());
1270 }