Rewriting and cleaning up
authorvestbo <vestbo@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 21 Sep 2001 14:01:08 +0000 (14:01 +0000)
committervestbo <vestbo@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 21 Sep 2001 14:01:08 +0000 (14:01 +0000)
HLT/hough/AliL3Hough.cxx
HLT/hough/AliL3Hough.h
HLT/hough/AliL3HoughEval.cxx
HLT/hough/AliL3HoughEval.h
HLT/hough/AliL3HoughMaxFinder.cxx
HLT/hough/AliL3HoughMaxFinder.h
HLT/hough/AliL3HoughTransformer.cxx
HLT/hough/AliL3HoughTransformer.h

index c4e2dda..f3c37e7 100644 (file)
@@ -22,132 +22,118 @@ ClassImp(AliL3Hough)
 AliL3Hough::AliL3Hough()
 {
 
-
-
+  
 }
 
 
-AliL3Hough::AliL3Hough(Int_t n_eta_segments,Int_t xbin,Double_t *xrange,Int_t ybin,Double_t *yrange)
+AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments)
 {
-  
+  fBinary = binary;
+  strcpy(fPath,path);
   fNEtaSegments = n_eta_segments;
-  fNxbin = xbin;
-  fNybin = ybin;
-  fXmin = xrange[0];
-  fXmax = xrange[1];
-  fYmin = yrange[0];
-  fYmax = yrange[1];
-
-  fMemHandler = new AliL3FileHandler();
-  fMaxFinder = new AliL3HoughMaxFinder("KappaPhi");
-  fEval = new AliL3HoughEval();
-  fTransform = new AliL3Transform();
-  fDeleteTrack = kTRUE;
-  fTracks = new AliL3TrackArray("AliL3HoughTrack");
+  Init();
 }
 
 
 AliL3Hough::~AliL3Hough()
 {
-  
-  if(fHoughTransformer)
-    delete fHoughTransformer;
   if(fMemHandler)
-    delete fMemHandler;
-  if(fMaxFinder)
-    delete fMaxFinder;
-  if(fEval)
-    delete fEval;
-  if(fTransform)
-    delete fTransform;
+    DeleteMemory();
+  if(fHoughTransformer)
+    DeleteTransformers();
+  if(fRootFile)
+    {
+      fRootFile->Close();
+      delete fRootFile;
+    }
 }
 
-void AliL3Hough::SetInput(Char_t *input,Bool_t binary)
+void AliL3Hough::DeleteTransformers()
 {
-  if(binary)
+  for(Int_t i=0; i<NPatches; i++)
     {
-      strcpy(fPath,input);
-      fUseBinary = kTRUE;
+      if(!fHoughTransformer[i]) continue;
+      delete fHoughTransformer;
     }
-  else
-    {
-      TFile *file = new TFile(input);
-      fMemHandler->SetAliInput(file);
-      fUseBinary = kFALSE;
-    }
-  
+  delete [] fHoughTransformer;
 }
 
-void AliL3Hough::ProcessSlice(Int_t slice)
+void AliL3Hough::DeleteMemory()
 {
-  
+  for(Int_t i=0; i<NPatches; i++)
+    {
+      if(!fMemHandler[i]) continue;
+      delete fMemHandler[i];
+    }
+  delete [] fMemHandler;
 }
 
-void AliL3Hough::ProcessPatch(Int_t slice,Int_t patch)
+void AliL3Hough::Init()
 {
-  
-  Char_t histname[50];
-  Int_t i;
-  
-  if(fHoughTransformer)
-    delete fHoughTransformer;
-  fHoughTransformer = new AliL3HoughTransformer(slice,patch);//,0,fNEtaSegments);
-  
-  fHistos = new AliL3Histogram*[fNEtaSegments];
-  printf("Allocating %d bytes to histograms\n",fNEtaSegments*sizeof(AliL3Histogram));
-  for(i=0; i<fNEtaSegments; i++)
+  fHoughTransformer = new AliL3HoughTransformer*[NPatches];
+  fMemHandler = new AliL3FileHandler*[NPatches];
+  for(Int_t i=0; i<NPatches; i++)
     {
-      sprintf(histname,"hist%d",i);
-      fHistos[i] = new AliL3Histogram(histname,"",fNxbin,fXmin,fXmax,fNybin,fYmin,fYmax);
+      fHoughTransformer[i] = new AliL3HoughTransformer(1,i,fNEtaSegments);
+      fMemHandler[i] = new AliL3FileHandler();
     }
+  if(!fBinary)
+    fRootFile = new TFile(fPath);
+}
 
