d20df6fc |
1 | /**************************************************************************\r |
3a039a31 |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r |
3 | * *\r |
4 | * Author: The ALICE Off-line Project. *\r |
5 | * Contributors are mentioned in the code where appropriate. *\r |
6 | * *\r |
7 | * Permission to use, copy, modify and distribute this software and its *\r |
8 | * documentation strictly for non-commercial purposes is hereby granted *\r |
9 | * without fee, provided that the above copyright notice appears in all *\r |
10 | * copies and that both the copyright notice and this permission notice *\r |
11 | * appear in the supporting documentation. The authors make no claims *\r |
12 | * about the suitability of this software for any purpose. It is *\r |
13 | * provided "as is" without express or implied warranty. *\r |
14 | **************************************************************************/\r |
d20df6fc |
15 | \r |
16 | /* $Id: AliTRDchamberTimeBin.cxx 23313 2008-01-11 14:56:43Z cblume $ */\r |
17 | \r |
18 | ///////////////////////////////////////////////////////////////////////////////\r |
19 | // //\r |
20 | // Organization of clusters at the level of 1 TRD chamber. //\r |
21 | // The data structure is used for tracking at the stack level. //\r |
22 | // //\r |
23 | // Functionalities: //\r |
24 | // 1. cluster organization and sorting //\r |
25 | // 2. fast data navigation //\r |
26 | // //\r |
27 | // Authors: //\r |
28 | // Alex Bercuci <A.Bercuci@gsi.de> //\r |
29 | // Markus Fasel <M.Fasel@gsi.de> //\r |
30 | // //\r |
31 | ///////////////////////////////////////////////////////////////////////////////\r |
32 | \r |
33 | #include <TObject.h>\r |
34 | #include <TROOT.h>\r |
35 | #include <TMath.h>\r |
36 | #include <TStopwatch.h>\r |
37 | #include <TTreeStream.h>\r |
38 | \r |
39 | #include "AliLog.h"\r |
40 | \r |
41 | #include "AliTRDcluster.h"\r |
42 | #include "AliTRDchamberTimeBin.h"\r |
43 | #include "AliTRDrecoParam.h"\r |
44 | #include "AliTRDReconstructor.h"\r |
45 | #include "AliTRDtrackerV1.h"\r |
46 | \r |
47 | \r |
48 | ClassImp(AliTRDchamberTimeBin)\r |
49 | \r |
50 | //_____________________________________________________________________________\r |
51 | AliTRDchamberTimeBin::AliTRDchamberTimeBin(Int_t plane, Int_t stack, Int_t sector, Double_t z0, Double_t zLength)\r |
52 | :TObject()\r |
3a039a31 |
53 | ,fReconstructor(0x0)\r |
d20df6fc |
54 | ,fOwner(kFALSE)\r |
55 | ,fPlane(plane)\r |
56 | ,fStack(stack)\r |
57 | ,fSector(sector)\r |
58 | ,fNRows(kMaxRows)\r |
59 | ,fN(0)\r |
60 | ,fX(0.)\r |
61 | ,fZ0(z0)\r |
62 | ,fZLength(zLength)\r |
63 | {\r |
64 | //\r |
65 | // Default constructor (Only provided to use AliTRDchamberTimeBin with arrays)\r |
66 | //\r |
67 | \r |
3a039a31 |
68 | for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r |
69 | for(int ic=0; ic<kMaxClustersLayer; ic++){\r |
70 | fClusters[ic] = 0x0;\r |
71 | fIndex[ic] = 0xffff;\r |
72 | }\r |
d20df6fc |
73 | }\r |
74 | \r |
75 | // //_____________________________________________________________________________\r |
76 | // AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer, Double_t\r |
77 | // z0, Double_t zLength, UChar_t stackNr):\r |
78 | // TObject()\r |
79 | // ,fOwner(kFALSE)\r |
80 | // ,fPlane(0xff)\r |
81 | // ,fStack(0xff)\r |
82 | // ,fSector(0xff)\r |
83 | // ,fNRows(kMaxRows)\r |
84 | // ,fN(0)\r |
85 | // ,fX(0.)\r |
86 | // ,fZ0(z0)\r |
87 | // ,fZLength(zLength)\r |
88 | // {\r |
89 | // // Standard constructor.\r |
90 | // // Initialize also the underlying AliTRDpropagationLayer using the copy constructor.\r |
91 | // \r |
92 | // SetT0(layer.IsT0());\r |
93 | // for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r |
94 | // for(int ic=0; ic<kMaxClustersLayer; ic++){\r |
95 | // fClusters[ic] = 0x0;\r |
96 | // fIndex[ic] = 0xffff;\r |
97 | // }\r |
98 | // }\r |
99 | // \r |
100 | // //_____________________________________________________________________________\r |
101 | // AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDpropagationLayer &layer):\r |
102 | // TObject()\r |
103 | // ,fOwner(kFALSE)\r |
104 | // ,fPlane(0xff)\r |
105 | // ,fStack(0xff)\r |
106 | // ,fSector(0xff)\r |
107 | // ,fNRows(kMaxRows)\r |
108 | // ,fN(0)\r |
109 | // ,fX(0.)\r |
110 | // ,fZ0(0)\r |
111 | // ,fZLength(0)\r |
112 | // {\r |
113 | // // Standard constructor using only AliTRDpropagationLayer.\r |
114 | // \r |
115 | // SetT0(layer.IsT0());\r |
116 | // for(int i=0; i<kMaxRows; i++) fPositions[i] = 0xff;\r |
117 | // for(int ic=0; ic<kMaxClustersLayer; ic++){\r |
118 | // fClusters[ic] = 0x0;\r |
119 | // fIndex[ic] = 0xffff;\r |
120 | // }\r |
121 | // }\r |
122 | // //_____________________________________________________________________________\r |
123 | // AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDpropagationLayer &layer)\r |
124 | // {\r |
125 | // // Assignment operator from an AliTRDpropagationLayer\r |
126 | // \r |
127 | // if (this != &layer) layer.Copy(*this);\r |
128 | // return *this;\r |
129 | // }\r |
130 | // \r |
131 | \r |
132 | //_____________________________________________________________________________\r |
133 | AliTRDchamberTimeBin::AliTRDchamberTimeBin(const AliTRDchamberTimeBin &layer):\r |
3a039a31 |
134 | TObject()\r |
135 | ,fReconstructor(layer.fReconstructor)\r |
136 | ,fOwner(layer.fOwner)\r |
d20df6fc |
137 | ,fPlane(layer.fPlane)\r |
138 | ,fStack(layer.fStack)\r |
139 | ,fSector(layer.fSector)\r |
3a039a31 |
140 | ,fNRows(layer.fNRows)\r |
d20df6fc |
141 | ,fN(layer.fN)\r |
142 | ,fX(layer.fX)\r |
3a039a31 |
143 | ,fZ0(layer.fZ0)\r |
144 | ,fZLength(layer.fZLength)\r |
d20df6fc |
145 | {\r |
146 | // Copy Constructor (performs a deep copy)\r |
3a039a31 |
147 | \r |
148 | SetT0(layer.IsT0());\r |
149 | for(int i=0; i<kMaxRows; i++) fPositions[i] = layer.fPositions[i];\r |
150 | memcpy(&fClusters[0], &layer.fClusters[0], kMaxClustersLayer*sizeof(UChar_t));\r |
151 | memcpy(&fIndex[0], &layer.fIndex[0], kMaxClustersLayer*sizeof(UInt_t));\r |
d20df6fc |
152 | \r |
153 | \r |
154 | // BuildIndices();\r |
155 | }\r |
156 | \r |
157 | //_____________________________________________________________________________\r |
158 | AliTRDchamberTimeBin &AliTRDchamberTimeBin::operator=(const AliTRDchamberTimeBin &layer)\r |
159 | {\r |
160 | // Assignment operator\r |
161 | \r |
3a039a31 |
162 | if (this != &layer) layer.Copy(*this);\r |
d20df6fc |
163 | return *this;\r |
164 | }\r |
165 | \r |
d611c74f |
166 | //_____________________________________________________________________________\r |
167 | void AliTRDchamberTimeBin::Clear(const Option_t *) \r |
168 | { \r |
3a039a31 |
169 | for (Int_t i = 0; i < fN; i++){ \r |
d611c74f |
170 | if(!fClusters[i]) continue;\r |
171 | if(fOwner) delete fClusters[i];\r |
172 | fClusters[i] = NULL;\r |
173 | }\r |
3a039a31 |
174 | fN = 0; \r |
d611c74f |
175 | }\r |
176 | \r |
d20df6fc |
177 | //_____________________________________________________________________________\r |
178 | void AliTRDchamberTimeBin::Copy(TObject &o) const\r |
179 | {\r |
180 | // Copy method. Performs a deep copy of all data from this object to object o.\r |
3a039a31 |
181 | \r |
182 | AliTRDchamberTimeBin &layer = (AliTRDchamberTimeBin &)o;\r |
183 | layer.fReconstructor = fReconstructor;\r |
184 | layer.fOwner = kFALSE;\r |
185 | layer.fPlane = fPlane;\r |
186 | layer.fStack = fStack;\r |
187 | layer.fSector = fSector;\r |
188 | layer.fNRows = fNRows;\r |
189 | layer.fN = fN;\r |
190 | layer.fX = fX;\r |
191 | layer.fZ0 = fZ0;\r |
192 | layer.fZLength = fZLength;\r |
193 | layer.SetT0(IsT0());\r |
194 | \r |
195 | for(int i = 0; i < kMaxRows; i++) layer.fPositions[i] = 0;\r |
196 | \r |
197 | for(int i=0; i<kMaxRows; i++) layer.fPositions[i] = fPositions[i];\r |
198 | memcpy(&layer.fClusters[0], &fClusters[0], kMaxClustersLayer*sizeof(UChar_t));\r |
199 | memcpy(&layer.fIndex[0], &fIndex[0], kMaxClustersLayer*sizeof(UInt_t));\r |
200 | \r |
201 | TObject::Copy(layer); // copies everything into layer\r |
202 | \r |
d20df6fc |
203 | // layer.BuildIndices();\r |
204 | }\r |
205 | \r |
206 | //_____________________________________________________________________________\r |
207 | AliTRDchamberTimeBin::~AliTRDchamberTimeBin()\r |
208 | {\r |
209 | // Destructor\r |
3a039a31 |
210 | if(fOwner) for(int ic=0; ic<fN; ic++) delete fClusters[ic];\r |
d20df6fc |
211 | }\r |
212 | \r |
213 | //_____________________________________________________________________________\r |
214 | void AliTRDchamberTimeBin::SetRange(const Float_t z0, const Float_t zLength)\r |
215 | {\r |
216 | // Sets the range in z-direction\r |
217 | //\r |
218 | // Parameters:\r |
219 | // z0 : starting position of layer in the z direction\r |
220 | // zLength : length of layer in the z direction \r |
221 | \r |
3a039a31 |
222 | fZ0 = (z0 <= z0 + zLength) ? z0 : z0 + zLength;\r |
223 | fZLength = TMath::Abs(zLength);\r |
d20df6fc |
224 | }\r |
225 | \r |
226 | //_____________________________________________________________________________\r |
227 | void AliTRDchamberTimeBin::InsertCluster(AliTRDcluster *c, UInt_t index) \r |
228 | {\r |
229 | //\r |
230 | // Insert cluster in cluster array.\r |
231 | // Clusters are sorted according to Y coordinate. \r |
232 | //\r |
233 | \r |
234 | //if (fTimeBinIndex < 0) { \r |
235 | //AliWarning("Attempt to insert cluster into non-sensitive time bin!\n");\r |
236 | //return;\r |
237 | //}\r |
238 | \r |
239 | if (fN == (Int_t) kMaxClustersLayer) {\r |
240 | //AliWarning("Too many clusters !\n"); \r |
241 | return;\r |
242 | }\r |
243 | \r |
244 | if (fN == 0) {\r |
245 | fIndex[0] = index; \r |
246 | fClusters[fN++] = c; \r |
247 | return;\r |
248 | }\r |
249 | \r |
250 | Int_t i = Find(c->GetY());\r |
251 | memmove(fClusters+i+1,fClusters+i,(fN-i)*sizeof(AliTRDcluster*));\r |
252 | memmove(fIndex +i+1,fIndex +i,(fN-i)*sizeof(UInt_t)); \r |
253 | fIndex[i] = index; \r |
254 | fClusters[i] = c; \r |
255 | fN++;\r |
256 | \r |
257 | } \r |
258 | \r |
259 | \r |
260 | //_____________________________________________________________________________\r |
261 | void AliTRDchamberTimeBin::BuildIndices(Int_t iter)\r |
262 | {\r |
263 | // Rearrangement of the clusters belonging to the propagation layer for the stack.\r |
264 | //\r |
265 | // Detailed description\r |
266 | //\r |
267 | // The array indices of all clusters in one PropagationLayer are stored in\r |
268 | // array. The array is divided into several bins.\r |
269 | // The clusters are sorted in increasing order of their y coordinate.\r |
270 | //\r |
271 | // Sorting algorithm: TreeSearch\r |
272 | //\r |
273 | \r |
3a039a31 |
274 | if(!fN) return;\r |
275 | \r |
276 | // Select clusters that belong to the Stack\r |
277 | Int_t nClStack = 0; // Internal counter\r |
278 | for(Int_t i = 0; i < fN; i++){\r |
279 | if(fClusters[i]->IsUsed()){\r |
280 | fClusters[i] = 0x0;\r |
281 | fIndex[i] = 0xffff;\r |
282 | } else nClStack++;\r |
283 | }\r |
284 | if(nClStack > kMaxClustersLayer) AliWarning(Form("Number of clusters in stack %d exceed buffer size %d. Truncating.", nClStack, kMaxClustersLayer));\r |
285 | \r |
286 | // Nothing in this time bin. Reset indexes \r |
287 | if(!nClStack){\r |
288 | fN = 0;\r |
289 | memset(&fPositions[0], 0xff, sizeof(UChar_t) * kMaxRows);\r |
290 | memset(&fClusters[0], 0x0, sizeof(AliTRDcluster*) * kMaxClustersLayer);\r |
291 | memset(&fIndex[0], 0xffff, sizeof(UInt_t) * kMaxClustersLayer);\r |
292 | return;\r |
293 | }\r |
294 | \r |
295 | // Make a copy\r |
296 | AliTRDcluster *helpCL[kMaxClustersLayer];\r |
297 | Int_t helpInd[kMaxClustersLayer];\r |
298 | nClStack = 0;\r |
299 | for(Int_t i = 0; i < TMath::Min(fN, kMaxClustersLayer); i++){\r |
300 | if(!fClusters[i]) continue;\r |
301 | helpCL[nClStack] = fClusters[i];\r |
302 | helpInd[nClStack] = fIndex[i];\r |
303 | fClusters[i] = 0x0;\r |
304 | fIndex[i] = 0xffff;\r |
305 | nClStack++;\r |
306 | }\r |
307 | \r |
308 | // do clusters arrangement\r |
309 | fX = 0.;\r |
310 | fN = nClStack;\r |
311 | nClStack = 0;\r |
312 | // Reset Positions array\r |
313 | memset(fPositions, 0, sizeof(UChar_t)*kMaxRows);\r |
314 | for(Int_t i = 0; i < fN; i++){\r |
315 | // boundary check\r |
316 | AliTRDcluster *cl = helpCL[i];\r |
317 | UChar_t rowIndex = cl->GetPadRow();\r |
318 | // Insert Leaf\r |
319 | Int_t pos = FindYPosition(cl->GetY(), rowIndex, i);\r |
320 | if(pos == -1){ // zbin is empty;\r |
321 | Int_t upper = (rowIndex == fNRows - 1) ? nClStack : fPositions[rowIndex + 1];\r |
322 | memmove(fClusters + upper + 1, fClusters + upper, (sizeof(AliTRDcluster *))*(nClStack-upper));\r |
323 | memmove(fIndex + upper + 1, fIndex + upper, (sizeof(UInt_t))*(nClStack-upper));\r |
324 | fClusters[upper] = cl;\r |
325 | fIndex[upper] = helpInd[i]; \r |
326 | // Move All pointer one position back\r |
327 | for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++;\r |
328 | nClStack++;\r |
329 | } else { // zbin not empty\r |
330 | memmove(fClusters + pos + 2, fClusters + pos+1, (sizeof(AliTRDcluster *))*(nClStack-(pos+1)));\r |
331 | memmove(fIndex + pos + 2, fIndex + pos+1, (sizeof(UInt_t))*(nClStack-(pos+1)));\r |
332 | fClusters[pos + 1] = cl; //fIndex[i];\r |
333 | fIndex[pos + 1] = helpInd[i];\r |
334 | // Move All pointer one position back\r |
335 | for(UChar_t j = rowIndex + 1; j < fNRows; j++) fPositions[j]++; \r |
336 | nClStack++;\r |
337 | }\r |
338 | \r |
339 | // calculate mean x\r |
340 | fX += cl->GetX(); \r |
341 | \r |
342 | // Debug Streaming\r |
343 | if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker) >= 3){\r |
344 | TTreeSRedirector &cstream = *AliTRDtrackerV1::DebugStreamer();\r |
345 | cstream << "BuildIndices"\r |
346 | << "Plane=" << fPlane\r |
347 | << "Stack=" << fStack\r |
348 | << "Sector=" << fSector\r |
349 | << "Iter=" << iter\r |
350 | << "C.=" << cl\r |
351 | << "rowIndex=" << rowIndex\r |
352 | << "\n";\r |
353 | }\r |
354 | }\r |
d20df6fc |
355 | \r |
356 | // AliInfo("Positions");\r |
357 | // for(int ir=0; ir<fNRows; ir++) printf("pos[%d] %d\n", ir, fPositions[ir]);\r |
358 | \r |
3a039a31 |
359 | fX /= fN;\r |
d20df6fc |
360 | }\r |
361 | \r |
362 | //_____________________________________________________________________________\r |
363 | Int_t AliTRDchamberTimeBin::Find(Float_t y) const\r |
364 | {\r |
365 | //\r |
366 | // Returns index of the cluster nearest in Y \r |
367 | //\r |
368 | \r |
369 | if (fN <= 0) return 0;\r |
370 | \r |
371 | if (y <= fClusters[0]->GetY()) return 0;\r |
372 | \r |
373 | if (y > fClusters[fN-1]->GetY()) return fN;\r |
374 | \r |
375 | \r |
376 | Int_t b = 0;\r |
377 | Int_t e = fN - 1;\r |
378 | Int_t m = (b + e) / 2;\r |
379 | \r |
380 | for ( ; b < e; m = (b + e) / 2) {\r |
381 | if (y > fClusters[m]->GetY()) b = m + 1;\r |
382 | else e = m;\r |
383 | }\r |
384 | \r |
385 | return m;\r |
386 | } \r |
387 | \r |
388 | //_____________________________________________________________________________\r |
389 | Int_t AliTRDchamberTimeBin::FindYPosition(Double_t y, UChar_t z, Int_t nClusters) const\r |
390 | {\r |
391 | //\r |
392 | // Tree search Algorithm to find the nearest left cluster for a given\r |
393 | // y-position in a certain z-bin (in fact AVL-tree). \r |
394 | // Making use of the fact that clusters are sorted in y-direction.\r |
395 | //\r |
396 | // Parameters:\r |
397 | // y : y position of the reference point in tracking coordinates\r |
398 | // z : z reference bin.\r |
399 | // nClusters : \r |
400 | //\r |
401 | // Output :\r |
402 | // Index of the nearest left cluster in the StackLayer indexing (-1 if no clusters are found)\r |
403 | //\r |
404 | \r |
3a039a31 |
405 | Int_t start = fPositions[z]; // starting Position of the bin\r |
406 | Int_t upper = (Int_t)((z != fNRows - 1) ? fPositions[z+1] : nClusters); // ending Position of the bin \r |
407 | Int_t end = upper - 1; // ending Position of the bin \r |
408 | if(end < start) return -1; // Bin is empty\r |
409 | Int_t middle = static_cast<Int_t>((start + end)/2);\r |
410 | // 1st Part: climb down the tree: get the next cluster BEFORE ypos\r |
411 | while(start + 1 < end){\r |
412 | if(y >= fClusters[middle]->GetY()) start = middle;\r |
413 | else end = middle;\r |
414 | middle = static_cast<Int_t>((start + end)/2);\r |
415 | }\r |
416 | if(y > fClusters[end]->GetY()) return end;\r |
417 | return start;\r |
d20df6fc |
418 | }\r |
419 | \r |
420 | //_____________________________________________________________________________\r |
421 | Int_t AliTRDchamberTimeBin::FindNearestYCluster(Double_t y, UChar_t z) const\r |
422 | {\r |
423 | //\r |
424 | // Tree search Algorithm to find the nearest cluster for a given\r |
425 | // y-position in a certain z-bin (in fact AVL-tree). \r |
426 | // Making use of the fact that clusters are sorted in y-direction.\r |
427 | //\r |
428 | // Parameters:\r |
429 | // y : y position of the reference point in tracking coordinates\r |
430 | // z : z reference bin.\r |
431 | //\r |
432 | // Output \r |
433 | // Index of the nearest cluster in the StackLayer indexing (-1 if no clusters are found)\r |
434 | //\r |
435 | \r |
3a039a31 |
436 | Int_t position = FindYPosition(y, z, fN);\r |
437 | if(position == -1) return position; // bin empty\r |
438 | // FindYPosition always returns the left Neighbor. We don't know if the left or the right Neighbor is nearest\r |
439 | // to the Reference y-position, so test both\r |
440 | Int_t upper = (Int_t)((z < fNRows-1) ? fPositions[z+1] : fN); // ending Position of the bin\r |
441 | if((position + 1) < (upper)){\r |
442 | if(TMath::Abs(y - fClusters[position + 1]->GetY()) < TMath::Abs(y - fClusters[position]->GetY())) return position + 1;\r |
443 | else return position;\r |
444 | }\r |
445 | return position;\r |
d20df6fc |
446 | }\r |
447 | \r |
448 | //_____________________________________________________________________________\r |
449 | Int_t AliTRDchamberTimeBin::SearchNearestCluster(Double_t y, Double_t z, Double_t maxroady, Double_t maxroadz) const\r |
450 | {\r |
451 | //\r |
452 | // Finds the nearest cluster from a given point in a defined range.\r |
453 | // Distance is determined in a 2D space by the 2-Norm.\r |
454 | //\r |
455 | // Parameters :\r |
456 | // y : y position of the reference point in tracking coordinates\r |
457 | // z : z reference bin.\r |
458 | // maxroady : maximum searching distance in y direction\r |
459 | // maxroadz : maximum searching distance in z direction\r |
460 | //\r |
461 | // Output \r |
462 | // Index of the nearest cluster in the StackLayer indexing (-1 if no cluster is found).\r |
463 | // Cluster can be accessed with the operator[] or GetCluster(Int_t index)\r |
464 | //\r |
465 | // Detail description\r |
466 | //\r |
467 | // The following steps are perfomed:\r |
468 | // 1. Get the expected z bins inside maxroadz.\r |
469 | // 2. For each z bin find nearest y cluster.\r |
470 | // 3. Select best candidate\r |
471 | //\r |
3a039a31 |
472 | Int_t index = -1;\r |
473 | // initial minimal distance will be represented as ellipse: semi-major = z-direction\r |
474 | // later 2-Norm will be used \r |
d20df6fc |
475 | // Float_t nExcentricity = TMath::Sqrt(maxroadz*maxroadz - maxroad*maxroad)/maxroadz;\r |
3a039a31 |
476 | Float_t mindist = maxroadz;\r |
477 | \r |
478 | // not very nice but unfortunately neccessarry: we have ho check the neighbors in both directions (+ and -) too. How \r |
479 | // much neighbors depends on the Quotient maxroadz/fZLength \r |
480 | UChar_t maxRows = 3;\r |
481 | UChar_t zpos[kMaxRows];\r |
d20df6fc |
482 | // Float_t mindist = TMath::Sqrt(maxroad*maxroad + maxroadz*maxroadz);\r |
483 | // UChar_t myZbin = FindTreePosition(z, fZ0 + fZLength/2, fZLength/4, 8, 8, kFALSE);\r |
3a039a31 |
484 | UChar_t myZbin = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - z)/fZLength * fNRows);\r |
485 | if(z < fZ0) myZbin = fNRows - 1;\r |
486 | if(z > fZ0 + fZLength) myZbin = 0;\r |
487 | //printf("\n%f < %f < %f [%d]\n", fZ0, z, fZ0 + fZLength, myZbin);\r |
488 | //for(int ic=0; ic<fN; ic++) printf("%d z = %f row %d\n", ic, fClusters[ic]->GetZ(), fClusters[ic]->GetPadRow());\r |
489 | \r |
490 | UChar_t nNeighbors = 0;\r |
491 | for(UChar_t i = 0; i < maxRows; i++){\r |
492 | if((myZbin - 1 + i) < 0) continue;\r |
493 | if((myZbin - 1 + i) > fNRows - 1) break;\r |
494 | zpos[nNeighbors] = myZbin - 1 + i;\r |
495 | nNeighbors++;\r |
496 | }\r |
497 | Float_t ycl = 0, zcl = 0;\r |
498 | for(UChar_t neighbor = 0; neighbor < nNeighbors; neighbor++){ // Always test the neighbors too\r |
499 | Int_t pos = FindNearestYCluster(y, zpos[neighbor]);\r |
500 | if(pos == -1) continue; // No cluster in bin\r |
501 | AliTRDcluster *c = (AliTRDcluster *) (fClusters[pos]);\r |
502 | if(c->IsUsed()) continue; // we are only interested in unused clusters\r |
503 | ycl = c->GetY();\r |
504 | // Too far away in y-direction (Prearrangement)\r |
505 | if (TMath::Abs(ycl - y) > maxroady){ \r |
506 | //printf("y[%f] ycl[%f] roady[%f]\n", y, ycl, maxroady);\r |
507 | continue;\r |
508 | }\r |
509 | zcl = c->GetZ();\r |
510 | // Too far away in z-Direction\r |
511 | // (Prearrangement since we have not so many bins to test)\r |
512 | if (TMath::Abs(zcl - z) > maxroadz) continue;\r |
513 | \r |
514 | Float_t dist; // distance defined as 2-Norm \r |
515 | // if we havent found a Particle that is in the ellipse around (y,z) with maxroad as semi-minor and\r |
516 | // maxroadz as semi-major, we take the radius of the ellipse concerning the cluster as mindist, later we \r |
517 | // take the 2-Norm when we found a cluster inside the ellipse (The value 10000 is taken because it is surely\r |
518 | // large enough to be usable as an indicator whether we have found a nearer cluster or not)\r |
d20df6fc |
519 | // if(mindist > 10000.){\r |
520 | // Float_t phi = ((zcl - z) == 0) ? TMath::Pi()/2 : TMath::ATan((ycl - y)/(zcl - z));\r |
521 | // mindist = maxroad/TMath::Sqrt(1 - nExcentricity*nExcentricity * (TMath::Cos(phi))*(TMath::Cos(phi)));\r |
522 | // }\r |
3a039a31 |
523 | dist = TMath::Max(TMath::Abs(y-ycl),TMath::Abs(z-zcl)); // infinity Norm\r |
d20df6fc |
524 | // dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z));\r |
3a039a31 |
525 | if((Int_t)(dist * 100000) < (Int_t)(mindist * 100000)){\r |
526 | //if((dist = TMath::Sqrt((ycl - y)*(ycl - y) + (zcl - z)*(zcl - z))) < mindist){\r |
527 | mindist = dist;\r |
528 | index = pos;\r |
529 | } \r |
530 | }\r |
531 | // This is the Array Position in fIndex2D of the Nearest cluster: if a\r |
532 | // cluster is called, then the function has to retrieve the Information\r |
533 | // which is Stored in the Array called, the function\r |
534 | return index;\r |
d20df6fc |
535 | }\r |
536 | \r |
537 | //_____________________________________________________________________________\r |
538 | void AliTRDchamberTimeBin::BuildCond(AliTRDcluster *cl, Double_t *cond, UChar_t Layer, Double_t theta, Double_t phi)\r |
539 | {\r |
540 | // Helper function to calculate the area where to expect a cluster in THIS\r |
541 | // layer. \r |
542 | //\r |
543 | // Parameters :\r |
544 | // cl : \r |
545 | // cond :\r |
546 | // Layer : \r |
547 | // theta : \r |
548 | // phi : \r |
549 | //\r |
550 | // Detail description\r |
551 | //\r |
552 | // Helper function to calculate the area where to expect a cluster in THIS\r |
553 | // layer. by using the information of a former cluster in another layer\r |
554 | // and the angle in theta- and phi-direction between layer 0 and layer 3.\r |
555 | // If the layer is zero, initial conditions are calculated. Otherwise a\r |
556 | // linear interpolation is performed. \r |
557 | //Begin_Html\r |
558 | //<img src="gif/build_cond.gif">\r |
559 | //End_Html\r |
560 | //\r |
561 | \r |
3a039a31 |
562 | if(!fReconstructor){\r |
563 | AliError("Reconstructor not set.");\r |
564 | return;\r |
565 | }\r |
566 | \r |
567 | if(Layer == 0){\r |
568 | cond[0] = cl->GetY(); // center: y-Direction\r |
569 | cond[1] = cl->GetZ(); // center: z-Direction\r |
570 | cond[2] = fReconstructor->GetRecoParam()->GetMaxPhi() * (cl->GetX() - GetX()) + 1.0; // deviation: y-Direction\r |
571 | cond[3] = fReconstructor->GetRecoParam()->GetMaxTheta() * (cl->GetX() - GetX()) + 1.0; // deviation: z-Direction\r |
572 | } else {\r |
573 | cond[0] = cl->GetY() + phi * (GetX() - cl->GetX()); \r |
574 | cond[1] = cl->GetZ() + theta * (GetX() - cl->GetX());\r |
575 | cond[2] = fReconstructor->GetRecoParam()->GetRoad0y() + phi;\r |
576 | cond[3] = fReconstructor->GetRecoParam()->GetRoad0z();\r |
577 | }\r |
d20df6fc |
578 | }\r |
579 | \r |
580 | //_____________________________________________________________________________\r |
581 | void AliTRDchamberTimeBin::GetClusters(Double_t *cond, Int_t *index, Int_t& ncl, Int_t BufferSize)\r |
582 | {\r |
583 | // Finds all clusters situated in this layer inside a rectangle given by the center an ranges.\r |
584 | //\r |
585 | // Parameters :\r |
586 | // cond :\r |
587 | // index : \r |
588 | // ncl :\r |
589 | // BufferSize :\r |
590 | //\r |
591 | // Output :\r |
592 | //\r |
593 | // Detail description\r |
594 | //\r |
595 | // Function returs an array containing the indices in the stacklayer of\r |
596 | // the clusters found an the number of found clusters in the stacklayer\r |
597 | \r |
3a039a31 |
598 | ncl = 0;\r |
599 | memset(index, 0, BufferSize*sizeof(Int_t));\r |
600 | if(fN == 0) return;\r |
601 | \r |
602 | //Boundary checks\r |
603 | Double_t zvals[2];\r |
604 | if(((cond[1] - cond[3]) >= (fZ0 + fZLength)) || (cond[1] + cond[3]) <= fZ0) return; // We are outside of the chamvber\r |
605 | zvals[0] = ((cond[1] - cond[3]) < fZ0) ? fZ0 : (cond[1] - cond[3]);\r |
606 | zvals[1] = ((cond[1] + cond[3]) < fZ0 + fZLength) ? (cond[1] + cond[3]) : fZ0 + fZLength - 1.E-3;\r |
d20df6fc |
607 | \r |
3a039a31 |
608 | UChar_t zhi = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[0])/fZLength * fNRows);\r |
609 | UChar_t zlo = fNRows - 1 - (UChar_t)(TMath::Abs(fZ0 - zvals[1])/fZLength * fNRows);\r |
8c3cbdc1 |
610 | \r |
611 | /* AliInfo(Form("yc[%f] zc[%f] dy[%f] dz[%f]", cond[0], cond[1], cond[2], cond[3]));\r |
3a039a31 |
612 | PrintClusters();\r |
613 | AliInfo(Form("zlo[%f] zhi[%f]", zvals[0], zvals[1]));\r |
614 | AliInfo(Form("zlo[%d] zhi[%d]", zlo, zhi));*/\r |
615 | \r |
616 | //Preordering in Direction z saves a lot of loops (boundary checked)\r |
617 | for(UChar_t z = zlo; z <= zhi; z++){\r |
618 | UInt_t upper = (z < fNRows-1) ? fPositions[z+1] : fN;\r |
619 | //AliInfo(Form("z[%d] y [%d %d]", z, fPositions[z], upper));\r |
620 | for(Int_t y = fPositions[z]; y < (Int_t)upper; y++){\r |
621 | if(ncl == BufferSize){\r |
622 | AliWarning("Buffer size riched. Some clusters may be lost.");\r |
623 | return; //Buffer filled\r |
624 | }\r |
625 | \r |
626 | if(fClusters[y]->GetY() > (cond[0] + cond[2])) break; // Abbortion conditions!!!\r |
627 | if(fClusters[y]->GetY() < (cond[0] - cond[2])) continue; // Too small\r |
628 | if(((Int_t)((fClusters[y]->GetZ())*1000) < (Int_t)(zvals[0]*1000)) || ((Int_t)((fClusters[y]->GetZ())*1000) > (Int_t)(zvals[1]*1000))){/*printf("exit z\n"); TODO*/ continue;}\r |
629 | index[ncl] = y;\r |
630 | ncl++;\r |
631 | }\r |
632 | }\r |
633 | if(ncl>fN) AliError(Form("Clusters found %d > %d (clusters in layer)", ncl, fN));\r |
d20df6fc |
634 | }\r |
635 | \r |
636 | //_____________________________________________________________________________\r |
637 | AliTRDcluster *AliTRDchamberTimeBin::GetNearestCluster(Double_t *cond)\r |
638 | {\r |
639 | // Function returning a pointer to the nearest cluster (nullpointer if not successfull).\r |
640 | //\r |
641 | // Parameters :\r |
642 | // cond :\r |
643 | //\r |
644 | // Output :\r |
645 | // pointer to the nearest cluster (nullpointer if not successfull).\r |
646 | // \r |
647 | // Detail description\r |
648 | //\r |
649 | // returns a pointer to the nearest cluster (nullpointer if not\r |
650 | // successfull) by the help of the method FindNearestCluster\r |
3a039a31 |
651 | \r |
652 | \r |
653 | Double_t maxroad = fReconstructor->GetRecoParam()->GetRoad2y();\r |
654 | Double_t maxroadz = fReconstructor->GetRecoParam()->GetRoad2z();\r |
655 | \r |
656 | Int_t index = SearchNearestCluster(cond[0],cond[1],maxroad,maxroadz);\r |
657 | AliTRDcluster *returnCluster = 0x0;\r |
658 | if(index != -1) returnCluster = (AliTRDcluster *) fClusters[index];\r |
659 | return returnCluster;\r |
d20df6fc |
660 | }\r |
661 | \r |
662 | //_____________________________________________________________________________\r |
663 | void AliTRDchamberTimeBin::PrintClusters() const\r |
664 | {\r |
665 | // Prints the position of each cluster in the stacklayer on the stdout\r |
666 | //\r |
3a039a31 |
667 | printf("\nnRows = %d\n", fNRows);\r |
668 | printf("Z0 = %f\n", fZ0);\r |
669 | printf("Z1 = %f\n", fZ0+fZLength);\r |
670 | printf("clusters in AliTRDchamberTimeBin %d\n", fN);\r |
671 | for(Int_t i = 0; i < fN; i++){\r |
672 | printf("AliTRDchamberTimeBin: index=%i, Cluster: X = %3.3f [%d] Y = %3.3f [%d] Z = %3.3f [%d]\n", i, fClusters[i]->GetX(), fClusters[i]->GetLocalTimeBin(), fClusters[i]->GetY(), fClusters[i]->GetPadCol(), fClusters[i]->GetZ(), fClusters[i]->GetPadRow());\r |
673 | if(fClusters[i]->IsUsed()) printf("cluster allready used. rejected in search algorithm\n");\r |
674 | }\r |
d20df6fc |
675 | }\r |