[u/mrichter/AliRoot.git] / HLT / TPCLib / tracking-ca / AliHLTTPCCASliceData.cxx

// **************************************************************************
// * This file is property of and copyright by the ALICE HLT Project        *
// * All rights reserved.                                                   *
// *                                                                        *
// * Primary Authors:                                                       *
// *     Copyright 2009       Matthias Kretz <kretz@kde.org>                *
// *                                                                        *
// * Permission to use, copy, modify and distribute this software and its   *
// * documentation strictly for non-commercial purposes is hereby granted   *
// * without fee, provided that the above copyright notice appears in all   *
// * copies and that both the copyright notice and this permission notice   *
// * appear in the supporting documentation. The authors make no claims     *
// * about the suitability of this software for any purpose. It is          *
// * provided "as is" without express or implied warranty.                  *
// **************************************************************************

#include "AliHLTTPCCASliceData.h"
#include "AliHLTTPCCAClusterData.h"
#include "AliHLTTPCCAMath.h"
#include "AliHLTArray.h"
#include "AliHLTTPCCAHit.h"
#include "AliHLTTPCCAParam.h"
#include "MemoryAssignmentHelpers.h"
#include "AliHLTTPCCAGPUConfig.h"
#include "AliHLTTPCCAGPUTracker.h"
#include <iostream>

// calculates an approximation for 1/sqrt(x)
// Google for 0x5f3759df :)
static inline float fastInvSqrt( float _x )
{
  // the function calculates fast inverse sqrt

  union { float f; int i; } x = { _x };
  const float xhalf = 0.5f * x.f;
  x.i = 0x5f3759df - ( x.i >> 1 );
  x.f = x.f * ( 1.5f - xhalf * x.f * x.f );
  return x.f;
}

inline void AliHLTTPCCASliceData::CreateGrid( AliHLTTPCCARow *row, const AliHLTTPCCAClusterData &data, int ClusterDataHitNumberOffset )
{
  // grid creation

  if ( row->NHits() <= 0 ) { // no hits or invalid data
    // grid coordinates don't matter, since there are no hits
    row->fGrid.CreateEmpty();
    return;
  }

  float yMin =  1.e3f;
  float yMax = -1.e3f;
  float zMin =  1.e3f;
  float zMax = -1.e3f;
  for ( int i = ClusterDataHitNumberOffset; i < ClusterDataHitNumberOffset + row->fNHits; ++i ) {
    const float y = data.Y( i );
    const float z = data.Z( i );
    if ( yMax < y ) yMax = y;
    if ( yMin > y ) yMin = y;
    if ( zMax < z ) zMax = z;
    if ( zMin > z ) zMin = z;
  }

  const float norm = fastInvSqrt( row->fNHits );
  row->fGrid.Create( yMin, yMax, zMin, zMax,
                     CAMath::Max( ( yMax - yMin ) * norm, 2.f ),
                     CAMath::Max( ( zMax - zMin ) * norm, 2.f ) );
}

inline void AliHLTTPCCASliceData::PackHitData( AliHLTTPCCARow* const row, const AliHLTArray<AliHLTTPCCAHit> &binSortedHits )
{
  // hit data packing

  static const float shortPackingConstant = 1.f / 65535.f;
  const float y0 = row->fGrid.YMin();
  const float z0 = row->fGrid.ZMin();
  const float stepY = ( row->fGrid.YMax() - y0 ) * shortPackingConstant;
  const float stepZ = ( row->fGrid.ZMax() - z0 ) * shortPackingConstant;
  const float stepYi = 1.f / stepY;
  const float stepZi = 1.f / stepZ;

  row->fHy0 = y0;
  row->fHz0 = z0;
  row->fHstepY = stepY;
  row->fHstepZ = stepZ;
  row->fHstepYi = stepYi;
  row->fHstepZi = stepZi;

  for ( int hitIndex = 0; hitIndex < row->fNHits; ++hitIndex ) {
    // bin sorted index!
    const int globalHitIndex = row->fHitNumberOffset + hitIndex;
    const AliHLTTPCCAHit &hh = binSortedHits[globalHitIndex];
    const float xx = ( ( hh.Y() - y0 ) * stepYi ) + .5 ;
    const float yy = ( ( hh.Z() - z0 ) * stepZi ) + .5 ;
    if ( xx < 0 || yy < 0 || xx >= 65536  || yy >= 65536 ) {
      std::cout << "!!!! hit packing error!!! " << xx << " " << yy << " " << std::endl;
    }
    // HitData is bin sorted
    fHitData[row->fHitNumberOffset + hitIndex].x = xx;
    fHitData[row->fHitNumberOffset + hitIndex].y = yy;
  }
}

void AliHLTTPCCASliceData::Clear()
{
  fNumberOfHits = 0;
}

void AliHLTTPCCASliceData::InitializeRows( const AliHLTTPCCAParam &p )
{
  // initialisation of rows
	if (!fRows) fRows = new AliHLTTPCCARow[HLTCA_ROW_COUNT + 1];
  for ( int i = 0; i < p.NRows(); ++i ) {
    fRows[i].fX = p.RowX( i );
    fRows[i].fMaxY = CAMath::Tan( p.DAlpha() / 2. ) * fRows[i].fX;
  }
}

#ifndef HLTCA_GPUCODE
	AliHLTTPCCASliceData::~AliHLTTPCCASliceData()
	{
		//Standard Destrcutor
		if (fRows)
		{
			if (!fIsGpuSliceData) delete[] fRows;
			fRows = NULL;
		}
		if (fMemory)
		{
			if (!fIsGpuSliceData) delete[] fMemory;
			fMemory = NULL;
		}

	}
#endif

GPUh() void AliHLTTPCCASliceData::SetGPUSliceDataMemory(void* const pSliceMemory, void* const pRowMemory)
{
	//Set Pointer to slice data memory to external memory
	fMemory = (char*) pSliceMemory;
	fRows = (AliHLTTPCCARow*) pRowMemory;
}

size_t AliHLTTPCCASliceData::SetPointers(const AliHLTTPCCAClusterData *data, bool allocate)
{
	//Set slice data internal pointers
  int hitMemCount = 0;
  for ( int rowIndex = data->FirstRow(); rowIndex <= data->LastRow(); ++rowIndex )
  {
	hitMemCount += NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(data->NumberOfClusters( rowIndex ));
  }
	//Calculate Memory needed to store hits in rows

  const int numberOfRows = data->LastRow() - data->FirstRow() + 1;
  const unsigned int kVectorAlignment = 256 /*sizeof( uint4 )*/ ;
  fNumberOfHitsPlusAlign = NextMultipleOf < ( kVectorAlignment > sizeof(HLTCA_GPU_ROWALIGNMENT) ? kVectorAlignment : sizeof(HLTCA_GPU_ROWALIGNMENT)) / sizeof( int ) > ( hitMemCount );
  fNumberOfHits = data->NumberOfClusters();
  const int firstHitInBinSize = (23 + sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(int)) * numberOfRows + 4 * fNumberOfHits + 3;
  //FIXME: sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(int) * numberOfRows is way to big and only to ensure to reserve enough memory for GPU Alignment.
  //Might be replaced by correct value

  const int memorySize =
    // LinkData, HitData
    fNumberOfHitsPlusAlign * 4 * sizeof( short ) +
    // FirstHitInBin
    NextMultipleOf<kVectorAlignment>( ( firstHitInBinSize ) * sizeof( int ) ) +
    // HitWeights, ClusterDataIndex
    fNumberOfHitsPlusAlign * 2 * sizeof( int );

  if ( 1||fMemorySize < memorySize ) { // release the memory on CPU
	fMemorySize = memorySize;
	if (allocate && !fIsGpuSliceData)
	{
		if (fMemory)
		{
			delete[] fMemory;
		}
	  fMemory = new char[fMemorySize + 4];// kVectorAlignment];
	}
  }

  char *mem = fMemory;
  AssignMemory( fLinkUpData,   mem, fNumberOfHitsPlusAlign );
  AssignMemory( fLinkDownData, mem, fNumberOfHitsPlusAlign );
  AssignMemory( fHitData,     mem, fNumberOfHitsPlusAlign );
  AssignMemory( fFirstHitInBin,  mem, firstHitInBinSize );
  fGpuMemorySize = mem - fMemory;

  //Memory Allocated below will not be copied to GPU but instead be initialized on the gpu itself. Therefore it must not be copied to GPU!
  AssignMemory( fHitWeights,   mem, fNumberOfHitsPlusAlign );
  AssignMemory( fClusterDataIndex, mem, fNumberOfHitsPlusAlign );
  return(mem - fMemory);
}

