HLT/TPCLib/AliHLTTPCClusterFinder.cxx

   1 // @(#) $Id$
   2 // Original: AliL3ClustFinderNew.cxx,v 1.29 2005/06/14 10:55:21 cvetan Exp
   3
   4 // Author: Anders Vestbo <mailto:vestbo@fi.uib.no>,
   5 //         Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de>
   6 //         Jochen Thaeder <mailto:thaeder@kip.uni-heidelberg.de>
   7
   8 //*-- Copyright &copy ALICE HLT Group
   9
  10
  11 #include "AliHLTTPCDigitReader.h"
  12 #include "AliHLTTPCRootTypes.h"
  13 #include "AliHLTTPCLogging.h"
  14 #include "AliHLTTPCClusterFinder.h"
  15 #include "AliHLTTPCDigitData.h"
  16 #include "AliHLTTPCTransform.h"
  17 #include "AliHLTTPCSpacePointData.h"
  18 #include "AliHLTTPCMemHandler.h"
  19 #include "AliHLTTPCPad.h"
  20
  21 #if __GNUC__ >= 3
  22 using namespace std;
  23 #endif
  24
  25 /** \class AliHLTTPCClusterFinder
  26 <pre>
  27 //_____________________________________________________________
  28 // AliHLTTPCClusterFinder
  29 //
  30 // The current cluster finder for HLT
  31 // (Based on STAR L3)
  32 //
  33 // The cluster finder is initialized with the Init function,
  34 // providing the slice and patch information to work on.
  35 //
  36 // The input is a provided by the AliHLTTPCDigitReader class,
  37 // using the init() funktion, and the next() funktion in order
  38 // to get the next bin. Either packed or unpacked data can be
  39 // processed, dependent if one uses AliHLTTPCDigitReaderPacked
  40 // class or AliHLTTPCDigitReaderUnpacked class in the
  41 // Clusterfinder Component.
  42 // The resulting space points will be in the
  43 // array given by the SetOutputArray function.
  44 //
  45 // There are several setters which control the behaviour:
  46 //
  47 // - SetXYError(Float_t):   set fixed error in XY direction
  48 // - SetZError(Float_t):    set fixed error in Z  direction
  49 //                            (used if errors are not calculated)
  50 // - SetDeconv(Bool_t):     switch on/off deconvolution
  51 // - SetThreshold(UInt_t):  set charge threshold for cluster
  52 // - SetMatchWidth(UInt_t): set the match distance in
  53 //                            time for sequences to be merged
  54 // - SetSTDOutput(Bool_t):  switch on/off output about found clusters
  55 // - SetCalcErr(Bool_t):    switch on/off calculation of
  56 //                          space point errors (or widths in raw system)
  57 // - SetRawSP(Bool_t):      switch on/off convertion to raw system
  58 //
  59 //
  60 // Example Usage:
  61 //
  62 // AliHLTTPCFileHandler *file = new AliHLTTPCFileHandler();
  63 // file->SetAliInput(digitfile); //give some input file
  64 // for(int slice=0; slice<=35; slice++){
  65 //   for(int patch=0; pat<6; pat++){
  66 //     file->Init(slice,patch);
  67 //     UInt_t ndigits=0;
  68 //     UInt_t maxclusters=100000;
  69 //     UInt_t pointsize = maxclusters*sizeof(AliHLTTPCSpacePointData);
  70 //     AliHLTTPCSpacePointData *points = (AliHLTTPCSpacePointData*)memory->Allocate(pointsize);
  71 //     AliHLTTPCDigitRowData *digits = (AliHLTTPCDigitRowData*)file->AliAltroDigits2Memory(ndigits,event);
  72 //     AliHLTTPCClusterFinder *cf = new AliHLTTPCClusterFinder();
  73 //     cf->SetMatchWidth(2);
  74 //     cf->InitSlice( slice, patch, row[0], row[1], maxPoints );
  75 //     cf->SetSTDOutput(kTRUE);    //Some output to standard IO
  76 //     cf->SetRawSP(kFALSE);       //Convert space points to local system
  77 //     cf->SetThreshold(5);        //Threshold of cluster charge
  78 //     cf->SetDeconv(kTRUE);       //Deconv in pad and time direction
  79 //     cf->SetCalcErr(kTRUE);      //Calculate the errors of the spacepoints
  80 //     cf->SetOutputArray(points); //Move the spacepoints to the array
  81 //     cf->Read(iter->fPtr, iter->fSize ); //give the data to the cf
  82 //     cf->ProcessDigits();        //process the rows given by init
  83 //     Int_t npoints = cf->GetNumberOfClusters();
  84 //     AliHLTTPCMemHandler *out= new AliHLTTPCMemHandler();
  85 //     out->SetBinaryOutput(fname);
  86 //     out->Memory2Binary(npoints,points); //store the spacepoints
  87 //     out->CloseBinaryOutput();
  88 //     delete out;
  89 //     file->free();
  90 //     delete cf;
  91 //   }
  92 // }
  93 </pre>
  94 */
  95
  96 ClassImp(AliHLTTPCClusterFinder)
  97
  98 AliHLTTPCClusterFinder::AliHLTTPCClusterFinder()
  99   :
 100   fMatch(1),
 101   fThreshold(10),
 102   fXYErr(0.2),
 103   fZErr(0.3),
 104   fDeconvPad(kTRUE),
 105   fDeconvTime(kTRUE),
 106   fStdout(kFALSE),
 107   fCalcerr(kTRUE),
 108   fRawSP(kFALSE),
 109   fFirstRow(0),
 110   fLastRow(0),
 111   fDigitReader(NULL)
 112 {
 113   //constructor
 114 }
 115
 116 AliHLTTPCClusterFinder::AliHLTTPCClusterFinder(const AliHLTTPCClusterFinder& src)
 117   :
 118   fMatch(src.fMatch),
 119   fThreshold(src.fThreshold),
 120   fXYErr(src.fXYErr),
 121   fZErr(src.fZErr),
 122   fDeconvPad(src.fDeconvPad),
 123   fDeconvTime(src.fDeconvTime),
 124   fStdout(src.fStdout),
 125   fCalcerr(src.fCalcerr),
 126   fRawSP(src.fRawSP),
 127   fFirstRow(src.fFirstRow),
 128   fLastRow(src.fLastRow),
 129   fDigitReader(src.fDigitReader)
 130 {
 131 }
 132
 133 AliHLTTPCClusterFinder& AliHLTTPCClusterFinder::operator=(const AliHLTTPCClusterFinder& src)
 134 {
 135   fMatch=src.fMatch;
 136   fThreshold=src.fThreshold;
 137   fXYErr=src.fXYErr;
 138   fZErr=src.fZErr;
 139   fDeconvPad=src.fDeconvPad;
 140   fDeconvTime=src.fDeconvTime;
 141   fStdout=src.fStdout;
 142   fCalcerr=src.fCalcerr;
 143   fRawSP=src.fRawSP;
 144   fFirstRow=src.fFirstRow;
 145   fLastRow=src.fLastRow;
 146   fDigitReader=src.fDigitReader;
 147   return (*this);
 148 }
 149
 150 AliHLTTPCClusterFinder::~AliHLTTPCClusterFinder()
 151 {
 152   //destructor
 153 }
 154
 155 void AliHLTTPCClusterFinder::InitSlice(Int_t slice,Int_t patch,Int_t firstrow, Int_t lastrow,Int_t nmaxpoints)
 156 {
 157   //init slice
 158   fNClusters = 0;
 159   fMaxNClusters = nmaxpoints;
 160   fCurrentSlice = slice;
 161   fCurrentPatch = patch;
 162   fFirstRow = firstrow;
 163   fLastRow = lastrow;
 164 }
 165
 166 void AliHLTTPCClusterFinder::InitSlice(Int_t slice,Int_t patch,Int_t nmaxpoints)
 167 {
 168   //init slice
 169   fNClusters = 0;
 170   fMaxNClusters = nmaxpoints;
 171   fCurrentSlice = slice;
 172   fCurrentPatch = patch;
 173   fFirstRow=AliHLTTPCTransform::GetFirstRow(patch);
 174   fLastRow=AliHLTTPCTransform::GetLastRow(patch);
 175 }
 176
 177 void AliHLTTPCClusterFinder::SetOutputArray(AliHLTTPCSpacePointData *pt)
 178 {
 179   //set pointer to output
 180   fSpacePointData = pt;
 181 }
 182
 183 void AliHLTTPCClusterFinder::Read(void* ptr,unsigned long size){
 184   //set input pointer
 185   fPtr = (UChar_t*)ptr;
 186   fSize = size;
 187 }
 188
 189 void AliHLTTPCClusterFinder::ProcessDigits()
 190 {
 191   bool readValue = true;
 192   Int_t newRow = 0;
 193   Int_t rowOffset = 0;
 194   UShort_t time=0,newTime=0;
 195   UInt_t pad=0,newPad=0;
 196   AliHLTTPCSignal_t charge=0;
 197
 198   fNClusters = 0;
 199
 200   // initialize block for reading packed data
 201   fDigitReader->InitBlock(fPtr,fSize,fFirstRow,fLastRow,fCurrentPatch,fCurrentSlice);
 202   readValue = fDigitReader->Next();
 203
 204   if (!readValue)return;
 205
 206   pad = fDigitReader->GetPad();
 207   time = fDigitReader->GetTime();
 208   fCurrentRow = fDigitReader->GetRow();
 209
 210   if ( fCurrentPatch >= 2 ) // Outer sector, patches 2, 3, 4, 5
 211     rowOffset = AliHLTTPCTransform::GetFirstRow( 2 );
 212
 213   fCurrentRow += rowOffset;
 214
 215   UInt_t lastpad = 123456789;
 216   AliClusterData *pad1[5000]; //2 lists for internal memory=2pads
 217   AliClusterData *pad2[5000]; //2 lists for internal memory=2pads
 218   AliClusterData clusterlist[10000]; //Clusterlist
 219
 220   AliClusterData **currentPt;  //List of pointers to the current pad
 221   AliClusterData **previousPt; //List of pointers to the previous pad
 222   currentPt = pad2;
 223   previousPt = pad1;
 224   UInt_t nprevious=0,ncurrent=0,ntotal=0;
 225
 226   /* quick implementation of baseline calculation and zero suppression
 227      open a pad object for each pad and delete it after processing.
 228      later a list of pad objects with base line history can be used
 229      The whole thing only works if we really get unprocessed raw data, if
 230      the data is already zero suppressed, there might be gaps in the time
 231      bins.
 232    */
 233   Int_t gatingGridOffset=50;
 234   AliHLTTPCPad baseline(gatingGridOffset, AliHLTTPCTransform::GetNTimeBins());
 235   // just to make later conversion to a list of objects easier
 236   AliHLTTPCPad* pCurrentPad=&baseline;
 237
 238   while ( readValue ){   // Reads through all digits in block
 239
 240     if(pad != lastpad){
 241       //This is a new pad
 242
 243       //Switch the lists:
 244       if(currentPt == pad2){
 245         currentPt = pad1;
 246         previousPt = pad2;
 247       }
 248       else {
 249         currentPt = pad2;
 250         previousPt = pad1;
 251       }
 252       nprevious = ncurrent;
 253       ncurrent = 0;
 254       if(pad != lastpad+1){
 255         //this happens if there is a pad with no signal.
 256         nprevious = ncurrent = 0;
 257       }
 258       lastpad = pad;
 259     }
 260
 261     Bool_t newcluster = kTRUE;
 262     UInt_t seqcharge=0,seqaverage=0,seqerror=0;
 263     UInt_t lastcharge=0,lastwas_falling=0;
 264     Int_t newbin=-1;
 265
 266
 267     if(fDeconvTime){
 268       redo: //This is a goto.
 269
 270       if(newbin > -1){
 271         //bin = newbin;
 272         newbin = -1;
 273       }
 274
 275       lastcharge=0;
 276       lastwas_falling = 0;
 277     }
 278
 279     while(1){ //Loop over time bins of current pad
 280       // read all the values for one pad at once to calculate the base line
 281       if (pCurrentPad) {
 282         if (!pCurrentPad->IsStarted()) {
 283           //HLTDebug("reading data for pad %d, padrow %d", fDigitReader->GetPad(), fDigitReader->GetRow()+rowOffset);
 284           pCurrentPad->SetID(fDigitReader->GetRow()+rowOffset,fDigitReader->GetPad());
 285           if ((pCurrentPad->StartEvent())>=0) {
 286             do {
 287               if ((fDigitReader->GetRow()+rowOffset)!=pCurrentPad->GetRowNumber()) break;
 288               if (fDigitReader->GetPad()!=pCurrentPad->GetPadNumber()) break;
 289               pCurrentPad->SetRawData(fDigitReader->GetTime(), fDigitReader->GetSignal());
 290               //HLTDebug("set raw data to pad: bin %d charge %d", fDigitReader->GetTime(), fDigitReader->GetSignal());
 291             } while ((readValue = fDigitReader->Next())!=0);
 292           }
 293           pCurrentPad->CalculateBaseLine(AliHLTTPCTransform::GetNTimeBins()/2);
 294           if (pCurrentPad->Next(kTRUE/*do zero suppression*/)==0) {
 295             HLTDebug("no data available after zero suppression");
 296             pCurrentPad->StopEvent();
 297             pCurrentPad->ResetHistory();
 298             break;
 299           }
 300           time=pCurrentPad->GetCurrentPosition();
 301           if (time>pCurrentPad->GetSize()) {
 302             HLTError("invalid time bin for pad");
 303             break;
 304           }
 305         }
 306       }
 307
 308       if (pCurrentPad) {
 309         charge = pCurrentPad->GetCorrectedData();
 310       } else {
 311       charge = fDigitReader->GetSignal();
 312       }
 313       //HLTDebug("get next charge value: position %d charge %d", time, charge);
 314
 315
 316       // CHARGE DEBUG
 317       if (fDigitReader->GetRow() == 90){
 318 /////     LOG(AliHLTTPCLog::kFatal,"AliHLTTPCClusterFinder::Row","row90")  << "PAD=" <<  fDigitReader->GetPad() << "  TIME=" <<  fDigitReader->GetTime()
 319           //                                       << "  SIGNAL=" <<  fDigitReader->GetSignal() << ENDLOG;
 320
 321       }
 322
 323       if(time >= AliHLTTPCTransform::GetNTimeBins()){
 324         HLTWarning("Timebin (%d) out of range (%d)", time, AliHLTTPCTransform::GetNTimeBins());
 325         break;
 326       }
 327
 328
 329       //Get the current ADC-value
 330       if(fDeconvTime){
 331
 332         //Check if the last pixel in the sequence is smaller than this
 333         if(charge > lastcharge){
 334           if(lastwas_falling){
 335             newbin = 1;
 336             break;
 337           }
 338         }
 339         else lastwas_falling = 1; //last pixel was larger than this
 340         lastcharge = charge;
 341       }
 342
 343       //Sum the total charge of this sequence
 344       seqcharge += charge;
 345       seqaverage += time*charge;
 346       seqerror += time*time*charge;
 347
 348       if (pCurrentPad) {
 349
 350         if((pCurrentPad->Next(kTRUE/*do zero suppression*/))==0) {
 351           pCurrentPad->StopEvent();
 352           pCurrentPad->ResetHistory();
 353           if(readValue) {
 354             newPad = fDigitReader->GetPad();
 355             newTime = fDigitReader->GetTime();
 356             newRow = fDigitReader->GetRow() + rowOffset;
 357           }
 358           break;
 359         }
 360
 361         newPad=pCurrentPad->GetPadNumber();
 362         newTime=pCurrentPad->GetCurrentPosition();
 363         newRow=pCurrentPad->GetRowNumber();
 364       } else {
 365       readValue = fDigitReader->Next();
 366       //Check where to stop:
 367       if(!readValue) break; //No more value
 368
 369       newPad = fDigitReader->GetPad();
 370       newTime = fDigitReader->GetTime();
 371       newRow = fDigitReader->GetRow() + rowOffset;
 372       }
 373
 374       if(newPad != pad)break; //new pad
 375       if(newTime != time+1) break; //end of sequence
 376
 377       // pad = newpad;    is equal
 378       time = newTime;
 379
 380     }//end loop over sequence
 381
 382     //HLTDebug("ended time bin sequence loop: seqcharge=%d readValue=%d", seqcharge, readValue);
 383     //HLTDebug("pad=%d newpad=%d current row=%d newrow=%d", pad, newPad, fCurrentRow, newRow);
 384     if (seqcharge<=0) {
 385       // with active zero suppression zero values are possible
 386       continue;
 387     }
 388
 389     //Calculate mean of sequence:
 390     Int_t seqmean=0;
 391     if(seqcharge)
 392       seqmean = seqaverage/seqcharge;
 393     else{
 394       LOG(AliHLTTPCLog::kFatal,"AliHLTTPCClusterFinder::ProcessRow","Data")
 395         <<"Error in data given to the cluster finder"<<ENDLOG;
 396       seqmean = 1;
 397       seqcharge = 1;
 398     }
 399
 400     //Calculate mean in pad direction:
 401     Int_t padmean = seqcharge*pad;
 402     Int_t paderror = pad*padmean;
 403
 404
 405     //Compare with results on previous pad:
 406     for(UInt_t p=0; p<nprevious; p++){
 407
 408       //dont merge sequences on the same pad twice
 409       if(previousPt[p]->fLastMergedPad==pad) continue;
 410
 411       Int_t difference = seqmean - previousPt[p]->fMean;
 412       if(difference < -fMatch) break;
 413
 414       if(difference <= fMatch){ //There is a match here!!
 415         AliClusterData *local = previousPt[p];
 416
 417         if(fDeconvPad){
 418           if(seqcharge > local->fLastCharge){
 419             if(local->fChargeFalling){ //The previous pad was falling
 420               break; //create a new cluster
 421             }
 422           }
 423           else local->fChargeFalling = 1;
 424           local->fLastCharge = seqcharge;
 425         }
 426
 427         //Don't create a new cluster, because we found a match
 428         newcluster = kFALSE;
 429
 430         //Update cluster on current pad with the matching one:
 431         local->fTotalCharge += seqcharge;
 432         local->fPad += padmean;
 433         local->fPad2 += paderror;
 434         local->fTime += seqaverage;
 435         local->fTime2 += seqerror;
 436         local->fMean = seqmean;
 437         local->fFlags++; //means we have more than one pad
 438         local->fLastMergedPad = pad;
 439
 440         currentPt[ncurrent] = local;
 441         ncurrent++;
 442
 443         break;
 444       } //Checking for match at previous pad
 445     } //Loop over results on previous pad.
 446
 447
 448     if(newcluster){
 449       //Start a new cluster. Add it to the clusterlist, and update
 450       //the list of pointers to clusters in current pad.
 451       //current pad will be previous pad on next pad.
 452
 453       //Add to the clusterlist:
 454       AliClusterData *tmp = &clusterlist[ntotal];
 455       tmp->fTotalCharge = seqcharge;
 456       tmp->fPad = padmean;
 457       tmp->fPad2 = paderror;
 458       tmp->fTime = seqaverage;
 459       tmp->fTime2 = seqerror;
 460       tmp->fMean = seqmean;
 461       tmp->fFlags = 0;  //flags for single pad clusters
 462       tmp->fLastMergedPad = pad;
 463
 464       if(fDeconvPad){
 465         tmp->fChargeFalling = 0;
 466         tmp->fLastCharge = seqcharge;
 467       }
 468
 469       //Update list of pointers to previous pad:
 470       currentPt[ncurrent] = &clusterlist[ntotal];
 471       ntotal++;
 472       ncurrent++;
 473     }
 474
 475     if(fDeconvTime)
 476       if(newbin >= 0) goto redo;
 477
 478     // to prevent endless loop
 479     if(time >= AliHLTTPCTransform::GetNTimeBins()){
 480       HLTWarning("Timebin (%d) out of range (%d)", time, AliHLTTPCTransform::GetNTimeBins());
 481       break;
 482     }
 483
 484
 485     if(!readValue) break; //No more value
 486
 487     if(fCurrentRow != newRow){
 488       WriteClusters(ntotal,clusterlist);
 489
 490       lastpad = 123456789;
 491
 492       currentPt = pad2;
 493       previousPt = pad1;
 494       nprevious=0;
 495       ncurrent=0;
 496       ntotal=0;
 497
 498       fCurrentRow = newRow;
 499     }
 500
 501     pad = newPad;
 502     time = newTime;
 503
 504   } // END while(readValue)
 505
 506   if (pCurrentPad) pCurrentPad->StopEvent();
 507
 508   WriteClusters(ntotal,clusterlist);
 509
 510   HLTInfo("ClusterFinder found %d clusters in slice %d patch %d", fNClusters, fCurrentSlice, fCurrentPatch);
 511
 512 } // ENDEND
 513
 514 void AliHLTTPCClusterFinder::WriteClusters(Int_t nclusters,AliClusterData *list)
 515 {
 516   //write cluster to output pointer
 517   Int_t thisrow,thissector;
 518   UInt_t counter = fNClusters;
 519
 520   for(int j=0; j<nclusters; j++)
 521     {
 522       if(!list[j].fFlags) continue; //discard single pad clusters
 523       if(list[j].fTotalCharge < fThreshold) continue; //noise cluster
 524
 525       Float_t xyz[3];
 526       Float_t fpad =(Float_t)list[j].fPad / list[j].fTotalCharge;
 527       Float_t fpad2=fXYErr*fXYErr; //fixed given error
 528       Float_t ftime =(Float_t)list[j].fTime / list[j].fTotalCharge;
 529       Float_t ftime2=fZErr*fZErr;  //fixed given error
 530
 531
 532
 533
 534
 535
 536       if(fCalcerr) { //calc the errors, otherwice take the fixed error
 537         Int_t patch = AliHLTTPCTransform::GetPatch(fCurrentRow);
 538         UInt_t q2=list[j].fTotalCharge*list[j].fTotalCharge;
 539         Float_t sy2=list[j].fPad2 * list[j].fTotalCharge - list[j].fPad * list[j].fPad;
 540         sy2/=q2;
 541         if(sy2 < 0) {
 542             LOG(AliHLTTPCLog::kError,"AliHLTTPCClusterFinder::WriteClusters","Cluster width")
 543               <<"SigmaY2 negative "<<sy2<<" on row "<<fCurrentRow<<" "<<fpad<<" "<<ftime<<ENDLOG;
 544             continue;
 545         } else {
 546           if(!fRawSP){
 547             fpad2 = (sy2 + 1./12)*AliHLTTPCTransform::GetPadPitchWidth(patch)*AliHLTTPCTransform::GetPadPitchWidth(patch);
 548             if(sy2 != 0){
 549               fpad2*=0.108; //constants are from offline studies
 550               if(patch<2)
 551                 fpad2*=2.07;
 552             }
 553           } else fpad2=sy2; //take the width not the error
 554         }
 555         Float_t sz2=list[j].fTime2*list[j].fTotalCharge - list[j].fTime*list[j].fTime;
 556         sz2/=q2;
 557         if(sz2 < 0){
 558           LOG(AliHLTTPCLog::kError,"AliHLTTPCClusterFinder::WriteClusters","Cluster width")
 559             <<"SigmaZ2 negative "<<sz2<<" on row "<<fCurrentRow<<" "<<fpad<<" "<<ftime<<ENDLOG;
 560           continue;
 561         } else {
 562           if(!fRawSP){
 563             ftime2 = (sz2 + 1./12)*AliHLTTPCTransform::GetZWidth()*AliHLTTPCTransform::GetZWidth();
 564             if(sz2 != 0) {
 565               ftime2 *= 0.169; //constants are from offline studies
 566               if(patch<2)
 567                 ftime2 *= 1.77;
 568             }
 569           } else ftime2=sz2; //take the width, not the error
 570         }
 571       }
 572       if(fStdout==kTRUE)
 573         cout<<"WriteCluster: padrow "<<fCurrentRow<<" pad "<<fpad << " +- "<<fpad2<<" time "<<ftime<<" +- "<<ftime2<<" charge "<<list[j].fTotalCharge<<endl;
 574
 575       if(!fRawSP){
 576         AliHLTTPCTransform::Slice2Sector(fCurrentSlice,fCurrentRow,thissector,thisrow);
 577         AliHLTTPCTransform::Raw2Local(xyz,thissector,thisrow,fpad,ftime);
 578
 579         if(xyz[0]==0) LOG(AliHLTTPCLog::kError,"AliHLTTPCClustFinder","Cluster Finder")
 580           <<AliHLTTPCLog::kDec<<"Zero cluster"<<ENDLOG;
 581         if(fNClusters >= fMaxNClusters)
 582           {
 583             LOG(AliHLTTPCLog::kError,"AliHLTTPCClustFinder::WriteClusters","Cluster Finder")
 584               <<AliHLTTPCLog::kDec<<"Too many clusters "<<fNClusters<<ENDLOG;
 585             return;
 586           }
 587
 588         fSpacePointData[counter].fX = xyz[0];
 589         fSpacePointData[counter].fY = xyz[1];
 590         fSpacePointData[counter].fZ = xyz[2];
 591
 592       } else {
 593         fSpacePointData[counter].fX = fCurrentRow;
 594         fSpacePointData[counter].fY = fpad;
 595         fSpacePointData[counter].fZ = ftime;
 596       }
 597
 598       fSpacePointData[counter].fCharge = list[j].fTotalCharge;
 599       fSpacePointData[counter].fPadRow = fCurrentRow;
 600       fSpacePointData[counter].fSigmaY2 = fpad2;
 601       fSpacePointData[counter].fSigmaZ2  = ftime2;
 602
 603       fSpacePointData[counter].fUsed = kFALSE;         // only used / set in AliHLTTPCDisplay
 604       fSpacePointData[counter].fTrackN = -1;           // only used / set in AliHLTTPCDisplay
 605
 606       Int_t patch=fCurrentPatch;
 607       if(patch==-1) patch=0; //never store negative patch number
 608       fSpacePointData[counter].fID = counter
 609         +((fCurrentSlice&0x7f)<<25)+((patch&0x7)<<22);//Uli
 610
 611 #ifdef do_mc
 612       Int_t trackID[3];
 613       GetTrackID((Int_t)rint(fpad),(Int_t)rint(ftime),trackID);
 614
 615       fSpacePointData[counter].fTrackID[0] = trackID[0];
 616       fSpacePointData[counter].fTrackID[1] = trackID[1];
 617       fSpacePointData[counter].fTrackID[2] = trackID[2];
 618
 619       //cout<<"padrow "<<fCurrentRow<<" pad "<<(Int_t)rint(fpad)<<" time "<<(Int_t)rint(ftime)<<" Trackid "<<trackID[0]<<endl;
 620 #endif
 621
 622       fNClusters++;
 623       counter++;
 624     }
 625 }
 626
 627 // STILL TO FIX  ----------------------------------------------------------------------------
 628
 629 #ifdef do_mc
 630 void AliHLTTPCClusterFinder::GetTrackID(Int_t pad,Int_t time,Int_t *trackID)
 631 {
 632   //get mc id
 633   AliHLTTPCDigitRowData *rowPt = (AliHLTTPCDigitRowData*)fDigitRowData;
 634
 635   trackID[0]=trackID[1]=trackID[2]=-2;
 636   //cout<<"Looking for pad "<<pad<<" time "<<time<<endl;
 637   for(Int_t i=fFirstRow; i<=fLastRow; i++){
 638     if(rowPt->fRow < (UInt_t)fCurrentRow){
 639       AliHLTTPCMemHandler::UpdateRowPointer(rowPt);
 640       continue;
 641     }
 642     AliHLTTPCDigitData *digPt = (AliHLTTPCDigitData*)rowPt->fDigitData;
 643     for(UInt_t j=0; j<rowPt->fNDigit; j++){
 644       Int_t cpad = digPt[j].fPad;
 645       Int_t ctime = digPt[j].fTime;
 646       if(cpad != pad) continue;
 647       if(ctime != time) continue;
 648
 649       trackID[0] = digPt[j].fTrackID[0];
 650       trackID[1] = digPt[j].fTrackID[1];
 651       trackID[2] = digPt[j].fTrackID[2];
 652
 653       //cout<<"Reading row "<<fCurrentRow<<" pad "<<cpad<<" time "<<ctime<<" trackID "<<digPt[j].fTrackID[0]<<endl;
 654       break;
 655     }
 656     break;
 657   }
 658 }
 659 #endif