StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
AliHLTTPCCAMerger.cxx
1 // $Id: AliHLTTPCCAMerger.cxx,v 1.3 2016/06/21 03:39:54 smirnovd Exp $
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 // Developed by: Igor Kulakov <I.Kulakov@gsi.de> *
11 // Maksym Zyzak <M.Zyzak@gsi.de> *
12 // *
13 // Permission to use, copy, modify and distribute this software and its *
14 // documentation strictly for non-commercial purposes is hereby granted *
15 // without fee, provided that the above copyright notice appears in all *
16 // copies and that both the copyright notice and this permission notice *
17 // appear in the supporting documentation. The authors make no claims *
18 // about the suitability of this software for any purpose. It is *
19 // provided "as is" without express or implied warranty. *
20 // *
21 //***************************************************************************
22 
23 
24 #include "AliHLTTPCCASliceTrack.h"
25 #include "AliHLTTPCCATracker.h"
26 #include "AliHLTTPCCAGBTrack.h"
27 #include "AliHLTTPCCATrackParam.h"
28 
29 #include "AliHLTTPCCAMerger.h"
30 
31 #include "AliHLTTPCCAMath.h"
32 #include "Stopwatch.h"
33 
34 #include "AliHLTTPCCATrackParam.h"
35 #include "AliHLTTPCCASliceTrack.h"
36 #include "AliHLTTPCCASliceOutput.h"
37 #include "AliHLTTPCCAMergedTrack.h"
38 #include "AliHLTTPCCAMergerOutput.h"
39 #include "AliHLTTPCCADataCompressor.h"
40 #include "AliHLTTPCCAParam.h"
41 #include "AliHLTTPCCATrackLinearisation.h"
42 
43 #include <iostream>
44 using std::cout;
45 using std::endl;
46 
47 #ifdef MAIN_DRAW
48 #include "AliHLTTPCCADisplay.h"
49 #endif // MAIN_DRAW
50 
51 #ifdef DO_TPCCATRACKER_EFF_PERFORMANCE
52 //#define DO_MERGER_PERF // TODO don't work for Standalone version
53 #else
54 // #include "AliHLTTPCCAMergerPerformance.h"
55 // #include "AliHLTTPCCAPerformance.h"
56 #endif //DO_TPCCATRACKER_EFF_PERFORMANCE
57 
58 #include "AliHLTTPCCATrackParamVector.h"
59 #include "AliHLTTPCCATrackLinearisationVector.h"
60 
61 /*struct AliHLTTPCCAMerger::AliHLTTPCCATrackMemory
62 {
63  short_v fHitIndex; // index of the current hit
64  ushort_v fNHits; // n track hits
65  TrackParamVector fInnerParam;
66  TrackParamVector fOuterParam;
67 };
68 
69 struct AliHLTTPCCAMerger::AliHLTTPCCAHitMemory
70 {
71  sfloat_v *XData; // packed x coordinate of the given (global) hit index
72  sfloat_v *YData; // packed y coordinate of the given (global) hit index
73  sfloat_v *ZData; // packed z coordinate of the given (global) hit index
74 
75  ushort_v *ISector; // sector number
76  ushort_v *IRow; // row number
77 
78  ushort_v *IClu;
79 
80  unsigned short FirstSectorHit[fgkNSlices];
81 };*/
82 
83 int AliHLTTPCCAMerger::fgDoNotMergeBorders = 0;
84 
86 {
87  public:
88  const AliHLTTPCCATrackParam &InnerParam() const { return fInnerParam; }
89  const AliHLTTPCCATrackParam &OuterParam() const { return fOuterParam; }
90  float InnerAlpha() const { return fInnerAlpha; }
91  float OuterAlpha() const { return fOuterAlpha; }
92  unsigned short NClusters() const { return fNClusters; }
93  int FirstClusterRef() const { return fFirstClusterRef; }
94  int PrevNeighbour() const { return fPrevNeighbour; }
95  int NextNeighbour() const { return fNextNeighbour; }
96  unsigned char SlicePrevNeighbour() const { return fSlicePrevNeighbour; }
97  unsigned char SliceNextNeighbour() const { return fSliceNextNeighbour; }
98  bool Used() const { return fUsed; }
99 
100  void SetInnerParam( const AliHLTTPCCATrackParam &v ) { fInnerParam = v; }
101  void SetOuterParam( const AliHLTTPCCATrackParam &v ) { fOuterParam = v; }
102  void SetInnerAlpha( float v ) { fInnerAlpha = v; }
103  void SetOuterAlpha( float v ) { fOuterAlpha = v; }
104  void SetNClusters ( unsigned short v ) { fNClusters = v; }
105  void SetFirstClusterRef( int v ) { fFirstClusterRef = v; }
106  void SetPrevNeighbour( int v ) { fPrevNeighbour = v; }
107  void SetNextNeighbour( int v ) { fNextNeighbour = v; }
108  void SetSlicePrevNeighbour( unsigned char v ) { fSlicePrevNeighbour = v; }
109  void SetSliceNextNeighbour( unsigned char v ) { fSliceNextNeighbour = v; }
110  void SetUsed( bool v ) { fUsed = v; }
111 #ifdef DO_TPCCATRACKER_EFF_PERFORMANCE
112  int orig_track_id;
113  unsigned char fSlice;
114  int number;
115 #endif // DO_TPCCATRACKER_EFF_PERFORMANCE
116 
117  public:
118 
119  float ChiPrev; //characteristic of the link quality to the inner neighbour (for overlaped tracks it is chi2, for not overlaped - distance between tracks)
120  float ChiNext; //characteristic of the link quality to the outer neighbour
121  unsigned char fInnerRow; // number of the inner row of the track
122  unsigned char fOuterRow; // number of the outer row of the track
123 
124  private:
125 
126  AliHLTTPCCATrackParam fInnerParam; // Parameters of the track at the inner point
127  AliHLTTPCCATrackParam fOuterParam; // Parameters of the track at the outer point
128  float fInnerAlpha; // The angle of the sector, where inner parameters are
129  float fOuterAlpha; // The angle of the sector, where outers parameters are
130  unsigned short fNClusters; //Clusters number of the track
131  int fFirstClusterRef; // index of the first track cluster in the global cluster array (AliHLTTPCCAMerger::fClusterInfos)
132  int fPrevNeighbour; // The number of the inner (previous) neighbour in the tracks array (AliHLTTPCCAMerger::fTrackInfos)
133  int fNextNeighbour; // The number of the outer (previous) neighbour in the tracks array (AliHLTTPCCAMerger::fTrackInfos)
134  unsigned char fSlicePrevNeighbour; // The number of the sector, which contains inner neighbour
135  unsigned char fSliceNextNeighbour; // The number of the sector, which contains outer neighbour
136  bool fUsed; // is the slice track already merged (=1 -> merged, 0 -> not merged)
137 };
138 
139 
140 AliHLTTPCCAMerger::AliHLTTPCCAMerger()
141  :
142  fSliceParam(),
143 
144  fMaxClusterInfos( 0 ),
145  fClusterInfos( 0 ),
146 
147  fMaxTrackInfos( 0 ),
148  fTrackInfos( 0 ),
149  fOutput( 0 )
150 
151 {
152  //* constructor
153  Clear();
154 
155  for(int iSlice=0; iSlice<fgkNSlices; iSlice++)
156  slices[iSlice] = 0;
157 }
158 
159 /*
160 AliHLTTPCCAMerger::AliHLTTPCCAMerger(const AliHLTTPCCAMerger&)
161  :
162  fSliceParam(),
163  fkSlices(0),
164  fOutput(0),
165  fTrackInfos(0),
166  fMaxTrackInfos(0),
167  fClusterInfos(0),
168  fMaxClusterInfos(0)
169 {
170 }
171 
172 const AliHLTTPCCAMerger &AliHLTTPCCAMerger::operator=(const AliHLTTPCCAMerger&) const
173 {
174  return *this;
175 }
176 */
177 
178 AliHLTTPCCAMerger::~AliHLTTPCCAMerger()
179 {
180  //* destructor
181  if ( fTrackInfos ) delete[] fTrackInfos;
182  if ( fClusterInfos ) delete[] fClusterInfos;
183  if ( fOutput ) delete[] ( ( char* )( fOutput ) );
184 }
185 
186 void AliHLTTPCCAMerger::Clear()
187 {
188  for ( int i = 0; i < fgkNSlices; ++i ) {
189  fkSlices[i] = 0;
190  }
191  for ( int i = 0; i < fNTimers; i++) fTimers[i] = 0;
192 }
193 
194 void AliHLTTPCCAMerger::SetSlices (int i, AliHLTTPCCATracker *sl )
195 {
196  //copy sector parameters information
197  slices[i] = sl;
198 }
199 
200 void AliHLTTPCCAMerger::SetSliceData( int index, const AliHLTTPCCASliceOutput *sliceData )
201 {
202  //copy sector output information: hits, reconstructed tracks...
203  fkSlices[index] = sliceData;
204 }
205 
206 void AliHLTTPCCAMerger::Reconstruct()
207 {
208  //* Main merging routine. Consist of 3 steps:
209 #ifdef USE_TIMERS
210  Stopwatch timer;
211  timer.Start();
212 #endif // USE_TIMERS
213 
214 // 1) copy information from the sector tracker
215  UnpackSlices();
216 #ifdef USE_TIMERS
217  timer.Stop();
218  fTimers[0] = timer.RealTime();
219 #endif // USE_TIMERS
220 
221 #ifdef DO_MERGER_PERF
222  AliHLTTPCCAPerformance::Instance().CreateHistos("Merger");
223  ((AliHLTTPCCAMergerPerformance*)(AliHLTTPCCAPerformance::Instance().GetSubPerformance("Merger")))->SetNewEvent(
224  AliHLTTPCCAPerformance::Instance().GetTracker(),
225  AliHLTTPCCAPerformance::Instance().GetHitLabels(),
226  AliHLTTPCCAPerformance::Instance().GetMCTracks(),
227  AliHLTTPCCAPerformance::Instance().GetMCPoints());
228  ((AliHLTTPCCAMergerPerformance*)(AliHLTTPCCAPerformance::Instance().GetSubPerformance("Merger")))->FillMC();
229 #endif //DO_MERGER_PERF
230 
231 #ifdef USE_TIMERS
232  timer.Start();
233 #endif // USE_TIMERS
234 
235 // 2) merge nonoverlaping tracks
236  Merging(1);
237 
238 #ifdef USE_TIMERS
239  timer.Stop();
240  fTimers[1] = timer.RealTime();
241 
242 // 3) merge overlaping tracks, store the tracks to the global tracker
243  timer.Start();
244 #endif // USE_TIMERS
245 
246  Merging(0);
247 
248 #ifdef USE_TIMERS
249  timer.Stop();
250  fTimers[2] = timer.RealTime();
251 #endif // USE_TIMERS
252 }
253 
254 void AliHLTTPCCAMerger::UnpackSlices()
255 {
256  //* unpack the cluster information from the slice tracks and initialize track info array
257 
258  // get N tracks and N clusters in event
259 #ifdef DO_TPCCATRACKER_EFF_PERFORMANCE
260  int NTracksPrev=0;
261 #endif // DO_TPCCATRACKER_EFF_PERFORMANCE
262  int nTracksTotal = 0;
263  int nTrackClustersTotal = 0;
264  for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
265  if ( !fkSlices[iSlice] ) continue;
266  nTracksTotal += fkSlices[iSlice]->NTracks();
267  nTrackClustersTotal += fkSlices[iSlice]->NTrackClusters();
268 #ifdef DO_TPCCATRACKER_EFF_PERFORMANCE
269  slices[iSlice]->fNOutTracks1 = 0;
270 #endif //DO_TPCCATRACKER_EFF_PERFORMANCE
271  }
272 
273  // book/clean memory if necessary
274  {
275  // if ( nTracksTotal > fMaxTrackInfos || ( fMaxTrackInfos > 100 && nTracksTotal < 0.5*fMaxTrackInfos ) )
276  {
277  if ( fTrackInfos ) delete[] fTrackInfos;
278  fMaxTrackInfos = ( int ) ( nTracksTotal );
279  fTrackInfos = new AliHLTTPCCASliceTrackInfo[fMaxTrackInfos];
280  }
281 
282 // if ( nTrackClustersTotal > fMaxClusterInfos || ( fMaxClusterInfos > 1000 && nTrackClustersTotal < 0.5*fMaxClusterInfos ) )
283  {
284  if ( fClusterInfos ) delete[] fClusterInfos;
285  fMaxClusterInfos = ( int ) ( nTrackClustersTotal );
286  fClusterInfos = new AliHLTTPCCAClusterInfo [fMaxClusterInfos];
287  }
288 
289  if ( fOutput ) delete[] ( ( char* )( fOutput ) );
290  int size = AliHLTTPCCAMergerOutput::EstimateSize( nTracksTotal, nTrackClustersTotal );
291  fOutput = ( AliHLTTPCCAMergerOutput* )( new float2[size/sizeof( float2 )+1] );
292  }
293 
295 /* if ( fHitMemory.XData ) delete[] fHitMemory.XData;
296  fHitMemory.XData = new float[nTracksTotal];
297  if ( fHitMemory.YData ) delete[] fHitMemory.YData;
298  fHitMemory.YData = new float[nTracksTotal];
299  if ( fHitMemory.ZData ) delete[] fHitMemory.ZData;
300  fHitMemory.ZData = new float[nTracksTotal];
301 
302  if ( fHitMemory.ISector ) delete[] fHitMemory.ISector;
303  fHitMemory.ISector = new unsigned char[nTracksTotal];
304  if ( fHitMemory.IRow ) delete[] fHitMemory.IRow;
305  fHitMemory.IRow = new unsigned char[nTracksTotal];
306  if ( fHitMemory.IClu ) delete[] fHitMemory.IClu;
307  fHitMemory.IClu = new short[nTracksTotal];
308 
309  for ( int iSec = 0; iSec < fgkNSlices; iSec++ ) {
310  if(iSec == 0) fHitMemory.FirstSectorHit[iSec] = 0;
311  if(iSec > 0) fHitMemory.FirstSectorHit[iSec] = fkSlices[iSec-1]->NTrackClusters();
312  const AliHLTTPCCASliceOutput &slice = *( fkSlices[iSlice] );
313 
314  for(int iHit=0; iHit < fkSlices[iSec-1]->NTrackClusters(); iHit++)
315  {
316  fHitMemory.ISector = iSec;
317  fHitMemory.IRow = slice.ClusterIDrc( iHit ).Row();
318  fHitMemory.IClu = slice.ClusterIDrc( iHit ).Cluster();
319  const float2 &yz = slice.ClusterUnpackedYZ(iHit);
320  fHitMemory.XData = slice.ClusterUnpackedX(iHit);
321  fHitMemory.YData = yz.x;
322  fHitMemory.ZData = yz.y;
323  }
324  }*/
326 
327  // unpack track and cluster information
328 
329  int nTracksCurrent = 0;
330  int nClustersCurrent = 0;
331  // TODO !t0.NDF()
332  for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
333 
334  fSliceTrackInfoStart[ iSlice ] = nTracksCurrent;
335  fSliceNTrackInfos[ iSlice ] = 0;
336 
337  if ( !fkSlices[iSlice] ) continue;
338 
339  const AliHLTTPCCASliceOutput &slice = *( fkSlices[iSlice] );
340 
341  for ( int itr = 0; itr < slice.NTracks(); itr += ushort_v::Size ) {
342 
343  int nTracksVector = ushort_v::Size;
344  if(slice.NTracks() - itr < ushort_v::Size )
345  nTracksVector = slice.NTracks() - itr;
346 
347  sfloat_v::Memory startAlpha;
348  sfloat_v::Memory endAlpha;
349 
350  int hits[1000][ushort_v::Size];
351  ushort_v::Memory nHits;
352  AliHLTTPCCATrackParam startPoint[ushort_v::Size];
353  AliHLTTPCCATrackParam endPoint[ushort_v::Size];
354  unsigned short nCluNew = 0;
355 
356  for(int iV=0; iV < nTracksVector; iV++)
357  {
358  const AliHLTTPCCASliceTrack &sTrack = slice.Track( itr + iV );
359 
360  nHits[iV] = 0;
361  for ( int iTrClu = 0; iTrClu < sTrack.NClusters(); iTrClu++ ) {
362 // unpack cluster information
363  AliHLTTPCCAClusterInfo &clu = fClusterInfos[nClustersCurrent + nCluNew + nHits[iV]];
364  int ic = sTrack.FirstClusterRef() + iTrClu;
365 
366  clu.SetISlice( iSlice );
367  clu.SetIRow( slice.ClusterIDrc( ic ).Row() );
368  clu.SetIClu( slice.ClusterIDrc( ic ).Cluster() );
369  float2 yz = slice.ClusterUnpackedYZ( ic );
370  clu.SetX( slice.ClusterUnpackedX( ic ) );
371  clu.SetY( yz.x );
372  clu.SetZ( yz.y );
373  hits[iTrClu][iV] = nClustersCurrent + nCluNew + iTrClu;
374  nHits[iV]++;
375  }
376  nCluNew += nHits[iV];
377 
378  startPoint[iV] = sTrack.Param();
379  endPoint[iV] = startPoint[iV];
380  startAlpha[iV] = slices[iSlice]->Param().Alpha();
381  endAlpha[iV] = startAlpha[iV];
382  }
383 
384 //when we turn off the extrapolation step in the tracklet constructor, we have parameters in the last point, not in the first!
385 //that's why the fitting direction should be changed
386  sfloat_m fitted = sfloat_m(true);
387  fitted &= static_cast<sfloat_m>(static_cast<ushort_v>(nHits) >= 3);
388  //? assert( fitted );
389 //#define DISABLE_HIT_SEARCH
390 #ifdef DISABLE_HIT_SEARCH
391 //TODO correct code
392  fitted &= FitTrack( startPoint, startAlpha, endPoint, endAlpha, hits, nHits, nTracksVector, 1 );
393  for(int iV=0; iV < nTracksVector; iV++)
394  {
395  endPoint[iV] = startPoint[iV];
396  endAlpha[iV] = startAlpha[iV];
397  }
398  fitted &= FitTrack( endPoint, endAlpha, startPoint, startAlpha, hits, nHits, nTracksVector, 0 );
399 #else
400 // refit the track
401  AliHLTTPCCATrackParamVector vStartPoint;
402  AliHLTTPCCATrackParamVector vEndPoint;
403  ConvertTrackParamToVector(startPoint,vStartPoint,nTracksVector);
404  sfloat_v vStartAlpha(Vc::Zero);
405  vStartAlpha.load(startAlpha);
406  sfloat_v vEndAlpha(Vc::Zero);
407 
408  fitted &= static_cast<sfloat_m>( ushort_v( Vc::IndexesFromZero ) < nTracksVector );
409 
410  ushort_v firstHits(Vc::Zero);
411 
412  vEndPoint = vStartPoint;
413  vEndAlpha = vStartAlpha;
414 //refit in the forward direction: going from the first hit to the last, mask "fitted" marks with 0 tracks, which are not fitted correctly
415  fitted &= FitTrack( vEndPoint, vEndAlpha, hits, firstHits, nHits, nTracksVector,fitted, 0 );
416 // if chi2 per degree of freedom > 3. sigma - mark track with 0
417  fitted &= vEndPoint.Chi2()/static_cast<sfloat_v>(vEndPoint.NDF()) < 9.f;
418 
419  vStartPoint = vEndPoint;
420  vStartAlpha = vEndAlpha;
421 //refit in the backward direction: going from the last hit to the first
422  fitted &= FitTrack( vStartPoint, vStartAlpha, hits, firstHits, nHits, nTracksVector,fitted, 1 );
423 // if chi2 per degree of freedom > 3. sigma - mark track with 0
424  fitted &= vStartPoint.Chi2()/static_cast<sfloat_v>(vStartPoint.NDF()) < 9.f;
425 
426  for(int iV=0; iV < nTracksVector; iV++)
427  {
428  if ( !fitted[iV] ) continue;
429  startPoint[iV] = AliHLTTPCCATrackParam( vStartPoint, iV );
430  endPoint[iV] = AliHLTTPCCATrackParam( vEndPoint, iV );
431  }
432  vStartAlpha.store(&(startAlpha[0])); // work around .store(startAlpha)
433  vEndAlpha.store(&(endAlpha[0]));
434 #endif
435  for(int iV=0; iV < nTracksVector; iV++)
436  {
437  const AliHLTTPCCASliceTrack &sTrack = slice.Track( itr + iV );
438  if ( !(fitted[iV]) ) continue;
439  if ( nHits[iV] < AliHLTTPCCAParameters::MinTrackPurity*sTrack.NClusters() ) continue;
440 // if the track fitted correctly store the track
441  AliHLTTPCCASliceTrackInfo &track = fTrackInfos[nTracksCurrent];
442 
443  track.SetInnerParam( startPoint[iV] );
444  track.SetInnerAlpha( startAlpha[iV] );
445  track.SetOuterParam( endPoint[iV] );
446  track.SetOuterAlpha( endAlpha[iV] );
447  track.SetFirstClusterRef( nClustersCurrent );
448  track.SetNClusters( nHits[iV] );
449  track.SetPrevNeighbour( -1 );
450  track.SetNextNeighbour( -1 );
451  track.SetSlicePrevNeighbour( -1 );
452  track.SetSliceNextNeighbour( -1 );
453  track.ChiPrev = 10000000;
454  track.ChiNext = 10000000;
455  track.SetUsed( 0 );
457 #ifdef DO_TPCCATRACKER_EFF_PERFORMANCE
458  track.orig_track_id = itr+iV;
459  track.fSlice = iSlice;
460  track.number = nTracksCurrent - fSliceTrackInfoStart[iSlice];
461 #endif // DO_TPCCATRACKER_EFF_PERFORMANCE
462  track.fInnerRow = (fClusterInfos[hits[0][iV]]).IRow();
463  track.fOuterRow = (fClusterInfos[hits[nHits[iV]-1][iV]]).IRow();
465  for ( int i = 0; i < nHits[iV]; i++ )
466  fClusterInfos[nClustersCurrent + i] = fClusterInfos[hits[i][iV]];
467  nTracksCurrent++;
468  fSliceNTrackInfos[ iSlice ]++;
469  nClustersCurrent += nHits[iV];
470  }
471  }
472 
473 #ifdef DO_TPCCATRACKER_EFF_PERFORMANCE
474  if (slices[iSlice]->fOutTracks1) delete[] slices[iSlice]->fOutTracks1;
475  slices[iSlice]->fOutTracks1 = new AliHLTTPCCAOutTrack [nTracksCurrent-NTracksPrev];
476  for (int i=0; i<nTracksCurrent-NTracksPrev; i++)
477  {
478  slices[iSlice]->fOutTracks1[i].SetStartPoint(fTrackInfos[i+NTracksPrev].InnerParam());
479  slices[iSlice]->fOutTracks1[i].SetEndPoint(fTrackInfos[i+NTracksPrev].OuterParam());
480  slices[iSlice]->fOutTracks1[i].SetOrigTrackID(fTrackInfos[i+NTracksPrev].orig_track_id);
481  slices[iSlice]->fNOutTracks1++;
482  }
483  NTracksPrev = nTracksCurrent;
484 #endif // DO_TPCCATRACKER_EFF_PERFORMANCE
485  //std::cout<<"Unpack slice "<<iSlice<<": ntracks "<<slice.NTracks()<<"/"<<fSliceNTrackInfos[iSlice]<<std::endl;
486  }
487 }
488 
489 sfloat_m AliHLTTPCCAMerger::FitTrack( AliHLTTPCCATrackParamVector &t, sfloat_v &Alpha0V,
490  int hits[2000][ushort_v::Size], ushort_v &firstHits,ushort_v::Memory &NTrackHits,
491  int &nTracksV, sfloat_m active0, bool dir )
492 {
493  // Fit the track
494 /*#ifdef MAIN_DRAW
495  AliHLTTPCCADisplay::Instance().ClearView();
496  AliHLTTPCCADisplay::Instance().SetTPCView();
497  AliHLTTPCCADisplay::Instance().DrawTPC();
498 #endif*/
500 
502 
503  bool first = 1;
504 
505  t.CalculateFitParameters( fitPar );
506 
507  int hitsNew[AliHLTTPCCAParameters::MaxNumberOfRows8*4][ushort_v::Size]; // 4 - reserve for several turn
508  ushort_v::Memory nHitsNew;
509  for(int iV=0; iV<ushort_v::Size; iV++)
510  nHitsNew[iV] = 0;
511 
512  ushort_v nHits(NTrackHits);
513  nHits.setZero(static_cast<ushort_m>(!active0));
514 
515  int nHitsMax = nHits.max();
516 
517  sfloat_m active = active0; // stop on a hit if this hit can't be added
518  for ( unsigned short ihit = 0; ihit < nHitsMax; ihit++ ) {
519 
520  active &= static_cast<sfloat_m>( ihit < nHits );
521  if(active.isEmpty()) continue;
522  sfloat_v::Memory xH, yH, zH, sliceAlpha;
523  ushort_v::Memory Row;
524 
525  for(int iV=0; iV < nTracksV; iV++)
526  {
527  if( !(active[iV]) ) continue;
528  if( ihit > NTrackHits[iV] - 1) continue;
529  const int jhit = dir ? ( NTrackHits[iV] - 1 - ihit ) : ihit;
530  AliHLTTPCCAClusterInfo &h = fClusterInfos[hits[static_cast <unsigned short>(firstHits[iV]) + jhit][iV]];
531  sliceAlpha[iV] = slices[h.ISlice()]->Param().Alpha();
532  xH[iV] = h.X();
533  yH[iV] = h.Y();
534  zH[iV] = h.Z();
535  Row[iV] = h.IRow();
536  }
537 
538  const sfloat_v xV(xH);
539  const sfloat_v yV(yH);
540  const sfloat_v zV(zH);
541  const sfloat_v sliceAlphaV(sliceAlpha);
542  const ushort_v RowV(Row);
543 
544  const sfloat_m savedActive = active;
545  const sfloat_v rotateA = sliceAlphaV - Alpha0V;
546  const sfloat_m &rotated = t.Rotate( rotateA, l, .999f, active);
547  active &= rotated;
548  const sfloat_v xLast = t.X();
549 
550  Alpha0V(active) = sliceAlphaV;
551 
552  const sfloat_m &transported = t.TransportToXWithMaterial( xV, l, fitPar, fSliceParam.cBz( ), 0.999f, active);
553  active &= transported;
554 
555  if ( first ) {
556  t.SetCov( 0, 10.f, active );
557  t.SetCov( 1, 0.f, active );
558  t.SetCov( 2, 10.f, active );
559  t.SetCov( 3, 0.f, active );
560  t.SetCov( 4, 0.f, active );
561  t.SetCov( 5, 1.f, active );
562  t.SetCov( 6, 0.f, active );
563  t.SetCov( 7, 0.f, active );
564  t.SetCov( 8, 0.f, active );
565  t.SetCov( 9, 1.f, active );
566  t.SetCov( 10, 0.f, active );
567  t.SetCov( 11, 0.f, active );
568  t.SetCov( 12, 0.f, active );
569  t.SetCov( 13, 0.f, active );
570  t.SetCov( 14, 10.f, active );
571  t.SetChi2( 0.f, active);
572  t.SetNDF( short_v(-5), static_cast<short_m>(active) );
573  t.CalculateFitParameters( fitPar );
574  }
575 
576  sfloat_v err2Y, err2Z;
577 
578  fSliceParam.GetClusterErrors2( RowV, t, &err2Y, &err2Z );
579  const sfloat_m &filtered = t.FilterWithMaterial(yV, zV, err2Y, err2Z, 0.999f,active);
580 
581  const sfloat_m broken = savedActive && (!rotated || !transported || !filtered);
582  if ( !broken.isEmpty() ) {
583  t.TransportToXWithMaterial( xLast, l, fitPar, fSliceParam.cBz( ), 0.999f, transported && !filtered ); // transport back if hit can't be added. TODO with out material
584  t.Rotate( -rotateA, l, .999f, rotated && !transported );
585  }
586 
587  active &= filtered;
588  if(active.isEmpty()) continue;
589 
590  first = 0;
591  for(int iV=0; iV < nTracksV; iV++)
592  {
593  if( !(active[iV]) ) continue;
594  const int jhit = dir ? ( NTrackHits[iV] - 1 - ihit ) : ihit;
595  hitsNew[nHitsNew[iV]][iV] = hits[static_cast <unsigned short>(firstHits[iV]) + jhit][iV];
596  nHitsNew[iV]++;
597  }
598  }
599  t.SetQPt( sfloat_v(1.e-8f), CAMath::Abs( t.QPt() ) < 1.e-8f );
600 
601  sfloat_m ok = active0 && ushort_v(nHitsNew) >= 3;
602 
603 // if track has infinite partameters or covariances, or negative diagonal elements of Cov. matrix, mark it as fitted uncorrectly
604  const sfloat_v *c = t.Cov();
605  for ( unsigned char i = 0; i < 15; i++ ) ok &= CAMath::Finite( c[i] );
606  for ( unsigned char i = 0; i < 5; i++ ) ok &= CAMath::Finite( t.Par()[i] );
608  ok &= (c[0] > Vc::Zero) && (c[2] > Vc::Zero) && (c[5] > Vc::Zero) && (c[9] > Vc::Zero) && (c[14] > Vc::Zero);
609 // ok &= (c[0] < 5.) && (c[2] < 5.) && (c[5] < 2.) && (c[9] < 2.) && (c[14] < 2.);
610  ok &= CAMath::Abs( t.SinPhi() ) < .999f;
611 
612  t.SetSignCosPhi( sfloat_v(1.f), ok && static_cast<sfloat_m>(l.CosPhi() >= Vc::Zero) );
613  t.SetSignCosPhi( sfloat_v(-1.f), ok && static_cast<sfloat_m>(l.CosPhi() < Vc::Zero) );
614 
615  for(int iV=0; iV < nTracksV; iV++)
616  {
617  if ( !ok[iV] ) continue;
618  NTrackHits[iV] = nHitsNew[iV];
619  for ( unsigned char i = 0; i < NTrackHits[iV]; i++ ) {
620  hits[static_cast <unsigned short>(firstHits[iV]) + (dir ?( NTrackHits[iV]-1-i ) :i)][iV] = hitsNew[i][iV];
621  }
622  }
623 
624  return ok;
625 }
626 
627 void AliHLTTPCCAMerger::InvertCholetsky(float a[15])
628 {
629  float d[5], uud, u[5][5];
630  for(int i=0; i<5; i++)
631  {
632  d[i]=0.f;
633  for(int j=0; j<5; j++)
634  u[i][j]=0.;
635  }
636 
637  for(int i=0; i<5; i++)
638  {
639  uud=0.;
640  for(int j=0; j<i; j++)
641  uud += u[j][i]*u[j][i]*d[j];
642  uud = a[i*(i+3)/2] - uud;
643 
644  if(fabs(uud)<1.e-12) uud = 1.e-12;
645  d[i] = uud/fabs(uud);
646  u[i][i] = sqrt(fabs(uud));
647 
648  for(int j=i+1; j<5; j++)
649  {
650  uud = 0.;
651  for(int k=0; k<i; k++)
652  uud += u[k][i]*u[k][j]*d[k];
653  uud = a[j*(j+1)/2+i] - uud;
654  u[i][j] = d[i]/u[i][i]*uud;
655  }
656  }
657 
658  float u1[5];
659 
660  for(int i=0; i<5; i++)
661  {
662  u1[i] = u[i][i];
663  u[i][i] = 1.f/u[i][i];
664  }
665  for(int i=0; i<4; i++)
666  {
667  u[i][i+1] = - u[i][i+1]*u[i][i]*u[i+1][i+1];
668  }
669  for(int i=0; i<3; i++)
670  {
671  u[i][i+2] = u[i][i+1]*u1[i+1]*u[i+1][i+2]-u[i][i+2]*u[i][i]*u[i+2][i+2];
672  }
673  for(int i=0; i<2; i++)
674  {
675  u[i][i+3] = u[i][i+2]*u1[i+2]*u[i+2][i+3] - u[i][i+3]*u[i][i]*u[i+3][i+3];
676  u[i][i+3] -= u[i][i+1]*u1[i+1]*(u[i+1][i+2]*u1[i+2]*u[i+2][i+3] - u[i+1][i+3]);
677  }
678  u[0][4] = u[0][2]*u1[2]*u[2][4] - u[0][4]*u[0][0]*u[4][4];
679  u[0][4] += u[0][1]*u1[1]*(u[1][4] - u[1][3]*u1[3]*u[3][4] - u[1][2]*u1[2]*u[2][4]);
680  u[0][4] += u[3][4]*u1[3]*(u[0][3] - u1[2]*u[2][3]*(u[0][2] - u[0][1]*u1[1]*u[1][2]));
681 
682  for(int i=0; i<5; i++)
683  a[i+10] = u[i][4]*d[4]*u[4][4];
684  for(int i=0; i<4; i++)
685  a[i+6] = u[i][3]*u[3][3]*d[3] + u[i][4]*u[3][4]*d[4];
686  for(int i=0; i<3; i++)
687  a[i+3] = u[i][2]*u[2][2]*d[2] + u[i][3]*u[2][3]*d[3] + u[i][4]*u[2][4]*d[4];
688  for(int i=0; i<2; i++)
689  a[i+1] = u[i][1]*u[1][1]*d[1] + u[i][2]*u[1][2]*d[2] + u[i][3]*u[1][3]*d[3] + u[i][4]*u[1][4]*d[4];
690  a[0] = u[0][0]*u[0][0]*d[0] + u[0][1]*u[0][1]*d[1] + u[0][2]*u[0][2]*d[2] + u[0][3]*u[0][3]*d[3] + u[0][4]*u[0][4]*d[4];
691 }
692 
693 void AliHLTTPCCAMerger::MultiplySS(float const C[15], float const V[15], float K[5][5])
694 {
695 //multiply 2 symmetric matricies
696  K[0][0] = C[0]*V[ 0] + C[1]*V[ 1] + C[3]*V[ 3] + C[6]*V[ 6] + C[10]*V[10];
697  K[0][1] = C[0]*V[ 1] + C[1]*V[ 2] + C[3]*V[ 4] + C[6]*V[ 7] + C[10]*V[11];
698  K[0][2] = C[0]*V[ 3] + C[1]*V[ 4] + C[3]*V[ 5] + C[6]*V[ 8] + C[10]*V[12];
699  K[0][3] = C[0]*V[ 6] + C[1]*V[ 7] + C[3]*V[ 8] + C[6]*V[ 9] + C[10]*V[13];
700  K[0][4] = C[0]*V[10] + C[1]*V[11] + C[3]*V[12] + C[6]*V[13] + C[10]*V[14];
701 
702  K[1][0] = C[1]*V[ 0] + C[2]*V[ 1] + C[4]*V[ 3] + C[7]*V[ 6] + C[11]*V[10];
703  K[1][1] = C[1]*V[ 1] + C[2]*V[ 2] + C[4]*V[ 4] + C[7]*V[ 7] + C[11]*V[11];
704  K[1][2] = C[1]*V[ 3] + C[2]*V[ 4] + C[4]*V[ 5] + C[7]*V[ 8] + C[11]*V[12];
705  K[1][3] = C[1]*V[ 6] + C[2]*V[ 7] + C[4]*V[ 8] + C[7]*V[ 9] + C[11]*V[13];
706  K[1][4] = C[1]*V[10] + C[2]*V[11] + C[4]*V[12] + C[7]*V[13] + C[11]*V[14];
707 
708  K[2][0] = C[3]*V[ 0] + C[4]*V[ 1] + C[5]*V[ 3] + C[8]*V[ 6] + C[12]*V[10];
709  K[2][1] = C[3]*V[ 1] + C[4]*V[ 2] + C[5]*V[ 4] + C[8]*V[ 7] + C[12]*V[11];
710  K[2][2] = C[3]*V[ 3] + C[4]*V[ 4] + C[5]*V[ 5] + C[8]*V[ 8] + C[12]*V[12];
711  K[2][3] = C[3]*V[ 6] + C[4]*V[ 7] + C[5]*V[ 8] + C[8]*V[ 9] + C[12]*V[13];
712  K[2][4] = C[3]*V[10] + C[4]*V[11] + C[5]*V[12] + C[8]*V[13] + C[12]*V[14];
713 
714  K[3][0] = C[6]*V[ 0] + C[7]*V[ 1] + C[8]*V[ 3] + C[9]*V[ 6] + C[13]*V[10];
715  K[3][1] = C[6]*V[ 1] + C[7]*V[ 2] + C[8]*V[ 4] + C[9]*V[ 7] + C[13]*V[11];
716  K[3][2] = C[6]*V[ 3] + C[7]*V[ 4] + C[8]*V[ 5] + C[9]*V[ 8] + C[13]*V[12];
717  K[3][3] = C[6]*V[ 6] + C[7]*V[ 7] + C[8]*V[ 8] + C[9]*V[ 9] + C[13]*V[13];
718  K[3][4] = C[6]*V[10] + C[7]*V[11] + C[8]*V[12] + C[9]*V[13] + C[13]*V[14];
719 
720  K[4][0] = C[10]*V[ 0] + C[11]*V[ 1] + C[12]*V[ 3] + C[13]*V[ 6] + C[14]*V[10];
721  K[4][1] = C[10]*V[ 1] + C[11]*V[ 2] + C[12]*V[ 4] + C[13]*V[ 7] + C[14]*V[11];
722  K[4][2] = C[10]*V[ 3] + C[11]*V[ 4] + C[12]*V[ 5] + C[13]*V[ 8] + C[14]*V[12];
723  K[4][3] = C[10]*V[ 6] + C[11]*V[ 7] + C[12]*V[ 8] + C[13]*V[ 9] + C[14]*V[13];
724  K[4][4] = C[10]*V[10] + C[11]*V[11] + C[12]*V[12] + C[13]*V[13] + C[14]*V[14];
725 }
726 
727 void AliHLTTPCCAMerger::MultiplyMS(float const C[5][5], float const V[15], float K[15])
728 {
729 //multiply symmetric and nonsymmetric matricies
730  K[0] = C[0][0]*V[0] + C[0][1]*V[1] + C[0][2]*V[3] + C[0][3]*V[6] + C[0][4]*V[10];
731 
732  K[1] = C[1][0]*V[0] + C[1][1]*V[1] + C[1][2]*V[3] + C[1][3]*V[6] + C[1][4]*V[10];
733  K[2] = C[1][0]*V[1] + C[1][1]*V[2] + C[1][2]*V[4] + C[1][3]*V[7] + C[1][4]*V[11];
734 
735  K[3] = C[2][0]*V[0] + C[2][1]*V[1] + C[2][2]*V[3] + C[2][3]*V[6] + C[2][4]*V[10];
736  K[4] = C[2][0]*V[1] + C[2][1]*V[2] + C[2][2]*V[4] + C[2][3]*V[7] + C[2][4]*V[11];
737  K[5] = C[2][0]*V[3] + C[2][1]*V[4] + C[2][2]*V[5] + C[2][3]*V[8] + C[2][4]*V[12];
738 
739  K[6] = C[3][0]*V[0] + C[3][1]*V[1] + C[3][2]*V[3] + C[3][3]*V[6] + C[3][4]*V[10];
740  K[7] = C[3][0]*V[1] + C[3][1]*V[2] + C[3][2]*V[4] + C[3][3]*V[7] + C[3][4]*V[11];
741  K[8] = C[3][0]*V[3] + C[3][1]*V[4] + C[3][2]*V[5] + C[3][3]*V[8] + C[3][4]*V[12];
742  K[9] = C[3][0]*V[6] + C[3][1]*V[7] + C[3][2]*V[8] + C[3][3]*V[9] + C[3][4]*V[13];
743 
744  K[10] = C[4][0]*V[ 0] + C[4][1]*V[ 1] + C[4][2]*V[ 3] + C[4][3]*V[ 6] + C[4][4]*V[10];
745  K[11] = C[4][0]*V[ 1] + C[4][1]*V[ 2] + C[4][2]*V[ 4] + C[4][3]*V[ 7] + C[4][4]*V[11];
746  K[12] = C[4][0]*V[ 3] + C[4][1]*V[ 4] + C[4][2]*V[ 5] + C[4][3]*V[ 8] + C[4][4]*V[12];
747  K[13] = C[4][0]*V[ 6] + C[4][1]*V[ 7] + C[4][2]*V[ 8] + C[4][3]*V[ 9] + C[4][4]*V[13];
748  K[14] = C[4][0]*V[10] + C[4][1]*V[11] + C[4][2]*V[12] + C[4][3]*V[13] + C[4][4]*V[14];
749 }
750 
751 void AliHLTTPCCAMerger::MultiplySR(float const C[15], float const r_in[5], float r_out[5])
752 {
753 //multiply vector and symmetric matrix
754  r_out[0] = r_in[0]*C[ 0] + r_in[1]*C[ 1] + r_in[2]*C[ 3] +r_in[3]*C[ 6] + r_in[4]*C[10];
755  r_out[1] = r_in[0]*C[ 1] + r_in[1]*C[ 2] + r_in[2]*C[ 4] +r_in[3]*C[ 7] + r_in[4]*C[11];
756  r_out[2] = r_in[0]*C[ 3] + r_in[1]*C[ 4] + r_in[2]*C[ 5] +r_in[3]*C[ 8] + r_in[4]*C[12];
757  r_out[3] = r_in[0]*C[ 6] + r_in[1]*C[ 7] + r_in[2]*C[ 8] +r_in[3]*C[ 9] + r_in[4]*C[13];
758  r_out[4] = r_in[0]*C[10] + r_in[1]*C[11] + r_in[2]*C[12] +r_in[3]*C[13] + r_in[4]*C[14];
759 }
760 
761 void AliHLTTPCCAMerger::FilterTracks(float const r[5], float const C[15], float const m[5], float const V[15], float R[5], float W[15], float &chi2)
762 {
763  float S[15];
764  for(int i=0; i<15; i++)
765  {
766  W[i] = C[i];
767  S[i] = C[i] + V[i];
768  }
769  for(int i=0; i<5; i++)
770  R[i] = r[i];
771 
772  InvertCholetsky(S);
773 
774  float K[5][5];
775  MultiplySS(C,S,K);
776  float dzeta[5];
777  for(int i=0; i<5; i++) dzeta[i] = m[i] - r[i];
778  float KC[15];
779  MultiplyMS(K,C,KC);
780  for(int i=0; i< 15; i++)
781  W[i] -= KC[i];
782 
783  float kd;
784  for(int i=0; i<5; i++)
785  {
786  kd = 0.f;
787  for(int j=0; j<5; j++)
788  kd += K[i][j]*dzeta[j];
789  R[i] += kd;
790  }
791  float S_dzeta[5];
792  MultiplySR(S, dzeta, S_dzeta);
793  chi2 = dzeta[0]*S_dzeta[0] + dzeta[1]*S_dzeta[1] + dzeta[2]*S_dzeta[2] + dzeta[3]*S_dzeta[3] + dzeta[4]*S_dzeta[4];
794 }
795 
796 void AliHLTTPCCAMerger::MakeBorderTracks(AliHLTTPCCABorderTrack B[], unsigned short &nB, unsigned char &iSlice)
797 {
798  //* prepare slice tracks for merging with next/previous/same sector
799 
800  for ( unsigned short itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ )
801  {
802  const AliHLTTPCCASliceTrackInfo &track = fTrackInfos[ fSliceTrackInfoStart[iSlice] + itr ];
803 
804 // if(track.NClusters() > 41) continue;
805 
806  AliHLTTPCCABorderTrack &bTr = B[nB];
807 
808  bTr.SetTrackID( itr );
809  bTr.Setb(track.InnerParam().DzDs());
810  bTr.SetbErr2(track.InnerParam().Err2DzDs());
811  bTr.Setp(track.InnerParam().QPt());
812  bTr.SetpErr2(track.InnerParam().Err2QPt());
813 
814  bTr.SetInnerRow( track.fInnerRow );
815  bTr.SetOuterRow( track.fOuterRow );
816 
817  nB++;
818  }
819 }
821 // #define BACK_ORDER_FOR_0 // little bit faster without it. (but why?)
822 void AliHLTTPCCAMerger::MergeBorderTracks( AliHLTTPCCABorderTrack B1[], int N1, int iSlice1, AliHLTTPCCABorderTrack B2[], int N2, int iSlice2, int number,unsigned short FirstTrIR[],unsigned short LastTrIR[])
823 {
824 // The function creates links to the inner and outer neighbours
825 
826 //if(number == 0) return;
827  const float factor2k = 64.f;
828 
829  sfloat_v dr_min2_localv(10000000.f);
830 // Find firts and last row for the nighbours finding. All tracks are sorted by the inner row.
831  for ( int i1 = 0; i1 < N1; i1++ ) {
832 
833  AliHLTTPCCABorderTrack &b1 = B1[i1];
834 
835  int lastIRow = fSliceParam.NRows()-1; // row to end finding of a clone
836  int firstIRow = 0;
837  int dir = -1; // direction from a start to an end row.
838  int ifirst2 = -1;
839  int ilast2 = N2;
840 
841  if(number == 1) { // find tracks with <= 1 common rows.
842  dir = -1; // Tracks2 sorted in order of decrease innerRow. For number==1 we'll go in the upper dir, so indices will decrease
843 
844  firstIRow = b1.OuterRow() + 0; // row to begin finding of clone
845  lastIRow = b1.OuterRow() + 11; // TODO choose parameter // todo parameters
846  lastIRow = (lastIRow < fSliceParam.NRows()) ? lastIRow : fSliceParam.NRows()-1;
847 
848  // find start index
849  ifirst2 = -1;
850  for(; firstIRow < fSliceParam.NRows(); firstIRow++)
851  if (LastTrIR[firstIRow] != 50000) {
852  ifirst2 = LastTrIR[firstIRow];
853  break;
854  }
855  // find end index
856  ilast2 = N2;
857  for(; lastIRow >= firstIRow; lastIRow--)
858  if(FirstTrIR[lastIRow] != 50000) {
859  ilast2 = FirstTrIR[lastIRow];
860  break;
861  }
862 
863  // std::cout << b1.OuterRow() << " " << firstIRow << " " << lastIRow << " " << ifirst2 << " " << ilast2 << std::endl; // dbg
864  }
865  else if(number == 0) { // find tracks with >= 2 common rows.
866 #ifdef BACK_ORDER_FOR_0
867  dir = 1; // Tracks2 sorted in order of decrease innerRow. For number==0 we'll go in the down dir, so indices will increase
868 
869  firstIRow = b1.OuterRow() - 1; // row to begin finding of a clone
870  lastIRow = b1.OuterRow() - 7; // TODO choose parameter // todo parameters
871  if ( firstIRow < 0 ) {
872  firstIRow = 0;
873  lastIRow = 0;
874  }
875  else if ( lastIRow < 0 ){
876  lastIRow = 0;
877  }
878 
879  // find start index
880  ifirst2 = N2;
881  for(; firstIRow >= 0; firstIRow--)
882  if (FirstTrIR[firstIRow] != 50000) {
883  ifirst2 = FirstTrIR[firstIRow];
884  break;
885  }
886  // find end index
887  ilast2 = -1;
888  for(; lastIRow <= firstIRow; lastIRow++)
889  if(LastTrIR[lastIRow] != 50000) {
890  ilast2 = LastTrIR[lastIRow];
891  break;
892  }
893 #else // BACK_ORDER_FOR_0
894  dir = -1;
895 
896  firstIRow = b1.OuterRow() - 7; // row to begin finding of a clone
897  lastIRow = b1.OuterRow() - 1; // TODO choose parameter // todo parameters
898  if ( lastIRow < 0 ) {
899  firstIRow = 0;
900  lastIRow = 0;
901  }
902  else if ( firstIRow < 0 ){
903  firstIRow = 0;
904  }
905 
906  // find start index
907  ifirst2 = -1;
908  for(; firstIRow < fSliceParam.NRows(); firstIRow++)
909  if (LastTrIR[firstIRow] != 50000) {
910  ifirst2 = LastTrIR[firstIRow];
911  break;
912  }
913  // find end index
914  ilast2 = N2;
915  for(; lastIRow >= firstIRow; lastIRow--)
916  if(FirstTrIR[lastIRow] != 50000) {
917  ilast2 = FirstTrIR[lastIRow];
918  break;
919  }
920 #endif // BACK_ORDER_FOR_0
921  }
922 
923 //std::cout << ifirst2 << " " << N2 << " " << N1 << std::endl;
924 // ifirst2 = (iSlice1 != iSlice2) ? 0 : i1+1;
925 
926  // rarely changed parameters
927  const AliHLTTPCCASliceTrackInfo *Tt1 = &fTrackInfos[fSliceTrackInfoStart[iSlice1] + b1.TrackID() ];
928  const float k1 = Tt1->InnerParam().QPt() * fSliceParam.cBz();
929 
930  sfloat_v maxLv(1e10f);
931  sfloat_v maxL2v(1e10f);
932 
933  float Tt2OuterAlpha=0;
934  float Tt2InnerAlpha=0;
935  float xE1_E2=0, xS1_S2=0, xS1_E2=0;
936  float yE1_E2=0, yS1_S2=0, yS1_E2=0;
937  float sinE1_E2=0, sinS1_S2=0, sinS1_E2=0;
938 
939  if (dir*ifirst2 <= dir*ilast2) {
940  Tt2OuterAlpha = fTrackInfos[fSliceTrackInfoStart[iSlice2] + B2[ifirst2].TrackID() ].OuterAlpha();
941  Tt2InnerAlpha = fTrackInfos[fSliceTrackInfoStart[iSlice2] + B2[ifirst2].TrackID() ].InnerAlpha();
942 // rotate coordinates to the same coordinate system
943  Tt1->OuterParam().RotateXY(Tt2OuterAlpha - Tt1->OuterAlpha(),xE1_E2,yE1_E2,sinE1_E2);
944  Tt1->InnerParam().RotateXY(Tt2InnerAlpha - Tt1->InnerAlpha(),xS1_S2,yS1_S2,sinS1_S2);
945  Tt1->InnerParam().RotateXY(Tt2OuterAlpha - Tt1->InnerAlpha(),xS1_E2,yS1_E2,sinS1_E2);
946  }
947 // Convert data of the first track to the SIMD vectors.
948  sfloat_v xE1_E2v(xE1_E2), xS1_S2v(xS1_S2), xS1_E2v(xS1_E2);
949  sfloat_v yE1_E2v(yE1_E2), yS1_S2v(yS1_S2), yS1_E2v(yS1_E2);
950  sfloat_v sinE1_E2v(sinE1_E2), sinS1_S2v(sinS1_S2), sinS1_E2v(sinS1_E2);
951  sfloat_v vTt2OuterAlpha(Tt2OuterAlpha), vTt2InnerAlpha(Tt2InnerAlpha);
952  int nVecElements = 0;
953 
954  ushort_v::Memory b2index;
955 
956  AliHLTTPCCATrackParamVector OutParT1,OutParT2,InParT1,InParT2;
957  sfloat_v OutAlphaT1( Tt1->OuterAlpha() ), OutAlphaT2, InAlphaT1( Tt1->InnerAlpha() ),InAlphaT2;
958  AliHLTTPCCATrackParam T2OuterParamMemory[ushort_v::Size],T2InnerParamMemory[ushort_v::Size];
959  sfloat_v::Memory T2InnerAlphaMemory,T2OuterAlphaMemory;
960 
961  OutParT1.SetX( sfloat_v( Tt1->OuterParam().X() ) );
962  OutParT1.SetSignCosPhi(sfloat_v( Tt1->OuterParam().SignCosPhi() ));
963  for(int iP=0; iP<5; iP++) { OutParT1.SetPar(iP,sfloat_v( Tt1->OuterParam().Par()[iP] )); }
964  for(int iC=0; iC<15; iC++){ OutParT1.SetCov(iC,sfloat_v( Tt1->OuterParam().Cov()[iC] )); }
965  OutParT1.SetChi2(sfloat_v( Tt1->OuterParam().Chi2() ));
966  OutParT1.SetNDF( Tt1->OuterParam().NDF() );
967 
968  InParT1.SetX(sfloat_v( Tt1->InnerParam().X() ));
969  InParT1.SetSignCosPhi(sfloat_v( Tt1->InnerParam().SignCosPhi() ));
970  for(int iP=0; iP<5; iP++) { InParT1.SetPar(iP,sfloat_v( Tt1->InnerParam().Par()[iP] )); }
971  for(int iC=0; iC<15; iC++){ InParT1.SetCov(iC,sfloat_v( Tt1->InnerParam().Cov()[iC] )); }
972  InParT1.SetChi2(sfloat_v( Tt1->InnerParam().Chi2() ));
973  InParT1.SetNDF(Tt1->InnerParam().NDF());
974 
975  const sfloat_v &yE1v = OutParT1.Y();
976  const sfloat_v &xE1v = OutParT1.X();
977 
978  for ( int i2 = ifirst2; dir*i2 <= dir*ilast2; i2 += dir ) {
979 
980  AliHLTTPCCABorderTrack &b2 = B2[i2];
981  bool Continue = 0;
982 // if dz/ds of the tracks differs more, than by several sigmas (now 8 is used) - they are not neighbours
983  float db2 = (b1.b()) - (b2.b());
984  db2 *= db2;
985  const float ddb2 = b1.bErr2() + b2.bErr2();
986  if( db2 > factor2k * ddb2 )
987  {
988  if(dir*i2 < dir*ilast2) continue;
989  else Continue = 1;
990  }
991 // if q/pt of the tracks differs more, than by several sigmas (now 8 is used) - they are not neighbours
992  float dp2 = (b1.p()) - (b2.p());
993  dp2 *= dp2;
994  const float ddp2 = b1.pErr2() + b2.pErr2();
995  if( dp2 > factor2k * ddp2 )
996  {
997  if(dir*i2 < dir*ilast2) continue;
998  else Continue = 1;
999  }
1000 // the track could not be itselfs neighbour
1001  if( ISUNLIKELY( iSlice1 == iSlice2 && b1.TrackID() == b2.TrackID() ) )
1002  {
1003  if(dir*i2 < dir*ilast2) continue;
1004  else Continue = 1;
1005  }
1006 
1007  AliHLTTPCCASliceTrackInfo *Tt2 = &fTrackInfos[fSliceTrackInfoStart[iSlice2] + b2.TrackID() ];
1008 
1009  if(!Continue){
1010 // store tracks, which passed previous cuts
1011  T2InnerParamMemory[nVecElements] = Tt2->InnerParam();
1012  T2OuterParamMemory[nVecElements] = Tt2->OuterParam();
1013  T2InnerAlphaMemory[nVecElements] = Tt2->InnerAlpha();
1014  T2OuterAlphaMemory[nVecElements] = Tt2->OuterAlpha();
1015  b2index[nVecElements] = i2;
1016  nVecElements++;
1017  }
1018  if(nVecElements == ushort_v::Size || (dir*i2 == dir*ilast2 && nVecElements > 0))
1019  {
1020 // convert parameters to SIMD vectors for the array of second tracks
1021  sfloat_m active = static_cast<sfloat_m>( ushort_v( Vc::IndexesFromZero ) < nVecElements );
1022 //std::cout << "1 "<<active<<std::endl;
1023  ConvertTrackParamToVector(T2InnerParamMemory,InParT2,nVecElements);
1024  ConvertTrackParamToVector(T2OuterParamMemory,OutParT2,nVecElements);
1025 
1026  OutAlphaT2.load(T2OuterAlphaMemory);
1027  InAlphaT2.load(T2InnerAlphaMemory);
1028 
1029  const sfloat_v &yE2v = OutParT2.Y();
1030  const sfloat_v &yS2v = InParT2.Y();
1031 
1032  const sfloat_v &xE2v = OutParT2.X();
1033  const sfloat_v &xS2v = InParT2.X();
1034 
1035  sfloat_v xS2_E1v(Vc::Zero);
1036  sfloat_v yS2_E1v(Vc::Zero);
1037  sfloat_v sinS2_E1v(Vc::Zero);
1038  InParT2.RotateXY(OutAlphaT1 - InAlphaT2,xS2_E1v,yS2_E1v,sinS2_E1v);
1039 
1040  const sfloat_m &rotate = active && (CAMath::Abs(vTt2OuterAlpha - OutAlphaT2) > 0.01f || CAMath::Abs(vTt2InnerAlpha - InAlphaT2)>0.01f);
1041 
1042  if(!(rotate.isEmpty()))
1043  {// recalculate values if they could change
1044  vTt2OuterAlpha(rotate) = OutAlphaT2;
1045  vTt2InnerAlpha(rotate) = InAlphaT2;
1046 
1047  OutParT1.RotateXY(vTt2OuterAlpha - OutAlphaT1,xE1_E2v,yE1_E2v,sinE1_E2v,rotate);
1048  InParT1.RotateXY(vTt2InnerAlpha - InAlphaT1,xS1_S2v,yS1_S2v,sinS1_S2v,rotate);
1049  InParT1.RotateXY(vTt2OuterAlpha - InAlphaT1,xS1_E2v,yS1_E2v,sinS1_E2v,rotate);
1050  }
1051 
1052  const sfloat_v dxArr[4] = {(xE1_E2v - xE2v),
1053  (xS1_S2v - xS2v),
1054  (xS1_E2v - xE2v),
1055  (xS2_E1v - xE1v)};
1056  const sfloat_v dyArr[4] = {(yE1_E2v - yE2v),
1057  (yS1_S2v - yS2v),
1058  (yS1_E2v - yE2v),
1059  (yS2_E1v - yE1v)};
1060 
1061  for(int iCycl=0; iCycl<1; iCycl++)
1062  {
1063  const sfloat_v dx_1 = xE2v - xS1_E2v;
1064  const sfloat_v t_sin1 = k1 * dx_1 + sinS1_E2v;
1065 // check, whether inner parameters of both tracks in the pair could be transported to the outer parameters of the other track in the pair
1066  const sfloat_v k2 = InParT2.QPt() * fSliceParam.cBz();
1067  const sfloat_v dx_2 = xE1v - xS2_E1v;
1068  const sfloat_v t_sin2 = k2 * dx_2 + sinS2_E1v;
1069  active &= CAMath::Abs( t_sin1 ) <= 0.999f && CAMath::Abs(t_sin2) <= 0.999f;
1070  if(active.isEmpty()) continue;
1071 // find the closest pair of the tracks parameters and store them to Thelp and Thelp1
1072  const sfloat_v dxyArr[4] = {dxArr[0]*dxArr[0] + dyArr[0]*dyArr[0],
1073  dxArr[1]*dxArr[1] + dyArr[1]*dyArr[1],
1074  dxArr[2]*dxArr[2] + dyArr[2]*dyArr[2],
1075  dxArr[3]*dxArr[3] + dyArr[3]*dyArr[3]};
1076 
1077  short_v idx(Vc::Zero), idx_temp1(Vc::Zero), idx_temp2((short int)2);
1078  sfloat_v dx = CAMath::Abs(dxyArr[0]);
1079  sfloat_v dx_temp1 = CAMath::Abs(dxyArr[0]), dx_temp2 = CAMath::Abs(dxyArr[2]);
1080  const sfloat_m &tmp1mask = dxyArr[1] < dx_temp1;
1081  const sfloat_m &tmp2mask = dxyArr[3] < dx_temp2;
1082  dx_temp1(tmp1mask) = dxyArr[1];
1083  idx_temp1(static_cast<short_m>(tmp1mask)) = 1;
1084  dx_temp2(tmp2mask) = dxyArr[3];
1085  idx_temp2(static_cast<short_m>(tmp2mask)) = 3;
1086 
1087  if(number == 1)
1088  { dx = dx_temp2; idx = idx_temp2; }
1089  else if (number == 0)
1090  {
1091  const sfloat_m &Isdx2 = dx_temp2 < dx_temp1;
1092  const short_m &Isdx2_short = static_cast<short_m>(Isdx2);
1093 
1094  dx(Isdx2) = dx_temp2;
1095  idx(Isdx2_short) = idx_temp2;
1096  dx(!Isdx2) = dx_temp1;
1097  idx(!Isdx2_short) = idx_temp1;
1098  }
1099 
1100  AliHLTTPCCATrackParamVector Thelp, Thelp1;
1101  sfloat_v a1(Vc::Zero);
1102  sfloat_v a2(Vc::Zero);
1103 
1104  if(number == 0)
1105  {
1106  const sfloat_m &IsIdx0 = active && CAMath::Abs(static_cast<sfloat_v>(idx) - 0.f) < 0.1f;
1107  const sfloat_m &IsIdx1 = active && CAMath::Abs(static_cast<sfloat_v>(idx) - 1.f) < 0.1f;
1108  Thelp.SetTrackParam(OutParT1, IsIdx0 );
1109  Thelp1.SetTrackParam(OutParT2, IsIdx0 );
1110  a1( IsIdx0 ) = OutAlphaT1;
1111  a2( IsIdx0 ) = OutAlphaT2;
1112  dx( IsIdx0 ) = dxArr[0];
1113  Thelp.SetTrackParam(InParT1, IsIdx1 );
1114  Thelp1.SetTrackParam(InParT2, IsIdx1 );
1115  a1( IsIdx1 ) = InAlphaT1;
1116  a2( IsIdx1 ) = InAlphaT2;
1117  dx( IsIdx1 ) = dxArr[1];
1118  }
1119  const sfloat_m &IsIdx2 = active && CAMath::Abs(static_cast<sfloat_v>(idx) - 2.f) < 0.1f;
1120  const sfloat_m &IsIdx3 = active && CAMath::Abs(static_cast<sfloat_v>(idx) - 3.f) < 0.1f;
1121  Thelp.SetTrackParam(InParT1, IsIdx2 );
1122  Thelp1.SetTrackParam(OutParT2, IsIdx2 );
1123  a1( IsIdx2 ) = InAlphaT1;
1124  a2( IsIdx2 ) = OutAlphaT2;
1125  dx( IsIdx2 ) = dxArr[2];
1126 
1127  Thelp.SetTrackParam(OutParT1, IsIdx3 );
1128  Thelp1.SetTrackParam(InParT2, IsIdx3 );
1129  a1( IsIdx3 ) = OutAlphaT1;
1130  a2( IsIdx3 ) = InAlphaT2;
1131  dx( IsIdx3 ) = dxArr[3];
1132  // on the first iteration we process only overlaped tracks
1133  if(number == 1) active &= (CAMath::Abs(dx)<0.1f || dx > 0.f);
1134 // on second - only nonoverlaped
1135  if(number == 0) active &= (CAMath::Abs(dx)<0.1f || dx < 0.f);
1136  if(active.isEmpty()) continue;
1137 // rotate track parameters to the same coordinate system
1138  active &= Thelp.Rotate( a2 - a1 , 0.999f, active);
1139  if(active.isEmpty()) continue;
1140 // calculate distance between tracks
1141  dr_min2_localv = (Thelp.Y()-Thelp1.Y())*(Thelp.Y()-Thelp1.Y()) + (Thelp.Z()-Thelp1.Z())*(Thelp.Z()-Thelp1.Z()) + dx*dx;
1142 
1143  if ( number == 1 ) active &= (dr_min2_localv <= maxLv*maxLv);
1144  if(active.isEmpty()) continue;
1145  maxL2v(active) = dr_min2_localv; // save dist^2 for lastIRow calculation
1146 // transport tracks to the middle point by X coordinate
1147  const sfloat_v &ToX = 0.5f * (Thelp.X() + Thelp1.X());
1148  active &= Thelp.TransportToX( ToX, sfloat_v(fSliceParam.cBz()), 0.999f, active);
1149  if(active.isEmpty()) continue;
1150  active &= Thelp1.TransportToX( ToX, sfloat_v(fSliceParam.cBz()), 0.999f, active);
1151  if(active.isEmpty()) continue;
1152 // if tracks are too far one from each other after propagation - they could not be neighbours
1153  active &= CAMath::Abs(Thelp1.Y() - Thelp.Y()) <= 10.f;
1154  active &= CAMath::Abs(Thelp1.Z() - Thelp.Z()) <= 10.f;
1155  active &= CAMath::Abs(Thelp1.SinPhi() - Thelp.SinPhi()) <= 0.15f;
1156  if(active.isEmpty()) continue;
1157 // merge parameters and calculate chi2
1158  sfloat_v C[15], r[5], chi2(Vc::Zero);
1159  FilterTracks(Thelp.GetPar(), Thelp.GetCov(), Thelp1.GetPar(), Thelp1.GetCov(), r, C, chi2, active);
1160 // if after merging parameters are infinite, or diagonal elements of cov. matrix are negative - tracks could not be neighbours
1161  for ( int i = 0; i < 15; i++ ) active &= CAMath::Finite( C[i] );
1162  for ( int i = 0; i < 5; i++ ) active &= CAMath::Finite( r[i] );
1163  active &= ( C[0] > 0 ) && ( C[2] > 0 ) && ( C[5] > 0 ) && ( C[9] > 0 ) && ( C[14] > 0 );
1164 // if obtained chi2 value is too high - tracks could not be neighbours
1165  if(number == 1) active &= (chi2<=100.f);
1166  if(number == 0) {
1167  active &= (chi2<=300.f);
1168  dr_min2_localv(active) = chi2;
1169  }
1170  if(active.isEmpty()) continue;
1171 
1172  for(int iV=0; iV < nVecElements; iV++)
1173  {
1174 // determine, whether neighbour is inner or outer
1175  if(!active[iV]) continue;
1176  AliHLTTPCCASliceTrackInfo *T1, *T2;
1177  int NextNew, SliceNextNew, PrevNew, SlicePrevNew;
1178 
1179  bool IsNext = (b1.InnerRow() < B2[b2index[iV]].InnerRow()) ||
1180  (b1.InnerRow() == B2[b2index[iV]].InnerRow() &&
1181  b1.OuterRow() < B2[b2index[iV]].OuterRow()) ||
1182  (b1.InnerRow() == B2[b2index[iV]].InnerRow() &&
1183  b1.OuterRow() == B2[b2index[iV]].OuterRow() &&
1184  b1.TrackID() < B2[b2index[iV]].TrackID());
1185  if(IsNext)
1186  {
1187  T1 = &fTrackInfos[fSliceTrackInfoStart[iSlice1] + b1.TrackID() ];
1188  T2 = &fTrackInfos[fSliceTrackInfoStart[iSlice2] + B2[b2index[iV]].TrackID() ];
1189 
1190  NextNew = B2[b2index[iV]].TrackID();
1191  SliceNextNew = iSlice2;
1192  PrevNew = b1.TrackID();
1193  SlicePrevNew = iSlice1;
1194  }
1195  else
1196  {
1197  T2 = &fTrackInfos[fSliceTrackInfoStart[iSlice1] + b1.TrackID() ];
1198  T1 = &fTrackInfos[fSliceTrackInfoStart[iSlice2] + B2[b2index[iV]].TrackID() ];
1199 
1200  NextNew = b1.TrackID();
1201  SliceNextNew = iSlice1;
1202  PrevNew = B2[b2index[iV]].TrackID();
1203  SlicePrevNew = iSlice2;
1204  }
1205 // if current neighbour is better, then previus one - make a reference to it
1206  if(T1->ChiNext > dr_min2_localv[iV] && T2->ChiPrev > dr_min2_localv[iV]) // clone was found
1207  {
1208  // reestimate end row and index
1209  maxLv = sqrt(maxL2v[iV]);
1210  if (number == 1) {
1211  // find end row
1212  const float x0 = fSliceParam.RowX(b1.OuterRow());
1213  for (; fSliceParam.RowX(lastIRow) - x0 > maxLv[iV]; lastIRow--);
1214  // find end index
1215  for(; lastIRow >= firstIRow; lastIRow--)
1216  if(FirstTrIR[lastIRow] != 50000) {
1217  ilast2 = FirstTrIR[lastIRow];
1218  break;
1219  }
1220  }
1221 #ifdef BACK_ORDER_FOR_0
1222  else if (number == 0) {
1223  // find end row
1224  const float x0 = fSliceParam.RowX(b1.OuterRow());
1225  for (; x0 - fSliceParam.RowX(lastIRow) > maxL; lastIRow++);
1226  // find end index
1227  for(; lastIRow <= firstIRow; lastIRow++)
1228  if(LastTrIR[lastIRow] != 50000) {
1229  // ilast2 = LastTrIR[lastIRow];
1230  break;
1231  }
1232  }
1233  // std::cout << b1.OuterRow() << " " << lastIRow << " " << sqrt(maxL) << " " << ifirst2 << " " << ilast2 << std::endl; // dbg
1234 #endif
1235  if(T1->NextNeighbour() > -1)
1236  {
1237  AliHLTTPCCASliceTrackInfo *T3 = &fTrackInfos[fSliceTrackInfoStart[T1->SliceNextNeighbour()] + T1->NextNeighbour() ];
1238  T3->SetPrevNeighbour(-2); // mark that T3 lose best neighbour (so there can be other neighbour for T3)
1239  T3->ChiPrev = 10000000;
1240  }
1241  if(T2->PrevNeighbour() > -1)
1242  {
1243  AliHLTTPCCASliceTrackInfo *T3 = &fTrackInfos[fSliceTrackInfoStart[T2->SlicePrevNeighbour()] + T2->PrevNeighbour() ];
1244  T3->SetNextNeighbour(-2);
1245  T3->ChiNext = 10000000;
1246  }
1247  T1->SetNextNeighbour( NextNew );
1248  T1->SetSliceNextNeighbour( SliceNextNew );
1249  T1->ChiNext = dr_min2_localv[iV];
1250 
1251  T2->SetPrevNeighbour( PrevNew );
1252  T2->SetSlicePrevNeighbour( SlicePrevNew );
1253  T2->ChiPrev = dr_min2_localv[iV];
1254  }
1255  }
1256  }
1257  nVecElements = 0;
1258  }
1259  }
1260  }
1261 
1262 }
1263 
1264 void AliHLTTPCCAMerger::Merging(int number)
1265 {
1266 #ifdef USE_TIMERS
1267  Stopwatch timer;
1268  timer.Start();
1269 #endif // USE_TIMERS
1270 
1271  //* track merging between slices
1272 
1273  fOutput->SetNTracks( 0 );
1274  fOutput->SetNTrackClusters( 0 );
1275  fOutput->SetPointers();
1276 
1277  // for each slice set number of the next neighbouring slice
1278  int nextSlice[fgkNSlices], oppSlice[fgkNSlices/2];
1279  // int prevSlice[fgkNSlices];
1280 
1281  int mid = fgkNSlices / 2 - 1 ;
1282  int last = fgkNSlices - 1 ;
1283 
1284  for ( int iSlice = 0; iSlice < fgkNSlices/2; iSlice++ ) {
1285  nextSlice[iSlice] = iSlice + 1;
1286 // prevSlice[iSlice] = iSlice - 1;
1287  oppSlice [iSlice] = 22 - iSlice;
1288  if(iSlice == 11) oppSlice[iSlice]=23;
1289  }
1290 
1291  for ( int iSlice = fgkNSlices/2; iSlice < fgkNSlices; iSlice++ ) {
1292  nextSlice[iSlice] = iSlice - 1;
1293 // prevSlice[iSlice] = iSlice + 1;
1294  }
1295 
1296  if ( mid < 0 ) mid = 0; // to avoid compiler warning
1297  if ( last < 0 ) last = 0; //
1298  nextSlice[ mid ] = 0;
1299 // prevSlice[ 0 ] = mid;
1300 // prevSlice[ last ] = fgkNSlices / 2;
1301  nextSlice[ fgkNSlices/2 ] = last;
1302 
1303  int maxNSliceTracks = 0;
1304  for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
1305  if ( maxNSliceTracks < fSliceNTrackInfos[iSlice] ) maxNSliceTracks = fSliceNTrackInfos[iSlice];
1306  }
1307 
1308  AliHLTTPCCABorderTrack *bCurrIR = new AliHLTTPCCABorderTrack[maxNSliceTracks*fgkNSlices];
1309  AliHLTTPCCABorderTrack *bCurrOR = new AliHLTTPCCABorderTrack[maxNSliceTracks*fgkNSlices];
1310  unsigned short FirstTrIR[fgkNSlices][AliHLTTPCCAParameters::MaxNumberOfRows8]; // index of the first track on row
1311  unsigned short LastTrIR[fgkNSlices][AliHLTTPCCAParameters::MaxNumberOfRows8];
1312 
1313  unsigned short nCurr[fgkNSlices];
1314  for(unsigned char iSl=0; iSl<fgkNSlices; iSl++)
1315  {
1316  for(unsigned char iRow=0; iRow<fSliceParam.NRows(); iRow++)
1317  {
1318  FirstTrIR[iSl][iRow] = 50000;
1319  LastTrIR[iSl][iRow] = 50000;
1320  }
1321 
1322  nCurr[iSl] = 0;
1323 // make border tracks for each sector, sort them by inner row
1324  MakeBorderTracks(bCurrIR + maxNSliceTracks*iSl, nCurr[iSl], iSl);
1325  std::sort(bCurrIR+maxNSliceTracks*iSl, bCurrIR+maxNSliceTracks*iSl+nCurr[iSl], CompareInnerRow);
1326 
1327  for(unsigned short itr=0; itr < nCurr[iSl]; itr++)
1328  bCurrOR[itr + maxNSliceTracks*iSl] = bCurrIR[itr + maxNSliceTracks*iSl];
1329  // std::sort(bCurrOR+maxNSliceTracks*iSl, bCurrOR+maxNSliceTracks*iSl+nCurr[iSl], CompareOuterRow);
1330 
1331  // save track indices range for each row
1332  if(nCurr[iSl] > 0)
1333  {
1334  unsigned char curRow = bCurrIR[maxNSliceTracks*iSl].InnerRow();
1335  FirstTrIR[iSl][curRow] = 0;
1336  for(unsigned short itr = 1; itr < nCurr[iSl]; itr++)
1337  {
1338  if( bCurrIR[maxNSliceTracks*iSl+itr].InnerRow() < curRow )
1339  {
1340  // for (int iR = bCurrIR[maxNSliceTracks*iSl+itr].InnerRow() + 1; iR < curRow; iR++) // check for a jump
1341  // FirstTrIR[iSl][iR] = FirstTrIR[iSl][curRow];
1342  LastTrIR[iSl][curRow] = itr - 1;
1343  curRow = bCurrIR[maxNSliceTracks*iSl + itr].InnerRow();
1344  FirstTrIR[iSl][curRow] = itr;
1345  }
1346  }
1347  LastTrIR[iSl][curRow] = nCurr[iSl] - 1;
1348  // for(unsigned char iRow=0; iRow<AliHLTTPCCAParameters::NumberOfRows; iRow++) { // dbg
1349  // std::cout << int(iSl) << " " << int(iRow) << " " << FirstTrIR[iSl][iRow] << " " << LastTrIR[iSl][iRow] << std::endl;
1350  // if (FirstTrIR[iSl][iRow] != 50000 ) std::cout << bCurrIR[maxNSliceTracks*iSl+FirstTrIR[iSl][iRow]].InnerRow() << std::endl;
1351  // if (LastTrIR[iSl][iRow] != 50000 ) std::cout << bCurrIR[maxNSliceTracks*iSl+LastTrIR[iSl][iRow]].InnerRow() << std::endl;
1352  // }
1353  }
1354 // create links to neighbour tracks clones in the same sector
1355  MergeBorderTracks(bCurrOR + maxNSliceTracks*iSl, nCurr[iSl], iSl,
1356  bCurrIR + maxNSliceTracks*iSl, nCurr[iSl], iSl,
1357  number, FirstTrIR[iSl], LastTrIR[iSl]);
1358  }
1359 
1360  if (! fgDoNotMergeBorders )
1361  for(int iSl=0; iSl<fgkNSlices; iSl++)
1362  {
1363  // create links to neighbour tracks in the neighbour sector in the same xy-plane
1364  MergeBorderTracks( bCurrOR+maxNSliceTracks*iSl, nCurr[iSl], iSl,
1365  bCurrIR+maxNSliceTracks*nextSlice[iSl], nCurr[nextSlice[iSl]], nextSlice[iSl],
1366  number, FirstTrIR[nextSlice[iSl]], LastTrIR[nextSlice[iSl]] ); // merge upper edges
1367  MergeBorderTracks( bCurrOR+maxNSliceTracks*nextSlice[iSl], nCurr[nextSlice[iSl]], nextSlice[iSl],
1368  bCurrIR+maxNSliceTracks*iSl, nCurr[iSl], iSl,
1369  number, FirstTrIR[iSl], LastTrIR[iSl] ); // merge lower edges
1370 
1371  if(iSl < fgkNSlices / 2)
1372  {
1373  // create links to neighbour tracks with the oposit sector (in z direction)
1374  MergeBorderTracks(bCurrOR+maxNSliceTracks*iSl, nCurr[iSl], iSl,
1375  bCurrIR+maxNSliceTracks*oppSlice[iSl], nCurr[oppSlice[iSl]], oppSlice[iSl],
1376  number, FirstTrIR[oppSlice[iSl]], LastTrIR[oppSlice[iSl]]);
1377  MergeBorderTracks(bCurrOR+maxNSliceTracks*oppSlice[iSl], nCurr[oppSlice[iSl]], oppSlice[iSl],
1378  bCurrIR+maxNSliceTracks*iSl, nCurr[iSl], iSl,
1379  number, FirstTrIR[iSl], LastTrIR[iSl]);
1380  }
1381  }
1382 
1383  if ( bCurrIR ) delete[] bCurrIR;
1384  if ( bCurrOR ) delete[] bCurrOR;
1385 
1386 #ifdef USE_TIMERS
1387  timer.Stop();
1388  fTimers[3+(1-number)] = timer.RealTime();
1389  timer.Start();
1390 #endif // USE_TIMERS
1391 
1392  int nOutTracks = 0;
1393  int nOutTrackClusters = 0;
1394 
1395  AliHLTTPCCAMergedTrack *outTracks = 0;
1396  DataCompressor::SliceRowCluster *outClusterIDsrc = 0;
1397  UChar_t *outClusterPackedAmp = 0;
1398 
1399  if(number == 0)
1400  {
1401  outTracks = new AliHLTTPCCAMergedTrack[fMaxTrackInfos];
1402  outClusterIDsrc = new DataCompressor::SliceRowCluster[fMaxClusterInfos];
1403  outClusterPackedAmp = new UChar_t [fMaxClusterInfos];
1404  }
1405 
1406  AliHLTTPCCAClusterInfo *tmpH = new AliHLTTPCCAClusterInfo[fMaxClusterInfos];
1407  AliHLTTPCCASliceTrackInfo *tmpT = new AliHLTTPCCASliceTrackInfo[fMaxTrackInfos];
1408  int nEndTracks = 0; // n tracks after merging.
1409  int tmpSliceTrackInfoStart[fgkNSlices];
1410 
1411  int nTrNew[fgkNSlices];
1412 
1413  for(int iSlice=0; iSlice < fgkNSlices; iSlice++) nTrNew[iSlice]=0;
1414 
1415  int nH = 0;
1416 
1418 
1419 // merge tracks, using obtained links to neighbours
1420  for ( int iSlice = 0; iSlice < fgkNSlices; iSlice++ ) {
1421 
1422  tmpSliceTrackInfoStart[iSlice] = nEndTracks;
1423  assert( iSlice == 0 || nEndTracks == tmpSliceTrackInfoStart[iSlice-1] + nTrNew[iSlice-1] );
1424 
1425  AliHLTTPCCASliceTrackInfo *vTrackOld[short_v::Size];
1426  AliHLTTPCCATrackParam StartPoint[short_v::Size];
1427  AliHLTTPCCATrackParam EndPoint[short_v::Size];
1428  sfloat_v::Memory StartAlpha, EndAlpha;
1429  int nVecElements = 0;
1430 
1431  sfloat_v vStartAlpha(Vc::Zero);
1432  sfloat_v vEndAlpha(Vc::Zero);
1433  AliHLTTPCCATrackParamVector vStartPoint;
1434  AliHLTTPCCATrackParamVector vEndPoint;
1435 
1436  AliHLTTPCCASliceTrackInfo *TInfos[10000];
1437  bool IsSorted[10000];
1438 // Resort tracks.
1439  int nSorted = 0;
1440 // First write tracks, which has inner (previous) neighbour
1441  for(int i=0; i< fSliceNTrackInfos[iSlice]; i++)
1442  {
1443  IsSorted[i] = 0;
1444  TInfos[i] = 0;
1445  AliHLTTPCCASliceTrackInfo *tr = fTrackInfos + fSliceTrackInfoStart[iSlice] + i;
1446  if(tr->PrevNeighbour() > -1)
1447  {
1448  TInfos[nSorted] = tr;
1449  IsSorted[i] = 1;
1450  nSorted++;
1451  }
1452  }
1453 // Then - tracks without any neighbours
1454  for(int i=0; i< fSliceNTrackInfos[iSlice]; i++)
1455  {
1456  AliHLTTPCCASliceTrackInfo *tr = fTrackInfos + fSliceTrackInfoStart[iSlice] + i;
1457  if(tr->PrevNeighbour() < 0 && tr->NextNeighbour() < 0 && !IsSorted[i])
1458  {
1459  TInfos[nSorted] = tr;
1460  IsSorted[i] = 1;
1461  nSorted++;
1462  }
1463  }
1464 // Then - tracks with outetr (next) neighbours only
1465  for(int i=0; i< fSliceNTrackInfos[iSlice]; i++)
1466  {
1467  AliHLTTPCCASliceTrackInfo *tr = fTrackInfos + fSliceTrackInfoStart[iSlice] + i;
1468  if(!IsSorted[i])
1469  {
1470  TInfos[nSorted] = tr;
1471  IsSorted[i] = 1;
1472  nSorted++;
1473  }
1474  }
1475 
1476  for ( int itr = 0; itr < fSliceNTrackInfos[iSlice]; itr++ ) {
1477 
1478  AliHLTTPCCASliceTrackInfo &trackOld = *TInfos[itr];
1479 
1480  bool Continue = 0;
1481 // If track is already used or has previous neighbours - do not use it
1482  if ( trackOld.Used() ) Continue = 1;
1483  if ( trackOld.PrevNeighbour() > -1 ) Continue = 1;
1484 
1485  if(!Continue){
1486 // Store selected tracks
1487  vTrackOld[nVecElements] = &trackOld;
1488  StartPoint[nVecElements] = trackOld.InnerParam();
1489  EndPoint[nVecElements] = trackOld.OuterParam();
1490  StartAlpha[nVecElements] = trackOld.InnerAlpha();
1491  EndAlpha[nVecElements] = trackOld.OuterAlpha();
1492  nVecElements++;
1493  }
1494  if((nVecElements == short_v::Size) || (itr == fSliceNTrackInfos[iSlice]-1))
1495  {
1496 // Convert their parameters to SIMD vectors
1497  ushort_m active_short = ushort_v( Vc::IndexesFromZero ) < nVecElements;
1498  sfloat_m active = static_cast<sfloat_m>( active_short );
1499 
1500  int hits[2000][ushort_v::Size];
1501  ushort_v::Memory nHits;
1502  ushort_v firstHit = (unsigned short int)1000;
1503  ushort_v::Memory SliceNextNeighbour;
1504  short_v::Memory NextNeighbour;
1505  ushort_v jSlice(iSlice);
1506  short_v jtr = (short int) -1;
1507  ushort_v vNHits(Vc::Zero);
1508 
1509  short_v vNextNeighbour;
1510  ushort_v vSliceNextNeighbour;
1511 
1512  ConvertTrackParamToVector(StartPoint,vStartPoint,nVecElements);
1513  ConvertTrackParamToVector(EndPoint,vEndPoint,nVecElements);
1514  vStartAlpha.load(StartAlpha);
1515  vEndAlpha.load(EndAlpha);
1516 
1517  const sfloat_m &invert1 = active && (vEndPoint.X() < vStartPoint.X());
1518  if( !(invert1.isEmpty()))
1519  {
1520  AliHLTTPCCATrackParamVector helpPoint = vEndPoint;
1521  vEndPoint.SetTrackParam(vStartPoint, invert1);
1522  vStartPoint.SetTrackParam(helpPoint, invert1);
1523  sfloat_v helpAlpha = vEndAlpha;
1524  vEndAlpha(invert1) = vStartAlpha;
1525  vStartAlpha(invert1) = helpAlpha;
1526  }
1527 
1528  for(int iV=0; iV<ushort_v::Size; iV++)
1529  {
1530  if(!active[iV]) continue;
1531  vTrackOld[iV]->SetUsed( 1 ); // ste track as used
1532 
1533  for ( int jhit = 0; jhit < vTrackOld[iV]->NClusters(); jhit++ ) {
1534  int id = vTrackOld[iV]->FirstClusterRef() + jhit;
1535  hits[static_cast <unsigned short>(firstHit[iV])+jhit][iV] = id;
1536  }
1537  nHits[iV] = vTrackOld[iV]->NClusters();
1538  NextNeighbour[iV] = vTrackOld[iV]->NextNeighbour(); // obtaine the outer neighbour index
1539  SliceNextNeighbour[iV] = vTrackOld[iV]->SliceNextNeighbour();
1540  }
1541  vNHits.load(nHits);
1542  vNextNeighbour.load(NextNeighbour);
1543  vSliceNextNeighbour.load(SliceNextNeighbour);
1544 
1545  const sfloat_m &isMerged = active && (static_cast<sfloat_v>(vNextNeighbour) > -1);
1546  sfloat_m IsNeighbour = isMerged;
1547 
1548  if(!(isMerged.isEmpty()))
1549  {
1550  jtr(static_cast<short_m>(isMerged)) = vNextNeighbour;
1551  jSlice(static_cast<ushort_m>(isMerged)) = vSliceNextNeighbour;
1552  }
1553 
1554  while ( !(IsNeighbour.isEmpty()) ) // while there are still outer neighbours
1555  {
1556  AliHLTTPCCATrackParam InnerParam[short_v::Size];
1557  AliHLTTPCCATrackParam OuterParam[short_v::Size];
1558  sfloat_v::Memory InnerAlpha;
1559  sfloat_v::Memory OuterAlpha;
1560  AliHLTTPCCATrackParamVector vInnerParam;
1561  AliHLTTPCCATrackParamVector vOuterParam;
1562  sfloat_v vInnerAlpha(Vc::Zero);
1563  sfloat_v vOuterAlpha(Vc::Zero);
1564 
1565  ushort_v::Memory NClusters;
1566  ushort_v vNClusters;
1567 
1568  ushort_v::Memory SegmSliceNextNeighbour;
1569  short_v::Memory SegmNextNeighbour;
1570  short_v vSegmNextNeighbour;
1571  ushort_v vSegmSliceNextNeighbour;
1572 
1573  AliHLTTPCCASliceTrackInfo *segment[ushort_v::Size];
1574  short_v::Memory Used;
1575  short_v vsUsed;
1576  short_m vUsed(false);
1577 // take the next neighbour
1578  for(int iV=0; iV<ushort_v::Size; iV++)
1579  {
1580  if(!IsNeighbour[iV]) continue;
1581  segment[iV] = &fTrackInfos[fSliceTrackInfoStart[jSlice[iV]] + jtr[iV]];
1582  InnerParam[iV] = segment[iV]->InnerParam();
1583  OuterParam[iV] = segment[iV]->OuterParam();
1584  InnerAlpha[iV] = segment[iV]->InnerAlpha();
1585  OuterAlpha[iV] = segment[iV]->OuterAlpha();
1586  Used[iV] = segment[iV]->Used();
1587  NClusters[iV] = segment[iV]->NClusters();
1588  SegmNextNeighbour[iV] = segment[iV]->NextNeighbour();
1589  SegmSliceNextNeighbour[iV] = segment[iV]->SliceNextNeighbour();
1590  }
1591  vsUsed.load(Used);
1592  vNClusters.load(NClusters);
1593  vUsed = vsUsed > 0;
1594  vSegmNextNeighbour.load(SegmNextNeighbour);
1595  vSegmSliceNextNeighbour.load(SegmSliceNextNeighbour);
1596 
1597  IsNeighbour &= !(static_cast<sfloat_m>(vUsed));
1598  if ( IsNeighbour.isEmpty() ) break;
1599 
1600  ConvertTrackParamToVector(InnerParam,vInnerParam,nVecElements);
1601  ConvertTrackParamToVector(OuterParam,vOuterParam,nVecElements);
1602  vInnerAlpha.load(InnerAlpha);
1603  vOuterAlpha.load(OuterAlpha);
1604 
1605  for(int iV=0; iV<ushort_v::Size; iV++)
1606  {
1607  if(!IsNeighbour[iV]) continue;
1608  segment[iV]->SetUsed( 1 );
1609  }
1610  sfloat_m dir(false);
1611  ushort_v startHit = firstHit + vNHits;
1612  const sfloat_v d00 = vStartPoint.GetDistXZ2( vInnerParam );
1613  const sfloat_v d01 = vStartPoint.GetDistXZ2( vOuterParam );
1614  const sfloat_v d10 = vEndPoint.GetDistXZ2( vInnerParam );
1615  const sfloat_v d11 = vEndPoint.GetDistXZ2( vOuterParam );
1616  const sfloat_v dz00 = abs( vStartPoint.Z() - vInnerParam.Z() );
1617  const sfloat_v dz01 = abs( vStartPoint.Z() - vOuterParam.Z() );
1618  const sfloat_v dz10 = abs( vEndPoint.Z() - vInnerParam.Z() );
1619  const sfloat_v dz11 = abs( vEndPoint.Z() - vOuterParam.Z() );
1620 
1621  const sfloat_m &case1 = IsNeighbour &&
1622  (d00 <= d01 && d00 <= d10 && d00 <= d11) &&
1623  (dz11 >= dz00 && dz11 >= dz01 && dz11 >= dz10 );
1624  /* 0----1
1625  * connection like : \
1626  * 0------1
1627  */
1628  const sfloat_m &case2 = IsNeighbour &&
1629  (d01 < d00 && d01 <= d10 && d01 <= d11) &&
1630  (dz10 > dz00 && dz10 >= dz01 && dz10 >= dz11 );
1631  /* 0----1
1632  * connection like : \
1633  * 0----1
1634  */
1635  const sfloat_m &case3 = IsNeighbour &&
1636  (d10 < d00 && d10 <= d01 && d10 <= d11) &&
1637  (dz01 > dz00 && dz01 >= dz10 && dz01 >= dz11 );
1638  /* 0---1
1639  * connection like : /
1640  * 0-------1
1641  */
1642  const sfloat_m &case4 = IsNeighbour &&
1643  (d11 < d00 && d11 <= d10 && d10 <= d01) &&
1644  (dz00 > dz01 && dz00 >= dz10 && dz00 > dz11);
1645  /* 0--1
1646  * connection like : \
1647  * 0------1
1648  */
1649  IsNeighbour &= case1 || case2 || case3 || case4;
1650 
1651  if(!(case1.isEmpty()))
1652  {
1653  vStartPoint.SetTrackParam( vOuterParam, case1);
1654  vStartAlpha(case1) = vOuterAlpha;
1655  firstHit(static_cast<ushort_m>(case1)) -= vNClusters;
1656  startHit(static_cast<ushort_m>(case1)) = firstHit;
1657  }
1658  if(!(case2.isEmpty()))
1659  {
1660  vStartPoint.SetTrackParam( vInnerParam, case2);
1661  vStartAlpha(case2) = vInnerAlpha;
1662  firstHit(static_cast<ushort_m>(case2)) -= vNClusters;
1663  startHit(static_cast<ushort_m>(case2)) = firstHit;
1664  }
1665  if(!(case3.isEmpty()))
1666  {
1667  vEndPoint.SetTrackParam( vOuterParam, case3);
1668  vEndAlpha(case3) = vOuterAlpha;
1669  }
1670  if(!(case4.isEmpty()))
1671  {
1672  vEndPoint.SetTrackParam( vInnerParam, case4);
1673  vEndAlpha(case4) = vInnerAlpha;
1674  }
1675 
1676  const sfloat_m &m1 = IsNeighbour && (vEndPoint.X() < vOuterParam.X());
1677  vEndPoint.SetTrackParam( vOuterParam, m1);
1678  vEndAlpha(m1) = vOuterAlpha;
1679 
1680  const sfloat_m &m2 = IsNeighbour && (vEndPoint.X() < vInnerParam.X());
1681  vEndPoint.SetTrackParam( vInnerParam, m2);
1682  vEndAlpha(m2) = vInnerAlpha;
1683 
1684  const sfloat_m &m3 = IsNeighbour && (vStartPoint.X() > vInnerParam.X());
1685  vStartPoint.SetTrackParam( vInnerParam, m3 );
1686  vStartAlpha(m3) = vInnerAlpha;
1687 
1688  const sfloat_m &m4 = IsNeighbour && (vStartPoint.X() > vOuterParam.X());
1689  vStartPoint.SetTrackParam( vOuterParam, m4 );
1690  vStartAlpha(m4) = vOuterAlpha;
1691 
1692  const sfloat_m dir_float = (case1 || case4);
1693  dir = dir_float;
1694 // add hits of the neighbour
1695  for(int iV=0; iV<ushort_v::Size; iV++)
1696  {
1697  if(!IsNeighbour[iV]) continue;
1698  for ( int jhit = 0; jhit < segment[iV]->NClusters(); jhit++ ) {
1699  int id = segment[iV]->FirstClusterRef() + jhit;
1700  hits[static_cast <unsigned short>(startHit[iV])+( dir[iV] ?( static_cast <unsigned short>(vNClusters[iV])-1-jhit ) :jhit )][iV] = id;
1701  }
1702  }
1703  vNHits(static_cast<ushort_m>(IsNeighbour)) += vNClusters;
1704 
1705  jtr = -1;
1706  IsNeighbour &= (static_cast<sfloat_v>(vSegmNextNeighbour) > -1);
1707  if(!(IsNeighbour.isEmpty()))
1708  {
1709  jtr(static_cast<ushort_m>(IsNeighbour)) = vSegmNextNeighbour;
1710  jSlice(static_cast<short_m>(IsNeighbour)) = vSegmSliceNextNeighbour;
1711  }
1712  }
1713 
1714  for(int iV=0; iV<ushort_v::Size; iV++)
1715  {
1716  if(!active[iV]) continue;
1717  if ( vEndPoint.X()[iV] < vStartPoint.X()[iV] ) { // swap
1718  for ( int i = 0; i < vNHits[iV]; i++ ) hits[i][iV] = hits[firstHit[iV]+vNHits[iV]-1-i][iV];
1719  firstHit[iV] = 0;
1720  }
1721  }
1722 
1723  const sfloat_m &invert2 = isMerged && (vEndPoint.X() < vStartPoint.X());
1724  if( !(invert2.isEmpty()))
1725  {
1726  AliHLTTPCCATrackParamVector helpPoint = vEndPoint;
1727  vEndPoint.SetTrackParam(vStartPoint, invert2);
1728  vStartPoint.SetTrackParam(helpPoint, invert2);
1729  sfloat_v helpAlpha = vEndAlpha;
1730  vEndAlpha(invert2) = vStartAlpha;
1731  vStartAlpha(invert2) = helpAlpha;
1732  }
1733 
1734  sfloat_m fitted = isMerged;
1735 // Refit tracks, which have been merged.
1736  vNHits.store(nHits);
1737  AliHLTTPCCATrackParamVector vHelpEndPoint = vStartPoint;
1738  sfloat_v vHelpEndAlpha = vStartAlpha;
1739  fitted &= FitTrack( vHelpEndPoint, vHelpEndAlpha, hits, firstHit, nHits, nVecElements,fitted, 0 );
1740  AliHLTTPCCATrackParamVector vHelpStartPoint = vHelpEndPoint;
1741  sfloat_v vHelpStartAlpha = vHelpEndAlpha;
1742  fitted &= FitTrack( vHelpStartPoint, vHelpStartAlpha, hits, firstHit, nHits, nVecElements,fitted, 1 );
1743  vNHits.load(nHits);
1744  active &= (!isMerged);
1745  active |= fitted;
1746 
1747  for(int iV=0; iV < ushort_v::Size; iV++)
1748  {
1749  if ( !active[iV] ) continue;
1750  if ( fitted[iV] ) { // merged
1751  StartPoint[iV] = AliHLTTPCCATrackParam( vHelpStartPoint, iV );
1752  EndPoint[iV] = AliHLTTPCCATrackParam( vHelpEndPoint, iV );
1753  }
1754  else { // not merged
1755  StartPoint[iV] = AliHLTTPCCATrackParam( vStartPoint, iV );
1756  EndPoint[iV] = AliHLTTPCCATrackParam( vEndPoint, iV );
1757  }
1758  }
1759  vHelpStartAlpha.store(&(StartAlpha[0]), active && fitted ); // work around for Vc-0.99.71
1760  vHelpEndAlpha.store(&(EndAlpha[0]), active && fitted );
1761  vStartAlpha.store(&(StartAlpha[0]), active && !fitted );
1762  vEndAlpha.store(&(EndAlpha[0]), active && !fitted );
1763 
1764  for(int iV=0; iV<ushort_v::Size; iV++)
1765  {
1766  if ( !active[iV] ) continue;
1767  int h[1000];
1768  for(int iClu = 0; iClu < vNHits[iV]; iClu++) h[iClu] = hits[iClu + firstHit[iV]][iV];
1769  int *usedHits = h; // get begin of array
1770 // If track has been merged, resort hits, rid of double hits.
1771  if (isMerged[iV]) {
1772  //std::sort( usedHits, usedHits + vNHits[iV], TrackHitsCompare(fClusterInfos) ); // sort hits by X (iRow) // TODO normal sort
1773  // rid of double hits
1774  int ihit2 = 0; // ihit in the output array
1775  for(int ihit = 1; ihit < vNHits[iV]; ihit++)
1776  {
1777  if( ISUNLIKELY(usedHits[ihit2] == usedHits[ihit]) ) continue;
1778  ++ihit2;
1779  if( ISUNLIKELY( ihit2 != ihit ) )
1780  usedHits[ihit2] = usedHits[ihit];
1781  }
1782  nHits[iV] = ihit2+1;
1783  }
1784 // on the final stage stor data to the global tracker
1785  if(number == 0)
1786  {
1787  AliHLTTPCCAMergedTrack &mergedTrack = outTracks[nOutTracks];
1788  mergedTrack.SetNClusters( nHits[iV] );
1789  mergedTrack.SetFirstClusterRef( nOutTrackClusters );
1790  mergedTrack.SetInnerParam( StartPoint[iV] );
1791  mergedTrack.SetInnerAlpha( StartAlpha[iV] );
1792  mergedTrack.SetOuterParam( EndPoint[iV] );
1793  mergedTrack.SetOuterAlpha( EndAlpha[iV] );
1794 
1795  for ( int i = 0; i < nHits[iV]; i++ ) {
1796  AliHLTTPCCAClusterInfo &clu = fClusterInfos[usedHits[i]];
1797  outClusterIDsrc[nOutTrackClusters+i] =
1798  DataCompressor::SliceRowCluster( clu.ISlice(), clu.IRow(), clu.IClu() );
1799  }
1800  nOutTracks++;
1801  nOutTrackClusters += nHits[iV];
1802  }
1803 // else restore tracks, obtained after merging
1804  if(number == 1)
1805  {
1806  AliHLTTPCCASliceTrackInfo &track = tmpT[nEndTracks];
1807 
1808  track = *vTrackOld[iV];
1809 
1810  track.SetFirstClusterRef(nH);
1811  track.SetNClusters(nHits[iV]);
1812  track.SetUsed(0);
1813  track.SetPrevNeighbour(-1);
1814  track.SetNextNeighbour(-1);
1815  track.SetSlicePrevNeighbour(-1);
1816  track.SetSliceNextNeighbour(-1);
1817  track.ChiPrev = 10000000;
1818  track.ChiNext = 10000000;
1819  track.SetInnerParam( StartPoint[iV] );
1820  track.SetInnerAlpha( StartAlpha[iV] );
1821  track.SetOuterParam( EndPoint[iV] );
1822  track.SetOuterAlpha( EndAlpha[iV] );
1823 
1824  track.fInnerRow = (fClusterInfos[usedHits[0]]).IRow();
1825  track.fOuterRow = (fClusterInfos[usedHits[nHits[iV]-1]]).IRow();
1826 
1827  for(int iClu=0; iClu < nHits[iV]; iClu++) tmpH[nH + iClu] = fClusterInfos[usedHits[iClu]];
1828  nH += nHits[iV];
1829  }
1830 
1831  nTrNew[iSlice]++;
1832  nEndTracks++;
1833  }
1834 
1835  nVecElements = 0;
1836  }
1837  }
1838  }
1839 
1840  if (fClusterInfos) delete[] fClusterInfos;
1841  fClusterInfos = tmpH;
1842 
1843  if (fClusterInfos) delete[] fTrackInfos;
1844  fTrackInfos = tmpT;
1845  for(int iSlice=0; iSlice < fgkNSlices; iSlice++ )
1846  {
1847  fSliceNTrackInfos[iSlice] = nTrNew[iSlice];
1848  fSliceTrackInfoStart[iSlice] = tmpSliceTrackInfoStart[iSlice];
1849  }
1850 
1851  if(number == 0)
1852  {
1853  fOutput->SetNTracks( nOutTracks );
1854  fOutput->SetNTrackClusters( nOutTrackClusters );
1855  fOutput->SetPointers();
1856 
1857  for ( int itr = 0; itr < nOutTracks; itr++ ) fOutput->SetTrack( itr, outTracks[itr] );
1858 
1859  for ( int ic = 0; ic < nOutTrackClusters; ic++ ) {
1860  fOutput->SetClusterIDsrc( ic, outClusterIDsrc[ic] );
1861  fOutput->SetClusterPackedAmp( ic, outClusterPackedAmp[ic] );
1862  }
1863 
1864  if (outTracks) delete[] outTracks;
1865  if (outClusterIDsrc) delete[] outClusterIDsrc;
1866  if (outClusterPackedAmp) delete[] outClusterPackedAmp;
1867  }
1868 
1869 #ifdef USE_TIMERS
1870  timer.Stop();
1871  fTimers[5+(1-number)] = timer.RealTime();
1872 #endif // USE_TIMERS
1873 }
void GetClusterErrors2(int iRow, const AliHLTTPCCATrackParam &t, float &Err2Y, float &Err2Z) const
mvz start 20.01.2010