[u/mrichter/AliRoot.git] / HLT / hough / AliL3HoughTransformerVhdl.cxx

//$Id$

// Author: Constantin Loizides <mailto:loizides@fi.uib.no>
//*-- Copyright&Copy CL

#include "AliL3StandardIncludes.h"

#include "AliL3Logging.h"
#include "AliL3MemHandler.h"
#include "AliL3Transform.h"
#include "AliL3DigitData.h"
#include "AliL3Histogram.h"

//try to be close to VHDL version (eg. LUTs)
//VESTBO: if you switch that of, you should get your version
//of the transformer!!!
#define VHDLVERSION

//dont use dynamic objects on heap
//VESTBO: this a desperately switched on in order 
//to find the bug, which i thought hat something
//to do with a wrong heap assignment! But the
//data still gets overwritten!
#define VHDLSTATIC

#include "AliL3HoughTransformerVhdl.h"

#if GCCVERSION == 3
using namespace std;
#endif

/** \class AliL3HoughTransformerVhdl
// <pre>
//_____________________________________________________________
// AliL3HoughTransformerVhdl
//
// Hough transformation class for VHDL comparism.
//
//</pre>
*/

ClassImp(AliL3HoughTransformerVhdl)

#ifdef VHDLVERSION

AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl()
  : AliL3HoughBaseTransformer()//, /*fMinRow(0),*/ fMaxRow(0)

{
  cout << "hallo" << endl;
  fParamSpace=0;

#ifndef VHDLSTATIC
  fLUTX=0;
  fLUTY=0;
  fLUTEta=0;
  fLUTphi0=0;
  fLUT2sinphi0=0;
  fLUT2cosphi0=0;
#endif
  
  fMinRow=0;
  fMaxRow=0;

  /*
  fNRows=0;
  fNEtas=0;
  fNPhi0=0;
  fSector=0;
  fSectorRow=0;
  fZSign=0;
  fZLengthPlusOff=0;
  fTimeWidth=0;
  fPadPitch=0;
  fEtaSlice=0;
  */
}

AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl(Int_t slice,Int_t patch,Int_t n_eta_segments) 
  : AliL3HoughBaseTransformer(slice,patch,n_eta_segments)
{
  AliL3HoughTransformerVhdl();
  cout << "wiesow hier??" << endl;
  Init(slice,patch,n_eta_segments);
}

AliL3HoughTransformerVhdl::~AliL3HoughTransformerVhdl()
{
  DeleteHistograms();

#ifndef VHDLSTATIC
  if(fNRows){
    delete[] fLUTX;
    delete[] fLUTY;
  }
  if(fNEtas) delete[] fLUTEta;
#endif
}

void AliL3HoughTransformerVhdl::Init(Int_t slice,Int_t patch,Int_t n_eta_segments)
{
  cout << "InitVhdl " << slice << " " << patch << " " << n_eta_segments << endl;
  AliL3HoughBaseTransformer::Init(slice,patch,n_eta_segments);

#ifndef VHDLSTATIC
  //delete old LUT tables
  if(fNRows){
    fNRows=0;
    delete[] fLUTX;
    delete[] fLUTY;
  }
  if(fNEtas){
    delete[] fLUTEta;
    fNEtas=0;
  }
#endif

  //set when histogram is filled
  //fLUTphi0=0;
  //fLUT2sinphi0=0;
  //fLUT2cosphi0=0;

  Int_t minrow_=AliL3Transform::GetFirstRow(patch);
  Int_t maxrow_=AliL3Transform::GetLastRow(patch);
  Int_t n_=AliL3Transform::GetNRows(patch);
  Int_t sector_=0,sectorrow_=0;
  AliL3Transform::Slice2Sector(slice,minrow_,sector_,sectorrow_);
  Float_t padpitch_=0.;
  if(sector_<AliL3Transform::GetNSectorLow())
    padpitch_=AliL3Transform::GetPadPitchWidthLow();
  else
    padpitch_=AliL3Transform::GetPadPitchWidthUp();  

  Float_t etamax_=GetEtaMax();
  Float_t etamin_=GetEtaMin();
  Float_t etaslice_=(etamax_-etamin_)/n_eta_segments;
#if 0
  //lookup tables for X and Y
#ifndef VHDLSTATIC
  fLUTX=new AliL3FFloat[n_];
  fLUTY=new AliL3FFloat[n_]; 
#endif
  for(Int_t rr=0;rr<n_;rr++){
    fLUTX[rr]=Float_t(AliL3Transform::Row2X(rr+minrow_));
    fLUTY[rr]=Float_t(0.5*(AliL3Transform::GetNPads(rr+minrow_)-1)*padpitch_);
    //VESTBO: uncomment to see values (and compare with Print function)
    //cout << rr << ": " << (Float_t)fLUTX[rr] << " " << (Float_t)fLUTY[rr] << endl;
  }

  //lookup tables for rz2s <=> etas
#ifndef VHDLSTATIC
  fLUTEta=new AliL3FFloat[n_eta_segments];
#endif
  for(Int_t rr=0;rr<n_eta_segments;rr++){
    fLUTEta[rr]=CalcRoverZ2(etamin_+(rr+1)*etaslice_);
    //VESTBO: uncomment to see values (and compare with Print function)
    //cout << rr << ": " << fLUTEta[rr] << endl;
  }
#endif /*uncomment arrays*/

  //member values
  fMinRow=minrow_;
  fMaxRow=maxrow_;
  fNRows=n_;
  fNEtas=n_eta_segments;
  fNPhi0=0;
  fSector=sector_;
  fSectorRow=sectorrow_;
  fZSign = slice < 18 ? 1:-1;
  fZLengthPlusOff=AliL3Transform::GetZLength()+AliL3Transform::GetZOffset();
  fTimeWidth=AliL3Transform::GetZWidth();
  fPadPitch=padpitch_;
  fEtaSlice=etaslice_;
}

#if 0
Float_t AliL3HoughTransformerVhdl::CalcRoverZ2(Float_t eta)
{
  Float_t e=exp(2*eta);
  Float_t ret=(e+1)/(e-1);
  ret*=ret;
  return ret;
}

Float_t AliL3HoughTransformerVhdl::CalcEta(Float_t roverz2)
{
  Float_t rz=sqrt(roverz2);
  if(fZSign>0) rz=-rz;
  Float_t ret=(1+rz)/(rz-1);
  ret=0.5*log(ret);
  return ret;
}

Float_t AliL3HoughTransformerVhdl::CalcX(Int_t row)
{
  return fLUTX[row];
}

Float_t AliL3HoughTransformerVhdl::CalcY(Int_t pad,Int_t row)
{
  return pad*fPadPitch-fLUTY[row];
}

Float_t AliL3HoughTransformerVhdl::CalcZ(Int_t time)
{
  Float_t ret=time*fTimeWidth;
  if(fZSign>0) ret=fZLengthPlusOff-ret;
  else ret=ret-fZLengthPlusOff;
  return ret;
}
#endif

void AliL3HoughTransformerVhdl::DeleteHistograms()
{
  if(!fParamSpace) return;
  for(Int_t i=0; i<GetNEtaSegments(); i++){
    if(!fParamSpace[i]) continue;
    delete fParamSpace[i];
  }
  delete[] fParamSpace;

#ifndef VHDLSTATIC
  if(fNPhi0){
    fNPhi0=0;
    delete[] fLUT2sinphi0;
    delete[] fLUT2cosphi0;
    delete[] fLUTphi0;
  }
#endif
}

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)
    
  AliL3FFloat x = AliL3Transform::GetBFieldValue()/pt_min;
  AliL3FFloat torad = AliL3Transform::Pi()/180;
  AliL3FFloat phimin_=phimin*torad;
  AliL3FFloat phimax_=phimax*torad;

  CreateHistograms(nxbin,-x,x,nybin,phimin_,phimax_);
}

void AliL3HoughTransformerVhdl::CreateHistograms(Int_t nxbin,Double_t xmin,Double_t xmax,
					         Int_t nybin,Double_t ymin,Double_t ymax)
{

  AliL3FFloat xminf=xmin;
  AliL3FFloat xmaxf=xmax;
  AliL3FFloat yminf=ymin;
  AliL3FFloat ymaxf=ymax;
  /*
  cout << xminf << endl;
  cout << xmaxf << " " << ((xmaxf-xminf)/nxbin) << endl;
  cout << yminf << endl;
  cout << ymaxf << " " << ((ymaxf-yminf)/nybin) << endl;
  */


  fParamSpace = new AliL3Histogram*[GetNEtaSegments()];

  Char_t histname[256];
  for(Int_t i=0; i<GetNEtaSegments(); i++)
    {
      sprintf(histname,"paramspace_vhdl_%d",i);
      fParamSpace[i] = new AliL3Histogram(histname,"id",nxbin,xminf,xmaxf,nybin,yminf,ymaxf);
    }
  
  //create lookup table for sin and cos
  fNPhi0=nybin+1;
#if 0
#ifndef VHDLSTATIC
  fLUTphi0=new AliL3FFloat[fNPhi0];
  fLUT2sinphi0=new AliL3FFloat[fNPhi0];
  fLUT2cosphi0=new AliL3FFloat[fNPhi0];
#endif
  AliL3FFloat diff=(ymax-ymin)/nybin;
  AliL3FFloat phi0=ymin-0.5*diff;
  for(Int_t i=0; i<fNPhi0; i++){
    phi0+=diff;
    fLUTphi0[i]=phi0;
    fLUT2sinphi0[i]=Float_t(2*sin(phi0));
    fLUT2cosphi0[i]=Float_t(2*cos(phi0));
    //VESTBO: uncomment to see values (and compare with Print function)
    //cout << i << ": " << fLUTphi0[i] << " " << fLUT2sinphi0[i] << " " << fLUT2cosphi0[i] << endl;
  }  
#endif /*arrays*/
}

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();
}

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];
}

#if 0
Int_t AliL3HoughTransformerVhdl::FindIndex(Double_t rz2)
{
  Int_t index=0;
  while((index<fNEtas)&&(rz2<=fLUTEta[index])){
    index++;
    //cout << index << ": " << rz2 << " " << fLUTEta[index] << endl;
  }
  return index;
}
#endif

Int_t AliL3HoughTransformerVhdl::GetEtaIndex(Double_t eta)
{
  //AliL3FFloat rz2=CalcRoverZ2(eta);
  //return FindIndex(rz2);
  return 0;
}

Double_t AliL3HoughTransformerVhdl::GetEta(Int_t eta_index,Int_t slice){
  if(eta_index >= fNEtas || eta_index < 0){
    LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::GetEta","Index")
      <<"Index out of range."<<ENDLOG;
    return 0.;
  }
  return 0;//(CalcEta(fLUTEta[eta_index])-0.5*fEtaSlice);
}

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 = 1/x * ( y/x*2*cos(phi0) - 2*sin(phi0) )

#if 0
  AliL3DigitRowData *tempPt = GetDataPointer();
  if(!tempPt)
    {
      LOG(AliL3Log::kError,"AliL3HoughTransformerVhdl::TransformCircle","Data")
	<<"No input data "<<ENDLOG;
      return;
    }

  AliL3FFloat x=0.,y=0.,z=0.,rz2=0.;
  //Loop over the padrows:
   for(Int_t i=fMinRow, row=0; i<=fMaxRow; i++, row++){

    //Get the data on this padrow:
    AliL3DigitData *digPt = tempPt->fDigitData;
    if(i != (Int_t)tempPt->fRow){
      LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Padrows")
	<<"Mismatching padrow numbering "<<i<<"!="<<(Int_t)tempPt->fRow<<ENDLOG;
      continue;
    }
      
    //Loop over the data on this padrow:
    for(UInt_t j=0; j<tempPt->fNDigit; j++){
      Int_t charge = digPt[j].fCharge;
      if((Int_t)charge <= GetLowerThreshold()) continue;
      
      Int_t pad = digPt[j].fPad;
      Int_t time = digPt[j].fTime;

      x=CalcX(row);
      y=CalcY(pad,row);
      z=CalcZ(time);

      //find eta slice
      rz2=1+(x*x+y*y)/z/z;
      Int_t eta_index = FindIndex(rz2);
      if(eta_index < 0 || eta_index >= GetNEtaSegments()){
	//LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"No histograms corresponding to eta index value of "<<eta_index<<"."<<ENDLOG;
	continue;
      }	  
      //Get the correct histogrampointer:
      AliL3Histogram *hist = fParamSpace[eta_index];
      if(!hist){
	//LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"Error getting histogram in index "<<eta_index<<"."<<ENDLOG;
	continue;
      }

      //Fill the histogram along the phirange
      AliL3FFloat kappa,ydx=y/x;
      for(Int_t b=0; b<fNPhi0; b++){
	kappa=ydx*fLUT2cosphi0[b]-fLUT2sinphi0[b];
	kappa/=x;
	hist->Fill(kappa,fLUTphi0[b],charge);
	//cout << kappa << " " << fLUTphi0[b] << " " << charge << endl;
      }

    }
    //Move the data pointer to the next padrow:
    AliL3MemHandler::UpdateRowPointer(tempPt);
  }
