ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitCommon/MathTools/src/Helix.cc
(Generate patch)

Comparing UserCode/MitCommon/MathTools/src/Helix.cc (file contents):
Revision 1.1 by loizides, Wed Sep 17 04:01:50 2008 UTC vs.
Revision 1.2 by loizides, Fri Mar 20 13:33:19 2009 UTC

# Line 1 | Line 1
1   // $Id$
2  
3   #include "MitCommon/MathTools/interface/Helix.h"
4 + #include <TSystem.h>
5  
6   using namespace mithep;
7  
8 + //--------------------------------------------------------------------------------------------------
9   Helix::Helix(const TVector3 &MomentumGev, const TVector3 &PositionCm,
10               double q, double BFieldTesla)
11   {
# Line 35 | Line 37 | Helix::Helix(const TVector3 &MomentumGev
37    }
38    // For the special case that we have no Helix, never happens for us .. right?
39    else {
40 <    std::cout << "MAJOR ERROR in HELIX !!!! STOP !!!!\n";
41 <  //  Line             auxLine(PositionCm,MomentumGev.unit());  
40 <  //  Z0             = auxLine.Z0();
41 <  //  Phi0           = auxLine.Phi0();
42 <  //  CotTheta       = auxLine.CotTheta();
43 <  //  D0             = auxLine.D0();
44 <  //  W              = 0.0;
45 <  //  
46 <  //  fCotTheta      = CotTheta;
47 <  //  fCurvature     = W/2;
48 <  //  fZ0            = Z0;
49 <  //  fD0            = D0;
50 <  //  fPhi0          = Phi0;
51 <  //  fIsStale       = 1;
52 <  //  fS             = -999.999;
53 <  //  fAa            = -999.999;
54 <  //  fSs            = -999.999;
55 <  //  fCc            = -999.999;
56 <  //  fSinPhi0       = 1.0;
57 <  //  fCosPhi0       = 1.0;
58 <  //  fSinTheta      = 1.0;
59 <  //  fCosTheta      = 1.0;
60 <  //  fVParameters   = 0;
61 <  //  fCenterIsValid = false;
62 <  //  fMx            = 0.0;
63 <  //  fMy            = 0.0;
40 >    std::cout << "MAJOR ERROR in HELIX!!!! Should not happen. STOP!!!!\n";
41 >    gSystem->Exit(123);
42    }
43   }
44  
45 + //--------------------------------------------------------------------------------------------------
46   bool Helix::operator == (const Helix &right) const
47   {
48    return
# Line 74 | Line 53 | bool Helix::operator == (const Helix &ri
53      fPhi0      == right.fPhi0;
54   }
55  
56 + //--------------------------------------------------------------------------------------------------
57   bool Helix::operator != (const Helix & right) const
58   {
59    return !((*this)==right);
60   }
61  
62 + //--------------------------------------------------------------------------------------------------
63   Helix::~Helix()
64   {
65    delete fVParameters;
66   }
67  
68 <
88 <
68 > //--------------------------------------------------------------------------------------------------
69   TVector3 Helix::Position(double s) const
70   {
71    fCacheSinesAndCosines(s);
# Line 102 | Line 82 | TVector3 Helix::Position(double s) const
82    }
83   }
84  
85 + //--------------------------------------------------------------------------------------------------
86   TVector3 Helix::Direction(double s) const
87   {
88    fCacheSinesAndCosines(s);
# Line 114 | Line 95 | TVector3 Helix::Direction(double s) cons
95    return TVector3(xtan,ytan,ztan);
96   }
97  
98 + //--------------------------------------------------------------------------------------------------
99   double Helix::PathLengthAtRhoEquals(double rho) const
100   {
101    return (SinTheta()?(L2DAtR(rho)/SinTheta()):0.0);
102   }
103  
104 + //--------------------------------------------------------------------------------------------------
105   double Helix::InverseRadius() const
106   {
107    return fCurvature*2.0;
108   }
109  
110 + //--------------------------------------------------------------------------------------------------
111   double Helix::Radius() const
112   {
113    return fabs(1.0/InverseRadius());
114   }
115  
116 + //--------------------------------------------------------------------------------------------------
117   SignedAngle Helix::TurningAngle(double s) const
118   {
119    return s/Radius();
120   }
121  
122 + //--------------------------------------------------------------------------------------------------
123   double Helix::Curvature() const
124   {
125    return fCurvature;
126   }
127  
128 + //--------------------------------------------------------------------------------------------------
129   double Helix::Helicity() const
130   {
131    return fCurvature>0 ? 1.0 : -1.0 ;
132   }
133  
134 + //--------------------------------------------------------------------------------------------------
135   double Helix::CotTheta() const
136   {
137    return fCotTheta;
138   }
139  
140 + //--------------------------------------------------------------------------------------------------
141   Angle Helix::Phi0() const
142   {
143    return fPhi0;
144   }
145  
146 + //--------------------------------------------------------------------------------------------------
147   double Helix::D0() const
148   {
149    return fD0;
150   }
151  
152 + //--------------------------------------------------------------------------------------------------
153   double Helix::SignLz() const
154   {
155    return (fD0>0) ? -1.0 : 1.0;
156   }
157  
158 + //--------------------------------------------------------------------------------------------------
159   double Helix::Z0() const
160   {
161    return fZ0;
162   }
163  
164 <
164 > //--------------------------------------------------------------------------------------------------
165   TVector3 Helix::SecondDerivative(double s) const
166   {
167    double phi1    = fPhi0+s*2.0*fCurvature*fSinTheta;
# Line 179 | Line 171 | TVector3 Helix::SecondDerivative(double
171    return TVector3(xsecond,ysecond,0.0);
172   }
173  
174 + //--------------------------------------------------------------------------------------------------
175   double Helix::SinPhi0() const
176   {
177    fRefreshCache();
178    return fSinPhi0;
179   }
180 +
181 + //--------------------------------------------------------------------------------------------------
182   double Helix::CosPhi0() const
183   {
184    fRefreshCache();
185    return fCosPhi0;
186   }
187 +
188 + //--------------------------------------------------------------------------------------------------
189   double Helix::SinTheta() const
190   {
191    fRefreshCache();
192    return fSinTheta;
193   }
194 +
195 + //--------------------------------------------------------------------------------------------------
196   double Helix::CosTheta() const
197   {
198    fRefreshCache();
199    return fCosTheta;
200   }
201  
202 + //--------------------------------------------------------------------------------------------------
203   Angle Helix::PhiAtR(double rho) const
204   {
205    double c = Curvature();
# Line 214 | Line 214 | Angle Helix::PhiAtR(double rho) const
214    return phi;
215   }
216  
217 + //--------------------------------------------------------------------------------------------------
218   double Helix::ZAtR(double rho) const
219   {
220    return fZ0 + CotTheta()*L2DAtR(rho);
221   }
222  
223 + //--------------------------------------------------------------------------------------------------
224   double Helix::L2DAtR(double rho) const {
225    double L2D;
226  
# Line 246 | Line 248 | double Helix::L2DAtR(double rho) const {
248    return L2D;
249   }
250  
251 + //--------------------------------------------------------------------------------------------------
252   double Helix::CosAlphaAtR(double rho) const {
253    double c   = Curvature();
254    double d   = D0();
# Line 265 | Line 268 | double Helix::CosAlphaAtR(double rho) co
268    double cosalpha = (x*dir_x + y*dir_y)/rho;
269    return cosalpha;
270   }
268 –
269 – //void Helix::Location(Trajectory::Location &loc, double s) const
270 – //{
271 – //  fCacheSinesAndCosines(s);
272 – //  double cP0sT = fCosPhi0*fSinTheta, sP0sT = fSinPhi0*fSinTheta;
273 – //  if (s && fCurvature) {
274 – //    loc.SetLocation(s,
275 – //                    (fCosPhi0*fSs-
276 – //                     fSinPhi0*(2.0*fCurvature*fD0+1.0-fCc))
277 – //                    /(2.0*fCurvature),
278 – //                    (fSinPhi0*fSs+
279 – //                     fCosPhi0*(2.0*fCurvature*fD0+1.0-fCc))
280 – //                    /(2.0*fCurvature),
281 – //                    s*fCosTheta + fZ0,
282 – //                    cP0sT*fCc-sP0sT*fSs,
283 – //                    cP0sT*fSs+sP0sT*fCc,
284 – //                    fCosTheta
285 – //                    );
286 – //  }
287 – //  else {
288 – //    loc.SetLocation(s,
289 – //                    -fD0*fSinPhi0+s*cP0sT,
290 – //                    fD0*fCosPhi0+s*sP0sT,
291 – //                    fZ0+s*fCosTheta,
292 – //                    cP0sT,sP0sT,fCosTheta
293 – //                    );
294 – //  }
295 – //}
296 –
297 – //Trajectory::Location* Helix::newIntersectionWith (const HepPlane3D& plane) const
298 – //{
299 – //  if (!SinTheta())                  // fastest way out of a screwy situation.
300 – //    return 0;
301 – //
302 – //  double deltaS=1.0,s=0.0;
303 – //  s = fabs(plane.d())/SinTheta();
304 – //  HepVector3D normal = plane.normal();
305 – //  Trajectory::Location *ploc = new Trajectory::Location ;
306 – //  for (int iteration=0;iteration<100;iteration++) {
307 – //    if (fabs(deltaS)>0.0001){
308 – //      Location(*ploc,s);
309 – //      deltaS = ((-(plane.distance(ploc->position())))
310 – //              /(normal.dot(ploc->direction())));
311 – //      s+=deltaS;
312 – //    }
313 – //    else {
314 – //      if (s>0)
315 – //      return ploc;
316 – //    }
317 – //  }
318 – //  delete ploc;
319 – //  return 0;
320 – //}

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines