Changes done to include new ALiL3HoughTransformerVhdl.
authorloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 29 May 2002 11:37:40 +0000 (11:37 +0000)
committerloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 29 May 2002 11:37:40 +0000 (11:37 +0000)
HLT/hough/AliL3Histogram.cxx
HLT/hough/AliL3Hough.cxx
HLT/hough/AliL3Hough.h
HLT/hough/AliL3HoughEval.cxx
HLT/hough/AliL3HoughLinkDef.h
HLT/hough/AliL3HoughTrack.cxx
HLT/hough/AliL3HoughTrack.h
HLT/hough/AliL3HoughTransformerVhdl.cxx [new file with mode: 0644]
HLT/hough/AliL3HoughTransformerVhdl.h [new file with mode: 0644]
HLT/hough/Makefile

index 2b182d1..9c0fbf5 100644 (file)
@@ -251,10 +251,10 @@ void AliL3Histogram::Draw(Char_t *option)
       fRootHisto->AddBinContent(bin,GetBinContent(bin));
     }
   
+  fRootHisto->SetStats(kFALSE);
   fRootHisto->Draw(option);
   return;
 #endif
   cerr<<"AliL3Histogram::Draw : You need to compile with ROOT in order to draw histogram"<<endl;
   
 }
-
index 1573bfb..acd6602 100644 (file)
@@ -13,6 +13,7 @@
 #include "AliL3Histogram.h"
 #include "AliL3Hough.h"
 #include "AliL3HoughTransformer.h"
+#include "AliL3HoughTransformerVhdl.h"
 #include "AliL3HoughMaxFinder.h"
 #ifdef use_aliroot
 #include "AliL3FileHandler.h"
@@ -61,8 +62,21 @@ AliL3Hough::AliL3Hough()
   fMerger = 0;
   fInterMerger = 0;
   fGlobalMerger = 0;
+  fversion = 0;
 }
 
+AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments,Int_t tversion=0)
+{
+  //Default ctor.
+
+  fBinary = binary;
+  strcpy(fPath,path);
+  fNEtaSegments = n_eta_segments;
+  fAddHistograms = kFALSE;
+  fDoIterative = kFALSE; 
+  fWriteDigits = kFALSE;
+  fversion = tversion;
+}
 
 AliL3Hough::~AliL3Hough()
 {
@@ -110,6 +124,7 @@ void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit
   fUse8bits = bit8;
   
   AliL3Transform::Init(fPath);
+
   fNPatches = AliL3Transform::GetNPatches();
   fHoughTransformer = new AliL3HoughBaseTransformer*[fNPatches];
   fMemHandler = new AliL3MemHandler*[fNPatches];
@@ -117,9 +132,18 @@ void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit
   fEval = new AliL3HoughEval*[fNPatches];
   for(Int_t i=0; i<fNPatches; i++)
     {
-      fHoughTransformer[i] = new AliL3HoughTransformer(0,i,fNEtaSegments);
-      fHoughTransformer[i]->CreateHistograms(64,0.1,64,-30,30);
-      fHoughTransformer[i]->SetThreshold(3);
+      switch (fversion){
+      case 1: 
+       fHoughTransformer[i] = new AliL3HoughTransformerVhdl(1,i,fNEtaSegments);
+       fHoughTransformer[i]->CreateHistograms(180,0.1,180,-90,90);
+       fHoughTransformer[i]->SetThreshold(3);
+       break;
+      default:
+       fHoughTransformer[i] = new AliL3HoughTransformer(1,i,fNEtaSegments);
+       fHoughTransformer[i]->CreateHistograms(64,0.1,64,-30,30);
+       fHoughTransformer[i]->SetThreshold(3);
+      }
+
       fEval[i] = new AliL3HoughEval();
       fTracks[i] = new AliL3TrackArray("AliL3HoughTrack");
       if(fUse8bits)
@@ -156,8 +180,6 @@ void AliL3Hough::Process(Int_t minslice,Int_t maxslice)
       Evaluate();
       fGlobalMerger->FillTracks(fTracks[0],i);
     }
-  
-  
 }
 
 void AliL3Hough::ReadData(Int_t slice,Int_t eventnr=0)
