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