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 "AliL3HoughMaxFinder.h"
  20 #include "AliL3Benchmark.h"
  21 #ifdef use_aliroot
  22 #include "AliL3FileHandler.h"
  23 #else
  24 #include "AliL3MemHandler.h"
  25 #endif
  26 #include "AliL3DataHandler.h"
  27 #include "AliL3DigitData.h"
  28 #include "AliL3HoughEval.h"
  29 #include "AliL3Transform.h"
  30 #include "AliL3TrackArray.h"
  31 #include "AliL3HoughTrack.h"
  32 #include "AliL3DDLDataFileHandler.h"
  33
  34 #if GCCVERSION == 3
  35 using namespace std;
  36 #endif
  37
  38 /** /class AliL3Hough
  39 //<pre>
  40 //_____________________________________________________________
  41 // AliL3Hough
  42 //
  43 // Interface class for the Hough transform
  44 //
  45 // Example how to use:
  46 //
  47 // AliL3Hough *hough = new AliL3Hough(path,kTRUE,NumberOfEtaSegments);
  48 // hough->ReadData(slice);
  49 // hough->Transform();
  50 // hough->FindTrackCandidates();
  51 //
  52 // AliL3TrackArray *tracks = hough->GetTracks(patch);
  53 //</pre>
  54 */
  55
  56 ClassImp(AliL3Hough)
  57
  58 AliL3Hough::AliL3Hough()
  59 {
  60   //Constructor
  61
  62   fBinary        = kFALSE;
  63   fAddHistograms = kFALSE;
  64   fDoIterative   = kFALSE;
  65   fWriteDigits   = kFALSE;
  66   fUse8bits      = kFALSE;
  67
  68   fMemHandler       = 0;
  69   fHoughTransformer = 0;
  70   fEval             = 0;
  71   fPeakFinder       = 0;
  72   fTracks           = 0;
  73   fGlobalTracks     = 0;
  74   fMerger           = 0;
  75   fInterMerger      = 0;
  76   fGlobalMerger     = 0;
  77   fBenchmark        = 0;
  78
  79   fNEtaSegments     = 0;
  80   fNPatches         = 0;
  81   fVersion          = 0;
  82   fCurrentSlice     = 0;
  83   fEvent            = 0;
  84
  85   SetTransformerParams();
  86   SetThreshold();
  87   SetNSaveIterations();
  88   SetPeakThreshold();
  89 }
  90
  91 AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile)
  92 {
  93   //Default ctor.
  94
  95   fBinary = binary;
  96   strcpy(fPath,path);
  97   fNEtaSegments  = n_eta_segments;
  98   fAddHistograms = kFALSE;
  99   fDoIterative   = kFALSE;
 100   fWriteDigits   = kFALSE;
 101   fUse8bits      = bit8;
 102   fVersion       = tv;
 103   if(!fBinary)
 104     fInputFile = infile;
 105   else
 106     fInputFile = 0;
 107 }
 108
 109 AliL3Hough::~AliL3Hough()
 110 {
 111   //dtor
 112
 113   CleanUp();
 114   if(fMerger)
 115     delete fMerger;
 116   if(fInterMerger)
 117     delete fInterMerger;
 118   if(fPeakFinder)
 119     delete fPeakFinder;
 120   if(fGlobalMerger)
 121     delete fGlobalMerger;
 122   if(fBenchmark)
 123     delete fBenchmark;
 124   if(fGlobalTracks)
 125     delete fGlobalTracks;
 126 }
 127
 128 void AliL3Hough::CleanUp()
 129 {
 130   //Cleanup memory
 131
 132   for(Int_t i=0; i<fNPatches; i++)
 133     {
 134       if(fTracks[i]) delete fTracks[i];
 135       if(fEval[i]) delete fEval[i];
 136       if(fHoughTransformer[i]) delete fHoughTransformer[i];
 137       if(fMemHandler[i]) delete fMemHandler[i];
 138     }
 139
 140   if(fTracks) delete [] fTracks;
 141   if(fEval) delete [] fEval;
 142   if(fHoughTransformer) delete [] fHoughTransformer;
 143   if(fMemHandler) delete [] fMemHandler;
 144
 145 }
 146
 147 void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile)
 148 {
 149   fBinary = binary;
 150   strcpy(fPath,path);
 151   fNEtaSegments = n_eta_segments;
 152   fWriteDigits  = kFALSE;
 153   fUse8bits     = bit8;
 154   fVersion      = tv;
 155   if(!fBinary)
 156     fInputFile = infile;
 157   else
 158     fInputFile = 0;
 159
 160   Init(); //do the rest
 161 }
 162
 163 void AliL3Hough::Init(Bool_t doit, Bool_t addhists)
 164 {
 165   fDoIterative   = doit;
 166   fAddHistograms = addhists;
 167
 168   fNPatches = AliL3Transform::GetNPatches();
 169
 170   fHoughTransformer = new AliL3HoughBaseTransformer*[fNPatches];
 171   fMemHandler = new AliL3MemHandler*[fNPatches];
 172
 173   fTracks = new AliL3TrackArray*[fNPatches];
 174   fEval = new AliL3HoughEval*[fNPatches];
 175
 176   fGlobalTracks = new AliL3TrackArray("AliL3HoughTrack");
 177
 178   for(Int_t i=0; i<fNPatches; i++)
 179     {
 180       switch (fVersion){ //choose Transformer
 181       case 1:
 182         fHoughTransformer[i] = new AliL3HoughTransformerLUT(0,i,fNEtaSegments);
 183         break;
 184       case 2:
 185         fHoughTransformer[i] = new AliL3HoughClusterTransformer(0,i,fNEtaSegments);
 186         break;
 187       case 3:
 188         fHoughTransformer[i] = new AliL3HoughTransformerVhdl(0,i,fNEtaSegments,fNSaveIterations);
 189         break;
 190       default:
 191         fHoughTransformer[i] = new AliL3HoughTransformer(0,i,fNEtaSegments,kFALSE,kFALSE);
 192       }
 193
 194       fHoughTransformer[i]->CreateHistograms(fNBinX[i],fLowPt[i],fNBinY[i],-fPhi[i],fPhi[i]);
 195       //fHoughTransformer[i]->CreateHistograms(fLowPt[i],fUpperPt[i],fPtRes[i],fNBinY[i],fPhi[i]);
 196       fHoughTransformer[i]->SetLowerThreshold(fThreshold[i]);
 197       fHoughTransformer[i]->SetUpperThreshold(100);
 198
 199       LOG(AliL3Log::kInformational,"AliL3Hough::Init","Version")
 200         <<"Initializing Hough transformer version "<<fVersion<<ENDLOG;
 201
 202       fEval[i] = new AliL3HoughEval();
 203       fTracks[i] = new AliL3TrackArray("AliL3HoughTrack");
 204       if(fUse8bits)
 205         fMemHandler[i] = new AliL3DataHandler();
 206       else
 207 #ifdef use_aliroot
 208         {
 209           if(!fInputFile) {
 210             /* In case of reading digits file */
 211             fMemHandler[i] = new AliL3FileHandler();
 212             if(!fBinary) {
 213               Char_t filename[1024];
 214               sprintf(filename,"%s/digitfile.root",fPath);
 215               fMemHandler[i]->SetAliInput(filename);
 216             }
 217           }
 218           else {
 219             /* In case of reading rawdata from ROOT file */
 220             fMemHandler[i] = new AliL3DDLDataFileHandler();
 221             fMemHandler[i]->SetReaderInput(fInputFile);
 222           }
 223         }
 224 #else
 225       fMemHandler[i] = new AliL3MemHandler();
 226 #endif
 227     }
 228
 229   fPeakFinder = new AliL3HoughMaxFinder("KappaPhi",1000);
 230   fMerger = new AliL3HoughMerger(fNPatches);
 231   fInterMerger = new AliL3HoughIntMerger();
 232   fGlobalMerger = 0;
 233   fBenchmark = new AliL3Benchmark();
 234 }
 235
 236 void AliL3Hough::SetTransformerParams(Float_t ptres,Float_t ptmin,Float_t ptmax,Int_t ny,Int_t patch)
 237 {
 238
 239   Int_t mrow;
 240   Float_t psi=0;
 241   if(patch==-1)
 242     mrow = 80;
 243   else
 244     mrow = AliL3Transform::GetLastRow(patch);
 245   if(ptmin)
 246     {
 247       Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2));
 248       Double_t kappa = -1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/ptmin;
 249       psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2);
 250       cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl;
 251     }
 252
 253   if(patch==-1)
 254     {
 255       Int_t i=0;
 256       while(i < 6)
 257         {
 258           fPtRes[i] = ptres;
 259           fLowPt[i] = ptmin;
 260           fUpperPt[i] = ptmax;
 261           fNBinY[i] = ny;
 262           fPhi[i] = psi;
 263           fNBinX[i]=0;
 264           i++;
 265         }
 266       return;
 267     }
 268
 269   fPtRes[patch] = ptres;
 270   fLowPt[patch] = ptmin;
 271   fUpperPt[patch] = ptmax;
 272   fNBinY[patch] = ny;
 273   fPhi[patch] = psi;
 274 }
 275
 276 void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t patch)
 277 {
 278
 279   Int_t mrow=80;
 280   Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2));
 281   Double_t kappa = -1*AliL3Transform::GetBField()*AliL3Transform::GetBFact()/ptmin;
 282   Double_t psi = AliL3Transform::Deg2Rad(10) - asin(lineradius*kappa/2);
 283   cout<<"Calculated psi range "<<psi<<" in patch "<<patch<<endl;
 284
 285   Int_t i=0;
 286   while(i < 6)
 287     {
 288       fLowPt[i] = ptmin;
 289       fNBinY[i] = ny;
 290       fNBinX[i] = nx;
 291       fPhi[i] = psi;
 292       i++;
 293     }
 294 }
 295
 296 void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t lpt,Float_t phi)
 297 {
 298   Int_t i=0;
 299   while(i < 6)
 300     {
 301       fLowPt[i] = lpt;
 302       fNBinY[i] = ny;
 303       fNBinX[i] = nx;
 304       fPhi[i] = phi;
 305       i++;
 306     }
 307 }
 308
 309 void AliL3Hough::SetThreshold(Int_t t3,Int_t patch)
 310 {
 311   if(patch==-1)
 312     {
 313       Int_t i=0;
 314       while(i < 6)
 315         fThreshold[i++]=t3;
 316       return;
 317     }
 318   fThreshold[patch]=t3;
 319 }
 320
 321 void AliL3Hough::SetPeakThreshold(Int_t threshold,Int_t patch)
 322 {
 323   if(patch==-1)
 324     {
 325       Int_t i=0;
 326       while(i < 6)
 327         fPeakThreshold[i++]=threshold;
 328       return;
 329     }
 330   fPeakThreshold[patch]=threshold;
 331 }
 332
 333 void AliL3Hough::DoBench(Char_t *name)
 334 {
 335   fBenchmark->Analyze(name);
 336 }
 337
 338 void AliL3Hough::Process(Int_t minslice,Int_t maxslice)
 339 {
 340   //Process all slices [minslice,maxslice].
 341   fGlobalMerger = new AliL3HoughGlobalMerger(minslice,maxslice);
 342
 343   for(Int_t i=minslice; i<=maxslice; i++)
 344     {
 345       ReadData(i);
 346       Transform();
 347       if(fAddHistograms)
 348         AddAllHistograms();
 349       FindTrackCandidates();
 350       //Evaluate();
 351       //fGlobalMerger->FillTracks(fTracks[0],i);
 352     }
 353 }
 354
 355 void AliL3Hough::ReadData(Int_t slice,Int_t eventnr)
 356 {
 357   //Read data from files, binary or root.
 358
 359   fEvent=eventnr;
 360   fCurrentSlice = slice;
 361   for(Int_t i=0; i<fNPatches; i++)
 362     {
 363       fMemHandler[i]->Free();
 364       UInt_t ndigits=0;
 365       AliL3DigitRowData *digits =0;
 366       Char_t name[256];
 367       fMemHandler[i]->Init(slice,i);
 368       if(fBinary)//take input data from binary files
 369         {
 370           if(fUse8bits)
 371             sprintf(name,"%s/binaries/digits_c8_%d_%d_%d.raw",fPath,eventnr,slice,i);
 372           else
 373             sprintf(name,"%s/binaries/digits_%d_%d_%d.raw",fPath,eventnr,slice,i);
 374
 375           fMemHandler[i]->SetBinaryInput(name);
 376           digits = (AliL3DigitRowData *)fMemHandler[i]->CompBinary2Memory(ndigits);
 377           fMemHandler[i]->CloseBinaryInput();
 378         }
 379       else //read data from root file
 380         {
 381 #ifdef use_aliroot
 382           digits=(AliL3DigitRowData *)fMemHandler[i]->AliAltroDigits2Memory(ndigits,eventnr);
 383 #else
 384           cerr<<"You cannot read from rootfile now"<<endl;
 385 #endif
 386         }
 387
 388       //set input data and init transformer
 389       fHoughTransformer[i]->SetInputData(ndigits,digits);
 390       fHoughTransformer[i]->Init(slice,i,fNEtaSegments);
 391     }
 392 }
 393
 394 void AliL3Hough::Transform(Int_t *row_range)
 395 {
 396   //Transform all data given to the transformer within the given slice
 397   //(after ReadData(slice))
 398
 399   Double_t initTime,cpuTime;
 400   initTime = GetCpuTime();
 401   for(Int_t i=0; i<fNPatches; i++)
 402     {
 403       fHoughTransformer[i]->Reset();//Reset the histograms
 404       fBenchmark->Start("Hough Transform");
 405       if(!row_range)
 406         fHoughTransformer[i]->TransformCircle();
 407       else
 408         fHoughTransformer[i]->TransformCircleC(row_range,1);
 409       fBenchmark->Stop("Hough Transform");
 410     }
 411   cpuTime = GetCpuTime() - initTime;
 412   LOG(AliL3Log::kInformational,"AliL3Hough::Transform()","Timing")
 413     <<"Transform done in average per patch of "<<cpuTime*1000/fNPatches<<" ms"<<ENDLOG;
 414 }
 415
 416 void AliL3Hough::MergePatches()
 417 {
 418   if(fAddHistograms) //Nothing to merge here
 419     return;
 420   fMerger->MergePatches(kTRUE);
 421 }
 422
 423 void AliL3Hough::MergeInternally()
 424 {
 425   if(fAddHistograms)
 426     fInterMerger->FillTracks(fTracks[0]);
 427   else
 428     fInterMerger->FillTracks(fMerger->GetOutTracks());
 429
 430   fInterMerger->MMerge();
 431 }
 432
 433 void AliL3Hough::ProcessSliceIter()
 434 {
 435   //Process current slice (after ReadData(slice)) iteratively.
 436
 437   if(!fAddHistograms)
 438     {
 439       for(Int_t i=0; i<fNPatches; i++)
 440         {
 441           ProcessPatchIter(i);
 442           fMerger->FillTracks(fTracks[i],i); //Copy tracks to merger
 443         }
 444     }
 445   else
 446     {
 447       for(Int_t i=0; i<10; i++)
 448         {
 449           Transform();
 450           AddAllHistograms();
 451           InitEvaluate();
 452           AliL3HoughBaseTransformer *tr = fHoughTransformer[0];
 453           for(Int_t j=0; j<fNEtaSegments; j++)
 454             {
 455               AliL3Histogram *hist = tr->GetHistogram(j);
 456               if(hist->GetNEntries()==0) continue;
 457               fPeakFinder->Reset();
 458               fPeakFinder->SetHistogram(hist);
 459               fPeakFinder->FindAbsMaxima();
 460               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[0]->NextTrack();
 461               track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0));
 462               track->SetEtaIndex(j);
 463               track->SetEta(tr->GetEta(j,fCurrentSlice));
 464               for(Int_t k=0; k<fNPatches; k++)
 465                 {
 466                   fEval[i]->SetNumOfPadsToLook(2);
 467                   fEval[i]->SetNumOfRowsToMiss(2);
 468                   fEval[i]->RemoveFoundTracks();
 469                   /*
 470                   Int_t nrows=0;
 471                   if(!fEval[i]->LookInsideRoad(track,nrows))
 472                     {
 473                       fTracks[0]->Remove(fTracks[0]->GetNTracks()-1);
 474                       fTracks[0]->Compress();
 475                     }
 476                   */
 477                 }
 478             }
 479
 480         }
 481
 482     }
 483 }
 484
 485 void AliL3Hough::ProcessPatchIter(Int_t patch)
 486 {
 487   //Process patch in a iterative way.
 488   //transform + peakfinding + evaluation + transform +...
 489
 490   Int_t num_of_tries = 5;
 491   AliL3HoughBaseTransformer *tr = fHoughTransformer[patch];
 492   AliL3TrackArray *tracks = fTracks[patch];
 493   tracks->Reset();
 494   AliL3HoughEval *ev = fEval[patch];
 495   ev->InitTransformer(tr);
 496   //ev->RemoveFoundTracks();
 497   ev->SetNumOfRowsToMiss(3);
 498   ev->SetNumOfPadsToLook(2);
 499   AliL3Histogram *hist;
 500   for(Int_t t=0; t<num_of_tries; t++)
 501     {
 502       tr->Reset();
 503       tr->TransformCircle();
 504       for(Int_t i=0; i<fNEtaSegments; i++)
 505         {
 506           hist = tr->GetHistogram(i);
 507           if(hist->GetNEntries()==0) continue;
 508           fPeakFinder->Reset();
 509           fPeakFinder->SetHistogram(hist);
 510           fPeakFinder->FindAbsMaxima();
 511           //fPeakFinder->FindPeak1();
 512           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
 513           track->SetTrackParameters(fPeakFinder->GetXPeak(0),fPeakFinder->GetYPeak(0),fPeakFinder->GetWeight(0));
 514           track->SetEtaIndex(i);
 515           track->SetEta(tr->GetEta(i,fCurrentSlice));
 516           /*
 517           Int_t nrows=0;
 518           if(!ev->LookInsideRoad(track,nrows))
 519             {
 520               tracks->Remove(tracks->GetNTracks()-1);
 521               tracks->Compress();
 522             }
 523           */
 524         }
 525     }
 526   fTracks[0]->QSort();
 527   LOG(AliL3Log::kInformational,"AliL3Hough::ProcessPatch","NTracks")
 528     <<AliL3Log::kDec<<"Found "<<tracks->GetNTracks()<<" tracks in patch "<<patch<<ENDLOG;
 529 }
 530
 531 void AliL3Hough::AddAllHistograms()
 532 {
 533   //Add the histograms within one etaslice.
 534   //Resulting histogram are in patch=0.
 535
 536   Double_t initTime,cpuTime;
 537   initTime = GetCpuTime();
 538   fBenchmark->Start("Add Histograms");
 539   for(Int_t i=0; i<fNEtaSegments; i++)
 540     {
 541       AliL3Histogram *hist0 = fHoughTransformer[0]->GetHistogram(i);
 542       for(Int_t j=1; j<fNPatches; j++)
 543         {
 544           AliL3Histogram *hist = fHoughTransformer[j]->GetHistogram(i);
 545           hist0->Add(hist);
 546         }
 547     }
 548   fBenchmark->Stop("Add Histograms");
 549   fAddHistograms = kTRUE;
 550   cpuTime = GetCpuTime() - initTime;
 551   LOG(AliL3Log::kInformational,"AliL3Hough::AddAllHistograms()","Timing")
 552     <<"Adding histograms in "<<cpuTime*1000<<" ms"<<ENDLOG;
 553 }
 554
 555 void AliL3Hough::AddTracks()
 556 {
 557   if(!fTracks[0])
 558     {
 559       cerr<<"AliL3Hough::AddTracks : No tracks"<<endl;
 560       return;
 561     }
 562   AliL3TrackArray *tracks = fTracks[0];
 563   for(Int_t i=0; i<tracks->GetNTracks(); i++)
 564     {
 565       AliL3Track *track = tracks->GetCheckedTrack(i);
 566       if(!track) continue;
 567       if(track->GetNHits()!=1) cerr<<"NHITS "<<track->GetNHits()<<endl;
 568       UInt_t *ids = track->GetHitNumbers();
 569       ids[0] = (fCurrentSlice&0x7f)<<25;
 570     }
 571
 572   fGlobalTracks->AddTracks(fTracks[0],0,fCurrentSlice);
 573 }
 574
 575 void AliL3Hough::FindTrackCandidates()
 576 {
 577   //Look for peaks in histograms, and find the track candidates
 578
 579   Int_t n_patches;
 580   if(fAddHistograms)
 581     n_patches = 1; //Histograms have been added.
 582   else
 583     n_patches = fNPatches;
 584
 585   Double_t initTime,cpuTime;
 586   initTime = GetCpuTime();
 587   fBenchmark->Start("Find Maxima");
 588   for(Int_t i=0; i<n_patches; i++)
 589     {
 590       AliL3HoughBaseTransformer *tr = fHoughTransformer[i];
 591       fTracks[i]->Reset();
 592
 593       for(Int_t j=0; j<fNEtaSegments; j++)
 594         {
 595           AliL3Histogram *hist = tr->GetHistogram(j);
 596           if(hist->GetNEntries()==0) continue;
 597           fPeakFinder->Reset();
 598           fPeakFinder->SetHistogram(hist);
 599
 600           fPeakFinder->SetThreshold(fPeakThreshold[i]);
 601           fPeakFinder->FindAdaptedPeaks(6,0.5);
 602
 603           //fPeakFinder->FindMaxima(fPeakThreshold[i]); //Simple maxima finder
 604
 605           for(Int_t k=0; k<fPeakFinder->GetEntries(); k++)
 606             {
 607               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack();
 608               track->SetTrackParameters(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetWeight(k));
 609               track->SetEtaIndex(j);
 610               track->SetEta(tr->GetEta(j,fCurrentSlice));
 611               track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));
 612             }
 613         }
 614       cout<<"Found "<<fTracks[i]->GetNTracks()<<" tracks in patch "<<i<<endl;
 615       fTracks[i]->QSort();
 616     }
 617   fBenchmark->Stop("Find Maxima");
 618   cpuTime = GetCpuTime() - initTime;
 619   LOG(AliL3Log::kInformational,"AliL3Hough::FindTrackCandidates()","Timing")
 620     <<"Maxima finding done in "<<cpuTime*1000<<" ms"<<ENDLOG;
 621 }
 622
 623 void AliL3Hough::InitEvaluate()
 624 {
 625   //Pass the transformer objects to the AliL3HoughEval objects:
 626   //This will provide the evaluation objects with all the necessary
 627   //data and parameters it needs.
 628
 629   for(Int_t i=0; i<fNPatches; i++)
 630     fEval[i]->InitTransformer(fHoughTransformer[i]);
 631 }
 632
 633 Int_t AliL3Hough::Evaluate(Int_t road_width,Int_t nrowstomiss)
 634 {
 635   //Evaluate the tracks, by looking along the road in the raw data.
 636   //If track does not cross all padrows - rows2miss, it is removed from the arrray.
 637   //If histograms were not added, the check is done locally in patch,
 638   //meaning that nrowstomiss is the number of padrows the road can miss with respect
 639   //to the number of rows in the patch.
 640   //If the histograms were added, the comparison is done globally in the _slice_,
 641   //meaing that nrowstomiss is the number of padrows the road can miss with
 642   //respect to the total number of padrows in the slice.
 643   //
 644   //Return value = number of tracks which were removed (only in case of fAddHistograms)
 645
 646   if(!fTracks[0])
 647     {
 648       LOG(AliL3Log::kError,"AliL3Hough::Evaluate","Track Array")
 649         <<"No tracks to work with..."<<ENDLOG;
 650       return 0;
 651     }
 652
 653   Int_t removed_tracks=0;
 654   AliL3TrackArray *tracks=0;
 655
 656   if(fAddHistograms)
 657     {
 658       tracks = fTracks[0];
 659       for(Int_t i=0; i<tracks->GetNTracks(); i++)
 660         {
 661           AliL3Track *track = tracks->GetCheckedTrack(i);
 662           if(!track) continue;
 663           track->SetNHits(0);
 664         }
 665     }
 666
 667   for(Int_t i=0; i<fNPatches; i++)
 668     EvaluatePatch(i,road_width,nrowstomiss);
 669
 670   //Here we check the tracks globally;
 671   //how many good rows (padrows with signal)
 672   //did it cross in the slice
 673   if(fAddHistograms)
 674     {
 675       for(Int_t j=0; j<tracks->GetNTracks(); j++)
 676         {
 677           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
 678
 679           if(track->GetNHits() < AliL3Transform::GetNRows() - nrowstomiss)
 680             {
 681               tracks->Remove(j);
 682               removed_tracks++;
 683             }
 684         }
 685       tracks->Compress();
 686       tracks->QSort();
 687     }
 688
 689   return removed_tracks;
 690 }
 691
 692 void AliL3Hough::EvaluatePatch(Int_t i,Int_t road_width,Int_t nrowstomiss)
 693 {
 694   //Evaluate patch i.
 695
 696   fEval[i]->InitTransformer(fHoughTransformer[i]);
 697   fEval[i]->SetNumOfPadsToLook(road_width);
 698   fEval[i]->SetNumOfRowsToMiss(nrowstomiss);
 699   //fEval[i]->RemoveFoundTracks();
 700
 701   AliL3TrackArray *tracks=0;
 702
 703   if(!fAddHistograms)
 704     tracks = fTracks[i];
 705   else
 706     tracks = fTracks[0];
 707
 708   Int_t nrows=0;
 709   for(Int_t j=0; j<tracks->GetNTracks(); j++)
 710     {
 711       AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
 712       if(!track)
 713         {
 714           LOG(AliL3Log::kWarning,"AliL3Hough::EvaluatePatch","Track array")
 715             <<"Track object missing!"<<ENDLOG;
 716           continue;
 717         }
 718       nrows=0;
 719       Int_t rowrange[2] = {AliL3Transform::GetFirstRow(i),AliL3Transform::GetLastRow(i)};
 720       Bool_t result = fEval[i]->LookInsideRoad(track,nrows,rowrange);
 721       if(fAddHistograms)
 722         {
 723           Int_t pre=track->GetNHits();
 724           track->SetNHits(pre+nrows);
 725         }
 726       else//the track crossed too few good padrows (padrows with signal) in the patch, so remove it
 727         {
 728           if(result == kFALSE)
 729             tracks->Remove(j);
 730         }
 731     }
 732
 733   tracks->Compress();
 734
 735 }
 736
 737 void AliL3Hough::MergeEtaSlices()
 738 {
 739   //Merge tracks found in neighbouring eta slices.
 740   //Removes the track with the lower weight.
 741
 742   fBenchmark->Start("Merge Eta-slices");
 743   AliL3TrackArray *tracks = fTracks[0];
 744   if(!tracks)
 745     {
 746       cerr<<"AliL3Hough::MergeEtaSlices : No tracks "<<endl;
 747       return;
 748     }
 749   for(Int_t j=0; j<tracks->GetNTracks(); j++)
 750     {
 751       AliL3HoughTrack *track1 = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
 752       if(!track1) continue;
 753       for(Int_t k=j+1; k<tracks->GetNTracks(); k++)
 754         {
 755           AliL3HoughTrack *track2 = (AliL3HoughTrack*)tracks->GetCheckedTrack(k);
 756           if(!track2) continue;
 757           if(abs(track1->GetEtaIndex() - track2->GetEtaIndex()) != 1) continue;
 758           if(fabs(track1->GetKappa()-track2->GetKappa()) < 0.006 &&
 759              fabs(track1->GetPsi()- track2->GetPsi()) < 0.1)
 760             {
 761               //cout<<"Merging track in slices "<<track1->GetEtaIndex()<<" "<<track2->GetEtaIndex()<<endl;
 762               if(track1->GetWeight() > track2->GetWeight())
 763                 tracks->Remove(k);
 764               else
 765                 tracks->Remove(j);
 766             }
 767         }
 768     }
 769   fBenchmark->Stop("Merge Eta-slices");
 770   tracks->Compress();
 771 }
 772
 773 void AliL3Hough::WriteTracks(Char_t *path)
 774 {
 775   //cout<<"AliL3Hough::WriteTracks : Sorting the tracsk"<<endl;
 776   //fGlobalTracks->QSort();
 777
 778   Char_t filename[1024];
 779   sprintf(filename,"%s/tracks_%d.raw",path,fEvent);
 780   AliL3MemHandler mem;
 781   mem.SetBinaryOutput(filename);
 782   mem.TrackArray2Binary(fGlobalTracks);
 783   mem.CloseBinaryOutput();
 784   fGlobalTracks->Reset();
 785 }
 786
 787 void AliL3Hough::WriteTracks(Int_t slice,Char_t *path)
 788 {
 789
 790   AliL3MemHandler mem;
 791   Char_t fname[100];
 792   if(fAddHistograms)
 793     {
 794       sprintf(fname,"%s/tracks_ho_%d_%d.raw",path,fEvent,slice);
 795       mem.SetBinaryOutput(fname);
 796       mem.TrackArray2Binary(fTracks[0]);
 797       mem.CloseBinaryOutput();
 798     }
 799   else
 800     {
 801       for(Int_t i=0; i<fNPatches; i++)
 802         {
 803           sprintf(fname,"%s/tracks_ho_%d_%d_%d.raw",path,fEvent,slice,i);
 804           mem.SetBinaryOutput(fname);
 805           mem.TrackArray2Binary(fTracks[i]);
 806           mem.CloseBinaryOutput();
 807         }
 808     }
 809 }
 810
 811 void AliL3Hough::WriteDigits(Char_t *outfile)
 812 {
 813 #ifdef use_aliroot
 814   //Write the current data to a new rootfile.
 815
 816   for(Int_t i=0; i<fNPatches; i++)
 817     {
 818       AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer[i]->GetDataPointer();
 819       fMemHandler[i]->AliDigits2RootFile(tempPt,outfile);
 820     }
 821 #else
 822   cerr<<"AliL3Hough::WriteDigits : You need to compile with AliROOT!"<<endl;
 823   return;
 824 #endif
 825 }
 826
 827 Double_t AliL3Hough::GetCpuTime()
 828 {
 829   //Return the Cputime in seconds.
 830  struct timeval tv;
 831  gettimeofday( &tv, NULL );
 832  return tv.tv_sec+(((Double_t)tv.tv_usec)/1000000.);
 833 }
 834