-  Char_t name[256];
+void AliL3Hough::TransformSlice(Int_t slice)
+{
   
-  UInt_t ndigits=0;
-  AliL3DigitRowData *digits =0;
-  // fMemHandler->Init(slice,patch,NRows[patch]);
-  //fMemHandler->Init(fTransform);
-  if(fUseBinary)
+  for(Int_t i=0; i<NPatches; i++)
     {
-      fMemHandler->Free();
-      sprintf(name,"%sdigits_%d_%d.raw",fPath,slice,patch);
-      fMemHandler->SetBinaryInput(name);
-      digits = (AliL3DigitRowData *)fMemHandler->CompBinary2Memory(ndigits);
-      fMemHandler->CloseBinaryInput();
+      fHoughTransformer[i]->CreateHistograms(64,-0.006,0.006,64,-0.26,0.26);
+      fHoughTransformer[i]->SetThreshold(3);
+      fMemHandler[i]->Free();
+      UInt_t ndigits=0;
+      AliL3DigitRowData *digits =0;
+      Char_t name[256];
+      if(fBinary)
+       {
+         sprintf(name,"%sdigits_%d_%d.raw",fPath,slice,i);
+         fMemHandler[i]->SetBinaryInput(name);
+         digits = (AliL3DigitRowData *)fMemHandler[i]->CompBinary2Memory(ndigits);
+         fMemHandler[i]->CloseBinaryInput();
+       }
+      else //read data from root file
+       {
+         fMemHandler[i]->SetAliInput(fRootFile);
+         fMemHandler[i]->Init(slice,i,NRows[i]);
+         digits=(AliL3DigitRowData *)fMemHandler[i]->AliDigits2Memory(ndigits); 
+       }
+      fHoughTransformer[i]->SetInputData(ndigits,digits);
+      fHoughTransformer[i]->TransformCircle();
     }
-  else
+  
+}
+
+AliL3Histogram *AliL3Hough::AddHistograms()
+{
+  AliL3Histogram *hist0 = fHoughTransformer[0]->GetHistogram(0);
+  for(Int_t i=1; i<NPatches; i++)
     {
-      digits=(AliL3DigitRowData *)fMemHandler->AliDigits2Memory(ndigits); 
+      AliL3Histogram *hist = fHoughTransformer[i]->GetHistogram(0);
+      hist0->Add(hist);
     }
-  printf("Setting up tables\n");
-  fHoughTransformer->SetHistogram(fHistos[0]);
-  fHoughTransformer->InitTables();
-  fHoughTransformer->SetInputData(ndigits,digits);
-  fEval->SetTransformer(fHoughTransformer);
   
-  AliL3HoughTrack *track;
-  Int_t good_count;
-  while(1)
+  return hist0;
+}
+
+void AliL3Hough::Evaluate(AliL3Histogram *hist)
+{
+  
+  AliL3HoughEval **eval = new AliL3HoughEval*[NPatches];
+  for(Int_t i=0; i<NPatches; i++)
     {
-      fHoughTransformer->TransformTables(fHistos);
-      
-      good_count=0;
-      for(Int_t e=0; e<fNEtaSegments; e++)
-       {
-         fMaxFinder->SetHistogram(fHistos[e]);
-         track = (AliL3HoughTrack*)fMaxFinder->FindPeak(3,0.95,5);
-         if(fEval->LookInsideRoad(track,e,fDeleteTrack))
-           {
-             //Found a good track here
-             fTracks->AddLast(track);
-             good_count++;
-           }
-       }
-      break;
-      if(good_count==0)
-       break;
+      eval[i] = new AliL3HoughEval(fHoughTransformer[i]);
+      eval[i]->DisplayEtaSlice(0,hist);
+      delete eval[i];
     }
-  printf("good_count %d\n",good_count);
   
-}
+  delete [] eval;
 
+}
index 5e66559..86d7982 100644 (file)
@@ -10,42 +10,37 @@ class AliL3FileHandler;
 class AliL3HoughEval;
 class AliL3Transform;
 class AliL3TrackArray;
