HLT/TPCLib/tracking-ca/AliHLTTPCCATracker.cxx

   1 // @(#) $Id$
   2 // **************************************************************************
   3 // This file is property of and copyright by the ALICE HLT Project          *
   4 // ALICE Experiment at CERN, All rights reserved.                           *
   5 //                                                                          *
   6 // Primary Authors: Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
   7 //                  Ivan Kisel <kisel@kip.uni-heidelberg.de>                *
   8 //                  for The ALICE HLT Project.                              *
   9 //                                                                          *
  10 // Permission to use, copy, modify and distribute this software and its     *
  11 // documentation strictly for non-commercial purposes is hereby granted     *
  12 // without fee, provided that the above copyright notice appears in all     *
  13 // copies and that both the copyright notice and this permission notice     *
  14 // appear in the supporting documentation. The authors make no claims       *
  15 // about the suitability of this software for any purpose. It is            *
  16 // provided "as is" without express or implied warranty.                    *
  17 //                                                                          *
  18 //***************************************************************************
  19
  20 #include "AliHLTTPCCATracker.h"
  21 #include "AliHLTTPCCAOutTrack.h"
  22 #include "AliHLTTPCCARow.h"
  23 #include "AliHLTTPCCATrack.h"
  24 #include "AliHLTTPCCATracklet.h"
  25 #include "AliHLTTPCCAMath.h"
  26 #include "MemoryAssignmentHelpers.h"
  27
  28 #include "TStopwatch.h"
  29 #include "AliHLTTPCCAHitArea.h"
  30 #include "AliHLTTPCCANeighboursFinder.h"
  31 #include "AliHLTTPCCANeighboursCleaner.h"
  32 #include "AliHLTTPCCAStartHitsFinder.h"
  33 #include "AliHLTTPCCATrackletConstructor.h"
  34 #include "AliHLTTPCCATrackletSelector.h"
  35 #include "AliHLTTPCCAProcess.h"
  36 #include "AliHLTTPCCASliceTrack.h"
  37 #include "AliHLTTPCCASliceOutput.h"
  38 #include "AliHLTTPCCAClusterData.h"
  39
  40 #include "AliHLTTPCCATrackParam.h"
  41
  42 #include "AliHLTTPCCAGPUConfig.h"
  43
  44 #if !defined(HLTCA_GPUCODE)
  45 #include <iostream>
  46 #include <iomanip>
  47 #include <string.h>
  48 #endif
  49
  50 //#define DRAW1
  51
  52 #ifdef DRAW1
  53 #include "AliHLTTPCCADisplay.h"
  54 #endif //DRAW
  55
  56 #ifdef HLTCA_INTERNAL_PERFORMANCE
  57 //#include "AliHLTTPCCAPerformance.h"
  58 #endif
  59
  60 #ifdef HLTCA_STANDALONE
  61 #include "AliHLTTPCCAStandaloneFramework.h"
  62 #endif
  63
  64 ClassImp( AliHLTTPCCATracker )
  65
  66
  67 #if !defined(HLTCA_GPUCODE)
  68
  69 AliHLTTPCCATracker::~AliHLTTPCCATracker()
  70 {
  71   // destructor
  72         if (!fIsGPUTracker)
  73         {
  74                 if (fCommonMem) delete fCommonMem;
  75                 if (fHitMemory) delete[] fHitMemory;
  76                 if (fTrackletMemory) delete[] fTrackletMemory;
  77                 if (fTrackMemory) delete[] fTrackMemory;
  78                 fCommonMem = NULL;
  79                 fHitMemory = fTrackMemory = NULL;
  80         }
  81 }
  82
  83 // ----------------------------------------------------------------------------------
  84 void AliHLTTPCCATracker::Initialize( const AliHLTTPCCAParam &param )
  85 {
  86   // initialisation
  87   fParam = param;
  88   fParam.Update();
  89   fData.InitializeRows( fParam );
  90
  91   StartEvent();
  92 }
  93
  94 void AliHLTTPCCATracker::StartEvent()
  95 {
  96   // start new event and fresh the memory
  97
  98   SetupCommonMemory();
  99 }
 100
 101 void AliHLTTPCCATracker::SetGPUTracker()
 102 {
 103         fIsGPUTracker = true;
 104         fData.SetGpuSliceData();
 105 }
 106
 107 char* AliHLTTPCCATracker::SetGPUTrackerCommonMemory(char* pGPUMemory)
 108 {
 109   fCommonMem = (commonMemoryStruct*) pGPUMemory;
 110   return(pGPUMemory + sizeof(commonMemoryStruct));
 111 }
 112
 113
 114 char* AliHLTTPCCATracker::SetGPUTrackerHitsMemory(char* pGPUMemory, int MaxNHits)
 115 {
 116         fHitMemory = (char*) pGPUMemory;
 117         SetPointersHits(MaxNHits);
 118         pGPUMemory += fHitMemorySize;
 119         AssignMemory(fTrackletTmpStartHits, pGPUMemory, NHitsTotal());
 120         AssignMemory(fRowStartHitCountOffset, pGPUMemory, Param().NRows());
 121
 122         return(pGPUMemory);
 123 }
 124
 125 char* AliHLTTPCCATracker::SetGPUTrackerTrackletsMemory(char* pGPUMemory, int MaxNTracks)
 126 {
 127         fTrackletMemory = (char*) pGPUMemory;
 128         SetPointersTracklets(MaxNTracks);
 129         pGPUMemory += fTrackletMemorySize;
 130         AssignMemory(fGPUTrackletTemp, pGPUMemory, MaxNTracks);
 131         AssignMemory(fRowBlockTracklets, pGPUMemory, MaxNTracks * 2 * (Param().NRows() / HLTCA_GPU_SCHED_ROW_STEP + 1));
 132         AssignMemory(fRowBlockPos, pGPUMemory, 2 * (Param().NRows() / HLTCA_GPU_SCHED_ROW_STEP + 1));
 133         AssignMemory(fBlockStartingTracklet, pGPUMemory, HLTCA_GPU_BLOCK_COUNT);
 134
 135         return(pGPUMemory);
 136 }
 137
 138 char* AliHLTTPCCATracker::SetGPUTrackerTracksMemory(char* pGPUMemory, int MaxNTracks, int MaxNHits )
 139 {
 140         fTrackMemory = (char*) pGPUMemory;
 141         SetPointersTracks(MaxNTracks, MaxNHits);
 142         pGPUMemory += fTrackMemorySize;
 143
 144         return(pGPUMemory);
 145 }
 146
 147 void AliHLTTPCCATracker::DumpSliceData(std::ostream &out)
 148 {
 149         out << "Slice Data:" << endl;
 150         for (int i = 0;i < Param().NRows();i++)
 151         {
 152                 out << "Row: " << i << endl;
 153                 for (int j = 0;j < Row(i).NHits();j++)
 154                 {
 155                         if (j && j % 16 == 0) out << std::endl;
 156                         out << j << '-' << Data().HitDataY(Row(i), j) << '-' << Data().HitDataZ(Row(i), j) << ", ";
 157                 }
 158                 out << endl;
 159         }
 160 }
 161
 162 void AliHLTTPCCATracker::DumpLinks(std::ostream &out)
 163 {
 164         out << "Hit Links:" << endl;
 165         for (int i = 0;i < Param().NRows();i++)
 166         {
 167                 out << "Row: " << i << endl;
 168                 for (int j = 0;j < Row(i).NHits();j++)
 169                 {
 170                         if (j && j % 32 == 0) out << endl;
 171                         out << HitLinkUpData(Row(i), j) << "/" << HitLinkDownData(Row(i), j) << ", ";
 172                 }
 173                 out << endl;
 174         }
 175 }
 176
 177 void AliHLTTPCCATracker::DumpHitWeights(std::ostream &out)
 178 {
 179     out << "Hit Weights:" << endl;
 180     for (int i = 0;i < Param().NRows();i++)
 181     {
 182         out << "Row: " << i << ":" << endl;
 183         for (int j = 0;j < Row(i).NHits();j++)
 184         {
 185             if (j && j % 32 == 0) out << endl;
 186             out << HitWeight(Row(i), j) << ", ";
 187         }
 188         out << endl;
 189     }
 190 }
 191
 192 int AliHLTTPCCATracker::StarthitSortComparison(const void*a, const void* b)
 193 {
 194         AliHLTTPCCAHitId* aa = (AliHLTTPCCAHitId*) a;
 195         AliHLTTPCCAHitId* bb = (AliHLTTPCCAHitId*) b;
 196
 197         if (aa->RowIndex() != bb->RowIndex()) return(aa->RowIndex() - bb->RowIndex());
 198         return(aa->HitIndex() - bb->HitIndex());
 199 }
 200
 201 void AliHLTTPCCATracker::DumpStartHits(std::ostream &out)
 202 {
 203         out << "Start Hits: (" << *NTracklets() << ")" << endl;
 204 #ifdef HLTCA_GPU_SORT_DUMPDATA
 205         qsort(TrackletStartHits(), *NTracklets(), sizeof(AliHLTTPCCAHitId), StarthitSortComparison);
 206 #endif
 207         for (int i = 0;i < *NTracklets();i++)
 208         {
 209                 out << TrackletStartHit(i).RowIndex() << "-" << TrackletStartHit(i).HitIndex() << endl;
 210         }
 211         out << endl;
 212 }
 213
 214 void AliHLTTPCCATracker::DumpTrackHits(std::ostream &out)
 215 {
 216         out << "Tracks: (" << *NTracks() << ")" << endl;
 217 #ifdef HLTCA_GPU_SORT_DUMPDATA
 218         for (int k = 0;k < Param().NRows();k++)
 219         {
 220                 for (int l = 0;l < Row(k).NHits();l++)
 221                 {
 222 #endif
 223                         for (int j = 0;j < *NTracks();j++)
 224                         {
 225                                 if (Tracks()[j].NHits() == 0 || !Tracks()[j].Alive()) continue;
 226 #ifdef HLTCA_GPU_SORT_DUMPDATA
 227                                 if (TrackHits()[Tracks()[j].FirstHitID()].RowIndex() == k && TrackHits()[Tracks()[j].FirstHitID()].HitIndex() == l)
 228                                 {
 229 #endif
 230                                         for (int i = 0;i < Tracks()[j].NHits();i++)
 231                                         {
 232                                                 out << TrackHits()[Tracks()[j].FirstHitID() + i].RowIndex() << "-" << TrackHits()[Tracks()[j].FirstHitID() + i].HitIndex() << ", ";
 233                                         }
 234                                         out << "(Track: " << j << ")" << endl;
 235 #ifdef HLTCA_GPU_SORT_DUMPDATA
 236                                 }
 237                         }
 238 #endif
 239                 }
 240 #ifdef HLTCA_GPU_SORT_DUMPDATA
 241         }
 242 #endif
 243 }
 244
 245 void AliHLTTPCCATracker::DumpTrackletHits(std::ostream &out)
 246 {
 247         out << "Tracklets: (" << *NTracklets() << ")" << endl;
 248 #ifdef HLTCA_GPU_SORT_DUMPDATA
 249         AliHLTTPCCAHitId* tmpIds = new AliHLTTPCCAHitId[*NTracklets()];
 250         AliHLTTPCCATracklet* tmpTracklets = new AliHLTTPCCATracklet[*NTracklets()];
 251         memcpy(tmpIds, TrackletStartHits(), *NTracklets() * sizeof(AliHLTTPCCAHitId));
 252         memcpy(tmpTracklets, Tracklets(), *NTracklets() * sizeof(AliHLTTPCCATracklet));
 253 #ifdef EXTERN_ROW_HITS
 254         int* tmpHits = new int[*NTracklets() * Param().NRows()];
 255         memcpy(tmpHits, TrackletRowHits(), *NTracklets() * Param().NRows() * sizeof(int));
 256 #endif
 257         qsort(TrackletStartHits(), *NTracklets(), sizeof(AliHLTTPCCAHitId), StarthitSortComparison);
 258         for (int i = 0;i < *NTracklets();i++)
 259         {
 260                 for (int j = 0;j < *NTracklets();j++)
 261                 {
 262                         if (tmpIds[i].RowIndex() == TrackletStartHit(j).RowIndex() && tmpIds[i].HitIndex() == TrackletStartHit(j).HitIndex())
 263                         {
 264                                 memcpy(&Tracklets()[j], &tmpTracklets[i], sizeof(AliHLTTPCCATracklet));
 265 #ifdef EXTERN_ROW_HITS
 266                                 if (tmpTracklets[i].NHits())
 267                                 {
 268                                         for (int k = tmpTracklets[i].FirstRow();k <= tmpTracklets[i].LastRow();k++)
 269                                         {
 270                                                 fTrackletRowHits[k * *NTracklets() + j] = tmpHits[k * *NTracklets() + i];
 271                                         }
 272                                 }
 273 #endif
 274                                 break;
 275                         }
 276                 }
 277         }
 278         delete[] tmpIds;
 279         delete[] tmpTracklets;
 280 #ifdef EXTERN_ROW_HITS
 281         delete[] tmpHits;
 282 #endif
 283 #endif
 284         for (int j = 0;j < *NTracklets();j++)
 285         {
 286                 out << "Tracklet " << j << " (Hits: " << setw(3) << Tracklets()[j].NHits() << ", Start: " << setw(3) << TrackletStartHit(j).RowIndex() << "-" << setw(3) << TrackletStartHit(j).HitIndex() << ") ";
 287                 if (Tracklets()[j].NHits() == 0);
 288                 else if (Tracklets()[j].LastRow() > Tracklets()[j].FirstRow() && (Tracklets()[j].FirstRow() >= Param().NRows() || Tracklets()[j].LastRow() >= Param().NRows()))
 289                 {
 290 #ifdef HLTCA_STANDALONE
 291                         printf("\nError: First %d Last %d Num %d", Tracklets()[j].FirstRow(), Tracklets()[j].LastRow(), Tracklets()[j].NHits());
 292 #endif
 293                 }
 294                 else if (Tracklets()[j].NHits() && Tracklets()[j].LastRow() > Tracklets()[j].FirstRow())
 295                 {
 296                         for (int i = Tracklets()[j].FirstRow();i <= Tracklets()[j].LastRow();i++)
 297                         {
 298                                 //if (Tracklets()[j].RowHit(i) != -1)
 299 #ifdef EXTERN_ROW_HITS
 300                                         out << i << "-" << fTrackletRowHits[i * fCommonMem->fNTracklets + j] << ", ";
 301 #else
 302                                         out << i << "-" << Tracklets()[j].RowHit(i) << ", ";
 303 #endif
 304                         }
 305                 }
 306                 out << endl;
 307         }
 308 }
 309
 310
 311 void AliHLTTPCCATracker::SetupCommonMemory()
 312 {
 313   // set up common memory
 314
 315   if (!fIsGPUTracker)
 316   {
 317     if ( !fCommonMem ) {
 318       // the 1600 extra bytes are not used unless fCommonMemorySize increases with a later event
 319       //fCommonMemory = reinterpret_cast<char*> ( new uint4 [ fCommonMemorySize/sizeof( uint4 ) + 100] );
 320           fCommonMem = new commonMemoryStruct;
 321     }
 322
 323         if (fHitMemory) delete[] fHitMemory;
 324         if (fTrackletMemory) delete[] fTrackletMemory;
 325         if (fTrackMemory) delete[] fTrackMemory;
 326   }
 327
 328   fHitMemory = fTrackletMemory = fTrackMemory = 0;
 329
 330   fData.Clear();
 331   fCommonMem->fNTracklets = 0;
 332   fCommonMem->fNTracks = 0 ;
 333   fCommonMem->fNTrackHits = 0;
 334 }
 335
 336 void AliHLTTPCCATracker::ReadEvent( AliHLTTPCCAClusterData *clusterData )
 337 {
 338   // read event
 339
 340   fClusterData = clusterData;
 341
 342   StartEvent();
 343
 344   //* Convert input hits, create grids, etc.
 345   fData.InitFromClusterData( *clusterData );
 346   {
 347     if (!fIsGPUTracker)
 348         {
 349                 SetPointersHits( fData.NumberOfHits() ); // to calculate the size
 350                 fHitMemory = reinterpret_cast<char*> ( new uint4 [ fHitMemorySize/sizeof( uint4 ) + 100] );
 351         }
 352     SetPointersHits( fData.NumberOfHits() ); // set pointers for hits
 353   }
 354 }
 355
 356 GPUhd() void  AliHLTTPCCATracker::SetPointersHits( int MaxNHits )
 357 {
 358   // set all pointers to the event memory
 359
 360   char *mem = fHitMemory;
 361
 362   // extra arrays for tpc clusters
 363
 364 #ifdef HLTCA_GPU_SORT_STARTHITS_2
 365   AssignMemory( fTrackletStartHits, mem, MaxNHits + 32);
 366 #else
 367   AssignMemory( fTrackletStartHits, mem, MaxNHits);
 368 #endif
 369
 370   // calculate the size
 371
 372   fHitMemorySize = mem - fHitMemory;
 373 }
 374
 375 GPUhd() void  AliHLTTPCCATracker::SetPointersTracklets( int MaxNTracklets )
 376 {
 377   // set all pointers to the tracklets memory
 378   char *mem = fTrackletMemory;
 379
 380   // memory for tracklets
 381
 382   AssignMemory( fTracklets, mem, MaxNTracklets );
 383 #ifdef EXTERN_ROW_HITS
 384   AssignMemory( fTrackletRowHits, mem, MaxNTracklets * Param().NRows());
 385 #endif
 386
 387   fTrackletMemorySize = mem - fTrackletMemory;
 388 }
 389
 390
 391 GPUhd() void  AliHLTTPCCATracker::SetPointersTracks( int MaxNTracks, int MaxNHits )
 392 {
 393   // set all pointers to the tracks memory
 394   char *mem = fTrackMemory;
 395
 396   // memory for selected tracks
 397
 398   AssignMemory( fTracks, mem, MaxNTracks );
 399   AssignMemory( fTrackHits, mem, 2 * MaxNHits );
 400
 401   // calculate the size
 402
 403   fTrackMemorySize = mem - fTrackMemory;
 404 }
 405
 406 GPUh() int AliHLTTPCCATracker::CheckEmptySlice() const
 407 {
 408   if ( NHitsTotal() < 1 ) {
 409     {
 410       AliHLTTPCCASliceOutput::Allocate(*fOutput, 0, 0, fOutputControl);
 411       AliHLTTPCCASliceOutput* useOutput = *fOutput;
 412       if (fOutput == NULL) return(1);
 413       useOutput->SetNTracks( 0 );
 414       useOutput->SetNTrackClusters( 0 );
 415       useOutput->SetNOutTracks(0);
 416       useOutput->SetNOutTrackHits(0);
 417     }
 418
 419     return 1;
 420   }
 421   return 0;
 422 }
 423
 424 void AliHLTTPCCATracker::RunNeighboursFinder()
 425 {
 426         AliHLTTPCCAProcess<AliHLTTPCCANeighboursFinder>( Param().NRows(), 1, *this );
 427 }
 428
 429 void AliHLTTPCCATracker::RunNeighboursCleaner()
 430 {
 431         AliHLTTPCCAProcess<AliHLTTPCCANeighboursCleaner>( Param().NRows() - 2, 1, *this );
 432 }
 433
 434 void AliHLTTPCCATracker::RunStartHitsFinder()
 435 {
 436         AliHLTTPCCAProcess<AliHLTTPCCAStartHitsFinder>( Param().NRows() - 4, 1, *this );
 437 }
 438
 439 void AliHLTTPCCATracker::RunTrackletConstructor()
 440 {
 441   AliHLTTPCCATrackletConstructor::AliHLTTPCCATrackletConstructorNewCPU(*this);
 442 }
 443
 444 void AliHLTTPCCATracker::RunTrackletSelector()
 445 {
 446   AliHLTTPCCAProcess<AliHLTTPCCATrackletSelector>( 1, fCommonMem->fNTracklets, *this );
 447 }
 448
 449 #ifdef HLTCA_STANDALONE
 450 void AliHLTTPCCATracker::StandalonePerfTime(int i)
 451 {
 452   if (fGPUDebugLevel >= 1)
 453   {
 454     AliHLTTPCCAStandaloneFramework::StandaloneQueryTime(&fPerfTimers[i]);
 455   }
 456 }
 457 #else
 458 void AliHLTTPCCATracker::StandalonePerfTime(int /*i*/) {}
 459 #endif
 460
 461 GPUh() void AliHLTTPCCATracker::Reconstruct()
 462 {
 463   //* reconstruction of event
 464   //std::cout<<"Reconstruct slice "<<fParam.ISlice()<<", nHits="<<NHitsTotal()<<std::endl;
 465
 466   fTimers[0] = 0; // find neighbours
 467   fTimers[1] = 0; // construct tracklets
 468   fTimers[2] = 0; // fit tracklets
 469   fTimers[3] = 0; // prolongation of tracklets
 470   fTimers[4] = 0; // selection
 471   fTimers[5] = 0; // write output
 472   fTimers[6] = 0;
 473   fTimers[7] = 0;
 474
 475   //if( fParam.ISlice()<1 ) return; //SG!!!
 476
 477   TStopwatch timer0;
 478
 479   StandalonePerfTime(0);
 480
 481   if (CheckEmptySlice()) return;
 482
 483 #ifdef DRAW1
 484   //if( fParam.ISlice()==2 || fParam.ISlice()==3)
 485     {
 486   AliHLTTPCCADisplay::Instance().ClearView();
 487   AliHLTTPCCADisplay::Instance().SetSliceView();
 488   AliHLTTPCCADisplay::Instance().SetCurrentSlice( this );
 489   AliHLTTPCCADisplay::Instance().DrawSlice( this, 1 );
 490   if ( NHitsTotal() > 0 ) {
 491     AliHLTTPCCADisplay::Instance().DrawSliceHits( kRed, .5 );
 492     AliHLTTPCCADisplay::Instance().Ask();
 493   }
 494   }
 495 #endif
 496
 497         fCommonMem->fNTracklets = fCommonMem->fNTracks = fCommonMem->fNTrackHits = 0;
 498
 499 #if !defined(HLTCA_GPUCODE)
 500
 501   /*if (fGPUDebugLevel >= 6)
 502   {
 503           *fGPUDebugOut << endl << endl << "Slice: " << Param().ISlice() << endl;
 504           *fGPUDebugOut << "Slice Data:" << endl;
 505           DumpSliceData(*fGPUDebugOut);
 506   }*/
 507
 508   StandalonePerfTime(1);
 509
 510   RunNeighboursFinder();
 511
 512   StandalonePerfTime(2);
 513
 514   if (fGPUDebugLevel >= 6) DumpLinks(*fGPUDebugOut);
 515
 516 #ifdef HLTCA_INTERNAL_PERFORMANCE
 517   //if( Param().ISlice()<=2 )
 518   //AliHLTTPCCAPerformance::Instance().LinkPerformance( Param().ISlice() );
 519 #endif
 520
 521
 522 #ifdef DRAW1
 523   if ( NHitsTotal() > 0 ) {
 524     AliHLTTPCCADisplay::Instance().DrawSliceLinks( -1, -1, 1 );
 525     AliHLTTPCCADisplay::Instance().Ask();
 526   }
 527 #endif
 528
 529   RunNeighboursCleaner();
 530
 531   StandalonePerfTime(3);
 532
 533   if (fGPUDebugLevel >= 6) DumpLinks(*fGPUDebugOut);
 534
 535   RunStartHitsFinder();
 536
 537   StandalonePerfTime(4);
 538   StandalonePerfTime(5);
 539
 540   if (fGPUDebugLevel >= 6) DumpStartHits(*fGPUDebugOut);
 541
 542   fData.ClearHitWeights();
 543
 544   SetPointersTracklets( fCommonMem->fNTracklets * 2 ); // to calculate the size
 545   fTrackletMemory = reinterpret_cast<char*> ( new uint4 [ fTrackletMemorySize/sizeof( uint4 ) + 100] );
 546   SetPointersTracklets( fCommonMem->fNTracklets * 2 ); // set pointers for hits
 547
 548   SetPointersTracks( fCommonMem->fNTracklets * 2, NHitsTotal() ); // to calculate the size
 549   fTrackMemory = reinterpret_cast<char*> ( new uint4 [ fTrackMemorySize/sizeof( uint4 ) + 100] );
 550   SetPointersTracks( fCommonMem->fNTracklets * 2, NHitsTotal() ); // set pointers for hits
 551
 552   StandalonePerfTime(6);
 553   StandalonePerfTime(7);
 554
 555   RunTrackletConstructor();
 556
 557   StandalonePerfTime(8);
 558
 559   if (fGPUDebugLevel >= 6) DumpTrackletHits(*fGPUDebugOut);
 560   if (fGPUDebugLevel >= 6) DumpHitWeights(*fGPUDebugOut);
 561
 562   //std::cout<<"Slice "<<Param().ISlice()<<": NHits="<<NHitsTotal()<<", NTracklets="<<*NTracklets()<<std::endl;
 563
 564   RunTrackletSelector();
 565
 566   StandalonePerfTime(9);
 567
 568   //std::cout<<"Slice "<<Param().ISlice()<<": N start hits/tracklets/tracks = "<<nStartHits<<" "<<nStartHits<<" "<<*fNTracks<<std::endl;
 569
 570   if (fGPUDebugLevel >= 6) DumpTrackHits(*fGPUDebugOut);
 571
 572   //std::cout<<"Memory used for slice "<<fParam.ISlice()<<" : "<<fCommonMemorySize/1024./1024.<<" + "<<fHitMemorySize/1024./1024.<<" + "<<fTrackMemorySize/1024./1024.<<" = "<<( fCommonMemorySize+fHitMemorySize+fTrackMemorySize )/1024./1024.<<" Mb "<<std::endl;
 573
 574   WriteOutput();
 575
 576   StandalonePerfTime(10);
 577
 578 #endif
 579
 580 #ifdef DRAW1
 581   {
 582     AliHLTTPCCADisplay &disp = AliHLTTPCCADisplay::Instance();
 583     AliHLTTPCCATracker &slice = *this;
 584     std::cout << "N out tracks = " << slice.NOutTracks() << std::endl;
 585     AliHLTTPCCADisplay::Instance().SetSliceView();
 586     AliHLTTPCCADisplay::Instance().SetCurrentSlice( this );
 587     AliHLTTPCCADisplay::Instance().DrawSlice( this, 1 );
 588     disp.DrawSliceHits( kRed, .5 );
 589     disp.Ask();
 590     for ( int itr = 0; itr < slice.NOutTracks(); itr++ ) {
 591       std::cout << "track N " << itr << ", nhits=" << slice.OutTracks()[itr].NHits() << std::endl;
 592       disp.DrawSliceOutTrack( itr, kBlue );
 593       //disp.Ask();
 594       //int id = slice.OutTracks()[itr].OrigTrackID();
 595       //AliHLTTPCCATrack &tr = Tracks()[id];
 596       //for( int ih=0; ih<tr.NHits(); ih++ ){
 597       //int ic = (fTrackHits[tr.FirstHitID()+ih]);
 598       //std::cout<<ih<<" "<<ID2IRow(ic)<<" "<<ID2IHit(ic)<<std::endl;
 599       //}
 600       //disp.DrawSliceTrack( id, kBlue );
 601       //disp.Ask();
 602     }
 603     disp.Ask();
 604   }
 605 #endif
 606
 607   timer0.Stop();
 608   fTimers[0] = timer0.CpuTime() / 100.;
 609
 610 }
 611
 612 GPUh() void AliHLTTPCCATracker::WriteOutput()
 613 {
 614   // write output
 615
 616   TStopwatch timer;
 617
 618   //cout<<"output: nTracks = "<<*fNTracks<<", nHitsTotal="<<NHitsTotal()<<std::endl;
 619
 620   if (fOutputControl == NULL) fOutputControl = new AliHLTTPCCASliceOutput::outputControlStruct;
 621   AliHLTTPCCASliceOutput::Allocate(*fOutput, fCommonMem->fNTracks, fCommonMem->fNTrackHits, fOutputControl);
 622   AliHLTTPCCASliceOutput* useOutput = *fOutput;
 623   if (useOutput == NULL) return;
 624
 625   if (fOutputControl->fDefaultOutput)
 626   {
 627           useOutput->SetNTracks( fCommonMem->fNTracks );
 628           useOutput->SetNTrackClusters( fCommonMem->fNTrackHits );
 629
 630           int nStoredHits = 0;
 631
 632           for ( int iTr = 0; iTr < fCommonMem->fNTracks; iTr++ ) {
 633                 AliHLTTPCCATrack &iTrack = fTracks[iTr];
 634
 635                 AliHLTTPCCASliceTrack out;
 636                 out.SetFirstClusterRef( nStoredHits );
 637                 out.SetNClusters( iTrack.NHits() );
 638                 out.SetParam( iTrack.Param() );
 639
 640                 useOutput->SetTrack( iTr, out );
 641
 642                 int iID = iTrack.FirstHitID();
 643                 for ( int ith = 0; ith < iTrack.NHits(); ith++ ) {
 644                   const AliHLTTPCCAHitId &ic = fTrackHits[iID + ith];
 645                   int iRow = ic.RowIndex();
 646                   int ih = ic.HitIndex();
 647
 648                   const AliHLTTPCCARow &row = fData.Row( iRow );
 649
 650                   //float y0 = row.Grid().YMin();
 651                   //float z0 = row.Grid().ZMin();
 652                   //float stepY = row.HstepY();
 653                   //float stepZ = row.HstepZ();
 654                   //float x = row.X();
 655
 656                   //const uint4 *tmpint4 = RowData() + row.FullOffset();
 657                   //const ushort2 *hits = reinterpret_cast<const ushort2*>(tmpint4);
 658                   //ushort2 hh = hits[ih];
 659                   //float y = y0 + hh.x*stepY;
 660                   //float z = z0 + hh.y*stepZ;
 661
 662                   int clusterIndex = fData.ClusterDataIndex( row, ih );
 663                   int clusterRowIndex = clusterIndex - fClusterData->RowOffset( iRow );
 664
 665                   if ( clusterIndex < 0 || clusterIndex >= fClusterData->NumberOfClusters() ) {
 666                         //std::cout << inpIDtot << ", " << fClusterData->NumberOfClusters()
 667                         //<< "; " << inpID << ", " << fClusterData->NumberOfClusters( iRow ) << std::endl;
 668                         //abort();
 669                         continue;
 670                   }
 671                   if ( clusterRowIndex < 0 || clusterRowIndex >= fClusterData->NumberOfClusters( iRow ) ) {
 672                         //std::cout << inpIDtot << ", " << fClusterData->NumberOfClusters()
 673                         //<< "; " << inpID << ", " << fClusterData->NumberOfClusters( iRow ) << std::endl;
 674                         //abort();
 675                         continue;
 676                   }
 677
 678                   float origX = fClusterData->X( clusterIndex );
 679                   float origY = fClusterData->Y( clusterIndex );
 680                   float origZ = fClusterData->Z( clusterIndex );
 681
 682
 683                   int id = fClusterData->Id( clusterIndex );
 684
 685                   unsigned short hPackedYZ = 0;
 686                   UChar_t hPackedAmp = 0;
 687                   float2 hUnpackedYZ;
 688                   hUnpackedYZ.x = origY;
 689                   hUnpackedYZ.y = origZ;
 690                   float hUnpackedX = origX;
 691
 692                   useOutput->SetClusterId( nStoredHits, id  );
 693                   useOutput->SetClusterRow( nStoredHits, ( unsigned char ) iRow  );
 694                   useOutput->SetClusterPackedYZ( nStoredHits, hPackedYZ );
 695                   useOutput->SetClusterPackedAmp( nStoredHits, hPackedAmp );
 696                   useOutput->SetClusterUnpackedYZ( nStoredHits, hUnpackedYZ );
 697                   useOutput->SetClusterUnpackedX( nStoredHits, hUnpackedX );
 698                   nStoredHits++;
 699                 }
 700           }
 701   }
 702
 703
 704   // old stuff
 705   if (fOutputControl->fObsoleteOutput)
 706   {
 707           useOutput->SetNOutTrackHits(0);
 708           useOutput->SetNOutTracks(0);
 709
 710
 711           for ( int iTr = 0; iTr < fCommonMem->fNTracks; iTr++ ) {
 712
 713                 const AliHLTTPCCATrack &iTrack = fTracks[iTr];
 714
 715                 //std::cout<<"iTr = "<<iTr<<", nHits="<<iTrack.NHits()<<std::endl;
 716
 717                 //if( !iTrack.Alive() ) continue;
 718                 if ( iTrack.NHits() < 3 ) continue;
 719                 AliHLTTPCCAOutTrack &out = useOutput->OutTracks()[useOutput->NOutTracks()];
 720                 out.SetFirstHitRef( useOutput->NOutTrackHits() );
 721                 out.SetNHits( 0 );
 722                 out.SetOrigTrackID( iTr );
 723                 out.SetStartPoint( iTrack.Param() );
 724                 out.SetEndPoint( iTrack.Param() );
 725
 726                 int iID = iTrack.FirstHitID();
 727                 int nOutTrackHitsOld = useOutput->NOutTrackHits();
 728
 729                 for ( int ith = 0; ith < iTrack.NHits(); ith++ ) {
 730                   const AliHLTTPCCAHitId &ic = fTrackHits[iID + ith];
 731                   const AliHLTTPCCARow &row = Row( ic );
 732                   int ih = ic.HitIndex();
 733                   useOutput->SetOutTrackHit(useOutput->NOutTrackHits(), HitInputID( row, ih ));
 734                   useOutput->SetNOutTrackHits(useOutput->NOutTrackHits() + 1 );
 735                   //std::cout<<"write i,row,hit,id="<<ith<<", "<<ID2IRow(ic)<<", "<<ih<<", "<<HitInputID( row, ih )<<std::endl;
 736                   if ( useOutput->NOutTrackHits() >= 10*NHitsTotal() ) {
 737                         std::cout << "fNOutTrackHits>NHitsTotal()" << std::endl;
 738                         //exit(0);
 739                         return;//SG!!!
 740                   }
 741                   out.SetNHits( out.NHits() + 1 );
 742                 }
 743                 if ( out.NHits() >= 2 ) {
 744                   useOutput->SetNOutTracks(useOutput->NOutTracks() + 1);
 745                 } else {
 746                   useOutput->SetNOutTrackHits(nOutTrackHitsOld);
 747                 }
 748           }
 749   }
 750
 751   timer.Stop();
 752   fTimers[5] += timer.CpuTime();
 753 }
 754
 755 #endif
 756
 757 GPUh() void AliHLTTPCCATracker::FitTrackFull( const AliHLTTPCCATrack &/**/, float * /**/ ) const
 758 {
 759   // fit track with material
 760 #ifdef XXX
 761   //* Fit the track
 762   FitTrack( iTrack, tt0 );
 763   if ( iTrack.NHits() <= 3 ) return;
 764
 765   AliHLTTPCCATrackParam &t = iTrack.Param();
 766   AliHLTTPCCATrackParam t0 = t;
 767
 768   t.Chi2() = 0;
 769   t.NDF() = -5;
 770   bool first = 1;
 771
 772   int iID = iTrack.FirstHitID();
 773   for ( int ih = 0; ih < iTrack.NHits(); ih++, iID++ ) {
 774     const AliHLTTPCCAHitId &ic = fTrackHits[iID];
 775     int iRow = ic.rowIndex();
 776     const AliHLTTPCCARow &row = fData.Row( iRow );
 777     if ( !t0.TransportToX( row.X() ) ) continue;
 778     float dy, dz;
 779     const AliHLTTPCCAHit &h = ic.hitIndex();
 780
 781     // check for wrong hits
 782     if ( 0 ) {
 783       dy = t0.GetY() - h.Y();
 784       dz = t0.GetZ() - h.Z();
 785
 786       //if( dy*dy > 3.5*3.5*(/*t0.GetErr2Y() + */h.ErrY()*h.ErrY() ) ) continue;//SG!!!
 787       //if( dz*dz > 3.5*3.5*(/*t0.GetErr2Z() + */h.ErrZ()*h.ErrZ() ) ) continue;
 788     }
 789
 790     if ( !t.TransportToX( row.X() ) ) continue;
 791
 792     //* Update the track
 793
 794     if ( first ) {
 795       t.Cov()[ 0] = .5 * .5;
 796       t.Cov()[ 1] = 0;
 797       t.Cov()[ 2] = .5 * .5;
 798       t.Cov()[ 3] = 0;
 799       t.Cov()[ 4] = 0;
 800       t.Cov()[ 5] = .2 * .2;
 801       t.Cov()[ 6] = 0;
 802       t.Cov()[ 7] = 0;
 803       t.Cov()[ 8] = 0;
 804       t.Cov()[ 9] = .2 * .2;
 805       t.Cov()[10] = 0;
 806       t.Cov()[11] = 0;
 807       t.Cov()[12] = 0;
 808       t.Cov()[13] = 0;
 809       t.Cov()[14] = .2 * .2;
 810       t.Chi2() = 0;
 811       t.NDF() = -5;
 812     }
 813     float err2Y, err2Z;
 814     GetErrors2( iRow, t, err2Y, err2Z );
 815
 816     if ( !t.Filter2( h.Y(), h.Z(), err2Y, err2Z ) ) continue;
 817
 818     first = 0;
 819   }
 820   /*
 821   float cosPhi = iTrack.Param().GetCosPhi();
 822   p0.Param().TransportToX(ID2Row( iTrack.PointID()[0] ).X());
 823   p2.Param().TransportToX(ID2Row( iTrack.PointID()[1] ).X());
 824   if( p0.Param().GetCosPhi()*cosPhi<0 ){ // change direction
 825   float *par = p0.Param().Par();
 826   float *cov = p0.Param().Cov();
 827   par[2] = -par[2]; // sin phi
 828   par[3] = -par[3]; // DzDs
 829     par[4] = -par[4]; // kappa
 830     cov[3] = -cov[3];
 831     cov[4] = -cov[4];
 832     cov[6] = -cov[6];
 833     cov[7] = -cov[7];
 834     cov[10] = -cov[10];
 835     cov[11] = -cov[11];
 836     p0.Param().CosPhi() = -p0.Param().GetCosPhi();
 837   }
 838   */
 839 #endif
 840 }
 841
 842 GPUh() void AliHLTTPCCATracker::FitTrack( const AliHLTTPCCATrack &/*track*/, float * /*t0[]*/ ) const
 843 {
 844   //* Fit the track
 845 #ifdef XXX
 846   AliHLTTPCCAEndPoint &p2 = ID2Point( track.PointID()[1] );
 847   const AliHLTTPCCAHit &c0 = ID2Hit( fTrackHits[p0.TrackHitID()].HitID() );
 848   const AliHLTTPCCAHit &c1 = ID2Hit( fTrackHits[track.HitID()[1]].HitID() );
 849   const AliHLTTPCCAHit &c2 = ID2Hit( fTrackHits[p2.TrackHitID()].HitID() );
 850   const AliHLTTPCCARow &row0 = ID2Row( fTrackHits[p0.TrackHitID()].HitID() );
 851   const AliHLTTPCCARow &row1 = ID2Row( fTrackHits[track.HitID()[1]].HitID() );
 852   const AliHLTTPCCARow &row2 = ID2Row( fTrackHits[p2.TrackHitID()].HitID() );
 853   float sp0[5] = {row0.X(), c0.Y(), c0.Z(), c0.ErrY(), c0.ErrZ() };
 854   float sp1[5] = {row1.X(), c1.Y(), c1.Z(), c1.ErrY(), c1.ErrZ() };
 855   float sp2[5] = {row2.X(), c2.Y(), c2.Z(), c2.ErrY(), c2.ErrZ() };
 856   //std::cout<<"Fit track, points ="<<sp0[0]<<" "<<sp0[1]<<" / "<<sp1[0]<<" "<<sp1[1]<<" / "<<sp2[0]<<" "<<sp2[1]<<std::endl;
 857   if ( track.NHits() >= 3 ) {
 858     p0.Param().ConstructXYZ3( sp0, sp1, sp2, p0.Param().CosPhi(), t0 );
 859     p2.Param().ConstructXYZ3( sp2, sp1, sp0, p2.Param().CosPhi(), t0 );
 860     //p2.Param() = p0.Param();
 861     //p2.Param().TransportToX(row2.X());
 862     //p2.Param().Par()[1] = -p2.Param().Par()[1];
 863     //p2.Param().Par()[4] = -p2.Param().Par()[4];
 864   } else {
 865     p0.Param().X() = row0.X();
 866     p0.Param().Y() = c0.Y();
 867     p0.Param().Z() = c0.Z();
 868     p0.Param().Err2Y() = c0.ErrY() * c0.ErrY();
 869     p0.Param().Err2Z() = c0.ErrZ() * c0.ErrZ();
 870     p2.Param().X() = row2.X();
 871     p2.Param().Y() = c2.Y();
 872     p2.Param().Z() = c2.Z();
 873     p2.Param().Err2Y() = c2.ErrY() * c2.ErrY();
 874     p2.Param().Err2Z() = c2.ErrZ() * c2.ErrZ();
 875   }
 876 #endif
 877 }
 878
 879
 880 GPUd() void AliHLTTPCCATracker::GetErrors2( int iRow, float z, float sinPhi, float cosPhi, float DzDs, float &Err2Y, float &Err2Z ) const
 881 {
 882   //
 883   // Use calibrated cluster error from OCDB
 884   //
 885
 886   fParam.GetClusterErrors2( iRow, z, sinPhi, cosPhi, DzDs, Err2Y, Err2Z );
 887   Err2Y*=fParam.ClusterError2CorrectionY();
 888   Err2Z*=fParam.ClusterError2CorrectionZ();
 889 }
 890
 891 GPUd() void AliHLTTPCCATracker::GetErrors2( int iRow, const AliHLTTPCCATrackParam &t, float &Err2Y, float &Err2Z ) const
 892 {
 893   //
 894   // Use calibrated cluster error from OCDB
 895   //
 896
 897   fParam.GetClusterErrors2( iRow, t.GetZ(), t.SinPhi(), t.GetCosPhi(), t.DzDs(), Err2Y, Err2Z );
 898 }
 899
 900
 901 #if !defined(HLTCA_GPUCODE)
 902
 903 GPUh() void AliHLTTPCCATracker::WriteEvent( std::ostream &out )
 904 {
 905   // write event to the file
 906   for ( int iRow = 0; iRow < fParam.NRows(); iRow++ ) {
 907     out << fData.Row( iRow ).HitNumberOffset() << " " << fData.Row( iRow ).NHits() << std::endl;
 908   }
 909   out << NHitsTotal() << std::endl;
 910
 911   AliHLTResizableArray<float> y( NHitsTotal() ), z( NHitsTotal() );
 912
 913   for ( int iRow = 0; iRow < fParam.NRows(); iRow++ ) {
 914     const AliHLTTPCCARow &row = Row( iRow );
 915     float y0 = row.Grid().YMin();
 916     float z0 = row.Grid().ZMin();
 917     float stepY = row.HstepY();
 918     float stepZ = row.HstepZ();
 919     for ( int ih = 0; ih < fData.Row( iRow ).NHits(); ih++ ) {
 920       int id = HitInputID( row, ih );
 921       y[id] = y0 + HitDataY( row, ih ) * stepY;
 922       z[id] = z0 + HitDataZ( row, ih ) * stepZ;
 923     }
 924   }
 925   for ( int ih = 0; ih < NHitsTotal(); ih++ ) {
 926     out << y[ih] << " " << z[ih] << std::endl;
 927   }
 928 }
 929
 930 GPUh() void AliHLTTPCCATracker::WriteTracks( std::ostream &out )
 931 {
 932   //* Write tracks to file
 933   AliHLTTPCCASliceOutput* useOutput = *fOutput;
 934
 935   out << fTimers[0] << std::endl;
 936   out << useOutput->NOutTrackHits() << std::endl;
 937   for ( int ih = 0; ih < useOutput->NOutTrackHits(); ih++ ) {
 938     out << useOutput->OutTrackHit(ih) << " ";
 939   }
 940   out << std::endl;
 941
 942   out << useOutput->NOutTracks() << std::endl;
 943
 944   for ( int itr = 0; itr < useOutput->NOutTracks(); itr++ ) {
 945     const AliHLTTPCCAOutTrack &t = useOutput->OutTrack(itr);
 946     AliHLTTPCCATrackParam p1 = t.StartPoint();
 947     AliHLTTPCCATrackParam p2 = t.EndPoint();
 948     out << t.NHits() << " ";
 949     out << t.FirstHitRef() << " ";
 950     out << t.OrigTrackID() << " ";
 951     out << std::endl;
 952     out << p1.X() << " ";
 953     out << p1.SignCosPhi() << " ";
 954     out << p1.Chi2() << " ";
 955     out << p1.NDF() << std::endl;
 956     for ( int i = 0; i < 5; i++ ) out << p1.Par()[i] << " ";
 957     out << std::endl;
 958     for ( int i = 0; i < 15; i++ ) out << p1.Cov()[i] << " ";
 959     out << std::endl;
 960     out << p2.X() << " ";
 961     out << p2.SignCosPhi() << " ";
 962     out << p2.Chi2() << " ";
 963     out << p2.NDF() << std::endl;
 964     for ( int i = 0; i < 5; i++ ) out << p2.Par()[i] << " ";
 965     out << std::endl;
 966     for ( int i = 0; i < 15; i++ ) out << p2.Cov()[i] << " ";
 967     out << std::endl;
 968   }
 969 }
 970
 971 GPUh() void AliHLTTPCCATracker::ReadTracks( std::istream &in )
 972 {
 973   //* Read tracks  from file
 974   AliHLTTPCCASliceOutput::Allocate(*fOutput, 4096, 16384, fOutputControl);//Just some max values
 975   AliHLTTPCCASliceOutput* useOutput = *fOutput;
 976
 977   int tmpval;
 978   in >> fTimers[0];
 979   in >> tmpval;
 980   useOutput->SetNOutTrackHits(tmpval);
 981
 982   for ( int ih = 0; ih < useOutput->NOutTrackHits(); ih++ ) {
 983     in >> tmpval;
 984         useOutput->SetOutTrackHit(ih, tmpval);
 985   }
 986   in >> tmpval;
 987   useOutput->SetNOutTracks(tmpval);
 988
 989   for ( int itr = 0; itr < useOutput->NOutTracks(); itr++ ) {
 990     AliHLTTPCCAOutTrack &t = useOutput->OutTracks()[itr];
 991     AliHLTTPCCATrackParam p1, p2;
 992     int i;
 993     float f;
 994     in >> i; t.SetNHits( i );
 995     in >> i; t.SetFirstHitRef( i );
 996     in >> i; t.SetOrigTrackID( i );
 997     in >> f; p1.SetX( f );
 998     in >> f; p1.SetSignCosPhi( f );
 999     in >> f; p1.SetChi2( f );
1000     in >> i; p1.SetNDF( i );
1001     for ( int j = 0; j < 5; j++ ) { in >> f; p1.SetPar( j, f ); }
1002     for ( int j = 0; j < 15; j++ ) { in >> f; p1.SetCov( j, f ); }
1003     in >> f; p2.SetX( f );
1004     in >> f; p2.SetSignCosPhi( f );
1005     in >> f; p2.SetChi2( f );
1006     in >> i; p2.SetNDF( i );
1007     for ( int j = 0; j < 5; j++ ) { in >> f; p2.SetPar( j, f ); }
1008     for ( int j = 0; j < 15; j++ ) { in >> f; p2.SetCov( j, f ); }
1009     t.SetStartPoint( p1 );
1010     t.SetEndPoint( p2 );
1011   }
1012 }
1013 #endif