4ec96209 |
1 | //$Id$ |
2 | |
b661165c |
3 | // Author: Anders Vestbo <mailto:vestbo@fi.uib.no> |
4 | //*-- Copyright © ASV |
46fbeb8a |
5 | |
6 | #include <iostream.h> |
10f815d9 |
7 | #include <math.h> |
46fbeb8a |
8 | #include "AliL3Logging.h" |
9 | #include "AliL3ClustFinderNew.h" |
10 | #include "AliL3DigitData.h" |
11 | #include "AliL3Transform.h" |
12 | #include "AliL3SpacePointData.h" |
10f815d9 |
13 | #include "AliL3MemHandler.h" |
46fbeb8a |
14 | |
b661165c |
15 | //_____________________________________________________________ |
16 | // AliL3ClustFinderNew |
17 | // |
18 | // The current cluster finder for HLT |
19 | // Based on STAR L3 |
20 | |
46fbeb8a |
21 | ClassImp(AliL3ClustFinderNew) |
22 | |
23 | AliL3ClustFinderNew::AliL3ClustFinderNew() |
24 | { |
25 | fMatch = 4; |
26 | fThreshold =10; |
99a6d5c9 |
27 | fDeconvPad = kTRUE; |
2bbf57d8 |
28 | fDeconvTime = kTRUE; |
46fbeb8a |
29 | } |
30 | |
46fbeb8a |
31 | AliL3ClustFinderNew::~AliL3ClustFinderNew() |
32 | { |
33 | |
34 | |
35 | } |
36 | |
37 | void AliL3ClustFinderNew::InitSlice(Int_t slice,Int_t patch,Int_t firstrow, Int_t lastrow,Int_t nmaxpoints) |
38 | { |
39 | fNClusters = 0; |
40 | fMaxNClusters = nmaxpoints; |
41 | fCurrentSlice = slice; |
42 | fCurrentPatch = patch; |
43 | fFirstRow = firstrow; |
44 | fLastRow = lastrow; |
99a6d5c9 |
45 | fDeconvTime = kTRUE; |
46 | fDeconvPad = kTRUE; |
46fbeb8a |
47 | } |
48 | |
49 | void AliL3ClustFinderNew::InitSlice(Int_t slice,Int_t patch,Int_t nmaxpoints) |
50 | { |
51 | fNClusters = 0; |
52 | fMaxNClusters = nmaxpoints; |
53 | fCurrentSlice = slice; |
54 | fCurrentPatch = patch; |
55 | } |
56 | |
57 | void AliL3ClustFinderNew::SetOutputArray(AliL3SpacePointData *pt) |
58 | { |
59 | |
60 | fSpacePointData = pt; |
61 | } |
62 | |
63 | |
64 | void AliL3ClustFinderNew::Read(UInt_t ndigits,AliL3DigitRowData *ptr) |
65 | { |
66 | fNDigitRowData = ndigits; |
67 | fDigitRowData = ptr; |
68 | |
69 | } |
70 | |
71 | void AliL3ClustFinderNew::ProcessDigits() |
72 | { |
73 | //Loop over rows, and call processrow |
74 | |
75 | |
76 | AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fDigitRowData; |
77 | |
78 | for(Int_t i=fFirstRow; i<=fLastRow; i++) |
79 | { |
80 | fCurrentRow = i; |
81 | ProcessRow(tempPt); |
82 | Byte_t *tmp = (Byte_t*)tempPt; |
83 | Int_t size = sizeof(AliL3DigitRowData) + tempPt->fNDigit*sizeof(AliL3DigitData); |
84 | tmp += size; |
85 | tempPt = (AliL3DigitRowData*)tmp; |
86 | } |
355debe1 |
87 | LOG(AliL3Log::kInformational,"AliL3ClustFinderNew::WriteClusters","Space points") |
88 | <<"Cluster finder found "<<fNClusters<<" clusters in slice "<<fCurrentSlice |
89 | <<" patch "<<fCurrentPatch<<ENDLOG; |
46fbeb8a |
90 | } |
91 | |
46fbeb8a |
92 | void AliL3ClustFinderNew::ProcessRow(AliL3DigitRowData *tempPt) |
93 | { |
94 | |
95 | UInt_t last_pad = 123456789; |
96 | |
97 | ClusterData *pad1[2500]; //2 lists for internal memory=2pads |
98 | ClusterData *pad2[2500]; //2 lists for internal memory=2pads |
99 | ClusterData clusterlist[5000]; //Clusterlist |
100 | |
101 | ClusterData **currentPt; //List of pointers to the current pad |
102 | ClusterData **previousPt; //List of pointers to the previous pad |
103 | currentPt = pad2; |
104 | previousPt = pad1; |
105 | UInt_t n_previous=0,n_current=0,n_total=0; |
106 | |
107 | //Loop over sequences of this row: |
108 | for(UInt_t bin=0; bin<tempPt->fNDigit; bin++) |
109 | { |
110 | |
111 | AliL3DigitData *data = tempPt->fDigitData; |
112 | if(data[bin].fPad != last_pad) |
113 | { |
114 | //This is a new pad |
46fbeb8a |
115 | |
116 | //Switch: |
117 | if(currentPt == pad2) |
118 | { |
119 | currentPt = pad1; |
120 | previousPt = pad2; |
121 | |
122 | } |
123 | else |
124 | { |
125 | currentPt = pad2; |
126 | previousPt = pad1; |
127 | } |
128 | n_previous = n_current; |
129 | n_current = 0; |
99a6d5c9 |
130 | if(bin[data].fPad != last_pad+1) |
131 | { |
132 | //this happens if there is a pad with no signal. |
133 | n_previous = n_current = 0; |
134 | } |
135 | last_pad = data[bin].fPad; |
46fbeb8a |
136 | } |
137 | |
138 | Bool_t new_cluster = kTRUE; |
139 | |
140 | UInt_t seq_charge=0,seq_average=0; |
141 | |
99a6d5c9 |
142 | UInt_t last_charge=0,last_was_falling=0; |
143 | Int_t new_bin=-1; |
144 | if(fDeconvTime) |
145 | { |
146 | redo: //This is a goto. |
147 | if(new_bin > -1) |
148 | { |
149 | bin = new_bin; |
150 | new_bin = -1; |
151 | } |
152 | |
153 | last_charge=0; |
154 | last_was_falling = 0; |
155 | } |
156 | |
157 | |
46fbeb8a |
158 | while(1) //Loop over current sequence |
159 | { |
494fad94 |
160 | if(data[bin].fTime >= AliL3Transform::GetNTimeBins()) |
adb0e845 |
161 | { |
162 | LOG(AliL3Log::kFatal,"AliL3ClustFinderNew::ProcessRow","Digits") |
163 | <<"Timebin out of range "<<(Int_t)data[bin].fTime<<ENDLOG; |
164 | break; |
165 | } |
166 | |
46fbeb8a |
167 | //Get the current ADC-value |
168 | UInt_t charge = data[bin].fCharge; |
99a6d5c9 |
169 | |
170 | if(fDeconvTime) |
171 | { |
172 | //Check if the last pixel in the sequence is smaller than this |
173 | if(charge > last_charge) |
174 | { |
175 | if(last_was_falling) |
176 | { |
177 | new_bin = bin; |
178 | break; |
179 | } |
180 | } |
181 | else last_was_falling = 1; //last pixel was larger than this |
182 | last_charge = charge; |
183 | } |
184 | |
46fbeb8a |
185 | //Sum the total charge of this sequence |
186 | seq_charge += charge; |
46fbeb8a |
187 | seq_average += data[bin].fTime*charge; |
188 | |
189 | //Check where to stop: |
190 | if(data[bin+1].fPad != data[bin].fPad) //new pad |
191 | break; |
192 | if(data[bin+1].fTime != data[bin].fTime+1) //end of sequence |
193 | break; |
194 | |
195 | bin++; |
196 | }//end loop over sequence |
197 | |
198 | //Calculate mean of sequence: |
199 | Int_t seq_mean=0; |
200 | if(seq_charge) |
201 | seq_mean = seq_average/seq_charge; |
202 | else |
99a6d5c9 |
203 | { |
adb0e845 |
204 | LOG(AliL3Log::kFatal,"AliL3ClustFinderNew::ProcessRow","Data") |
205 | <<"Error in data given to the cluster finder"<<ENDLOG; |
99a6d5c9 |
206 | seq_mean = 1; |
207 | seq_charge = 1; |
208 | } |
46fbeb8a |
209 | |
210 | //Calculate mean in pad direction: |
211 | Int_t pad_mean = seq_charge*data[bin].fPad; |
212 | |
213 | //Compare with results on previous pad: |
214 | for(UInt_t p=0; p<n_previous; p++) |
215 | { |
216 | Int_t difference = seq_mean - previousPt[p]->fMean; |
217 | if(difference < -fMatch) |
218 | break; |
219 | |
220 | if(difference <= fMatch)//There is a match here!! |
221 | { |
99a6d5c9 |
222 | |
223 | ClusterData *local = previousPt[p]; |
224 | if(fDeconvPad) |
225 | { |
226 | if(seq_charge > local->fLastCharge) |
227 | { |
228 | if(local->fChargeFalling)//The previous pad was falling |
229 | { |
230 | break;//create a new cluster |
231 | } |
232 | } |
233 | else |
234 | local->fChargeFalling = 1; |
235 | local->fLastCharge = seq_charge; |
236 | } |
237 | |
46fbeb8a |
238 | //Don't create a new cluster, because we found a match |
239 | new_cluster = kFALSE; |
240 | |
241 | //Update cluster on current pad with the matching one: |
99a6d5c9 |
242 | //currentPt[n_current] = previousPt[p]; |
243 | |
244 | local->fTotalCharge += seq_charge; |
245 | local->fPad += pad_mean; |
246 | local->fTime += seq_average; |
247 | local->fMean = seq_mean; |
248 | local->fFlags++; //means we have more than one pad |
249 | |
250 | currentPt[n_current] = local; |
46fbeb8a |
251 | n_current++; |
252 | |
46fbeb8a |
253 | break; |
254 | }//Checking for match at previous pad |
255 | }//Loop over results on previous pad. |
256 | |
257 | if(new_cluster) |
258 | { |
259 | //Start a new cluster. Add it to the clusterlist, and update |
260 | //the list of pointers to clusters in current pad. |
99a6d5c9 |
261 | //current pad will be previous pad on next pad. |
46fbeb8a |
262 | |
263 | //Add to the clusterlist: |
264 | ClusterData *tmp = &clusterlist[n_total]; |
265 | tmp->fTotalCharge = seq_charge; |
266 | tmp->fPad = pad_mean; |
267 | tmp->fTime = seq_average; |
268 | tmp->fMean = seq_mean; |
269 | tmp->fFlags = 0; //flags for 1 pad clusters |
99a6d5c9 |
270 | if(fDeconvPad) |
271 | { |
272 | tmp->fChargeFalling = 0; |
273 | tmp->fLastCharge = seq_charge; |
274 | } |
275 | |
46fbeb8a |
276 | //Update list of pointers to previous pad: |
277 | currentPt[n_current] = &clusterlist[n_total]; |
278 | n_total++; |
279 | n_current++; |
280 | } |
99a6d5c9 |
281 | if(fDeconvTime) |
282 | if(new_bin >= 0) goto redo; |
283 | |
46fbeb8a |
284 | }//Loop over digits on this padrow |
285 | |
286 | WriteClusters(n_total,clusterlist); |
46fbeb8a |
287 | } |
288 | |
289 | void AliL3ClustFinderNew::WriteClusters(Int_t n_clusters,ClusterData *list) |
290 | { |
291 | Int_t thisrow,thissector; |
292 | UInt_t counter = fNClusters; |
293 | |
294 | for(int j=0; j<n_clusters; j++) |
295 | { |
296 | if(!list[j].fFlags) continue; //discard 1 pad clusters |
297 | if(list[j].fTotalCharge < fThreshold) continue; //noise cluster |
d3b44290 |
298 | |
299 | Float_t xyz[3]; |
46fbeb8a |
300 | Float_t fpad=(Float_t)list[j].fPad/(Float_t)list[j].fTotalCharge; |
301 | Float_t ftime=(Float_t)list[j].fTime/(Float_t)list[j].fTotalCharge; |
eeddc64d |
302 | |
d3b44290 |
303 | //cout<<"WriteCluster: padrow "<<fCurrentRow<<" pad "<<fpad<<" time "<<ftime<<" charge "<<list[j].fTotalCharge<<endl; |
304 | |
494fad94 |
305 | AliL3Transform::Slice2Sector(fCurrentSlice,fCurrentRow,thissector,thisrow); |
306 | AliL3Transform::Raw2Local(xyz,thissector,thisrow,fpad,ftime); |
46fbeb8a |
307 | if(xyz[0]==0) LOG(AliL3Log::kError,"AliL3ClustFinder","Cluster Finder") |
308 | <<AliL3Log::kDec<<"Zero cluster"<<ENDLOG; |
309 | if(fNClusters >= fMaxNClusters) |
310 | { |
311 | LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Finder") |
312 | <<AliL3Log::kDec<<"Too many clusters"<<ENDLOG; |
313 | return; |
d3b44290 |
314 | } |
95a00d93 |
315 | fSpacePointData[counter].fCharge = list[j].fTotalCharge; |
46fbeb8a |
316 | fSpacePointData[counter].fX = xyz[0]; |
317 | fSpacePointData[counter].fY = xyz[1]; |
318 | fSpacePointData[counter].fZ = xyz[2]; |
319 | fSpacePointData[counter].fPadRow = fCurrentRow; |
320 | fSpacePointData[counter].fXYErr = fXYErr; |
321 | fSpacePointData[counter].fZErr = fZErr; |
322 | fSpacePointData[counter].fID = counter |
d3b44290 |
323 | +((fCurrentSlice&0x7f)<<25)+((fCurrentPatch&0x7)<<22);//Uli |
10f815d9 |
324 | #ifdef do_mc |
325 | Int_t trackID[3]; |
326 | GetTrackID((Int_t)rint(fpad),(Int_t)rint(ftime),trackID); |
6ab8a12a |
327 | //cout<<"padrow "<<fCurrentRow<<" pad "<<(Int_t)rint(fpad)<<" time "<<(Int_t)rint(ftime)<<" Trackid "<<trackID[0]<<endl; |
10f815d9 |
328 | fSpacePointData[counter].fTrackID[0] = trackID[0]; |
329 | fSpacePointData[counter].fTrackID[1] = trackID[1]; |
330 | fSpacePointData[counter].fTrackID[2] = trackID[2]; |
331 | #endif |
46fbeb8a |
332 | fNClusters++; |
333 | counter++; |
46fbeb8a |
334 | } |
46fbeb8a |
335 | } |
10f815d9 |
336 | |
337 | #ifdef do_mc |
338 | void AliL3ClustFinderNew::GetTrackID(Int_t pad,Int_t time,Int_t *trackID) |
339 | { |
340 | AliL3DigitRowData *rowPt = (AliL3DigitRowData*)fDigitRowData; |
341 | |
355debe1 |
342 | trackID[0]=trackID[1]=trackID[2]=-2; |
10f815d9 |
343 | //cout<<"Looking for pad "<<pad<<" time "<<time<<endl; |
344 | for(Int_t i=fFirstRow; i<=fLastRow; i++) |
345 | { |
346 | if(rowPt->fRow < (UInt_t)fCurrentRow) |
347 | { |
348 | AliL3MemHandler::UpdateRowPointer(rowPt); |
349 | continue; |
350 | } |
351 | AliL3DigitData *digPt = (AliL3DigitData*)rowPt->fDigitData; |
352 | for(UInt_t j=0; j<rowPt->fNDigit; j++) |
353 | { |
354 | Int_t cpad = digPt[j].fPad; |
355 | Int_t ctime = digPt[j].fTime; |
355debe1 |
356 | if(cpad != pad) continue; |
357 | if(ctime != time) continue; |
358 | //if(cpad != pad && ctime != ctime) continue; |
10f815d9 |
359 | //cout<<"Reading row "<<fCurrentRow<<" pad "<<cpad<<" time "<<ctime<<" trackID "<<digPt[j].fTrackID[0]<<endl; |
360 | trackID[0] = digPt[j].fTrackID[0]; |
361 | trackID[1] = digPt[j].fTrackID[1]; |
362 | trackID[2] = digPt[j].fTrackID[2]; |
355debe1 |
363 | break; |
10f815d9 |
364 | //cout<<"Reading trackID "<<trackID[0]<<endl; |
365 | } |
366 | break; |
367 | } |
368 | |
369 | } |
370 | #endif |