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 ////////////////////////////////////////////////////
30 #include <TObjArray.h>
34 #include "AliESDVertex.h"
35 #include "AliESDtrack.h"
36 #include "AliITSRiemannFit.h"
37 #include "AliITSVertexer.h"
38 #include "AliITSclusterTable.h"
39 #include "AliITSclusterV2.h"
40 #include "AliITSgeom.h"
41 #include "AliITStrackSA.h"
42 #include "AliITStrackerSA.h"
45 ClassImp(AliITStrackerSA)
47 //____________________________________________________________________________
48 AliITStrackerSA::AliITStrackerSA():AliITStrackerMI(){
49 // Default constructor
53 //____________________________________________________________________________
54 AliITStrackerSA::AliITStrackerSA(AliITSgeom *geom):AliITStrackerMI(geom)
56 // Standard constructor (Vertex is known and passed to this obj.)
63 //____________________________________________________________________________
64 AliITStrackerSA::AliITStrackerSA(AliITSgeom *geom, AliESDVertex *vert):AliITStrackerMI(geom)
66 // Standard constructor (Vertex is known and passed to this obj.)
73 //____________________________________________________________________________
74 AliITStrackerSA::AliITStrackerSA(AliITSgeom *geom, AliITSVertexer *vertexer):AliITStrackerMI(geom)
76 // Standard constructor (Vertex is unknown - vertexer is passed to this obj)
83 //____________________________________________________________________________
84 AliITStrackerSA::AliITStrackerSA(AliITStrackerSA& tracker):AliITStrackerMI(){
86 fPhiEstimate = tracker.fPhiEstimate;
87 for(Int_t i=0;i<2;i++){
88 fPoint1[i]=tracker.fPoint1[i];
89 fPoint2[i]=tracker.fPoint2[i];
90 fPoint3[i]=tracker.fPoint3[i];
91 fPointc[i]=tracker.fPointc[i];
93 fLambdac = tracker.fLambdac;
94 fPhic = tracker.fPhic;
95 fCoef1 = tracker.fCoef1;
96 fCoef2 = tracker.fCoef2;
97 fCoef3 = tracker.fCoef3;
98 fNloop = tracker.fNloop;
99 fPhiWin = tracker.fPhiWin;
100 fLambdaWin = tracker.fLambdaWin;
101 if(tracker.fVertexer && tracker.fVert){
102 fVert = new AliESDVertex(*tracker.fVert);
105 fVert = tracker.fVert;
107 fVertexer = tracker.fVertexer;
108 fGeom = tracker.fGeom;
109 fTable = tracker.fTable;
110 fListOfTracks = tracker.fListOfTracks;
113 //____________________________________________________________________________
114 AliITStrackerSA::~AliITStrackerSA(){
116 // if fVertexer is not null, the AliESDVertex obj. is owned by this class
117 // and is deleted here
119 if(fVert)delete fVert;
124 if(fPhiWin)delete []fPhiWin;
125 if(fLambdaWin)delete []fLambdaWin;
127 fListOfTracks->Delete();
130 //____________________________________________________________________________
131 void AliITStrackerSA::Init(){
132 // Reset all data members
134 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
149 fListOfTracks=new TObjArray(0,0);
151 //_______________________________________________________________________
152 void AliITStrackerSA::ResetForFinding(){
153 // Reset data members used in all loops during track finding
155 for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
163 fListOfTracks->Delete();
165 //____________________________________________________________________________
166 void AliITStrackerSA::FindTracks(TTree *out,Int_t evnumber){
168 /**************************************************************************
169 * This function finds primary tracks.
172 * Example: to execute function with only the ITS (no combined tracking *
173 * with TPC+ITS) and requiring 5/6 points to define a good track *
174 * call SetSixPoinbts(kFALSE) in advance and then *
175 * use: FindTracks(treein,treeout,evnumber) *
176 * to execute combined tracking, before using FindTracks, use *
178 *************************************************************************/
181 Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
186 if(fVert)delete fVert;
187 fVert = fVertexer->FindVertexForCurrentEvent(evnumber);
190 gAlice->GetEvent(evnumber);
192 Fatal("FindTracks","Vertex is missing\n");
196 Double_t primaryVertex[3];
197 Double_t errorsprimvert[3];
198 fVert->GetXYZ(primaryVertex);
199 fVert->GetSigmaXYZ(errorsprimvert);
200 if(errorsprimvert[0]==0 || errorsprimvert[1]==0){
201 Warning("FindTracks","Set errors on vertex positions x and y at 0.0001");
202 errorsprimvert[0]=0.0001;
203 errorsprimvert[1]=0.0001;
205 fVert->PrintStatus();
208 //Fill array with cluster indices for each module
210 fTable = new AliITSclusterTable(fGeom,this,primaryVertex);
211 fTable->FillArray(fITSclusters);
212 fTable->FillArrayCoorAngles();
216 //Fill tree for found tracks
217 AliITStrackV2* outrack=0;
218 TBranch* branch=out->Branch("tracks","AliITStrackV2",&outrack,32000,0);
219 if (!branch) out->Branch("tracks","AliITStrackV2",&outrack,32000,3);
220 else branch->SetAddress(&outrack);
223 Int_t * firstmod = new Int_t[fGeom->GetNlayers()];
224 for(Int_t i=0;i<fGeom->GetNlayers();i++){
225 firstmod[i]=fGeom->GetModuleIndex(i+1,1,1);
227 // firstmod [i] number of the first module in the ITS layer i.
229 AliITSlayer &layer=fgLayers[0]; // first layer
231 Int_t dim=layer.GetNumberOfClusters();
232 //loop on the different windows
233 for(Int_t nloop=0;nloop<fNloop;nloop++){
234 for(Int_t ncl=0;ncl<dim;ncl++){
235 //loop starting from layer 0
238 AliITSclusterV2* cl = layer.GetCluster(ncl);
239 if(cl->IsUsed()==1) continue;
240 if(cl->TestBit(kSAflag)==kTRUE) continue;
242 fPhic = fTable->GetPhiCluster(0,ncl);
243 fLambdac = fTable->GetLambdaCluster(0,ncl);
244 fPhiEstimate = fPhic;
245 AliITStrackSA* trs = new AliITStrackSA();
246 fPoint1[0]=primaryVertex[0];
247 fPoint1[1]=primaryVertex[1];
248 fPoint2[0]=fTable->GetXCluster(0,ncl);
249 fPoint2[1]=fTable->GetYCluster(0,ncl);
251 Int_t * nn = new Int_t[fGeom->GetNlayers()];//counter for clusters on each layer
252 for(Int_t i=0;i<fGeom->GetNlayers();i++){ nn[i]=0;}
253 nn[0] = SearchClusters(0,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
254 nn[1] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
258 fPoint3[0] = fPointc[0];
259 fPoint3[1] = fPointc[1];
261 nn[2] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
262 if(nn[1]==0 && nn[2]==0) pflag=0;
263 if(nn[2]!=0 && nn[1]!=0){ pflag=1; UpdatePoints();}
264 if(nn[2]!=0 && nn[1]==0){
266 fPoint3[0]=fPointc[0];
267 fPoint3[1]=fPointc[1];
270 nn[3] = SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
272 if(nn[3]!=0) UpdatePoints();
273 nn[4] = SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
275 if(nn[4]!=0) UpdatePoints();
276 nn[5] = SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
280 Int_t numberofpoints;
281 if(fSixPoints) numberofpoints=6; //check of the candidate track
282 else numberofpoints=5; //if track is good (with the required number
283 for(Int_t nnp=0;nnp<fGeom->GetNlayers();nnp++){ //of points) it is written on file
284 if(nn[nnp]!=0) layOK+=1;
286 if(layOK>=numberofpoints){
287 AliITStrackV2* tr2 = FitTrack(trs,primaryVertex,errorsprimvert);
289 Int_t nct = trs->GetNumberOfClustersSA();
291 Int_t index = trs->GetClusterIndexSA(nct);
292 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
293 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
302 Int_t nct = tr2->GetNumberOfClusters();
305 Int_t index = tr2->GetClusterIndex(nct);
306 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
312 Int_t nct = trs->GetNumberOfClustersSA();
314 Int_t index = trs->GetClusterIndexSA(nct);
315 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
316 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
322 }//end loop on clusters of layer1
326 //if 5/6 points are required, second loop starting
327 //from second layer, to find tracks with point of
331 // counter for clusters on each layer
332 Int_t * nn = new Int_t[fGeom->GetNlayers()-1];
333 for(Int_t nloop=0;nloop<fNloop;nloop++){
334 AliITSlayer &layer2=fgLayers[1]; //loop on layer 2
335 Int_t ncl2=layer2.GetNumberOfClusters();
336 while(ncl2--){ //loop starting from layer 2
339 AliITSclusterV2* cl = layer2.GetCluster(ncl2);
340 if(cl->IsUsed()==1) continue;
341 if(cl->TestBit(kSAflag)==kTRUE) continue;
342 fPhic = fTable->GetPhiCluster(1,ncl2);
343 fLambdac = fTable->GetLambdaCluster(1,ncl2);
344 fPhiEstimate = fPhic;
345 AliITStrackSA* trs = new AliITStrackSA();
347 fPoint1[0]=primaryVertex[0];
348 fPoint1[1]=primaryVertex[1];
349 fPoint2[0]=fTable->GetXCluster(1,ncl2);;
350 fPoint2[1]=fTable->GetYCluster(1,ncl2);;
352 for(Int_t kk=0;kk<fGeom->GetNlayers()-1;kk++)nn[kk] = 0;
353 for(Int_t kk=0;kk<fGeom->GetNlayers()-1;kk++){
354 nn[kk] = SearchClusters(kk+1,fPhiWin[nloop],fLambdaWin[nloop],
355 trs,primaryVertex[2],pflag);
363 for(Int_t nnp=0;nnp<fGeom->GetNlayers()-1;nnp++){
364 if(nn[nnp]!=0) fl+=1;
367 AliITStrackV2* tr2 = FitTrack(trs,primaryVertex,errorsprimvert);
369 Int_t nct = trs->GetNumberOfClustersSA();
371 Int_t index = trs->GetClusterIndexSA(nct);
372 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
373 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
380 Int_t nct = tr2->GetNumberOfClusters();
382 Int_t index = tr2->GetClusterIndex(nct);
383 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
389 Int_t nct = trs->GetNumberOfClustersSA();
391 Int_t index = trs->GetClusterIndexSA(nct);
392 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
394 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
398 }//end loop on clusters of layer2
401 } // if(!fSixPoints....
404 delete fTable; fTable=0;
408 //______________________________________________________________________
409 Int_t AliITStrackerSA::FindTracks(AliESD* event){
411 // Track finder using the ESD object
415 Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
420 Double_t errorsprimvert[3];
421 Double_t primaryVertex[3];
422 event->GetVertex()->GetXYZ(primaryVertex);
423 event->GetVertex()->GetSigmaXYZ(errorsprimvert);
425 if(errorsprimvert[0]==0 || errorsprimvert[1]==0){
426 // Warning("FindTracks","Set errors on vertex positions x and y at 0.005");
427 errorsprimvert[0]=0.005;
428 errorsprimvert[1]=0.005;
431 //Fill array with cluster indices for each module
433 fTable = new AliITSclusterTable(fGeom,this,primaryVertex);
434 fTable->FillArray(fITSclusters);
435 fTable->FillArrayCoorAngles();
438 Int_t * firstmod = new Int_t[fGeom->GetNlayers()];
439 for(Int_t i=0;i<fGeom->GetNlayers();i++){
440 firstmod[i]=fGeom->GetModuleIndex(i+1,1,1);
442 // firstmod [i] number of the first module in the ITS layer i.
445 AliITSlayer &layer=fgLayers[0];
447 Int_t dim=layer.GetNumberOfClusters();
448 //loop on the different windows
449 for(Int_t nloop=0;nloop<fNloop;nloop++){
450 for(Int_t ncl=0;ncl<dim;ncl++){ //loop starting from layer 0
454 AliITSclusterV2* cl = (AliITSclusterV2*)layer.GetCluster(ncl);
456 if(cl->IsUsed()==1) continue;
457 if(cl->TestBit(kSAflag)==kTRUE) continue;
458 if (cl->GetQ()<=0) continue;
460 fPhic = fTable->GetPhiCluster(0,ncl);
461 fLambdac = fTable->GetLambdaCluster(0,ncl);
463 if (TMath::Abs(fLambdac)>0.26*TMath::Pi()) continue;
465 fPhiEstimate = fPhic;
466 AliITStrackSA* trs = new AliITStrackSA();
467 fPoint1[0]=primaryVertex[0];
468 fPoint1[1]=primaryVertex[1];
469 fPoint2[0]=fTable->GetXCluster(0,ncl);
470 fPoint2[1]=fTable->GetYCluster(0,ncl);
471 Int_t * nn = new Int_t[fGeom->GetNlayers()];//counter for clusters on each layer
472 for(Int_t i=0;i<fGeom->GetNlayers();i++){ nn[i]=0;}
473 nn[0] = SearchClusters(0,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
475 nn[1] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
478 fPoint3[0] = fPointc[0];
479 fPoint3[1] = fPointc[1];
481 nn[2] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
482 if(nn[1]==0 && nn[2]==0) pflag=0;
483 if(nn[2]!=0 && nn[1]!=0){ pflag=1; UpdatePoints();}
484 if(nn[2]!=0 && nn[1]==0){
486 fPoint3[0]=fPointc[0];
487 fPoint3[1]=fPointc[1];
490 nn[3] = SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
492 if(nn[3]!=0) UpdatePoints();
493 nn[4] = SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
495 if(nn[4]!=0) UpdatePoints();
496 nn[5] = SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
500 Int_t numberofpoints;
501 if(fSixPoints) numberofpoints=6; //check of the candidate track
502 else numberofpoints=5; //if track is good (with the required number
503 for(Int_t nnp=0;nnp<fGeom->GetNlayers();nnp++){ //of points) it is written on file
504 if(nn[nnp]!=0) layOK+=1;
506 if(layOK>=numberofpoints){
507 AliITStrackV2* tr2 = FitTrack(trs,primaryVertex,errorsprimvert);
509 Int_t nct = trs->GetNumberOfClustersSA();
511 Int_t index = trs->GetClusterIndexSA(nct);
512 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
513 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
519 AliESDtrack outtrack;
520 outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
521 event->AddTrack(&outtrack);
523 Int_t nct = tr2->GetNumberOfClusters();
525 Int_t index = tr2->GetClusterIndex(nct);
526 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
532 Int_t nct = trs->GetNumberOfClustersSA();
534 Int_t index = trs->GetClusterIndexSA(nct);
535 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
536 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
543 }//end loop on clusters of layer1
549 //if 5/6 points are required, second loop starting
550 //from second layer, to find tracks with point of
554 // counter for clusters on each layer
555 Int_t * nn = new Int_t[fGeom->GetNlayers()-1];
556 for(Int_t nloop=0;nloop<fNloop;nloop++){
557 AliITSlayer &layer2=fgLayers[1];
558 Int_t ncl2=layer2.GetNumberOfClusters();
559 while(ncl2--){ //loop starting from layer 2
562 AliITSclusterV2* cl = layer2.GetCluster(ncl2);
563 if(cl->IsUsed()==1) continue;
564 if(cl->TestBit(kSAflag)==kTRUE) continue;
565 fPhic = fTable->GetPhiCluster(1,ncl2);
566 fLambdac = fTable->GetLambdaCluster(1,ncl2);
567 fPhiEstimate = fPhic;
568 AliITStrackSA* trs = new AliITStrackSA();
569 fPoint1[0]=primaryVertex[0];
570 fPoint1[1]=primaryVertex[1];
571 fPoint2[0]=fTable->GetXCluster(1,ncl2);
572 fPoint2[1]=fTable->GetYCluster(1,ncl2);
574 for(Int_t kk=0;kk<fGeom->GetNlayers()-1;kk++)nn[kk] = 0;
575 for(Int_t kk=0;kk<fGeom->GetNlayers()-1;kk++){
576 nn[kk] = SearchClusters(kk+1,fPhiWin[nloop],fLambdaWin[nloop],
577 trs,primaryVertex[2],pflag);
585 for(Int_t nnp=0;nnp<fGeom->GetNlayers()-1;nnp++){
586 if(nn[nnp]!=0) fl+=1;
589 AliITStrackV2* tr2 = FitTrack(trs,primaryVertex,errorsprimvert);
591 Int_t nct = trs->GetNumberOfClustersSA();
593 Int_t index = trs->GetClusterIndexSA(nct);
594 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
595 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
602 AliESDtrack outtrack;
603 outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
604 event->AddTrack(&outtrack);
606 Int_t nct = tr2->GetNumberOfClusters();
608 Int_t index = tr2->GetClusterIndex(nct);
609 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
615 Int_t nct = trs->GetNumberOfClustersSA();
617 Int_t index = trs->GetClusterIndexSA(nct);
618 AliITSclusterV2* kl = (AliITSclusterV2*)GetCluster(index);
620 if(kl->TestBit(kSAflag)==kTRUE) kl->ResetBit(kSAflag);
625 }//end loop on clusters of layer2
631 delete fTable;fTable=0;
632 Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
638 //________________________________________________________________________
640 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex,Double_t *errorsprimvert){
641 //fit of the found track
644 Int_t * firstmod = new Int_t[fGeom->GetNlayers()];
645 for(Int_t i=0;i<fGeom->GetNlayers();i++){
646 firstmod[i]=fGeom->GetModuleIndex(i+1,1,1);
649 Int_t nclusters = tr->GetNumberOfClustersSA();
650 TObjArray** listlayer = new TObjArray*[fGeom->GetNlayers()];
651 for(Int_t i=0;i<fGeom->GetNlayers();i++){
652 listlayer[i] = new TObjArray(0,0);
662 Int_t * nnn = new Int_t[fGeom->GetNlayers()];
663 for(Int_t i=0;i<fGeom->GetNlayers();i++)nnn[i]=0;
665 for(Int_t ncl=0;ncl<nclusters;ncl++){
666 Int_t index = tr->GetClusterIndexSA(ncl);
667 AliITSclusterV2* cl = (AliITSclusterV2*)GetCluster(index);
669 if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
670 Int_t lay = (index & 0xf0000000) >> 28;
671 if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
672 if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
673 if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
674 if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
675 if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
676 if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}
681 Int_t * end = new Int_t[fGeom->GetNlayers()];
682 for(Int_t i=0;i<fGeom->GetNlayers();i++){
683 if(listlayer[i]->GetEntries()==0) end[i]=1;
684 else end[i]=listlayer[i]->GetEntries();
687 for(Int_t l1=0;l1<end[0];l1++){//loop on layer 1
688 AliITSclusterV2* cl0 = (AliITSclusterV2*)listlayer[0]->At(l1);
689 TVector3** recp = new TVector3*[3];
690 TVector3** errs = new TVector3*[3];
691 recp[0] = new TVector3(primaryVertex[0],primaryVertex[1],primaryVertex[2]);
692 errs[0] = new TVector3(errorsprimvert[0],errorsprimvert[1],errorsprimvert[2]);
693 Double_t x1,y1,z1,sx1,sy1,sz1;
694 Double_t x2,y2,z2,sx2,sy2,sz2;
695 AliITSclusterV2* p1=0;
696 AliITSclusterV2* p2=0;
697 Int_t index1=clind0[l1];
699 for(Int_t l2=0;l2<end[1];l2++){//loop on layer 2
700 AliITSclusterV2* cl1 = (AliITSclusterV2*)listlayer[1]->At(l2);
702 for(Int_t l3=0;l3<end[2];l3++){ //loop on layer 3
703 AliITSclusterV2* cl2 = (AliITSclusterV2*)listlayer[2]->At(l3);
705 if(cl0==0 && cl1!=0) {
706 p2 = cl2;index1=clind2[l3];
710 if(cl0!=0 && cl1==0){
712 p2=cl2;index2=clind2[l3];
714 if(cl0!=0 && cl1!=0){
718 Int_t lay1=(index1 & 0xf0000000) >> 28;
719 Int_t cln1=(index1 & 0x0fffffff) >> 00;
720 Int_t lay2=(index2 & 0xf0000000) >> 28;
721 Int_t cln2=(index2 & 0x0fffffff) >> 00;
722 x1 = fTable->GetXCluster(lay1,cln1);
723 x2 = fTable->GetXCluster(lay2,cln2);
724 y1 = fTable->GetYCluster(lay1,cln1);
725 y2 = fTable->GetYCluster(lay2,cln2);
726 z1 = fTable->GetZCluster(lay1,cln1);
727 z2 = fTable->GetZCluster(lay2,cln2);
728 sx1 = fTable->GetXClusterError(lay1,cln1);
729 sx2 = fTable->GetXClusterError(lay2,cln2);
730 sy1 = fTable->GetYClusterError(lay1,cln1);
731 sy2 = fTable->GetYClusterError(lay2,cln2);
732 sz1 = fTable->GetZClusterError(lay1,cln1);
733 sz2 = fTable->GetZClusterError(lay2,cln2);
734 Double_t phi1 = fTable->GetPhiCluster(lay1,cln1);
735 Int_t module1 = p1->GetDetectorIndex()+firstmod[0];
736 recp[1] = new TVector3(x1,y1,z1);
737 errs[1] = new TVector3(sx1,sy1,sz1);
738 recp[2] = new TVector3(x2,y2,z2);
739 errs[2] = new TVector3(sx2,sy2,sz2);
741 //fit on the Riemann sphere
742 Float_t seed1,seed2,seed3;
743 AliITSRiemannFit fit;
744 Int_t rf = fit.FitHelix(3,recp,errs,seed1,seed2,seed3); //this gives phi,tgl,curvature to start Kalman Filter
746 for(Int_t i=1;i<3;i++){
756 if(seed1>-TMath::Pi() && seed1<-0.5*TMath::Pi()){
757 phi=seed1+1.5*TMath::Pi();
760 if(seed1>-0.5*TMath::Pi() && seed1<0.5*TMath::Pi()){
761 phi=seed1+0.5*TMath::Pi();
764 if(seed1>0.5*TMath::Pi() && seed1<TMath::Pi()){
765 phi=seed1-0.5*TMath::Pi();
770 if(seed1>-TMath::Pi() && seed1<-0.5*TMath::Pi()){
771 phi=seed1+0.5*TMath::Pi();
774 if(seed1>-0.5*TMath::Pi() && seed1<0.5*TMath::Pi()){
775 phi=seed1-0.5*TMath::Pi();
778 if(seed1>0.5*TMath::Pi() && seed1<TMath::Pi()){
779 phi=seed1-1.5*TMath::Pi();
784 Int_t layer,ladder,detector;
785 fGeom->GetModuleId(module1,layer,ladder,detector);
786 Float_t yclu1 = p1->GetY();
787 Float_t zclu1 = p1->GetZ();
788 Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
790 for(Int_t l4=0;l4<end[3];l4++){ //loop on layer 4
791 AliITSclusterV2* cl3 = (AliITSclusterV2*)listlayer[3]->At(l4);
792 for(Int_t l5=0;l5<end[4];l5++){ //loop on layer 5
793 AliITSclusterV2* cl4 = (AliITSclusterV2*)listlayer[4]->At(l5);
794 for(Int_t l6=0;l6<end[5];l6++){ //loop on layer 6
795 AliITSclusterV2* cl5 = (AliITSclusterV2*)listlayer[5]->At(l6);
796 AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi,tgl,cv,1);
798 if(cl5!=0) trac->AddClusterV2(5,(clind5[l6] & 0x0fffffff)>>0);
799 if(cl4!=0) trac->AddClusterV2(4,(clind4[l5] & 0x0fffffff)>>0);
800 if(cl3!=0) trac->AddClusterV2(3,(clind3[l4] & 0x0fffffff)>>0);
801 if(cl2!=0) trac->AddClusterV2(2,(clind2[l3] & 0x0fffffff)>>0);
802 if(cl1!=0) trac->AddClusterV2(1,(clind1[l2] & 0x0fffffff)>>0);
803 if(cl0!=0) trac->AddClusterV2(0,(clind0[l1] & 0x0fffffff)>>0);
805 //fit with Kalman filter using AliITStrackerMI::RefitAt()
807 AliITStrackMI* ot = new AliITStrackSA(*trac);
809 ot->ResetCovariance();
812 if(RefitAt(49.,ot,trac)){ //fit from layer 1 to layer 6
814 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
815 otrack2->ResetCovariance();
816 otrack2->ResetClusters();
817 //fit from layer 6 to layer 1
818 if(RefitAt(3.7,otrack2,ot)) {
819 fListOfTracks->AddLast(otrack2);
832 for(Int_t i=1;i<3;i++){
846 Int_t dim=fListOfTracks->GetEntries();
848 for(Int_t i=0;i<fGeom->GetNlayers();i++){
856 AliITStrackV2* otrack =(AliITStrackV2*)FindTrackLowChiSquare(fListOfTracks,dim);
859 for(Int_t i=0;i<fGeom->GetNlayers();i++){
866 Int_t * indexc = new Int_t[fGeom->GetNlayers()];
867 for(Int_t i=0;i<fGeom->GetNlayers();i++) indexc[i]=0;
868 for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
869 indexc[nind] = otrack->GetClusterIndex(nind);
871 AliITSclusterV2* cl0 = (AliITSclusterV2*)GetCluster(indexc[0]);
872 AliITSclusterV2* cl1 = (AliITSclusterV2*)GetCluster(indexc[1]);
873 AliITSclusterV2* cl2 = (AliITSclusterV2*)GetCluster(indexc[2]);
874 AliITSclusterV2* cl3 = (AliITSclusterV2*)GetCluster(indexc[3]);
875 AliITSclusterV2* cl4 = (AliITSclusterV2*)GetCluster(indexc[4]);
876 Int_t labl[3]={-1,-1,-1};
877 if(otrack->GetNumberOfClusters()==fGeom->GetNlayers()){
878 AliITSclusterV2* cl5 = (AliITSclusterV2*)GetCluster(indexc[5]);
879 labl[0]=cl5->GetLabel(0);
880 labl[1]=cl5->GetLabel(1);
881 labl[2]=cl5->GetLabel(2);
884 if(otrack->GetNumberOfClusters()==(fGeom->GetNlayers()-1)){
889 Int_t numberofpoints;
890 if(fSixPoints) numberofpoints=6;
891 else numberofpoints=5;
892 CookLabel(otrack,0.); //MI change - to see fake ratio
893 Int_t label = Label(cl0->GetLabel(0),cl1->GetLabel(0),
894 cl2->GetLabel(0),cl3->GetLabel(0),
895 cl4->GetLabel(0),labl[0],
896 cl0->GetLabel(1),cl1->GetLabel(1),
897 cl2->GetLabel(1),cl3->GetLabel(1),
898 cl4->GetLabel(1),labl[1],
899 cl0->GetLabel(2),cl1->GetLabel(2),
900 cl2->GetLabel(2),cl3->GetLabel(2),
901 cl4->GetLabel(2),labl[2],numberofpoints);
903 otrack->SetLabel(label);
904 for(Int_t i=0;i<fGeom->GetNlayers();i++){
913 //_______________________________________________________________________
914 void AliITStrackerSA::UseFoundTracksV2(Int_t evnum,TTree* treev2){
915 // Marks as used clusters belonging to tracks found with V2 TPC+ITS tracking
916 //(or AliITStrackV2 tracks found with function FindTracks of this class)
921 if(fVert)delete fVert;
922 fVert = fVertexer->FindVertexForCurrentEvent(evnum);
925 gAlice->GetEvent(evnum);
927 Fatal("FindTracks","Vertex is missing\n");
931 Double_t primaryVertex[3];
932 fVert->GetXYZ(primaryVertex);
935 fTable = new AliITSclusterTable(fGeom,this,primaryVertex);
936 fTable->FillArray(fITSclusters);
937 fTable->FillArrayCoorAngles();
940 TBranch* bra = (TBranch*)treev2->GetBranch("tracks");
941 if(!bra) Warning("UseFoundTracksV2","No branch for track tree");
942 AliITStrackV2* ttrrt = new AliITStrackV2;
943 bra->SetAddress(&ttrrt);
945 for(Int_t nj=0;nj<treev2->GetEntries();nj++){
946 treev2->GetEvent(nj);
947 Int_t ncl = ttrrt->GetNumberOfClusters();
948 for(Int_t k=0;k<ncl;k++){
949 Int_t index = ttrrt->GetClusterIndex(k);
950 AliITSclusterV2* clui = (AliITSclusterV2*)GetCluster(index);
951 if(clui->IsUsed()==0) clui->Use();
959 //_______________________________________________________________________
960 void AliITStrackerSA::UseFoundTracksV2(AliESD *event){
961 // Marks as used clusters belonging to tracks found with V2 TPC+ITS tracking
965 Double_t primaryVertex[3];
966 event->GetVertex()->GetXYZ(primaryVertex);
969 fTable = new AliITSclusterTable(fGeom,this,primaryVertex);
970 fTable->FillArray(fITSclusters);
971 fTable->FillArrayCoorAngles();
974 Int_t ntracks = event->GetNumberOfTracks();
976 AliESDtrack *esd=event->GetTrack(ntracks);
977 if ((esd->GetStatus()&
978 AliESDtrack::kITSin|AliESDtrack::kTPCin)==0) continue;
980 Int_t ncl = esd->GetITSclusters(idx);
981 for(Int_t clu=0; clu<ncl; clu++){
982 AliITSclusterV2* cl = (AliITSclusterV2*)GetCluster(idx[clu]);
983 if(cl->IsUsed()==0) cl->Use();
987 Info("UseFoundTracksV2","Clusters of tracks prolonged from TPC deleted");
993 //_______________________________________________________
994 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t zvertex,Int_t pflag){
995 //function used to to find the clusters associated to the track
997 AliITSlayer &lay = fgLayers[layer];
998 Int_t * firstmod = new Int_t[fGeom->GetNlayers()];
999 for(Int_t i=0;i<fGeom->GetNlayers();i++){
1000 firstmod[i]=fGeom->GetModuleIndex(i+1,1,1);
1004 Float_t cx1,cx2,cy1,cy2;
1005 FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
1006 Int_t fun = FindIntersection(fCoef1,fCoef2,fCoef3,-(lay.GetR()*lay.GetR()),cx1,cy1,cx2,cy2);
1012 Double_t fi1 =TMath::ATan2(cy1,cx1);
1013 Double_t fi2 =TMath::ATan2(cy2,cx2);
1014 fPhiEstimate = ChoosePoint(fi1,fi2,fPhic);
1017 Double_t zed = TMath::Tan(fLambdac)*lay.GetR()+zvertex;
1018 Double_t zed1 = TMath::Tan(fLambdac+lambdawindow)*lay.GetR()+zvertex;
1019 Double_t zed2 = TMath::Tan(fLambdac-lambdawindow)*lay.GetR()+zvertex;
1021 Double_t fi = fPhiEstimate;
1022 Int_t nmod = lay.FindDetectorIndex(fi,zed);
1027 nmod += firstmod[layer];
1029 Int_t nm[8]={0,0,0,0,0,0,0,0};
1030 nm[0] = lay.FindDetectorIndex(fi+phiwindow,zed);
1031 nm[1] = lay.FindDetectorIndex(fi-phiwindow,zed);
1032 nm[2] = lay.FindDetectorIndex(fi,zed1);
1033 nm[3] = lay.FindDetectorIndex(fi,zed2);
1034 nm[4] = lay.FindDetectorIndex(fi+phiwindow,zed1);
1035 nm[5] = lay.FindDetectorIndex(fi-phiwindow,zed1);
1036 nm[6] = lay.FindDetectorIndex(fi+phiwindow,zed2);
1037 nm[7] = lay.FindDetectorIndex(fi-phiwindow,zed2);
1041 TArrayI* array =(TArrayI*)fTable->GetListOfClusters(nmod);
1042 TArrayI* listc = new TArrayI(array->GetSize());
1043 for(Int_t i=0;i<array->GetSize();i++){
1044 Int_t in=(Int_t)array->At(i);
1045 listc->AddAt(in,nn);
1052 for(Int_t h=k+1;h<8;h++){
1065 for(Int_t ii=0;ii<8;ii++){
1066 if(nm[ii]!=value && nm[ii]!=nmod && nm[ii]>=0){
1067 TArrayI* ar =(TArrayI*)fTable->GetListOfClusters(nm[ii]+firstmod[layer]);
1068 listc->Set(listc->GetSize()+ar->GetSize());
1069 for(Int_t j=0;j<ar->GetSize();j++){
1070 Int_t in=(Int_t)ar->At(j);
1071 listc->AddAt(in,nn);
1079 for(Int_t i=0;i<listc->GetSize();i++){
1080 Int_t index = (Int_t)listc->At(i);
1081 AliITSclusterV2* cllay = lay.GetCluster(index);
1082 if(cllay==0) continue;
1083 if(cllay->IsUsed()==1) continue;
1084 if(cllay->TestBit(kSAflag)==kTRUE) continue;
1085 Double_t phi = fTable->GetPhiCluster(layer,index);
1086 Double_t lambda= fTable->GetLambdaCluster(layer,index);
1088 if(TMath::Abs(fLambdac-lambda)<lambdawindow &&
1089 TMath::Abs(fPhiEstimate-phi)<phiwindow){
1092 if(trs->GetNumberOfClustersSA()==15){
1097 trs->AddClusterSA(layer,index);
1098 cllay->SetBit(kSAflag);
1100 fPointc[0]=fTable->GetXCluster(layer,index);
1101 fPointc[1]=fTable->GetYCluster(layer,index);
1113 //_______________________________________________________
1114 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t zvertex,Int_t pflag){
1115 //function used to to find the clusters associated to the track
1117 AliITSlayer &lay = fgLayers[layer];
1118 Double_t r=lay.GetR(),tgl=TMath::Tan(fLambdac);
1121 Float_t cx1,cx2,cy1,cy2;
1122 FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
1123 if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
1125 Double_t fi1=TMath::ATan2(cy1,cx1);
1126 Double_t fi2=TMath::ATan2(cy2,cx2);
1127 fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
1130 Double_t dz=r*lambdawindow*TMath::Sqrt(1+tgl*tgl) + 0.3*TMath::Abs(tgl);
1131 Double_t zmax=r*tgl + zvertex + dz;
1132 Double_t zmin=r*tgl + zvertex - dz;
1134 Int_t indmin = lay.FindClusterIndex(zmin);
1135 Int_t indmax = lay.FindClusterIndex(zmax);
1136 for (Int_t index=indmin; index<indmax; index++) {
1137 AliITSclusterV2 *c=lay.GetCluster(index);
1138 if (c->IsUsed()) continue;
1139 if (c->GetQ()<=0) continue;
1140 if (c->TestBit(kSAflag)==kTRUE) continue;
1142 Double_t phi =fTable->GetPhiCluster(layer,index);
1143 if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
1145 Double_t lambda=fTable->GetLambdaCluster(layer,index);
1146 if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
1148 if(trs->GetNumberOfClustersSA()==15) return 0;
1150 trs->AddClusterSA(layer,index);
1154 fPointc[0]=fTable->GetXCluster(layer,index);
1155 fPointc[1]=fTable->GetYCluster(layer,index);
1163 //________________________________________________________________
1164 void AliITStrackerSA::UpdatePoints(){
1165 //update of points for the estimation of the curvature
1167 //fPoint1[0]=fPoint2[0];
1168 //fPoint1[1]=fPoint2[1];
1169 fPoint2[0]=fPoint3[0];
1170 fPoint2[1]=fPoint3[1];
1171 fPoint3[0]=fPointc[0];
1172 fPoint3[1]=fPointc[1];
1177 //___________________________________________________________________
1178 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){
1180 //given (x,y) of three recpoints (in global coordinates)
1181 //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
1183 Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1184 if(den==0) return 0;
1185 a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1186 b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1187 c = -x1*x1-y1*y1-a*x1-b*y1;
1190 //__________________________________________________________________________
1191 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){
1193 //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
1194 //c2 is -rlayer*rlayer
1198 Double_t aA = (b1*b1)/(a1*a1)+1;
1199 Double_t bB = (-2*m*b1/(a1*a1));
1200 Double_t cC = c2+(m*m)/(a1*a1);
1201 Double_t dD = bB*bB-4*aA*cC;
1204 y1 = (-bB+TMath::Sqrt(dD))/(2*aA);
1205 y2 = (-bB-TMath::Sqrt(dD))/(2*aA);
1206 x1 = (c2-c1-b1*y1)/a1;
1207 x2 = (c2-c1-b1*y2)/a1;
1211 //____________________________________________________________________
1212 Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t
1213 x2,Double_t y2,Double_t x3,Double_t y3){
1215 //calculates the curvature of track
1216 Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1217 if(den==0) return 0;
1218 Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1219 Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1220 Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1223 if((a*a+b*b-4*c)<0) return 0;
1224 Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1225 if(rad==0) return 0;
1227 if((x1>0 && y1>0 && x1<xc)) rad*=-1;
1228 if((x1<0 && y1>0 && x1<xc)) rad*=-1;
1229 // if((x1<0 && y1<0 && x1<xc)) rad*=-1;
1230 // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1235 //____________________________________________________________________
1236 Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1238 //Returns the point closest to pp
1240 Double_t diff1 = p1-pp;
1241 Double_t diff2 = p2-pp;
1243 if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1244 else fPhiEstimate=p2;
1245 return fPhiEstimate;
1250 //_________________________________________________________________
1251 AliITStrackV2* AliITStrackerSA::FindTrackLowChiSquare(TObjArray* tracklist, Int_t dim) const {
1252 // returns track with lowes chi square
1254 AliITStrackV2* trk = (AliITStrackV2*)tracklist->At(0);
1257 if(dim==0) return 0;
1258 Double_t * chi2 = new Double_t[dim];
1259 Int_t * index = new Int_t[dim];
1260 for(Int_t i=0;i<dim;i++){
1261 AliITStrackV2* trk = (AliITStrackV2*)tracklist->At(i);
1262 chi2[i]=trk->GetChi2();
1266 Int_t w=0;Double_t value;
1269 for(Int_t j=w+1;j<dim;j++){
1270 if(chi2[w]<chi2[j]){
1282 AliITStrackV2* trk = (AliITStrackV2*)tracklist->At(index[dim-1]);
1289 //__________________________________________________________
1290 Int_t AliITStrackerSA::FindLabel(Int_t l1, Int_t l2, Int_t l3, Int_t l4, Int_t l5, Int_t l6){
1292 //function used to determine the track label
1294 Int_t lb[6] = {l1,l2,l3,l4,l5,l6};
1295 Int_t aa[6]={1,1,1,1,1,1};
1298 Int_t k=0;Int_t w=0;Int_t num=6;
1299 if(lb[5]==-1) num=5;
1303 for(Int_t i=k+1;i<num;i++){
1305 if(lb[k]==lb[i] && aa[k]!=0){
1316 for(Int_t j=0;j<6;j++){
1328 if(num==6) return lb[5];
1332 //_____________________________________________________________________________
1333 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,
1334 Int_t gl12, Int_t gl13, Int_t gl14,Int_t gl15, Int_t gl16, Int_t gl17, Int_t gl18, Int_t numberofpoints){
1337 //function used to assign label to the found track. If track is fake, the label is negative
1339 Int_t lb0[6] = {gl1,gl2,gl3,gl4,gl5,gl6};
1340 Int_t lb1[6] = {gl7,gl8,gl9,gl10,gl11,gl12};
1341 Int_t lb2[6] = {gl13,gl14,gl15,gl16,gl17,gl18};
1342 Int_t ll=FindLabel(lb0[0],lb0[1],lb0[2],lb0[3],lb0[4],lb0[5]);
1343 Int_t lflag=0;Int_t num=6;
1344 if(lb0[5]==-1 && lb1[5]==-1 && lb2[5]==-1) num=5;
1346 for(Int_t i=0;i<num;i++){
1347 if(lb0[i]==ll || lb1[i]==ll || lb2[i]==ll) lflag+=1;
1350 if(lflag>=numberofpoints) return ll;
1356 //_____________________________________________________________________________
1357 void AliITStrackerSA::SetWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
1358 // Set sizes of the phi and lambda windows used for track finding
1360 if(phi){ // user defined values
1361 fPhiWin = new Double_t[fNloop];
1362 fLambdaWin = new Double_t[fNloop];
1363 for(Int_t k=0;k<fNloop;k++){
1365 fLambdaWin[k]=lam[k];
1368 else { // default values
1370 Double_t phid[33] = {0.002,0.003,0.004,0.0045,0.0047,
1371 0.005,0.0053,0.0055,
1372 0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1373 0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
1374 0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1375 Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1377 0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1378 0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1379 0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
1386 fPhiWin = new Double_t[fNloop];
1387 fLambdaWin = new Double_t[fNloop];
1389 for(Int_t k=0;k<fNloop;k++){
1391 fLambdaWin[k]=lambdad[k];