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"
45 ClassImp(AliITStrackerSA)
47 //____________________________________________________________________________
48 AliITStrackerSA::AliITStrackerSA():AliITStrackerMI(),
71 // Default constructor
75 //____________________________________________________________________________
76 AliITStrackerSA::AliITStrackerSA(const Char_t *geom):AliITStrackerMI(0),
100 // Standard constructor (Vertex is known and passed to this obj.)
102 AliWarning("\"geom\" is actually a dummy argument !");
110 //____________________________________________________________________________
111 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliESDVertex *vert):AliITStrackerMI(0),
135 // Standard constructor (Vertex is known and passed to this obj.)
137 AliWarning("\"geom\" is actually a dummy argument !");
143 //____________________________________________________________________________
144 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliITSVertexer *vertexer):AliITStrackerMI(0),
168 // Standard constructor (Vertex is unknown - vertexer is passed to this obj)
170 AliWarning("\"geom\" is actually a dummy argument !");
173 fVertexer = vertexer;
177 //____________________________________________________________________________
178 AliITStrackerSA::AliITStrackerSA(const AliITStrackerSA& tracker):AliITStrackerMI(),
179 fPhiEstimate(tracker.fPhiEstimate),
180 fITSStandAlone(tracker.fITSStandAlone),
181 fLambdac(tracker.fLambdac),
182 fPhic(tracker.fPhic),
183 fCoef1(tracker.fCoef1),
184 fCoef2(tracker.fCoef2),
185 fCoef3(tracker.fCoef3),
186 fNloop(tracker.fNloop),
187 fPhiWin(tracker.fPhiWin),
188 fLambdaWin(tracker.fLambdaWin),
189 fVert(tracker.fVert),
190 fVertexer(tracker.fVertexer),
191 fListOfTracks(tracker.fListOfTracks),
192 fListOfSATracks(tracker.fListOfSATracks),
193 fITSclusters(tracker.fITSclusters),
194 fInwardFlag(tracker.fInwardFlag),
195 fOuterStartLayer(tracker.fOuterStartLayer),
196 fInnerStartLayer(tracker.fInnerStartLayer),
197 fMinNPoints(tracker.fMinNPoints),
198 fMinQ(tracker.fMinQ),
199 fCluLayer(tracker.fCluLayer),
200 fCluCoord(tracker.fCluCoord) {
202 for(Int_t i=0;i<2;i++){
203 fPoint1[i]=tracker.fPoint1[i];
204 fPoint2[i]=tracker.fPoint2[i];
205 fPoint3[i]=tracker.fPoint3[i];
206 fPointc[i]=tracker.fPointc[i];
208 if(tracker.fVertexer && tracker.fVert){
209 fVert = new AliESDVertex(*tracker.fVert);
212 fVert = tracker.fVert;
214 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
215 fCluLayer[i] = tracker.fCluLayer[i];
216 fCluCoord[i] = tracker.fCluCoord[i];
219 //______________________________________________________________________
220 AliITStrackerSA& AliITStrackerSA::operator=(const AliITStrackerSA& source){
221 // Assignment operator.
222 this->~AliITStrackerSA();
223 new(this) AliITStrackerSA(source);
228 //____________________________________________________________________________
229 AliITStrackerSA::~AliITStrackerSA(){
231 // if fVertexer is not null, the AliESDVertex obj. is owned by this class
232 // and is deleted here
234 if(fVert)delete fVert;
239 if(fPhiWin)delete []fPhiWin;
240 if(fLambdaWin)delete []fLambdaWin;
241 fListOfTracks->Delete();
242 delete fListOfTracks;
243 fListOfSATracks->Delete();
244 delete fListOfSATracks;
246 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
248 fCluLayer[i]->Delete();
255 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
257 fCluCoord[i]->Delete();
266 //____________________________________________________________________________
267 Int_t AliITStrackerSA::Clusters2Tracks(AliESDEvent *event){
268 // This method is used to find and fit the tracks. By default the corresponding
269 // method in the parent class is invoked. In this way a combined tracking
270 // TPC+ITS is performed. If the flag fITSStandAlone is true, the tracking
271 // is done in the ITS only. In the standard reconstruction chain this option
272 // can be set via AliReconstruction::SetOption("ITS","onlyITS")
275 rc=AliITStrackerMI::Clusters2Tracks(event);
278 AliDebug(1,"Stand Alone flag set: doing tracking in ITS alone\n");
280 if(!rc) rc=FindTracks(event);
284 //____________________________________________________________________________
285 void AliITStrackerSA::Init(){
286 // Reset all data members
288 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
298 Int_t nLoops=AliITSReconstructor::GetRecoParam()->GetNLoopsSA();
300 SetFixedWindowSizes();
302 Double_t phimin=AliITSReconstructor::GetRecoParam()->GetMinPhiSA();
303 Double_t phimax=AliITSReconstructor::GetRecoParam()->GetMaxPhiSA();
304 Double_t lambmin=AliITSReconstructor::GetRecoParam()->GetMinLambdaSA();
305 Double_t lambmax=AliITSReconstructor::GetRecoParam()->GetMaxLambdaSA();
306 SetCalculatedWindowSizes(nLoops,phimin,phimax,lambmin,lambmax);
308 fMinQ=AliITSReconstructor::GetRecoParam()->GetSAMinClusterCharge();
310 SetOuterStartLayer(1);
312 fListOfTracks=new TObjArray(0,0);
313 fListOfSATracks=new TObjArray(0,0);
318 //_______________________________________________________________________
319 void AliITStrackerSA::ResetForFinding(){
320 // Reset data members used in all loops during track finding
322 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
330 fListOfTracks->Delete();
331 fListOfSATracks->Delete();
336 //______________________________________________________________________
337 Int_t AliITStrackerSA::FindTracks(AliESDEvent* event){
339 // Track finder using the ESD object
341 AliDebug(2,Form(" field is %f",event->GetMagneticField()));
342 AliDebug(2,Form("SKIPPING %d %d %d %d %d %d",ForceSkippingOfLayer(0),ForceSkippingOfLayer(1),ForceSkippingOfLayer(2),ForceSkippingOfLayer(3),ForceSkippingOfLayer(4),ForceSkippingOfLayer(5)));
345 Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
348 //Reads event and mark clusters of traks already found, with flag kITSin
349 Int_t nentr=event->GetNumberOfTracks();
350 if(AliITSReconstructor::GetRecoParam()->GetSAUseAllClusters()==kFALSE) {
352 AliESDtrack *track=event->GetTrack(nentr);
353 if ((track->GetStatus()&AliESDtrack::kITSin) == AliESDtrack::kITSin){
355 Int_t ncl = track->GetITSclusters(idx);
356 for(Int_t k=0;k<ncl;k++){
357 AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
358 cll->SetBit(kSAflag);
364 Double_t primaryVertex[3];
365 event->GetVertex()->GetXYZ(primaryVertex);
366 //Creates TClonesArray with clusters for each layer. The clusters already used
367 //by AliITStrackerMI are not considered
368 Int_t nclusters[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
369 Int_t dmar[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
370 if (fCluLayer == 0) {
371 fCluLayer = new TClonesArray*[AliITSgeomTGeo::kNLayers];
372 fCluCoord = new TClonesArray*[AliITSgeomTGeo::kNLayers];
373 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
378 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
379 AliITSlayer &layer=fgLayers[i];
380 if (!ForceSkippingOfLayer(i)) {
381 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
382 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
383 if(cls->TestBit(kSAflag)==kTRUE) continue; //clusters used by TPC prol.
384 if(cls->GetQ()==0) continue; //fake clusters dead zones
385 if(i>1 && cls->GetQ()<=fMinQ) continue; // cut on SDD and SSD cluster charge
391 fCluLayer[i] = new TClonesArray("AliITSRecPoint",nclusters[i]);
393 fCluCoord[i] = new TClonesArray("AliITSclusterTable",nclusters[i]);
396 for(Int_t ilay=0;ilay<AliITSgeomTGeo::GetNLayers();ilay++){
397 TClonesArray &clulay = *fCluLayer[ilay];
398 TClonesArray &clucoo = *fCluCoord[ilay];
399 AliITSlayer &layer=fgLayers[ilay];
400 if (!ForceSkippingOfLayer(ilay)) {
401 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
402 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
403 if(cls->TestBit(kSAflag)==kTRUE) continue;
404 if(cls->GetQ()==0) continue;
405 if(ilay>1 && cls->GetQ()<=fMinQ) continue;
406 Double_t phi=0;Double_t lambda=0;
407 Float_t x=0;Float_t y=0;Float_t z=0;
408 Float_t sx=0;Float_t sy=0;Float_t sz=0;
409 GetCoorAngles(cls,phi,lambda,x,y,z,primaryVertex);
410 GetCoorErrors(cls,sx,sy,sz);
411 new (clulay[dmar[ilay]]) AliITSRecPoint(*cls);
412 new (clucoo[dmar[ilay]]) AliITSclusterTable(x,y,z,sx,sy,sz,phi,lambda,cli);
421 static Int_t nClusLay[AliITSgeomTGeo::kNLayers];//counter for clusters on each layer
423 if(AliITSReconstructor::GetRecoParam()->GetAllowProlongationWithEmptyRoad()) fMinNPoints=2;
426 // loop on minimum number of points
427 for(Int_t iMinNPoints=AliITSgeomTGeo::GetNLayers(); iMinNPoints>=fMinNPoints; iMinNPoints--) {
429 if(!fInwardFlag){ // Tracking outwards from the inner layers
430 // loop on starting layer for track finding
431 for(Int_t innLay=0; innLay<=fOuterStartLayer; innLay++) {
432 if(ForceSkippingOfLayer(innLay)) continue;
433 Int_t minNPoints=iMinNPoints-innLay;
434 for(Int_t i=innLay+1;i<AliITSgeomTGeo::GetNLayers();i++)
435 if(ForceSkippingOfLayer(i))
437 if(minNPoints<fMinNPoints) continue;
439 // loop on phi and lambda window size
440 for(Int_t nloop=0;nloop<fNloop;nloop++){
441 Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
444 Bool_t useRP=SetFirstPoint(innLay,nclInnLay,primaryVertex);
446 AliITStrackSA* trs = new AliITStrackSA();
450 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
453 nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
454 trs,primaryVertex[2],pflag);
455 for(Int_t nextLay=innLay+1; nextLay<AliITSgeomTGeo::GetNLayers(); nextLay++) {
457 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
458 trs,primaryVertex[2],pflag);
462 fPoint3[0]=fPointc[0];
463 fPoint3[1]=fPointc[1];
471 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
472 if(nClusLay[nnp]!=0) layOK+=1;
474 if(layOK>=minNPoints){
475 AliDebug(2,Form("---NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]));
476 AliITStrackV2* tr2 = 0;
477 tr2 = FitTrack(trs,primaryVertex);
479 AliDebug(2,Form("---NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
481 StoreTrack(tr2,event);
487 }//end loop on clusters of innLay
488 } //end loop on window sizes
489 } //end loop on innLay
490 }else{// Tracking inwards from the outer layers
491 // loop on starting layer for track finding
492 for(Int_t outLay=AliITSgeomTGeo::GetNLayers()-1; outLay>=fInnerStartLayer; outLay--) {
494 if(ForceSkippingOfLayer(outLay)) continue;
495 Int_t minNPoints=iMinNPoints-(AliITSgeomTGeo::GetNLayers()-1-outLay);
496 for(Int_t i=0;i<outLay-1;i++)
497 if(ForceSkippingOfLayer(i))
499 if(minNPoints<fMinNPoints) continue;
501 // loop on phi and lambda window size
502 for(Int_t nloop=0;nloop<fNloop;nloop++){
503 Int_t nclOutLay=fCluLayer[outLay]->GetEntries();
506 Bool_t useRP=SetFirstPoint(outLay,nclOutLay,primaryVertex);
508 AliITStrackSA* trs = new AliITStrackSA();
512 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
515 nClusLay[kk] = SearchClusters(outLay,fPhiWin[nloop],fLambdaWin[nloop],
516 trs,primaryVertex[2],pflag);
517 for(Int_t nextLay=outLay-1; nextLay>=0; nextLay--) {
519 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
520 trs,primaryVertex[2],pflag);
524 fPoint3[0]=fPointc[0];
525 fPoint3[1]=fPointc[1];
533 for(Int_t nnp=outLay; nnp>=0; nnp--){
534 if(nClusLay[nnp]!=0) layOK+=1;
536 if(layOK>=minNPoints){
537 AliDebug(2,Form("---NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]));
538 AliITStrackV2* tr2 = 0;
539 tr2 = FitTrack(trs,primaryVertex);
541 AliDebug(2,Form("---NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
543 StoreTrack(tr2,event);
549 }//end loop on clusters of innLay
550 } //end loop on window sizes
551 } //end loop on innLay
552 } // end if (fInwardFlag)
553 }//end loop on min points
555 // search for 1-point tracks in SPD, only for cosmics
556 // (A.Dainese 21.03.08)
557 if(AliITSReconstructor::GetRecoParam()->GetSAOnePointTracks() &&
558 TMath::Abs(event->GetMagneticField())<0.01) {
559 Int_t outerLayer=1; // only SPD
560 for(Int_t innLay=0; innLay<=TMath::Min(1,fOuterStartLayer); innLay++) {
561 // counter for clusters on each layer
563 for(Int_t nloop=0;nloop<fNloop;nloop++){
564 Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
565 while(nclInnLay--){ //loop starting from layer innLay
567 Bool_t useRP=SetFirstPoint(innLay,nclInnLay,primaryVertex);
569 AliITStrackSA* trs = new AliITStrackSA();
573 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
576 nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
577 trs,primaryVertex[2],pflag);
578 for(Int_t nextLay=innLay+1; nextLay<=outerLayer; nextLay++) {
580 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
581 trs,primaryVertex[2],pflag);
585 fPoint3[0]=fPointc[0];
586 fPoint3[1]=fPointc[1];
594 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
595 if(nClusLay[nnp]!=0) layOK+=1;
598 AliDebug(2,Form("----NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]));
599 AliITStrackV2* tr2 = 0;
600 Bool_t onePoint = kTRUE;
601 tr2 = FitTrack(trs,primaryVertex,onePoint);
603 AliDebug(2,Form("----NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
605 StoreTrack(tr2,event);
611 }//end loop on clusters of innLay
612 } //end loop on window sizes
614 } //end loop on innLay
615 } // end search 1-point tracks
617 Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
622 //________________________________________________________________________
624 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex,Bool_t onePoint) {
625 //fit of the found track (most general case, also <6 points, layers missing)
626 // A.Dainese 16.11.07
629 const Int_t kMaxClu=AliITStrackSA::kMaxNumberOfClusters;
631 static Int_t firstmod[AliITSgeomTGeo::kNLayers];
633 static Int_t clind[AliITSgeomTGeo::kNLayers][kMaxClu];
634 static Int_t clmark[AliITSgeomTGeo::kNLayers][kMaxClu];
635 static Int_t nnn[AliITSgeomTGeo::kNLayers];
636 static Int_t kkk[AliITSgeomTGeo::kNLayers];
637 static Int_t end[AliITSgeomTGeo::kNLayers];
638 static AliITSRecPoint *listlayer[AliITSgeomTGeo::kNLayers][kMaxClu];
640 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
641 firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
645 for(Int_t j=0;j<kMaxClu; j++){
653 Int_t nclusters = tr->GetNumberOfClustersSA();
654 for(Int_t ncl=0;ncl<nclusters;ncl++){
655 Int_t index = tr->GetClusterIndexSA(ncl);
656 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
657 if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
658 Int_t lay = (index & 0xf0000000) >> 28;
659 Int_t nInLay=end[lay];
660 listlayer[lay][nInLay]=cl;
662 clind[lay][ind]=index;
666 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
667 for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
668 Int_t mark = tr->GetClusterMark(nlay,ncl);
670 clmark[nlay][ind]=mark;
675 Int_t firstLay=-1,secondLay=-1;
676 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
682 } else if(secondLay==-1) {
688 if(firstLay==-1 || (secondLay==-1 && !onePoint)) return 0;
690 for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
691 AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0][l0];
692 for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
693 AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1][l1];
694 for(Int_t l2=0;l2<end[2];l2++){ //loop on layer 3
695 AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2][l2];
696 for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4
697 AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3][l3];
698 for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
699 AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4][l4];
700 for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6
701 AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5][l5];
704 Double_t x1,y1,z1,sx1,sy1,sz1;
705 Double_t x2,y2,z2,sx2,sy2,sz2;
706 AliITSRecPoint* p1=0;
707 AliITSRecPoint* p2=0;
708 Int_t index1=0,index2=0;
714 index1=clind[0][l0];mrk1=clmark[0][l0];
718 index1=clind[1][l1];mrk1=clmark[1][l1];
722 index1=clind[2][l2];mrk1=clmark[2][l2];
726 index1=clind[3][l3];mrk1=clmark[3][l3];
730 index1=clind[4][l4];mrk1=clmark[4][l4];
737 index2=clind[1][l1];mrk2=clmark[1][l1];
741 index2=clind[2][l2];mrk2=clmark[2][l2];
745 index2=clind[3][l3];mrk2=clmark[3][l3];
749 index2=clind[4][l4];mrk2=clmark[4][l4];
753 index2=clind[5][l5];mrk2=clmark[5][l5];
761 Int_t module1 = p1->GetDetectorIndex()+firstmod[firstLay];
762 Int_t layer,ladder,detector;
763 AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
764 Float_t yclu1 = p1->GetY();
765 Float_t zclu1 = p1->GetZ();
766 Double_t cv=0,tgl2=0,phi2=0;
769 AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(firstLay,cln1);
779 AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(secondLay,cln2);
786 cv = Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
787 tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
788 phi2 = TMath::ATan2((y2-y1),(x2-x1));
789 } else { // special case of 1-point tracks, only for cosmics (B=0)
790 x2 = primaryVertex[0];
791 y2 = primaryVertex[1];
792 z2 = primaryVertex[2];
794 tgl2 = (z1-z2)/TMath::Sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
795 phi2 = TMath::ATan2((y1-y2),(x1-x2));
799 AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
803 trac->AddClusterV2(5,(clind[5][l5] & 0x0fffffff)>>0);
804 trac->AddClusterMark(5,clmark[5][l5]);
807 trac->AddClusterV2(4,(clind[4][l4] & 0x0fffffff)>>0);
808 trac->AddClusterMark(4,clmark[4][l4]);
811 trac->AddClusterV2(3,(clind[3][l3] & 0x0fffffff)>>0);
812 trac->AddClusterMark(3,clmark[3][l3]);
815 trac->AddClusterV2(2,(clind[2][l2] & 0x0fffffff)>>0);
816 trac->AddClusterMark(2,clmark[2][l2]);
819 trac->AddClusterV2(1,(clind[1][l1] & 0x0fffffff)>>0);
820 trac->AddClusterMark(1,clmark[1][l1]);
823 trac->AddClusterV2(0,(clind[0][l0] & 0x0fffffff)>>0);
824 trac->AddClusterMark(0,clmark[0][l0]);
827 //fit with Kalman filter using AliITStrackerMI::RefitAt()
828 AliITStrackMI* ot = new AliITStrackSA(*trac);
830 ot->ResetCovariance(10.);
833 // Propagate inside the innermost layer with a cluster
834 if(ot->Propagate(ot->GetX()-0.1*ot->GetX())) {
836 if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
837 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
838 otrack2->ResetCovariance(10.);
839 otrack2->ResetClusters();
840 //fit from layer 6 to layer 1
841 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
842 fListOfTracks->AddLast(otrack2);
843 fListOfSATracks->AddLast(trac);
862 if(fListOfTracks->GetEntries()==0) return 0;
864 Int_t lowchi2 = FindTrackLowChiSquare();
865 AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
866 AliITStrackSA* trsa = (AliITStrackSA*)fListOfSATracks->At(lowchi2);
868 if(otrack==0) return 0;
870 Int_t indexc[AliITSgeomTGeo::kNLayers];
871 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
872 for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
873 indexc[nind] = otrack->GetClusterIndex(nind);
876 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
877 if(i<otrack->GetNumberOfClusters()) {
878 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc[i]);
879 labl[i][0]=cl->GetLabel(0);
880 labl[i][1]=cl->GetLabel(1);
881 labl[i][2]=cl->GetLabel(2);
889 CookLabel(otrack,0.); //MI change - to see fake ratio
891 Int_t label=FindLabel(labl[0][0],labl[1][0],labl[2][0],labl[3][0],labl[4][0],labl[5][0]);
893 for(Int_t i=0;i<otrack->GetNumberOfClusters();i++)
894 if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
896 if(lflag<otrack->GetNumberOfClusters()) label = -label;
897 otrack->SetLabel(label);
899 //remove clusters of found track
900 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
901 for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
902 Int_t index = trsa->GetClusterMark(nlay,cln);
903 fCluLayer[nlay]->RemoveAt(index);
904 RemoveClusterCoord(nlay,index);
905 fCluLayer[nlay]->Compress();
913 //_______________________________________________________
914 void AliITStrackerSA::StoreTrack(AliITStrackV2 *t,AliESDEvent *event) const
917 // Add new track to the ESD
919 AliESDtrack outtrack;
920 outtrack.UpdateTrackParams(t,AliESDtrack::kITSin);
921 for(Int_t i=0;i<12;i++) {
922 outtrack.SetITSModuleIndex(i,t->GetModuleIndex(i));
924 Double_t sdedx[4]={0.,0.,0.,0.};
925 for(Int_t i=0; i<4; i++) sdedx[i]=t->GetSampledEdx(i);
926 outtrack.SetITSdEdxSamples(sdedx);
927 event->AddTrack(&outtrack);
933 //_______________________________________________________
934 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t /*zvertex*/,Int_t pflag){
935 //function used to to find the clusters associated to the track
937 if(ForceSkippingOfLayer(layer)) return 0;
940 AliITSlayer &lay = fgLayers[layer];
941 Double_t r=lay.GetR();
943 Float_t cx1,cx2,cy1,cy2;
944 FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
945 if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
947 Double_t fi1=TMath::ATan2(cy1-fPoint1[1],cx1-fPoint1[0]);
948 Double_t fi2=TMath::ATan2(cy2-fPoint1[1],cx2-fPoint1[0]);
949 fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
953 Int_t ncl = fCluLayer[layer]->GetEntries();
954 for (Int_t index=0; index<ncl; index++) {
955 AliITSRecPoint *c = (AliITSRecPoint*)fCluLayer[layer]->At(index);
957 if (c->GetQ()<=0) continue;
958 if(layer>1 && c->GetQ()<=fMinQ) continue;
960 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(layer,index);
961 Double_t phi = arr->GetPhi();
962 if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
964 Double_t lambda = arr->GetLambda();
965 if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
967 if(trs->GetNumberOfClustersSA()==trs->GetMaxNumberOfClusters()) return 0;
968 if(trs->GetNumberOfMarked(layer)==trs->GetMaxNMarkedPerLayer()) return 0;
969 Int_t orind = arr->GetOrInd();
970 trs->AddClusterSA(layer,orind);
971 trs->AddClusterMark(layer,index);
977 fPointc[0]=arr->GetX();
978 fPointc[1]=arr->GetY();
984 //________________________________________________________________
985 Bool_t AliITStrackerSA::SetFirstPoint(Int_t lay, Int_t clu, Double_t* primaryVertex){
986 // Sets the first point (seed) for tracking
988 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[lay]->At(clu);
989 if(!cl) return kFALSE;
990 if (cl->GetQ()<=0) return kFALSE;
991 if(lay>1 && cl->GetQ()<=fMinQ) return kFALSE;
993 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(lay,clu);
994 fPhic = arr->GetPhi();
995 fLambdac = arr->GetLambda();
996 fPhiEstimate = fPhic;
997 fPoint1[0]=primaryVertex[0];
998 fPoint1[1]=primaryVertex[1];
999 fPoint2[0]=arr->GetX();
1000 fPoint2[1]=arr->GetY();
1004 //________________________________________________________________
1005 void AliITStrackerSA::UpdatePoints(){
1006 //update of points for the estimation of the curvature
1008 fPoint2[0]=fPoint3[0];
1009 fPoint2[1]=fPoint3[1];
1010 fPoint3[0]=fPointc[0];
1011 fPoint3[1]=fPointc[1];
1016 //___________________________________________________________________
1017 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){
1019 //given (x,y) of three recpoints (in global coordinates)
1020 //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
1022 Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1023 if(den==0) return 0;
1024 a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1025 b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1026 c = -x1*x1-y1*y1-a*x1-b*y1;
1029 //__________________________________________________________________________
1030 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){
1032 //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
1033 //c2 is -rlayer*rlayer
1037 Double_t aA = (b1*b1)/(a1*a1)+1;
1038 Double_t bB = (-2*m*b1/(a1*a1));
1039 Double_t cC = c2+(m*m)/(a1*a1);
1040 Double_t dD = bB*bB-4*aA*cC;
1043 y1 = (-bB+TMath::Sqrt(dD))/(2*aA);
1044 y2 = (-bB-TMath::Sqrt(dD))/(2*aA);
1045 x1 = (c2-c1-b1*y1)/a1;
1046 x2 = (c2-c1-b1*y2)/a1;
1050 //____________________________________________________________________
1051 Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t
1052 x2,Double_t y2,Double_t x3,Double_t y3){
1054 //calculates the curvature of track
1055 Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1056 if(den==0) return 0;
1057 Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1058 Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1059 Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1062 if((a*a+b*b-4*c)<0) return 0;
1063 Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1064 if(rad==0) return 0;
1066 if((x1>0 && y1>0 && x1<xc)) rad*=-1;
1067 if((x1<0 && y1>0 && x1<xc)) rad*=-1;
1068 // if((x1<0 && y1<0 && x1<xc)) rad*=-1;
1069 // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1076 //____________________________________________________________________
1077 Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1079 //Returns the point closest to pp
1081 Double_t diff1 = p1-pp;
1082 Double_t diff2 = p2-pp;
1084 if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1085 else fPhiEstimate=p2;
1086 return fPhiEstimate;
1091 //_________________________________________________________________
1092 Int_t AliITStrackerSA::FindTrackLowChiSquare() const {
1093 // returns track with lowest chi square
1094 Int_t dim=fListOfTracks->GetEntries();
1095 if(dim<=1) return 0;
1096 AliITStrackV2* trk = (AliITStrackV2*)fListOfTracks->At(0);
1097 Double_t minChi2=trk->GetChi2();
1099 for(Int_t i=1;i<dim;i++){
1100 trk = (AliITStrackV2*)fListOfTracks->At(i);
1101 Double_t chi2=trk->GetChi2();
1110 //__________________________________________________________
1111 Int_t AliITStrackerSA::FindLabel(Int_t l0, Int_t l1, Int_t l2, Int_t l3, Int_t l4, Int_t l5){
1113 //function used to determine the track label
1115 Int_t lb[6] = {l0,l1,l2,l3,l4,l5};
1116 Int_t aa[6]={1,1,1,1,1,1};
1119 Int_t k=0;Int_t w=0;Int_t num=6;
1120 for(Int_t i=5;i>=0;i--) if(lb[i]==-1) num=i;
1124 for(Int_t i=k+1;i<num;i++){
1126 if(lb[k]==lb[i] && aa[k]!=0){
1137 for(Int_t j=0;j<6;j++){
1150 if(num<1) return -1;
1154 //_____________________________________________________________________________
1155 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,
1156 Int_t gl12, Int_t gl13, Int_t gl14,Int_t gl15, Int_t gl16, Int_t gl17, Int_t gl18, Int_t minNPoints){
1159 //function used to assign label to the found track. If track is fake, the label is negative
1161 Int_t lb0[6] = {gl1,gl2,gl3,gl4,gl5,gl6};
1162 Int_t lb1[6] = {gl7,gl8,gl9,gl10,gl11,gl12};
1163 Int_t lb2[6] = {gl13,gl14,gl15,gl16,gl17,gl18};
1164 Int_t ll=FindLabel(lb0[0],lb0[1],lb0[2],lb0[3],lb0[4],lb0[5]);
1165 Int_t lflag=0;Int_t num=6;
1166 if(lb0[5]==-1 && lb1[5]==-1 && lb2[5]==-1) num=5;
1168 for(Int_t i=0;i<num;i++){
1169 if(lb0[i]==ll || lb1[i]==ll || lb2[i]==ll) lflag+=1;
1172 if(lflag>=minNPoints) return ll;
1177 //_____________________________________________________________________________
1178 void AliITStrackerSA::SetCalculatedWindowSizes(Int_t n, Float_t phimin, Float_t phimax, Float_t lambdamin, Float_t lambdamax){
1179 // Set sizes of the phi and lambda windows used for track finding
1181 if(fPhiWin) delete [] fPhiWin;
1182 if(fLambdaWin) delete [] fLambdaWin;
1183 fPhiWin = new Double_t[fNloop];
1184 fLambdaWin = new Double_t[fNloop];
1185 Float_t stepPhi=(phimax-phimin)/(Float_t)(fNloop-1);
1186 Float_t stepLambda=(lambdamax-lambdamin)/(Float_t)(fNloop-1);
1187 for(Int_t k=0;k<fNloop;k++){
1188 Float_t phi=phimin+k*stepPhi;
1189 Float_t lam=lambdamin+k*stepLambda;
1194 //_____________________________________________________________________________
1195 void AliITStrackerSA::SetFixedWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
1196 // Set sizes of the phi and lambda windows used for track finding
1198 if(phi){ // user defined values
1199 fPhiWin = new Double_t[fNloop];
1200 fLambdaWin = new Double_t[fNloop];
1201 for(Int_t k=0;k<fNloop;k++){
1203 fLambdaWin[k]=lam[k];
1206 else { // default values
1208 Double_t phid[33] = {0.002,0.003,0.004,0.0045,0.0047,
1209 0.005,0.0053,0.0055,
1210 0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1211 0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
1212 0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1213 Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1215 0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1216 0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1217 0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
1224 fPhiWin = new Double_t[fNloop];
1225 fLambdaWin = new Double_t[fNloop];
1227 Double_t factor=AliITSReconstructor::GetRecoParam()->GetFactorSAWindowSizes(); // possibility to enlarge windows for cosmics reco with large misalignments (A.Dainese)
1229 for(Int_t k=0;k<fNloop;k++){
1230 fPhiWin[k]=phid[k]*factor;
1231 fLambdaWin[k]=lambdad[k]*factor;
1237 //_______________________________________________________________________
1238 void AliITStrackerSA::GetCoorAngles(AliITSRecPoint* cl,Double_t &phi,Double_t &lambda, Float_t &x, Float_t &y,Float_t &z,Double_t* vertex){
1239 //Returns values of phi (azimuthal) and lambda angles for a given cluster
1241 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
1242 Int_t lay,lad,det; fGeom->GetModuleId(module,lay,lad,det);
1243 Float_t tx,ty,tz; fGeom->GetTrans(lay,lad,det,tx,ty,tz);
1245 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1246 Double_t phi1=TMath::Pi()/2+alpha;
1247 if (lay==1) phi1+=TMath::Pi();
1249 Float_t cp=TMath::Cos(phi1), sp=TMath::Sin(phi1);
1250 Float_t r=tx*cp+ty*sp;
1252 xyz= r*cp - cl->GetY()*sp;
1253 y= r*sp + cl->GetY()*cp;
1257 cl->GetGlobalXYZ(xyz);
1262 phi=TMath::ATan2(y-vertex[1],x-vertex[0]);
1263 lambda=TMath::ATan2(z-vertex[2],TMath::Sqrt((x-vertex[0])*(x-vertex[0])+(y-vertex[1])*(y-vertex[1])));
1266 //________________________________________________________________________
1267 void AliITStrackerSA::GetCoorErrors(AliITSRecPoint* cl,Float_t &sx,Float_t &sy, Float_t &sz){
1269 //returns sigmax, y, z of cluster in global coordinates
1271 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
1273 AliITSgeomTGeo::GetModuleId(module,lay,lad,det);
1275 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1276 Double_t phi=TMath::Pi()/2+alpha;
1277 if (lay==1) phi+=TMath::Pi();
1279 Float_t cp=TMath::Cos(phi), sp=TMath::Sin(phi);
1282 cl->GetGlobalCov(covm);
1283 sx=TMath::Sqrt(covm[0]);
1284 sy=TMath::Sqrt(covm[3]);
1285 sz=TMath::Sqrt(covm[5]);
1287 sx = TMath::Sqrt(sp*sp*cl->GetSigmaY2());
1288 sy = TMath::Sqrt(cp*cp*cl->GetSigmaY2());
1289 sz = TMath::Sqrt(cl->GetSigmaZ2());