]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITSVertexerTracks.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Removing obsolete macros
[u/mrichter/AliRoot.git] / ITS / AliITSVertexerTracks.cxx
CommitLineData
cab6ff9b 1/**************************************************************************
2 * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-----------------------------------------------------------------
17// Implementation of the vertexer from tracks
18//
19// Origin: A.Dainese, Padova, andrea.dainese@pd.infn.it
20// M.Masera, Torino, massimo.masera@to.infn.it
21//-----------------------------------------------------------------
22
23//---- standard headers ----
24#include <Riostream.h>
25//---- Root headers --------
26#include <TFile.h>
27#include <TTree.h>
28#include <TVector3.h>
29#include <TMatrixD.h>
30#include <TObjArray.h>
31#include <TRandom.h> // TEMPORARY !!!!!!!
32//---- AliRoot headers -----
33#include <AliRun.h>
34#include "AliKalmanTrack.h"
35#include "AliITSStrLine.h"
36#include "AliITStrackV2.h"
37#include "AliITSVertex.h"
38#include "AliITSVertexerTracks.h"
39
40
41ClassImp(AliITSVertexerTracks)
42
43
44//----------------------------------------------------------------------------
45AliITSVertexerTracks::AliITSVertexerTracks():AliITSVertexer() {
46//
47// Default constructor
48//
49
50 SetVtxStart();
51 SetMinTracks();
52 SetUseThrustFrame();
53 SetPhiThrust();
54 fTrksToSkip = 0;
55 fNTrksToSkip = 0;
56 for(Int_t i=0; i<3; i++)fInitPos[i] = 0.;
57}
58
59//----------------------------------------------------------------------------
60AliITSVertexerTracks::AliITSVertexerTracks(TFile *inFile,TFile *outFile,
61 Double_t field,
62 Double_t xStart,Double_t yStart,
63 Int_t useThFr)
64 :AliITSVertexer(inFile,outFile) {
65//
66// Standard constructor
67//
68
69 SetField(field);
70 SetVtxStart(xStart,yStart);
71 SetMinTracks();
72 SetUseThrustFrame(useThFr);
73 SetPhiThrust();
74 fTrksToSkip = 0;
75 fNTrksToSkip = 0;
76 for(Int_t i=0; i<3; i++)fInitPos[i] = 0.;
77}
78//----------------------------------------------------------------------------
79Bool_t AliITSVertexerTracks::CheckField() const {
80//
81// Check if the conv. const. has been set
82//
83 AliITStrackV2 t;
84 Double_t cc = t.GetConvConst();
85 Double_t field = 100./0.299792458/cc;
86
87 if(field<0.1 || field>0.6) {
88 printf("AliITSVertexerTracks::CheckField():\n ERROR: AliKalmanTrack::fConvConst not set\n Use AliKalmanTrack::SetConvConst() or AliITSVertexerTracks::SetField()\n");
89 return kFALSE;
90 }
91 printf("AliITSVertexerTracks::CheckField(): Using B = %3.1f T\n",field);
92 return kTRUE;
93}
94//---------------------------------------------------------------------------
95void AliITSVertexerTracks::ComputeMaxChi2PerTrack(Int_t nTracks) {
96//
97// Max. contr. to the chi2 has been tuned as a function of multiplicity
98//
99 if(nTracks < 7) { fMaxChi2PerTrack = 1.e6;
100 } else { fMaxChi2PerTrack = 100.; }
101
102 return;
103}
104//---------------------------------------------------------------------------
105void AliITSVertexerTracks::FindVertices() {
106//
107// Vertices for all events from fFirstEvent to fLastEvent
108//
109
110 // Check if the conv. const. has been set
111 if(!CheckField()) return;
112
113
114 // loop over events
115 for(Int_t ev=fFirstEvent; ev<=fLastEvent; ev++) {
116 if(ev % 100 == 0 || fDebug) printf("--- Processing event %d of %d ---\n",ev,fLastEvent);
117
118 FindVertexForCurrentEvent(ev);
119
120 if(!fCurrentVertex) {
121 printf("AliITSVertexerTracks::FindVertixes(): no tracks tree for event %d\n",ev);
122 continue;
123 }
124
125 if(fDebug) fCurrentVertex->PrintStatus();
d4221964 126 TString vtxName = "Vertex_";
127 vtxName += ev;
128 fCurrentVertex->SetName(vtxName.Data());
129 fCurrentVertex->SetTitle("VertexerTracks");
cab6ff9b 130 WriteCurrentVertex();
131 } // loop over events
132
133 return;
134}
135//----------------------------------------------------------------------------
136Int_t AliITSVertexerTracks::PrepareTracks(TTree &trkTree) {
137//
138// Propagate tracks to initial vertex position and store them in a TObjArray
139//
140 Double_t maxd0rphi = 3.;
141 Double_t alpha,xlStart,d0rphi;
142 Int_t nTrks = 0;
143 Bool_t skipThis;
144
145 Int_t nEntries = (Int_t)trkTree.GetEntries();
146
147 if(!fTrkArray.IsEmpty()) fTrkArray.Clear();
148 fTrkArray.Expand(nEntries);
149
150 if(fDebug) {
151 printf(" PrepareTracks()\n");
152 trkTree.Print();
153 }
154
155 for(Int_t i=0; i<nEntries; i++) {
156 // check tracks to skip
157 skipThis = kFALSE;
158 for(Int_t j=0; j<fNTrksToSkip; j++) {
159 if(i==fTrksToSkip[j]) {
160 if(fDebug) printf("skipping track: %d\n",i);
161 skipThis = kTRUE;
162 }
163 }
164 if(skipThis) continue;
165
166 AliITStrackV2 *itstrack = new AliITStrackV2;
167 trkTree.SetBranchAddress("tracks",&itstrack);
168 trkTree.GetEvent(i);
169
170
171 // propagate track to vtxSeed
172 alpha = itstrack->GetAlpha();
173 xlStart = fNominalPos[0]*TMath::Cos(alpha)+fNominalPos[1]*TMath::Sin(alpha);
174 itstrack->PropagateTo(3.,0.0023,65.19); // to beam pipe (0.8 mm of Be)
175 itstrack->PropagateTo(xlStart,0.,0.); // to vtxSeed
176
177 // select tracks with d0rphi < maxd0rphi
178 d0rphi = TMath::Abs(itstrack->GetD(fNominalPos[0],fNominalPos[1]));
179 if(d0rphi > maxd0rphi) { delete itstrack; continue; }
180
181 fTrkArray.AddLast(itstrack);
182
183 nTrks++;
184 }
185
186 if(fTrksToSkip) delete [] fTrksToSkip;
187
188 return nTrks;
189}
190//----------------------------------------------------------------------------
191void AliITSVertexerTracks::PrintStatus() const {
192//
193// Print status
194//
195 printf(" Initial position (%f,%f)\n",fNominalPos[0],fNominalPos[1]);
196 printf(" Vertex position after vertex finder (%f, %f, %f)\n",fInitPos[0],fInitPos[1],fInitPos[2]);
197 printf(" Number of tracks in array: %d\n",(Int_t)fTrkArray.GetEntriesFast());
198 printf(" Minimum # tracks required in fit: %d\n",fMinTracks);
199 printf(" Using Thrust Frame: %d fPhiThrust = %f\n",fUseThrustFrame,fPhiThrust);
200
201 return;
202}
203//----------------------------------------------------------------------------
204AliITSVertex* AliITSVertexerTracks::FindVertexForCurrentEvent(Int_t evnumb) {
205//
206// Vertex for current event
207//
208 fCurrentVertex = 0;
209
210 // get tree with tracks from input file
211 TString treeName = "TreeT_ITS_";
212 treeName += evnumb;
213 TTree *trkTree=(TTree*)fInFile->Get(treeName.Data());
214 if(!trkTree) return fCurrentVertex;
215
216
217 // get tracks and propagate them to initial vertex position
218 Int_t nTrks = PrepareTracks(*trkTree);
219 delete trkTree;
220 if(fDebug) printf(" tracks prepared: %d\n",nTrks);
221 if(nTrks < fMinTracks) { TooFewTracks(); return fCurrentVertex; }
222
223 // VERTEX FINDER
224 VertexFinder();
225
226 // VERTEX FITTER
227 ComputeMaxChi2PerTrack(nTrks);
228 if(fUseThrustFrame) ThrustFinderXY();
229 if(fDebug) printf(" thrust found: phi = %f\n",fPhiThrust);
230 VertexFitter();
231 if(fDebug) printf(" vertex fit completed\n");
232
233 TString vtxName;
d4221964 234 vtxName = "Vertex_";
cab6ff9b 235 vtxName += evnumb;
236 fCurrentVertex->SetName(vtxName.Data());
237 return fCurrentVertex;
238}
239//---------------------------------------------------------------------------
240void AliITSVertexerTracks::SetSkipTracks(Int_t n,Int_t *skipped) {
241//
242// Mark the tracks not ot be used in the vertex finding
243//
244 fNTrksToSkip = n;
245 fTrksToSkip = new Int_t[n];
246 for(Int_t i=0;i<n;i++) fTrksToSkip[i] = skipped[i];
247 return;
248}
249//----------------------------------------------------------------------------
250Double_t AliITSVertexerTracks::SumPl(TTree &momTree,Double_t phi) const {
251//
252// Function to be maximized for thrust determination
253//
254 TVector3 u(1.,1.,0);
255 u.SetMag(1.);
256 u.SetPhi(phi);
257 Double_t pl;
258
259 Double_t sum = 0.;
260
261 TVector3* mom=0;
262 momTree.SetBranchAddress("momenta",&mom);
263 Int_t entries = (Int_t)momTree.GetEntries();
264
265 for(Int_t i=0; i<entries; i++) {
266 momTree.GetEvent(i);
267 pl = mom->Dot(u);
268 if(pl>0.) sum += pl;
269 }
270
271 delete mom;
272
273 return sum;
274}
275//---------------------------------------------------------------------------
276void AliITSVertexerTracks::ThrustFinderXY() {
277//
278// This function looks for the thrust direction, \vec{u}, in the (x,y) plane.
279// The following function is maximized:
280// \Sum_{\vec{p}\cdot\vec{u}} \vec{p}\cdot\vec{u} / \Sum |\vec{p}|
281// where \vec{p} = (p_x,p_y)
282//
283 Double_t pt,alpha,phi;
284 Double_t totPt;
285
286 // tree for thrust determination
287 TVector3 *ioMom = new TVector3;
288 TTree *t = new TTree("Tree_Momenta","Tree with momenta");
289 t->Branch("momenta","TVector3",&ioMom);
290 totPt = 0.;
291
292 AliITStrackV2 *itstrack = 0;
293 Int_t arrEntries = (Int_t)fTrkArray.GetEntries();
294
295 // loop on tracks
296 for(Int_t i=0; i<arrEntries; i++) {
297 itstrack = (AliITStrackV2*)fTrkArray.At(i);
298 // momentum of the track at the vertex
299 pt = 1./TMath::Abs(itstrack->Get1Pt());
300 alpha = itstrack->GetAlpha();
301 phi = alpha+TMath::ASin(itstrack->GetSnp());
302 ioMom->SetX(pt*TMath::Cos(phi));
303 ioMom->SetY(pt*TMath::Sin(phi));
304 ioMom->SetZ(0.);
305
306 totPt += ioMom->Pt();
307 t->Fill();
308 } // end loop on tracks
309
310 Double_t tValue=0.,tPhi=0.;
311 Double_t maxSumPl = 0.;
312 Double_t thisSumPl;
313 Double_t dPhi;
314 Int_t nSteps,iStep;
315
316 phi = 0.;
317 nSteps = 100;
318 dPhi = 2.*TMath::Pi()/(Double_t)nSteps;
319
320 for(iStep=0; iStep<nSteps; iStep++) {
321 phi += dPhi;
322 thisSumPl = SumPl(*t,phi);
323 if(thisSumPl > maxSumPl) {
324 maxSumPl = thisSumPl;
325 tPhi = phi;
326 }
327 }
328
329 phi = tPhi-dPhi;
330 nSteps = 10;
331 dPhi /= 5.;
332
333 for(iStep=0; iStep<nSteps; iStep++) {
334 phi += dPhi;
335 thisSumPl = SumPl(*t,phi);
336 if(thisSumPl > maxSumPl) {
337 maxSumPl = thisSumPl;
338 tPhi = phi;
339 }
340 }
341
342 tValue = 2.*maxSumPl/totPt;
343 if(tPhi<0.) tPhi += 2.*TMath::Pi();
344 if(tPhi>2.*TMath::Pi()) tPhi -= 2.*TMath::Pi();
345
346 SetPhiThrust(tPhi);
347
348 delete t;
349 delete ioMom;
350
351 return;
352}
353//---------------------------------------------------------------------------
354void AliITSVertexerTracks::TooFewTracks() {
355//
356// When the number of tracks is < fMinTracks the vertex is set to (0,0,0)
357// and the number of tracks to -1
358//
359 fCurrentVertex = new AliITSVertex(0.,0.,-1);
360 return;
361}
362//---------------------------------------------------------------------------
363void AliITSVertexerTracks::VertexFinder() {
364
365 // Get estimate of vertex position in (x,y) from tracks DCA
366 // Then this estimate is stored to the data member fInitPos
367 // (previous values are overwritten)
368
369
370
371 /*
372******* TEMPORARY!!! FOR TEST ONLY!!! **********************************
373
374fInitPos[0] = fNominalPos[0]+gRandom->Gaus(0.,0.0100); // 100 micron gaussian smearing
375fInitPos[1] = fNominalPos[1]+gRandom->Gaus(0.,0.0100); // 100 micron gaussian smearing
376 */
377
378 fInitPos[2] = 0.;
379 for(Int_t i=0;i<2;i++)fInitPos[i]=fNominalPos[i];
380
381 Int_t nacc = (Int_t)fTrkArray.GetEntriesFast();
382
383 Double_t aver[3]={0.,0.,0.};
384 Int_t ncombi = 0;
385 AliITStrackV2 *track1;
386 AliITStrackV2 *track2;
387 for(Int_t i=0; i<nacc; i++){
388 track1 = (AliITStrackV2*)fTrkArray.At(i);
389 if(fDebug>5){
390 Double_t xv,par[5];
391 track1->GetExternalParameters(xv,par);
392 cout<<"Track in position "<<i<<" xr= "<<xv<<endl;
393 for(Int_t ii=0;ii<5;ii++)cout<<par[ii]<<" ";
394 cout<<endl;
395 }
396 Double_t mom1[3];
397 Double_t alpha = track1->GetAlpha();
398 Double_t azim = TMath::ASin(track1->GetSnp())+alpha;
399 Double_t theta = TMath::Pi()/2. - TMath::ATan(track1->GetTgl());
400 mom1[0] = TMath::Sin(theta)*TMath::Cos(azim);
401 mom1[1] = TMath::Sin(theta)*TMath::Sin(azim);
402 mom1[2] = TMath::Cos(theta);
403
404 Double_t pos1[3];
405 Double_t mindist = TMath::Cos(alpha)*fNominalPos[0]+TMath::Sin(alpha)*fNominalPos[1];
406 track1->GetGlobalXYZat(mindist,pos1[0],pos1[1],pos1[2]);
407 AliITSStrLine *line1 = new AliITSStrLine(pos1,mom1);
408 for(Int_t j=i+1; j<nacc; j++){
409 track2 = (AliITStrackV2*)fTrkArray.At(j);
410 Double_t mom2[3];
411 alpha = track2->GetAlpha();
412 azim = TMath::ASin(track2->GetSnp())+alpha;
413 theta = TMath::Pi()/2. - TMath::ATan(track2->GetTgl());
414 mom2[0] = TMath::Sin(theta)*TMath::Cos(azim);
415 mom2[1] = TMath::Sin(theta)*TMath::Sin(azim);
416 mom2[2] = TMath::Cos(theta);
417 Double_t pos2[3];
418 mindist = TMath::Cos(alpha)*fNominalPos[0]+TMath::Sin(alpha)*fNominalPos[1];
419 track2->GetGlobalXYZat(mindist,pos2[0],pos2[1],pos2[2]);
420 AliITSStrLine *line2 = new AliITSStrLine(pos2,mom2);
421 Double_t crosspoint[3];
422 Int_t retcode = line2->Cross(line1,crosspoint);
423 if(retcode<0){
424 if(fDebug>10)cout<<" i= "<<i<<", j= "<<j<<endl;
425 if(fDebug>10)cout<<"bad intersection\n";
426 line1->PrintStatus();
427 line2->PrintStatus();
428 }
429 else {
430 ncombi++;
431 for(Int_t jj=0;jj<3;jj++)aver[jj]+=crosspoint[jj];
432 if(fDebug>10)cout<<" i= "<<i<<", j= "<<j<<endl;
433 if(fDebug>10)cout<<"\n Cross point: ";
434 if(fDebug>10)cout<<crosspoint[0]<<" "<<crosspoint[1]<<" "<<crosspoint[2]<<endl;
435 }
436 delete line2;
437 }
438 delete line1;
439 }
440 if(ncombi>0){
441 for(Int_t jj=0;jj<3;jj++)fInitPos[jj] = aver[jj]/ncombi;
442 }
443 else {
444 Warning("VertexFinder","Finder did not succed");
445 }
446
447
448 //************************************************************************
449 return;
450
451}
452
453//---------------------------------------------------------------------------
454void AliITSVertexerTracks::VertexFitter() {
455//
456// The optimal estimate of the vertex position is given by a "weighted
457// average of tracks positions"
458// Original method: CMS Note 97/0051
459//
460 if(fDebug) {
461 printf(" VertexFitter(): start\n");
462 PrintStatus();
463 }
464
465
466 Int_t i,j,k,step=0;
467 TMatrixD rv(3,1);
468 TMatrixD V(3,3);
469 rv(0,0) = fInitPos[0];
470 rv(1,0) = fInitPos[1];
471 rv(2,0) = 0.;
472 Double_t xlStart,alpha;
473 Double_t rotAngle;
474 Double_t cosRot,sinRot;
475 Double_t cc[15];
476 Int_t nUsedTrks;
477 Double_t chi2,chi2i;
478 Int_t arrEntries = (Int_t)fTrkArray.GetEntries();
479 AliITStrackV2 *t = 0;
480 Int_t failed = 0;
481
482 Int_t *skipTrack = new Int_t[arrEntries];
483 for(i=0; i<arrEntries; i++) skipTrack[i]=0;
484
485 // 3 steps:
486 // 1st - first estimate of vtx using all tracks
487 // 2nd - apply cut on chi2 max per track
488 // 3rd - estimate of global chi2
489 for(step=0; step<3; step++) {
490 if(fDebug) printf(" step = %d\n",step);
491 chi2 = 0.;
492 nUsedTrks = 0;
493
494 TMatrixD SumWiri(3,1);
495 TMatrixD SumWi(3,3);
496 for(i=0; i<3; i++) {
497 SumWiri(i,0) = 0.;
498 for(j=0; j<3; j++) SumWi(j,i) = 0.;
499 }
500
501 // loop on tracks
502 for(k=0; k<arrEntries; k++) {
503 if(skipTrack[k]) continue;
504 // get track from track array
505 t = (AliITStrackV2*)fTrkArray.At(k);
506 alpha = t->GetAlpha();
507 xlStart = fInitPos[0]*TMath::Cos(alpha)+fInitPos[1]*TMath::Sin(alpha);
508 t->PropagateTo(xlStart,0.,0.); // to vtxSeed
509 rotAngle = alpha-fPhiThrust;
510 if(alpha<0.) rotAngle += 2.*TMath::Pi();
511 cosRot = TMath::Cos(rotAngle);
512 sinRot = TMath::Sin(rotAngle);
513
514 // vector of track global coordinates
515 TMatrixD ri(3,1);
516 ri(0,0) = t->GetX()*cosRot-t->GetY()*sinRot;
517 ri(1,0) = t->GetX()*sinRot+t->GetY()*cosRot;
518 ri(2,0) = t->GetZ();
519
520 // matrix to go from global (x,y,z) to local (y,z);
521 TMatrixD Qi(2,3);
522 Qi(0,0) = -sinRot;
523 Qi(0,1) = cosRot;
524 Qi(0,2) = 0.;
525 Qi(1,0) = 0.;
526 Qi(1,1) = 0.;
527 Qi(1,2) = 1.;
528
529 // covariance matrix of local (y,z) - inverted
530 TMatrixD Ui(2,2);
531 t->GetExternalCovariance(cc);
532 Ui(0,0) = cc[0];
533 Ui(0,1) = cc[1];
534 Ui(1,0) = cc[1];
535 Ui(1,1) = cc[2];
536
537 // weights matrix: Wi = QiT * UiInv * Qi
538 if(Ui.Determinant() <= 0.) continue;
539 TMatrixD UiInv(TMatrixD::kInverted,Ui);
540 TMatrixD UiInvQi(UiInv,TMatrixD::kMult,Qi);
541 TMatrixD Wi(Qi,TMatrixD::kTransposeMult,UiInvQi);
542
543 // track chi2
544 TMatrixD deltar = rv; deltar -= ri;
545 TMatrixD Wideltar(Wi,TMatrixD::kMult,deltar);
546 chi2i = deltar(0,0)*Wideltar(0,0)+
547 deltar(1,0)*Wideltar(1,0)+
548 deltar(2,0)*Wideltar(2,0);
549
550
551 if(step==1 && chi2i > fMaxChi2PerTrack) {
552 skipTrack[k] = 1;
553 continue;
554 }
555
556 // add to total chi2
557 chi2 += chi2i;
558
559 TMatrixD Wiri(Wi,TMatrixD::kMult,ri);
560
561 SumWiri += Wiri;
562 SumWi += Wi;
563
564 nUsedTrks++;
565 } // end loop on tracks
566
567 if(nUsedTrks < fMinTracks) {
568 failed=1;
569 continue;
570 }
571
572 Double_t determinant = SumWi.Determinant();
573 //cerr<<" determinant: "<<determinant<<endl;
574 if(determinant < 100.) {
575 printf("det(V) = 0\n");
576 failed=1;
577 continue;
578 }
579
580 // inverted of weights matrix
581 TMatrixD InvSumWi(TMatrixD::kInverted,SumWi);
582 V = InvSumWi;
583
584 // position of primary vertex
585 rv.Mult(V,SumWiri);
586
587 } // end loop on the 3 steps
588
589 delete [] skipTrack;
590 delete t;
591
592 if(failed) {
593 TooFewTracks();
594 return;
595 }
596
597 Double_t position[3];
598 position[0] = rv(0,0);
599 position[1] = rv(1,0);
600 position[2] = rv(2,0);
601 Double_t covmatrix[6];
602 covmatrix[0] = V(0,0);
603 covmatrix[1] = V(0,1);
604 covmatrix[2] = V(1,1);
605 covmatrix[3] = V(0,2);
606 covmatrix[4] = V(1,2);
607 covmatrix[5] = V(2,2);
608
609 // store data in the vertex object
610 fCurrentVertex = new AliITSVertex(fPhiThrust,position,covmatrix,chi2,nUsedTrks);
611
612 if(fDebug) {
613 printf(" VertexFitter(): finish\n");
614 printf(" rv = ( %f , %f , %f )\n\n",rv(0,0),rv(1,0),rv(2,0));
615 fCurrentVertex->PrintStatus();
616 }
617
618 return;
619}
620//----------------------------------------------------------------------------
621AliITSVertex *AliITSVertexerTracks::VertexOnTheFly(TTree &trkTree) {
622//
623// Return vertex from tracks in trkTree
624//
625 if(fCurrentVertex) fCurrentVertex = 0;
626
627 // get tracks and propagate them to initial vertex position
628 Int_t nTrks = PrepareTracks(*(&trkTree));
629 if(fDebug) printf(" tracks prepared: %d\n",nTrks);
630 if(nTrks < fMinTracks) { TooFewTracks(); return fCurrentVertex; }
631
632 // VERTEX FINDER
633 VertexFinder();
634
635 // VERTEX FITTER
636 ComputeMaxChi2PerTrack(nTrks);
637 if(fUseThrustFrame) ThrustFinderXY();
638 if(fDebug) printf(" thrust found: phi = %f\n",fPhiThrust);
639 VertexFitter();
640 if(fDebug) printf(" vertex fit completed\n");
641
642 return fCurrentVertex;
643}
644//----------------------------------------------------------------------------
645
646
647