HLT/TPCLib/tracking/AliHLTTPCHoughTransformerRow.cxx

   1 // @(#) $Id$
   2 // origin hough/AliL3HoughTransformerRow.cxx,v 1.21 Tue Mar 28 18:05:12 2006 UTC by alibrary
   3
   4 /** @file   AliHLTTPCHoughTransformerRow.cxx
   5     @author Cvetan Cheshkov <mailto:cvetan.cheshkov@cern.ch>
   6     @date
   7     @brief  Implementation of fast HLT TPC hough transform tracking. */
   8
   9 //_____________________________________________________________
  10 // AliHLTTPCHoughTransformerRow
  11 //
  12 // Impelementation of the so called "TPC rows Hough transformation" class
  13 //
  14 // Transforms the TPC data into the hough space and counts the missed TPC
  15 // rows corresponding to each track cnadidate - hough space bin
  16
  17 #include "AliHLTStdIncludes.h"
  18
  19 #include "AliHLTTPCLogging.h"
  20 #include "AliHLTTPCMemHandler.h"
  21 #include "AliHLTTPCTransform.h"
  22 #include "AliHLTTPCDigitData.h"
  23 #include "AliHLTTPCHistogramAdaptive.h"
  24 #include "AliHLTTPCHoughTrack.h"
  25 #include "AliHLTTPCHoughTransformerRow.h"
  26 #include "AliTPCRawStream.h"
  27
  28 #if __GNUC__ >= 3
  29 using namespace std;
  30 #endif
  31
  32 /** ROOT macro for the implementation of ROOT specific class methods */
  33 ClassImp(AliHLTTPCHoughTransformerRow);
  34
  35 Float_t AliHLTTPCHoughTransformerRow::fgBeta1 = 1.0/AliHLTTPCTransform::Row2X(84);
  36 Float_t AliHLTTPCHoughTransformerRow::fgBeta2 = 1.0/(AliHLTTPCTransform::Row2X(158)*(1.0+tan(AliHLTTPCTransform::Pi()*10/180)*tan(AliHLTTPCTransform::Pi()*10/180)));
  37 Float_t AliHLTTPCHoughTransformerRow::fgDAlpha = 0.22;
  38 Float_t AliHLTTPCHoughTransformerRow::fgDEta   = 0.40;
  39 Double_t AliHLTTPCHoughTransformerRow::fgEtaCalcParam1 = 1.0289;
  40 Double_t AliHLTTPCHoughTransformerRow::fgEtaCalcParam2 = 0.15192;
  41 Double_t AliHLTTPCHoughTransformerRow::fgEtaCalcParam3 = 1./(32.*600.*600.);
  42
  43 AliHLTTPCHoughTransformerRow::AliHLTTPCHoughTransformerRow()
  44 {
  45   //Default constructor
  46   fParamSpace = 0;
  47
  48   fGapCount = 0;
  49   fCurrentRowCount = 0;
  50 #ifdef do_mc
  51   fTrackID = 0;
  52 #endif
  53   fTrackNRows = 0;
  54   fTrackFirstRow = 0;
  55   fTrackLastRow = 0;
  56   fInitialGapCount = 0;
  57
  58   fPrevBin = 0;
  59   fNextBin = 0;
  60   fNextRow = 0;
  61
  62   fStartPadParams = 0;
  63   fEndPadParams = 0;
  64   fLUTr = 0;
  65   fLUTforwardZ = 0;
  66   fLUTbackwardZ = 0;
  67
  68 }
  69
  70 AliHLTTPCHoughTransformerRow::AliHLTTPCHoughTransformerRow(Int_t slice,Int_t patch,Int_t netasegments,Bool_t /*DoMC*/,Float_t zvertex) : AliHLTTPCHoughTransformer(slice,patch,netasegments,zvertex)
  71 {
  72   //Normal constructor
  73   fParamSpace = 0;
  74
  75   fGapCount = 0;
  76   fCurrentRowCount = 0;
  77 #ifdef do_mc
  78   fTrackID = 0;
  79 #endif
  80
  81   fTrackNRows = 0;
  82   fTrackFirstRow = 0;
  83   fTrackLastRow = 0;
  84   fInitialGapCount = 0;
  85
  86   fPrevBin = 0;
  87   fNextBin = 0;
  88   fNextRow = 0;
  89
  90   fStartPadParams = 0;
  91   fEndPadParams = 0;
  92   fLUTr = 0;
  93   fLUTforwardZ = 0;
  94   fLUTbackwardZ = 0;
  95
  96 }
  97
  98 AliHLTTPCHoughTransformerRow::~AliHLTTPCHoughTransformerRow()
  99 {
 100   //Destructor
 101   if(fLastTransformer) return;
 102   DeleteHistograms();
 103 #ifdef do_mc
 104   if(fTrackID)
 105     {
 106       for(Int_t i=0; i<GetNEtaSegments(); i++)
 107         {
 108           if(!fTrackID[i]) continue;
 109           delete fTrackID[i];
 110         }
 111       delete [] fTrackID;
 112       fTrackID = 0;
 113     }
 114 #endif
 115
 116   if(fGapCount)
 117     {
 118       for(Int_t i=0; i<GetNEtaSegments(); i++)
 119         {
 120           if(!fGapCount[i]) continue;
 121           delete [] fGapCount[i];
 122         }
 123       delete [] fGapCount;
 124       fGapCount = 0;
 125     }
 126   if(fCurrentRowCount)
 127     {
 128       for(Int_t i=0; i<GetNEtaSegments(); i++)
 129         {
 130           if(fCurrentRowCount[i])
 131             delete [] fCurrentRowCount[i];
 132         }
 133       delete [] fCurrentRowCount;
 134       fCurrentRowCount = 0;
 135     }
 136   if(fTrackNRows)
 137     {
 138       delete [] fTrackNRows;
 139       fTrackNRows = 0;
 140     }
 141   if(fTrackFirstRow)
 142     {
 143       delete [] fTrackFirstRow;
 144       fTrackFirstRow = 0;
 145     }
 146   if(fTrackLastRow)
 147     {
 148       delete [] fTrackLastRow;
 149       fTrackLastRow = 0;
 150     }
 151   if(fInitialGapCount)
 152     {
 153       delete [] fInitialGapCount;
 154       fInitialGapCount = 0;
 155     }
 156   if(fPrevBin)
 157     {
 158       for(Int_t i=0; i<GetNEtaSegments(); i++)
 159         {
 160           if(!fPrevBin[i]) continue;
 161           delete [] fPrevBin[i];
 162         }
 163       delete [] fPrevBin;
 164       fPrevBin = 0;
 165     }
 166   if(fNextBin)
 167     {
 168       for(Int_t i=0; i<GetNEtaSegments(); i++)
 169         {
 170           if(!fNextBin[i]) continue;
 171           delete [] fNextBin[i];
 172         }
 173       delete [] fNextBin;
 174       fNextBin = 0;
 175     }
 176   if(fNextRow)
 177     {
 178       for(Int_t i=0; i<GetNEtaSegments(); i++)
 179         {
 180           if(!fNextRow[i]) continue;
 181           delete [] fNextRow[i];
 182         }
 183       delete [] fNextRow;
 184       fNextRow = 0;
 185     }
 186   if(fStartPadParams)
 187     {
 188       for(Int_t i = AliHLTTPCTransform::GetFirstRow(0); i<=AliHLTTPCTransform::GetLastRow(5); i++)
 189         {
 190           if(!fStartPadParams[i]) continue;
 191           delete [] fStartPadParams[i];
 192         }
 193       delete [] fStartPadParams;
 194       fStartPadParams = 0;
 195     }
 196   if(fEndPadParams)
 197     {
 198       for(Int_t i = AliHLTTPCTransform::GetFirstRow(0); i<=AliHLTTPCTransform::GetLastRow(5); i++)
 199         {
 200           if(!fEndPadParams[i]) continue;
 201           delete [] fEndPadParams[i];
 202         }
 203       delete [] fEndPadParams;
 204       fEndPadParams = 0;
 205     }
 206   if(fLUTr)
 207     {
 208       for(Int_t i = AliHLTTPCTransform::GetFirstRow(0); i<=AliHLTTPCTransform::GetLastRow(5); i++)
 209         {
 210           if(!fLUTr[i]) continue;
 211           delete [] fLUTr[i];
 212         }
 213       delete [] fLUTr;
 214       fLUTr = 0;
 215     }
 216   if(fLUTforwardZ)
 217     {
 218       delete[] fLUTforwardZ;
 219       fLUTforwardZ=0;
 220     }
 221   if(fLUTbackwardZ)
 222     {
 223       delete[] fLUTbackwardZ;
 224       fLUTbackwardZ=0;
 225     }
 226 }
 227
 228 void AliHLTTPCHoughTransformerRow::DeleteHistograms()
 229 {
 230   // Clean up
 231   if(!fParamSpace)
 232     return;
 233   for(Int_t i=0; i<GetNEtaSegments(); i++)
 234     {
 235       if(!fParamSpace[i]) continue;
 236       delete fParamSpace[i];
 237     }
 238   delete [] fParamSpace;
 239 }
 240
 241 void AliHLTTPCHoughTransformerRow::CreateHistograms(Int_t nxbin,Float_t xmin,Float_t xmax,
 242                                                 Int_t nybin,Float_t ymin,Float_t ymax)
 243 {
 244   //Create the histograms (parameter space)
 245   //nxbin = #bins in X
 246   //nybin = #bins in Y
 247   //xmin xmax ymin ymax = histogram limits in X and Y
 248   if(fLastTransformer) {
 249     SetTransformerArrays((AliHLTTPCHoughTransformerRow *)fLastTransformer);
 250     return;
 251   }
 252   fParamSpace = new AliHLTTPCHistogram*[GetNEtaSegments()];
 253
 254   Char_t histname[256];
 255   for(Int_t i=0; i<GetNEtaSegments(); i++)
 256     {
 257       sprintf(histname,"paramspace_%d",i);
 258       fParamSpace[i] = new AliHLTTPCHistogram(histname,"",nxbin,xmin,xmax,nybin,ymin,ymax);
 259     }
 260 #ifdef do_mc
 261   {
 262     AliHLTTPCHistogram *hist = fParamSpace[0];
 263     Int_t ncellsx = (hist->GetNbinsX()+3)/2;
 264     Int_t ncellsy = (hist->GetNbinsY()+3)/2;
 265     Int_t ncells = ncellsx*ncellsy;
 266     if(!fTrackID)
 267       {
 268         LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 269           <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(AliHLTTrackIndex)<<" bytes to fTrackID"<<ENDLOG;
 270         fTrackID = new AliHLTTrackIndex*[GetNEtaSegments()];
 271         for(Int_t i=0; i<GetNEtaSegments(); i++)
 272           fTrackID[i] = new AliHLTTrackIndex[ncells];
 273       }
 274   }
 275 #endif
 276   AliHLTTPCHistogram *hist = fParamSpace[0];
 277   Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
 278   if(!fGapCount)
 279     {
 280       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 281         <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fGapCount"<<ENDLOG;
 282       fGapCount = new UChar_t*[GetNEtaSegments()];
 283       for(Int_t i=0; i<GetNEtaSegments(); i++)
 284         fGapCount[i] = new UChar_t[ncells];
 285     }
 286   if(!fCurrentRowCount)
 287     {
 288       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 289         <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fCurrentRowCount"<<ENDLOG;
 290       fCurrentRowCount = new UChar_t*[GetNEtaSegments()];
 291       for(Int_t i=0; i<GetNEtaSegments(); i++)
 292         fCurrentRowCount[i] = new UChar_t[ncells];
 293     }
 294   if(!fPrevBin)
 295     {
 296       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 297         <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fPrevBin"<<ENDLOG;
 298       fPrevBin = new UChar_t*[GetNEtaSegments()];
 299       for(Int_t i=0; i<GetNEtaSegments(); i++)
 300         fPrevBin[i] = new UChar_t[ncells];
 301     }
 302   if(!fNextBin)
 303     {
 304       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 305         <<"Transformer: Allocating "<<GetNEtaSegments()*ncells*sizeof(UChar_t)<<" bytes to fNextBin"<<ENDLOG;
 306       fNextBin = new UChar_t*[GetNEtaSegments()];
 307       for(Int_t i=0; i<GetNEtaSegments(); i++)
 308         fNextBin[i] = new UChar_t[ncells];
 309     }
 310   Int_t ncellsy = hist->GetNbinsY()+2;
 311   if(!fNextRow)
 312     {
 313       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 314         <<"Transformer: Allocating "<<GetNEtaSegments()*ncellsy*sizeof(UChar_t)<<" bytes to fNextRow"<<ENDLOG;
 315       fNextRow = new UChar_t*[GetNEtaSegments()];
 316       for(Int_t i=0; i<GetNEtaSegments(); i++)
 317         fNextRow[i] = new UChar_t[ncellsy];
 318     }
 319
 320   if(!fTrackNRows)
 321     {
 322       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 323         <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fTrackNRows"<<ENDLOG;
 324       fTrackNRows = new UChar_t[ncells];
 325       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 326         <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fTrackFirstRow"<<ENDLOG;
 327       fTrackFirstRow = new UChar_t[ncells];
 328       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 329         <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fTrackLastRow"<<ENDLOG;
 330       fTrackLastRow = new UChar_t[ncells];
 331       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 332         <<"Transformer: Allocating "<<ncells*sizeof(UChar_t)<<" bytes to fInitialGapCount"<<ENDLOG;
 333       fInitialGapCount = new UChar_t[ncells];
 334
 335
 336       AliHLTTPCHoughTrack track;
 337       Int_t xmin = hist->GetFirstXbin();
 338       Int_t xmax = hist->GetLastXbin();
 339       Int_t xmiddle = (hist->GetNbinsX()+1)/2;
 340       Int_t ymin = hist->GetFirstYbin();
 341       Int_t ymax = hist->GetLastYbin();
 342       Int_t nxbins = hist->GetNbinsX()+2;
 343       Int_t nxgrid = (hist->GetNbinsX()+3)/2+1;
 344       Int_t nygrid = hist->GetNbinsY()+3;
 345
 346       AliHLTTrackLength *tracklength = new AliHLTTrackLength[nxgrid*nygrid];
 347       memset(tracklength,0,nxgrid*nygrid*sizeof(AliHLTTrackLength));
 348
 349       for(Int_t ybin=ymin-1; ybin<=(ymax+1); ybin++)
 350         {
 351           for(Int_t xbin=xmin-1; xbin<=xmiddle; xbin++)
 352             {
 353               fTrackNRows[xbin + ybin*nxbins] = 255;
 354               for(Int_t deltay = 0; deltay <= 1; deltay++) {
 355                 for(Int_t deltax = 0; deltax <= 1; deltax++) {
 356
 357                   AliHLTTrackLength *curtracklength = &tracklength[(xbin + deltax) + (ybin + deltay)*nxgrid];
 358                   UInt_t maxfirstrow = 0;
 359                   UInt_t maxlastrow = 0;
 360                   Float_t maxtrackpt = 0;
 361                   if(curtracklength->fIsFilled) {
 362                     maxfirstrow = curtracklength->fFirstRow;
 363                     maxlastrow = curtracklength->fLastRow;
 364                     maxtrackpt = curtracklength->fTrackPt;
 365                   }
 366                   else {
 367                     Float_t xtrack = hist->GetPreciseBinCenterX((Float_t)xbin+0.5*(Float_t)(2*deltax-1));
 368                     Float_t ytrack = hist->GetPreciseBinCenterY((Float_t)ybin+0.5*(Float_t)(2*deltay-1));
 369
 370                     Float_t psi = atan((xtrack-ytrack)/(fgBeta1-fgBeta2));
 371                     Float_t kappa = 2.0*(xtrack*cos(psi)-fgBeta1*sin(psi));
 372                     track.SetTrackParameters(kappa,psi,1);
 373                     maxtrackpt = track.GetPt();
 374                     if(maxtrackpt < 0.9*0.1*AliHLTTPCTransform::GetSolenoidField())
 375                       {
 376                         maxfirstrow = maxlastrow = 0;
 377                         curtracklength->fIsFilled = kTRUE;
 378                         curtracklength->fFirstRow = maxfirstrow;
 379                         curtracklength->fLastRow = maxlastrow;
 380                         curtracklength->fTrackPt = maxtrackpt;
 381                       }
 382                     else
 383                       {
 384                         Bool_t firstrow = kFALSE;
 385                         UInt_t curfirstrow = 0;
 386                         UInt_t curlastrow = 0;
 387
 388                         Double_t centerx = track.GetCenterX();
 389                         Double_t centery = track.GetCenterY();
 390                         Double_t radius = track.GetRadius();
 391
 392                         for(Int_t j=AliHLTTPCTransform::GetFirstRow(0); j<=AliHLTTPCTransform::GetLastRow(5); j++)
 393                           {
 394                             Float_t hit[3];
 395                             //                if(!track.GetCrossingPoint(j,hit)) continue;
 396                             hit[0] = AliHLTTPCTransform::Row2X(j);
 397                             Double_t aa = (hit[0] - centerx)*(hit[0] - centerx);
 398                             Double_t r2 = radius*radius;
 399                             if(aa > r2)
 400                               continue;
 401
 402                             Double_t aa2 = sqrt(r2 - aa);
 403                             Double_t y1 = centery + aa2;
 404                             Double_t y2 = centery - aa2;
 405                             hit[1] = y1;
 406                             if(fabs(y2) < fabs(y1)) hit[1] = y2;
 407
 408                             hit[2] = 0;
 409
 410                             AliHLTTPCTransform::LocHLT2Raw(hit,0,j);
 411                             hit[1] += 0.5;
 412                             if(hit[1]>=0 && hit[1]<AliHLTTPCTransform::GetNPads(j))
 413                               {
 414                                 if(!firstrow) {
 415                                   curfirstrow = j;
 416                                   firstrow = kTRUE;
 417                                 }
 418                                 curlastrow = j;
 419                               }
 420                             else {
 421                               if(firstrow) {
 422                                 firstrow = kFALSE;
 423                                 if((curlastrow-curfirstrow) >= (maxlastrow-maxfirstrow)) {
 424                                   maxfirstrow = curfirstrow;
 425                                   maxlastrow = curlastrow;
 426                                 }
 427                               }
 428                             }
 429                           }
 430                         if((curlastrow-curfirstrow) >= (maxlastrow-maxfirstrow)) {
 431                           maxfirstrow = curfirstrow;
 432                           maxlastrow = curlastrow;
 433                         }
 434
 435                         curtracklength->fIsFilled = kTRUE;
 436                         curtracklength->fFirstRow = maxfirstrow;
 437                         curtracklength->fLastRow = maxlastrow;
 438                         curtracklength->fTrackPt = maxtrackpt;
 439                       }
 440                   }
 441                   if((maxlastrow-maxfirstrow) < fTrackNRows[xbin + ybin*nxbins]) {
 442                     fTrackNRows[xbin + ybin*nxbins] = maxlastrow-maxfirstrow;
 443                     fInitialGapCount[xbin + ybin*nxbins] = 1;
 444                     if((maxlastrow-maxfirstrow+1)<=MIN_TRACK_LENGTH)
 445                       fInitialGapCount[xbin + ybin*nxbins] = MAX_N_GAPS+1;
 446                     if(maxtrackpt < 0.9*0.1*AliHLTTPCTransform::GetSolenoidField())
 447                       fInitialGapCount[xbin + ybin*nxbins] = MAX_N_GAPS;
 448                     fTrackFirstRow[xbin + ybin*nxbins] = maxfirstrow;
 449                     fTrackLastRow[xbin + ybin*nxbins] = maxlastrow;
 450
 451                     Int_t xmirror = xmax - xbin + 1;
 452                     Int_t ymirror = ymax - ybin + 1;
 453                     fTrackNRows[xmirror + ymirror*nxbins] = fTrackNRows[xbin + ybin*nxbins];
 454                     fInitialGapCount[xmirror + ymirror*nxbins] = fInitialGapCount[xbin + ybin*nxbins];
 455                     fTrackFirstRow[xmirror + ymirror*nxbins] = fTrackFirstRow[xbin + ybin*nxbins];
 456                     fTrackLastRow[xmirror + ymirror*nxbins] = fTrackLastRow[xbin + ybin*nxbins];
 457                   }
 458                 }
 459               }
 460               //              cout<<" fTrackNRows "<<(Int_t)fInitialGapCount[xbin + ybin*nxbins]<<" "<<xbin<<" "<<ybin<<" "<<(Int_t)fTrackNRows[xbin + ybin*nxbins]<<" "<<(Int_t)fTrackFirstRow[xbin + ybin*nxbins]<<" "<<(Int_t)fTrackLastRow[xbin + ybin*nxbins]<<" "<<endl;
 461             }
 462         }
 463       delete [] tracklength;
 464     }
 465
 466   if(!fStartPadParams)
 467     {
 468       Int_t nrows = AliHLTTPCTransform::GetLastRow(5) - AliHLTTPCTransform::GetFirstRow(0) + 1;
 469       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 470         <<"Transformer: Allocating about "<<nrows*100*sizeof(AliHLTPadHoughParams)<<" bytes to fStartPadParams"<<ENDLOG;
 471       fStartPadParams = new AliHLTPadHoughParams*[nrows];
 472       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 473         <<"Transformer: Allocating about "<<nrows*100*sizeof(AliHLTPadHoughParams)<<" bytes to fEndPadParams"<<ENDLOG;
 474       fEndPadParams = new AliHLTPadHoughParams*[nrows];
 475       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 476         <<"Transformer: Allocating about "<<nrows*100*sizeof(Float_t)<<" bytes to fLUTr"<<ENDLOG;
 477       fLUTr = new Float_t*[nrows];
 478
 479       Float_t beta1 = fgBeta1;
 480       Float_t beta2 = fgBeta2;
 481       Float_t beta1minusbeta2 = fgBeta1 - fgBeta2;
 482       Float_t ymin = hist->GetYmin();
 483       Float_t histbin = hist->GetBinWidthY();
 484       Float_t xmin = hist->GetXmin();
 485       Float_t xmax = hist->GetXmax();
 486       Float_t xbin = (xmax-xmin)/hist->GetNbinsX();
 487       Int_t firstbinx = hist->GetFirstXbin();
 488       Int_t lastbinx = hist->GetLastXbin();
 489       Int_t nbinx = hist->GetNbinsX()+2;
 490       Int_t firstbin = hist->GetFirstYbin();
 491       Int_t lastbin = hist->GetLastYbin();
 492       for(Int_t i=AliHLTTPCTransform::GetFirstRow(0); i<=AliHLTTPCTransform::GetLastRow(5); i++)
 493         {
 494           Int_t npads = AliHLTTPCTransform::GetNPads(i);
 495           Int_t ipatch = AliHLTTPCTransform::GetPatch(i);
 496           Double_t padpitch = AliHLTTPCTransform::GetPadPitchWidth(ipatch);
 497           Float_t x = AliHLTTPCTransform::Row2X(i);
 498           Float_t x2 = x*x;
 499
 500           fStartPadParams[i] = new AliHLTPadHoughParams[npads];
 501           fEndPadParams[i] = new AliHLTPadHoughParams[npads];
 502           fLUTr[i] = new Float_t[npads];
 503           for(Int_t pad=0; pad<npads; pad++)
 504             {
 505               Float_t y = (pad-0.5*(npads-1))*padpitch;
 506               fLUTr[i][pad] = sqrt(x2 + y*y);
 507               Float_t starty = (pad-0.5*npads)*padpitch;
 508               Float_t r1 = x2 + starty*starty;
 509               Float_t xoverr1 = x/r1;
 510               Float_t startyoverr1 = starty/r1;
 511               Float_t endy = (pad-0.5*(npads-2))*padpitch;
 512               Float_t r2 = x2 + endy*endy;
 513               Float_t xoverr2 = x/r2;
 514               Float_t endyoverr2 = endy/r2;
 515               Float_t a1 = beta1minusbeta2/(xoverr1-beta2);
 516               Float_t b1 = (xoverr1-beta1)/(xoverr1-beta2);
 517               Float_t a2 = beta1minusbeta2/(xoverr2-beta2);
 518               Float_t b2 = (xoverr2-beta1)/(xoverr2-beta2);
 519
 520               Float_t alpha1 = (a1*startyoverr1+b1*ymin-xmin)/xbin;
 521               Float_t deltaalpha1 = b1*histbin/xbin;
 522               if(b1<0)
 523                 alpha1 += deltaalpha1;
 524               Float_t alpha2 = (a2*endyoverr2+b2*ymin-xmin)/xbin;
 525               Float_t deltaalpha2 = b2*histbin/xbin;
 526               if(b2>=0)
 527                 alpha2 += deltaalpha2;
 528
 529               fStartPadParams[i][pad].fAlpha = alpha1;
 530               fStartPadParams[i][pad].fDeltaAlpha = deltaalpha1;
 531               fEndPadParams[i][pad].fAlpha = alpha2;
 532               fEndPadParams[i][pad].fDeltaAlpha = deltaalpha2;
 533
 534               //Find the first and last bin rows to be filled
 535               Bool_t binfound1 = kFALSE;
 536               Bool_t binfound2 = kFALSE;
 537               Int_t firstbin1 = lastbin;
 538               Int_t firstbin2 = lastbin;
 539               Int_t lastbin1 = firstbin;
 540               Int_t lastbin2 = firstbin;
 541               for(Int_t b=firstbin; b<=lastbin; b++, alpha1 += deltaalpha1, alpha2 += deltaalpha2)
 542                 {
 543                   Int_t binx1 = 1 + (Int_t)alpha1;
 544                   if(binx1<=lastbinx) {
 545                     UChar_t initialgapcount;
 546                     if(binx1>=firstbinx)
 547                       initialgapcount = fInitialGapCount[binx1 + b*nbinx];
 548                     else
 549                       initialgapcount = fInitialGapCount[firstbinx + b*nbinx];
 550                     if(initialgapcount != MAX_N_GAPS) {
 551                       if(!binfound1) {
 552                         firstbin1 = b;
 553                         binfound1 = kTRUE;
 554                       }
 555                       lastbin1 = b;
 556                     }
 557                   }
 558                   Int_t binx2 = 1 + (Int_t)alpha2;
 559                   if(binx2>=firstbin) {
 560                     UChar_t initialgapcount;
 561                     if(binx2<=lastbinx)
 562                       initialgapcount = fInitialGapCount[binx2 + b*nbinx];
 563                     else
 564                       initialgapcount = fInitialGapCount[lastbinx + b*nbinx];
 565                     if(initialgapcount != MAX_N_GAPS) {
 566                       if(!binfound2) {
 567                         firstbin2 = b;
 568                         binfound2 = kTRUE;
 569                       }
 570                       lastbin2 = b;
 571                     }
 572                   }
 573                 }
 574               fStartPadParams[i][pad].fFirstBin=firstbin1;
 575               fStartPadParams[i][pad].fLastBin=lastbin1;
 576               fEndPadParams[i][pad].fFirstBin=firstbin2;
 577               fEndPadParams[i][pad].fLastBin=lastbin2;
 578             }
 579         }
 580     }
 581
 582   //create lookup table for z of the digits
 583   if(!fLUTforwardZ)
 584     {
 585       Int_t ntimebins = AliHLTTPCTransform::GetNTimeBins();
 586       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 587         <<"Transformer: Allocating "<<ntimebins*sizeof(Float_t)<<" bytes to fLUTforwardZ"<<ENDLOG;
 588       fLUTforwardZ = new Float_t[ntimebins];
 589       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCHoughTransformerRow::CreateHistograms()","")
 590         <<"Transformer: Allocating "<<ntimebins*sizeof(Float_t)<<" bytes to fLUTbackwardZ"<<ENDLOG;
 591       fLUTbackwardZ = new Float_t[ntimebins];
 592       for(Int_t i=0; i<ntimebins; i++){
 593         Float_t z;
 594         z=AliHLTTPCTransform::GetZFast(0,i,GetZVertex());
 595         fLUTforwardZ[i]=z;
 596         z=AliHLTTPCTransform::GetZFast(18,i,GetZVertex());
 597         fLUTbackwardZ[i]=z;
 598       }
 599     }
 600 }
 601
 602 void AliHLTTPCHoughTransformerRow::Reset()
 603 {
 604   //Reset all the histograms. Should be done when processing new slice
 605   if(fLastTransformer) return;
 606   if(!fParamSpace)
 607     {
 608       LOG(AliHLTTPCLog::kWarning,"AliHLTTPCHoughTransformer::Reset","Histograms")
 609         <<"No histograms to reset"<<ENDLOG;
 610       return;
 611     }
 612
 613   for(Int_t i=0; i<GetNEtaSegments(); i++)
 614     fParamSpace[i]->Reset();
 615
 616 #ifdef do_mc
 617   {
 618     AliHLTTPCHistogram *hist = fParamSpace[0];
 619     Int_t ncellsx = (hist->GetNbinsX()+3)/2;
 620     Int_t ncellsy = (hist->GetNbinsY()+3)/2;
 621     Int_t ncells = ncellsx*ncellsy;
 622     for(Int_t i=0; i<GetNEtaSegments(); i++)
 623       memset(fTrackID[i],0,ncells*sizeof(AliHLTTrackIndex));
 624   }
 625 #endif
 626   AliHLTTPCHistogram *hist = fParamSpace[0];
 627   Int_t ncells = (hist->GetNbinsX()+2)*(hist->GetNbinsY()+2);
 628   for(Int_t i=0; i<GetNEtaSegments(); i++)
 629     {
 630       memcpy(fGapCount[i],fInitialGapCount,ncells*sizeof(UChar_t));
 631       memcpy(fCurrentRowCount[i],fTrackFirstRow,ncells*sizeof(UChar_t));
 632     }
 633 }
 634
 635 Int_t AliHLTTPCHoughTransformerRow::GetEtaIndex(Double_t eta) const
 636 {
 637   //Return the histogram index of the corresponding eta.
 638
 639   Double_t etaslice = (GetEtaMax() - GetEtaMin())/GetNEtaSegments();
 640   Double_t index = (eta-GetEtaMin())/etaslice;
 641   return (Int_t)index;
 642 }
 643
 644 inline AliHLTTPCHistogram *AliHLTTPCHoughTransformerRow::GetHistogram(Int_t etaindex)
 645 {
 646   // Return a pointer to the histogram which contains etaindex eta slice
 647   if(!fParamSpace || etaindex >= GetNEtaSegments() || etaindex < 0)
 648     return 0;
 649   if(!fParamSpace[etaindex])
 650     return 0;
 651   return fParamSpace[etaindex];
 652 }
 653
 654 Double_t AliHLTTPCHoughTransformerRow::GetEta(Int_t etaindex,Int_t /*slice*/) const
 655 {
 656   // Return eta calculated in the middle of the eta slice
 657   Double_t etaslice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments();
 658   Double_t eta=0;
 659   eta=(Double_t)((etaindex+0.5)*etaslice);
 660   return eta;
 661 }
 662
 663 void AliHLTTPCHoughTransformerRow::TransformCircle()
 664 {
 665   // This method contains the hough transformation
 666   // Depending on the option selected, it reads as an input
 667   // either the preloaded array with digits or directly raw
 668   // data stream (option fast_raw)
 669   if(GetDataPointer())
 670     TransformCircleFromDigitArray();
 671   else if(fTPCRawStream)
 672     TransformCircleFromRawStream();
 673 }
 674 void AliHLTTPCHoughTransformerRow::TransformCircleFromDigitArray()
 675 {
 676   //Do the Hough Transform
 677
 678   //Load the parameters used by the fast calculation of eta
 679   Double_t etaparam1 = GetEtaCalcParam1();
 680   Double_t etaparam2 = GetEtaCalcParam2();
 681   Double_t etaparam3 = GetEtaCalcParam3();
 682
 683   Int_t netasegments = GetNEtaSegments();
 684   Double_t etamin = GetEtaMin();
 685   Double_t etaslice = (GetEtaMax() - etamin)/netasegments;
 686
 687   Int_t lowerthreshold = GetLowerThreshold();
 688
 689   //Assumes that all the etaslice histos are the same!
 690   AliHLTTPCHistogram *h = fParamSpace[0];
 691   Int_t firstbiny = h->GetFirstYbin();
 692   Int_t firstbinx = h->GetFirstXbin();
 693   Int_t lastbinx = h->GetLastXbin();
 694   Int_t nbinx = h->GetNbinsX()+2;
 695
 696   UChar_t lastpad;
 697   Int_t lastetaindex=-1;
 698   AliHLTEtaRow *etaclust = new AliHLTEtaRow[netasegments];
 699
 700   AliHLTTPCDigitRowData *tempPt = GetDataPointer();
 701   if(!tempPt)
 702     {
 703       LOG(AliHLTTPCLog::kError,"AliHLTTPCHoughTransformer::TransformCircle","Data")
 704         <<"No input data "<<ENDLOG;
 705       return;
 706     }
 707
 708   Int_t ipatch = GetPatch();
 709   Int_t ilastpatch = GetLastPatch();
 710   Int_t islice = GetSlice();
 711   Float_t *lutz;
 712   if(islice < 18)
 713     lutz = fLUTforwardZ;
 714   else
 715     lutz = fLUTbackwardZ;
 716
 717   //Loop over the padrows:
 718   for(UChar_t i=AliHLTTPCTransform::GetFirstRow(ipatch); i<=AliHLTTPCTransform::GetLastRow(ipatch); i++)
 719     {
 720       lastpad = 255;
 721       //Flush eta clusters array
 722       memset(etaclust,0,netasegments*sizeof(AliHLTEtaRow));
 723
 724       Float_t radius=0;
 725
 726       //Get the data on this padrow:
 727       AliHLTTPCDigitData *digPt = tempPt->fDigitData;
 728       if((Int_t)i != (Int_t)tempPt->fRow)
 729         {
 730           LOG(AliHLTTPCLog::kError,"AliHLTTPCHoughTransformerRow::TransformCircle","Data")
 731             <<"AliHLTTPCHoughTransform::TransformCircle : Mismatching padrow numbering "<<(Int_t)i<<" "<<(Int_t)tempPt->fRow<<ENDLOG;
 732           continue;
 733         }
 734       //      cout<<" Starting row "<<i<<endl;
 735       //Loop over the data on this padrow:
 736       for(UInt_t j=0; j<tempPt->fNDigit; j++)
 737         {
 738           UShort_t charge = digPt[j].fCharge;
 739           if((Int_t)charge <= lowerthreshold)
 740             continue;
 741           UChar_t pad = digPt[j].fPad;
 742           UShort_t time = digPt[j].fTime;
 743
 744           if(pad != lastpad)
 745             {
 746               radius = fLUTr[(Int_t)i][(Int_t)pad];
 747               lastetaindex = -1;
 748             }
 749
 750           Float_t z = lutz[(Int_t)time];
 751           Double_t radiuscorr = radius*(1.+etaparam3*radius*radius);
 752           Double_t zovr = z/radiuscorr;
 753           Double_t eta = (etaparam1-etaparam2*fabs(zovr))*zovr;
 754           //Get the corresponding index, which determines which histogram to fill:
 755           Int_t etaindex = (Int_t)((eta-etamin)/etaslice);
 756
 757 #ifndef do_mc
 758           if(etaindex == lastetaindex) continue;
 759 #endif
 760           //      cout<<" Digit at patch "<<ipatch<<" row "<<i<<" pad "<<(Int_t)pad<<" time "<<time<<" etaslice "<<etaindex<<endl;
 761
 762           if(etaindex < 0 || etaindex >= netasegments)
 763             continue;
 764
 765           if(!etaclust[etaindex].fIsFound)
 766             {
 767               etaclust[etaindex].fStartPad = pad;
 768               etaclust[etaindex].fEndPad = pad;
 769               etaclust[etaindex].fIsFound = 1;
 770 #ifdef do_mc
 771               FillClusterMCLabels(digPt[j],&etaclust[etaindex]);
 772 #endif
 773               continue;
 774             }
 775           else
 776             {
 777               if(pad <= (etaclust[etaindex].fEndPad + 1))
 778                 {
 779                   etaclust[etaindex].fEndPad = pad;
 780 #ifdef do_mc
 781                   FillClusterMCLabels(digPt[j],&etaclust[etaindex]);
 782 #endif
 783                 }
 784               else
 785                 {
 786                   FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
 787
 788                   etaclust[etaindex].fStartPad = pad;
 789                   etaclust[etaindex].fEndPad = pad;
 790
 791 #ifdef do_mc
 792                   memset(etaclust[etaindex].fMcLabels,0,MaxTrack);
 793                   FillClusterMCLabels(digPt[j],&etaclust[etaindex]);
 794 #endif
 795                 }
 796             }
 797           lastpad = pad;
 798           lastetaindex = etaindex;
 799         }
 800       //Write remaining clusters
 801       for(Int_t etaindex = 0;etaindex < netasegments;etaindex++)
 802         {
 803           //Check for empty row
 804           if((etaclust[etaindex].fStartPad == 0) && (etaclust[etaindex].fEndPad == 0)) continue;
 805
 806           FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
 807         }
 808
 809       //Move the data pointer to the next padrow:
 810       AliHLTTPCMemHandler::UpdateRowPointer(tempPt);
 811     }
 812
 813   delete [] etaclust;
 814 }
 815
 816 void AliHLTTPCHoughTransformerRow::TransformCircleFromRawStream()
 817 {
 818   //Do the Hough Transform
 819
 820   //Load the parameters used by the fast calculation of eta
 821   Double_t etaparam1 = GetEtaCalcParam1();
 822   Double_t etaparam2 = GetEtaCalcParam2();
 823   Double_t etaparam3 = GetEtaCalcParam3();
 824
 825   Int_t netasegments = GetNEtaSegments();
 826   Double_t etamin = GetEtaMin();
 827   Double_t etaslice = (GetEtaMax() - etamin)/netasegments;
 828
 829   Int_t lowerthreshold = GetLowerThreshold();
 830
 831   //Assumes that all the etaslice histos are the same!
 832   AliHLTTPCHistogram *h = fParamSpace[0];
 833   Int_t firstbiny = h->GetFirstYbin();
 834   Int_t firstbinx = h->GetFirstXbin();
 835   Int_t lastbinx = h->GetLastXbin();
 836   Int_t nbinx = h->GetNbinsX()+2;
 837
 838   Int_t lastetaindex = -1;
 839   AliHLTEtaRow *etaclust = new AliHLTEtaRow[netasegments];
 840
 841   if(!fTPCRawStream)
 842     {
 843       LOG(AliHLTTPCLog::kError,"AliHLTTPCHoughTransformer::TransformCircle","Data")
 844         <<"No input data "<<ENDLOG;
 845       return;
 846     }
 847
 848   Int_t ipatch = GetPatch();
 849   UChar_t rowmin = AliHLTTPCTransform::GetFirstRowOnDDL(ipatch);
 850   UChar_t rowmax = AliHLTTPCTransform::GetLastRowOnDDL(ipatch);
 851   //  Int_t ntimebins = AliHLTTPCTransform::GetNTimeBins();
 852   Int_t ilastpatch = GetLastPatch();
 853   Int_t islice = GetSlice();
 854   Float_t *lutz;
 855   if(islice < 18)
 856     lutz = fLUTforwardZ;
 857   else
 858     lutz = fLUTbackwardZ;
 859
 860   //Flush eta clusters array
 861   memset(etaclust,0,netasegments*sizeof(AliHLTEtaRow));
 862
 863   UChar_t i=0;
 864   Int_t npads=0;
 865   Float_t radius=0;
 866   UChar_t pad=0;
 867
 868   //Loop over the rawdata:
 869   while (fTPCRawStream->Next()) {
 870
 871     if(fTPCRawStream->IsNewSector() || fTPCRawStream->IsNewRow()) {
 872
 873       //Write remaining clusters
 874       for(Int_t etaindex = 0;etaindex < netasegments;etaindex++)
 875         {
 876           //Check for empty row
 877           if((etaclust[etaindex].fStartPad == 0) && (etaclust[etaindex].fEndPad == 0)) continue;
 878
 879           FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
 880         }
 881
 882       Int_t sector=fTPCRawStream->GetSector();
 883       Int_t row=fTPCRawStream->GetRow();
 884       Int_t slice,srow;
 885       AliHLTTPCTransform::Sector2Slice(slice,srow,sector,row);
 886       if(slice!=islice){
 887         LOG(AliHLTTPCLog::kError,"AliHLTDDLDataFileHandler::DDLDigits2Memory","Slice")
 888           <<AliHLTTPCLog::kDec<<"Found slice "<<slice<<", expected "<<islice<<ENDLOG;
 889         continue;
 890       }
 891
 892       i=(UChar_t)srow;
 893       npads = AliHLTTPCTransform::GetNPads(srow)-1;
 894
 895       //Flush eta clusters array
 896       memset(etaclust,0,netasegments*sizeof(AliHLTEtaRow));
 897
 898       radius=0;
 899
 900     }
 901
 902     if((i<rowmin)||(i>rowmax))continue;
 903
 904     //    cout<<" Starting row "<<(UInt_t)i<<endl;
 905     //Loop over the data on this padrow:
 906     if(fTPCRawStream->IsNewRow() || fTPCRawStream->IsNewPad()) {
 907       pad=fTPCRawStream->GetPad();
 908       /*
 909       if((pad<0)||(pad>=(npads+1))){
 910         LOG(AliHLTTPCLog::kError,"AliHLTDDLDataFileHandler::DDLDigits2Memory","Pad")
 911           <<AliHLTTPCLog::kDec<<"Pad value out of bounds "<<pad<<" "
 912           <<npads+1<<ENDLOG;
 913         continue;
 914       }
 915       */
 916       radius = fLUTr[(Int_t)i][(Int_t)pad];
 917       lastetaindex = -1;
 918     }
 919
 920     UShort_t time=fTPCRawStream->GetTime();
 921     /*
 922     if((time<0)||(time>=ntimebins)){
 923       LOG(AliHLTTPCLog::kError,"AliHLTDDLDataFileHandler::DDLDigits2Memory","Time")
 924         <<AliHLTTPCLog::kDec<<"Time out of bounds "<<time<<" "
 925         <<AliHLTTPCTransform::GetNTimeBins()<<ENDLOG;
 926       continue;
 927     }
 928     */
 929
 930     if(fTPCRawStream->GetSignal() <= lowerthreshold)
 931       continue;
 932
 933     Float_t z = lutz[(Int_t)time];
 934     Double_t radiuscorr = radius*(1.+etaparam3*radius*radius);
 935     Double_t zovr = z/radiuscorr;
 936     Double_t eta = (etaparam1-etaparam2*fabs(zovr))*zovr;
 937     //Get the corresponding index, which determines which histogram to fill:
 938     Int_t etaindex = (Int_t)((eta-etamin)/etaslice);
 939
 940 #ifndef do_mc
 941     if(etaindex == lastetaindex) continue;
 942 #endif
 943     //    cout<<" Digit at patch "<<ipatch<<" row "<<i<<" pad "<<(Int_t)pad<<" time "<<time<<" etaslice "<<etaindex<<endl;
 944
 945     if(etaindex < 0 || etaindex >= netasegments)
 946       continue;
 947
 948     if(!etaclust[etaindex].fIsFound)
 949       {
 950         etaclust[etaindex].fStartPad = pad;
 951         etaclust[etaindex].fEndPad = pad;
 952         etaclust[etaindex].fIsFound = 1;
 953         continue;
 954       }
 955     else
 956       {
 957         if(pad <= (etaclust[etaindex].fEndPad + 1))
 958           {
 959             etaclust[etaindex].fEndPad = pad;
 960           }
 961         else
 962           {
 963             FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
 964
 965             etaclust[etaindex].fStartPad = pad;
 966             etaclust[etaindex].fEndPad = pad;
 967
 968           }
 969       }
 970     lastetaindex = etaindex;
 971   }
 972
 973   //Write remaining clusters
 974   for(Int_t etaindex = 0;etaindex < netasegments;etaindex++)
 975     {
 976       //Check for empty row
 977       if((etaclust[etaindex].fStartPad == 0) && (etaclust[etaindex].fEndPad == 0)) continue;
 978
 979       FillCluster(i,etaindex,etaclust,ilastpatch,firstbinx,lastbinx,nbinx,firstbiny);
 980     }
 981
 982   delete [] etaclust;
 983 }
 984
 985 #ifndef do_mc
 986 Int_t AliHLTTPCHoughTransformerRow::GetTrackID(Int_t /*etaindex*/,Double_t /*alpha1*/,Double_t /*alpha2*/) const
 987 {
 988   // Does nothing if do_mc undefined
 989   LOG(AliHLTTPCLog::kWarning,"AliHLTTPCHoughTransformerRow::GetTrackID","Data")
 990     <<"Flag switched off"<<ENDLOG;
 991   return -1;
 992 #else
 993 Int_t AliHLTTPCHoughTransformerRow::GetTrackID(Int_t etaindex,Double_t alpha1,Double_t alpha2) const
 994 {
 995   // Returns the MC label for a given peak found in the Hough space
 996   if(etaindex < 0 || etaindex > GetNEtaSegments())
 997     {
 998       LOG(AliHLTTPCLog::kWarning,"AliHLTTPCHoughTransformerRow::GetTrackID","Data")
 999         <<"Wrong etaindex "<<etaindex<<ENDLOG;
1000       return -1;
1001     }
1002   AliHLTTPCHistogram *hist = fParamSpace[etaindex];
1003   Int_t bin = hist->FindLabelBin(alpha1,alpha2);
1004   if(bin==-1) {
1005     LOG(AliHLTTPCLog::kWarning,"AliHLTTPCHoughTransformerRow::GetTrackID()","")
1006       <<"Track candidate outside Hough space boundaries: "<<alpha1<<" "<<alpha2<<ENDLOG;
1007     return -1;
1008   }
1009   Int_t label=-1;
1010   Int_t max=0;
1011   for(UInt_t i=0; i<(MaxTrack-1); i++)
1012     {
1013       Int_t nhits=fTrackID[etaindex][bin].fNHits[i];
1014       if(nhits == 0) break;
1015       if(nhits > max)
1016         {
1017           max = nhits;
1018           label = fTrackID[etaindex][bin].fLabel[i];
1019         }
1020     }
1021   Int_t label2=-1;
1022   Int_t max2=0;
1023   for(UInt_t i=0; i<(MaxTrack-1); i++)
1024     {
1025       Int_t nhits=fTrackID[etaindex][bin].fNHits[i];
1026       if(nhits == 0) break;
1027       if(nhits > max2)
1028         {
1029           if(fTrackID[etaindex][bin].fLabel[i]!=label) {
1030             max2 = nhits;
1031             label2 = fTrackID[etaindex][bin].fLabel[i];
1032           }
1033         }
1034     }
1035   if(max2 !=0 ) {
1036     LOG(AliHLTTPCLog::kDebug,"AliHLTTPCHoughTransformerRow::GetTrackID()","")
1037       <<" TrackID"<<" label "<<label<<" max "<<max<<" label2 "<<label2<<" max2 "<<max2<<" "<<(Float_t)max2/(Float_t)max<<" "<<fTrackID[etaindex][bin].fLabel[MaxTrack-1]<<" "<<(Int_t)fTrackID[etaindex][bin].fNHits[MaxTrack-1]<<ENDLOG;
1038   }
1039   return label;
1040 #endif
1041 }
1042
1043 Int_t AliHLTTPCHoughTransformerRow::GetTrackLength(Double_t alpha1,Double_t alpha2,Int_t *rows) const
1044 {
1045   // Returns the track length for a given peak found in the Hough space
1046
1047   AliHLTTPCHistogram *hist = fParamSpace[0];
1048   Int_t bin = hist->FindBin(alpha1,alpha2);
1049   if(bin==-1) {
1050     LOG(AliHLTTPCLog::kWarning,"AliHLTTPCHoughTransformerRow::GetTrackLength()","")
1051       <<"Track candidate outside Hough space boundaries: "<<alpha1<<" "<<alpha2<<ENDLOG;
1052     return -1;
1053   }
1054   rows[0] = fTrackFirstRow[bin];
1055   rows[1] = fTrackLastRow[bin];
1056
1057   return 0;
1058 }
1059
1060 inline void AliHLTTPCHoughTransformerRow::FillClusterRow(UChar_t i,Int_t binx1,Int_t binx2,UChar_t *ngaps2,UChar_t *currentrow2,UChar_t *lastrow2
1061 #ifdef do_mc
1062                            ,AliHLTEtaRow etaclust,AliHLTTrackIndex *trackid
1063 #endif
1064                            )
1065 {
1066   // The method is a part of the fast hough transform.
1067   // It fills one row of the hough space.
1068   // It is called by FillCluster() method inside the
1069   // loop over alpha2 bins
1070
1071   for(Int_t bin=binx1;bin<=binx2;bin++)
1072     {
1073       if(ngaps2[bin] < MAX_N_GAPS) {
1074         if(i < lastrow2[bin] && i > currentrow2[bin]) {
1075           ngaps2[bin] += (i-currentrow2[bin]-1);
1076           currentrow2[bin]=i;
1077         }
1078 #ifdef do_mc
1079         if(i < lastrow2[bin] && i >= currentrow2[bin]) {
1080           for(UInt_t t=0;t<(MaxTrack-1); t++)
1081             {
1082               Int_t label = etaclust.fMcLabels[t];
1083               if(label == 0) break;
1084               UInt_t c;
1085               Int_t tempbin2 = (Int_t)(bin/2);
1086               for(c=0; c<(MaxTrack-1); c++)
1087                 if(trackid[tempbin2].fLabel[c] == label || trackid[tempbin2].fNHits[c] == 0)
1088                   break;
1089               trackid[tempbin2].fLabel[c] = label;
1090               if(trackid[tempbin2].fCurrentRow[c] != i) {
1091                 trackid[tempbin2].fNHits[c]++;
1092                 trackid[tempbin2].fCurrentRow[c] = i;
1093               }
1094             }
1095         }
1096 #endif
1097       }
1098     }
1099
1100 }
1101
1102 inline void AliHLTTPCHoughTransformerRow::FillCluster(UChar_t i,Int_t etaindex,AliHLTEtaRow *etaclust,Int_t ilastpatch,Int_t firstbinx,Int_t lastbinx,Int_t nbinx,Int_t firstbiny)
1103 {
1104   // The method is a part of the fast hough transform.
1105   // It fills a TPC cluster into the hough space.
1106
1107   UChar_t *ngaps = fGapCount[etaindex];
1108   UChar_t *currentrow = fCurrentRowCount[etaindex];
1109   UChar_t *lastrow = fTrackLastRow;
1110   UChar_t *prevbin = fPrevBin[etaindex];
1111   UChar_t *nextbin = fNextBin[etaindex];
1112   UChar_t *nextrow = fNextRow[etaindex];
1113 #ifdef do_mc
1114   AliHLTTrackIndex *trackid = fTrackID[etaindex];
1115 #endif
1116
1117   //Do the transformation:
1118   AliHLTPadHoughParams *startparams = &fStartPadParams[(Int_t)i][etaclust[etaindex].fStartPad];
1119   AliHLTPadHoughParams *endparams = &fEndPadParams[(Int_t)i][etaclust[etaindex].fEndPad];
1120
1121   Float_t alpha1 = startparams->fAlpha;
1122   Float_t deltaalpha1 = startparams->fDeltaAlpha;
1123   Float_t alpha2 = endparams->fAlpha;
1124   Float_t deltaalpha2 = endparams->fDeltaAlpha;
1125
1126   Int_t firstbin1 = startparams->fFirstBin;
1127   Int_t firstbin2 = endparams->fFirstBin;
1128   Int_t firstbin = firstbin1;
1129   if(firstbin>firstbin2) firstbin = firstbin2;
1130
1131   Int_t lastbin1 = startparams->fLastBin;
1132   Int_t lastbin2 = endparams->fLastBin;
1133   Int_t lastbin = lastbin1;
1134   if(lastbin<lastbin2) lastbin = lastbin2;
1135
1136   alpha1 += (firstbin-firstbiny)*deltaalpha1;
1137   alpha2 += (firstbin-firstbiny)*deltaalpha2;
1138
1139   //Fill the histogram along the alpha2 range
1140   if(ilastpatch == -1) {
1141     for(Int_t b=firstbin; b<=lastbin; b++, alpha1 += deltaalpha1, alpha2 += deltaalpha2)
1142       {
1143         Int_t binx1 = 1 + (Int_t)alpha1;
1144         if(binx1>lastbinx) continue;
1145         if(binx1<firstbinx) binx1 = firstbinx;
1146         Int_t binx2 = 1 + (Int_t)alpha2;
1147         if(binx2<firstbinx) continue;
1148         if(binx2>lastbinx) binx2 = lastbinx;
1149 #ifdef do_mc
1150         if(binx2<binx1) {
1151           LOG(AliHLTTPCLog::kWarning,"AliHLTTPCHoughTransformerRow::TransformCircle()","")
1152             <<"Wrong filling "<<binx1<<" "<<binx2<<" "<<i<<" "<<etaclust[etaindex].fStartPad<<" "<<etaclust[etaindex].fEndPad<<ENDLOG;
1153         }
1154 #endif
1155         Int_t tempbin = b*nbinx;
1156         UChar_t *ngaps2 = ngaps + tempbin;
1157         UChar_t *currentrow2 = currentrow + tempbin;
1158         UChar_t *lastrow2 = lastrow + tempbin;
1159 #ifdef do_mc
1160         Int_t tempbin2 = ((Int_t)(b/2))*((Int_t)((nbinx+1)/2));
1161         AliHLTTrackIndex *trackid2 = trackid + tempbin2;
1162 #endif
1163         FillClusterRow(i,binx1,binx2,ngaps2,currentrow2,lastrow2
1164 #ifdef do_mc
1165                        ,etaclust[etaindex],trackid2
1166 #endif
1167                        );
1168       }
1169   }
1170   else {
1171     for(Int_t b=firstbin; b<=lastbin; b++, alpha1 += deltaalpha1, alpha2 += deltaalpha2)
1172       {
1173         Int_t binx1 = 1 + (Int_t)alpha1;
1174         if(binx1>lastbinx) continue;
1175         if(binx1<firstbinx) binx1 = firstbinx;
1176         Int_t binx2 = 1 + (Int_t)alpha2;
1177         if(binx2<firstbinx) continue;
1178         if(binx2>lastbinx) binx2 = lastbinx;
1179 #ifdef do_mc
1180         if(binx2<binx1) {
1181           LOG(AliHLTTPCLog::kWarning,"AliHLTTPCHoughTransformerRow::TransformCircle()","")
1182             <<"Wrong filling "<<binx1<<" "<<binx2<<" "<<i<<" "<<etaclust[etaindex].fStartPad<<" "<<etaclust[etaindex].fEndPad<<ENDLOG;
1183         }
1184 #endif
1185         if(nextrow[b] > b) {
1186           Int_t deltab = (nextrow[b] - b - 1);
1187           b += deltab;
1188           alpha1 += deltaalpha1*deltab;
1189           alpha2 += deltaalpha2*deltab;
1190           continue;
1191         }
1192         Int_t tempbin = b*nbinx;
1193         binx1 = (UInt_t)nextbin[binx1 + tempbin];
1194         binx2 = (UInt_t)prevbin[binx2 + tempbin];
1195         if(binx2<binx1) continue;
1196         UChar_t *ngaps2 = ngaps + tempbin;
1197         UChar_t *currentrow2 = currentrow + tempbin;
1198         UChar_t *lastrow2 = lastrow + tempbin;
1199 #ifdef do_mc
1200         Int_t tempbin2 = ((Int_t)(b/2))*((Int_t)((nbinx+1)/2));
1201         AliHLTTrackIndex *trackid2 = trackid + tempbin2;
1202 #endif
1203         FillClusterRow(i,binx1,binx2,ngaps2,currentrow2,lastrow2
1204 #ifdef do_mc
1205                        ,etaclust[etaindex],trackid2
1206 #endif
1207                        );
1208       }
1209   }
1210 }
1211
1212 #ifdef do_mc
1213 inline void AliHLTTPCHoughTransformerRow::FillClusterMCLabels(AliHLTTPCDigitData digpt,AliHLTEtaRow *etaclust)
1214 {
1215   // The method is a part of the fast hough transform.
1216   // It fills the MC labels of a TPC cluster into a
1217   // special hough space array.
1218   for(Int_t t=0; t<3; t++)
1219     {
1220       Int_t label = digpt.fTrackID[t];
1221       if(label < 0) break;
1222       UInt_t c;
1223       for(c=0; c<(MaxTrack-1); c++)
1224         if(etaclust->fMcLabels[c] == label || etaclust->fMcLabels[c] == 0)
1225           break;
1226
1227       etaclust->fMcLabels[c] = label;
1228     }
1229 }
1230 #endif
1231
1232 void AliHLTTPCHoughTransformerRow::SetTransformerArrays(AliHLTTPCHoughTransformerRow *tr)
1233 {
1234   // In case of sequential filling of the hough space, the method is used to
1235   // transmit the pointers to the hough arrays from one transformer to the
1236   // following one.
1237
1238     fGapCount = tr->fGapCount;
1239     fCurrentRowCount = tr->fCurrentRowCount;
1240 #ifdef do_mc
1241     fTrackID = tr->fTrackID;
1242 #endif
1243     fTrackNRows = tr->fTrackNRows;
1244     fTrackFirstRow = tr->fTrackFirstRow;
1245     fTrackLastRow = tr->fTrackLastRow;
1246     fInitialGapCount = tr->fInitialGapCount;
1247
1248     fPrevBin = tr->fPrevBin;
1249     fNextBin = tr->fNextBin;
1250     fNextRow = tr->fNextRow;
1251
1252     fStartPadParams = tr->fStartPadParams;
1253     fEndPadParams = tr->fEndPadParams;
1254     fLUTr = tr->fLUTr;
1255     fLUTforwardZ = tr->fLUTforwardZ;
1256     fLUTbackwardZ = tr->fLUTbackwardZ;
1257
1258     fParamSpace = tr->fParamSpace;
1259
1260     return;
1261 }