void AliHLTTPCCASliceData::InitFromClusterData( const AliHLTTPCCAClusterData &data )
{
  // initialisation from cluster data

  ////////////////////////////////////
  // 1. prepare arrays
  ////////////////////////////////////

  const int numberOfRows = data.LastRow() - data.FirstRow() + 1;
  fNumberOfHits = data.NumberOfClusters();

  /* TODO Vectorization
  for ( int rowIndex = data.FirstRow(); rowIndex <= data.LastRow(); ++rowIndex ) {
    int NumberOfClusters( int rowIndex ) const;
  }
  const int memorySize = fNumberOfHits * sizeof( short_v::Type )
  */
  SetPointers(&data, true);

  ////////////////////////////////////
  // 2. fill HitData and FirstHitInBin
  ////////////////////////////////////

  for ( int rowIndex = 0; rowIndex < data.FirstRow(); ++rowIndex ) {
    AliHLTTPCCARow &row = fRows[rowIndex];
    row.fGrid.CreateEmpty();
    row.fNHits = 0;
    row.fFullSize = 0;
    row.fHitNumberOffset = 0;
    row.fFirstHitInBinOffset = 0;

    row.fHy0 = 0.f;
    row.fHz0 = 0.f;
    row.fHstepY = 1.f;
    row.fHstepZ = 1.f;
    row.fHstepYi = 1.f;
    row.fHstepZi = 1.f;
  }
  for ( int rowIndex = data.LastRow() + 1; rowIndex < HLTCA_ROW_COUNT + 1; ++rowIndex ) {
    AliHLTTPCCARow &row = fRows[rowIndex];
    row.fGrid.CreateEmpty();
    row.fNHits = 0;
    row.fFullSize = 0;
    row.fHitNumberOffset = 0;
    row.fFirstHitInBinOffset = 0;

    row.fHy0 = 0.f;
    row.fHz0 = 0.f;
    row.fHstepY = 1.f;
    row.fHstepZ = 1.f;
    row.fHstepYi = 1.f;
    row.fHstepZi = 1.f;
  }


  AliHLTResizableArray<AliHLTTPCCAHit> binSortedHits( fNumberOfHits + sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v) * numberOfRows + 1 );

  int gridContentOffset = 0;
  int hitOffset = 0;

  int binCreationMemorySize = 103 * 2 + fNumberOfHits;
  AliHLTResizableArray<unsigned short> binCreationMemory( binCreationMemorySize );

  fGPUSharedDataReq = 0;

  for ( int rowIndex = data.FirstRow(); rowIndex <= data.LastRow(); ++rowIndex ) {
    AliHLTTPCCARow &row = fRows[rowIndex];
    row.fNHits = data.NumberOfClusters( rowIndex );
    assert( row.fNHits < ( 1 << sizeof( unsigned short ) * 8 ) );
	row.fHitNumberOffset = hitOffset;
	hitOffset += NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(data.NumberOfClusters( rowIndex ));

    row.fFirstHitInBinOffset = gridContentOffset;

    CreateGrid( &row, data, data.RowOffset( rowIndex ) );
    const AliHLTTPCCAGrid &grid = row.fGrid;
    const int numberOfBins = grid.N();

    int binCreationMemorySizeNew;
    if ( ( binCreationMemorySizeNew = numberOfBins * 2 + 6 + row.fNHits + sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(unsigned short) * numberOfRows + 1) > binCreationMemorySize ) {
      binCreationMemorySize = binCreationMemorySizeNew;
      binCreationMemory.Resize( binCreationMemorySize );
    }

    AliHLTArray<unsigned short> c = binCreationMemory;           // number of hits in all previous bins
    AliHLTArray<unsigned short> bins = c + ( numberOfBins + 3 ); // cache for the bin index for every hit in this row
    AliHLTArray<unsigned short> filled = bins + row.fNHits;      // counts how many hits there are per bin

    for ( unsigned int bin = 0; bin < row.fGrid.N() + 3; ++bin ) {
      filled[bin] = 0; // initialize filled[] to 0
    }

    for ( int hitIndex = 0; hitIndex < row.fNHits; ++hitIndex ) {
      const int globalHitIndex = data.RowOffset( rowIndex ) + hitIndex;
      const unsigned short bin = row.fGrid.GetBin( data.Y( globalHitIndex ), data.Z( globalHitIndex ) );

      bins[hitIndex] = bin;
      ++filled[bin];
    }

    unsigned short n = 0;
    for ( int bin = 0; bin < numberOfBins + 3; ++bin ) {
      c[bin] = n;
      n += filled[bin];
    }

    for ( int hitIndex = 0; hitIndex < row.fNHits; ++hitIndex ) {
      const unsigned short bin = bins[hitIndex];
      --filled[bin];
      const unsigned short ind = c[bin] + filled[bin]; // generate an index for this hit that is >= c[bin] and < c[bin + 1]
      const int globalBinsortedIndex = row.fHitNumberOffset + ind;
      const int globalHitIndex = data.RowOffset( rowIndex ) + hitIndex;

      // allows to find the global hit index / coordinates from a global bin sorted hit index
      fClusterDataIndex[globalBinsortedIndex] = globalHitIndex;
      binSortedHits[globalBinsortedIndex].SetY( data.Y( globalHitIndex ) );
      binSortedHits[globalBinsortedIndex].SetZ( data.Z( globalHitIndex ) );
    }

    PackHitData( &row, binSortedHits );

    for ( int i = 0; i < numberOfBins; ++i ) {
      fFirstHitInBin[row.fFirstHitInBinOffset + i] = c[i]; // global bin-sorted hit index
    }
    const unsigned short a = c[numberOfBins];
    // grid.N is <= row.fNHits
    const int nn = numberOfBins + grid.Ny() + 3;
    for ( int i = numberOfBins; i < nn; ++i ) {
      assert( (signed) row.fFirstHitInBinOffset + i < 23 * numberOfRows + 4 * fNumberOfHits + 3 );
      fFirstHitInBin[row.fFirstHitInBinOffset + i] = a;
    }

    row.fFullSize = nn;
    gridContentOffset += nn;

	if (NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(row.fNHits) + nn > (unsigned) fGPUSharedDataReq)
		fGPUSharedDataReq = NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(row.fNHits) + nn;

	//Make pointer aligned
	gridContentOffset = NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(gridContentOffset);
  }

