]>
Commit | Line | Data |
---|---|---|
6de2bc40 | 1 | // ************************************************************************** |
2 | // * This file is property of and copyright by the ALICE HLT Project * | |
3 | // * All rights reserved. * | |
4 | // * * | |
5 | // * Primary Authors: * | |
6 | // * Copyright 2009 Matthias Kretz <kretz@kde.org> * | |
7 | // * * | |
8 | // * Permission to use, copy, modify and distribute this software and its * | |
9 | // * documentation strictly for non-commercial purposes is hereby granted * | |
10 | // * without fee, provided that the above copyright notice appears in all * | |
11 | // * copies and that both the copyright notice and this permission notice * | |
12 | // * appear in the supporting documentation. The authors make no claims * | |
13 | // * about the suitability of this software for any purpose. It is * | |
14 | // * provided "as is" without express or implied warranty. * | |
15 | // ************************************************************************** | |
4acc2401 | 16 | |
17 | #include "AliHLTTPCCASliceData.h" | |
18 | #include "AliHLTTPCCAClusterData.h" | |
19 | #include "AliHLTTPCCAMath.h" | |
20 | #include "AliHLTArray.h" | |
21 | #include "AliHLTTPCCAHit.h" | |
22 | #include "AliHLTTPCCAParam.h" | |
23 | #include "MemoryAssignmentHelpers.h" | |
24 | #include <iostream> | |
25 | ||
26 | // calculates an approximation for 1/sqrt(x) | |
27 | // Google for 0x5f3759df :) | |
28 | static inline float fastInvSqrt( float _x ) | |
29 | { | |
6de2bc40 | 30 | // the function calculates fast inverse sqrt |
31 | ||
4acc2401 | 32 | union { float f; int i; } x = { _x }; |
33 | const float xhalf = 0.5f * x.f; | |
34 | x.i = 0x5f3759df - ( x.i >> 1 ); | |
35 | x.f = x.f * ( 1.5f - xhalf * x.f * x.f ); | |
36 | return x.f; | |
37 | } | |
38 | ||
6de2bc40 | 39 | inline void AliHLTTPCCASliceData::CreateGrid( AliHLTTPCCARow *row, const AliHLTTPCCAClusterData &data ) |
4acc2401 | 40 | { |
6de2bc40 | 41 | // grid creation |
42 | ||
4acc2401 | 43 | if ( row->NHits() <= 0 ) { // no hits or invalid data |
44 | // grid coordinates don't matter, since there are no hits | |
45 | row->fGrid.CreateEmpty(); | |
46 | return; | |
47 | } | |
48 | ||
49 | float yMin = 1.e3f; | |
50 | float yMax = -1.e3f; | |
51 | float zMin = 1.e3f; | |
52 | float zMax = -1.e3f; | |
53 | for ( int i = row->fHitNumberOffset; i < row->fHitNumberOffset + row->fNHits; ++i ) { | |
54 | const float y = data.Y( i ); | |
55 | const float z = data.Z( i ); | |
56 | if ( yMax < y ) yMax = y; | |
57 | if ( yMin > y ) yMin = y; | |
58 | if ( zMax < z ) zMax = z; | |
59 | if ( zMin > z ) zMin = z; | |
60 | } | |
61 | ||
62 | const float norm = fastInvSqrt( row->fNHits ); | |
63 | row->fGrid.Create( yMin, yMax, zMin, zMax, | |
64 | CAMath::Max( ( yMax - yMin ) * norm, 2.f ), | |
65 | CAMath::Max( ( zMax - zMin ) * norm, 2.f ) ); | |
66 | } | |
67 | ||
68 | inline void AliHLTTPCCASliceData::PackHitData( AliHLTTPCCARow *row, const AliHLTArray<AliHLTTPCCAHit> &binSortedHits ) | |
69 | { | |
6de2bc40 | 70 | // hit data packing |
71 | ||
4acc2401 | 72 | static const float shortPackingConstant = 1.f / 65535.f; |
73 | const float y0 = row->fGrid.YMin(); | |
74 | const float z0 = row->fGrid.ZMin(); | |
75 | const float stepY = ( row->fGrid.YMax() - y0 ) * shortPackingConstant; | |
76 | const float stepZ = ( row->fGrid.ZMax() - z0 ) * shortPackingConstant; | |
77 | const float stepYi = 1.f / stepY; | |
78 | const float stepZi = 1.f / stepZ; | |
79 | ||
80 | row->fHy0 = y0; | |
81 | row->fHz0 = z0; | |
82 | row->fHstepY = stepY; | |
83 | row->fHstepZ = stepZ; | |
84 | row->fHstepYi = stepYi; | |
85 | row->fHstepZi = stepZi; | |
86 | ||
87 | for ( int hitIndex = 0; hitIndex < row->fNHits; ++hitIndex ) { | |
88 | // bin sorted index! | |
89 | const int globalHitIndex = row->fHitNumberOffset + hitIndex; | |
90 | const AliHLTTPCCAHit &hh = binSortedHits[globalHitIndex]; | |
91 | const float xx = ( ( hh.Y() - y0 ) * stepYi ) + .5 ; | |
92 | const float yy = ( ( hh.Z() - z0 ) * stepZi ) + .5 ; | |
93 | if ( xx < 0 || yy < 0 || xx >= 65536 || yy >= 65536 ) { | |
94 | std::cout << "!!!! hit packing error!!! " << xx << " " << yy << " " << std::endl; | |
95 | } | |
96 | // HitData is bin sorted | |
97 | fHitDataY[row->fHitNumberOffset + hitIndex] = xx; | |
98 | fHitDataZ[row->fHitNumberOffset + hitIndex] = yy; | |
99 | } | |
100 | } | |
101 | ||
102 | void AliHLTTPCCASliceData::Clear() | |
103 | { | |
104 | fNumberOfHits = 0; | |
105 | } | |
106 | ||
107 | void AliHLTTPCCASliceData::InitializeRows( const AliHLTTPCCAParam &p ) | |
108 | { | |
6de2bc40 | 109 | // initialisation of rows |
110 | ||
4acc2401 | 111 | for ( int i = 0; i < p.NRows(); ++i ) { |
112 | fRows[i].fX = p.RowX( i ); | |
113 | fRows[i].fMaxY = CAMath::Tan( p.DAlpha() / 2. ) * fRows[i].fX; | |
114 | } | |
115 | } | |
116 | ||
7be9b0d7 | 117 | GPUh() char* AliHLTTPCCASliceData::SetGPUSliceDataMemory(char* pGPUMemory, const AliHLTTPCCAClusterData *data) |
4acc2401 | 118 | { |
7be9b0d7 | 119 | fMemory = (char*) pGPUMemory; |
120 | return(pGPUMemory + SetPointers(data, false)); | |
121 | } | |
4acc2401 | 122 | |
7be9b0d7 | 123 | size_t AliHLTTPCCASliceData::SetPointers(const AliHLTTPCCAClusterData *data, bool allocate) |
124 | { | |
125 | const int numberOfRows = data->LastRow() - data->FirstRow(); | |
6de2bc40 | 126 | enum { kVectorAlignment = sizeof( int ) }; |
127 | const int numberOfHitsPlusAlignment = NextMultipleOf < kVectorAlignment / sizeof( int ) > ( fNumberOfHits ); | |
4acc2401 | 128 | const int memorySize = |
129 | // LinkData, HitData | |
130 | numberOfHitsPlusAlignment * 4 * sizeof( short ) + | |
131 | // FirstHitInBin | |
b8139972 | 132 | NextMultipleOf<kVectorAlignment>( ( 23 * numberOfRows + 4 * fNumberOfHits + 3 ) * sizeof( int ) ) + |
4acc2401 | 133 | // HitWeights, ClusterDataIndex |
134 | numberOfHitsPlusAlignment * 2 * sizeof( int ); | |
135 | ||
136 | if ( fMemorySize < memorySize ) { | |
7be9b0d7 | 137 | fMemorySize = memorySize; |
138 | if (allocate) | |
139 | { | |
140 | delete[] fMemory; | |
141 | fMemory = new char[fMemorySize + 4];// kVectorAlignment]; | |
142 | } | |
4acc2401 | 143 | } |
144 | ||
145 | char *mem = fMemory; | |
146 | AssignMemory( fLinkUpData, mem, numberOfHitsPlusAlignment ); | |
147 | AssignMemory( fLinkDownData, mem, numberOfHitsPlusAlignment ); | |
148 | AssignMemory( fHitDataY, mem, numberOfHitsPlusAlignment ); | |
149 | AssignMemory( fHitDataZ, mem, numberOfHitsPlusAlignment ); | |
6de2bc40 | 150 | AssignMemory( fFirstHitInBin, mem, 23 * numberOfRows + 4 * fNumberOfHits + 3 ); |
4acc2401 | 151 | AssignMemory( fHitWeights, mem, numberOfHitsPlusAlignment ); |
152 | AssignMemory( fClusterDataIndex, mem, numberOfHitsPlusAlignment ); | |
7be9b0d7 | 153 | return(mem - fMemory); |
154 | } | |
155 | ||
156 | void AliHLTTPCCASliceData::InitFromClusterData( const AliHLTTPCCAClusterData &data ) | |
157 | { | |
158 | // initialisation from cluster data | |
159 | ||
160 | //////////////////////////////////// | |
161 | // 1. prepare arrays | |
162 | //////////////////////////////////// | |
163 | ||
f0fb467d | 164 | //const int numberOfRows = data.LastRow() - data.FirstRow(); |
7be9b0d7 | 165 | fNumberOfHits = data.NumberOfClusters(); |
166 | ||
167 | /* TODO Vectorization | |
168 | for ( int rowIndex = data.FirstRow(); rowIndex <= data.LastRow(); ++rowIndex ) { | |
169 | int NumberOfClusters( int rowIndex ) const; | |
170 | } | |
171 | const int memorySize = fNumberOfHits * sizeof( short_v::Type ) | |
172 | */ | |
173 | SetPointers(&data, true); | |
4acc2401 | 174 | |
175 | //////////////////////////////////// | |
176 | // 2. fill HitData and FirstHitInBin | |
177 | //////////////////////////////////// | |
b8139972 | 178 | |
4acc2401 | 179 | for ( int rowIndex = 0; rowIndex < data.FirstRow(); ++rowIndex ) { |
180 | AliHLTTPCCARow &row = fRows[rowIndex]; | |
181 | row.fGrid.CreateEmpty(); | |
182 | row.fNHits = 0; | |
183 | row.fFullSize = 0; | |
184 | row.fHitNumberOffset = 0; | |
185 | row.fFirstHitInBinOffset = 0; | |
186 | ||
187 | row.fHy0 = 0.f; | |
188 | row.fHz0 = 0.f; | |
189 | row.fHstepY = 1.f; | |
190 | row.fHstepZ = 1.f; | |
191 | row.fHstepYi = 1.f; | |
192 | row.fHstepZi = 1.f; | |
193 | } | |
194 | for ( int rowIndex = data.LastRow() + 1; rowIndex < 160; ++rowIndex ) { | |
195 | AliHLTTPCCARow &row = fRows[rowIndex]; | |
196 | row.fGrid.CreateEmpty(); | |
197 | row.fNHits = 0; | |
198 | row.fFullSize = 0; | |
199 | row.fHitNumberOffset = 0; | |
200 | row.fFirstHitInBinOffset = 0; | |
201 | ||
202 | row.fHy0 = 0.f; | |
203 | row.fHz0 = 0.f; | |
204 | row.fHstepY = 1.f; | |
205 | row.fHstepZ = 1.f; | |
206 | row.fHstepYi = 1.f; | |
207 | row.fHstepZi = 1.f; | |
208 | } | |
209 | ||
6de2bc40 | 210 | |
4acc2401 | 211 | AliHLTResizableArray<AliHLTTPCCAHit> binSortedHits( fNumberOfHits ); |
212 | ||
213 | int gridContentOffset = 0; | |
214 | ||
215 | int binCreationMemorySize = 103 * 2 + fNumberOfHits; | |
216 | AliHLTResizableArray<unsigned short> binCreationMemory( binCreationMemorySize ); | |
217 | ||
218 | for ( int rowIndex = data.FirstRow(); rowIndex <= data.LastRow(); ++rowIndex ) { | |
219 | AliHLTTPCCARow &row = fRows[rowIndex]; | |
4acc2401 | 220 | row.fNHits = data.NumberOfClusters( rowIndex ); |
221 | assert( row.fNHits < ( 1 << sizeof( unsigned short ) * 8 ) ); | |
222 | row.fHitNumberOffset = data.RowOffset( rowIndex ); | |
223 | row.fFirstHitInBinOffset = gridContentOffset; | |
224 | ||
6de2bc40 | 225 | CreateGrid( &row, data ); |
4acc2401 | 226 | const AliHLTTPCCAGrid &grid = row.fGrid; |
227 | const int numberOfBins = grid.N(); | |
228 | ||
229 | int binCreationMemorySizeNew; | |
230 | if ( ( binCreationMemorySizeNew = numberOfBins * 2 + 6 + row.fNHits ) > binCreationMemorySize ) { | |
231 | binCreationMemorySize = binCreationMemorySizeNew; | |
232 | binCreationMemory.Resize( binCreationMemorySize ); | |
233 | } | |
b8139972 | 234 | |
4acc2401 | 235 | AliHLTArray<unsigned short> c = binCreationMemory; // number of hits in all previous bins |
236 | AliHLTArray<unsigned short> bins = c + ( numberOfBins + 3 ); // cache for the bin index for every hit in this row | |
237 | AliHLTArray<unsigned short> filled = bins + row.fNHits; // counts how many hits there are per bin | |
238 | ||
239 | for ( unsigned int bin = 0; bin < row.fGrid.N() + 3; ++bin ) { | |
240 | filled[bin] = 0; // initialize filled[] to 0 | |
241 | } | |
b8139972 | 242 | |
4acc2401 | 243 | for ( int hitIndex = 0; hitIndex < row.fNHits; ++hitIndex ) { |
244 | const int globalHitIndex = row.fHitNumberOffset + hitIndex; | |
245 | const unsigned short bin = row.fGrid.GetBin( data.Y( globalHitIndex ), data.Z( globalHitIndex ) ); | |
246 | bins[hitIndex] = bin; | |
247 | ++filled[bin]; | |
248 | } | |
249 | ||
250 | unsigned short n = 0; | |
251 | for ( int bin = 0; bin < numberOfBins + 3; ++bin ) { | |
252 | c[bin] = n; | |
253 | n += filled[bin]; | |
254 | } | |
255 | ||
256 | for ( int hitIndex = 0; hitIndex < row.fNHits; ++hitIndex ) { | |
257 | const unsigned short bin = bins[hitIndex]; | |
258 | --filled[bin]; | |
259 | const unsigned short ind = c[bin] + filled[bin]; // generate an index for this hit that is >= c[bin] and < c[bin + 1] | |
260 | const int globalBinsortedIndex = row.fHitNumberOffset + ind; | |
261 | const int globalHitIndex = row.fHitNumberOffset + hitIndex; | |
262 | ||
263 | // allows to find the global hit index / coordinates from a global bin sorted hit index | |
264 | fClusterDataIndex[globalBinsortedIndex] = globalHitIndex; | |
265 | binSortedHits[globalBinsortedIndex].SetY( data.Y( globalHitIndex ) ); | |
266 | binSortedHits[globalBinsortedIndex].SetZ( data.Z( globalHitIndex ) ); | |
267 | } | |
268 | ||
269 | PackHitData( &row, binSortedHits ); | |
270 | ||
271 | for ( int i = 0; i < numberOfBins; ++i ) { | |
272 | fFirstHitInBin[row.fFirstHitInBinOffset + i] = c[i]; // global bin-sorted hit index | |
273 | } | |
274 | const unsigned short a = c[numberOfBins]; | |
275 | // grid.N is <= row.fNHits | |
b8139972 | 276 | const int nn = numberOfBins + grid.Ny() + 3; |
4acc2401 | 277 | for ( int i = numberOfBins; i < nn; ++i ) { |
6de2bc40 | 278 | assert( row.fFirstHitInBinOffset + i < 23 * numberOfRows + 4 * fNumberOfHits + 3 ); |
4acc2401 | 279 | fFirstHitInBin[row.fFirstHitInBinOffset + i] = a; |
280 | } | |
281 | ||
282 | row.fFullSize = nn; | |
283 | gridContentOffset += nn; | |
284 | } | |
285 | ||
286 | #if 0 | |
287 | //SG cell finder - test code | |
288 | ||
289 | if ( fTmpHitInputIDs ) delete[] fTmpHitInputIDs; | |
290 | fTmpHitInputIDs = new int [NHits]; | |
291 | const float areaY = .5; | |
292 | const float areaZ = .5; | |
293 | int newRowNHitsTotal = 0; | |
294 | bool *usedHits = new bool [NHits]; | |
295 | for ( int iHit = 0; iHit < NHits; iHit++ ) usedHits[iHit] = 0; | |
296 | for ( int iRow = 0; iRow < fParam.NRows(); iRow++ ) { | |
297 | rowHeaders[iRow*2 ] = newRowNHitsTotal; // new first hit | |
298 | rowHeaders[iRow*2+1] = 0; // new N hits | |
299 | int newRowNHits = 0; | |
300 | int oldRowFirstHit = RowFirstHit[iRow]; | |
301 | int oldRowLastHit = oldRowFirstHit + RowNHits[iRow]; | |
302 | for ( int iHit = oldRowFirstHit; iHit < oldRowLastHit; iHit++ ) { | |
303 | if ( usedHits[iHit] ) continue; | |
304 | float x0 = X[iHit]; | |
305 | float y0 = Y[iHit]; | |
306 | float z0 = Z[iHit]; | |
307 | float cx = x0; | |
308 | float cy = y0; | |
309 | float cz = z0; | |
310 | int nclu = 1; | |
311 | usedHits[iHit] = 1; | |
312 | if ( 0 ) for ( int jHit = iHit + 1; jHit < oldRowLastHit; jHit++ ) {//SG!!! | |
313 | //if( usedHits[jHit] ) continue; | |
314 | float dy = Y[jHit] - y0; | |
315 | float dz = Z[jHit] - z0; | |
316 | if ( CAMath::Abs( dy ) < areaY && CAMath::Abs( dz ) < areaZ ) { | |
317 | cx += X[jHit]; | |
318 | cy += Y[jHit]; | |
319 | cz += Z[jHit]; | |
320 | nclu++; | |
321 | usedHits[jHit] = 1; | |
322 | } | |
323 | } | |
324 | int id = newRowNHitsTotal + newRowNHits; | |
325 | hitsXYZ[id*3+0 ] = cx / nclu; | |
326 | hitsXYZ[id*3+1 ] = cy / nclu; | |
327 | hitsXYZ[id*3+2 ] = cz / nclu; | |
328 | fTmpHitInputIDs[id] = iHit; | |
329 | newRowNHits++; | |
330 | } | |
331 | rowHeaders[iRow*2+1] = newRowNHits; | |
332 | newRowNHitsTotal += newRowNHits; | |
333 | } | |
334 | NHitsTotal() = newRowNHitsTotal; | |
335 | reinterpret_cast<int*>( fInputEvent )[1+fParam.NRows()*2] = newRowNHitsTotal; | |
336 | ||
337 | delete[] usedHits; | |
338 | #endif | |
339 | } | |
340 | ||
341 | void AliHLTTPCCASliceData::ClearHitWeights() | |
342 | { | |
6de2bc40 | 343 | // clear hit weights |
344 | ||
4acc2401 | 345 | #ifdef ENABLE_VECTORIZATION |
346 | const int_v v0( Zero ); | |
347 | const int *const end = fHitWeights + fNumberOfHits; | |
348 | for ( int *mem = fHitWeights; mem < end; mem += v0.Size ) { | |
349 | v0.store( mem ); | |
350 | } | |
351 | #else | |
352 | for ( int i = 0; i < fNumberOfHits; ++i ) { | |
353 | fHitWeights[i] = 0; | |
354 | } | |
355 | #endif | |
356 | } | |
357 | ||
358 | void AliHLTTPCCASliceData::ClearLinks() | |
359 | { | |
6de2bc40 | 360 | // link cleaning |
361 | ||
4acc2401 | 362 | #ifdef ENABLE_VECTORIZATION |
363 | const short_v v0( -1 ); | |
364 | const short *const end1 = fLinkUpData + fNumberOfHits; | |
365 | for ( short *mem = fLinkUpData; mem < end; mem += v0.Size ) { | |
366 | v0.store( mem ); | |
367 | } | |
368 | const short *const end2 = fLinkDownData + fNumberOfHits; | |
369 | for ( short *mem = fLinkDownData; mem < end; mem += v0.Size ) { | |
370 | v0.store( mem ); | |
371 | } | |
372 | #else | |
373 | for ( int i = 0; i < fNumberOfHits; ++i ) { | |
374 | fLinkUpData[i] = -1; | |
375 | } | |
376 | for ( int i = 0; i < fNumberOfHits; ++i ) { | |
377 | fLinkDownData[i] = -1; | |
378 | } | |
379 | #endif | |
380 | } | |
381 |