Moved all the constants from AliL3Defs.h to static AliL3Transform. No
[u/mrichter/AliRoot.git] / HLT / hough / AliL3HoughEval.cxx
CommitLineData
95a00d93 1//$Id$
2
b1886074 3// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
4//*-- Copyright &copy ASV
f80b98cb 5
162e2e5a 6#include <math.h>
95a00d93 7#include <iostream.h>
8#ifdef use_root
b46b21a5 9#include <TH1.h>
162e2e5a 10#include <TFile.h>
95a00d93 11#endif
b46b21a5 12#include "AliL3MemHandler.h"
162e2e5a 13#include "GetGoodParticles.h"
14#include "AliL3TrackArray.h"
15#include "AliL3Logging.h"
99e7186b 16#include "AliL3HoughEval.h"
237d3f5c 17#include "AliL3HoughBaseTransformer.h"
99e7186b 18#include "AliL3DigitData.h"
4de874d1 19#include "AliL3HoughTrack.h"
99e7186b 20#include "AliL3Transform.h"
21#include "AliL3Histogram.h"
b46b21a5 22#include "AliL3Histogram1D.h"
99e7186b 23#include "AliL3Defs.h"
4de874d1 24
b1886074 25//_____________________________________________________________
26// AliL3HoughEval
27//
28// Evaluation class for tracklets produced by the Hough transform.
29
4de874d1 30ClassImp(AliL3HoughEval)
31
4de874d1 32AliL3HoughEval::AliL3HoughEval()
4de874d1 33{
4fc9a6a4 34
99e7186b 35 fRemoveFoundTracks = kFALSE;
f80b98cb 36 fNumOfPadsToLook = 1;
99e7186b 37 fNumOfRowsToMiss = 1;
b46b21a5 38 fEtaHistos=0;
b1886074 39 fRowPointers = 0;
4de874d1 40}
41
b46b21a5 42
4de874d1 43AliL3HoughEval::~AliL3HoughEval()
44{
b46b21a5 45 fHoughTransformer = 0;
99e7186b 46 if(fRowPointers)
b1886074 47 {
48 for(Int_t i=0; i<fNrows; i++)
49 fRowPointers[i] = 0;
50 delete [] fRowPointers;
51 }
4de874d1 52}
53
237d3f5c 54void AliL3HoughEval::InitTransformer(AliL3HoughBaseTransformer *transformer)
b46b21a5 55{
56 fHoughTransformer = transformer;
57 fSlice = fHoughTransformer->GetSlice();
58 fPatch = fHoughTransformer->GetPatch();
59 fNrows = NRows[fPatch][1] - NRows[fPatch][0] + 1;
60 fNEtaSegments = fHoughTransformer->GetNEtaSegments();
61 fEtaMin = fHoughTransformer->GetEtaMin();
62 fEtaMax = fHoughTransformer->GetEtaMax();
63 GenerateLUT();
64}
65
99e7186b 66void AliL3HoughEval::GenerateLUT()
4de874d1 67{
b46b21a5 68 //Generate a Look-up table, to limit the access to raw data
99e7186b 69
b1886074 70 if(!fRowPointers)
71 fRowPointers = new AliL3DigitRowData*[fNrows];
99e7186b 72
73 AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fHoughTransformer->GetDataPointer();
74 if(!tempPt)
b46b21a5 75 printf("\nAliL3HoughEval::GenerateLUT : Zero data pointer\n");
f80b98cb 76
99e7186b 77 for(Int_t i=NRows[fPatch][0]; i<=NRows[fPatch][1]; i++)
f80b98cb 78 {
99e7186b 79 Int_t prow = i - NRows[fPatch][0];
80 fRowPointers[prow] = tempPt;
b46b21a5 81 AliL3MemHandler::UpdateRowPointer(tempPt);
f80b98cb 82 }
83
99e7186b 84}
f80b98cb 85
99e7186b 86Bool_t AliL3HoughEval::LookInsideRoad(AliL3HoughTrack *track,Int_t eta_index,Bool_t remove)
87{
88 //Look at rawdata along the road specified by the track candidates.
89 //If track is good, return true, if not return false.
162e2e5a 90
f80b98cb 91 Int_t sector,row;
99e7186b 92
5de3a7f2 93 Int_t nrow=0,npixs=0,rows_crossed=0;
f80b98cb 94 Float_t xyz[3];
99e7186b 95
b46b21a5 96 Int_t total_charge=0;//total charge along the road
99e7186b 97 Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments;
95a00d93 98
99
100 //Check if the track is leaving the sector at some point
5de3a7f2 101 Float_t maxrow=300;
102 Double_t angle=Pi/18;
103 track->CalculateEdgePoint(angle);
104 if(!track->IsPoint())
105 {
106 track->CalculateEdgePoint(-1.*angle);
107 if(track->IsPoint())
108 maxrow = track->GetPointX();
109 }
110 else
111 maxrow = track->GetPointX();
112
99e7186b 113 for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
f80b98cb 114 {
4ab9f8f0 115 if(AliL3Transform::Row2X(padrow) > maxrow) break;//The track has left this slice
5de3a7f2 116 rows_crossed++;
99e7186b 117 Int_t prow = padrow - NRows[fPatch][0];
f80b98cb 118 if(!track->GetCrossingPoint(padrow,xyz))
119 {
5de3a7f2 120 //printf("AliL3HoughEval::LookInsideRoad : Track does not cross line!!; pt %f phi0 %f\n",track->GetPt(),track->GetPhi0());
f80b98cb 121 continue;
122 }
123
4ab9f8f0 124 AliL3Transform::Slice2Sector(fSlice,padrow,sector,row);
125 AliL3Transform::Local2Raw(xyz,sector,row);
f80b98cb 126 npixs=0;
99e7186b 127
95a00d93 128
99e7186b 129 //Get the timebins for this pad
130 AliL3DigitRowData *tempPt = fRowPointers[prow];
131 if(!tempPt)
f80b98cb 132 {
99e7186b 133 printf("AliL3HoughEval::LookInsideRoad : Zero data pointer\n");
134 continue;
f80b98cb 135 }
136
99e7186b 137 //Look at both sides of the pad:
162e2e5a 138 for(Int_t p=(Int_t)rint(xyz[1])-fNumOfPadsToLook; p<=(Int_t)rint(xyz[1])+fNumOfPadsToLook; p++)
99e7186b 139 {
140 AliL3DigitData *digPt = tempPt->fDigitData;
141 for(UInt_t j=0; j<tempPt->fNDigit; j++)
142 {
b1886074 143 //if(digPt->fCharge <= fHoughTransformer->GetThreshold()) continue;
99e7186b 144 UChar_t pad = digPt[j].fPad;
145 if(pad < p) continue;
146 if(pad > p) break;
147 UShort_t time = digPt[j].fTime;
4ab9f8f0 148 Double_t eta = AliL3Transform::GetEta(padrow,pad,time);
99e7186b 149 Int_t pixel_index = (Int_t)(eta/etaslice);
150 if(pixel_index > eta_index) continue;
151 if(pixel_index != eta_index) break;
b46b21a5 152 total_charge += digPt[j].fCharge;
99e7186b 153 if(remove)
b46b21a5 154 digPt[j].fCharge = 0; //Erease the track from image
99e7186b 155 npixs++;
156 }
157 }
158
162e2e5a 159 if(npixs > 1)//At least 2 digits on this padrow
f80b98cb 160 {
161 nrow++;
162 }
163 }
99e7186b 164 if(remove)
165 return kTRUE;
166
5de3a7f2 167 if(nrow >= rows_crossed - fNumOfRowsToMiss)//this was a good track
f80b98cb 168 {
b46b21a5 169 Double_t eta_track = (Double_t)eta_index*etaslice;
f80b98cb 170 track->SetEtaIndex(eta_index);
b46b21a5 171 track->SetWeight(total_charge,kTRUE);
172 track->SetEta(eta_track);
b1886074 173 track->SetRowRange(NRows[fPatch][0],NRows[fPatch][1]);
174 track->SetSlice(fSlice);
99e7186b 175 if(fRemoveFoundTracks)
176 LookInsideRoad(track,eta_index,kTRUE);
f80b98cb 177 return kTRUE;
178 }
179 else
180 return kFALSE;
181}
182
b46b21a5 183void AliL3HoughEval::FindEta(AliL3TrackArray *tracks)
184{
185
186 Int_t sector,row;
187 Float_t xyz[3];
188
189 Int_t ntracks = tracks->GetNTracks();
190 fEtaHistos = new AliL3Histogram1D*[ntracks];
191
192 Char_t hname[100];
193 for(Int_t i=0; i<ntracks; i++)
194 {
195 sprintf(hname,"etahist_%d",i);
196 fEtaHistos[i] = new AliL3Histogram1D(hname,hname,100,0,1);
197 }
198 Double_t etaslice = (fEtaMax - fEtaMin)/fNEtaSegments;
199
200 for(Int_t ntr=0; ntr<ntracks; ntr++)
201 {
202 AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(ntr);
203 if(!track) continue;
204 for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
205 {
206 Int_t prow = padrow - NRows[fPatch][0];
207
208 if(!track->GetCrossingPoint(padrow,xyz))
209 {
210 printf("AliL3HoughEval::LookInsideRoad : Track does not cross line!!\n");
211 continue;
212 }
213
4ab9f8f0 214 AliL3Transform::Slice2Sector(fSlice,padrow,sector,row);
215 AliL3Transform::Local2Raw(xyz,sector,row);
b46b21a5 216
217 //Get the timebins for this pad
218 AliL3DigitRowData *tempPt = fRowPointers[prow];
219 if(!tempPt)
220 {
221 printf("AliL3HoughEval::LookInsideRoad : Zero data pointer\n");
222 continue;
223 }
224
225 //Look at both sides of the pad:
226 for(Int_t p=(Int_t)rint(xyz[1])-fNumOfPadsToLook; p<=(Int_t)rint(xyz[1])+fNumOfPadsToLook; p++)
227 {
228 AliL3DigitData *digPt = tempPt->fDigitData;
229 for(UInt_t j=0; j<tempPt->fNDigit; j++)
230 {
231 UChar_t pad = digPt[j].fPad;
232
233 if(pad < p) continue;
234 if(pad > p) break;
235 UShort_t time = digPt[j].fTime;
4ab9f8f0 236 Double_t eta = AliL3Transform::GetEta(padrow,pad,time);
b46b21a5 237 Int_t pixel_index = (Int_t)(eta/etaslice);
238 if(pixel_index > track->GetEtaIndex()+1) continue;
239 if(pixel_index < track->GetEtaIndex()-1) break;
240 fEtaHistos[ntr]->Fill(eta,digPt[j].fCharge);
241 }
242 }
243 }
244 }
245
246 for(Int_t i=0; i<ntracks; i++)
247 {
248 AliL3Histogram1D *hist = fEtaHistos[i];
249 Int_t max_bin = hist->GetMaximumBin();
250 Double_t max_value = hist->GetBinContent(max_bin);
251 AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(i);
252 if(!track) continue;
253 if(hist->GetBinContent(max_bin-1)<max_value && hist->GetBinContent(max_bin+1)<max_value)
254 {
255 track->SetWeight((Int_t)max_value,kTRUE);
256 track->SetEta(hist->GetBinCenter(max_bin));
257 track->SetNHits(track->GetWeight());
258 }
259 else
260 {
261 track->SetWeight(0);
262 tracks->Remove(i); //remove this track, because it was not a peak
263 }
264 }
265 tracks->Compress();
266
267 //for(Int_t i=0; i<ntracks; i++)
268 //delete fEtaHistos[i];
269