]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITStrackerSA.cxx
Modifications to ITS standalone tracker: 1) improve efficiency at low pt for p-p...
[u/mrichter/AliRoot.git] / ITS / AliITStrackerSA.cxx
CommitLineData
13918578 1/**************************************************************************
2 * Copyright(c) 1998-2003, 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
bf6adc12 16/* $Id$ */
17
13918578 18////////////////////////////////////////////////////
19// Stand alone tracker class //
20// Origin: Elisabetta Crescio //
21// e-mail: crescio@to.infn.it //
22// tracks are saved as AliITStrackV2 objects //
23////////////////////////////////////////////////////
24
aa4b78e1 25#include <stdlib.h>
bf6adc12 26
27#include <TArrayI.h>
13918578 28#include <TBranch.h>
13918578 29#include <TObjArray.h>
30#include <TTree.h>
bf6adc12 31
af885e0f 32#include "AliESDEvent.h"
bf6adc12 33#include "AliESDVertex.h"
34#include "AliESDtrack.h"
bf6adc12 35#include "AliITSVertexer.h"
13918578 36#include "AliITSclusterTable.h"
00a7cc50 37#include "AliITSRecPoint.h"
cc088660 38#include "AliITSgeomTGeo.h"
13918578 39#include "AliITStrackSA.h"
bf6adc12 40#include "AliITStrackerSA.h"
e340bb86 41#include "AliITSReconstructor.h"
bf6adc12 42#include "AliRun.h"
13918578 43
44ClassImp(AliITStrackerSA)
45
46//____________________________________________________________________________
e399f961 47AliITStrackerSA::AliITStrackerSA():AliITStrackerMI(),
48fPhiEstimate(0),
1966f03c 49fITSStandAlone(0),
e399f961 50fLambdac(0),
51fPhic(0),
52fCoef1(0),
53fCoef2(0),
54fCoef3(0),
55fNloop(0),
56fPhiWin(0),
57fLambdaWin(0),
58fVert(0),
59fVertexer(0),
e399f961 60fListOfTracks(0),
deae0246 61fListOfSATracks(0),
e399f961 62fITSclusters(0),
63fSixPoints(0),
b8ed1a92 64fOuterStartLayer(0),
3733ccd2 65fMinQ(0.),
e399f961 66fCluLayer(0),
67fCluCoord(0){
13918578 68 // Default constructor
69 Init();
627b6db6 70
13918578 71}
2257f27e 72//____________________________________________________________________________
cc088660 73AliITStrackerSA::AliITStrackerSA(const Char_t *geom):AliITStrackerMI(0),
e399f961 74fPhiEstimate(0),
1966f03c 75fITSStandAlone(0),
e399f961 76fLambdac(0),
77fPhic(0),
78fCoef1(0),
79fCoef2(0),
80fCoef3(0),
81fNloop(0),
82fPhiWin(0),
83fLambdaWin(0),
84fVert(0),
85fVertexer(0),
e399f961 86fListOfTracks(0),
deae0246 87fListOfSATracks(0),
e399f961 88fITSclusters(0),
89fSixPoints(0),
b8ed1a92 90fOuterStartLayer(0),
3733ccd2 91fMinQ(0.),
e399f961 92fCluLayer(0),
93fCluCoord(0)
2257f27e 94{
95 // Standard constructor (Vertex is known and passed to this obj.)
cc088660 96 if (geom) {
97 AliWarning("\"geom\" is actually a dummy argument !");
98 }
99
2257f27e 100 Init();
101 fVert = 0;
627b6db6 102
2257f27e 103}
104
13918578 105//____________________________________________________________________________
cc088660 106AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliESDVertex *vert):AliITStrackerMI(0),
e399f961 107fPhiEstimate(0),
1966f03c 108fITSStandAlone(0),
e399f961 109fLambdac(0),
110fPhic(0),
111fCoef1(0),
112fCoef2(0),
113fCoef3(0),
114fNloop(0),
115fPhiWin(0),
116fLambdaWin(0),
117fVert(vert),
118fVertexer(0),
e399f961 119fListOfTracks(0),
deae0246 120fListOfSATracks(0),
e399f961 121fITSclusters(0),
122fSixPoints(0),
b8ed1a92 123fOuterStartLayer(0),
3733ccd2 124fMinQ(0.),
e399f961 125fCluLayer(0),
1966f03c 126fCluCoord(0)
13918578 127{
128 // Standard constructor (Vertex is known and passed to this obj.)
cc088660 129 if (geom) {
130 AliWarning("\"geom\" is actually a dummy argument !");
131 }
13918578 132 Init();
627b6db6 133
13918578 134}
135
13918578 136//____________________________________________________________________________
cc088660 137AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliITSVertexer *vertexer):AliITStrackerMI(0),
e399f961 138fPhiEstimate(0),
1966f03c 139fITSStandAlone(0),
e399f961 140fLambdac(0),
141fPhic(0),
142fCoef1(0),
143fCoef2(0),
144fCoef3(0),
145fNloop(0),
146fPhiWin(0),
147fLambdaWin(0),
148fVert(),
149fVertexer(vertexer),
e399f961 150fListOfTracks(0),
deae0246 151fListOfSATracks(0),
e399f961 152fITSclusters(0),
153fSixPoints(0),
b8ed1a92 154fOuterStartLayer(0),
3733ccd2 155fMinQ(0.),
e399f961 156fCluLayer(0),
1966f03c 157fCluCoord(0)
13918578 158{
159 // Standard constructor (Vertex is unknown - vertexer is passed to this obj)
cc088660 160 if (geom) {
161 AliWarning("\"geom\" is actually a dummy argument !");
162 }
13918578 163 Init();
164 fVertexer = vertexer;
627b6db6 165
13918578 166}
167
168//____________________________________________________________________________
e399f961 169AliITStrackerSA::AliITStrackerSA(const AliITStrackerSA& tracker):AliITStrackerMI(),
170fPhiEstimate(tracker.fPhiEstimate),
1966f03c 171fITSStandAlone(tracker.fITSStandAlone),
e399f961 172fLambdac(tracker.fLambdac),
173fPhic(tracker.fPhic),
174fCoef1(tracker.fCoef1),
175fCoef2(tracker.fCoef2),
176fCoef3(tracker.fCoef3),
177fNloop(tracker.fNloop),
178fPhiWin(tracker.fPhiWin),
179fLambdaWin(tracker.fLambdaWin),
180fVert(tracker.fVert),
181fVertexer(tracker.fVertexer),
e399f961 182fListOfTracks(tracker.fListOfTracks),
deae0246 183fListOfSATracks(tracker.fListOfSATracks),
e399f961 184fITSclusters(tracker.fITSclusters),
185fSixPoints(tracker.fSixPoints),
b8ed1a92 186fOuterStartLayer(tracker.fOuterStartLayer),
3733ccd2 187fMinQ(tracker.fMinQ),
e399f961 188fCluLayer(tracker.fCluLayer),
189fCluCoord(tracker.fCluCoord) {
13918578 190 // Copy constructor
13918578 191 for(Int_t i=0;i<2;i++){
192 fPoint1[i]=tracker.fPoint1[i];
193 fPoint2[i]=tracker.fPoint2[i];
194 fPoint3[i]=tracker.fPoint3[i];
195 fPointc[i]=tracker.fPointc[i];
196 }
13918578 197 if(tracker.fVertexer && tracker.fVert){
d681bb2d 198 fVert = new AliESDVertex(*tracker.fVert);
13918578 199 }
200 else {
201 fVert = tracker.fVert;
202 }
cc088660 203 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
7e5bf5af 204 fCluLayer[i] = tracker.fCluLayer[i];
205 fCluCoord[i] = tracker.fCluCoord[i];
206 }
13918578 207}
7d62fb64 208//______________________________________________________________________
e399f961 209AliITStrackerSA& AliITStrackerSA::operator=(const AliITStrackerSA& source){
210 // Assignment operator.
211 this->~AliITStrackerSA();
212 new(this) AliITStrackerSA(source);
213 return *this;
214
7d62fb64 215}
13918578 216
217//____________________________________________________________________________
218AliITStrackerSA::~AliITStrackerSA(){
219 // destructor
d681bb2d 220 // if fVertexer is not null, the AliESDVertex obj. is owned by this class
13918578 221 // and is deleted here
222 if(fVertexer){
223 if(fVert)delete fVert;
224 }
225 fVert = 0;
226 fVertexer = 0;
7e5bf5af 227
13918578 228 if(fPhiWin)delete []fPhiWin;
229 if(fLambdaWin)delete []fLambdaWin;
4e05ab9a 230 fListOfTracks->Delete();
adda6c2b 231 delete fListOfTracks;
deae0246 232 fListOfSATracks->Delete();
233 delete fListOfSATracks;
7e5bf5af 234 if(fCluLayer){
cc088660 235 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
7e5bf5af 236 if(fCluLayer[i]){
237 fCluLayer[i]->Delete();
238 delete fCluLayer[i];
239 }
240 }
241 delete [] fCluLayer;
242 }
243 if(fCluCoord){
cc088660 244 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
7e5bf5af 245 if(fCluCoord[i]){
246 fCluCoord[i]->Delete();
247 delete fCluCoord[i];
248 }
249 }
250 delete [] fCluCoord;
627b6db6 251 }
7e5bf5af 252
7e5bf5af 253}
13918578 254
1966f03c 255//____________________________________________________________________________
af885e0f 256Int_t AliITStrackerSA::Clusters2Tracks(AliESDEvent *event){
1966f03c 257// This method is used to find and fit the tracks. By default the corresponding
258// method in the parent class is invoked. In this way a combined tracking
259// TPC+ITS is performed. If the flag fITSStandAlone is true, the tracking
260// is done in the ITS only. In the standard reconstruction chain this option
261// can be set via AliReconstruction::SetOption("ITS","onlyITS")
262 Int_t rc=0;
263 if(!fITSStandAlone){
264 rc=AliITStrackerMI::Clusters2Tracks(event);
265 }
266 else {
267 AliDebug(1,"Stand Alone flag set: doing tracking in ITS alone\n");
268 }
269 if(!rc) rc=FindTracks(event);
270 return rc;
271}
272
13918578 273//____________________________________________________________________________
274void AliITStrackerSA::Init(){
275 // Reset all data members
276 fPhiEstimate=0;
277 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
278 fLambdac=0;
279 fPhic=0;
280 fCoef1=0;
281 fCoef2=0;
282 fCoef3=0;
283 fPointc[0]=0;
284 fPointc[1]=0;
285 fVert = 0;
286 fVertexer = 0;
5b4dedce 287 Int_t nLoops=AliITSReconstructor::GetRecoParam()->GetNLoopsSA();
288 if(nLoops==33){
289 SetFixedWindowSizes();
290 }else{
291 Double_t phimin=AliITSReconstructor::GetRecoParam()->GetMinPhiSA();
292 Double_t phimax=AliITSReconstructor::GetRecoParam()->GetMaxPhiSA();
293 Double_t lambmin=AliITSReconstructor::GetRecoParam()->GetMinLambdaSA();
294 Double_t lambmax=AliITSReconstructor::GetRecoParam()->GetMaxLambdaSA();
295 SetCalculatedWindowSizes(nLoops,phimin,phimax,lambmin,lambmax);
296 }
3733ccd2 297 fMinQ=AliITSReconstructor::GetRecoParam()->GetSAMinClusterCharge();
2257f27e 298 fITSclusters = 0;
299 SetSixPoints();
b8ed1a92 300 SetOuterStartLayer(0);
1966f03c 301 SetSAFlag(kFALSE);
4e05ab9a 302 fListOfTracks=new TObjArray(0,0);
deae0246 303 fListOfSATracks=new TObjArray(0,0);
7e5bf5af 304 fCluLayer = 0;
305 fCluCoord = 0;
627b6db6 306 }
13918578 307//_______________________________________________________________________
308void AliITStrackerSA::ResetForFinding(){
309 // Reset data members used in all loops during track finding
310 fPhiEstimate=0;
311 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
312 fLambdac=0;
313 fPhic=0;
314 fCoef1=0;
315 fCoef2=0;
316 fCoef3=0;
317 fPointc[0]=0;
318 fPointc[1]=0;
4e05ab9a 319 fListOfTracks->Delete();
deae0246 320 fListOfSATracks->Delete();
13918578 321}
bef31448 322
7e5bf5af 323
2257f27e 324
bef31448 325//______________________________________________________________________
af885e0f 326Int_t AliITStrackerSA::FindTracks(AliESDEvent* event){
bef31448 327
7e5bf5af 328// Track finder using the ESD object
6f178e14 329
330
deae0246 331
2257f27e 332 if(!fITSclusters){
333 Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
334 return -1;
bef31448 335 }
7e5bf5af 336 //Reads event and mark clusters of traks already found, with flag kITSin
ebda3140 337 Int_t nentr=event->GetNumberOfTracks();
5a03f353 338 if(AliITSReconstructor::GetRecoParam()->GetSAUseAllClusters()==kFALSE) {
339 while (nentr--) {
340 AliESDtrack *track=event->GetTrack(nentr);
5af4a2d0 341 if ((track->GetStatus()&AliESDtrack::kITSin) == AliESDtrack::kITSin){
5a03f353 342 Int_t idx[12];
343 Int_t ncl = track->GetITSclusters(idx);
344 for(Int_t k=0;k<ncl;k++){
345 AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
346 cll->SetBit(kSAflag);
347 }
ebda3140 348 }
349 }
350 }
351 //Get primary vertex
352 Double_t primaryVertex[3];
353 event->GetVertex()->GetXYZ(primaryVertex);
354 //Creates TClonesArray with clusters for each layer. The clusters already used
355 //by AliITStrackerMI are not considered
deae0246 356 Int_t nclusters[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
357 Int_t dmar[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
ebda3140 358 if (fCluLayer == 0) {
deae0246 359 fCluLayer = new TClonesArray*[AliITSgeomTGeo::kNLayers];
360 fCluCoord = new TClonesArray*[AliITSgeomTGeo::kNLayers];
ebda3140 361 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
362 fCluLayer[i]=0;
363 fCluCoord[i]=0;
364 }
365 }
ebda3140 366 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
ebda3140 367 AliITSlayer &layer=fgLayers[i];
25be2b21 368 if (!AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) {
369 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
370 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
371 if(cls->TestBit(kSAflag)==kTRUE) continue; //clusters used by TPC prol.
372 if(cls->GetQ()==0) continue; //fake clusters dead zones
3733ccd2 373 if(i>1 && cls->GetQ()<=fMinQ) continue; // cut on SDD and SSD cluster charge
25be2b21 374 nclusters[i]++;
375 }
ebda3140 376 }
377 dmar[i]=0;
378 delete fCluLayer[i];
379 fCluLayer[i] = new TClonesArray("AliITSRecPoint",nclusters[i]);
380 delete fCluCoord[i];
381 fCluCoord[i] = new TClonesArray("AliITSclusterTable",nclusters[i]);
382 }
383
384 for(Int_t ilay=0;ilay<AliITSgeomTGeo::GetNLayers();ilay++){
385 TClonesArray &clulay = *fCluLayer[ilay];
386 TClonesArray &clucoo = *fCluCoord[ilay];
387 AliITSlayer &layer=fgLayers[ilay];
25be2b21 388 if (!AliITSReconstructor::GetRecoParam()->GetLayersToSkip(ilay)) {
389 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
390 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
391 if(cls->TestBit(kSAflag)==kTRUE) continue;
392 if(cls->GetQ()==0) continue;
3733ccd2 393 if(ilay>1 && cls->GetQ()<=fMinQ) continue;
25be2b21 394 Double_t phi=0;Double_t lambda=0;
395 Float_t x=0;Float_t y=0;Float_t z=0;
396 Float_t sx=0;Float_t sy=0;Float_t sz=0;
397 GetCoorAngles(cls,phi,lambda,x,y,z,primaryVertex);
398 GetCoorErrors(cls,sx,sy,sz);
399 new (clulay[dmar[ilay]]) AliITSRecPoint(*cls);
400 new (clucoo[dmar[ilay]]) AliITSclusterTable(x,y,z,sx,sy,sz,phi,lambda,cli);
401 dmar[ilay]++;
402 }
ebda3140 403 }
404 }
7e5bf5af 405
6f178e14 406 // track counter
407 Int_t ntrack=0;
408
ebda3140 409 Int_t minNPoints = (fSixPoints ? AliITSgeomTGeo::GetNLayers() : AliITSgeomTGeo::GetNLayers()-1);
410 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
411 if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) {
412 minNPoints--;
413 }
414 }
7e5bf5af 415
deae0246 416 static Int_t nClusLay[AliITSgeomTGeo::kNLayers];//counter for clusters on each layer
6f178e14 417
418 //loop on different minNPoints
419 Int_t minMinNPoints=minNPoints;
9870291d 420 if(AliITSReconstructor::GetRecoParam()->GetAllowProlongationWithEmptyRoad()) minMinNPoints=2;
6f178e14 421 for(Int_t iMinNPoints=minNPoints; iMinNPoints>=minMinNPoints; iMinNPoints--) {
422 //loop on the different windows
423 for(Int_t nloop=0;nloop<fNloop;nloop++){
424 for(Int_t ncl=0;ncl<fCluLayer[0]->GetEntries();ncl++){ //loop starting from layer 0
ebda3140 425
6f178e14 426 ResetForFinding();
427 Int_t pflag=0;
ebda3140 428
6f178e14 429 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[0]->At(ncl);
ebda3140 430
6f178e14 431 if(!cl) continue;
432 if (cl->GetQ()<=0) continue;
3733ccd2 433
6f178e14 434 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(0,ncl);
435 fPhic = arr->GetPhi();
436 fLambdac = arr->GetLambda();
6f178e14 437 fPhiEstimate = fPhic;
438 AliITStrackSA* trs = new AliITStrackSA();
439 fPoint1[0]=primaryVertex[0];
440 fPoint1[1]=primaryVertex[1];
441
442
443 fPoint2[0]=arr->GetX();
444 fPoint2[1]=arr->GetY();
445 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) nClusLay[i]=0;
446 nClusLay[0] = SearchClusters(0,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
447 nClusLay[1] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
448 if(nClusLay[1]>0){
449 pflag=1;
450 fPoint3[0] = fPointc[0];
451 fPoint3[1] = fPointc[1];
452 }
453 nClusLay[2] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
454 if(nClusLay[1]==0 && nClusLay[2]==0) pflag=0;
455 if(nClusLay[2]!=0 && nClusLay[1]!=0){ pflag=1; UpdatePoints();}
456 if(nClusLay[2]!=0 && nClusLay[1]==0){
457 pflag=1;
458 fPoint3[0]=fPointc[0];
459 fPoint3[1]=fPointc[1];
460 }
461
462 nClusLay[3] = SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
463 pflag=1;
464 if(nClusLay[3]!=0) UpdatePoints();
465 nClusLay[4] = SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
466 pflag=1;
467 if(nClusLay[4]!=0) UpdatePoints();
468 nClusLay[5] = SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
469
470
471 Int_t layOK=0;
472 //check of the candidate track
473 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers();nnp++) {
474 if(nClusLay[nnp]!=0) layOK+=1;
475 }
476
477 if(layOK>=iMinNPoints) {
8403ece4 478 //printf("-NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]);
6f178e14 479 AliITStrackV2* tr2 = 0;
480 tr2 = FitTrack(trs,primaryVertex);
481 if(!tr2) continue;
8403ece4 482 //printf("-NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
483
c421499f 484 StoreTrack(tr2,event);
6f178e14 485 ntrack++;
486 }
487 delete trs;
488 }//end loop on clusters of layer1
489 }//end loop on windows
490 }//end loop on min points
2755f080 491
7e5bf5af 492
ebda3140 493 minNPoints--;
7e5bf5af 494
ebda3140 495 //if 5/6 points are required, second loop starting
496 //from second layer (SPD2), to find tracks with point of
8403ece4 497 //layer 1 missing. Not done for cosmics.
498 if(!fSixPoints && fOuterStartLayer==0) {
ebda3140 499 //printf("looking from SPD2\n");
500 // counter for clusters on each layer
ebda3140 501 for(Int_t nloop=0;nloop<fNloop;nloop++){
502 Int_t ncl2=fCluLayer[1]->GetEntries();
503 while(ncl2--){ //loop starting from layer 2
504 ResetForFinding();
505 Int_t pflag=0;
506 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[1]->At(ncl2);
507
508 if(!cl) continue;
509 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(1,ncl2);
510 fPhic = arr->GetPhi();
511 fLambdac = arr->GetLambda();
512 fPhiEstimate = fPhic;
513
514 AliITStrackSA* trs = new AliITStrackSA();
515 fPoint1[0]=primaryVertex[0];
516 fPoint1[1]=primaryVertex[1];
517
518 fPoint2[0]=arr->GetX();
519 fPoint2[1]=arr->GetY();
deae0246 520 for(Int_t kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
521 nClusLay[0] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],
ebda3140 522 trs,primaryVertex[2],pflag);
deae0246 523 nClusLay[1] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],
ebda3140 524 trs,primaryVertex[2],pflag);
deae0246 525 if(nClusLay[1]!=0){
ebda3140 526 pflag=1;
527 fPoint3[0]=fPointc[0];
528 fPoint3[1]=fPointc[1];
529 }
deae0246 530 nClusLay[2]= SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],
ebda3140 531 trs,primaryVertex[2],pflag);
deae0246 532 if(nClusLay[2]!=0){
ebda3140 533 pflag=1;
534 UpdatePoints();
535 }
deae0246 536 nClusLay[3]= SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],
ebda3140 537 trs,primaryVertex[2],pflag);
deae0246 538 if(nClusLay[3]!=0){
ebda3140 539 pflag=1;
540 UpdatePoints();
541 }
deae0246 542 nClusLay[4]=SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],
ebda3140 543 trs,primaryVertex[2],pflag);
544
5a03f353 545 Int_t layOK=0;
ebda3140 546 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-1;nnp++){
deae0246 547 if(nClusLay[nnp]!=0) layOK+=1;
ebda3140 548 }
8403ece4 549 if(layOK>=minNPoints){ // 5/6
550 //printf("--NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]);
551
ebda3140 552 AliITStrackV2* tr2 = 0;
553 tr2 = FitTrack(trs,primaryVertex);
554 if(!tr2) continue;
8403ece4 555 //printf("--NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
ebda3140 556
c421499f 557 StoreTrack(tr2,event);
ebda3140 558 ntrack++;
559
560 }
561
bef31448 562 delete trs;
563 }//end loop on clusters of layer2
ebda3140 564 }
565
ebda3140 566 } //end opt="5/6"
567
568 // search for tracks starting from SPD2, SDD1, SDD2, SSD2
569 // for cosmics (A. Dainese 31.07.07)
9870291d 570 if(fOuterStartLayer>0 && !AliITSReconstructor::GetRecoParam()->GetAllowProlongationWithEmptyRoad()) {
ebda3140 571 for(Int_t innLay=1; innLay<=fOuterStartLayer; innLay++) {
572 //printf("Searching from layer %d outward\n",innLay);
573 minNPoints=AliITSgeomTGeo::GetNLayers()-innLay;
574 for(Int_t i=innLay;i<AliITSgeomTGeo::GetNLayers();i++)
575 if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i))
576 minNPoints--;
577 // counter for clusters on each layer
deae0246 578
ebda3140 579 for(Int_t nloop=0;nloop<fNloop;nloop++){
580 Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
581 while(nclInnLay--){ //loop starting from layer innLay
582 ResetForFinding();
583 Int_t pflag=0;
584 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[innLay]->At(nclInnLay);
585
586 if(!cl) continue;
587 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(innLay,nclInnLay);
588 fPhic = arr->GetPhi();
589 fLambdac = arr->GetLambda();
590 fPhiEstimate = fPhic;
591
592 AliITStrackSA* trs = new AliITStrackSA();
593 fPoint1[0]=primaryVertex[0];
594 fPoint1[1]=primaryVertex[1];
595 fPoint2[0]=arr->GetX();
596 fPoint2[1]=arr->GetY();
597
598 Int_t kk;
deae0246 599 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
ebda3140 600
601 kk=0;
deae0246 602 nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
b8ed1a92 603 trs,primaryVertex[2],pflag);
ebda3140 604 for(Int_t nextLay=innLay+1; nextLay<AliITSgeomTGeo::GetNLayers(); nextLay++) {
605 kk++;
deae0246 606 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
ebda3140 607 trs,primaryVertex[2],pflag);
deae0246 608 if(nClusLay[kk]!=0){
ebda3140 609 pflag=1;
610 if(kk==1) {
611 fPoint3[0]=fPointc[0];
612 fPoint3[1]=fPointc[1];
613 } else {
614 UpdatePoints();
615 }
616 }
617 }
618
5a03f353 619 Int_t layOK=0;
ebda3140 620 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
deae0246 621 if(nClusLay[nnp]!=0) layOK+=1;
ebda3140 622 }
5a03f353 623 if(layOK>=minNPoints){
8403ece4 624 //printf("---NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]);
ebda3140 625 AliITStrackV2* tr2 = 0;
626 tr2 = FitTrack(trs,primaryVertex);
627 if(!tr2) continue;
8403ece4 628 //printf("---NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
ebda3140 629
c421499f 630 StoreTrack(tr2,event);
ebda3140 631 ntrack++;
632
633 }
634
635 delete trs;
636 }//end loop on clusters of innLay
637 } //end loop on window sizes
638
ebda3140 639 } //end loop on innLay
640 } //end if(fOuterStartLayer>0)
641
5a03f353 642
8403ece4 643 // search for 1-point tracks in SPD, only for cosmics
5a03f353 644 // (A.Dainese 21.03.08)
645 if(AliITSReconstructor::GetRecoParam()->GetSAOnePointTracks() &&
646 TMath::Abs(event->GetMagneticField())<0.01) {
8403ece4 647 Int_t outerLayer=1; // only SPD
648 for(Int_t innLay=0; innLay<=TMath::Min(1,fOuterStartLayer); innLay++) {
5a03f353 649 // counter for clusters on each layer
deae0246 650
5a03f353 651 for(Int_t nloop=0;nloop<fNloop;nloop++){
652 Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
653 while(nclInnLay--){ //loop starting from layer innLay
654 ResetForFinding();
655 Int_t pflag=0;
656 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[innLay]->At(nclInnLay);
657
658 if(!cl) continue;
659 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(innLay,nclInnLay);
660 fPhic = arr->GetPhi();
661 fLambdac = arr->GetLambda();
662 fPhiEstimate = fPhic;
663
664 AliITStrackSA* trs = new AliITStrackSA();
665 fPoint1[0]=primaryVertex[0];
666 fPoint1[1]=primaryVertex[1];
667 fPoint2[0]=arr->GetX();
668 fPoint2[1]=arr->GetY();
669
670 Int_t kk;
deae0246 671 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
5a03f353 672
673 kk=0;
deae0246 674 nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
5a03f353 675 trs,primaryVertex[2],pflag);
8403ece4 676 for(Int_t nextLay=innLay+1; nextLay<=outerLayer; nextLay++) {
5a03f353 677 kk++;
deae0246 678 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
5a03f353 679 trs,primaryVertex[2],pflag);
deae0246 680 if(nClusLay[kk]!=0){
5a03f353 681 pflag=1;
682 if(kk==1) {
683 fPoint3[0]=fPointc[0];
684 fPoint3[1]=fPointc[1];
685 } else {
686 UpdatePoints();
687 }
688 }
689 }
690
691 Int_t layOK=0;
692 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
deae0246 693 if(nClusLay[nnp]!=0) layOK+=1;
5a03f353 694 }
695 if(layOK==1) {
8403ece4 696 //printf("----NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]);
5a03f353 697 AliITStrackV2* tr2 = 0;
698 Bool_t onePoint = kTRUE;
699 tr2 = FitTrack(trs,primaryVertex,onePoint);
700 if(!tr2) continue;
8403ece4 701 //printf("----NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
5a03f353 702
c421499f 703 StoreTrack(tr2,event);
5a03f353 704 ntrack++;
705
706 }
707
708 delete trs;
709 }//end loop on clusters of innLay
710 } //end loop on window sizes
711
5a03f353 712 } //end loop on innLay
713 } // end search 1-point tracks
ebda3140 714
bef31448 715 Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
2257f27e 716 return 0;
bef31448 717
718}
7e5bf5af 719
13918578 720//________________________________________________________________________
721
5a03f353 722AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex,Bool_t onePoint) {
ebda3140 723 //fit of the found track (most general case, also <6 points, layers missing)
724 // A.Dainese 16.11.07
13918578 725
13918578 726
deae0246 727 const Int_t kMaxClu=AliITStrackSA::kMaxNumberOfClusters;
13918578 728
deae0246 729 static Int_t firstmod[AliITSgeomTGeo::kNLayers];
730
731 static Int_t clind[AliITSgeomTGeo::kNLayers][kMaxClu];
732 static Int_t clmark[AliITSgeomTGeo::kNLayers][kMaxClu];
733 static Int_t nnn[AliITSgeomTGeo::kNLayers];
734 static Int_t kkk[AliITSgeomTGeo::kNLayers];
735 static Int_t end[AliITSgeomTGeo::kNLayers];
736 static AliITSRecPoint *listlayer[AliITSgeomTGeo::kNLayers][kMaxClu];
13918578 737
ebda3140 738 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
b8ed1a92 739 firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
ebda3140 740 nnn[i]=0;
741 kkk[i]=0;
deae0246 742 end[i]=0;
743 for(Int_t j=0;j<kMaxClu; j++){
ebda3140 744 clind[i][j]=0;
745 clmark[i][j]=0;
deae0246 746 listlayer[i][j]=0;
747 }
b8ed1a92 748 }
2755f080 749
2755f080 750
751 Int_t nclusters = tr->GetNumberOfClustersSA();
2755f080 752 for(Int_t ncl=0;ncl<nclusters;ncl++){
753 Int_t index = tr->GetClusterIndexSA(ncl);
754 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
755 if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
756 Int_t lay = (index & 0xf0000000) >> 28;
deae0246 757 Int_t nInLay=end[lay];
758 listlayer[lay][nInLay]=cl;
ebda3140 759 Int_t ind=nnn[lay];
760 clind[lay][ind]=index;
deae0246 761 end[lay]++;
2755f080 762 }
2755f080 763
764 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
765 for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
766 Int_t mark = tr->GetClusterMark(nlay,ncl);
ebda3140 767 Int_t ind=kkk[nlay];
768 clmark[nlay][ind]=mark;
2755f080 769 }
770 }
771
2755f080 772
773 Int_t firstLay=-1,secondLay=-1;
2755f080 774 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
deae0246 775 if(end[i]==0) {
2755f080 776 end[i]=1;
deae0246 777 }else{
2755f080 778 if(firstLay==-1) {
779 firstLay=i;
780 } else if(secondLay==-1) {
781 secondLay=i;
782 }
783 }
784 }
785
deae0246 786 if(firstLay==-1 || (secondLay==-1 && !onePoint)) return 0;
787
2755f080 788 for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
deae0246 789 AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0][l0];
2755f080 790 for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
deae0246 791 AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1][l1];
2755f080 792 for(Int_t l2=0;l2<end[2];l2++){ //loop on layer 3
deae0246 793 AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2][l2];
2755f080 794 for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4
deae0246 795 AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3][l3];
2755f080 796 for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
deae0246 797 AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4][l4];
2755f080 798 for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6
deae0246 799 AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5][l5];
2755f080 800
801
802 Double_t x1,y1,z1,sx1,sy1,sz1;
803 Double_t x2,y2,z2,sx2,sy2,sz2;
804 AliITSRecPoint* p1=0;
805 AliITSRecPoint* p2=0;
806 Int_t index1=0,index2=0;
807 Int_t mrk1=0,mrk2=0;
808
809 switch(firstLay) {
810 case 0:
811 p1=cl0;
ebda3140 812 index1=clind[0][l0];mrk1=clmark[0][l0];
2755f080 813 break;
814 case 1:
815 p1=cl1;
ebda3140 816 index1=clind[1][l1];mrk1=clmark[1][l1];
2755f080 817 break;
818 case 2:
819 p1=cl2;
ebda3140 820 index1=clind[2][l2];mrk1=clmark[2][l2];
2755f080 821 break;
822 case 3:
823 p1=cl3;
ebda3140 824 index1=clind[3][l3];mrk1=clmark[3][l3];
2755f080 825 break;
826 case 4:
827 p1=cl4;
ebda3140 828 index1=clind[4][l4];mrk1=clmark[4][l4];
2755f080 829 break;
830 }
831
832 switch(secondLay) {
833 case 1:
834 p2=cl1;
ebda3140 835 index2=clind[1][l1];mrk2=clmark[1][l1];
2755f080 836 break;
837 case 2:
838 p2=cl2;
ebda3140 839 index2=clind[2][l2];mrk2=clmark[2][l2];
2755f080 840 break;
841 case 3:
842 p2=cl3;
ebda3140 843 index2=clind[3][l3];mrk2=clmark[3][l3];
2755f080 844 break;
845 case 4:
846 p2=cl4;
ebda3140 847 index2=clind[4][l4];mrk2=clmark[4][l4];
2755f080 848 break;
849 case 5:
850 p2=cl5;
ebda3140 851 index2=clind[5][l5];mrk2=clmark[5][l5];
2755f080 852 break;
5a03f353 853 default:
854 p2=0;
855 index2=-1;mrk2=-1;
856 break;
2755f080 857 }
858
859 Int_t module1 = p1->GetDetectorIndex()+firstmod[firstLay];
5a03f353 860 Int_t layer,ladder,detector;
861 AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
862 Float_t yclu1 = p1->GetY();
863 Float_t zclu1 = p1->GetZ();
864 Double_t cv=0,tgl2=0,phi2=0;
2755f080 865
866 Int_t cln1=mrk1;
2755f080 867 AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(firstLay,cln1);
2755f080 868 x1 = arr1->GetX();
2755f080 869 y1 = arr1->GetY();
2755f080 870 z1 = arr1->GetZ();
2755f080 871 sx1 = arr1->GetSx();
2755f080 872 sy1 = arr1->GetSy();
2755f080 873 sz1 = arr1->GetSz();
225a692b 874
5a03f353 875 if(secondLay>0) {
876 Int_t cln2=mrk2;
877 AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(secondLay,cln2);
878 x2 = arr2->GetX();
879 y2 = arr2->GetY();
880 z2 = arr2->GetZ();
881 sx2 = arr2->GetSx();
882 sy2 = arr2->GetSy();
883 sz2 = arr2->GetSz();
884 cv = Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
885 tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
886 phi2 = TMath::ATan2((y2-y1),(x2-x1));
887 } else { // special case of 1-point tracks, only for cosmics (B=0)
888 x2 = primaryVertex[0];
889 y2 = primaryVertex[1];
890 z2 = primaryVertex[2];
891 cv = 0;
892 tgl2 = (z1-z2)/TMath::Sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
893 phi2 = TMath::ATan2((y1-y2),(x1-x2));
894 }
895
896
2755f080 897 AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
225a692b 898
899
2755f080 900 if(cl5!=0) {
ebda3140 901 trac->AddClusterV2(5,(clind[5][l5] & 0x0fffffff)>>0);
902 trac->AddClusterMark(5,clmark[5][l5]);
b8ed1a92 903 }
904 if(cl4!=0){
ebda3140 905 trac->AddClusterV2(4,(clind[4][l4] & 0x0fffffff)>>0);
906 trac->AddClusterMark(4,clmark[4][l4]);
b8ed1a92 907 }
908 if(cl3!=0){
ebda3140 909 trac->AddClusterV2(3,(clind[3][l3] & 0x0fffffff)>>0);
910 trac->AddClusterMark(3,clmark[3][l3]);
b8ed1a92 911 }
912 if(cl2!=0){
ebda3140 913 trac->AddClusterV2(2,(clind[2][l2] & 0x0fffffff)>>0);
914 trac->AddClusterMark(2,clmark[2][l2]);
b8ed1a92 915 }
916 if(cl1!=0){
ebda3140 917 trac->AddClusterV2(1,(clind[1][l1] & 0x0fffffff)>>0);
918 trac->AddClusterMark(1,clmark[1][l1]);
b8ed1a92 919 }
920 if(cl0!=0){
ebda3140 921 trac->AddClusterV2(0,(clind[0][l0] & 0x0fffffff)>>0);
922 trac->AddClusterMark(0,clmark[0][l0]);
b8ed1a92 923 }
924
b8ed1a92 925 //fit with Kalman filter using AliITStrackerMI::RefitAt()
926 AliITStrackMI* ot = new AliITStrackSA(*trac);
927
928 ot->ResetCovariance(10.);
929 ot->ResetClusters();
930
e50912db 931 if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
b8ed1a92 932 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
933 otrack2->ResetCovariance(10.);
934 otrack2->ResetClusters();
935 //fit from layer 6 to layer 1
e50912db 936 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
b8ed1a92 937 fListOfTracks->AddLast(otrack2);
deae0246 938 fListOfSATracks->AddLast(trac);
b8ed1a92 939 } else {
940 delete otrack2;
deae0246 941 delete trac;
b8ed1a92 942 }
943
944 }
945
946 delete ot;
b8ed1a92 947 }//end loop layer 6
948 }//end loop layer 5
949 }//end loop layer 4
950 }//end loop layer 3
951 }//end loop layer 2
952 }//end loop layer 1
953
b8ed1a92 954
225a692b 955
b8ed1a92 956
deae0246 957 if(fListOfTracks->GetEntries()==0) return 0;
b8ed1a92 958
deae0246 959 Int_t lowchi2 = FindTrackLowChiSquare();
b8ed1a92 960 AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
deae0246 961 AliITStrackSA* trsa = (AliITStrackSA*)fListOfSATracks->At(lowchi2);
b8ed1a92 962
deae0246 963 if(otrack==0) return 0;
964
965 Int_t indexc[AliITSgeomTGeo::kNLayers];
b8ed1a92 966 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
967 for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
968 indexc[nind] = otrack->GetClusterIndex(nind);
969 }
970 Int_t labl[6][3];
971 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
972 if(i<otrack->GetNumberOfClusters()) {
973 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc[i]);
974 labl[i][0]=cl->GetLabel(0);
975 labl[i][1]=cl->GetLabel(1);
976 labl[i][2]=cl->GetLabel(2);
977 } else {
978 labl[i][0]=-1;
979 labl[i][1]=-1;
980 labl[i][2]=-1;
981 }
982 }
2755f080 983
b8ed1a92 984 CookLabel(otrack,0.); //MI change - to see fake ratio
985
986 Int_t label=FindLabel(labl[0][0],labl[1][0],labl[2][0],labl[3][0],labl[4][0],labl[5][0]);
987 Int_t lflag=0;
2755f080 988 for(Int_t i=0;i<otrack->GetNumberOfClusters();i++)
b8ed1a92 989 if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
990
2755f080 991 if(lflag<otrack->GetNumberOfClusters()) label = -label;
b8ed1a92 992 otrack->SetLabel(label);
993
994 //remove clusters of found track
995 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
996 for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
997 Int_t index = trsa->GetClusterMark(nlay,cln);
998 fCluLayer[nlay]->RemoveAt(index);
999 RemoveClusterCoord(nlay,index);
1000 fCluLayer[nlay]->Compress();
7e5bf5af 1001 }
1002 }
7e5bf5af 1003
13918578 1004 return otrack;
1005
1006}
1007
c421499f 1008//_______________________________________________________
1009void AliITStrackerSA::StoreTrack(AliITStrackV2 *t,AliESDEvent *event) const
1010{
1011 //
1012 // Add new track to the ESD
1013 //
1014 AliESDtrack outtrack;
1015 outtrack.UpdateTrackParams(t,AliESDtrack::kITSin);
1016 for(Int_t i=0;i<12;i++) {
1017 outtrack.SetITSModuleIndex(i,t->GetModuleIndex(i));
1018 }
1019 event->AddTrack(&outtrack);
1020
1021 return;
1022}
bef31448 1023
4e05ab9a 1024
bef31448 1025//_______________________________________________________
7e5bf5af 1026Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t /*zvertex*/,Int_t pflag){
13918578 1027 //function used to to find the clusters associated to the track
2755f080 1028
1029 if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(layer)) return 0;
1030
13918578 1031 Int_t nc=0;
1032 AliITSlayer &lay = fgLayers[layer];
7e5bf5af 1033 Double_t r=lay.GetR();
2755f080 1034 if(pflag==1){
cfe729e0 1035 Float_t cx1,cx2,cy1,cy2;
1036 FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
1037 if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
1038 return 0;
1966f03c 1039 Double_t fi1=TMath::ATan2(cy1-fPoint1[1],cx1-fPoint1[0]);
1040 Double_t fi2=TMath::ATan2(cy2-fPoint1[1],cx2-fPoint1[0]);
cfe729e0 1041 fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
1042 }
13918578 1043
7e5bf5af 1044
1045 Int_t ncl = fCluLayer[layer]->GetEntries();
1046 for (Int_t index=0; index<ncl; index++) {
00a7cc50 1047 AliITSRecPoint *c = (AliITSRecPoint*)fCluLayer[layer]->At(index);
7e5bf5af 1048 if (!c) continue;
1049 if (c->GetQ()<=0) continue;
3733ccd2 1050 if(layer>1 && c->GetQ()<=fMinQ) continue;
7e5bf5af 1051
1052 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(layer,index);
1053 Double_t phi = arr->GetPhi();
cfe729e0 1054 if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
494fda7c 1055
7e5bf5af 1056 Double_t lambda = arr->GetLambda();
cfe729e0 1057 if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
1058
7e5bf5af 1059 if(trs->GetNumberOfClustersSA()==trs->GetMaxNumberOfClusters()) return 0;
1060 if(trs->GetNumberOfMarked(layer)==trs->GetMaxNMarkedPerLayer()) return 0;
1061 Int_t orind = arr->GetOrInd();
1062 trs->AddClusterSA(layer,orind);
1063 trs->AddClusterMark(layer,index);
1064
cfe729e0 1065 nc++;
1066 fLambdac=lambda;
1067 fPhiEstimate=phi;
cfe729e0 1068
7e5bf5af 1069 fPointc[0]=arr->GetX();
1070 fPointc[1]=arr->GetY();
1071
cfe729e0 1072 }
1073 return nc;
1074}
13918578 1075
bef31448 1076//________________________________________________________________
1077void AliITStrackerSA::UpdatePoints(){
1078 //update of points for the estimation of the curvature
13918578 1079
bef31448 1080 fPoint2[0]=fPoint3[0];
1081 fPoint2[1]=fPoint3[1];
1082 fPoint3[0]=fPointc[0];
1083 fPoint3[1]=fPointc[1];
13918578 1084
1085
13918578 1086}
13918578 1087
13918578 1088//___________________________________________________________________
1089Int_t AliITStrackerSA::FindEquation(Float_t x1, Float_t y1, Float_t x2, Float_t y2, Float_t x3, Float_t y3,Float_t& a, Float_t& b, Float_t& c){
1090
1091 //given (x,y) of three recpoints (in global coordinates)
1092 //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
1093
1094 Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1095 if(den==0) return 0;
1096 a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1097 b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1098 c = -x1*x1-y1*y1-a*x1-b*y1;
1099 return 1;
1100 }
1101//__________________________________________________________________________
1102 Int_t AliITStrackerSA::FindIntersection(Float_t a1, Float_t b1, Float_t c1, Float_t c2,Float_t& x1,Float_t& y1, Float_t& x2, Float_t& y2){
1103
1104 //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
1105 //c2 is -rlayer*rlayer
1106
1107 if(a1==0) return 0;
3a2f227d 1108 Double_t m = c2-c1;
1109 Double_t aA = (b1*b1)/(a1*a1)+1;
1110 Double_t bB = (-2*m*b1/(a1*a1));
1111 Double_t cC = c2+(m*m)/(a1*a1);
1112 Double_t dD = bB*bB-4*aA*cC;
1113 if(dD<0) return 0;
13918578 1114
3a2f227d 1115 y1 = (-bB+TMath::Sqrt(dD))/(2*aA);
1116 y2 = (-bB-TMath::Sqrt(dD))/(2*aA);
13918578 1117 x1 = (c2-c1-b1*y1)/a1;
1118 x2 = (c2-c1-b1*y2)/a1;
1119
1120 return 1;
1121}
1122//____________________________________________________________________
1123Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t
1124x2,Double_t y2,Double_t x3,Double_t y3){
1125
1126 //calculates the curvature of track
1127 Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1128 if(den==0) return 0;
1129 Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1130 Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1131 Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1132 Double_t xc=-a/2.;
1133
1134 if((a*a+b*b-4*c)<0) return 0;
1135 Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1136 if(rad==0) return 0;
1137
1138 if((x1>0 && y1>0 && x1<xc)) rad*=-1;
1139 if((x1<0 && y1>0 && x1<xc)) rad*=-1;
1140 // if((x1<0 && y1<0 && x1<xc)) rad*=-1;
1141 // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1142
1143 return 1/rad;
1144
1145}
7e5bf5af 1146
1147
13918578 1148//____________________________________________________________________
1149Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1150
1151 //Returns the point closest to pp
1152
1153 Double_t diff1 = p1-pp;
1154 Double_t diff2 = p2-pp;
1155
1156 if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1157 else fPhiEstimate=p2;
1158 return fPhiEstimate;
1159
1160}
1161
1162
1163//_________________________________________________________________
deae0246 1164Int_t AliITStrackerSA::FindTrackLowChiSquare() const {
1165 // returns track with lowest chi square
1166 Int_t dim=fListOfTracks->GetEntries();
1167 if(dim<=1) return 0;
1168 AliITStrackV2* trk = (AliITStrackV2*)fListOfTracks->At(0);
1169 Double_t minChi2=trk->GetChi2();
1170 Int_t index=0;
1171 for(Int_t i=1;i<dim;i++){
4adcf390 1172 trk = (AliITStrackV2*)fListOfTracks->At(i);
deae0246 1173 Double_t chi2=trk->GetChi2();
1174 if(chi2<minChi2){
1175 minChi2=chi2;
1176 index=i;
13918578 1177 }
13918578 1178 }
deae0246 1179 return index;
13918578 1180}
1181
1182//__________________________________________________________
2755f080 1183Int_t AliITStrackerSA::FindLabel(Int_t l0, Int_t l1, Int_t l2, Int_t l3, Int_t l4, Int_t l5){
13918578 1184
1185 //function used to determine the track label
1186
2755f080 1187 Int_t lb[6] = {l0,l1,l2,l3,l4,l5};
13918578 1188 Int_t aa[6]={1,1,1,1,1,1};
1189 Int_t ff=0;
1190 Int_t ll=0;
1191 Int_t k=0;Int_t w=0;Int_t num=6;
8b3779ef 1192 for(Int_t i=5;i>=0;i--) if(lb[i]==-1) num=i;
13918578 1193
1194 while(k<num){
1195
1196 for(Int_t i=k+1;i<num;i++){
1197
1198 if(lb[k]==lb[i] && aa[k]!=0){
1199
1200 aa[k]+=1;
1201 aa[i]=0;
1202 }
1203 }
8b3779ef 1204 k++;
13918578 1205 }
1206
1207 while(w<num){
8b3779ef 1208
13918578 1209 for(Int_t j=0;j<6;j++){
1210 if(aa[w]<aa[j]){
1211 ff=aa[w];
1212 aa[w]=aa[j];
1213 aa[j]=ff;
1214 ll=lb[w];
1215 lb[w]=lb[j];
1216 lb[j]=ll;
1217 }
1218 }
8b3779ef 1219 w++;
13918578 1220 }
8b3779ef 1221
1222 if(num<1) return -1;
1223 return lb[num-1];
13918578 1224}
1225
1226//_____________________________________________________________________________
1227Int_t AliITStrackerSA::Label(Int_t gl1, Int_t gl2, Int_t gl3, Int_t gl4, Int_t gl5, Int_t gl6,Int_t gl7, Int_t gl8, Int_t gl9, Int_t gl10,Int_t gl11,
2755f080 1228Int_t gl12, Int_t gl13, Int_t gl14,Int_t gl15, Int_t gl16, Int_t gl17, Int_t gl18, Int_t minNPoints){
13918578 1229
1230
1231 //function used to assign label to the found track. If track is fake, the label is negative
1232
1233 Int_t lb0[6] = {gl1,gl2,gl3,gl4,gl5,gl6};
1234 Int_t lb1[6] = {gl7,gl8,gl9,gl10,gl11,gl12};
1235 Int_t lb2[6] = {gl13,gl14,gl15,gl16,gl17,gl18};
1236 Int_t ll=FindLabel(lb0[0],lb0[1],lb0[2],lb0[3],lb0[4],lb0[5]);
1237 Int_t lflag=0;Int_t num=6;
1238 if(lb0[5]==-1 && lb1[5]==-1 && lb2[5]==-1) num=5;
1239
1240 for(Int_t i=0;i<num;i++){
1241 if(lb0[i]==ll || lb1[i]==ll || lb2[i]==ll) lflag+=1;
1242 }
1243
2755f080 1244 if(lflag>=minNPoints) return ll;
13918578 1245 else return -ll;
1246
1247
1248}
13918578 1249//_____________________________________________________________________________
5b4dedce 1250void AliITStrackerSA::SetCalculatedWindowSizes(Int_t n, Float_t phimin, Float_t phimax, Float_t lambdamin, Float_t lambdamax){
1251 // Set sizes of the phi and lambda windows used for track finding
1252 fNloop = n;
1253 if(fPhiWin) delete [] fPhiWin;
1254 if(fLambdaWin) delete [] fLambdaWin;
1255 fPhiWin = new Double_t[fNloop];
1256 fLambdaWin = new Double_t[fNloop];
1257 Float_t stepPhi=(phimax-phimin)/(Float_t)(fNloop-1);
1258 Float_t stepLambda=(lambdamax-lambdamin)/(Float_t)(fNloop-1);
1259 for(Int_t k=0;k<fNloop;k++){
1260 Float_t phi=phimin+k*stepPhi;
1261 Float_t lam=lambdamin+k*stepLambda;
1262 fPhiWin[k]=phi;
1263 fLambdaWin[k]=lam;
1264 }
1265}
1266//_____________________________________________________________________________
1267void AliITStrackerSA::SetFixedWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
13918578 1268 // Set sizes of the phi and lambda windows used for track finding
1269 fNloop = n;
1270 if(phi){ // user defined values
1271 fPhiWin = new Double_t[fNloop];
1272 fLambdaWin = new Double_t[fNloop];
1273 for(Int_t k=0;k<fNloop;k++){
1274 fPhiWin[k]=phi[k];
1275 fLambdaWin[k]=lam[k];
1276 }
1277 }
1278 else { // default values
627b6db6 1279
1280 Double_t phid[33] = {0.002,0.003,0.004,0.0045,0.0047,
bef31448 1281 0.005,0.0053,0.0055,
1282 0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1283 0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
627b6db6 1284 0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1285 Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1286 0.005,0.005,0.006,
1287 0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1288 0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1289 0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
bef31448 1290
627b6db6 1291 if(fNloop!=33){
1292 fNloop = 33;
13918578 1293 }
627b6db6 1294
1295
13918578 1296 fPhiWin = new Double_t[fNloop];
1297 fLambdaWin = new Double_t[fNloop];
e340bb86 1298
e340bb86 1299 Double_t factor=AliITSReconstructor::GetRecoParam()->GetFactorSAWindowSizes(); // possibility to enlarge windows for cosmics reco with large misalignments (A.Dainese)
2755f080 1300
1301 for(Int_t k=0;k<fNloop;k++){
e340bb86 1302 fPhiWin[k]=phid[k]*factor;
1303 fLambdaWin[k]=lambdad[k]*factor;
13918578 1304 }
1305
1306 }
1307
1308}
7e5bf5af 1309//_______________________________________________________________________
cc088660 1310void AliITStrackerSA::GetCoorAngles(AliITSRecPoint* cl,Double_t &phi,Double_t &lambda, Float_t &x, Float_t &y,Float_t &z,Double_t* vertex){
7e5bf5af 1311 //Returns values of phi (azimuthal) and lambda angles for a given cluster
cc088660 1312/*
7e5bf5af 1313 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
1314 Int_t lay,lad,det; fGeom->GetModuleId(module,lay,lad,det);
1315 Float_t tx,ty,tz; fGeom->GetTrans(lay,lad,det,tx,ty,tz);
1316
1317 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1318 Double_t phi1=TMath::Pi()/2+alpha;
1319 if (lay==1) phi1+=TMath::Pi();
1320
1321 Float_t cp=TMath::Cos(phi1), sp=TMath::Sin(phi1);
1322 Float_t r=tx*cp+ty*sp;
1323
cc088660 1324 xyz= r*cp - cl->GetY()*sp;
7e5bf5af 1325 y= r*sp + cl->GetY()*cp;
1326 z=cl->GetZ();
cc088660 1327*/
1328 Float_t xyz[3];
1329 cl->GetGlobalXYZ(xyz);
1330 x=xyz[0];
1331 y=xyz[1];
1332 z=xyz[2];
1333
1966f03c 1334 phi=TMath::ATan2(y-vertex[1],x-vertex[0]);
7e5bf5af 1335 lambda=TMath::ATan2(z-vertex[2],TMath::Sqrt((x-vertex[0])*(x-vertex[0])+(y-vertex[1])*(y-vertex[1])));
1336}
13918578 1337
7e5bf5af 1338//________________________________________________________________________
cc088660 1339void AliITStrackerSA::GetCoorErrors(AliITSRecPoint* cl,Float_t &sx,Float_t &sy, Float_t &sz){
7e5bf5af 1340
cc088660 1341 //returns sigmax, y, z of cluster in global coordinates
1342/*
7e5bf5af 1343 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
cc088660 1344 Int_t lay,lad,det;
1345 AliITSgeomTGeo::GetModuleId(module,lay,lad,det);
7e5bf5af 1346
1347 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1348 Double_t phi=TMath::Pi()/2+alpha;
1349 if (lay==1) phi+=TMath::Pi();
1350
1351 Float_t cp=TMath::Cos(phi), sp=TMath::Sin(phi);
cc088660 1352*/
1353 Float_t covm[6];
1354 cl->GetGlobalCov(covm);
1355 sx=TMath::Sqrt(covm[0]);
1356 sy=TMath::Sqrt(covm[3]);
1357 sz=TMath::Sqrt(covm[5]);
1358/*
7e5bf5af 1359 sx = TMath::Sqrt(sp*sp*cl->GetSigmaY2());
1360 sy = TMath::Sqrt(cp*cp*cl->GetSigmaY2());
1361 sz = TMath::Sqrt(cl->GetSigmaZ2());
cc088660 1362*/
7e5bf5af 1363}