2 // Main authors: Matevz Tadel & Alja Mrak-Tadel & Bogdan Vulpescu: 2006, 2007
4 /**************************************************************************
5 * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. *
6 * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for *
7 * full copyright notice. *
8 **************************************************************************/
10 #include "AliEveMUONTrack.h"
12 #include <EveBase/AliEveEventManager.h>
16 #include <AliESDMuonTrack.h>
17 #include <AliESDEvent.h>
18 #include <AliESDVertex.h>
19 #include <AliRunLoader.h>
22 #include <AliMUONTrack.h>
23 #include <AliMUONTriggerTrack.h>
24 #include <AliMUONTrackParam.h>
25 #include <AliMUONConstants.h>
26 #include <AliMUONESDInterface.h>
27 #include <AliMUONVCluster.h>
29 #include <TClonesArray.h>
30 #include <TGeoGlobalMagField.h>
31 #include <TParticle.h>
32 #include <TParticlePDG.h>
36 #include <Riostream.h>
38 //==============================================================================
40 //==============================================================================
41 //==============================================================================
43 //______________________________________________________________________________
44 // Produce TEveUtil:TEveTrack from AliMUONTrack with dipole field model
46 ClassImp(AliEveMUONTrack)
49 //______________________________________________________________________________
50 AliEveMUONTrack::AliEveMUONTrack(TEveRecTrack* t, TEveTrackPropagator* rs) :
56 fIsMUONTriggerTrack(kFALSE),
66 //______________________________________________________________________________
67 AliEveMUONTrack::~AliEveMUONTrack()
73 if (fIsRefTrack || fIsESDTrack) delete fTrack;
74 if (fIsMCTrack) delete fPart;
78 //______________________________________________________________________________
79 void AliEveMUONTrack::PrintMCTrackInfo()
82 // information about the MC particle
88 cout << " ! no particle ..." << endl;
93 cout << " MC track parameters at vertex" << endl;
94 cout << " -------------------------------------------------------------------------------------" << endl;
95 cout << " PDG code Vx Vy Vz Px Py Pz " << endl;
98 setw(8) << setprecision(0) <<
99 fPart->GetPdgCode() << " " <<
100 setw(8) << setprecision(3) <<
101 fPart->Vx() << " " <<
102 setw(8) << setprecision(3) <<
103 fPart->Vy() << " " <<
104 setw(8) << setprecision(3) <<
105 fPart->Vz() << " " <<
106 setw(8) << setprecision(3) <<
107 fPart->Px() << " " <<
108 setw(8) << setprecision(3) <<
109 fPart->Py() << " " <<
110 setw(8) << setprecision(4) <<
111 fPart->Pz() << " " <<
115 pt = TMath::Sqrt(fPart->Px()*fPart->Px()+fPart->Py()*fPart->Py());
116 p = TMath::Sqrt(fPart->Px()*fPart->Px()+fPart->Py()*fPart->Py()+fPart->Pz()*fPart->Pz());
120 setw(8) << setprecision(3) <<
121 pt << " GeV/c" << endl;
124 setw(8) << setprecision(4) <<
125 p << " GeV/c" << endl;
129 //______________________________________________________________________________
130 void AliEveMUONTrack::PrintMUONTrackInfo()
133 // information about the reconstructed/reference track; at hits and at vertex
136 Double_t radDeg = 180.0/TMath::Pi();
139 Float_t pt, bc, nbc, zc;
140 AliMUONTrackParam *mtp;
141 TClonesArray *trackParamAtCluster;
144 cout << " ! no reconstructed track ..." << endl;
150 cout << " TEveTrack number " << fLabel << endl;
151 cout << " ---------------------------------------------------------------------------------------------------------------------------------" << endl;
153 cout << " Number of clusters " << fTrack->GetNClusters() << endl;
154 cout << " Match to trigger " << fTrack->GetMatchTrigger() << endl;
155 if (fTrack->GetMatchTrigger()) {
156 cout << " Chi2 tracking-trigger " << fTrack->GetChi2MatchTrigger() << endl;
157 cout << " Local trigger number " << fTrack->GetLoTrgNum() << endl;
163 cout << " TEveTrack reference number " << fLabel << endl;
164 cout << " ---------------------------------------------------------------------------------------------------------------------------------" << endl;
166 cout << " Number of clusters " << fTrack->GetNClusters() << endl;
169 trackParamAtCluster = fTrack->GetTrackParamAtCluster();
170 nparam = trackParamAtCluster->GetEntries();
173 cout << " trackParamAtCluster entries " << nparam << "" << endl;
174 cout << " ---------------------------------------------------------------------------------------------------------------------------------" << endl;
175 cout << " Number InvBendMom BendSlope NonBendSlope BendCoord NonBendCoord Z Px Py Pz P" << endl;
177 for (Int_t i = 0; i < nparam; i++) {
179 mtp = (AliMUONTrackParam*)trackParamAtCluster->At(i);
182 setw(9)<< setprecision(3) <<
185 setw(8) << setprecision(3) <<
186 mtp->GetInverseBendingMomentum() << " " <<
188 setw(8) << setprecision(3) <<
189 mtp->GetBendingSlope()*radDeg << " " <<
191 setw(8) << setprecision(3) <<
192 mtp->GetNonBendingSlope()*radDeg << " " <<
194 setw(8) << setprecision(4) <<
195 mtp->GetBendingCoor() << " " <<
197 setw(8) << setprecision(4) <<
198 mtp->GetNonBendingCoor() << " " <<
200 setw(10) << setprecision(6) <<
201 mtp->GetZ() << " " <<
203 setw(8) << setprecision(4) <<
206 setw(8) << setprecision(4) <<
209 setw(8) << setprecision(4) <<
212 setw(8) << setprecision(4) <<
220 cout << " TEveTrack parameters at vertex" << endl;
221 cout << " --------------------------------------------------------------------------------------------------------------------" << endl;
222 cout << " InvBendMom BendSlope NonBendSlope BendCoord NonBendCoord Z Px Py Pz P" << endl;
224 mtp = (AliMUONTrackParam*)fTrack->GetTrackParamAtVertex();
226 bc = mtp->GetBendingCoor();
227 nbc = mtp->GetNonBendingCoor();
229 if (bc < 0.001) bc = 0.0;
230 if (nbc < 0.001) nbc = 0.0;
231 if (zc < 0.001) zc = 0.0;
234 setw(8) << setprecision(3) <<
235 mtp->GetInverseBendingMomentum() << " " <<
237 setw(8) << setprecision(3) <<
238 mtp->GetBendingSlope()*radDeg << " " <<
240 setw(8) << setprecision(3) <<
241 mtp->GetNonBendingSlope()*radDeg << " " <<
243 setw(8) << setprecision(4) <<
246 setw(8) << setprecision(4) <<
249 setw(10) << setprecision(6) <<
252 setw(8) << setprecision(4) <<
255 setw(8) << setprecision(4) <<
258 setw(8) << setprecision(4) <<
261 setw(8) << setprecision(4) <<
266 pt = TMath::Sqrt(mtp->Px()*mtp->Px()+mtp->Py()*mtp->Py());
270 setw(8) << setprecision(3) <<
271 pt << " GeV/c" << endl;
275 //______________________________________________________________________________
276 void AliEveMUONTrack::PrintMUONTriggerTrackInfo()
279 // information about the trigger track
282 // Double_t radDeg = 180.0/TMath::Pi();
286 //______________________________________________________________________________
287 void AliEveMUONTrack::PrintESDTrackInfo()
290 // information about the reconstructed ESD track at vertex
293 Double_t radDeg = 180.0/TMath::Pi();
296 AliMUONTrackParam *mtp = (AliMUONTrackParam*)fTrack->GetTrackParamAtVertex();
299 cout << " ESD muon track " << endl;
300 cout << " -----------------------------------------------------------------------------------------------------------" << endl;
301 cout << " InvBendMom BendSlope NonBendSlope BendCoord NonBendCoord Z Px Py Pz" << endl;
305 setw(8) << setprecision(4) <<
306 mtp->GetInverseBendingMomentum() << " " <<
308 setw(8) << setprecision(3) <<
309 mtp->GetBendingSlope()*radDeg << " " <<
311 setw(8) << setprecision(3) <<
312 mtp->GetNonBendingSlope()*radDeg << " " <<
314 setw(8) << setprecision(4) <<
315 mtp->GetBendingCoor() << " " <<
317 setw(8) << setprecision(4) <<
318 mtp->GetNonBendingCoor() << " " <<
320 setw(10) << setprecision(6) <<
321 mtp->GetZ() << " " <<
323 setw(8) << setprecision(3) <<
326 setw(8) << setprecision(3) <<
329 setw(8) << setprecision(3) <<
334 pt = TMath::Sqrt(mtp->Px()*mtp->Px()+mtp->Py()*mtp->Py());
338 setw(8) << setprecision(3) <<
339 pt << " GeV/c" << endl;
342 setw(8) << setprecision(4) <<
343 mtp->P() << " GeV/c" << endl;
345 AliESDEvent* esd = AliEveEventManager::AssertESD();
347 Double_t spdVertexX = 0;
348 Double_t spdVertexY = 0;
349 Double_t spdVertexZ = 0;
350 Double_t esdVertexX = 0;
351 Double_t esdVertexY = 0;
352 Double_t esdVertexZ = 0;
354 AliESDVertex* spdVertex = (AliESDVertex*) esd->GetVertex();
355 if (spdVertex->GetNContributors()) {
356 spdVertexZ = spdVertex->GetZv();
357 spdVertexY = spdVertex->GetYv();
358 spdVertexX = spdVertex->GetXv();
361 AliESDVertex* esdVertex = (AliESDVertex*) esd->GetPrimaryVertex();
362 if (esdVertex->GetNContributors()) {
363 esdVertexZ = esdVertex->GetZv();
364 esdVertexY = esdVertex->GetYv();
365 esdVertexX = esdVertex->GetXv();
368 Float_t t0v = esd->GetT0zVertex();
372 cout << "External vertex SPD: " <<
374 spdVertex->GetNContributors() << " " <<
375 setw(8) << setprecision(3) <<
378 spdVertexZ << " " << endl;
379 cout << "External vertex ESD: " <<
381 esdVertex->GetNContributors() << " " <<
382 setw(8) << setprecision(3) <<
385 esdVertexZ << " " << endl;
386 cout << "External vertex T0: " <<
387 setw(8) << setprecision(3) <<
392 //______________________________________________________________________________
393 void AliEveMUONTrack::MUONTrackInfo()
403 if (fIsMUONTrack || fIsRefTrack) {
404 PrintMUONTrackInfo();
411 if (fIsMUONTriggerTrack) {
412 PrintMUONTriggerTrackInfo();
418 cout << " (slopes [deg], coord [cm], p [GeV/c])" << endl;
422 //______________________________________________________________________________
423 void AliEveMUONTrack::MUONTriggerInfo()
430 TEveUtil::TEveUtil::LoadMacro("MUON_trigger_info.C");
431 gROOT->ProcessLine(Form("MUON_trigger_info(%d);", fLabel));
434 cout << "This is a reference track!" << endl;
437 cout << "This is a Monte-Carlo track!" << endl;
441 AliESDEvent* esd = AliEveEventManager::AssertESD();
442 ULong64_t triggerMask = esd->GetTriggerMask();
445 cout << ">>>>>#########################################################################################################################" << endl;
448 cout << " ESD track trigger info" << endl;
449 cout << " -----------------------------------------------------" << endl;
452 cout << " Match to trigger " << fTrack->GetMatchTrigger() << endl;
454 cout << " ESD trigger mask = " << triggerMask << endl;
457 cout << "#########################################################################################################################<<<<<" << endl;
464 //______________________________________________________________________________
465 void AliEveMUONTrack::MakeMUONTrack(AliMUONTrack *mtrack)
468 // builds the track with dipole field
473 fIsMUONTrack = kTRUE;
476 fTrack = new AliMUONTrack(*mtrack);
481 Float_t ax, bx, ay, by;
482 Float_t xr[28], yr[28], zr[28];
483 Float_t xrc[28], yrc[28], zrc[28];
487 TMatrixD smatrix(2,2);
491 // middle z between the two detector planes of the trigger chambers
492 Float_t zg[4] = { -1603.5, -1620.5, -1703.5, -1720.5 };
495 Float_t pv[3] = {0., 0., 0.};
498 if (mtrack->GetMatchTrigger()) {
499 snprintf(form,1000,"MUONTrack %2d (MT)", fLabel);
501 snprintf(form,1000,"MUONTrack %2d ", fLabel);
507 AliMUONTrackParam *trackParam = 0x0;
508 if (fIsMUONTrack || fIsESDTrack) {
509 trackParam = mtrack->GetTrackParamAtVertex();
510 SetPoint(fCount,trackParam->GetNonBendingCoor(),trackParam->GetBendingCoor(),trackParam->GetZ());
514 for (Int_t i = 0; i < 28; i++) {
515 xr[i]=yr[i]=zr[i]=0.0;
519 Int_t nTrackHits = mtrack->GetNClusters();
520 TClonesArray* trackParamAtCluster = mtrack->GetTrackParamAtCluster();
521 for (Int_t iHit = 0; iHit < nTrackHits; iHit++){
523 trackParam = (AliMUONTrackParam*) trackParamAtCluster->At(iHit);
526 if (IsMUONTrack() || IsESDTrack()) {
527 pt = TMath::Sqrt(trackParam->Px()*trackParam->Px()+trackParam->Py()*trackParam->Py());
528 SetLineColor(ColorIndex(pt));
530 pv[0] = trackParam->Px();
531 pv[1] = trackParam->Py();
532 pv[2] = trackParam->Pz();
536 xr[iHit] = trackParam->GetNonBendingCoor();
537 yr[iHit] = trackParam->GetBendingCoor();
538 zr[iHit] = trackParam->GetZ();
539 chr[iHit] = trackParam->GetClusterPtr()->GetChamberId();
543 SetPoint(fCount,xr[0],yr[0],zr[0]);
545 for (Int_t iHit = 1; iHit < nTrackHits; iHit++) {
546 if (chr[iHit] > 3 && chr[iHit-1] < 6) Propagate(xr,yr,zr,iHit-1,iHit);
547 SetPoint(fCount,xr[iHit],yr[iHit],zr[iHit]);
551 if (!fIsMUONTrack && !fIsESDTrack) return;
554 if (mtrack->GetMatchTrigger()) {
556 for (Int_t i = 0; i < nTrackHits; i++) {
557 if (TMath::Abs(zr[i]) > 1000.0) {
558 //printf("TEveHit %d x %f y %f z %f \n",iHit,xr[i],yr[i],zr[i]);
571 for (Int_t i = 0; i < nrc; i++) {
572 xv = (Double_t)zrc[i];
573 yv = (Double_t)xrc[i];
574 //printf("x-z: xv %f yv %f \n",xv,yv);
576 smatrix(1,1) += xv*xv;
582 res = smatrix.Invert() * sums;
589 for (Int_t i = 0; i < nrc; i++) {
590 xv = (Double_t)zrc[i];
591 yv = (Double_t)yrc[i];
592 //printf("y-z: xv %f yv %f \n",xv,yv);
594 smatrix(1,1) += xv*xv;
600 res = smatrix.Invert() * sums;
604 Float_t xtc, ytc, ztc;
605 for (Int_t ii = 0; ii < 4; ii++) {
611 //printf("tc: x %f y %f z %f \n",xtc,ytc,ztc);
613 SetPoint(fCount,xtc,ytc,ztc);
618 } // end match trigger
622 //______________________________________________________________________________
623 void AliEveMUONTrack::MakeMUONTriggerTrack(AliMUONTriggerTrack *mtrack)
626 // builds the trigger track from one point and direction
629 Float_t x1 = mtrack->GetX11();
630 Float_t y1 = mtrack->GetY11();
631 Float_t thex = mtrack->GetThetax();
632 Float_t they = mtrack->GetThetay();
634 Float_t z11 = -1600.0;
635 Float_t z22 = -1724.0;
636 Float_t dz = z22-z11;
638 Float_t x2 = x1 + dz*TMath::Tan(thex);
639 Float_t y2 = y1 + dz*TMath::Tan(they);
641 SetPoint(fCount,x1,y1,z11); fCount++;
642 SetPoint(fCount,x2,y2,z22); fCount++;
646 snprintf(form,1000,"MUONTriggerTrack %2d",mtrack->GetLoTrgNum());
652 //______________________________________________________________________________
653 void AliEveMUONTrack::MakeESDTrack(AliESDMuonTrack *mtrack)
656 // builds the track with dipole field starting from the TParticle
662 if (mtrack->GetMatchTrigger()) {
663 snprintf(form,1000,"ESDTrack %2d (MT)", fLabel);
665 snprintf(form,1000,"ESDTrack %2d ", fLabel);
671 fTrack = new AliMUONTrack();
673 // create a simple track from the ESD track
674 AliMUONESDInterface::ESDToMUON(*mtrack,*fTrack);
676 // reset track parameters at vertex to the ones at DCA
677 AliMUONTrackParam paramAtDCA;
678 AliMUONESDInterface::GetParamAtDCA(*mtrack, paramAtDCA);
679 fTrack->SetTrackParamAtVertex(¶mAtDCA);
681 MakeMUONTrack(fTrack);
685 //______________________________________________________________________________
686 void AliEveMUONTrack::MakeMCTrack(TParticle *part)
689 // builds the track with dipole field starting from the TParticle
694 fPart = new TParticle(*part);
697 snprintf(form,1000,"TEveMCTrack %2d ", fLabel);
702 Double_t vect[7], vout[7];
711 vect[0] = fPart->Vx();
712 vect[1] = fPart->Vy();
713 vect[2] = fPart->Vz();
714 vect[3] = fPart->Px()/fPart->P();
715 vect[4] = fPart->Py()/fPart->P();
716 vect[5] = fPart->Pz()/fPart->P();
717 vect[6] = fPart->P();
719 TParticlePDG *ppdg = fPart->GetPDG(1);
720 Int_t charge = (Int_t)(ppdg->Charge()/3.0);
722 Double_t zMax = -1750.0;
723 Double_t rMax = 350.0;
727 while ((vect[2] > zMax) && (nSteps < 10000) && (r < rMax)) {
729 OneStepRungekutta(charge, step, vect, vout);
730 SetPoint(fCount,vout[0],vout[1],vout[2]);
732 for (Int_t i = 0; i < 7; i++) {
735 r = TMath::Sqrt(vect[0]*vect[0]+vect[1]*vect[1]);
740 //______________________________________________________________________________
741 void AliEveMUONTrack::MakeRefTrack(AliMUONTrack *mtrack)
744 // builds the track with dipole field starting from the TParticle
750 snprintf(form,1000,"RefTrack %2d ", fLabel);
755 MakeMUONTrack(mtrack);
759 //______________________________________________________________________________
760 void AliEveMUONTrack::Propagate(Float_t *xr, Float_t *yr, Float_t *zr, Int_t i1, Int_t i2)
763 // propagate in magnetic field between hits of indices i1 and i2
766 Double_t vect[7], vout[7];
770 AliMUONTrackParam *trackParam = 0;
771 TClonesArray *trackParamAtCluster = 0;
774 zMax = zr[i1]+1.5*step;
776 zMax = zr[i2]+1.5*step;
779 trackParamAtCluster = fTrack->GetTrackParamAtCluster();
781 if (IsMUONTrack() || IsESDTrack() || IsRefTrack()) {
782 trackParam = (AliMUONTrackParam*)trackParamAtCluster->At(i1);
783 charge = (Int_t)trackParam->GetCharge();
791 vect[6] = trackParam->P();
792 vect[3] = trackParam->Px()/vect[6];
793 vect[4] = trackParam->Py()/vect[6];
794 vect[5] = trackParam->Pz()/vect[6];
797 while ((vect[2] > zMax) && (nSteps < 10000)) {
799 OneStepRungekutta(charge, step, vect, vout);
800 SetPoint(fCount,vout[0],vout[1],vout[2]);
802 for (Int_t i = 0; i < 7; i++) {
809 //______________________________________________________________________________
810 void AliEveMUONTrack::OneStepRungekutta(Double_t charge, Double_t step,
811 Double_t* vect, Double_t* vout)
813 /// ******************************************************************
815 /// * Runge-Kutta method for tracking a particle through a magnetic *
816 /// * field. Uses Nystroem algorithm (See Handbook Nat. Bur. of *
817 /// * Standards, procedure 25.5.20) *
819 /// * Input parameters *
820 /// * CHARGE Particle charge *
821 /// * STEP Step size *
822 /// * VECT Initial co-ords,direction cosines,momentum *
823 /// * Output parameters *
824 /// * VOUT Output co-ords,direction cosines,momentum *
825 /// * User routine called *
826 /// * CALL GUFLD(X,F) *
828 /// * ==>Called by : <USER>, GUSWIM *
829 /// * Authors R.Brun, M.Hansroul ********* *
830 /// * V.Perevoztchikov (CUT STEP implementation) *
833 /// ******************************************************************
835 Double_t h2, h4, f[4];
836 Double_t xyzt[3], a, b, c, ph,ph2;
837 Double_t secxs[4],secys[4],seczs[4],hxp[3];
838 Double_t g1, g2, g3, g4, g5, g6, ang2, dxt, dyt, dzt;
839 Double_t est, at, bt, ct, cba;
840 Double_t f1, f2, f3, f4, rho, tet, hnorm, hp, rho1, sint, cost;
850 Double_t maxit = 1992;
851 Double_t maxcut = 11;
853 const Double_t kdlt = 1e-4;
854 const Double_t kdlt32 = kdlt/32.;
855 const Double_t kthird = 1./3.;
856 const Double_t khalf = 0.5;
857 const Double_t kec = 2.9979251e-4;
859 const Double_t kpisqua = 9.86960440109;
863 const Int_t kipx = 3;
864 const Int_t kipy = 4;
865 const Int_t kipz = 5;
868 // *. ------------------------------------------------------------------
870 // * this constant is for units cm,gev/c and kgauss
874 for(Int_t j = 0; j < 7; j++)
877 Double_t pinv = kec * charge / vect[6];
885 if (TMath::Abs(h) > TMath::Abs(rest)) h = rest;
886 //cmodif: call gufld(vout,f) changed into:
887 TGeoGlobalMagField::Instance()->Field(vout,f);
890 // * start of integration
903 secxs[0] = (b * f[2] - c * f[1]) * ph2;
904 secys[0] = (c * f[0] - a * f[2]) * ph2;
905 seczs[0] = (a * f[1] - b * f[0]) * ph2;
906 ang2 = (secxs[0]*secxs[0] + secys[0]*secys[0] + seczs[0]*seczs[0]);
907 if (ang2 > kpisqua) break;
909 dxt = h2 * a + h4 * secxs[0];
910 dyt = h2 * b + h4 * secys[0];
911 dzt = h2 * c + h4 * seczs[0];
916 // * second intermediate point
919 est = TMath::Abs(dxt) + TMath::Abs(dyt) + TMath::Abs(dzt);
921 if (ncut++ > maxcut) break;
930 //cmodif: call gufld(xyzt,f) changed into:
931 TGeoGlobalMagField::Instance()->Field(xyzt,f);
937 secxs[1] = (bt * f[2] - ct * f[1]) * ph2;
938 secys[1] = (ct * f[0] - at * f[2]) * ph2;
939 seczs[1] = (at * f[1] - bt * f[0]) * ph2;
943 secxs[2] = (bt * f[2] - ct * f[1]) * ph2;
944 secys[2] = (ct * f[0] - at * f[2]) * ph2;
945 seczs[2] = (at * f[1] - bt * f[0]) * ph2;
946 dxt = h * (a + secxs[2]);
947 dyt = h * (b + secys[2]);
948 dzt = h * (c + seczs[2]);
952 at = a + 2.*secxs[2];
953 bt = b + 2.*secys[2];
954 ct = c + 2.*seczs[2];
956 est = TMath::Abs(dxt)+TMath::Abs(dyt)+TMath::Abs(dzt);
957 if (est > 2.*TMath::Abs(h)) {
958 if (ncut++ > maxcut) break;
967 //cmodif: call gufld(xyzt,f) changed into:
968 TGeoGlobalMagField::Instance()->Field(xyzt,f);
970 z = z + (c + (seczs[0] + seczs[1] + seczs[2]) * kthird) * h;
971 y = y + (b + (secys[0] + secys[1] + secys[2]) * kthird) * h;
972 x = x + (a + (secxs[0] + secxs[1] + secxs[2]) * kthird) * h;
974 secxs[3] = (bt*f[2] - ct*f[1])* ph2;
975 secys[3] = (ct*f[0] - at*f[2])* ph2;
976 seczs[3] = (at*f[1] - bt*f[0])* ph2;
977 a = a+(secxs[0]+secxs[3]+2. * (secxs[1]+secxs[2])) * kthird;
978 b = b+(secys[0]+secys[3]+2. * (secys[1]+secys[2])) * kthird;
979 c = c+(seczs[0]+seczs[3]+2. * (seczs[1]+seczs[2])) * kthird;
981 est = TMath::Abs(secxs[0]+secxs[3] - (secxs[1]+secxs[2]))
982 + TMath::Abs(secys[0]+secys[3] - (secys[1]+secys[2]))
983 + TMath::Abs(seczs[0]+seczs[3] - (seczs[1]+seczs[2]));
985 if (est > kdlt && TMath::Abs(h) > 1.e-4) {
986 if (ncut++ > maxcut) break;
992 // * if too many iterations, go to helix
993 if (iter++ > maxit) break;
998 cba = 1./ TMath::Sqrt(a*a + b*b + c*c);
1006 if (step < 0.) rest = -rest;
1007 if (rest < 1.e-5*TMath::Abs(step)) return;
1011 // angle too big, use helix
1016 f4 = TMath::Sqrt(f1*f1+f2*f2+f3*f3);
1025 hxp[0] = f2*vect[kipz] - f3*vect[kipy];
1026 hxp[1] = f3*vect[kipx] - f1*vect[kipz];
1027 hxp[2] = f1*vect[kipy] - f2*vect[kipx];
1029 hp = f1*vect[kipx] + f2*vect[kipy] + f3*vect[kipz];
1032 sint = TMath::Sin(tet);
1033 cost = 2.*TMath::Sin(khalf*tet)*TMath::Sin(khalf*tet);
1037 g3 = (tet-sint) * hp*rho1;
1042 vout[kix] = vect[kix] + g1*vect[kipx] + g2*hxp[0] + g3*f1;
1043 vout[kiy] = vect[kiy] + g1*vect[kipy] + g2*hxp[1] + g3*f2;
1044 vout[kiz] = vect[kiz] + g1*vect[kipz] + g2*hxp[2] + g3*f3;
1046 vout[kipx] = vect[kipx] + g4*vect[kipx] + g5*hxp[0] + g6*f1;
1047 vout[kipy] = vect[kipy] + g4*vect[kipy] + g5*hxp[1] + g6*f2;
1048 vout[kipz] = vect[kipz] + g4*vect[kipz] + g5*hxp[2] + g6*f3;
1053 //______________________________________________________________________________
1054 Int_t AliEveMUONTrack::ColorIndex(Float_t val)
1057 // returns color index in the palette for a give value
1060 Float_t threshold = 0.0;
1061 Float_t maxVal = 10.0;
1063 Float_t div = TMath::Max(1, (Int_t)(maxVal - threshold));
1064 Int_t nCol = gStyle->GetNumberOfColors();
1065 Int_t cBin = (Int_t) TMath::Nint(nCol*(val - threshold)/div);
1067 return gStyle->GetColorPalette(TMath::Min(nCol - 1, cBin));
1071 //==============================================================================
1072 // Temporary AliEveMUONTrackList
1073 //==============================================================================
1075 //______________________________________________________________________________
1076 void AliEveMUONTrackList::HackMomentumLimits(Bool_t recurse)
1078 // Find momentum limits from included tracks.
1082 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
1084 TEveTrack* track = dynamic_cast<TEveTrack*>(*i);
1087 fLimPt = TMath::Max(fLimPt, track->GetMomentum().Perp());
1088 fLimP = TMath::Max(fLimP, track->GetMomentum().Mag());
1091 FindMomentumLimits(*i, recurse);
1094 fLimPt = RoundMomentumLimit(fLimPt);
1095 fLimP = RoundMomentumLimit(fLimP);
1096 if (fMaxPt == 0) fMaxPt = fLimPt;
1097 if (fMaxP == 0) fMaxP = fLimP;