]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - TRD/AliTRDseed.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
coverity fixes
[u/mrichter/AliRoot.git] / TRD / AliTRDseed.cxx
... / ...
CommitLineData
1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
20// The seed of a local TRD track //
21// //
22///////////////////////////////////////////////////////////////////////////////
23
24#include "TMath.h"
25#include "TLinearFitter.h"
26
27#include "AliMathBase.h"
28
29#include "AliTRDseed.h"
30#include "AliTRDcluster.h"
31#include "AliTRDtracker.h"
32#include "AliTRDtrackerV1.h"
33
34ClassImp(AliTRDseed)
35
36//_____________________________________________________________________________
37AliTRDseed::AliTRDseed()
38 :TObject()
39 ,fTilt(0)
40 ,fPadLength(0)
41 ,fX0(0)
42 ,fSigmaY(0)
43 ,fSigmaY2(0)
44 ,fMeanz(0)
45 ,fZProb(0)
46 ,fN(0)
47 ,fN2(0)
48 ,fNUsed(0)
49 ,fFreq(0)
50 ,fNChange(0)
51 ,fMPads(0)
52 ,fC(0)
53 ,fCC(0)
54 ,fChi2(1.0e10)
55 ,fChi2Z(1.0e10)
56{
57 //
58 // Default constructor
59 //
60
61 for (Int_t i = 0; i < knTimebins; i++) {
62 fX[i] = 0; // x position
63 fY[i] = 0; // y position
64 fZ[i] = 0; // z position
65 fIndexes[i] = 0; // Indexes
66 fClusters[i] = NULL; // Clusters
67 fUsable[i] = 0; // Indication - usable cluster
68 }
69
70 for (Int_t i = 0; i < 2; i++) {
71 fYref[i] = 0; // Reference y
72 fZref[i] = 0; // Reference z
73 fYfit[i] = 0; // Y fit position +derivation
74 fYfitR[i] = 0; // Y fit position +derivation
75 fZfit[i] = 0; // Z fit position
76 fZfitR[i] = 0; // Z fit position
77 fLabels[i] = 0; // Labels
78 }
79
80}
81
82//_____________________________________________________________________________
83AliTRDseed::AliTRDseed(const AliTRDseed &s)
84 :TObject(s)
85 ,fTilt(s.fTilt)
86 ,fPadLength(s.fPadLength)
87 ,fX0(s.fX0)
88 ,fSigmaY(s.fSigmaY)
89 ,fSigmaY2(s.fSigmaY2)
90 ,fMeanz(s.fMeanz)
91 ,fZProb(s.fZProb)
92 ,fN(s.fN)
93 ,fN2(s.fN2)
94 ,fNUsed(s.fNUsed)
95 ,fFreq(s.fFreq)
96 ,fNChange(s.fNChange)
97 ,fMPads(s.fMPads)
98 ,fC(s.fC)
99 ,fCC(s.fCC)
100 ,fChi2(s.fChi2)
101 ,fChi2Z(s.fChi2Z)
102{
103 //
104 // Copy constructor
105 //
106
107 for (Int_t i = 0; i < knTimebins; i++) {
108 fX[i] = s.fX[i]; // x position
109 fY[i] = s.fY[i]; // y position
110 fZ[i] = s.fZ[i]; // z position
111 fIndexes[i] = s.fIndexes[i]; // Indexes
112 fClusters[i] = s.fClusters[i]; // Clusters
113 fUsable[i] = s.fUsable[i]; // Indication - usable cluster
114 }
115
116 for (Int_t i = 0; i < 2; i++) {
117 fYref[i] = s.fYref[i]; // Reference y
118 fZref[i] = s.fZref[i]; // Reference z
119 fYfit[i] = s.fYfit[i]; // Y fit position +derivation
120 fYfitR[i] = s.fYfitR[i]; // Y fit position +derivation
121 fZfit[i] = s.fZfit[i]; // Z fit position
122 fZfitR[i] = s.fZfitR[i]; // Z fit position
123 fLabels[i] = s.fLabels[i]; // Labels
124 }
125
126}
127
128//_____________________________________________________________________________
129AliTRDseed &AliTRDseed::operator=(const AliTRDseed &s)
130{
131 //
132 // Assignment operator
133 //
134
135 if (this != &s) {
136 ((AliTRDseed &) s).Copy(*this);
137 }
138
139 return *this;
140
141}
142
143//_____________________________________________________________________________
144void AliTRDseed::Copy(TObject &o) const
145{
146 //printf("AliTRDseed::Copy()\n");
147
148 AliTRDseed &seed = (AliTRDseed &)o;
149
150 seed.fTilt = fTilt;
151 seed.fPadLength = fPadLength;
152 seed.fX0 = fX0;
153 seed.fSigmaY = fSigmaY;
154 seed.fSigmaY2 = fSigmaY2;
155 seed.fMeanz = fMeanz;
156 seed.fZProb = fZProb;
157 seed.fN = fN;
158 seed.fN2 = fN2;
159 seed.fNUsed = fNUsed;
160 seed.fFreq = fFreq;
161 seed.fNChange = fNChange;
162 seed.fMPads = fMPads;
163 seed.fC = fC;
164 seed.fCC = fCC;
165 seed.fChi2 = fChi2;
166 seed.fChi2Z = fChi2Z;
167 for (Int_t i = 0; i < knTimebins; i++) {
168 seed.fX[i] = fX[i];
169 seed.fY[i] = fY[i];
170 seed.fZ[i] = fZ[i];
171 seed.fIndexes[i] = fIndexes[i];
172 seed.fClusters[i] = fClusters[i];
173 seed.fUsable[i] = fUsable[i];
174 }
175
176 for (Int_t i = 0; i < 2; i++) {
177 seed.fYref[i] = fYref[i];
178 seed.fZref[i] = fZref[i];
179 seed.fYfit[i] = fYfit[i];
180 seed.fYfitR[i] = fYfitR[i];
181 seed.fZfit[i] = fZfit[i];
182 seed.fZfitR[i] = fZfitR[i];
183 seed.fLabels[i] = fLabels[i];
184 }
185
186 TObject::Copy(seed);
187
188}
189
190//_____________________________________________________________________________
191void AliTRDseed::Reset()
192{
193 //
194 // Reset seed
195 //
196
197 for (Int_t i = 0; i < knTimebins; i++) {
198 fX[i] = 0; // X position
199 fY[i] = 0; // Y position
200 fZ[i] = 0; // Z position
201 fIndexes[i] = 0; // Indexes
202 fClusters[i] = NULL; // Clusters
203 fUsable[i] = kFALSE;
204 }
205
206 for (Int_t i = 0; i < 2; i++) {
207 fYref[i] = 0; // Reference y
208 fZref[i] = 0; // Reference z
209 fYfit[i] = 0; // Y fit position +derivation
210 fYfitR[i] = 0; // Y fit position +derivation
211 fZfit[i] = 0; // Z fit position
212 fZfitR[i] = 0; // Z fit position
213 fLabels[i] = -1; // Labels
214 }
215 fSigmaY = 0; // "Robust" sigma in y
216 fSigmaY2 = 0; // "Robust" sigma in y
217 fMeanz = 0; // Mean vaue of z
218 fZProb = 0; // Max probbable z
219 fMPads = 0;
220 fN = 0; // Number of associated clusters
221 fN2 = 0; // Number of not crossed
222 fNUsed = 0; // Number of used clusters
223 fNChange = 0; // Change z counter
224
225}
226
227//_____________________________________________________________________________
228void AliTRDseed::CookLabels()
229{
230 //
231 // Cook 2 labels for seed
232 //
233
234 Int_t labels[200];
235 Int_t out[200];
236 Int_t nlab = 0;
237
238 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
239 if (!fClusters[i]) continue;
240 for (Int_t ilab = 0; ilab < 3; ilab++) {
241 if (fClusters[i]->GetLabel(ilab) >= 0) {
242 labels[nlab] = fClusters[i]->GetLabel(ilab);
243 nlab++;
244 }
245 }
246 }
247
248 Int_t nlab2 = AliTRDtracker::Freq(nlab,labels,out,kTRUE);
249 fLabels[0] = out[0];
250 if ((nlab2 > 1) &&
251 (out[3] > 1)) {
252 fLabels[1] = out[2];
253 }
254
255}
256
257//_____________________________________________________________________________
258void AliTRDseed::UseClusters()
259{
260 //
261 // Use clusters
262 //
263
264 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
265 if (!fClusters[i]) continue;
266 if (!(fClusters[i]->IsUsed())) fClusters[i]->Use();
267 }
268
269}
270
271//_____________________________________________________________________________
272void AliTRDseed::Update()
273{
274 //
275 // Update the seed.
276 //
277
278
279
280 // linear fit on the y direction with respect to the reference direction.
281 // The residuals for each x (x = xc - x0) are deduced from:
282 // dy = y - yt (1)
283 // the tilting correction is written :
284 // y = yc + h*(zc-zt) (2)
285 // yt = y0+dy/dx*x (3)
286 // zt = z0+dz/dx*x (4)
287 // from (1),(2),(3) and (4)
288 // dy = yc - y0 - (dy/dx + h*dz/dx)*x + h*(zc-z0)
289 // the last term introduces the correction on y direction due to tilting pads. There are 2 ways to account for this:
290 // 1. use tilting correction for calculating the y
291 // 2. neglect tilting correction here and account for it in the error parametrization of the tracklet.
292
293 const Float_t kRatio = 0.8;
294 const Int_t kClmin = 5;
295 const Float_t kmaxtan = 2;
296
297 if (TMath::Abs(fYref[1]) > kmaxtan){
298 //printf("Exit: Abs(fYref[1]) = %3.3f, kmaxtan = %3.3f\n", TMath::Abs(fYref[1]), kmaxtan);
299 return; // Track inclined too much
300 }
301
302 Float_t sigmaexp = 0.05 + TMath::Abs(fYref[1] * 0.25); // Expected r.m.s in y direction
303 Float_t ycrosscor = fPadLength * fTilt * 0.5; // Y correction for crossing
304 fNChange = 0;
305
306 Double_t sumw;
307 Double_t sumwx;
308 Double_t sumwx2;
309 Double_t sumwy;
310 Double_t sumwxy;
311 Double_t sumwz;
312 Double_t sumwxz;
313
314 // Buffering: Leave it constant fot Performance issues
315 Int_t zints[knTimebins]; // Histograming of the z coordinate
316 // Get 1 and second max probable coodinates in z
317 Int_t zouts[2*knTimebins];
318 Float_t allowedz[knTimebins]; // Allowed z for given time bin
319 memset(allowedz, 0, knTimebins*sizeof(Float_t));
320 Float_t yres[knTimebins]; // Residuals from reference
321 //Float_t anglecor = fTilt * fZref[1]; // Correction to the angle
322
323
324 fN = 0;
325 fN2 = 0;
326 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins(); i++) {
327 yres[i] = 10000.0;
328 if (!fClusters[i]) continue;
329 if(!fClusters[i]->IsInChamber()) continue;
330 // Residual y
331 //yres[i] = fY[i] - fYref[0] - (fYref[1] + anglecor) * fX[i] + fTilt*(fZ[i] - fZref[0]);
332 yres[i] = fY[i] - fTilt*(fZ[i] - (fZref[0] - fX[i]*fZref[1]));
333 zints[fN] = Int_t(fZ[i]);
334 fN++;
335 }
336
337 if (fN < kClmin){
338 //printf("Exit fN < kClmin: fN = %d\n", fN);
339 return;
340 }
341 Int_t nz = AliTRDtracker::Freq(fN, zints, zouts, kFALSE);
342 fZProb = zouts[0];
343 if (nz <= 1) zouts[3] = 0;
344 if (zouts[1] + zouts[3] < kClmin) {
345 //printf("Exit zouts[1] = %d, zouts[3] = %d\n",zouts[1],zouts[3]);
346 return;
347 }
348
349 // Z distance bigger than pad - length
350 if (TMath::Abs(zouts[0]-zouts[2]) > 12.0) zouts[3] = 0;
351
352 Int_t breaktime = -1;
353 Bool_t mbefore = kFALSE;
354 Int_t cumul[knTimebins][2];
355 Int_t counts[2] = { 0, 0 };
356
357 if (zouts[3] >= 3) {
358
359 //
360 // Find the break time allowing one chage on pad-rows
361 // with maximal number of accepted clusters
362 //
363 fNChange = 1;
364 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins(); i++) {
365 cumul[i][0] = counts[0];
366 cumul[i][1] = counts[1];
367 if (TMath::Abs(fZ[i]-zouts[0]) < 2) counts[0]++;
368 if (TMath::Abs(fZ[i]-zouts[2]) < 2) counts[1]++;
369 }
370 Int_t maxcount = 0;
371 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins(); i++) {
372 Int_t after = cumul[AliTRDtrackerV1::GetNTimeBins()][0] - cumul[i][0];
373 Int_t before = cumul[i][1];
374 if (after + before > maxcount) {
375 maxcount = after + before;
376 breaktime = i;
377 mbefore = kFALSE;
378 }
379 after = cumul[AliTRDtrackerV1::GetNTimeBins()-1][1] - cumul[i][1];
380 before = cumul[i][0];
381 if (after + before > maxcount) {
382 maxcount = after + before;
383 breaktime = i;
384 mbefore = kTRUE;
385 }
386 }
387
388 breaktime -= 1;
389
390 }
391
392 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
393 if (i > breaktime) allowedz[i] = mbefore ? zouts[2] : zouts[0];
394 if (i <= breaktime) allowedz[i] = (!mbefore) ? zouts[2] : zouts[0];
395 }
396
397 if (((allowedz[0] > allowedz[AliTRDtrackerV1::GetNTimeBins()]) && (fZref[1] < 0)) ||
398 ((allowedz[0] < allowedz[AliTRDtrackerV1::GetNTimeBins()]) && (fZref[1] > 0))) {
399 //
400 // Tracklet z-direction not in correspondance with track z direction
401 //
402 fNChange = 0;
403 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
404 allowedz[i] = zouts[0]; // Only longest taken
405 }
406 }
407
408 if (fNChange > 0) {
409 //
410 // Cross pad -row tracklet - take the step change into account
411 //
412 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
413 if (!fClusters[i]) continue;
414 if(!fClusters[i]->IsInChamber()) continue;
415 if (TMath::Abs(fZ[i] - allowedz[i]) > 2) continue;
416 // Residual y
417 //yres[i] = fY[i] - fYref[0] - (fYref[1] + anglecor) * fX[i] /*+ fTilt*(fZ[i] - fZref[0])*/;
418 yres[i] = fY[i] - fTilt*(fZ[i] - (fZref[0] - fX[i]*fZref[1]));
419/* if (TMath::Abs(fZ[i] - fZProb) > 2) {
420 if (fZ[i] > fZProb) yres[i] += fTilt * fPadLength;
421 if (fZ[i] < fZProb) yres[i] -= fTilt * fPadLength;
422 }*/
423 }
424 }
425
426 Double_t yres2[knTimebins];
427 Double_t mean;
428 Double_t sigma;
429 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
430 if (!fClusters[i]) continue;
431 if(!fClusters[i]->IsInChamber()) continue;
432 if (TMath::Abs(fZ[i] - allowedz[i]) > 2) continue;
433 yres2[fN2] = yres[i];
434 fN2++;
435 }
436 if (fN2 < kClmin) {
437 //printf("Exit fN2 < kClmin: fN2 = %d\n", fN2);
438 fN2 = 0;
439 return;
440 }
441 AliMathBase::EvaluateUni(fN2,yres2,mean,sigma, Int_t(fN2*kRatio-2.));
442 if (sigma < sigmaexp * 0.8) {
443 sigma = sigmaexp;
444 }
445 fSigmaY = sigma;
446
447 // Reset sums
448 sumw = 0;
449 sumwx = 0;
450 sumwx2 = 0;
451 sumwy = 0;
452 sumwxy = 0;
453 sumwz = 0;
454 sumwxz = 0;
455
456 fN2 = 0;
457 fMeanz = 0;
458 fMPads = 0;
459
460 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
461
462 fUsable[i] = kFALSE;
463 if (!fClusters[i]) continue;
464 if (!fClusters[i]->IsInChamber()) continue;
465 if (TMath::Abs(fZ[i] - allowedz[i]) > 2){fClusters[i] = NULL; continue;}
466 if (TMath::Abs(yres[i] - mean) > 4.0 * sigma){fClusters[i] = NULL; continue;}
467 fUsable[i] = kTRUE;
468 fN2++;
469 fMPads += fClusters[i]->GetNPads();
470 Float_t weight = 1.0;
471 if (fClusters[i]->GetNPads() > 4) weight = 0.5;
472 if (fClusters[i]->GetNPads() > 5) weight = 0.2;
473
474
475 Double_t x = fX[i];
476 //printf("x = %7.3f dy = %7.3f fit %7.3f\n", x, yres[i], fY[i]-yres[i]);
477
478 sumw += weight;
479 sumwx += x * weight;
480 sumwx2 += x*x * weight;
481 sumwy += weight * yres[i];
482 sumwxy += weight * (yres[i]) * x;
483 sumwz += weight * fZ[i];
484 sumwxz += weight * fZ[i] * x;
485
486 }
487
488 if (fN2 < kClmin){
489 //printf("Exit fN2 < kClmin(2): fN2 = %d\n",fN2);
490 fN2 = 0;
491 return;
492 }
493 fMeanz = sumwz / sumw;
494 Float_t correction = 0;
495 if (fNChange > 0) {
496 // Tracklet on boundary
497 if (fMeanz < fZProb) correction = ycrosscor;
498 if (fMeanz > fZProb) correction = -ycrosscor;
499 }
500
501 Double_t det = sumw * sumwx2 - sumwx * sumwx;
502 fYfitR[0] = (sumwx2 * sumwy - sumwx * sumwxy) / det;
503 fYfitR[1] = (sumw * sumwxy - sumwx * sumwy) / det;
504
505 fSigmaY2 = 0;
506 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
507 if (!fUsable[i]) continue;
508 Float_t delta = yres[i] - fYfitR[0] - fYfitR[1] * fX[i];
509 fSigmaY2 += delta*delta;
510 }
511 fSigmaY2 = TMath::Sqrt(fSigmaY2 / Float_t(fN2-2));
512 // TEMPORARY UNTIL covariance properly calculated
513 fSigmaY2 = TMath::Max(fSigmaY2, Float_t(.1));
514
515 fZfitR[0] = (sumwx2 * sumwz - sumwx * sumwxz) / det;
516 fZfitR[1] = (sumw * sumwxz - sumwx * sumwz) / det;
517 fZfit[0] = (sumwx2 * sumwz - sumwx * sumwxz) / det;
518 fZfit[1] = (sumw * sumwxz - sumwx * sumwz) / det;
519// fYfitR[0] += fYref[0] + correction;
520// fYfitR[1] += fYref[1];
521 fYfit[0] = fYfitR[0];
522 fYfit[1] = -fYfitR[1];
523
524 //printf("y0 = %7.3f tgy = %7.3f z0 = %7.3f tgz = %7.3f \n", fYfitR[0], fYfitR[1], fZfitR[0], fZfitR[1]);
525
526 UpdateUsed();
527
528}
529
530//_____________________________________________________________________________
531void AliTRDseed::UpdateUsed()
532{
533 //
534 // Update used seed
535 //
536
537 fNUsed = 0;
538 for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins(); i++) {
539 if (!fClusters[i]) continue;
540 if(!fUsable[i]) continue;
541 if ((fClusters[i]->IsUsed())) fNUsed++;
542 }
543
544}
545
546//_____________________________________________________________________________
547Float_t AliTRDseed::FitRiemanTilt(AliTRDseed * cseed, Bool_t terror)
548{
549 //
550 // Fit the Rieman tilt
551 //
552
553 // Fitting with tilting pads - kz not fixed
554 TLinearFitter fitterT2(4,"hyp4");
555 fitterT2.StoreData(kTRUE);
556
557 Float_t xref2 = (cseed[2].fX0 + cseed[3].fX0) * 0.5; // Reference x0 for z
558
559 Int_t npointsT = 0;
560 fitterT2.ClearPoints();
561
562 for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
563
564 if (!cseed[iLayer].IsOK()) continue;
565 Double_t tilt = cseed[iLayer].fTilt;
566
567 for (Int_t itime = 0; itime < AliTRDtrackerV1::GetNTimeBins()+1; itime++) {
568
569 if (!cseed[iLayer].fUsable[itime]) continue;
570 // x relative to the midle chamber
571 Double_t x = cseed[iLayer].fX[itime] + cseed[iLayer].fX0 - xref2;
572 Double_t y = cseed[iLayer].fY[itime];
573 Double_t z = cseed[iLayer].fZ[itime];
574
575 //
576 // Tilted rieman
577 //
578 Double_t uvt[6];
579 Double_t x2 = cseed[iLayer].fX[itime] + cseed[iLayer].fX0; // Global x
580 Double_t t = 1.0 / (x2*x2 + y*y);
581 uvt[1] = t;
582 uvt[0] = 2.0 * x2 * uvt[1];
583 uvt[2] = 2.0 * tilt * uvt[1];
584 uvt[3] = 2.0 * tilt *uvt[1] * x;
585 uvt[4] = 2.0 * (y + tilt * z) * uvt[1];
586
587 Double_t error = 2.0 * uvt[1];
588 if (terror) {
589 error *= cseed[iLayer].fSigmaY;
590 }
591 else {
592 error *= 0.2; //Default error
593 }
594 fitterT2.AddPoint(uvt,uvt[4],error);
595 npointsT++;
596
597 }
598
599 }
600
601 fitterT2.Eval();
602 Double_t rpolz0 = fitterT2.GetParameter(3);
603 Double_t rpolz1 = fitterT2.GetParameter(4);
604
605 //
606 // Linear fitter - not possible to make boundaries
607 // non accept non possible z and dzdx combination
608 //
609 Bool_t acceptablez = kTRUE;
610 for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
611 if (!cseed[iLayer].IsOK()) continue;
612 Double_t zT2 = rpolz0 + rpolz1 * (cseed[iLayer].fX0 - xref2);
613 if (TMath::Abs(cseed[iLayer].fZProb - zT2) > cseed[iLayer].fPadLength * 0.5 + 1.0) acceptablez = kFALSE;
614 }
615 if (!acceptablez) {
616 Double_t zmf = cseed[2].fZref[0] + cseed[2].fZref[1] * (xref2 - cseed[2].fX0);
617 Double_t dzmf = (cseed[2].fZref[1] + cseed[3].fZref[1]) * 0.5;
618 fitterT2.FixParameter(3,zmf);
619 fitterT2.FixParameter(4,dzmf);
620 fitterT2.Eval();
621 fitterT2.ReleaseParameter(3);
622 fitterT2.ReleaseParameter(4);
623 rpolz0 = fitterT2.GetParameter(3);
624 rpolz1 = fitterT2.GetParameter(4);
625 }
626
627 Double_t chi2TR = fitterT2.GetChisquare() / Float_t(npointsT);
628 Double_t params[3];
629 params[0] = fitterT2.GetParameter(0);
630 params[1] = fitterT2.GetParameter(1);
631 params[2] = fitterT2.GetParameter(2);
632 Double_t curvature = 1.0 + params[1] * params[1] - params[2] * params[0];
633
634
635 for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
636 Double_t x = cseed[iLayer].fX0;
637 Double_t y = 0;
638 Double_t dy = 0;
639 Double_t z = 0;
640 Double_t dz = 0;
641
642 // y
643 Double_t res2 = (x * params[0] + params[1]);
644 res2 *= res2;
645 res2 = 1.0 - params[2]*params[0] + params[1]*params[1] - res2;
646 if (res2 >= 0) {
647 res2 = TMath::Sqrt(res2);
648 y = (1.0 - res2) / params[0];
649 }
650
651 //dy
652 Double_t x0 = -params[1] / params[0];
653 if (-params[2]*params[0] + params[1]*params[1] + 1 > 0) {
654 Double_t rm1 = params[0] / TMath::Sqrt(-params[2]*params[0] + params[1]*params[1] + 1);
655 if (1.0/(rm1*rm1) - (x-x0) * (x-x0) > 0.0) {
656 Double_t res = (x - x0) / TMath::Sqrt((1./rm1-(x-x0))*(1./rm1+(x-x0)));
657 if (params[0] < 0) res *= -1.0;
658 dy = res;
659 }
660 }
661 z = rpolz0 + rpolz1 * (x - xref2);
662 dz = rpolz1;
663 cseed[iLayer].fYref[0] = y;
664 cseed[iLayer].fYref[1] = dy;
665 cseed[iLayer].fZref[0] = z;
666 cseed[iLayer].fZref[1] = dz;
667 cseed[iLayer].fC = curvature;
668
669 }
670
671 return chi2TR;
672
673}