HLT/hough/AliL3Hough.cxx

   1 //$Id$
   2
   3 // Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
   4 //*-- Copyright &copy ASV
   5
   6
   7 #include <string.h>
   8 #include <TCanvas.h>
   9 #include <TFile.h>
  10
  11 #include "AliL3HoughMerger.h"
  12 #include "AliL3HoughIntMerger.h"
  13 #include "AliL3HoughGlobalMerger.h"
  14 #include "AliL3Logging.h"
  15 #include "AliL3Histogram.h"
  16 #include "AliL3Hough.h"
  17 #include "AliL3HoughTransformer.h"
  18 #include "AliL3HoughMaxFinder.h"
  19 #include "AliL3FileHandler.h"
  20 #include "AliL3DigitData.h"
  21 #include "AliL3HoughEval.h"
  22 #include "AliL3Transform.h"
  23 #include "AliL3Defs.h"
  24 #include "AliL3TrackArray.h"
  25 #include "AliL3HoughTrack.h"
  26 #include "AliL3Benchmark.h"
  27
  28 //_____________________________________________________________
  29 // AliL3Hough
  30 //
  31 // Base class for the Hough transform
  32 //
  33
  34
  35 ClassImp(AliL3Hough)
  36
  37 AliL3Hough::AliL3Hough()
  38 {
  39   //Constructor
  40
  41   fBinary = kFALSE;
  42   fNEtaSegments = 0;
  43   fAddHistograms = kFALSE;
  44   fDoIterative = kFALSE;
  45   fWriteDigits=kFALSE;
  46   fNPatches=0;
  47 }
  48
  49
  50 AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments)
  51 {
  52   fBinary = binary;
  53   strcpy(fPath,path);
  54   fNEtaSegments = n_eta_segments;
  55   fAddHistograms = kFALSE;
  56   fDoIterative = kFALSE;
  57   fWriteDigits = kFALSE;
  58   Init();
  59 }
  60
  61
  62 AliL3Hough::~AliL3Hough()
  63 {
  64   CleanUp();
  65   if(fMerger)
  66     delete fMerger;
  67   if(fInterMerger)
  68     delete fInterMerger;
  69   if(fPeakFinder)
  70     delete fPeakFinder;
  71   if(fGlobalMerger)
  72     delete fGlobalMerger;
  73 }
  74
  75 void AliL3Hough::CleanUp()
  76 {
  77   //Cleanup memory
  78
  79   for(Int_t i=0; i<fNPatches; i++)
  80     {
  81       if(fTracks[i]) delete fTracks[i];
  82       if(fEval[i]) delete fEval[i];
  83       if(fHoughTransformer[i]) delete fHoughTransformer[i];
  84       if(fMemHandler[i]) delete fMemHandler[i];
  85     }
  86
  87   /*Shitty compiler doesn't allow this:
  88
  89   if(fTracks) delete [] fTracks;
  90   if(fEval) delete [] fEval;
  91   if(fHoughTransformer) delete [] fHoughTransformer;
  92   if(fMemHandler) delete [] fMemHandler;
  93   */
  94 }
  95
  96 void AliL3Hough::Init()
  97 {
  98   fNPatches = NPatches;
  99   fHoughTransformer = new AliL3HoughTransformer*[fNPatches];
 100   fMemHandler = new AliL3FileHandler*[fNPatches];
 101   fTracks = new AliL3TrackArray*[fNPatches];
 102   fEval = new AliL3HoughEval*[fNPatches];
 103   for(Int_t i=0; i<fNPatches; i++)
 104     {
 105       fHoughTransformer[i] = new AliL3HoughTransformer(1,i,fNEtaSegments);
 106       fHoughTransformer[i]->CreateHistograms(64,-0.003,0.003,64,-0.26,0.26);
 107       fHoughTransformer[i]->SetThreshold(3);
 108       fEval[i] = new AliL3HoughEval();
 109       fTracks[i] = new AliL3TrackArray("AliL3HoughTrack");
 110       fMemHandler[i] = new AliL3FileHandler();
 111       if(!fBinary)
 112         fMemHandler[i]->SetAliInput(fPath);
 113     }
 114   fPeakFinder = new AliL3HoughMaxFinder("KappaPhi");
 115   fMerger = new AliL3HoughMerger(fNPatches);
 116   fInterMerger = new AliL3HoughIntMerger();
 117 }
 118
 119 void AliL3Hough::Process(Int_t minslice,Int_t maxslice)
 120 {
 121   //Process all slices [minslice,maxslice].
 122   fGlobalMerger = new AliL3HoughGlobalMerger(minslice,maxslice);
 123
 124   for(Int_t i=minslice; i<=maxslice; i++)
 125     {
 126       ReadData(i);
 127       Transform();
 128       if(fAddHistograms)
 129         AddAllHistograms();
 130       FindTrackCandidates();
 131       Evaluate();
 132       if(fWriteDigits)
 133         WriteDigits();
 134       fGlobalMerger->FillTracks(fTracks[0],i);
 135     }
 136
 137
 138 }
 139
 140 void AliL3Hough::ReadData(Int_t slice)
 141 {
 142   //Read data from files, binary or root.
 143
 144   for(Int_t i=0; i<fNPatches; i++)
 145     {
 146       fMemHandler[i]->Free();
 147       UInt_t ndigits=0;
 148       AliL3DigitRowData *digits =0;
 149       Char_t name[256];
 150       if(fBinary)//take input data from binary files
 151         {
 152           sprintf(name,"%sdigits_%d_%d.raw",fPath,slice,i);
 153           fMemHandler[i]->SetBinaryInput(name);
 154           digits = (AliL3DigitRowData *)fMemHandler[i]->CompBinary2Memory(ndigits);
 155           fMemHandler[i]->CloseBinaryInput();
 156         }
 157       else //read data from root file
 158         {
 159           fMemHandler[i]->Init(slice,i,NRows[i]);
 160           digits=(AliL3DigitRowData *)fMemHandler[i]->AliDigits2Memory(ndigits);
 161         }
 162       fHoughTransformer[i]->SetInputData(ndigits,digits);
 163     }
 164 }
 165
 166 void AliL3Hough::Transform()
 167 {
 168   //Transform all data given to the transformer within the given slice
 169   //(after ReadData(slice))
 170
 171   Double_t initTime,finalTime;
 172   for(Int_t i=0; i<fNPatches; i++)
 173     {
 174       fHoughTransformer[i]->Reset();//Reset the histograms
 175       initTime = AliL3Benchmark::GetCpuTime();
 176       fHoughTransformer[i]->TransformCircle();
 177       finalTime = AliL3Benchmark::GetCpuTime();
 178       LOG(AliL3Log::kInformational,"AliL3Hough::Transform","Timing")
 179         <<AliL3Log::kDec<<"Transform finished in "<<(finalTime-initTime)*1000<<"ms"<<ENDLOG;
 180     }
 181 }
 182
 183 void AliL3Hough::MergePatches()
 184 {
 185   if(fAddHistograms) //Nothing to merge here
 186     return;
 187   AliL3Transform *tr = new AliL3Transform();
 188   fMerger->SetTransformer(tr);
 189   fMerger->MergePatches(kTRUE);
 190   delete tr;
 191 }
 192
 193 void AliL3Hough::MergeInternally()
 194 {
 195   if(fAddHistograms)
 196     fInterMerger->FillTracks(fTracks[0]);
 197   else
 198     fInterMerger->FillTracks(fMerger->GetOutTracks());
 199
 200   fInterMerger->MMerge();
 201 }
 202
 203 void AliL3Hough::ProcessSliceIter()
 204 {
 205   //Process current slice (after ReadData(slice)) iteratively.
 206
 207   for(Int_t i=0; i<fNPatches; i++)
 208     {
 209       ProcessPatchIter(i);
 210       fMerger->FillTracks(fTracks[i],i); //Copy tracks to merger
 211     }
 212
 213 }
 214
 215 void AliL3Hough::ProcessPatchIter(Int_t patch)
 216 {
 217   //Process patch in a iterative way.
 218   //transform + peakfinding + evaluation + transform +...
 219
 220   Int_t num_of_tries = 10;
 221   AliL3HoughTransformer *tr = fHoughTransformer[patch];
 222   AliL3TrackArray *tracks = fTracks[patch];
 223   tracks->Reset();
 224   AliL3HoughEval *ev = fEval[patch];
 225   ev->InitTransformer(tr);
 226   ev->RemoveFoundTracks();
 227   ev->SetNumOfRowsToMiss(2);
 228   AliL3Histogram *hist;
 229   for(Int_t t=0; t<num_of_tries; t++)
 230     {
 231       tr->Reset();
 232       tr->TransformCircle();
 233       for(Int_t i=0; i<fNEtaSegments; i++)
 234         {
 235           hist = tr->GetHistogram(i);
 236           if(hist->GetNEntries()==0) continue;
 237           fPeakFinder->SetHistogram(hist);
 238           Int_t n=1;
 239           Int_t x[n],y[n];
 240           fPeakFinder->FindAbsMaxima(*x,*y);
 241           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
 242           track->SetTrackParameters(hist->GetBinCenterX(*x),hist->GetBinCenterY(*y),1);
 243           if(!ev->LookInsideRoad(track,i))
 244             {
 245               tracks->Remove(tracks->GetNTracks()-1);
 246               tracks->Compress();
 247             }
 248         }
 249     }
 250   LOG(AliL3Log::kInformational,"AliL3Hough::ProcessPatch","NTracks")
 251     <<AliL3Log::kDec<<"Found "<<tracks->GetNTracks()<<" tracks in patch "<<patch<<ENDLOG;
 252 }
 253
 254
 255 void AliL3Hough::AddAllHistograms()
 256 {
 257   //Add the histograms within one etaslice.
 258   //Resulting histogram are in patch=0.
 259
 260   for(Int_t i=0; i<fNEtaSegments; i++)
 261     {
 262       AliL3Histogram *hist0 = fHoughTransformer[0]->GetHistogram(i);
 263       for(Int_t j=1; j<fNPatches; j++)
 264         {
 265           AliL3Histogram *hist = fHoughTransformer[j]->GetHistogram(i);
 266           hist0->Add(hist);
 267         }
 268     }
 269   fAddHistograms = kTRUE;
 270 }
 271
 272 void AliL3Hough::FindTrackCandidates()
 273 {
 274   //Look for peaks in histograms, and find the track candidates
 275
 276   Int_t n_patches;
 277   if(fAddHistograms)
 278     n_patches = 1; //Histograms has been added.
 279   else
 280     n_patches = fNPatches;
 281
 282   for(Int_t i=0; i<n_patches; i++)
 283     {
 284       AliL3HoughTransformer *tr = fHoughTransformer[i];
 285       fTracks[i]->Reset();
 286       for(Int_t j=0; j<fNEtaSegments; j++)
 287         {
 288           AliL3Histogram *hist = tr->GetHistogram(j);
 289           if(hist->GetNEntries()==0) continue;
 290           fPeakFinder->SetHistogram(hist);
 291           Int_t n=10;
 292           Float_t x[10];
 293           Float_t y[10];
 294           Int_t weight[10];
 295           fPeakFinder->FindPeak1(x,y,weight,n,1);
 296           for(Int_t k=0; k<n; k++)
 297             {
 298               if(weight[k] == 0) continue;
 299               AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack();
 300               track->SetTrackParameters(x[k],y[k],weight[k]);
 301               track->SetEtaIndex(j);
 302               track->SetEta((Double_t)(j*tr->GetEtaSlice()));
 303               track->SetRowRange(NRows[0][0],NRows[5][1]);
 304             }
 305         }
 306       fTracks[i]->QSort();
 307     }
 308 }
 309
 310 void AliL3Hough::Evaluate(Int_t road_width)
 311 {
 312   //Evaluate the tracks, by looking along the road in the raw data.
 313
 314
 315   if(!fTracks[0])
 316     {
 317       LOG(AliL3Log::kError,"AliL3Hough::Evaluate","Track Array")
 318         <<"No tracks to work with..."<<ENDLOG;
 319       return;
 320     }
 321
 322   printf("Number of tracks before evaluation %d\n",fTracks[0]->GetNTracks());
 323   AliL3TrackArray *tracks;
 324   for(Int_t i=0; i<fNPatches; i++)
 325     {
 326       fEval[i]->InitTransformer(fHoughTransformer[i]);
 327       fEval[i]->SetNumOfRowsToMiss(2);
 328       fEval[i]->SetNumOfPadsToLook(road_width);
 329       if(fAddHistograms)
 330         tracks = fTracks[0];
 331       else
 332         tracks = fTracks[i];
 333       for(Int_t j=0; j<tracks->GetNTracks(); j++)
 334         {
 335           AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
 336           if(!track)
 337             {
 338               LOG(AliL3Log::kWarning,"AliL3Hough::Evaluate","Track array")
 339                 <<"Track object missing!"<<ENDLOG;
 340               continue;
 341             }
 342
 343           if(!fEval[i]->LookInsideRoad(track,track->GetEtaIndex()))
 344             tracks->Remove(j);
 345           if(fAddHistograms)
 346             track->SetRowRange(NRows[0][0],NRows[5][1]);//All rows included
 347         }
 348       tracks->Compress();
 349       tracks->QSort(); //Sort the tracks according to weight
 350
 351       if(!fAddHistograms)
 352         fMerger->FillTracks(tracks,i); //Copy tracks to the track merger
 353     }
 354
 355 }
 356
 357 void AliL3Hough::EvaluateWithEta()
 358 {
 359   if(!fTracks[0])
 360     {
 361       printf("AliL3Hough::EvaluateWithEta: NO TRACKS\n");
 362       return;
 363     }
 364   printf("Number of tracks before evaluation %d\n",fTracks[0]->GetNTracks());
 365
 366   for(Int_t i=0; i<fNPatches; i++)
 367     {
 368       fEval[i]->InitTransformer(fHoughTransformer[i]);
 369       fEval[i]->FindEta(fTracks[0]);
 370     }
 371   fMerger->FillTracks(fTracks[0],0);
 372   fMerger->MergeEtaSlices(0);
 373 }
 374
 375 void AliL3Hough::WriteTracks()
 376 {
 377   AliL3MemHandler *mem = new AliL3MemHandler();
 378   mem->SetBinaryOutput("tracks.raw");
 379   mem->TrackArray2Binary(fTracks[0]);
 380   mem->CloseBinaryOutput();
 381   delete mem;
 382
 383 }
 384
 385 void AliL3Hough::WriteDigits(Char_t *outfile)
 386 {
 387   //Write the current data to a new rootfile.
 388
 389   for(Int_t i=0; i<fNPatches; i++)
 390     {
 391       AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer[i]->GetDataPointer();
 392       fMemHandler[i]->AliDigits2RootFile(tempPt,outfile);
 393     }
 394
 395 }