1 /**************************************************************************
2 * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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 ////////////////////////////////////////////////////
29 #include <TObjArray.h>
32 #include "AliESDEvent.h"
33 #include "AliESDVertex.h"
34 #include "AliESDtrack.h"
35 #include "AliITSVertexer.h"
36 #include "AliITSclusterTable.h"
37 #include "AliITSRecPoint.h"
38 #include "AliITSgeomTGeo.h"
39 #include "AliITStrackSA.h"
40 #include "AliITStrackerSA.h"
41 #include "AliITSReconstructor.h"
44 ClassImp(AliITStrackerSA)
46 //____________________________________________________________________________
47 AliITStrackerSA::AliITStrackerSA():AliITStrackerMI(),
67 // Default constructor
71 //____________________________________________________________________________
72 AliITStrackerSA::AliITStrackerSA(const Char_t *geom):AliITStrackerMI(0),
93 // Standard constructor (Vertex is known and passed to this obj.)
95 AliWarning("\"geom\" is actually a dummy argument !");
103 //____________________________________________________________________________
104 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliESDVertex *vert):AliITStrackerMI(0),
125 // Standard constructor (Vertex is known and passed to this obj.)
127 AliWarning("\"geom\" is actually a dummy argument !");
133 //____________________________________________________________________________
134 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliITSVertexer *vertexer):AliITStrackerMI(0),
155 // Standard constructor (Vertex is unknown - vertexer is passed to this obj)
157 AliWarning("\"geom\" is actually a dummy argument !");
160 fVertexer = vertexer;
164 //____________________________________________________________________________
165 AliITStrackerSA::AliITStrackerSA(const AliITStrackerSA& tracker):AliITStrackerMI(),
166 fPhiEstimate(tracker.fPhiEstimate),
167 fITSStandAlone(tracker.fITSStandAlone),
168 fLambdac(tracker.fLambdac),
169 fPhic(tracker.fPhic),
170 fCoef1(tracker.fCoef1),
171 fCoef2(tracker.fCoef2),
172 fCoef3(tracker.fCoef3),
173 fNloop(tracker.fNloop),
174 fPhiWin(tracker.fPhiWin),
175 fLambdaWin(tracker.fLambdaWin),
176 fVert(tracker.fVert),
177 fVertexer(tracker.fVertexer),
178 fListOfTracks(tracker.fListOfTracks),
179 fListOfSATracks(tracker.fListOfSATracks),
180 fITSclusters(tracker.fITSclusters),
181 fSixPoints(tracker.fSixPoints),
182 fOuterStartLayer(tracker.fOuterStartLayer),
183 fCluLayer(tracker.fCluLayer),
184 fCluCoord(tracker.fCluCoord) {
186 for(Int_t i=0;i<2;i++){
187 fPoint1[i]=tracker.fPoint1[i];
188 fPoint2[i]=tracker.fPoint2[i];
189 fPoint3[i]=tracker.fPoint3[i];
190 fPointc[i]=tracker.fPointc[i];
192 if(tracker.fVertexer && tracker.fVert){
193 fVert = new AliESDVertex(*tracker.fVert);
196 fVert = tracker.fVert;
198 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
199 fCluLayer[i] = tracker.fCluLayer[i];
200 fCluCoord[i] = tracker.fCluCoord[i];
203 //______________________________________________________________________
204 AliITStrackerSA& AliITStrackerSA::operator=(const AliITStrackerSA& source){
205 // Assignment operator.
206 this->~AliITStrackerSA();
207 new(this) AliITStrackerSA(source);
212 //____________________________________________________________________________
213 AliITStrackerSA::~AliITStrackerSA(){
215 // if fVertexer is not null, the AliESDVertex obj. is owned by this class
216 // and is deleted here
218 if(fVert)delete fVert;
223 if(fPhiWin)delete []fPhiWin;
224 if(fLambdaWin)delete []fLambdaWin;
225 fListOfTracks->Delete();
226 delete fListOfTracks;
227 fListOfSATracks->Delete();
228 delete fListOfSATracks;
230 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
232 fCluLayer[i]->Delete();
239 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
241 fCluCoord[i]->Delete();
250 //____________________________________________________________________________
251 Int_t AliITStrackerSA::Clusters2Tracks(AliESDEvent *event){
252 // This method is used to find and fit the tracks. By default the corresponding
253 // method in the parent class is invoked. In this way a combined tracking
254 // TPC+ITS is performed. If the flag fITSStandAlone is true, the tracking
255 // is done in the ITS only. In the standard reconstruction chain this option
256 // can be set via AliReconstruction::SetOption("ITS","onlyITS")
259 rc=AliITStrackerMI::Clusters2Tracks(event);
262 AliDebug(1,"Stand Alone flag set: doing tracking in ITS alone\n");
264 if(!rc) rc=FindTracks(event);
268 //____________________________________________________________________________
269 void AliITStrackerSA::Init(){
270 // Reset all data members
272 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
282 Int_t nLoops=AliITSReconstructor::GetRecoParam()->GetNLoopsSA();
284 SetFixedWindowSizes();
286 Double_t phimin=AliITSReconstructor::GetRecoParam()->GetMinPhiSA();
287 Double_t phimax=AliITSReconstructor::GetRecoParam()->GetMaxPhiSA();
288 Double_t lambmin=AliITSReconstructor::GetRecoParam()->GetMinLambdaSA();
289 Double_t lambmax=AliITSReconstructor::GetRecoParam()->GetMaxLambdaSA();
290 SetCalculatedWindowSizes(nLoops,phimin,phimax,lambmin,lambmax);
294 SetOuterStartLayer(0);
296 fListOfTracks=new TObjArray(0,0);
297 fListOfSATracks=new TObjArray(0,0);
301 //_______________________________________________________________________
302 void AliITStrackerSA::ResetForFinding(){
303 // Reset data members used in all loops during track finding
305 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
313 fListOfTracks->Delete();
314 fListOfSATracks->Delete();
319 //______________________________________________________________________
320 Int_t AliITStrackerSA::FindTracks(AliESDEvent* event){
322 // Track finder using the ESD object
327 Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
330 //Reads event and mark clusters of traks already found, with flag kITSin
331 Int_t nentr=event->GetNumberOfTracks();
332 if(AliITSReconstructor::GetRecoParam()->GetSAUseAllClusters()==kFALSE) {
334 AliESDtrack *track=event->GetTrack(nentr);
335 if (track->GetStatus()&AliESDtrack::kITSin==AliESDtrack::kITSin){
337 Int_t ncl = track->GetITSclusters(idx);
338 for(Int_t k=0;k<ncl;k++){
339 AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
340 cll->SetBit(kSAflag);
346 Double_t primaryVertex[3];
347 event->GetVertex()->GetXYZ(primaryVertex);
348 //Creates TClonesArray with clusters for each layer. The clusters already used
349 //by AliITStrackerMI are not considered
350 Int_t nclusters[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
351 Int_t dmar[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
352 if (fCluLayer == 0) {
353 fCluLayer = new TClonesArray*[AliITSgeomTGeo::kNLayers];
354 fCluCoord = new TClonesArray*[AliITSgeomTGeo::kNLayers];
355 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
360 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
361 AliITSlayer &layer=fgLayers[i];
362 if (!AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) {
363 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
364 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
365 if(cls->TestBit(kSAflag)==kTRUE) continue; //clusters used by TPC prol.
366 if(cls->GetQ()==0) continue; //fake clusters dead zones
372 fCluLayer[i] = new TClonesArray("AliITSRecPoint",nclusters[i]);
374 fCluCoord[i] = new TClonesArray("AliITSclusterTable",nclusters[i]);
377 for(Int_t ilay=0;ilay<AliITSgeomTGeo::GetNLayers();ilay++){
378 TClonesArray &clulay = *fCluLayer[ilay];
379 TClonesArray &clucoo = *fCluCoord[ilay];
380 AliITSlayer &layer=fgLayers[ilay];
381 if (!AliITSReconstructor::GetRecoParam()->GetLayersToSkip(ilay)) {
382 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
383 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
384 if(cls->TestBit(kSAflag)==kTRUE) continue;
385 if(cls->GetQ()==0) continue;
386 Double_t phi=0;Double_t lambda=0;
387 Float_t x=0;Float_t y=0;Float_t z=0;
388 Float_t sx=0;Float_t sy=0;Float_t sz=0;
389 GetCoorAngles(cls,phi,lambda,x,y,z,primaryVertex);
390 GetCoorErrors(cls,sx,sy,sz);
391 new (clulay[dmar[ilay]]) AliITSRecPoint(*cls);
392 new (clucoo[dmar[ilay]]) AliITSclusterTable(x,y,z,sx,sy,sz,phi,lambda,cli);
401 Int_t minNPoints = (fSixPoints ? AliITSgeomTGeo::GetNLayers() : AliITSgeomTGeo::GetNLayers()-1);
402 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
403 if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) {
408 static Int_t nClusLay[AliITSgeomTGeo::kNLayers];//counter for clusters on each layer
410 //loop on different minNPoints
411 Int_t minMinNPoints=minNPoints;
412 if(AliITSReconstructor::GetRecoParam()->GetSAOnePointTracks()) minMinNPoints=2;
413 for(Int_t iMinNPoints=minNPoints; iMinNPoints>=minMinNPoints; iMinNPoints--) {
414 //loop on the different windows
415 for(Int_t nloop=0;nloop<fNloop;nloop++){
416 for(Int_t ncl=0;ncl<fCluLayer[0]->GetEntries();ncl++){ //loop starting from layer 0
421 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[0]->At(ncl);
424 if (cl->GetQ()<=0) continue;
426 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(0,ncl);
427 fPhic = arr->GetPhi();
428 fLambdac = arr->GetLambda();
429 if (TMath::Abs(fLambdac)>0.26*TMath::Pi()) continue;
430 fPhiEstimate = fPhic;
431 AliITStrackSA* trs = new AliITStrackSA();
432 fPoint1[0]=primaryVertex[0];
433 fPoint1[1]=primaryVertex[1];
436 fPoint2[0]=arr->GetX();
437 fPoint2[1]=arr->GetY();
438 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) nClusLay[i]=0;
439 nClusLay[0] = SearchClusters(0,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
440 nClusLay[1] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
443 fPoint3[0] = fPointc[0];
444 fPoint3[1] = fPointc[1];
446 nClusLay[2] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
447 if(nClusLay[1]==0 && nClusLay[2]==0) pflag=0;
448 if(nClusLay[2]!=0 && nClusLay[1]!=0){ pflag=1; UpdatePoints();}
449 if(nClusLay[2]!=0 && nClusLay[1]==0){
451 fPoint3[0]=fPointc[0];
452 fPoint3[1]=fPointc[1];
455 nClusLay[3] = SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
457 if(nClusLay[3]!=0) UpdatePoints();
458 nClusLay[4] = SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
460 if(nClusLay[4]!=0) UpdatePoints();
461 nClusLay[5] = SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
465 //check of the candidate track
466 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers();nnp++) {
467 if(nClusLay[nnp]!=0) layOK+=1;
470 if(layOK>=iMinNPoints) {
471 //printf("-NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]);
472 AliITStrackV2* tr2 = 0;
473 tr2 = FitTrack(trs,primaryVertex);
475 //printf("-NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
478 AliESDtrack outtrack;
479 outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
480 event->AddTrack(&outtrack);
484 }//end loop on clusters of layer1
485 }//end loop on windows
486 }//end loop on min points
491 //if 5/6 points are required, second loop starting
492 //from second layer (SPD2), to find tracks with point of
493 //layer 1 missing. Not done for cosmics.
494 if(!fSixPoints && fOuterStartLayer==0) {
495 //printf("looking from SPD2\n");
496 // counter for clusters on each layer
497 for(Int_t nloop=0;nloop<fNloop;nloop++){
498 Int_t ncl2=fCluLayer[1]->GetEntries();
499 while(ncl2--){ //loop starting from layer 2
502 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[1]->At(ncl2);
505 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(1,ncl2);
506 fPhic = arr->GetPhi();
507 fLambdac = arr->GetLambda();
508 fPhiEstimate = fPhic;
510 AliITStrackSA* trs = new AliITStrackSA();
511 fPoint1[0]=primaryVertex[0];
512 fPoint1[1]=primaryVertex[1];
514 fPoint2[0]=arr->GetX();
515 fPoint2[1]=arr->GetY();
516 for(Int_t kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
517 nClusLay[0] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],
518 trs,primaryVertex[2],pflag);
519 nClusLay[1] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],
520 trs,primaryVertex[2],pflag);
523 fPoint3[0]=fPointc[0];
524 fPoint3[1]=fPointc[1];
526 nClusLay[2]= SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],
527 trs,primaryVertex[2],pflag);
532 nClusLay[3]= SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],
533 trs,primaryVertex[2],pflag);
538 nClusLay[4]=SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],
539 trs,primaryVertex[2],pflag);
542 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-1;nnp++){
543 if(nClusLay[nnp]!=0) layOK+=1;
545 if(layOK>=minNPoints){ // 5/6
546 //printf("--NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]);
548 AliITStrackV2* tr2 = 0;
549 tr2 = FitTrack(trs,primaryVertex);
551 //printf("--NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
553 AliESDtrack outtrack;
554 outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
555 event->AddTrack(&outtrack);
561 }//end loop on clusters of layer2
566 // search for tracks starting from SPD2, SDD1, SDD2, SSD2
567 // for cosmics (A. Dainese 31.07.07)
568 if(fOuterStartLayer>0 && !AliITSReconstructor::GetRecoParam()->GetSAOnePointTracks()) {
569 for(Int_t innLay=1; innLay<=fOuterStartLayer; innLay++) {
570 //printf("Searching from layer %d outward\n",innLay);
571 minNPoints=AliITSgeomTGeo::GetNLayers()-innLay;
572 for(Int_t i=innLay;i<AliITSgeomTGeo::GetNLayers();i++)
573 if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i))
575 // counter for clusters on each layer
577 for(Int_t nloop=0;nloop<fNloop;nloop++){
578 Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
579 while(nclInnLay--){ //loop starting from layer innLay
582 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[innLay]->At(nclInnLay);
585 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(innLay,nclInnLay);
586 fPhic = arr->GetPhi();
587 fLambdac = arr->GetLambda();
588 fPhiEstimate = fPhic;
590 AliITStrackSA* trs = new AliITStrackSA();
591 fPoint1[0]=primaryVertex[0];
592 fPoint1[1]=primaryVertex[1];
593 fPoint2[0]=arr->GetX();
594 fPoint2[1]=arr->GetY();
597 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
600 nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
601 trs,primaryVertex[2],pflag);
602 for(Int_t nextLay=innLay+1; nextLay<AliITSgeomTGeo::GetNLayers(); nextLay++) {
604 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
605 trs,primaryVertex[2],pflag);
609 fPoint3[0]=fPointc[0];
610 fPoint3[1]=fPointc[1];
618 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
619 if(nClusLay[nnp]!=0) layOK+=1;
621 if(layOK>=minNPoints){
622 //printf("---NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]);
623 AliITStrackV2* tr2 = 0;
624 tr2 = FitTrack(trs,primaryVertex);
626 //printf("---NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
628 AliESDtrack outtrack;
629 outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
630 event->AddTrack(&outtrack);
636 }//end loop on clusters of innLay
637 } //end loop on window sizes
639 } //end loop on innLay
640 } //end if(fOuterStartLayer>0)
643 // search for 1-point tracks in SPD, only for cosmics
644 // (A.Dainese 21.03.08)
645 if(AliITSReconstructor::GetRecoParam()->GetSAOnePointTracks() &&
646 TMath::Abs(event->GetMagneticField())<0.01) {
647 Int_t outerLayer=1; // only SPD
648 for(Int_t innLay=0; innLay<=TMath::Min(1,fOuterStartLayer); innLay++) {
649 // counter for clusters on each layer
651 for(Int_t nloop=0;nloop<fNloop;nloop++){
652 Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
653 while(nclInnLay--){ //loop starting from layer innLay
656 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[innLay]->At(nclInnLay);
659 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(innLay,nclInnLay);
660 fPhic = arr->GetPhi();
661 fLambdac = arr->GetLambda();
662 fPhiEstimate = fPhic;
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();
671 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
674 nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
675 trs,primaryVertex[2],pflag);
676 for(Int_t nextLay=innLay+1; nextLay<=outerLayer; nextLay++) {
678 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
679 trs,primaryVertex[2],pflag);
683 fPoint3[0]=fPointc[0];
684 fPoint3[1]=fPointc[1];
692 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
693 if(nClusLay[nnp]!=0) layOK+=1;
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]);
697 AliITStrackV2* tr2 = 0;
698 Bool_t onePoint = kTRUE;
699 tr2 = FitTrack(trs,primaryVertex,onePoint);
701 //printf("----NPOINTS fit: %d\n",tr2->GetNumberOfClusters());
703 AliESDtrack outtrack;
704 outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
705 event->AddTrack(&outtrack);
711 }//end loop on clusters of innLay
712 } //end loop on window sizes
714 } //end loop on innLay
715 } // end search 1-point tracks
717 Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
722 //________________________________________________________________________
724 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex,Bool_t onePoint) {
725 //fit of the found track (most general case, also <6 points, layers missing)
726 // A.Dainese 16.11.07
729 const Int_t kMaxClu=AliITStrackSA::kMaxNumberOfClusters;
731 static Int_t firstmod[AliITSgeomTGeo::kNLayers];
733 static Int_t clind[AliITSgeomTGeo::kNLayers][kMaxClu];
734 static Int_t clmark[AliITSgeomTGeo::kNLayers][kMaxClu];
735 static Int_t nnn[AliITSgeomTGeo::kNLayers];
736 static Int_t kkk[AliITSgeomTGeo::kNLayers];
737 static Int_t end[AliITSgeomTGeo::kNLayers];
738 static AliITSRecPoint *listlayer[AliITSgeomTGeo::kNLayers][kMaxClu];
740 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
741 firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
745 for(Int_t j=0;j<kMaxClu; j++){
753 Int_t nclusters = tr->GetNumberOfClustersSA();
754 for(Int_t ncl=0;ncl<nclusters;ncl++){
755 Int_t index = tr->GetClusterIndexSA(ncl);
756 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
757 if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
758 Int_t lay = (index & 0xf0000000) >> 28;
759 Int_t nInLay=end[lay];
760 listlayer[lay][nInLay]=cl;
762 clind[lay][ind]=index;
766 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
767 for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
768 Int_t mark = tr->GetClusterMark(nlay,ncl);
770 clmark[nlay][ind]=mark;
775 Int_t firstLay=-1,secondLay=-1;
776 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
782 } else if(secondLay==-1) {
788 if(firstLay==-1 || (secondLay==-1 && !onePoint)) return 0;
790 for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
791 AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0][l0];
792 for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
793 AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1][l1];
794 for(Int_t l2=0;l2<end[2];l2++){ //loop on layer 3
795 AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2][l2];
796 for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4
797 AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3][l3];
798 for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
799 AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4][l4];
800 for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6
801 AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5][l5];
804 Double_t x1,y1,z1,sx1,sy1,sz1;
805 Double_t x2,y2,z2,sx2,sy2,sz2;
806 AliITSRecPoint* p1=0;
807 AliITSRecPoint* p2=0;
808 Int_t index1=0,index2=0;
814 index1=clind[0][l0];mrk1=clmark[0][l0];
818 index1=clind[1][l1];mrk1=clmark[1][l1];
822 index1=clind[2][l2];mrk1=clmark[2][l2];
826 index1=clind[3][l3];mrk1=clmark[3][l3];
830 index1=clind[4][l4];mrk1=clmark[4][l4];
837 index2=clind[1][l1];mrk2=clmark[1][l1];
841 index2=clind[2][l2];mrk2=clmark[2][l2];
845 index2=clind[3][l3];mrk2=clmark[3][l3];
849 index2=clind[4][l4];mrk2=clmark[4][l4];
853 index2=clind[5][l5];mrk2=clmark[5][l5];
861 Int_t module1 = p1->GetDetectorIndex()+firstmod[firstLay];
862 Int_t layer,ladder,detector;
863 AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
864 Float_t yclu1 = p1->GetY();
865 Float_t zclu1 = p1->GetZ();
866 Double_t cv=0,tgl2=0,phi2=0;
869 AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(firstLay,cln1);
879 AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(secondLay,cln2);
886 cv = Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
887 tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
888 phi2 = TMath::ATan2((y2-y1),(x2-x1));
889 } else { // special case of 1-point tracks, only for cosmics (B=0)
890 x2 = primaryVertex[0];
891 y2 = primaryVertex[1];
892 z2 = primaryVertex[2];
894 tgl2 = (z1-z2)/TMath::Sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
895 phi2 = TMath::ATan2((y1-y2),(x1-x2));
899 AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
903 trac->AddClusterV2(5,(clind[5][l5] & 0x0fffffff)>>0);
904 trac->AddClusterMark(5,clmark[5][l5]);
907 trac->AddClusterV2(4,(clind[4][l4] & 0x0fffffff)>>0);
908 trac->AddClusterMark(4,clmark[4][l4]);
911 trac->AddClusterV2(3,(clind[3][l3] & 0x0fffffff)>>0);
912 trac->AddClusterMark(3,clmark[3][l3]);
915 trac->AddClusterV2(2,(clind[2][l2] & 0x0fffffff)>>0);
916 trac->AddClusterMark(2,clmark[2][l2]);
919 trac->AddClusterV2(1,(clind[1][l1] & 0x0fffffff)>>0);
920 trac->AddClusterMark(1,clmark[1][l1]);
923 trac->AddClusterV2(0,(clind[0][l0] & 0x0fffffff)>>0);
924 trac->AddClusterMark(0,clmark[0][l0]);
927 //fit with Kalman filter using AliITStrackerMI::RefitAt()
928 AliITStrackMI* ot = new AliITStrackSA(*trac);
930 ot->ResetCovariance(10.);
933 if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
934 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
935 otrack2->ResetCovariance(10.);
936 otrack2->ResetClusters();
937 //fit from layer 6 to layer 1
938 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
939 fListOfTracks->AddLast(otrack2);
940 fListOfSATracks->AddLast(trac);
959 if(fListOfTracks->GetEntries()==0) return 0;
961 Int_t lowchi2 = FindTrackLowChiSquare();
962 AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
963 AliITStrackSA* trsa = (AliITStrackSA*)fListOfSATracks->At(lowchi2);
965 if(otrack==0) return 0;
967 Int_t indexc[AliITSgeomTGeo::kNLayers];
968 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
969 for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
970 indexc[nind] = otrack->GetClusterIndex(nind);
973 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
974 if(i<otrack->GetNumberOfClusters()) {
975 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc[i]);
976 labl[i][0]=cl->GetLabel(0);
977 labl[i][1]=cl->GetLabel(1);
978 labl[i][2]=cl->GetLabel(2);
986 CookLabel(otrack,0.); //MI change - to see fake ratio
988 Int_t label=FindLabel(labl[0][0],labl[1][0],labl[2][0],labl[3][0],labl[4][0],labl[5][0]);
990 for(Int_t i=0;i<otrack->GetNumberOfClusters();i++)
991 if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
993 if(lflag<otrack->GetNumberOfClusters()) label = -label;
994 otrack->SetLabel(label);
996 //remove clusters of found track
997 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
998 for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
999 Int_t index = trsa->GetClusterMark(nlay,cln);
1000 fCluLayer[nlay]->RemoveAt(index);
1001 RemoveClusterCoord(nlay,index);
1002 fCluLayer[nlay]->Compress();
1012 //_______________________________________________________
1013 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t /*zvertex*/,Int_t pflag){
1014 //function used to to find the clusters associated to the track
1016 if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(layer)) return 0;
1019 AliITSlayer &lay = fgLayers[layer];
1020 Double_t r=lay.GetR();
1022 Float_t cx1,cx2,cy1,cy2;
1023 FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
1024 if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
1026 Double_t fi1=TMath::ATan2(cy1-fPoint1[1],cx1-fPoint1[0]);
1027 Double_t fi2=TMath::ATan2(cy2-fPoint1[1],cx2-fPoint1[0]);
1028 fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
1032 Int_t ncl = fCluLayer[layer]->GetEntries();
1033 for (Int_t index=0; index<ncl; index++) {
1034 AliITSRecPoint *c = (AliITSRecPoint*)fCluLayer[layer]->At(index);
1036 if (c->GetQ()<=0) continue;
1038 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(layer,index);
1039 Double_t phi = arr->GetPhi();
1040 if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
1042 Double_t lambda = arr->GetLambda();
1043 if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
1045 if(trs->GetNumberOfClustersSA()==trs->GetMaxNumberOfClusters()) return 0;
1046 if(trs->GetNumberOfMarked(layer)==trs->GetMaxNMarkedPerLayer()) return 0;
1047 Int_t orind = arr->GetOrInd();
1048 trs->AddClusterSA(layer,orind);
1049 trs->AddClusterMark(layer,index);
1055 fPointc[0]=arr->GetX();
1056 fPointc[1]=arr->GetY();
1062 //________________________________________________________________
1063 void AliITStrackerSA::UpdatePoints(){
1064 //update of points for the estimation of the curvature
1066 fPoint2[0]=fPoint3[0];
1067 fPoint2[1]=fPoint3[1];
1068 fPoint3[0]=fPointc[0];
1069 fPoint3[1]=fPointc[1];
1074 //___________________________________________________________________
1075 Int_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){
1077 //given (x,y) of three recpoints (in global coordinates)
1078 //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
1080 Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1081 if(den==0) return 0;
1082 a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1083 b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1084 c = -x1*x1-y1*y1-a*x1-b*y1;
1087 //__________________________________________________________________________
1088 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){
1090 //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
1091 //c2 is -rlayer*rlayer
1095 Double_t aA = (b1*b1)/(a1*a1)+1;
1096 Double_t bB = (-2*m*b1/(a1*a1));
1097 Double_t cC = c2+(m*m)/(a1*a1);
1098 Double_t dD = bB*bB-4*aA*cC;
1101 y1 = (-bB+TMath::Sqrt(dD))/(2*aA);
1102 y2 = (-bB-TMath::Sqrt(dD))/(2*aA);
1103 x1 = (c2-c1-b1*y1)/a1;
1104 x2 = (c2-c1-b1*y2)/a1;
1108 //____________________________________________________________________
1109 Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t
1110 x2,Double_t y2,Double_t x3,Double_t y3){
1112 //calculates the curvature of track
1113 Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1114 if(den==0) return 0;
1115 Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1116 Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1117 Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1120 if((a*a+b*b-4*c)<0) return 0;
1121 Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1122 if(rad==0) return 0;
1124 if((x1>0 && y1>0 && x1<xc)) rad*=-1;
1125 if((x1<0 && y1>0 && x1<xc)) rad*=-1;
1126 // if((x1<0 && y1<0 && x1<xc)) rad*=-1;
1127 // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1134 //____________________________________________________________________
1135 Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1137 //Returns the point closest to pp
1139 Double_t diff1 = p1-pp;
1140 Double_t diff2 = p2-pp;
1142 if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1143 else fPhiEstimate=p2;
1144 return fPhiEstimate;
1149 //_________________________________________________________________
1150 Int_t AliITStrackerSA::FindTrackLowChiSquare() const {
1151 // returns track with lowest chi square
1152 Int_t dim=fListOfTracks->GetEntries();
1153 if(dim<=1) return 0;
1154 AliITStrackV2* trk = (AliITStrackV2*)fListOfTracks->At(0);
1155 Double_t minChi2=trk->GetChi2();
1157 for(Int_t i=1;i<dim;i++){
1158 trk = (AliITStrackV2*)fListOfTracks->At(i);
1159 Double_t chi2=trk->GetChi2();
1168 //__________________________________________________________
1169 Int_t AliITStrackerSA::FindLabel(Int_t l0, Int_t l1, Int_t l2, Int_t l3, Int_t l4, Int_t l5){
1171 //function used to determine the track label
1173 Int_t lb[6] = {l0,l1,l2,l3,l4,l5};
1174 Int_t aa[6]={1,1,1,1,1,1};
1177 Int_t k=0;Int_t w=0;Int_t num=6;
1178 for(Int_t i=5;i>=0;i--) if(lb[i]==-1) num=i;
1182 for(Int_t i=k+1;i<num;i++){
1184 if(lb[k]==lb[i] && aa[k]!=0){
1195 for(Int_t j=0;j<6;j++){
1208 if(num<1) return -1;
1212 //_____________________________________________________________________________
1213 Int_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,
1214 Int_t gl12, Int_t gl13, Int_t gl14,Int_t gl15, Int_t gl16, Int_t gl17, Int_t gl18, Int_t minNPoints){
1217 //function used to assign label to the found track. If track is fake, the label is negative
1219 Int_t lb0[6] = {gl1,gl2,gl3,gl4,gl5,gl6};
1220 Int_t lb1[6] = {gl7,gl8,gl9,gl10,gl11,gl12};
1221 Int_t lb2[6] = {gl13,gl14,gl15,gl16,gl17,gl18};
1222 Int_t ll=FindLabel(lb0[0],lb0[1],lb0[2],lb0[3],lb0[4],lb0[5]);
1223 Int_t lflag=0;Int_t num=6;
1224 if(lb0[5]==-1 && lb1[5]==-1 && lb2[5]==-1) num=5;
1226 for(Int_t i=0;i<num;i++){
1227 if(lb0[i]==ll || lb1[i]==ll || lb2[i]==ll) lflag+=1;
1230 if(lflag>=minNPoints) return ll;
1235 //_____________________________________________________________________________
1236 void AliITStrackerSA::SetCalculatedWindowSizes(Int_t n, Float_t phimin, Float_t phimax, Float_t lambdamin, Float_t lambdamax){
1237 // Set sizes of the phi and lambda windows used for track finding
1239 if(fPhiWin) delete [] fPhiWin;
1240 if(fLambdaWin) delete [] fLambdaWin;
1241 fPhiWin = new Double_t[fNloop];
1242 fLambdaWin = new Double_t[fNloop];
1243 Float_t stepPhi=(phimax-phimin)/(Float_t)(fNloop-1);
1244 Float_t stepLambda=(lambdamax-lambdamin)/(Float_t)(fNloop-1);
1245 for(Int_t k=0;k<fNloop;k++){
1246 Float_t phi=phimin+k*stepPhi;
1247 Float_t lam=lambdamin+k*stepLambda;
1252 //_____________________________________________________________________________
1253 void AliITStrackerSA::SetFixedWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
1254 // Set sizes of the phi and lambda windows used for track finding
1256 if(phi){ // user defined values
1257 fPhiWin = new Double_t[fNloop];
1258 fLambdaWin = new Double_t[fNloop];
1259 for(Int_t k=0;k<fNloop;k++){
1261 fLambdaWin[k]=lam[k];
1264 else { // default values
1266 Double_t phid[33] = {0.002,0.003,0.004,0.0045,0.0047,
1267 0.005,0.0053,0.0055,
1268 0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1269 0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
1270 0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1271 Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1273 0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1274 0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1275 0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
1282 fPhiWin = new Double_t[fNloop];
1283 fLambdaWin = new Double_t[fNloop];
1285 Double_t factor=AliITSReconstructor::GetRecoParam()->GetFactorSAWindowSizes(); // possibility to enlarge windows for cosmics reco with large misalignments (A.Dainese)
1287 for(Int_t k=0;k<fNloop;k++){
1288 fPhiWin[k]=phid[k]*factor;
1289 fLambdaWin[k]=lambdad[k]*factor;
1295 //_______________________________________________________________________
1296 void AliITStrackerSA::GetCoorAngles(AliITSRecPoint* cl,Double_t &phi,Double_t &lambda, Float_t &x, Float_t &y,Float_t &z,Double_t* vertex){
1297 //Returns values of phi (azimuthal) and lambda angles for a given cluster
1299 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
1300 Int_t lay,lad,det; fGeom->GetModuleId(module,lay,lad,det);
1301 Float_t tx,ty,tz; fGeom->GetTrans(lay,lad,det,tx,ty,tz);
1303 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1304 Double_t phi1=TMath::Pi()/2+alpha;
1305 if (lay==1) phi1+=TMath::Pi();
1307 Float_t cp=TMath::Cos(phi1), sp=TMath::Sin(phi1);
1308 Float_t r=tx*cp+ty*sp;
1310 xyz= r*cp - cl->GetY()*sp;
1311 y= r*sp + cl->GetY()*cp;
1315 cl->GetGlobalXYZ(xyz);
1320 phi=TMath::ATan2(y-vertex[1],x-vertex[0]);
1321 lambda=TMath::ATan2(z-vertex[2],TMath::Sqrt((x-vertex[0])*(x-vertex[0])+(y-vertex[1])*(y-vertex[1])));
1324 //________________________________________________________________________
1325 void AliITStrackerSA::GetCoorErrors(AliITSRecPoint* cl,Float_t &sx,Float_t &sy, Float_t &sz){
1327 //returns sigmax, y, z of cluster in global coordinates
1329 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
1331 AliITSgeomTGeo::GetModuleId(module,lay,lad,det);
1333 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1334 Double_t phi=TMath::Pi()/2+alpha;
1335 if (lay==1) phi+=TMath::Pi();
1337 Float_t cp=TMath::Cos(phi), sp=TMath::Sin(phi);
1340 cl->GetGlobalCov(covm);
1341 sx=TMath::Sqrt(covm[0]);
1342 sy=TMath::Sqrt(covm[3]);
1343 sz=TMath::Sqrt(covm[5]);
1345 sx = TMath::Sqrt(sp*sp*cl->GetSigmaY2());
1346 sy = TMath::Sqrt(cp*cp*cl->GetSigmaY2());
1347 sz = TMath::Sqrt(cl->GetSigmaZ2());