BOSS 7.0.5
BESIII Offline Software System
Loading...
Searching...
No Matches
Public Member Functions | Friends | List of all members
TRunge Class Reference

A class to represent a track in tracking. More...

#include <TRunge.h>

+ Inheritance diagram for TRunge:

Public Member Functions

 TRunge ()
 Constructors.
 
 TRunge (const TTrack &)
 
 TRunge (const Helix &)
 
 TRunge (const TRunge &)
 
 TRunge (const RecMdcTrack &)
 
 ~TRunge ()
 Destructor.
 
unsigned objectType (void) const
 returns object type
 
double dr (void) const
 Track parameters (at pivot)
 
double phi0 (void) const
 
double kappa (void) const
 
double dz (void) const
 
double tanl (void) const
 
const HepPoint3D & pivot (void) const
 pivot position
 
const Vector & a (void) const
 returns helix parameter
 
const SymMatrix & Ea (void) const
 returns error matrix
 
Helix helix (void) const
 returns helix class
 
unsigned ndf (void) const
 returns NDF
 
double chi2 (void) const
 returns chi2.
 
double reducedchi2 (void) const
 returns reduced-chi2
 
int BfieldID (void) const
 returns B field ID
 
double StepSize (void) const
 returns step size
 
const double * Yscal (void) const
 return error parameters for fitting with step size control
 
double Eps (void) const
 
double StepSizeMax (void) const
 
double StepSizeMin (void) const
 
float Mass (void) const
 return mass
 
double MaxFlightLength (void) const
 return max flight length
 
int approach (TMLink &, const RkFitMaterial, bool sagCorrection=true)
 
int approach (TMLink &, float &tof, HepVector3D &p, const RkFitMaterial, bool sagCorrection=true)
 
int approach_line (TMLink &, const HepPoint3D &, const HepVector3D &, HepPoint3D &onLine, HepPoint3D &onTrack, const RkFitMaterial material)
 caluculate closest points between a line and this track
 
int approach_line (TMLink &, const HepPoint3D &, const HepVector3D &, HepPoint3D &onLine, HepPoint3D &onTrack, float &tof, HepVector3D &p, const RkFitMaterial material)
 
int approach_line (TMLink &, const HepPoint3D &, const HepVector3D &, HepPoint3D &onLine, HepPoint3D &onTrack, float &tof, HepVector3D &p, const RkFitMaterial material, unsigned &stepNum)
 
int approach_point (TMLink &, const HepPoint3D &, HepPoint3D &onTrack, const RkFitMaterial material) const
 caluculate closest point between a point and this track
 
const HepPoint3D & pivot (const HepPoint3D &)
 set new pivot
 
const Vector & a (const Vector &)
 set helix parameter
 
const SymMatrix & Ea (const SymMatrix &)
 set helix error matrix
 
int BfieldID (int)
 set B field map ID
 
double StepSize (double)
 set step size to calc. trajectory
 
const double * Yscal (const double *)
 set error parameters for fitting with step size control
 
double Eps (double)
 
double StepSizeMax (double)
 
double StepSizeMin (double)
 
float Mass (float)
 set mass for tof calc.
 
double MaxFlightLength (double)
 
double DisToWire (TMLink &, double, HepPoint3D, HepPoint3D &)
 
double dEpath (double, double, double) const
 
void eloss (double path, const RkFitMaterial *material, double mass, double y[6], int index) const
 
unsigned Fly (const RkFitMaterial material) const
 make the trajectory in cache, return the number of step
 
unsigned Fly_SC (void) const
 
double intersect_cylinder (double r) const
 
double intersect_yz_plane (const HepTransform3D &plane, double x) const
 
double intersect_zx_plane (const HepTransform3D &plane, double y) const
 
double intersect_xy_plane (double z) const
 
void Propagate (double y[6], const double &step, const RkFitMaterial material) const
 propagate the track using 4th-order Runge-Kutta method
 
void Function (const double y[6], double f[6]) const
 
void Propagate1 (const double y[6], const double dydx[6], const double &step, double yout[6]) const
 
void Propagate_QC (double y[6], double dydx[6], const double &steptry, const double &eps, const double yscal[6], double &stepdid, double &stepnext) const
 
void SetFirst (double y[6]) const
 set first point (position, momentum) s=0, phi=0
 
unsigned Nstep (void) const
 access to trajectory cache
 
int GetXP (unsigned stepNum, double y[6]) const
 
int GetStep (unsigned stepNum, double &step) const
 
double SetFlightLength (void)
 set flight length using wire hits
 
 TRunge ()
 Constructors.
 
 TRunge (const TTrack &)
 
 TRunge (const Helix &)
 
 TRunge (const TRunge &)
 
 TRunge (const RecMdcTrack &)
 
 ~TRunge ()
 Destructor.
 
unsigned objectType (void) const
 returns object type
 
double dr (void) const
 Track parameters (at pivot)
 
double phi0 (void) const
 
double kappa (void) const
 
double dz (void) const
 
double tanl (void) const
 
const HepPoint3D & pivot (void) const
 pivot position
 
const Vector & a (void) const
 returns helix parameter
 
const SymMatrix & Ea (void) const
 returns error matrix
 
Helix helix (void) const
 returns helix class
 
unsigned ndf (void) const
 returns NDF
 
double chi2 (void) const
 returns chi2.
 
double reducedchi2 (void) const
 returns reduced-chi2
 
int BfieldID (void) const
 returns B field ID
 
double StepSize (void) const
 returns step size
 
const double * Yscal (void) const
 return error parameters for fitting with step size control
 
double Eps (void) const
 
double StepSizeMax (void) const
 
double StepSizeMin (void) const
 
float Mass (void) const
 return mass
 
double MaxFlightLength (void) const
 return max flight length
 
int approach (TMLink &, const RkFitMaterial, bool sagCorrection=true)
 
int approach (TMLink &, float &tof, HepVector3D &p, const RkFitMaterial, bool sagCorrection=true)
 
int approach_line (TMLink &, const HepPoint3D &, const HepVector3D &, HepPoint3D &onLine, HepPoint3D &onTrack, const RkFitMaterial material)
 caluculate closest points between a line and this track
 
int approach_line (TMLink &, const HepPoint3D &, const HepVector3D &, HepPoint3D &onLine, HepPoint3D &onTrack, float &tof, HepVector3D &p, const RkFitMaterial material)
 
int approach_line (TMLink &, const HepPoint3D &, const HepVector3D &, HepPoint3D &onLine, HepPoint3D &onTrack, float &tof, HepVector3D &p, const RkFitMaterial material, unsigned &stepNum)
 
int approach_point (TMLink &, const HepPoint3D &, HepPoint3D &onTrack, const RkFitMaterial material) const
 caluculate closest point between a point and this track
 
const HepPoint3D & pivot (const HepPoint3D &)
 set new pivot
 
const Vector & a (const Vector &)
 set helix parameter
 
const SymMatrix & Ea (const SymMatrix &)
 set helix error matrix
 
int BfieldID (int)
 set B field map ID
 
double StepSize (double)
 set step size to calc. trajectory
 
const double * Yscal (const double *)
 set error parameters for fitting with step size control
 
double Eps (double)
 
double StepSizeMax (double)
 
double StepSizeMin (double)
 
float Mass (float)
 set mass for tof calc.
 
double MaxFlightLength (double)
 
double DisToWire (TMLink &, double, HepPoint3D, HepPoint3D &)
 
double dEpath (double, double, double) const
 
void eloss (double path, const RkFitMaterial *material, double mass, double y[6], int index) const
 
unsigned Fly (const RkFitMaterial material) const
 make the trajectory in cache, return the number of step
 
unsigned Fly_SC (void) const
 
double intersect_cylinder (double r) const
 
double intersect_yz_plane (const HepTransform3D &plane, double x) const
 
double intersect_zx_plane (const HepTransform3D &plane, double y) const
 
double intersect_xy_plane (double z) const
 
void Propagate (double y[6], const double &step, const RkFitMaterial material) const
 propagate the track using 4th-order Runge-Kutta method
 
void Function (const double y[6], double f[6]) const
 
void Propagate1 (const double y[6], const double dydx[6], const double &step, double yout[6]) const
 
void Propagate_QC (double y[6], double dydx[6], const double &steptry, const double &eps, const double yscal[6], double &stepdid, double &stepnext) const
 
void SetFirst (double y[6]) const
 set first point (position, momentum) s=0, phi=0
 
unsigned Nstep (void) const
 access to trajectory cache
 
int GetXP (unsigned stepNum, double y[6]) const
 
int GetStep (unsigned stepNum, double &step) const
 
double SetFlightLength (void)
 set flight length using wire hits
 
- Public Member Functions inherited from TTrackBase
 TTrackBase ()
 Constructor.
 
 TTrackBase (const AList< TMLink > &links)
 Constructor.
 
virtual ~TTrackBase ()
 Destructor.
 
virtual unsigned objectType (void) const
 returns object type.
 
virtual unsigned type (void) const
 returns type. Definition is depending on an object class.
 
