Fixed coding violations
[u/mrichter/AliRoot.git] / TRD / AliTRDseedV1.cxx
CommitLineData
e4f2f73d 1/**************************************************************************
29b87567 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**************************************************************************/
e4f2f73d 15
16/* $Id$ */
17
18////////////////////////////////////////////////////////////////////////////
19// //
20// The TRD track seed //
21// //
22// Authors: //
23// Alex Bercuci <A.Bercuci@gsi.de> //
24// Markus Fasel <M.Fasel@gsi.de> //
25// //
26////////////////////////////////////////////////////////////////////////////
27
28#include "TMath.h"
29#include "TLinearFitter.h"
eb38ed55 30#include "TClonesArray.h" // tmp
31#include <TTreeStream.h>
e4f2f73d 32
33#include "AliLog.h"
34#include "AliMathBase.h"
d937ad7a 35#include "AliCDBManager.h"
36#include "AliTracker.h"
e4f2f73d 37
03cef9b2 38#include "AliTRDpadPlane.h"
e4f2f73d 39#include "AliTRDcluster.h"
f3d3af1b 40#include "AliTRDseedV1.h"
41#include "AliTRDtrackV1.h"
e4f2f73d 42#include "AliTRDcalibDB.h"
eb38ed55 43#include "AliTRDchamberTimeBin.h"
44#include "AliTRDtrackingChamber.h"
45#include "AliTRDtrackerV1.h"
46#include "AliTRDReconstructor.h"
e4f2f73d 47#include "AliTRDrecoParam.h"
d937ad7a 48
0906e73e 49#include "Cal/AliTRDCalPID.h"
d937ad7a 50#include "Cal/AliTRDCalROC.h"
51#include "Cal/AliTRDCalDet.h"
e4f2f73d 52
e4f2f73d 53ClassImp(AliTRDseedV1)
54
55//____________________________________________________________________
ae4e8b84 56AliTRDseedV1::AliTRDseedV1(Int_t det)
e4f2f73d 57 :AliTRDseed()
3a039a31 58 ,fReconstructor(0x0)
ae4e8b84 59 ,fClusterIter(0x0)
60 ,fClusterIdx(0)
61 ,fDet(det)
0906e73e 62 ,fMom(0.)
bcb6fb78 63 ,fSnp(0.)
64 ,fTgl(0.)
65 ,fdX(0.)
6e4d4425 66 ,fXref(0.)
d937ad7a 67 ,fExB(0.)
e4f2f73d 68{
69 //
70 // Constructor
71 //
29b87567 72 for(int islice=0; islice < knSlices; islice++) fdEdx[islice] = 0.;
73 for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) fProb[ispec] = -1.;
6e4d4425 74 fRefCov[0] = 1.; fRefCov[1] = 0.; fRefCov[2] = 1.;
d937ad7a 75 // covariance matrix [diagonal]
76 // default sy = 200um and sz = 2.3 cm
77 fCov[0] = 4.e-4; fCov[1] = 0.; fCov[2] = 5.3;
e4f2f73d 78}
79
80//____________________________________________________________________
0906e73e 81AliTRDseedV1::AliTRDseedV1(const AliTRDseedV1 &ref)
e4f2f73d 82 :AliTRDseed((AliTRDseed&)ref)
43d6ad34 83 ,fReconstructor(ref.fReconstructor)
ae4e8b84 84 ,fClusterIter(0x0)
85 ,fClusterIdx(0)
86 ,fDet(ref.fDet)
0906e73e 87 ,fMom(ref.fMom)
bcb6fb78 88 ,fSnp(ref.fSnp)
89 ,fTgl(ref.fTgl)
90 ,fdX(ref.fdX)
6e4d4425 91 ,fXref(ref.fXref)
d937ad7a 92 ,fExB(ref.fExB)
e4f2f73d 93{
94 //
95 // Copy Constructor performing a deep copy
96 //
97
29b87567 98 //printf("AliTRDseedV1::AliTRDseedV1(const AliTRDseedV1 &)\n");
99 SetBit(kOwner, kFALSE);
100 for(int islice=0; islice < knSlices; islice++) fdEdx[islice] = ref.fdEdx[islice];
101 for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) fProb[ispec] = ref.fProb[ispec];
6e4d4425 102 memcpy(fRefCov, ref.fRefCov, 3*sizeof(Double_t));
d937ad7a 103 memcpy(fCov, ref.fCov, 3*sizeof(Double_t));
fbb2ea06 104}
d9950a5a 105
0906e73e 106
e4f2f73d 107//____________________________________________________________________
108AliTRDseedV1& AliTRDseedV1::operator=(const AliTRDseedV1 &ref)
109{
110 //
111 // Assignment Operator using the copy function
112 //
113
29b87567 114 if(this != &ref){
115 ref.Copy(*this);
116 }
221ab7e0 117 SetBit(kOwner, kFALSE);
118
29b87567 119 return *this;
e4f2f73d 120
121}
122
123//____________________________________________________________________
124AliTRDseedV1::~AliTRDseedV1()
125{
126 //
127 // Destructor. The RecoParam object belongs to the underlying tracker.
128 //
129
29b87567 130 //printf("I-AliTRDseedV1::~AliTRDseedV1() : Owner[%s]\n", IsOwner()?"YES":"NO");
e4f2f73d 131
29b87567 132 if(IsOwner())
133 for(int itb=0; itb<knTimebins; itb++){
134 if(!fClusters[itb]) continue;
135 //AliInfo(Form("deleting c %p @ %d", fClusters[itb], itb));
136 delete fClusters[itb];
137 fClusters[itb] = 0x0;
138 }
e4f2f73d 139}
140
141//____________________________________________________________________
142void AliTRDseedV1::Copy(TObject &ref) const
143{
144 //
145 // Copy function
146 //
147
29b87567 148 //AliInfo("");
149 AliTRDseedV1 &target = (AliTRDseedV1 &)ref;
150
ae4e8b84 151 target.fClusterIter = 0x0;
152 target.fClusterIdx = 0;
153 target.fDet = fDet;
29b87567 154 target.fMom = fMom;
155 target.fSnp = fSnp;
156 target.fTgl = fTgl;
157 target.fdX = fdX;
6e4d4425 158 target.fXref = fXref;
d937ad7a 159 target.fExB = fExB;
29b87567 160 target.fReconstructor = fReconstructor;
161
162 for(int islice=0; islice < knSlices; islice++) target.fdEdx[islice] = fdEdx[islice];
163 for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) target.fProb[ispec] = fProb[ispec];
6e4d4425 164 memcpy(target.fRefCov, fRefCov, 3*sizeof(Double_t));
d937ad7a 165 memcpy(target.fCov, fCov, 3*sizeof(Double_t));
29b87567 166
167 AliTRDseed::Copy(target);
e4f2f73d 168}
169
0906e73e 170
171//____________________________________________________________
f3d3af1b 172Bool_t AliTRDseedV1::Init(AliTRDtrackV1 *track)
0906e73e 173{
174// Initialize this tracklet using the track information
175//
176// Parameters:
177// track - the TRD track used to initialize the tracklet
178//
179// Detailed description
180// The function sets the starting point and direction of the
181// tracklet according to the information from the TRD track.
182//
183// Caution
184// The TRD track has to be propagated to the beginning of the
185// chamber where the tracklet will be constructed
186//
187
29b87567 188 Double_t y, z;
189 if(!track->GetProlongation(fX0, y, z)) return kFALSE;
b1957d3c 190 UpDate(track);
29b87567 191 return kTRUE;
0906e73e 192}
193
bcb6fb78 194
195//____________________________________________________________________
b1957d3c 196void AliTRDseedV1::UpDate(const AliTRDtrackV1 *trk)
197{
198 // update tracklet reference position from the TRD track
199 // Funny name to avoid the clash with the function AliTRDseed::Update() (has to be made obselete)
200
201 fSnp = trk->GetSnp();
202 fTgl = trk->GetTgl();
203 fMom = trk->GetP();
204 fYref[1] = fSnp/(1. - fSnp*fSnp);
205 fZref[1] = fTgl;
206 SetCovRef(trk->GetCovariance());
207
208 Double_t dx = trk->GetX() - fX0;
209 fYref[0] = trk->GetY() - dx*fYref[1];
210 fZref[0] = trk->GetZ() - dx*fZref[1];
211}
212
213//____________________________________________________________________
bcb6fb78 214void AliTRDseedV1::CookdEdx(Int_t nslices)
215{
216// Calculates average dE/dx for all slices and store them in the internal array fdEdx.
217//
218// Parameters:
219// nslices : number of slices for which dE/dx should be calculated
220// Output:
221// store results in the internal array fdEdx. This can be accessed with the method
222// AliTRDseedV1::GetdEdx()
223//
224// Detailed description
225// Calculates average dE/dx for all slices. Depending on the PID methode
226// the number of slices can be 3 (LQ) or 8(NN).
3ee48d6e 227// The calculation of dQ/dl are done using the tracklet fit results (see AliTRDseedV1::GetdQdl(Int_t))
bcb6fb78 228//
229// The following effects are included in the calculation:
230// 1. calibration values for t0 and vdrift (using x coordinate to calculate slice)
231// 2. cluster sharing (optional see AliTRDrecoParam::SetClusterSharing())
232// 3. cluster size
233//
234
29b87567 235 Int_t nclusters[knSlices];
236 for(int i=0; i<knSlices; i++){
237 fdEdx[i] = 0.;
238 nclusters[i] = 0;
239 }
3ee48d6e 240 Float_t pathLength = (.5 * AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
29b87567 241
3ee48d6e 242 AliTRDcluster *c = 0x0;
29b87567 243 for(int ic=0; ic<AliTRDtrackerV1::GetNTimeBins(); ic++){
3ee48d6e 244 if(!(c = fClusters[ic]) && !(c = fClusters[knTimebins+ic])) continue;
245 Float_t x = c->GetX();
29b87567 246
247 // Filter clusters for dE/dx calculation
248
249 // 1.consider calibration effects for slice determination
250 Int_t slice;
3ee48d6e 251 if(c->IsInChamber()) slice = Int_t(TMath::Abs(fX0 - x) * nslices / pathLength);
29b87567 252 else slice = x < fX0 ? 0 : nslices-1;
253
254 // 2. take sharing into account
3ee48d6e 255 Float_t w = c->IsShared() ? .5 : 1.;
29b87567 256
257 // 3. take into account large clusters TODO
258 //w *= c->GetNPads() > 3 ? .8 : 1.;
259
260 //CHECK !!!
261 fdEdx[slice] += w * GetdQdl(ic); //fdQdl[ic];
262 nclusters[slice]++;
263 } // End of loop over clusters
264
cd40b287 265 //if(fReconstructor->GetPIDMethod() == AliTRDReconstructor::kLQPID){
266 if(nslices == AliTRDReconstructor::kLQslices){
29b87567 267 // calculate mean charge per slice (only LQ PID)
268 for(int is=0; is<nslices; is++){
269 if(nclusters[is]) fdEdx[is] /= nclusters[is];
270 }
271 }
bcb6fb78 272}
273
d937ad7a 274//____________________________________________________________________
275void AliTRDseedV1::GetClusterXY(const AliTRDcluster *c, Double_t &x, Double_t &y)
276{
277// Return corrected position of the cluster taking into
278// account variation of the drift velocity with drift length.
279
280
281 // drift velocity correction TODO to be moved to the clusterizer
282 const Float_t cx[] = {
283 -9.6280e-02, 1.3091e-01,-1.7415e-02,-9.9221e-02,-1.2040e-01,-9.5493e-02,
284 -5.0041e-02,-1.6726e-02, 3.5756e-03, 1.8611e-02, 2.6378e-02, 3.3823e-02,
285 3.4811e-02, 3.5282e-02, 3.5386e-02, 3.6047e-02, 3.5201e-02, 3.4384e-02,
286 3.2864e-02, 3.1932e-02, 3.2051e-02, 2.2539e-02,-2.5154e-02,-1.2050e-01,
287 -1.2050e-01
288 };
289
290 // PRF correction TODO to be replaced by the gaussian
291 // approximation with full error parametrization and // moved to the clusterizer
292 const Float_t cy[AliTRDgeometry::kNlayer][3] = {
293 { 4.014e-04, 8.605e-03, -6.880e+00},
294 {-3.061e-04, 9.663e-03, -6.789e+00},
295 { 1.124e-03, 1.105e-02, -6.825e+00},
296 {-1.527e-03, 1.231e-02, -6.777e+00},
297 { 2.150e-03, 1.387e-02, -6.783e+00},
298 {-1.296e-03, 1.486e-02, -6.825e+00}
299 };
300
301 Int_t ily = AliTRDgeometry::GetLayer(c->GetDetector());
302 x = c->GetX() - cx[c->GetLocalTimeBin()];
303 y = c->GetY() + cy[ily][0] + cy[ily][1] * TMath::Sin(cy[ily][2] * c->GetCenter());
304 return;
305}
b83573da 306
bcb6fb78 307//____________________________________________________________________
308Float_t AliTRDseedV1::GetdQdl(Int_t ic) const
309{
3ee48d6e 310// Using the linear approximation of the track inside one TRD chamber (TRD tracklet)
311// the charge per unit length can be written as:
312// BEGIN_LATEX
313// #frac{dq}{dl} = #frac{q_{c}}{dx * #sqrt{1 + #(){#frac{dy}{dx}}^{2}_{fit} + #(){#frac{dy}{dx}}^{2}_{ref}}}
314// END_LATEX
315// where qc is the total charge collected in the current time bin and dx is the length
316// of the time bin. For the moment (Jan 20 2009) only pad row cross corrections are
317// considered for the charge but none are applied for drift velocity variations along
318// the drift region or assymetry of the TRF
319//
320// Author : Alex Bercuci <A.Bercuci@gsi.de>
321//
322 Float_t dq = 0.;
323 if(fClusters[ic]) dq += TMath::Abs(fClusters[ic]->GetQ());
324 if(fClusters[knTimebins+ic]) dq += TMath::Abs(fClusters[knTimebins+ic]->GetQ());
325
326 return dq/fdX/TMath::Sqrt(1. + fYfit[1]*fYfit[1] + fZref[1]*fZref[1]);
bcb6fb78 327}
328
0906e73e 329//____________________________________________________________________
330Double_t* AliTRDseedV1::GetProbability()
331{
332// Fill probability array for tracklet from the DB.
333//
334// Parameters
335//
336// Output
337// returns pointer to the probability array and 0x0 if missing DB access
338//
339// Detailed description
340
29b87567 341
342 // retrive calibration db
0906e73e 343 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
344 if (!calibration) {
345 AliError("No access to calibration data");
346 return 0x0;
347 }
348
3a039a31 349 if (!fReconstructor) {
350 AliError("Reconstructor not set.");
4ba1d6ae 351 return 0x0;
352 }
353
0906e73e 354 // Retrieve the CDB container class with the parametric detector response
3a039a31 355 const AliTRDCalPID *pd = calibration->GetPIDObject(fReconstructor->GetPIDMethod());
0906e73e 356 if (!pd) {
357 AliError("No access to AliTRDCalPID object");
358 return 0x0;
359 }
29b87567 360 //AliInfo(Form("Method[%d] : %s", fReconstructor->GetRecoParam() ->GetPIDMethod(), pd->IsA()->GetName()));
10f75631 361
29b87567 362 // calculate tracklet length TO DO
0906e73e 363 Float_t length = (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
364 /// TMath::Sqrt((1.0 - fSnp[iPlane]*fSnp[iPlane]) / (1.0 + fTgl[iPlane]*fTgl[iPlane]));
365
366 //calculate dE/dx
3a039a31 367 CookdEdx(fReconstructor->GetNdEdxSlices());
0906e73e 368
369 // Sets the a priori probabilities
370 for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) {
ae4e8b84 371 fProb[ispec] = pd->GetProbability(ispec, fMom, &fdEdx[0], length, GetPlane());
0906e73e 372 }
373
29b87567 374 return &fProb[0];
0906e73e 375}
376
377//____________________________________________________________________
e4f2f73d 378Float_t AliTRDseedV1::GetQuality(Bool_t kZcorr) const
379{
380 //
381 // Returns a quality measurement of the current seed
382 //
383
29b87567 384 Float_t zcorr = kZcorr ? fTilt * (fZProb - fZref[0]) : 0.;
385 return
386 .5 * TMath::Abs(18.0 - fN2)
387 + 10.* TMath::Abs(fYfit[1] - fYref[1])
388 + 5. * TMath::Abs(fYfit[0] - fYref[0] + zcorr)
389 + 2. * TMath::Abs(fMeanz - fZref[0]) / fPadLength;
e4f2f73d 390}
391
392//____________________________________________________________________
d937ad7a 393void AliTRDseedV1::GetCovAt(Double_t x, Double_t *cov) const
0906e73e 394{
d937ad7a 395// Computes covariance in the y-z plane at radial point x (in tracking coordinates)
396// and returns the results in the preallocated array cov[3] as :
397// cov[0] = Var(y)
398// cov[1] = Cov(yz)
399// cov[2] = Var(z)
400//
401// Details
402//
403// For the linear transformation
404// BEGIN_LATEX
405// Y = T_{x} X^{T}
406// END_LATEX
407// The error propagation has the general form
408// BEGIN_LATEX
409// C_{Y} = T_{x} C_{X} T_{x}^{T}
410// END_LATEX
411// We apply this formula 2 times. First to calculate the covariance of the tracklet
412// at point x we consider:
413// BEGIN_LATEX
414// T_{x} = (1 x); X=(y0 dy/dx); C_{X}=#(){#splitline{Var(y0) Cov(y0, dy/dx)}{Cov(y0, dy/dx) Var(dy/dx)}}
415// END_LATEX
416// and secondly to take into account the tilt angle
417// BEGIN_LATEX
418// T_{#alpha} = #(){#splitline{cos(#alpha) __ sin(#alpha)}{-sin(#alpha) __ cos(#alpha)}}; X=(y z); C_{X}=#(){#splitline{Var(y) 0}{0 Var(z)}}
419// END_LATEX
420//
421// using simple trigonometrics one can write for this last case
422// BEGIN_LATEX
423// C_{Y}=#frac{1}{1+tg^{2}#alpha} #(){#splitline{(#sigma_{y}^{2}+tg^{2}#alpha#sigma_{z}^{2}) __ tg#alpha(#sigma_{z}^{2}-#sigma_{y}^{2})}{tg#alpha(#sigma_{z}^{2}-#sigma_{y}^{2}) __ (#sigma_{z}^{2}+tg^{2}#alpha#sigma_{y}^{2})}}
424// END_LATEX
425// which can be aproximated for small alphas (2 deg) with
426// BEGIN_LATEX
427// C_{Y}=#(){#splitline{#sigma_{y}^{2} __ (#sigma_{z}^{2}-#sigma_{y}^{2})tg#alpha}{((#sigma_{z}^{2}-#sigma_{y}^{2})tg#alpha __ #sigma_{z}^{2}}}
428// END_LATEX
429//
430// before applying the tilt rotation we also apply systematic uncertainties to the tracklet
431// position which can be tunned from outside via the AliTRDrecoParam::SetSysCovMatrix(). They might
432// account for extra misalignment/miscalibration uncertainties.
433//
434// Author :
435// Alex Bercuci <A.Bercuci@gsi.de>
436// Date : Jan 8th 2009
437//
b1957d3c 438
439
d937ad7a 440 Double_t xr = fX0-x;
441 Double_t sy2 = fCov[0] +2.*xr*fCov[1] + xr*xr*fCov[2];
442 Double_t sz2 = fPadLength*fPadLength/12.;
0906e73e 443
d937ad7a 444 // insert systematic uncertainties
445 Double_t sys[15];
446 fReconstructor->GetRecoParam()->GetSysCovMatrix(sys);
447 sy2 += sys[0];
448 sz2 += sys[1];
449
450 // rotate covariance matrix
451 Double_t t2 = fTilt*fTilt;
452 Double_t correction = 1./(1. + t2);
453 cov[0] = (sy2+t2*sz2)*correction;
454 cov[1] = fTilt*(sz2 - sy2)*correction;
455 cov[2] = (t2*sy2+sz2)*correction;
456}
eb38ed55 457
0906e73e 458
d937ad7a 459//____________________________________________________________________
460void AliTRDseedV1::SetExB()
461{
462// Retrive the tg(a_L) from OCDB. The following information are used
463// - detector index
464// - column and row position of first attached cluster.
465//
466// If no clusters are attached to the tracklet a random central chamber
467// position (c=70, r=7) will be used to retrieve the drift velocity.
468//
469// Author :
470// Alex Bercuci <A.Bercuci@gsi.de>
471// Date : Jan 8th 2009
472//
eb38ed55 473
d937ad7a 474 AliCDBManager *cdb = AliCDBManager::Instance();
475 if(cdb->GetRun() < 0){
476 AliError("OCDB manager not properly initialized");
477 return;
478 }
0906e73e 479
d937ad7a 480 AliTRDcalibDB *fCalib = AliTRDcalibDB::Instance();
481 AliTRDCalROC *fCalVdriftROC = fCalib->GetVdriftROC(fDet);
482 const AliTRDCalDet *fCalVdriftDet = fCalib->GetVdriftDet();
483
484 Int_t col = 70, row = 7;
485 AliTRDcluster **c = &fClusters[0];
486 if(fN){
487 Int_t ic = 0;
488 while (ic<AliTRDseed::knTimebins && !(*c)){ic++; c++;}
489 if(*c){
490 col = (*c)->GetPadCol();
491 row = (*c)->GetPadRow();
492 }
493 }
3a039a31 494
d937ad7a 495 Double_t vd = fCalVdriftDet->GetValue(fDet) * fCalVdriftROC->GetValue(col, row);
496 fExB = fCalib->GetOmegaTau(vd, -0.1*AliTracker::GetBz());
0906e73e 497}
498
0906e73e 499//____________________________________________________________________
29b87567 500void AliTRDseedV1::SetOwner()
0906e73e 501{
29b87567 502 //AliInfo(Form("own [%s] fOwner[%s]", own?"YES":"NO", fOwner?"YES":"NO"));
503
504 if(TestBit(kOwner)) return;
505 for(int ic=0; ic<knTimebins; ic++){
506 if(!fClusters[ic]) continue;
507 fClusters[ic] = new AliTRDcluster(*fClusters[ic]);
508 }
509 SetBit(kOwner);
0906e73e 510}
511
512//____________________________________________________________________
eb38ed55 513Bool_t AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t quality, Bool_t kZcorr, AliTRDcluster *c)
e4f2f73d 514{
515 //
516 // Iterative process to register clusters to the seed.
517 // In iteration 0 we try only one pad-row and if quality not
518 // sufficient we try 2 pad-rows (about 5% of tracks cross 2 pad-rows)
519 //
29b87567 520 // debug level 7
521 //
522
523 if(!fReconstructor->GetRecoParam() ){
524 AliError("Seed can not be used without a valid RecoParam.");
525 return kFALSE;
526 }
527
528 AliTRDchamberTimeBin *layer = 0x0;
a8276d32 529 if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7){
e8037fda 530 AliTRDtrackingChamber ch(*chamber);
531 ch.SetOwner();
29f95561 532 TTreeSRedirector &cstreamer = *fReconstructor->GetDebugStream(AliTRDReconstructor::kTracker);
533 cstreamer << "AttachClustersIter"
e8037fda 534 << "chamber.=" << &ch
29b87567 535 << "tracklet.=" << this
29b87567 536 << "\n";
537 }
538
35c24814 539 Float_t tquality;
29b87567 540 Double_t kroady = fReconstructor->GetRecoParam() ->GetRoad1y();
541 Double_t kroadz = fPadLength * .5 + 1.;
35c24814 542
543 // initialize configuration parameters
544 Float_t zcorr = kZcorr ? fTilt * (fZProb - fZref[0]) : 0.;
545 Int_t niter = kZcorr ? 1 : 2;
546
29b87567 547 Double_t yexp, zexp;
548 Int_t ncl = 0;
35c24814 549 // start seed update
550 for (Int_t iter = 0; iter < niter; iter++) {
29b87567 551 ncl = 0;
552 for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
553 if(!(layer = chamber->GetTB(iTime))) continue;
554 if(!Int_t(*layer)) continue;
555
556 // define searching configuration
557 Double_t dxlayer = layer->GetX() - fX0;
558 if(c){
559 zexp = c->GetZ();
560 //Try 2 pad-rows in second iteration
561 if (iter > 0) {
562 zexp = fZref[0] + fZref[1] * dxlayer - zcorr;
563 if (zexp > c->GetZ()) zexp = c->GetZ() + fPadLength*0.5;
564 if (zexp < c->GetZ()) zexp = c->GetZ() - fPadLength*0.5;
565 }
566 } else zexp = fZref[0] + (kZcorr ? fZref[1] * dxlayer : 0.);
35c24814 567 yexp = fYref[0] + fYref[1] * dxlayer - zcorr;
29b87567 568
569 // Get and register cluster
570 Int_t index = layer->SearchNearestCluster(yexp, zexp, kroady, kroadz);
571 if (index < 0) continue;
572 AliTRDcluster *cl = (*layer)[index];
35c24814 573
29b87567 574 fIndexes[iTime] = layer->GetGlobalIndex(index);
575 fClusters[iTime] = cl;
576 fY[iTime] = cl->GetY();
577 fZ[iTime] = cl->GetZ();
578 ncl++;
579 }
35c24814 580 if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7) AliInfo(Form("iter = %d ncl [%d] = %d", iter, fDet, ncl));
29b87567 581
582 if(ncl>1){
583 // calculate length of the time bin (calibration aware)
584 Int_t irp = 0; Float_t x[2]; Int_t tb[2];
585 for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
586 if(!fClusters[iTime]) continue;
587 x[irp] = fClusters[iTime]->GetX();
588 tb[irp] = iTime;
589 irp++;
590 if(irp==2) break;
591 }
592 fdX = (x[1] - x[0]) / (tb[0] - tb[1]);
593
594 // update X0 from the clusters (calibration/alignment aware)
595 for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
596 if(!(layer = chamber->GetTB(iTime))) continue;
597 if(!layer->IsT0()) continue;
598 if(fClusters[iTime]){
599 fX0 = fClusters[iTime]->GetX();
600 break;
601 } else { // we have to infere the position of the anode wire from the other clusters
602 for (Int_t jTime = iTime+1; jTime < AliTRDtrackerV1::GetNTimeBins(); jTime++) {
603 if(!fClusters[jTime]) continue;
604 fX0 = fClusters[jTime]->GetX() + fdX * (jTime - iTime);
f660dce9 605 break;
29b87567 606 }
29b87567 607 }
608 }
609
610 // update YZ reference point
611 // TODO
612
613 // update x reference positions (calibration/alignment aware)
614 for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
615 if(!fClusters[iTime]) continue;
0849f128 616 fX[iTime] = fX0 - fClusters[iTime]->GetX();
29b87567 617 }
618
619 AliTRDseed::Update();
620 }
35c24814 621 if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7) AliInfo(Form("iter = %d nclFit [%d] = %d", iter, fDet, fN2));
29b87567 622
623 if(IsOK()){
624 tquality = GetQuality(kZcorr);
625 if(tquality < quality) break;
626 else quality = tquality;
627 }
628 kroadz *= 2.;
629 } // Loop: iter
630 if (!IsOK()) return kFALSE;
631
804bb02e 632 if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=1) CookLabels();
d937ad7a 633
b1957d3c 634 // set ExB angle
635 SetExB();
29b87567 636 UpdateUsed();
637 return kTRUE;
e4f2f73d 638}
639
640//____________________________________________________________________
b1957d3c 641Bool_t AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *chamber, Bool_t tilt)
e4f2f73d 642{
643 //
644 // Projective algorithm to attach clusters to seeding tracklets
645 //
646 // Parameters
647 //
648 // Output
649 //
650 // Detailed description
651 // 1. Collapse x coordinate for the full detector plane
652 // 2. truncated mean on y (r-phi) direction
653 // 3. purge clusters
654 // 4. truncated mean on z direction
655 // 5. purge clusters
656 // 6. fit tracklet
657 //
b1957d3c 658 Bool_t kPRINT = kFALSE;
29b87567 659 if(!fReconstructor->GetRecoParam() ){
660 AliError("Seed can not be used without a valid RecoParam.");
661 return kFALSE;
662 }
b1957d3c 663 // Initialize reco params for this tracklet
664 // 1. first time bin in the drift region
665 Int_t t0 = 4;
666 Int_t kClmin = Int_t(fReconstructor->GetRecoParam() ->GetFindableClusters()*AliTRDtrackerV1::GetNTimeBins());
29b87567 667
b1957d3c 668 Double_t syRef = TMath::Sqrt(fRefCov[0]);
29b87567 669 //define roads
b1957d3c 670 Double_t kroady = 1.;
671 //fReconstructor->GetRecoParam() ->GetRoad1y();
29b87567 672 Double_t kroadz = fPadLength * 1.5 + 1.;
b1957d3c 673 if(kPRINT) printf("AttachClusters() sy[%f] road[%f]\n", syRef, kroady);
29b87567 674
675 // working variables
b1957d3c 676 const Int_t kNrows = 16;
677 AliTRDcluster *clst[kNrows][knTimebins];
678 Double_t cond[4], dx, dy, yt, zt,
679 yres[kNrows][knTimebins];
680 Int_t idxs[kNrows][knTimebins], ncl[kNrows], ncls = 0;
681 memset(ncl, 0, kNrows*sizeof(Int_t));
682 memset(clst, 0, kNrows*knTimebins*sizeof(AliTRDcluster*));
683
29b87567 684 // Do cluster projection
b1957d3c 685 AliTRDcluster *c = 0x0;
29b87567 686 AliTRDchamberTimeBin *layer = 0x0;
b1957d3c 687 Bool_t kBUFFER = kFALSE;
688 for (Int_t it = 0; it < AliTRDtrackerV1::GetNTimeBins(); it++) {
689 if(!(layer = chamber->GetTB(it))) continue;
29b87567 690 if(!Int_t(*layer)) continue;
691
b1957d3c 692 dx = fX0 - layer->GetX();
693 yt = fYref[0] - fYref[1] * dx;
694 zt = fZref[0] - fZref[1] * dx;
695 if(kPRINT) printf("\t%2d dx[%f] yt[%f] zt[%f]\n", it, dx, yt, zt);
696
697 // select clusters on a 5 sigmaKalman level
698 cond[0] = yt; cond[2] = kroady;
699 cond[1] = zt; cond[3] = kroadz;
700 Int_t n=0, idx[6];
701 layer->GetClusters(cond, idx, n, 6);
702 for(Int_t ic = n; ic--;){
703 c = (*layer)[idx[ic]];
704 dy = yt - c->GetY();
705 dy += tilt ? fTilt * (c->GetZ() - zt) : 0.;
706 // select clusters on a 3 sigmaKalman level
707/* if(tilt && TMath::Abs(dy) > 3.*syRef){
708 printf("too large !!!\n");
709 continue;
710 }*/
711 Int_t r = c->GetPadRow();
712 if(kPRINT) printf("\t\t%d dy[%f] yc[%f] r[%d]\n", ic, TMath::Abs(dy), c->GetY(), r);
713 clst[r][ncl[r]] = c;
714 idxs[r][ncl[r]] = idx[ic];
715 yres[r][ncl[r]] = dy;
716 ncl[r]++; ncls++;
717
718 if(ncl[r] >= knTimebins) {
719 AliWarning(Form("Cluster candidates reached limit %d. Some may be lost.", knTimebins));
720 kBUFFER = kTRUE;
29b87567 721 break;
722 }
723 }
b1957d3c 724 if(kBUFFER) break;
29b87567 725 }
b1957d3c 726 if(kPRINT) printf("Found %d clusters\n", ncls);
727 if(ncls<kClmin) return kFALSE;
728
729 // analyze each row individualy
730 Double_t mean, syDis;
731 Int_t nrow[] = {0, 0, 0}, nr = 0, lr=-1;
732 for(Int_t ir=kNrows; ir--;){
733 if(!(ncl[ir])) continue;
734 if(lr>0 && lr-ir != 1){
735 if(kPRINT) printf("W - gap in rows attached !!\n");
29b87567 736 }
b1957d3c 737 if(kPRINT) printf("\tir[%d] lr[%d] n[%d]\n", ir, lr, ncl[ir]);
738 // Evaluate truncated mean on the y direction
739 if(ncl[ir] > 3) AliMathBase::EvaluateUni(ncl[ir], yres[ir], mean, syDis, Int_t(ncl[ir]*.8));
740 else {
741 mean = 0.; syDis = 0.;
742 }
743
744 // TODO check mean and sigma agains cluster resolution !!
745 if(kPRINT) printf("\tr[%2d] m[%f %5.3fsigma] s[%f]\n", ir, mean, TMath::Abs(mean/syRef), syDis);
746 // select clusters on a 3 sigmaDistr level
747 Bool_t kFOUND = kFALSE;
748 for(Int_t ic = ncl[ir]; ic--;){
749 if(yres[ir][ic] - mean > 3. * syDis){
750 clst[ir][ic] = 0x0; continue;
751 }
752 nrow[nr]++; kFOUND = kTRUE;
753 }
754 // exit loop
755 if(kFOUND) nr++;
756 lr = ir; if(nr>=3) break;
29b87567 757 }
b1957d3c 758 if(kPRINT) printf("lr[%d] nr[%d] nrow[0]=%d nrow[1]=%d nrow[2]=%d\n", lr, nr, nrow[0], nrow[1], nrow[2]);
759
760 // classify cluster rows
761 Int_t row = -1;
762 switch(nr){
763 case 1:
764 row = lr;
765 break;
766 case 2:
767 SetBit(kRowCross, kTRUE); // mark pad row crossing
768 if(nrow[0] > nrow[1]){ row = lr+1; lr = -1;}
769 else{
770 row = lr; lr = 1;
771 nrow[2] = nrow[1];
772 nrow[1] = nrow[0];
773 nrow[0] = nrow[2];
29b87567 774 }
b1957d3c 775 break;
776 case 3:
777 SetBit(kRowCross, kTRUE); // mark pad row crossing
778 break;
29b87567 779 }
b1957d3c 780 if(kPRINT) printf("\trow[%d] n[%d]\n\n", row, nrow[0]);
781 if(row<0) return kFALSE;
29b87567 782
b1957d3c 783 // Select and store clusters
784 // We should consider here :
785 // 1. How far is the chamber boundary
786 // 2. How big is the mean
29b87567 787 fN2 = 0;
b1957d3c 788 for (Int_t ir = 0; ir < nr; ir++) {
789 Int_t jr = row + ir*lr;
790 if(kPRINT) printf("\tattach %d clusters for row %d\n", ncl[jr], jr);
791 for (Int_t ic = 0; ic < ncl[jr]; ic++) {
792 if(!(c = clst[jr][ic])) continue;
793 Int_t it = c->GetPadTime();
794 // TODO proper indexing of clusters !!
795 fIndexes[it+35*ir] = chamber->GetTB(it)->GetGlobalIndex(idxs[jr][ic]);
796 fClusters[it+35*ir] = c;
29b87567 797
b1957d3c 798 //printf("\tid[%2d] it[%d] idx[%d]\n", ic, it, fIndexes[it]);
799
800 fN2++;
801 }
802 }
803
29b87567 804 // number of minimum numbers of clusters expected for the tracklet
e4f2f73d 805 if (fN2 < kClmin){
29b87567 806 AliWarning(Form("Not enough clusters to fit the tracklet %d [%d].", fN2, kClmin));
e4f2f73d 807 fN2 = 0;
808 return kFALSE;
809 }
0906e73e 810
b1957d3c 811 // update used clusters and select
29b87567 812 fNUsed = 0;
b1957d3c 813 for (Int_t it = 0; it < AliTRDtrackerV1::GetNTimeBins(); it++) {
814 if(fClusters[it] && fClusters[it]->IsUsed()) fNUsed++;
815 if(fClusters[it+35] && fClusters[it+35]->IsUsed()) fNUsed++;
29b87567 816 }
0906e73e 817 if (fN2-fNUsed < kClmin){
b1957d3c 818 //AliWarning(Form("Too many clusters already in use %d (from %d).", fNUsed, fN2));
0906e73e 819 fN2 = 0;
820 return kFALSE;
821 }
b1957d3c 822
823 // set the Lorentz angle for this tracklet
824 SetExB();
825
826 // calculate dx for time bins in the drift region (calibration aware)
827 Int_t irp = 0; Float_t x[2]; Int_t tb[2];
828 for (Int_t it = t0; it < AliTRDtrackerV1::GetNTimeBins(); it++) {
829 if(!fClusters[it]) continue;
830 x[irp] = fClusters[it]->GetX();
831 tb[irp] = it;
832 irp++;
833 if(irp==2) break;
834 }
835 fdX = (x[1] - x[0]) / (tb[0] - tb[1]);
836
837 // update X0 from the clusters (calibration/alignment aware) TODO remove dependence on x0 !!
838 for (Int_t it = 0; it < AliTRDtrackerV1::GetNTimeBins(); it++) {
839 if(!(layer = chamber->GetTB(it))) continue;
840 if(!layer->IsT0()) continue;
841 if(fClusters[it]){
842 fX0 = fClusters[it]->GetX();
843 break;
844 } else { // we have to infere the position of the anode wire from the other clusters
845 for (Int_t jt = it+1; jt < AliTRDtrackerV1::GetNTimeBins(); jt++) {
846 if(!fClusters[jt]) continue;
847 fX0 = fClusters[jt]->GetX() + fdX * (jt - it);
848 break;
849 }
850 }
851 }
852
29b87567 853 return kTRUE;
e4f2f73d 854}
855
03cef9b2 856//____________________________________________________________
857void AliTRDseedV1::Bootstrap(const AliTRDReconstructor *rec)
858{
859// Fill in all derived information. It has to be called after recovery from file or HLT.
860// The primitive data are
861// - list of clusters
862// - detector (as the detector will be removed from clusters)
863// - position of anode wire (fX0) - temporary
864// - track reference position and direction
865// - momentum of the track
866// - time bin length [cm]
867//
868// A.Bercuci <A.Bercuci@gsi.de> Oct 30th 2008
869//
870 fReconstructor = rec;
871 AliTRDgeometry g;
872 AliTRDpadPlane *pp = g.GetPadPlane(fDet);
873 fTilt = TMath::Tan(TMath::DegToRad()*pp->GetTiltingAngle());
874 fPadLength = pp->GetLengthIPad();
875 fSnp = fYref[1]/TMath::Sqrt(1+fYref[1]*fYref[1]);
876 fTgl = fZref[1];
877 fN = 0; fN2 = 0; fMPads = 0.;
878 AliTRDcluster **cit = &fClusters[0];
879 for(Int_t ic = knTimebins; ic--; cit++){
880 if(!(*cit)) return;
881 fN++; fN2++;
882 fX[ic] = (*cit)->GetX() - fX0;
883 fY[ic] = (*cit)->GetY();
884 fZ[ic] = (*cit)->GetZ();
885 }
886 Update(); // Fit();
887 CookLabels();
888 GetProbability();
889}
890
891
e4f2f73d 892//____________________________________________________________________
d937ad7a 893Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
e4f2f73d 894{
895 //
896 // Linear fit of the tracklet
897 //
898 // Parameters :
899 //
900 // Output :
901 // True if successful
902 //
903 // Detailed description
904 // 2. Check if tracklet crosses pad row boundary
905 // 1. Calculate residuals in the y (r-phi) direction
906 // 3. Do a Least Square Fit to the data
907 //
908
29b87567 909 const Int_t kClmin = 8;
010d62b0 910
9462866a 911
912 // cluster error parametrization parameters
010d62b0 913 // 1. sy total charge
9462866a 914 const Float_t sq0inv = 0.019962; // [1/q0]
915 const Float_t sqb = 1.0281564; //[cm]
010d62b0 916 // 2. sy for the PRF
917 const Float_t scy[AliTRDgeometry::kNlayer][4] = {
d937ad7a 918 {2.827e-02, 9.600e-04, 4.296e-01, 2.271e-02},
919 {2.952e-02,-2.198e-04, 4.146e-01, 2.339e-02},
920 {3.090e-02, 1.514e-03, 4.020e-01, 2.402e-02},
921 {3.260e-02,-2.037e-03, 3.946e-01, 2.509e-02},
922 {3.439e-02,-3.601e-04, 3.883e-01, 2.623e-02},
923 {3.510e-02, 2.066e-03, 3.651e-01, 2.588e-02},
010d62b0 924 };
925 // 3. sy parallel to the track
d937ad7a 926 const Float_t sy0 = 2.649e-02; // [cm]
927 const Float_t sya = -8.864e-04; // [cm]
928 const Float_t syb = -2.435e-01; // [cm]
929
010d62b0 930 // 4. sx parallel to the track
d937ad7a 931 const Float_t sxgc = 5.427e-02;
932 const Float_t sxgm = 7.783e-01;
933 const Float_t sxgs = 2.743e-01;
934 const Float_t sxe0 =-2.065e+00;
935 const Float_t sxe1 =-2.978e-02;
936
010d62b0 937 // 5. sx perpendicular to the track
d937ad7a 938// const Float_t sxd0 = 1.881e-02;
939// const Float_t sxd1 =-4.101e-01;
940// const Float_t sxd2 = 1.572e+00;
941
2f7d6ac8 942 // get track direction
943 Double_t y0 = fYref[0];
944 Double_t dydx = fYref[1];
945 Double_t z0 = fZref[0];
946 Double_t dzdx = fZref[1];
947 Double_t yt, zt;
ae4e8b84 948
29b87567 949 const Int_t kNtb = AliTRDtrackerV1::GetNTimeBins();
b1957d3c 950 //AliTRDtrackerV1::AliTRDLeastSquare fitterZ;
24d8660e 951 TLinearFitter fitterY(1, "pol1");
29b87567 952 // convertion factor from square to gauss distribution for sigma
b1957d3c 953 //Double_t convert = 1./TMath::Sqrt(12.);
ae4e8b84 954
29b87567 955 // book cluster information
b1957d3c 956 Double_t q, xc[knTimebins], yc[knTimebins], zc[knTimebins], sy[knTimebins]/*, sz[knTimebins]*/;
957// Int_t zRow[knTimebins];
9462866a 958
010d62b0 959 Int_t ily = AliTRDgeometry::GetLayer(fDet);
b1957d3c 960 fN = 0; //fXref = 0.; Double_t ssx = 0.;
9eb2d46c 961 AliTRDcluster *c=0x0, **jc = &fClusters[0];
9eb2d46c 962 for (Int_t ic=0; ic<kNtb; ic++, ++jc) {
b1957d3c 963 //zRow[ic] = -1;
29b87567 964 xc[ic] = -1.;
965 yc[ic] = 999.;
966 zc[ic] = 999.;
967 sy[ic] = 0.;
b1957d3c 968 //sz[ic] = 0.;
9eb2d46c 969 if(!(c = (*jc))) continue;
29b87567 970 if(!c->IsInChamber()) continue;
9462866a 971
29b87567 972 Float_t w = 1.;
973 if(c->GetNPads()>4) w = .5;
974 if(c->GetNPads()>5) w = .2;
010d62b0 975
b1957d3c 976 //zRow[fN] = c->GetPadRow();
d937ad7a 977 // correct cluster position for PRF and v drift
978 Double_t x, y; GetClusterXY(c, x, y);
979 xc[fN] = fX0 - x;
980 yc[fN] = y;
2f7d6ac8 981 zc[fN] = c->GetZ();
982
983 // extrapolated y value for the track
984 yt = y0 - xc[fN]*dydx;
985 // extrapolated z value for the track
986 zt = z0 - xc[fN]*dzdx;
987 // tilt correction
988 if(tilt) yc[fN] -= fTilt*(zc[fN] - zt);
989
010d62b0 990 // ELABORATE CLUSTER ERROR
991 // TODO to be moved to AliTRDcluster
9462866a 992 q = TMath::Abs(c->GetQ());
d937ad7a 993 Double_t tgg = (dydx-fExB)/(1.+dydx*fExB); tgg *= tgg;
010d62b0 994 // basic y error (|| to track).
d937ad7a 995 sy[fN] = xc[fN] < AliTRDgeometry::CamHght() ? 2. : sy0 + sya*TMath::Exp(1./(xc[fN]+syb));
996 //printf("cluster[%d]\n\tsy[0] = %5.3e [um]\n", fN, sy[fN]*1.e4);
010d62b0 997 // y error due to total charge
998 sy[fN] += sqb*(1./q - sq0inv);
d937ad7a 999 //printf("\tsy[1] = %5.3e [um]\n", sy[fN]*1.e4);
010d62b0 1000 // y error due to PRF
1001 sy[fN] += scy[ily][0]*TMath::Gaus(c->GetCenter(), scy[ily][1], scy[ily][2]) - scy[ily][3];
d937ad7a 1002 //printf("\tsy[2] = %5.3e [um]\n", sy[fN]*1.e4);
1003
010d62b0 1004 sy[fN] *= sy[fN];
1005
1006 // ADD ERROR ON x
9462866a 1007 // error of drift length parallel to the track
d937ad7a 1008 Double_t sx = sxgc*TMath::Gaus(xc[fN], sxgm, sxgs) + TMath::Exp(sxe0+sxe1*xc[fN]); // [cm]
1009 //printf("\tsx[0] = %5.3e [um]\n", sx*1.e4);
9462866a 1010 // error of drift length perpendicular to the track
1011 //sx += sxd0 + sxd1*d + sxd2*d*d;
d937ad7a 1012 sx *= sx; // square sx
1013 // update xref
b1957d3c 1014 //fXref += xc[fN]/sx; ssx+=1./sx;
d937ad7a 1015
9462866a 1016 // add error from ExB
d937ad7a 1017 if(errors>0) sy[fN] += fExB*fExB*sx;
1018 //printf("\tsy[3] = %5.3e [um^2]\n", sy[fN]*1.e8);
1019
1020 // global radial error due to misalignment/miscalibration
1021 Double_t sx0 = 0.; sx0 *= sx0;
1022 // add sx contribution to sy due to track angle
1023 if(errors>1) sy[fN] += tgg*(sx+sx0);
1024 // TODO we should add tilt pad correction here
1025 //printf("\tsy[4] = %5.3e [um^2]\n", sy[fN]*1.e8);
1026 c->SetSigmaY2(sy[fN]);
1027
9462866a 1028 sy[fN] = TMath::Sqrt(sy[fN]);
24d8660e 1029 fitterY.AddPoint(&xc[fN], yc[fN]/*-yt*/, sy[fN]);
2f7d6ac8 1030 fN++;
29b87567 1031 }
47d5d320 1032 // to few clusters
2f7d6ac8 1033 if (fN < kClmin) return kFALSE;
1034
d937ad7a 1035 // fit XY
2f7d6ac8 1036 fitterY.Eval();
d937ad7a 1037 fYfit[0] = fitterY.GetParameter(0);
1038 fYfit[1] = -fitterY.GetParameter(1);
1039 // store covariance
1040 Double_t *p = fitterY.GetCovarianceMatrix();
1041 fCov[0] = p[0]; // variance of y0
1042 fCov[1] = p[1]; // covariance of y0, dydx
1043 fCov[2] = p[3]; // variance of dydx
b1957d3c 1044 // the ref radial position is set at the minimum of
1045 // the y variance of the tracklet
1046 fXref = -fCov[1]/fCov[2]; //fXref = fX0 - fXref;
1047
1048 // fit XZ
1049 if(IsRowCross()){
1050 // TODO pad row cross position estimation !!!
1051 //AliInfo(Form("Padrow cross in detector %d", fDet));
1052 fZfit[0] = .5*(zc[0]+zc[fN-1]); fZfit[1] = 0.;
1053 } else {
1054 fZfit[0] = zc[0]; fZfit[1] = 0.;
29b87567 1055 }
1056
29b87567 1057
b1957d3c 1058// // determine z offset of the fit
1059// Float_t zslope = 0.;
1060// Int_t nchanges = 0, nCross = 0;
1061// if(nz==2){ // tracklet is crossing pad row
1062// // Find the break time allowing one chage on pad-rows
1063// // with maximal number of accepted clusters
1064// Int_t padRef = zRow[0];
1065// for (Int_t ic=1; ic<fN; ic++) {
1066// if(zRow[ic] == padRef) continue;
1067//
1068// // debug
1069// if(zRow[ic-1] == zRow[ic]){
1070// printf("ERROR in pad row change!!!\n");
1071// }
1072//
1073// // evaluate parameters of the crossing point
1074// Float_t sx = (xc[ic-1] - xc[ic])*convert;
1075// fCross[0] = .5 * (xc[ic-1] + xc[ic]);
1076// fCross[2] = .5 * (zc[ic-1] + zc[ic]);
1077// fCross[3] = TMath::Max(dzdx * sx, .01);
1078// zslope = zc[ic-1] > zc[ic] ? 1. : -1.;
1079// padRef = zRow[ic];
1080// nCross = ic;
1081// nchanges++;
1082// }
1083// }
1084//
1085// // condition on nCross and reset nchanges TODO
1086//
1087// if(nchanges==1){
1088// if(dzdx * zslope < 0.){
1089// AliInfo("Tracklet-Track mismatch in dzdx. TODO.");
1090// }
1091//
1092//
1093// //zc[nc] = fitterZ.GetFunctionParameter(0);
1094// fCross[1] = fYfit[0] - fCross[0] * fYfit[1];
1095// fCross[0] = fX0 - fCross[0];
1096// }
29b87567 1097
2389e96f 1098 UpdateUsed();
29b87567 1099 return kTRUE;
e4f2f73d 1100}
1101
e4f2f73d 1102
1103//___________________________________________________________________
203967fc 1104void AliTRDseedV1::Print(Option_t *o) const
e4f2f73d 1105{
1106 //
1107 // Printing the seedstatus
1108 //
1109
203967fc 1110 AliInfo(Form("Det[%3d] Tilt[%+6.2f] Pad[%5.2f]", fDet, fTilt, fPadLength));
1111 AliInfo(Form("Nattach[%2d] Nfit[%2d] Nuse[%2d] pads[%f]", fN, fN2, fNUsed, fMPads));
1112 AliInfo(Form("x[%7.2f] y[%7.2f] z[%7.2f] dydx[%5.2f] dzdx[%5.2f]", fX0, fYfit[0], fZfit[0], fYfit[1], fZfit[1]));
1113 AliInfo(Form("Ref y[%7.2f] z[%7.2f] dydx[%5.2f] dzdx[%5.2f]", fYref[0], fZref[0], fYref[1], fZref[1]))
1114
1115
1116 if(strcmp(o, "a")!=0) return;
1117
4dc4dc2e 1118 AliTRDcluster* const* jc = &fClusters[0];
1119 for(int ic=0; ic<AliTRDtrackerV1::GetNTimeBins(); ic++, jc++) {
1120 if(!(*jc)) continue;
203967fc 1121 (*jc)->Print(o);
4dc4dc2e 1122 }
e4f2f73d 1123}
47d5d320 1124
203967fc 1125
1126//___________________________________________________________________
1127Bool_t AliTRDseedV1::IsEqual(const TObject *o) const
1128{
1129 // Checks if current instance of the class has the same essential members
1130 // as the given one
1131
1132 if(!o) return kFALSE;
1133 const AliTRDseedV1 *inTracklet = dynamic_cast<const AliTRDseedV1*>(o);
1134 if(!inTracklet) return kFALSE;
1135
1136 for (Int_t i = 0; i < 2; i++){
1137 if ( fYref[i] != inTracklet->GetYref(i) ) return kFALSE;
1138 if ( fZref[i] != inTracklet->GetZref(i) ) return kFALSE;
1139 }
1140
1141 if ( fSigmaY != inTracklet->GetSigmaY() ) return kFALSE;
1142 if ( fSigmaY2 != inTracklet->GetSigmaY2() ) return kFALSE;
1143 if ( fTilt != inTracklet->GetTilt() ) return kFALSE;
1144 if ( fPadLength != inTracklet->GetPadLength() ) return kFALSE;
1145
1146 for (Int_t i = 0; i < knTimebins; i++){
1147 if ( fX[i] != inTracklet->GetX(i) ) return kFALSE;
1148 if ( fY[i] != inTracklet->GetY(i) ) return kFALSE;
1149 if ( fZ[i] != inTracklet->GetZ(i) ) return kFALSE;
1150 if ( fIndexes[i] != inTracklet->GetIndexes(i) ) return kFALSE;
1151 if ( fUsable[i] != inTracklet->IsUsable(i) ) return kFALSE;
1152 }
1153
1154 for (Int_t i=0; i < 2; i++){
1155 if ( fYfit[i] != inTracklet->GetYfit(i) ) return kFALSE;
1156 if ( fZfit[i] != inTracklet->GetZfit(i) ) return kFALSE;
1157 if ( fYfitR[i] != inTracklet->GetYfitR(i) ) return kFALSE;
1158 if ( fZfitR[i] != inTracklet->GetZfitR(i) ) return kFALSE;
1159 if ( fLabels[i] != inTracklet->GetLabels(i) ) return kFALSE;
1160 }
1161
1162 if ( fMeanz != inTracklet->GetMeanz() ) return kFALSE;
1163 if ( fZProb != inTracklet->GetZProb() ) return kFALSE;
1164 if ( fN2 != inTracklet->GetN2() ) return kFALSE;
1165 if ( fNUsed != inTracklet->GetNUsed() ) return kFALSE;
1166 if ( fFreq != inTracklet->GetFreq() ) return kFALSE;
1167 if ( fNChange != inTracklet->GetNChange() ) return kFALSE;
1168 if ( fNChange != inTracklet->GetNChange() ) return kFALSE;
1169
1170 if ( fC != inTracklet->GetC() ) return kFALSE;
1171 if ( fCC != inTracklet->GetCC() ) return kFALSE;
1172 if ( fChi2 != inTracklet->GetChi2() ) return kFALSE;
1173 // if ( fChi2Z != inTracklet->GetChi2Z() ) return kFALSE;
1174
1175 if ( fDet != inTracklet->GetDetector() ) return kFALSE;
1176 if ( fMom != inTracklet->GetMomentum() ) return kFALSE;
1177 if ( fdX != inTracklet->GetdX() ) return kFALSE;
1178
1179 for (Int_t iCluster = 0; iCluster < knTimebins; iCluster++){
1180 AliTRDcluster *curCluster = fClusters[iCluster];
1181 AliTRDcluster *inCluster = inTracklet->GetClusters(iCluster);
1182 if (curCluster && inCluster){
1183 if (! curCluster->IsEqual(inCluster) ) {
1184 curCluster->Print();
1185 inCluster->Print();
1186 return kFALSE;
1187 }
1188 } else {
1189 // if one cluster exists, and corresponding
1190 // in other tracklet doesn't - return kFALSE
1191 if(curCluster || inCluster) return kFALSE;
1192 }
1193 }
1194 return kTRUE;
1195}