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