#if 0
  //SG cell finder - test code

  if ( fTmpHitInputIDs ) delete[] fTmpHitInputIDs;
  fTmpHitInputIDs = new int [NHits];
  const float areaY = .5;
  const float areaZ = .5;
  int newRowNHitsTotal = 0;
  bool *usedHits = new bool [NHits];
  for ( int iHit = 0; iHit < NHits; iHit++ ) usedHits[iHit] = 0;
  for ( int iRow = 0; iRow < fParam.NRows(); iRow++ ) {
    rowHeaders[iRow*2  ] = newRowNHitsTotal; // new first hit
    rowHeaders[iRow*2+1] = 0; // new N hits
    int newRowNHits = 0;
    int oldRowFirstHit = RowFirstHit[iRow];
    int oldRowLastHit = oldRowFirstHit + RowNHits[iRow];
    for ( int iHit = oldRowFirstHit; iHit < oldRowLastHit; iHit++ ) {
      if ( usedHits[iHit] ) continue;
      float x0 = X[iHit];
      float y0 = Y[iHit];
      float z0 = Z[iHit];
      float cx = x0;
      float cy = y0;
      float cz = z0;
      int nclu = 1;
      usedHits[iHit] = 1;
      if ( 0 ) for ( int jHit = iHit + 1; jHit < oldRowLastHit; jHit++ ) {//SG!!!
          //if( usedHits[jHit] ) continue;
          float dy = Y[jHit] - y0;
          float dz = Z[jHit] - z0;
          if ( CAMath::Abs( dy ) < areaY && CAMath::Abs( dz ) < areaZ ) {
            cx += X[jHit];
            cy += Y[jHit];
            cz += Z[jHit];
            nclu++;
            usedHits[jHit] = 1;
          }
        }
      int id = newRowNHitsTotal + newRowNHits;
      hitsXYZ[id*3+0 ] = cx / nclu;
      hitsXYZ[id*3+1 ] = cy / nclu;
      hitsXYZ[id*3+2 ] = cz / nclu;
      fTmpHitInputIDs[id] = iHit;
      newRowNHits++;
    }
    rowHeaders[iRow*2+1] = newRowNHits;
    newRowNHitsTotal += newRowNHits;
  }
  NHitsTotal() = newRowNHitsTotal;
  reinterpret_cast<int*>( fInputEvent )[1+fParam.NRows()*2] = newRowNHitsTotal;

  delete[] usedHits;
#endif
}

void AliHLTTPCCASliceData::ClearHitWeights()
{
  // clear hit weights

#ifdef ENABLE_VECTORIZATION
  const int_v v0( Zero );
  const int *const end = fHitWeights + fNumberOfHits;
  for ( int *mem = fHitWeights; mem < end; mem += v0.Size ) {
    v0.store( mem );
  }
#else
  for ( int i = 0; i < fNumberOfHitsPlusAlign; ++i ) {
    fHitWeights[i] = 0;
  }
#endif
}

void AliHLTTPCCASliceData::ClearLinks()
{
  // link cleaning

#ifdef ENABLE_VECTORIZATION
  const short_v v0( -1 );
  const short *const end1 = fLinkUpData + fNumberOfHits;
  for ( short *mem = fLinkUpData; mem < end; mem += v0.Size ) {
    v0.store( mem );
  }
  const short *const end2 = fLinkDownData + fNumberOfHits;
  for ( short *mem = fLinkDownData; mem < end; mem += v0.Size ) {
    v0.store( mem );
  }
#else
  for ( int i = 0; i < fNumberOfHits; ++i ) {
    fLinkUpData[i] = -1;
  }
  for ( int i = 0; i < fNumberOfHits; ++i ) {
    fLinkDownData[i] = -1;
  }
#endif
}
Commit	Line	Data
6de2bc40	1	// **************************************************************************
	2	// * This file is property of and copyright by the ALICE HLT Project *
	3	// * All rights reserved. *
	4	// * *
	5	// * Primary Authors: *
	6	// * Copyright 2009 Matthias Kretz <kretz@kde.org> *
	7	// * *
	8	// * Permission to use, copy, modify and distribute this software and its *
	9	// * documentation strictly for non-commercial purposes is hereby granted *
	10	// * without fee, provided that the above copyright notice appears in all *
	11	// * copies and that both the copyright notice and this permission notice *
	12	// * appear in the supporting documentation. The authors make no claims *
	13	// * about the suitability of this software for any purpose. It is *
	14	// * provided "as is" without express or implied warranty. *
	15	// **************************************************************************
4acc2401	16
	17	#include "AliHLTTPCCASliceData.h"
	18	#include "AliHLTTPCCAClusterData.h"
	19	#include "AliHLTTPCCAMath.h"
	20	#include "AliHLTArray.h"
	21	#include "AliHLTTPCCAHit.h"
	22	#include "AliHLTTPCCAParam.h"
	23	#include "MemoryAssignmentHelpers.h"
b22af1bf	24	#include "AliHLTTPCCAGPUConfig.h"
b22af1bf	25	#include "AliHLTTPCCAGPUTracker.h"
4acc2401	26	#include <iostream>
	27
	28	// calculates an approximation for 1/sqrt(x)
	29	// Google for 0x5f3759df :)
	30	static inline float fastInvSqrt( float _x )
	31	{
6de2bc40	32	// the function calculates fast inverse sqrt
6de2bc40	33
4acc2401	34	union { float f; int i; } x = { _x };
	35	const float xhalf = 0.5f * x.f;
	36	x.i = 0x5f3759df - ( x.i >> 1 );
	37	x.f = x.f * ( 1.5f - xhalf * x.f * x.f );
	38	return x.f;
	39	}
	40
