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