]>
Commit | Line | Data |
---|---|---|
3e87ef69 | 1 | // @(#) $Id$ |
1e843b89 | 2 | |
3e87ef69 | 3 | // Author: Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de> |
4 | //*-- Copyright & copy ALICE HLT Group | |
8e1f5064 | 5 | |
24dbb695 | 6 | #include "AliL3StandardIncludes.h" |
3e87ef69 | 7 | #include "AliL3RootTypes.h" |
8 | #include "AliL3Logging.h" | |
9 | #include "AliL3AltroMemHandler.h" | |
8e1f5064 | 10 | |
3e87ef69 | 11 | //#define VHDLDEBUG |
8e1f5064 | 12 | #include "AliL3VHDLClusterFinder.h" |
13 | ||
0bd0c1ef | 14 | #if __GNUC__ == 3 |
24dbb695 | 15 | using namespace std; |
16 | #endif | |
17 | ||
8e1f5064 | 18 | /** \class AliL3VHDLClusterFinder |
3e87ef69 | 19 | <pre> |
8e1f5064 | 20 | //____________________________________________________ |
21 | // AliL3VHDLClusterFinder | |
22 | // | |
23 | // The current VHDL cluster finder for HLT | |
24 | // Based on STAR L3 | |
25 | // | |
26 | // Most important parameters: | |
27 | // fThreshold - threshold for noise clusters | |
28 | // fMatch - length in time for overlapping sequences | |
3e87ef69 | 29 | </pre> |
6aec85a4 | 30 | */ |
8e1f5064 | 31 | |
c8e2510a | 32 | |
8e1f5064 | 33 | ClassImp(AliL3VHDLClusterFinder) |
34 | ||
35 | AliL3VHDLClusterFinder::AliL3VHDLClusterFinder() | |
36 | { | |
37 | fMatch = 4; | |
38 | fThreshold = 10; | |
39 | fMinMerge = 1; | |
40 | fNClusters=0; | |
41 | ||
42 | fXYErr = 0.2; | |
43 | fZErr = 0.3; | |
44 | fDeconvPad = kTRUE; | |
45 | fDeconvTime = kTRUE; | |
46 | fstdout = kFALSE; | |
47 | fcalcerr = kTRUE; | |
48 | ||
49 | Clear(); | |
3e87ef69 | 50 | #ifdef VHDLDEBUG |
8e1f5064 | 51 | fdeb=fopen("vhdlclusterfinder.debug","w"); |
52 | //fdeb=stderr; | |
53 | #endif | |
54 | } | |
55 | ||
56 | AliL3VHDLClusterFinder::~AliL3VHDLClusterFinder() | |
57 | { | |
3e87ef69 | 58 | #ifdef VHDLDEBUG |
8e1f5064 | 59 | fclose(fdeb); |
60 | #endif | |
61 | } | |
62 | ||
63 | void AliL3VHDLClusterFinder::ProcessDigits() | |
64 | { | |
65 | //Loop over data like the analyzer of the VHDL code | |
66 | const UChar_t n=255; | |
67 | UShort_t rrow=0,rtime=0; | |
68 | UChar_t rpad=0,i=n; | |
69 | UShort_t *charges=new UShort_t[n]; | |
8e1f5064 | 70 | |
71 | fNClusters=0; | |
72 | fRow=0; | |
73 | fPad=0; | |
74 | Clear(); | |
75 | ||
76 | //loop over input data | |
77 | while(fAltromem.ReadSequence(rrow,rpad,rtime,i,&charges)){ | |
8e1f5064 | 78 | #if 0 |
79 | cout << "Padrow " << (int)rrow << " pad " << (int)rpad << " time " <<(int)rtime << " charges "; | |
80 | for(UChar_t ii=0;ii<i;ii++) cout << (int)charges[ii] << " "; | |
81 | cout << endl; | |
82 | #endif | |
3e87ef69 | 83 | #ifdef VHDLDEBUG |
8e1f5064 | 84 | fprintf(fdeb,"ProcessDigits: Input Data: %d %d %d charges:",(int)rrow,(int)rpad,(int)rtime); |
85 | for(UChar_t ii=0;ii<i;ii++) fprintf(fdeb," %d",(int)charges[ii]); | |
86 | fprintf(fdeb,"\n"); | |
87 | #endif | |
88 | ||
89 | fNRow=rrow; | |
90 | fNPad=rpad; | |
c8e2510a | 91 | fTC=0,fMT=0,fST=0; |
8e1f5064 | 92 | |
93 | //calculate sequence values | |
c8e2510a | 94 | if(fDeconvTime){ |
95 | UChar_t ii=0; | |
96 | Int_t sl=0; | |
97 | Int_t charge=0,lcharge=0; | |
98 | Bool_t falling=kFALSE; | |
99 | while(ii<i){ | |
100 | charge=charges[ii]; | |
101 | if((falling)&&(charge>lcharge)){ | |
102 | fSM=rtime-sl/2; | |
103 | MakeSequence(); | |
104 | ProcessSequence(); | |
105 | OutputMemory(); | |
106 | ||
107 | rtime=rtime-sl; | |
108 | falling=kFALSE; | |
109 | sl=0; | |
110 | fTC=0,fMT=0,fST=0; | |
111 | } else if(charge<lcharge) falling=kTRUE; | |
112 | ||
113 | fTC+=charge; | |
114 | fMT+=(rtime-sl)*charge; | |
115 | fST+=(rtime-sl)*(rtime-sl)*charge; | |
116 | sl++; | |
117 | ii++; | |
118 | lcharge=charge; | |
119 | } //end loop over sequence | |
120 | fSM=rtime-sl/2; | |
121 | } else { /* no deconvolution */ | |
122 | for(UChar_t ii=0;ii<i;ii++){ | |
123 | fTC+=charges[ii]; | |
124 | fMT+=(rtime-ii)*charges[ii]; | |
125 | fST+=(rtime-ii)*(rtime-ii)*charges[ii]; | |
126 | } | |
127 | fSM=rtime-i/2; | |
8e1f5064 | 128 | } |
c8e2510a | 129 | |
130 | MakeSequence(); | |
8e1f5064 | 131 | ProcessSequence(); |
c8e2510a | 132 | OutputMemory(); |
8e1f5064 | 133 | |
3e87ef69 | 134 | #ifdef VHDLDEBUG |
8e1f5064 | 135 | fflush(fdeb); |
136 | #endif | |
137 | i=n; //store size of charges array | |
138 | } //loop over data | |
139 | ||
140 | //flush everything left | |
141 | FlushMemory(); | |
142 | while(fOP!=fFP) OutputMemory(); | |
143 | } | |
144 | ||
c8e2510a | 145 | void AliL3VHDLClusterFinder::MakeSequence(){ |
146 | if(!fTC) return; | |
147 | ||
148 | Int_t mp=fNPad*fTC; | |
149 | Int_t sp=fNPad*fNPad*fTC; | |
150 | ||
151 | fSeq.fTotalCharge=fTC; | |
152 | fSeq.fPad=mp; | |
153 | fSeq.fTime=fMT; | |
154 | fSeq.fPad2=sp; | |
155 | fSeq.fTime2=fST; | |
156 | fSeq.fMean=0; | |
157 | fSeq.fMean=fSM; | |
158 | fSeq.fMerge=0; | |
159 | fSeq.fRow=fNRow; | |
160 | fSeq.fLastPad=fNPad; | |
161 | fSeq.fChargeFalling=0; | |
162 | if(fDeconvPad) fSeq.fLastCharge=fTC; | |
163 | else fSeq.fLastCharge=0; | |
164 | } | |
165 | ||
8e1f5064 | 166 | void AliL3VHDLClusterFinder::ProcessSequence() |
167 | { | |
168 | if(fNRow!=fRow) FlushMemory(); | |
169 | else if(fNPad==fPad+1) PrepareMemory(); | |
170 | else if(fNPad!=fPad) FlushMemory(); | |
171 | ||
172 | //store new row and pad values | |
173 | fRow=fNRow; | |
174 | fPad=fNPad; | |
175 | ||
3e87ef69 | 176 | #ifdef VHDLDEBUG |
8e1f5064 | 177 | fprintf(fdeb,"ProcessSequence: Mean=%d TC=%d ",fSeq.fMean,fSeq.fTotalCharge); |
178 | #endif | |
179 | ||
180 | CompareSeq(); //merge or insert | |
181 | } | |
182 | ||
183 | void AliL3VHDLClusterFinder::PrepareMemory() | |
184 | { | |
3e87ef69 | 185 | #ifdef VHDLDEBUG |
8e1f5064 | 186 | fprintf(fdeb,"PrepareMemory %d %d %d\n",fRP,fEP,fWP); |
187 | #endif | |
188 | fRP=fEP; | |
189 | fEP=fWP; | |
190 | } | |
191 | ||
192 | void AliL3VHDLClusterFinder::FlushMemory() | |
193 | { | |
3e87ef69 | 194 | #ifdef VHDLDEBUG |
8e1f5064 | 195 | fprintf(fdeb,"FlushMemory %d %d %d %d\n",fFP,fRP,fEP,fWP); |
196 | #endif | |
197 | fFP=fWP; | |
198 | fRP=fWP; | |
199 | fEP=fWP; | |
200 | } | |
201 | ||
202 | void AliL3VHDLClusterFinder::CompareSeq() | |
203 | { | |
8e1f5064 | 204 | while(fRP!=fEP){ |
205 | Int_t diff=fSeqs[fPList[fRP]].fMean-fSeq.fMean; | |
206 | ||
207 | if(diff>fMatch){ //no match | |
208 | IncRPointer(); //cluster finished | |
209 | continue; | |
210 | } else if(diff<-fMatch){ //no match | |
211 | break; //insert new cluster | |
c8e2510a | 212 | } else { //match found, merge it |
8e1f5064 | 213 | MergeSeq(); |
214 | return; | |
215 | } | |
216 | } | |
217 | ||
218 | InsertSeq(); //start new cluster | |
219 | } | |
220 | ||
221 | void AliL3VHDLClusterFinder::MergeSeq() | |
222 | { | |
3e87ef69 | 223 | #ifdef VHDLDEBUG |
6aec85a4 | 224 | fprintf(fdeb,"merged with Mean=%d TC=%d (new Merge=%d) %d %d\n",fSeqs[fPList[fRP]].fMean,fSeqs[fPList[fRP]].fTotalCharge,fSeqs[fPList[fRP]].fMerge+1,fRow,fPad); |
8e1f5064 | 225 | #endif |
fbd800c2 | 226 | if(fSeqs[fPList[fRP]].fRow!=fSeq.fRow){ |
8e1f5064 | 227 | LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::","Memory Check") |
228 | <<"Sequences can be merged on the same rows only."<<ENDLOG; | |
229 | } | |
230 | if(fSeqs[fPList[fRP]].fLastPad+1!=fSeq.fLastPad){ | |
231 | LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::","Memory Check") | |
232 | <<"Sequences can be merged on consecutive pads only."<<ENDLOG; | |
233 | } | |
c8e2510a | 234 | |
235 | if(fDeconvPad){ | |
236 | if(fSeq.fTotalCharge > fSeqs[fPList[fRP]].fLastCharge){ | |
237 | if(fSeqs[fPList[fRP]].fChargeFalling){ //The previous pad was falling | |
238 | IncRPointer(); | |
239 | InsertSeq(); //start a new cluster | |
240 | return; | |
241 | } | |
242 | } | |
243 | else fSeqs[fPList[fRP]].fChargeFalling = 1; | |
8e1f5064 | 244 | |
c8e2510a | 245 | fSeqs[fPList[fRP]].fLastCharge = fSeq.fTotalCharge; |
246 | } | |
247 | ||
8e1f5064 | 248 | fSeqs[fPList[fRP]].fMean=fSeq.fMean; //take the new mean |
249 | fSeqs[fPList[fRP]].fLastPad=fSeq.fLastPad; | |
250 | fSeqs[fPList[fRP]].fTotalCharge+=fSeq.fTotalCharge; | |
251 | fSeqs[fPList[fRP]].fPad+=fSeq.fPad; | |
252 | fSeqs[fPList[fRP]].fTime+=fSeq.fTime; | |
253 | fSeqs[fPList[fRP]].fPad2+=fSeq.fPad2; | |
254 | fSeqs[fPList[fRP]].fTime2+=fSeq.fTime2; | |
255 | fSeqs[fPList[fRP]].fMerge+=1; | |
256 | ||
257 | fPList[fWP]=fPList[fRP]; | |
258 | fPList[fRP]=N_clmem; //mark it to be free | |
259 | ||
260 | IncRPointer(); | |
261 | IncWPointer(); | |
262 | } | |
263 | ||
264 | void AliL3VHDLClusterFinder::InsertSeq() | |
265 | { | |
3e87ef69 | 266 | #ifdef VHDLDEBUG |
6aec85a4 | 267 | fprintf(fdeb,"inserted %d %d\n",fRow,fPad); |
8e1f5064 | 268 | #endif |
269 | NextFreeIndex(); //get next index | |
270 | fSeqs[fFirst]=fSeq; //store data | |
271 | fPList[fWP]=fFirst; //store index | |
272 | IncWPointer(); | |
273 | } | |
274 | ||
275 | void AliL3VHDLClusterFinder::OutputMemory() | |
276 | { | |
277 | Float_t mtime=0,mpad=0; | |
278 | Float_t mtime2=0,mpad2=0; | |
279 | UInt_t tc,row,mno; | |
280 | ||
281 | if(fOP!=fFP){ | |
282 | //while(fOP!=fFP){ | |
283 | ||
284 | UInt_t index=fPList[fOP]; | |
285 | fPList[fOP]=N_clmem; | |
286 | //cout << fOP << " - " << fFP << " - " << index << endl; | |
287 | IncPointer(fOP); | |
288 | if(index>=N_clmem) return; //nothing to do | |
289 | //if(index>=N_clmem) continue; //nothing to do | |
290 | ||
291 | tc=fSeqs[index].fTotalCharge; | |
292 | mno=fSeqs[index].fMerge; | |
293 | row=fSeqs[index].fRow; | |
294 | if(tc!=0){ | |
295 | mtime=(Float_t)(fSeqs[index].fTime)/tc; | |
296 | mtime2=sqrt((Float_t)(fSeqs[index].fTime2)/tc-mtime*mtime); | |
297 | } | |
298 | if(tc!=0){ | |
299 | mpad=(Float_t)(fSeqs[index].fPad)/tc; | |
300 | mpad2=sqrt((Float_t)(fSeqs[index].fPad2)/tc-mpad*mpad); | |
301 | } | |
302 | ||
303 | if(mno<fMinMerge){ //noise cluster | |
3e87ef69 | 304 | #ifdef VHDLDEBUG |
8e1f5064 | 305 | fprintf(fdeb,"OutputMemory %d %d: noise cluster (merge): Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fOP,fFP,row,mpad,mtime,tc,mno); |
306 | #endif | |
307 | FreeSeq(index); | |
308 | return; | |
309 | //continue; | |
310 | } | |
311 | ||
312 | if(tc<fThreshold){ //total charge below threshold | |
3e87ef69 | 313 | #ifdef VHDLDEBUG |
8e1f5064 | 314 | fprintf(fdeb,"OutputMemory %d %d: noise cluster (threshold): Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fOP,fFP,row,mpad,mtime,tc,mno); |
315 | #endif | |
316 | FreeSeq(index); | |
317 | return; | |
318 | //continue; | |
319 | } | |
320 | ||
3e87ef69 | 321 | #ifdef VHDLDEBUG |
8e1f5064 | 322 | fprintf(fdeb,"OutputMemory %d: Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fNClusters,row,mpad,mtime,tc,mno); |
323 | #endif | |
324 | ||
325 | if(stdout) | |
326 | cout<<"WriteCluster: padrow "<<row<<" pad "<<mpad<< " +- "<<mpad2<<" time "<<mtime<<" +- "<<mtime2<<" charge "<<tc<<endl; | |
327 | ||
328 | fNClusters++; | |
329 | FreeSeq(index); | |
330 | } | |
331 | } | |
332 | ||
333 | inline void AliL3VHDLClusterFinder::FreeSeq(UShort_t i) | |
334 | { | |
335 | ClearSeq(i); | |
336 | IncPointer(fLast,1,N_clmem); | |
337 | } | |
338 | ||
339 | void AliL3VHDLClusterFinder::Clear() | |
340 | { | |
341 | fFirst=0; //first list pointer | |
342 | fLast=0; //last list pointer | |
343 | ||
344 | fWP=0; //write pointer | |
345 | fRP=0; //read pointer | |
346 | fEP=0; //end read pointer | |
347 | fFP=0; //end flush pointer | |
348 | fOP=0; //output pointer | |
349 | ||
350 | fSeq.fRow=0; | |
351 | fSeq.fLastPad=0; | |
352 | fSeq.fTotalCharge=0; | |
353 | fSeq.fPad=0; | |
354 | fSeq.fTime=0; | |
355 | fSeq.fPad2=0; | |
356 | fSeq.fTime2=0; | |
357 | fSeq.fMean=0; | |
358 | fSeq.fMerge=0; | |
359 | ||
360 | for(UInt_t i=0;i<N_mem;i++) fPList[i]=N_clmem; | |
361 | for(UInt_t i=0;i<N_clmem;i++) ClearSeq(i); | |
362 | } | |
363 | ||
364 | void AliL3VHDLClusterFinder::ClearSeq(UShort_t i){ | |
365 | fSeqs[i].fRow=0; | |
366 | fSeqs[i].fLastPad=0; | |
367 | fSeqs[i].fTotalCharge=0; | |
368 | fSeqs[i].fPad=0; | |
369 | fSeqs[i].fTime=0; | |
370 | fSeqs[i].fPad2=0; | |
371 | fSeqs[i].fTime2=0; | |
372 | fSeqs[i].fMean=0; | |
373 | fSeqs[i].fMerge=0; | |
374 | } | |
375 | ||
376 | void AliL3VHDLClusterFinder::IncPointer(UShort_t &p, Short_t add, UShort_t N){ | |
377 | Short_t pp=p; | |
378 | pp+=add; | |
379 | if(pp>=N) p=UShort_t(pp-N); | |
380 | else if(pp<0) p=UShort_t(pp+N); | |
381 | else p=UShort_t(pp); | |
382 | } | |
383 | ||
384 | inline void AliL3VHDLClusterFinder::IncRPointer(){ | |
385 | IncPointer(fRP); | |
386 | } | |
387 | ||
388 | inline void AliL3VHDLClusterFinder::IncWPointer(){ | |
389 | IncPointer(fWP); | |
390 | ||
391 | if(fWP==fOP){ | |
392 | LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::IncWPointer","Memory Check") | |
393 | <<"Write pointer overwrites output pointer."<<ENDLOG; | |
394 | } | |
395 | } | |
396 | ||
397 | inline void AliL3VHDLClusterFinder::NextFreeIndex(){ | |
398 | IncPointer(fFirst,1,N_clmem); | |
399 | if(fFirst==fLast) { | |
400 | LOG(AliL3Log::kFatal,"AliL3VHDLClusterFinder::GetFreeIndex","Memory Check") | |
401 | <<"No space left in sequence list: "<<fFirst<<"=="<<fLast<<ENDLOG; | |
402 | } | |
403 | } | |
404 | ||
405 | #if 0 | |
406 | void AliL3ClustFinderNew::WriteClusters(Int_t n_clusters,ClusterData *list) | |
407 | { | |
408 | Int_t thisrow,thissector; | |
409 | UInt_t counter = fNClusters; | |
410 | ||
411 | for(int j=0; j<n_clusters; j++) | |
412 | { | |
413 | if(!list[j].fFlags) continue; //discard 1 pad clusters | |
414 | if(list[j].fTotalCharge < fThreshold) continue; //noise cluster | |
415 | ||
416 | Float_t xyz[3]; | |
417 | Float_t fpad =(Float_t)list[j].fPad /(Float_t)list[j].fTotalCharge; | |
418 | Float_t fpad2=fXYErr; | |
419 | Float_t ftime =(Float_t)list[j].fTime /(Float_t)list[j].fTotalCharge; | |
420 | Float_t ftime2=fZErr; | |
421 | ||
422 | if(fcalcerr) { | |
423 | fpad2=(Float_t)list[j].fPad2/(Float_t)list[j].fTotalCharge - fpad*fpad; | |
424 | fpad2 = sqrt(fpad2); | |
425 | ftime2=(Float_t)list[j].fTime2/(Float_t)list[j].fTotalCharge - ftime*ftime; | |
426 | ftime2 = sqrt(ftime2); | |
427 | } | |
428 | ||
429 | if(fstdout==kTRUE) | |
430 | cout<<"WriteCluster: padrow "<<fCurrentRow<<" pad "<<fpad << "+-"<<fpad2<<" time "<<ftime<<"+-"<<ftime2<<" charge "<<list[j].fTotalCharge<<endl; | |
431 | ||
432 | AliL3Transform::Slice2Sector(fCurrentSlice,fCurrentRow,thissector,thisrow); | |
433 | AliL3Transform::Raw2Local(xyz,thissector,thisrow,fpad,ftime); | |
434 | if(xyz[0]==0) LOG(AliL3Log::kError,"AliL3ClustFinder","Cluster Finder") | |
435 | <<AliL3Log::kDec<<"Zero cluster"<<ENDLOG; | |
436 | if(fNClusters >= fMaxNClusters) | |
437 | { | |
438 | LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Finder") | |
439 | <<AliL3Log::kDec<<"Too many clusters"<<ENDLOG; | |
440 | return; | |
441 | } | |
442 | fSpacePointData[counter].fCharge = list[j].fTotalCharge; | |
443 | fSpacePointData[counter].fX = xyz[0]; | |
444 | fSpacePointData[counter].fY = xyz[1]; | |
445 | fSpacePointData[counter].fZ = xyz[2]; | |
446 | fSpacePointData[counter].fPadRow = fCurrentRow; | |
447 | fSpacePointData[counter].fXYErr = fpad2; | |
448 | fSpacePointData[counter].fZErr = ftime2; | |
449 | fSpacePointData[counter].fID = counter | |
450 | +((fCurrentSlice&0x7f)<<25)+((fCurrentPatch&0x7)<<22);//Uli | |
451 | ||
452 | fNClusters++; | |
453 | counter++; | |
454 | } | |
455 | } | |
456 | #endif | |
457 |