#endif 
}

void AliL3HoughTransformerVhdl::Print()
{
  cout << "fSlice: " << GetSlice() << endl;
  cout << "fPatch: " << GetPatch() << endl;
  cout << "fSector: " << fSector << endl;
  cout << "fSectorRow: " << fSectorRow << endl;
  cout << "fMinRow: " << fMinRow << endl;
  cout << "fMaxRow: " << fMaxRow << endl;
  cout << "fNRows: " << fNRows << endl;
  cout << "fNEtas: " << fNEtas << endl;
  cout << "fNPhi0: " << fNPhi0 << endl;
  cout << "fZSign: " << fZSign << endl;
  cout << "fZLengthPlusOff: " << fZLengthPlusOff << endl;
  cout << "fPadPitch: " << fPadPitch << endl;
  cout << "fTimeWidth: " << fTimeWidth << endl;
#if 0
  if(!fNRows) return;
  cout << "fLUTX " << fNRows << endl;
  for(Int_t i=0;i<fNRows;i++) cout << "fLUTX[" << i << "]=" << (Float_t)fLUTX[i] << endl;
  cout << "fLUTY " << fNRows << endl;
  for(Int_t i=0;i<fNRows;i++) cout << "fLUTY[" << i << "]=" << fLUTY[i] << endl;
  if(!fNEtas) return;
  cout << "fLUTEta " << fNEtas << endl;
  for(Int_t i=0;i<fNEtas;i++) cout << "fLUTEta[" << i << "]=" << fLUTEta[i] << endl;
  if(!fNPhi0) return;
  cout << "fLUTphi0 " << fNPhi0 << endl;
  for(Int_t i=0;i<fNPhi0;i++) cout << "fLUTPhi0[" << i << "]=" << fLUTphi0[i] << endl;
  cout << "fLUT2sinphi0 " << fNPhi0 << endl;
  for(Int_t i=0;i<fNPhi0;i++) cout << "fLUT2sinphi0[" << i << "]=" << fLUT2sinphi0[i] << endl;
  cout << "fLUT2cosphi0 " << fNPhi0 << endl;
  for(Int_t i=0;i<fNPhi0;i++) cout << "fLUT2cosphi0[" << i << "]=" << fLUT2cosphi0[i] << endl;
#endif 
}

//end vhdl version
#else 
//standard version

AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl()
{
  fParamSpace=0;
}

AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl(Int_t slice,Int_t patch,Int_t n_eta_segments) 
            : AliL3HoughBaseTransformer(slice,patch,n_eta_segments)
{
  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 x = AliL3Transform::GetBFieldValue()/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)
{

  AliL3FFloat xminf(xmin);
  AliL3FFloat xmaxf(xmax);
  AliL3FFloat yminf(ymin);
  AliL3FFloat ymaxf(ymax);

  cout << xminf << endl;
  cout << xmaxf << " " << ((xmaxf-xminf)/nxbin) << endl;
  cout << yminf << endl;
  cout << ymaxf << " " << ((ymaxf-yminf)/nybin) << endl;

  
  fParamSpace = new AliL3Histogram*[GetNEtaSegments()];
  
  Char_t histname[256];
  for(Int_t i=0; i<GetNEtaSegments(); i++)
    {
      sprintf(histname,"paramspace_vhdl_%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();
}

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];
}

Int_t AliL3HoughTransformerVhdl::GetEtaIndex(Double_t eta)
{
  //Return the histogram index of the corresponding eta. 
  AliL3FFloat etaslice = (GetEtaMax() - GetEtaMin())/GetNEtaSegments();
  AliL3FFloat index = (eta-GetEtaMin())/etaslice;

  //if((Int_t)index==0) cout << index << " " << (Int_t) index << endl;
  return (Int_t)index;
}

Double_t AliL3HoughTransformerVhdl::GetEta(Int_t eta_index)
{
  if(eta_index >= GetNEtaSegments() || eta_index < 0){
    return 0;
  }
  Double_t eta_slice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments();
  Double_t eta=(Double_t)((eta_index+0.5)*eta_slice);
  if(GetSlice()>17) eta*=-1;
  return eta;
}

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)
	{
      LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Padrows")
	<<"Mismatching padrow numbering."<<ENDLOG;
	  continue;
	}
      

      //Loop over the data on this padrow:
      for(UInt_t j=0; j<tempPt->fNDigit; j++)
	{
	  UShort_t charge = digPt[j].fCharge;
	  if((Int_t)charge <= GetLowerThreshold()) continue;

	  UChar_t pad = digPt[j].fPad;
	  UShort_t time = digPt[j].fTime;

	  Int_t sector,row;
	  Float_t xyz[3];
	  AliL3FFloat fxyz[3];

	  //Transform data to local cartesian coordinates:
	  AliL3Transform::Slice2Sector(GetSlice(),i,sector,row);
	  AliL3Transform::Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);

	  //trunc to fixed format.
	  for(int i=0;i<3;i++) {
	    fxyz[i]=xyz[i];
	    //cout << fxyz[i];
	    xyz[i]=fxyz[i];
	  }
	  //cout << endl;

	  //Calculate the eta:
	  AliL3FFloat 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()){
	    //LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"No histograms corresponding to eta index value of "<<eta_index<<"."<<ENDLOG;
	    continue;
	  }	  
	  //Get the correct histogrampointer:
	  AliL3Histogram *hist = fParamSpace[eta_index];
	  if(!hist)
	    {
	      //LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"Error getting histogram in index "<<eta_index<<"."<<ENDLOG;
	      continue;
	    }

	  //Do the transformation:
	  AliL3FFloat R = sqrt(fxyz[0]*fxyz[0] + fxyz[1]*fxyz[1]); 
	  AliL3FFloat phi = AliL3Transform::GetPhi(xyz);
	  
	  //Fill the histogram along the phirange
	  for(Int_t b=hist->GetFirstYbin(); b<=hist->GetLastYbin(); b++)
	    {
	      AliL3FFloat phi0 = hist->GetBinCenterY(b);
	      AliL3FFloat kappa = 2*sin(phi - phi0)/R;
	      //hist->Fill(kappa.GetExactVal(),phi0.GetExactVal(),charge);
	      hist->Fill(kappa,phi0,charge);
	    }
	}
      
      //Move the data pointer to the next padrow:
      AliL3MemHandler::UpdateRowPointer(tempPt);
    }
}

