fNEtas=0;
fNPhi0=0;
fSector=0;
+ fSlice=0;
fSectorRow=0;
fZSign=0;
fZLengthPlusOff=0;
AliL3HoughTransformerLUT::AliL3HoughTransformerLUT(Int_t slice,Int_t patch,Int_t n_eta_segments) : AliL3HoughBaseTransformer(slice,patch,n_eta_segments)
{
- AliL3HoughTransformerLUT();
+ fMinRow=0;
+ fMaxRow=0;
+
+ fNRows=0;
+ fNEtas=0;
+ fNPhi0=0;
+ fSector=0;
+ fSectorRow=0;
+ fZSign=0;
+ fZLengthPlusOff=0;
+ fTimeWidth=0;
+ fPadPitch=0;
+ fEtaSlice=0;
+
+ fLUTX=0;
+ fLUTY=0;
+ fLUTEta=0;
+ fLUTphi0=0;
+ fLUT2sinphi0=0;
+ fLUT2cosphi0=0;
Init(slice,patch,n_eta_segments);
}
#ifdef do_mc
if(fTrackID)
{
- for(Int_t i=0; i<GetNEtaSegments(); i++)
+ for(Int_t i=0; i<fNEtas; i++)
{
if(!fTrackID[i]) continue;
delete fTrackID[i];
delete[] fLUTY;
fNRows=0;
}
+
if(fNEtas){
delete[] fLUTEta;
fNEtas=0;
}
if(!fParamSpace){
- for(Int_t i=0; i<GetNEtaSegments(); i++)
+ for(Int_t i=0; i<fNEtas; i++)
{
if(!fParamSpace[i]) continue;
delete fParamSpace[i];
fMaxRow=AliL3Transform::GetLastRow(patch);;
fNRows=AliL3Transform::GetNRows(patch);;
fNEtas=n_eta_segments;
+ if(fNEtas!=GetNEtaSegments()) {
+ cerr << "AliL3HoughTransformerLUT::Init -> Error: Number of Etas must be equal!" << endl;
+ exit(1);
+ }
AliL3Transform::Slice2Sector(slice,fMinRow,fSector,fSectorRow);
fZSign = slice < 18 ? 1:-1;
+ fSlice = slice;
fZLengthPlusOff=AliL3Transform::GetZLength()+AliL3Transform::GetZOffset();
fTimeWidth=AliL3Transform::GetZWidth();
//lookup tables for rz2s <=> etas
fLUTEta=new Float_t[fNEtas];
for(Int_t rr=0;rr<fNEtas;rr++){
- fLUTEta[rr]=CalcRoverZ2(etamin_+(rr+1)*fEtaSlice);
- //cout << rr << ": " << fLUTEta[rr] << endl;
+ Float_t eta=etamin_+(rr+1)*fEtaSlice;
+ fLUTEta[rr]=CalcRoverZ2(eta);
+ //cout << rr << ": " << eta << " " << fLUTEta[rr] << " " << GetEta(rr,fSlice) << endl;
}
}
{
AliL3Histogram *hist = fParamSpace[0];
Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
- cout<<"Allocating "<<GetNEtaSegments()*ncells*sizeof(TrackIndex)<<" bytes to fTrackID"<<endl;
- fTrackID = new TrackIndex*[GetNEtaSegments()];
- for(Int_t i=0; i<GetNEtaSegments(); i++)
+ cout<<"Allocating "<<fNEtas*ncells*sizeof(TrackIndex)<<" bytes to fTrackID"<<endl;
+ fTrackID = new TrackIndex*[fNEtas];
+ for(Int_t i=0; i<fNEtas; i++)
fTrackID[i] = new TrackIndex[ncells];
}
#endif
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));
- cout << i << ": " << fLUTphi0[i] << " " << fLUT2sinphi0[i] << " " << fLUT2cosphi0[i] << endl;
+ fLUT2sinphi0[i]=2.*sin(phi0);
+ fLUT2cosphi0[i]=2.*cos(phi0);
+ //cout << i << ": " << fLUTphi0[i] << " " << fLUT2sinphi0[i] << " " << fLUT2cosphi0[i] << endl;
}
}
return;
}
- for(Int_t i=0; i<GetNEtaSegments(); i++)
+ for(Int_t i=0; i<fNEtas; i++)
fParamSpace[i]->Reset();
#ifdef do_mc
if(fDoMC)
{
AliL3Histogram *hist = fParamSpace[0];
Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
- for(Int_t i=0; i<GetNEtaSegments(); i++)
+ for(Int_t i=0; i<fNEtas; i++)
memset(fTrackID[i],0,ncells*sizeof(TrackIndex));
}
#endif
}
-
Int_t AliL3HoughTransformerLUT::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;
+
+ Float_t rz2=CalcRoverZ2(eta);
+ return FindIndex(rz2);
+}
+
+Int_t AliL3HoughTransformerLUT::FindIndex(Float_t rz2)
+{
+ Int_t index=0; //could improve search through devide and conquere strategy
+ while((index<fNEtas)&&(rz2<=fLUTEta[index])){
+ index++;
+ //cout << index << ": " << rz2 << " " << fLUTEta[index] << endl;
+ }
+ return index;
}
inline AliL3Histogram *AliL3HoughTransformerLUT::GetHistogram(Int_t eta_index)
{
- if(!fParamSpace || eta_index >= GetNEtaSegments() || eta_index < 0)
+ if(!fParamSpace || eta_index >= fNEtas || eta_index < 0)
return 0;
if(!fParamSpace[eta_index])
return 0;
Float_t AliL3HoughTransformerLUT::CalcEta(Float_t roverz2)
{
Float_t rz=sqrt(roverz2);
- if(fZSign>0) rz=-rz;
+ if(fZSign<0) rz=-rz;
Float_t ret=(1+rz)/(rz-1);
ret=0.5*log(ret);
return ret;
Double_t AliL3HoughTransformerLUT::GetEta(Int_t eta_index,Int_t slice)
{
- Double_t eta_slice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments();
- Double_t eta=(Double_t)((eta_index+0.5)*eta_slice);
- if(slice>17) eta*=-1;
- return eta;
+ if(eta_index >= fNEtas || eta_index < 0){
+ //LOG(AliL3Log::kWarning,"AliL3HoughTransformerLUT::GetEta","Index") << "Index out of range."<<ENDLOG;
+ return 0.;
+ }
+ if(slice != fSlice){
+ //LOG(AliL3Log::kWarning,"AliL3HoughTransformerLUT::GetEta","Index") << "Given slice does not match internal slice."<<ENDLOG;
+ return 0.;
+ }
+
+ return (CalcEta(fLUTEta[eta_index])-0.5*fEtaSlice);
}
void AliL3HoughTransformerLUT::TransformCircle()
{
-#if 0
//Transform the input data with a circle HT.
- //The function loops over all the data, and transforms each pixel with the equations:
+ //The function loops over all the data, and transforms each pixel with the equation:
//
- //kappa = 2/R*sin(phi - phi0)
+ //kappa = 2/R^2 * (y*cos(phi0) - x*sin(phi0) )
//
- //where R = sqrt(x*x +y*y), and phi = arctan(y/x)
+ //where R^2 = x^2 + y^2
//
//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,"AliL3HoughTransformer::TransformCircle","Data")
+ LOG(AliL3Log::kError,"AliL3HoughTransformerLUT::TransformCircle","Data")
<<"No input data "<<ENDLOG;
return;
}
//Loop over the padrows:
- for(Int_t i=AliL3Transform::GetFirstRow(GetPatch()); i<=AliL3Transform::GetLastRow(GetPatch()); i++)
+ 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)
{
- cerr<<"AliL3HoughTransform::TransformCircle : Mismatching padrow numbering "<<i<<" "<<(Int_t)tempPt->fRow<<endl;
+ LOG(AliL3Log::kError,"AliL3HoughTransformerLUT::TransformCircle","Data")
+ <<"AliL3HoughTransformerLUT::TransformCircle : Mismatching padrow numbering "<<i<<" != "<<(Int_t)tempPt->fRow<<endl;
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;
+
+ //check threshold
if((Int_t)charge <= GetLowerThreshold() || (Int_t)charge > GetUpperThreshold())
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())
+
+ UChar_t pad = digPt[j].fPad;
+ UShort_t time = digPt[j].fTime;
+
+ Float_t x = CalcX(row);
+ Float_t y = CalcY(pad,row);
+ Float_t z = CalcZ(time);
+
+ //find eta slice
+ Float_t rz2 = 1 + (x*x+y*y)/(z*z);
+ Int_t eta_index = FindIndex(rz2);
+
+ if(eta_index < 0 || eta_index >= fNEtas){
+ //LOG(AliL3Log::kWarning,"AliL3HoughTransformerLUT::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)
- {
- printf("AliL3HoughTransformer::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);
-
+ if(!hist){
+ //LOG(AliL3Log::kWarning,"AliL3HoughTransformerLUT::TransformCircle","Histograms")<<"Error getting histogram in index "<<eta_index<<"."<<ENDLOG;
+ continue;
+ }
+
//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);
+ Float_t R2=1/(x*x+y*y);
+ for(Int_t b=0; b<fNPhi0; b++){
+ Float_t kappa=R2*(y*fLUT2cosphi0[b]-x*fLUT2sinphi0[b]);
+
+ hist->Fill(kappa,fLUTphi0[b],charge);
+ //cout << kappa << " " << fLUTphi0[b] << " " << charge << endl;
+
#ifdef do_mc
- if(fDoMC)
- {
- Int_t bin = hist->FindBin(kappa,phi0);
- for(Int_t t=0; t<3; t++)
- {
- Int_t label = digPt[j].fTrackID[t];
- if(label < 0) break;
- UInt_t c;
- for(c=0; c<MaxTrack; c++)
- if(fTrackID[eta_index][bin].fLabel[c] == label || fTrackID[eta_index][bin].fNHits[c] == 0)
- break;
- if(c == MaxTrack-1) cerr<<"AliL3HoughTransformer::TransformCircle : Array reached maximum!! "<<c<<endl;
- fTrackID[eta_index][bin].fLabel[c] = label;
- fTrackID[eta_index][bin].fNHits[c]++;
- }
- }
+ if(fDoMC)
+ {
+ Int_t bin = hist->FindBin(kappa,phi0);
+ for(Int_t t=0; t<3; t++)
+ {
+ Int_t label = digPt[j].fTrackID[t];
+ if(label < 0) break;
+ UInt_t c;
+ for(c=0; c<MaxTrack; c++)
+ if(fTrackID[eta_index][bin].fLabel[c] == label || fTrackID[eta_index][bin].fNHits[c] == 0)
+ break;
+ if(c == MaxTrack-1) LOG(AliL3Log::kError,"AliL3HoughTransformerLUT::TransformCircle","MCData") <<"Array reached maximum!! "<<c<<endl;
+ fTrackID[eta_index][bin].fLabel[c] = label;
+ fTrackID[eta_index][bin].fNHits[c]++;
+ }
+ }
#endif
- }
+ }
}
-
//Move the data pointer to the next padrow:
AliL3MemHandler::UpdateRowPointer(tempPt);
}
-#endif //to do later
}
Int_t AliL3HoughTransformerLUT::GetTrackID(Int_t eta_index,Double_t kappa,Double_t psi)
}
#ifdef do_mc
- if(eta_index < 0 || eta_index > GetNEtaSegments())
+ if(eta_index < 0 || eta_index > fNEtas)
{
cerr<<"AliL3HoughTransformer::GetTrackID : Wrong etaindex "<<eta_index<<endl;
return -1;
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff