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