#endif
Commit	Line	Data
b46b53c1	1	//$Id$
	2
	3	// Author: Constantin Loizides <mailto:loizides@fi.uib.no>
6173606e	4	//*-- Copyright&Copy CL
6173606e	5
e6cdf93b	6	#include "AliL3StandardIncludes.h"
b46b53c1	7
b46b53c1	8	#include "AliL3Logging.h"
e6cdf93b	9	#include "AliL3MemHandler.h"
b46b53c1	10	#include "AliL3Transform.h"
	11	#include "AliL3DigitData.h"
	12	#include "AliL3Histogram.h"
6173606e	13
6173606e	14	//try to be close to VHDL version (eg. LUTs)
	15	//VESTBO: if you switch that of, you should get your version
	16	//of the transformer!!!
	17	#define VHDLVERSION
	18
	19	//dont use dynamic objects on heap
	20	//VESTBO: this a desperately switched on in order
	21	//to find the bug, which i thought hat something
	22	//to do with a wrong heap assignment! But the
	23	//data still gets overwritten!
	24	#define VHDLSTATIC
	25
b46b53c1	26	#include "AliL3HoughTransformerVhdl.h"
b46b53c1	27
e6cdf93b	28	#if GCCVERSION == 3
	29	using namespace std;
	30	#endif
	31
6173606e	32	/** \class AliL3HoughTransformerVhdl
6173606e	33	// <pre>
b46b53c1	34	//_____________________________________________________________
	35	// AliL3HoughTransformerVhdl
	36	//
6173606e	37	// Hough transformation class for VHDL comparism.
b46b53c1	38	//
6173606e	39	//</pre>
6173606e	40	*/
b46b53c1	41
	42	ClassImp(AliL3HoughTransformerVhdl)
	43
6173606e	44	#ifdef VHDLVERSION
6173606e	45
b46b53c1	46	AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl()
e6cdf93b	47	: AliL3HoughBaseTransformer()//, /fMinRow(0),/ fMaxRow(0)
e6cdf93b	48
b46b53c1	49	{
e6cdf93b	50	cout << "hallo" << endl;
6173606e	51	fParamSpace=0;
e6cdf93b	52
6173606e	53	#ifndef VHDLSTATIC
	54	fLUTX=0;
	55	fLUTY=0;
	56	fLUTEta=0;
	57	fLUTphi0=0;
	58	fLUT2sinphi0=0;
	59	fLUT2cosphi0=0;
	60	#endif
e6cdf93b	61
6173606e	62	fMinRow=0;
6173606e	63	fMaxRow=0;
e6cdf93b	64
e6cdf93b	65	/*
6173606e	66	fNRows=0;
	67	fNEtas=0;
	68	fNPhi0=0;
	69	fSector=0;
	70	fSectorRow=0;
	71	fZSign=0;
e6cdf93b	72	fZLengthPlusOff=0;
	73	fTimeWidth=0;
	74	fPadPitch=0;
	75	fEtaSlice=0;
	76	*/
b46b53c1	77	}
b46b53c1	78
6173606e	79	AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl(Int_t slice,Int_t patch,Int_t n_eta_segments)
6173606e	80	: AliL3HoughBaseTransformer(slice,patch,n_eta_segments)
b46b53c1	81	{
6173606e	82	AliL3HoughTransformerVhdl();
e6cdf93b	83	cout << "wiesow hier??" << endl;
6173606e	84	Init(slice,patch,n_eta_segments);
	85	}
	86
	87	AliL3HoughTransformerVhdl::~AliL3HoughTransformerVhdl()
	88	{
	89	DeleteHistograms();
	90
	91	#ifndef VHDLSTATIC
	92	if(fNRows){
	93	delete[] fLUTX;
	94	delete[] fLUTY;
	95	}
	96	if(fNEtas) delete[] fLUTEta;
	97	#endif
	98	}
	99
