CGEM BOSS 6.6.5.f
BESIII Offline Software System
Loading...
Searching...
No Matches
Public Member Functions | List of all members
TrkPocaXY Class Reference

#include <TrkPocaXY.h>

+ Inheritance diagram for TrkPocaXY:

Public Member Functions

 TrkPocaXY (const Trajectory &traj, double flt, const HepPoint3D &pt, double precision=1.0e-4)
 
 TrkPocaXY (const Trajectory &traj1, double flt1, const Trajectory &traj2, double flt2, double precision=1.0e-4)
 
 ~TrkPocaXY ()
 
double docaXY () const
 
double fltl1 () const
 
double fltl2 () const
 
 TrkPocaXY (const Trajectory &traj, double flt, const HepPoint3D &pt, double precision=1.0e-4)
 
 TrkPocaXY (const Trajectory &traj1, double flt1, const Trajectory &traj2, double flt2, double precision=1.0e-4)
 
 ~TrkPocaXY ()
 
double docaXY () const
 
double fltl1 () const
 
double fltl2 () const
 
- Public Member Functions inherited from TrkPocaBase
const TrkErrCode & status () const
 
double flt1 () const
 
double flt2 () const
 
double precision ()
 
const TrkErrCode & status () const
 
double flt1 () const
 
double flt2 () const
 
double precision ()
 

Additional Inherited Members

- Protected Member Functions inherited from TrkPocaBase
 TrkPocaBase (double flt1, double flt2, double precision)
 
 TrkPocaBase (double flt1, double precision)
 
virtual ~TrkPocaBase ()
 
void minimize (const Trajectory &traj1, double f1, const Trajectory &traj2, double f2)
 
void minimize (const Trajectory &traj1, double f1, const HepPoint3D &pt)
 
void stepTowardPoca (const Trajectory &traj1, const Trajectory &traj2)
 
void stepToPointPoca (const Trajectory &traj, const HepPoint3D &pt)
 
 TrkPocaBase (double flt1, double flt2, double precision)
 
 TrkPocaBase (double flt1, double precision)
 
virtual ~TrkPocaBase ()
 
void minimize (const Trajectory &traj1, double f1, const Trajectory &traj2, double f2)
 
void minimize (const Trajectory &traj1, double f1, const HepPoint3D &pt)
 
void stepTowardPoca (const Trajectory &traj1, const Trajectory &traj2)
 
void stepToPointPoca (const Trajectory &traj, const HepPoint3D &pt)
 
- Protected Attributes inherited from TrkPocaBase
double _precision
 
double _flt1
 
double _flt2
 
TrkErrCode _status
 
- Static Protected Attributes inherited from TrkPocaBase
static double _maxDist = 1.e7
 
static int _maxTry = 500
 
static double _extrapToler = 2.
 

Detailed Description

Definition at line 27 of file InstallArea/include/TrkBase/TrkBase/TrkPocaXY.h.

Constructor & Destructor Documentation

◆ TrkPocaXY() [1/4]

TrkPocaXY::TrkPocaXY ( const Trajectory &  traj,
double  flt,
const HepPoint3D &  pt,
double  precision = 1.0e-4 
)

Definition at line 23 of file TrkPocaXY.cxx.

25 : TrkPocaBase(flt,prec) , _docaxy(-9999.0) {
26// construct a line trajectory through the given point
27 Hep3Vector zaxis(0.0,0.0,1.0);
28// length is set to 200cm (shouldn't matter)
29// TrkLineTraj zline(pt, zaxis, 200.0);//yzhang del
30 TrkLineTraj zline(pt, zaxis, 2000.0);//yzhang change
31// create a 2-traj poca with this
32 double zlen(pt.z());
33 TrkPoca zlinepoca(traj,flt,zline,zlen,prec);
34// transfer over the data
35 _status = zlinepoca.status();
36 if(status().success()){
37 _flt1 = zlinepoca.flt1();
38 _flt2 = zlinepoca.flt2();
39 _docaxy = zlinepoca.doca(); // doca should be perpendicular to the zline so 2d by construction
40 }
41}

