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

//Author:        Anders Strand Vestbo
//Last Modified: 14.09.01

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

ClassImp(AliL3HoughTransformer)

AliL3HoughTransformer::AliL3HoughTransformer()
{
  //Default constructor
  
}

AliL3HoughTransformer::AliL3HoughTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments)
{
  fSlice = slice;
  fPatch = patch;
  fNEtaSegments = n_eta_segments;
  fEtaMin = 0;
  fEtaMax = fSlice < 18 ? 0.9 : -0.9;
  fTransform = new AliL3Transform();
  fThreshold = 0;
  fNDigitRowData=0;
  fDigitRowData=0;
}

AliL3HoughTransformer::~AliL3HoughTransformer()
{
  if(fTransform)
    delete fTransform;
  DeleteHistograms();
}

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

void AliL3HoughTransformer::CreateHistograms(Int_t nxbin,Double_t pt_min,
					     Int_t nybin,Double_t phimin,Double_t phimax)
{
  //Set the minimum absolute pt value, and phi0 angles given in degrees.

  Double_t torad = 3.1415/180;
  Double_t bfact = 0.0029980;
  Double_t bfield = 0.2;
  Double_t x = bfact*bfield/pt_min;
  CreateHistograms(nxbin,-1.*x,x,nybin,phimin*torad,phimax*torad);
}

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];
  for(Int_t i=0; i<fNEtaSegments; i++)
    {
      sprintf(histname,"paramspace_%d",i);
      fParamSpace[i] = new AliL3Histogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
    }
}

void AliL3HoughTransformer::Reset()
{
  //Reset all the histograms

  if(!fParamSpace)
    {
      LOG(AliL3Log::kWarning,"AliL3HoughTransformer::Reset","Histograms")
	<<"No histograms to reset"<<ENDLOG;
      return;
    }
  
  for(Int_t i=0; i<fNEtaSegments; i++)
    fParamSpace[i]->Reset();
}

void AliL3HoughTransformer::SetInputData(UInt_t ndigits,AliL3DigitRowData *ptr)
{
  fNDigitRowData = ndigits;
  fDigitRowData = ptr;
}