e6cdf93b	100	void AliL3HoughTransformerVhdl::Init(Int_t slice,Int_t patch,Int_t n_eta_segments)
6173606e	101	{
	102	cout << "InitVhdl " << slice << " " << patch << " " << n_eta_segments << endl;
	103	AliL3HoughBaseTransformer::Init(slice,patch,n_eta_segments);
	104
	105	#ifndef VHDLSTATIC
	106	//delete old LUT tables
	107	if(fNRows){
	108	fNRows=0;
	109	delete[] fLUTX;
	110	delete[] fLUTY;
	111	}
	112	if(fNEtas){
	113	delete[] fLUTEta;
	114	fNEtas=0;
	115	}
	116	#endif
	117
	118	//set when histogram is filled
	119	//fLUTphi0=0;
	120	//fLUT2sinphi0=0;
	121	//fLUT2cosphi0=0;
	122
	123	Int_t minrow_=AliL3Transform::GetFirstRow(patch);
	124	Int_t maxrow_=AliL3Transform::GetLastRow(patch);
	125	Int_t n_=AliL3Transform::GetNRows(patch);
	126	Int_t sector_=0,sectorrow_=0;
	127	AliL3Transform::Slice2Sector(slice,minrow_,sector_,sectorrow_);
	128	Float_t padpitch_=0.;
	129	if(sector_<AliL3Transform::GetNSectorLow())
	130	padpitch_=AliL3Transform::GetPadPitchWidthLow();
	131	else
	132	padpitch_=AliL3Transform::GetPadPitchWidthUp();
	133
	134	Float_t etamax_=GetEtaMax();
	135	Float_t etamin_=GetEtaMin();
	136	Float_t etaslice_=(etamax_-etamin_)/n_eta_segments;
e6cdf93b	137	#if 0
6173606e	138	//lookup tables for X and Y
	139	#ifndef VHDLSTATIC
	140	fLUTX=new AliL3FFloat[n_];
	141	fLUTY=new AliL3FFloat[n_];
	142	#endif
	143	for(Int_t rr=0;rr<n_;rr++){
	144	fLUTX[rr]=Float_t(AliL3Transform::Row2X(rr+minrow_));
	145	fLUTY[rr]=Float_t(0.5(AliL3Transform::GetNPads(rr+minrow_)-1)padpitch_);
	146	//VESTBO: uncomment to see values (and compare with Print function)
	147	//cout << rr << ": " << (Float_t)fLUTX[rr] << " " << (Float_t)fLUTY[rr] << endl;
	148	}
	149
	150	//lookup tables for rz2s <=> etas
	151	#ifndef VHDLSTATIC
	152	fLUTEta=new AliL3FFloat[n_eta_segments];
	153	#endif
	154	for(Int_t rr=0;rr<n_eta_segments;rr++){
	155	fLUTEta[rr]=CalcRoverZ2(etamin_+(rr+1)*etaslice_);
	156	//VESTBO: uncomment to see values (and compare with Print function)
	157	//cout << rr << ": " << fLUTEta[rr] << endl;
	158	}
e6cdf93b	159	#endif /uncomment arrays/
6173606e	160
	161	//member values
	162	fMinRow=minrow_;
	163	fMaxRow=maxrow_;
	164	fNRows=n_;
	165	fNEtas=n_eta_segments;
	166	fNPhi0=0;
	167	fSector=sector_;
	168	fSectorRow=sectorrow_;
	169	fZSign = slice < 18 ? 1:-1;
	170	fZLengthPlusOff=AliL3Transform::GetZLength()+AliL3Transform::GetZOffset();
	171	fTimeWidth=AliL3Transform::GetZWidth();
	172	fPadPitch=padpitch_;
	173	fEtaSlice=etaslice_;
	174	}
	175
e6cdf93b	176	#if 0
6173606e	177	Float_t AliL3HoughTransformerVhdl::CalcRoverZ2(Float_t eta)
	178	{
	179	Float_t e=exp(2*eta);
	180	Float_t ret=(e+1)/(e-1);
	181	ret*=ret;
	182	return ret;
	183	}
	184
	185	Float_t AliL3HoughTransformerVhdl::CalcEta(Float_t roverz2)
	186	{
	187	Float_t rz=sqrt(roverz2);
	188	if(fZSign>0) rz=-rz;
	189	Float_t ret=(1+rz)/(rz-1);
	190	ret=0.5*log(ret);
	191	return ret;
	192	}
	193
e6cdf93b	194	Float_t AliL3HoughTransformerVhdl::CalcX(Int_t row)
6173606e	195	{
	196	return fLUTX[row];
	197	}
	198
e6cdf93b	199	Float_t AliL3HoughTransformerVhdl::CalcY(Int_t pad,Int_t row)
6173606e	200	{
	201	return pad*fPadPitch-fLUTY[row];
	202	}
	203
e6cdf93b	204	Float_t AliL3HoughTransformerVhdl::CalcZ(Int_t time)
6173606e	205	{
	206	Float_t ret=time*fTimeWidth;
	207	if(fZSign>0) ret=fZLengthPlusOff-ret;
	208	else ret=ret-fZLengthPlusOff;
	209	return ret;
	210	}
e6cdf93b	211	#endif
6173606e	212
	213	void AliL3HoughTransformerVhdl::DeleteHistograms()
	214	{
	215	if(!fParamSpace) return;
	216	for(Int_t i=0; i<GetNEtaSegments(); i++){
	217	if(!fParamSpace[i]) continue;
	218	delete fParamSpace[i];
	219	}
	220	delete[] fParamSpace;
	221
	222	#ifndef VHDLSTATIC
	223	if(fNPhi0){
	224	fNPhi0=0;
	225	delete[] fLUT2sinphi0;
	226	delete[] fLUT2cosphi0;
	227	delete[] fLUTphi0;
	228	}
	229	#endif
	230	}
	231
	232	void AliL3HoughTransformerVhdl::CreateHistograms(Int_t nxbin,Double_t pt_min,Int_t nybin,Double_t phimin,Double_t phimax)
	233	{
	234	//Create the histograms (parameter space).
	235	//These are 2D histograms, span by kappa (curvature of track) and phi0 (emission angle with x-axis).
	236	//The arguments give the range and binning;
	237	//nxbin = #bins in kappa
	238	//nybin = #bins in phi0
	239	//pt_min = mimium Pt of track (corresponding to maximum kappa)
	240	//phi_min = mimimum phi0 (degrees)
	241	//phi_max = maximum phi0 (degrees)
	242
	243	AliL3FFloat x = AliL3Transform::GetBFieldValue()/pt_min;
	244	AliL3FFloat torad = AliL3Transform::Pi()/180;
	245	AliL3FFloat phimin_=phimin*torad;
	246	AliL3FFloat phimax_=phimax*torad;
	247
	248	CreateHistograms(nxbin,-x,x,nybin,phimin_,phimax_);
	249	}
	250
	251	void AliL3HoughTransformerVhdl::CreateHistograms(Int_t nxbin,Double_t xmin,Double_t xmax,
	252	Int_t nybin,Double_t ymin,Double_t ymax)
	253	{
	254
	255	AliL3FFloat xminf=xmin;
	256	AliL3FFloat xmaxf=xmax;
	257	AliL3FFloat yminf=ymin;
	258	AliL3FFloat ymaxf=ymax;
	259	/*
	260	cout << xminf << endl;
	261	cout << xmaxf << " " << ((xmaxf-xminf)/nxbin) << endl;
	262	cout << yminf << endl;
	263	cout << ymaxf << " " << ((ymaxf-yminf)/nybin) << endl;
	264	*/
	265
e6cdf93b	266
6173606e	267	fParamSpace = new AliL3Histogram*[GetNEtaSegments()];
e6cdf93b	268
6173606e	269	Char_t histname[256];
	270	for(Int_t i=0; i<GetNEtaSegments(); i++)
	271	{
	272	sprintf(histname,"paramspace_vhdl_%d",i);
e6cdf93b	273	fParamSpace[i] = new AliL3Histogram(histname,"id",nxbin,xminf,xmaxf,nybin,yminf,ymaxf);
6173606e	274	}
	275
	276	//create lookup table for sin and cos
	277	fNPhi0=nybin+1;
e6cdf93b	278	#if 0
6173606e	279	#ifndef VHDLSTATIC
	280	fLUTphi0=new AliL3FFloat[fNPhi0];
	281	fLUT2sinphi0=new AliL3FFloat[fNPhi0];
	282	fLUT2cosphi0=new AliL3FFloat[fNPhi0];
	283	#endif
	284	AliL3FFloat diff=(ymax-ymin)/nybin;
	285	AliL3FFloat phi0=ymin-0.5*diff;
	286	for(Int_t i=0; i<fNPhi0; i++){
	287	phi0+=diff;
	288	fLUTphi0[i]=phi0;
	289	fLUT2sinphi0[i]=Float_t(2*sin(phi0));
	290	fLUT2cosphi0[i]=Float_t(2*cos(phi0));
	291	//VESTBO: uncomment to see values (and compare with Print function)
	292	//cout << i << ": " << fLUTphi0[i] << " " << fLUT2sinphi0[i] << " " << fLUT2cosphi0[i] << endl;
	293	}
e6cdf93b	294	#endif /arrays/
6173606e	295	}
	296
	297	void AliL3HoughTransformerVhdl::Reset()
	298	{
	299	//Reset all the histograms
	300
	301	if(!fParamSpace)
	302	{
	303	LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::Reset","Histograms")
	304	<<"No histograms to reset"<<ENDLOG;
	305	return;
	306	}
	307
	308	for(Int_t i=0; i<GetNEtaSegments(); i++)
	309	fParamSpace[i]->Reset();
	310	}
	311
