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");
281 rc=FindTracks(event,kFALSE);
282 if(AliITSReconstructor::GetRecoParam()->GetSAUseAllClusters()==kTRUE) {
283 rc=FindTracks(event,kTRUE);
290 //____________________________________________________________________________
291 void AliITStrackerSA::Init(){
292 // Reset all data members
294 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
304 Int_t nLoops=AliITSReconstructor::GetRecoParam()->GetNLoopsSA();
306 SetFixedWindowSizes();
308 Double_t phimin=AliITSReconstructor::GetRecoParam()->GetMinPhiSA();
309 Double_t phimax=AliITSReconstructor::GetRecoParam()->GetMaxPhiSA();
310 Double_t lambmin=AliITSReconstructor::GetRecoParam()->GetMinLambdaSA();
311 Double_t lambmax=AliITSReconstructor::GetRecoParam()->GetMaxLambdaSA();
312 SetCalculatedWindowSizes(nLoops,phimin,phimax,lambmin,lambmax);
314 fMinQ=AliITSReconstructor::GetRecoParam()->GetSAMinClusterCharge();
316 SetOuterStartLayer(1);
318 fListOfTracks=new TClonesArray("AliITStrackMI",100);
319 fListOfSATracks=new TClonesArray("AliITStrackSA",100);
324 //_______________________________________________________________________
325 void AliITStrackerSA::ResetForFinding(){
326 // Reset data members used in all loops during track finding
328 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
336 fListOfTracks->Clear();
337 fListOfSATracks->Clear();
342 //______________________________________________________________________
343 Int_t AliITStrackerSA::FindTracks(AliESDEvent* event, Bool_t useAllClusters){
345 // Track finder using the ESD object
347 AliDebug(2,Form(" field is %f",event->GetMagneticField()));
348 AliDebug(2,Form("SKIPPING %d %d %d %d %d %d",ForceSkippingOfLayer(0),ForceSkippingOfLayer(1),ForceSkippingOfLayer(2),ForceSkippingOfLayer(3),ForceSkippingOfLayer(4),ForceSkippingOfLayer(5)));
351 Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
354 //Reads event and mark clusters of traks already found, with flag kITSin
355 Int_t nentr=event->GetNumberOfTracks();
356 if(!useAllClusters) {
358 AliESDtrack *track=event->GetTrack(nentr);
359 if ((track->GetStatus()&AliESDtrack::kITSin) == AliESDtrack::kITSin){
361 Int_t ncl = track->GetITSclusters(idx);
362 for(Int_t k=0;k<ncl;k++){
363 AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
364 cll->SetBit(kSAflag);
370 AliESDtrack *track=event->GetTrack(nentr);
371 if ((track->GetStatus()&AliESDtrack::kITSin) == AliESDtrack::kITSin){
373 Int_t ncl = track->GetITSclusters(idx);
374 for(Int_t k=0;k<ncl;k++){
375 AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
376 cll->ResetBit(kSAflag);
382 Double_t primaryVertex[3];
383 event->GetVertex()->GetXYZ(primaryVertex);
384 //Creates TClonesArray with clusters for each layer. The clusters already used
385 //by AliITStrackerMI are not considered
386 Int_t nclusters[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
387 Int_t dmar[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
388 if (fCluLayer == 0) {
389 fCluLayer = new TClonesArray*[AliITSgeomTGeo::kNLayers];
390 fCluCoord = new TClonesArray*[AliITSgeomTGeo::kNLayers];
391 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
396 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
397 AliITSlayer &layer=fgLayers[i];
398 if (!ForceSkippingOfLayer(i)) {
399 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
400 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
401 if(cls->TestBit(kSAflag)==kTRUE) continue; //clusters used by TPC prol.
402 if(cls->GetQ()==0) continue; //fake clusters dead zones
403 if(i>1 && cls->GetQ()<=fMinQ) continue; // cut on SDD and SSD cluster charge
409 fCluLayer[i] = new TClonesArray("AliITSRecPoint",nclusters[i]);
411 fCluLayer[i]->Delete();
412 fCluLayer[i]->Expand(nclusters[i]);
415 fCluCoord[i] = new TClonesArray("AliITSclusterTable",nclusters[i]);
417 fCluCoord[i]->Delete();
418 fCluCoord[i]->Expand(nclusters[i]);
422 for(Int_t ilay=0;ilay<AliITSgeomTGeo::GetNLayers();ilay++){
423 TClonesArray &clulay = *fCluLayer[ilay];
424 TClonesArray &clucoo = *fCluCoord[ilay];
425 AliITSlayer &layer=fgLayers[ilay];
426 if (!ForceSkippingOfLayer(ilay)) {
427 for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
428 AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
429 if(cls->TestBit(kSAflag)==kTRUE) continue;
430 if(cls->GetQ()==0) continue;
431 if(ilay>1 && cls->GetQ()<=fMinQ) continue;
432 Double_t phi=0;Double_t lambda=0;
433 Float_t x=0;Float_t y=0;Float_t z=0;
434 Float_t sx=0;Float_t sy=0;Float_t sz=0;
435 GetCoorAngles(cls,phi,lambda,x,y,z,primaryVertex);
436 GetCoorErrors(cls,sx,sy,sz);
437 new (clulay[dmar[ilay]]) AliITSRecPoint(*cls);
438 new (clucoo[dmar[ilay]]) AliITSclusterTable(x,y,z,sx,sy,sz,phi,lambda,cli);
447 static Int_t nClusLay[AliITSgeomTGeo::kNLayers];//counter for clusters on each layer
448 Int_t startLayForSeed=0;
449 Int_t lastLayForSeed=fOuterStartLayer;
450 Int_t nSeedSteps=lastLayForSeed-startLayForSeed;
453 startLayForSeed=AliITSgeomTGeo::GetNLayers()-1;
454 lastLayForSeed=fInnerStartLayer;
455 nSeedSteps=startLayForSeed-lastLayForSeed;
459 // loop on minimum number of points
460 for(Int_t iMinNPoints=AliITSgeomTGeo::GetNLayers(); iMinNPoints>=fMinNPoints; iMinNPoints--) {
462 // loop on starting layer for track finding
463 for(Int_t iSeedLay=0; iSeedLay<=nSeedSteps; iSeedLay++) {
464 Int_t theLay=startLayForSeed+iSeedLay*seedStep;
465 if(ForceSkippingOfLayer(theLay)) continue;
466 Int_t minNPoints=iMinNPoints-theLay;
467 if(fInwardFlag) minNPoints=iMinNPoints-(AliITSgeomTGeo::GetNLayers()-1-theLay);
468 for(Int_t i=theLay+1;i<AliITSgeomTGeo::GetNLayers();i++)
469 if(ForceSkippingOfLayer(i))
471 if(minNPoints<fMinNPoints) continue;
473 // loop on phi and lambda window size
474 for(Int_t nloop=0;nloop<fNloop;nloop++){
475 Int_t nclTheLay=fCluLayer[theLay]->GetEntries();
478 Bool_t useRP=SetFirstPoint(theLay,nclTheLay,primaryVertex);
484 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
487 nClusLay[kk] = SearchClusters(theLay,fPhiWin[nloop],fLambdaWin[nloop],
488 &trs,primaryVertex[2],pflag);
489 Int_t nextLay=theLay+seedStep;
493 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
494 &trs,primaryVertex[2],pflag);
498 fPoint3[0]=fPointc[0];
499 fPoint3[1]=fPointc[1];
505 if(nextLay<0 || nextLay==6) goon=kFALSE;
511 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-theLay;nnp++){
512 if(nClusLay[nnp]!=0) layOK+=1;
515 for(Int_t nnp=theLay; nnp>=0; nnp--){
516 if(nClusLay[nnp]!=0) layOK+=1;
519 if(layOK>=minNPoints){
520 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]));
521 AliITStrackV2* tr2 = 0;
522 tr2 = FitTrack(&trs,primaryVertex);
526 AliDebug(2,Form("---NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
528 StoreTrack(tr2,event,useAllClusters);
533 }//end loop on clusters of theLay
534 } //end loop on window sizes
535 } //end loop on theLay
536 }//end loop on min points
538 // search for 1-point tracks in SPD, only for cosmics
539 // (A.Dainese 21.03.08)
540 if(AliITSReconstructor::GetRecoParam()->GetSAOnePointTracks() &&
541 TMath::Abs(event->GetMagneticField())<0.01) {
542 Int_t outerLayer=1; // only SPD
543 for(Int_t innLay=0; innLay<=TMath::Min(1,fOuterStartLayer); innLay++) {
544 // counter for clusters on each layer
546 for(Int_t nloop=0;nloop<fNloop;nloop++){
547 Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
548 while(nclInnLay--){ //loop starting from layer innLay
550 Bool_t useRP=SetFirstPoint(innLay,nclInnLay,primaryVertex);
556 for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
559 nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
560 &trs,primaryVertex[2],pflag);
561 for(Int_t nextLay=innLay+1; nextLay<=outerLayer; nextLay++) {
563 nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
564 &trs,primaryVertex[2],pflag);
568 fPoint3[0]=fPointc[0];
569 fPoint3[1]=fPointc[1];
577 for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
578 if(nClusLay[nnp]!=0) layOK+=1;
581 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]));
582 AliITStrackV2* tr2 = 0;
583 Bool_t onePoint = kTRUE;
584 tr2 = FitTrack(&trs,primaryVertex,onePoint);
588 AliDebug(2,Form("----NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
590 StoreTrack(tr2,event,useAllClusters);
595 }//end loop on clusters of innLay
596 } //end loop on window sizes
598 } //end loop on innLay
599 } // end search 1-point tracks
601 if(!useAllClusters) AliInfo(Form("Number of found tracks: %d",event->GetNumberOfTracks()));
607 //________________________________________________________________________
609 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex,Bool_t onePoint) {
610 //fit of the found track (most general case, also <6 points, layers missing)
611 // A.Dainese 16.11.07
614 const Int_t kMaxClu=AliITStrackSA::kMaxNumberOfClusters;
616 static Int_t firstmod[AliITSgeomTGeo::kNLayers];
618 static Int_t clind[AliITSgeomTGeo::kNLayers][kMaxClu];
619 static Int_t clmark[AliITSgeomTGeo::kNLayers][kMaxClu];
620 static Int_t end[AliITSgeomTGeo::kNLayers];
621 static AliITSRecPoint *listlayer[AliITSgeomTGeo::kNLayers][kMaxClu];
623 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
624 firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
626 for(Int_t j=0;j<kMaxClu; j++){
634 Int_t nclusters = tr->GetNumberOfClustersSA();
635 for(Int_t ncl=0;ncl<nclusters;ncl++){
636 Int_t index = tr->GetClusterIndexSA(ncl);
637 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
638 Int_t lay = (index & 0xf0000000) >> 28;
639 Int_t nInLay=end[lay];
640 listlayer[lay][nInLay]=cl;
641 clind[lay][nInLay]=index;
645 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
646 for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
647 Int_t mark = tr->GetClusterMark(nlay,ncl);
648 clmark[nlay][ncl]=mark;
653 Int_t firstLay=-1,secondLay=-1;
654 for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
660 } else if(secondLay==-1) {
666 if(firstLay==-1 || (secondLay==-1 && !onePoint)) return 0;
667 TClonesArray &arrMI= *fListOfTracks;
668 TClonesArray &arrSA= *fListOfSATracks;
669 Int_t nFoundTracks=0;
671 for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
672 AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0][l0];
673 for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
674 AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1][l1];
675 for(Int_t l2=0;l2<end[2];l2++){ //loop on layer 3
676 AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2][l2];
677 for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4
678 AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3][l3];
679 for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
680 AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4][l4];
681 for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6
682 AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5][l5];
685 Double_t x1,y1,z1,sx1,sy1,sz1;
686 Double_t x2,y2,z2,sx2,sy2,sz2;
687 AliITSRecPoint* p1=0;
688 AliITSRecPoint* p2=0;
689 Int_t index1=0,index2=0;
695 index1=clind[0][l0];mrk1=clmark[0][l0];
699 index1=clind[1][l1];mrk1=clmark[1][l1];
703 index1=clind[2][l2];mrk1=clmark[2][l2];
707 index1=clind[3][l3];mrk1=clmark[3][l3];
711 index1=clind[4][l4];mrk1=clmark[4][l4];
718 index2=clind[1][l1];mrk2=clmark[1][l1];
722 index2=clind[2][l2];mrk2=clmark[2][l2];
726 index2=clind[3][l3];mrk2=clmark[3][l3];
730 index2=clind[4][l4];mrk2=clmark[4][l4];
734 index2=clind[5][l5];mrk2=clmark[5][l5];
742 Int_t module1 = p1->GetDetectorIndex()+firstmod[firstLay];
743 Int_t layer,ladder,detector;
744 AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
745 Float_t yclu1 = p1->GetY();
746 Float_t zclu1 = p1->GetZ();
747 Double_t cv=0,tgl2=0,phi2=0;
750 AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(firstLay,cln1);
760 AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(secondLay,cln2);
767 cv = Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
768 tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
769 phi2 = TMath::ATan2((y2-y1),(x2-x1));
770 } else { // special case of 1-point tracks, only for cosmics (B=0)
771 x2 = primaryVertex[0];
772 y2 = primaryVertex[1];
773 z2 = primaryVertex[2];
775 tgl2 = (z1-z2)/TMath::Sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
776 phi2 = TMath::ATan2((y1-y2),(x1-x2));
780 AliITStrackSA trac(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
784 trac.AddClusterV2(5,(clind[5][l5] & 0x0fffffff)>>0);
785 trac.AddClusterMark(5,clmark[5][l5]);
788 trac.AddClusterV2(4,(clind[4][l4] & 0x0fffffff)>>0);
789 trac.AddClusterMark(4,clmark[4][l4]);
792 trac.AddClusterV2(3,(clind[3][l3] & 0x0fffffff)>>0);
793 trac.AddClusterMark(3,clmark[3][l3]);
796 trac.AddClusterV2(2,(clind[2][l2] & 0x0fffffff)>>0);
797 trac.AddClusterMark(2,clmark[2][l2]);
800 trac.AddClusterV2(1,(clind[1][l1] & 0x0fffffff)>>0);
801 trac.AddClusterMark(1,clmark[1][l1]);
804 trac.AddClusterV2(0,(clind[0][l0] & 0x0fffffff)>>0);
805 trac.AddClusterMark(0,clmark[0][l0]);
808 //fit with Kalman filter using AliITStrackerMI::RefitAt()
809 AliITStrackSA ot(trac);
811 ot.ResetCovariance(10.);
814 // Propagate inside the innermost layer with a cluster
815 if(ot.Propagate(ot.GetX()-0.1*ot.GetX())) {
817 if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),&ot,&trac)){ //fit from layer 1 to layer 6
818 AliITStrackMI otrack2(ot);
819 otrack2.ResetCovariance(10.);
820 otrack2.ResetClusters();
821 //fit from layer 6 to layer 1
822 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),&otrack2,&ot)) {
823 new(arrMI[nFoundTracks]) AliITStrackMI(otrack2);
824 new(arrSA[nFoundTracks]) AliITStrackSA(trac);
840 if(fListOfTracks->GetEntries()==0) return 0;
842 Int_t lowchi2 = FindTrackLowChiSquare();
843 AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
844 AliITStrackSA* trsa = (AliITStrackSA*)fListOfSATracks->At(lowchi2);
846 if(otrack==0) return 0;
848 CookLabel(otrack,0.); //MI change - to see fake ratio
849 Int_t label=FindLabel(otrack);
850 otrack->SetLabel(label);
853 otrack->CookdEdx(low,up);
855 //remove clusters of found track
856 for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
857 for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
858 Int_t index = trsa->GetClusterMark(nlay,cln);
859 fCluLayer[nlay]->RemoveAt(index);
860 RemoveClusterCoord(nlay,index);
861 fCluLayer[nlay]->Compress();
869 //_______________________________________________________
870 void AliITStrackerSA::StoreTrack(AliITStrackV2 *t,AliESDEvent *event, Bool_t pureSA) const
873 // Add new track to the ESD
875 AliESDtrack outtrack;
876 outtrack.UpdateTrackParams(t,AliESDtrack::kITSin);
877 if(pureSA) outtrack.SetStatus(AliESDtrack::kITSpureSA);
878 for(Int_t i=0;i<12;i++) {
879 outtrack.SetITSModuleIndex(i,t->GetModuleIndex(i));
881 Double_t sdedx[4]={0.,0.,0.,0.};
882 for(Int_t i=0; i<4; i++) sdedx[i]=t->GetSampledEdx(i);
883 outtrack.SetITSdEdxSamples(sdedx);
884 event->AddTrack(&outtrack);
890 //_______________________________________________________
891 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t /*zvertex*/,Int_t pflag){
892 //function used to to find the clusters associated to the track
894 if(ForceSkippingOfLayer(layer)) return 0;
897 AliITSlayer &lay = fgLayers[layer];
898 Double_t r=lay.GetR();
900 Float_t cx1,cx2,cy1,cy2;
901 FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
902 if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
904 Double_t fi1=TMath::ATan2(cy1-fPoint1[1],cx1-fPoint1[0]);
905 Double_t fi2=TMath::ATan2(cy2-fPoint1[1],cx2-fPoint1[0]);
906 fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
910 Int_t ncl = fCluLayer[layer]->GetEntries();
911 for (Int_t index=0; index<ncl; index++) {
912 AliITSRecPoint *c = (AliITSRecPoint*)fCluLayer[layer]->At(index);
914 if (c->GetQ()<=0) continue;
915 if(layer>1 && c->GetQ()<=fMinQ) continue;
917 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(layer,index);
918 Double_t phi = arr->GetPhi();
919 if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
921 Double_t lambda = arr->GetLambda();
922 if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
924 if(trs->GetNumberOfClustersSA()==trs->GetMaxNumberOfClusters()) return 0;
925 if(trs->GetNumberOfMarked(layer)==trs->GetMaxNMarkedPerLayer()) return 0;
926 Int_t orind = arr->GetOrInd();
927 trs->AddClusterSA(layer,orind);
928 trs->AddClusterMark(layer,index);
934 fPointc[0]=arr->GetX();
935 fPointc[1]=arr->GetY();
941 //________________________________________________________________
942 Bool_t AliITStrackerSA::SetFirstPoint(Int_t lay, Int_t clu, Double_t* primaryVertex){
943 // Sets the first point (seed) for tracking
945 AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[lay]->At(clu);
946 if(!cl) return kFALSE;
947 if (cl->GetQ()<=0) return kFALSE;
948 if(lay>1 && cl->GetQ()<=fMinQ) return kFALSE;
950 AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(lay,clu);
951 fPhic = arr->GetPhi();
952 fLambdac = arr->GetLambda();
953 fPhiEstimate = fPhic;
954 fPoint1[0]=primaryVertex[0];
955 fPoint1[1]=primaryVertex[1];
956 fPoint2[0]=arr->GetX();
957 fPoint2[1]=arr->GetY();
961 //________________________________________________________________
962 void AliITStrackerSA::UpdatePoints(){
963 //update of points for the estimation of the curvature
965 fPoint2[0]=fPoint3[0];
966 fPoint2[1]=fPoint3[1];
967 fPoint3[0]=fPointc[0];
968 fPoint3[1]=fPointc[1];
973 //___________________________________________________________________
974 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){
976 //given (x,y) of three recpoints (in global coordinates)
977 //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
979 Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
981 a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
982 b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
983 c = -x1*x1-y1*y1-a*x1-b*y1;
986 //__________________________________________________________________________
987 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){
989 //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
990 //c2 is -rlayer*rlayer
994 Double_t aA = (b1*b1)/(a1*a1)+1;
995 Double_t bB = (-2*m*b1/(a1*a1));
996 Double_t cC = c2+(m*m)/(a1*a1);
997 Double_t dD = bB*bB-4*aA*cC;
1000 y1 = (-bB+TMath::Sqrt(dD))/(2*aA);
1001 y2 = (-bB-TMath::Sqrt(dD))/(2*aA);
1002 x1 = (c2-c1-b1*y1)/a1;
1003 x2 = (c2-c1-b1*y2)/a1;
1007 //____________________________________________________________________
1008 Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t
1009 x2,Double_t y2,Double_t x3,Double_t y3){
1011 //calculates the curvature of track
1012 Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1013 if(den==0) return 0;
1014 Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1015 Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1016 Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1019 if((a*a+b*b-4*c)<0) return 0;
1020 Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1021 if(rad==0) return 0;
1023 if((x1>0 && y1>0 && x1<xc)) rad*=-1;
1024 if((x1<0 && y1>0 && x1<xc)) rad*=-1;
1025 // if((x1<0 && y1<0 && x1<xc)) rad*=-1;
1026 // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1033 //____________________________________________________________________
1034 Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1036 //Returns the point closest to pp
1038 Double_t diff1 = p1-pp;
1039 Double_t diff2 = p2-pp;
1041 if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1042 else fPhiEstimate=p2;
1043 return fPhiEstimate;
1048 //_________________________________________________________________
1049 Int_t AliITStrackerSA::FindTrackLowChiSquare() const {
1050 // returns track with lowest chi square
1051 Int_t dim=fListOfTracks->GetEntries();
1052 if(dim<=1) return 0;
1053 AliITStrackV2* trk = (AliITStrackV2*)fListOfTracks->At(0);
1054 Double_t minChi2=trk->GetChi2();
1056 for(Int_t i=1;i<dim;i++){
1057 trk = (AliITStrackV2*)fListOfTracks->At(i);
1058 Double_t chi2=trk->GetChi2();
1067 //__________________________________________________________
1068 Int_t AliITStrackerSA::FindLabel(AliITStrackV2* track){
1071 Int_t labl[AliITSgeomTGeo::kNLayers][3];
1072 Int_t cnts[AliITSgeomTGeo::kNLayers][3];
1073 for(Int_t j=0;j<AliITSgeomTGeo::GetNLayers();j++){
1074 for(Int_t k=0;k<3;k++){
1080 for(Int_t i=0;i<track->GetNumberOfClusters(); i++) {
1081 Int_t indexc = track->GetClusterIndex(i);
1082 AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc);
1083 Int_t iLayer=cl->GetLayer();
1084 for(Int_t k=0;k<3;k++){
1085 labl[iLayer][k]=cl->GetLabel(k);
1086 if(labl[iLayer][k]<0) iNotLabel++;
1089 if(iNotLabel==3*track->GetNumberOfClusters()) return -2;
1091 for(Int_t j1=0;j1<AliITSgeomTGeo::kNLayers; j1++) {
1092 for(Int_t j2=0; j2<j1; j2++){
1093 for(Int_t k1=0; k1<3; k1++){
1094 for(Int_t k2=0; k2<3; k2++){
1095 if(labl[j1][k1]>=0 && labl[j1][k1]==labl[j2][k2] && cnts[j2][k2]>0){
1107 for(Int_t j=0;j<AliITSgeomTGeo::kNLayers;j++){
1108 for(Int_t k=0;k<3;k++){
1109 if(cnts[j][k]>cntMax && labl[j][k]>=0){
1117 for(Int_t i=0;i<AliITSgeomTGeo::kNLayers;i++)
1118 if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
1120 if(lflag<track->GetNumberOfClusters()) label = -label;
1123 //_____________________________________________________________________________
1124 void AliITStrackerSA::SetCalculatedWindowSizes(Int_t n, Float_t phimin, Float_t phimax, Float_t lambdamin, Float_t lambdamax){
1125 // Set sizes of the phi and lambda windows used for track finding
1127 if(fPhiWin) delete [] fPhiWin;
1128 if(fLambdaWin) delete [] fLambdaWin;
1129 fPhiWin = new Double_t[fNloop];
1130 fLambdaWin = new Double_t[fNloop];
1131 Float_t stepPhi=(phimax-phimin)/(Float_t)(fNloop-1);
1132 Float_t stepLambda=(lambdamax-lambdamin)/(Float_t)(fNloop-1);
1133 for(Int_t k=0;k<fNloop;k++){
1134 Float_t phi=phimin+k*stepPhi;
1135 Float_t lam=lambdamin+k*stepLambda;
1140 //_____________________________________________________________________________
1141 void AliITStrackerSA::SetFixedWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
1142 // Set sizes of the phi and lambda windows used for track finding
1144 if(phi){ // user defined values
1145 fPhiWin = new Double_t[fNloop];
1146 fLambdaWin = new Double_t[fNloop];
1147 for(Int_t k=0;k<fNloop;k++){
1149 fLambdaWin[k]=lam[k];
1152 else { // default values
1154 Double_t phid[33] = {0.002,0.003,0.004,0.0045,0.0047,
1155 0.005,0.0053,0.0055,
1156 0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1157 0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
1158 0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1159 Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1161 0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1162 0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1163 0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
1170 fPhiWin = new Double_t[fNloop];
1171 fLambdaWin = new Double_t[fNloop];
1173 Double_t factor=AliITSReconstructor::GetRecoParam()->GetFactorSAWindowSizes(); // possibility to enlarge windows for cosmics reco with large misalignments (A.Dainese)
1175 for(Int_t k=0;k<fNloop;k++){
1176 fPhiWin[k]=phid[k]*factor;
1177 fLambdaWin[k]=lambdad[k]*factor;
1183 //_______________________________________________________________________
1184 void AliITStrackerSA::GetCoorAngles(AliITSRecPoint* cl,Double_t &phi,Double_t &lambda, Float_t &x, Float_t &y,Float_t &z,Double_t* vertex){
1185 //Returns values of phi (azimuthal) and lambda angles for a given cluster
1187 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
1188 Int_t lay,lad,det; fGeom->GetModuleId(module,lay,lad,det);
1189 Float_t tx,ty,tz; fGeom->GetTrans(lay,lad,det,tx,ty,tz);
1191 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1192 Double_t phi1=TMath::Pi()/2+alpha;
1193 if (lay==1) phi1+=TMath::Pi();
1195 Float_t cp=TMath::Cos(phi1), sp=TMath::Sin(phi1);
1196 Float_t r=tx*cp+ty*sp;
1198 xyz= r*cp - cl->GetY()*sp;
1199 y= r*sp + cl->GetY()*cp;
1203 cl->GetGlobalXYZ(xyz);
1208 phi=TMath::ATan2(y-vertex[1],x-vertex[0]);
1209 lambda=TMath::ATan2(z-vertex[2],TMath::Sqrt((x-vertex[0])*(x-vertex[0])+(y-vertex[1])*(y-vertex[1])));
1212 //________________________________________________________________________
1213 void AliITStrackerSA::GetCoorErrors(AliITSRecPoint* cl,Float_t &sx,Float_t &sy, Float_t &sz){
1215 //returns sigmax, y, z of cluster in global coordinates
1217 Double_t rot[9]; fGeom->GetRotMatrix(module,rot);
1219 AliITSgeomTGeo::GetModuleId(module,lay,lad,det);
1221 Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1222 Double_t phi=TMath::Pi()/2+alpha;
1223 if (lay==1) phi+=TMath::Pi();
1225 Float_t cp=TMath::Cos(phi), sp=TMath::Sin(phi);
1228 cl->GetGlobalCov(covm);
1229 sx=TMath::Sqrt(covm[0]);
1230 sy=TMath::Sqrt(covm[3]);
1231 sz=TMath::Sqrt(covm[5]);
1233 sx = TMath::Sqrt(sp*sp*cl->GetSigmaY2());
1234 sy = TMath::Sqrt(cp*cp*cl->GetSigmaY2());
1235 sz = TMath::Sqrt(cl->GetSigmaZ2());