+class TFile;
 
 class AliL3Hough : public TObject {
   
  private:
-
   Char_t fPath[256];
+  Bool_t fBinary;
   Int_t fNEtaSegments;
-  AliL3Histogram **fHistos; //!
-  AliL3FileHandler *fMemHandler; //!
-  AliL3HoughMaxFinder *fMaxFinder; 
-  AliL3HoughEval *fEval;
-  AliL3HoughTransformer *fHoughTransformer;
-  AliL3Transform *fTransform; //!
-  Bool_t fUseBinary;
-  Bool_t fDeleteTrack;
-  AliL3TrackArray *fTracks; //!
-
-  Int_t fNxbin;
-  Int_t fNybin;
-  Double_t fXmin;
-  Double_t fXmax;
-  Double_t fYmin;
-  Double_t fYmax;
+  AliL3FileHandler **fMemHandler; //!
+  AliL3HoughTransformer **fHoughTransformer; //!
+  TFile *fRootFile; //!
 
+  void DeleteTransformers();
+  void DeleteMemory();
+  void Init();
+  
  public:
-
+  
   AliL3Hough(); 
-  AliL3Hough(Int_t n_eta_segments,Int_t xbin,Double_t *xrange,Int_t ybin,Double_t *yrange);
+  AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments=100);
   virtual ~AliL3Hough();
   
-  void SetInput(Char_t *input,Bool_t binary);
-  void ProcessSlice(Int_t slice);
-  void ProcessPatch(Int_t slice,Int_t patch);
-  void SetDeleteTrack(Bool_t f) {fDeleteTrack = (Bool_t)f;}
+  void TransformSlice(Int_t slice);
+  AliL3Histogram *AddHistograms();
+  void Evaluate(AliL3Histogram *hist);
+
+  //Setters
+  void SetNEtaSegments(Int_t i) {fNEtaSegments = i;}
   
-  AliL3TrackArray *GetTracks() {return fTracks;}
+  //Getters
+  AliL3HoughTransformer *GetTransformer(Int_t i) {if(!fHoughTransformer[i]) return 0; return fHoughTransformer[i];}
 
   ClassDef(AliL3Hough,1)
 
index 8982be4..4d50d2e 100644 (file)
@@ -13,12 +13,12 @@ ClassImp(AliL3HoughEval)
 
 AliL3HoughEval::AliL3HoughEval()
 {
-
+    
 }
 
 AliL3HoughEval::AliL3HoughEval(AliL3HoughTransformer *transformer)
 {
-
+  
   fHoughTransformer = transformer;
   fTransform = new AliL3Transform();
   
@@ -56,7 +56,7 @@ void AliL3HoughEval::GenerateLUT()
     {
       Int_t prow = i - NRows[fPatch][0];
       fRowPointers[prow] = tempPt;
-      tempPt = fHoughTransformer->UpdateDataPointer(tempPt);
+      fHoughTransformer->UpdateDataPointer(tempPt);
     }
   
 }
@@ -161,6 +161,7 @@ void AliL3HoughEval::DisplayEtaSlice(Int_t eta_index,AliL3Histogram *hist)
          UChar_t pad = digPt[j].fPad;
          UChar_t charge = digPt[j].fCharge;
          UShort_t time = digPt[j].fTime;
+         if(charge < fHoughTransformer->GetThreshold()) continue;
          Float_t xyz[3];
          Int_t sector,row;
          fTransform->Slice2Sector(fSlice,padrow,sector,row);
index 8dfaf51..4d2cf98 100644 (file)
@@ -40,7 +40,9 @@ class AliL3HoughEval : public TObject {
   
   //Setters:
   void RemoveFoundTracks() {fRemoveFoundTracks = kTRUE;}
-  
+  void SetNumOfRowsToMiss(Int_t i) {fNumOfRowsToMiss = i;}
+  void SetNumOfPadsToLook(Int_t i) {fNumOfPadsToLook = i;}
+
   ClassDef(AliL3HoughEval,1)
 
 };
index 441ff92..a5cda4d 100644 (file)
@@ -2,6 +2,8 @@
 //Last Modified: 28.6.01
 
 #include <string.h>
+#include <math.h>
+#include <stdlib.h>
 
 #include "AliL3Histogram.h"
 #include "AliL3TrackArray.h"
@@ -39,21 +41,20 @@ AliL3HoughMaxFinder::~AliL3HoughMaxFinder()
 
 }
 
