HLT/TPCLib/AliHLTTPCConfMapper.cxx

   1 // @(#) $Id$
   2 // Original: AliHLTConfMapper.cxx,v 1.26 2005/06/14 10:55:21 cvetan Exp $
   3
   4 //**************************************************************************
   5 //* This file is property of and copyright by the ALICE HLT Project        *
   6 //* ALICE Experiment at CERN, All rights reserved.                         *
   7 //*                                                                        *
   8 //* Primary Authors: Anders Vestbo, maintained by
   9 //*                  Matthias Richter <Matthias.Richter@ift.uib.no>        *
  10 //*                  for The ALICE HLT Project.                            *
  11 //*                                                                        *
  12 //* Permission to use, copy, modify and distribute this software and its   *
  13 //* documentation strictly for non-commercial purposes is hereby granted   *
  14 //* without fee, provided that the above copyright notice appears in all   *
  15 //* copies and that both the copyright notice and this permission notice   *
  16 //* appear in the supporting documentation. The authors make no claims     *
  17 //* about the suitability of this software for any purpose. It is          *
  18 //* provided "as is" without express or implied warranty.                  *
  19 //**************************************************************************
  20
  21 /** @file   AliHLTTPCConfMapper.cxx
  22     @author Anders Vestbo, Matthias Richter
  23     @date   Conformal mapping base class.
  24     @brief
  25 */
  26
  27 #include <cassert>
  28 #include <sys/time.h>
  29
  30 #include "AliHLTTPCRootTypes.h"
  31 #include "AliHLTTPCSpacePointData.h"
  32 #include "AliHLTTPCLogging.h"
  33 #include "AliHLTTPCVertex.h"
  34 #include "AliHLTTPCConfMapTrack.h"
  35 #include "AliHLTTPCConfMapPoint.h"
  36 #include "AliHLTTPCTrackArray.h"
  37 #include "AliHLTTPCTransform.h"
  38 #include "AliHLTTPCConfMapper.h"
  39
  40 #if __GNUC__ >= 3
  41 using namespace std;
  42 #endif
  43
  44 ClassImp(AliHLTTPCConfMapper)
  45
  46 AliHLTTPCConfMapper::AliHLTTPCConfMapper()
  47   :
  48   fBench(kTRUE),
  49   fNTracks(0),
  50   fVertex(NULL),
  51   fVertexFinder(kFALSE),
  52   fHit(),
  53   fTrack(NULL),
  54   fMaxDca(0.0), // no clue whether this is reasonable, but at least better than without initialization
  55   fVolume(NULL),
  56   fRow(NULL),
  57   fNumRowSegment(0),
  58   fNumPhiSegment(0),
  59   fNumEtaSegment(0),
  60   fNumRowSegmentPlusOne(0),
  61   fNumPhiSegmentPlusOne(0),
  62   fNumEtaSegmentPlusOne(0),
  63   fNumPhiEtaSegmentPlusOne(0),
  64   fBounds(0),
  65   fPhiHitsOutOfRange(0),
  66   fEtaHitsOutOfRange(0),
  67   fPhiMin(0),
  68   fPhiMax(0),
  69   fEtaMin(0),
  70   fEtaMax(0),
  71   fRowMin(0),
  72   fRowMax(0),
  73   fVertexConstraint(kTRUE),
  74   fGoodDist(0.0),
  75   fMaxPhi(0.0),
  76   fMaxEta(0.0),
  77   fMainVertexTracks(0),
  78   fClustersUnused(0),
  79   fClusterCutZ(-1)
  80 {
  81   //Default constructor
  82   fParamSet[0]=0;
  83   fParamSet[1]=0;
  84 }
  85
  86 AliHLTTPCConfMapper::~AliHLTTPCConfMapper()
  87 {
  88   // Destructor.
  89   if(fRow) {
  90     delete [] fRow;
  91   }
  92   if(fTrack) {
  93     delete fTrack;
  94   }
  95 }
  96
  97 void AliHLTTPCConfMapper::InitVolumes()
  98 {
  99   //Data organization.
 100   //Allocate volumes, set conformal coordinates and pointers.
 101
 102   //Should be done after setting the track parameters
 103
 104   fNumRowSegmentPlusOne = AliHLTTPCTransform::GetNRows();//NumRows[0]; //Maximum 32.
 105   fNumPhiSegmentPlusOne = fNumPhiSegment+1;
 106   fNumEtaSegmentPlusOne = fNumEtaSegment+1;
 107   fNumPhiEtaSegmentPlusOne = fNumPhiSegmentPlusOne*fNumEtaSegmentPlusOne;
 108   fBounds = fNumRowSegmentPlusOne * fNumPhiSegmentPlusOne * fNumEtaSegmentPlusOne;
 109
 110   Reset();
 111
 112   fTrack = new AliHLTTPCTrackArray("AliHLTTPCConfMapTrack",10);
 113 }
 114
 115 void AliHLTTPCConfMapper::Reset()
 116 {
 117   if(fVolume) delete [] fVolume;
 118   fVolume=NULL;
 119   if(fRow) delete [] fRow;
 120   fRow=NULL;
 121
 122   fClustersUnused=0;
 123   fHit.clear();
 124 }
 125
 126 void AliHLTTPCConfMapper::InitSector(Int_t sector,Int_t *rowrange,Float_t *etarange)
 127 { //sector means slice here
 128   //Initialize tracker for tracking in a given sector.
 129   //Resets track and hit arrays.
 130   //Here it is also possible to specify a subsector, by defining
 131   //rowrange[0]=innermost row;
 132   //rowrange[1]=outermostrow;
 133   //Finally you can specify etaslices to save time (assuming a good seed from TRD...)
 134
 135   //Define tracking area:
 136   if(rowrange)
 137     {
 138       fRowMin = rowrange[0];
 139       fRowMax = rowrange[1];
 140     }
 141   else //complete sector
 142     {
 143       fRowMin = 0;
 144       fRowMax = AliHLTTPCTransform::GetNRows() - 1;
 145     }
 146   if(etarange)
 147     {
 148       fEtaMin = etarange[0];
 149       fEtaMax = sector < 18 ? etarange[1] : -etarange[1];
 150     }
 151   else
 152     {
 153       fEtaMin = 0;
 154       fEtaMax = sector < 18 ? 0.9 : -0.9;
 155     }
 156
 157   //Set the angles to sector 2:
 158   fPhiMin = -10*AliHLTTPCTransform::ToRad();//fParam->GetAngle(sector) - 10/todeg;
 159   fPhiMax =  10*AliHLTTPCTransform::ToRad();//fParam->GetAngle(sector) + 10/todeg;
 160
 161   fNTracks=0;
 162   fMainVertexTracks = 0;
 163   fClustersUnused = 0;
 164   fEtaHitsOutOfRange=0;
 165   fPhiHitsOutOfRange=0;
 166
 167   fNumRowSegment = fRowMax - fRowMin; //number of rows to be considered by tracker
 168   LOG(AliHLTTPCLog::kInformational,"AliHLTTPCConfMapper::InitSector","B-field")
 169     <<"Tracker initializing with a magnetic field of "<<AliHLTTPCTransform::GetBField()<<ENDLOG;
 170
 171   fTrack->Reset();
 172 }
 173
 174 Bool_t AliHLTTPCConfMapper::ReadHits(UInt_t count, AliHLTTPCSpacePointData* hits )
 175 {
 176   //read hits
 177   if(fClusterCutZ == -1){
 178     if (fHit.size()<fClustersUnused+count) fHit.resize(fClustersUnused+count);
 179     assert(fHit.size()>=fClustersUnused+count);
 180     for (Int_t i=0;(UInt_t)i<count;i++)
 181       {
 182         fHit[i+fClustersUnused].Reset();
 183         fHit[i+fClustersUnused].Read(hits[i]);
 184       }
 185     fClustersUnused += count;
 186   }
 187   else{
 188     //Skipping clusters with high Z.
 189     UInt_t skipped=0;
 190
 191     if (fHit.size()<fClustersUnused+count) fHit.resize(fClustersUnused+count);
 192     assert(fHit.size()>=fClustersUnused+count);
 193     for (Int_t i=0;(UInt_t)i<count;i++)
 194       {
 195         if(hits[i].fZ > fClusterCutZ || hits[i].fZ < -1*fClusterCutZ){
 196           ++skipped;
 197           continue;
 198         }
 199         fHit[i+fClustersUnused-skipped].Reset();
 200         fHit[i+fClustersUnused-skipped].Read(hits[i]);
 201       }
 202     fClustersUnused += count - skipped;
 203     fHit.resize(fClustersUnused);
 204   }
 205
 206   LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapper::ReadHits","#hits")
 207     <<AliHLTTPCLog::kDec<<"#hits: "<<count<<" total: "<<fClustersUnused<<ENDLOG;
 208
 209   return true;
 210 }
 211
 212 void AliHLTTPCConfMapper::SetPointers()
 213 {
 214   //Check if there are not enough clusters to make a track in this sector
 215   //Can happen in pp events.
 216
 217   if(fClustersUnused < fMinPoints[fVertexConstraint])
 218     return;
 219
 220   //Allocate detector volumes
 221   if (fVolume==NULL) {
 222     LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapper::InitVolumes","Memory")<<AliHLTTPCLog::kDec<<
 223       "Allocating "<<fBounds*sizeof(AliHLTTPCConfMapContainer)<<" Bytes to fVolume"<<ENDLOG;
 224     fVolume = new AliHLTTPCConfMapContainer[fBounds];
 225   }
 226
 227   if (fRow==NULL) {
 228     LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapper::InitVolumes","Memory")<<AliHLTTPCLog::kDec<<
 229       "Allocating "<<fNumRowSegmentPlusOne*sizeof(AliHLTTPCConfMapContainer)<<" Bytes to fRow"<<ENDLOG;
 230     fRow = new AliHLTTPCConfMapContainer[fNumRowSegmentPlusOne];
 231   }
 232
 233   memset(fVolume,0,fBounds*sizeof(AliHLTTPCConfMapContainer));
 234   memset(fRow,0,fNumRowSegmentPlusOne*sizeof(AliHLTTPCConfMapContainer));
 235
 236   Float_t phiSlice = (fPhiMax-fPhiMin)/fNumPhiSegment;
 237   Float_t etaSlice = (fEtaMax-fEtaMin)/fNumEtaSegment;
 238
 239   Int_t volumeIndex;
 240   Int_t localcounter=0;
 241   assert((int)fHit.size()>=fClustersUnused);
 242   for(Int_t j=0; j<fClustersUnused; j++)
 243     {
 244       AliHLTTPCConfMapPoint *thisHit = &(fHit[j]);
 245
 246       thisHit->Setup(fVertex);
 247
 248       Int_t localrow = thisHit->GetPadRow();
 249
 250       if(localrow < fRowMin || localrow > fRowMax)
 251         continue;
 252
 253       //Get indexes:
 254       thisHit->SetPhiIndex((Int_t)((thisHit->GetPhi()-fPhiMin)/phiSlice +1));
 255
 256       if(thisHit->GetPhiIndex()<1 || thisHit->GetPhiIndex()>fNumPhiSegment)
 257         {
 258           //cout << "Phiindex: " << thisHit->phiIndex << " " << thisHit->GetPhi() << endl;
 259           fPhiHitsOutOfRange++;
 260           continue;
 261         }
 262
 263       thisHit->SetEtaIndex((Int_t)((thisHit->GetEta()-fEtaMin)/etaSlice + 1));
 264       if(thisHit->GetEtaIndex()<1 || thisHit->GetEtaIndex()>fNumEtaSegment)
 265         {
 266           //cout << "Etaindex: " << thisHit->etaIndex << " " << thisHit->GetEta() << endl;
 267           fEtaHitsOutOfRange++;
 268           continue;
 269         }
 270       localcounter++;
 271
 272       volumeIndex = (localrow-fRowMin)*fNumPhiEtaSegmentPlusOne +
 273                     thisHit->GetPhiIndex()*fNumEtaSegmentPlusOne+thisHit->GetEtaIndex();
 274
 275       if(fVolume[volumeIndex].first == NULL)
 276         fVolume[volumeIndex].first = (void *)thisHit;
 277       else
 278         ((AliHLTTPCConfMapPoint *)fVolume[volumeIndex].last)->SetNextVolumeHit(thisHit);
 279       fVolume[volumeIndex].last = (void *)thisHit;
 280
 281
 282       //set row pointers
 283       if(fRow[(localrow-fRowMin)].first == NULL)
 284         fRow[(localrow-fRowMin)].first = (void *)thisHit;
 285       else
 286         ((AliHLTTPCConfMapPoint *)(fRow[(localrow-fRowMin)].last))->SetNextRowHit(thisHit);
 287         fRow[(localrow-fRowMin)].last = (void *)thisHit;
 288     }
 289
 290   //If a cluster has an Eta outside the Eta or Phi range set in the Tracker, it will go in to
 291   //the if here. This has been seen for high Eta clusters most likely from signal from the gating grid.
 292   //These clusters are read in, but not used in the Tracking.
 293 #ifdef PACKAGE_STRING
 294   if(fClustersUnused>0 && localcounter==0)
 295     LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapper::SetPointers","Parameters")
 296       <<AliHLTTPCLog::kDec<<"No points passed to track finder, hits out of range: "
 297       <<fEtaHitsOutOfRange+fPhiHitsOutOfRange<<ENDLOG;
 298
 299   Int_t hits_accepted=fClustersUnused-(fEtaHitsOutOfRange+fPhiHitsOutOfRange);
 300   LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapper::SetPointers","Setup")
 301     <<"Setup finished, hits out of range: "<<fEtaHitsOutOfRange+fPhiHitsOutOfRange
 302     <<" hits accepted "<<hits_accepted<<ENDLOG;
 303 #endif //PACKAGE_STRING
 304 }
 305
 306 void AliHLTTPCConfMapper::MainVertexTrackingA()
 307 {
 308   //Tracking with vertex constraint.
 309
 310   if(!fParamSet[(Int_t)kTRUE])
 311     {
 312       LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapper::MainVertexTracking","Parameters")<<AliHLTTPCLog::kDec<<
 313         "Tracking parameters not set!"<<ENDLOG;
 314       return;
 315     }
 316
 317   Double_t initCpuTime,cpuTime;
 318   initCpuTime = CpuTime();
 319
 320   SetPointers();
 321   SetVertexConstraint(true);
 322   cpuTime = CpuTime() - initCpuTime;
 323   if(fBench)
 324     LOG(AliHLTTPCLog::kBenchmark,"AliHLTTPCConfMapper::MainVertexTrackingA","Timing")
 325       <<AliHLTTPCLog::kDec<<"Setup finished in "<<cpuTime*1000<<" ms"<<ENDLOG;
 326
 327 }
 328
 329 void AliHLTTPCConfMapper::MainVertexTrackingB()
 330 {
 331   //Tracking with vertex constraint.
 332
 333   if(!fParamSet[(Int_t)kTRUE])
 334     {
 335       LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapper::MainVertexTracking","Parameters")<<AliHLTTPCLog::kDec<<
 336         "Tracking parameters not set!"<<ENDLOG;
 337       return;
 338     }
 339   Double_t initCpuTime,cpuTime;
 340   initCpuTime = CpuTime();
 341
 342   ClusterLoop();
 343
 344   cpuTime = CpuTime() - initCpuTime;
 345   if(fBench)
 346     LOG(AliHLTTPCLog::kBenchmark,"AliHLTTPCConfMapper::MainVertexTrackingB","Timing")
 347       <<AliHLTTPCLog::kDec<<"Main Tracking finished in "<<cpuTime*1000<<" ms"<<ENDLOG;
 348 }
 349
 350 void AliHLTTPCConfMapper::MainVertexTracking()
 351 {
 352   //Tracking with vertex constraint.
 353
 354   if(!fParamSet[(Int_t)kTRUE])
 355     {
 356       LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapper::MainVertexTracking","Parameters")<<AliHLTTPCLog::kDec<<
 357         "Tracking parameters not set!"<<ENDLOG;
 358       return;
 359     }
 360
 361   Double_t initCpuTime,cpuTime;
 362   initCpuTime = CpuTime();
 363
 364   SetPointers();
 365
 366   SetVertexConstraint(true);
 367
 368   ClusterLoop();
 369
 370   cpuTime = CpuTime() - initCpuTime;
 371   if(fBench)
 372     LOG(AliHLTTPCLog::kBenchmark,"AliHLTTPCConfMapper::MainVertexTracking","Timing")<<AliHLTTPCLog::kDec<<
 373       "Tracking finished in "<<cpuTime*1000<<" ms"<<ENDLOG;
 374
 375   return;
 376 }
 377
 378 void AliHLTTPCConfMapper::NonVertexTracking()
 379 {
 380   //Tracking with no vertex constraint. This should be called after doing MainVertexTracking,
 381   //in order to do tracking on the remaining clusters.
 382   //The conformal mapping is now done with respect to the first cluster
 383   //assosciated with this track.
 384
 385   if(!fParamSet[(Int_t)kFALSE])
 386     {
 387       LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapper::NonVertexTracking","Parameters")<<AliHLTTPCLog::kDec<<
 388         "Tracking parameters not set!"<<ENDLOG;
 389       return;
 390     }
 391
 392   SetVertexConstraint(false);
 393
 394   SetPointers(); //To be able to do only nonvertextracking (more testing)
 395
 396   ClusterLoop();
 397   LOG(AliHLTTPCLog::kInformational,"AliHLTTPCConfMapper::NonVertexTracking","ntracks")<<AliHLTTPCLog::kDec<<
 398     "Number of nonvertex tracks found: "<<(fNTracks-fMainVertexTracks)<<ENDLOG;
 399   return;
 400 }
 401
 402 void AliHLTTPCConfMapper::MainVertexSettings(Int_t trackletlength, Int_t tracklength,
 403                                          Int_t rowscopetracklet, Int_t rowscopetrack,
 404                                          Double_t maxphi,Double_t maxeta)
 405 {
 406   //Settings for main vertex tracking. The cuts are:
 407   //TrackletLength:      #hits on segment, before trying to build a track
 408   //TrackLength:         Minimum hits on a track
 409   //RowScopeTracklet:    Row search range for segments
 410   //RowScopeTrack:       Row search range for tracks
 411
 412   SetTrackletLength(trackletlength,(Bool_t)true);
 413   SetRowScopeTracklet(rowscopetracklet, (Bool_t) true);
 414   SetRowScopeTrack(rowscopetrack, (Bool_t) true);
 415   SetMinPoints(tracklength,(Bool_t)true);
 416   fMaxPhi=maxphi;
 417   fMaxEta=maxeta;
 418   SetParamDone(kTRUE);
 419 }
 420
 421 void AliHLTTPCConfMapper::NonVertexSettings(Int_t trackletlength, Int_t tracklength,
 422                                         Int_t rowscopetracklet, Int_t rowscopetrack)
 423 {
 424   //set parameters for non-vertex tracking
 425   SetTrackletLength(trackletlength,(Bool_t)false);
 426   SetRowScopeTracklet(rowscopetracklet, (Bool_t)false);
 427   SetRowScopeTrack(rowscopetrack, (Bool_t)false);
 428   SetMinPoints(tracklength,(Bool_t)false);
 429   SetParamDone(kFALSE);
 430 }
 431
 432 void AliHLTTPCConfMapper::SetTrackCuts(Double_t hitChi2Cut, Double_t goodHitChi2, Double_t trackChi2Cut,Int_t maxdist,Bool_t vertexconstraint)
 433 {
 434   //Settings for tracks. The cuts are:
 435   //HitChi2Cut:     Maximum hit chi2
 436   //goodHitChi2:    Chi2 to stop look for next hit
 437   //trackChi2Cut:   Maximum track chi2
 438   //maxdist:        Maximum distance between two clusters when forming segments
 439
 440   SetHitChi2Cut(hitChi2Cut,vertexconstraint);
 441   SetGoodHitChi2(goodHitChi2,vertexconstraint);
 442   SetTrackChi2Cut(trackChi2Cut,vertexconstraint);
 443   SetMaxDist(maxdist,vertexconstraint);
 444 }
 445
 446 void AliHLTTPCConfMapper::SetTrackletCuts(Double_t maxangle,Double_t goodDist, Bool_t vc)
 447 {
 448   //Sets cuts of tracklets. Right now this is only:
 449   //maxangle:  Maximum angle when forming segments (if trackletlength > 2)
 450
 451   fGoodDist=goodDist;
 452   SetMaxAngleTracklet(maxangle, vc);
 453 }
 454
 455 void AliHLTTPCConfMapper::ClusterLoop()
 456 {
 457   //Loop over hits, starting at outermost padrow, and trying to build segments.
 458
 459   //Check if there are not enough clusters to make a track in this sector
 460   //Can happen in pp events.
 461   if(fClustersUnused < fMinPoints[fVertexConstraint])
 462     return;
 463
 464   Int_t rowsegm,lastrow = fRowMin + fMinPoints[fVertexConstraint];
 465   AliHLTTPCConfMapPoint *hit;
 466
 467   //Loop over rows, and try to create tracks from the hits.
 468   //Starts at the outermost row, and loops as long as a track can be build, due to length.
 469
 470   for(rowsegm = fRowMax; rowsegm >= lastrow; rowsegm--)
 471     {
 472       if(fRow[(rowsegm-fRowMin)].first && ((AliHLTTPCConfMapPoint*)fRow[(rowsegm-fRowMin)].first)->GetPadRow() < fRowMin + 1)
 473         break;
 474
 475       for(hit = (AliHLTTPCConfMapPoint*)fRow[(rowsegm-fRowMin)].first; hit!=0; hit=hit->GetNextRowHit())
 476         {
 477           if(hit->GetUsage() == true)
 478             continue;
 479           else
 480             CreateTrack(hit);
 481         }
 482     }
 483
 484   return;
 485 }
 486
 487
 488 void AliHLTTPCConfMapper::CreateTrack(AliHLTTPCConfMapPoint *hit)
 489 {
 490   //Tries to create a track from the initial hit given by ClusterLoop()
 491
 492   AliHLTTPCConfMapPoint *closesthit = NULL;
 493   AliHLTTPCConfMapTrack *track = NULL;
 494
 495   Int_t point;
 496   Int_t tracks = fNTracks;
 497   fNTracks++;
 498
 499   track = (AliHLTTPCConfMapTrack*)fTrack->NextTrack();
 500
 501   //reset hit parameters:
 502   track->Reset();
 503
 504   UInt_t *trackhitnumber = track->GetHitNumbers();
 505
 506   //set conformal coordinates if we are looking for non vertex tracks
 507   if(!fVertexConstraint)
 508     {
 509       hit->SetAllCoord(hit);
 510     }
 511
 512   //fill fit parameters of initial track:
 513   track->UpdateParam(hit); //here the number of hits is incremented.
 514   trackhitnumber[track->GetNumberOfPoints()-1] = hit->GetHitNumber();
 515
 516   Double_t dx,dy;
 517   //create tracklets:
 518
 519   for(point=1; point<fTrackletLength[fVertexConstraint]; point++)
 520     {
 521       if((closesthit = GetNextNeighbor(hit)))
 522         {//closest hit exist
 523
 524           //   Calculate track length in sz plane
 525           dx = ((AliHLTTPCConfMapPoint*)closesthit)->GetX() - ((AliHLTTPCConfMapPoint*)hit)->GetX();
 526           dy = ((AliHLTTPCConfMapPoint*)closesthit)->GetY() - ((AliHLTTPCConfMapPoint*)hit)->GetY();
 527           //track->fLength += (Double_t)sqrt ( dx * dx + dy * dy ) ;
 528           Double_t length = track->GetLength()+(Double_t)sqrt ( dx * dx + dy * dy );
 529           track->SetLength(length);
 530
 531           //closesthit->SetS(track->fLength);
 532           closesthit->SetS(track->GetLength());
 533
 534           //update fit parameters
 535           track->UpdateParam(closesthit);
 536           trackhitnumber[track->GetNumberOfPoints()-1] = closesthit->GetHitNumber();
 537
 538           hit = closesthit;
 539         }
 540       else
 541         {
 542           //closest hit does not exist:
 543           track->DeleteCandidate();
 544           fTrack->RemoveLast();
 545           fNTracks--;
 546           point = fTrackletLength[fVertexConstraint];
 547         }
 548     }
 549
 550   //tracklet is long enough to be extended to a track
 551   if(track->GetNumberOfPoints() == fTrackletLength[fVertexConstraint])
 552     {
 553
 554       track->SetProperties(true);
 555
 556       if(TrackletAngle(track) > fMaxAngleTracklet[fVertexConstraint])
 557         {//proof if the first points seem to be a beginning of a track
 558           track->SetProperties(false);
 559           track->DeleteCandidate();
 560           fTrack->RemoveLast();
 561           fNTracks--;
 562         }
 563
 564       else//good tracklet ->proceed, follow the trackfit
 565         {
 566           tracks++;
 567
 568           //define variables to keep the total chi:
 569           Double_t xyChi2 = track->GetChiSq1();
 570           Double_t szChi2 = track->GetChiSq2();
 571
 572           for(point = fTrackletLength[fVertexConstraint]; point <= fNumRowSegment; point++)
 573             {
 574               track->SetChiSq1(fHitChi2Cut[fVertexConstraint]);
 575               closesthit = GetNextNeighbor((AliHLTTPCConfMapPoint*)track->GetLastHit(),track);
 576
 577               if(closesthit)
 578                 {
 579                   //keep total chi:
 580                   Double_t lxyChi2 = track->GetChiSq1()-track->GetChiSq2();
 581                   xyChi2 += lxyChi2;
 582                   closesthit->SetXYChi2(lxyChi2);
 583
 584                   //update track length:
 585                   track->SetLength(closesthit->GetS());
 586                   szChi2 += track->GetChiSq2();
 587                   closesthit->SetSZChi2(track->GetChiSq2());
 588
 589                   track->UpdateParam(closesthit);
 590                   trackhitnumber[track->GetNumberOfPoints()-1] = closesthit->GetHitNumber();
 591
 592                   //add closest hit to track
 593                   closesthit->SetUsage(true);
 594                   closesthit->SetTrackNumber(tracks-1);
 595
 596                 }//closesthit
 597
 598               else
 599                 {
 600                   //closest hit does not exist
 601                   point = fNumRowSegment; //continue with next hit in segment
 602                 }//else
 603
 604             }//create tracks
 605
 606           //store track chi2:
 607           track->SetChiSq1(xyChi2);
 608           track->SetChiSq2(szChi2);
 609           Double_t normalizedchi2 = (track->GetChiSq1()+track->GetChiSq2())/track->GetNumberOfPoints();
 610
 611           //remove tracks with not enough points already now
 612           if(track->GetNumberOfPoints() < fMinPoints[fVertexConstraint] || normalizedchi2 > fTrackChi2Cut[fVertexConstraint])
 613             {
 614               track->SetProperties(false);
 615               fNTracks--;
 616               track->DeleteCandidate();
 617               fTrack->RemoveLast();
 618               tracks--;
 619             }
 620
 621           else
 622             {
 623               fClustersUnused -= track->GetNumberOfPoints();
 624               track->ComesFromMainVertex(fVertexConstraint);
 625               //mark track as main vertex track or not
 626               track->SetSector(2); //only needed for testing purposes.
 627               track->SetRowRange(fRowMin,fRowMax);
 628
 629               if(fVertexConstraint)
 630                 fMainVertexTracks++;
 631             }
 632
 633         }//good tracklet
 634
 635     }
 636
 637   return;
 638 }
 639
 640 AliHLTTPCConfMapPoint *AliHLTTPCConfMapper::GetNextNeighbor(AliHLTTPCConfMapPoint *starthit,
 641                                           AliHLTTPCConfMapTrack *track)
 642 {
 643   //When forming segments: Finds closest hit to input hit
 644   //When forming tracks: Find closest hit to track fit.
 645
 646   Double_t dist,closestdist = fMaxDist[fVertexConstraint];
 647
 648   AliHLTTPCConfMapPoint *hit = NULL;
 649   AliHLTTPCConfMapPoint *closesthit = NULL;
 650
 651   Int_t subrowsegm;
 652   Int_t subphisegm;
 653   Int_t subetasegm;
 654   Int_t volumeIndex;
 655   Int_t testhit;
 656
 657   Int_t maxrow = starthit->GetPadRow()-1;
 658   Int_t minrow;
 659
 660   if(track) //finding hit close to trackfit
 661     {
 662       minrow = starthit->GetPadRow()-fRowScopeTrack[fVertexConstraint];
 663     }
 664   else
 665     {
 666       minrow = starthit->GetPadRow()-fRowScopeTracklet[fVertexConstraint];
 667     }
 668
 669   //make a smart loop
 670   Int_t loopeta[25] = {0,0,0,-1,-1,-1,1,1,1, 0,0,-1,-1,1,1,-2,-2,-2,-2,-2,2,2,2,2,2};
 671   Int_t loopphi[25] = {0,-1,1,0,-1,1,0,-1,1, -2,2,-2,2,-2,2,-2,-1,0,1,2,-2,-1,0,1,2};
 672
 673   if(minrow < fRowMin)
 674     minrow = fRowMin;
 675   if(maxrow < fRowMin)
 676     return 0;  //reached the last padrow under consideration
 677
 678   else
 679     {
 680       //loop over sub rows
 681       for(subrowsegm=maxrow; subrowsegm>=minrow; subrowsegm--)
 682         {
 683           //loop over subsegments, in the order defined above.
 684           for(Int_t i=0; i<9; i++)
 685             {
 686               subphisegm = starthit->GetPhiIndex() + loopphi[i];
 687
 688               if(subphisegm < 0 || subphisegm >= fNumPhiSegment)
 689                 continue;
 690               /*
 691                 if(subphisegm<0)
 692                 subphisegm += fNumPhiSegment;
 693
 694                 else if(subphisegm >=fNumPhiSegment)
 695                 subphisegm -= fNumPhiSegment;
 696               */
 697               //loop over sub eta segments
 698
 699               subetasegm = starthit->GetEtaIndex() + loopeta[i];
 700
 701               if(subetasegm < 0 || subetasegm >=fNumEtaSegment)
 702                 continue;//segment exceeds bounds->skip it
 703
 704               //loop over hits in this sub segment:
 705               volumeIndex=(subrowsegm-fRowMin)*fNumPhiEtaSegmentPlusOne +
 706                 subphisegm*fNumEtaSegmentPlusOne + subetasegm;
 707
 708               if(volumeIndex<0)
 709                 {//debugging
 710                   LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapper::GetNextNeighbor","Memory")<<AliHLTTPCLog::kDec<<
 711                     "VolumeIndex error "<<volumeIndex<<ENDLOG;
 712                 }
 713
 714               assert(fVolume!=NULL);
 715               for(hit = (AliHLTTPCConfMapPoint*)fVolume[volumeIndex].first;
 716                   hit!=0; hit = hit->GetNextVolumeHit())
 717                 {
 718
 719                   if(!hit->GetUsage())
 720                     {//hit was not used before
 721
 722                       //set conformal mapping if looking for nonvertex tracks:
 723                       if(!fVertexConstraint)
 724                         {
 725                           hit->SetAllCoord(starthit);
 726                         }
 727
 728                       if(track)//track search - look for nearest neighbor to extrapolated track
 729                         {
 730                           if (fVertexConstraint) {
 731                             if(!VerifyRange(starthit,hit))
 732                               continue;
 733                           }
 734                           testhit = EvaluateHit(starthit,hit,track);
 735
 736                           if(testhit == 0)//chi2 not good enough, keep looking
 737                             continue;
 738                           else if(testhit==2)//chi2 good enough, return it
 739                             return hit;
 740                           else
 741                             closesthit = hit;//chi2 acceptable, but keep looking
 742
 743                         }//track search
 744
 745                       else //tracklet search, look for nearest neighbor
 746                         {
 747
 748                           if((dist=CalcDistance(starthit,hit)) < closestdist)
 749                             {
 750                               if (fVertexConstraint) {
 751                                 if(!VerifyRange(starthit,hit))
 752                                   continue;
 753                               }
 754                               closestdist = dist;
 755                               closesthit = hit;
 756
 757                               //if this hit is good enough, return it:
 758                               if(closestdist < fGoodDist)
 759                                 return closesthit;
 760                             }
 761                           else
 762                             continue;//sub hit was farther away than a hit before
 763
 764                         }//tracklet search
 765
 766                     }//hit not used before
 767
 768                   else continue; //sub hit was used before
 769
 770                 }//loop over hits in sub segment
 771
 772             }//loop over sub segments
 773
 774         }//loop over subrows
 775
 776     }//else
 777
 778   //closest hit found:
 779   if(closesthit)// && closestdist < mMaxDist)
 780     return closesthit;
 781   else
 782     return 0;
 783 }
 784
 785 Int_t AliHLTTPCConfMapper::EvaluateHit(AliHLTTPCConfMapPoint *starthit,AliHLTTPCConfMapPoint *hit,AliHLTTPCConfMapTrack *track)
 786 {
 787   //Check if space point gives a fit with acceptable chi2.
 788
 789   Double_t temp,dxy,lchi2,dx,dy,slocal,dsz,lszChi2;
 790   temp = (track->GetA2Xy()*hit->GetXprime()-hit->GetYprime()+track->GetA1Xy());
 791   dxy = temp*temp/(track->GetA2Xy()*track->GetA2Xy() + 1.);
 792
 793   //Calculate chi2
 794   lchi2 = (dxy*hit->GetXYWeight());
 795
 796   if(lchi2 > track->GetChiSq1())//chi2 was worse than before.
 797     return 0;
 798
 799   //calculate s and the distance hit-line
 800   dx = starthit->GetX()-hit->GetX();
 801   dy = starthit->GetY()-hit->GetY();
 802   //slocal = track->fLength+sqrt(dx*dx+dy*dy);
 803   slocal = track->GetLength()+sqrt(dx*dx+dy*dy);
 804
 805   temp = (track->GetA2Sz()*slocal-hit->GetZ()+track->GetA1Sz());
 806   dsz = temp*temp/(track->GetA2Sz()*track->GetA2Sz()+1);
 807
 808   //calculate chi2
 809   lszChi2 = dsz*hit->GetZWeight();
 810   lchi2 += lszChi2;
 811
 812
 813   //check whether chi2 is better than previous one:
 814   if(lchi2 < track->GetChiSq1())
 815     {
 816       track->SetChiSq1(lchi2);
 817       track->SetChiSq2(lszChi2);
 818
 819       hit->SetS(slocal);
 820
 821       //if chi2 good enough, stop here:
 822       if(lchi2 < fGoodHitChi2[fVertexConstraint])
 823         return 2;
 824
 825       return 1;
 826     }
 827
 828   return 0;
 829
 830 }
 831
 832 Double_t AliHLTTPCConfMapper::CalcDistance(const AliHLTTPCConfMapPoint *hit1,const AliHLTTPCConfMapPoint *hit2) const
 833 {
 834   //Return distance between two clusters, defined by Pablo
 835
 836   Double_t phidiff = fabs( hit1->GetPhi() - hit2->GetPhi() );
 837   if (phidiff > AliHLTTPCTransform::Pi()) phidiff = AliHLTTPCTransform::TwoPi() - phidiff;
 838
 839   return AliHLTTPCTransform::ToDeg()*fabs((Float_t)((hit1->GetPadRow() - hit2->GetPadRow()) *
 840          (phidiff + fabs( hit1->GetEta() - hit2->GetEta()))));
 841 }
 842
 843 Bool_t AliHLTTPCConfMapper::VerifyRange(const AliHLTTPCConfMapPoint *hit1,const AliHLTTPCConfMapPoint *hit2) const
 844 {
 845   //Check if the hit are within reasonable range in phi and eta
 846   Double_t dphi,deta;//maxphi=0.1,maxeta=0.1;
 847   dphi = fabs(hit1->GetPhi() - hit2->GetPhi());
 848   if(dphi > AliHLTTPCTransform::Pi()) dphi = fabs(AliHLTTPCTransform::TwoPi() - dphi);
 849   if(dphi > fMaxPhi) return false;
 850
 851   deta = fabs(hit1->GetEta() - hit2->GetEta());
 852   if(deta > fMaxEta) return false;
 853
 854   return true;
 855
 856 }
 857
 858 Double_t AliHLTTPCConfMapper::TrackletAngle(AliHLTTPCConfMapTrack *track,Int_t n) const
 859 {
 860   // Returns the angle 'between' the last three points (started at point number n) on this track.
 861
 862   if(n > track->GetNumberOfPoints())
 863     n = track->GetNumberOfPoints();
 864
 865   if(n<3)
 866     return 0;
 867
 868   Double_t x1[2]={0,0};
 869   Double_t x2[2]={0,0};
 870   Double_t x3[2]={0,0};
 871   Double_t angle1,angle2;
 872   Int_t counter=0;
 873   for(track->StartLoop(); track->LoopDone(); track->GetNextHit())
 874     {
 875       AliHLTTPCConfMapPoint *p = (AliHLTTPCConfMapPoint*)track->GetCurrentHit();
 876       if( (n-1) == counter)
 877         {
 878           x1[0] = p->GetX();
 879           x1[1] = p->GetY();
 880         }
 881       else if( (n-2) == counter)
 882         {
 883           x2[0] = p->GetX();
 884           x2[1] = p->GetY();
 885         }
 886       else if( (n-3) == counter)
 887         {
 888           x3[0] = p->GetX();
 889           x3[1] = p->GetY();
 890         }
 891       counter++;
 892     }
 893
 894   angle1 = atan2(x2[1]-x3[1],x2[0]-x3[0]);
 895   angle2 = atan2(x1[1]-x2[1],x1[0]-x2[0]);
 896
 897   return fabs(angle1-angle2);
 898
 899   /*
 900     Double_t x1[2];
 901   Double_t x2[2];
 902   Double_t angle1,angle2;
 903   TObjArray *hits = track->GetHits();
 904
 905   if (n > track->GetNumberOfPoints()) {
 906     n = track->GetNumberOfPoints();
 907   }
 908
 909   if (n<3)
 910     return 0;
 911
 912
 913   x1[0] = ((AliHLTTPCConfMapPoint *)hits->At(n-2))->GetX() - ((AliHLTTPCConfMapPoint *)hits->At(n-3))->GetX();
 914   x1[1] = ((AliHLTTPCConfMapPoint *)hits->At(n-2))->GetY() - ((AliHLTTPCConfMapPoint *)hits->At(n-3))->GetY();
 915
 916   x2[0] = ((AliHLTTPCConfMapPoint *)hits->At(n-1))->GetX() - ((AliHLTTPCConfMapPoint *)hits->At(n-2))->GetX();
 917   x2[1] = ((AliHLTTPCConfMapPoint *)hits->At(n-1))->GetY() - ((AliHLTTPCConfMapPoint *)hits->At(n-2))->GetY();
 918
 919   angle1 = atan2(x1[1],x1[0]);
 920   angle2 = atan2(x2[1],x1[0]);
 921   return fabs(angle1-angle2);
 922   */
 923 }
 924
 925 Int_t AliHLTTPCConfMapper::FillTracks()
 926 {
 927   //Fill track parameters. Which basically means do a fit of helix in real space,
 928   //which should be done in order to get nice tracks.
 929
 930   Int_t numoftracks = fNTracks;
 931   if(fNTracks == 0)
 932     {
 933       LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapper::FillTracks","fNTracks")<<AliHLTTPCLog::kDec<<
 934         "No tracks found!!"<<ENDLOG;
 935       return 0;
 936     }
 937
 938   LOG(AliHLTTPCLog::kInformational,"AliHLTTPCConfMapper::FillTracks","fNTracks")<<AliHLTTPCLog::kDec<<
 939     "Number of found tracks: "<<fNTracks<<ENDLOG;
 940
 941   //  fTrack->Sort();
 942   for(Int_t i=0; i<numoftracks; i++)
 943     {
 944       AliHLTTPCConfMapTrack *track = (AliHLTTPCConfMapTrack*)fTrack->GetTrack(i);
 945       track->Fill(fVertex,fMaxDca);
 946     }
 947   return 1;
 948 }
 949
 950 Double_t AliHLTTPCConfMapper::CpuTime()
 951 {
 952   //Return the Cputime in seconds.
 953  struct timeval tv;
 954  gettimeofday( &tv, NULL );
 955  return tv.tv_sec+(((Double_t)tv.tv_usec)/1000000.);
 956 }