b22af1bf	41	inline void AliHLTTPCCASliceData::CreateGrid( AliHLTTPCCARow *row, const AliHLTTPCCAClusterData &data, int ClusterDataHitNumberOffset )
4acc2401	42	{
6de2bc40	43	// grid creation
6de2bc40	44
4acc2401	45	if ( row->NHits() <= 0 ) { // no hits or invalid data
	46	// grid coordinates don't matter, since there are no hits
	47	row->fGrid.CreateEmpty();
	48	return;
	49	}
	50
	51	float yMin = 1.e3f;
	52	float yMax = -1.e3f;
	53	float zMin = 1.e3f;
	54	float zMax = -1.e3f;
b22af1bf	55	for ( int i = ClusterDataHitNumberOffset; i < ClusterDataHitNumberOffset + row->fNHits; ++i ) {
4acc2401	56	const float y = data.Y( i );
	57	const float z = data.Z( i );
	58	if ( yMax < y ) yMax = y;
	59	if ( yMin > y ) yMin = y;
	60	if ( zMax < z ) zMax = z;
	61	if ( zMin > z ) zMin = z;
	62	}
	63
	64	const float norm = fastInvSqrt( row->fNHits );
	65	row->fGrid.Create( yMin, yMax, zMin, zMax,
	66	CAMath::Max( ( yMax - yMin ) * norm, 2.f ),
	67	CAMath::Max( ( zMax - zMin ) * norm, 2.f ) );
	68	}
	69
444e5682	70	inline void AliHLTTPCCASliceData::PackHitData( AliHLTTPCCARow* const row, const AliHLTArray<AliHLTTPCCAHit> &binSortedHits )
4acc2401	71	{
6de2bc40	72	// hit data packing
6de2bc40	73
4acc2401	74	static const float shortPackingConstant = 1.f / 65535.f;
	75	const float y0 = row->fGrid.YMin();
	76	const float z0 = row->fGrid.ZMin();
	77	const float stepY = ( row->fGrid.YMax() - y0 ) * shortPackingConstant;
	78	const float stepZ = ( row->fGrid.ZMax() - z0 ) * shortPackingConstant;
	79	const float stepYi = 1.f / stepY;
	80	const float stepZi = 1.f / stepZ;
	81
	82	row->fHy0 = y0;
	83	row->fHz0 = z0;
	84	row->fHstepY = stepY;
	85	row->fHstepZ = stepZ;
	86	row->fHstepYi = stepYi;
	87	row->fHstepZi = stepZi;
	88
	89	for ( int hitIndex = 0; hitIndex < row->fNHits; ++hitIndex ) {
	90	// bin sorted index!
	91	const int globalHitIndex = row->fHitNumberOffset + hitIndex;
	92	const AliHLTTPCCAHit &hh = binSortedHits[globalHitIndex];
	93	const float xx = ( ( hh.Y() - y0 ) * stepYi ) + .5 ;
	94	const float yy = ( ( hh.Z() - z0 ) * stepZi ) + .5 ;
	95	if ( xx < 0 \|\| yy < 0 \|\| xx >= 65536 \|\| yy >= 65536 ) {
	96	std::cout << "!!!! hit packing error!!! " << xx << " " << yy << " " << std::endl;
	97	}
	98	// HitData is bin sorted
b22af1bf	99	fHitData[row->fHitNumberOffset + hitIndex].x = xx;
b22af1bf	100	fHitData[row->fHitNumberOffset + hitIndex].y = yy;
4acc2401	101	}
	102	}
	103
	104	void AliHLTTPCCASliceData::Clear()
	105	{
	106	fNumberOfHits = 0;
	107	}
	108
	109	void AliHLTTPCCASliceData::InitializeRows( const AliHLTTPCCAParam &p )
	110	{
6de2bc40	111	// initialisation of rows
b22af1bf	112	if (!fRows) fRows = new AliHLTTPCCARow[HLTCA_ROW_COUNT + 1];
4acc2401	113	for ( int i = 0; i < p.NRows(); ++i ) {
	114	fRows[i].fX = p.RowX( i );
	115	fRows[i].fMaxY = CAMath::Tan( p.DAlpha() / 2. ) * fRows[i].fX;
	116	}
	117	}
	118
b22af1bf	119	#ifndef HLTCA_GPUCODE
	120	AliHLTTPCCASliceData::~AliHLTTPCCASliceData()
	121	{
	122	//Standard Destrcutor
	123	if (fRows)
	124	{
	125	if (!fIsGpuSliceData) delete[] fRows;
	126	fRows = NULL;
	127	}
	128	if (fMemory)
	129	{
	130	if (!fIsGpuSliceData) delete[] fMemory;
	131	fMemory = NULL;
	132	}
	133
	134	}
	135	#endif
	136
	137	GPUh() void AliHLTTPCCASliceData::SetGPUSliceDataMemory(void* const pSliceMemory, void* const pRowMemory)