void AliL3HoughTransformer::TransformCircle()
{
  //Transform the input data with a circle HT.
  

  //Set pointer to the data
  AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fDigitRowData;
  if(!tempPt || fNDigitRowData==0)
    {
      printf("\nAliL3HoughTransformer::TransformCircle : 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::TransformCircle : 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);
	  Double_t eta = fTransform->GetEta(xyz);
	  Int_t eta_index = (Int_t)((eta-fEtaMin)/etaslice);
	  if(eta_index < 0 || eta_index >= fNEtaSegments)
	    continue;
	  
	  //Get the correct histogrampointer:
	  AliL3Histogram *hist = fParamSpace[eta_index];
	  if(!hist)
	    {
	      printf("AliL3HoughTransformer::TransformCircle : Error getting histogram in index %d\n",eta_index);
	      continue;
	    }

	  //Start transformation
	  Float_t R = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]); // + xyz[2]*xyz[2]);
	  Float_t phi = fTransform->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);
	    }
	}
      AliL3MemHandler::UpdateRowPointer(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);
	    }
	}
      AliL3MemHandler::UpdateRowPointer(tempPt);
    }
  
}
Commit	Line	Data
	1	//Author: Anders Strand Vestbo
	2	//Last Modified: 14.09.01
	3
	4	#include "AliL3MemHandler.h"
	5	#include "AliL3Logging.h"
	6	#include "AliL3HoughTransformer.h"
	7	#include "AliL3Defs.h"
	8	#include "AliL3Transform.h"
	9	#include "AliL3DigitData.h"
	10	#include "AliL3Histogram.h"
	11
	12	ClassImp(AliL3HoughTransformer)
	13
	14	AliL3HoughTransformer::AliL3HoughTransformer()
	15	{
	16	//Default constructor
	17
	18	}
	19
	20	AliL3HoughTransformer::AliL3HoughTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments)
	21	{
	22	fSlice = slice;
	23	fPatch = patch;
	24	fNEtaSegments = n_eta_segments;
	25	fEtaMin = 0;
	26	fEtaMax = fSlice < 18 ? 0.9 : -0.9;
	27	fTransform = new AliL3Transform();
	28	fThreshold = 0;
	29	fNDigitRowData=0;
	30	fDigitRowData=0;
	31	}
	32
	33	AliL3HoughTransformer::~AliL3HoughTransformer()
	34	{
	35	if(fTransform)
	36	delete fTransform;
	37	DeleteHistograms();
	38	}
	39
	40	void AliL3HoughTransformer::DeleteHistograms()
	41	{
	42	if(!fParamSpace)
	43	return;
	44	for(Int_t i=0; i<fNEtaSegments; i++)
	45	{
	46	if(!fParamSpace[i]) continue;
	47	delete fParamSpace[i];
	48	}
	49	delete [] fParamSpace;
	50	}
	51
	52	void AliL3HoughTransformer::CreateHistograms(Int_t nxbin,Double_t pt_min,
	53	Int_t nybin,Double_t phimin,Double_t phimax)
	54	{
	55	//Set the minimum absolute pt value, and phi0 angles given in degrees.
	56
	57	Double_t torad = 3.1415/180;
	58	Double_t bfact = 0.0029980;
	59	Double_t bfield = 0.2;
	60	Double_t x = bfact*bfield/pt_min;
	61	CreateHistograms(nxbin,-1.x,x,nybin,phimintorad,phimax*torad);
	62	}
	63
	64	void AliL3HoughTransformer::CreateHistograms(Int_t nxbin,Double_t xmin,Double_t xmax,
	65	Int_t nybin,Double_t ymin,Double_t ymax)
	66	{
	67
	68	fParamSpace = new AliL3Histogram*[fNEtaSegments];
	69
	70	Char_t histname[256];
	71	for(Int_t i=0; i<fNEtaSegments; i++)
	72	{
	73	sprintf(histname,"paramspace_%d",i);
	74	fParamSpace[i] = new AliL3Histogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
	75	}
	76	}
	77
	78	void AliL3HoughTransformer::Reset()
	79	{
	80	//Reset all the histograms
	81
	82	if(!fParamSpace)
	83	{
	84	LOG(AliL3Log::kWarning,"AliL3HoughTransformer::Reset","Histograms")
	85	<<"No histograms to reset"<<ENDLOG;
	86	return;
	87	}
	88
	89	for(Int_t i=0; i<fNEtaSegments; i++)
	90	fParamSpace[i]->Reset();
	91	}
	92
	93	void AliL3HoughTransformer::SetInputData(UInt_t ndigits,AliL3DigitRowData *ptr)
	94	{
	95	fNDigitRowData = ndigits;
	96	fDigitRowData = ptr;
	97	}
	98
	99	void AliL3HoughTransformer::TransformCircle()
	100	{
	101	//Transform the input data with a circle HT.
	102
	103
	104	//Set pointer to the data
	105	AliL3DigitRowData tempPt = (AliL3DigitRowData)fDigitRowData;
	106	if(!tempPt \|\| fNDigitRowData==0)
	107	{
	108	printf("\nAliL3HoughTransformer::TransformCircle : No input data!!!\n\n");
	109	return;
	110	}
	111
	112	Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments;
	113	for(Int_t i=NRows[fPatch][0]; i<=NRows[fPatch][1]; i++)
	114	{
	115	AliL3DigitData *digPt = tempPt->fDigitData;
	116	if(i != (Int_t)tempPt->fRow)
	117	{
	118	printf("AliL3HoughTransform::TransformCircle : Mismatching padrow numbering\n");
	119	continue;
	120	}
	121	for(UInt_t j=0; j<tempPt->fNDigit; j++)
	122	{
	123	UShort_t charge = digPt[j].fCharge;
	124	UChar_t pad = digPt[j].fPad;
	125	UShort_t time = digPt[j].fTime;
	126	if(charge <= fThreshold)
	127	continue;
	128	Int_t sector,row;
	129	Float_t xyz[3];
	130	fTransform->Slice2Sector(fSlice,i,sector,row);
	131	fTransform->Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
	132	Double_t eta = fTransform->GetEta(xyz);
	133	Int_t eta_index = (Int_t)((eta-fEtaMin)/etaslice);
	134	if(eta_index < 0 \|\| eta_index >= fNEtaSegments)
	135	continue;
	136
	137	//Get the correct histogrampointer:
	138	AliL3Histogram *hist = fParamSpace[eta_index];
	139	if(!hist)
	140	{
	141	printf("AliL3HoughTransformer::TransformCircle : Error getting histogram in index %d\n",eta_index);
	142	continue;
	143	}
	144
	145	//Start transformation
	146	Float_t R = sqrt(xyz[0]xyz[0] + xyz[1]xyz[1]); // + xyz[2]*xyz[2]);
	147	Float_t phi = fTransform->GetPhi(xyz);
	148
	149	//Fill the histogram along the phirange
	150	for(Int_t b=hist->GetFirstYbin(); b<=hist->GetLastYbin(); b++)
	151	{
	152	Float_t phi0 = hist->GetBinCenterY(b);
	153	Float_t kappa = 2*sin(phi - phi0)/R;
	154	hist->Fill(kappa,phi0,charge);
	155	}
	156	}
	157	AliL3MemHandler::UpdateRowPointer(tempPt);
	158	}
	159	}
	160
	161	void AliL3HoughTransformer::TransformLine()
	162	{
	163	//Do a line transform on the data.
	164
	165
	166	AliL3DigitRowData tempPt = (AliL3DigitRowData)fDigitRowData;
	167	if(!tempPt \|\| fNDigitRowData==0)
	168	{
	169	printf("\nAliL3HoughTransformer::TransformLine : No input data!!!\n\n");
	170	return;
	171	}
	172
	173	Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments;
	174	for(Int_t i=NRows[fPatch][0]; i<=NRows[fPatch][1]; i++)
	175	{
	176	AliL3DigitData *digPt = tempPt->fDigitData;
	177	if(i != (Int_t)tempPt->fRow)
	178	{
	179	printf("AliL3HoughTransform::TransformLine : Mismatching padrow numbering\n");
	180	continue;
	181	}
	182	for(UInt_t j=0; j<tempPt->fNDigit; j++)
	183	{
	184	UShort_t charge = digPt[j].fCharge;
	185	UChar_t pad = digPt[j].fPad;
	186	UShort_t time = digPt[j].fTime;
	187	if(charge < fThreshold)
	188	continue;
	189	Int_t sector,row;
	190	Float_t xyz[3];
	191	fTransform->Slice2Sector(fSlice,i,sector,row);
	192	fTransform->Raw2Local(xyz,sector,row,(Int_t)pad,(Int_t)time);
	193	Float_t eta = fTransform->GetEta(xyz);
	194	Int_t eta_index = (Int_t)(eta/etaslice);
	195	if(eta_index < 0 \|\| eta_index >= fNEtaSegments)
	196	continue;
	197
	198	//Get the correct histogram:
	199	AliL3Histogram *hist = fParamSpace[eta_index];
	200	if(!hist)
	201	{
	202	printf("AliL3HoughTransformer::TransformLine : Error getting histogram in index %d\n",eta_index);
	203	continue;
	204	}
	205	for(Int_t xbin=hist->GetFirstXbin(); xbin<hist->GetLastXbin(); xbin++)
	206	{
	207	Double_t theta = hist->GetBinCenterX(xbin);
	208	Double_t rho = xyz[0]cos(theta) + xyz[1]sin(theta);
	209	hist->Fill(theta,rho,charge);
	210	}
	211	}
	212	AliL3MemHandler::UpdateRowPointer(tempPt);
	213	}
	214
	215	}