#include <Helix.h>
Inheritance diagram for KalmanFit::Helix:
Public Member Functions | |
| const HepVector & | a (const HepVector &newA) |
| sets helix parameters. | |
| const HepVector & | a (void) const |
| returns helix parameters. | |
| const HepVector & | a (const HepVector &newA) |
| sets helix parameters. | |
| const HepVector & | a (void) const |
| returns helix parameters. | |
| double | alpha (void) const |
| double | alpha (void) const |
| double | approach (HepPoint3D pfwd, HepPoint3D pbwd, bool doSagCorrection) const |
| double | approach (KalFitHitMdc &hit, bool doSagCorrection) const |
| double | approach (HepPoint3D pfwd, HepPoint3D pbwd, bool doSagCorrection) const |
| double | approach (KalFitHitMdc &hit, bool doSagCorrection) const |
| double | bFieldZ (void) const |
| double | bFieldZ (double) |
| sets/returns z componet of the magnetic field. | |
| double | bFieldZ (void) const |
| double | bFieldZ (double) |
| sets/returns z componet of the magnetic field. | |
| const HepPoint3D & | center (void) const |
| returns position of helix center(z = 0.); | |
| const HepPoint3D & | center (void) const |
| returns position of helix center(z = 0.); | |
| double | cosPhi0 (void) const |
| double | cosPhi0 (void) const |
| double | curv (void) const |
| double | curv (void) const |
| HepMatrix | del4MDelA (double phi, double mass) const |
| HepMatrix | del4MDelA (double phi, double mass) const |
| HepMatrix | del4MXDelA (double phi, double mass) const |
| HepMatrix | del4MXDelA (double phi, double mass) const |
| HepMatrix | delApDelA (const HepVector &ap) const |
| HepMatrix | delApDelA (const HepVector &ap) const |
| HepMatrix | delMDelA (double phi) const |
| HepMatrix | delMDelA (double phi) const |
| HepMatrix | delXDelA (double phi) const |
| HepMatrix | delXDelA (double phi) const |
| Hep3Vector | direction (double dPhi=0.) const |
| returns direction vector after rotating angle dPhi in phi direction. | |
| Hep3Vector | direction (double dPhi=0.) const |
| returns direction vector after rotating angle dPhi in phi direction. | |
| double | dr (void) const |
| returns an element of parameters. | |
| double | dr (void) const |
| returns an element of parameters. | |
| double | dz (void) const |
| double | dz (void) const |
| const HepSymMatrix & | Ea (const HepSymMatrix &newdA) |
| sets helix paramters and error matrix. | |
| const HepSymMatrix & | Ea (void) const |
| returns error matrix. | |
| const HepSymMatrix & | Ea (const HepSymMatrix &newdA) |
| sets helix paramters and error matrix. | |
| const HepSymMatrix & | Ea (void) const |
| returns error matrix. | |
| Helix (const HepPoint3D &position, const Hep3Vector &momentum, double charge) | |
| Constructor with position, momentum, and charge. | |
| Helix (const HepPoint3D &pivot, const HepVector &a) | |
| Constructor without error matrix. | |
| Helix (const HepPoint3D &pivot, const HepVector &a, const HepSymMatrix &Ea) | |
| Constructor with pivot, helix parameter a, and its error matrix. | |
| Helix (const HepPoint3D &position, const Hep3Vector &momentum, double charge) | |
| Constructor with position, momentum, and charge. | |
| Helix (const HepPoint3D &pivot, const HepVector &a) | |
| Constructor without error matrix. | |
| Helix (const HepPoint3D &pivot, const HepVector &a, const HepSymMatrix &Ea) | |
| Constructor with pivot, helix parameter a, and its error matrix. | |
| void | ignoreErrorMatrix (void) |
| unsets error matrix. Error calculations will be ignored after this function call until an error matrix be set again. 0 matrix will be return as a return value for error matrix when you call functions which returns an error matrix. | |
| void | ignoreErrorMatrix (void) |
| unsets error matrix. Error calculations will be ignored after this function call until an error matrix be set again. 0 matrix will be return as a return value for error matrix when you call functions which returns an error matrix. | |
| double | kappa (void) const |
| double | kappa (void) const |
| HepLorentzVector | momentum (double dPhi, double mass, HepPoint3D &x, HepSymMatrix &Emx) const |
| returns 4momentum vector after rotating angle dPhi in phi direction. | |
| HepLorentzVector | momentum (double dPhi, double mass, HepSymMatrix &Em) const |
| returns 4momentum vector after rotating angle dPhi in phi direction. | |
| HepLorentzVector | momentum (double dPhi, double mass) const |
| returns 4momentum vector after rotating angle dPhi in phi direction. | |
| Hep3Vector | momentum (double dPhi, HepSymMatrix &Em) const |
| returns momentum vector after rotating angle dPhi in phi direction. | |
| Hep3Vector | momentum (double dPhi=0.) const |
| returns momentum vector after rotating angle dPhi in phi direction. | |
| HepLorentzVector | momentum (double dPhi, double mass, HepPoint3D &x, HepSymMatrix &Emx) const |
| returns 4momentum vector after rotating angle dPhi in phi direction. | |
| HepLorentzVector | momentum (double dPhi, double mass, HepSymMatrix &Em) const |
| returns 4momentum vector after rotating angle dPhi in phi direction. | |
| HepLorentzVector | momentum (double dPhi, double mass) const |
| returns 4momentum vector after rotating angle dPhi in phi direction. | |
| Hep3Vector | momentum (double dPhi, HepSymMatrix &Em) const |
| returns momentum vector after rotating angle dPhi in phi direction. | |
| Hep3Vector | momentum (double dPhi=0.) const |
| returns momentum vector after rotating angle dPhi in phi direction. | |
| Helix & | operator= (const Helix &) |
| Copy operator. | |
| Helix & | operator= (const Helix &) |
| Copy operator. | |
| double | phi0 (void) const |
| double | phi0 (void) const |
| const HepPoint3D & | pivot (const HepPoint3D &newPivot) |
| sets pivot position. | |
| const HepPoint3D & | pivot (void) const |
| returns pivot position. | |
| const HepPoint3D & | pivot (const HepPoint3D &newPivot) |
| sets pivot position. | |
| const HepPoint3D & | pivot (void) const |
| returns pivot position. | |
| double | radius (void) const |
| returns radious of helix. | |
| double | radius (void) const |
| returns radious of helix. | |
| void | set (const HepPoint3D &pivot, const HepVector &a, const HepSymMatrix &Ea) |
| sets helix pivot position, parameters, and error matrix. | |
| void | set (const HepPoint3D &pivot, const HepVector &a, const HepSymMatrix &Ea) |
| sets helix pivot position, parameters, and error matrix. | |
| double | sinPhi0 (void) const |
| double | sinPhi0 (void) const |
| double | tanl (void) const |
| double | tanl (void) const |
| HepPoint3D | x (double dPhi, HepSymMatrix &Ex) const |
| returns position and convariance matrix(Ex) after rotation. | |
| double * | x (double dPhi, double p[3]) const |
| HepPoint3D | x (double dPhi=0.) const |
| returns position after rotating angle dPhi in phi direction. | |
| HepPoint3D | x (double dPhi, HepSymMatrix &Ex) const |
| returns position and convariance matrix(Ex) after rotation. | |
| double * | x (double dPhi, double p[3]) const |
| HepPoint3D | x (double dPhi=0.) const |
| returns position after rotating angle dPhi in phi direction. | |
| virtual | ~Helix () |
| Destructor. | |
| virtual | ~Helix () |
| Destructor. | |
Static Public Attributes | |
| const double | ConstantAlpha = 333.564095 |
| Constant alpha for uniform field. | |
Private Member Functions | |
| void | updateCache (void) |
| void | updateCache (void) |
Private Attributes | |
| HepVector | m_a |
| double | m_ac [5] |
| double | m_alpha |
| double | m_bField |
| HepPoint3D | m_center |
| double | m_cp |
| HepSymMatrix | m_Ea |
| bool | m_matrixValid |
| HepPoint3D | m_pivot |
| IMagneticFieldSvc * | m_pmgnIMF |
| IMagneticFieldSvc * | m_pmgnIMF |
| double | m_pt |
| double | m_r |
| double | m_sp |
|
||||||||||||||||
|
Constructor with pivot, helix parameter a, and its error matrix.
00052 : //m_bField(-10.0), 00053 //m_alpha(-333.564095), 00054 m_pivot(pivot), 00055 m_a(a), 00056 m_matrixValid(true), 00057 m_Ea(Ea) { 00058 StatusCode scmgn = Gaudi::svcLocator()->service("MagneticFieldSvc",m_pmgnIMF); 00059 if(scmgn!=StatusCode::SUCCESS) { 00060 std::cout<< "Unable to open Magnetic field service"<<std::endl; 00061 } 00062 m_bField = -10000*(m_pmgnIMF->getReferField()); 00063 m_alpha = 10000. / 2.99792458 / m_bField; 00064 // m_alpha = 10000. / 2.99792458 / m_bField; 00065 // m_alpha = 333.564095; 00066 updateCache(); 00067 }
|
|
||||||||||||
|
Constructor without error matrix.
00071 : //m_bField(-10.0), 00072 //m_alpha(-333.564095), 00073 m_pivot(pivot), 00074 m_a(a), 00075 m_matrixValid(false), 00076 m_Ea(HepSymMatrix(5,0)) { 00077 StatusCode scmgn = Gaudi::svcLocator()->service("MagneticFieldSvc",m_pmgnIMF); 00078 if(scmgn!=StatusCode::SUCCESS) { 00079 // log << MSG::ERROR << "Unable to open Magnetic field service"<<endreq; 00080 std::cout<< "Unable to open Magnetic field service"<<std::endl; 00081 } 00082 m_bField = -10000*(m_pmgnIMF->getReferField()); 00083 m_alpha = 10000. / 2.99792458 / m_bField; 00084 // m_alpha = 333.564095; 00085 //cout<<"MdcFastTrakAlg:: bField,alpha: "<<m_bField<<" , "<<m_alpha<<endl; 00086 updateCache(); 00087 }
|
|
||||||||||||||||
|
Constructor with position, momentum, and charge.
00092 : //m_bField(-10.0), 00093 //m_alpha(-333.564095), 00094 m_pivot(position), 00095 m_a(HepVector(5,0)), 00096 m_matrixValid(false), 00097 m_Ea(HepSymMatrix(5,0)) { 00098 StatusCode scmgn = Gaudi::svcLocator()->service("MagneticFieldSvc",m_pmgnIMF); 00099 if(scmgn!=StatusCode::SUCCESS) { 00100 // log << MSG::ERROR << "Unable to open Magnetic field service"<<endreq; 00101 std::cout<< "Unable to open Magnetic field service"<<std::endl; 00102 } 00103 m_bField = -10000*(m_pmgnIMF->getReferField()); 00104 m_alpha = 10000. / 2.99792458 / m_bField; 00105 00106 m_a[0] = 0.; 00107 m_a[1] = fmod(atan2(- momentum.x(), momentum.y()) 00108 + M_PI4, M_PI2); 00109 m_a[3] = 0.; 00110 double perp(momentum.perp()); 00111 if (perp != 0.0) { 00112 m_a[2] = charge / perp; 00113 m_a[4] = momentum.z() / perp; 00114 } 00115 else { 00116 m_a[2] = charge * (DBL_MAX); 00117 if (momentum.z() >= 0) { 00118 m_a[4] = (DBL_MAX); 00119 } else { 00120 m_a[4] = -(DBL_MAX); 00121 } 00122 } 00123 // m_alpha = 333.564095; 00124 updateCache(); 00125 }
|
|
|
Destructor.
00127 {
00128 }
|
|
||||||||||||||||
|
Constructor with pivot, helix parameter a, and its error matrix.
|
|
||||||||||||
|
Constructor without error matrix.
|
|
||||||||||||||||
|
Constructor with position, momentum, and charge.
|
|
|
Destructor.
|
|
|
sets helix parameters.
|
|
|
returns helix parameters.
|
|
|
sets helix parameters.
00266 {
00267 m_a = i;
00268 updateCache();
00269 return m_a;
00270 }
|
|
|
returns helix parameters.
00254 {
00255 return m_a;
00256 }
|
|
|
|
|
|
00292 {
00293
00294 return m_alpha;
00295 }
|
|
||||||||||||||||
|
|
|
||||||||||||
|
|
|
||||||||||||||||
|
00209 {
00210 // ...Cal. dPhi to rotate...
00211 // const TMDCWire & w = * link.wire();
00212 HepPoint3D positionOnWire, positionOnTrack;
00213 HepPoint3D pv = pivot();
00214 HepVector Va = a();
00215 HepSymMatrix Ma = Ea();
00216
00217 Helix _helix(pv, Va ,Ma);
00218 Hep3Vector Wire;
00219 Wire[0] = (pfwd - pbwd).x();
00220 Wire[1] = (pfwd - pbwd).y();
00221 Wire[2] = (pfwd - pbwd).z();
00222 // xyPosition(), returns middle position of a wire. z componet is 0.
00223 // w.xyPosition(wp);
00224 double wp[3];
00225 wp[0] = 0.5*(pfwd + pbwd).x();
00226 wp[1] = 0.5*(pfwd + pbwd).y();
00227 wp[2] = 0.;
00228 double wb[3];
00229 // w.backwardPosition(wb);
00230 wb[0] = pbwd.x();
00231 wb[1] = pbwd.y();
00232 wb[2] = pbwd.z();
00233 double v[3];
00234 v[0] = Wire.unit().x();
00235 v[1] = Wire.unit().y();
00236 v[2] = Wire.unit().z();
00237 // std::cout<<"Wire.unit() is "<<Wire.unit()<<std::endl;
00238
00239 // ...Sag correction...
00240 /* if (doSagCorrection) {
00241 HepVector3D dir = w.direction();
00242 HepPoint3D xw(wp[0], wp[1], wp[2]);
00243 HepPoint3D wireBackwardPosition(wb[0], wb[1], wb[2]);
00244 w.wirePosition(link.positionOnTrack().z(),
00245 xw,
00246 wireBackwardPosition,
00247 dir);
00248 v[0] = dir.x();
00249 v[1] = dir.y();
00250 v[2] = dir.z();
00251 wp[0] = xw.x();
00252 wp[1] = xw.y();
00253 wp[2] = xw.z();
00254 wb[0] = wireBackwardPosition.x();
00255 wb[1] = wireBackwardPosition.y();
00256 wb[2] = wireBackwardPosition.z();
00257 }
00258 */
00259 // ...Cal. dPhi to rotate...
00260 const HepPoint3D & xc = _helix.center();
00261 double xt[3];
00262 _helix.x(0., xt);
00263 double x0 = - xc.x();
00264 double y0 = - xc.y();
00265 double x1 = wp[0] + x0;
00266 double y1 = wp[1] + y0;
00267 x0 += xt[0];
00268 y0 += xt[1];
00269 //std::cout<<" x0 is: "<<x0<<" y0 is: "<<y0<<std::endl;
00270 //std::cout<<" x1 is: "<<x1<<" y1 is: "<<y1<<std::endl;
00271 //std::cout<<" xt[0] is: "<<xt[0]<<" xt[1] is: "<<xt[1]<<std::endl;
00272
00273 double dPhi = atan2(x0 * y1 - y0 * x1, x0 * x1 + y0 * y1);
00274 //std::cout<<" x0 * y1 - y0 * x1 is: "<<(x0 * y1 - y0 * x1)<<std::endl;
00275 //std::cout<<" x0 * x1 + y0 * y1 is: "<<(x0 * x1 + y0 * y1)<<std::endl;
00276 //std::cout<<" before loop dPhi is "<<dPhi<<std::endl;
00277 //...Setup...
00278 double kappa = _helix.kappa();
00279 double phi0 = _helix.phi0();
00280
00281 //...Axial case...
00282 /* if (!w.stereo()) {
00283 positionOnTrack = _helix.x(dPhi);
00284 HepPoint3D x(wp[0], wp[1], wp[2]);
00285 x.setZ(positionOnTrack.z());
00286 positionOnWire = x;
00287 //link.dPhi(dPhi);
00288 std::cout<<" in axial wire : positionOnTrack is "<<positionOnTrack
00289 <<" positionOnWire is "<<positionOnWire<<std::endl;
00290 return (positionOnTrack - positionOnWire).mag();
00291 }
00292 */
00293 double firstdfdphi = 0.;
00294 static bool first = true;
00295 if (first) {
00296 // extern BelleTupleManager * BASF_Histogram;
00297 // BelleTupleManager * m = BASF_Histogram;
00298 // h_nTrial = m->histogram("TTrack::approach nTrial", 100, 0., 100.);
00299 }
00300 //#endif
00301
00302 //...Stereo case...
00303 double rho = Helix::ConstantAlpha / kappa;
00304 double tanLambda = _helix.tanl();
00305 static HepPoint3D x;
00306 double t_x[3];
00307 double t_dXdPhi[3];
00308 const double convergence = 1.0e-5;
00309 double l;
00310 unsigned nTrial = 0;
00311 while (nTrial < 100) {
00312
00313 // x = link.positionOnTrack(_helix->x(dPhi));
00314 positionOnTrack = _helix.x(dPhi);
00315 x = _helix.x(dPhi);
00316 t_x[0] = x[0];
00317 t_x[1] = x[1];
00318 t_x[2] = x[2];
00319
00320 l = v[0] * t_x[0] + v[1] * t_x[1] + v[2] * t_x[2]
00321 - v[0] * wb[0] - v[1] * wb[1] - v[2] * wb[2];
00322
00323 double rcosPhi = rho * cos(phi0 + dPhi);
00324 double rsinPhi = rho * sin(phi0 + dPhi);
00325 t_dXdPhi[0] = rsinPhi;
00326 t_dXdPhi[1] = - rcosPhi;
00327 t_dXdPhi[2] = - rho * tanLambda;
00328
00329 //...f = d(Distance) / d phi...
00330 double t_d2Xd2Phi[2];
00331 t_d2Xd2Phi[0] = rcosPhi;
00332 t_d2Xd2Phi[1] = rsinPhi;
00333
00334 //...iw new...
00335 double n[3];
00336 n[0] = t_x[0] - wb[0];
00337 n[1] = t_x[1] - wb[1];
00338 n[2] = t_x[2] - wb[2];
00339
00340 double a[3];
00341 a[0] = n[0] - l * v[0];
00342 a[1] = n[1] - l * v[1];
00343 a[2] = n[2] - l * v[2];
00344 double dfdphi = a[0] * t_dXdPhi[0]
00345 + a[1] * t_dXdPhi[1]
00346 + a[2] * t_dXdPhi[2];
00347
00348 if (nTrial == 0) {
00349 // break;
00350 firstdfdphi = dfdphi;
00351 }
00352
00353 //...Check bad case...
00354 if (nTrial > 3) {
00355 // std::cout<<" BAD CASE!!, calculate approach ntrial = "<<nTrial<< endl;
00356 }
00357 //...Is it converged?...
00358 if (fabs(dfdphi) < convergence)
00359 break;
00360
00361 double dv = v[0] * t_dXdPhi[0]
00362 + v[1] * t_dXdPhi[1]
00363 + v[2] * t_dXdPhi[2];
00364 double t0 = t_dXdPhi[0] * t_dXdPhi[0]
00365 + t_dXdPhi[1] * t_dXdPhi[1]
00366 + t_dXdPhi[2] * t_dXdPhi[2];
00367 double d2fd2phi = t0 - dv * dv
00368 + a[0] * t_d2Xd2Phi[0]
00369 + a[1] * t_d2Xd2Phi[1];
00370 // + a[2] * t_d2Xd2Phi[2];
00371
00372 dPhi -= dfdphi / d2fd2phi;
00373 ++nTrial;
00374 }
00375 //std::cout<<" dPhi is: "<<dPhi<<std::endl;
00376 //...Cal. positions...
00377 positionOnWire[0] = wb[0] + l * v[0];
00378 positionOnWire[1] = wb[1] + l * v[1];
00379 positionOnWire[2] = wb[2] + l * v[2];
00380
00381 //std::cout<<"wb[0] is: "<<wb[0]<<" l is: "<<l<<" v[0] is: "<<v[0]<<std::endl;
00382 //std::cout<<"wb[1] is: "<<wb[1]<<" v[1] is: "<<v[1]<<std::endl;
00383 //std::cout<<"wb[2] is: "<<wb[2]<<" v[2] is: "<<v[2]<<std::endl;
00384
00385 //std::cout<<" positionOnTrack is "<<positionOnTrack
00386 // <<" positionOnWire is "<<positionOnWire<<std::endl;
00387
00388 return (positionOnTrack - positionOnWire).mag();
00389 // link.dPhi(dPhi);
00390 // return nTrial;
00391 }
|
|
||||||||||||
|
00016 {
00017 //...Cal. dPhi to rotate...
00018 //const TMDCWire & w = * link.wire();
00019 HepPoint3D positionOnWire, positionOnTrack;
00020 HepPoint3D pv = pivot();
00021 //std::cout<<"the track pivot in approach is "<<pv<<std::endl;
00022 HepVector Va = a();
00023 //std::cout<<"the track parameters is "<<Va<<std::endl;
00024 HepSymMatrix Ma = Ea();
00025 //std::cout<<"the error matrix is "<<Ma<<std::endl;
00026
00027 Helix _helix(pv, Va ,Ma);
00028 //_helix.pivot(IP);
00029 const KalFitWire& w = hit.wire();
00030 Hep3Vector Wire = w.fwd() - w.bck();
00031 //xyPosition(), returns middle position of a wire. z componet is 0.
00032 //w.xyPosition(wp);
00033 double wp[3];
00034 wp[0] = 0.5*(w.fwd() + w.bck()).x();
00035 wp[1] = 0.5*(w.fwd() + w.bck()).y();
00036 wp[2] = 0.;
00037 double wb[3];
00038 //w.backwardPosition(wb);
00039 wb[0] = w.bck().x();
00040 wb[1] = w.bck().y();
00041 wb[2] = w.bck().z();
00042 double v[3];
00043 v[0] = Wire.unit().x();
00044 v[1] = Wire.unit().y();
00045 v[2] = Wire.unit().z();
00046 //std::cout<<"Wire.unit() is "<<Wire.unit()<<std::endl;
00047
00048 //...Sag correction...
00049 /* if (doSagCorrection) {
00050 HepVector3D dir = w.direction();
00051 HepPoint3D xw(wp[0], wp[1], wp[2]);
00052 HepPoint3D wireBackwardPosition(wb[0], wb[1], wb[2]);
00053 w.wirePosition(link.positionOnTrack().z(),
00054 xw,
00055 wireBackwardPosition,
00056 dir);
00057 v[0] = dir.x();
00058 v[1] = dir.y();
00059 v[2] = dir.z();
00060 wp[0] = xw.x();
00061 wp[1] = xw.y();
00062 wp[2] = xw.z();
00063 wb[0] = wireBackwardPosition.x();
00064 wb[1] = wireBackwardPosition.y();
00065 wb[2] = wireBackwardPosition.z();
00066 }
00067 */
00068 //...Cal. dPhi to rotate...
00069 const HepPoint3D & xc = _helix.center();
00070
00071 //std::cout<<" helix center: "<<xc<<std::endl;
00072
00073 double xt[3]; _helix.x(0., xt);
00074 double x0 = - xc.x();
00075 double y0 = - xc.y();
00076 double x1 = wp[0] + x0;
00077 double y1 = wp[1] + y0;
00078 x0 += xt[0];
00079 y0 += xt[1];
00080 double dPhi = atan2(x0 * y1 - y0 * x1, x0 * x1 + y0 * y1);
00081
00082 //std::cout<<"dPhi is "<<dPhi<<std::endl;
00083
00084 //...Setup...
00085 double kappa = _helix.kappa();
00086 double phi0 = _helix.phi0();
00087
00088 //...Axial case...
00089 /* if (!w.stereo()) {
00090 positionOnTrack = _helix.x(dPhi);
00091 HepPoint3D x(wp[0], wp[1], wp[2]);
00092 x.setZ(positionOnTrack.z());
00093 positionOnWire = x;
00094 //link.dPhi(dPhi);
00095 std::cout<<" in axial wire : positionOnTrack is "<<positionOnTrack
00096 <<" positionOnWire is "<<positionOnWire<<std::endl;
00097 return (positionOnTrack - positionOnWire).mag();
00098 }
00099 */
00100 double firstdfdphi = 0.;
00101 static bool first = true;
00102 if (first) {
00103 // extern BelleTupleManager * BASF_Histogram;
00104 // BelleTupleManager * m = BASF_Histogram;
00105 // h_nTrial = m->histogram("TTrack::approach nTrial", 100, 0., 100.);
00106 }
00107 //#endif
00108
00109 //...Stereo case...
00110 double rho = Helix::ConstantAlpha / kappa;
00111 double tanLambda = _helix.tanl();
00112 static HepPoint3D x;
00113 double t_x[3];
00114 double t_dXdPhi[3];
00115 const double convergence = 1.0e-5;
00116 double l;
00117 unsigned nTrial = 0;
00118 while (nTrial < 100) {
00119
00120 // x = link.positionOnTrack(_helix->x(dPhi));
00121 positionOnTrack = _helix.x(dPhi);
00122 x = _helix.x(dPhi);
00123 t_x[0] = x[0];
00124 t_x[1] = x[1];
00125 t_x[2] = x[2];
00126
00127 l = v[0] * t_x[0] + v[1] * t_x[1] + v[2] * t_x[2]
00128 - v[0] * wb[0] - v[1] * wb[1] - v[2] * wb[2];
00129
00130 double rcosPhi = rho * cos(phi0 + dPhi);
00131 double rsinPhi = rho * sin(phi0 + dPhi);
00132 t_dXdPhi[0] = rsinPhi;
00133 t_dXdPhi[1] = - rcosPhi;
00134 t_dXdPhi[2] = - rho * tanLambda;
00135
00136 //...f = d(Distance) / d phi...
00137 double t_d2Xd2Phi[2];
00138 t_d2Xd2Phi[0] = rcosPhi;
00139 t_d2Xd2Phi[1] = rsinPhi;
00140
00141 //...iw new...
00142 double n[3];
00143 n[0] = t_x[0] - wb[0];
00144 n[1] = t_x[1] - wb[1];
00145 n[2] = t_x[2] - wb[2];
00146
00147 double a[3];
00148 a[0] = n[0] - l * v[0];
00149 a[1] = n[1] - l * v[1];
00150 a[2] = n[2] - l * v[2];
00151 double dfdphi = a[0] * t_dXdPhi[0]
00152 + a[1] * t_dXdPhi[1]
00153 + a[2] * t_dXdPhi[2];
00154
00155 //#ifdef TRKRECO_DEBUG
00156 if (nTrial == 0) {
00157 // break;
00158 firstdfdphi = dfdphi;
00159 }
00160
00161 //...Check bad case...
00162 if (nTrial > 3) {
00163 std::cout<<" bad case, calculate approach ntrial = "<<nTrial<< endl;
00164 }
00165 //#endif
00166
00167 //...Is it converged?...
00168 if (fabs(dfdphi) < convergence)
00169 break;
00170
00171 double dv = v[0] * t_dXdPhi[0]
00172 + v[1] * t_dXdPhi[1]
00173 + v[2] * t_dXdPhi[2];
00174 double t0 = t_dXdPhi[0] * t_dXdPhi[0]
00175 + t_dXdPhi[1] * t_dXdPhi[1]
00176 + t_dXdPhi[2] * t_dXdPhi[2];
00177 double d2fd2phi = t0 - dv * dv
00178 + a[0] * t_d2Xd2Phi[0]
00179 + a[1] * t_d2Xd2Phi[1];
00180 // + a[2] * t_d2Xd2Phi[2];
00181
00182 dPhi -= dfdphi / d2fd2phi;
00183
00184 // cout << "nTrial=" << nTrial << endl;
00185 // cout << "iw f,df,dphi=" << dfdphi << "," << d2fd2phi << "," << dPhi << endl;
00186
00187 ++nTrial;
00188 }
00189 //...Cal. positions...
00190 positionOnWire[0] = wb[0] + l * v[0];
00191 positionOnWire[1] = wb[1] + l * v[1];
00192 positionOnWire[2] = wb[2] + l * v[2];
00193
00194 //std::cout<<" positionOnTrack is "<<positionOnTrack
00195 // <<" positionOnWire is "<<positionOnWire<<std::endl;
00196 return (positionOnTrack - positionOnWire).mag();
00197
00198 //link.dPhi(dPhi);
00199
00200 // #ifdef TRKRECO_DEBUG
00201 // h_nTrial->accumulate((float) nTrial + .5);
00202 // #endif
00203 // return nTrial;
00204 }
|
|
|
|
|
|
sets/returns z componet of the magnetic field.
|
|
|
00300 {
00301 return m_bField;
00302 }
|
|
|
sets/returns z componet of the magnetic field.
00280 {
00281 m_bField = a;
00282 m_alpha = 10000. / 2.99792458 / m_bField;
00283
00284 //std::cout<<"m_alpha: "<<m_alpha<<std::endl;
00285 updateCache();
00286 return m_bField;
00287 }
|
|
|
returns position of helix center(z = 0.);
|
|
|
returns position of helix center(z = 0.);
00194 {
00195 return m_center;
00196 }
|
|
|
|
|
|
00312 {
00313 return m_cp;
00314 }
|
|
|
|
|
|
00248 {
00249 return m_r;
00250 }
|
|
||||||||||||
|
|
|
||||||||||||
|
00605 {
00606 //
00607 // Calculate Jacobian (@4m/@a)
00608 // Vector a is helix parameters and phi is internal parameter.
00609 // Vector 4m is 4 momentum.
00610 //
00611
00612 HepMatrix d4MDA(4,5,0);
00613
00614 double phi0 = m_ac[1];
00615 double cpa = m_ac[2];
00616 double tnl = m_ac[4];
00617
00618 double cosf0phi = cos(phi0+phi);
00619 double sinf0phi = sin(phi0+phi);
00620
00621 double rho;
00622 if(cpa != 0.)rho = 1./cpa;
00623 else rho = (DBL_MAX);
00624
00625 double charge = 1.;
00626 if(cpa < 0.)charge = -1.;
00627
00628 double E = sqrt(rho*rho*(1.+tnl*tnl)+mass*mass);
00629
00630 d4MDA[0][1] = -fabs(rho)*cosf0phi;
00631 d4MDA[0][2] = charge*rho*rho*sinf0phi;
00632
00633 d4MDA[1][1] = -fabs(rho)*sinf0phi;
00634 d4MDA[1][2] = -charge*rho*rho*cosf0phi;
00635
00636 d4MDA[2][2] = -charge*rho*rho*tnl;
00637 d4MDA[2][4] = fabs(rho);
00638
00639 if (cpa != 0.0 && E != 0.0) {
00640 d4MDA[3][2] = (-1.-tnl*tnl)/(cpa*cpa*cpa*E);
00641 d4MDA[3][4] = tnl/(cpa*cpa*E);
00642 } else {
00643 d4MDA[3][2] = (DBL_MAX);
00644 d4MDA[3][4] = (DBL_MAX);
00645 }
00646 return d4MDA;
00647 }
|
|
||||||||||||
|
|
|
||||||||||||
|
00651 {
00652 //
00653 // Calculate Jacobian (@4mx/@a)
00654 // Vector a is helix parameters and phi is internal parameter.
00655 // Vector 4xm is 4 momentum and position.
00656 //
00657
00658 HepMatrix d4MXDA(7,5,0);
00659
00660 const double & dr = m_ac[0];
00661 const double & phi0 = m_ac[1];
00662 const double & cpa = m_ac[2];
00663 const double & dz = m_ac[3];
00664 const double & tnl = m_ac[4];
00665
00666 double cosf0phi = cos(phi0+phi);
00667 double sinf0phi = sin(phi0+phi);
00668
00669 double rho;
00670 if(cpa != 0.)rho = 1./cpa;
00671 else rho = (DBL_MAX);
00672
00673 double charge = 1.;
00674 if(cpa < 0.)charge = -1.;
00675
00676 double E = sqrt(rho * rho * (1. + tnl * tnl) + mass * mass);
00677
00678 d4MXDA[0][1] = - fabs(rho) * cosf0phi;
00679 d4MXDA[0][2] = charge * rho * rho * sinf0phi;
00680
00681 d4MXDA[1][1] = - fabs(rho) * sinf0phi;
00682 d4MXDA[1][2] = - charge * rho * rho * cosf0phi;
00683
00684 d4MXDA[2][2] = - charge * rho * rho * tnl;
00685 d4MXDA[2][4] = fabs(rho);
00686
00687 if (cpa != 0.0 && E != 0.0) {
00688 d4MXDA[3][2] = (- 1. - tnl * tnl) / (cpa * cpa * cpa * E);
00689 d4MXDA[3][4] = tnl / (cpa * cpa * E);
00690 } else {
00691 d4MXDA[3][2] = (DBL_MAX);
00692 d4MXDA[3][4] = (DBL_MAX);
00693 }
00694
00695 d4MXDA[4][0] = m_cp;
00696 d4MXDA[4][1] = - dr * m_sp + m_r * (- m_sp + sinf0phi);
00697 if (cpa != 0.0) {
00698 d4MXDA[4][2] = - (m_r / cpa) * (m_cp - cosf0phi);
00699 } else {
00700 d4MXDA[4][2] = (DBL_MAX);
00701 }
00702
00703 d4MXDA[5][0] = m_sp;
00704 d4MXDA[5][1] = dr * m_cp + m_r * (m_cp - cosf0phi);
00705 if (cpa != 0.0) {
00706 d4MXDA[5][2] = - (m_r / cpa) * (m_sp - sinf0phi);
00707
00708 d4MXDA[6][2] = (m_r / cpa) * tnl * phi;
00709 } else {
00710 d4MXDA[5][2] = (DBL_MAX);
00711
00712 d4MXDA[6][2] = (DBL_MAX);
00713 }
00714
00715 d4MXDA[6][3] = 1.;
00716 d4MXDA[6][4] = - m_r * phi;
00717
00718 return d4MXDA;
00719 }
|
|
|
|
|
|
00446 {
00447 //
00448 // Calculate Jacobian (@ap/@a)
00449 // Vector ap is new helix parameters and a is old helix parameters.
00450 //
00451
00452 HepMatrix dApDA(5,5,0);
00453
00454 const double & dr = m_ac[0];
00455 const double & phi0 = m_ac[1];
00456 const double & cpa = m_ac[2];
00457 const double & dz = m_ac[3];
00458 const double & tnl = m_ac[4];
00459
00460 double drp = ap[0];
00461 double phi0p = ap[1];
00462 double cpap = ap[2];
00463 double dzp = ap[3];
00464 double tnlp = ap[4];
00465
00466 double rdr = m_r + dr;
00467 double rdrpr;
00468 if ((m_r + drp) != 0.0) {
00469 rdrpr = 1. / (m_r + drp);
00470 } else {
00471 rdrpr = (DBL_MAX);
00472 }
00473 // double csfd = cos(phi0)*cos(phi0p) + sin(phi0)*sin(phi0p);
00474 // double snfd = cos(phi0)*sin(phi0p) - sin(phi0)*cos(phi0p);
00475 double csfd = cos(phi0p - phi0);
00476 double snfd = sin(phi0p - phi0);
00477 double phid = fmod(phi0p - phi0 + M_PI8, M_PI2);
00478 if (phid > M_PI) phid = phid - M_PI2;
00479
00480 dApDA[0][0] = csfd;
00481 dApDA[0][1] = rdr*snfd;
00482 if(cpa!=0.0) {
00483 dApDA[0][2] = (m_r/cpa)*( 1.0 - csfd );
00484 } else {
00485 dApDA[0][2] = (DBL_MAX);
00486 }
00487
00488 dApDA[1][0] = - rdrpr*snfd;
00489 dApDA[1][1] = rdr*rdrpr*csfd;
00490 if(cpa!=0.0) {
00491 dApDA[1][2] = (m_r/cpa)*rdrpr*snfd;
00492 } else {
00493 dApDA[1][2] = (DBL_MAX);
00494 }
00495
00496 dApDA[2][2] = 1.0;
00497
00498 dApDA[3][0] = m_r*rdrpr*tnl*snfd;
00499 dApDA[3][1] = m_r*tnl*(1.0 - rdr*rdrpr*csfd);
00500 if(cpa!=0.0) {
00501 dApDA[3][2] = (m_r/cpa)*tnl*(phid - m_r*rdrpr*snfd);
00502 } else {
00503 dApDA[3][2] = (DBL_MAX);
00504 }
00505 dApDA[3][3] = 1.0;
00506 dApDA[3][4] = - m_r*phid;
00507
00508 dApDA[4][4] = 1.0;
00509
00510 return dApDA;
00511 }
|
|
|
|
|
|
00568 {
00569 //
00570 // Calculate Jacobian (@m/@a)
00571 // Vector a is helix parameters and phi is internal parameter.
00572 // Vector m is momentum.
00573 //
00574
00575 HepMatrix dMDA(3,5,0);
00576
00577 const double & phi0 = m_ac[1];
00578 const double & cpa = m_ac[2];
00579 const double & tnl = m_ac[4];
00580
00581 double cosf0phi = cos(phi0+phi);
00582 double sinf0phi = sin(phi0+phi);
00583
00584 double rho;
00585 if(cpa != 0.)rho = 1./cpa;
00586 else rho = (DBL_MAX);
00587
00588 double charge = 1.;
00589 if(cpa < 0.)charge = -1.;
00590
00591 dMDA[0][1] = -fabs(rho)*cosf0phi;
00592 dMDA[0][2] = charge*rho*rho*sinf0phi;
00593
00594 dMDA[1][1] = -fabs(rho)*sinf0phi;
00595 dMDA[1][2] = -charge*rho*rho*cosf0phi;
00596
00597 dMDA[2][2] = -charge*rho*rho*tnl;
00598 dMDA[2][4] = fabs(rho);
00599
00600 return dMDA;
00601 }
|
|
|
|
|
|
00514 {
00515 //
00516 // Calculate Jacobian (@x/@a)
00517 // Vector a is helix parameters and phi is internal parameter
00518 // which specifys the point to be calculated for Ex(phi).
00519 //
00520
00521 HepMatrix dXDA(3,5,0);
00522
00523 const double & dr = m_ac[0];
00524 const double & phi0 = m_ac[1];
00525 const double & cpa = m_ac[2];
00526 const double & dz = m_ac[3];
00527 const double & tnl = m_ac[4];
00528
00529 double cosf0phi = cos(phi0 + phi);
00530 double sinf0phi = sin(phi0 + phi);
00531
00532 dXDA[0][0] = m_cp;
00533 dXDA[0][1] = - dr * m_sp + m_r * (- m_sp + sinf0phi);
00534 if(cpa!=0.0) {
00535 dXDA[0][2] = - (m_r / cpa) * (m_cp - cosf0phi);
00536 } else {
00537 dXDA[0][2] = (DBL_MAX);
00538 }
00539 // dXDA[0][3] = 0.0;
00540 // dXDA[0][4] = 0.0;
00541
00542 dXDA[1][0] = m_sp;
00543 dXDA[1][1] = dr * m_cp + m_r * (m_cp - cosf0phi);
00544 if(cpa!=0.0) {
00545 dXDA[1][2] = - (m_r / cpa) * (m_sp - sinf0phi);
00546 } else {
00547 dXDA[1][2] = (DBL_MAX);
00548 }
00549 // dXDA[1][3] = 0.0;
00550 // dXDA[1][4] = 0.0;
00551
00552 // dXDA[2][0] = 0.0;
00553 // dXDA[2][1] = 0.0;
00554 if(cpa!=0.0) {
00555 dXDA[2][2] = (m_r / cpa) * tnl * phi;
00556 } else {
00557 dXDA[2][2] = (DBL_MAX);
00558 }
00559 dXDA[2][3] = 1.0;
00560 dXDA[2][4] = - m_r * phi;
00561
00562 return dXDA;
00563 }
|
|
|
returns direction vector after rotating angle dPhi in phi direction.
|
|
|
returns direction vector after rotating angle dPhi in phi direction.
00212 {
00213 return momentum(phi).unit();
00214 }
|
|
|
returns an element of parameters.
|
|
|
returns an element of parameters.
00218 {
00219 return m_ac[0];
00220 }
|
|
|
|
|
|
00236 {
00237 return m_ac[3];
00238 }
|
|
|
sets helix paramters and error matrix.
|
|
|
returns error matrix.
|
|
|
sets helix paramters and error matrix.
00274 {
00275 return m_Ea = i;
00276 }
|
|
|
returns error matrix.
00260 {
00261 return m_Ea;
00262 }
|
|
|
unsets error matrix. Error calculations will be ignored after this function call until an error matrix be set again. 0 matrix will be return as a return value for error matrix when you call functions which returns an error matrix.
|
|
|
unsets error matrix. Error calculations will be ignored after this function call until an error matrix be set again. 0 matrix will be return as a return value for error matrix when you call functions which returns an error matrix.
00722 {
00723 m_matrixValid = false;
00724 m_Ea *= 0.;
00725 }
|
|
|
|
|
|
00230 {
00231 return m_ac[2];
00232 }
|
|
||||||||||||||||||||
|
returns 4momentum vector after rotating angle dPhi in phi direction.
|
|
||||||||||||||||
|
returns 4momentum vector after rotating angle dPhi in phi direction.
|
|
||||||||||||
|
returns 4momentum vector after rotating angle dPhi in phi direction.
|
|
||||||||||||
|
returns momentum vector after rotating angle dPhi in phi direction.
|
|
|
returns momentum vector after rotating angle dPhi in phi direction.
|
|
||||||||||||||||||||
|
returns 4momentum vector after rotating angle dPhi in phi direction.
00280 {
00281 //
00282 // Calculate momentum.
00283 //
00284 // Pt = | 1/kappa | (GeV/c)
00285 //
00286 // Px = -Pt * sin(phi0 + phi)
00287 // Py = Pt * cos(phi0 + phi)
00288 // Pz = Pt * tan(lambda)
00289 //
00290 // E = sqrt( 1/kappa/kappa * (1+tan(lambda)*tan(lambda)) + mass*mass )
00291
00292 double pt = fabs(m_pt);
00293 double px = - pt * sin(m_ac[1] + phi);
00294 double py = pt * cos(m_ac[1] + phi);
00295 double pz = pt * m_ac[4];
00296 double E = sqrt(pt * pt * (1. + m_ac[4] * m_ac[4]) + mass * mass);
00297
00298 x.setX(m_pivot.x() + m_ac[0] * m_cp + m_r * (m_cp - cos(m_ac[1] + phi)));
00299 x.setY(m_pivot.y() + m_ac[0] * m_sp + m_r * (m_sp - sin(m_ac[1] + phi)));
00300 x.setZ(m_pivot.z() + m_ac[3] - m_r * m_ac[4] * phi);
00301
00302 if (m_matrixValid) Emx = m_Ea.similarity(del4MXDelA(phi,mass));
00303 else Emx = m_Ea;
00304
00305 return HepLorentzVector(px, py, pz, E);
00306 }
|
|
||||||||||||||||
|
returns 4momentum vector after rotating angle dPhi in phi direction.
00252 {
00253 //
00254 // Calculate momentum.
00255 //
00256 // Pt = | 1/kappa | (GeV/c)
00257 //
00258 // Px = -Pt * sin(phi0 + phi)
00259 // Py = Pt * cos(phi0 + phi)
00260 // Pz = Pt * tan(lambda)
00261 //
00262 // E = sqrt( 1/kappa/kappa * (1+tan(lambda)*tan(lambda)) + mass*mass )
00263
00264 double pt = fabs(m_pt);
00265 double px = - pt * sin(m_ac[1] + phi);
00266 double py = pt * cos(m_ac[1] + phi);
00267 double pz = pt * m_ac[4];
00268 double E = sqrt(pt*pt*(1.+m_ac[4]*m_ac[4])+mass*mass);
00269
00270 if (m_matrixValid) Em = m_Ea.similarity(del4MDelA(phi,mass));
00271 else Em = m_Ea;
00272
00273 return HepLorentzVector(px, py, pz, E);
00274 }
|
|
||||||||||||
|
returns 4momentum vector after rotating angle dPhi in phi direction.
00229 {
00230 //
00231 // Calculate momentum.
00232 //
00233 // Pt = | 1/kappa | (GeV/c)
00234 //
00235 // Px = -Pt * sin(phi0 + phi)
00236 // Py = Pt * cos(phi0 + phi)
00237 // Pz = Pt * tan(lambda)
00238 //
00239 // E = sqrt( 1/kappa/kappa * (1+tan(lambda)*tan(lambda)) + mass*mass )
00240
00241 double pt = fabs(m_pt);
00242 double px = - pt * sin(m_ac[1] + phi);
00243 double py = pt * cos(m_ac[1] + phi);
00244 double pz = pt * m_ac[4];
00245 double E = sqrt(pt*pt*(1.+m_ac[4]*m_ac[4])+mass*mass);
00246
00247 return HepLorentzVector(px, py, pz, E);
00248 }
|
|
||||||||||||
|
returns momentum vector after rotating angle dPhi in phi direction.
00206 {
00207 //
00208 // Calculate momentum.
00209 //
00210 // Pt = | 1/kappa | (GeV/c)
00211 //
00212 // Px = -Pt * sin(phi0 + phi)
00213 // Py = Pt * cos(phi0 + phi)
00214 // Pz = Pt * tan(lambda)
00215 //
00216
00217 double pt = fabs(m_pt);
00218 double px = - pt * sin(m_ac[1] + phi);
00219 double py = pt * cos(m_ac[1] + phi);
00220 double pz = pt * m_ac[4];
00221
00222 if (m_matrixValid) Em = m_Ea.similarity(delMDelA(phi));
00223 else Em = m_Ea;
00224
00225 return Hep3Vector(px, py, pz);
00226 }
|
|
|
returns momentum vector after rotating angle dPhi in phi direction.
00186 {
00187 //
00188 // Calculate momentum.
00189 //
00190 // Pt = | 1/kappa | (GeV/c)
00191 //
00192 // Px = -Pt * sin(phi0 + phi)
00193 // Py = Pt * cos(phi0 + phi)
00194 // Pz = Pt * tan(lambda)
00195 //
00196
00197 double pt = fabs(m_pt);
00198 double px = - pt * sin(m_ac[1] + phi);
00199 double py = pt * cos(m_ac[1] + phi);
00200 double pz = pt * m_ac[4];
00201
00202 return Hep3Vector(px, py, pz);
00203 }
|
|
|
Copy operator.
|
|
|
Copy operator.
00386 {
00387 if (this == & i) return * this;
00388
00389 m_bField = i.m_bField;
00390 m_alpha = i.m_alpha;
00391 m_pivot = i.m_pivot;
00392 m_a = i.m_a;
00393 m_Ea = i.m_Ea;
00394 m_matrixValid = i.m_matrixValid;
00395
00396 m_center = i.m_center;
00397 m_cp = i.m_cp;
00398 m_sp = i.m_sp;
00399 m_pt = i.m_pt;
00400 m_r = i.m_r;
00401 m_ac[0] = i.m_ac[0];
00402 m_ac[1] = i.m_ac[1];
00403 m_ac[2] = i.m_ac[2];
00404 m_ac[3] = i.m_ac[3];
00405 m_ac[4] = i.m_ac[4];
00406
00407 return * this;
00408 }
|
|
|
|
|
|
00224 {
00225 return m_ac[1];
00226 }
|
|
|
sets pivot position.
|
|
|
returns pivot position.
|
|
|
sets pivot position.
00310 {
00311
00312 //std::cout<<" in Helix::pivot:"<<std::endl;
00313 //std::cout<<" m_alpha: "<<m_alpha<<std::endl;
00314
00315 const double & dr = m_ac[0];
00316 const double & phi0 = m_ac[1];
00317 const double & kappa = m_ac[2];
00318 const double & dz = m_ac[3];
00319 const double & tanl = m_ac[4];
00320
00321 double rdr = dr + m_r;
00322 double phi = fmod(phi0 + M_PI4, M_PI2);
00323 double csf0 = cos(phi);
00324 double snf0 = (1. - csf0) * (1. + csf0);
00325 snf0 = sqrt((snf0 > 0.) ? snf0 : 0.);
00326 if(phi > M_PI) snf0 = - snf0;
00327
00328 double xc = m_pivot.x() + rdr * csf0;
00329 double yc = m_pivot.y() + rdr * snf0;
00330 double csf, snf;
00331 if(m_r != 0.0) {
00332 csf = (xc - newPivot.x()) / m_r;
00333 snf = (yc - newPivot.y()) / m_r;
00334 double anrm = sqrt(csf * csf + snf * snf);
00335 if(anrm != 0.0) {
00336 csf /= anrm;
00337 snf /= anrm;
00338 phi = atan2(snf, csf);
00339 } else {
00340 csf = 1.0;
00341 snf = 0.0;
00342 phi = 0.0;
00343 }
00344 } else {
00345 csf = 1.0;
00346 snf = 0.0;
00347 phi = 0.0;
00348 }
00349 double phid = fmod(phi - phi0 + M_PI8, M_PI2);
00350 if(phid > M_PI) phid = phid - M_PI2;
00351 double drp = (m_pivot.x() + dr * csf0 + m_r * (csf0 - csf) - newPivot.x())
00352 * csf
00353 + (m_pivot.y() + dr * snf0 + m_r * (snf0 - snf) - newPivot.y()) * snf;
00354 double dzp = m_pivot.z() + dz - m_r * tanl * phid - newPivot.z();
00355
00356 HepVector ap(5);
00357 ap[0] = drp;
00358 ap[1] = fmod(phi + M_PI4, M_PI2);
00359 ap[2] = kappa;
00360 ap[3] = dzp;
00361 ap[4] = tanl;
00362
00363 // if (m_matrixValid) m_Ea.assign(delApDelA(ap) * m_Ea * delApDelA(ap).T());
00364 if (m_matrixValid) m_Ea = m_Ea.similarity(delApDelA(ap));
00365
00366 m_a = ap;
00367 m_pivot = newPivot;
00368
00369 //...Are these needed?...iw...
00370 updateCache();
00371 return m_pivot;
00372 }
|
|
|
returns pivot position.
00200 {
00201 return m_pivot;
00202 }
|
|
|
returns radious of helix.
|
|
|
returns radious of helix.
00206 {
00207 return m_r;
00208 }
|
|
||||||||||||||||
|
sets helix pivot position, parameters, and error matrix.
|
|
||||||||||||||||
|
sets helix pivot position, parameters, and error matrix.
00377 {
00378 m_pivot = pivot;
00379 m_a = a;
00380 m_Ea = Ea;
00381 m_matrixValid = true;
00382 updateCache();
00383 }
|
|
|
|
|
|
00306 {
00307 return m_sp;
00308 }
|
|
|
|
|
|
00242 {
00243 return m_ac[4];
00244 }
|
|
|
|
|
|
00411 {
00412 //
00413 // Calculate Helix center( xc, yc ).
00414 //
00415 // xc = x0 + (dr + (alpha / kappa)) * cos(phi0) (cm)
00416 // yc = y0 + (dr + (alpha / kappa)) * sin(phi0) (cm)
00417 //
00418
00419 //std::cout<<" in updateCache, m_alpha: "<<m_alpha<<std::endl;
00420
00421 m_ac[0] = m_a[0];
00422 m_ac[1] = m_a[1];
00423 m_ac[2] = m_a[2];
00424 m_ac[3] = m_a[3];
00425 m_ac[4] = m_a[4];
00426
00427 m_cp = cos(m_ac[1]);
00428 m_sp = sin(m_ac[1]);
00429 if (m_ac[2] != 0.0) {
00430 m_pt = 1. / m_ac[2];
00431 m_r = m_alpha / m_ac[2];
00432 }
00433 else {
00434 m_pt = (DBL_MAX);
00435 m_r = (DBL_MAX);
00436 }
00437
00438 double x = m_pivot.x() + (m_ac[0] + m_r) * m_cp;
00439 double y = m_pivot.y() + (m_ac[0] + m_r) * m_sp;
00440 m_center.setX(x);
00441 m_center.setY(y);
00442 m_center.setZ(0.);
00443 }
|
|
||||||||||||
|
returns position and convariance matrix(Ex) after rotation.
|
|
||||||||||||
|
|
|
|
returns position after rotating angle dPhi in phi direction.
|
|
||||||||||||
|
returns position and convariance matrix(Ex) after rotation.
00165 {
00166 double x = m_pivot.x() + m_ac[0] * m_cp + m_r * (m_cp - cos(m_ac[1] +phi));
00167 double y = m_pivot.y() + m_ac[0] * m_sp + m_r * (m_sp - sin(m_ac[1] +phi));
00168 double z = m_pivot.z() + m_ac[3] - m_r * m_ac[4] * phi;
00169
00170 //
00171 // Calculate position error matrix.
00172 // Ex(phi) = (@x/@a)(Ea)(@x/@a)^T, phi is deflection angle to specify the
00173 // point to be calcualted.
00174 //
00175 // HepMatrix dXDA(3, 5, 0);
00176 // dXDA = delXDelA(phi);
00177 // Ex.assign(dXDA * m_Ea * dXDA.T());
00178
00179 if (m_matrixValid) Ex = m_Ea.similarity(delXDelA(phi));
00180 else Ex = m_Ea;
00181
00182 return HepPoint3D(x, y, z);
00183 }
|
|
||||||||||||
|
00148 {
00149 //
00150 // Calculate position (x,y,z) along helix.
00151 //
00152 // x = x0 + dr * cos(phi0) + (alpha / kappa) * (cos(phi0) - cos(phi0+phi))
00153 // y = y0 + dr * sin(phi0) + (alpha / kappa) * (sin(phi0) - sin(phi0+phi))
00154 // z = z0 + dz - (alpha / kappa) * tan(lambda) * phi
00155 //
00156
00157 p[0] = m_pivot.x() + m_ac[0] * m_cp + m_r * (m_cp - cos(m_ac[1] + phi));
00158 p[1] = m_pivot.y() + m_ac[0] * m_sp + m_r * (m_sp - sin(m_ac[1] + phi));
00159 p[2] = m_pivot.z() + m_ac[3] - m_r * m_ac[4] * phi;
00160
00161 return p;
00162 }
|
|
|
returns position after rotating angle dPhi in phi direction.
00131 {
00132 //
00133 // Calculate position (x,y,z) along helix.
00134 //
00135 // x = x0 + dr * cos(phi0) + (alpha / kappa) * (cos(phi0) - cos(phi0+phi))
00136 // y = y0 + dr * sin(phi0) + (alpha / kappa) * (sin(phi0) - sin(phi0+phi))
00137 // z = z0 + dz - (alpha / kappa) * tan(lambda) * phi
00138 //
00139
00140 double x = m_pivot.x() + m_ac[0] * m_cp + m_r * (m_cp - cos(m_ac[1] +phi));
00141 double y = m_pivot.y() + m_ac[0] * m_sp + m_r * (m_sp - sin(m_ac[1] +phi));
00142 double z = m_pivot.z() + m_ac[3] - m_r * m_ac[4] * phi;
00143
00144 return HepPoint3D(x, y, z);
00145 }
|
|
|
Constant alpha for uniform field.
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
1.3.9.1