HLT/TPCLib/tracking/AliHLTTPCHoughTrack.cxx

   1 // @(#) $Id$
   2 // origin: hough/AliL3HoughTrack.cxx,v 1.23 ue Mar 28 18:05:12 2006 UTC by alibrary
   3
   4 // Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
   5 //*-- Copyright &copy ALICE HLT Group
   6
   7 #include "AliHLTStandardIncludes.h"
   8
   9 #include "AliHLTLogging.h"
  10 #include "AliHLTTrack.h"
  11 #include "AliHLTTPCHoughTrack.h"
  12 #include "AliHLTTransform.h"
  13 #include "AliHLTHoughTransformerRow.h"
  14
  15 #if __GNUC__ >= 3
  16 using namespace std;
  17 #endif
  18
  19 /** \class AliHLTTPCHoughTrack
  20 <pre>
  21 //_____________________________________________________________
  22 // AliHLTTPCHoughTrack
  23 //
  24 // Track class for Hough tracklets
  25 //
  26 </pre>
  27 */
  28
  29 ClassImp(AliHLTTPCHoughTrack)
  30
  31
  32     AliHLTTPCHoughTrack::AliHLTTPCHoughTrack() : AliHLTTrack()
  33 {
  34   //Constructor
  35
  36   fWeight = 0;
  37   fMinDist=0;
  38   fDLine = 0;
  39   fPsiLine = 0;
  40   fIsHelix = true;
  41   fEtaIndex = -1;
  42   fEta = 0;
  43   ComesFromMainVertex(kTRUE);
  44 }
  45
  46 AliHLTTPCHoughTrack::~AliHLTTPCHoughTrack()
  47 {
  48   //dtor
  49 }
  50
  51 void AliHLTTPCHoughTrack::Set(AliHLTTrack *track)
  52 {
  53   //Basically copy constructor
  54   AliHLTTPCHoughTrack *tpt = (AliHLTTPCHoughTrack*)track;
  55   SetTrackParameters(tpt->GetKappa(),tpt->GetPsi(),tpt->GetWeight());
  56   SetEtaIndex(tpt->GetEtaIndex());
  57   SetEta(tpt->GetEta());
  58   SetTgl(tpt->GetTgl());
  59   SetPsi(tpt->GetPsi());
  60   SetPterr(tpt->GetPterr());
  61   SetTglerr(tpt->GetTglerr());
  62   SetPsierr(tpt->GetPsierr());
  63   SetCenterX(tpt->GetCenterX());
  64   SetCenterY(tpt->GetCenterY());
  65   SetFirstPoint(tpt->GetFirstPointX(),tpt->GetFirstPointY(),tpt->GetFirstPointZ());
  66   SetLastPoint(tpt->GetLastPointX(),tpt->GetLastPointY(),tpt->GetLastPointZ());
  67   SetCharge(tpt->GetCharge());
  68   SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
  69   SetSlice(tpt->GetSlice());
  70   SetHits(tpt->GetNHits(),(UInt_t *)tpt->GetHitNumbers());
  71   SetMCid(tpt->GetMCid());
  72   SetBinXY(tpt->GetBinX(),tpt->GetBinY(),tpt->GetSizeX(),tpt->GetSizeY());
  73   SetSector(tpt->GetSector());
  74   return;
  75
  76 //    fWeight = tpt->GetWeight();
  77 //    fDLine = tpt->GetDLine();
  78 //    fPsiLine = tpt->GetPsiLine();
  79 //    SetNHits(tpt->GetWeight());
  80 //    SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
  81 //    fIsHelix = false;
  82 }
  83
  84 Int_t AliHLTTPCHoughTrack::Compare(const AliHLTTrack *tpt) const
  85 {
  86   //Compare 2 hough tracks according to their weight
  87   AliHLTTPCHoughTrack *track = (AliHLTTPCHoughTrack*)tpt;
  88   if(track->GetWeight() < GetWeight()) return 1;
  89   if(track->GetWeight() > GetWeight()) return -1;
  90   return 0;
  91 }
  92
  93 void AliHLTTPCHoughTrack::SetEta(Double_t f)
  94 {
  95   //Set eta, and calculate fTanl, which is the tan of dipangle
  96
  97   fEta = f;
  98   Double_t theta = 2*atan(exp(-1.*fEta));
  99   Double_t dipangle = AliHLTTransform::PiHalf() - theta;
 100   Double_t tgl = tan(dipangle);
 101   SetTgl(tgl);
 102 }
 103
 104 void AliHLTTPCHoughTrack::UpdateToFirstRow()
 105 {
 106   //Update the track parameters to the point where track cross
 107   //its first padrow.`
 108
 109   //Get the crossing point with the first padrow:
 110   Float_t xyz[3];
 111   if(!GetCrossingPoint(GetFirstRow(),xyz))
 112     LOG(AliHLTLog::kWarning,"AliHLTTPCHoughTrack::UpdateToFirstRow()","Track parameters")
 113       <<AliHLTLog::kDec<<"Track does not cross padrow "<<GetFirstRow()<<" centerx "
 114       <<GetCenterX()<<" centery "<<GetCenterY()<<" Radius "<<GetRadius()<<" tgl "<<GetTgl()<<ENDLOG;
 115
 116   //printf("Track with eta %f tgl %f crosses at x %f y %f z %f on padrow %d\n",GetEta(),GetTgl(),xyz[0],xyz[1],xyz[2],GetFirstRow());
 117   //printf("Before: first %f %f %f tgl %f center %f %f charge %d\n",GetFirstPointX(),GetFirstPointY(),GetFirstPointZ(),GetTgl(),GetCenterX(),GetCenterY(),GetCharge());
 118
 119   Double_t radius = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]);
 120
 121   //Get the track parameters
 122
 123   /*
 124     Double_t x0    = GetR0() * cos(GetPhi0()) ;
 125     Double_t y0    = GetR0() * sin(GetPhi0()) ;
 126   */
 127   Double_t rc    = GetRadius();//fabs(GetPt()) / AliHLTTransform::GetBFieldValue();
 128   Double_t tPhi0 = GetPsi() + GetCharge() * AliHLTTransform::PiHalf() / abs(GetCharge()) ;
 129   Double_t xc    = GetCenterX();//x0 - rc * cos(tPhi0) ;
 130   Double_t yc    = GetCenterY();//y0 - rc * sin(tPhi0) ;
 131
 132   //Check helix and cylinder intersect
 133   Double_t fac1 = xc*xc + yc*yc ;
 134   Double_t sfac = sqrt( fac1 ) ;
 135
 136   if ( fabs(sfac-rc) > radius || fabs(sfac+rc) < radius ) {
 137     LOG(AliHLTLog::kError,"AliHLTTPCHoughTrack::UpdateToFirstRow","Tracks")<<AliHLTLog::kDec<<
 138       "Track does not intersect"<<ENDLOG;
 139     return;
 140   }
 141
 142   //Find intersection
 143   Double_t fac2   = (radius*radius + fac1 - rc*rc) / (2.00 * radius * sfac ) ;
 144   Double_t phi    = atan2(yc,xc) + GetCharge()*acos(fac2) ;
 145   Double_t td     = atan2(radius*sin(phi) - yc,radius*cos(phi) - xc) ;
 146
 147   //Intersection in z
 148   if ( td < 0 ) td = td + AliHLTTransform::TwoPi();
 149   Double_t deltat = fmod((-GetCharge()*td + GetCharge()*tPhi0),AliHLTTransform::TwoPi());
 150   if ( deltat < 0. ) deltat += AliHLTTransform::TwoPi();
 151   else if ( deltat > AliHLTTransform::TwoPi() ) deltat -= AliHLTTransform::TwoPi();
 152   Double_t z = GetZ0() + rc * GetTgl() * deltat ;
 153
 154   Double_t xExtra = radius * cos(phi) ;
 155   Double_t yExtra = radius * sin(phi) ;
 156
 157   Double_t tPhi = atan2(yExtra-yc,xExtra-xc);
 158
 159   //if ( tPhi < 0 ) tPhi += 2. * M_PI ;
 160   Double_t tPsi = tPhi - GetCharge() * AliHLTTransform::PiHalf() / abs(GetCharge()) ;
 161   if ( tPsi > AliHLTTransform::TwoPi() ) tPsi -= AliHLTTransform::TwoPi() ;
 162   else if ( tPsi < 0. ) tPsi += AliHLTTransform::TwoPi();
 163
 164   //And finally, update the track parameters
 165   SetR0(radius);
 166   SetPhi0(phi);
 167   SetZ0(z);
 168   SetPsi(tPsi);
 169   SetFirstPoint(xyz[0],xyz[1],z);
 170   //printf("After: first %f %f %f tgl %f center %f %f charge %d\n",GetFirstPointX(),GetFirstPointY(),GetFirstPointZ(),GetTgl(),GetCenterX(),GetCenterY(),GetCharge());
 171
 172   //printf("First point set %f %f %f\n",xyz[0],xyz[1],z);
 173
 174   //Also, set the coordinates of the point where track crosses last padrow:
 175   GetCrossingPoint(GetLastRow(),xyz);
 176   SetLastPoint(xyz[0],xyz[1],xyz[2]);
 177   //printf("last point %f %f %f\n",xyz[0],xyz[1],xyz[2]);
 178 }
 179
 180 void AliHLTTPCHoughTrack::SetTrackParameters(Double_t kappa,Double_t eangle,Int_t weight)
 181 {
 182   //Set track parameters - sort of ctor
 183   fWeight = weight;
 184   fMinDist = 100000;
 185   SetKappa(kappa);
 186   Double_t pt = fabs(AliHLTTransform::GetBFieldValue()/kappa);
 187   SetPt(pt);
 188   Double_t radius = 1/fabs(kappa);
 189   SetRadius(radius);
 190   SetFirstPoint(0,0,0);
 191   SetPsi(eangle); //Psi = emission angle when first point is vertex
 192   SetPhi0(0);     //not defined for vertex reference point
 193   SetR0(0);
 194   Double_t charge = -1.*kappa;
 195   SetCharge((Int_t)copysign(1.,charge));
 196   Double_t trackPhi0 = GetPsi() + charge*0.5*AliHLTTransform::Pi()/fabs(charge);
 197   Double_t xc = GetFirstPointX() - GetRadius() * cos(trackPhi0) ;
 198   Double_t yc = GetFirstPointY() - GetRadius() * sin(trackPhi0) ;
 199   SetCenterX(xc);
 200   SetCenterY(yc);
 201   SetNHits(1); //just for the trackarray IO
 202   fIsHelix = true;
 203 }
 204
 205 void AliHLTTPCHoughTrack::SetTrackParametersRow(Double_t alpha1,Double_t alpha2,Double_t eta,Int_t weight)
 206 {
 207   //Set track parameters for HoughTransformerRow
 208   //This includes curvature,emission angle and eta
 209   Double_t psi = atan((alpha1-alpha2)/(AliHLTHoughTransformerRow::GetBeta1()-AliHLTHoughTransformerRow::GetBeta2()));
 210   Double_t kappa = 2.0*(alpha1*cos(psi)-AliHLTHoughTransformerRow::GetBeta1()*sin(psi));
 211   SetTrackParameters(kappa,psi,weight);
 212
 213   Double_t zovr;
 214   Double_t etaparam1 = AliHLTHoughTransformerRow::GetEtaCalcParam1();
 215   Double_t etaparam2 = AliHLTHoughTransformerRow::GetEtaCalcParam2();
 216   if(eta>0)
 217     zovr = (etaparam1 - sqrt(etaparam1*etaparam1 - 4.*etaparam2*eta))/(2.*etaparam2);
 218   else
 219     zovr = -1.*(etaparam1 - sqrt(etaparam1*etaparam1 + 4.*etaparam2*eta))/(2.*etaparam2);
 220   Double_t r = sqrt(1.+zovr*zovr);
 221   Double_t exacteta = 0.5*log((1+zovr/r)/(1-zovr/r));
 222   SetEta(exacteta);
 223 }
 224
 225 void AliHLTTPCHoughTrack::SetLineParameters(Double_t psi,Double_t D,Int_t weight,Int_t *rowrange,Int_t /*ref_row*/)
 226 {
 227   //Initialize a track piece, not yet a track
 228   //Used in case of straight line transformation
 229
 230   //Transform line parameters to coordinate system of slice:
 231
 232   //  D = D + fTransform->Row2X(ref_row)*cos(psi);
 233
 234   fDLine = D;
 235   fPsiLine = psi;
 236   fWeight = weight;
 237   SetNHits(1);
 238   SetRowRange(rowrange[0],rowrange[1]);
 239   fIsHelix = false;
 240 }
 241
 242 void AliHLTTPCHoughTrack::SetBestMCid(Int_t mcid,Double_t mindist)
 243 {
 244   //Finds and set the closest mc label
 245   if(mindist < fMinDist)
 246     {
 247       fMinDist = mindist;
 248       SetMCid(mcid);
 249     }
 250 }
 251
 252 void AliHLTTPCHoughTrack::GetLineCrossingPoint(Int_t padrow,Float_t *xy)
 253 {
 254   //Returns the crossing point of the track with a given padrow
 255   if(fIsHelix)
 256     {
 257       printf("AliHLTTPCHoughTrack::GetLineCrossingPoint : Track is not a line\n");
 258       return;
 259     }
 260
 261   Float_t xhit = AliHLTTransform::Row2X(padrow) - AliHLTTransform::Row2X(GetFirstRow());
 262   Float_t a = -1/tan(fPsiLine);
 263   Float_t b = fDLine/sin(fPsiLine);
 264   Float_t yhit = a*xhit + b;
 265   xy[0] = xhit;
 266   xy[1] = yhit;
 267 }