4acc2401	138	{
b22af1bf	139	//Set Pointer to slice data memory to external memory
	140	fMemory = (char*) pSliceMemory;
	141	fRows = (AliHLTTPCCARow*) pRowMemory;
7be9b0d7	142	}
4acc2401	143
7be9b0d7	144	size_t AliHLTTPCCASliceData::SetPointers(const AliHLTTPCCAClusterData *data, bool allocate)
7be9b0d7	145	{
b22af1bf	146	//Set slice data internal pointers
	147	int hitMemCount = 0;
	148	for ( int rowIndex = data->FirstRow(); rowIndex <= data->LastRow(); ++rowIndex )
	149	{
	150	hitMemCount += NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(data->NumberOfClusters( rowIndex ));
	151	}
	152	//Calculate Memory needed to store hits in rows
	153
	154	const int numberOfRows = data->LastRow() - data->FirstRow() + 1;
a59a784e	155	const unsigned int kVectorAlignment = 256 /sizeof( uint4 )/ ;
a59a784e	156	fNumberOfHitsPlusAlign = NextMultipleOf < ( kVectorAlignment > sizeof(HLTCA_GPU_ROWALIGNMENT) ? kVectorAlignment : sizeof(HLTCA_GPU_ROWALIGNMENT)) / sizeof( int ) > ( hitMemCount );
b22af1bf	157	fNumberOfHits = data->NumberOfClusters();
	158	const int firstHitInBinSize = (23 + sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(int)) * numberOfRows + 4 * fNumberOfHits + 3;
	159	//FIXME: sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(int) * numberOfRows is way to big and only to ensure to reserve enough memory for GPU Alignment.
	160	//Might be replaced by correct value
	161
4acc2401	162	const int memorySize =
4acc2401	163	// LinkData, HitData
b22af1bf	164	fNumberOfHitsPlusAlign * 4 * sizeof( short ) +
4acc2401	165	// FirstHitInBin
b22af1bf	166	NextMultipleOf<kVectorAlignment>( ( firstHitInBinSize ) * sizeof( int ) ) +
4acc2401	167	// HitWeights, ClusterDataIndex
b22af1bf	168	fNumberOfHitsPlusAlign * 2 * sizeof( int );
4acc2401	169
103bbd80	170	if ( 1\|\|fMemorySize < memorySize ) { // release the memory on CPU
7be9b0d7	171	fMemorySize = memorySize;
b22af1bf	172	if (allocate && !fIsGpuSliceData)
7be9b0d7	173	{
b22af1bf	174	if (fMemory)
	175	{
	176	delete[] fMemory;
	177	}
7be9b0d7	178	fMemory = new char[fMemorySize + 4];// kVectorAlignment];
7be9b0d7	179	}
4acc2401	180	}
	181
	182	char *mem = fMemory;
b22af1bf	183	AssignMemory( fLinkUpData, mem, fNumberOfHitsPlusAlign );
	184	AssignMemory( fLinkDownData, mem, fNumberOfHitsPlusAlign );
	185	AssignMemory( fHitData, mem, fNumberOfHitsPlusAlign );
	186	AssignMemory( fFirstHitInBin, mem, firstHitInBinSize );
	187	fGpuMemorySize = mem - fMemory;
	188
	189	//Memory Allocated below will not be copied to GPU but instead be initialized on the gpu itself. Therefore it must not be copied to GPU!
	190	AssignMemory( fHitWeights, mem, fNumberOfHitsPlusAlign );
	191	AssignMemory( fClusterDataIndex, mem, fNumberOfHitsPlusAlign );
7be9b0d7	192	return(mem - fMemory);
	193	}
	194
	195	void AliHLTTPCCASliceData::InitFromClusterData( const AliHLTTPCCAClusterData &data )
	196	{
	197	// initialisation from cluster data
	198
	199	////////////////////////////////////
	200	// 1. prepare arrays
	201	////////////////////////////////////
	202
b22af1bf	203	const int numberOfRows = data.LastRow() - data.FirstRow() + 1;
7be9b0d7	204	fNumberOfHits = data.NumberOfClusters();
	205
	206	/* TODO Vectorization
	207	for ( int rowIndex = data.FirstRow(); rowIndex <= data.LastRow(); ++rowIndex ) {
	208	int NumberOfClusters( int rowIndex ) const;
	209	}
	210	const int memorySize = fNumberOfHits * sizeof( short_v::Type )
	211	*/
	212	SetPointers(&data, true);
4acc2401	213
	214	////////////////////////////////////
	215	// 2. fill HitData and FirstHitInBin
	216	////////////////////////////////////
b8139972	217
4acc2401	218	for ( int rowIndex = 0; rowIndex < data.FirstRow(); ++rowIndex ) {
	219	AliHLTTPCCARow &row = fRows[rowIndex];
	220	row.fGrid.CreateEmpty();
	221	row.fNHits = 0;
	222	row.fFullSize = 0;
	223	row.fHitNumberOffset = 0;
	224	row.fFirstHitInBinOffset = 0;
	225
	226	row.fHy0 = 0.f;
	227	row.fHz0 = 0.f;
	228	row.fHstepY = 1.f;
	229	row.fHstepZ = 1.f;
	230	row.fHstepYi = 1.f;
	231	row.fHstepZi = 1.f;
	232	}
b22af1bf	233	for ( int rowIndex = data.LastRow() + 1; rowIndex < HLTCA_ROW_COUNT + 1; ++rowIndex ) {
4acc2401	234	AliHLTTPCCARow &row = fRows[rowIndex];
	235	row.fGrid.CreateEmpty();
	236	row.fNHits = 0;
	237	row.fFullSize = 0;
	238	row.fHitNumberOffset = 0;
	239	row.fFirstHitInBinOffset = 0;
	240
	241	row.fHy0 = 0.f;
	242	row.fHz0 = 0.f;
	243	row.fHstepY = 1.f;
	244	row.fHstepZ = 1.f;
	245	row.fHstepYi = 1.f;
	246	row.fHstepZi = 1.f;
	247	}
	248
6de2bc40	249
b22af1bf	250	AliHLTResizableArray<AliHLTTPCCAHit> binSortedHits( fNumberOfHits + sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v) * numberOfRows + 1 );
4acc2401	251
4acc2401	252	int gridContentOffset = 0;
b22af1bf	253	int hitOffset = 0;
4acc2401	254
	255	int binCreationMemorySize = 103 * 2 + fNumberOfHits;
	256	AliHLTResizableArray<unsigned short> binCreationMemory( binCreationMemorySize );
	257
