HLT/TPCLib/tracking-ca/AliHLTTPCCAMerger.cxx

   1 // $Id: AliHLTTPCCAMerger.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
  21 #include "AliHLTTPCCASliceTrack.h"
  22 #include "AliHLTTPCCATracker.h"
  23 #include "AliHLTTPCCAGBTrack.h"
  24 #include "AliHLTTPCCATrackParam.h"
  25
  26 #include "AliHLTTPCCAMerger.h"
  27
  28 #include "AliHLTTPCCAMath.h"
  29 #include "TStopwatch.h"
  30
  31 #include "AliHLTTPCCATrackParam.h"
  32 #include "AliHLTTPCCASliceTrack.h"
  33 #include "AliHLTTPCCASliceOutput.h"
  34 #include "AliHLTTPCCAMergedTrack.h"
  35 #include "AliHLTTPCCAMergerOutput.h"
  36 #include "AliHLTTPCCADataCompressor.h"
  37 #include "AliHLTTPCCAParam.h"
  38
  39
  40
  41 class AliHLTTPCCAMerger::AliHLTTPCCAClusterInfo
  42 {
  43
  44 public:
  45
  46   UInt_t  ISlice()    const { return fISlice;    }
  47   UInt_t  IRow()      const { return fIRow;      }
  48   UInt_t  IClu()      const { return fIClu;      }
  49   UChar_t PackedAmp() const { return fPackedAmp; }
  50   Float_t Y()         const { return fY;         }
  51   Float_t Z()         const { return fZ;         }
  52   Float_t Err2Y()     const { return fErr2Y;     }
  53   Float_t Err2Z()     const { return fErr2Z;     }
  54
  55   void SetISlice    ( UInt_t v  ) { fISlice    = v; }
  56   void SetIRow      ( UInt_t v  ) { fIRow      = v; }
  57   void SetIClu      ( UInt_t v  ) { fIClu      = v; }
  58   void SetPackedAmp ( UChar_t v ) { fPackedAmp = v; }
  59   void SetY         ( Float_t v ) { fY         = v; }
  60   void SetZ         ( Float_t v ) { fZ         = v; }
  61   void SetErr2Y     ( Float_t v ) { fErr2Y     = v; }
  62   void SetErr2Z     ( Float_t v ) { fErr2Z     = v; }
  63
  64 private:
  65
  66   UInt_t fISlice;            // slice number
  67   UInt_t fIRow;              // row number
  68   UInt_t fIClu;              // cluster number
  69   UChar_t fPackedAmp; // packed cluster amplitude
  70   Float_t fY;                // y position (slice coord.system)
  71   Float_t fZ;                // z position (slice coord.system)
  72   Float_t fErr2Y;            // Squared measurement error of y position
  73   Float_t fErr2Z;            // Squared measurement error of z position
  74 };
  75
  76
  77 class AliHLTTPCCAMerger::AliHLTTPCCASliceTrackInfo
  78 {
  79
  80 public:
  81
  82   const AliHLTTPCCATrackParam &InnerParam() const { return fInnerParam;      }
  83   const AliHLTTPCCATrackParam &OuterParam() const { return fOuterParam;      }
  84   Float_t InnerAlpha()                      const { return fInnerAlpha;      }
  85   Float_t OuterAlpha()                      const { return fOuterAlpha;      }
  86   Int_t   NClusters()                       const { return fNClusters;       }
  87   Int_t   FirstClusterRef()                 const { return fFirstClusterRef; }
  88   Int_t   PrevNeighbour()                   const { return fPrevNeighbour;   }
  89   Int_t   NextNeighbour()                   const { return fNextNeighbour;   }
  90   Bool_t  Used()                            const { return fUsed;            }
  91
  92   void SetInnerParam( const AliHLTTPCCATrackParam &v ) { fInnerParam = v;      }
  93   void SetOuterParam( const AliHLTTPCCATrackParam &v ) { fOuterParam = v;      }
  94   void SetInnerAlpha( Float_t v )                      { fInnerAlpha = v;      }
  95   void SetOuterAlpha( Float_t v )                      { fOuterAlpha = v;      }
  96   void SetNClusters ( Int_t v )                        { fNClusters = v;       }
  97   void SetFirstClusterRef( Int_t v )                   { fFirstClusterRef = v; }
  98   void SetPrevNeighbour( Int_t v )                     { fPrevNeighbour = v;   }
  99   void SetNextNeighbour( Int_t v )                     { fNextNeighbour = v;   }
 100   void SetUsed( Bool_t v )                             { fUsed = v;            }
 101
 102 private:
 103
 104   AliHLTTPCCATrackParam fInnerParam; // inner parameters
 105   AliHLTTPCCATrackParam fOuterParam; // outer parameters
 106   Float_t fInnerAlpha;               // alpha angle for inner parameters
 107   Float_t fOuterAlpha;               // alpha angle for outer parameters
 108   Int_t fNClusters;                  // N clusters
 109   Int_t fFirstClusterRef;  // index of the first track cluster in the global cluster array
 110   Int_t fPrevNeighbour;    // neighbour in the previous slise
 111   Int_t fNextNeighbour;    // neighbour in the next slise
 112   Bool_t fUsed;            // is the slice track already merged
 113
 114 };
 115
 116
 117 class AliHLTTPCCAMerger::AliHLTTPCCABorderTrack
 118 {
 119
 120 public:
 121
 122   const AliHLTTPCCATrackParam &Param() const { return fParam;     }
 123   Int_t   TrackID()                    const { return fTrackID;   }
 124   Int_t   NClusters()                  const { return fNClusters; }
 125   Int_t   IRow()                       const { return fIRow;      }
 126   Float_t X()                          const { return fX;         }
 127   Bool_t  OK()                         const { return fOK;        }
 128
 129   void SetParam     ( const AliHLTTPCCATrackParam &v ) { fParam     = v; }
 130   void SetTrackID   ( Int_t v )                        { fTrackID   = v; }
 131   void SetNClusters ( Int_t v )                        { fNClusters = v; }
 132   void SetIRow      ( Int_t v )                        { fIRow      = v; }
 133   void SetX         ( Float_t v )                      { fX         = v; }
 134   void SetOK        ( Bool_t v )                       { fOK        = v; }
 135
 136 private:
 137
 138   AliHLTTPCCATrackParam fParam;  // track parameters at the border
 139   Int_t   fTrackID;              // track index
 140   Int_t   fNClusters;            // n clusters
 141   Int_t   fIRow;                 // row number of the closest cluster
 142   Float_t fX;                    // X coordinate of the closest cluster
 143   Bool_t  fOK;                   // is the track rotated and extrapolated correctly
 144
 145 };
 146
 147
 148
 149 AliHLTTPCCAMerger::AliHLTTPCCAMerger()
 150   :
 151   fSliceParam(),
 152   fOutput(0),
 153   fTrackInfos(0),
 154   fMaxTrackInfos(0),
 155   fClusterInfos(0),
 156   fMaxClusterInfos(0)
 157 {
 158   //* constructor
 159   Clear();
 160 }
 161
 162 /*
 163 AliHLTTPCCAMerger::AliHLTTPCCAMerger(const AliHLTTPCCAMerger&)
 164   :
 165   fSliceParam(),
 166   fkSlices(0),
 167   fOutput(0),
 168   fTrackInfos(0),
 169   fMaxTrackInfos(0),
 170   fClusterInfos(0),
 171   fMaxClusterInfos(0)
 172 {
 173 }
 174
 175 const AliHLTTPCCAMerger &AliHLTTPCCAMerger::operator=(const AliHLTTPCCAMerger&) const
 176 {
 177   return *this;
 178 }
 179 */
 180
 181 AliHLTTPCCAMerger::~AliHLTTPCCAMerger()
 182 {
 183   //* destructor
 184   if( fTrackInfos ) delete[] fTrackInfos;
 185   if( fClusterInfos ) delete[] fClusterInfos;
 186   if( fOutput ) delete[] ((char*)(fOutput));
 187 }
 188
 189 void AliHLTTPCCAMerger::Clear()
 190 {
 191   for ( int i = 0; i < fgkNSlices; ++i ) {
 192     fkSlices[i] = 0;
 193   }
 194 }
 195
 196
 197 void AliHLTTPCCAMerger::SetSliceData( int index, const AliHLTTPCCASliceOutput *SliceData )
 198 {
 199   fkSlices[index] = SliceData;
 200 }
 201
 202 void AliHLTTPCCAMerger::Reconstruct()
 203 {
 204   //* main merging routine
 205
 206   UnpackSlices();
 207   Merging();
 208 }
 209
 210 void AliHLTTPCCAMerger::UnpackSlices()
 211 {
 212   //* unpack the cluster information from the slice tracks and initialize track info array
 213
 214   // get N tracks and N clusters in event
 215
 216   Int_t nTracksTotal=0;
 217   Int_t nTrackClustersTotal=0;
 218   for( Int_t iSlice=0; iSlice<fgkNSlices; iSlice++ ){
 219     if( !fkSlices[iSlice] ) continue;
 220     nTracksTotal+=fkSlices[iSlice]->NTracks();
 221     nTrackClustersTotal+=fkSlices[iSlice]->NTrackClusters();
 222   }
 223
 224   // book/clean memory if necessary
 225   {
 226     if( nTracksTotal>fMaxTrackInfos || ( fMaxTrackInfos>100 && nTracksTotal< 0.5*fMaxTrackInfos ) ){
 227       if( fTrackInfos ) delete[] fTrackInfos;
 228       fMaxTrackInfos = (Int_t ) (nTracksTotal*1.2);
 229       fTrackInfos = new AliHLTTPCCASliceTrackInfo[fMaxTrackInfos];
 230     }
 231
 232     if( nTrackClustersTotal>fMaxClusterInfos || ( fMaxClusterInfos>1000 && nTrackClustersTotal< 0.5*fMaxClusterInfos ) ){
 233       if( fClusterInfos ) delete[] fClusterInfos;
 234       fMaxClusterInfos = (Int_t ) (nTrackClustersTotal*1.2);
 235       fClusterInfos = new AliHLTTPCCAClusterInfo [fMaxClusterInfos];
 236     }
 237
 238     if( fOutput ) delete[] ( (char*)(fOutput));
 239     Int_t size = fOutput->EstimateSize(nTracksTotal, nTrackClustersTotal);
 240     fOutput = (AliHLTTPCCAMergerOutput*)(new float2[size/sizeof(float2)+1]);
 241   }
 242
 243   // unpack track and cluster information
 244
 245   Int_t nTracksCurrent = 0;
 246   Int_t nClustersCurrent = 0;
 247
 248   for( Int_t iSlice=0; iSlice<fgkNSlices; iSlice++ ){
 249
 250     fSliceTrackInfoStart[ iSlice ] = nTracksCurrent;
 251     fSliceNTrackInfos[ iSlice ] = 0;
 252
 253     if( !fkSlices[iSlice] ) continue;
 254
 255     const AliHLTTPCCASliceOutput &slice = *(fkSlices[iSlice]);
 256
 257     for( Int_t itr=0; itr<slice.NTracks(); itr++ ){
 258
 259       const AliHLTTPCCASliceTrack &sTrack = slice.Track( itr );
 260       AliHLTTPCCATrackParam t0 = sTrack.Param();
 261       Int_t nCluNew = 0;
 262
 263       for( Int_t iTrClu=0; iTrClu<sTrack.NClusters(); iTrClu++ ){
 264
 265         // unpack cluster information
 266
 267         AliHLTTPCCAClusterInfo &clu = fClusterInfos[nClustersCurrent + nCluNew];
 268         Int_t ic = sTrack.FirstClusterRef() + iTrClu;
 269
 270         clu.SetISlice( iSlice );
 271         clu.SetIRow( AliHLTTPCCADataCompressor::IDrc2IRow( slice.ClusterIDrc( ic ) ) );
 272         clu.SetIClu( AliHLTTPCCADataCompressor::IDrc2IClu( slice.ClusterIDrc( ic ) ) );
 273         clu.SetPackedAmp( slice.ClusterPackedAmp( ic ) );
 274         float2 yz = slice.ClusterUnpackedYZ( ic );
 275         clu.SetY( yz.x );
 276         clu.SetZ( yz.y );
 277
 278         if( !t0.TransportToX( fSliceParam.RowX( clu.IRow() ), .999 ) ) continue;
 279
 280         Float_t err2Y, err2Z;
 281         fSliceParam.GetClusterErrors2( clu.IRow(), clu.Z(), t0.SinPhi(), t0.CosPhi(), t0.DzDs(), err2Y, err2Z );
 282
 283         clu.SetErr2Y( err2Y );
 284         clu.SetErr2Z( err2Z );
 285         nCluNew++ ;
 286       }
 287
 288       if( nCluNew<.8*sTrack.NClusters() ) continue;
 289
 290       // refit the track
 291
 292       Int_t hits[1000];
 293       Int_t nHits = nCluNew;
 294       for( Int_t i=0; i<nHits; i++ ) hits[i] = nClustersCurrent + i;
 295
 296       AliHLTTPCCATrackParam startPoint = sTrack.Param();
 297       AliHLTTPCCATrackParam endPoint = startPoint;
 298       Float_t startAlpha = fSliceParam.Alpha( iSlice );
 299       Float_t endAlpha = startAlpha;
 300
 301       if( !FitTrack( endPoint, endAlpha, startPoint, startAlpha, hits, nHits, 0 ) ) continue;
 302
 303       startPoint = endPoint;
 304       startAlpha = endAlpha;
 305       if( !FitTrack( startPoint, startAlpha, endPoint, endAlpha, hits, nHits, 1 ) ) continue;
 306
 307       if( nHits<.8*sTrack.NClusters() ) continue;
 308
 309       // store the track
 310
 311       AliHLTTPCCASliceTrackInfo &track = fTrackInfos[nTracksCurrent];
 312
 313       track.SetInnerParam( startPoint );
 314       track.SetInnerAlpha( startAlpha );
 315       track.SetOuterParam( endPoint );
 316       track.SetOuterAlpha( endAlpha );
 317       track.SetFirstClusterRef( nClustersCurrent );
 318       track.SetNClusters( nHits );
 319       track.SetPrevNeighbour( -1 );
 320       track.SetNextNeighbour( -1 );
 321       track.SetUsed( 0 );
 322
 323       for( Int_t i=0; i<nHits; i++ )
 324         fClusterInfos[nClustersCurrent + i] = fClusterInfos[hits[i]];
 325       nTracksCurrent++;
 326       fSliceNTrackInfos[ iSlice ]++;
 327       nClustersCurrent+=nHits;
 328     }
 329   }
 330 }
 331
 332
 333
 334 Bool_t AliHLTTPCCAMerger::FitTrack( AliHLTTPCCATrackParam &T, Float_t &Alpha,
 335                                     AliHLTTPCCATrackParam t0, Float_t Alpha0,
 336                                     Int_t hits[], Int_t &NTrackHits, Bool_t dir )
 337 {
 338   // Fit the track
 339
 340   AliHLTTPCCATrackParam::AliHLTTPCCATrackFitParam fitPar;
 341   AliHLTTPCCATrackParam t = t0;
 342
 343   Bool_t first = 1;
 344
 345   t0.CalculateFitParameters( fitPar, fSliceParam.Bz() );
 346
 347   Int_t hitsNew[1000];
 348   Int_t nHitsNew = 0;
 349
 350   for( Int_t ihit=0; ihit<NTrackHits; ihit++){
 351
 352     Int_t jhit = dir ?(NTrackHits-1-ihit) :ihit;
 353     AliHLTTPCCAClusterInfo &h = fClusterInfos[hits[jhit]];
 354
 355     Int_t iSlice = h.ISlice();
 356
 357     Float_t sliceAlpha =  fSliceParam.Alpha( iSlice );
 358
 359     if( CAMath::Abs( sliceAlpha - Alpha0)>1.e-4 ){
 360       if( ! t.RotateNoCos(  sliceAlpha - Alpha0, t0, .999 ) ) continue;
 361       Alpha0 = sliceAlpha;
 362     }
 363
 364     Float_t x = fSliceParam.RowX( h.IRow() );
 365
 366     if( !t.TransportToXWithMaterial( x, t0, fitPar ) ) continue;
 367
 368     if( first ){
 369       t.SetCov( 0, 10 );
 370       t.SetCov( 1,  0 );
 371       t.SetCov( 2, 10 );
 372       t.SetCov( 3,  0 );
 373       t.SetCov( 4,  0 );
 374       t.SetCov( 5,  1 );
 375       t.SetCov( 6,  0 );
 376       t.SetCov( 7,  0 );
 377       t.SetCov( 8,  0 );
 378       t.SetCov( 9,  1 );
 379       t.SetCov(10,  0 );
 380       t.SetCov(11,  0 );
 381       t.SetCov(12,  0 );
 382       t.SetCov(13,  0 );
 383       t.SetCov(14,  1 );
 384       t.SetChi2( 0 );
 385       t.SetNDF( -5 );
 386       t0.CalculateFitParameters( fitPar, fSliceParam.Bz() );
 387     }
 388
 389     if( !t.Filter2NoCos( h.Y(), h.Z(), h.Err2Y(), h.Err2Z() ) ) continue;
 390
 391     first = 0;
 392
 393     hitsNew[nHitsNew++] = hits[jhit];
 394   }
 395
 396   if( CAMath::Abs(t.Kappa())<1.e-8 ) t.SetKappa( 1.e-8 );
 397
 398   Bool_t ok=1;
 399
 400   const Float_t *c = t.Cov();
 401   for( Int_t i=0; i<15; i++ ) ok = ok && finite(c[i]);
 402   for( Int_t i=0; i<5; i++ ) ok = ok && finite(t.Par()[i]);
 403   ok = ok && (t.GetX()>50);
 404
 405   if( c[0]<=0 || c[2]<=0 || c[5]<=0 || c[9]<=0 || c[14]<=0 ) ok = 0;
 406   if( c[0]>5. || c[2]>5. || c[5]>2. || c[9]>2 || c[14]>2 ) ok = 0;
 407
 408   if( CAMath::Abs(t.SinPhi())>.99 ) ok = 0;
 409   else if( t0.CosPhi()>=0 ) t.SetCosPhi( CAMath::Sqrt(1.-t.SinPhi()*t.SinPhi()) );
 410   else t.SetCosPhi( -CAMath::Sqrt(1.-t.SinPhi()*t.SinPhi()) );
 411
 412   if( ok ){
 413     T = t;
 414     Alpha = Alpha0;
 415     NTrackHits = nHitsNew;
 416     for( Int_t i=0; i<NTrackHits; i++ ){
 417       hits[dir ?(NTrackHits-1-i) :i] = hitsNew[i];
 418     }
 419   }
 420   return ok;
 421 }
 422
 423
 424 Float_t AliHLTTPCCAMerger::GetChi2( Float_t x1, Float_t y1, Float_t a00, Float_t a10, Float_t a11,
 425                                     Float_t x2, Float_t y2, Float_t b00, Float_t b10, Float_t b11  )
 426 {
 427   //* Calculate Chi2/ndf deviation
 428
 429   Float_t d[2]={ x1-x2, y1-y2 };
 430
 431   Float_t mSi[3] = { a00 + b00, a10 + b10, a11 + b11 };
 432
 433   Float_t s = ( mSi[0]*mSi[2] - mSi[1]*mSi[1] );
 434
 435   if( s < 1.E-10 ) return 10000.;
 436
 437   Float_t mS[3] = { mSi[2], -mSi[1], mSi[0] };
 438
 439   return TMath::Abs( ( ( mS[0]*d[0] + mS[1]*d[1] )*d[0]
 440                        +(mS[1]*d[0] + mS[2]*d[1] )*d[1] )/s/2);
 441
 442 }
 443
 444
 445
 446 void AliHLTTPCCAMerger::MakeBorderTracks( Int_t iSlice, Int_t iBorder, AliHLTTPCCABorderTrack B[], Int_t &nB )
 447 {
 448   //* prepare slice tracks for merging with next/previous/same sector
 449   //* each track transported to the border line,
 450   //* in some cases both inner and outer parameters of the track are transported
 451
 452   static int statAll=0, statOK=0;
 453   nB = 0;
 454   Float_t dAlpha = fSliceParam.DAlpha() /2;
 455   Float_t x0 = 0;
 456
 457   if( iBorder==0 ){ // transport to the left age of the sector and rotate horisontally
 458     dAlpha = dAlpha - CAMath::Pi()/2 ;
 459   } else if( iBorder==1 ){ //  transport to the right age of the sector and rotate horisontally
 460     dAlpha = -dAlpha - CAMath::Pi()/2 ;
 461   } else if( iBorder==2 ){ // transport to the left age of the sector and rotate vertically
 462     dAlpha = dAlpha;
 463     x0 = fSliceParam.RowX( 63 );
 464   }else if( iBorder==3 ){ // transport to the right age of the sector and rotate vertically
 465     dAlpha = -dAlpha;
 466     x0 =  fSliceParam.RowX( 63 );
 467   } else if( iBorder==4 ){ // transport to the middle of the sector, w/o rotation
 468     dAlpha = 0;
 469     x0 = fSliceParam.RowX(63);
 470   }
 471
 472   for (Int_t itr=0; itr<fSliceNTrackInfos[iSlice]; itr++) {
 473
 474     const AliHLTTPCCASliceTrackInfo &track = fTrackInfos[ fSliceTrackInfoStart[iSlice] + itr ];
 475
 476     AliHLTTPCCATrackParam t0 = track.InnerParam();
 477     AliHLTTPCCATrackParam t1 = track.OuterParam();
 478
 479     const Float_t maxSin = CAMath::Sin(60./180.*CAMath::Pi());
 480
 481     Bool_t ok0 = t0.Rotate( dAlpha, maxSin );
 482     Bool_t ok1 = t1.Rotate( dAlpha, maxSin );
 483
 484     Bool_t do0 = ok0;
 485     Bool_t do1 = ok1 && ( !ok0 || t1.CosPhi()*t0.CosPhi()<0 );
 486
 487     if( ok0 && !do1 && ok1 && (t1.X() < t0.X()) ){
 488       do0 = 0;
 489       do1 = 1;
 490     }
 491
 492     if( do0 ){
 493       AliHLTTPCCABorderTrack &b = B[nB];
 494       b.SetX( t0.GetX() );
 495       if( t0.TransportToX( x0, maxSin ) ){
 496         b.SetOK( 1 );
 497         b.SetTrackID( itr );
 498         b.SetNClusters( track.NClusters() );
 499         b.SetIRow( fClusterInfos[ track.FirstClusterRef() + 0 ].IRow() );
 500         b.SetParam( t0 );
 501         nB++;
 502       }
 503     }
 504     if( do1 ){
 505       AliHLTTPCCABorderTrack &b = B[nB];
 506       b.SetX( t1.GetX() );
 507       if( t1.TransportToX( x0, maxSin ) ){
 508         b.SetOK( 1 );
 509         b.SetTrackID( itr );
 510         b.SetNClusters( track.NClusters() );
 511         b.SetIRow( fClusterInfos[ track.FirstClusterRef() + track.NClusters()-1 ].IRow() );
 512         b.SetParam( t1 );
 513         nB++;
 514       }
 515     }
 516     if( do0 || do1 ) statOK++;
 517     statAll++;
 518   }
 519 }
 520
 521
 522
 523 void AliHLTTPCCAMerger::SplitBorderTracks( Int_t iSlice1, AliHLTTPCCABorderTrack B1[], Int_t N1,
 524                                            Int_t iSlice2, AliHLTTPCCABorderTrack B2[], Int_t N2
 525                                            )
 526 {
 527   //* split two sets of tracks
 528
 529   Float_t factor2ys = 1.;//1.5;//SG!!!
 530   Float_t factor2zt = 1.;//1.5;//SG!!!
 531   Float_t factor2k = 2.0;//2.2;
 532
 533   factor2k  = 3.5*3.5*factor2k*factor2k;
 534   factor2ys = 3.5*3.5*factor2ys*factor2ys;
 535   factor2zt = 3.5*3.5*factor2zt*factor2zt;
 536
 537   Int_t minNPartHits = 10;//SG!!!
 538   Int_t minNTotalHits = 20;
 539
 540   //Float_t maxDX = fSliceParam.RowX(40) -  fSliceParam.RowX(0);
 541
 542   for (Int_t i1=0; i1<N1; i1++) {
 543     AliHLTTPCCABorderTrack &b1 = B1[i1];
 544     if( !b1.OK() ) continue;
 545     if( b1.NClusters() < minNPartHits ) continue;
 546     const AliHLTTPCCATrackParam &t1 = b1.Param();
 547     Int_t iBest2 = -1;
 548     Int_t lBest2 = 0;
 549     Int_t start2 = (iSlice1!=iSlice2) ?0 :i1+1;
 550     for (Int_t i2=start2; i2<N2; i2++) {
 551       AliHLTTPCCABorderTrack &b2 = B2[i2];
 552       if( !b2.OK() ) continue;
 553       if( b2.NClusters() < minNPartHits ) continue;
 554       if( b2.NClusters() < lBest2 ) continue;
 555       if( b1.NClusters() + b2.NClusters() < minNTotalHits ) continue;
 556
 557       //if( TMath::Abs(b1.fX - b2.fX)>maxDX ) continue;
 558
 559       const AliHLTTPCCATrackParam &t2 = b2.Param();
 560
 561       Float_t c= t2.CosPhi()*t1.CosPhi()>=0 ?1 :-1;
 562       Float_t dk = t2.Kappa() - c*t1.Kappa();
 563       Float_t s2k = t2.Err2Kappa() + t1.Err2Kappa();
 564
 565       if( dk*dk>factor2k*s2k ) continue;
 566
 567
 568       Float_t chi2ys = GetChi2( t1.Y(),c*t1.SinPhi(),t1.Cov()[0],c*t1.Cov()[3],t1.Cov()[5],
 569                                 t2.Y(),  t2.SinPhi(),t2.Cov()[0],  t2.Cov()[3],t2.Cov()[5] );
 570
 571       if( chi2ys>factor2ys ) continue;
 572
 573       Float_t chi2zt = GetChi2( t1.Z(),c*t1.DzDs(),t1.Cov()[2],c*t1.Cov()[7],t1.Cov()[9],
 574                                 t2.Z(),  t2.DzDs(),t2.Cov()[2],  t2.Cov()[7],t2.Cov()[9] );
 575
 576       if( chi2zt>factor2zt ) continue;
 577
 578       lBest2 = b2.NClusters();
 579       iBest2 = b2.TrackID();
 580     }
 581
 582     if( iBest2 <0 ) continue;
 583
 584     AliHLTTPCCASliceTrackInfo &newTrack1 = fTrackInfos[fSliceTrackInfoStart[iSlice1]+b1.TrackID() ];
 585     AliHLTTPCCASliceTrackInfo &newTrack2 = fTrackInfos[fSliceTrackInfoStart[iSlice2]+iBest2 ];
 586
 587     Int_t old1 = newTrack2.PrevNeighbour();
 588
 589     if( old1 >= 0 ){
 590       AliHLTTPCCASliceTrackInfo &oldTrack1 = fTrackInfos[fSliceTrackInfoStart[iSlice1]+old1];
 591       if( oldTrack1.NClusters()  < newTrack1.NClusters() ){
 592         newTrack2.SetPrevNeighbour(-1);
 593         oldTrack1.SetNextNeighbour(-1);
 594       } else continue;
 595     }
 596     Int_t old2 = newTrack1.NextNeighbour();
 597     if( old2 >= 0 ){
 598       AliHLTTPCCASliceTrackInfo &oldTrack2 = fTrackInfos[fSliceTrackInfoStart[iSlice2]+old2];
 599       if( oldTrack2.NClusters() < newTrack2.NClusters() ){
 600         oldTrack2.SetPrevNeighbour(-1);
 601       } else continue;
 602     }
 603     newTrack1.SetNextNeighbour(iBest2);
 604     newTrack2.SetPrevNeighbour(b1.TrackID());
 605   }
 606
 607 }
 608
 609
 610 void AliHLTTPCCAMerger::Merging()
 611 {
 612   //* track merging between slices
 613
 614   fOutput->SetNTracks( 0 );
 615   fOutput->SetNTrackClusters( 0 );
 616   fOutput->SetPointers();
 617
 618
 619   // for each slice set number of the next neighbouring slice
 620
 621   Int_t nextSlice[100], prevSlice[100];
 622
 623   for( Int_t iSlice=0; iSlice<fgkNSlices; iSlice++ ){
 624     nextSlice[iSlice] = iSlice + 1;
 625     prevSlice[iSlice] = iSlice - 1;
 626   }
 627   Int_t mid = fgkNSlices/2 - 1 ;
 628   Int_t last = fgkNSlices - 1 ;
 629   if( mid<0 ) mid = 0; // to avoid compiler warning
 630   if( last<0 ) last = 0; //
 631   nextSlice[ mid ] = 0;
 632   prevSlice[ 0 ] = mid;
 633   nextSlice[ last ] = fgkNSlices/2;
 634   prevSlice[ fgkNSlices/2 ] = last;
 635
 636   Int_t maxNSliceTracks = 0;
 637   for( Int_t iSlice=0; iSlice<fgkNSlices; iSlice++ ){
 638     if( maxNSliceTracks < fSliceNTrackInfos[iSlice] ) maxNSliceTracks = fSliceNTrackInfos[iSlice];
 639   }
 640
 641   if(1){// merging track segments withing one slice
 642
 643     AliHLTTPCCABorderTrack bord[maxNSliceTracks*10];
 644
 645     AliHLTTPCCASliceTrackInfo *tmpT = new AliHLTTPCCASliceTrackInfo[maxNSliceTracks];
 646     AliHLTTPCCAClusterInfo *tmpH = new AliHLTTPCCAClusterInfo[fMaxClusterInfos];
 647
 648     for( Int_t iSlice=0; iSlice<fgkNSlices; iSlice++ ){
 649
 650       Int_t nBord=0;
 651       MakeBorderTracks( iSlice, 4, bord, nBord );
 652       SplitBorderTracks( iSlice, bord, nBord, iSlice, bord, nBord );
 653
 654       Int_t nTr=0, nH=0;
 655       Int_t sliceFirstClusterRef = 0;
 656       for( Int_t itr=0; itr<fSliceNTrackInfos[iSlice]; itr++ ){
 657         AliHLTTPCCASliceTrackInfo &track = fTrackInfos[ fSliceTrackInfoStart[iSlice]+itr];
 658         if( itr==0 ) sliceFirstClusterRef = track.FirstClusterRef();
 659         track.SetPrevNeighbour( -1 );
 660         if( track.NextNeighbour() == -2 ){
 661           track.SetNextNeighbour( -1 );
 662           continue;
 663         }
 664         AliHLTTPCCASliceTrackInfo &trackNew = tmpT[nTr];
 665         trackNew = track;
 666         trackNew.SetFirstClusterRef( sliceFirstClusterRef + nH);
 667
 668         for( Int_t ih=0; ih<track.NClusters(); ih++ ) tmpH[nH+ih] = fClusterInfos[track.FirstClusterRef()+ih];
 669         nTr++;
 670         nH+=track.NClusters();
 671
 672         int jtr =  track.NextNeighbour();
 673
 674         if( jtr<0 ) continue;
 675         AliHLTTPCCASliceTrackInfo &neighTrack = fTrackInfos[ fSliceTrackInfoStart[iSlice]+jtr];
 676
 677         track.SetNextNeighbour(-1);
 678         neighTrack.SetNextNeighbour(-2);
 679
 680         for( Int_t ih=0; ih<neighTrack.NClusters(); ih++ )
 681           tmpH[nH+ih] = fClusterInfos[neighTrack.FirstClusterRef()+ih];
 682
 683         trackNew.SetNClusters( trackNew.NClusters() + neighTrack.NClusters() );
 684         trackNew.SetNextNeighbour( -1 );
 685         nH+=neighTrack.NClusters();
 686         if( neighTrack.InnerParam().X() < track.InnerParam().X() ) trackNew.SetInnerParam( neighTrack.InnerParam());
 687         if( neighTrack.OuterParam().X() > track.OuterParam().X() ) trackNew.SetOuterParam( neighTrack.OuterParam());
 688       }
 689
 690       fSliceNTrackInfos[iSlice] = nTr;
 691       for( Int_t itr=0; itr<nTr; itr++ ) fTrackInfos[ fSliceTrackInfoStart[iSlice]+itr] = tmpT[itr];
 692       for( Int_t ih=0; ih<nH; ih++ ) fClusterInfos[sliceFirstClusterRef + ih] = tmpH[ih];
 693
 694     }
 695     delete[] tmpT;
 696     delete[] tmpH;
 697   }
 698
 699
 700   //* merging tracks between slices
 701
 702
 703   // arrays for the rotated track parameters
 704
 705   AliHLTTPCCABorderTrack
 706     *bCurr0 = new AliHLTTPCCABorderTrack[maxNSliceTracks*10],
 707     *bNext0 = new AliHLTTPCCABorderTrack[maxNSliceTracks*10],
 708     *bCurr = new AliHLTTPCCABorderTrack[maxNSliceTracks*10],
 709     *bNext = new AliHLTTPCCABorderTrack[maxNSliceTracks*10];
 710
 711   for( Int_t iSlice=0; iSlice<fgkNSlices; iSlice++ ){
 712
 713     Int_t jSlice = nextSlice[iSlice];
 714
 715     Int_t nCurr0 = 0, nNext0 = 0;
 716     Int_t nCurr = 0, nNext = 0;
 717
 718     MakeBorderTracks( iSlice, 0, bCurr, nCurr );
 719     MakeBorderTracks( jSlice, 1, bNext, nNext );
 720     MakeBorderTracks( iSlice, 2, bCurr0, nCurr0 );
 721     MakeBorderTracks( jSlice, 3, bNext0, nNext0 );
 722
 723     SplitBorderTracks( iSlice, bCurr0, nCurr0, jSlice, bNext0, nNext0 );
 724     SplitBorderTracks( iSlice, bCurr, nCurr, jSlice, bNext, nNext );
 725   }
 726
 727   if( bCurr0 ) delete[] bCurr0;
 728   if( bNext0 ) delete[] bNext0;
 729   if( bCurr  ) delete[] bCurr;
 730   if( bNext  ) delete[] bNext;
 731
 732
 733   //static Int_t nRejected = 0;
 734
 735   Int_t nOutTracks = 0;
 736   Int_t nOutTrackClusters = 0;
 737
 738   AliHLTTPCCAMergedTrack *outTracks = new AliHLTTPCCAMergedTrack[fMaxTrackInfos];
 739   UInt_t   *outClusterIDsrc = new UInt_t [fMaxClusterInfos];
 740   UChar_t  *outClusterPackedAmp = new UChar_t [fMaxClusterInfos];
 741
 742   for( Int_t iSlice = 0; iSlice<fgkNSlices; iSlice++ ){
 743
 744     for( Int_t itr=0; itr<fSliceNTrackInfos[iSlice]; itr++ ){
 745
 746       AliHLTTPCCASliceTrackInfo &track = fTrackInfos[fSliceTrackInfoStart[iSlice]+itr];
 747
 748       if( track.Used() ) continue;
 749       if( track.PrevNeighbour()>=0 ) continue;
 750
 751       AliHLTTPCCATrackParam startPoint = track.InnerParam(), endPoint = track.OuterParam();
 752       Float_t startAlpha = track.InnerAlpha(), endAlpha = track.OuterAlpha();
 753
 754       Int_t hits[2000];
 755       Int_t firstHit = 1000;
 756       Int_t nHits = 0;
 757       Int_t jSlice = iSlice;
 758       Int_t jtr = itr;
 759
 760       {
 761         track.SetUsed( 1 );
 762         for( Int_t jhit=0; jhit<track.NClusters(); jhit++){
 763           Int_t id = track.FirstClusterRef() + jhit;
 764           hits[firstHit+jhit] = id;
 765         }
 766         nHits=track.NClusters();
 767         jtr = track.NextNeighbour();
 768         jSlice = nextSlice[iSlice];
 769       }
 770
 771       while( jtr >=0 ){
 772         AliHLTTPCCASliceTrackInfo &segment = fTrackInfos[fSliceTrackInfoStart[jSlice]+jtr];
 773         if( segment.Used() ) break;
 774         segment.SetUsed( 1 );
 775         Bool_t dir = 0;
 776         Int_t startHit = firstHit+ nHits;
 777         Float_t d00 = startPoint.GetDistXZ2(segment.InnerParam() );
 778         Float_t d01 = startPoint.GetDistXZ2(segment.OuterParam() );
 779         Float_t d10 = endPoint.GetDistXZ2(segment.InnerParam() );
 780         Float_t d11 = endPoint.GetDistXZ2(segment.OuterParam() );
 781         if( d00<=d01 && d00<=d10 && d00<=d11 ){
 782           startPoint = segment.OuterParam();
 783           startAlpha = segment.OuterAlpha();
 784           dir = 1;
 785           firstHit -= segment.NClusters();
 786           startHit = firstHit;
 787         }else if( d01<=d10 && d01<=d11 ){
 788           startPoint = segment.InnerParam();
 789           startAlpha = segment.InnerAlpha();
 790           dir = 0;
 791           firstHit -= segment.NClusters();
 792           startHit = firstHit;
 793         }else if( d10<=d11 ){
 794           endPoint = segment.OuterParam();
 795           endAlpha = segment.OuterAlpha();
 796           dir = 0;
 797         }else{
 798           endPoint = segment.InnerParam();
 799           endAlpha = segment.InnerAlpha();
 800           dir = 1;
 801         }
 802
 803         for( Int_t jhit=0; jhit<segment.NClusters(); jhit++){
 804           Int_t id = segment.FirstClusterRef() + jhit;
 805           hits[startHit+(dir ?(segment.NClusters()-1-jhit) :jhit)] = id;
 806         }
 807         nHits+=segment.NClusters();
 808         jtr = segment.NextNeighbour();
 809         jSlice = nextSlice[jSlice];
 810       }
 811
 812       if( endPoint.X() < startPoint.X() ){ // swap
 813         for( Int_t i=0; i<nHits; i++ ) hits[i] = hits[firstHit+nHits-1-i];
 814         firstHit = 0;
 815       }
 816
 817       if( nHits < 30 ) continue;     //SG!!!
 818
 819       // refit
 820
 821       // need best t0!!!SG
 822
 823       endPoint = startPoint;
 824       if( !FitTrack( endPoint, endAlpha, startPoint, startAlpha, hits+firstHit, nHits, 0 ) ) continue;
 825       if( !FitTrack( startPoint, startAlpha, endPoint, endAlpha, hits+firstHit, nHits, 1 ) ) continue;
 826
 827       if( nHits < 30 ) continue;     //SG!!!
 828
 829       AliHLTTPCCATrackParam &p = startPoint;
 830       AliHLTTPCCATrackParam::AliHLTTPCCATrackFitParam fitPar;
 831       p.CalculateFitParameters( fitPar, fSliceParam.Bz() );
 832
 833       {
 834         Double_t xTPC=83.65; //SG!!!
 835         Double_t dAlpha = 0.00609235;
 836
 837         if( p.TransportToXWithMaterial( xTPC, fitPar ) ){
 838           Double_t y=p.GetY();
 839           Double_t ymax=xTPC*CAMath::Tan(dAlpha/2.);
 840           if (y > ymax) {
 841             if( p.Rotate( dAlpha ) ){ startAlpha+=dAlpha;  p.TransportToXWithMaterial( xTPC, fitPar ); }
 842           } else if (y <-ymax) {
 843             if( p.Rotate( -dAlpha ) ){  startAlpha-=dAlpha; p.TransportToXWithMaterial( xTPC, fitPar );}
 844           }
 845         }
 846       }
 847
 848       {
 849         Bool_t ok=1;
 850
 851         const Float_t *c = p.Cov();
 852         for( Int_t i=0; i<15; i++ ) ok = ok && finite(c[i]);
 853         for( Int_t i=0; i<5; i++ ) ok = ok && finite(p.Par()[i]);
 854         ok = ok && (p.GetX()>50);
 855
 856         if( c[0]<=0 || c[2]<=0 || c[5]<=0 || c[9]<=0 || c[14]<=0 ) ok = 0;
 857         if( c[0]>5. || c[2]>5. || c[5]>2. || c[9]>2 || c[14]>2 ) ok = 0;
 858         if(!ok) continue;
 859       }
 860
 861       AliHLTTPCCAMergedTrack &mergedTrack = outTracks[nOutTracks];
 862       mergedTrack.SetNClusters(nHits);
 863       mergedTrack.SetFirstClusterRef(nOutTrackClusters);
 864       mergedTrack.SetInnerParam(startPoint);
 865       mergedTrack.SetInnerAlpha(startAlpha);
 866       mergedTrack.SetOuterParam(endPoint);
 867       mergedTrack.SetOuterAlpha(endAlpha);
 868
 869       for( Int_t i = 0; i<nHits; i++ ){
 870         AliHLTTPCCAClusterInfo &clu = fClusterInfos[hits[firstHit+i]];
 871         outClusterIDsrc[nOutTrackClusters+i] =
 872           AliHLTTPCCADataCompressor::ISliceIRowIClu2IDsrc(clu.ISlice(), clu.IRow(), clu.IClu());
 873         outClusterPackedAmp[nOutTrackClusters+i]=clu.PackedAmp();
 874       }
 875
 876       nOutTracks++;
 877       nOutTrackClusters+= nHits;
 878     }
 879   }
 880
 881   fOutput->SetNTracks( nOutTracks );
 882   fOutput->SetNTrackClusters( nOutTrackClusters );
 883   fOutput->SetPointers();
 884
 885   for( Int_t itr=0; itr<nOutTracks; itr++ ) fOutput->SetTrack( itr, outTracks[itr] );
 886
 887   for( Int_t ic=0; ic<nOutTrackClusters; ic++ ){
 888     fOutput->SetClusterIDsrc( ic, outClusterIDsrc[ic] );
 889     fOutput->SetClusterPackedAmp( ic, outClusterPackedAmp[ic] );
 890   }
 891
 892   delete[] outTracks;
 893   delete[] outClusterIDsrc;
 894   delete[] outClusterPackedAmp;
 895 }