◆ TrkPocaXY() [2/4]

TrkPocaXY::TrkPocaXY ( const Trajectory &  traj1,
double  flt1,
const Trajectory &  traj2,
double  flt2,
double  precision = 1.0e-4 
)

Definition at line 43 of file TrkPocaXY.cxx.

46 : TrkPocaBase(flt1,flt2,prec) , _docaxy(-9999.0) {
47
48 // this traj-traj version starts by approximating the trejectories as
49 // either a line or a circle and treats separately the line-line,
50 // circle-circle, and line-circle cases.
51
52 _flt1 = flt1;
53 _flt2 = flt2;
54 double delta1(10*prec);
55 double delta2(10*prec);
56 unsigned niter(0);
57 static const unsigned maxiter(20);
58 while(niter < maxiter &&
59 ( delta1 > prec || delta2 > prec ) ){
60 // get positions and directions
61 HepPoint3D pos1 = traj1.position(_flt1);
62 HepPoint3D pos1b;
63 if(niter == 0) pos1b = pos1;
64 HepPoint3D pos2 = traj2.position(_flt2);
65 Hep3Vector dir1 = traj1.direction(_flt1);
66 Hep3Vector dir2 = traj2.direction(_flt2);
67 Hep3Vector dd1 = traj1.delDirect(_flt1);
68 Hep3Vector dd2 = traj2.delDirect(_flt2);
69 double curv1 = traj1.curvature(_flt1);
70 double curv2 = traj2.curvature(_flt2);
71 double m1, m2, q1, q2;
72 double r1, r2, xc1,xc2,yc1,yc2,sr1,sr2;
73 if(curv1 == 0) {
74 // convert to m*x+q
75 if(dir1[0] == 0){
76 m1 = 1.0e12;
77 }else{
78 m1 = dir1[1]/dir1[0];
79 }
80 q1 = pos1.y()-pos1.x()*m1;
81 }else{
82 // get circle parameters
83 r1 = (1- dir1[2]*dir1[2])/curv1;
84 sr1=r1;
85 if(dir1[0]*dd1[1] - dir1[1]*dd1[0] < 0) sr1 = -r1;
86 double cosphi1 = dir1[0]/sqrt(dir1[0]*dir1[0]+dir1[1]*dir1[1]);
87 double sinphi1 = cosphi1 * dir1[1]/dir1[0];
88 xc1 = pos1.x() - sr1 * sinphi1;
89 yc1 = pos1.y() + sr1 * cosphi1;
90 }
91 if(curv2 == 0) {
92 if(dir2[0] == 0){
93 m2 = 1.0e12;
94 }else{
95 m2 = dir2[1]/dir2[0];
96 }
97 q2 = pos2.y()-pos2.x()*m2;
98 }else{
99 r2 = (1-dir2[2]*dir2[2])/curv2;
100 sr2 = r2;
101 if(dir2[0]*dd2[1] - dir2[1]*dd2[0] < 0) sr2 = -r2;
102 double cosphi2 = dir2[0]/sqrt(dir2[0]*dir2[0]+dir2[1]*dir2[1]);
103 double sinphi2 = cosphi2 * dir2[1]/dir2[0];
104 xc2 = pos2.x() - sr2 * sinphi2;
105 yc2 = pos2.y() + sr2 * cosphi2;
106 }
107 double xint, yint, xint1, xint2, yint1, yint2, absdoca;
108 _docaxy = 0;
109 _status.setSuccess(3);
110
111 // First the line-line case
112
113 if(curv1==0 && curv2==0){
114 // do the intersection in 2d
115 interTwoLines(m1, q1, m2, q2, xint, yint);
116 }
117
118 // next do the two circle case
119
120 if(curv1 !=0 && curv2 !=0){
121 //There are four cases
122 double cdist = sqrt((xc1-xc2)*(xc1-xc2)+(yc1-yc2)*(yc1-yc2));
123 // a - coincident centers
124 if (fabs(cdist) < 1.e-12 ) {
125 //the algorithm fails because the points have all the same distance
126 _status.setFailure(12,
127 "TrkPocaXY:: the two circles have the same center...");
128 return;
129 }
130 // b - Actual intersection
131 if ( (fabs(r1-r2) < cdist) && (cdist < r1+r2 ) ) {
132 interTwoCircles(xc1,yc1,r1,xc2,yc2,r2,xint1,yint1,xint2,yint2);
133 double dist1 = (xint1-pos1b.x())*(xint1-pos1b.x()) +
134 (yint1-pos1b.y())*(yint1-pos1b.y());
135 double dist2 = (xint2-pos1b.x())*(xint2-pos1b.x()) +
136 (yint2-pos1b.y())*(yint2-pos1b.y());
137
138 //choose closest to begining of track1
139 if(dist1<dist2){
140 xint = xint1;
141 yint = yint1;
142 } else {
143 xint = xint2;
144 yint = yint2;
145 }
146 }
147
148 // c - nested circles and d - separated circles
149
150 if ( (fabs(r1-r2) > cdist) || // nested circles
151 ( cdist > (r1+r2) )) { // separated circles
152 // line going through the centers of the two circles
153
154 double m = (yc1-yc2)/(xc1-xc2); // y = m*x+q
155 double q = yc1 - m*xc1;
156
157 // intersection points between the line and the two circles
158
159 double xint1, yint1, xint2, yint2, zOfDCrossT1;
160
161 interLineCircle(m, q, xc1, yc1, r1, xint1, yint1, xint2, yint2);
162
163 double xint3, yint3, xint4, yint4 ;
164 interLineCircle(m, q, xc2, yc2, r2, xint3, yint3, xint4, yint4);
165 if (fabs(r1-r2) > cdist) { // nested circles
166 absdoca = fabs(r1-r2)-cdist;
167#ifdef MDCPATREC_DEBUG
168 cout << "ErrMsg(debugging) doing nested circles in TrkPocaXY " << endl;
169#endif
170
171 double dist1_3 = pow((xint1-xint3),2.) + pow((yint1-yint3),2.);
172 double dist2_4 = pow((xint2-xint4),2.) + pow((yint2-yint4),2.);
173
174 if(dist1_3<dist2_4) {
175 xint = 0.5*(xint1+xint3);
176 yint = 0.5*(yint1+yint3);
177 zOfDCrossT1 = (xint3-xint1)*dir1[1]-(yint3-yint1)*dir1[0];
178 } else {
179 xint = 0.5*(xint2+xint4);
180 yint = 0.5*(yint2+yint4);
181 zOfDCrossT1 = (xint4-xint2)*dir1[1]-(yint4-yint2)*dir1[0];
182 }
183 _docaxy = absdoca;
184 if( zOfDCrossT1 > 0) _docaxy = -absdoca;
185 }
186 if( cdist > (r1+r2) ) { //separated circles
187 absdoca = cdist - (r1+r2);
188#ifdef MDCPATREC_DEBUG
189 cout << "ErrMsg(debugging) doing separated circles in TrkPocaXY"<< endl;
190#endif
191
192 double dist2_3 = pow((xint2-xint3),2.) + pow((yint2-yint3),2.);
193 double dist1_4 = pow((xint1-xint4),2.) + pow((yint1-yint4),2.);
194 if (dist2_3<dist1_4){
195 xint = 0.5*(xint2+xint3);
196 yint = 0.5*(yint2+yint3);
197 zOfDCrossT1 = (xint3-xint1)*dir1[1]-(yint3-yint1)*dir1[0];
198 } else {
199 xint = 0.5*(xint1+xint4);
200 yint = 0.5*(yint1+yint4);
201 zOfDCrossT1 = (xint4-xint2)*dir1[1]-(yint4-yint2)*dir1[0];
202 }
203 _docaxy = absdoca;
204 if( zOfDCrossT1 > 0) _docaxy = -absdoca;
205 }
206 }
207 }
208
209 // Now do the line-circle case
210
211 if((curv1 == 0 && curv2 !=0) || (curv1 != 0 && curv2 == 0)) {
212 // distance between the line and the circle center
213 HepVector dirT;
214 double m,q,r,xc,yc, zOfDCrossT1;
215 if(curv1==0) {m=m1; q=q1; r=r2; xc=xc2; yc=yc2; dirT=dir2;}
216 else {m=m2; q=q2; r=r1; xc=xc1; yc=yc1; dirT=dir1;}
217
218 double dist = fabs((m*xc-yc+q)/sqrt(1+m*m));
219 if(dist <= r) {
220
221 // the intersection points
222
223 interLineCircle(m,q,xc,yc,r, xint1, yint1, xint2, yint2);
224 double dist1 = (xint1-pos1b.x())*(xint1-pos1b.x()) +
225 (yint1-pos1b.y())*(yint1-pos1b.y());
226 double dist2 = (xint2-pos1b.x())*(xint2-pos1b.x()) +
227 (yint2-pos1b.y())*(yint2-pos1b.y());
228 //choose closest to the beginning of track 1
229 if(dist1<dist2){
230 xint = xint1;
231 yint = yint1;
232 } else {
233 xint = xint2;
234 yint = yint2;
235 }
236 } else { // no intersection points
237#ifdef MDCPATREC_DEBUG
238 cout << "ErrMsg(debugging) doing separated line-circle in TrkPocaXY"<<endl;
239#endif
240
241 // line going through the circle center and perpendicular to traj1
242
243 double mperp = -1./m;
244 double qperp = yc - mperp*xc;
245
246 // intersection between this line and the two trajectories
247
248 interLineCircle(mperp, qperp, xc, yc, r, xint1,yint1,xint2,yint2);
249
250 double xint3,yint3;
251 interTwoLines(m, q, mperp, qperp, xint3, yint3);
252
253 double dist1_3 = pow((xint1-xint3),2.) + pow((yint1-yint3),2.);
254 double dist2_3 = pow((xint2-xint3),2.) + pow((yint2-yint3),2.);
255 if (dist1_3<dist2_3) {
256 xint = 0.5*(xint3 + xint1);
257 yint = 0.5*(yint3 + yint1);
258 absdoca = sqrt((xint3-xint1)*(xint3-xint1)+
259 (yint3-yint1)*(yint3-yint1));
260 zOfDCrossT1 = (xint3-xint1)*dirT[1]-(yint3-yint1)*dirT[0];
261 }
262 else {
263 xint = 0.5*(xint3 + xint2);
264 yint = 0.5*(yint3 + yint2);
265 absdoca = sqrt((xint3-xint2)*(xint3-xint2)+
266 (yint3-yint2)*(yint3-yint2));
267 zOfDCrossT1 = (xint3-xint2)*dirT[1]-(yint3-yint2)*dirT[0];
268 }
269 _docaxy = absdoca;
270 if( zOfDCrossT1 > 0) _docaxy = -absdoca;
271 }
272 }
273
274 // get the flight lengths for the intersection point
275
276 HepPoint3D intpt( xint, yint, 0.0);
277 TrkPocaXY poca1(traj1,_flt1,intpt,prec);
278 TrkPocaXY poca2(traj2,_flt2,intpt,prec);
279 _status = poca2.status();
280 if(poca1.status().success() && poca2.status().success()) {
281 delta1 = fabs(_flt1 - poca1.flt1());
282 delta2 = fabs(_flt2 - poca2.flt1());
283 _flt1 = poca1.flt1();
284 _flt2 = poca2.flt1();
285 }else{
286#ifdef MDCPATREC_WARNING
287 cout << "ErrMsg(warning)" << " poca fialure " << endl;
288#endif
289 if(poca1.status().failure()) _status=poca1.status();
290 break;
291 }
292 niter++;
293 }
294#ifdef MDCPATREC_DEBUG
295 cout <<"ErrMsg(debugging)" <<"TrkPocaXY iterations = " << niter-1
296 << " flight lengths " << _flt1 <<" " << _flt2 << endl
297 << " deltas " << delta1 <<" " << delta2 << endl;
298#endif
299 if(niter == maxiter){
300#ifdef MDCPATREC_ROUTINE
301 cout << "ErrMsg(routine)" << "TrkPocaXY:: did not converge" << endl;
302#endif
303 _status.setFailure(13, "TrkPocaXY:: did not converge");
304 }
305}
q
****INTEGER imax DOUBLE PRECISION m_pi *DOUBLE PRECISION m_amfin DOUBLE PRECISION m_Chfin DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_sinw2 DOUBLE PRECISION m_GFermi DOUBLE PRECISION m_MfinMin DOUBLE PRECISION m_ta2 INTEGER m_out INTEGER m_KeyFSR INTEGER m_KeyQCD *COMMON c_Semalib $ !copy of input $ !CMS energy $ !beam mass $ !final mass $ !beam charge $ !final charge $ !smallest final mass $ !Z mass $ !Z width $ !EW mixing angle $ !Gmu Fermi $ alphaQED at q
Definition: KKsem.h:33
Trajectory::position
virtual HepPoint3D position(double) const =0
Trajectory::delDirect
virtual Hep3Vector delDirect(double) const =0
Trajectory::direction
virtual Hep3Vector direction(double) const =0
Trajectory::curvature
virtual double curvature(double) const =0
TrkErrCode::setSuccess
void setSuccess(int i, const char *str=0)
Definition: InstallArea/include/TrkBase/TrkBase/TrkErrCode.h:84
TrkErrCode::setFailure
void setFailure(int i, const char *str=0)
Definition: InstallArea/include/TrkBase/TrkBase/TrkErrCode.h:79

◆ ~TrkPocaXY() [1/2]

TrkPocaXY::~TrkPocaXY ( )
inline

Definition at line 43 of file InstallArea/include/TrkBase/TrkBase/TrkPocaXY.h.

43{}

◆ TrkPocaXY() [3/4]

TrkPocaXY::TrkPocaXY ( const Trajectory &  traj,
double  flt,
const HepPoint3D &  pt,
double  precision = 1.0e-4 
)

◆ TrkPocaXY() [4/4]

TrkPocaXY::TrkPocaXY ( const Trajectory &  traj1,
double  flt1,
const Trajectory &  traj2,
double  flt2,
double  precision = 1.0e-4 
)

◆ ~TrkPocaXY() [2/2]

TrkPocaXY::~TrkPocaXY ( )
inline

Definition at line 43 of file Reconstruction/MdcPatRec/TrkBase/TrkBase-00-01-12/TrkBase/TrkPocaXY.h.

43{}

Member Function Documentation

◆ docaXY() [1/2]

double TrkPocaXY::docaXY ( ) const
inline

Definition at line 73 of file InstallArea/include/TrkBase/TrkBase/TrkPocaXY.h.

73{return _docaxy;}

◆ docaXY() [2/2]

double TrkPocaXY::docaXY ( ) const
inline

◆ fltl1() [1/2]

double TrkPocaXY::fltl1 ( ) const
inline

Definition at line 49 of file InstallArea/include/TrkBase/TrkBase/TrkPocaXY.h.

49{ return flt1(); }

◆ fltl1() [2/2]

double TrkPocaXY::fltl1 ( ) const
inline

Definition at line 49 of file Reconstruction/MdcPatRec/TrkBase/TrkBase-00-01-12/TrkBase/TrkPocaXY.h.

49{ return flt1(); }

◆ fltl2() [1/2]

double TrkPocaXY::fltl2 ( ) const
inline

Definition at line 50 of file InstallArea/include/TrkBase/TrkBase/TrkPocaXY.h.

50{ return flt2(); }

◆ fltl2() [2/2]

double TrkPocaXY::fltl2 ( ) const
inline

Definition at line 50 of file Reconstruction/MdcPatRec/TrkBase/TrkBase-00-01-12/TrkBase/TrkPocaXY.h.

50{ return flt2(); }
The documentation for this class was generated from the following files: