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