Adding debug stream for first track reconstrunction (Marian)
[u/mrichter/AliRoot.git] / TPC / AliTPCtrackerMI.cxx
CommitLineData
1c53abe2 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
1c53abe2 16
17//-------------------------------------------------------
18// Implementation of the TPC tracker
19//
20// Origin: Marian Ivanov Marian.Ivanov@cern.ch
21//
34acb742 22// AliTPC parallel tracker
6d493ea0 23//
24// The track fitting is based on Kalaman filtering approach
25
26// The track finding steps:
27// 1. Seeding - with and without vertex constraint
28// - seeding with vertex constain done at first n^2 proble
29// - seeding without vertex constraint n^3 problem
30// 2. Tracking - follow prolongation road - find cluster - update kalman track
31
32// The seeding and tracking is repeated several times, in different seeding region.
33// This approach enables to find the track which cannot be seeded in some region of TPC
34// This can happen because of low momenta (track do not reach outer radius), or track is currently in the ded region between sectors, or the track is for the moment overlapped with other track (seed quality is poor) ...
35
36// With this approach we reach almost 100 % efficiency also for high occupancy events.
37// (If the seeding efficiency in a region is about 90 % than with logical or of several
38// regions we will reach 100% (in theory - supposing independence)
39
40// Repeating several seeding - tracking procedures some of the tracks can be find
41// several times.
42
43// The procedures to remove multi find tacks are impremented:
44// RemoveUsed2 - fast procedure n problem -
45// Algorithm - Sorting tracks according quality
46// remove tracks with some shared fraction
47// Sharing in respect to all tacks
48// Signing clusters in gold region
49// FindSplitted - slower algorithm n^2
50// Sort the tracks according quality
51// Loop over pair of tracks
52// If overlap with other track bigger than threshold - remove track
53//
54// FindCurling - Finds the pair of tracks which are curling
55// - About 10% of tracks can be find with this procedure
56// The combinatorial background is too big to be used in High
57// multiplicity environment
58// - n^2 problem - Slow procedure - currently it is disabled because of
59// low efficiency
60//
61// The number of splitted tracks can be reduced disabling the sharing of the cluster.
62// tpcRecoParam-> SetClusterSharing(kFALSE);
63// IT IS HIGHLY non recomended to use it in high flux enviroonment
64// Even using this switch some tracks can be found more than once
65// (because of multiple seeding and low quality tracks which will not cross full chamber)
66//
67//
68// The tracker itself can be debugged - the information about tracks can be stored in several // phases of the reconstruction
69// To enable storage of the TPC tracks in the ESD friend track
70// use AliTPCReconstructor::SetStreamLevel(n); where nis bigger 0
71//
72// The debug level - different procedure produce tree for numerical debugging
73// To enable them set AliTPCReconstructor::SetStreamLevel(n); where nis bigger 1
74//
92f513f5 75
76//
77// Adding systematic errors to the covariance:
78//
79// The systematic errors due to the misalignment and miscalibration are added to the covariance matrix
80// of the tracks (not to the clusters as they are dependent):
81// The parameters form AliTPCRecoParam are used AliTPCRecoParam::GetSystematicError
82// The systematic errors are expressed there in RMS - position (cm), angle (rad), curvature (1/cm)
83// The default values are 0.
84//
85// The sytematic errors are added to the covariance matrix in following places:
86//
87// 1. During fisrt itteration - AliTPCtrackerMI::FillESD
88// 2. Second iteration -
89// 2.a ITS->TPC - AliTPCtrackerMI::ReadSeeds
90// 2.b TPC->TRD - AliTPCtrackerMI::PropagateBack
91// 3. Third iteration -
92// 3.a TRD->TPC - AliTPCtrackerMI::ReadSeeds
93// 3.b TPC->ITS - AliTPCtrackerMI::RefitInward
94//
fd065ea2 95// There are several places in the code which can be numerically debuged
96// This code is keeped in order to enable code development and to check the calibration implementtion
97//
98// 1. ErrParam stream (Log level 9) - dump information about
99// 1.a) cluster
100// 2.a) cluster error estimate
101// 3.a) cluster shape estimate
102//
103//
1c53abe2 104//-------------------------------------------------------
47966a6d 105
106
107/* $Id$ */
108
cc5e9db0 109#include "Riostream.h"
6d171107 110#include <TClonesArray.h>
111#include <TFile.h>
112#include <TObjArray.h>
113#include <TTree.h>
a3232aae 114#include "AliLog.h"
47966a6d 115#include "AliComplexCluster.h"
af885e0f 116#include "AliESDEvent.h"
aad72f45 117#include "AliESDtrack.h"
118#include "AliESDVertex.h"
6c94f330 119#include "AliKink.h"
120#include "AliV0.h"
91162307 121#include "AliHelix.h"
91162307 122#include "AliRunLoader.h"
6d171107 123#include "AliTPCClustersRow.h"
124#include "AliTPCParam.h"
9996a03b 125#include "AliTPCReconstructor.h"
6d171107 126#include "AliTPCpolyTrack.h"
81e97e0d 127#include "AliTPCreco.h"
9350f379 128#include "AliTPCseed.h"
129
130#include "AliTPCtrackerSector.h"
b67e07dc 131#include "AliTPCtrackerMI.h"
6d171107 132#include "TStopwatch.h"
81e97e0d 133#include "AliTPCReconstructor.h"
81e97e0d 134#include "AliPID.h"
eea478d3 135#include "TTreeStream.h"
5d837844 136#include "AliAlignObj.h"
137#include "AliTrackPointArray.h"
6d493ea0 138#include "TRandom.h"
24db6af7 139#include "AliTPCcalibDB.h"
140#include "AliTPCTransform.h"
fd065ea2 141#include "AliTPCClusterParam.h"
19b00333 142
6d171107 143//
c9427e08 144
91162307 145ClassImp(AliTPCtrackerMI)
c9427e08 146
147
002af263 148
b67e07dc 149class AliTPCFastMath {
91162307 150public:
b67e07dc 151 AliTPCFastMath();
91162307 152 static Double_t FastAsin(Double_t x);
153 private:
b67e07dc 154 static Double_t fgFastAsin[20000]; //lookup table for fast asin computation
91162307 155};
c9427e08 156
b67e07dc 157Double_t AliTPCFastMath::fgFastAsin[20000];
2274b54b 158AliTPCFastMath gAliTPCFastMath; // needed to fill the LUT
c9427e08 159
b67e07dc 160AliTPCFastMath::AliTPCFastMath(){
161 //
162 // initialized lookup table;
91162307 163 for (Int_t i=0;i<10000;i++){
164 fgFastAsin[2*i] = TMath::ASin(i/10000.);
165 fgFastAsin[2*i+1] = (TMath::ASin((i+1)/10000.)-fgFastAsin[2*i]);
166 }
c9427e08 167}
168
b67e07dc 169Double_t AliTPCFastMath::FastAsin(Double_t x){
170 //
171 // return asin using lookup table
91162307 172 if (x>0){
173 Int_t index = int(x*10000);
174 return fgFastAsin[2*index]+(x*10000.-index)*fgFastAsin[2*index+1];
175 }
176 x*=-1;
177 Int_t index = int(x*10000);
178 return -(fgFastAsin[2*index]+(x*10000.-index)*fgFastAsin[2*index+1]);
1c53abe2 179}
e046d791 180//__________________________________________________________________
181AliTPCtrackerMI::AliTPCtrackerMI()
182 :AliTracker(),
183 fkNIS(0),
184 fInnerSec(0),
185 fkNOS(0),
186 fOuterSec(0),
187 fN(0),
188 fSectors(0),
189 fInput(0),
190 fOutput(0),
191 fSeedTree(0),
192 fTreeDebug(0),
193 fEvent(0),
194 fDebug(0),
195 fNewIO(kFALSE),
196 fNtracks(0),
197 fSeeds(0),
198 fIteration(0),
199 fParam(0),
200 fDebugStreamer(0)
201{
202 //
203 // default constructor
204 //
205}
206//_____________________________________________________________________
1c53abe2 207
208
209
91162307 210Int_t AliTPCtrackerMI::UpdateTrack(AliTPCseed * track, Int_t accept){
b67e07dc 211 //
212 //update track information using current cluster - track->fCurrentCluster
213
1c53abe2 214
b9671574 215 AliTPCclusterMI* c =track->GetCurrentCluster();
216 if (accept>0) track->SetCurrentClusterIndex1(track->GetCurrentClusterIndex1() | 0x8000); //sign not accepted clusters
c9427e08 217
b9671574 218 UInt_t i = track->GetCurrentClusterIndex1();
1c53abe2 219
220 Int_t sec=(i&0xff000000)>>24;
91162307 221 //Int_t row = (i&0x00ff0000)>>16;
b9671574 222 track->SetRow((i&0x00ff0000)>>16);
223 track->SetSector(sec);
1c53abe2 224 // Int_t index = i&0xFFFF;
b9671574 225 if (sec>=fParam->GetNInnerSector()) track->SetRow(track->GetRow()+fParam->GetNRowLow());
226 track->SetClusterIndex2(track->GetRow(), i);
91162307 227 //track->fFirstPoint = row;
228 //if ( track->fLastPoint<row) track->fLastPoint =row;
229 // if (track->fRow<0 || track->fRow>160) {
230 // printf("problem\n");
231 //}
b9671574 232 if (track->GetFirstPoint()>track->GetRow())
233 track->SetFirstPoint(track->GetRow());
234 if (track->GetLastPoint()<track->GetRow())
235 track->SetLastPoint(track->GetRow());
91162307 236
237
b9671574 238 track->SetClusterPointer(track->GetRow(),c);
1c53abe2 239 //
240
3f82c4f2 241 Double_t angle2 = track->GetSnp()*track->GetSnp();
1c53abe2 242 //
243 //SET NEW Track Point
244 //
85e2b57d 245 if (angle2<1) //PH sometimes angle2 is very big. To be investigated...
91162307 246 {
85e2b57d 247 angle2 = TMath::Sqrt(angle2/(1-angle2));
b9671574 248 AliTPCTrackerPoint &point =*(track->GetTrackPoint(track->GetRow()));
1c53abe2 249 //
b9671574 250 point.SetSigmaY(c->GetSigmaY2()/track->GetCurrentSigmaY2());
251 point.SetSigmaZ(c->GetSigmaZ2()/track->GetCurrentSigmaZ2());
252 point.SetErrY(sqrt(track->GetErrorY2()));
253 point.SetErrZ(sqrt(track->GetErrorZ2()));
1c53abe2 254 //
91162307 255 point.SetX(track->GetX());
256 point.SetY(track->GetY());
257 point.SetZ(track->GetZ());
258 point.SetAngleY(angle2);
259 point.SetAngleZ(track->GetTgl());
b9671574 260 if (point.IsShared()){
261 track->SetErrorY2(track->GetErrorY2()*4);
262 track->SetErrorZ2(track->GetErrorZ2()*4);
91162307 263 }
264 }
265
b9671574 266 Double_t chi2 = track->GetPredictedChi2(track->GetCurrentCluster());
91162307 267 //
1d91d749 268// track->SetErrorY2(track->GetErrorY2()*1.3);
269// track->SetErrorY2(track->GetErrorY2()+0.01);
270// track->SetErrorZ2(track->GetErrorZ2()*1.3);
271// track->SetErrorZ2(track->GetErrorZ2()+0.005);
91162307 272 //}
273 if (accept>0) return 0;
274 if (track->GetNumberOfClusters()%20==0){
275 // if (track->fHelixIn){
276 // TClonesArray & larr = *(track->fHelixIn);
277 // Int_t ihelix = larr.GetEntriesFast();
278 // new(larr[ihelix]) AliHelix(*track) ;
279 //}
1c53abe2 280 }
b9671574 281 track->SetNoCluster(0);
91162307 282 return track->Update(c,chi2,i);
283}
284
285
286
00055a22 287Int_t AliTPCtrackerMI::AcceptCluster(AliTPCseed * seed, AliTPCclusterMI * cluster)
91162307 288{
1c53abe2 289 //
91162307 290 // decide according desired precision to accept given
291 // cluster for tracking
00055a22 292 Double_t sy2=ErrY2(seed,cluster);
293 Double_t sz2=ErrZ2(seed,cluster);
1c53abe2 294
91162307 295 Double_t sdistancey2 = sy2+seed->GetSigmaY2();
296 Double_t sdistancez2 = sz2+seed->GetSigmaZ2();
297
b9671574 298 Double_t rdistancey2 = (seed->GetCurrentCluster()->GetY()-seed->GetY())*
299 (seed->GetCurrentCluster()->GetY()-seed->GetY())/sdistancey2;
300 Double_t rdistancez2 = (seed->GetCurrentCluster()->GetZ()-seed->GetZ())*
301 (seed->GetCurrentCluster()->GetZ()-seed->GetZ())/sdistancez2;
91162307 302
303 Double_t rdistance2 = rdistancey2+rdistancez2;
304 //Int_t accept =0;
1c53abe2 305
0de3e3c0 306 if (AliTPCReconstructor::StreamLevel()>5 && seed->GetNumberOfClusters()>20) {
fd065ea2 307 Float_t rmsy2 = seed->GetCurrentSigmaY2();
308 Float_t rmsz2 = seed->GetCurrentSigmaZ2();
e0e13b88 309 Float_t rmsy2p30 = seed->GetCMeanSigmaY2p30();
310 Float_t rmsz2p30 = seed->GetCMeanSigmaZ2p30();
42aec1f1 311 Float_t rmsy2p30R = seed->GetCMeanSigmaY2p30R();
312 Float_t rmsz2p30R = seed->GetCMeanSigmaZ2p30R();
e0e13b88 313
fd065ea2 314 AliExternalTrackParam param(*seed);
eea6b724 315 static TVectorD gcl(3),gtr(3);
316 Float_t gclf[3];
317 param.GetXYZ(gcl.GetMatrixArray());
318 cluster->GetGlobalXYZ(gclf);
319 gcl[0]=gclf[0]; gcl[1]=gclf[1]; gcl[2]=gclf[2];
fd065ea2 320 (*fDebugStreamer)<<"ErrParam"<<
e0e13b88 321 "Cl.="<<cluster<<
322 "T.="<<&param<<
eea6b724 323 "gcl.="<<&gcl<<
324 "gtr.="<<&gtr<<
e0e13b88 325 "erry2="<<sy2<<
326 "errz2="<<sz2<<
327 "rmsy2="<<rmsy2<<
328 "rmsz2="<<rmsz2<<
329 "rmsy2p30="<<rmsy2p30<<
330 "rmsz2p30="<<rmsz2p30<<
42aec1f1 331 "rmsy2p30R="<<rmsy2p30R<<
332 "rmsz2p30R="<<rmsz2p30R<<
e0e13b88 333 "\n";
fd065ea2 334 }
e0e13b88 335
91162307 336 if (rdistance2>16) return 3;
337
338
00055a22 339 if ((rdistancey2>9. || rdistancez2>9.) && cluster->GetType()==0)
91162307 340 return 2; //suspisiouce - will be changed
341
00055a22 342 if ((rdistancey2>6.25 || rdistancez2>6.25) && cluster->GetType()>0)
91162307 343 // strict cut on overlaped cluster
344 return 2; //suspisiouce - will be changed
345
00055a22 346 if ( (rdistancey2>1. || rdistancez2>6.25 )
91162307 347 && cluster->GetType()<0){
b9671574 348 seed->SetNFoundable(seed->GetNFoundable()-1);
91162307 349 return 2;
1c53abe2 350 }
91162307 351 return 0;
352}
353
354
1c53abe2 355
1c53abe2 356
1c53abe2 357//_____________________________________________________________________________
f8aae377 358AliTPCtrackerMI::AliTPCtrackerMI(const AliTPCParam *par):
e046d791 359AliTracker(),
360 fkNIS(par->GetNInnerSector()/2),
361 fInnerSec(0),
362 fkNOS(par->GetNOuterSector()/2),
363 fOuterSec(0),
364 fN(0),
365 fSectors(0),
366 fInput(0),
367 fOutput(0),
368 fSeedTree(0),
369 fTreeDebug(0),
370 fEvent(0),
371 fDebug(0),
372 fNewIO(0),
373 fNtracks(0),
374 fSeeds(0),
375 fIteration(0),
376 fParam(0),
377 fDebugStreamer(0)
1c53abe2 378{
379 //---------------------------------------------------------------------
380 // The main TPC tracker constructor
381 //---------------------------------------------------------------------
9350f379 382 fInnerSec=new AliTPCtrackerSector[fkNIS];
383 fOuterSec=new AliTPCtrackerSector[fkNOS];
91162307 384
1c53abe2 385 Int_t i;
386 for (i=0; i<fkNIS; i++) fInnerSec[i].Setup(par,0);
387 for (i=0; i<fkNOS; i++) fOuterSec[i].Setup(par,1);
388
91162307 389 fParam = par;
390 Int_t nrowlow = par->GetNRowLow();
391 Int_t nrowup = par->GetNRowUp();
392
393
394 for (Int_t i=0;i<nrowlow;i++){
395 fXRow[i] = par->GetPadRowRadiiLow(i);
396 fPadLength[i]= par->GetPadPitchLength(0,i);
397 fYMax[i] = fXRow[i]*TMath::Tan(0.5*par->GetInnerAngle());
398 }
399
400
401 for (Int_t i=0;i<nrowup;i++){
402 fXRow[i+nrowlow] = par->GetPadRowRadiiUp(i);
403 fPadLength[i+nrowlow] = par->GetPadPitchLength(60,i);
404 fYMax[i+nrowlow] = fXRow[i+nrowlow]*TMath::Tan(0.5*par->GetOuterAngle());
405 }
e046d791 406
81e97e0d 407 fDebugStreamer = new TTreeSRedirector("TPCdebug.root");
1c53abe2 408}
2fc0c115 409//________________________________________________________________________
58251ea0 410AliTPCtrackerMI::AliTPCtrackerMI(const AliTPCtrackerMI &t):
411 AliTracker(t),
e046d791 412 fkNIS(t.fkNIS),
413 fInnerSec(0),
414 fkNOS(t.fkNOS),
415 fOuterSec(0),
416 fN(0),
417 fSectors(0),
418 fInput(0),
419 fOutput(0),
420 fSeedTree(0),
421 fTreeDebug(0),
422 fEvent(0),
423 fDebug(0),
424 fNewIO(kFALSE),
425 fNtracks(0),
426 fSeeds(0),
427 fIteration(0),
428 fParam(0),
429 fDebugStreamer(0)
58251ea0 430{
2fc0c115 431 //------------------------------------
432 // dummy copy constructor
433 //------------------------------------------------------------------
e046d791 434 fOutput=t.fOutput;
2fc0c115 435}
436AliTPCtrackerMI & AliTPCtrackerMI::operator=(const AliTPCtrackerMI& /*r*/){
437 //------------------------------
438 // dummy
439 //--------------------------------------------------------------
440 return *this;
441}
1c53abe2 442//_____________________________________________________________________________
443AliTPCtrackerMI::~AliTPCtrackerMI() {
444 //------------------------------------------------------------------
445 // TPC tracker destructor
446 //------------------------------------------------------------------
447 delete[] fInnerSec;
448 delete[] fOuterSec;
449 if (fSeeds) {
450 fSeeds->Delete();
451 delete fSeeds;
452 }
81e97e0d 453 if (fDebugStreamer) delete fDebugStreamer;
1c53abe2 454}
455
1c53abe2 456
d26d9159 457void AliTPCtrackerMI::FillESD(TObjArray* arr)
91162307 458{
47966a6d 459 //
460 //
461 //fill esds using updated tracks
91162307 462 if (fEvent){
463 // write tracks to the event
464 // store index of the track
d26d9159 465 Int_t nseed=arr->GetEntriesFast();
51ad6848 466 //FindKinks(arr,fEvent);
91162307 467 for (Int_t i=0; i<nseed; i++) {
d26d9159 468 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
91162307 469 if (!pt) continue;
51ad6848 470 pt->UpdatePoints();
92f513f5 471 AddCovariance(pt);
8cecaa87 472 if (AliTPCReconstructor::StreamLevel()>1) {
473 (*fDebugStreamer)<<"Track0"<<
474 "Tr0.="<<pt<<
475 "\n";
476 }
eea478d3 477 // pt->PropagateTo(fParam->GetInnerRadiusLow());
478 if (pt->GetKinkIndex(0)<=0){ //don't propagate daughter tracks
479 pt->PropagateTo(fParam->GetInnerRadiusLow());
480 }
51ad6848 481
482 if (( pt->GetPoints()[2]- pt->GetPoints()[0])>5 && pt->GetPoints()[3]>0.8){
483 AliESDtrack iotrack;
484 iotrack.UpdateTrackParams(pt,AliESDtrack::kTPCin);
485 iotrack.SetTPCPoints(pt->GetPoints());
486 iotrack.SetKinkIndexes(pt->GetKinkIndexes());
81e97e0d 487 iotrack.SetV0Indexes(pt->GetV0Indexes());
488 // iotrack.SetTPCpid(pt->fTPCr);
51ad6848 489 //iotrack.SetTPCindex(i);
490 fEvent->AddTrack(&iotrack);
491 continue;
492 }
493
b9671574 494 if ( (pt->GetNumberOfClusters()>70)&& (Float_t(pt->GetNumberOfClusters())/Float_t(pt->GetNFoundable()))>0.55) {
d26d9159 495 AliESDtrack iotrack;
51ad6848 496 iotrack.UpdateTrackParams(pt,AliESDtrack::kTPCin);
497 iotrack.SetTPCPoints(pt->GetPoints());
d26d9159 498 //iotrack.SetTPCindex(i);
51ad6848 499 iotrack.SetKinkIndexes(pt->GetKinkIndexes());
81e97e0d 500 iotrack.SetV0Indexes(pt->GetV0Indexes());
501 // iotrack.SetTPCpid(pt->fTPCr);
d26d9159 502 fEvent->AddTrack(&iotrack);
a42a6bae 503 continue;
504 }
51ad6848 505 //
506 // short tracks - maybe decays
507
b9671574 508 if ( (pt->GetNumberOfClusters()>30) && (Float_t(pt->GetNumberOfClusters())/Float_t(pt->GetNFoundable()))>0.70) {
a42a6bae 509 Int_t found,foundable,shared;
510 pt->GetClusterStatistic(0,60,found, foundable,shared,kFALSE);
b9671574 511 if ( (found>20) && (pt->GetNShared()/float(pt->GetNumberOfClusters())<0.2)){
a42a6bae 512 AliESDtrack iotrack;
513 iotrack.UpdateTrackParams(pt,AliESDtrack::kTPCin);
514 //iotrack.SetTPCindex(i);
51ad6848 515 iotrack.SetTPCPoints(pt->GetPoints());
516 iotrack.SetKinkIndexes(pt->GetKinkIndexes());
81e97e0d 517 iotrack.SetV0Indexes(pt->GetV0Indexes());
518 //iotrack.SetTPCpid(pt->fTPCr);
a42a6bae 519 fEvent->AddTrack(&iotrack);
520 continue;
521 }
522 }
523
b9671574 524 if ( (pt->GetNumberOfClusters()>20) && (Float_t(pt->GetNumberOfClusters())/Float_t(pt->GetNFoundable()))>0.8) {
a42a6bae 525 Int_t found,foundable,shared;
526 pt->GetClusterStatistic(0,60,found, foundable,shared,kFALSE);
527 if (found<20) continue;
b9671574 528 if (pt->GetNShared()/float(pt->GetNumberOfClusters())>0.2) continue;
a42a6bae 529 //
530 AliESDtrack iotrack;
531 iotrack.UpdateTrackParams(pt,AliESDtrack::kTPCin);
51ad6848 532 iotrack.SetTPCPoints(pt->GetPoints());
533 iotrack.SetKinkIndexes(pt->GetKinkIndexes());
81e97e0d 534 iotrack.SetV0Indexes(pt->GetV0Indexes());
535 //iotrack.SetTPCpid(pt->fTPCr);
51ad6848 536 //iotrack.SetTPCindex(i);
537 fEvent->AddTrack(&iotrack);
538 continue;
539 }
540 // short tracks - secondaties
541 //
542 if ( (pt->GetNumberOfClusters()>30) ) {
543 Int_t found,foundable,shared;
544 pt->GetClusterStatistic(128,158,found, foundable,shared,kFALSE);
b9671574 545 if ( (found>20) && (pt->GetNShared()/float(pt->GetNumberOfClusters())<0.2) &&float(found)/float(foundable)>0.8){
51ad6848 546 AliESDtrack iotrack;
547 iotrack.UpdateTrackParams(pt,AliESDtrack::kTPCin);
548 iotrack.SetTPCPoints(pt->GetPoints());
549 iotrack.SetKinkIndexes(pt->GetKinkIndexes());
81e97e0d 550 iotrack.SetV0Indexes(pt->GetV0Indexes());
551 //iotrack.SetTPCpid(pt->fTPCr);
51ad6848 552 //iotrack.SetTPCindex(i);
553 fEvent->AddTrack(&iotrack);
554 continue;
555 }
556 }
557
558 if ( (pt->GetNumberOfClusters()>15)) {
559 Int_t found,foundable,shared;
560 pt->GetClusterStatistic(138,158,found, foundable,shared,kFALSE);
561 if (found<15) continue;
e7eb17e4 562 if (foundable<=0) continue;
b9671574 563 if (pt->GetNShared()/float(pt->GetNumberOfClusters())>0.2) continue;
51ad6848 564 if (float(found)/float(foundable)<0.8) continue;
565 //
566 AliESDtrack iotrack;
567 iotrack.UpdateTrackParams(pt,AliESDtrack::kTPCin);
568 iotrack.SetTPCPoints(pt->GetPoints());
569 iotrack.SetKinkIndexes(pt->GetKinkIndexes());
81e97e0d 570 iotrack.SetV0Indexes(pt->GetV0Indexes());
571 // iotrack.SetTPCpid(pt->fTPCr);
a42a6bae 572 //iotrack.SetTPCindex(i);
573 fEvent->AddTrack(&iotrack);
574 continue;
575 }
91162307 576 }
1c53abe2 577 }
51ad6848 578 printf("Number of filled ESDs-\t%d\n",fEvent->GetNumberOfTracks());
d26d9159 579}
580
d26d9159 581
1c53abe2 582
1c53abe2 583
584
91162307 585Double_t AliTPCtrackerMI::ErrY2(AliTPCseed* seed, AliTPCclusterMI * cl){
586 //
587 //
fd065ea2 588 // Use calibrated cluster error from OCDB
91162307 589 //
fd065ea2 590 AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
91162307 591 //
a1ec4d07 592 Float_t z = TMath::Abs(fParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
91162307 593 Int_t ctype = cl->GetType();
fd065ea2 594 Int_t type = (cl->GetRow()<63) ? 0: (cl->GetRow()>126) ? 1:2;
595 Double_t angle = seed->GetSnp()*seed->GetSnp();
e0e13b88 596 angle = TMath::Sqrt(TMath::Abs(angle/(1.-angle)));
fd065ea2 597 Double_t erry2 = clparam->GetError0Par(0,type, z,angle);
598 if (ctype<0) {
599 erry2+=0.5; // edge cluster
600 }
601 erry2*=erry2;
602 seed->SetErrorY2(erry2);
603 //
604 return erry2;
605
606//calculate look-up table at the beginning
607// static Bool_t ginit = kFALSE;
608// static Float_t gnoise1,gnoise2,gnoise3;
609// static Float_t ggg1[10000];
610// static Float_t ggg2[10000];
611// static Float_t ggg3[10000];
612// static Float_t glandau1[10000];
613// static Float_t glandau2[10000];
614// static Float_t glandau3[10000];
615// //
616// static Float_t gcor01[500];
617// static Float_t gcor02[500];
618// static Float_t gcorp[500];
619// //
1c53abe2 620
fd065ea2 621// //
622// if (ginit==kFALSE){
623// for (Int_t i=1;i<500;i++){
624// Float_t rsigma = float(i)/100.;
625// gcor02[i] = TMath::Max(0.78 +TMath::Exp(7.4*(rsigma-1.2)),0.6);
626// gcor01[i] = TMath::Max(0.72 +TMath::Exp(3.36*(rsigma-1.2)),0.6);
627// gcorp[i] = TMath::Max(TMath::Power((rsigma+0.5),1.5),1.2);
628// }
629
630// //
631// for (Int_t i=3;i<10000;i++){
632// //
633// //
634// // inner sector
635// Float_t amp = float(i);
636// Float_t padlength =0.75;
637// gnoise1 = 0.0004/padlength;
638// Float_t nel = 0.268*amp;
639// Float_t nprim = 0.155*amp;
640// ggg1[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.001*nel/(padlength*padlength))/nel;
641// glandau1[i] = (2.+0.12*nprim)*0.5* (2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
642// if (glandau1[i]>1) glandau1[i]=1;
643// glandau1[i]*=padlength*padlength/12.;
644// //
645// // outer short
646// padlength =1.;
647// gnoise2 = 0.0004/padlength;
648// nel = 0.3*amp;
649// nprim = 0.133*amp;
650// ggg2[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
651// glandau2[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
652// if (glandau2[i]>1) glandau2[i]=1;
653// glandau2[i]*=padlength*padlength/12.;
654// //
655// //
656// // outer long
657// padlength =1.5;
658// gnoise3 = 0.0004/padlength;
659// nel = 0.3*amp;
660// nprim = 0.133*amp;
661// ggg3[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
662// glandau3[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
663// if (glandau3[i]>1) glandau3[i]=1;
664// glandau3[i]*=padlength*padlength/12.;
665// //
666// }
667// ginit = kTRUE;
668// }
669// //
670// //
671// //
672// Int_t amp = int(TMath::Abs(cl->GetQ()));
673// if (amp>9999) {
674// seed->SetErrorY2(1.);
675// return 1.;
676// }
677// Float_t snoise2;
678// Float_t z = TMath::Abs(fParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
679// Int_t ctype = cl->GetType();
680// Float_t padlength= GetPadPitchLength(seed->GetRow());
681// Double_t angle2 = seed->GetSnp()*seed->GetSnp();
682// angle2 = angle2/(1-angle2);
683// //
684// //cluster "quality"
685// Int_t rsigmay = int(100.*cl->GetSigmaY2()/(seed->GetCurrentSigmaY2()));
686// Float_t res;
687// //
688// if (fSectors==fInnerSec){
689// snoise2 = gnoise1;
690// res = ggg1[amp]*z+glandau1[amp]*angle2;
691// if (ctype==0) res *= gcor01[rsigmay];
692// if ((ctype>0)){
693// res+=0.002;
694// res*= gcorp[rsigmay];
695// }
696// }
697// else {
698// if (padlength<1.1){
699// snoise2 = gnoise2;
700// res = ggg2[amp]*z+glandau2[amp]*angle2;
701// if (ctype==0) res *= gcor02[rsigmay];
702// if ((ctype>0)){
703// res+=0.002;
704// res*= gcorp[rsigmay];
705// }
706// }
707// else{
708// snoise2 = gnoise3;
709// res = ggg3[amp]*z+glandau3[amp]*angle2;
710// if (ctype==0) res *= gcor02[rsigmay];
711// if ((ctype>0)){
712// res+=0.002;
713// res*= gcorp[rsigmay];
714// }
715// }
716// }
717
718// if (ctype<0){
719// res+=0.005;
720// res*=2.4; // overestimate error 2 times
721// }
722// res+= snoise2;
91162307 723
fd065ea2 724// if (res<2*snoise2)
725// res = 2*snoise2;
91162307 726
fd065ea2 727// seed->SetErrorY2(res);
728// return res;
1c53abe2 729
730
91162307 731}
c9427e08 732
733
734
91162307 735Double_t AliTPCtrackerMI::ErrZ2(AliTPCseed* seed, AliTPCclusterMI * cl){
736 //
737 //
fd065ea2 738 // Use calibrated cluster error from OCDB
91162307 739 //
fd065ea2 740 AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
91162307 741 //
a1ec4d07 742 Float_t z = TMath::Abs(fParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
91162307 743 Int_t ctype = cl->GetType();
fd065ea2 744 Int_t type = (cl->GetRow()<63) ? 0: (cl->GetRow()>126) ? 1:2;
91162307 745 //
3f82c4f2 746 Double_t angle2 = seed->GetSnp()*seed->GetSnp();
91162307 747 angle2 = seed->GetTgl()*seed->GetTgl()*(1+angle2/(1-angle2));
e0e13b88 748 Double_t angle = TMath::Sqrt(TMath::Abs(angle2));
fd065ea2 749 Double_t errz2 = clparam->GetError0Par(1,type, z,angle);
750 if (ctype<0) {
751 errz2+=0.5; // edge cluster
91162307 752 }
fd065ea2 753 errz2*=errz2;
754 seed->SetErrorZ2(errz2);
755 //
756 return errz2;
91162307 757
fd065ea2 758
759
760// //seed->SetErrorY2(0.1);
761// //return 0.1;
762// //calculate look-up table at the beginning
763// static Bool_t ginit = kFALSE;
764// static Float_t gnoise1,gnoise2,gnoise3;
765// static Float_t ggg1[10000];
766// static Float_t ggg2[10000];
767// static Float_t ggg3[10000];
768// static Float_t glandau1[10000];
769// static Float_t glandau2[10000];
770// static Float_t glandau3[10000];
771// //
772// static Float_t gcor01[1000];
773// static Float_t gcor02[1000];
774// static Float_t gcorp[1000];
775// //
776
777// //
778// if (ginit==kFALSE){
779// for (Int_t i=1;i<1000;i++){
780// Float_t rsigma = float(i)/100.;
781// gcor02[i] = TMath::Max(0.81 +TMath::Exp(6.8*(rsigma-1.2)),0.6);
782// gcor01[i] = TMath::Max(0.72 +TMath::Exp(2.04*(rsigma-1.2)),0.6);
783// gcorp[i] = TMath::Max(TMath::Power((rsigma+0.5),1.5),1.2);
784// }
785
786// //
787// for (Int_t i=3;i<10000;i++){
788// //
789// //
790// // inner sector
791// Float_t amp = float(i);
792// Float_t padlength =0.75;
793// gnoise1 = 0.0004/padlength;
794// Float_t nel = 0.268*amp;
795// Float_t nprim = 0.155*amp;
796// ggg1[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.001*nel/(padlength*padlength))/nel;
797// glandau1[i] = (2.+0.12*nprim)*0.5* (2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
798// if (glandau1[i]>1) glandau1[i]=1;
799// glandau1[i]*=padlength*padlength/12.;
800// //
801// // outer short
802// padlength =1.;
803// gnoise2 = 0.0004/padlength;
804// nel = 0.3*amp;
805// nprim = 0.133*amp;
806// ggg2[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
807// glandau2[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
808// if (glandau2[i]>1) glandau2[i]=1;
809// glandau2[i]*=padlength*padlength/12.;
810// //
811// //
812// // outer long
813// padlength =1.5;
814// gnoise3 = 0.0004/padlength;
815// nel = 0.3*amp;
816// nprim = 0.133*amp;
817// ggg3[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
818// glandau3[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
819// if (glandau3[i]>1) glandau3[i]=1;
820// glandau3[i]*=padlength*padlength/12.;
821// //
822// }
823// ginit = kTRUE;
824// }
825// //
826// //
827// //
828// Int_t amp = int(TMath::Abs(cl->GetQ()));
829// if (amp>9999) {
830// seed->SetErrorY2(1.);
831// return 1.;
832// }
833// Float_t snoise2;
834// Float_t z = TMath::Abs(fParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
835// Int_t ctype = cl->GetType();
836// Float_t padlength= GetPadPitchLength(seed->GetRow());
837// //
838// Double_t angle2 = seed->GetSnp()*seed->GetSnp();
839// // if (angle2<0.6) angle2 = 0.6;
840// angle2 = seed->GetTgl()*seed->GetTgl()*(1+angle2/(1-angle2));
841// //
842// //cluster "quality"
843// Int_t rsigmaz = int(100.*cl->GetSigmaZ2()/(seed->GetCurrentSigmaZ2()));
844// Float_t res;
845// //
846// if (fSectors==fInnerSec){
847// snoise2 = gnoise1;
848// res = ggg1[amp]*z+glandau1[amp]*angle2;
849// if (ctype==0) res *= gcor01[rsigmaz];
850// if ((ctype>0)){
851// res+=0.002;
852// res*= gcorp[rsigmaz];
853// }
854// }
855// else {
856// if (padlength<1.1){
857// snoise2 = gnoise2;
858// res = ggg2[amp]*z+glandau2[amp]*angle2;
859// if (ctype==0) res *= gcor02[rsigmaz];
860// if ((ctype>0)){
861// res+=0.002;
862// res*= gcorp[rsigmaz];
863// }
864// }
865// else{
866// snoise2 = gnoise3;
867// res = ggg3[amp]*z+glandau3[amp]*angle2;
868// if (ctype==0) res *= gcor02[rsigmaz];
869// if ((ctype>0)){
870// res+=0.002;
871// res*= gcorp[rsigmaz];
872// }
873// }
874// }
875
876// if (ctype<0){
877// res+=0.002;
878// res*=1.3;
879// }
880// if ((ctype<0) &&amp<70){
881// res+=0.002;
882// res*=1.3;
883// }
884// res += snoise2;
885// if (res<2*snoise2)
886// res = 2*snoise2;
887// if (res>3) res =3;
888// seed->SetErrorZ2(res);
889// return res;
91162307 890}
891
892
893
91162307 894
1c53abe2 895
c9427e08 896void AliTPCtrackerMI::RotateToLocal(AliTPCseed *seed)
897{
898 //rotate to track "local coordinata
899 Float_t x = seed->GetX();
900 Float_t y = seed->GetY();
901 Float_t ymax = x*TMath::Tan(0.5*fSectors->GetAlpha());
91162307 902
c9427e08 903 if (y > ymax) {
b9671574 904 seed->SetRelativeSector((seed->GetRelativeSector()+1) % fN);
c9427e08 905 if (!seed->Rotate(fSectors->GetAlpha()))
906 return;
907 } else if (y <-ymax) {
b9671574 908 seed->SetRelativeSector((seed->GetRelativeSector()-1+fN) % fN);
c9427e08 909 if (!seed->Rotate(-fSectors->GetAlpha()))
910 return;
911 }
1c53abe2 912
c9427e08 913}
1c53abe2 914
915
916
1c53abe2 917//_____________________________________________________________________________
b67e07dc 918Double_t AliTPCtrackerMI::F1old(Double_t x1,Double_t y1,
1c53abe2 919 Double_t x2,Double_t y2,
920 Double_t x3,Double_t y3)
921{
922 //-----------------------------------------------------------------
923 // Initial approximation of the track curvature
924 //-----------------------------------------------------------------
925 Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
926 Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
927 (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
928 Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
929 (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
930
931 Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
91162307 932 if ( xr*xr+yr*yr<=0.00000000000001) return 100;
1c53abe2 933 return -xr*yr/sqrt(xr*xr+yr*yr);
934}
935
936
91162307 937
1c53abe2 938//_____________________________________________________________________________
b67e07dc 939Double_t AliTPCtrackerMI::F1(Double_t x1,Double_t y1,
91162307 940 Double_t x2,Double_t y2,
941 Double_t x3,Double_t y3)
942{
943 //-----------------------------------------------------------------
944 // Initial approximation of the track curvature
945 //-----------------------------------------------------------------
946 x3 -=x1;
947 x2 -=x1;
948 y3 -=y1;
949 y2 -=y1;
950 //
951 Double_t det = x3*y2-x2*y3;
952 if (det==0) {
953 return 100;
954 }
955 //
956 Double_t u = 0.5* (x2*(x2-x3)+y2*(y2-y3))/det;
957 Double_t x0 = x3*0.5-y3*u;
958 Double_t y0 = y3*0.5+x3*u;
959 Double_t c2 = 1/TMath::Sqrt(x0*x0+y0*y0);
960 if (det<0) c2*=-1;
961 return c2;
962}
963
964
b67e07dc 965Double_t AliTPCtrackerMI::F2(Double_t x1,Double_t y1,
1c53abe2 966 Double_t x2,Double_t y2,
967 Double_t x3,Double_t y3)
91162307 968{
969 //-----------------------------------------------------------------
970 // Initial approximation of the track curvature
971 //-----------------------------------------------------------------
972 x3 -=x1;
973 x2 -=x1;
974 y3 -=y1;
975 y2 -=y1;
976 //
977 Double_t det = x3*y2-x2*y3;
978 if (det==0) {
979 return 100;
980 }
981 //
982 Double_t u = 0.5* (x2*(x2-x3)+y2*(y2-y3))/det;
983 Double_t x0 = x3*0.5-y3*u;
984 Double_t y0 = y3*0.5+x3*u;
985 Double_t c2 = 1/TMath::Sqrt(x0*x0+y0*y0);
986 if (det<0) c2*=-1;
987 x0+=x1;
988 x0*=c2;
989 return x0;
990}
991
992
993
994//_____________________________________________________________________________
b67e07dc 995Double_t AliTPCtrackerMI::F2old(Double_t x1,Double_t y1,
91162307 996 Double_t x2,Double_t y2,
997 Double_t x3,Double_t y3)
1c53abe2 998{
999 //-----------------------------------------------------------------
1000 // Initial approximation of the track curvature times center of curvature
1001 //-----------------------------------------------------------------
1002 Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
1003 Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
1004 (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
1005 Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
1006 (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
1007
1008 Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
1009
1010 return -a/(d*y1-b)*xr/sqrt(xr*xr+yr*yr);
1011}
1012
1013//_____________________________________________________________________________
b67e07dc 1014Double_t AliTPCtrackerMI::F3(Double_t x1,Double_t y1,
1c53abe2 1015 Double_t x2,Double_t y2,
1016 Double_t z1,Double_t z2)
1017{
1018 //-----------------------------------------------------------------
1019 // Initial approximation of the tangent of the track dip angle
1020 //-----------------------------------------------------------------
1021 return (z1 - z2)/sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
1022}
1023
1024
b67e07dc 1025Double_t AliTPCtrackerMI::F3n(Double_t x1,Double_t y1,
91162307 1026 Double_t x2,Double_t y2,
1027 Double_t z1,Double_t z2, Double_t c)
1c53abe2 1028{
91162307 1029 //-----------------------------------------------------------------
1030 // Initial approximation of the tangent of the track dip angle
1031 //-----------------------------------------------------------------
1032
1033 // Double_t angle1;
1034
1035 //angle1 = (z1-z2)*c/(TMath::ASin(c*x1-ni)-TMath::ASin(c*x2-ni));
1c53abe2 1036 //
91162307 1037 Double_t d = TMath::Sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
1038 if (TMath::Abs(d*c*0.5)>1) return 0;
1039 // Double_t angle2 = TMath::ASin(d*c*0.5);
b67e07dc 1040 // Double_t angle2 = AliTPCFastMath::FastAsin(d*c*0.5);
1041 Double_t angle2 = (d*c*0.5>0.1)? TMath::ASin(d*c*0.5): AliTPCFastMath::FastAsin(d*c*0.5);
91162307 1042
1043 angle2 = (z1-z2)*c/(angle2*2.);
1044 return angle2;
1045}
1046
1047Bool_t AliTPCtrackerMI::GetProlongation(Double_t x1, Double_t x2, Double_t x[5], Double_t &y, Double_t &z)
1048{//-----------------------------------------------------------------
1049 // This function find proloncation of a track to a reference plane x=x2.
1050 //-----------------------------------------------------------------
1051
1052 Double_t dx=x2-x1;
1053
1054 if (TMath::Abs(x[4]*x1 - x[2]) >= 0.999) {
1055 return kFALSE;
1c53abe2 1056 }
f8aae377 1057
91162307 1058 Double_t c1=x[4]*x1 - x[2], r1=sqrt(1.- c1*c1);
1059 Double_t c2=x[4]*x2 - x[2], r2=sqrt(1.- c2*c2);
1060 y = x[0];
1061 z = x[1];
1062
1063 Double_t dy = dx*(c1+c2)/(r1+r2);
1064 Double_t dz = 0;
1065 //
1066 Double_t delta = x[4]*dx*(c1+c2)/(c1*r2 + c2*r1);
1067
1068 if (TMath::Abs(delta)>0.01){
1069 dz = x[3]*TMath::ASin(delta)/x[4];
1070 }else{
b67e07dc 1071 dz = x[3]*AliTPCFastMath::FastAsin(delta)/x[4];
91162307 1072 }
1073
b67e07dc 1074 //dz = x[3]*AliTPCFastMath::FastAsin(delta)/x[4];
f8aae377 1075
91162307 1076 y+=dy;
1077 z+=dz;
1078
1079 return kTRUE;
1c53abe2 1080}
1081
d26d9159 1082Int_t AliTPCtrackerMI::LoadClusters (TTree *tree)
1083{
1084 //
1085 //
1086 fInput = tree;
1087 return LoadClusters();
1088}
91162307 1089
af86c1fd 1090
1091Int_t AliTPCtrackerMI::LoadClusters(TObjArray *arr)
1092{
1093 //
1094 // load clusters to the memory
1095 AliTPCClustersRow *clrow = 0x0;
1096 Int_t lower = arr->LowerBound();
1097 Int_t entries = arr->GetEntriesFast();
1098
1099 for (Int_t i=lower; i<entries; i++) {
1100 clrow = (AliTPCClustersRow*) arr->At(i);
1101 if(!clrow) continue;
1102 if(!clrow->GetArray()) continue;
1103
1104 //
1105 Int_t sec,row;
1106 fParam->AdjustSectorRow(clrow->GetID(),sec,row);
1107
1108 for (Int_t icl=0; icl<clrow->GetArray()->GetEntriesFast(); icl++){
1109 Transform((AliTPCclusterMI*)(clrow->GetArray()->At(icl)));
1110 }
1111 //
1112 if (clrow->GetArray()->GetEntriesFast()<=0) continue;
bd26fa83 1113 AliTPCtrackerRow * tpcrow=0;
af86c1fd 1114 Int_t left=0;
1115 if (sec<fkNIS*2){
1116 tpcrow = &(fInnerSec[sec%fkNIS][row]);
1117 left = sec/fkNIS;
1118 }
1119 else{
1120 tpcrow = &(fOuterSec[(sec-fkNIS*2)%fkNOS][row]);
1121 left = (sec-fkNIS*2)/fkNOS;
1122 }
1123 if (left ==0){
1124 tpcrow->SetN1(clrow->GetArray()->GetEntriesFast());
1125 tpcrow->SetClusters1(new AliTPCclusterMI[tpcrow->GetN1()]);
1126 for (Int_t i=0;i<tpcrow->GetN1();i++)
1127 tpcrow->SetCluster1(i, *(AliTPCclusterMI*)(clrow->GetArray()->At(i)));
1128 }
1129 if (left ==1){
1130 tpcrow->SetN2(clrow->GetArray()->GetEntriesFast());
1131 tpcrow->SetClusters2(new AliTPCclusterMI[tpcrow->GetN2()]);
1132 for (Int_t i=0;i<tpcrow->GetN2();i++)
1133 tpcrow->SetCluster2(i,*(AliTPCclusterMI*)(clrow->GetArray()->At(i)));
1134 }
1135 }
1136 //
1137 delete clrow;
1138 LoadOuterSectors();
1139 LoadInnerSectors();
1140 return 0;
1141}
1142
1143
91162307 1144Int_t AliTPCtrackerMI::LoadClusters()
1c53abe2 1145{
1146 //
1147 // load clusters to the memory
91162307 1148 AliTPCClustersRow *clrow= new AliTPCClustersRow;
1149 clrow->SetClass("AliTPCclusterMI");
1150 clrow->SetArray(0);
1151 clrow->GetArray()->ExpandCreateFast(10000);
1152 //
1153 // TTree * tree = fClustersArray.GetTree();
1154
1155 TTree * tree = fInput;
1156 TBranch * br = tree->GetBranch("Segment");
1157 br->SetAddress(&clrow);
1158 //
1159 Int_t j=Int_t(tree->GetEntries());
1c53abe2 1160 for (Int_t i=0; i<j; i++) {
91162307 1161 br->GetEntry(i);
1162 //
1163 Int_t sec,row;
1164 fParam->AdjustSectorRow(clrow->GetID(),sec,row);
0201b65c 1165 for (Int_t icl=0; icl<clrow->GetArray()->GetEntriesFast(); icl++){
75fb37cc 1166 Transform((AliTPCclusterMI*)(clrow->GetArray()->At(icl)));
0201b65c 1167 }
91162307 1168 //
bd26fa83 1169 AliTPCtrackerRow * tpcrow=0;
91162307 1170 Int_t left=0;
1171 if (sec<fkNIS*2){
1172 tpcrow = &(fInnerSec[sec%fkNIS][row]);
1173 left = sec/fkNIS;
1174 }
1175 else{
1176 tpcrow = &(fOuterSec[(sec-fkNIS*2)%fkNOS][row]);
1177 left = (sec-fkNIS*2)/fkNOS;
1178 }
1179 if (left ==0){
b9671574 1180 tpcrow->SetN1(clrow->GetArray()->GetEntriesFast());
1181 tpcrow->SetClusters1(new AliTPCclusterMI[tpcrow->GetN1()]);
1182 for (Int_t i=0;i<tpcrow->GetN1();i++)
1183 tpcrow->SetCluster1(i, *(AliTPCclusterMI*)(clrow->GetArray()->At(i)));
91162307 1184 }
1185 if (left ==1){
b9671574 1186 tpcrow->SetN2(clrow->GetArray()->GetEntriesFast());
1187 tpcrow->SetClusters2(new AliTPCclusterMI[tpcrow->GetN2()]);
1188 for (Int_t i=0;i<tpcrow->GetN2();i++)
1189 tpcrow->SetCluster2(i,*(AliTPCclusterMI*)(clrow->GetArray()->At(i)));
91162307 1190 }
1c53abe2 1191 }
91162307 1192 //
1193 delete clrow;
1194 LoadOuterSectors();
1195 LoadInnerSectors();
1196 return 0;
1c53abe2 1197}
1198
1199
91162307 1200void AliTPCtrackerMI::UnloadClusters()
1201{
1202 //
1203 // unload clusters from the memory
1204 //
1205 Int_t nrows = fOuterSec->GetNRows();
1206 for (Int_t sec = 0;sec<fkNOS;sec++)
1207 for (Int_t row = 0;row<nrows;row++){
bd26fa83 1208 AliTPCtrackerRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
982aff31 1209 // if (tpcrow){
1210 // if (tpcrow->fClusters1) delete []tpcrow->fClusters1;
1211 // if (tpcrow->fClusters2) delete []tpcrow->fClusters2;
1212 //}
1213 tpcrow->ResetClusters();
1c53abe2 1214 }
91162307 1215 //
1216 nrows = fInnerSec->GetNRows();
1217 for (Int_t sec = 0;sec<fkNIS;sec++)
1218 for (Int_t row = 0;row<nrows;row++){
bd26fa83 1219 AliTPCtrackerRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
982aff31 1220 //if (tpcrow){
1221 // if (tpcrow->fClusters1) delete []tpcrow->fClusters1;
1222 //if (tpcrow->fClusters2) delete []tpcrow->fClusters2;
1223 //}
1224 tpcrow->ResetClusters();
91162307 1225 }
1226
1227 return ;
1c53abe2 1228}
1229
002af263 1230void AliTPCtrackerMI::FillClusterArray(TObjArray* array) const{
1231 //
1232 // Filling cluster to the array - For visualization purposes
1233 //
1234 Int_t nrows=0;
1235 nrows = fOuterSec->GetNRows();
1236 for (Int_t sec = 0;sec<fkNOS;sec++)
1237 for (Int_t row = 0;row<nrows;row++){
1238 AliTPCtrackerRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
1239 if (!tpcrow) continue;
1240 for (Int_t icl = 0;icl<tpcrow->GetN();icl++){
1241 array->AddLast((TObject*)((*tpcrow)[icl]));
1242 }
1243 }
1244 nrows = fInnerSec->GetNRows();
1245 for (Int_t sec = 0;sec<fkNIS;sec++)
1246 for (Int_t row = 0;row<nrows;row++){
1247 AliTPCtrackerRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
1248 if (!tpcrow) continue;
1249 for (Int_t icl = 0;icl<tpcrow->GetN();icl++){
1250 array->AddLast((TObject*)(*tpcrow)[icl]);
1251 }
1252 }
1253}
1254
1255
75fb37cc 1256void AliTPCtrackerMI::Transform(AliTPCclusterMI * cluster){
893468d9 1257 //
1258 //
1259 //
24db6af7 1260 AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
1261 if (!transform) {
1262 AliFatal("Tranformations not in calibDB");
1263 }
1264 Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()};
1265 Int_t i[1]={cluster->GetDetector()};
022ee144 1266 transform->Transform(x,i,0,1);
1267 if (!AliTPCReconstructor::GetRecoParam()->GetBYMirror()){
1268 if (cluster->GetDetector()%36>17){
1269 x[1]*=-1;
1270 }
1271 }
1272
24db6af7 1273 //
1274 // in debug mode check the transformation
1275 //
9e83cf6d 1276 if (AliTPCReconstructor::StreamLevel()>1) {
af86c1fd 1277 Float_t gx[3];
1278 cluster->GetGlobalXYZ(gx);
002af263 1279 Int_t event = (fEvent==NULL)? 0: fEvent->GetEventNumberInFile();
24db6af7 1280 TTreeSRedirector &cstream = *fDebugStreamer;
1281 cstream<<"Transform"<<
002af263 1282 "event="<<event<<
24db6af7 1283 "x0="<<x[0]<<
1284 "x1="<<x[1]<<
1285 "x2="<<x[2]<<
af86c1fd 1286 "gx0="<<gx[0]<<
1287 "gx1="<<gx[1]<<
1288 "gx2="<<gx[2]<<
24db6af7 1289 "Cl.="<<cluster<<
1290 "\n";
1291 }
1292 cluster->SetX(x[0]);
1293 cluster->SetY(x[1]);
1294 cluster->SetZ(x[2]);
19b00333 1295 // The old stuff:
1296
0201b65c 1297 //
1298 //
1299 //
1300 //if (!fParam->IsGeoRead()) fParam->ReadGeoMatrices();
1301 TGeoHMatrix *mat = fParam->GetClusterMatrix(cluster->GetDetector());
e7eb17e4 1302 //TGeoHMatrix mat;
0201b65c 1303 Double_t pos[3]= {cluster->GetX(),cluster->GetY(),cluster->GetZ()};
a7760332 1304 Double_t posC[3]={cluster->GetX(),cluster->GetY(),cluster->GetZ()};
1305 if (mat) mat->LocalToMaster(pos,posC);
1306 else{
1307 // chack Loading of Geo matrices from GeoManager - TEMPORARY FIX
1308 }
24db6af7 1309 cluster->SetX(posC[0]);
1310 cluster->SetY(posC[1]);
1311 cluster->SetZ(posC[2]);
0201b65c 1312}
1c53abe2 1313
1314//_____________________________________________________________________________
91162307 1315Int_t AliTPCtrackerMI::LoadOuterSectors() {
1c53abe2 1316 //-----------------------------------------------------------------
91162307 1317 // This function fills outer TPC sectors with clusters.
1c53abe2 1318 //-----------------------------------------------------------------
91162307 1319 Int_t nrows = fOuterSec->GetNRows();
1320 UInt_t index=0;
1321 for (Int_t sec = 0;sec<fkNOS;sec++)
1322 for (Int_t row = 0;row<nrows;row++){
bd26fa83 1323 AliTPCtrackerRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
91162307 1324 Int_t sec2 = sec+2*fkNIS;
1325 //left
b9671574 1326 Int_t ncl = tpcrow->GetN1();
91162307 1327 while (ncl--) {
b9671574 1328 AliTPCclusterMI *c= (tpcrow->GetCluster1(ncl));
91162307 1329 index=(((sec2<<8)+row)<<16)+ncl;
1330 tpcrow->InsertCluster(c,index);
1331 }
1332 //right
b9671574 1333 ncl = tpcrow->GetN2();
91162307 1334 while (ncl--) {
b9671574 1335 AliTPCclusterMI *c= (tpcrow->GetCluster2(ncl));
91162307 1336 index=((((sec2+fkNOS)<<8)+row)<<16)+ncl;
1337 tpcrow->InsertCluster(c,index);
1338 }
1339 //
1340 // write indexes for fast acces
1341 //
1342 for (Int_t i=0;i<510;i++)
b9671574 1343 tpcrow->SetFastCluster(i,-1);
91162307 1344 for (Int_t i=0;i<tpcrow->GetN();i++){
1345 Int_t zi = Int_t((*tpcrow)[i]->GetZ()+255.);
b9671574 1346 tpcrow->SetFastCluster(zi,i); // write index
91162307 1347 }
1348 Int_t last = 0;
1349 for (Int_t i=0;i<510;i++){
b9671574 1350 if (tpcrow->GetFastCluster(i)<0)
1351 tpcrow->SetFastCluster(i,last);
91162307 1352 else
b9671574 1353 last = tpcrow->GetFastCluster(i);
91162307 1354 }
1355 }
1356 fN=fkNOS;
1357 fSectors=fOuterSec;
1358 return 0;
1359}
1360
1361
1362//_____________________________________________________________________________
1363Int_t AliTPCtrackerMI::LoadInnerSectors() {
1364 //-----------------------------------------------------------------
1365 // This function fills inner TPC sectors with clusters.
1366 //-----------------------------------------------------------------
1367 Int_t nrows = fInnerSec->GetNRows();
1368 UInt_t index=0;
1369 for (Int_t sec = 0;sec<fkNIS;sec++)
1370 for (Int_t row = 0;row<nrows;row++){
bd26fa83 1371 AliTPCtrackerRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
91162307 1372 //
1373 //left
b9671574 1374 Int_t ncl = tpcrow->GetN1();
91162307 1375 while (ncl--) {
b9671574 1376 AliTPCclusterMI *c= (tpcrow->GetCluster1(ncl));
91162307 1377 index=(((sec<<8)+row)<<16)+ncl;
1378 tpcrow->InsertCluster(c,index);
1379 }
1380 //right
b9671574 1381 ncl = tpcrow->GetN2();
91162307 1382 while (ncl--) {
b9671574 1383 AliTPCclusterMI *c= (tpcrow->GetCluster2(ncl));
91162307 1384 index=((((sec+fkNIS)<<8)+row)<<16)+ncl;
1385 tpcrow->InsertCluster(c,index);
1386 }
1387 //
1388 // write indexes for fast acces
1389 //
1390 for (Int_t i=0;i<510;i++)
b9671574 1391 tpcrow->SetFastCluster(i,-1);
91162307 1392 for (Int_t i=0;i<tpcrow->GetN();i++){
1393 Int_t zi = Int_t((*tpcrow)[i]->GetZ()+255.);
b9671574 1394 tpcrow->SetFastCluster(zi,i); // write index
91162307 1395 }
1396 Int_t last = 0;
1397 for (Int_t i=0;i<510;i++){
b9671574 1398 if (tpcrow->GetFastCluster(i)<0)
1399 tpcrow->SetFastCluster(i,last);
91162307 1400 else
b9671574 1401 last = tpcrow->GetFastCluster(i);
91162307 1402 }
1c53abe2 1403
91162307 1404 }
1405
1c53abe2 1406 fN=fkNIS;
1407 fSectors=fInnerSec;
91162307 1408 return 0;
1409}
1410
1411
1412
1413//_________________________________________________________________________
1414AliTPCclusterMI *AliTPCtrackerMI::GetClusterMI(Int_t index) const {
1415 //--------------------------------------------------------------------
1416 // Return pointer to a given cluster
1417 //--------------------------------------------------------------------
e7eb17e4 1418 if (index<0) return 0; // no cluster
91162307 1419 Int_t sec=(index&0xff000000)>>24;
1420 Int_t row=(index&0x00ff0000)>>16;
d26d9159 1421 Int_t ncl=(index&0x00007fff)>>00;
91162307 1422
bd26fa83 1423 const AliTPCtrackerRow * tpcrow=0;
91162307 1424 AliTPCclusterMI * clrow =0;
3da363a1 1425
ad23441a 1426 if (sec<0 || sec>=fkNIS*4) {
1427 AliWarning(Form("Wrong sector %d",sec));
1428 return 0x0;
1429 }
1430
91162307 1431 if (sec<fkNIS*2){
1432 tpcrow = &(fInnerSec[sec%fkNIS][row]);
3da363a1 1433 if (tpcrow==0) return 0;
1434
1435 if (sec<fkNIS) {
b9671574 1436 if (tpcrow->GetN1()<=ncl) return 0;
1437 clrow = tpcrow->GetClusters1();
3da363a1 1438 }
1439 else {
b9671574 1440 if (tpcrow->GetN2()<=ncl) return 0;
1441 clrow = tpcrow->GetClusters2();
3da363a1 1442 }
91162307 1443 }
3da363a1 1444 else {
91162307 1445 tpcrow = &(fOuterSec[(sec-fkNIS*2)%fkNOS][row]);
3da363a1 1446 if (tpcrow==0) return 0;
1447
1448 if (sec-2*fkNIS<fkNOS) {
b9671574 1449 if (tpcrow->GetN1()<=ncl) return 0;
1450 clrow = tpcrow->GetClusters1();
3da363a1 1451 }
1452 else {
b9671574 1453 if (tpcrow->GetN2()<=ncl) return 0;
1454 clrow = tpcrow->GetClusters2();
3da363a1 1455 }
91162307 1456 }
3da363a1 1457
91162307 1458 return &(clrow[ncl]);
1459
1c53abe2 1460}
1461
91162307 1462
1463
1c53abe2 1464Int_t AliTPCtrackerMI::FollowToNext(AliTPCseed& t, Int_t nr) {
1465 //-----------------------------------------------------------------
1466 // This function tries to find a track prolongation to next pad row
1467 //-----------------------------------------------------------------
1c53abe2 1468 //
91162307 1469 Double_t x= GetXrow(nr), ymax=GetMaxY(nr);
4d158c36 1470 AliTPCclusterMI *cl=0;
1471 Int_t tpcindex= t.GetClusterIndex2(nr);
1472 //
1473 // update current shape info every 5 pad-row
1474 // if ( (nr%5==0) || t.GetNumberOfClusters()<2 || (t.fCurrentSigmaY2<0.0001) ){
1475 GetShape(&t,nr);
1476 //}
1477 //
1478 if (fIteration>0 && tpcindex>=-1){ //if we have already clusters
1479 //
1480 if (tpcindex==-1) return 0; //track in dead zone
1481 if (tpcindex>0){ //
b9671574 1482 cl = t.GetClusterPointer(nr);
e0656451 1483 if ( (cl==0) ) cl = GetClusterMI(tpcindex);
b9671574 1484 t.SetCurrentClusterIndex1(tpcindex);
4d158c36 1485 }
1486 if (cl){
1487 Int_t relativesector = ((tpcindex&0xff000000)>>24)%18; // if previously accepted cluster in different sector
1488 Float_t angle = relativesector*fSectors->GetAlpha()+fSectors->GetAlphaShift();
1489 //
1490 if (angle<-TMath::Pi()) angle += 2*TMath::Pi();
1491 if (angle>=TMath::Pi()) angle -= 2*TMath::Pi();
1492
1493 if (TMath::Abs(angle-t.GetAlpha())>0.001){
1494 Double_t rotation = angle-t.GetAlpha();
b9671574 1495 t.SetRelativeSector(relativesector);
3f82c4f2 1496 if (!t.Rotate(rotation)) return 0;
4d158c36 1497 }
3f82c4f2 1498 if (!t.PropagateTo(x)) return 0;
4d158c36 1499 //
b9671574 1500 t.SetCurrentCluster(cl);
1501 t.SetRow(nr);
00055a22 1502 Int_t accept = AcceptCluster(&t,t.GetCurrentCluster());
4d158c36 1503 if ((tpcindex&0x8000)==0) accept =0;
1504 if (accept<3) {
1505 //if founded cluster is acceptible
1506 if (cl->IsUsed(11)) { // id cluster is shared inrease uncertainty
b9671574 1507 t.SetErrorY2(t.GetErrorY2()+0.03);
1508 t.SetErrorZ2(t.GetErrorZ2()+0.03);
1509 t.SetErrorY2(t.GetErrorY2()*3);
1510 t.SetErrorZ2(t.GetErrorZ2()*3);
4d158c36 1511 }
b9671574 1512 t.SetNFoundable(t.GetNFoundable()+1);
4d158c36 1513 UpdateTrack(&t,accept);
1514 return 1;
1515 }
1516 }
1627d1c4 1517 }
3f82c4f2 1518 if (TMath::Abs(t.GetSnp())>AliTPCReconstructor::GetMaxSnpTracker()) return 0; // cut on angle
1519 if (fIteration>1){
1520 // not look for new cluster during refitting
b9671574 1521 t.SetNFoundable(t.GetNFoundable()+1);
3f82c4f2 1522 return 0;
1523 }
91162307 1524 //
4d158c36 1525 UInt_t index=0;
ca142b1f 1526 // if (TMath::Abs(t.GetSnp())>0.95 || TMath::Abs(x*t.GetC()-t.GetEta())>0.95) return 0;// patch 28 fev 06
4d158c36 1527 Double_t y=t.GetYat(x);
91162307 1528 if (TMath::Abs(y)>ymax){
1529 if (y > ymax) {
b9671574 1530 t.SetRelativeSector((t.GetRelativeSector()+1) % fN);
91162307 1531 if (!t.Rotate(fSectors->GetAlpha()))
1532 return 0;
1533 } else if (y <-ymax) {
b9671574 1534 t.SetRelativeSector((t.GetRelativeSector()-1+fN) % fN);
91162307 1535 if (!t.Rotate(-fSectors->GetAlpha()))
1536 return 0;
1537 }
4d158c36 1538 //return 1;
91162307 1539 }
1540 //
4d158c36 1541 if (!t.PropagateTo(x)) {
b9671574 1542 if (fIteration==0) t.SetRemoval(10);
4d158c36 1543 return 0;
91162307 1544 }
4d158c36 1545 y=t.GetY();
1546 Double_t z=t.GetZ();
1547 //
a3232aae 1548
b9671574 1549 if (!IsActive(t.GetRelativeSector(),nr)) {
1550 t.SetInDead(kTRUE);
a3232aae 1551 t.SetClusterIndex2(nr,-1);
1552 return 0;
1553 }
1554 //AliInfo(Form("A - Sector%d phi %f - alpha %f", t.fRelativeSector,y/x, t.GetAlpha()));
b9671574 1555 Bool_t isActive = IsActive(t.GetRelativeSector(),nr);
1556 Bool_t isActive2 = (nr>=fInnerSec->GetNRows()) ? fOuterSec[t.GetRelativeSector()][nr-fInnerSec->GetNRows()].GetN()>0:fInnerSec[t.GetRelativeSector()][nr].GetN()>0;
a3232aae 1557
1558 if (!isActive || !isActive2) return 0;
1559
bd26fa83 1560 const AliTPCtrackerRow &krow=GetRow(t.GetRelativeSector(),nr);
91162307 1561 if ( (t.GetSigmaY2()<0) || t.GetSigmaZ2()<0) return 0;
1562 Double_t roady =1.;
1563 Double_t roadz = 1.;
1564 //
b9671574 1565 if (TMath::Abs(TMath::Abs(y)-ymax)<krow.GetDeadZone()){
1566 t.SetInDead(kTRUE);
91162307 1567 t.SetClusterIndex2(nr,-1);
1c53abe2 1568 return 0;
1569 }
1570 else
1571 {
a1ec4d07 1572 if (TMath::Abs(z)<(AliTPCReconstructor::GetCtgRange()*x+10) && TMath::Abs(z)<fParam->GetZLength(0) && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker()))
b9671574 1573 t.SetNFoundable(t.GetNFoundable()+1);
1627d1c4 1574 else
1575 return 0;
1c53abe2 1576 }
1577 //calculate
91162307 1578 if (krow) {
1579 // cl = krow.FindNearest2(y+10.,z,roady,roadz,index);
1580 cl = krow.FindNearest2(y,z,roady,roadz,index);
b9671574 1581 if (cl) t.SetCurrentClusterIndex1(krow.GetIndex(index));
91162307 1582 }
91162307 1583 if (cl) {
b9671574 1584 t.SetCurrentCluster(cl);
1585 t.SetRow(nr);
4d158c36 1586 if (fIteration==2&&cl->IsUsed(10)) return 0;
00055a22 1587 Int_t accept = AcceptCluster(&t,t.GetCurrentCluster());
4d158c36 1588 if (fIteration==2&&cl->IsUsed(11)) {
b9671574 1589 t.SetErrorY2(t.GetErrorY2()+0.03);
1590 t.SetErrorZ2(t.GetErrorZ2()+0.03);
1591 t.SetErrorY2(t.GetErrorY2()*3);
1592 t.SetErrorZ2(t.GetErrorZ2()*3);
4d158c36 1593 }
d26d9159 1594 /*
91162307 1595 if (t.fCurrentCluster->IsUsed(10)){
1596 //
1597 //
c9427e08 1598
91162307 1599 t.fNShared++;
1600 if (t.fNShared>0.7*t.GetNumberOfClusters()) {
1601 t.fRemoval =10;
1602 return 0;
1603 }
1604 }
d26d9159 1605 */
91162307 1606 if (accept<3) UpdateTrack(&t,accept);
c9427e08 1607
91162307 1608 } else {
b9671574 1609 if ( fIteration==0 && t.GetNFoundable()*0.5 > t.GetNumberOfClusters()) t.SetRemoval(10);
91162307 1610
1611 }
1612 return 1;
1613}
c9427e08 1614
91162307 1615Int_t AliTPCtrackerMI::FollowToNextFast(AliTPCseed& t, Int_t nr) {
1616 //-----------------------------------------------------------------
1617 // This function tries to find a track prolongation to next pad row
1618 //-----------------------------------------------------------------
1619 //
1620 Double_t x= GetXrow(nr), ymax=GetMaxY(nr);
1621 Double_t y,z;
1622 if (!t.GetProlongation(x,y,z)) {
b9671574 1623 t.SetRemoval(10);
91162307 1624 return 0;
1625 }
1626 //
1627 //
1628 if (TMath::Abs(y)>ymax){
91162307 1629
1c53abe2 1630 if (y > ymax) {
b9671574 1631 t.SetRelativeSector((t.GetRelativeSector()+1) % fN);
1c53abe2 1632 if (!t.Rotate(fSectors->GetAlpha()))
1633 return 0;
1634 } else if (y <-ymax) {
b9671574 1635 t.SetRelativeSector((t.GetRelativeSector()-1+fN) % fN);
1c53abe2 1636 if (!t.Rotate(-fSectors->GetAlpha()))
1637 return 0;
91162307 1638 }
1639 if (!t.PropagateTo(x)) {
1640 return 0;
1641 }
1642 t.GetProlongation(x,y,z);
1643 }
1644 //
42aec1f1 1645 // update current shape info every 2 pad-row
1646 if ( (nr%2==0) || t.GetNumberOfClusters()<2 || (t.GetCurrentSigmaY2()<0.0001) ){
91162307 1647 // t.fCurrentSigmaY = GetSigmaY(&t);
1648 //t.fCurrentSigmaZ = GetSigmaZ(&t);
1649 GetShape(&t,nr);
1650 }
1651 //
1652 AliTPCclusterMI *cl=0;
1653 UInt_t index=0;
1654
1655
1656 //Int_t nr2 = nr;
bd26fa83 1657 const AliTPCtrackerRow &krow=GetRow(t.GetRelativeSector(),nr);
91162307 1658 if ( (t.GetSigmaY2()<0) || t.GetSigmaZ2()<0) return 0;
1659 Double_t roady =1.;
1660 Double_t roadz = 1.;
1661 //
1662 Int_t row = nr;
b9671574 1663 if (TMath::Abs(TMath::Abs(y)-ymax)<krow.GetDeadZone()){
1664 t.SetInDead(kTRUE);
91162307 1665 t.SetClusterIndex2(row,-1);
1666 return 0;
1667 }
1668 else
1669 {
d7a11555 1670 if (TMath::Abs(z)>(AliTPCReconstructor::GetCtgRange()*x+10)) t.SetClusterIndex2(row,-1);
1c53abe2 1671 }
91162307 1672 //calculate
1673
1674 if ((cl==0)&&(krow)) {
1675 // cl = krow.FindNearest2(y+10,z,roady,roadz,index);
1676 cl = krow.FindNearest2(y,z,roady,roadz,index);
1677
b9671574 1678 if (cl) t.SetCurrentClusterIndex1(krow.GetIndex(index));
91162307 1679 }
1680
1681 if (cl) {
b9671574 1682 t.SetCurrentCluster(cl);
00055a22 1683 // Int_t accept = AcceptCluster(&t,t.fCurrentCluster);
91162307 1684 //if (accept<3){
1685 t.SetClusterIndex2(row,index);
b9671574 1686 t.SetClusterPointer(row, cl);
91162307 1687 //}
1c53abe2 1688 }
1689 return 1;
1690}
1691
1692
91162307 1693
5d837844 1694//_________________________________________________________________________
1695Bool_t AliTPCtrackerMI::GetTrackPoint(Int_t index, AliTrackPoint &p ) const
1696{
1697 // Get track space point by index
1698 // return false in case the cluster doesn't exist
1699 AliTPCclusterMI *cl = GetClusterMI(index);
1700 if (!cl) return kFALSE;
1701 Int_t sector = (index&0xff000000)>>24;
0201b65c 1702 // Int_t row = (index&0x00ff0000)>>16;
5d837844 1703 Float_t xyz[3];
0201b65c 1704 // xyz[0] = fParam->GetPadRowRadii(sector,row);
1705 xyz[0] = cl->GetX();
5d837844 1706 xyz[1] = cl->GetY();
1707 xyz[2] = cl->GetZ();
1708 Float_t sin,cos;
1709 fParam->AdjustCosSin(sector,cos,sin);
1710 Float_t x = cos*xyz[0]-sin*xyz[1];
1711 Float_t y = cos*xyz[1]+sin*xyz[0];
1712 Float_t cov[6];
1713 Float_t sigmaY2 = 0.027*cl->GetSigmaY2();
1714 if (sector < fParam->GetNInnerSector()) sigmaY2 *= 2.07;
1715 Float_t sigmaZ2 = 0.066*cl->GetSigmaZ2();
1716 if (sector < fParam->GetNInnerSector()) sigmaZ2 *= 1.77;
1717 cov[0] = sin*sin*sigmaY2;
1718 cov[1] = -sin*cos*sigmaY2;
1719 cov[2] = 0.;
1720 cov[3] = cos*cos*sigmaY2;
1721 cov[4] = 0.;
1722 cov[5] = sigmaZ2;
1723 p.SetXYZ(x,y,xyz[2],cov);
ae079791 1724 AliGeomManager::ELayerID iLayer;
5d837844 1725 Int_t idet;
1726 if (sector < fParam->GetNInnerSector()) {
ae079791 1727 iLayer = AliGeomManager::kTPC1;
5d837844 1728 idet = sector;
1729 }
1730 else {
ae079791 1731 iLayer = AliGeomManager::kTPC2;
5d837844 1732 idet = sector - fParam->GetNInnerSector();
1733 }
ae079791 1734 UShort_t volid = AliGeomManager::LayerToVolUID(iLayer,idet);
5d837844 1735 p.SetVolumeID(volid);
1736 return kTRUE;
1737}
1738
1739
1740
91162307 1741Int_t AliTPCtrackerMI::UpdateClusters(AliTPCseed& t, Int_t nr) {
1c53abe2 1742 //-----------------------------------------------------------------
1743 // This function tries to find a track prolongation to next pad row
1744 //-----------------------------------------------------------------
b9671574 1745 t.SetCurrentCluster(0);
1746 t.SetCurrentClusterIndex1(0);
1c53abe2 1747
1748 Double_t xt=t.GetX();
91162307 1749 Int_t row = GetRowNumber(xt)-1;
1750 Double_t ymax= GetMaxY(nr);
1751
1c53abe2 1752 if (row < nr) return 1; // don't prolongate if not information until now -
ca142b1f 1753// if (TMath::Abs(t.GetSnp())>0.9 && t.GetNumberOfClusters()>40. && fIteration!=2) {
1754// t.fRemoval =10;
1755// return 0; // not prolongate strongly inclined tracks
1756// }
1757// if (TMath::Abs(t.GetSnp())>0.95) {
1758// t.fRemoval =10;
1759// return 0; // not prolongate strongly inclined tracks
1760// }// patch 28 fev 06
91162307 1761
1762 Double_t x= GetXrow(nr);
1763 Double_t y,z;
1764 //t.PropagateTo(x+0.02);
1765 //t.PropagateTo(x+0.01);
1627d1c4 1766 if (!t.PropagateTo(x)){
1627d1c4 1767 return 0;
1768 }
1c53abe2 1769 //
91162307 1770 y=t.GetY();
1771 z=t.GetZ();
1c53abe2 1772
91162307 1773 if (TMath::Abs(y)>ymax){
1774 if (y > ymax) {
b9671574 1775 t.SetRelativeSector((t.GetRelativeSector()+1) % fN);
91162307 1776 if (!t.Rotate(fSectors->GetAlpha()))
1777 return 0;
1778 } else if (y <-ymax) {
b9671574 1779 t.SetRelativeSector((t.GetRelativeSector()-1+fN) % fN);
91162307 1780 if (!t.Rotate(-fSectors->GetAlpha()))
1781 return 0;
1782 }
982aff31 1783 // if (!t.PropagateTo(x)){
1784 // return 0;
1785 //}
1786 return 1;
1787 //y = t.GetY();
1c53abe2 1788 }
91162307 1789 //
3f82c4f2 1790 if (TMath::Abs(t.GetSnp())>AliTPCReconstructor::GetMaxSnpTracker()) return 0;
a3232aae 1791
b9671574 1792 if (!IsActive(t.GetRelativeSector(),nr)) {
1793 t.SetInDead(kTRUE);
a3232aae 1794 t.SetClusterIndex2(nr,-1);
1795 return 0;
1796 }
1797 //AliInfo(Form("A - Sector%d phi %f - alpha %f", t.fRelativeSector,y/x, t.GetAlpha()));
1798
bd26fa83 1799 AliTPCtrackerRow &krow=GetRow(t.GetRelativeSector(),nr);
1c53abe2 1800
b9671574 1801 if (TMath::Abs(TMath::Abs(y)-ymax)<krow.GetDeadZone()){
1802 t.SetInDead(kTRUE);
91162307 1803 t.SetClusterIndex2(nr,-1);
1c53abe2 1804 return 0;
1805 }
1806 else
1807 {
b9671574 1808 if (TMath::Abs(t.GetZ())<(AliTPCReconstructor::GetCtgRange()*t.GetX()+10) && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker())) t.SetNFoundable(t.GetNFoundable()+1);
1627d1c4 1809 else
1810 return 0;
1c53abe2 1811 }
1c53abe2 1812
91162307 1813 // update current
42aec1f1 1814 if ( (nr%2==0) || t.GetNumberOfClusters()<2){
91162307 1815 // t.fCurrentSigmaY = GetSigmaY(&t);
1816 //t.fCurrentSigmaZ = GetSigmaZ(&t);
1817 GetShape(&t,nr);
1818 }
1c53abe2 1819
91162307 1820 AliTPCclusterMI *cl=0;
d184f295 1821 Int_t index=0;
91162307 1822 //
1823 Double_t roady = 1.;
1824 Double_t roadz = 1.;
1825 //
d26d9159 1826
1827 if (!cl){
1828 index = t.GetClusterIndex2(nr);
1829 if ( (index>0) && (index&0x8000)==0){
b9671574 1830 cl = t.GetClusterPointer(nr);
d26d9159 1831 if ( (cl==0) && (index>0)) cl = GetClusterMI(index);
b9671574 1832 t.SetCurrentClusterIndex1(index);
d26d9159 1833 if (cl) {
b9671574 1834 t.SetCurrentCluster(cl);
d26d9159 1835 return 1;
1836 }
1837 }
1838 }
1839
3d172d79 1840 // if (index<0) return 0;
1841 UInt_t uindex = TMath::Abs(index);
d184f295 1842
91162307 1843 if (krow) {
d184f295 1844 //cl = krow.FindNearest2(y+10,z,roady,roadz,uindex);
1845 cl = krow.FindNearest2(y,z,roady,roadz,uindex);
91162307 1846 }
d26d9159 1847
b9671574 1848 if (cl) t.SetCurrentClusterIndex1(krow.GetIndex(uindex));
1849 t.SetCurrentCluster(cl);
d26d9159 1850
91162307 1851 return 1;
1852}
1c53abe2 1853
1c53abe2 1854
91162307 1855Int_t AliTPCtrackerMI::FollowToNextCluster(AliTPCseed & t, Int_t nr) {
1856 //-----------------------------------------------------------------
1857 // This function tries to find a track prolongation to next pad row
1858 //-----------------------------------------------------------------
1c53abe2 1859
91162307 1860 //update error according neighborhoud
1c53abe2 1861
b9671574 1862 if (t.GetCurrentCluster()) {
1863 t.SetRow(nr);
00055a22 1864 Int_t accept = AcceptCluster(&t,t.GetCurrentCluster());
91162307 1865
b9671574 1866 if (t.GetCurrentCluster()->IsUsed(10)){
91162307 1867 //
1868 //
1869 // t.fErrorZ2*=2;
1870 // t.fErrorY2*=2;
b9671574 1871 t.SetNShared(t.GetNShared()+1);
1872 if (t.GetNShared()>0.7*t.GetNumberOfClusters()) {
1873 t.SetRemoval(10);
91162307 1874 return 0;
1875 }
b364ca79 1876 }
d26d9159 1877 if (fIteration>0) accept = 0;
b364ca79 1878 if (accept<3) UpdateTrack(&t,accept);
1879
1c53abe2 1880 } else {
91162307 1881 if (fIteration==0){
b9671574 1882 if ( ( (t.GetSigmaY2()+t.GetSigmaZ2())>0.16)&& t.GetNumberOfClusters()>18) t.SetRemoval(10);
1883 if ( t.GetChi2()/t.GetNumberOfClusters()>6 &&t.GetNumberOfClusters()>18) t.SetRemoval(10);
91162307 1884
b9671574 1885 if (( (t.GetNFoundable()*0.5 > t.GetNumberOfClusters()) || t.GetNoCluster()>15)) t.SetRemoval(10);
1c53abe2 1886 }
1887 }
1888 return 1;
1889}
1890
1891
1892
91162307 1893//_____________________________________________________________________________
1894Int_t AliTPCtrackerMI::FollowProlongation(AliTPCseed& t, Int_t rf, Int_t step) {
1c53abe2 1895 //-----------------------------------------------------------------
91162307 1896 // This function tries to find a track prolongation.
1c53abe2 1897 //-----------------------------------------------------------------
1898 Double_t xt=t.GetX();
1899 //
91162307 1900 Double_t alpha=t.GetAlpha() - fSectors->GetAlphaShift();
1c53abe2 1901 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
1902 if (alpha < 0. ) alpha += 2.*TMath::Pi();
91162307 1903 //
b9671574 1904 t.SetRelativeSector(Int_t(alpha/fSectors->GetAlpha()+0.0001)%fN);
1c53abe2 1905
91162307 1906 Int_t first = GetRowNumber(xt)-1;
51ad6848 1907 for (Int_t nr= first; nr>=rf; nr-=step) {
1908 // update kink info
1909 if (t.GetKinkIndexes()[0]>0){
1910 for (Int_t i=0;i<3;i++){
1911 Int_t index = t.GetKinkIndexes()[i];
1912 if (index==0) break;
1913 if (index<0) continue;
1914 //
6c94f330 1915 AliKink * kink = (AliKink*)fEvent->GetKink(index-1);
51ad6848 1916 if (!kink){
1917 printf("PROBLEM\n");
1918 }
1919 else{
eea478d3 1920 Int_t kinkrow = kink->GetTPCRow0()+2+Int_t(0.5/(0.05+kink->GetAngle(2)));
51ad6848 1921 if (kinkrow==nr){
1922 AliExternalTrackParam paramd(t);
1923 kink->SetDaughter(paramd);
eea478d3 1924 kink->SetStatus(2,5);
51ad6848 1925 kink->Update();
51ad6848 1926 }
1927 }
1928 }
1929 }
1930
1931 if (nr==80) t.UpdateReference();
982aff31 1932 if (nr<fInnerSec->GetNRows())
1933 fSectors = fInnerSec;
1934 else
1935 fSectors = fOuterSec;
91162307 1936 if (FollowToNext(t,nr)==0)
4d158c36 1937 if (!t.IsActive())
1938 return 0;
91162307 1939
1940 }
1941 return 1;
1942}
1943
1c53abe2 1944
1945//_____________________________________________________________________________
91162307 1946Int_t AliTPCtrackerMI::FollowProlongationFast(AliTPCseed& t, Int_t rf, Int_t step) {
1c53abe2 1947 //-----------------------------------------------------------------
1948 // This function tries to find a track prolongation.
1949 //-----------------------------------------------------------------
1950 Double_t xt=t.GetX();
1951 //
1952 Double_t alpha=t.GetAlpha() - fSectors->GetAlphaShift();
1953 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
1954 if (alpha < 0. ) alpha += 2.*TMath::Pi();
b9671574 1955 t.SetRelativeSector(Int_t(alpha/fSectors->GetAlpha()+0.0001)%fN);
91162307 1956
1957 for (Int_t nr=GetRowNumber(xt)-1; nr>=rf; nr-=step) {
1958
1959 if (FollowToNextFast(t,nr)==0)
1960 if (!t.IsActive()) return 0;
1c53abe2 1961
1c53abe2 1962 }
1963 return 1;
1964}
1965
1966
91162307 1967
1968
1969
1c53abe2 1970Int_t AliTPCtrackerMI::FollowBackProlongation(AliTPCseed& t, Int_t rf) {
1971 //-----------------------------------------------------------------
1972 // This function tries to find a track prolongation.
1973 //-----------------------------------------------------------------
1c53abe2 1974 //
eea478d3 1975 Double_t xt=t.GetX();
1c53abe2 1976 Double_t alpha=t.GetAlpha() - fSectors->GetAlphaShift();
1977 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
1978 if (alpha < 0. ) alpha += 2.*TMath::Pi();
b9671574 1979 t.SetRelativeSector(Int_t(alpha/fSectors->GetAlpha()+0.0001)%fN);
1c53abe2 1980
b9671574 1981 Int_t first = t.GetFirstPoint();
eea478d3 1982 if (first<GetRowNumber(xt)+1) first = GetRowNumber(xt)+1;
91162307 1983 //
1984 if (first<0) first=0;
4d158c36 1985 for (Int_t nr=first; nr<=rf; nr++) {
ca142b1f 1986 // if ( (TMath::Abs(t.GetSnp())>0.95)) break;//patch 28 fev 06
51ad6848 1987 if (t.GetKinkIndexes()[0]<0){
1988 for (Int_t i=0;i<3;i++){
1989 Int_t index = t.GetKinkIndexes()[i];
1990 if (index==0) break;
1991 if (index>0) continue;
1992 index = TMath::Abs(index);
6c94f330 1993 AliKink * kink = (AliKink*)fEvent->GetKink(index-1);
51ad6848 1994 if (!kink){
1995 printf("PROBLEM\n");
1996 }
1997 else{
eea478d3 1998 Int_t kinkrow = kink->GetTPCRow0()-2-Int_t(0.5/(0.05+kink->GetAngle(2)));
51ad6848 1999 if (kinkrow==nr){
2000 AliExternalTrackParam paramm(t);
2001 kink->SetMother(paramm);
eea478d3 2002 kink->SetStatus(2,1);
51ad6848 2003 kink->Update();
51ad6848 2004 }
2005 }
eea478d3 2006 }
51ad6848 2007 }
eea478d3 2008 //
d26d9159 2009 if (nr<fInnerSec->GetNRows())
2010 fSectors = fInnerSec;
2011 else
2012 fSectors = fOuterSec;
c274e255 2013
d26d9159 2014 FollowToNext(t,nr);
1c53abe2 2015 }
2016 return 1;
2017}
2018
2019
2020
2021
2022
2023Float_t AliTPCtrackerMI::OverlapFactor(AliTPCseed * s1, AliTPCseed * s2, Int_t &sum1, Int_t & sum2)
2024{
2025 //
2026 //
2027 sum1=0;
2028 sum2=0;
2029 Int_t sum=0;
1c53abe2 2030 //
2031 Float_t dz2 =(s1->GetZ() - s2->GetZ());
c9427e08 2032 dz2*=dz2;
91162307 2033
c9427e08 2034 Float_t dy2 =TMath::Abs((s1->GetY() - s2->GetY()));
1c53abe2 2035 dy2*=dy2;
2036 Float_t distance = TMath::Sqrt(dz2+dy2);
c9427e08 2037 if (distance>4.) return 0; // if there are far away - not overlap - to reduce combinatorics
1c53abe2 2038
91162307 2039 // Int_t offset =0;
b9671574 2040 Int_t firstpoint = TMath::Min(s1->GetFirstPoint(),s2->GetFirstPoint());
2041 Int_t lastpoint = TMath::Max(s1->GetLastPoint(),s2->GetLastPoint());
c9427e08 2042 if (lastpoint>160)
2043 lastpoint =160;
2044 if (firstpoint<0)
2045 firstpoint = 0;
91162307 2046 if (firstpoint>lastpoint) {
2047 firstpoint =lastpoint;
2048 // lastpoint =160;
c9427e08 2049 }
2050
2051
91162307 2052 for (Int_t i=firstpoint-1;i<lastpoint+1;i++){
2053 if (s1->GetClusterIndex2(i)>0) sum1++;
2054 if (s2->GetClusterIndex2(i)>0) sum2++;
2055 if (s1->GetClusterIndex2(i)==s2->GetClusterIndex2(i) && s1->GetClusterIndex2(i)>0) {
1c53abe2 2056 sum++;
2057 }
2058 }
91162307 2059 if (sum<5) return 0;
2060
1627d1c4 2061 Float_t summin = TMath::Min(sum1+1,sum2+1);
2062 Float_t ratio = (sum+1)/Float_t(summin);
1c53abe2 2063 return ratio;
2064}
2065
2066void AliTPCtrackerMI::SignShared(AliTPCseed * s1, AliTPCseed * s2)
2067{
2068 //
2069 //
a0f4d6a6 2070 Float_t thetaCut = 0.2;//+10.*TMath::Sqrt(s1->GetSigmaTglZ()+ s2->GetSigmaTglZ());
2071 if (TMath::Abs(s1->GetTgl()-s2->GetTgl())>thetaCut) return;
2072 Float_t minCl = TMath::Min(s1->GetNumberOfClusters(),s2->GetNumberOfClusters());
2073 Int_t cutN0 = TMath::Max(5,TMath::Nint(0.1*minCl));
91162307 2074
1c53abe2 2075 //
91162307 2076 Int_t sumshared=0;
2077 //
a0f4d6a6 2078 //Int_t firstpoint = TMath::Max(s1->GetFirstPoint(),s2->GetFirstPoint());
2079 //Int_t lastpoint = TMath::Min(s1->GetLastPoint(),s2->GetLastPoint());
2080 Int_t firstpoint = 0;
2081 Int_t lastpoint = 160;
91162307 2082 //
a0f4d6a6 2083 // if (firstpoint>=lastpoint-5) return;;
1c53abe2 2084
91162307 2085 for (Int_t i=firstpoint;i<lastpoint;i++){
2086 // if ( (s1->GetClusterIndex2(i)&0xFFFF8FFF)==(s2->GetClusterIndex2(i)&0xFFFF8FFF) && s1->GetClusterIndex2(i)>0) {
2087 if ( (s1->GetClusterIndex2(i))==(s2->GetClusterIndex2(i)) && s1->GetClusterIndex2(i)>0) {
2088 sumshared++;
19b00333 2089 s1->SetSharedMapBit(i, kTRUE);
2090 s2->SetSharedMapBit(i, kTRUE);
91162307 2091 }
19b00333 2092 if (s1->GetClusterIndex2(i)>0)
2093 s1->SetClusterMapBit(i, kTRUE);
91162307 2094 }
a0f4d6a6 2095 if (sumshared>cutN0){
91162307 2096 // sign clusters
2097 //
2098 for (Int_t i=firstpoint;i<lastpoint;i++){
2099 // if ( (s1->GetClusterIndex2(i)&0xFFFF8FFF)==(s2->GetClusterIndex2(i)&0xFFFF8FFF) && s1->GetClusterIndex2(i)>0) {
2100 if ( (s1->GetClusterIndex2(i))==(s2->GetClusterIndex2(i)) && s1->GetClusterIndex2(i)>0) {
2101 AliTPCTrackerPoint *p1 = s1->GetTrackPoint(i);
2102 AliTPCTrackerPoint *p2 = s2->GetTrackPoint(i);;
2103 if (s1->IsActive()&&s2->IsActive()){
b9671574 2104 p1->SetShared(kTRUE);
2105 p2->SetShared(kTRUE);
91162307 2106 }
2107 }
2108 }
2109 }
2110 //
a0f4d6a6 2111 if (sumshared>cutN0){
91162307 2112 for (Int_t i=0;i<4;i++){
b9671574 2113 if (s1->GetOverlapLabel(3*i)==0){
2114 s1->SetOverlapLabel(3*i, s2->GetLabel());
2115 s1->SetOverlapLabel(3*i+1,sumshared);
2116 s1->SetOverlapLabel(3*i+2,s2->GetUniqueID());
91162307 2117 break;
2118 }
2119 }
2120 for (Int_t i=0;i<4;i++){
b9671574 2121 if (s2->GetOverlapLabel(3*i)==0){
2122 s2->SetOverlapLabel(3*i, s1->GetLabel());
2123 s2->SetOverlapLabel(3*i+1,sumshared);
2124 s2->SetOverlapLabel(3*i+2,s1->GetUniqueID());
91162307 2125 break;
2126 }
2127 }
2128 }
91162307 2129}
1c53abe2 2130
91162307 2131void AliTPCtrackerMI::SignShared(TObjArray * arr)
2132{
1c53abe2 2133 //
91162307 2134 //sort trackss according sectors
2135 //
c9427e08 2136 for (Int_t i=0; i<arr->GetEntriesFast(); i++) {
2137 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
91162307 2138 if (!pt) continue;
2139 //if (pt) RotateToLocal(pt);
b9671574 2140 pt->SetSort(0);
c9427e08 2141 }
91162307 2142 arr->UnSort();
6d493ea0 2143 arr->Sort(); // sorting according relative sectors
1c53abe2 2144 arr->Expand(arr->GetEntries());
91162307 2145 //
2146 //
1c53abe2 2147 Int_t nseed=arr->GetEntriesFast();
1c53abe2 2148 for (Int_t i=0; i<nseed; i++) {
2149 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
91162307 2150 if (!pt) continue;
2151 for (Int_t j=0;j<=12;j++){
b9671574 2152 pt->SetOverlapLabel(j,0);
1c53abe2 2153 }
91162307 2154 }
2155 for (Int_t i=0; i<nseed; i++) {
2156 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
2157 if (!pt) continue;
b9671574 2158 if (pt->GetRemoval()>10) continue;
1c53abe2 2159 for (Int_t j=i+1; j<nseed; j++){
2160 AliTPCseed *pt2=(AliTPCseed*)arr->UncheckedAt(j);
6d493ea0 2161 if (TMath::Abs(pt->GetRelativeSector()-pt2->GetRelativeSector())>1) continue;
91162307 2162 // if (pt2){
b9671574 2163 if (pt2->GetRemoval()<=10) {
6d493ea0 2164 //if ( TMath::Abs(pt->GetRelativeSector()-pt2->GetRelativeSector())>0) break;
91162307 2165 SignShared(pt,pt2);
c9427e08 2166 }
91162307 2167 }
2168 }
2169}
2170
91162307 2171
47966a6d 2172void AliTPCtrackerMI::SortTracks(TObjArray * arr, Int_t mode) const
91162307 2173{
2174 //
2175 //sort tracks in array according mode criteria
2176 Int_t nseed = arr->GetEntriesFast();
2177 for (Int_t i=0; i<nseed; i++) {
2178 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
2179 if (!pt) {
2180 continue;
2181 }
b9671574 2182 pt->SetSort(mode);
91162307 2183 }
2184 arr->UnSort();
2185 arr->Sort();
1c53abe2 2186}
c9427e08 2187
51ad6848 2188
2189void AliTPCtrackerMI::RemoveUsed2(TObjArray * arr, Float_t factor1, Float_t factor2, Int_t minimal)
2190{
51ad6848 2191 //
6d493ea0 2192 // Loop over all tracks and remove overlaped tracks (with lower quality)
2193 // Algorithm:
2194 // 1. Unsign clusters
2195 // 2. Sort tracks according quality
2196 // Quality is defined by the number of cluster between first and last points
2197 //
2198 // 3. Loop over tracks - decreasing quality order
2199 // a.) remove - If the fraction of shared cluster less than factor (1- n or 2)
2200 // b.) remove - If the minimal number of clusters less than minimal and not ITS
2201 // c.) if track accepted - sign clusters
51ad6848 2202 //
6d493ea0 2203 //Called in - AliTPCtrackerMI::Clusters2Tracks()
2204 // - AliTPCtrackerMI::PropagateBack()
2205 // - AliTPCtrackerMI::RefitInward()
2206 //
2207
2208
51ad6848 2209 UnsignClusters();
2210 //
2211 Int_t nseed = arr->GetEntriesFast();
2212 Float_t * quality = new Float_t[nseed];
2213 Int_t * indexes = new Int_t[nseed];
2214 Int_t good =0;
2215 //
2216 //
2217 for (Int_t i=0; i<nseed; i++) {
2218 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
2219 if (!pt){
2220 quality[i]=-1;
2221 continue;
2222 }
2223 pt->UpdatePoints(); //select first last max dens points
2224 Float_t * points = pt->GetPoints();
2225 if (points[3]<0.8) quality[i] =-1;
51ad6848 2226 quality[i] = (points[2]-points[0])+pt->GetNumberOfClusters();
6d493ea0 2227 //prefer high momenta tracks if overlaps
2228 quality[i] *= TMath::Sqrt(TMath::Abs(pt->Pt())+0.5);
51ad6848 2229 }
2230 TMath::Sort(nseed,quality,indexes);
2231 //
2232 //
2233 for (Int_t itrack=0; itrack<nseed; itrack++) {
2234 Int_t trackindex = indexes[itrack];
6d493ea0 2235 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(trackindex);
2236 if (!pt) continue;
2237 //
51ad6848 2238 if (quality[trackindex]<0){
2239 if (pt) {
2240 delete arr->RemoveAt(trackindex);
2241 }
2242 else{
2243 arr->RemoveAt(trackindex);
2244 }
2245 continue;
2246 }
2247 //
6d493ea0 2248 //
51ad6848 2249 Int_t first = Int_t(pt->GetPoints()[0]);
2250 Int_t last = Int_t(pt->GetPoints()[2]);
b9671574 2251 Double_t factor = (pt->GetBConstrain()) ? factor1: factor2;
51ad6848 2252 //
2253 Int_t found,foundable,shared;
6d493ea0 2254 pt->GetClusterStatistic(first,last, found, foundable,shared,kFALSE); // better to get statistic in "high-dens" region do't use full track as in line bellow
2255 // pt->GetClusterStatistic(0,160, found, foundable,shared,kFALSE);
b406fdb0 2256 Bool_t itsgold =kFALSE;
b9671574 2257 if (pt->GetESD()){
ef7253ac 2258 Int_t dummy[12];
b9671574 2259 if (pt->GetESD()->GetITSclusters(dummy)>4) itsgold= kTRUE;
51ad6848 2260 }
b406fdb0 2261 if (!itsgold){
2262 //
2263 if (Float_t(shared+1)/Float_t(found+1)>factor){
2264 if (pt->GetKinkIndexes()[0]!=0) continue; //don't remove tracks - part of the kinks
2265 delete arr->RemoveAt(trackindex);
2266 continue;
2267 }
2268 if (pt->GetNumberOfClusters()<50&&(found-0.5*shared)<minimal){ //remove short tracks
2269 if (pt->GetKinkIndexes()[0]!=0) continue; //don't remove tracks - part of the kinks
2270 delete arr->RemoveAt(trackindex);
2271 continue;
2272 }
51ad6848 2273 }
2274
2275 good++;
6d493ea0 2276 //if (sharedfactor>0.4) continue;
b406fdb0 2277 if (pt->GetKinkIndexes()[0]>0) continue;
6d493ea0 2278 //Remove tracks with undefined properties - seems
2279 if (pt->GetSigmaY2()<kAlmost0) continue; // ? what is the origin ?
2280 //
51ad6848 2281 for (Int_t i=first; i<last; i++) {
2282 Int_t index=pt->GetClusterIndex2(i);
2283 // if (index<0 || index&0x8000 ) continue;
2284 if (index<0 || index&0x8000 ) continue;
b9671574 2285 AliTPCclusterMI *c= pt->GetClusterPointer(i);
51ad6848 2286 if (!c) continue;
2287 c->Use(10);
2288 }
2289 }
2290 fNtracks = good;
2291 if (fDebug>0){
2292 Info("RemoveUsed","\n*****\nNumber of good tracks after shared removal\t%d\n",fNtracks);
2293 }
2294 delete []quality;
2295 delete []indexes;
2296}
2297
91162307 2298void AliTPCtrackerMI::UnsignClusters()
1c53abe2 2299{
91162307 2300 //
2301 // loop over all clusters and unsign them
2302 //
2303
2304 for (Int_t sec=0;sec<fkNIS;sec++){
2305 for (Int_t row=0;row<fInnerSec->GetNRows();row++){
b9671574 2306 AliTPCclusterMI *cl = fInnerSec[sec][row].GetClusters1();
2307 for (Int_t icl =0;icl< fInnerSec[sec][row].GetN1();icl++)
91162307 2308 // if (cl[icl].IsUsed(10))
2309 cl[icl].Use(-1);
b9671574 2310 cl = fInnerSec[sec][row].GetClusters2();
2311 for (Int_t icl =0;icl< fInnerSec[sec][row].GetN2();icl++)
91162307 2312 //if (cl[icl].IsUsed(10))
2313 cl[icl].Use(-1);
2314 }
2315 }
2316
2317 for (Int_t sec=0;sec<fkNOS;sec++){
2318 for (Int_t row=0;row<fOuterSec->GetNRows();row++){
b9671574 2319 AliTPCclusterMI *cl = fOuterSec[sec][row].GetClusters1();
2320 for (Int_t icl =0;icl< fOuterSec[sec][row].GetN1();icl++)
91162307 2321 //if (cl[icl].IsUsed(10))
2322 cl[icl].Use(-1);
b9671574 2323 cl = fOuterSec[sec][row].GetClusters2();
2324 for (Int_t icl =0;icl< fOuterSec[sec][row].GetN2();icl++)
91162307 2325 //if (cl[icl].IsUsed(10))
2326 cl[icl].Use(-1);
2327 }
2328 }
2329
1c53abe2 2330}
2331
91162307 2332
2333
732b9e78 2334void AliTPCtrackerMI::SignClusters(TObjArray * arr, Float_t fnumber, Float_t fdensity)
1c53abe2 2335{
2336 //
91162307 2337 //sign clusters to be "used"
2338 //
2339 // snumber and sdensity sign number of sigmas - bellow mean value to be accepted
2340 // loop over "primaries"
2341
2342 Float_t sumdens=0;
2343 Float_t sumdens2=0;
2344 Float_t sumn =0;
2345 Float_t sumn2 =0;
2346 Float_t sumchi =0;
2347 Float_t sumchi2 =0;
2348
2349 Float_t sum =0;
2350
1c53abe2 2351 TStopwatch timer;
91162307 2352 timer.Start();
1c53abe2 2353
91162307 2354 Int_t nseed = arr->GetEntriesFast();
2355 for (Int_t i=0; i<nseed; i++) {
2356 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
2357 if (!pt) {
2358 continue;
2359 }
2360 if (!(pt->IsActive())) continue;
b9671574 2361 Float_t dens = pt->GetNumberOfClusters()/Float_t(pt->GetNFoundable());
91162307 2362 if ( (dens>0.7) && (pt->GetNumberOfClusters()>70)){
2363 sumdens += dens;
2364 sumdens2+= dens*dens;
2365 sumn += pt->GetNumberOfClusters();
2366 sumn2 += pt->GetNumberOfClusters()*pt->GetNumberOfClusters();
2367 Float_t chi2 = pt->GetChi2()/pt->GetNumberOfClusters();
2368 if (chi2>5) chi2=5;
2369 sumchi +=chi2;
2370 sumchi2 +=chi2*chi2;
2371 sum++;
2372 }
1627d1c4 2373 }
91162307 2374
2375 Float_t mdensity = 0.9;
2376 Float_t meann = 130;
2377 Float_t meanchi = 1;
2378 Float_t sdensity = 0.1;
2379 Float_t smeann = 10;
2380 Float_t smeanchi =0.4;
1627d1c4 2381
91162307 2382
2383 if (sum>20){
2384 mdensity = sumdens/sum;
2385 meann = sumn/sum;
2386 meanchi = sumchi/sum;
2387 //
2388 sdensity = sumdens2/sum-mdensity*mdensity;
c1ea348f 2389 if (sdensity >= 0)
2390 sdensity = TMath::Sqrt(sdensity);
2391 else
2392 sdensity = 0.1;
91162307 2393 //
2394 smeann = sumn2/sum-meann*meann;
c1ea348f 2395 if (smeann >= 0)
2396 smeann = TMath::Sqrt(smeann);
2397 else
2398 smeann = 10;
91162307 2399 //
2400 smeanchi = sumchi2/sum - meanchi*meanchi;
c1ea348f 2401 if (smeanchi >= 0)
2402 smeanchi = TMath::Sqrt(smeanchi);
2403 else
2404 smeanchi = 0.4;
91162307 2405 }
2406
2407
2408 //REMOVE SHORT DELTAS or tracks going out of sensitive volume of TPC
2409 //
2410 for (Int_t i=0; i<nseed; i++) {
2411 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
2412 if (!pt) {
2413 continue;
1c53abe2 2414 }
b9671574 2415 if (pt->GetBSigned()) continue;
2416 if (pt->GetBConstrain()) continue;
91162307 2417 //if (!(pt->IsActive())) continue;
2418 /*
2419 Int_t found,foundable,shared;
2420 pt->GetClusterStatistic(0,160,found, foundable,shared);
2421 if (shared/float(found)>0.3) {
2422 if (shared/float(found)>0.9 ){
2423 //delete arr->RemoveAt(i);
2424 }
2425 continue;
c9427e08 2426 }
91162307 2427 */
2428 Bool_t isok =kFALSE;
b9671574 2429 if ( (pt->GetNShared()/pt->GetNumberOfClusters()<0.5) &&pt->GetNumberOfClusters()>60)
91162307 2430 isok = kTRUE;
b9671574 2431 if ((TMath::Abs(1/pt->GetC())<100.) && (pt->GetNShared()/pt->GetNumberOfClusters()<0.7))
91162307 2432 isok =kTRUE;
2433 if (TMath::Abs(pt->GetZ()/pt->GetX())>1.1)
2434 isok =kTRUE;
2435 if ( (TMath::Abs(pt->GetSnp()>0.7) && pt->GetD(0,0)>60.))
2436 isok =kTRUE;
2437
2438 if (isok)
2439 for (Int_t i=0; i<160; i++) {
2440 Int_t index=pt->GetClusterIndex2(i);
2441 if (index<0) continue;
b9671574 2442 AliTPCclusterMI *c= pt->GetClusterPointer(i);
91162307 2443 if (!c) continue;
2444 //if (!(c->IsUsed(10))) c->Use();
2445 c->Use(10);
2446 }
2447 }
2448
c9427e08 2449
1c53abe2 2450 //
91162307 2451 Double_t maxchi = meanchi+2.*smeanchi;
2452
2453 for (Int_t i=0; i<nseed; i++) {
2454 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
2455 if (!pt) {
1c53abe2 2456 continue;
91162307 2457 }
2458 //if (!(pt->IsActive())) continue;
b9671574 2459 if (pt->GetBSigned()) continue;
91162307 2460 Double_t chi = pt->GetChi2()/pt->GetNumberOfClusters();
2461 if (chi>maxchi) continue;
2462
2463 Float_t bfactor=1;
b9671574 2464 Float_t dens = pt->GetNumberOfClusters()/Float_t(pt->GetNFoundable());
91162307 2465
2466 //sign only tracks with enoug big density at the beginning
2467
2468 if ((pt->GetDensityFirst(40)<0.75) && pt->GetNumberOfClusters()<meann) continue;
2469
2470
3e5d0aa2 2471 Double_t mindens = TMath::Max(double(mdensity-sdensity*fdensity*bfactor),0.65);
91162307 2472 Double_t minn = TMath::Max(Int_t(meann-fnumber*smeann*bfactor),50);
2473
2474 // if (pt->fBConstrain) mindens = TMath::Max(mdensity-sdensity*fdensity*bfactor,0.65);
b9671574 2475 if ( (pt->GetRemoval()==10) && (pt->GetSnp()>0.8)&&(dens>mindens))
91162307 2476 minn=0;
2477
2478 if ((dens>mindens && pt->GetNumberOfClusters()>minn) && chi<maxchi ){
2479 //Int_t noc=pt->GetNumberOfClusters();
b9671574 2480 pt->SetBSigned(kTRUE);
91162307 2481 for (Int_t i=0; i<160; i++) {
2482
2483 Int_t index=pt->GetClusterIndex2(i);
2484 if (index<0) continue;
b9671574 2485 AliTPCclusterMI *c= pt->GetClusterPointer(i);
91162307 2486 if (!c) continue;
2487 // if (!(c->IsUsed(10))) c->Use();
2488 c->Use(10);
2489 }
1c53abe2 2490 }
91162307 2491 }
2492 // gLastCheck = nseed;
2493 // arr->Compress();
2494 if (fDebug>0){
2495 timer.Print();
2496 }
1c53abe2 2497}
2498
2499
47966a6d 2500void AliTPCtrackerMI::StopNotActive(TObjArray * arr, Int_t row0, Float_t th0, Float_t th1, Float_t th2) const
91162307 2501{
2502 // stop not active tracks
2503 // take th1 as threshold for number of founded to number of foundable on last 10 active rows
2504 // take th2 as threshold for number of founded to number of foundable on last 20 active rows
2505 Int_t nseed = arr->GetEntriesFast();
2506 //
2507 for (Int_t i=0; i<nseed; i++) {
2508 AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
2509 if (!pt) {
2510 continue;
2511 }
2512 if (!(pt->IsActive())) continue;
2513 StopNotActive(pt,row0,th0, th1,th2);
2514 }
2515}
1c53abe2 2516
1c53abe2 2517
1c53abe2 2518
91162307 2519void AliTPCtrackerMI::StopNotActive(AliTPCseed * seed, Int_t row0, Float_t th0, Float_t th1,
47966a6d 2520 Float_t th2) const
91162307 2521{
2522 // stop not active tracks
2523 // take th1 as threshold for number of founded to number of foundable on last 10 active rows
2524 // take th2 as threshold for number of founded to number of foundable on last 20 active rows
2525 Int_t sumgood1 = 0;
2526 Int_t sumgood2 = 0;
2527 Int_t foundable = 0;
b9671574 2528 Int_t maxindex = seed->GetLastPoint(); //last foundable row
2529 if (seed->GetNFoundable()*th0 > seed->GetNumberOfClusters()) {
91162307 2530 seed->Desactivate(10) ;
2531 return;
2532 }
1c53abe2 2533
3e5d0aa2 2534 for (Int_t i=row0; i<maxindex; i++){
91162307 2535 Int_t index = seed->GetClusterIndex2(i);
2536 if (index!=-1) foundable++;
2537 //if (!c) continue;
2538 if (foundable<=30) sumgood1++;
2539 if (foundable<=50) {
2540 sumgood2++;
2541 }
2542 else{
2543 break;
2544 }
2545 }
2546 if (foundable>=30.){
2547 if (sumgood1<(th1*30.)) seed->Desactivate(10);
2548 }
2549 if (foundable>=50)
2550 if (sumgood2<(th2*50.)) seed->Desactivate(10);
2551}
1c53abe2 2552
1c53abe2 2553
af885e0f 2554Int_t AliTPCtrackerMI::RefitInward(AliESDEvent *event)
d26d9159 2555{
2556 //
2557 // back propagation of ESD tracks
2558 //
2559 //return 0;
6d493ea0 2560 const Int_t kMaxFriendTracks=2000;
d26d9159 2561 fEvent = event;
2562 ReadSeeds(event,2);
2563 fIteration=2;
982aff31 2564 //PrepareForProlongation(fSeeds,1);
2565 PropagateForward2(fSeeds);
6d493ea0 2566 RemoveUsed2(fSeeds,0.4,0.4,20);
2567
a0f4d6a6 2568 TObjArray arraySeed(fSeeds->GetEntries());
2569 for (Int_t i=0;i<fSeeds->GetEntries();i++) {
2570 arraySeed.AddAt(fSeeds->At(i),i);
2571 }
2572 SignShared(&arraySeed);
6d493ea0 2573 // FindCurling(fSeeds, event,2); // find multi found tracks
2574 FindSplitted(fSeeds, event,2); // find multi found tracks
af32720d 2575
4d158c36 2576 Int_t ntracks=0;
d26d9159 2577 Int_t nseed = fSeeds->GetEntriesFast();
2578 for (Int_t i=0;i<nseed;i++){
2579 AliTPCseed * seed = (AliTPCseed*) fSeeds->UncheckedAt(i);
4d158c36 2580 if (!seed) continue;
eea478d3 2581 if (seed->GetKinkIndex(0)>0) UpdateKinkQualityD(seed); // update quality informations for kinks
2582
d26d9159 2583 seed->PropagateTo(fParam->GetInnerRadiusLow());
51ad6848 2584 seed->UpdatePoints();
92f513f5 2585 AddCovariance(seed);
19b00333 2586 MakeBitmaps(seed);
d26d9159 2587 AliESDtrack *esd=event->GetTrack(i);
2588 seed->CookdEdx(0.02,0.6);
2589 CookLabel(seed,0.1); //For comparison only
19b00333 2590 esd->SetTPCClusterMap(seed->GetClusterMap());
2591 esd->SetTPCSharedMap(seed->GetSharedMap());
af32720d 2592 //
6d493ea0 2593 if (AliTPCReconstructor::StreamLevel()>1 && seed!=0&&esd!=0) {
af32720d 2594 TTreeSRedirector &cstream = *fDebugStreamer;
2595 cstream<<"Crefit"<<
2596 "Esd.="<<esd<<
2597 "Track.="<<seed<<
2598 "\n";
2599 }
c274e255 2600
51ad6848 2601 if (seed->GetNumberOfClusters()>15){
4d158c36 2602 esd->UpdateTrackParams(seed,AliESDtrack::kTPCrefit);
51ad6848 2603 esd->SetTPCPoints(seed->GetPoints());
b9671574 2604 esd->SetTPCPointsF(seed->GetNFoundable());
2605 Int_t ndedx = seed->GetNCDEDX(0)+seed->GetNCDEDX(1)+seed->GetNCDEDX(2)+seed->GetNCDEDX(3);
2606 Float_t sdedx = (seed->GetSDEDX(0)+seed->GetSDEDX(1)+seed->GetSDEDX(2)+seed->GetSDEDX(3))*0.25;
b406fdb0 2607 Float_t dedx = seed->GetdEdx();
2608 esd->SetTPCsignal(dedx, sdedx, ndedx);
fdedfdec 2609 //
2610 // add seed to the esd track in Calib level
2611 //
6d493ea0 2612 Bool_t storeFriend = gRandom->Rndm()<(kMaxFriendTracks)/Float_t(nseed);
2613 if (AliTPCReconstructor::StreamLevel()>0 &&storeFriend){
fdedfdec 2614 AliTPCseed * seedCopy = new AliTPCseed(*seed, kTRUE);
2615 esd->AddCalibObject(seedCopy);
2616 }
4d158c36 2617 ntracks++;
d26d9159 2618 }
2619 else{
2620 //printf("problem\n");
2621 }
2622 }
51ad6848 2623 //FindKinks(fSeeds,event);
4d158c36 2624 Info("RefitInward","Number of refitted tracks %d",ntracks);
d26d9159 2625 return 0;
2626}
2627
1c53abe2 2628
af885e0f 2629Int_t AliTPCtrackerMI::PropagateBack(AliESDEvent *event)
91162307 2630{
2631 //
2632 // back propagation of ESD tracks
2633 //
1c53abe2 2634
91162307 2635 fEvent = event;
d26d9159 2636 fIteration = 1;
5d837844 2637 ReadSeeds(event,1);
b406fdb0 2638 PropagateBack(fSeeds);
2639 RemoveUsed2(fSeeds,0.4,0.4,20);
6d493ea0 2640 //FindCurling(fSeeds, fEvent,1);
2641 FindSplitted(fSeeds, event,1); // find multi found tracks
2642
b406fdb0 2643 //
91162307 2644 Int_t nseed = fSeeds->GetEntriesFast();
4d158c36 2645 Int_t ntracks=0;
91162307 2646 for (Int_t i=0;i<nseed;i++){
2647 AliTPCseed * seed = (AliTPCseed*) fSeeds->UncheckedAt(i);
d9e62e7e 2648 if (!seed) continue;
51ad6848 2649 if (seed->GetKinkIndex(0)<0) UpdateKinkQualityM(seed); // update quality informations for kinks
2650 seed->UpdatePoints();
92f513f5 2651 AddCovariance(seed);
91162307 2652 AliESDtrack *esd=event->GetTrack(i);
d26d9159 2653 seed->CookdEdx(0.02,0.6);
91162307 2654 CookLabel(seed,0.1); //For comparison only
51ad6848 2655 if (seed->GetNumberOfClusters()>15){
4d158c36 2656 esd->UpdateTrackParams(seed,AliESDtrack::kTPCout);
51ad6848 2657 esd->SetTPCPoints(seed->GetPoints());
b9671574 2658 esd->SetTPCPointsF(seed->GetNFoundable());
2659 Int_t ndedx = seed->GetNCDEDX(0)+seed->GetNCDEDX(1)+seed->GetNCDEDX(2)+seed->GetNCDEDX(3);
2660 Float_t sdedx = (seed->GetSDEDX(0)+seed->GetSDEDX(1)+seed->GetSDEDX(2)+seed->GetSDEDX(3))*0.25;
167c41ab 2661 Float_t dedx = seed->GetdEdx();
2662 esd->SetTPCsignal(dedx, sdedx, ndedx);
4d158c36 2663 ntracks++;
31fd97b2 2664 Int_t eventnumber = event->GetEventNumberInFile();// patch 28 fev 06
2665 // This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number
8cecaa87 2666 if (AliTPCReconstructor::StreamLevel()>1 && esd) {
6d493ea0 2667 (*fDebugStreamer)<<"Cback"<<
2668 "Tr0.="<<seed<<
8cecaa87 2669 "esd.="<<esd<<
6d493ea0 2670 "EventNrInFile="<<eventnumber<<
8cecaa87 2671 "\n";
6d493ea0 2672 }
4d158c36 2673 }
91162307 2674 }
51ad6848 2675 //FindKinks(fSeeds,event);
4d158c36 2676 Info("PropagateBack","Number of back propagated tracks %d",ntracks);
91162307 2677 fEvent =0;
6d493ea0 2678
91162307 2679 return 0;
2680}
2681
2682
2683void AliTPCtrackerMI::DeleteSeeds()
2684{
b67e07dc 2685 //
2686 //delete Seeds
6c94f330 2687
91162307 2688 Int_t nseed = fSeeds->GetEntriesFast();
2689 for (Int_t i=0;i<nseed;i++){
2690 AliTPCseed * seed = (AliTPCseed*)fSeeds->At(i);
2691 if (seed) delete fSeeds->RemoveAt(i);
2692 }
2693 delete fSeeds;
6c94f330 2694
91162307 2695 fSeeds =0;
2696}
2697
af885e0f 2698void AliTPCtrackerMI::ReadSeeds(AliESDEvent *event, Int_t direction)
91162307 2699{
2700 //
2701 //read seeds from the event
2702
2703 Int_t nentr=event->GetNumberOfTracks();
6bdc18d6 2704 if (fDebug>0){
2705 Info("ReadSeeds", "Number of ESD tracks: %d\n", nentr);
2706 }
91162307 2707 if (fSeeds)
2708 DeleteSeeds();
2709 if (!fSeeds){
4d158c36 2710 fSeeds = new TObjArray(nentr);
91162307 2711 }
4d158c36 2712 UnsignClusters();
2713 // Int_t ntrk=0;
91162307 2714 for (Int_t i=0; i<nentr; i++) {
2715 AliESDtrack *esd=event->GetTrack(i);
51ad6848 2716 ULong_t status=esd->GetStatus();
2717 if (!(status&AliESDtrack::kTPCin)) continue;
b67e07dc 2718 AliTPCtrack t(*esd);
5d837844 2719 t.SetNumberOfClusters(0);
eea478d3 2720 // AliTPCseed *seed = new AliTPCseed(t,t.GetAlpha());
2721 AliTPCseed *seed = new AliTPCseed(t/*,t.GetAlpha()*/);
a0f4d6a6 2722 seed->SetUniqueID(esd->GetID());
92f513f5 2723 AddCovariance(seed); //add systematic ucertainty
eea478d3 2724 for (Int_t ikink=0;ikink<3;ikink++) {
2725 Int_t index = esd->GetKinkIndex(ikink);
2726 seed->GetKinkIndexes()[ikink] = index;
2727 if (index==0) continue;
2728 index = TMath::Abs(index);
2729 AliESDkink * kink = fEvent->GetKink(index-1);
2730 if (kink&&esd->GetKinkIndex(ikink)<0){
2731 if ((status & AliESDtrack::kTRDrefit) != 0) kink->SetStatus(1,2);
2732 if ((status & AliESDtrack::kITSout) != 0) kink->SetStatus(1,0);
2733 }
2734 if (kink&&esd->GetKinkIndex(ikink)>0){
2735 if ((status & AliESDtrack::kTRDrefit) != 0) kink->SetStatus(1,6);
2736 if ((status & AliESDtrack::kITSout) != 0) kink->SetStatus(1,4);
2737 }
2738
2739 }
6c94f330 2740 if (((status&AliESDtrack::kITSout)==0)&&(direction==1)) seed->ResetCovariance(10.);
2741 if ( direction ==2 &&(status & AliESDtrack::kTRDrefit) == 0 ) seed->ResetCovariance(10.);
4d158c36 2742 if ( direction ==2 && ((status & AliESDtrack::kTPCout) == 0) ) {
2743 fSeeds->AddAt(0,i);
2744 delete seed;
2745 continue;
2746 }
2747 if ( direction ==2 &&(status & AliESDtrack::kTRDrefit) > 0 ) {
c0b978f0 2748 Double_t par0[5],par1[5],alpha,x;
2749 esd->GetInnerExternalParameters(alpha,x,par0);
4d158c36 2750 esd->GetExternalParameters(x,par1);
2751 Double_t delta1 = TMath::Abs(par0[4]-par1[4])/(0.000000001+TMath::Abs(par0[4]+par1[4]));
2752 Double_t delta2 = TMath::Abs(par0[3]-par1[3]);
51ad6848 2753 Double_t trdchi2=0;
2754 if (esd->GetTRDncls()>0) trdchi2 = esd->GetTRDchi2()/esd->GetTRDncls();
4d158c36 2755 //reset covariance if suspicious
51ad6848 2756 if ( (delta1>0.1) || (delta2>0.006) ||trdchi2>7.)
6c94f330 2757 seed->ResetCovariance(10.);
4d158c36 2758 }
982aff31 2759
91162307 2760 //
2761 //
2762 // rotate to the local coordinate system
eea478d3 2763 //
2764 fSectors=fInnerSec; fN=fkNIS;
91162307 2765 Double_t alpha=seed->GetAlpha() - fSectors->GetAlphaShift();
2766 if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
2767 if (alpha < 0. ) alpha += 2.*TMath::Pi();
2768 Int_t ns=Int_t(alpha/fSectors->GetAlpha())%fN;
2769 alpha =ns*fSectors->GetAlpha() + fSectors->GetAlphaShift();
4d158c36 2770 if (alpha<-TMath::Pi()) alpha += 2*TMath::Pi();
2771 if (alpha>=TMath::Pi()) alpha -= 2*TMath::Pi();
91162307 2772 alpha-=seed->GetAlpha();
d9b8978b 2773 if (!seed->Rotate(alpha)) {
2774 delete seed;
2775 continue;
2776 }
b9671574 2777 seed->SetESD(esd);
4d158c36 2778 // sign clusters
b406fdb0 2779 if (esd->GetKinkIndex(0)<=0){
2780 for (Int_t irow=0;irow<160;irow++){
2781 Int_t index = seed->GetClusterIndex2(irow);
2782 if (index>0){
2783 //
2784 AliTPCclusterMI * cl = GetClusterMI(index);
b9671574 2785 seed->SetClusterPointer(irow,cl);
b406fdb0 2786 if (cl){
2787 if ((index & 0x8000)==0){
2788 cl->Use(10); // accepted cluster
2789 }else{
2790 cl->Use(6); // close cluster not accepted
2791 }
4d158c36 2792 }else{
b406fdb0 2793 Info("ReadSeeds","Not found cluster");
2794 }
4d158c36 2795 }
2796 }
2797 }
2798 fSeeds->AddAt(seed,i);
91162307 2799 }
2800}
2801
2802
2803
2804//_____________________________________________________________________________
2805void AliTPCtrackerMI::MakeSeeds3(TObjArray * arr, Int_t sec, Int_t i1, Int_t i2, Float_t cuts[4],
2806 Float_t deltay, Int_t ddsec) {
2807 //-----------------------------------------------------------------
2808 // This function creates track seeds.
2809 // SEEDING WITH VERTEX CONSTRAIN
2810 //-----------------------------------------------------------------
2811 // cuts[0] - fP4 cut
2812 // cuts[1] - tan(phi) cut
2813 // cuts[2] - zvertex cut
2814 // cuts[3] - fP3 cut
2815 Int_t nin0 = 0;
2816 Int_t nin1 = 0;
2817 Int_t nin2 = 0;
2818 Int_t nin = 0;
2819 Int_t nout1 = 0;
2820 Int_t nout2 = 0;
2821
2822 Double_t x[5], c[15];
2823 // Int_t di = i1-i2;
2824 //
316c6cd9 2825 AliTPCseed * seed = new AliTPCseed();
91162307 2826 Double_t alpha=fSectors->GetAlpha(), shift=fSectors->GetAlphaShift();
2827 Double_t cs=cos(alpha), sn=sin(alpha);
2828 //
2829 // Double_t x1 =fOuterSec->GetX(i1);
2830 //Double_t xx2=fOuterSec->GetX(i2);
2831
2832 Double_t x1 =GetXrow(i1);
2833 Double_t xx2=GetXrow(i2);
2834
2835 Double_t x3=GetX(), y3=GetY(), z3=GetZ();
2836
2837 Int_t imiddle = (i2+i1)/2; //middle pad row index
2838 Double_t xm = GetXrow(imiddle); // radius of middle pad-row
bd26fa83 2839 const AliTPCtrackerRow& krm=GetRow(sec,imiddle); //middle pad -row
91162307 2840 //
2841 Int_t ns =sec;
2842
bd26fa83 2843 const AliTPCtrackerRow& kr1=GetRow(ns,i1);
b9671574 2844 Double_t ymax = GetMaxY(i1)-kr1.GetDeadZone()-1.5;
2845 Double_t ymaxm = GetMaxY(imiddle)-kr1.GetDeadZone()-1.5;
91162307 2846
2847 //
2848 // change cut on curvature if it can't reach this layer
2849 // maximal curvature set to reach it
2850 Double_t dvertexmax = TMath::Sqrt((x1-x3)*(x1-x3)+(ymax+5-y3)*(ymax+5-y3));
2851 if (dvertexmax*0.5*cuts[0]>0.85){
2852 cuts[0] = 0.85/(dvertexmax*0.5+1.);
2853 }
2854 Double_t r2min = 1/(cuts[0]*cuts[0]); //minimal square of radius given by cut
2855
2856 // Int_t ddsec = 1;
2857 if (deltay>0) ddsec = 0;
2858 // loop over clusters
2859 for (Int_t is=0; is < kr1; is++) {
2860 //
2861 if (kr1[is]->IsUsed(10)) continue;
2862 Double_t y1=kr1[is]->GetY(), z1=kr1[is]->GetZ();
2863 //if (TMath::Abs(y1)>ymax) continue;
2864
2865 if (deltay>0 && TMath::Abs(ymax-TMath::Abs(y1))> deltay ) continue; // seed only at the edge
2866
2867 // find possible directions
2868 Float_t anglez = (z1-z3)/(x1-x3);
2869 Float_t extraz = z1 - anglez*(x1-xx2); // extrapolated z
2870 //
2871 //
2872 //find rotation angles relative to line given by vertex and point 1
2873 Double_t dvertex2 = (x1-x3)*(x1-x3)+(y1-y3)*(y1-y3);
2874 Double_t dvertex = TMath::Sqrt(dvertex2);
2875 Double_t angle13 = TMath::ATan((y1-y3)/(x1-x3));
2876 Double_t cs13 = cos(-angle13), sn13 = sin(-angle13);
2877
2878 //
2879 // loop over 2 sectors
2880 Int_t dsec1=-ddsec;
2881 Int_t dsec2= ddsec;
2882 if (y1<0) dsec2= 0;
2883 if (y1>0) dsec1= 0;
2884
2885 Double_t dddz1=0; // direction of delta inclination in z axis
2886 Double_t dddz2=0;
2887 if ( (z1-z3)>0)
2888 dddz1 =1;
2889 else
2890 dddz2 =1;
2891 //
2892 for (Int_t dsec = dsec1; dsec<=dsec2;dsec++){
2893 Int_t sec2 = sec + dsec;
2894 //
bd26fa83 2895 // AliTPCtrackerRow& kr2 = fOuterSec[(sec2+fkNOS)%fkNOS][i2];
2896 //AliTPCtrackerRow& kr2m = fOuterSec[(sec2+fkNOS)%fkNOS][imiddle];
2897 AliTPCtrackerRow& kr2 = GetRow((sec2+fkNOS)%fkNOS,i2);
2898 AliTPCtrackerRow& kr2m = GetRow((sec2+fkNOS)%fkNOS,imiddle);
91162307 2899 Int_t index1 = TMath::Max(kr2.Find(extraz-0.6-dddz1*TMath::Abs(z1)*0.05)-1,0);
2900 Int_t index2 = TMath::Min(kr2.Find(extraz+0.6+dddz2*TMath::Abs(z1)*0.05)+1,kr2);
2901
2902 // rotation angles to p1-p3
2903 Double_t cs13r = cos(-angle13+dsec*alpha)/dvertex, sn13r = sin(-angle13+dsec*alpha)/dvertex;
2904 Double_t x2, y2, z2;
2905 //
2906 // Double_t dymax = maxangle*TMath::Abs(x1-xx2);
2907
2908 //
2909 Double_t dxx0 = (xx2-x3)*cs13r;
2910 Double_t dyy0 = (xx2-x3)*sn13r;
2911 for (Int_t js=index1; js < index2; js++) {
2912 const AliTPCclusterMI *kcl = kr2[js];
2913 if (kcl->IsUsed(10)) continue;
2914 //
2915 //calcutate parameters
2916 //
2917 Double_t yy0 = dyy0 +(kcl->GetY()-y3)*cs13r;
2918 // stright track
2919 if (TMath::Abs(yy0)<0.000001) continue;
2920 Double_t xx0 = dxx0 -(kcl->GetY()-y3)*sn13r;
2921 Double_t y0 = 0.5*(xx0*xx0+yy0*yy0-xx0)/yy0;
2922 Double_t r02 = (0.25+y0*y0)*dvertex2;
2923 //curvature (radius) cut
2924 if (r02<r2min) continue;
2925
2926 nin0++;
2927 //
2928 Double_t c0 = 1/TMath::Sqrt(r02);
2929 if (yy0>0) c0*=-1.;
2930
2931
2932 //Double_t dfi0 = 2.*TMath::ASin(dvertex*c0*0.5);
2933 //Double_t dfi1 = 2.*TMath::ASin(TMath::Sqrt(yy0*yy0+(1-xx0)*(1-xx0))*dvertex*c0*0.5);
b67e07dc 2934 Double_t dfi0 = 2.*AliTPCFastMath::FastAsin(dvertex*c0*0.5);
2935 Double_t dfi1 = 2.*AliTPCFastMath::FastAsin(TMath::Sqrt(yy0*yy0+(1-xx0)*(1-xx0))*dvertex*c0*0.5);
91162307 2936 //
2937 //
2938 Double_t z0 = kcl->GetZ();
2939 Double_t zzzz2 = z1-(z1-z3)*dfi1/dfi0;
2940 if (TMath::Abs(zzzz2-z0)>0.5) continue;
2941 nin1++;
2942 //
2943 Double_t dip = (z1-z0)*c0/dfi1;
2944 Double_t x0 = (0.5*cs13+y0*sn13)*dvertex*c0;
2945 //
2946 y2 = kcl->GetY();
2947 if (dsec==0){
2948 x2 = xx2;
2949 z2 = kcl->GetZ();
2950 }
2951 else
2952 {
2953 // rotation
2954 z2 = kcl->GetZ();
2955 x2= xx2*cs-y2*sn*dsec;
2956 y2=+xx2*sn*dsec+y2*cs;
2957 }
2958
2959 x[0] = y1;
2960 x[1] = z1;
2961 x[2] = x0;
2962 x[3] = dip;
2963 x[4] = c0;
2964 //
2965 //
2966 // do we have cluster at the middle ?
2967 Double_t ym,zm;
2968 GetProlongation(x1,xm,x,ym,zm);
2969 UInt_t dummy;
2970 AliTPCclusterMI * cm=0;
2971 if (TMath::Abs(ym)-ymaxm<0){
2972 cm = krm.FindNearest2(ym,zm,1.0,0.6,dummy);
2973 if ((!cm) || (cm->IsUsed(10))) {
2974 continue;
2975 }
2976 }
2977 else{
2978 // rotate y1 to system 0
2979 // get state vector in rotated system
2980 Double_t yr1 = (-0.5*sn13+y0*cs13)*dvertex*c0;
2981 Double_t xr2 = x0*cs+yr1*sn*dsec;
2982 Double_t xr[5]={kcl->GetY(),kcl->GetZ(), xr2, dip, c0};
2983 //
2984 GetProlongation(xx2,xm,xr,ym,zm);
2985 if (TMath::Abs(ym)-ymaxm<0){
2986 cm = kr2m.FindNearest2(ym,zm,1.0,0.6,dummy);
2987 if ((!cm) || (cm->IsUsed(10))) {
2988 continue;
2989 }
2990 }
2991 }
2992
2993
2994 Double_t dym = 0;
2995 Double_t dzm = 0;
2996 if (cm){
2997 dym = ym - cm->GetY();
2998 dzm = zm - cm->GetZ();
2999 }
3000 nin2++;
3001
3002
3003 //
3004 //
3005 Double_t sy1=kr1[is]->GetSigmaY2()*2., sz1=kr1[is]->GetSigmaZ2()*2.;
3006 Double_t sy2=kcl->GetSigmaY2()*2., sz2=kcl->GetSigmaZ2()*2.;
3007 //Double_t sy3=400*3./12., sy=0.1, sz=0.1;
3008 Double_t sy3=25000*x[4]*x[4]+0.1, sy=0.1, sz=0.1;
3009 //Double_t sy3=25000*x[4]*x[4]*60+0.5, sy=0.1, sz=0.1;
3010
b67e07dc 3011 Double_t f40=(F1(x1,y1+sy,x2,y2,x3,y3)-x[4])/sy;
3012 Double_t f42=(F1(x1,y1,x2,y2+sy,x3,y3)-x[4])/sy;
3013 Double_t f43=(F1(x1,y1,x2,y2,x3,y3+sy)-x[4])/sy;
3014 Double_t f20=(F2(x1,y1+sy,x2,y2,x3,y3)-x[2])/sy;
3015 Double_t f22=(F2(x1,y1,x2,y2+sy,x3,y3)-x[2])/sy;
3016 Double_t f23=(F2(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy;
91162307 3017
b67e07dc 3018 Double_t f30=(F3(x1,y1+sy,x2,y2,z1,z2)-x[3])/sy;
3019 Double_t f31=(F3(x1,y1,x2,y2,z1+sz,z2)-x[3])/sz;
3020 Double_t f32=(F3(x1,y1,x2,y2+sy,z1,z2)-x[3])/sy;
3021 Double_t f34=(F3(x1,y1,x2,y2,z1,z2+sz)-x[3])/sz;
91162307 3022
1c53abe2 3023 c[0]=sy1;
3024 c[1]=0.; c[2]=sz1;
3025 c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f23*sy3*f23;
3026 c[6]=f30*sy1; c[7]=f31*sz1; c[8]=f30*sy1*f20+f32*sy2*f22;
3027 c[9]=f30*sy1*f30+f31*sz1*f31+f32*sy2*f32+f34*sz2*f34;
3028 c[10]=f40*sy1; c[11]=0.; c[12]=f40*sy1*f20+f42*sy2*f22+f43*sy3*f23;
3029 c[13]=f30*sy1*f40+f32*sy2*f42;
3030 c[14]=f40*sy1*f40+f42*sy2*f42+f43*sy3*f43;
91162307 3031
3032 // if (!BuildSeed(kr1[is],kcl,0,x1,x2,x3,x,c)) continue;
3033
1c53abe2 3034 UInt_t index=kr1.GetIndex(is);
316c6cd9 3035 seed->~AliTPCseed(); // this does not set the pointer to 0...
6c94f330 3036 AliTPCseed *track=new(seed) AliTPCseed(x1, ns*alpha+shift, x, c, index);
91162307 3037
b9671574 3038 track->SetIsSeeding(kTRUE);
3039 track->SetSeed1(i1);
3040 track->SetSeed2(i2);
3041 track->SetSeedType(3);
c9427e08 3042
91162307 3043
3044 //if (dsec==0) {
3045 FollowProlongation(*track, (i1+i2)/2,1);
3046 Int_t foundable,found,shared;
3047 track->GetClusterStatistic((i1+i2)/2,i1, found, foundable, shared, kTRUE);
3048 if ((found<0.55*foundable) || shared>0.5*found || (track->GetSigmaY2()+track->GetSigmaZ2())>0.5){
3049 seed->Reset();
3050 seed->~AliTPCseed();
3051 continue;
3052 }
3053 //}
3054
3055 nin++;
3056 FollowProlongation(*track, i2,1);
3057
3058
3059 //Int_t rc = 1;
b9671574 3060 track->SetBConstrain(1);
91162307 3061 // track->fLastPoint = i1+fInnerSec->GetNRows(); // first cluster in track position
b9671574 3062 track->SetLastPoint(i1); // first cluster in track position
3063 track->SetFirstPoint(track->GetLastPoint());
91162307 3064
3065 if (track->GetNumberOfClusters()<(i1-i2)*0.5 ||
b9671574 3066 track->GetNumberOfClusters() < track->GetNFoundable()*0.6 ||
3067 track->GetNShared()>0.4*track->GetNumberOfClusters() ) {
91162307 3068 seed->Reset();
3069 seed->~AliTPCseed();
c9427e08 3070 continue;
3071 }
91162307 3072 nout1++;
3073 // Z VERTEX CONDITION
c274e255 3074 Double_t zv, bz=GetBz();
3075 if ( !track->GetZAt(0.,bz,zv) ) continue;
91162307 3076 if (TMath::Abs(zv-z3)>cuts[2]) {
3077 FollowProlongation(*track, TMath::Max(i2-20,0));
c274e255 3078 if ( !track->GetZAt(0.,bz,zv) ) continue;
91162307 3079 if (TMath::Abs(zv-z3)>cuts[2]){
3080 FollowProlongation(*track, TMath::Max(i2-40,0));
c274e255 3081 if ( !track->GetZAt(0.,bz,zv) ) continue;
b9671574 3082 if (TMath::Abs(zv-z3)>cuts[2] &&(track->GetNumberOfClusters() > track->GetNFoundable()*0.7)){
91162307 3083 // make seed without constrain
3084 AliTPCseed * track2 = MakeSeed(track,0.2,0.5,1.);
3085 FollowProlongation(*track2, i2,1);
b9671574 3086 track2->SetBConstrain(kFALSE);
3087 track2->SetSeedType(1);
91162307 3088 arr->AddLast(track2);
3089 seed->Reset();
3090 seed->~AliTPCseed();
3091 continue;
3092 }
3093 else{
3094 seed->Reset();
3095 seed->~AliTPCseed();
3096 continue;
3097
3098 }
3099 }
c9427e08 3100 }
316c6cd9 3101
b9671574 3102 track->SetSeedType(0);
91162307 3103 arr->AddLast(track);
3104 seed = new AliTPCseed;
3105 nout2++;
3106 // don't consider other combinations
b9671574 3107 if (track->GetNumberOfClusters() > track->GetNFoundable()*0.8)
91162307 3108 break;
1c53abe2 3109 }
3110 }
3111 }
6bdc18d6 3112 if (fDebug>3){
3113 Info("MakeSeeds3","\nSeeding statistic:\t%d\t%d\t%d\t%d\t%d\t%d",nin0,nin1,nin2,nin,nout1,nout2);
91162307 3114 }
3115 delete seed;
1c53abe2 3116}
3117
1627d1c4 3118
91162307 3119void AliTPCtrackerMI::MakeSeeds5(TObjArray * arr, Int_t sec, Int_t i1, Int_t i2, Float_t cuts[4],
3120 Float_t deltay) {
3121
3122
3123
1627d1c4 3124 //-----------------------------------------------------------------
91162307 3125 // This function creates track seeds.
1627d1c4 3126 //-----------------------------------------------------------------
91162307 3127 // cuts[0] - fP4 cut
3128 // cuts[1] - tan(phi) cut
3129 // cuts[2] - zvertex cut
3130 // cuts[3] - fP3 cut
3131
3132
3133 Int_t nin0 = 0;
3134 Int_t nin1 = 0;
3135 Int_t nin2 = 0;
3136 Int_t nin = 0;
3137 Int_t nout1 = 0;
3138 Int_t nout2 = 0;
3139 Int_t nout3 =0;
3140 Double_t x[5], c[15];
3141 //
3142 // make temporary seed
3143 AliTPCseed * seed = new AliTPCseed;
1627d1c4 3144 Double_t alpha=fOuterSec->GetAlpha(), shift=fOuterSec->GetAlphaShift();
3145 // Double_t cs=cos(alpha), sn=sin(alpha);
91162307 3146 //
3147 //
1627d1c4 3148
91162307 3149 // first 3 padrows
3150 Double_t x1 = GetXrow(i1-1);
bd26fa83 3151 const AliTPCtrackerRow& kr1=GetRow(sec,i1-1);
b9671574 3152 Double_t y1max = GetMaxY(i1-1)-kr1.GetDeadZone()-1.5;
91162307 3153 //
3154 Double_t x1p = GetXrow(i1);
bd26fa83 3155 const AliTPCtrackerRow& kr1p=GetRow(sec,i1);
91162307 3156 //
3157 Double_t x1m = GetXrow(i1-2);
bd26fa83 3158 const AliTPCtrackerRow& kr1m=GetRow(sec,i1-2);
1627d1c4 3159
91162307 3160 //
3161 //last 3 padrow for seeding
bd26fa83 3162 AliTPCtrackerRow& kr3 = GetRow((sec+fkNOS)%fkNOS,i1-7);
91162307 3163 Double_t x3 = GetXrow(i1-7);
3164 // Double_t y3max= GetMaxY(i1-7)-kr3.fDeadZone-1.5;
3165 //
bd26fa83 3166 AliTPCtrackerRow& kr3p = GetRow((sec+fkNOS)%fkNOS,i1-6);
91162307 3167 Double_t x3p = GetXrow(i1-6);
3168 //
bd26fa83 3169 AliTPCtrackerRow& kr3m = GetRow((sec+fkNOS)%fkNOS,i1-8);
91162307 3170 Double_t x3m = GetXrow(i1-8);
1627d1c4 3171
91162307 3172 //
3173 //
3174 // middle padrow
3175 Int_t im = i1-4; //middle pad row index
3176 Double_t xm = GetXrow(im); // radius of middle pad-row
bd26fa83 3177 const AliTPCtrackerRow& krm=GetRow(sec,im); //middle pad -row
91162307 3178 // Double_t ymmax = GetMaxY(im)-kr1.fDeadZone-1.5;
3179 //
3180 //
3181 Double_t deltax = x1-x3;
3182 Double_t dymax = deltax*cuts[1];
3183 Double_t dzmax = deltax*cuts[3];
3184 //
3185 // loop over clusters
3186 for (Int_t is=0; is < kr1; is++) {
1627d1c4 3187 //
91162307 3188 if (kr1[is]->IsUsed(10)) continue;
3189 Double_t y1=kr1[is]->GetY(), z1=kr1[is]->GetZ();
1627d1c4 3190 //
91162307