-Int_t *AliL3HoughMaxFinder::FindAbsMaxima()
+void AliL3HoughMaxFinder::FindAbsMaxima(Int_t &max_xbin,Int_t &max_ybin)
 {
   if(!fCurrentHisto)
     {
       printf("AliL3HoughMaxFinder::FindAbsMaxima : No histogram!\n");
-      return 0;
+      return;
     }
   AliL3Histogram *hist = fCurrentHisto;
-
+  
   Int_t xmin = hist->GetFirstXbin();
   Int_t xmax = hist->GetLastXbin();
   Int_t ymin = hist->GetFirstYbin();
   Int_t ymax = hist->GetLastYbin();  
   Int_t bin;
-  Int_t *max_bin = new Int_t[2];
   Stat_t value,max_value=0;
 
   for(Int_t xbin=xmin; xbin<=xmax; xbin++)
@@ -65,12 +66,12 @@ Int_t *AliL3HoughMaxFinder::FindAbsMaxima()
          if(value>max_value)
            {
              max_value = value;
-             max_bin[0] = xbin;
-             max_bin[1] = ybin;
+             max_xbin = xbin;
+             max_ybin = ybin;
            }
        }
     }
-  return max_bin;
+  
 }
 
 AliL3TrackArray *AliL3HoughMaxFinder::FindBigMaxima(AliL3Histogram *hist)
@@ -115,7 +116,7 @@ AliL3TrackArray *AliL3HoughMaxFinder::FindBigMaxima(AliL3Histogram *hist)
              Double_t max_x = hist->GetBinCenterX(xbin);
              Double_t max_y = hist->GetBinCenterY(ybin);
              track = (AliL3HoughTrack*)tracks->NextTrack();
-             track->SetTrackParameters(max_x,max_y,value[12]);
+             track->SetTrackParameters(max_x,max_y,(Int_t)value[12]);
            }
        }
     }