b22af1bf	258	fGPUSharedDataReq = 0;
b22af1bf	259
4acc2401	260	for ( int rowIndex = data.FirstRow(); rowIndex <= data.LastRow(); ++rowIndex ) {
4acc2401	261	AliHLTTPCCARow &row = fRows[rowIndex];
4acc2401	262	row.fNHits = data.NumberOfClusters( rowIndex );
4acc2401	263	assert( row.fNHits < ( 1 << sizeof( unsigned short ) * 8 ) );
b22af1bf	264	row.fHitNumberOffset = hitOffset;
	265	hitOffset += NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(data.NumberOfClusters( rowIndex ));
	266
4acc2401	267	row.fFirstHitInBinOffset = gridContentOffset;
4acc2401	268
b22af1bf	269	CreateGrid( &row, data, data.RowOffset( rowIndex ) );
4acc2401	270	const AliHLTTPCCAGrid &grid = row.fGrid;
	271	const int numberOfBins = grid.N();
	272
	273	int binCreationMemorySizeNew;
b22af1bf	274	if ( ( binCreationMemorySizeNew = numberOfBins * 2 + 6 + row.fNHits + sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(unsigned short) * numberOfRows + 1) > binCreationMemorySize ) {
4acc2401	275	binCreationMemorySize = binCreationMemorySizeNew;
	276	binCreationMemory.Resize( binCreationMemorySize );
	277	}
b8139972	278
4acc2401	279	AliHLTArray<unsigned short> c = binCreationMemory; // number of hits in all previous bins
	280	AliHLTArray<unsigned short> bins = c + ( numberOfBins + 3 ); // cache for the bin index for every hit in this row
	281	AliHLTArray<unsigned short> filled = bins + row.fNHits; // counts how many hits there are per bin
	282
	283	for ( unsigned int bin = 0; bin < row.fGrid.N() + 3; ++bin ) {
	284	filled[bin] = 0; // initialize filled[] to 0
	285	}
b8139972	286
4acc2401	287	for ( int hitIndex = 0; hitIndex < row.fNHits; ++hitIndex ) {
b22af1bf	288	const int globalHitIndex = data.RowOffset( rowIndex ) + hitIndex;
4acc2401	289	const unsigned short bin = row.fGrid.GetBin( data.Y( globalHitIndex ), data.Z( globalHitIndex ) );
b22af1bf	290
4acc2401	291	bins[hitIndex] = bin;
	292	++filled[bin];
	293	}
	294
	295	unsigned short n = 0;
	296	for ( int bin = 0; bin < numberOfBins + 3; ++bin ) {
	297	c[bin] = n;
	298	n += filled[bin];
	299	}
	300
	301	for ( int hitIndex = 0; hitIndex < row.fNHits; ++hitIndex ) {
	302	const unsigned short bin = bins[hitIndex];
	303	--filled[bin];
	304	const unsigned short ind = c[bin] + filled[bin]; // generate an index for this hit that is >= c[bin] and < c[bin + 1]
	305	const int globalBinsortedIndex = row.fHitNumberOffset + ind;
b22af1bf	306	const int globalHitIndex = data.RowOffset( rowIndex ) + hitIndex;
4acc2401	307
	308	// allows to find the global hit index / coordinates from a global bin sorted hit index
	309	fClusterDataIndex[globalBinsortedIndex] = globalHitIndex;
	310	binSortedHits[globalBinsortedIndex].SetY( data.Y( globalHitIndex ) );
	311	binSortedHits[globalBinsortedIndex].SetZ( data.Z( globalHitIndex ) );
	312	}
	313
	314	PackHitData( &row, binSortedHits );
	315
	316	for ( int i = 0; i < numberOfBins; ++i ) {
	317	fFirstHitInBin[row.fFirstHitInBinOffset + i] = c[i]; // global bin-sorted hit index
	318	}
	319	const unsigned short a = c[numberOfBins];
	320	// grid.N is <= row.fNHits
b8139972	321	const int nn = numberOfBins + grid.Ny() + 3;
4acc2401	322	for ( int i = numberOfBins; i < nn; ++i ) {
b22af1bf	323	assert( (signed) row.fFirstHitInBinOffset + i < 23 * numberOfRows + 4 * fNumberOfHits + 3 );
4acc2401	324	fFirstHitInBin[row.fFirstHitInBinOffset + i] = a;
	325	}
	326
	327	row.fFullSize = nn;
	328	gridContentOffset += nn;
b22af1bf	329
	330	if (NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(row.fNHits) + nn > (unsigned) fGPUSharedDataReq)
	331	fGPUSharedDataReq = NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(row.fNHits) + nn;
	332
	333	//Make pointer aligned
	334	gridContentOffset = NextMultipleOf<sizeof(HLTCA_GPU_ROWALIGNMENT) / sizeof(ushort_v)>(gridContentOffset);
4acc2401	335	}
	336
	337	#if 0
	338	//SG cell finder - test code
	339
	340	if ( fTmpHitInputIDs ) delete[] fTmpHitInputIDs;
	341	fTmpHitInputIDs = new int [NHits];
	342	const float areaY = .5;
	343	const float areaZ = .5;
	344	int newRowNHitsTotal = 0;
	345	bool *usedHits = new bool [NHits];
	346	for ( int iHit = 0; iHit < NHits; iHit++ ) usedHits[iHit] = 0;
	347	for ( int iRow = 0; iRow < fParam.NRows(); iRow++ ) {
	348	rowHeaders[iRow*2 ] = newRowNHitsTotal; // new first hit
	349	rowHeaders[iRow*2+1] = 0; // new N hits
	350	int newRowNHits = 0;
	351	int oldRowFirstHit = RowFirstHit[iRow];
	352	int oldRowLastHit = oldRowFirstHit + RowNHits[iRow];
	353	for ( int iHit = oldRowFirstHit; iHit < oldRowLastHit; iHit++ ) {
	354	if ( usedHits[iHit] ) continue;
	355	float x0 = X[iHit];
	356	float y0 = Y[iHit];
	357	float z0 = Z[iHit];
	358	float cx = x0;
	359	float cy = y0;
	360	float cz = z0;
	361	int nclu = 1;
	362	usedHits[iHit] = 1;
	363	if ( 0 ) for ( int jHit = iHit + 1; jHit < oldRowLastHit; jHit++ ) {//SG!!!
	364	//if( usedHits[jHit] ) continue;
	365	float dy = Y[jHit] - y0;
	366	float dz = Z[jHit] - z0;
	367	if ( CAMath::Abs( dy ) < areaY && CAMath::Abs( dz ) < areaZ ) {
	368	cx += X[jHit];
	369	cy += Y[jHit];
	370	cz += Z[jHit];
	371	nclu++;
	372	usedHits[jHit] = 1;
	373	}
	374	}
	375	int id = newRowNHitsTotal + newRowNHits;
	376	hitsXYZ[id*3+0 ] = cx / nclu;
	377	hitsXYZ[id*3+1 ] = cy / nclu;
	378	hitsXYZ[id*3+2 ] = cz / nclu;
	379	fTmpHitInputIDs[id] = iHit;
	380	newRowNHits++;
	381	}
	382	rowHeaders[iRow*2+1] = newRowNHits;
	383	newRowNHitsTotal += newRowNHits;
	384	}
	385	NHitsTotal() = newRowNHitsTotal;
	386	reinterpret_cast<int>( fInputEvent )[1+fParam.NRows()2] = newRowNHitsTotal;
	387
	388	delete[] usedHits;
	389	#endif
	390	}
	391
	392	void AliHLTTPCCASliceData::ClearHitWeights()
	393	{
6de2bc40	394	// clear hit weights
6de2bc40	395
4acc2401	396	#ifdef ENABLE_VECTORIZATION
	397	const int_v v0( Zero );
	398	const int *const end = fHitWeights + fNumberOfHits;
	399	for ( int *mem = fHitWeights; mem < end; mem += v0.Size ) {
	400	v0.store( mem );
	401	}
	402	#else
b22af1bf	403	for ( int i = 0; i < fNumberOfHitsPlusAlign; ++i ) {
4acc2401	404	fHitWeights[i] = 0;
	405	}
	406	#endif
	407	}
	408
	409	void AliHLTTPCCASliceData::ClearLinks()
	410	{
6de2bc40	411	// link cleaning
6de2bc40	412
4acc2401	413	#ifdef ENABLE_VECTORIZATION
	414	const short_v v0( -1 );
	415	const short *const end1 = fLinkUpData + fNumberOfHits;
	416	for ( short *mem = fLinkUpData; mem < end; mem += v0.Size ) {
	417	v0.store( mem );
	418	}
	419	const short *const end2 = fLinkDownData + fNumberOfHits;
	420	for ( short *mem = fLinkDownData; mem < end; mem += v0.Size ) {
	421	v0.store( mem );
	422	}
	423	#else
	424	for ( int i = 0; i < fNumberOfHits; ++i ) {
	425	fLinkUpData[i] = -1;
	426	}
	427	for ( int i = 0; i < fNumberOfHits; ++i ) {
	428	fLinkDownData[i] = -1;
	429	}
	430	#endif
	431	}
	432