e6cdf93b	312	AliL3Histogram *AliL3HoughTransformerVhdl::GetHistogram(Int_t eta_index)
6173606e	313	{
	314	if(!fParamSpace \|\| eta_index >= GetNEtaSegments() \|\| eta_index < 0) return 0;
	315	if(!fParamSpace[eta_index]) return 0;
	316	return fParamSpace[eta_index];
	317	}
	318
e6cdf93b	319	#if 0
e6cdf93b	320	Int_t AliL3HoughTransformerVhdl::FindIndex(Double_t rz2)
6173606e	321	{
	322	Int_t index=0;
	323	while((index<fNEtas)&&(rz2<=fLUTEta[index])){
	324	index++;
	325	//cout << index << ": " << rz2 << " " << fLUTEta[index] << endl;
	326	}
	327	return index;
	328	}
e6cdf93b	329	#endif
6173606e	330
	331	Int_t AliL3HoughTransformerVhdl::GetEtaIndex(Double_t eta)
	332	{
e6cdf93b	333	//AliL3FFloat rz2=CalcRoverZ2(eta);
	334	//return FindIndex(rz2);
	335	return 0;
6173606e	336	}
6173606e	337
afd8fed4	338	Double_t AliL3HoughTransformerVhdl::GetEta(Int_t eta_index,Int_t slice){
6173606e	339	if(eta_index >= fNEtas \|\| eta_index < 0){
	340	LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::GetEta","Index")
	341	<<"Index out of range."<<ENDLOG;
	342	return 0.;
	343	}
e6cdf93b	344	return 0;//(CalcEta(fLUTEta[eta_index])-0.5*fEtaSlice);
6173606e	345	}
	346
	347	void AliL3HoughTransformerVhdl::TransformCircle()
	348	{
	349	//Transform the input data with a circle HT.
	350	//The function loops over all the data, and transforms each pixel with the equations:
	351	//
	352	//kappa = 1/x * ( y/x2cos(phi0) - 2*sin(phi0) )
	353
e6cdf93b	354	#if 0
6173606e	355	AliL3DigitRowData *tempPt = GetDataPointer();
	356	if(!tempPt)
	357	{
	358	LOG(AliL3Log::kError,"AliL3HoughTransformerVhdl::TransformCircle","Data")
	359	<<"No input data "<<ENDLOG;
	360	return;
	361	}
	362
	363	AliL3FFloat x=0.,y=0.,z=0.,rz2=0.;
	364	//Loop over the padrows:
	365	for(Int_t i=fMinRow, row=0; i<=fMaxRow; i++, row++){
	366
	367	//Get the data on this padrow:
	368	AliL3DigitData *digPt = tempPt->fDigitData;
	369	if(i != (Int_t)tempPt->fRow){
	370	LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Padrows")
	371	<<"Mismatching padrow numbering "<<i<<"!="<<(Int_t)tempPt->fRow<<ENDLOG;
	372	continue;
	373	}
	374
	375	//Loop over the data on this padrow:
	376	for(UInt_t j=0; j<tempPt->fNDigit; j++){
	377	Int_t charge = digPt[j].fCharge;
	378	if((Int_t)charge <= GetLowerThreshold()) continue;
	379
	380	Int_t pad = digPt[j].fPad;
	381	Int_t time = digPt[j].fTime;
	382
	383	x=CalcX(row);
	384	y=CalcY(pad,row);
	385	z=CalcZ(time);
	386
	387	//find eta slice
	388	rz2=1+(xx+yy)/z/z;
	389	Int_t eta_index = FindIndex(rz2);
	390	if(eta_index < 0 \|\| eta_index >= GetNEtaSegments()){
	391	//LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"No histograms corresponding to eta index value of "<<eta_index<<"."<<ENDLOG;
	392	continue;
	393	}
	394	//Get the correct histogrampointer:
	395	AliL3Histogram *hist = fParamSpace[eta_index];
	396	if(!hist){
	397	//LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"Error getting histogram in index "<<eta_index<<"."<<ENDLOG;
	398	continue;
	399	}
	400
	401	//Fill the histogram along the phirange
	402	AliL3FFloat kappa,ydx=y/x;
	403	for(Int_t b=0; b<fNPhi0; b++){
	404	kappa=ydx*fLUT2cosphi0[b]-fLUT2sinphi0[b];
	405	kappa/=x;
	406	hist->Fill(kappa,fLUTphi0[b],charge);
	407	//cout << kappa << " " << fLUTphi0[b] << " " << charge << endl;
	408	}
	409
	410	}
	411	//Move the data pointer to the next padrow:
	412	AliL3MemHandler::UpdateRowPointer(tempPt);
	413	}
e6cdf93b	414	#endif
6173606e	415	}
	416
	417	void AliL3HoughTransformerVhdl::Print()
	418	{
	419	cout << "fSlice: " << GetSlice() << endl;
	420	cout << "fPatch: " << GetPatch() << endl;
	421	cout << "fSector: " << fSector << endl;
	422	cout << "fSectorRow: " << fSectorRow << endl;
	423	cout << "fMinRow: " << fMinRow << endl;
	424	cout << "fMaxRow: " << fMaxRow << endl;
	425	cout << "fNRows: " << fNRows << endl;
	426	cout << "fNEtas: " << fNEtas << endl;
	427	cout << "fNPhi0: " << fNPhi0 << endl;
	428	cout << "fZSign: " << fZSign << endl;
	429	cout << "fZLengthPlusOff: " << fZLengthPlusOff << endl;
	430	cout << "fPadPitch: " << fPadPitch << endl;
	431	cout << "fTimeWidth: " << fTimeWidth << endl;
e6cdf93b	432	#if 0
6173606e	433	if(!fNRows) return;
	434	cout << "fLUTX " << fNRows << endl;
	435	for(Int_t i=0;i<fNRows;i++) cout << "fLUTX[" << i << "]=" << (Float_t)fLUTX[i] << endl;
	436	cout << "fLUTY " << fNRows << endl;
	437	for(Int_t i=0;i<fNRows;i++) cout << "fLUTY[" << i << "]=" << fLUTY[i] << endl;
	438	if(!fNEtas) return;
	439	cout << "fLUTEta " << fNEtas << endl;
	440	for(Int_t i=0;i<fNEtas;i++) cout << "fLUTEta[" << i << "]=" << fLUTEta[i] << endl;
	441	if(!fNPhi0) return;
	442	cout << "fLUTphi0 " << fNPhi0 << endl;
	443	for(Int_t i=0;i<fNPhi0;i++) cout << "fLUTPhi0[" << i << "]=" << fLUTphi0[i] << endl;
	444	cout << "fLUT2sinphi0 " << fNPhi0 << endl;
	445	for(Int_t i=0;i<fNPhi0;i++) cout << "fLUT2sinphi0[" << i << "]=" << fLUT2sinphi0[i] << endl;
	446	cout << "fLUT2cosphi0 " << fNPhi0 << endl;
	447	for(Int_t i=0;i<fNPhi0;i++) cout << "fLUT2cosphi0[" << i << "]=" << fLUT2cosphi0[i] << endl;
e6cdf93b	448	#endif
6173606e	449	}
	450
	451	//end vhdl version
	452	#else
	453	//standard version
	454
	455	AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl()
	456	{
	457	fParamSpace=0;
	458	}
	459
	460	AliL3HoughTransformerVhdl::AliL3HoughTransformerVhdl(Int_t slice,Int_t patch,Int_t n_eta_segments)
	461	: AliL3HoughBaseTransformer(slice,patch,n_eta_segments)
	462	{
	463	fParamSpace=0;
b46b53c1	464	}
	465
	466	AliL3HoughTransformerVhdl::~AliL3HoughTransformerVhdl()
	467	{
	468	DeleteHistograms();
	469	}
	470
	471	void AliL3HoughTransformerVhdl::DeleteHistograms()
	472	{
	473	if(!fParamSpace)
	474	return;
	475	for(Int_t i=0; i<GetNEtaSegments(); i++)
	476	{
	477	if(!fParamSpace[i]) continue;
	478	delete fParamSpace[i];
	479	}
6173606e	480	delete[] fParamSpace;
b46b53c1	481	}
b46b53c1	482
6173606e	483	void AliL3HoughTransformerVhdl::CreateHistograms(Int_t nxbin,Double_t pt_min,Int_t nybin,Double_t phimin,Double_t phimax)
b46b53c1	484	{
	485	//Create the histograms (parameter space).
	486	//These are 2D histograms, span by kappa (curvature of track) and phi0 (emission angle with x-axis).
	487	//The arguments give the range and binning;
	488	//nxbin = #bins in kappa
	489	//nybin = #bins in phi0
	490	//pt_min = mimium Pt of track (corresponding to maximum kappa)
	491	//phi_min = mimimum phi0 (degrees)
	492	//phi_max = maximum phi0 (degrees)
	493
6173606e	494	Double_t x = AliL3Transform::GetBFieldValue()/pt_min;
b46b53c1	495	Double_t torad = AliL3Transform::Pi()/180;
6173606e	496
b46b53c1	497	CreateHistograms(nxbin,-1.x,x,nybin,phimintorad,phimax*torad);
	498	}
	499
	500	void AliL3HoughTransformerVhdl::CreateHistograms(Int_t nxbin,Double_t xmin,Double_t xmax,
6173606e	501	Int_t nybin,Double_t ymin,Double_t ymax)
b46b53c1	502	{
6173606e	503
	504	AliL3FFloat xminf(xmin);
	505	AliL3FFloat xmaxf(xmax);
	506	AliL3FFloat yminf(ymin);
	507	AliL3FFloat ymaxf(ymax);
e6cdf93b	508
6173606e	509	cout << xminf << endl;
	510	cout << xmaxf << " " << ((xmaxf-xminf)/nxbin) << endl;
	511	cout << yminf << endl;
	512	cout << ymaxf << " " << ((ymaxf-yminf)/nybin) << endl;
e6cdf93b	513
b46b53c1	514
	515	fParamSpace = new AliL3Histogram*[GetNEtaSegments()];
	516
	517	Char_t histname[256];
	518	for(Int_t i=0; i<GetNEtaSegments(); i++)
	519	{
6173606e	520	sprintf(histname,"paramspace_vhdl_%d",i);
e6cdf93b	521	fParamSpace[i] = new AliL3Histogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
b46b53c1	522	}
	523	}
	524
	525	void AliL3HoughTransformerVhdl::Reset()
	526	{
	527	//Reset all the histograms
	528
	529	if(!fParamSpace)
	530	{
	531	LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::Reset","Histograms")
	532	<<"No histograms to reset"<<ENDLOG;
	533	return;
	534	}
	535
	536	for(Int_t i=0; i<GetNEtaSegments(); i++)
	537	fParamSpace[i]->Reset();
	538	}
	539
e6cdf93b	540	AliL3Histogram *AliL3HoughTransformerVhdl::GetHistogram(Int_t eta_index)
6173606e	541	{
	542	if(!fParamSpace \|\| eta_index >= GetNEtaSegments() \|\| eta_index < 0) return 0;
	543	if(!fParamSpace[eta_index]) return 0;
	544	return fParamSpace[eta_index];
	545	}
b46b53c1	546
	547	Int_t AliL3HoughTransformerVhdl::GetEtaIndex(Double_t eta)
	548	{
	549	//Return the histogram index of the corresponding eta.
6173606e	550	AliL3FFloat etaslice = (GetEtaMax() - GetEtaMin())/GetNEtaSegments();
6173606e	551	AliL3FFloat index = (eta-GetEtaMin())/etaslice;
b46b53c1	552
6173606e	553	//if((Int_t)index==0) cout << index << " " << (Int_t) index << endl;
b46b53c1	554	return (Int_t)index;
	555	}
	556
6173606e	557	Double_t AliL3HoughTransformerVhdl::GetEta(Int_t eta_index)
	558	{
	559	if(eta_index >= GetNEtaSegments() \|\| eta_index < 0){
	560	return 0;
	561	}
	562	Double_t eta_slice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments();
	563	Double_t eta=(Double_t)((eta_index+0.5)*eta_slice);
	564	if(GetSlice()>17) eta*=-1;
	565	return eta;
	566	}
	567
b46b53c1	568	void AliL3HoughTransformerVhdl::TransformCircle()
	569	{
	570	//Transform the input data with a circle HT.
	571	//The function loops over all the data, and transforms each pixel with the equations:
	572	//
	573	//kappa = 2/R*sin(phi - phi0)
	574	//
	575	//where R = sqrt(xx +yy), and phi = arctan(y/x)
	576	//
	577	//Each pixel then transforms into a curve in the (kappa,phi0)-space. In order to find
	578	//which histogram in which the pixel should be transformed, the eta-value is calcluated
	579	//and the proper histogram index is found by GetEtaIndex(eta).
	580
	581
	582	AliL3DigitRowData *tempPt = GetDataPointer();
	583	if(!tempPt)
	584	{
	585	LOG(AliL3Log::kError,"AliL3HoughTransformerVhdl::TransformCircle","Data")
	586	<<"No input data "<<ENDLOG;
	587	return;
	588	}
	589
	590	//Loop over the padrows:
	591	for(Int_t i=AliL3Transform::GetFirstRow(GetPatch()); i<=AliL3Transform::GetLastRow(GetPatch()); i++)
	592	{
	593	//Get the data on this padrow:
	594	AliL3DigitData *digPt = tempPt->fDigitData;
	595	if(i != (Int_t)tempPt->fRow)
	596	{
6173606e	597	LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Padrows")
6173606e	598	<<"Mismatching padrow numbering."<<ENDLOG;
b46b53c1	599	continue;
	600	}
	601
6173606e	602
b46b53c1	603	//Loop over the data on this padrow:
	604	for(UInt_t j=0; j<tempPt->fNDigit; j++)
	605	{
	606	UShort_t charge = digPt[j].fCharge;
6173606e	607	if((Int_t)charge <= GetLowerThreshold()) continue;
6173606e	608
b46b53c1	609	UChar_t pad = digPt[j].fPad;
b46b53c1	610	UShort_t time = digPt[j].fTime;
6173606e	611
b46b53c1	612	Int_t sector,row;
b46b53c1	613	Float_t xyz[3];
6173606e	614	AliL3FFloat fxyz[3];
6173606e	615
b46b53c1	616	//Transform data to local cartesian coordinates:
	617	AliL3Transform::Slice2Sector(GetSlice(),i,sector,row);
	618	AliL3Transform::Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
6173606e	619
	620	//trunc to fixed format.
	621	for(int i=0;i<3;i++) {
	622	fxyz[i]=xyz[i];
	623	//cout << fxyz[i];
	624	xyz[i]=fxyz[i];
	625	}
	626	//cout << endl;
	627
b46b53c1	628	//Calculate the eta:
6173606e	629	AliL3FFloat eta = AliL3Transform::GetEta(xyz);
b46b53c1	630
	631	//Get the corresponding index, which determines which histogram to fill:
	632	Int_t eta_index = GetEtaIndex(eta);
6173606e	633	if(eta_index < 0 \|\| eta_index >= GetNEtaSegments()){
6173606e	634	//LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"No histograms corresponding to eta index value of "<<eta_index<<"."<<ENDLOG;
b46b53c1	635	continue;
6173606e	636	}
b46b53c1	637	//Get the correct histogrampointer:
	638	AliL3Histogram *hist = fParamSpace[eta_index];
	639	if(!hist)
	640	{
6173606e	641	//LOG(AliL3Log::kWarning,"AliL3HoughTransformerVhdl::TransformCircle","Histograms")<<"Error getting histogram in index "<<eta_index<<"."<<ENDLOG;
b46b53c1	642	continue;
	643	}
	644
	645	//Do the transformation:
6173606e	646	AliL3FFloat R = sqrt(fxyz[0]fxyz[0] + fxyz[1]fxyz[1]);
6173606e	647	AliL3FFloat phi = AliL3Transform::GetPhi(xyz);
b46b53c1	648
	649	//Fill the histogram along the phirange
	650	for(Int_t b=hist->GetFirstYbin(); b<=hist->GetLastYbin(); b++)
	651	{
6173606e	652	AliL3FFloat phi0 = hist->GetBinCenterY(b);
	653	AliL3FFloat kappa = 2*sin(phi - phi0)/R;
	654	//hist->Fill(kappa.GetExactVal(),phi0.GetExactVal(),charge);
b46b53c1	655	hist->Fill(kappa,phi0,charge);
	656	}
	657	}
	658
	659	//Move the data pointer to the next padrow:
	660	AliL3MemHandler::UpdateRowPointer(tempPt);
	661	}
	662	}
	663
6173606e	664	#endif