View Javadoc

1   package org.lcsim.util.swim;
2   
3   import hep.physics.vec.Hep3Vector;
4   
5   import org.lcsim.event.LCIOParameters;
6   import org.lcsim.event.Track;
7   import org.lcsim.spacegeom.CartesianPoint;
8   import org.lcsim.spacegeom.CartesianVector;
9   import org.lcsim.spacegeom.SpacePoint;
10  import org.lcsim.spacegeom.SpaceVector;
11  
12  import static java.lang.Math.atan;
13  import static java.lang.Math.sqrt;
14  import static org.lcsim.constants.Constants.fieldConversion;
15  import static org.lcsim.event.LCIOParameters.ParameterName.phi0;
16  import static org.lcsim.event.LCIOParameters.ParameterName.omega;
17  import static org.lcsim.event.LCIOParameters.ParameterName.tanLambda;
18  
19  /**
20   * A simple helix swimmer for use in org.lcsim. Uses standard lcsim units Tesla,
21   * mm, GeV. This swimmer works for charged and neutral tracks.
22   * <p>
23   * For more info on swimming see <a href="doc-files/transport.pdf">this paper</a>
24   * by Paul Avery.
25   * 
26   * @author jstrube
27   * @version $Id: HelixSwimmer.java,v 1.22 2011/04/17 14:22:13 grefe Exp $
28   */
29  public class HelixSwimmer {
30      protected double field;
31      protected Trajectory _trajectory;
32      protected SpaceVector _momentum;
33      protected double _charge;
34  
35      /**
36       * Creates a new instance of HelixSwimmer
37       * 
38       * @param B
39       *                field strength in Tesla; uniform, solenoidal, directed
40       *                along z-axis
41       */
42      public HelixSwimmer(double B) {
43          field = B * fieldConversion;
44      }
45  
46      /**
47       * Sets parameters for helix swimmmer.
48       * 
49       * @param p
50       *                3-momentum (px,py,pz)
51       * @param r0
52       *                initial position(x0,y0,z0)
53       * @param iq
54       *                charge iq = q/|e| = +1/0/-1
55       */
56      public void setTrack(Hep3Vector p0, SpacePoint r0, int iq) {
57          SpaceVector p = new CartesianVector(p0.v());
58          double phi = Math.atan2(p.y(), p.x());
59          double lambda = Math.atan2(p.z(), p.rxy());
60  
61          if (iq != 0 && field != 0) {
62              double radius = p.rxy() / (iq * field);
63              _trajectory = new Helix(r0, radius, phi, lambda);
64          } else {
65              _trajectory = new Line(r0, phi, lambda);
66          }
67          _momentum = p;
68          _charge = iq;
69      }
70      
71      /**
72       * Sets parameters for helix swimmmer.
73       * 
74       * @param p
75       *                3-momentum (px,py,pz)
76       * @param r0
77       *                initial position(x0,y0,z0)
78       * @param q
79       *                charge q = q/|e| = +1/0/-1
80       */
81      public void setTrack(Hep3Vector p0, SpacePoint r0, double q) {
82          SpaceVector p = new CartesianVector(p0.v());
83          double phi = Math.atan2(p.y(), p.x());
84          double lambda = Math.atan2(p.z(), p.rxy());
85  
86          if (q != 0 && field != 0) {
87              double radius = p.rxy() / (q * field);
88              _trajectory = new Helix(r0, radius, phi, lambda);
89          } else {
90              _trajectory = new Line(r0, phi, lambda);
91          }
92          _momentum = p;
93          _charge = q;
94      }
95  
96      /**
97       * Sets parameters for helix swimmmer.
98       * 
99       * @param p
100      *                3-momentum (px,py,pz)
101      * @param r0
102      *                initial position(x0,y0,z0)
103      * @param iq
104      *                charge iq = q/|e| = +1/0/-1
105      * @deprecated in favor of
106      *             {@link setTrack(SpaceVector p, SpacePoint r0, int iq)}
107      *             because this method has an ambiguous signature
108      */
109     @Deprecated
110     public void setTrack(Hep3Vector p, Hep3Vector r0, int iq) {
111         double pt = sqrt(p.x() * p.x() + p.y() * p.y());
112         double phi = Math.atan2(p.y(), p.x());
113         double lambda = Math.atan2(p.z(), pt);
114 
115         if (iq != 0 && field != 0) {
116             double radius = pt / (iq * field);
117             _trajectory = new Helix(r0, radius, phi, lambda);
118         } else {
119             _trajectory = new Line(r0, phi, lambda);
120         }
121         _momentum = new CartesianVector(p.v());
122         _charge = iq;
123     }
124 
125     /**
126      * Sets the parameters for the helix swimmer. Uses the LCIOParameters class
127      * for conversion between track parameters and space and momentum
128      * representation
129      * 
130      * @param t
131      *                The track to approximate with a helix
132      */
133     public void setTrack(Track t) {
134         double pt = sqrt(t.getPX() * t.getPX() + t.getPY() * t.getPY());
135         LCIOParameters parameters = new LCIOParameters(t.getTrackParameters(), pt);
136 
137         SpacePoint ref = new CartesianPoint(t.getReferencePoint());
138         SpacePoint origin = LCIOParameters.Parameters2Position(parameters, ref);
139         _trajectory = new Helix(origin, 1 / parameters.get(omega), parameters
140                 .get(phi0), atan(parameters.get(tanLambda)));
141         // Hep3Vectors have too many shortcomings
142         _momentum = new CartesianVector(LCIOParameters.Parameters2Momentum(parameters).v());
143     }
144 
145     /**
146      * 
147      * @param alpha
148      * @return a {@link SpacePoint} at the length alpha from the origin along
149      *         the track
150      * @deprecated in favor of {@link getPointAtLength}
151      */
152     @Deprecated
153     public SpacePoint getPointAtDistance(double alpha) {
154         return getPointAtLength(alpha);
155     }
156 
157     /**
158      * Returns a SpacePoint at the length alpha from the origin along the track
159      * 
160      * @param alpha
161      * @return a {@link SpacePoint} at the length alpha from the origin along
162      *         the track
163      */
164     public SpacePoint getPointAtLength(double alpha) {
165         if (_trajectory == null) {
166             throw new RuntimeException("Trajectory not set");
167         }
168         return _trajectory.getPointAtDistance(alpha);
169     }
170 
171     public double getDistanceToRadius(double r) {
172         if (_trajectory == null) {
173             throw new RuntimeException("Trajectory not set");
174         }
175         return _trajectory.getDistanceToInfiniteCylinder(r);
176     }
177     public double getDistanceToPolyhedra(double r, int nsides)
178     {
179         if (_trajectory == null) {
180             throw new RuntimeException("Trajectory not set");
181         }
182         double s = _trajectory.getDistanceToInfiniteCylinder(r);
183         if(Double.isNaN(s))return s;
184         return ((Helix) (_trajectory)).getDistanceToPolyhedra(r,nsides);
185     }
186 
187     public double getDistanceToZ(double z) {
188         if (_trajectory == null) {
189             throw new RuntimeException("Trajectory not set");
190         }
191         double result = _trajectory.getDistanceToZPlane(z);
192         if (result < 0) {
193             result = _trajectory.getDistanceToZPlane(-z);
194         }
195         return result;
196     }
197 
198     public double getDistanceToCylinder(double r, double z) {
199         double x1 = getDistanceToRadius(r);
200         double x2 = getDistanceToZ(z);
201         return Double.isNaN(x1) ? x2 : Math.min(x1, x2);
202     }
203 
204     /**
205      * Returns the distance along the trajectory to get to the point of closest
206      * approach
207      * 
208      * @param point
209      *                The point to swim as close as possible to
210      * @return the length parameter by how much the trajectory has to be swum
211      */
212     public double getTrackLengthToPoint(Hep3Vector point) {
213         return _trajectory.getDistanceToPoint(point);
214     }
215 
216     public double getDistanceToPoint(Hep3Vector point) {
217         return getTrackLengthToPoint(point);
218     }
219     /**
220      * Returns the momentum on a point on the track at a distance from the
221      * origin
222      * 
223      * @param alpha
224      *                The 2D distance from the origin along the track
225      * @return The components of the momentum in a SpacePoint
226      */
227     public SpaceVector getMomentumAtLength(double alpha) {
228         // the trajectory can only return the unit direction of the momentum
229         SpaceVector unitDirection = _trajectory.getUnitTangentAtLength(alpha);
230         double magnitude = _momentum.rxy();
231         // System.out.println("HelixSwimmer: momentum.magnitude= "+magnitude);
232         return VectorArithmetic.multiply(unitDirection, magnitude);
233     }
234 
235     public Trajectory getTrajectory() {
236         return _trajectory;
237     }
238 }