@@ -194,7 +216,7 @@ void AliL3Hough::ReadData(Int_t slice,Int_t eventnr=0)
     }
 }
 
-void AliL3Hough::Transform(Int_t row_range)
+void AliL3Hough::Transform(Int_t row_range=-1)
 {
   //Transform all data given to the transformer within the given slice
   //(after ReadData(slice))
@@ -303,7 +325,7 @@ void AliL3Hough::FindTrackCandidates()
   
   Int_t n_patches;
   if(fAddHistograms)
-    n_patches = 1; //Histograms has been added.
+    n_patches = 1; //Histograms have been added.
   else
     n_patches = fNPatches;
 
@@ -320,15 +342,15 @@ void AliL3Hough::FindTrackCandidates()
          fPeakFinder->Reset();
          fPeakFinder->SetHistogram(hist);
          fPeakFinder->FindMaxima(0,0); //Simple maxima finder
-         //fPeakFinder->FindAbsMaxima();
          cout<<"Found "<<fPeakFinder->GetEntries()<<endl;
+
          for(Int_t k=0; k<fPeakFinder->GetEntries(); k++)
            {
              if(fPeakFinder->GetWeight(k) == 0) continue;
              AliL3HoughTrack *track = (AliL3HoughTrack*)fTracks[i]->NextTrack();
              track->SetTrackParameters(fPeakFinder->GetXPeak(k),fPeakFinder->GetYPeak(k),fPeakFinder->GetWeight(k));
              track->SetEtaIndex(j);
-             track->SetEta((Double_t)(j*eta_slice));
+             track->SetEta((Double_t)((j+0.5)*eta_slice));
              track->SetRowRange(AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5));
            }
        }
@@ -363,6 +385,7 @@ void AliL3Hough::Evaluate(Int_t road_width)
   AliL3TrackArray *tracks;
   for(Int_t i=0; i<fNPatches; i++)
     {
+      fEval[i]->InitTransformer(fHoughTransformer[i]);
       fEval[i]->SetNumOfRowsToMiss(2);
       fEval[i]->SetNumOfPadsToLook(road_width);
       if(fAddHistograms)
index 1c034a3..e6c7354 100644 (file)
@@ -25,6 +25,8 @@ class AliL3Hough {
   Bool_t fUse8bits;
   Int_t fNEtaSegments;
   Int_t fNPatches;
+  Int_t fversion;
+
   AliL3MemHandler **fMemHandler; //!
   AliL3HoughBaseTransformer **fHoughTransformer; //!
   AliL3HoughEval **fEval; //!
@@ -39,6 +41,7 @@ class AliL3Hough {
  public:
   
   AliL3Hough(); 
+  AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments=100,Int_t tversion=0);
   virtual ~AliL3Hough();
   
   void Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8=kFALSE);
@@ -74,7 +77,6 @@ class AliL3Hough {
   AliL3HoughMaxFinder *GetMaxFinder() {return fPeakFinder;}
 
   ClassDef(AliL3Hough,1) //Hough transform base class
-
 };
 
 #endif
index a26cc42..0053e04 100644 (file)
@@ -145,7 +145,7 @@ Bool_t AliL3HoughEval::LookInsideRoad(AliL3HoughTrack *track,Int_t eta_index,Boo
              if(pixel_index != eta_index) break;
              total_charge += digPt[j].fCharge;
              if(remove)
-               digPt[j].fCharge = 0; //Erease the track from image
+               digPt[j].fCharge = 0; //Erase the track from image
              npixs++;
            }
        }
index 0a1140b..f7c0ccc 100644 (file)
@@ -6,6 +6,7 @@
 
 #pragma link C++ class AliL3Hough; 
 #pragma link C++ class AliL3HoughTransformer; 
+#pragma link C++ class AliL3HoughTransformerVhdl; 
 #pragma link C++ class AliL3HoughBaseTransformer; 
 #pragma link C++ class AliL3HoughTrack; 
 #pragma link C++ class AliL3HoughMaxFinder;
index ba6c1a4..2978901 100644 (file)
@@ -41,7 +41,7 @@ void AliL3HoughTrack::Set(AliL3Track *track)
 {
   
   AliL3HoughTrack *tpt = (AliL3HoughTrack*)track;
-  SetTrackParameters(tpt->GetKappa(),tpt->GetPhi0(),tpt->GetWeight());
+  SetTrackParameters(tpt->GetKappa(),tpt->GetPsi(),tpt->GetWeight());
   SetEtaIndex(tpt->GetEtaIndex());
   SetEta(tpt->GetEta());
   SetTgl(tpt->GetTgl());
@@ -129,7 +129,7 @@ void AliL3HoughTrack::UpdateToFirstRow()
   }
   
   //Find intersection
-  Double_t fac2   = ( radius*radius + fac1 - rc*rc) / (2.00 * radius * sfac ) ;
+  Double_t fac2   = (radius*radius + fac1 - rc*rc) / (2.00 * radius * sfac ) ;
   Double_t phi    = atan2(yc,xc) + GetCharge()*acos(fac2) ;
   Double_t td     = atan2(radius*sin(phi) - yc,radius*cos(phi) - xc) ;
   
@@ -168,7 +168,7 @@ void AliL3HoughTrack::UpdateToFirstRow()
   //printf("last point %f %f %f\n",xyz[0],xyz[1],xyz[2]);
 }
 
