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