HLT/TPCLib/tracking-ca/AliHLTTPCCANeighboursFinder.cxx

   1 // @(#) $Id: AliHLTTPCCANeighboursFinder1.cxx 27042 2008-07-02 12:06:02Z richterm $
   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 "AliHLTTPCCANeighboursFinder.h"
  21 #include "AliHLTTPCCAMath.h"
  22 #include "AliHLTTPCCAHitArea.h"
  23 #include "AliHLTTPCCATracker.h"
  24
  25 //#define DRAW
  26
  27 #ifdef DRAW
  28 #include "AliHLTTPCCADisplay.h"
  29 #endif
  30
  31 GPUd() void AliHLTTPCCANeighboursFinder::Thread
  32 ( int /*nBlocks*/, int nThreads, int iBlock, int iThread, int iSync,
  33   AliHLTTPCCASharedMemory &s, AliHLTTPCCATracker &tracker )
  34 {
  35   //* find neighbours
  36
  37   if ( iSync == 0 ) {
  38     if ( iThread == 0 ) {
  39       s.fNRows = tracker.Param().NRows();
  40       s.fIRow = iBlock;
  41       if ( s.fIRow < s.fNRows ) {
  42         const AliHLTTPCCARow &row = tracker.Row( s.fIRow );
  43         s.fFirst = row.FirstHit();
  44         s.fNHits = row.NHits();
  45
  46         s.fHitLinkUp = ( ( short* )( tracker.RowData() + row.FullOffset() ) ) + row.FullLinkOffset();
  47         s.fHitLinkDn = s.fHitLinkUp + row.NHits();
  48
  49         if ( ( s.fIRow >= 2 ) && ( s.fIRow <= s.fNRows - 3 ) ) {
  50           s.fIRowUp = s.fIRow + 2;
  51           s.fIRowDn = s.fIRow - 2;
  52           s.fFirstDn = tracker.Row( s.fIRowDn ).FirstHit();
  53           s.fFirstUp = tracker.Row( s.fIRowUp ).FirstHit();
  54           float xDn = tracker.Row( s.fIRowDn ).X();
  55           float x = tracker.Row( s.fIRow ).X();
  56           float xUp = tracker.Row( s.fIRowUp ).X();
  57           s.fUpNHits = tracker.Row( s.fIRowUp ).NHits();
  58           s.fDnNHits = s.fFirst - s.fFirstDn;
  59           s.fUpDx = xUp - x;
  60           s.fDnDx = xDn - x;
  61           s.fUpTx = xUp / x;
  62           s.fDnTx = xDn / x;
  63           s.fGridUp = tracker.Row( s.fIRowUp ).Grid();
  64           s.fGridDn = tracker.Row( s.fIRowDn ).Grid();
  65         }
  66       }
  67     }
  68   } else if ( iSync == 1 ) {
  69     if ( s.fIRow < s.fNRows ) {
  70       if ( ( s.fIRow == 0 ) || ( s.fIRow == s.fNRows - 1 ) || ( s.fIRow == 1 ) || ( s.fIRow == s.fNRows - 2 ) ) {
  71         for ( int ih = iThread; ih < s.fNHits; ih += nThreads ) {
  72           s.fHitLinkUp[ih] = -1;
  73           s.fHitLinkDn[ih] = -1;
  74         }
  75       } else {
  76         const AliHLTTPCCARow &rowUp = tracker.Row( s.fIRowUp );
  77         const AliHLTTPCCARow &rowDn = tracker.Row( s.fIRowDn );
  78         const unsigned short *gContentUp = ( reinterpret_cast<const unsigned short*>( tracker.RowData() + rowUp.FullOffset() ) ) + rowUp.FullGridOffset();
  79         const unsigned short *gContentDn = ( reinterpret_cast<const unsigned short*>( tracker.RowData() + rowDn.FullOffset() ) ) + rowDn.FullGridOffset();
  80
  81         for ( unsigned int ih = iThread; ih < s.fGridUp.N() + s.fGridUp.Ny() + 2; ih += nThreads ) {
  82           s.fGridContentUp[ih] = gContentUp[ih];
  83         }
  84         for ( unsigned int ih = iThread; ih < s.fGridDn.N() + s.fGridDn.Ny() + 2; ih += nThreads ) {
  85           s.fGridContentDn[ih] = gContentDn[ih];
  86         }
  87       }
  88     }
  89   } else if ( iSync == 2 ) {
  90     if ( ( s.fIRow <= 1 ) || ( s.fIRow >= s.fNRows - 2 ) ) return;
  91
  92     //const float kAreaSize = 3;
  93     float chi2Cut = 3.*3.*4 * ( s.fUpDx * s.fUpDx + s.fDnDx * s.fDnDx );
  94     const float kAreaSize = 3;
  95     //float chi2Cut = 3.*3.*(s.fUpDx*s.fUpDx + s.fDnDx*s.fDnDx ); //SG
  96     const int kMaxN = 20;
  97
  98     const AliHLTTPCCARow &row = tracker.Row( s.fIRow );
  99     const AliHLTTPCCARow &rowUp = tracker.Row( s.fIRowUp );
 100     const AliHLTTPCCARow &rowDn = tracker.Row( s.fIRowDn );
 101     float y0 = row.Grid().YMin();
 102     float z0 = row.Grid().ZMin();
 103     float stepY = row.HstepY();
 104     float stepZ = row.HstepZ();
 105     const uint4* tmpint4 = tracker.RowData() + row.FullOffset();
 106     const ushort2 *hits = reinterpret_cast<const ushort2*>( tmpint4 );
 107
 108     for ( int ih = iThread; ih < s.fNHits; ih += nThreads ) {
 109
 110       unsigned short *neighUp = s.fB[iThread];
 111       float2 *yzUp = s.fA[iThread];
 112       //unsigned short neighUp[5];
 113       //float2 yzUp[5];
 114
 115       int linkUp = -1;
 116       int linkDn = -1;
 117
 118       if ( s.fDnNHits > 0 && s.fUpNHits > 0 ) {
 119
 120         int nNeighUp = 0;
 121         AliHLTTPCCAHit h0;
 122         {
 123           ushort2 hh = hits[ih];
 124           h0.SetY( y0 + hh.x*stepY );
 125           h0.SetZ( z0 + hh.y*stepZ );
 126         }
 127         //h0 = tracker.Hits()[ s.fFirst + ih ];
 128
 129         float y = h0.Y();
 130         float z = h0.Z();
 131
 132         AliHLTTPCCAHitArea areaDn, areaUp;
 133         areaUp.Init( s.fGridUp,  s.fGridContentUp, s.fFirstUp, y*s.fUpTx, z*s.fUpTx, kAreaSize, kAreaSize );
 134         areaDn.Init( s.fGridDn,  s.fGridContentDn, s.fFirstDn, y*s.fDnTx, z*s.fDnTx, kAreaSize, kAreaSize );
 135
 136         do {
 137           AliHLTTPCCAHit h;
 138           int i = areaUp.GetNext( tracker, rowUp, s.fGridContentUp, h );
 139           if ( i < 0 ) break;
 140           neighUp[nNeighUp] = ( unsigned short ) i;
 141           yzUp[nNeighUp] = CAMath::MakeFloat2( s.fDnDx * ( h.Y() - y ), s.fDnDx * ( h.Z() - z ) );
 142           if ( ++nNeighUp >= kMaxN ) break;
 143         } while ( 1 );
 144
 145         int nNeighDn = 0;
 146
 147         if ( nNeighUp > 0 ) {
 148
 149           int bestDn = -1, bestUp = -1;
 150           float bestD = 1.e10;
 151
 152           do {
 153             AliHLTTPCCAHit h;
 154             int i = areaDn.GetNext( tracker, rowDn, s.fGridContentDn, h );
 155             if ( i < 0 ) break;
 156
 157             nNeighDn++;
 158             float2 yzdn = CAMath::MakeFloat2( s.fUpDx * ( h.Y() - y ), s.fUpDx * ( h.Z() - z ) );
 159
 160             for ( int iUp = 0; iUp < nNeighUp; iUp++ ) {
 161               float2 yzup = yzUp[iUp];
 162               float dy = yzdn.x - yzup.x;
 163               float dz = yzdn.y - yzup.y;
 164               float d = dy * dy + dz * dz;
 165               if ( d < bestD ) {
 166                 bestD = d;
 167                 bestDn = i;
 168                 bestUp = iUp;
 169               }
 170             }
 171           } while ( 1 );
 172
 173           if ( bestD <= chi2Cut ) {
 174             linkUp = neighUp[bestUp];
 175             linkDn = bestDn;
 176           }
 177         }
 178 #ifdef DRAW
 179         std::cout << "n NeighUp = " << nNeighUp << ", n NeighDn = " << nNeighDn << std::endl;
 180 #endif
 181
 182       }
 183
 184       s.fHitLinkUp[ih] = linkUp;
 185       s.fHitLinkDn[ih] = linkDn;
 186 #ifdef DRAW
 187       std::cout << "Links for row " << s.fIRow << ", hit " << ih << ": " << linkUp << " " << linkDn << std::endl;
 188       if ( s.fIRow == 22 && ih == 5 ) {
 189         AliHLTTPCCADisplay::Instance().DrawSliceLink( s.fIRow, ih, -1, -1, 1 );
 190         AliHLTTPCCADisplay::Instance().DrawSliceHit( s.fIRow, ih, kBlue, 1. );
 191         AliHLTTPCCADisplay::Instance().Ask();
 192       }
 193 #endif
 194     }
 195   }
 196 }
 197