-void AliL3HoughTrack::SetTrackParameters(Double_t kappa,Double_t phi,Int_t weight)
+void AliL3HoughTrack::SetTrackParameters(Double_t kappa,Double_t eangle,Int_t weight)
 {
 
   fWeight = weight;
@@ -179,11 +179,11 @@ void AliL3HoughTrack::SetTrackParameters(Double_t kappa,Double_t phi,Int_t weigh
   Double_t radius = 1/fabs(kappa);
   SetRadius(radius);
   SetFirstPoint(0,0,0);
-  SetPsi(phi); 
+  SetPsi(eangle); //Psi = emission angle when first point is vertex
+  SetPhi0(0);     //not defined for vertex reference point
   SetR0(0);
   Double_t charge = -1.*kappa;
   SetCharge((Int_t)copysign(1.,charge));
-  
   Double_t trackPhi0 = GetPsi() + charge*0.5*AliL3Transform::Pi()/fabs(charge);
   Double_t xc = GetFirstPointX() - GetRadius() * cos(trackPhi0) ;
   Double_t yc = GetFirstPointY() - GetRadius() * sin(trackPhi0) ;
@@ -191,8 +191,6 @@ void AliL3HoughTrack::SetTrackParameters(Double_t kappa,Double_t phi,Int_t weigh
   SetCenterY(yc);
   SetNHits(1); //just for the trackarray IO
   fIsHelix = true;
-  
-  
 }
 
 void AliL3HoughTrack::SetLineParameters(Double_t psi,Double_t D,Int_t weight,Int_t *rowrange,Int_t ref_row)
index b687308..5a9a912 100644 (file)
@@ -27,7 +27,7 @@ class AliL3HoughTrack : public AliL3Track {
   virtual void CalculateHelix();
   
   void UpdateToFirstRow();
-  void SetTrackParameters(Double_t kappa,Double_t phi,Int_t weight);  
+  void SetTrackParameters(Double_t kappa,Double_t eangle,Int_t weight);  
   void SetLineParameters(Double_t psi,Double_t D,Int_t weight,Int_t *rowrange,Int_t ref_row);
 
   Int_t GetWeight()  const {return fWeight;}
diff --git a/HLT/hough/AliL3HoughTransformerVhdl.cxx b/HLT/hough/AliL3HoughTransformerVhdl.cxx
new file mode 100644 (file)
index 0000000..801e3ac
--- /dev/null
@@ -0,0 +1,333 @@
+//$Id$
+
+// Author: Constantin Loizides <mailto:loizides@fi.uib.no>
+//*-- Copyright &copy CL and ASV
+
+#include "AliL3MemHandler.h"
+#include "AliL3Logging.h"
+#include "AliL3Transform.h"
+#include "AliL3DigitData.h"
+#include "AliL3Histogram.h"
+#include "AliL3HoughTransformerVhdl.h"
+
+//_____________________________________________________________
+// AliL3HoughTransformerVhdl
+//
+// Hough transformation class
+//
+
+ClassImp(AliL3HoughTransformerVhdl)
+
+AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl()
+{
+  //Default constructor
+  fParamSpace = 0;
+}
+
+AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl(Int_t slice,Int_t patch,Int_t n_eta_segments) : AliL3HoughBaseTransformer(slice,patch,n_eta_segments)
+{
+  //Normal constructor
+  fParamSpace = 0;
+}
+
+AliL3HoughTransformerVhdl::~AliL3HoughTransformerVhdl()
+{
+  DeleteHistograms();
+}
+
+void AliL3HoughTransformerVhdl::DeleteHistograms()
+{
+  if(!fParamSpace)
+    return;
+  for(Int_t i=0; i<GetNEtaSegments(); i++)
+    {
+      if(!fParamSpace[i]) continue;
+      delete fParamSpace[i];
+    }
+  delete [] fParamSpace;
+}
+
+void AliL3HoughTransformerVhdl::CreateHistograms(Int_t nxbin,Double_t pt_min,
+                                            Int_t nybin,Double_t phimin,Double_t phimax)
+{
+  //Create the histograms (parameter space).
+  //These are 2D histograms, span by kappa (curvature of track) and phi0 (emission angle with x-axis).
+  //The arguments give the range and binning; 
+  //nxbin = #bins in kappa
+  //nybin = #bins in phi0
+  //pt_min = mimium Pt of track (corresponding to maximum kappa)
+  //phi_min = mimimum phi0 (degrees)
+  //phi_max = maximum phi0 (degrees)
+    
+  Double_t bfact = 0.0029980;
+  Double_t bfield = AliL3Transform::GetBField();
+  Double_t x = bfact*bfield/pt_min;
+  Double_t torad = AliL3Transform::Pi()/180;
+  CreateHistograms(nxbin,-1.*x,x,nybin,phimin*torad,phimax*torad);
+}
+
+void AliL3HoughTransformerVhdl::CreateHistograms(Int_t nxbin,Double_t xmin,Double_t xmax,
+                                            Int_t nybin,Double_t ymin,Double_t ymax)
+{
+  
+  fParamSpace = new AliL3Histogram*[GetNEtaSegments()];
+  
+  Char_t histname[256];
+  for(Int_t i=0; i<GetNEtaSegments(); i++)
+    {
+      sprintf(histname,"paramspace_%d",i);
+      fParamSpace[i] = new AliL3Histogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
+    }
+}
+
+void AliL3HoughTransformerVhdl::Reset()
+{
+  //Reset all the histograms
+
+  if(!fParamSpace)
+    {
+      LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::Reset","Histograms")
+       <<"No histograms to reset"<<ENDLOG;
+      return;
+    }
+  
+  for(Int_t i=0; i<GetNEtaSegments(); i++)
+    fParamSpace[i]->Reset();
+}
+
+
+Int_t AliL3HoughTransformerVhdl::GetEtaIndex(Double_t eta)
+{
+  //Return the histogram index of the corresponding eta. 
+
+  Double_t etaslice = (GetEtaMax() - GetEtaMin())/GetNEtaSegments();
+  Double_t index = (eta-GetEtaMin())/etaslice;
+  return (Int_t)index;
+}
+
+void AliL3HoughTransformerVhdl::TransformCircle()
+{
+  //Transform the input data with a circle HT.
+  //The function loops over all the data, and transforms each pixel with the equations:
+  // 
+  //kappa = 2/R*sin(phi - phi0)
+  //
+  //where R = sqrt(x*x +y*y), and phi = arctan(y/x)
+  //
+  //Each pixel then transforms into a curve in the (kappa,phi0)-space. In order to find
+  //which histogram in which the pixel should be transformed, the eta-value is calcluated
+  //and the proper histogram index is found by GetEtaIndex(eta).
+
+
+  AliL3DigitRowData *tempPt = GetDataPointer();
+  if(!tempPt)
+    {
+      LOG(AliL3Log::kError,"AliL3HoughTransformerVhdl::TransformCircle","Data")
+       <<"No input data "<<ENDLOG;
+      return;
+    }
+  
+  //Loop over the padrows:
+  for(Int_t i=AliL3Transform::GetFirstRow(GetPatch()); i<=AliL3Transform::GetLastRow(GetPatch()); i++)
+    {
+      //Get the data on this padrow:
+      AliL3DigitData *digPt = tempPt->fDigitData;
+      if(i != (Int_t)tempPt->fRow)
+       {
+         printf("AliL3HoughTransform::TransformCircle : Mismatching padrow numbering\n");
+         continue;
+       }
+      
+      //Loop over the data on this padrow:
+      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 <= GetThreshold())
+           continue;
+         Int_t sector,row;
+         Float_t xyz[3];
+         
+         //Transform data to local cartesian coordinates:
+         AliL3Transform::Slice2Sector(GetSlice(),i,sector,row);
+         AliL3Transform::Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
+         
+         //Calculate the eta:
+         Double_t eta = AliL3Transform::GetEta(xyz);
+         
+         //Get the corresponding index, which determines which histogram to fill:
+         Int_t eta_index = GetEtaIndex(eta);
+         if(eta_index < 0 || eta_index >= GetNEtaSegments())
+           continue;
+         
+         //Get the correct histogrampointer:
+         AliL3Histogram *hist = fParamSpace[eta_index];
+         if(!hist)
+           {
+             printf("AliL3HoughTransformerVhdl::TransformCircle : Error getting histogram in index %d\n",eta_index);
+             continue;
+           }
+
+         //Do the transformation:
+         Float_t R = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]); 
+         Float_t phi = AliL3Transform::GetPhi(xyz);
+         
+         //Fill the histogram along the phirange
+         for(Int_t b=hist->GetFirstYbin(); b<=hist->GetLastYbin(); b++)
+           {
+             Float_t phi0 = hist->GetBinCenterY(b);
+             Float_t kappa = 2*sin(phi - phi0)/R;
+             hist->Fill(kappa,phi0,charge);
+           }
+       }
+      
+      //Move the data pointer to the next padrow:
+      AliL3MemHandler::UpdateRowPointer(tempPt);
+    }
+}
+
+void AliL3HoughTransformerVhdl::TransformCircleC(Int_t row_range)
+{
+  //Circle transform, using combinations of every 2 points lying
+  //on different padrows and within the same etaslice.
+  
+  AliL3DigitRowData *tempPt = GetDataPointer();
+  if(!tempPt)
+    LOG(AliL3Log::kError,"AliL3HoughTransformerVhdl::TransformCircleC","Data")
+      <<"No input data "<<ENDLOG;
+  Int_t counter=0;
+  for(Int_t i=AliL3Transform::GetFirstRow(GetPatch()); i<=AliL3Transform::GetLastRow(GetPatch()); i++)
+    {
+      counter += tempPt->fNDigit;
+      AliL3MemHandler::UpdateRowPointer(tempPt);
+    }
+  
+  struct Digit {
+    Int_t row;
+    Double_t r;
+    Double_t phi;
+    Int_t eta_index;
+    Int_t charge;
+  };
+  
+  Digit *digits = new Digit[counter];
+  cout<<"Allocating "<<counter*sizeof(Digit)<<" bytes to digitsarray"<<endl;
+  
+  Int_t total_digits=counter;
+  Int_t sector,row,tot_charge,pad,time,charge;
+  Double_t r1,r2,phi1,phi2,eta,kappa,phi_0;
+  Float_t xyz[3];
+  
+  counter=0;
+  tempPt = GetDataPointer();
+  
+  for(Int_t i=AliL3Transform::GetFirstRow(GetPatch()); i<=AliL3Transform::GetLastRow(GetPatch()); i++)
+    {
+      AliL3DigitData *digPt = tempPt->fDigitData;
+      for(UInt_t di=0; di<tempPt->fNDigit; di++)
+       {
+         charge = digPt[di].fCharge;
+         pad = digPt[di].fPad;
+         time = digPt[di].fTime;
+         AliL3Transform::Slice2Sector(GetSlice(),i,sector,row);
+         AliL3Transform::Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
+         eta = AliL3Transform::GetEta(xyz);
+         digits[counter].row = i;
+         digits[counter].r = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]);
+         digits[counter].phi = atan2(xyz[1],xyz[0]);
+         digits[counter].eta_index = GetEtaIndex(eta);
+         digits[counter].charge = charge;
+         counter++;
+       }
+      AliL3MemHandler::UpdateRowPointer(tempPt);
+    }
+  
+  for(Int_t i=0; i<total_digits; i++)
+    {
+      if(digits[i].eta_index < 0 || digits[i].eta_index >= GetNEtaSegments()) continue;
+      Int_t ind = digits[i].eta_index;
+      
+      for(Int_t j=i+1; j<total_digits; j++)
+       {
+         if(digits[i].row == digits[j].row) continue;
+         if(digits[i].eta_index != digits[j].eta_index) continue;
+         if(digits[i].row + row_range < digits[j].row) break;
+         
+         //Get the correct histogrampointer:
+         AliL3Histogram *hist = fParamSpace[ind];
+         if(!hist)
+           {
+             printf("AliL3HoughTransformerVhdl::TransformCircleC() : No histogram at index %d\n",ind);
+             continue;
+           }
+         
+         r1 = digits[i].r;
+         phi1 = digits[i].phi;
+         r2 = digits[j].r;
+         phi2 = digits[j].phi;
+         phi_0 = atan( (r2*sin(phi1)-r1*sin(phi2))/(r2*cos(phi1)-r1*cos(phi2)) );
+         kappa = 2*sin(phi2-phi_0)/r2;
+         tot_charge = digits[i].charge + digits[j].charge;
+         hist->Fill(kappa,phi_0,tot_charge);
+       }
+    }
+  delete [] digits;
+}
+
+void AliL3HoughTransformerVhdl::TransformLine()
+{
+  //Do a line transform on the data.
+
+  
+  AliL3DigitRowData *tempPt = GetDataPointer();
+  if(!tempPt)
+    {
+      LOG(AliL3Log::kError,"AliL3HoughTransformerVhdl::TransformLine","Data")
+       <<"No input data "<<ENDLOG;
+      return;
+    }
+    
+  for(Int_t i=AliL3Transform::GetFirstRow(GetPatch()); i<=AliL3Transform::GetLastRow(GetPatch()); 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 < GetThreshold())
+           continue;
+         Int_t sector,row;
+         Float_t xyz[3];
+         AliL3Transform::Slice2Sector(GetSlice(),i,sector,row);
+         AliL3Transform::Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
+         Float_t eta = AliL3Transform::GetEta(xyz);
+         Int_t eta_index = GetEtaIndex(eta);//(Int_t)(eta/etaslice);
+         if(eta_index < 0 || eta_index >= GetNEtaSegments())
+           continue;
+         
+         //Get the correct histogram:
+         AliL3Histogram *hist = fParamSpace[eta_index];
+         if(!hist)
+           {
+             printf("AliL3HoughTransformerVhdl::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);
+           }
+       }
+      AliL3MemHandler::UpdateRowPointer(tempPt);
+    }
+  
+}
diff --git a/HLT/hough/AliL3HoughTransformerVhdl.h b/HLT/hough/AliL3HoughTransformerVhdl.h
new file mode 100644 (file)
index 0000000..db4d026
--- /dev/null
@@ -0,0 +1,50 @@
+#ifndef ALIL3_HOUGHTRANSFORMERVHDL
+#define ALIL3_HOUGHTRANSFORMERVDHL
+
+#include "AliL3RootTypes.h"
+#include "AliL3HoughBaseTransformer.h"
+
+class AliL3Histogram;
+
+class AliL3HoughTransformerVhdl : public AliL3HoughBaseTransformer {
+  
+ private:
+  
+  AliL3Histogram **fParamSpace; //!
+  void DeleteHistograms();
+
+ public:
+  AliL3HoughTransformerVhdl(); 
+  AliL3HoughTransformerVhdl(Int_t slice,Int_t patch,Int_t n_eta_segments);
+  virtual ~AliL3HoughTransformerVhdl();
+  
+  void CreateHistograms(Int_t nxbin,Double_t ptmin,Int_t nybin,Double_t phimin,Double_t phimax);
+  void CreateHistograms(Int_t nxbin,Double_t xmin,Double_t xmax,
+                       Int_t nybin,Double_t ymin,Double_t ymax);
+  void Reset();
+  void TransformCircle();
+  void TransformCircleC(Int_t row_range);
+  void TransformLine();
+
+  Int_t GetEtaIndex(Double_t eta);
+  AliL3Histogram *GetHistogram(Int_t eta_index);
+  
+
+  ClassDef(AliL3HoughTransformerVhdl,1) //Normal Hough transformation class
+
+};
+
+inline AliL3Histogram *AliL3HoughTransformerVhdl::GetHistogram(Int_t eta_index)
+{
+  if(!fParamSpace || eta_index >= GetNEtaSegments() || eta_index < 0)
+    return 0;
+  if(!fParamSpace[eta_index])
+    return 0;
+  return fParamSpace[eta_index];
+}
+
+#endif
+
+
+
+
index 49ad596..c343655 100644 (file)
@@ -40,6 +40,7 @@ endif
 
 #Save the particle id's
 #DEFSTR += -Ddo_mc
+DEFSTR += $(EXTRADEF)
 
 #Use logging classes
 ifndef NOLOGGING
@@ -55,7 +56,8 @@ endif
 SRCS   = AliL3HoughTransformer.cxx AliL3Hough.cxx AliL3HoughTrack.cxx\
          AliL3HoughMaxFinder.cxx AliL3HoughEval.cxx AliL3HoughMerger.cxx \
          AliL3Histogram.cxx AliL3Histogram1D.cxx AliL3HoughBaseTransformer.cxx \
-         AliL3HoughIntMerger.cxx AliL3HoughGlobalMerger.cxx 
+         AliL3HoughIntMerger.cxx AliL3HoughGlobalMerger.cxx \
+         AliL3HoughTransformerVhdl.cxx
 
 DICT = AliL3HoughCint.cxx
 DICTH = AliL3HoughCint.h
@@ -90,4 +92,4 @@ clean :
        rm -f $(DICT) $(DICTH) 
 
 so: 
-       rm -f $(LIBDIR)/libAliL3Hough.so
\ No newline at end of file
+       rm -f $(LIBDIR)/libAliL3Hough.so