TRD/AliTRDtrackOnline.cxx

   1 #include <limits>
   2
   3 #include "TObject.h"
   4 #include "TList.h"
   5 #include "TMath.h"
   6 #include "Math/Minimizer.h"
   7
   8 #include "AliLog.h"
   9 #include "AliVTrdTracklet.h"
  10 #include "AliTRDgeometry.h"
  11 #include "AliTRDpadPlane.h"
  12
  13 #include "AliTRDtrackOnline.h"
  14
  15 AliTRDgeometry *AliTRDtrackOnline::fgGeometry = new AliTRDgeometry();
  16
  17 AliTRDtrackOnline::AliTRDtrackOnline() :
  18   TObject(),
  19   fNTracklets(0),
  20   fTracklets(),
  21   fTrackParametrizations(),
  22   fErrorY(0.05)
  23 {
  24
  25 }
  26
  27
  28 AliTRDtrackOnline::~AliTRDtrackOnline()
  29 {
  30
  31 }
  32
  33
  34 void AliTRDtrackOnline::AddTracklet(AliVTrdTracklet *trkl)
  35 {
  36   if (fNTracklets == fgkMaxTracklets)
  37     return;
  38   else
  39     fTracklets[fNTracklets++] = trkl;
  40 }
  41
  42
  43 Bool_t AliTRDtrackOnline::Fit(ROOT::Math::Minimizer *minim)
  44 {
  45   // fit all attached parametrizations
  46
  47   Bool_t minSuccess = kFALSE;
  48
  49   if (minim) {
  50     minSuccess = kTRUE;
  51
  52     TIter param(&fTrackParametrizations);
  53
  54     while (AliTRDtrackParametrization *par = (AliTRDtrackParametrization*) param()) {
  55
  56       AliTRDtrackResiduals res(this, par);
  57       minim->Clear();
  58       minim->SetFunction(res);
  59       par->SetParams(minim);
  60       minSuccess &= minim->Minimize();
  61       par->GetParams(minim);
  62     }
  63   }
  64
  65   return minSuccess;
  66 }
  67
  68
  69 AliTRDtrackPosition AliTRDtrackOnline::ExtrapolateToLayer(Int_t /* layer */)
  70 {
  71   Int_t maxLayer = -1;
  72   AliVTrdTracklet *trklBest = 0x0;
  73   for (Int_t iTracklet = fNTracklets-1; iTracklet > -1; iTracklet--) {
  74     AliVTrdTracklet *trkl = (AliVTrdTracklet*) fTracklets[iTracklet];
  75     if (trkl->GetDetector() % 6 >= maxLayer) {
  76       maxLayer = trkl->GetDetector() % 6;
  77       trklBest = trkl;
  78     }
  79   }
  80   if (trklBest)
  81     return AliTRDtrackPosition(trklBest->GetLocalY(), GetZ(trklBest));
  82   else {
  83     AliFatal("No tracklet in this track");
  84     return AliTRDtrackPosition(std::numeric_limits<float>::quiet_NaN(), std::numeric_limits<float>::quiet_NaN());
  85   }
  86 }
  87
  88
  89 void AliTRDtrackOnline::Print(Option_t* /* option */) const
  90 {
  91   printf("track with %i tracklets:\n", GetNTracklets());
  92   for (Int_t iTracklet = 0; iTracklet < fNTracklets; iTracklet++) {
  93     printf("  0x%08x %i %4.1f %4.1f\n",
  94            ((AliVTrdTracklet*) fTracklets[iTracklet])->GetTrackletWord(),
  95            ((AliVTrdTracklet*) fTracklets[iTracklet])->GetDetector() % 6,
  96            ((AliVTrdTracklet*) fTracklets[iTracklet])->GetLocalY(),
  97            GetZ((AliVTrdTracklet*) fTracklets[iTracklet]));
  98   }
  99   TIter next(&fTrackParametrizations);
 100   while (AliTRDtrackParametrization *param = (AliTRDtrackParametrization*) next()) {
 101     param->Print();
 102   }
 103 }
 104
 105 AliTRDtrackPosition::AliTRDtrackPosition(Float_t y, Float_t z, Float_t dy) :
 106   TObject(),
 107   fY(y),
 108   fZ(z),
 109   fDy(dy)
 110 {
 111
 112 }
 113
 114 AliTRDtrackPosition::~AliTRDtrackPosition()
 115 {
 116
 117 }
 118
 119 Float_t AliTRDtrackPosition::Distance(AliVTrdTracklet *trkl) const
 120 {
 121   return TMath::Hypot(trkl->GetLocalY() - fY, AliTRDtrackOnline::GetZ(trkl) - fZ);
 122 }
 123
 124
 125 AliTRDtrackParametrization::AliTRDtrackParametrization(const char* name, const char* title) :
 126   TNamed(name, title),
 127   fFitGood(kFALSE)
 128 {
 129
 130 }
 131
 132 AliTRDtrackParametrizationStraightLine::AliTRDtrackParametrizationStraightLine() :
 133   AliTRDtrackParametrization("straight line", "straight line"),
 134   fOffsetY(0),
 135   fSlopeY(0),
 136   fOffsetZ(0),
 137   fSlopeZ(0)
 138 {
 139
 140 }
 141
 142 AliTRDtrackParametrizationStraightLine::AliTRDtrackParametrizationStraightLine(Double_t offsetY, Double_t slopeY,
 143                                                                                Double_t offsetZ, Double_t slopeZ) :
 144   AliTRDtrackParametrization("straight line", Form("straight line: y = %4.2f + %4.2f * x, z = %4.2f + %4.2f *x",
 145                                                    offsetY, slopeY, offsetZ, slopeZ)),
 146   fOffsetY(offsetY),
 147   fSlopeY(slopeY),
 148   fOffsetZ(offsetZ),
 149   fSlopeZ(slopeZ)
 150 {
 151
 152 }
 153
 154 void AliTRDtrackParametrizationStraightLine::SetParams(ROOT::Math::Minimizer * minim)
 155 {
 156   minim->SetVariable(0, "offsety", 0., 0.1);
 157   minim->SetVariable(1, "slopey", 0., 0.1);
 158   // minim->SetVariable(2, "offsetz", 0., 0.1);
 159   minim->SetFixedVariable(2, "offsetz", 0.);
 160   minim->SetVariable(3, "slopez", 0., 0.1);
 161 }
 162
 163 void AliTRDtrackParametrizationStraightLine::GetParams(ROOT::Math::Minimizer * minim)
 164 {
 165   fOffsetY = minim->X()[0];
 166   fSlopeY  = minim->X()[1];
 167   fOffsetZ = minim->X()[2];
 168   fSlopeZ  = minim->X()[3];
 169 }
 170
 171 void AliTRDtrackParametrizationStraightLine::SetValues(const Double_t *par)
 172 {
 173   fOffsetY = par[0];
 174   fSlopeY  = par[1];
 175   fOffsetZ = par[2];
 176   fSlopeZ  = par[3];
 177 }
 178
 179 AliTRDtrackPosition AliTRDtrackParametrizationStraightLine::ExtrapolateToLayer(Int_t layer)
 180 {
 181   Float_t y = fOffsetY + fSlopeY * AliTRDtrackOnline::fgGeometry->GetTime0(layer);
 182   Float_t z = fOffsetZ + fSlopeZ * AliTRDtrackOnline::fgGeometry->GetTime0(layer);
 183   return AliTRDtrackPosition(y, z, fSlopeY*3.);
 184 }
 185
 186 AliTRDtrackPosition AliTRDtrackParametrizationStraightLine::ExtrapolateToX(Float_t x)
 187 {
 188   Float_t y = fOffsetY + fSlopeY * x;
 189   Float_t z = fOffsetZ + fSlopeZ * x;
 190   return AliTRDtrackPosition(y, z, fSlopeY*3.);
 191 }
 192
 193 void AliTRDtrackParametrizationStraightLine::Print(Option_t * /* option */) const
 194 {
 195   printf("straight line: offsetY = %4.1f, slopeY = %4.1f; offsetZ = %4.1f, slopeZ = %4.1f\n",
 196          fOffsetY, fSlopeY, fOffsetZ, fSlopeZ);
 197 }
 198
 199 AliTRDtrackParametrizationCurved::AliTRDtrackParametrizationCurved() :
 200   AliTRDtrackParametrization("helix", "helix"),
 201   fRadiusInv(0.),
 202   fOffsetY(0.),
 203   fOffsetZ(0.),
 204   fSlopeZ(0.),
 205   fOffsetX(300.)
 206 {
 207
 208 }
 209
 210
 211 void AliTRDtrackParametrizationCurved::SetParams(ROOT::Math::Minimizer * minim)
 212 {
 213   minim->SetVariable(0, "offsety", 0., 0.1);
 214   minim->SetVariable(1, "invradius", 0., 0.1);
 215   // minim->SetVariable(2, "offsetz", 1., 0.1);
 216   minim->SetFixedVariable(2, "offsetz", 0.);
 217   minim->SetVariable(3, "slopez", 0., 0.1);
 218 }
 219
 220
 221 void AliTRDtrackParametrizationCurved::GetParams(ROOT::Math::Minimizer * minim)
 222 {
 223   this->SetValues(minim->X());
 224 }
 225
 226
 227 void AliTRDtrackParametrizationCurved::SetValues(const Double_t *par)
 228 {
 229   fOffsetY    = par[0];
 230   fRadiusInv  = par[1];
 231   fOffsetZ    = par[2];
 232   fSlopeZ     = par[3];
 233 }
 234
 235 AliTRDtrackPosition AliTRDtrackParametrizationCurved::ExtrapolateToLayer(Int_t layer)
 236 {
 237   return ExtrapolateToX(AliTRDtrackOnline::fgGeometry->GetTime0(layer));
 238 }
 239
 240 AliTRDtrackPosition AliTRDtrackParametrizationCurved::ExtrapolateToX(Float_t x)
 241 {
 242   Double_t yext1 = GetY(x);
 243   Double_t yext2 = GetY(x + 3.);
 244
 245   Double_t zext = fOffsetZ + fSlopeZ * x;
 246
 247   return AliTRDtrackPosition(yext1, zext, yext2-yext1);
 248 }
 249
 250 Float_t AliTRDtrackParametrizationCurved::GetY(Float_t x)
 251 {
 252  Double_t yext = 0.;
 253   // use Taylor expansion for small 1/R
 254   if (TMath::Abs(fRadiusInv) < 1.) {
 255     // offset
 256     yext  = fOffsetY * x/fOffsetX;
 257     // linear term
 258     yext += - (fOffsetX - x) * x * fRadiusInv /
 259       (2 * (fOffsetX*1./TMath::Sqrt(fOffsetX*fOffsetX + fOffsetY*fOffsetY)) *
 260        (fOffsetX*1./TMath::Sqrt(fOffsetX*fOffsetX + fOffsetY*fOffsetY)) *
 261        (fOffsetX*1./TMath::Sqrt(fOffsetX*fOffsetX + fOffsetY*fOffsetY)));
 262   }
 263   else {
 264     Double_t disc = 1./(fOffsetX*fOffsetX + fOffsetY*fOffsetY) - fRadiusInv*fRadiusInv/4.;
 265     if (disc < 0) {
 266       AliError("Discriminant < 0");
 267       return 1000.;
 268     }
 269     yext = TMath::Sqrt(disc) -
 270       TMath::Sqrt((fRadiusInv*fOffsetY/2. + fOffsetX * TMath::Sqrt(disc)) *
 271                   (fRadiusInv*fOffsetY/2. + fOffsetX * TMath::Sqrt(disc)) /
 272                   (fOffsetX*fOffsetX) -
 273                   fRadiusInv*fRadiusInv/(fOffsetX*fOffsetX)* x*x +
 274                   fRadiusInv*fRadiusInv/fOffsetX * x +
 275                   2 * fRadiusInv * fOffsetY * (x - fOffsetX)/(fOffsetX*fOffsetX) * TMath::Sqrt(disc));
 276     yext = fOffsetY/2. - fOffsetX * yext / fRadiusInv;
 277   }
 278
 279   return yext;
 280 }
 281
 282 void AliTRDtrackParametrizationCurved::Print(Option_t * /* option */) const
 283 {
 284   printf("helix curve: 1/R = %f, y = %4.1f\n", fRadiusInv, fOffsetY);
 285 }
 286
 287
 288 AliTRDtrackResiduals::AliTRDtrackResiduals(const AliTRDtrackOnline *track, AliTRDtrackParametrization *param) :
 289   ROOT::Math::IBaseFunctionMultiDim(),
 290   fTrack(track),
 291   fParam(param)
 292 {
 293
 294 }
 295
 296 AliTRDtrackResiduals::AliTRDtrackResiduals(const AliTRDtrackResiduals &rhs) :
 297   ROOT::Math::IBaseFunctionMultiDim(rhs),
 298   fTrack(rhs.fTrack),
 299   fParam(rhs.fParam)
 300 {
 301
 302 }
 303
 304 AliTRDtrackResiduals& AliTRDtrackResiduals::operator=(const AliTRDtrackResiduals &rhs)
 305 {
 306   if (&rhs != this) {
 307     ROOT::Math::IBaseFunctionMultiDim::operator=(rhs);
 308     fTrack = rhs.fTrack;
 309     fParam = rhs.fParam;
 310   }
 311
 312   return *this;
 313 }
 314
 315 AliTRDtrackResiduals* AliTRDtrackResiduals::Clone() const
 316 {
 317   return new AliTRDtrackResiduals(*this);
 318 }
 319
 320 Double_t AliTRDtrackResiduals::DoEval(const Double_t *par) const
 321 {
 322   // calculate chi2 for the given values for the parametrization
 323
 324   // initialisation
 325   Float_t deltaY = 0.;
 326   Float_t deltaZ = 0.;
 327   Float_t chi2 = 0.;
 328
 329   // actually set the values for the parametrization
 330   fParam->SetValues(par);
 331
 332   // loop over all contributing tracklets
 333   for (Int_t iTracklet = 0; iTracklet < fTrack->GetNTracklets(); iTracklet++) {
 334     AliVTrdTracklet *trkl = fTrack->GetTracklet(iTracklet);
 335
 336     // Int_t layer = trkl->GetDetector() % 6;
 337
 338     AliTRDtrackPosition pos = fParam->ExtrapolateToX(AliTRDtrackOnline::GetX(trkl));
 339     Float_t yext = pos.GetY();
 340     Float_t zext = pos.GetZ();
 341
 342     AliTRDpadPlane *pp = fgGeometry->GetPadPlane(trkl->GetDetector());
 343     Float_t zlen = 0.5 * pp->GetRowSize(trkl->GetBinZ());
 344     Float_t zpad = AliTRDtrackOnline::GetZ(trkl);
 345     Float_t zrel = zext - zpad;
 346     if (zrel > zlen)
 347       zrel = zlen;
 348     else if (zrel < -zlen)
 349       zrel = -zlen;
 350
 351     Float_t ycorr = trkl->GetLocalY() + TMath::Tan(TMath::Pi()/180.*pp->GetTiltingAngle()) * zrel;
 352
 353     deltaY = ycorr - yext;
 354     deltaZ = zpad  - zext;
 355     deltaY /= fTrack->GetErrorY();
 356     deltaZ /= pp->GetRowSize(trkl->GetBinZ()) / TMath::Sqrt(12.);
 357     // printf("for tracklet %i: deltaY = %f, deltaZ = %f\n", iTracklet, deltaY, deltaZ);
 358
 359     chi2 += deltaY*deltaY + deltaZ*deltaZ;
 360   }
 361
 362 //   printf("chi2 = %f\n", chi2);
 363   return chi2;
 364 }