This commit was generated by cvs2svn to compensate for changes in r3176,
[u/mrichter/AliRoot.git] / HLT / src / AliL3ClustFinder.cxx
CommitLineData
108615fc 1//Author: Anders Strand Vestbo
2//Last Modified: 14.12.2000
3
4#include <iostream.h>
5
6#include "AliL3Logging.h"
7#include "AliL3ClustFinder.h"
8#include "AliL3Transform.h"
9#include "AliL3DigitData.h"
10#include "AliL3SpacePointData.h"
11
12//
13//Simple clusterfinder. Initial version by Tonko.
14//
15
16ClassImp(AliL3ClustFinder)
17
18AliL3ClustFinder::AliL3ClustFinder()
19{
20
21 fDeconvTime=false;
22 fDeconvPad=false;
23 fXYErr = 0;
24 fZErr = 0;
25}
26
27AliL3ClustFinder::AliL3ClustFinder(AliL3Transform *transform)
28{
29
30 fTransform = transform;
31
32 fThreshold = 20;
33 fDeconvTime = true;
34 fDeconvPad = true;
35}
36
37
38AliL3ClustFinder::~AliL3ClustFinder()
39{
40 //destructor
41}
42
43
44 void AliL3ClustFinder::SetTransformer( AliL3Transform *transform )
45 {
46 fTransform = transform;
47 }
48
49
50void AliL3ClustFinder::InitSlice(Int_t slice,Int_t patch,Int_t firstrow, Int_t lastrow,Int_t nmaxpoints)
51{
52 fNClusters = 0;
53 fMaxNClusters = nmaxpoints;
54 fCurrentSlice = slice;
55 fCurrentPatch = patch;
56 fFirstRow = firstrow;
57 fLastRow = lastrow;
58}
59
60
61void AliL3ClustFinder::InitSlice(Int_t slice,Int_t patch,Int_t nmaxpoints)
62{
63 fNClusters = 0;
64 fMaxNClusters = nmaxpoints;
65 fCurrentSlice = slice;
66 fCurrentPatch = patch;
67}
68
69void AliL3ClustFinder::SetOutputArray(AliL3SpacePointData *pt)
70{
71
72 fSpacePointData = pt;
73}
74
75void AliL3ClustFinder::Read(UInt_t ndigits,AliL3DigitRowData *ptr)
76{
77 fNDigitRowData = ndigits;
78 fDigitRowData = ptr;
79
80}
81
82void AliL3ClustFinder::ProcessDigits()
83{
84 //Loop over rows, and call processrow
85
86
87 AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fDigitRowData;
88
89 for(Int_t i=fFirstRow; i<=fLastRow; i++)
90 {
91 fCurrentRow = i;
92 ProcessRow(tempPt);
93 Byte_t *tmp = (Byte_t*)tempPt;
94 Int_t size = sizeof(AliL3DigitRowData) + tempPt->fNDigit*sizeof(AliL3DigitData);
95 tmp += size;
96 tempPt = (AliL3DigitRowData*)tmp;
97 }
98
99 LOG(AliL3Log::kInformational,"AliL3ClustFinder","Cluster Finder")
100 <<AliL3Log::kDec<<"ClusterFinder found "<<fNClusters<<" clusters"<<ENDLOG;
101}
102
103void AliL3ClustFinder::ProcessRow(AliL3DigitRowData* tempPt)
104{
105
106 Int_t cl_found,trackID[3];
107 //UInt_t rows;
108
109 //Reserve array for found clusters:
110 resx aliresx[5000];
111 bzero((char *)aliresx,5000*sizeof(resx));
112
113 //local results
114 resx *r;
115// UInt_t pres1[2500],pres2[2500];
116 resx *pres1[2500];
117 resx *pres2[2500];
118 resx rr_local;
119 resx *rr = &rr_local;
120
121 //Set variables for each row:
122 Int_t pres_cou1,pres_cou2;
123// UInt_t *r1,*r2;
124 resx **r1,**r2;
125 pres_cou1=pres_cou2=0;
126 r1=pres2;
127 r2=pres1;
128 r=aliresx;
129
130 //variables for each pad:
131 UInt_t start;
132 Int_t start_new=-1;
133 Int_t go_on;
134 Int_t cl_counter=0,rows;
135 Int_t prev_start;
136
137 UInt_t last_pad = 123456789;
138
139 //Loop over digits in this row:
140 for(start=0; start<tempPt->fNDigit; start++)
141 {
142 AliL3DigitData *bins = tempPt->fDigitData;
143
144 if(bins[start].fPad!=last_pad)
145 {//new pad
146
147 last_pad=bins[start].fPad;
148
149 if(r2==pres2)
150 {
151 r1=pres2;
152 r2=pres1;
153 }
154 else
155 {
156 r1=pres1;
157 r2=pres2;
158 }
159 pres_cou1=pres_cou2;
160 pres_cou2=0;
161
162 start_new=-1;
163
164 cl_counter=0;
165 prev_start=0;
166
167
168 }//Looping over pads
169
170 UInt_t av,charge;
171 UInt_t mean;
172 UInt_t flags;
173 Int_t last_falling=0;
174// UInt_t *ri;
175 resx **ri;
176 UInt_t tmp_charge;
177
178 //Set flag:
179 flags=0;
180
181 cl_counter++;
182 if(cl_counter>MAX_C) {LOG(AliL3Log::kWarning,"AliL3ClustFinder::ProcessRow","Cluster Finder")
183 <<AliL3Log::kDec<<"Too many clusters"<<ENDLOG; return;} //too many clusters
184
185
186 if(fDeconvTime)
187 {
188
189 //this is a goto label:
190 redo: ;
191 //divide sequenz in two in case there is a change in slope:
192 if(start_new>-1)
193 {
194 start=start_new;
195 start_new=-1;
196 }
197 }
198 //average and charge
199 av=charge=0;
200
201 //ok, create a new cluster:
202 go_on=1;
203
204 if(fDeconvTime)
205 {
206
207 last_falling=flags=0;
208 }
209
210 if(fDeconvPad)
211 {
212 flags=0;
213 }
214
215
216 //block to introduce new variables:
217 {
218 UInt_t last_a;
219 last_a=0;
220
221 //Loop over this sequenz:
222 while(1)
223 {
224 UInt_t aa;
225 //Int_t start_temp=start;
226
227 //get the adc-value:
228 aa=bins[start].fCharge;
229
230
231 if(fDeconvTime)
232 {
233 //Check if the last pixel in this sequenz is bigger or smaller than this:
234 if(aa>last_a)
235 {
236 if(last_falling)
237 {
238 start_new=start;
239 break;
240 }
241 }
242 else last_falling=1;
243 last_a=aa;
244
245 }
246 //sum of total charge of this cluster on this pad
247 charge +=aa;
248 //this one is needed to determine the mean over all pads:
249 //av += start * (aa);
250 av+=bins[start].fTime*(aa);
251
252 if((bins[start+1].fTime)!=(bins[start].fTime+1) ||
253 ((bins[start+1].fPad)!=(bins[start].fPad))) break;
254
255 start++;
256 }//Loop over sequenz
257 }//block of new variables
258
259 if(fDeconvTime)
260 {
261 if(start_new>0) flags = FLAG_DOUBLE_T;
262 }
263 //calculate mean in time direction for this sequenz:
264 if(charge)
265 {
266 //mean=cog for this sequenz:
267 mean=av/charge;
268 }
269 else
270 {
271 charge=1;
272 mean=1;
273 }
274 //get the pointer to results on previous pad(s)
275 ri=r1;
276
277 //Center of gravity in pad direction
278 //tmp_charge=charge*bins[start-1].fPad;
279 tmp_charge=charge*bins[start].fPad;
280
281 //compare with prevously stored results (previous pads)
282 for(int k=0; k<pres_cou1; k++)
283 {
284 //variable to store difference:
285 Int_t v;
286
287 //structure that contains means (+more) of last pads(s)
288 resx *rr_tmp;
289
290 //go get it
291 rr_tmp = (resx *) *ri;
292
293 //difference between mean and mean on previous pads
294 v=mean-rr_tmp->mean;
295
296
297 if(v<-PARAM1) break;
298
299 //goto next mean on previous pad
300 ri++;
301
302 //is there a matching cluster on previous pad:
303 if(v<=PARAM1)
304 {
305 //ok, we got a new one, so we will add the new and the old
306 rr=&rr_local;
307
308 //Add this function somewhere
309 memcpy7((UInt_t *)rr,(UInt_t *)rr_tmp);
310
311 /* in case we found another sequenz-section on this pad
312 don't loop over all sequenzes on last pad but start
313 where we found a matching candidate for this sequenz-
314 section, adjust upper limit (=pres_cou1)*/
315 r1=ri;
316 pres_cou1-=1+k;
317
318 if(fDeconvPad)
319 {
320 if(charge>rr->scharge)
321 {
322 if(rr->falling)
323 {
324 //previous state was falling
325 flags |= FLAG_DOUBLE_PAD;
326 rr_tmp->flags |= flags;
327 //and create a new one
328 break;
329 }
330 }
331 else
332 {
333 rr->falling = 1;
334 }
335
336 }
337 //don't create a new one
338 go_on=0;
339
340 //store the pointer to the matched in "2"
341// r2[pres_cou2++] = (UInt_t) rr_tmp;
342 r2[pres_cou2++] = rr_tmp;
343
344 //calculate and fill the new means, charge etc. in rr_tmp
345 //charge on this pad to determine whether we cut in pad direction
346 rr->scharge = charge;
347 //unset FLAG_ONEPAD since we have at least 2 pads
348 rr->flags &= (~FLAG_ONEPAD);
349 //summ the charges for de/dx
350 rr->charge += charge;
351 //calculate mean in pad direction
352 rr->pad += tmp_charge;
353 //calculate mean in time direction
354 rr->t+=av;
355 //store the mean of this sequenz on this pad
356 rr->mean = mean;
357
358 //Copy it back to storage
359 memcpy7((UInt_t *)rr_tmp,(UInt_t *)rr);
360
361 //we found a matching one so stop looping over pads on previous pad
362 break;
363
364 }//matching cluster on previous pads
365 }//loop: means on previous pad
366
367 //here we create from scratch new clusters
368 if(go_on)
369 {
370 //store this one beacuse it is the first
371// r2[pres_cou2++]= (UInt_t )r;
372 r2[pres_cou2++]= r;
373
374 mstore((UInt_t *)r,av,tmp_charge,charge,flags|FLAG_ONEPAD,mean,trackID);
375 r++;
376 }
377
378 if(fDeconvTime)
379 {
380
381 //maybe we still have to do something
382 if(start_new>=0) goto redo;
383 }
384
385 }//end of loop over digits
386
387
388 //number of clusters found
389 cl_found = r-aliresx;
390
391 //tot_cl_found+=cl_found;
392
393 //there was something on this padrow
394 rows++;
395
396 WriteClusters(cl_found,aliresx);
397
398}
399
400void AliL3ClustFinder::memcpy7(UInt_t *dst, UInt_t *src)
401{
402 Int_t i;
403 for(i=0; i<7; i++)
404 {
405 *dst++ = *src++;
406 }
407 return;
408}
409
410void AliL3ClustFinder::mstore(UInt_t *r,UInt_t av,UInt_t pad,UInt_t ch,UInt_t flags,UInt_t mean,Int_t *trackID)
411{
412 resx *rr = (resx *) r;
413
414 rr->mean = mean;
415 rr->charge = ch;
416 rr->flags =flags;
417 rr->scharge = ch;
418 rr->pad = pad;
419 rr->t = av;
420 rr->falling=0;
421
422 rr->trackID[0]=trackID[0];
423 rr->trackID[1]=trackID[1];
424 rr->trackID[2]=trackID[2];
425 return;
426}
427
428void AliL3ClustFinder::WriteClusters(Int_t ncl,resx *r)
429{
430/*
431 if(fNClusters >= fMaxNClusters)
432 {
433 LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Finder")
434 <<AliL3Log::kDec<<"Too many clusters"<<ENDLOG;
435 return;
436 }
437*/
438 Int_t thisrow,thissector;
439 UInt_t counter = fNClusters;
440
441 if(!fXYErr || !fZErr)
442 {
443 LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Errors")
444 <<"No Errors"<<ENDLOG;
445 return;
446 }
447
448 for(int j=0; j<ncl; j++)
449 {
450 if(r[j].flags & FLAG_ONEPAD) continue; //discard 1 pad clusters
451 if(r[j].charge < fThreshold) continue; //noise cluster
452 Float_t xyz[3];
453
454 Float_t fpad=(Float_t)r[j].pad/(Float_t)r[j].charge;
455 Float_t ftime=(Float_t)r[j].t/(Float_t)r[j].charge;
456
457 if(fCurrentRow > 54) {thisrow = fCurrentRow-55; thissector = fCurrentSlice+36;}
458 else {thisrow = fCurrentRow; thissector = fCurrentSlice;}
459 fTransform->Raw2Local(xyz,thissector,thisrow,fpad,ftime);
460 if(xyz[0]==0) LOG(AliL3Log::kError,"AliL3ClustFinder","Cluster Finder")
461 <<AliL3Log::kDec<<"Zero cluster"<<ENDLOG;
462 if(fNClusters >= fMaxNClusters)
463 {
464 LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Finder")
465 <<AliL3Log::kDec<<"Too many clusters"<<ENDLOG;
466 return;
467 }
468 fSpacePointData[counter].fX = xyz[0];
469 fSpacePointData[counter].fY = xyz[1];
470 fSpacePointData[counter].fZ = xyz[2];
471 fSpacePointData[counter].fPadRow = fCurrentRow;
472 fSpacePointData[counter].fXYErr = fXYErr;
473 fSpacePointData[counter].fZErr = fZErr;
474 fSpacePointData[counter].fID = counter
475 +((fCurrentSlice&0x7f)<<25)+((fCurrentPatch&0x7)<<22);//uli
476
477
478 fNClusters++;
479 counter++;
480
481 }
482
483
484}