@@ -213,7 +214,6 @@ AliL3TrackArray *AliL3HoughMaxFinder::LookForPeaks(AliL3Histogram *hist,Int_t nb
   Int_t ymax = hist->GetLastYbin();
   
   Int_t weight_loc;
-  Int_t candidate=0;
   for(Int_t xbin=xmin+nbins; xbin <= xmax - nbins; xbin++)
     {
       for(Int_t ybin=ymin+nbins; ybin <= ymax - nbins; ybin++)
@@ -264,38 +264,6 @@ AliL3TrackArray *AliL3HoughMaxFinder::LookForPeaks(AliL3Histogram *hist,Int_t nb
   return tracks;
 }
 
-AliL3TrackArray *AliL3HoughMaxFinder::LookInWindows(AliL3Histogram *hist,Int_t nbins,Int_t t1,Double_t t2,Int_t t3)
-{
-
-  AliL3TrackArray *tracks = new AliL3TrackArray("AliL3HoughTrack");
-  AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
-
-  Int_t xmin = hist->GetFirstXbin();
-  Int_t xmax = hist->GetLastXbin();
-  Int_t ymin = hist->GetFirstYbin();
-  Int_t ymax = hist->GetLastYbin();
-  
-  for(Int_t xbin=xmin+nbins; xbin < xmax - nbins; xbin += nbins)
-    {
-      for(Int_t ybin=ymin+nbins; ybin<ymax-nbins; ybin += nbins)
-       {
-         //Int_t bin = hist->GetBin(xbin,ybin);
-         //if((Int_t)hist->GetBinContent(bin)==0) continue;
-         
-         Int_t xrange[2] = {xbin-nbins,xbin+nbins};
-         Int_t yrange[2] = {ybin-nbins,ybin+nbins};
-         if(LocatePeak(hist,track,xrange,yrange,t1,t2,t3))
-           track = (AliL3HoughTrack*)tracks->NextTrack();
-         else
-           continue;
-       }
-    }
-  tracks->QSort();
-  tracks->RemoveLast();
-  return tracks;
-  
-}
-
 AliL3HoughTrack *AliL3HoughMaxFinder::CalculatePeakInWindow(Int_t *maxbin,Int_t t0,Int_t t1,Double_t t2,Int_t t3)
 {
   //Try to expand the area around the maxbin +- t0
@@ -485,11 +453,223 @@ AliL3HoughTrack *AliL3HoughMaxFinder::CalculatePeakInWindow(Int_t *maxbin,Int_t
   
 }
 
+AliL3HoughTrack *AliL3HoughMaxFinder::FindPeakLine(Double_t rho,Double_t theta)
+{
+  //Peak finder based on a second line transformation on kappa-phi space, 
+  //to use as a baseline.
+
+  if(!fCurrentHisto)
+    {
+      printf("AliL3HoughTransformer::FindPeakLine : No input histogram\n");
+      return 0;
+    }
+  
+  //get the line parameters:
+  Double_t a = -1./tan(theta);
+  Double_t b = rho/sin(theta);
+  
+  printf("rho %f theta %f\n",rho,theta);
+  //now, start looking along the line.
+  Int_t xmin = fCurrentHisto->GetFirstXbin();
+  Int_t xmax = fCurrentHisto->GetLastXbin();
+    
+  Int_t max_weight=0;
+  Double_t max_bin[2];
+  for(Int_t xbin=xmin; xbin<=xmax; xbin++)
+    {
+      Double_t x = fCurrentHisto->GetBinCenterX(xbin);
+      Double_t y = a*x + b;
+      Int_t bin = fCurrentHisto->FindBin(x,y);
+      //printf("x %f y %f weight %d\n",x,y,fCurrentHisto->GetBinContent(bin));
+      if(fCurrentHisto->GetBinContent(bin) > max_weight)
+       {
+         max_weight = (Int_t)fCurrentHisto->GetBinContent(bin);
+         max_bin[0] = x;
+         max_bin[1] = y;
+       }
+    }
+  
+  AliL3HoughTrack *track = new AliL3HoughTrack();
+  track->SetTrackParameters(max_bin[0],max_bin[1],max_weight);
+  return track;
+}
+
+
+
+void AliL3HoughMaxFinder::FindPeak1(Float_t *xpeaks,Float_t *ypeaks,Int_t &n)
+{
+  //testing mutliple peakfinding
+  
+  Int_t max_entries = n;
+  n=0;
+  if(!fCurrentHisto)
+    {
+      printf("AliL3HoughMaxFinder::FindPeak1 : No input histogram\n");
+      return;
+    }  
+  Int_t y_window=2;
+  Int_t x_bin_sides=2;
+  Int_t max_sum=0;
+  
+  Int_t xmin = fCurrentHisto->GetFirstXbin();
+  Int_t xmax = fCurrentHisto->GetLastXbin();
+  Int_t ymin = fCurrentHisto->GetFirstYbin();
+  Int_t ymax = fCurrentHisto->GetLastYbin();
+  Int_t nbinsx = fCurrentHisto->GetNbinsX()+1;
+  
+  AxisWindow **windowPt = new AxisWindow*[nbinsx];
+  AxisWindow **anotherPt = new AxisWindow*[nbinsx];
+  
+  for(Int_t i=0; i<nbinsx; i++)
+    {
+      windowPt[i] = new AxisWindow;
+      bzero((void*)windowPt[i],sizeof(AxisWindow));
+      anotherPt[i] = windowPt[i];
+    }
+  
+  for(Int_t xbin=xmin; xbin<=xmax; xbin++)
+    {
+      max_sum = 0;
+      for(Int_t ybin=ymin; ybin<=ymax-y_window; ybin++)
+       {
+         Int_t sum_in_window=0;
+         for(Int_t b=ybin; b<ybin+y_window; b++)
+           {
+             //inside window
+             Int_t bin = fCurrentHisto->GetBin(xbin,b);
+             sum_in_window += (Int_t)fCurrentHisto->GetBinContent(bin);
+           }
+         
+         if(sum_in_window > max_sum)
+           {
+             max_sum = sum_in_window;
+             windowPt[xbin]->ymin = ybin;
+             windowPt[xbin]->ymax = ybin + y_window;
+             windowPt[xbin]->weight = sum_in_window;
+             windowPt[xbin]->xbin = xbin;
+           }
+       }
+    }
+
+  //Sort the windows according to the weight
+  SortPeaks(windowPt,0,nbinsx);
+  
+  //for(Int_t i=0; i<nbinsx; i++)
+  //printf("xbin %f weight %d\n",fCurrentHisto->GetBinCenterX(windowPt[i]->xbin),windowPt[i]->weight);
+  
+  Float_t top,butt;
+  for(Int_t i=0; i<nbinsx; i++)
+    {
+      top=butt=0;
+      Int_t xbin = windowPt[i]->xbin;
+      
+      if(xbin<xmin || xbin > xmax-1) continue;
+      
+      //Check if this is really a peak
+      if(anotherPt[xbin-1]->weight > windowPt[i]->weight ||
+        anotherPt[xbin+1]->weight > windowPt[i]->weight)
+       continue;
+
+      for(Int_t j=windowPt[i]->ymin; j<windowPt[i]->ymax; j++)
+       {
+         //Calculate the mean in y direction:
+         Int_t bin = fCurrentHisto->GetBin(windowPt[i]->xbin,j);
+         top += (fCurrentHisto->GetBinCenterY(j))*(fCurrentHisto->GetBinContent(bin));
+         butt += fCurrentHisto->GetBinContent(bin);
+       }
+      xpeaks[n] = fCurrentHisto->GetBinCenterX(windowPt[i]->xbin);
+      ypeaks[n] = top/butt;
+      n++;
+      if(n==max_entries) break;
+    }
+
+  
+  //Improve the peaks by including the region around in x.
+  Float_t ytop,ybutt;
+  for(Int_t i=0; i<n; i++)
+    {
+      Int_t xbin = fCurrentHisto->FindXbin(xpeaks[i]);
+      if(xbin - x_bin_sides < xmin || xbin + x_bin_sides > xmax) continue;
+      top=butt=0;
+      ytop=0,ybutt=0;    
+      for(Int_t j=xbin-x_bin_sides; j<=xbin+x_bin_sides; j++)
+       {
+         top += fCurrentHisto->GetBinCenterX(j)*anotherPt[j]->weight;
+         butt += anotherPt[j]->weight;
+         for(Int_t k=anotherPt[j]->ymin; k<anotherPt[j]->ymax; k++)
+           {
+             Int_t bin = fCurrentHisto->GetBin(j,k);
+             ytop += (fCurrentHisto->GetBinCenterY(k))*(fCurrentHisto->GetBinContent(bin));
+             ybutt += fCurrentHisto->GetBinContent(bin);
+           }
+       }
+      xpeaks[i] = top/butt;
+      ypeaks[i] = ytop/ybutt;
+    }
+  
+  for(Int_t i=0; i<nbinsx; i++)
+    delete windowPt[i];
+  delete [] windowPt;
+  delete [] anotherPt;
+  
+}
+
+void AliL3HoughMaxFinder::SortPeaks(struct AxisWindow **a,Int_t first,Int_t last)
+{
+  //General sorting routine
+  //Sort according to PeakCompare()
+
+  static struct AxisWindow *tmp;
+  static int i;           // "static" to save stack space
+  int j;
+  
+  while (last - first > 1) {
+    i = first;
+    j = last;
+    for (;;) {
+      while (++i < last && PeakCompare(a[i], a[first]) < 0)
+       ;
+      while (--j > first && PeakCompare(a[j], a[first]) > 0)
+       ;
+      if (i >= j)
+       break;
+      
+      tmp  = a[i];
+      a[i] = a[j];
+      a[j] = tmp;
+    }
+    if (j == first) {
+      ++first;
+      continue;
+    }
+    tmp = a[first];
+    a[first] = a[j];
+    a[j] = tmp;
+    if (j - first < last - (j + 1)) {
+      SortPeaks(a, first, j);
+      first = j + 1;   // QSort(j + 1, last);
+    } else {
+      SortPeaks(a, j + 1, last);
+      last = j;        // QSort(first, j);
+    }
+  }
+  
+}
+
+Int_t AliL3HoughMaxFinder::PeakCompare(struct AxisWindow *a,struct AxisWindow *b)
+{
+  if(a->weight < b->weight) return 1;
+  if(a->weight > b->weight) return -1;
+  return 0;
+
+}
 
 AliL3HoughTrack *AliL3HoughMaxFinder::FindPeak(Int_t t1,Double_t t2,Int_t t3)
 {
   //Attempt of a more sophisticated peak finder.
-  
+  //Finds the best peak in the histogram, and returns the corresponding
+  //track object.
+
   if(!fCurrentHisto)
     {
       printf("AliL3HoughMaxFinder::FindPeak : No histogram!!\n");
index 2776e4f..229f918 100644 (file)
@@ -7,13 +7,20 @@ class AliL3Histogram;
 class AliL3TrackArray;
 class AliL3HoughTrack;
 
+struct AxisWindow
+{
+  Int_t ymin;
+  Int_t ymax;
+  Int_t xbin;
+  Int_t weight;
+};
 
 class AliL3HoughMaxFinder : public TObject {
   
  private:
 
   Int_t fThreshold;
-  AliL3Histogram *fCurrentHisto;
+  AliL3Histogram *fCurrentHisto;  //!
 
   Char_t fHistoType;
 
@@ -22,14 +29,17 @@ class AliL3HoughMaxFinder : public TObject {
   AliL3HoughMaxFinder(Char_t *histotype,AliL3Histogram *hist=0);
   virtual ~AliL3HoughMaxFinder();
 
-  Int_t *FindAbsMaxima();
+  void FindAbsMaxima(Int_t &max_xbin,Int_t &max_ybin);
   AliL3TrackArray *FindBigMaxima(AliL3Histogram *hist);
   AliL3TrackArray *FindMaxima(AliL3Histogram *hist,Int_t *rowrange=0,Int_t ref_row=0);
   AliL3TrackArray *LookForPeaks(AliL3Histogram *hist,Int_t nbins);
-  AliL3TrackArray *LookInWindows(AliL3Histogram *hist,Int_t nbins,Int_t t1,Double_t t2,Int_t t3);
-  Bool_t LocatePeak(AliL3Histogram *hist,AliL3HoughTrack *track,Int_t *xrange,Int_t *yrange,Int_t t1,Double_t t2,Int_t t3);
+  
   AliL3HoughTrack *FindPeak(Int_t t1,Double_t t2,Int_t t3);
+  AliL3HoughTrack *FindPeakLine(Double_t rho,Double_t theta);
   AliL3HoughTrack *CalculatePeakInWindow(Int_t *maxbin,Int_t t0,Int_t t1,Double_t t2,Int_t t3);
+  void FindPeak1(Float_t *xpeaks,Float_t *ypeaks,Int_t &n);
+  void SortPeaks(struct AxisWindow **a,Int_t first,Int_t last);
+  Int_t PeakCompare(struct AxisWindow *a,struct AxisWindow *b);
   
   void SetThreshold(Int_t f) {fThreshold = f;}
   
index bf3136e..eb6cb9a 100644 (file)
@@ -39,14 +39,17 @@ void AliL3HoughTransformer::DeleteHistograms()
   if(!fParamSpace)
     return;
   for(Int_t i=0; i<fNEtaSegments; i++)
-    delete fParamSpace[i];
+    {
+      if(!fParamSpace[i]) continue;
+      delete fParamSpace[i];
+    }
   delete [] fParamSpace;
 }
 
 void AliL3HoughTransformer::CreateHistograms(Int_t nxbin,Double_t xmin,Double_t xmax,
                                             Int_t nybin,Double_t ymin,Double_t ymax)
 {
-      
+  
   fParamSpace = new AliL3Histogram*[fNEtaSegments];
   
   Char_t histname[256];
@@ -63,17 +66,17 @@ void AliL3HoughTransformer::SetInputData(UInt_t ndigits,AliL3DigitRowData *ptr)
   fDigitRowData = ptr;
 }
 
-AliL3DigitRowData *AliL3HoughTransformer::UpdateDataPointer(AliL3DigitRowData *tempPt)
+void AliL3HoughTransformer::UpdateDataPointer(AliL3DigitRowData *&tempPt)
 {
   //Update the data pointer to the next padrow
-
+  
   Byte_t *tmp = (Byte_t*)tempPt;
   Int_t size = sizeof(AliL3DigitRowData) + tempPt->fNDigit*sizeof(AliL3DigitData);
   tmp += size;
-  return (AliL3DigitRowData*)tmp;
+  tempPt = (AliL3DigitRowData*)tmp;
 }
 
-void AliL3HoughTransformer::Transform()
+void AliL3HoughTransformer::TransformCircle()
 {
   //Transform the input data with a circle HT.
   
@@ -82,7 +85,7 @@ void AliL3HoughTransformer::Transform()
   AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fDigitRowData;
   if(!tempPt || fNDigitRowData==0)
     {
-      printf("\nAliL3HoughTransformer::TransformTables : No input data!!!\n\n");
+      printf("\nAliL3HoughTransformer::TransformCircle : No input data!!!\n\n");
       return;
     }
   
@@ -92,7 +95,7 @@ void AliL3HoughTransformer::Transform()
       AliL3DigitData *digPt = tempPt->fDigitData;
       if(i != (Int_t)tempPt->fRow)
        {
-         printf("AliL3HoughTransform::Transform : Mismatching padrow numbering\n");
+         printf("AliL3HoughTransform::TransformCircle : Mismatching padrow numbering\n");
          continue;
        }
       for(UInt_t j=0; j<tempPt->fNDigit; j++)
@@ -106,7 +109,7 @@ void AliL3HoughTransformer::Transform()
          Float_t xyz[3];
          fTransform->Slice2Sector(fSlice,i,sector,row);
          fTransform->Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
-         Float_t eta = fTransform->GetEta(xyz);
+         Double_t eta = fTransform->GetEta(xyz);
          Int_t eta_index = (Int_t)(eta/etaslice);
          if(eta_index < 0 || eta_index >= fNEtaSegments)
            continue;
@@ -115,7 +118,7 @@ void AliL3HoughTransformer::Transform()
          AliL3Histogram *hist = fParamSpace[eta_index];
          if(!hist)
            {
-             printf("AliL3HoughTransformer::Transform : Error getting histogram in index %d\n",eta_index);
+             printf("AliL3HoughTransformer::TransformCircle : Error getting histogram in index %d\n",eta_index);
              continue;
            }
 
@@ -131,7 +134,62 @@ void AliL3HoughTransformer::Transform()
              hist->Fill(kappa,phi0,charge);
            }
        }
-      tempPt = UpdateDataPointer(tempPt);
+      UpdateDataPointer(tempPt);
     }
 }
 
+void AliL3HoughTransformer::TransformLine()
+{
+  //Do a line transform on the data.
+
+  
+  AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fDigitRowData;
+  if(!tempPt || fNDigitRowData==0)
+    {
+      printf("\nAliL3HoughTransformer::TransformLine : No input data!!!\n\n");
+      return;
+    }
+  
+  Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments;
+  for(Int_t i=NRows[fPatch][0]; i<=NRows[fPatch][1]; i++)
+    {
+      AliL3DigitData *digPt = tempPt->fDigitData;
+      if(i != (Int_t)tempPt->fRow)
+       {
+         printf("AliL3HoughTransform::TransformLine : Mismatching padrow numbering\n");
+         continue;
+       }
+      for(UInt_t j=0; j<tempPt->fNDigit; j++)
+       {
+         UShort_t charge = digPt[j].fCharge;
+         UChar_t pad = digPt[j].fPad;
+         UShort_t time = digPt[j].fTime;
+         if(charge < fThreshold)
+           continue;
+         Int_t sector,row;
+         Float_t xyz[3];
+         fTransform->Slice2Sector(fSlice,i,sector,row);
+         fTransform->Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
+         Float_t eta = fTransform->GetEta(xyz);
+         Int_t eta_index = (Int_t)(eta/etaslice);
+         if(eta_index < 0 || eta_index >= fNEtaSegments)
+           continue;
+         
+         //Get the correct histogram:
+         AliL3Histogram *hist = fParamSpace[eta_index];
+         if(!hist)
+           {
+             printf("AliL3HoughTransformer::TransformLine : Error getting histogram in index %d\n",eta_index);
+             continue;
+           }
+         for(Int_t xbin=hist->GetFirstXbin(); xbin<hist->GetLastXbin(); xbin++)
+           {
+             Double_t theta = hist->GetBinCenterX(xbin);
+             Double_t rho = xyz[0]*cos(theta) + xyz[1]*sin(theta);
+             hist->Fill(theta,rho,charge);
+           }
+       }
+      UpdateDataPointer(tempPt);
+    }
+  
+}
index 3c0b40e..aad4f7c 100644 (file)
@@ -19,7 +19,7 @@ class AliL3HoughTransformer : public TObject {
   Int_t fThreshold;
   AliL3Transform *fTransform; //!
   
-  //Pointer to histograms
+  //Pointers to histograms
   AliL3Histogram **fParamSpace; //!
 
   //Data pointers
@@ -34,18 +34,24 @@ class AliL3HoughTransformer : public TObject {
   virtual ~AliL3HoughTransformer();
   
   void SetInputData(UInt_t ndigits,AliL3DigitRowData *ptr);
-  AliL3DigitRowData *UpdateDataPointer(AliL3DigitRowData *tempPt);
-  void CreateHistograms(Int_t nxbin=70,Double_t xmin=-0.006,Double_t xmax=0.006,
-                       Int_t nybin=70,Double_t ymin=-0.26,Double_t ymax=0.26);
-  void Transform();
-  
+  void UpdateDataPointer(AliL3DigitRowData *&tempPt);
+  void CreateHistograms(Int_t nxbin=64,Double_t xmin=-0.006,Double_t xmax=0.006,
+                       Int_t nybin=64,Double_t ymin=-0.26,Double_t ymax=0.26);
+  void TransformCircle();
+  void TransformLine();
+
+  //Getters
   Int_t GetSlice() {return fSlice;}
   Int_t GetPatch() {return fPatch;}
   Int_t GetNEtaSegments() {return fNEtaSegments;}
+  Int_t GetThreshold() {return fThreshold;}
   Double_t GetEtaMin() {return fEtaMin;}
   Double_t GetEtaMax() {return fEtaMax;}
   void *GetDataPointer() {return (void*)fDigitRowData;}
   AliL3Histogram *GetHistogram(Int_t eta_index);
+  
+  //setters
+  void SetThreshold(Int_t i) {fThreshold = i;}
 
   ClassDef(AliL3HoughTransformer,1)