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