Line data Source code
1 : // $Id: AliHLTTPCGMMerger.cxx 30732 2009-01-22 23:02:02Z sgorbuno $
2 : // **************************************************************************
3 : // This file is property of and copyright by the ALICE HLT Project *
4 : // ALICE Experiment at CERN, All rights reserved. *
5 : // *
6 : // Primary Authors: Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
7 : // Ivan Kisel <kisel@kip.uni-heidelberg.de> *
8 : // for The ALICE HLT Project. *
9 : // *
10 : // Permission to use, copy, modify and distribute this software and its *
11 : // documentation strictly for non-commercial purposes is hereby granted *
12 : // without fee, provided that the above copyright notice appears in all *
13 : // copies and that both the copyright notice and this permission notice *
14 : // appear in the supporting documentation. The authors make no claims *
15 : // about the suitability of this software for any purpose. It is *
16 : // provided "as is" without express or implied warranty. *
17 : // *
18 : //***************************************************************************
19 :
20 : #include <stdio.h>
21 : #include "AliHLTTPCCASliceOutTrack.h"
22 : #include "AliHLTTPCCATracker.h"
23 : #include "AliHLTTPCCATrackParam.h"
24 : #include "AliHLTTPCGMCluster.h"
25 :
26 : #include "AliHLTTPCGMMerger.h"
27 :
28 : #include "AliHLTTPCCAMath.h"
29 : #include "TStopwatch.h"
30 :
31 : #include "AliHLTTPCCATrackParam.h"
32 : #include "AliHLTTPCCASliceOutput.h"
33 : #include "AliHLTTPCGMMergedTrack.h"
34 : #include "AliHLTTPCCADataCompressor.h"
35 : #include "AliHLTTPCCAParam.h"
36 : #include "AliHLTTPCCATrackLinearisation.h"
37 : #include "AliHLTTPCCADataCompressor.h"
38 :
39 : #include "AliHLTTPCGMTrackParam.h"
40 : #include "AliHLTTPCGMTrackLinearisation.h"
41 : #include "AliHLTTPCGMSliceTrack.h"
42 : #include "AliHLTTPCGMBorderTrack.h"
43 : #include <cmath>
44 :
45 : #include <algorithm>
46 :
47 : #include "AliHLTTPCCAGPUConfig.h"
48 : #include "MemoryAssignmentHelpers.h"
49 :
50 : #define GLOBAL_TRACKS_SPECIAL_TREATMENT
51 :
52 : AliHLTTPCGMMerger::AliHLTTPCGMMerger()
53 : :
54 0 : fSliceParam(),
55 0 : fNOutputTracks( 0 ),
56 0 : fNOutputTrackClusters( 0 ),
57 0 : fOutputTracks( 0 ),
58 0 : fOutputClusterIds(0),
59 0 : fSliceTrackInfos( 0 ),
60 0 : fMaxSliceTracks(0),
61 0 : fClusterX(0),
62 0 : fClusterY(0),
63 0 : fClusterZ(0),
64 0 : fClusterRowType(0),
65 0 : fClusterAngle(0),
66 0 : fBorderMemory(0),
67 0 : fBorderRangeMemory(0),
68 0 : fGPUTracker(NULL),
69 0 : fDebugLevel(0),
70 0 : fNClusters(0)
71 0 : {
72 : //* constructor
73 :
74 0 : for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
75 0 : fNextSliceInd[iSlice] = iSlice + 1;
76 0 : fPrevSliceInd[iSlice] = iSlice - 1;
77 : }
78 : int mid = fgkNSlices / 2 - 1 ;
79 : int last = fgkNSlices - 1 ;
80 0 : fNextSliceInd[ mid ] = 0;
81 0 : fPrevSliceInd[ 0 ] = mid; fNextSliceInd[ last ] = fgkNSlices / 2;
82 0 : fPrevSliceInd[ fgkNSlices/2 ] = last;
83 : {
84 : const double kCLight = 0.000299792458;
85 0 : double constBz = fSliceParam.BzkG() * kCLight;
86 :
87 0 : fPolinomialFieldBz[0] = constBz * ( 0.999286 );
88 0 : fPolinomialFieldBz[1] = constBz * ( -4.54386e-7 );
89 0 : fPolinomialFieldBz[2] = constBz * ( 2.32950e-5 );
90 0 : fPolinomialFieldBz[3] = constBz * ( -2.99912e-7 );
91 0 : fPolinomialFieldBz[4] = constBz * ( -2.03442e-8 );
92 0 : fPolinomialFieldBz[5] = constBz * ( 9.71402e-8 );
93 : }
94 :
95 0 : Clear();
96 0 : }
97 :
98 :
99 : AliHLTTPCGMMerger::AliHLTTPCGMMerger(const AliHLTTPCGMMerger&)
100 : :
101 0 : fSliceParam(),
102 0 : fNOutputTracks( 0 ),
103 0 : fNOutputTrackClusters( 0 ),
104 0 : fOutputTracks( 0 ),
105 0 : fOutputClusterIds(0),
106 0 : fSliceTrackInfos( 0 ),
107 0 : fMaxSliceTracks(0),
108 0 : fClusterX(0),
109 0 : fClusterY(0),
110 0 : fClusterZ(0),
111 0 : fClusterRowType(0),
112 0 : fClusterAngle(0),
113 0 : fBorderMemory(0),
114 0 : fBorderRangeMemory(0),
115 0 : fGPUTracker(NULL),
116 0 : fDebugLevel(0),
117 0 : fNClusters(0)
118 0 : {
119 : //* dummy
120 0 : for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
121 0 : fNextSliceInd[iSlice] = 0;
122 0 : fPrevSliceInd[iSlice] = 0;
123 : }
124 : {
125 : const double kCLight = 0.000299792458;
126 0 : double constBz = fSliceParam.BzkG() * kCLight;
127 :
128 0 : fPolinomialFieldBz[0] = constBz * ( 0.999286 );
129 0 : fPolinomialFieldBz[1] = constBz * ( -4.54386e-7 );
130 0 : fPolinomialFieldBz[2] = constBz * ( 2.32950e-5 );
131 0 : fPolinomialFieldBz[3] = constBz * ( -2.99912e-7 );
132 0 : fPolinomialFieldBz[4] = constBz * ( -2.03442e-8 );
133 0 : fPolinomialFieldBz[5] = constBz * ( 9.71402e-8 );
134 : }
135 0 : Clear();
136 0 : }
137 :
138 : const AliHLTTPCGMMerger &AliHLTTPCGMMerger::operator=(const AliHLTTPCGMMerger&) const
139 : {
140 : //* dummy
141 0 : return *this;
142 : }
143 :
144 :
145 : AliHLTTPCGMMerger::~AliHLTTPCGMMerger()
146 0 : {
147 : //* destructor
148 0 : ClearMemory();
149 0 : }
150 :
151 : void AliHLTTPCGMMerger::Clear()
152 : {
153 0 : for ( int i = 0; i < fgkNSlices; ++i ) {
154 0 : fkSlices[i] = 0;
155 0 : fSliceNTrackInfos[ i ] = 0;
156 0 : fSliceTrackInfoStart[ i ] = 0;
157 : }
158 0 : ClearMemory();
159 0 : }
160 :
161 :
162 : void AliHLTTPCGMMerger::ClearMemory()
163 : {
164 0 : if (fOutputClusterIds) delete[] fOutputClusterIds;
165 0 : if (fSliceTrackInfos) delete[] fSliceTrackInfos;
166 0 : if (!(fGPUTracker && fGPUTracker->IsInitialized()))
167 : {
168 0 : if (fOutputTracks) delete[] fOutputTracks;
169 0 : if (fClusterX) delete[] fClusterX;
170 0 : if (fClusterY) delete[] fClusterY;
171 0 : if (fClusterZ) delete[] fClusterZ;
172 0 : if (fClusterRowType) delete[] fClusterRowType;
173 0 : if (fClusterAngle) delete[] fClusterAngle;
174 : }
175 0 : if (fBorderMemory) delete[] fBorderMemory;
176 0 : if (fBorderRangeMemory) delete[] fBorderRangeMemory;
177 :
178 0 : fNOutputTracks = 0;
179 0 : fOutputTracks = 0;
180 0 : fOutputClusterIds = 0;
181 0 : fSliceTrackInfos = 0;
182 0 : fMaxSliceTracks = 0;
183 0 : fClusterX = 0;
184 0 : fClusterY = 0;
185 0 : fClusterZ = 0;
186 0 : fClusterRowType = 0;
187 0 : fClusterAngle = 0;
188 0 : fBorderMemory = 0;
189 0 : fBorderRangeMemory = 0;
190 0 : }
191 :
192 :
193 : void AliHLTTPCGMMerger::SetSliceData( int index, const AliHLTTPCCASliceOutput *sliceData )
194 : {
195 0 : fkSlices[index] = sliceData;
196 0 : }
197 :
198 :
199 : bool AliHLTTPCGMMerger::Reconstruct()
200 : {
201 : //* main merging routine
202 :
203 : {
204 : const double kCLight = 0.000299792458;
205 0 : double constBz = fSliceParam.BzkG() * kCLight;
206 :
207 0 : fPolinomialFieldBz[0] = constBz * ( 0.999286 );
208 0 : fPolinomialFieldBz[1] = constBz * ( -4.54386e-7 );
209 0 : fPolinomialFieldBz[2] = constBz * ( 2.32950e-5 );
210 0 : fPolinomialFieldBz[3] = constBz * ( -2.99912e-7 );
211 0 : fPolinomialFieldBz[4] = constBz * ( -2.03442e-8 );
212 0 : fPolinomialFieldBz[5] = constBz * ( 9.71402e-8 );
213 : }
214 :
215 : int nIter = 1;
216 0 : TStopwatch timer;
217 : #ifdef HLTCA_STANDALONE
218 : unsigned long long int a, b, c, d, e, f, g;
219 : AliHLTTPCCATracker::StandaloneQueryFreq(&g);
220 : #endif
221 : //cout<<"Merger..."<<endl;
222 0 : for( int iter=0; iter<nIter; iter++ ){
223 0 : if( !AllocateMemory() ) return 0;
224 : #ifdef HLTCA_STANDALONE
225 : AliHLTTPCCATracker::StandaloneQueryTime(&a);
226 : #endif
227 0 : UnpackSlices();
228 : #ifdef HLTCA_STANDALONE
229 : AliHLTTPCCATracker::StandaloneQueryTime(&b);
230 : #endif
231 0 : MergeWithingSlices();
232 : #ifdef HLTCA_STANDALONE
233 : AliHLTTPCCATracker::StandaloneQueryTime(&c);
234 : #endif
235 0 : MergeSlices();
236 : #ifdef HLTCA_STANDALONE
237 : AliHLTTPCCATracker::StandaloneQueryTime(&d);
238 : #endif
239 0 : CollectMergedTracks();
240 : #ifdef HLTCA_STANDALONE
241 : AliHLTTPCCATracker::StandaloneQueryTime(&e);
242 : #endif
243 0 : Refit();
244 : #ifdef HLTCA_STANDALONE
245 : AliHLTTPCCATracker::StandaloneQueryTime(&f);
246 : if (fDebugLevel > 0)
247 : {
248 : printf("Merge Time:\tUnpack Slices:\t%lld us\n", (b - a) * 1000000 / g);
249 : printf("\t\tMerge Within:\t%lld us\n", (c - b) * 1000000 / g);
250 : printf("\t\tMerge Slices:\t%lld us\n", (d - c) * 1000000 / g);
251 : printf("\t\tCollect:\t%lld us\n", (e - d) * 1000000 / g);
252 : printf("\t\tRefit:\t\t%lld us\n", (f - e) * 1000000 / g);
253 : }
254 : int newTracks = 0;
255 : for (int i = 0;i < fNOutputTracks;i++) if (fOutputTracks[i].OK()) newTracks++;
256 : printf("Output Tracks: %d\n", newTracks);
257 : #endif
258 : }
259 0 : timer.Stop();
260 : //cout<<"\nMerger time = "<<timer.CpuTime()*1.e3/nIter<<" ms\n"<<endl;
261 :
262 0 : return 1;
263 0 : }
264 :
265 :
266 :
267 : bool AliHLTTPCGMMerger::AllocateMemory()
268 : {
269 : //* memory allocation
270 :
271 0 : ClearMemory();
272 :
273 : int nTracks = 0;
274 0 : fNClusters = 0;
275 0 : fMaxSliceTracks = 0;
276 :
277 0 : for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
278 0 : if ( !fkSlices[iSlice] ) continue;
279 0 : nTracks += fkSlices[iSlice]->NTracks();
280 0 : fNClusters += fkSlices[iSlice]->NTrackClusters();
281 0 : if( fMaxSliceTracks < fkSlices[iSlice]->NTracks() ) fMaxSliceTracks = fkSlices[iSlice]->NTracks();
282 : }
283 :
284 : //cout<<"\nMerger: input "<<nTracks<<" tracks, "<<nClusters<<" clusters"<<endl;
285 :
286 0 : fOutputClusterIds = new UInt_t[fNClusters];
287 0 : fSliceTrackInfos = new AliHLTTPCGMSliceTrack[nTracks];
288 0 : if (fGPUTracker && fGPUTracker->IsInitialized())
289 : {
290 0 : char* basemem = fGPUTracker->MergerBaseMemory();
291 0 : AssignMemory(fClusterX, basemem, fNClusters);
292 0 : AssignMemory(fClusterY, basemem, fNClusters);
293 0 : AssignMemory(fClusterZ, basemem, fNClusters);
294 0 : AssignMemory(fClusterAngle, basemem, fNClusters);
295 0 : AssignMemory(fClusterRowType, basemem, fNClusters);
296 0 : AssignMemory(fOutputTracks, basemem, nTracks);
297 0 : }
298 : else
299 : {
300 0 : fOutputTracks = new AliHLTTPCGMMergedTrack[nTracks];
301 0 : fClusterX = new float[fNClusters];
302 0 : fClusterY = new float[fNClusters];
303 0 : fClusterZ = new float[fNClusters];
304 0 : fClusterRowType = new UInt_t[fNClusters];
305 0 : fClusterAngle = new float[fNClusters];
306 : }
307 0 : fBorderMemory = new AliHLTTPCGMBorderTrack[fMaxSliceTracks*2];
308 0 : fBorderRangeMemory = new AliHLTTPCGMBorderTrack::Range[fMaxSliceTracks*2];
309 :
310 0 : return ( ( fOutputTracks!=NULL )
311 0 : && ( fOutputClusterIds!=NULL )
312 0 : && ( fSliceTrackInfos!=NULL )
313 0 : && ( fClusterX!=NULL )
314 0 : && ( fClusterY!=NULL )
315 0 : && ( fClusterZ!=NULL )
316 0 : && ( fClusterRowType!=NULL )
317 0 : && ( fClusterAngle!=NULL )
318 0 : && ( fBorderMemory!=NULL )
319 0 : && ( fBorderRangeMemory!=NULL )
320 : );
321 : }
322 :
323 :
324 :
325 : void AliHLTTPCGMMerger::UnpackSlices()
326 : {
327 : //* unpack the cluster information from the slice tracks and initialize track info array
328 :
329 : int nTracksCurrent = 0;
330 :
331 0 : const AliHLTTPCCASliceOutTrack* firstGlobalTracks[fgkNSlices];
332 :
333 0 : for( int i=0; i<fgkNSlices; i++) firstGlobalTracks[i] = 0;
334 :
335 : #ifdef GLOBAL_TRACKS_SPECIAL_TREATMENT
336 0 : const int kMaxTrackIdInSlice = AliHLTTPCCASliceOutTrack::MaxTrackId();
337 0 : int* TrackIds = (int*) malloc(kMaxTrackIdInSlice * fgkNSlices * sizeof(int));
338 0 : for (int i = 0;i < kMaxTrackIdInSlice * fgkNSlices;i++) TrackIds[i] = -1;
339 : #endif
340 :
341 0 : for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
342 :
343 0 : fSliceTrackInfoStart[ iSlice ] = nTracksCurrent;
344 :
345 0 : fSliceNTrackInfos[ iSlice ] = 0;
346 :
347 0 : if ( !fkSlices[iSlice] ) continue;
348 :
349 0 : float alpha = fSliceParam.Alpha( iSlice );
350 :
351 0 : const AliHLTTPCCASliceOutput &slice = *( fkSlices[iSlice] );
352 0 : const AliHLTTPCCASliceOutTrack *sliceTr = slice.GetFirstTrack();
353 :
354 : #ifdef GLOBAL_TRACKS_SPECIAL_TREATMENT
355 0 : for ( int itr = 0; itr < slice.NLocalTracks(); itr++, sliceTr = sliceTr->GetNextTrack() ) {
356 : #else
357 : for ( int itr = 0; itr < slice.NTracks(); itr++, sliceTr = sliceTr->GetNextTrack() ) {
358 : #endif
359 0 : AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[nTracksCurrent];
360 0 : track.Set( sliceTr, alpha );
361 0 : if( !track.FilterErrors( fSliceParam, .999 ) ) continue;
362 0 : track.SetPrevNeighbour( -1 );
363 0 : track.SetNextNeighbour( -1 );
364 0 : track.SetSliceNeighbour( -1 );
365 0 : track.SetUsed( 0 );
366 : #ifdef GLOBAL_TRACKS_SPECIAL_TREATMENT
367 0 : track.SetGlobalTrackId(0, -1);
368 0 : track.SetGlobalTrackId(1, -1);
369 0 : TrackIds[iSlice * kMaxTrackIdInSlice + sliceTr->LocalTrackId()] = nTracksCurrent;
370 : #endif
371 0 : nTracksCurrent++;
372 0 : fSliceNTrackInfos[ iSlice ]++;
373 0 : }
374 0 : firstGlobalTracks[iSlice] = sliceTr;
375 :
376 : //std::cout<<"Unpack slice "<<iSlice<<": ntracks "<<slice.NTracks()<<"/"<<fSliceNTrackInfos[iSlice]<<std::endl;
377 0 : }
378 : #ifdef GLOBAL_TRACKS_SPECIAL_TREATMENT
379 0 : for (int iSlice = 0;iSlice < fgkNSlices;iSlice++)
380 : {
381 0 : fSliceTrackGlobalInfoStart[iSlice] = nTracksCurrent;
382 0 : fSliceNGlobalTrackInfos[iSlice] = 0;
383 :
384 0 : if ( !fkSlices[iSlice] ) continue;
385 0 : float alpha = fSliceParam.Alpha( iSlice );
386 :
387 0 : const AliHLTTPCCASliceOutput &slice = *( fkSlices[iSlice] );
388 0 : const AliHLTTPCCASliceOutTrack *sliceTr = firstGlobalTracks[iSlice];
389 0 : for (int itr = slice.NLocalTracks();itr < slice.NTracks();itr++, sliceTr = sliceTr->GetNextTrack())
390 : {
391 0 : if (TrackIds[sliceTr->LocalTrackId()] == -1) continue;
392 0 : AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[nTracksCurrent];
393 0 : track.Set( sliceTr, alpha );
394 : //if( !track.FilterErrors( fSliceParam, .999 ) ) continue;
395 0 : track.SetPrevNeighbour( -1 );
396 0 : track.SetNextNeighbour( -1 );
397 0 : track.SetSliceNeighbour( -1 );
398 0 : track.SetUsed( 0 );
399 0 : track.SetLocalTrackId(TrackIds[sliceTr->LocalTrackId()]);
400 0 : nTracksCurrent++;
401 0 : fSliceNGlobalTrackInfos[ iSlice ]++;
402 0 : }
403 0 : }
404 :
405 0 : free(TrackIds);
406 : #endif
407 0 : }
408 :
409 :
410 :
411 :
412 :
413 : void AliHLTTPCGMMerger::MakeBorderTracks( int iSlice, int iBorder, AliHLTTPCGMBorderTrack B[], int &nB )
414 : {
415 : //* prepare slice tracks for merging with next/previous/same sector
416 : //* each track transported to the border line
417 :
418 0 : float fieldBz = fSliceParam.ConstBz();
419 :
420 0 : nB = 0;
421 :
422 0 : float dAlpha = fSliceParam.DAlpha() / 2;
423 : float x0 = 0;
424 :
425 0 : if ( iBorder == 0 ) { // transport to the left age of the sector and rotate horisontally
426 0 : dAlpha = dAlpha - CAMath::Pi() / 2 ;
427 0 : } else if ( iBorder == 1 ) { // transport to the right age of the sector and rotate horisontally
428 0 : dAlpha = -dAlpha - CAMath::Pi() / 2 ;
429 0 : } else if ( iBorder == 2 ) { // transport to the left age of the sector and rotate vertically
430 : dAlpha = dAlpha;
431 0 : x0 = fSliceParam.RowX( 63 );
432 0 : } else if ( iBorder == 3 ) { // transport to the right age of the sector and rotate vertically
433 0 : dAlpha = -dAlpha;
434 0 : x0 = fSliceParam.RowX( 63 );
435 0 : } else if ( iBorder == 4 ) { // transport to the middle of the sector, w/o rotation
436 : dAlpha = 0;
437 0 : x0 = fSliceParam.RowX( 63 );
438 0 : }
439 :
440 0 : const float maxSin = CAMath::Sin( 60. / 180.*CAMath::Pi() );
441 0 : float cosAlpha = AliHLTTPCCAMath::Cos( dAlpha );
442 0 : float sinAlpha = AliHLTTPCCAMath::Sin( dAlpha );
443 :
444 0 : for ( int itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ ) {
445 :
446 0 : AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[ fSliceTrackInfoStart[iSlice] + itr ];
447 :
448 0 : if( track.Used() ) continue;
449 0 : AliHLTTPCGMBorderTrack &b = B[nB];
450 :
451 0 : if( track.TransportToXAlpha( x0, sinAlpha, cosAlpha, fieldBz, b, maxSin)){
452 0 : b.SetTrackID( itr );
453 0 : b.SetNClusters( track.NClusters() );
454 0 : for (int i = 0;i < 4;i++) if (fabs(b.Cov()[i]) >= 5.0) b.SetCov(i, 5.0);
455 0 : if (fabs(b.Cov()[4]) >= 0.5) b.SetCov(4, 0.5);
456 0 : nB++;
457 0 : }
458 0 : }
459 0 : }
460 :
461 :
462 : void AliHLTTPCGMMerger::MergeBorderTracks ( int iSlice1, AliHLTTPCGMBorderTrack B1[], int N1,
463 : int iSlice2, AliHLTTPCGMBorderTrack B2[], int N2 )
464 : {
465 : //* merge two sets of tracks
466 :
467 : //std::cout<<" Merge slices "<<iSlice1<<"+"<<iSlice2<<": tracks "<<N1<<"+"<<N2<<std::endl;
468 : int statAll=0, statMerged=0;
469 : float factor2ys = 1.5;//1.5;//SG!!!
470 : float factor2zt = 1.5;//1.5;//SG!!!
471 : float factor2k = 2.0;//2.2;
472 :
473 : factor2k = 3.5 * 3.5 * factor2k * factor2k;
474 : factor2ys = 3.5 * 3.5 * factor2ys * factor2ys;
475 : factor2zt = 3.5 * 3.5 * factor2zt * factor2zt;
476 :
477 : int minNPartHits = 10;//SG!!!
478 : int minNTotalHits = 20;
479 :
480 0 : AliHLTTPCGMBorderTrack::Range *range1 = fBorderRangeMemory;
481 0 : AliHLTTPCGMBorderTrack::Range *range2 = fBorderRangeMemory + N1;
482 :
483 0 : bool sameSlice = (iSlice1 == iSlice2);
484 : {
485 0 : for ( int itr = 0; itr < N1; itr++ ){
486 0 : AliHLTTPCGMBorderTrack &b = B1[itr];
487 : // if( iSlice1==7 && iSlice2==8 ){
488 : //cout<<b.TrackID()<<": "<<b.Cov()[0]<<" "<<b.Cov()[1]<<endl;
489 : //}
490 0 : float d = 3.5*sqrt(b.Cov()[1]);
491 0 : range1[itr].fId = itr;
492 0 : range1[itr].fMin = b.Par()[1] - d;
493 0 : range1[itr].fMax = b.Par()[1] + d;
494 : }
495 0 : std::sort(range1,range1+N1,AliHLTTPCGMBorderTrack::Range::CompMin);
496 0 : if( sameSlice ){
497 0 : for(int i=0; i<N1; i++) range2[i]= range1[i];
498 0 : std::sort(range2,range2+N1,AliHLTTPCGMBorderTrack::Range::CompMax);
499 : N2 = N1;
500 : B2 = B1;
501 0 : }else{
502 0 : for ( int itr = 0; itr < N2; itr++ ){
503 0 : AliHLTTPCGMBorderTrack &b = B2[itr];
504 0 : float d = 3.5*sqrt(b.Cov()[1]);
505 0 : range2[itr].fId = itr;
506 0 : range2[itr].fMin = b.Par()[1] - d;
507 0 : range2[itr].fMax = b.Par()[1] + d;
508 : }
509 0 : std::sort(range2,range2+N2,AliHLTTPCGMBorderTrack::Range::CompMax);
510 : }
511 : }
512 :
513 : int i2 = 0;
514 0 : for ( int i1 = 0; i1 < N1; i1++ ) {
515 :
516 0 : AliHLTTPCGMBorderTrack::Range r1 = range1[i1];
517 0 : while( i2<N2 && range2[i2].fMax< r1.fMin ) i2++;
518 :
519 0 : AliHLTTPCGMBorderTrack &b1 = B1[r1.fId];
520 0 : if ( b1.NClusters() < minNPartHits ) continue;
521 : int iBest2 = -1;
522 : int lBest2 = 0;
523 0 : statAll++;
524 0 : for( int k2 = i2; k2<N2; k2++){
525 :
526 0 : AliHLTTPCGMBorderTrack::Range r2 = range2[k2];
527 0 : if( r2.fMin > r1.fMax ) break;
528 0 : if( sameSlice && (r1.fId >= r2.fId) ) continue;
529 : // do check
530 0 : AliHLTTPCGMBorderTrack &b2 = B2[r2.fId];
531 0 : if ( b2.NClusters() < lBest2 ) continue;
532 :
533 0 : if( !b1.CheckChi2Y(b2, factor2ys ) ) continue;
534 : //if( !b1.CheckChi2Z(b2, factor2zt ) ) continue;
535 0 : if( !b1.CheckChi2QPt(b2, factor2k ) ) continue;
536 0 : if( !b1.CheckChi2YS(b2, factor2ys ) ) continue;
537 0 : if( !b1.CheckChi2ZT(b2, factor2zt ) ) continue;
538 0 : if ( b2.NClusters() < minNPartHits ) continue;
539 0 : if ( b1.NClusters() + b2.NClusters() < minNTotalHits ) continue;
540 :
541 0 : lBest2 = b2.NClusters();
542 0 : iBest2 = b2.TrackID();
543 0 : }
544 :
545 0 : if ( iBest2 < 0 ) continue;
546 0 : statMerged++;
547 0 : AliHLTTPCGMSliceTrack &newTrack1 = fSliceTrackInfos[fSliceTrackInfoStart[iSlice1] + b1.TrackID() ];
548 0 : AliHLTTPCGMSliceTrack &newTrack2 = fSliceTrackInfos[fSliceTrackInfoStart[iSlice2] + iBest2 ];
549 :
550 0 : int old1 = newTrack2.PrevNeighbour();
551 :
552 0 : if ( old1 >= 0 ) {
553 0 : AliHLTTPCGMSliceTrack &oldTrack1 = fSliceTrackInfos[fSliceTrackInfoStart[iSlice1] + old1];
554 0 : if ( oldTrack1.NClusters() < newTrack1.NClusters() ) {
555 0 : newTrack2.SetPrevNeighbour( -1 );
556 0 : oldTrack1.SetNextNeighbour( -1 );
557 0 : } else continue;
558 0 : }
559 0 : int old2 = newTrack1.NextNeighbour();
560 0 : if ( old2 >= 0 ) {
561 0 : AliHLTTPCGMSliceTrack &oldTrack2 = fSliceTrackInfos[fSliceTrackInfoStart[iSlice2] + old2];
562 0 : if ( oldTrack2.NClusters() < newTrack2.NClusters() ) {
563 0 : oldTrack2.SetPrevNeighbour( -1 );
564 0 : } else continue;
565 0 : }
566 0 : newTrack1.SetNextNeighbour( iBest2 );
567 0 : newTrack2.SetPrevNeighbour( b1.TrackID() );
568 0 : }
569 : //cout<<"slices "<<iSlice1<<","<<iSlice2<<": all "<<statAll<<" merged "<<statMerged<<endl;
570 0 : }
571 :
572 :
573 : void AliHLTTPCGMMerger::MergeWithingSlices()
574 : {
575 : //* merge track segments withing one slice
576 :
577 0 : float x0 = fSliceParam.RowX( 63 );
578 0 : const float maxSin = CAMath::Sin( 60. / 180.*CAMath::Pi() );
579 :
580 0 : for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
581 :
582 : int nBord = 0;
583 0 : for ( int itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ ) {
584 0 : AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[ fSliceTrackInfoStart[iSlice] + itr ];
585 : //track.SetPrevNeighbour( -1 );
586 : //track.SetNextNeighbour( -1 );
587 : //track.SetSliceNeighbour( -1 );
588 : //track.SetUsed(0);
589 :
590 0 : AliHLTTPCGMBorderTrack &b = fBorderMemory[nBord];
591 0 : if( track.TransportToX( x0, fSliceParam.ConstBz(), b, maxSin) ){
592 0 : b.SetTrackID( itr );
593 0 : b.SetNClusters( track.NClusters() );
594 0 : nBord++;
595 0 : }
596 : }
597 :
598 0 : MergeBorderTracks( iSlice, fBorderMemory, nBord, iSlice, fBorderMemory, nBord );
599 :
600 0 : for ( int itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ ) {
601 0 : AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[ fSliceTrackInfoStart[iSlice] + itr];
602 0 : if( track.PrevNeighbour()>=0 || track.Used() ) continue;
603 0 : int jtr = track.NextNeighbour();
604 0 : track.SetSliceNeighbour( jtr );
605 0 : track.SetNextNeighbour(-1);
606 0 : while( jtr>=0 ){
607 0 : AliHLTTPCGMSliceTrack &trackN = fSliceTrackInfos[ fSliceTrackInfoStart[iSlice] + jtr];
608 0 : if( trackN.NClusters()>track.NClusters() ) track.CopyParamFrom(trackN);
609 0 : trackN.SetUsed(2);
610 0 : jtr = trackN.NextNeighbour();
611 0 : trackN.SetSliceNeighbour( jtr );
612 0 : trackN.SetNextNeighbour(-1);
613 0 : trackN.SetPrevNeighbour(-1);
614 : }
615 0 : }
616 : }
617 0 : }
618 :
619 :
620 :
621 :
622 : void AliHLTTPCGMMerger::MergeSlices()
623 : {
624 : //* track merging between slices
625 :
626 : //for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
627 : //for ( int itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ ) {
628 : //AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[ fSliceTrackInfoStart[iSlice] + itr ];
629 : //track.SetPrevNeighbour( -1 );
630 : //track.SetNextNeighbour( -1 );
631 : //}
632 : //}
633 :
634 : AliHLTTPCGMBorderTrack
635 0 : *bCurr = fBorderMemory,
636 0 : *bNext = fBorderMemory + fMaxSliceTracks;
637 :
638 0 : for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
639 0 : int jSlice = fNextSliceInd[iSlice];
640 0 : int nCurr = 0, nNext = 0;
641 0 : MakeBorderTracks( iSlice, 2, bCurr, nCurr );
642 0 : MakeBorderTracks( jSlice, 3, bNext, nNext );
643 0 : MergeBorderTracks( iSlice, bCurr, nCurr, jSlice, bNext, nNext );
644 0 : MakeBorderTracks( iSlice, 0, bCurr, nCurr );
645 0 : MakeBorderTracks( jSlice, 1, bNext, nNext );
646 0 : MergeBorderTracks( iSlice, bCurr, nCurr, jSlice, bNext, nNext );
647 0 : }
648 0 : }
649 :
650 :
651 :
652 :
653 :
654 : void AliHLTTPCGMMerger::CollectMergedTracks()
655 : {
656 : //*
657 :
658 : //for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
659 : //for ( int itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ ) {
660 : //AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[ fSliceTrackInfoStart[iSlice] + itr ];
661 : //if( track.Used()!=2 ) track.SetUsed(0);
662 : //}
663 : //}
664 :
665 : //Resolve connections for global tracks first
666 : #ifdef GLOBAL_TRACKS_SPECIAL_TREATMENT
667 0 : for (int iSlice = 0;iSlice < fgkNSlices;iSlice++)
668 : {
669 0 : for (int itr = 0;itr < fSliceNGlobalTrackInfos[iSlice];itr++)
670 : {
671 0 : AliHLTTPCGMSliceTrack &globalTrack = fSliceTrackInfos[fSliceTrackGlobalInfoStart[iSlice] + itr];
672 0 : AliHLTTPCGMSliceTrack &localTrack = fSliceTrackInfos[globalTrack.LocalTrackId()];
673 0 : localTrack.SetGlobalTrackId(localTrack.GlobalTrackId(0) != -1, fSliceTrackGlobalInfoStart[iSlice] + itr);
674 : }
675 : }
676 : #endif
677 :
678 : //Now collect the merged tracks
679 0 : fNOutputTracks = 0;
680 : int nOutTrackClusters = 0;
681 : const int kMaxParts = 400;
682 : const int kMaxClusters = 1000;
683 :
684 0 : const AliHLTTPCGMSliceTrack *trackParts[kMaxParts];
685 :
686 0 : for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
687 :
688 0 : for ( int itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ ) {
689 :
690 0 : AliHLTTPCGMSliceTrack &track = fSliceTrackInfos[fSliceTrackInfoStart[iSlice] + itr];
691 :
692 0 : if ( track.Used() ) continue;
693 0 : if ( track.PrevNeighbour() >= 0 ) continue;
694 : int nParts = 0;
695 : int jSlice = iSlice;
696 : AliHLTTPCGMSliceTrack *trbase = &track, *tr = &track;
697 0 : tr->SetUsed( 1 );
698 : do{
699 0 : if( nParts >= kMaxParts ) break;
700 0 : trackParts[nParts++] = tr;
701 : #ifdef GLOBAL_TRACKS_SPECIAL_TREATMENT
702 0 : for (int i = 0;i < 2;i++) if (tr->GlobalTrackId(i) != -1) trackParts[nParts++] = &fSliceTrackInfos[tr->GlobalTrackId(i)];
703 : #endif
704 0 : int jtr = tr->SliceNeighbour();
705 0 : if( jtr >= 0 ) {
706 0 : tr = &(fSliceTrackInfos[fSliceTrackInfoStart[jSlice] + jtr]);
707 0 : tr->SetUsed( 2 );
708 0 : continue;
709 : }
710 0 : jtr = trbase->NextNeighbour();
711 0 : if( jtr>=0 ){
712 0 : jSlice = fNextSliceInd[jSlice];
713 0 : trbase = &(fSliceTrackInfos[fSliceTrackInfoStart[jSlice] + jtr]);
714 : tr = trbase;
715 0 : if( tr->Used() ) break;
716 0 : tr->SetUsed( 1 );
717 0 : continue;
718 : }
719 0 : break;
720 : }while(1);
721 :
722 : // unpack and sort clusters
723 :
724 0 : std::sort(trackParts, trackParts+nParts, CompareTrackParts );
725 :
726 0 : AliHLTTPCCASliceOutCluster tmp[kMaxClusters];
727 : int currCluster = nOutTrackClusters;
728 : int nHits = 0;
729 0 : for( int ipart=0; ipart<nParts; ipart++ ){
730 0 : const AliHLTTPCGMSliceTrack *t = trackParts[ipart];
731 0 : int nTrackHits = t->NClusters();
732 0 : if( nHits + nTrackHits >= kMaxClusters ) break;
733 0 : const AliHLTTPCCASliceOutCluster *c= t->OrigTrack()->Clusters();
734 0 : AliHLTTPCCASliceOutCluster *c2 = tmp+nHits + nTrackHits-1;
735 0 : for( int i=0; i<nTrackHits; i++, c++, c2-- ) *c2 = *c;
736 0 : float alpha = t->Alpha();
737 0 : for( int i=0; i<nTrackHits; i++) fClusterAngle[currCluster++] = alpha;
738 : nHits+=nTrackHits;
739 0 : }
740 :
741 0 : if ( nHits < 30 ) continue;
742 :
743 0 : float *clX = fClusterX + nOutTrackClusters;
744 0 : clX[0] = tmp[0].GetX();
745 : int ordered = 1;
746 0 : for( int i=1; i<nHits; i++ )
747 : {
748 0 : clX[i] = tmp[i].GetX();
749 0 : if (clX[i] > clX[i - 1]) ordered = 0;
750 : }
751 :
752 0 : if (ordered == 0)
753 : {
754 0 : int2 tmpFilter[kMaxClusters];
755 0 : for( int i = 0;i < nHits;i++)
756 : {
757 0 : tmpFilter[i].x = i;
758 0 : tmpFilter[i].y = tmp[i].GetId();
759 : }
760 0 : qsort(tmpFilter, nHits, sizeof(int2), CompareClusterIds);
761 0 : bool tmpHitsOk[kMaxClusters];
762 0 : tmpHitsOk[tmpFilter[0].x] = true;
763 0 : for (int i = 1;i < nHits;i++)
764 : {
765 0 : tmpHitsOk[tmpFilter[i].x] = (tmpFilter[i].y != tmpFilter[i - 1].y);
766 : }
767 : int nFilteredHits = 0;
768 0 : float *clA = fClusterAngle + nOutTrackClusters;
769 0 : for (int i = 0;i < nHits;i++)
770 : {
771 0 : if (tmpHitsOk[i])
772 : {
773 0 : if (nFilteredHits < i)
774 : {
775 0 : tmp[nFilteredHits] = tmp[i];
776 0 : clA[nFilteredHits] = clA[i];
777 0 : }
778 0 : nFilteredHits++;
779 0 : }
780 : }
781 :
782 : nHits = nFilteredHits;
783 0 : for( int i=0; i<nHits; i++ ) clX[i] = tmp[i].GetX();
784 0 : }
785 :
786 0 : UInt_t *clId = fOutputClusterIds + nOutTrackClusters;
787 0 : for( int i=0; i<nHits; i++ ) clId[i] = tmp[i].GetId();
788 :
789 0 : UInt_t *clT = fClusterRowType + nOutTrackClusters;
790 0 : for( int i=0; i<nHits; i++ ) clT[i] = tmp[i].GetRowType();
791 :
792 0 : float *clY = fClusterY + nOutTrackClusters;
793 0 : for( int i=0; i<nHits; i++ ) clY[i] = tmp[i].GetY();
794 :
795 0 : float *clZ = fClusterZ + nOutTrackClusters;
796 0 : for( int i=0; i<nHits; i++ ) clZ[i] = tmp[i].GetZ();
797 :
798 0 : AliHLTTPCGMMergedTrack &mergedTrack = fOutputTracks[fNOutputTracks];
799 0 : mergedTrack.SetOK(1);
800 0 : mergedTrack.SetNClusters( nHits );
801 0 : mergedTrack.SetFirstClusterRef( nOutTrackClusters );
802 0 : AliHLTTPCGMTrackParam &p1 = mergedTrack.Param();
803 0 : const AliHLTTPCGMSliceTrack &p2 = *(trackParts[0]);
804 :
805 0 : p1.X() = p2.X();
806 0 : p1.Y() = p2.Y();
807 0 : p1.Z() = p2.Z();
808 0 : p1.SinPhi() = p2.SinPhi();
809 0 : p1.DzDs() = p2.DzDs();
810 0 : p1.QPt() = p2.QPt();
811 0 : mergedTrack.SetAlpha( p2.Alpha() );
812 :
813 0 : fNOutputTracks++;
814 0 : nOutTrackClusters += nHits;
815 0 : }
816 : }
817 0 : fNOutputTrackClusters = nOutTrackClusters;
818 0 : }
819 :
820 : void AliHLTTPCGMMerger::Refit()
821 : {
822 : //* final refit
823 : #ifdef HLTCA_GPU_MERGER
824 : if (fGPUTracker && fGPUTracker->IsInitialized())
825 : {
826 : fGPUTracker->RefitMergedTracks(this);
827 : }
828 : else
829 : #endif
830 : {
831 : #ifdef HLTCA_STANDALONE
832 : #pragma omp parallel for
833 : #endif
834 0 : for ( int itr = 0; itr < fNOutputTracks; itr++ ) {
835 :
836 0 : AliHLTTPCGMMergedTrack &track = fOutputTracks[itr];
837 0 : if( !track.OK() ) continue;
838 :
839 0 : int nTrackHits = track.NClusters();
840 :
841 0 : AliHLTTPCGMTrackParam t = track.Param();
842 0 : float Alpha = track.Alpha();
843 :
844 0 : t.Fit( fPolinomialFieldBz,
845 0 : fClusterX+track.FirstClusterRef(),
846 0 : fClusterY+track.FirstClusterRef(),
847 0 : fClusterZ+track.FirstClusterRef(),
848 0 : fClusterRowType+track.FirstClusterRef(),
849 0 : fClusterAngle+track.FirstClusterRef(),
850 0 : fSliceParam, nTrackHits, Alpha, 0 );
851 :
852 0 : if ( fabs( t.QPt() ) < 1.e-4 ) t.QPt() = 1.e-4 ;
853 :
854 0 : bool ok = nTrackHits >= 30 && t.CheckNumericalQuality() && fabs( t.SinPhi() ) <= .999;
855 :
856 0 : track.SetOK(ok);
857 0 : if (!ok) continue;
858 :
859 : if( 1 ){//SG!!!
860 0 : track.SetNClusters( nTrackHits );
861 0 : track.Param() = t;
862 0 : track.Alpha() = Alpha;
863 : }
864 :
865 : {
866 0 : int ind = track.FirstClusterRef();
867 0 : float alpha = fClusterAngle[ind];
868 0 : float x = fClusterX[ind];
869 0 : float y = fClusterY[ind];
870 0 : float z = fClusterZ[ind];
871 0 : float sinA = AliHLTTPCCAMath::Sin( alpha - track.Alpha());
872 0 : float cosA = AliHLTTPCCAMath::Cos( alpha - track.Alpha());
873 0 : track.SetLastX( x*cosA - y*sinA );
874 0 : track.SetLastY( x*sinA + y*cosA );
875 0 : track.SetLastZ( z );
876 : }
877 0 : }
878 : }
879 0 : }
880 :
|