virtual void dump (const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
 dumps debug information.
 
const AList< TMLink > & links (unsigned mask=0) const
 returns a list of masked TMLinks assigned to this track. 'mask' will be applied if mask is not 0.
 
unsigned nLinks (unsigned mask=0) const
 returns # of masked TMLinks assigned to this track object.
 
const AList< TMLink > & cores (unsigned mask=0) const
 returns a list of masked TMLinks for fit. 'mask' will be applied if mask is not 0.
 
unsigned nCores (unsigned mask=0) const
 returns # of masked TMLinks for fit. 'mask' will be applied if mask is not 0.
 
void update (void) const
 update cache.
 
void append (TMLink &)
 appends a TMLink.
 
void append (const AList< TMLink > &)
 appends TMLinks.
 
void appendByApproach (AList< TMLink > &list, double maxSigma)
 appends TMLinks by approach. 'list' is an input. Unappended TMLinks will be removed from 'list' when returned.
 
void appendByDistance (AList< TMLink > &list, double maxDistance)
 appends TMLinks by distance. 'list' is an input. Unappended TMLinks will be removed from 'list' when returned.
 
void remove (TMLink &a)
 removes a TMLink.
 
void remove (const AList< TMLink > &)
 removes TMLinks.
 
virtual void refine (AList< TMLink > &list, double maxSigma)
 removes bad points by pull. The bad points are removed from the track, and are returned in 'list'.
 
virtual void refine (double maxSigma)
 removes bad points by pull. The bad points are masked not to be used in fit.
 
virtual int DropWorst ()
 
virtual void removeLinks (void)
 
virtual double distance (const TMLink &) const
 returns distance to a position of TMLink in TMLink space.
 
virtual int approach (TMLink &) const
 calculates the closest approach to a wire in real space. Results are stored in TMLink. Return value is negative if error happened.
 
unsigned testByApproach (const TMLink &list, double sigma) const
 returns # of good hits to be appended.
 
unsigned testByApproach (const AList< TMLink > &list, double sigma) const
 
virtual int fit (void)
 fits itself by a default fitter. Error was happened if return value is not zero.
 
const TMFitter *const fitter (void) const
 returns a pointer to a default fitter.
 
const TMFitter *const fitter (const TMFitter *)
 sets a default fitter.
 
void falseFit ()
 false Fit
 
TMLink * operator[] (unsigned i) const
 
const TTrackHEP *const hep (void) const
 returns TTrackHEP.
 
unsigned nHeps (void) const
 returns # of contributed TTrackHEP tracks.
 
const TTrackMC *const mc (void) const
 returns a pointer to TTrackMC.
 
bool fitted (void) const
 returns true if fitted.
 
bool fittedWithCathode (void) const
 returns true if fitted with cathode hits(TEMPORARY).
 
 TTrackBase ()
 Constructor.
 
 TTrackBase (const AList< TMLink > &links)
 Constructor.
 
virtual ~TTrackBase ()
 Destructor.
 
virtual unsigned objectType (void) const
 returns object type.
 
virtual unsigned type (void) const
 returns type. Definition is depending on an object class.
 
virtual void dump (const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
 dumps debug information.
 
const AList< TMLink > & links (unsigned mask=0) const
 returns a list of masked TMLinks assigned to this track. 'mask' will be applied if mask is not 0.
 
unsigned nLinks (unsigned mask=0) const
 returns # of masked TMLinks assigned to this track object.
 
const AList< TMLink > & cores (unsigned mask=0) const
 returns a list of masked TMLinks for fit. 'mask' will be applied if mask is not 0.
 
unsigned nCores (unsigned mask=0) const
 returns # of masked TMLinks for fit. 'mask' will be applied if mask is not 0.
 
void update (void) const
 update cache.
 
void append (TMLink &)
 appends a TMLink.
 
void append (const AList< TMLink > &)
 appends TMLinks.
 
void appendByApproach (AList< TMLink > &list, double maxSigma)
 appends TMLinks by approach. 'list' is an input. Unappended TMLinks will be removed from 'list' when returned.
 
void appendByDistance (AList< TMLink > &list, double maxDistance)
 appends TMLinks by distance. 'list' is an input. Unappended TMLinks will be removed from 'list' when returned.
 
void remove (TMLink &a)
 removes a TMLink.
 
void remove (const AList< TMLink > &)
 removes TMLinks.
 
virtual void refine (AList< TMLink > &list, double maxSigma)
 removes bad points by pull. The bad points are removed from the track, and are returned in 'list'.
 
virtual void refine (double maxSigma)
 removes bad points by pull. The bad points are masked not to be used in fit.
 
virtual int DropWorst ()
 
virtual void removeLinks (void)
 
virtual double distance (const TMLink &) const
 returns distance to a position of TMLink in TMLink space.
 
virtual int approach (TMLink &) const
 calculates the closest approach to a wire in real space. Results are stored in TMLink. Return value is negative if error happened.
 
unsigned testByApproach (const TMLink &list, double sigma) const
 returns # of good hits to be appended.
 
unsigned testByApproach (const AList< TMLink > &list, double sigma) const
 
virtual int fit (void)
 fits itself by a default fitter. Error was happened if return value is not zero.
 
const TMFitter *const fitter (void) const
 returns a pointer to a default fitter.
 
const TMFitter *const fitter (const TMFitter *)
 sets a default fitter.
 
void falseFit ()
 false Fit
 
TMLink * operator[] (unsigned i) const
 
const TTrackHEP *const hep (void) const
 returns TTrackHEP.
 
unsigned nHeps (void) const
 returns # of contributed TTrackHEP tracks.
 
const TTrackMC *const mc (void) const
 returns a pointer to TTrackMC.
 
bool fitted (void) const
 returns true if fitted.
 
bool fittedWithCathode (void) const
 returns true if fitted with cathode hits(TEMPORARY).
 

Friends

class TRungeFitter
 

Additional Inherited Members

- Protected Attributes inherited from TTrackBase
AList< TMLink > _links
 
bool _fitted
 
bool _fittedWithCathode
 
TTrackMC * _mc
 

Detailed Description

A class to represent a track in tracking.

Definition at line 65 of file InstallArea/include/TrkReco/TrkReco/TRunge.h.

Constructor & Destructor Documentation

◆ TRunge() [1/10]

TRunge::TRunge ( )

Constructors.

Definition at line 47 of file TRunge.cxx.

48 : TTrackBase(),
49 _pivot(ORIGIN),_a(Vector(5,0)),_Ea(SymMatrix(5,0)),
50 _chi2(0),_ndf(0),_charge(1),_bfieldID(21),_Nstep(0),
51 _stepSize(default_stepSize),_mass(0.140),_eps(EPS),
52 _stepSizeMax(default_stepSizeMax),_stepSizeMin(default_stepSizeMin),
53 _maxflightlength(default_maxflightlength) {
54
55 //...Set a default fitter...
56 fitter(& TRunge::_fitter);
57
58 _fitted = false;
59 _fittedWithCathode = false;
60
61 _bfield = Bfield::getBfield(_bfieldID);
62
63 _mass2=_mass*_mass;
64
65 _yscal[0]=0.1; //1mm -> error = 1mm*EPS
66 _yscal[1]=0.1; //1mm
67 _yscal[2]=0.1; //1mm
68 _yscal[3]=0.001; //1MeV
69 _yscal[4]=0.001; //1MeV
70 _yscal[5]=0.001; //1MeV
71 //jialk
72 StatusCode scmgn = Gaudi::svcLocator()->service("MagneticFieldSvc",m_pmgnIMF);
73 if(scmgn!=StatusCode::SUCCESS) {
74 std::cout<< __FILE__<<" Unable to open Magnetic field service"<<std::endl;
75 }
76
77}
SymMatrix
CLHEP::HepSymMatrix SymMatrix
Definition: InstallArea/include/TrackUtil/TrackUtil/Helix.h:40
ORIGIN
const HepPoint3D ORIGIN
Constants.
Definition: TMDCUtil.cxx:47
default_maxflightlength
const double default_maxflightlength
Definition: TRunge.cxx:45
default_stepSizeMin
const double default_stepSizeMin
Definition: TRunge.cxx:41
default_stepSize
const double default_stepSize
Definition: TRunge.cxx:38
default_stepSizeMax
const double default_stepSizeMax
Definition: TRunge.cxx:40
EPS
const double EPS
Definition: TRunge.cxx:43
Bfield::getBfield
static Bfield * getBfield(int)
returns Bfield object.
Definition: KalFitAlg/KalFitAlg-00-07-55-p03/src/coil/Bfield.cxx:52
TTrackBase::fitter
const TMFitter *const fitter(void) const
returns a pointer to a default fitter.
Definition: InstallArea/include/TrkReco/TrkReco/TTrackBase.h:255
TTrackBase::TTrackBase
TTrackBase()
Constructor.
Definition: TTrackBase.cxx:40

◆ TRunge() [2/10]

TRunge::TRunge ( const TTrack &  a)

Definition at line 114 of file TRunge.cxx.

115 : TTrackBase(a),
116 _pivot(a.helix().pivot()),_a(a.helix().a()),_Ea(a.helix().Ea()),
117 _chi2(0),_ndf(0),_bfieldID(21),_Nstep(0),
118 _stepSize(default_stepSize),_mass(0.140),_eps(EPS),
119 _stepSizeMax(default_stepSizeMax),_stepSizeMin(default_stepSizeMin),
120 _maxflightlength(default_maxflightlength) {
121// _a[2]=-_a[2];
122 //...Set a default fitter...
123 fitter(& TRunge::_fitter);
124
125 _fitted = false;
126 _fittedWithCathode = false;
127
128 if(_a[2]<0) _charge=-1;
129 else _charge=1;
130
131 _bfield = Bfield::getBfield(_bfieldID);
132
133 _mass2=_mass*_mass;
134
135 _yscal[0]=0.1; //1mm -> error = 1mm*EPS
136 _yscal[1]=0.1; //1mm
137 _yscal[2]=0.1; //1mm
138 _yscal[3]=0.001; //1MeV
139 _yscal[4]=0.001; //1MeV
140 _yscal[5]=0.001; //1MeV
141 StatusCode scmgn = Gaudi::svcLocator()->service("MagneticFieldSvc",m_pmgnIMF);
142 if(scmgn!=StatusCode::SUCCESS) {
143 std::cout<< "Unable to open Magnetic field service"<<std::endl;
144 }
145 //SetFlightLength();
146 //cout<<"TR:: _maxflightlength="<<_maxflightlength<<endl;
147}
TRunge::a
const Vector & a(void) const
returns helix parameter
Definition: TRunge.cxx:230

◆ TRunge() [3/10]

TRunge::TRunge ( const Helix &  h)

Definition at line 149 of file TRunge.cxx.

150 : TTrackBase(),
151 _pivot(h.pivot()),_a(h.a()),_Ea(h.Ea()),
152 _chi2(0),_ndf(0),_bfieldID(21),_Nstep(0),
153 _stepSize(default_stepSize),_mass(0.140),_eps(EPS),
154 _stepSizeMax(default_stepSizeMax),_stepSizeMin(default_stepSizeMin),
155 _maxflightlength(default_maxflightlength) {
156
157 //...Set a default fitter...
158 fitter(& TRunge::_fitter);
159
160 _fitted = false;
161 _fittedWithCathode = false;
162
163 if(_a[2]<0) _charge=-1;
164 else _charge=1;
165
166 _bfield = Bfield::getBfield(_bfieldID);
167
168 _mass2=_mass*_mass;
169
170 _yscal[0]=0.1; //1mm -> error = 1mm*EPS
171 _yscal[1]=0.1; //1mm
172 _yscal[2]=0.1; //1mm
173 _yscal[3]=0.001; //1MeV
174 _yscal[4]=0.001; //1MeV
175 _yscal[5]=0.001; //1MeV
176}
Helix::Ea
const HepSymMatrix & Ea(void) const
returns error matrix.
Definition: InstallArea/include/TrackUtil/TrackUtil/Helix.h:269
Helix::a
const HepVector & a(void) const
returns helix parameters.
Definition: InstallArea/include/TrackUtil/TrackUtil/Helix.h:263
Helix::pivot
const HepPoint3D & pivot(void) const
returns pivot position.
Definition: InstallArea/include/TrackUtil/TrackUtil/Helix.h:209

◆ TRunge() [4/10]

TRunge::TRunge ( const TRunge &  a)

Definition at line 178 of file TRunge.cxx.

179 : TTrackBase((TTrackBase &) a),
180 _pivot(a.pivot()),_a(a.a()),_Ea(a.Ea()),
181 _chi2(a.chi2()),_ndf(a.ndf()),_bfieldID(a.BfieldID()),_Nstep(0),
182 _stepSize(a.StepSize()),_mass(a.Mass()),_eps(a.Eps()),
183 _stepSizeMax(a.StepSizeMax()),_stepSizeMin(a.StepSizeMin()),
184 _maxflightlength(a.MaxFlightLength()) {
185
186 //...Set a default fitter...
187 fitter(& TRunge::_fitter);
188
189 _fitted = false;
190 _fittedWithCathode = false;
191
192 if(_a[2]<0) _charge=-1;
193 else _charge=1;
194
195 _bfield = Bfield::getBfield(_bfieldID);
196
197 _mass2=_mass*_mass;
198
199 for(unsigned i=0;i<6;i++) _yscal[i]=a.Yscal()[i];
200}
TTrackBase
A virtual class for a track class in tracking.
Definition: InstallArea/include/TrkReco/TrkReco/TTrackBase.h:46

◆ TRunge() [5/10]

TRunge::TRunge ( const RecMdcTrack &  a)

Definition at line 79 of file TRunge.cxx.

80 : TTrackBase(),
81 _pivot(a.getPivot()),_a(a.helix()),_Ea(a.err()),
82 _chi2(a.chi2()),_ndf(a.ndof()),_bfieldID(21),_Nstep(0),
83 _stepSize(default_stepSize),_mass(0.140),_eps(EPS),
84 _stepSizeMax(default_stepSizeMax),_stepSizeMin(default_stepSizeMin),
85 _maxflightlength(default_maxflightlength) {
86// _a[2]=-_a[2];
87 //...Set a default fitter...
88 fitter(& TRunge::_fitter);
89
90 _fitted = false;
91 _fittedWithCathode = false;
92
93 if(_a[2]<0) _charge=-1;
94 else _charge=1;
95
96 _bfield = Bfield::getBfield(_bfieldID);
97
98 _mass2=_mass*_mass;
99
100 _yscal[0]=0.1; //1mm -> error = 1mm*EPS
101 _yscal[1]=0.1; //1mm
102 _yscal[2]=0.1; //1mm
103 _yscal[3]=0.001; //1MeV
104 _yscal[4]=0.001; //1MeV
105 _yscal[5]=0.001; //1MeV
106 StatusCode scmgn = Gaudi::svcLocator()->service("MagneticFieldSvc",m_pmgnIMF);
107 if(scmgn!=StatusCode::SUCCESS) {
108 std::cout<< "Unable to open Magnetic field service"<<std::endl;
109 }
110 //SetFlightLength();
111 //cout<<"TR:: _maxflightlength="<<_maxflightlength<<endl;
112
113 }

◆ ~TRunge() [1/2]

TRunge::~TRunge ( )

Destructor.

Definition at line 203 of file TRunge.cxx.

203 {
204}

◆ TRunge() [6/10]

TRunge::TRunge ( )

Constructors.

◆ TRunge() [7/10]

TRunge::TRunge ( const TTrack &  )

◆ TRunge() [8/10]

TRunge::TRunge ( const Helix &  )

◆ TRunge() [9/10]

TRunge::TRunge ( const TRunge &  )

◆ TRunge() [10/10]

TRunge::TRunge ( const RecMdcTrack &  )

◆ ~TRunge() [2/2]

TRunge::~TRunge ( )

Destructor.

Member Function Documentation

◆ a() [1/4]

const Vector & TRunge::a ( const Vector &  ta)

set helix parameter

Definition at line 299 of file TRunge.cxx.

299 {
300 _a=ta;
301 if(_a[2]<0) _charge=-1;
302 else _charge=1;
303 _Nstep=0;
304 return(_a);
305}

◆ a() [2/4]

const Vector & TRunge::a ( const Vector &  )

set helix parameter

◆ a() [3/4]

const Vector & TRunge::a ( void  ) const

returns helix parameter

Definition at line 230 of file TRunge.cxx.

230 {
231 return(_a);
232}

Referenced by dEpath(), DisToWire(), and TRunge().

◆ a() [4/4]

const Vector & TRunge::a ( void  ) const

returns helix parameter

◆ approach() [1/4]

int TRunge::approach ( TMLink &  l,
const RkFitMaterial  material,
bool  sagCorrection = true 
)

calculates the closest approach to a wire in real space. Results are stored in TMLink. Return value is negative if error happened.

Definition at line 360 of file TRunge.cxx.

360 {
361 float tof;
362 HepVector3D p;
363 return(approach(l,tof,p,material,doSagCorrection));
364}
TRunge::approach
int approach(TMLink &, const RkFitMaterial, bool sagCorrection=true)
Definition: TRunge.cxx:360

Referenced by approach().

◆ approach() [2/4]

int TRunge::approach ( TMLink &  ,
const  RkFitMaterial,
bool  sagCorrection = true 
)

calculates the closest approach to a wire in real space. Results are stored in TMLink. Return value is negative if error happened.

◆ approach() [3/4]

int TRunge::approach ( TMLink &  l,
float &  tof,
HepVector3D &  p,
const RkFitMaterial  material,
bool  sagCorrection = true 
)

Definition at line 366 of file TRunge.cxx.

367 {
368 HepPoint3D xw;
369 HepPoint3D wireBackwardPosition;
370 HepVector3D v;
371 HepPoint3D onWire,onTrack;
372 double onWire_x,onWire_y, onWire_z, zwf, zwb;
373 const TMDCWire& w=*l.wire();
374 xw = w.xyPosition();
375 wireBackwardPosition = w.backwardPosition();
376 v = w.direction();
377
378 unsigned stepNum=0;
379 if(approach_line(l,wireBackwardPosition,v,onWire,onTrack,tof,p,material,stepNum)<0){
380 //cout<<"TR::error approach_line"<<endl;
381 return(-1);
382 }
383 zwf = w.forwardPosition().z();
384 zwb = w.backwardPosition().z();
385 // protection for strange wire hit info. (Only 'onWire' is corrected)
386 if(onWire.z() > zwf)
387 w.wirePosition(zwf,onWire,wireBackwardPosition,(HepVector3D&)v);
388 else if(onWire.z() < zwb)
389 w.wirePosition(zwb,onWire,wireBackwardPosition,(HepVector3D&)v);
390
391 // onWire,onTrack filled
392
393 if(!doSagCorrection){
394 l.positionOnWire(onWire);
395 l.positionOnTrack(onTrack);
396 double phipoint=atan((onTrack.y()-_pivot.y())/onTrack.x()-_pivot.x());
397// cout<<"dphi track:"<<l.dPhi()<<endl;
398 // cout<<"dphi rk:"<<phipoint-_a[0]<<endl;
399// l.dPhi(phipoint-_a[0]);
400 return(0); // no sag correction
401 }
402 // Sag correction
403 // loop for sag correction
404 onWire_y = onWire.y();
405 onWire_z = onWire.z();
406
407 unsigned nTrial = 1;
408 while(nTrial<100){
409 //cout<<"TR: nTrial "<<nTrial<<endl;
410 w.wirePosition(onWire_z,xw,wireBackwardPosition,(HepVector3D&)v);
411 if(approach_line(l,wireBackwardPosition,v,onWire,onTrack,tof,p,material,stepNum)<0)
412 return(-1);
413 if(fabs(onWire_y - onWire.y())<0.0001) break; // |dy|< 1 micron
414 onWire_y = onWire.y();
415 onWire_z += (onWire.z()-onWire_z)/2;
416 // protection for strange wire hit info.
417 if(onWire_z > zwf) onWire_z=zwf;
418 else if(onWire_z < zwb) onWire_z=zwb;
419
420 nTrial++;
421 }
422 // cout<<"TR nTrial="<<nTrial<<endl;
423 l.positionOnWire(onWire);
424 l.positionOnTrack(onTrack);
425 return(nTrial);
426}
w
double w
Definition: EvtPhokharaDef.hh:26
v
**********Class see also m_nmax DOUBLE PRECISION m_amel DOUBLE PRECISION m_x2 DOUBLE PRECISION m_alfinv DOUBLE PRECISION m_Xenph INTEGER m_KeyWtm INTEGER m_idyfs DOUBLE PRECISION m_zini DOUBLE PRECISION m_q2 DOUBLE PRECISION m_Wt_KF DOUBLE PRECISION m_WtCut INTEGER m_KFfin *COMMON c_KarLud $ !Input CMS energy[GeV] $ !CMS energy after beam spread beam strahlung[GeV] $ !Beam energy spread[GeV] $ !z boost due to beam spread $ !electron beam mass *ff pair spectrum $ !minimum v
Definition: KarLud.h:35
TMDCWire
A class to represent a wire in MDC.
Definition: InstallArea/include/TrkReco/TrkReco/TMDCWire.h:55
TRunge::approach_line
int approach_line(TMLink &, const HepPoint3D &, const HepVector3D &, HepPoint3D &onLine, HepPoint3D &onTrack, const RkFitMaterial material)
caluculate closest points between a line and this track
Definition: TRunge.cxx:428

◆ approach() [4/4]

int TRunge::approach ( TMLink &  ,
float &  tof,
HepVector3D &  p,
const  RkFitMaterial,
bool  sagCorrection = true 
)

◆ approach_line() [1/6]

int TRunge::approach_line ( TMLink &  l,
const HepPoint3D &  w0,
const HepVector3D &  v,
HepPoint3D &  onLine,
HepPoint3D &  onTrack,
const RkFitMaterial  material 
)

caluculate closest points between a line and this track

Definition at line 428 of file TRunge.cxx.

429 {
430 float tof;
431 HepVector3D p;
432 return(approach_line(l,w0,v,onLine,onTrack,tof,p,material));
433}

Referenced by approach(), and approach_line().

◆ approach_line() [2/6]

int TRunge::approach_line ( TMLink &  ,
const HepPoint3D &  ,
const HepVector3D &  ,
HepPoint3D &  onLine,
HepPoint3D &  onTrack,
const RkFitMaterial  material 
)

caluculate closest points between a line and this track

◆ approach_line() [3/6]

int TRunge::approach_line ( TMLink &  l,
const HepPoint3D &  w0,
const HepVector3D &  v,
HepPoint3D &  onLine,
HepPoint3D &  onTrack,
float &  tof,
HepVector3D &  p,
const RkFitMaterial  material 
)

Definition at line 435 of file TRunge.cxx.

437 {
438 unsigned stepNum=0;
439 return(approach_line(l,w0,v,onLine,onTrack,tof,p,material ,stepNum));
440}

◆ approach_line() [4/6]

int TRunge::approach_line ( TMLink &  ,
const HepPoint3D &  ,
const HepVector3D &  ,
HepPoint3D &  onLine,
HepPoint3D &  onTrack,
float &  tof,
HepVector3D &  p,
const RkFitMaterial  material 
)

◆ approach_line() [5/6]

int TRunge::approach_line ( TMLink &  l,
const HepPoint3D &  w0,
const HepVector3D &  v,
HepPoint3D &  onLine,
HepPoint3D &  onTrack,
float &  tof,
HepVector3D &  p,
const RkFitMaterial  material,
unsigned &  stepNum 
)

Definition at line 442 of file TRunge.cxx.

444 {
445 // line = [w0] + t * [v] -> [onLine]
446 if(_Nstep==0){
447 if(_stepSize==0) Fly_SC();
448 else Fly(material);
449
450 }
451 //cout<<"TR::approach_line stepNum="<<stepNum<<endl;
452
453 const double w0x = w0.x();
454 const double w0y = w0.y();
455 const double w0z = w0.z();
456 //cout<<"TR::line w0="<<w0x<<","<<w0y<<","<<w0z<<endl;
457 const double vx = v.x();
458 const double vy = v.y();
459 const double vz = v.z();
460 const double v_2 = vx*vx+vy*vy+vz*vz;
461
462 const float clight=29.9792458; //[cm/ns]
463 //const float M2=_mass*_mass;
464 //const float& M2=_mass2;
465 const float p2 = _y[0][3]*_y[0][3]+_y[0][4]*_y[0][4]+_y[0][5]*_y[0][5];
466 const float tof_factor=1./clight*sqrt(1+_mass2/p2);
467
468 // search for the closest point in cache (point - line)
469 // unsigned stepNum;
470 double l2_old= DBL_MAX;
471 if(stepNum>_Nstep) stepNum=0;
472 unsigned stepNumLo;
473 unsigned stepNumHi;
474 if(stepNum==0 && _stepSize!=0){ // skip
475 const double dx = _y[0][0]-w0x;
476 const double dy = _y[0][1]-w0y;
477 stepNum=(unsigned)
478 (sqrt( (dx*dx+dy*dy)*(1+_a[4]*_a[4]) )/_stepSize);
479 }
480 unsigned mergin;
481 if(_stepSize==0){
482 mergin=10; //10 step back
483 }else{
484 mergin=(unsigned)(1.0/_stepSize); //1mm back
485 }
486 if(stepNum>mergin) stepNum-=mergin;
487 else stepNum=0;
488 if(stepNum>=_Nstep) stepNum=_Nstep-1;
489 // hunt
490 // unsigned inc=1;
491 unsigned inc=(mergin>>1)+1;
492 stepNumLo=stepNum;
493 stepNumHi=stepNum;
494 for(;;){
495 const double dx = _y[stepNumHi][0]-w0x;
496 const double dy = _y[stepNumHi][1]-w0y;
497 const double dz = _y[stepNumHi][2]-w0z;
498 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
499// const double l2 = (dx*dx+dy*dy+dz*dz)-(t*t)/v_2;
500 HepVector3D zvec;
501 zvec.setX(0);
502 zvec.setY(0);
503 zvec.setZ(1);
504 HepVector3D hitv=t*v;
505 HepPoint3D pp,onw;
506 pp.setX(_y[stepNumHi][0]);
507 pp.setY(_y[stepNumHi][1]);
508 pp.setZ(_y[stepNumHi][2]);
509 double zwire=hitv.dot(zvec)+w0z;
510 double dtl=DisToWire(l,zwire,pp,onw);
511 const double l2=dtl*dtl;
512 if(l2 > l2_old) break;
513 l2_old=l2;
514 stepNumLo=stepNumHi;
515 // inc+=inc;
516 stepNumHi+=inc;
517 if(stepNumHi>=_Nstep){
518 stepNumHi=_Nstep;
519 break;
520 }
521 }
522 // locate (2-bun-hou, bisection method)
523 while(stepNumHi-stepNumLo>1){
524 unsigned j=(stepNumHi+stepNumLo)>>1;
525 const double dx = _y[j][0]-w0x;
526 const double dy = _y[j][1]-w0y;
527 const double dz = _y[j][2]-w0z;
528 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
529// const double l2 = (dx*dx+dy*dy+dz*dz)-(t*t)/v_2;
530 HepVector3D zv;
531 zv.setX(0);
532 zv.setY(0);
533 zv.setZ(1);
534 HepVector3D hitv=t*v;
535 HepPoint3D point,onwir;
536 point.setX(_y[j][0]);
537 point.setY(_y[j][1]);
538 point.setZ(_y[j][2]);
539 double zwir=hitv.dot(zv)+w0z;
540 double dtll=DisToWire(l,zwir,point,onwir);
541 const double l2 =dtll*dtll;
542 if(l2 > l2_old){
543 stepNumHi=j;
544 }else{
545 l2_old=l2;
546 stepNumLo=j;
547 }
548 }
549 //stepNum=stepNumHi;
550 stepNum=stepNumLo;
551 /*
552 for(;stepNum<_Nstep;stepNum++){
553 const double dx = _y[stepNum][0]-w0x;
554 const double dy = _y[stepNum][1]-w0y;
555 const double dz = _y[stepNum][2]-w0z;
556 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
557 const double l2 = (dx*dx+dy*dy+dz*dz)-(t*t)/v_2;
558 if(l2 > l2_old) break;
559 l2_old=l2;
560 // const float p2 = _y[stepNum][3]*_y[stepNum][3]+
561 // _y[stepNum][4]*_y[stepNum][4]+
562 // _y[stepNum][5]*_y[stepNum][5];
563 // tof+=_stepSize/clight*sqrt(1+M2/p2);
564 }
565 */
566 // cout<<"TR stepNum="<<stepNum<<endl;
567 //if(stepNum>=_Nstep) return(-1); // not found
568 //stepNum--;
569 if(_stepSize==0){
570 double dstep=0;
571 for(unsigned i=0;i<stepNum;i++) dstep+=_h[i];
572 tof=dstep*tof_factor;
573 }else{
574 tof=stepNum*_stepSize*tof_factor;
575 }
576
577 // propagate the track and search for the closest point
578 //2-bun-hou (bisection method)
579 /*
580 double y[6],y_old[6];
581 for(unsigned i=0;i<6;i++) y[i]=_y[stepNum][i];
582 double step=_stepSize;
583 double step2=0;
584 for(;;){
585 for(unsigned i=0;i<6;i++) y_old[i]=y[i];
586 Propagate(y,step);
587 const double dx = y[0]-w0x;
588 const double dy = y[1]-w0y;
589 const double dz = y[2]-w0z;
590 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
591 const double l2 = (dx*dx+dy*dy+dz*dz)-(t*t)/v_2;
592 if(l2 > l2_old){ // back
593 for(unsigned i=0;i<6;i++) y[i]=y_old[i];
594 }else{ // propagate
595 l2_old=l2;
596 // const float p2=y[3]*y[3]+y[4]*y[4]+y[5]*y[5];
597 // tof+=step/clight*sqrt(1+M2/p2);
598 step2+=step;
599 }
600 step/=2;
601 if(step < 0.0001) break; // step < 1 [um]
602 }
603 */
604 // Hasamiuchi-Hou (false position method)
605 /*
606 double y[6],y1[6],y2[6];
607 for(unsigned i=0;i<6;i++) y1[i]=_y[stepNum][i];
608 for(unsigned i=0;i<6;i++) y2[i]=_y[stepNum+2][i];
609 double minStep=0;
610 double maxStep=_stepSize*2;
611 double step2=0;
612 const double A[3][3]={{1-vx*vx/v_2, -vx*vy/v_2, -vx*vz/v_2},
613 { -vy*vx/v_2,1-vy*vy/v_2, -vy*vz/v_2},
614 { -vz*vx/v_2, -vz*vy/v_2,1-vz*vz/v_2}};
615
616 double Y1=0;
617 {
618 const double dx = y1[0]-w0x;
619 const double dy = y1[1]-w0y;
620 const double dz = y1[2]-w0z;
621 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
622 const double d[3]={dx-t*vx,dy-t*vy,dz-t*vz};
623 const double l = sqrt( (dx*dx+dy*dy+dz*dz)-(t*t)/v_2 );
624 const double pmag=sqrt(y1[3]*y1[3]+y1[4]*y1[4]+y1[5]*y1[5]);
625 for(int j=0;j<3;j++){
626 double g=0;
627 for(int k=0;k<3;k++) g+=A[j][k]*d[k];
628 Y1+=y1[j+3]*g;
629 }
630 Y1=Y1/pmag/l;
631 }
632 double Y2=0;
633 {
634 const double dx = y2[0]-w0x;
635 const double dy = y2[1]-w0y;
636 const double dz = y2[2]-w0z;
637 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
638 const double d[3]={dx-t*vx,dy-t*vy,dz-t*vz};
639 const double l = sqrt( (dx*dx+dy*dy+dz*dz)-(t*t)/v_2 );
640 const double pmag=sqrt(y2[3]*y2[3]+y2[4]*y2[4]+y2[5]*y2[5]);
641 for(int j=0;j<3;j++){
642 double g=0;
643 for(int k=0;k<3;k++) g+=A[j][k]*d[k];
644 Y2+=y2[j+3]*g;
645 }
646 Y2=Y2/pmag/l;
647 }
648 if(Y1*Y2>=0) cout<<"TR fatal error!"<<endl;
649 double step_old= DBL_MAX;
650 for(;;){
651 const double step=(minStep*Y2-maxStep*Y1)/(Y2-Y1);
652 if(fabs(step-step_old)<0.0001) break;
653 step_old=step;
654 for(unsigned i=0;i<6;i++) y[i]=y1[i];
655 Propagate(y,step-minStep);
656 double Y=0;
657 {
658 const double dx = y[0]-w0x;
659 const double dy = y[1]-w0y;
660 const double dz = y[2]-w0z;
661 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
662 const double d[3]={dx-t*vx,dy-t*vy,dz-t*vz};
663 const double l = sqrt( (dx*dx+dy*dy+dz*dz)-(t*t)/v_2 );
664 const double pmag=sqrt(y[3]*y[3]+y[4]*y[4]+y[5]*y[5]);
665 for(int j=0;j<3;j++){
666 double g=0;
667 for(int k=0;k<3;k++) g+=A[j][k]*d[k];
668 Y+=y[j+3]*g;
669 }
670 Y=Y/pmag/l;
671 }
672 if(Y1*Y<0){ //Y->Y2
673 Y2=Y;
674 for(unsigned i=0;i<6;i++) y2[i]=y[i];
675 maxStep=step;
676 }else{ //Y->Y1
677 Y1=Y;
678 for(unsigned i=0;i<6;i++) y1[i]=y[i];
679 minStep=step;
680 }
681 if(Y1==Y2) break;
682 }
683 step2=step_old;
684 */
685 // Newton method
686 double y[6];
687 for(unsigned i=0;i<6;i++) y[i]=_y[stepNum][i];
688 //memcpy(y,_y[stepNum],sizeof(double)*6);
689 double step2=0;
690 const double A[3][3]={{1-vx*vx/v_2, -vx*vy/v_2, -vx*vz/v_2},
691 { -vy*vx/v_2,1-vy*vy/v_2, -vy*vz/v_2},
692 { -vz*vx/v_2, -vz*vy/v_2,1-vz*vz/v_2}};
693 double factor=1;
694 double g[3];
695 double f[6];
696 double Af[3],Af2[3];
697 unsigned i,j,k;
698 HepPoint3D point,onwir;
699 for(i=0;i<10;i++){
700 //cout<<"TR::line "<<y[0]<<","<<y[1]<<","<<y[2]<<","<<y[3]<<","<<y[4]<<","<<y[5]<<endl;
701 const double dx = y[0]-w0x;
702 const double dy = y[1]-w0y;
703 const double dz = y[2]-w0z;
704 const double t = (dx*vx+dy*vy+dz*vz)/v_2;
705 const double d[3]={dx-t*vx,dy-t*vy,dz-t*vz};
706 // const double l2 = (dx*dx+dy*dy+dz*dz)-(t*t)/v_2;
707 HepVector3D zv;
708 zv.setX(0);
709 zv.setY(0);
710 zv.setZ(1);
711 HepVector3D hitv=t*v;
712 point.setX(y[0]);
713 point.setY(y[1]);
714 point.setZ(y[2]);
715 double zwir=hitv.dot(zv)+w0z;
716 double dtll=DisToWire(l,zwir,point,onwir);
717 const double l2=dtll*dtll;
718 for(j=0;j<3;j++){
719 g[j]=0;
720 for(k=0;k<3;k++) g[j]+=A[j][k]*d[k];
721 }
722 //cout<<"g="<<g[0]<<","<<g[1]<<","<<g[2]<<endl;
723 Function(y,f);
724 //cout<<"f="<<f[0]<<","<<f[1]<<","<<f[2]<<","<<f[3]<<","<<f[4]<<","<<f[5]<<endl;
725 //cout<<"A="<<A[0][0]<<","<<A[0][1]<<","<<A[0][2]<<endl;
726 //cout<<"A="<<A[1][0]<<","<<A[1][1]<<","<<A[1][2]<<endl;
727 //cout<<"A="<<A[2][0]<<","<<A[2][1]<<","<<A[2][2]<<endl;
728 double Y=0;
729 for(j=0;j<3;j++) Y+=y[j+3]*g[j];
730 double dYds=0;
731 for(j=0;j<3;j++) dYds+=f[j+3]*g[j];
732 //cout<<"dYds="<<dYds<<endl;
733 for(j=0;j<3;j++){
734 Af[j]=0;
735 Af2[j]=0;
736 }
737 for(j=0;j<3;j++){
738 //Af[j]=0;
739 //Af2[j]=0;
740 for(k=0;k<3;k++) Af[j]+=(A[j][k]*f[k]);
741 for(k=0;k<3;k++) Af2[j]+=(A[j][k]*Af[k]);
742 dYds+=y[j+3]*Af2[j];
743 //cout<<j<<" dYds="<<dYds<<endl;
744 }
745 //cout<<"dYds="<<dYds<<endl;
746 const double step=-Y/dYds*factor;
747 //cout<<"TR step="<<step<<" i="<<i<<endl;
748 if(fabs(step) < 0.00001) break; // step < 1 [um]
749 if(l2 > l2_old) factor/=2;
750 l2_old=l2;
751 Propagate(y,step,material);
752 step2+=step;
753 }
754
755 tof+=step2*tof_factor;
756
757 onTrack.setX(y[0]);
758 onTrack.setY(y[1]);
759 onTrack.setZ(y[2]);
760 p.setX(y[3]);
761 p.setY(y[4]);
762 p.setZ(y[5]);
763
764 const double dx = y[0]-w0x;
765 const double dy = y[1]-w0y;
766 const double dz = y[2]-w0z;
767 const double s = (dx*vx+dy*vy+dz*vz)/v_2;
768 // onLine.setX(w0x+s*vx);
769 // onLine.setY(w0y+s*vy);
770 // onLine.setZ(w0z+s*vz);
771 onLine.setX(onwir.x());
772 onLine.setY(onwir.y());
773 onLine.setZ(onwir.z());
774 return(0);
775}
s
XmlRpcServer s
Definition: HelloServer.cpp:11
t
TTree * t
Definition: binning.cxx:23
TRunge::Fly_SC
unsigned Fly_SC(void) const
Definition: TRunge.cxx:1097
TRunge::Fly
unsigned Fly(const RkFitMaterial material) const
make the trajectory in cache, return the number of step
Definition: TRunge.cxx:836
TRunge::DisToWire
double DisToWire(TMLink &, double, HepPoint3D, HepPoint3D &)
Definition: TRunge.cxx:1215
TRunge::dz
double dz(void) const
Definition: TRunge.cxx:218
TRunge::Function
void Function(const double y[6], double f[6]) const
Definition: TRunge.cxx:916
TRunge::Propagate
void Propagate(double y[6], const double &step, const RkFitMaterial material) const
propagate the track using 4th-order Runge-Kutta method
Definition: TRunge.cxx:862
y
double y[1000]
Definition: draw_charge_space_coarse_II.cxx:8
f
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")

◆ approach_line() [6/6]

int TRunge::approach_line ( TMLink &  ,
const HepPoint3D &  ,
const HepVector3D &  ,
HepPoint3D &  onLine,
HepPoint3D &  onTrack,
float &  tof,
HepVector3D &  p,
const RkFitMaterial  material,
unsigned &  stepNum 
)

◆ approach_point() [1/2]

int TRunge::approach_point ( TMLink &  l,
const HepPoint3D &  p0,
HepPoint3D &  onTrack,
const RkFitMaterial  material 
) const

caluculate closest point between a point and this track

Definition at line 777 of file TRunge.cxx.

777 {
778 if(_Nstep==0){
779 if(_stepSize==0) Fly_SC();
780 else Fly(material);
781 }
782
783 double x0=p0.x();
784 double y0=p0.y();
785 double z0=p0.z();
786
787 //tof=0;
788 //const float clight=29.9792458; //[cm/ns]
789 //const float M2=_mass*_mass;
790
791 // search for the closest point in cache
792 unsigned stepNum;
793 double l2_old= DBL_MAX;
794 for(stepNum=0;stepNum<_Nstep;stepNum++){
795 double l2=(_y[stepNum][0]-x0)*(_y[stepNum][0]-x0)+
796 (_y[stepNum][1]-y0)*(_y[stepNum][1]-y0)+
797 (_y[stepNum][2]-z0)*(_y[stepNum][2]-z0);
798 if(l2 > l2_old) break;
799 l2_old=l2;
800 //const double p2 = _y[stepNum][3]*_y[stepNum][3]+
801 // _y[stepNum][4]*_y[stepNum][4]+
802 // _y[stepNum][5]*_y[stepNum][5];
803 //tof+=_stepSize/clight*sqrt(1+M2/p2);
804 }
805 if(stepNum>=_Nstep) return(-1); // not found
806 stepNum--;
807
808 // propagate the track and search for the closest point
809 double y[6],y_old[6];
810 for(unsigned i=0;i<6;i++) y[i]=_y[stepNum][i];
811 double step=_stepSize;
812 for(;;){
813 for(unsigned i=0;i<6;i++) y_old[i]=y[i];
814 Propagate(y,step,material);
815 double l2=(y[0]-x0)*(y[0]-x0)+(y[1]-y0)*(y[1]-y0)+(y[2]-z0)*(y[2]-z0);
816 if(l2 > l2_old){ // back
817 for(unsigned i=0;i<6;i++) y[i]=y_old[i];
818 }else{ // propagate
819 l2_old=l2;
820 //const double p2=y[3]*y[3]+y[4]*y[4]+y[5]*y[5];
821 //tof+=step/clight*sqrt(1+M2/p2);
822 }
823 step/=2;
824 if(step < 0.0001) break; // step < 1 [um]
825 }
826
827 onTrack.setX(y[0]);
828 onTrack.setY(y[1]);
829 onTrack.setZ(y[2]);
830 //p.setX(y[3]);
831 //p.setY(y[4]);
832 //p.setZ(y[5]);
833 return(0);
834}

◆ approach_point() [2/2]

int TRunge::approach_point ( TMLink &  ,
const HepPoint3D &  ,
HepPoint3D &  onTrack,
const RkFitMaterial  material 
) const

caluculate closest point between a point and this track

◆ BfieldID() [1/4]

int TRunge::BfieldID ( int  id)

set B field map ID

Definition at line 313 of file TRunge.cxx.

313 {
314 _bfieldID=id;
315 _bfield = Bfield::getBfield(_bfieldID);
316 _Nstep=0;
317 return(_bfieldID);
318}

◆ BfieldID() [2/4]

int TRunge::BfieldID ( int  )

set B field map ID

◆ BfieldID() [3/4]

int TRunge::BfieldID ( void  ) const

returns B field ID

Definition at line 258 of file TRunge.cxx.

258 {
259 return(_bfieldID);
260}

◆ BfieldID() [4/4]

int TRunge::BfieldID ( void  ) const

returns B field ID

◆ chi2() [1/2]

double TRunge::chi2 ( void  ) const

returns chi2.

Definition at line 246 of file TRunge.cxx.

246 {
247 return _chi2;
248}

◆ chi2() [2/2]

double TRunge::chi2 ( void  ) const

returns chi2.

◆ dEpath() [1/2]

double TRunge::dEpath ( double  mass,
double  path,
double  p 
) const

Definition at line 1239 of file TRunge.cxx.

1239 {
1240 double z=4.72234;
1241 double a=9.29212;
1242 double Density=0.00085144;
1243 double coeff=Density*z/a;
1244 double i=51.4709;
1245 double isq=i * i * 1e-18;
1246// double sigma0_2 = 0.1569*coeff*path;
1247//
1248// if (sigma0_2<0) return 0;
1249//
1250// double psq = p * p;
1251// double bsq = psq / (psq + mass * mass);
1252//
1253// // Correction for relativistic particles :
1254// double sigma_2 = sigma0_2*(1-0.5*bsq)/(1-bsq);
1255//
1256// if (sigma_2<0) return 0;
1257//
1258// // Return sigma in GeV !!
1259// return sqrt(sigma_2)*0.001;
1260const double Me = 0.000510999;
1261 double psq = p * p;
1262 double bsq = psq / (psq + mass * mass);
1263 double esq = psq / (mass * mass);
1264
1265 double s = Me / mass;
1266 double w = (4 * Me * esq
1267 / (1 + 2 * s * sqrt(1 + esq)
1268 + s * s));
1269
1270 // Density correction :
1271 double cc, x0;
1272 cc = 1+2*log(sqrt(isq)/(28.8E-09*sqrt(coeff)));
1273 if (cc < 5.215)
1274 x0 = 0.2;
1275 else
1276 x0 = 0.326*cc-1.5;
1277 double x1(3), xa(cc/4.606), aa;
1278 aa = 4.606*(xa-x0)/((x1-x0)*(x1-x0)*(x1-x0));
1279 double delta(0);
1280double x(log10(sqrt(esq)));
1281 if (x > x0){
1282 delta = 4.606*x-cc;
1283 if (x < x1) delta=delta+aa*(x1-x)*(x1-x)*(x1-x);
1284 }
1285
1286 // Shell correction :
1287 float f1, f2, f3, f4, f5, ce;
1288 f1 = 1/esq;
1289 f2 = f1*f1;
1290 f3 = f1*f2;
1291 f4 = (f1*0.42237+f2*0.0304-f3*0.00038)*1E12;
1292 f5 = (f1*3.858-f2*0.1668+f3*0.00158)*1E18;
1293 ce = f4*isq+f5*isq*sqrt(isq);
1294
1295 return (0.0001535 * coeff / bsq
1296 * (log(Me * esq * w / isq)
1297 - 2 * bsq-delta-2.0*ce/z)) * path;
1298}
delta
const double delta
Definition: CheckRung/check/gainVsDoca.cxx:32
mass
double mass
Definition: CosmicGenerator.cxx:138
f1
TFile * f1
Definition: DataBase/tau_mode.c:4
x
Double_t x[10]
Definition: DataBase/tau_mode.c:57

◆ dEpath() [2/2]

double TRunge::dEpath ( double  ,
double  ,
double   
) const

◆ DisToWire() [1/2]

double TRunge::DisToWire ( TMLink &  l,
double  z,
HepPoint3D  a,
HepPoint3D &  b 
)

Definition at line 1215 of file TRunge.cxx.

1215 {// by liucy 2011/7/28
1216 double zinter[400];
1217 double ddist=999;
1218 double finalz=0;
1219 const TMDCWire& w=*l.wire();
1220 HepPoint3D onwire;
1221// for(int i=0;i<20;i++){
1222 // zinter[i]=z-0.0002*i+0.002;
1223 HepPoint3D point=w.xyPosition(z);
1224// HepPoint3D point=w.xyPosition(zinter[i]);
1225 onwire.setX(point.x());
1226 onwire.setY(point.y());
1227 // onwire.setZ(zinter[i]);
1228 onwire.setZ(z);
1229// double dist=sqrt((point.x()-a.x())*(point.x()-a.x())+(point.y()-a.y())*(point.y()-a.y())+(zinter[i]-a.z())*(zinter[i]-a.z()));
1230 double dist=sqrt((point.x()-a.x())*(point.x()-a.x())+(point.y()-a.y())*(point.y()-a.y())+(z-a.z())*(z-a.z()));
1231 if(dist<ddist){
1232 ddist=dist;
1233// finalz=zinter[i];
1234 b=onwire;
1235 }
1236// }
1237 return ddist;
1238}
b
const double b
Definition: slope.cxx:9

Referenced by approach_line().

◆ DisToWire() [2/2]

double TRunge::DisToWire ( TMLink &  ,
double  ,
HepPoint3D  ,
HepPoint3D &   
)

◆ dr() [1/2]

double TRunge::dr ( void  ) const

Track parameters (at pivot)

Definition at line 206 of file TRunge.cxx.

206 {
207 return(_a[0]);
208}

◆ dr() [2/2]

double TRunge::dr ( void  ) const

Track parameters (at pivot)

◆ dz() [1/2]

double TRunge::dz ( void  ) const

Definition at line 218 of file TRunge.cxx.

218 {
219 return(_a[3]);
220}

Referenced by approach_line(), and intersect_xy_plane().

◆ dz() [2/2]

double TRunge::dz ( void  ) const

◆ Ea() [1/4]

const SymMatrix & TRunge::Ea ( const SymMatrix &  tEa)

set helix error matrix

Definition at line 307 of file TRunge.cxx.

307 {
308 _Ea=tEa;
309 _Nstep=0;
310 return(_Ea);
311}

◆ Ea() [2/4]

const SymMatrix & TRunge::Ea ( const SymMatrix &  )

set helix error matrix

◆ Ea() [3/4]

const SymMatrix & TRunge::Ea ( void  ) const

returns error matrix

Definition at line 234 of file TRunge.cxx.

234 {
235 return(_Ea);
236}

◆ Ea() [4/4]

const SymMatrix & TRunge::Ea ( void  ) const

returns error matrix

◆ eloss() [1/2]

void TRunge::eloss ( double  path,
const RkFitMaterial *  material,
double  mass,
double  y[6],
int  index 
) const

Definition at line 1299 of file TRunge.cxx.

1299 {
1300double psq=y[5]*y[5]+y[3]*y[3]+y[4]*y[4];
1301
1302double p=sqrt(psq);
1303
1304double dE=materiral->dE(mass,path,p);
1305if (index > 0)
1306 psq += dE * (dE + 2 * sqrt(mass * mass + psq));
1307else
1308 psq += dE * (dE - 2 * sqrt(mass * mass + psq));
1309double pnew=sqrt(psq);
1310double coeff=pnew/p;
1311
1312y[3]=y[3]*coeff;
1313y[4]=y[4]*coeff;
1314y[5]=y[5]*coeff;
1315}

Referenced by RkFitCylinder::updateTrack().

◆ eloss() [2/2]

void TRunge::eloss ( double  path,
const RkFitMaterial *  material,
double  mass,
double  y[6],
int  index 
) const

◆ Eps() [1/4]

double TRunge::Eps ( double  eps)

Definition at line 331 of file TRunge.cxx.

331 {
332 _eps=eps;
333 _Nstep=0;
334 return(_eps);
335}
eps
EvtTensor3C eps(const EvtVector3R &v)
Definition: EvtTensor3C.cc:307

◆ Eps() [2/4]

double TRunge::Eps ( double  )

◆ Eps() [3/4]

double TRunge::Eps ( void  ) const

Definition at line 269 of file TRunge.cxx.

269 {
270 return(_eps);
271}

◆ Eps() [4/4]

double TRunge::Eps ( void  ) const

◆ Fly() [1/2]

unsigned TRunge::Fly ( const RkFitMaterial  material) const

make the trajectory in cache, return the number of step

Definition at line 836 of file TRunge.cxx.

836 {
837 double y[6];
838 unsigned Nstep;
839 double flightlength;
840
841 flightlength=0;
842 SetFirst(y);
843 for(Nstep=0;Nstep<TRunge_MAXstep;Nstep++){
844 for(unsigned j=0;j<6;j++) _y[Nstep][j]=y[j];
845 //memcpy(_y[Nstep],y,sizeof(double)*6);
846
847 Propagate(y,_stepSize,material);
848
849 if( y[2]>160 || y[2]<-80 || (y[0]*y[0]+y[1]*y[1])>8100 ) break;
850 //if position is out side of MDC, stop to fly
851 // R>90cm, z<-80cm,160cm<z
852
853 flightlength += _stepSize;
854 if(flightlength>_maxflightlength) break;
855
856 _h[Nstep]=_stepSize;
857 }
858 _Nstep=Nstep+1;
859 return(_Nstep);
860}
TRunge::SetFirst
void SetFirst(double y[6]) const
set first point (position, momentum) s=0, phi=0
Definition: TRunge.cxx:956
TRunge::Nstep
unsigned Nstep(void) const
access to trajectory cache
Definition: TRunge.cxx:975

Referenced by approach_line(), and approach_point().

◆ Fly() [2/2]

unsigned TRunge::Fly ( const RkFitMaterial  material) const

make the trajectory in cache, return the number of step

◆ Fly_SC() [1/2]

unsigned TRunge::Fly_SC ( void  ) const

Definition at line 1097 of file TRunge.cxx.

1097 {//fly the particle track with stepsize control
1098 double y[6],dydx[6],yscal[6];
1099 unsigned Nstep;
1100 double step,stepdid,stepnext;
1101 unsigned i;
1102 double flightlength;
1103
1104 //yscal[0]=0.1; //1mm -> error = 1mm*EPS
1105 //yscal[1]=0.1; //1mm
1106 //yscal[2]=0.1; //1mm
1107 //yscal[3]=0.001; //1MeV
1108 //yscal[4]=0.001; //1MeV
1109 //yscal[5]=0.001; //1MeV
1110
1111 step=default_stepSize0;
1112 flightlength=0;
1113 SetFirst(y);
1114 for(Nstep=0;Nstep<TRunge_MAXstep;Nstep++){
1115 for(i=0;i<6;i++) _y[Nstep][i]=y[i]; // save each step
1116 //memcpy(_y[Nstep],y,sizeof(double)*6);
1117
1118 Function(y,dydx);
1119 //for(i=0;i<6;i++) yscal[i]=abs(y[i])+abs(dydx[i]*step)+TINY;
1120
1121 Propagate_QC(y,dydx,step,_eps,_yscal,stepdid,stepnext);
1122
1123 if( y[2]>160 || y[2]<-80 || (y[0]*y[0]+y[1]*y[1])>8100 ) break;
1124 //if position is out side of MDC, stop to fly
1125 // R>90cm, z<-80cm,160cm<z
1126
1127 flightlength += stepdid;
1128 if(flightlength>_maxflightlength) break;
1129
1130 _h[Nstep]=stepdid;
1131 //cout<<"TR:"<<Nstep<<":step="<<stepdid<<endl;
1132 if(stepnext<_stepSizeMin) step=_stepSizeMin;
1133 else if(stepnext>_stepSizeMax) step=_stepSizeMax;
1134 else step=stepnext;
1135 }
1136 _Nstep=Nstep+1;
1137 //cout<<"TR: Nstep="<<_Nstep<<endl;
1138 return(_Nstep);
1139}
default_stepSize0
const double default_stepSize0
Definition: TRunge.cxx:39
TRunge::Propagate_QC
void Propagate_QC(double y[6], double dydx[6], const double &steptry, const double &eps, const double yscal[6], double &stepdid, double &stepnext) const
Definition: TRunge.cxx:1048

Referenced by approach_line(), and approach_point().

◆ Fly_SC() [2/2]

unsigned TRunge::Fly_SC ( void  ) const

◆ Function() [1/2]

void TRunge::Function ( const double  y[6],
double  f[6] 
) const

Definition at line 916 of file TRunge.cxx.

916 {
917 // return the value of formula
918 // y[6] = (x,y,z,px,py,pz)
919 // f[6] = ( dx[3]/ds, dp[3]/ds )
920 // dx/ds = p/|p|
921 // dp/ds = e/|e| 1/alpha (p x B)/|p||B|
922 // alpha = 1/cB = 333.5640952/B[Tesla] [cm/(GeV/c)]
923 // const float Bx = _bfield->bx(y[0],y[1],y[2]);
924 // const float By = _bfield->by(y[0],y[1],y[2]);
925 // const float Bz = _bfield->bz(y[0],y[1],y[2]); //[kGauss]
926 double pmag;
927 double factor;
928 HepVector3D B;
929 HepPoint3D pos;
930 pos.setX((float)y[0]);
931 pos.setY((float)y[1]);
932 pos.setZ((float)y[2]);
933 //cout<<"TR::pos="<<pos[0]<<","<<pos[1]<<","<<pos[2]<<endl;
934// _bfield->fieldMap(pos,B);
935 //cout<<"TR::B="<<B[0]<<","<<B[1]<<","<<B[2]<<endl;
936 // const double Bmag = sqrt(Bx*Bx+By*By+Bz*Bz);
937 m_pmgnIMF->fieldVector(10*pos,B);
938 pmag = sqrt(y[3]*y[3]+y[4]*y[4]+y[5]*y[5]);
939 f[0] = y[3]/pmag; // p/|p|
940 f[1] = y[4]/pmag;
941 f[2] = y[5]/pmag;
942
943 // const factor = _charge/(alpha2/Bmag)/pmag/Bmag;
944 // factor = ((double)_charge)/alpha2/10.; //[(GeV/c)/(cm kG)]
945 // factor = ((double)_charge)/(333.564095/(-1000*(m_pmgnIMF->getReferField())))/10.; //[(GeV/c)/(cm kG)]
946 double Bmag=sqrt(B.x()*B.x()+B.y()*B.y()+B.z()*B.z());
947 factor = ((double)_charge)/(0.333564095);
948 // f[3] = factor*(f[1]*Bz-f[2]*By);
949 // f[4] = factor*(f[2]*Bx-f[0]*Bz);
950 // f[5] = factor*(f[0]*By-f[1]*Bx);
951 f[3] = factor*(f[1]*B.z()-f[2]*B.y());
952 f[4] = factor*(f[2]*B.x()-f[0]*B.z());
953 f[5] = factor*(f[0]*B.y()-f[1]*B.x());
954}
IMagneticFieldSvc::fieldVector
virtual StatusCode fieldVector(const HepGeom::Point3D< double > &xyz, HepGeom::Vector3D< double > &fvec) const =0

Referenced by approach_line(), Fly_SC(), Propagate(), Propagate1(), and Propagate_QC().

◆ Function() [2/2]

void TRunge::Function ( const double  y[6],
double  f[6] 
) const

◆ GetStep() [1/2]

int TRunge::GetStep ( unsigned  stepNum,
double &  step 
) const

Definition at line 985 of file TRunge.cxx.

985 {
986 if(stepNum>=_Nstep || stepNum>=TRunge_MAXstep) return(-1);
987 step=_h[stepNum];
988 return(0);
989}

◆ GetStep() [2/2]

int TRunge::GetStep ( unsigned  stepNum,
double &  step 
) const

◆ GetXP() [1/2]

int TRunge::GetXP ( unsigned  stepNum,
double  y[6] 
) const

Definition at line 979 of file TRunge.cxx.

979 {
980 if(stepNum>=_Nstep || stepNum>=TRunge_MAXstep) return(-1);
981
982 for(unsigned i=0;i<6;i++) y[i]=_y[stepNum][i];
983 return(0);
984}

◆ GetXP() [2/2]

int TRunge::GetXP ( unsigned  stepNum,
double  y[6] 
) const

◆ helix() [1/2]

Helix TRunge::helix ( void  ) const

returns helix class

Definition at line 238 of file TRunge.cxx.

238 {
239 return(Helix(_pivot,_a,_Ea));
240}

Referenced by RkFitCylinder::intersect(), intersect_cylinder(), intersect_xy_plane(), intersect_yz_plane(), intersect_zx_plane(), pivot(), and SetFlightLength().

◆ helix() [2/2]

Helix TRunge::helix ( void  ) const

returns helix class

◆ intersect_cylinder() [1/2]

double TRunge::intersect_cylinder ( double  r) const

Definition at line 1316 of file TRunge.cxx.

1316 {
1317 double m_rad = helix().radius();
1318 double l = helix().center().perp();
1319 double cosPhi = (m_rad * m_rad + l * l - r * r) / (2 * m_rad * l);
1320
1321 if(cosPhi < -1 || cosPhi > 1) return 0.;
1322
1323 double dPhi = helix().center().phi() - acos(cosPhi) - helix().phi0();
1324 if(dPhi < -M_PI) dPhi += 2 * M_PI;
1325
1326 return dPhi;
1327}
M_PI
#define M_PI
Definition: TConstant.h:4
Helix::center
const HepPoint3D & center(void) const
returns position of helix center(z = 0.);
Definition: InstallArea/include/TrackUtil/TrackUtil/Helix.h:203
Helix::radius
double radius(void) const
returns radious of helix.
Definition: InstallArea/include/TrackUtil/TrackUtil/Helix.h:215
TRunge::helix
Helix helix(void) const
returns helix class
Definition: TRunge.cxx:238

Referenced by RkFitCylinder::intersect(), intersect_yz_plane(), and intersect_zx_plane().

◆ intersect_cylinder() [2/2]

double TRunge::intersect_cylinder ( double  r) const

◆ intersect_xy_plane() [1/2]

double TRunge::intersect_xy_plane ( double  z) const

Definition at line 1362 of file TRunge.cxx.

1362 {
1363 if (helix().tanl() != 0 && helix().radius() != 0)
1364 return (helix().pivot().z() + helix().dz() - z) / (helix().radius() * helix().tanl());
1365 else return 0;
1366}
TRunge::pivot
const HepPoint3D & pivot(void) const
pivot position
Definition: TRunge.cxx:226
TRunge::tanl
double tanl(void) const
Definition: TRunge.cxx:222

Referenced by RkFitCylinder::intersect().

◆ intersect_xy_plane() [2/2]

double TRunge::intersect_xy_plane ( double  z) const

◆ intersect_yz_plane() [1/2]

double TRunge::intersect_yz_plane ( const HepTransform3D &  plane,
double  x 
) const

Definition at line 1345 of file TRunge.cxx.

1345 {
1346 HepPoint3D xc = plane * helix().center();
1347 double r = helix().radius();
1348 double d = r * r - (x - xc.x()) * (x - xc.x());
1349 if(d < 0) return 0;
1350
1351 double yy = xc.y();
1352 if(yy > 0) yy -= sqrt(d);
1353 else yy += sqrt(d);
1354
1355 double l = (plane.inverse() *
1356 HepPoint3D(x, yy, 0)).perp();
1357
1358 return intersect_cylinder(l);
1359}
TRunge::intersect_cylinder
double intersect_cylinder(double r) const
Definition: TRunge.cxx:1316

◆ intersect_yz_plane() [2/2]

double TRunge::intersect_yz_plane ( const HepTransform3D &  plane,
double  x 
) const

◆ intersect_zx_plane() [1/2]

double TRunge::intersect_zx_plane ( const HepTransform3D &  plane,
double  y 
) const

Definition at line 1328 of file TRunge.cxx.

1328 {
1329 HepPoint3D xc = plane * helix().center();
1330 double r = helix().radius();
1331 double d = r * r - (y - xc.y()) * (y - xc.y());
1332 if(d < 0) return 0;
1333
1334 double xx = xc.x();
1335 if(xx > 0) xx -= sqrt(d);
1336 else xx += sqrt(d);
1337
1338 double l = (plane.inverse() *
1339 HepPoint3D(xx, y, 0)).perp();
1340
1341 return intersect_cylinder(l);
1342}

◆ intersect_zx_plane() [2/2]

double TRunge::intersect_zx_plane ( const HepTransform3D &  plane,
double  y 
) const

◆ kappa() [1/2]

double TRunge::kappa ( void  ) const

Definition at line 214 of file TRunge.cxx.

214 {
215 return(_a[2]);
216}

◆ kappa() [2/2]

double TRunge::kappa ( void  ) const

◆ Mass() [1/4]

float TRunge::Mass ( float  mass)

set mass for tof calc.

Definition at line 347 of file TRunge.cxx.

347 {
348 _mass=mass;
349 _mass2=_mass*_mass;
350 return(_mass);
351}

◆ Mass() [2/4]

float TRunge::Mass ( float  )

set mass for tof calc.

◆ Mass() [3/4]

float TRunge::Mass ( void  ) const

return mass

Definition at line 279 of file TRunge.cxx.

279 {
280 return(_mass);
281}

◆ Mass() [4/4]

float TRunge::Mass ( void  ) const

return mass

◆ MaxFlightLength() [1/4]

double TRunge::MaxFlightLength ( double  length)

Definition at line 353 of file TRunge.cxx.

353 {
354 if(length>0) _maxflightlength=length;
355 else _maxflightlength=default_maxflightlength;
356 _Nstep=0;
357 return(_maxflightlength);
358}

◆ MaxFlightLength() [2/4]

double TRunge::MaxFlightLength ( double  )

◆ MaxFlightLength() [3/4]

double TRunge::MaxFlightLength ( void  ) const

return max flight length

Definition at line 283 of file TRunge.cxx.

283 {
284 return(_maxflightlength);
285}

◆ MaxFlightLength() [4/4]

double TRunge::MaxFlightLength ( void  ) const

return max flight length

◆ ndf() [1/2]

unsigned TRunge::ndf ( void  ) const

returns NDF

Definition at line 242 of file TRunge.cxx.

242 {
243 return _ndf;
244}

◆ ndf() [2/2]

unsigned TRunge::ndf ( void  ) const

returns NDF

◆ Nstep() [1/2]

unsigned TRunge::Nstep ( void  ) const

access to trajectory cache

Definition at line 975 of file TRunge.cxx.

975 {
976 return(_Nstep);
977}

Referenced by Fly(), and Fly_SC().

◆ Nstep() [2/2]

unsigned TRunge::Nstep ( void  ) const

access to trajectory cache

◆ objectType() [1/2]

unsigned TRunge::objectType ( void  ) const
inlinevirtual

returns object type

Reimplemented from TTrackBase.

Definition at line 256 of file InstallArea/include/TrkReco/TrkReco/TRunge.h.

256 {
257 return Runge;
258}

◆ objectType() [2/2]

unsigned TRunge::objectType ( void  ) const
virtual

returns object type

Reimplemented from TTrackBase.

◆ phi0() [1/2]

double TRunge::phi0 ( void  ) const

Definition at line 210 of file TRunge.cxx.

210 {
211 return(_a[1]);
212}

◆ phi0() [2/2]

double TRunge::phi0 ( void  ) const

◆ pivot() [1/4]

const HepPoint3D & TRunge::pivot ( const HepPoint3D &  newpivot)

set new pivot

!!!!! under construction !!!!! track parameter should be extracted after track propagation.

Definition at line 287 of file TRunge.cxx.

287 {
288 /// !!!!! under construction !!!!!
289 /// track parameter should be extracted after track propagation.
290 Helix tHelix(helix());
291 tHelix.pivot(newpivot);
292 _pivot=newpivot;
293 _a=tHelix.a();
294 _Ea=tHelix.Ea();
295 _Nstep=0;
296 return(_pivot);
297}

◆ pivot() [2/4]

const HepPoint3D & TRunge::pivot ( const HepPoint3D &  )

set new pivot

◆ pivot() [3/4]

const HepPoint3D & TRunge::pivot ( void  ) const

pivot position

Definition at line 226 of file TRunge.cxx.

226 {
227 return(_pivot);
228}

Referenced by intersect_xy_plane().

◆ pivot() [4/4]

const HepPoint3D & TRunge::pivot ( void  ) const

pivot position

◆ Propagate() [1/2]

void TRunge::Propagate ( double  y[6],
const double &  step,
const RkFitMaterial  material 
) const

propagate the track using 4th-order Runge-Kutta method

Definition at line 862 of file TRunge.cxx.

862 {
863 // y[6] = (x,y,z,px,py,pz)
864 double f1[6],f2[6],f3[6],f4[6],yt[6];
865 double hh;
866 double h6;
867 unsigned i;
868 const double mpi=0.13957;
869 double me=0.000511;
870 //eloss(step,&material, me,y,0);
871 hh = step*0.5;
872 //cout<<"TR:Pro hh="<<hh<<endl;
873 Function(y,f1);
874 for(i=0;i<6;i++) yt[i]=y[i]+hh*f1[i];
875 //{
876 // register double *a=y;
877 // register double *b=f1;
878 // register double *t=yt;
879 // register double *e=y+6;
880 // for(;a<e;a++,b++,t++) (*t)=(*a)+hh*(*b);
881 //}
882 Function(yt,f2);
883 for(i=0;i<6;i++) yt[i]=y[i]+hh*f2[i];
884 //{
885 // register double *a=y;
886 // register double *b=f2;
887 // register double *t=yt;
888 // register double *e=y+6;
889 // for(;a<e;a++,b++,t++) (*t)=(*a)+hh*(*b);
890 //}
891 Function(yt,f3);
892 for(i=0;i<6;i++) yt[i]=y[i]+step*f3[i];
893 //{
894 // register double *a=y;
895 // register double *b=f3;
896 // register double *t=yt;
897 // register double *e=y+6;
898 // for(;a<e;a++,b++,t++) (*t)=(*a)+step*(*b);
899 //}
900 Function(yt,f4);
901
902 h6 = step/6;
903 for(i=0;i<6;i++) y[i]+=h6*(f1[i]+2*f2[i]+2*f3[i]+f4[i]);
904 //{
905 // register double *a=f1;
906 // register double *b=f2;
907 // register double *c=f3;
908 // register double *d=f4;
909 // register double *t=y;
910 // register double *e=y+6;
911 // for(;t<e;a++,b++,c++,d++,t++)
912 // (*t)+=h6*((*a)+2*(*b)+2*(*c)+(*d));
913 //}
914}
mpi
const double mpi
Definition: Gam4pikp.cxx:47
me
const double me
Definition: PipiJpsi.cxx:48

Referenced by approach_line(), approach_point(), and Fly().

◆ Propagate() [2/2]

void TRunge::Propagate ( double  y[6],
const double &  step,
const RkFitMaterial  material 
) const

propagate the track using 4th-order Runge-Kutta method

◆ Propagate1() [1/2]

void TRunge::Propagate1 ( const double  y[6],
const double  dydx[6],
const double &  step,
double  yout[6] 
) const

Definition at line 991 of file TRunge.cxx.

992 {
993 // y[6] = (x,y,z,px,py,pz)
994 double f2[6],f3[6],f4[6],yt[6];
995 double hh;
996 double h6;
997 unsigned i;
998
999 hh = step*0.5;
1000 for(i=0;i<6;i++) yt[i]=y[i]+hh*dydx[i]; // 1st step
1001 //{
1002 // const register double *a=y;
1003 // const register double *b=dydx;
1004 // register double *t=yt;
1005 // const register double *e=y+6;
1006 // for(;a<e;a++,b++,t++) (*t)=(*a)+hh*(*b);
1007 //}
1008 Function(yt,f2); // 2nd step
1009 for(i=0;i<6;i++) yt[i]=y[i]+hh*f2[i];
1010 //{
1011 // const register double *a=y;
1012 // const register double *b=f2;
1013 // register double *t=yt;
1014 // const register double *e=y+6;
1015 // for(;a<e;a++,b++,t++) (*t)=(*a)+hh*(*b);
1016 //}
1017 Function(yt,f3); // 3rd step
1018 for(i=0;i<6;i++) yt[i]=y[i]+step*f3[i];
1019 //{
1020 // const register double *a=y;
1021 // const register double *b=f3;
1022 // register double *t=yt;
1023 // const register double *e=y+6;
1024 // for(;a<e;a++,b++,t++) (*t)=(*a)+step*(*b);
1025 //}
1026 Function(yt,f4); // 4th step
1027
1028 h6 = step/6;
1029 for(i=0;i<6;i++) yout[i]=y[i]+h6*(dydx[i]+2*f2[i]+2*f3[i]+f4[i]);
1030 //{
1031 // const register double *a=dydx;
1032 // const register double *b=f2;
1033 // const register double *c=f3;
1034 // const register double *d=f4;
1035 // const register double *e=y;
1036 // register double *t=yout;
1037 // const register double *s=yout+6;
1038 // for(;t<s;a++,b++,c++,d++,e++,t++)
1039 // (*t)=(*e)+h6*((*a)+2*(*b)+2*(*c)+(*d));
1040 //}
1041}

Referenced by Propagate_QC().

◆ Propagate1() [2/2]

void TRunge::Propagate1 ( const double  y[6],
const double  dydx[6],
const double &  step,
double  yout[6] 
) const

◆ Propagate_QC() [1/2]

void TRunge::Propagate_QC ( double  y[6],
double  dydx[6],
const double &  steptry,
const double &  eps,
const double  yscal[6],
double &  stepdid,
double &  stepnext 
) const

Definition at line 1048 of file TRunge.cxx.

1050 {
1051 //propagate with quality check, step will be scaled automatically
1052 double ysav[6],dysav[6],ytemp[6];
1053 double h,hh,errmax,temp;
1054 unsigned i;
1055
1056 for(i=0;i<6;i++) ysav[i]=y[i];
1057 //memcpy(ysav,y,sizeof(double)*6);
1058 for(i=0;i<6;i++) dysav[i]=dydx[i];
1059 //memcpy(dysav,dydx,sizeof(double)*6);
1060
1061 h=steptry;
1062 for(;;){
1063 hh=h*0.5; // 2 half step
1064 Propagate1(ysav,dysav,hh,ytemp);
1065 Function(ytemp,dydx);
1066 Propagate1(ytemp,dydx,hh,y);
1067 //if check step size
1068 Propagate1(ysav,dysav,h,ytemp); // 1 full step
1069 // error calculation
1070 errmax=0.0;
1071 for(i=0;i<6;i++){
1072 ytemp[i]=y[i]-ytemp[i];
1073 temp=fabs(ytemp[i]/yscal[i]);
1074 if(errmax < temp) errmax=temp;
1075 }
1076 //cout<<"TR: errmax="<<errmax<<endl;
1077 errmax/=eps;
1078 if(errmax <= 1.0){ // step O.K. calc. for next step
1079 stepdid=h;
1080 stepnext=(errmax>ERRCON ? SAFETY*h*exp(PGROW*log(errmax)) : 4.0*h);
1081 break;
1082 }
1083 h=SAFETY*h*exp(PSHRNK*log(errmax));
1084 }
1085 for(i=0;i<6;i++) y[i]+=ytemp[i]*FCOR;
1086 //{
1087 // register double *a=ytemp;
1088 // register double *t=y;
1089 // register double *e=y+6;
1090 // for(;t<e;a++,t++) (*t)+=(*a)*FCOR;
1091 //}
1092
1093}
exp
EvtComplex exp(const EvtComplex &c)
Definition: EvtComplex.hh:252
FCOR
#define FCOR
Definition: TRunge.cxx:1045
SAFETY
#define SAFETY
Definition: TRunge.cxx:1046
ERRCON
#define ERRCON
Definition: TRunge.cxx:1047
PSHRNK
#define PSHRNK
Definition: TRunge.cxx:1044
PGROW
#define PGROW
Definition: TRunge.cxx:1043
TRunge::Propagate1
void Propagate1(const double y[6], const double dydx[6], const double &step, double yout[6]) const
Definition: TRunge.cxx:991

Referenced by Fly_SC().

◆ Propagate_QC() [2/2]

void TRunge::Propagate_QC ( double  y[6],
double  dydx[6],
const double &  steptry,
const double &  eps,
const double  yscal[6],
double &  stepdid,
double &  stepnext 
) const

◆ reducedchi2() [1/2]

double TRunge::reducedchi2 ( void  ) const

returns reduced-chi2

Definition at line 250 of file TRunge.cxx.

250 {
251 if(_ndf==0){
252 std::cout<<"error at TRunge::reducedchi2 ndf=0"<<std::endl;
253 return 0;
254 }
255 return (_chi2/_ndf);
256}

◆ reducedchi2() [2/2]

double TRunge::reducedchi2 ( void  ) const

returns reduced-chi2

◆ SetFirst() [1/2]

void TRunge::SetFirst ( double  y[6]) const

set first point (position, momentum) s=0, phi=0

Definition at line 956 of file TRunge.cxx.

956 {
957 // y[6] = (x,y,z,px,py,pz)
958 const double cosPhi0=cos(_a[1]);
959 const double sinPhi0=sin(_a[1]);
960 const double invKappa=1/abs(_a[2]);
961
962 y[0]= _pivot.x()+_a[0]*cosPhi0; //x
963 y[1]= _pivot.y()+_a[0]*sinPhi0; //y
964 y[2]= _pivot.z()+_a[3]; //z
965 y[3]= -sinPhi0*invKappa; //px
966 y[4]= cosPhi0*invKappa; //py
967 y[5]= _a[4]*invKappa; //pz
968}
sin
double sin(const BesAngle a)
Definition: InstallArea/include/MdcGeom/MdcGeom/BesAngle.h:210
cos
double cos(const BesAngle a)
Definition: InstallArea/include/MdcGeom/MdcGeom/BesAngle.h:213

Referenced by Fly(), and Fly_SC().

◆ SetFirst() [2/2]

void TRunge::SetFirst ( double  y[6]) const

set first point (position, momentum) s=0, phi=0

◆ SetFlightLength() [1/2]

double TRunge::SetFlightLength ( void  )

set flight length using wire hits

Definition at line 1141 of file TRunge.cxx.

1141 {
1142
1143 // get a list of links to a wire hit
1144 const AList<TMLink>& cores=this->cores();
1145 //cout<<"TR:: cores.length="<<cores.length()<<endl;
1146
1147 const Helix hel(this->helix());
1148 double tanl=hel.tanl();
1149 //double curv=hel.curv();
1150 //double rho = Helix::ConstantAlpha / hel.kappa();
1151 // double rho = 333.564095 / hel.kappa();
1152 const double Bz = 1000*m_pmgnIMF->getReferField();
1153 double rho = 333.564095/(hel.kappa()*Bz);
1154
1155 // search max phi
1156 double phi_max=- DBL_MAX;
1157 double phi_min= DBL_MAX;
1158 // loop over link
1159 for(int j=0;j<cores.length();j++){
1160 TMLink& l=*cores[j];
1161 // fitting valid?
1162 const TMDCWire& wire=*l.wire();
1163 double Wz=0;
1164 //double mindist;
1165 HepPoint3D onWire;
1166 HepPoint3D dummy_bP;
1167 HepVector3D dummy_dir;
1168 Helix th(hel);
1169 for(int i=0;i<10;i++){
1170 wire.wirePosition(Wz,onWire,dummy_bP,dummy_dir);
1171 th.pivot(onWire);
1172 if(abs(th.dz())<0.1){ //<1mm
1173 //mindist=abs(hel.dr());
1174 break;
1175 }
1176 Wz+=th.dz();
1177 }
1178 //onWire is closest point , th.x() is closest point on trk
1179 //Wz is z position of closest point
1180 //Z->dphi
1181 /*
1182 // Calculate position (x,y,z) along helix.
1183 // x = x0 + dr * cos(phi0) + (alpha / kappa) * (cos(phi0) - cos(phi0+phi))
1184 // y = y0 + dr * sin(phi0) + (alpha / kappa) * (sin(phi0) - sin(phi0+phi))
1185 // z = z0 + dz - (alpha / kappa) * tan(lambda) * phi
1186 cout<<"TR:: Wz="<<Wz<<" tanl="<<tanl<<endl;
1187 double phi=( hel.pivot().z()+hel.dz()-Wz )/( curv*tanl );
1188 */
1189 //...Cal. dPhi to rotate...
1190 const HepPoint3D & xc = hel.center();
1191 const HepPoint3D & xt = hel.x();
1192 HepVector3D v0, v1;
1193 v0 = xt - xc;
1194 v1 = th.x() - xc;
1195 const double vCrs = v0.x() * v1.y() - v0.y() * v1.x();
1196 const double vDot = v0.x() * v1.x() + v0.y() * v1.y();
1197 double phi = atan2(vCrs, vDot);
1198
1199 double tz=hel.x(phi).z();
1200 //cout<<"TR:: Wz="<<Wz<<" tz="<<tz<<endl;
1201
1202 if(phi>phi_max) phi_max=phi;
1203 if(phi<phi_min) phi_min=phi;
1204 //cout<<"TR:: phi="<<phi<<endl;
1205 }
1206 //cout<<"TR:: phi_max="<<phi_max<<" phi_min="<<phi_min
1207 // <<" rho="<<rho<<" tanl="<<tanl<<endl;
1208 //phi->length, set max filght length
1209 //tanl*=1.2; // for mergin
1210 _maxflightlength=
1211 abs( rho*(phi_max-phi_min)*sqrt(1+tanl*tanl) )*1.2; //x1.1 mergin
1212
1213 return(_maxflightlength);
1214}
IMagneticFieldSvc::getReferField
virtual double getReferField()=0
TMDCWire::wirePosition
void wirePosition(float zPosition, HepPoint3D &xyPosition, HepPoint3D &backwardPosition, HepVector3D &direction) const
calculates position and direction vector with sag correction.
Definition: TMDCWire.cxx:541
TTrackBase::cores
const AList< TMLink > & cores(unsigned mask=0) const
returns a list of masked TMLinks for fit. 'mask' will be applied if mask is not 0.
Definition: TTrackBase.cxx:317

◆ SetFlightLength() [2/2]

double TRunge::SetFlightLength ( void  )

set flight length using wire hits

◆ StepSize() [1/4]

double TRunge::StepSize ( double  step)

set step size to calc. trajectory

Definition at line 320 of file TRunge.cxx.

320 {
321 _stepSize=step;
322 _Nstep=0;
323 return(_stepSize);
324}

◆ StepSize() [2/4]

double TRunge::StepSize ( double  )

set step size to calc. trajectory

◆ StepSize() [3/4]

double TRunge::StepSize ( void  ) const

returns step size

Definition at line 262 of file TRunge.cxx.

262 {
263 return(_stepSize);
264}

◆ StepSize() [4/4]

double TRunge::StepSize ( void  ) const

returns step size

◆ StepSizeMax() [1/4]

double TRunge::StepSizeMax ( double  step)

Definition at line 336 of file TRunge.cxx.

336 {
337 _stepSizeMax=step;
338 _Nstep=0;
339 return(_stepSizeMax);
340}

◆ StepSizeMax() [2/4]

double TRunge::StepSizeMax ( double  )

◆ StepSizeMax() [3/4]

double TRunge::StepSizeMax ( void  ) const

Definition at line 272 of file TRunge.cxx.

272 {
273 return(_stepSizeMax);
274}

◆ StepSizeMax() [4/4]

double TRunge::StepSizeMax ( void  ) const

◆ StepSizeMin() [1/4]

double TRunge::StepSizeMin ( double  step)

Definition at line 341 of file TRunge.cxx.

341 {
342 _stepSizeMin=step;
343 _Nstep=0;
344 return(_stepSizeMin);
345}

◆ StepSizeMin() [2/4]

double TRunge::StepSizeMin ( double  )

◆ StepSizeMin() [3/4]

double TRunge::StepSizeMin ( void  ) const

Definition at line 275 of file TRunge.cxx.

275 {
276 return(_stepSizeMin);
277}

◆ StepSizeMin() [4/4]

double TRunge::StepSizeMin ( void  ) const

◆ tanl() [1/2]

double TRunge::tanl ( void  ) const

Definition at line 222 of file TRunge.cxx.

222 {
223 return(_a[4]);
224}

Referenced by intersect_xy_plane(), and SetFlightLength().

◆ tanl() [2/2]

double TRunge::tanl ( void  ) const

◆ Yscal() [1/4]

const double * TRunge::Yscal ( const double *  )

set error parameters for fitting with step size control

◆ Yscal() [2/4]

const double * TRunge::Yscal ( const double *  )

set error parameters for fitting with step size control

◆ Yscal() [3/4]

const double * TRunge::Yscal ( void  ) const

return error parameters for fitting with step size control

Definition at line 266 of file TRunge.cxx.

266 {
267 return(_yscal);
268}

◆ Yscal() [4/4]

const double * TRunge::Yscal ( void  ) const

return error parameters for fitting with step size control

Friends And Related Function Documentation

◆ TRungeFitter

TRungeFitter
friend

Definition at line 224 of file InstallArea/include/TrkReco/TrkReco/TRunge.h.

The documentation for this class was generated from the following files: