1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 ///////////////////////////////////////////////////////////////////////////
20 // Representation and extrapolation of AliTracks in a magnetic field.
22 // This class is meant to provide a means to display and extrapolate
23 // AliTrack objects in the presence of a constant homogeneous magnetic field.
25 // For track/event displays the line width, colour etc... can be set using the
26 // standard facilities (see TAttLine).
27 // By default the linewith is set to 2 and the colour set to -1 in the constructor.
28 // The latter results in an automatic colour coding according to the track charge
29 // with the convention positive=red neutral=green negative=blue.
31 // To indicate the track starting point, the memberfunction SetMarker()
33 // By default no marker will be displayed.
38 // Display and extrapolation of individual tracks
39 // ----------------------------------------------
48 // r1.SetVector(vec,"car");
53 // p.SetVector(vec,"car");
56 // t.SetBeginPoint(r1);
60 // // The magnetic field vector in Tesla
65 // b.SetVector(vec,"car");
67 // AliHelix* helix=new AliHelix();
69 // helix->SetTofmax(1e-7);
71 // TCanvas* c1=new TCanvas("c1","c1");
72 // TView* view=new TView(1);
73 // view->SetRange(-1000,-1000,-1000,1000,1000,1000);
77 // Double_t range[2]={0,600};
78 // helix->Display(&t,range,3);
79 // t.SetCharge(-t.GetCharge());
80 // helix->Display(&t);
82 // // Track extrapolation
83 // Double_t pars[3]={550,0.001,3};
84 // AliPosition* rext=helix->Extrapolate(&t,pars);
85 // if (rext) rext->Data();
86 // ======================================================================
88 // Online display of events generated via AliCollider
89 // --------------------------------------------------
90 // Int_t nevents=5; // Number of events to be generated
91 // Int_t jrun=1; // The run number of this batch of generated events
93 // cout << " ***" << endl;
94 // cout << " *** AliCollider run for " << nevents << " events." << endl;
95 // cout << " ***" << endl;
97 // AliCollider* gen=new AliCollider();
99 // gen->OpenFortranFile(6,"dump.log");
101 // gen->SetVertexMode(2);
102 // gen->SetResolution(1e-4);
104 // gen->SetRunNumber(jrun);
105 // gen->SetPrintFreq(1);
107 // gen->SetSpectatorPmin(0.01);
114 // gen->Init("fixt",zp,ap,zt,at,158);
116 // AliHelix* helix=new AliHelix();
117 // Float_t vec[3]={0,2,0};
119 // b.SetVector(vec,"car");
122 // helix->Refresh(-1); // Refresh display after each event
124 // TCanvas* c1=new TCanvas("c1","c1");
125 // TView* view=new TView(1);
126 // view->SetRange(-200,-200,-200,200,200,200);
129 // // Prepare random number sequence for this run
130 // // to obtain the number of participants for each event
131 // AliRandom rndm(abs(jrun));
132 // Float_t* rans=new Float_t[nevents];
133 // rndm.Uniform(rans,nevents,2,ap+at);
136 // for (Int_t i=0; i<nevents; i++)
139 // gen->MakeEvent(npart);
140 // AliEvent* evt=gen->GetEvent();
143 // helix->Display(evt);
145 // gSystem->Sleep(5000); // Some delay to keep the display on screen
148 // ======================================================================
150 //--- Author: Nick van Eijndhoven 17-jun-2004 Utrecht University
151 //- Modified: NvE $Date$ Utrecht University
152 ///////////////////////////////////////////////////////////////////////////
154 #include "AliHelix.h"
155 #include "Riostream.h"
157 ClassImp(AliHelix) // Class implementation to enable ROOT I/O
159 AliHelix::AliHelix() : THelix()
161 // Default constructor
174 ///////////////////////////////////////////////////////////////////////////
175 AliHelix::~AliHelix()
177 // Destructor to delete dynamically allocated memory.
189 ///////////////////////////////////////////////////////////////////////////
190 AliHelix::AliHelix(const AliHelix& h) : THelix(h)
201 ///////////////////////////////////////////////////////////////////////////
202 void AliHelix::SetB(Ali3Vector& b)
204 // Set the magnetic field vector in Tesla.
210 fB.GetVector(axis,"car");
214 ///////////////////////////////////////////////////////////////////////////
215 Ali3Vector& AliHelix::GetB()
217 // Provide the magnetic field vector in Tesla.
220 ///////////////////////////////////////////////////////////////////////////
221 void AliHelix::SetTofmax(Float_t tof)
223 // Set the maximum time of flight for straight tracks in seconds.
224 // This maximum tof will be used for drawing etc... in case no begin
225 // and endpoints can be determined from the track info.
228 // 1) In case the user specifies an explicit range, it will override
229 // the maximum tof limit.
230 // 2) By default the tofmax is set to 10 ns in the AliHelix constructor.
233 ///////////////////////////////////////////////////////////////////////////
234 Float_t AliHelix::GetTofmax() const
236 // Provide the maximum time of flight for straight tracks in seconds.
239 ///////////////////////////////////////////////////////////////////////////
240 void AliHelix::SetMarker(Int_t style,Float_t size,Int_t col)
242 // Specify the marker (style, size and colour) to indicate the starting point
243 // of a track in a display.
244 // In case col<0 the marker will have the same color as the track itself.
246 // Defaults are style=8, size=0.2 and col=-1.
252 ///////////////////////////////////////////////////////////////////////////
253 void AliHelix::UseEndPoint(Int_t mode)
255 // Select usage of track endpoint in drawing and extrapolation.
256 // This allows correct event displays even for very long tracks.
258 // mode = 0 : Do not use the track endpoint
259 // 1 : Use the track endpoint
261 // The default value is mode=1 (which is also set in the constructor).
263 if (mode==0 || mode==1) fEnduse=mode;
265 ///////////////////////////////////////////////////////////////////////////
266 void AliHelix::MakeCurve(AliTrack* t,Double_t* range,Int_t iaxis,Double_t scale)
268 // Make the helix curve for the specified AliTrack.
269 // Detailed information of all the helix points can be obtained via the
270 // GetN() and GetP() memberfunctions of TPolyLine3D.
271 // In case one wants to display or extrapolate an AliTrack it is preferable
272 // to use the Display() or Extrapolate() memberfunctions.
273 // It is assumed that the track charge is stored in elementary units
274 // (i.e. charge=1 for a proton).
275 // The input argument "scale" specifies the unit scale for the various
276 // locations where scale=0.01 indicates unit scales in cm etc...
277 // In case scale<=0, the unit scale for locations is determined from the
278 // begin, reference or endpoint of the track. If neither of these
279 // positions is present, all locations are assumed to be given in meter.
280 // The lower and upper bounds for the range are specified by range[0] and
281 // range[1] and the argument "iaxis" indicates along which axis this range
283 // The range can be specified either in the LAB frame or in the Helix frame.
284 // The latter is the frame in which the Z axis points in the B direction.
286 // The conventions for the "iaxis" argument are the following :
287 // iaxis = 1 ==> X axis in the LAB frame
288 // 2 ==> Y axis in the LAB frame
289 // 3 ==> Z axis in the LAB frame
290 // -1 ==> X axis in the Helix frame
291 // -2 ==> Y axis in the Helix frame
292 // -3 ==> Z axis in the Helix frame
294 // In case range=0 the begin/end/reference points of the AliTrack and the
295 // maximum time of flight (see the SetTofmax() memberfunction) will be used
296 // and an appropriate choice for the iaxis parameter will be made automatically
297 // based on the track kinematics.
298 // In case the reference point is not present, the begin or endpoint will be used
299 // as reference point for the 3-momentum specification. If neither of these positions
300 // is present, (0,0,0) will be taken as the reference point.
302 // The default values are range=0, iaxis=3 and scale=-1.
304 SetPolyLine(0); // Reset the polyline data points
306 if (!t || (range && !iaxis)) return;
308 Double_t energy=t->GetEnergy(1); // Track energy in GeV
309 Double_t betanorm=t->GetBeta();
311 if (energy<=0 || betanorm<=0) return;
313 AliPosition* rbeg=t->GetBeginPoint();
315 if (fEnduse) rend=t->GetEndPoint();
316 AliPosition* rref=t->GetReferencePoint();
318 // Magnetic field vector or default Z-direction
319 Double_t bvec[3]={0,0,1};
320 if (fB.GetNorm()>0) fB.GetVector(bvec,"car");
322 // The unit scale for locations if not specified by the user
325 scale=1; // Set default to meter
328 scale=rbeg->GetUnitScale();
332 scale=rend->GetUnitScale();
336 scale=rref->GetUnitScale();
340 Double_t c=2.99792458e8/scale; // Lightspeed in the selected unit scale
342 // The helix angular frequency
343 Double_t w=9e7*(t->GetCharge()*fB.GetNorm())/energy;
345 // The particle velocity in the LAB frame
346 Ali3Vector beta=t->GetBetaVector();
349 v.GetVector(vel,"car");
351 // The particle velocity in the Helix frame
352 Ali3Vector betaprim=beta.GetPrimed(fRotMat);
353 v=v.GetPrimed(fRotMat);
355 v.GetVector(velprim,"car");
357 // Check compatibility of velocity and range specification.
361 if (iaxis>0) beta.GetVector(betavec,"car");
362 if (iaxis<0) betaprim.GetVector(betavec,"car");
363 if (fabs(betavec[abs(iaxis)-1])/betanorm<1e-10) return;
366 // The LAB location in which the velocity of the particle is defined
367 Double_t loc[3]={0,0,0};
372 rx=(Ali3Vector)(*rref);
373 scalex=rref->GetUnitScale();
377 rx=(Ali3Vector)(*rbeg);
378 scalex=rbeg->GetUnitScale();
382 rx=(Ali3Vector)(*rend);
383 scalex=rend->GetUnitScale();
386 if (scalex>0 && (scalex/scale>1.1 || scale/scalex>1.1)) rx*=scalex/scale;
387 rx.GetVector(loc,"car");
389 // Initialisation of Helix kinematics
390 SetHelix(loc,vel,w,0,kUnchanged,bvec);
393 if (fabs(w)>0 && fabs(fVt)>0) bend=1;
395 // Flight time boundaries.
396 // The time origin t=0 is chosen to indicate the position in which
397 // the particle velocity was defined.
398 // The total flight time is initialised to the (user specified) tofmax.
399 Double_t tmin=0,tmax=0;
400 Double_t tof=fTofmax;
403 // The trajectory begin and end points
404 Double_t vec1[3]={0,0,0};
405 Double_t vec2[3]={0,0,0};
413 ////////////////////////////////////////
414 // Treatment of straight trajectories //
415 ////////////////////////////////////////
417 if (range) // Specified range allows for exact flight time boundaries
421 tmin=(range[0]-loc[iaxis-1])/vel[iaxis-1];
422 tmax=(range[1]-loc[iaxis-1])/vel[iaxis-1];
429 tmin=(range[0]-loc[abs(iaxis)-1])/velprim[abs(iaxis)-1];
430 tmax=(range[1]-loc[abs(iaxis)-1])/velprim[abs(iaxis)-1];
438 // Make the 'curve' in the LAB frame and exit.
439 // Use the parametrisation : r(t)=r0+t*v
440 // using the range based flight time boundaries.
441 // An additional point in the middle of the trajectory is
442 // generated in view of accuracy in the case of extrapolations.
448 r1.GetVector(vec1,"car");
449 SetNextPoint(float(vec1[0]),float(vec1[1]),float(vec1[2]));
452 r2.GetVector(vec2,"car");
453 SetNextPoint(float(vec2[0]),float(vec2[1]),float(vec2[2]));
456 r2.GetVector(vec2,"car");
457 SetNextPoint(float(vec2[0]),float(vec2[1]),float(vec2[2]));
459 else // Automatic range determination
461 // Initially the point with Z=0 in the Helix frame is taken as a starting point.
462 // In case this point can't be reached, the point in which the particle velocity
463 // was defined is taken as the starting point.
464 // The endpoint is initially obtained by applying the tofmax from the start point.
466 if (fabs(fVz)>0) tmin=-fZ0/fVz;
472 // Override the initial begin and endpoint settings by the track data
475 r1=(Ali3Vector)(*rbeg);
476 scale1=rbeg->GetUnitScale();
477 // All coordinates in the selected unit scale
478 if (scale1/scale>1.1 || scale/scale1>1.1) r1*=scale1/scale;
485 r2=(Ali3Vector)(*rend);
486 scale2=rend->GetUnitScale();
487 // All coordinates in the selected unit scale
488 if (scale2/scale>1.1 || scale/scale2>1.1) r2*=scale2/scale;
491 r1.GetVector(vec1,"car");
492 r2.GetVector(vec2,"car");
494 // Make the 'curve' in the LAB frame and exit.
495 SetNextPoint(float(vec1[0]),float(vec1[1]),float(vec1[2]));
496 SetNextPoint(float(vec2[0]),float(vec2[1]),float(vec2[2]));
501 //////////////////////////////////////
502 // Treatment of curved trajectories //
503 //////////////////////////////////////
505 // Initialisation of the flight time boundaries.
506 // Based on the constant motion of the particle along the Helix Z-axis,
507 // the parametrisation z(t)=z0+fVz*t in the Helix frame is used.
508 // If possible the point with Z=0 in the Helix frame is taken as a starting point.
509 // In case this point can't be reached, the point in which the particle velocity
510 // was defined is taken as the starting point.
512 if (fabs(fVz)>0) tmin=-fZ0/fVz;
522 // Determination of the range in the helix frame
524 if (!range) // Automatic range determination
530 r1=rbeg->GetPrimed(fRotMat);
531 scale1=rbeg->GetUnitScale();
532 // All coordinates in the selected unit scale
533 if (scale1/scale>1.1 || scale/scale1>1.1) r1*=scale1/scale;
534 // Re-calculate the tmin for this new starting point
535 r1.GetVector(vec1,"car");
536 if (fabs(fVz)>0) tmin=(vec1[2]-fZ0)/fVz;
541 r2=rend->GetPrimed(fRotMat);
542 scale2=rend->GetUnitScale();
543 // All coordinates in the selected unit scale
544 if (scale2/scale>1.1 || scale/scale2>1.1) r2*=scale2/scale;
545 r2.GetVector(vec2,"car");
546 if (fabs(fVz)>0) tmax=(vec2[2]-fZ0)/fVz;
548 // Make the curve on basis of the flight time boundaries and exit
555 SetRange(tmin,tmax,kHelixT);
557 else // User explicitly specified range
559 vec1[abs(iaxis)-1]=range[0];
560 vec2[abs(iaxis)-1]=range[1];
561 r1.SetVector(vec1,"car");
562 r2.SetVector(vec2,"car");
563 if (iaxis>0) // Range specified in LAB frame
565 r1=r1.GetPrimed(fRotMat);
566 r1.GetVector(vec1,"car");
567 r2=r2.GetPrimed(fRotMat);
568 r2.GetVector(vec2,"car");
570 // Determination of the axis component with the
571 // largest range difference
575 for (Int_t i=0; i<3; i++)
577 test=fabs(vec1[i]-vec2[i]);
585 Double_t rmin=vec1[imax];
586 Double_t rmax=vec2[imax];
594 // The kinematic range boundaries in the helix frame
595 Double_t xmin=fX0-fVt/fW;
596 Double_t xmax=fX0+fVt/fW;
597 Double_t ymin=fY0-fVt/fW;
598 Double_t ymax=fY0+fVt/fW;
613 // Set the range for the helix
614 if (imax==2 && dmax>0) SetRange(rmin,rmax,kHelixZ);
617 // Limit range to kinematic boundaries if needed
618 if (rmin<=ymin) rmin=ymin+1e-6*dmax;
619 if (rmax>=ymax) rmax=ymax-1e-6*dmax;
620 if (rmin<rmax) SetRange(rmin,rmax,kHelixY);
624 // Limit range to kinematic boundaries if needed
625 if (rmin<=xmin) rmin=xmin+1e-6*dmax;
626 if (rmax>=xmax) rmax=xmax-1e-6*dmax;
627 if (rmin<rmax) SetRange(rmin,rmax,kHelixX);
633 ///////////////////////////////////////////////////////////////////////////
634 void AliHelix::Display(AliTrack* t,Double_t* range,Int_t iaxis,Double_t scale)
636 // Display the helix curve of an AliTrack.
637 // Various curves can be displayed together or individually; please refer to
638 // the memberfunction Refresh() for further details.
639 // It is assumed that the track charge is stored in elementary units
640 // (i.e. charge=1 for a proton).
641 // The input argument "scale" specifies the unit scale for the various
642 // locations where scale=0.01 indicates unit scales in cm etc...
643 // In case scale<=0, the unit scale for locations is determined from the
644 // begin, reference or endpoint of the track. If neither of these
645 // positions is present, all locations are assumed to be given in meter.
646 // The lower and upper bounds for the range are specified by range[0] and
647 // range[1] and the argument "iaxis" indicates along which axis this range
649 // The range can be specified either in the LAB frame or in the Helix frame.
650 // The latter is the frame in which the Z axis points in the B direction.
652 // The conventions for the "iaxis" argument are the following :
653 // iaxis = 1 ==> X axis in the LAB frame
654 // 2 ==> Y axis in the LAB frame
655 // 3 ==> Z axis in the LAB frame
656 // -1 ==> X axis in the Helix frame
657 // -2 ==> Y axis in the Helix frame
658 // -3 ==> Z axis in the Helix frame
660 // In case range=0 the begin/end/reference points of the AliTrack and the
661 // maximum time of flight (see the SetTofmax() memberfunction) will be used
662 // and an appropriate choice for the iaxis parameter will be made automatically
663 // based on the track kinematics.
664 // In case the reference point is not present, the begin or endpoint will be used
665 // as reference point for the 3-momentum specification. If neither of these positions
666 // is present, (0,0,0) will be taken as the reference point.
668 // The default values are range=0, iaxis=3 and scale=-1.
672 // Before any display activity, a TCanvas and a TView have to be initiated
673 // first by the user like for instance
675 // TCanvas* c1=new TCanvas("c1","c1");
676 // TView* view=new TView(1);
677 // view->SetRange(-1000,-1000,-1000,1000,1000,1000);
680 // The user can also use the 3D viewing facilities from the TCanvas menu
681 // to open an appropriate view.
683 if (!t || (range && !iaxis)) return;
685 MakeCurve(t,range,iaxis,scale);
687 if (fRefresh>0) Refresh(fRefresh);
689 Int_t np=GetLastPoint()+1;
692 Float_t* points=GetP();
693 TPolyLine3D* curve=new TPolyLine3D(np,points);
695 curve->SetLineWidth(fLineWidth);
698 Float_t q=t->GetCharge();
699 curve->SetLineColor(kGreen);
700 if (q>0) curve->SetLineColor(kRed);
701 if (q<0) curve->SetLineColor(kBlue);
705 curve->SetLineColor(fLineColor);
711 fCurves=new TObjArray();
716 // Display the marker for the track starting point
719 TPolyMarker3D* m=new TPolyMarker3D();
720 m->SetPoint(0,points[0],points[1],points[2]);
721 m->SetMarkerStyle(fMstyle);
722 m->SetMarkerSize(fMsize);
723 Int_t col=curve->GetLineColor();
724 if (fMcol>0) col=fMcol;
725 m->SetMarkerColor(col);
730 ///////////////////////////////////////////////////////////////////////////
731 void AliHelix::Refresh(Int_t mode)
733 // Refresh the display screen before showing the next curve.
735 // mode = 0 : refreshing fully under user control.
736 // 1 : the display screen will be refreshed automatically
737 // at each individual track display.
738 // -1 : the display screen will be refreshed automatically
739 // at each event display.
741 // The default is mode=0.
743 if (abs(mode)<2) fRefresh=mode;
744 if (fCurves) fCurves->Delete();
746 ///////////////////////////////////////////////////////////////////////////
747 void AliHelix::Display(AliEvent* evt,Double_t* range,Int_t iaxis,Double_t scale)
749 // Display the helix curves of all tracks of the specified event.
750 // Various events can be displayed together or individually; please refer to
751 // the memberfunction Refresh() for further details.
752 // Please refer to the track display memberfunction for further details
753 // on the input arguments.
755 // The default values are range=0, iaxis=3 and scale=-1.
759 // Before any display activity, a TCanvas and a TView have to be initiated
760 // first by the user like for instance
762 // TCanvas* c1=new TCanvas("c1","c1");
763 // TView* view=new TView(1);
764 // view->SetRange(-1000,-1000,-1000,1000,1000,1000);
767 // The user can also use the 3D viewing facilities from the TCanvas menu
768 // to open an appropriate view.
772 if (fRefresh<0) Refresh(fRefresh);
774 Int_t ntk=evt->GetNtracks();
775 for (Int_t jtk=1; jtk<=ntk; jtk++)
777 AliTrack* tx=evt->GetTrack(jtk);
778 if (tx) Display(tx,range,iaxis,scale);
781 ///////////////////////////////////////////////////////////////////////////
782 void AliHelix::Display(TObjArray* arr,Double_t* range,Int_t iaxis,Double_t scale)
784 // Display the helix curves of all tracks in the specified array.
785 // A convenient way to obtain an array with selected tracks from e.g. an AliEvent
786 // is to make use of its GetTracks() selection facility.
787 // Various arrays can be displayed together or individually; please refer to
788 // the memberfunction Refresh() for further details.
789 // Please refer to the track display memberfunction for further details
790 // on the input arguments.
792 // The default values are range=0, iaxis=3 and scale=-1.
796 // Before any display activity, a TCanvas and a TView have to be initiated
797 // first by the user like for instance
799 // TCanvas* c1=new TCanvas("c1","c1");
800 // TView* view=new TView(1);
801 // view->SetRange(-1000,-1000,-1000,1000,1000,1000);
804 // The user can also use the 3D viewing facilities from the TCanvas menu
805 // to open an appropriate view.
809 Int_t ntk=arr->GetEntries();
810 for (Int_t jtk=0; jtk<ntk; jtk++)
812 TObject* obj=arr->At(jtk);
814 if (!(obj->InheritsFrom("AliTrack"))) continue;
815 AliTrack* tx=(AliTrack*)obj;
816 Display(tx,range,iaxis,scale);
819 ///////////////////////////////////////////////////////////////////////////
820 AliPosition* AliHelix::Extrapolate(AliTrack* t,Double_t* pars,Double_t scale)
822 // Extrapolate an AliTrack according to the corresponding helix curve
823 // and provide a pointer to the impact position w.r.t. a specified plane.
824 // In case the track can never reach the specified plane, the returned
825 // position pointer is zero.
826 // Detailed information of all the helix points used in the extrapolation
827 // can be obtained via the GetN() and GetP() memberfunctions of TPolyLine3D.
828 // It is assumed that the track charge is stored in elementary units
829 // (i.e. charge=1 for a proton).
830 // The input argument "scale" specifies the unit scale for the various
831 // locations where scale=0.01 indicates unit scales in cm etc...
832 // In case scale<=0, the unit scale for locations is determined from the
833 // begin, reference or endpoint of the track. If neither of these
834 // positions is present, all locations are assumed to be given in meter.
835 // The extrapolation parameters for the impact plane and required accuracy
836 // are specified by pars[0], pars[1] and pars[2], respectively.
837 // pars[0] = coordinate value of the plane for the impact point
838 // pars[1] = required accuracy on the specified impact plane coordinate
839 // pars[2] = the axis along which the value of par[0] is specified
841 // The parameters can be specified either w.r.t. the LAB frame or the Helix frame.
842 // The latter is the frame in which the Z axis points in the B direction.
844 // The conventions for the par[2] argument are the following :
845 // par[2] = 1 ==> X axis in the LAB frame
846 // 2 ==> Y axis in the LAB frame
847 // 3 ==> Z axis in the LAB frame
848 // -1 ==> X axis in the Helix frame
849 // -2 ==> Y axis in the Helix frame
850 // -3 ==> Z axis in the Helix frame
854 // To obtain an extrapolation to the plane Z=0 in the LAB frame
855 // with an accuracy of 0.001 cm the input arguments would be
856 // pars[0]=0 pars[1]=0.001 pars[2]=3 scale=0.01
858 // Note : The default value for the scale is -1.
866 if (!t || !pars) return fExt;
868 AliPosition* rbeg=t->GetBeginPoint();
869 AliPosition* rend=t->GetEndPoint();
870 AliPosition* rref=t->GetReferencePoint();
872 // The unit scale for locations if not specified by the user
875 scale=1; // Set default to meter
878 scale=rbeg->GetUnitScale();
882 scale=rend->GetUnitScale();
886 scale=rref->GetUnitScale();
891 range[0]=pars[0]-fabs(pars[1])/2.;
892 range[1]=pars[0]+fabs(pars[1])/2.;
894 Int_t iaxis=int(pars[2]);
896 MakeCurve(t,range,iaxis,scale);
898 Int_t np=GetLastPoint()+1;
899 if (!np) return fExt;
901 Float_t* points=GetP();
903 // First point of the curve around the impact
905 Float_t first[3]={points[3*ip],points[3*ip+1],points[3*ip+2]};
907 // Last point of the curve around the impact
909 Float_t last[3]={points[3*ip],points[3*ip+1],points[3*ip+2]};
911 // The accuracy on the impact point
913 err[0]=fabs(first[0]-last[0]);
914 err[1]=fabs(first[1]-last[1]);
915 err[2]=fabs(first[2]-last[2]);
917 // Take the middle point as impact location
919 Float_t imp[3]={points[3*ip],points[3*ip+1],points[3*ip+2]};
921 fExt=new AliPosition();
922 fExt->SetUnitScale(scale);
923 fExt->SetPosition(imp,"car");
924 fExt->SetPositionErrors(err,"car");
928 ///////////////////////////////////////////////////////////////////////////