5 modified files
lcsim/src/org/lcsim/util/swim
diff -u -r1.11 -r1.12
--- Helix.java 26 Aug 2005 17:37:04 -0000 1.11
+++ Helix.java 26 Aug 2005 21:02:26 -0000 1.12
@@ -1,16 +1,13 @@
package org.lcsim.util.swim;
-import hep.physics.vec.BasicHep3Vector;
-import hep.physics.vec.Hep3Vector;
-import hep.physics.vec.VecOp;
-
+import static java.lang.Math.asin;
+import static java.lang.Math.atan2;
import static java.lang.Math.cos;
import static java.lang.Math.sin;
import static java.lang.Math.sqrt;
-import static java.lang.Math.abs;
-import static java.lang.Math.atan2;
-import static java.lang.Math.asin;
-import static java.lang.Math.PI;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
@@ -19,7 +16,7 @@
* All quantities in this class are dimensionless. It has no dependencies
* except for Hep3Vector (which could easily be removed).
* @author tonyj
- * @version $Id: Helix.java,v 1.11 2005/08/26 17:37:04 jstrube Exp $
+ * @version $Id: Helix.java,v 1.12 2005/08/26 21:02:26 jstrube Exp $
*/
public class Helix implements Trajectory
{
@@ -80,12 +77,11 @@
/**
* Swims the helix along its trajectory to the
* point of closest approach to the given SpacePoint.
- * The equation to solve for s is an O(2) Taylor-expansion.
+ * @author jstrube
* @param point Point in Space to swim to.
* @return the length Parameter s
*/
- // FIXME works only in the first period apparently. Fix: translate along z into first period.
- public double getDistanceFromHelixToPoint(Hep3Vector point) {
+ public double getDistanceToPoint(Hep3Vector point) {
double tanLambda = sinLambda/cosLambda;
double px = radius*cosPhi;
double py = radius*sinPhi;
@@ -98,38 +94,58 @@
// first, the point needs to be translated into the first period
Hep3Vector xDiff = VecOp.sub(origin, point);
double zPos = xDiff.z();
- int addedPeriods = 0;
- // FIXME translate to first period
-// while (zPos > tanLambda*2*Math.PI) {
-// zPos -= tanLambda*2*Math.PI;
-// ++addedPeriods;
-// }
+ int addedQuarterPeriods = 0;
+ // these are two mutually exclusive cases and two while loops
+ // may not be the best way to express this
+ while (zPos > radius*tanLambda*Math.PI/2) {
+ zPos -= radius*tanLambda*Math.PI/2;
+ ++addedQuarterPeriods;
+ }
+ while (zPos < -radius*tanLambda*Math.PI/2) {
+ zPos += radius*tanLambda*Math.PI/2;
+ ++addedQuarterPeriods;
+ }
+ // Make sure the helix is in the right quadrant for the atan
+ if (addedQuarterPeriods % 2 != 0) ++addedQuarterPeriods;
+ xDiff = new BasicHep3Vector(xDiff.x(), xDiff.y(), zPos);
double factorA1 = pMag - pz*pz/pMag - (VecOp.dot(xDiff, pCrossB))*rho;
double factorA2 = (VecOp.dot(xDiff, p0) - xDiff.z()*pz)*rho;
-
- return Math.atan(factorA2/factorA1) / -rho;
+ return addedQuarterPeriods*(radius/cosLambda*Math.PI/2) + Math.atan(factorA2/factorA1) / -rho;
}
- /**
- * Calculates the <em>signed</em> distance in mm between the Helix and an arbitrary point in Space
- * @param point the point in space to calculate the distance to
- * @return the distance in mm between the point and the helix at the point of closest approach
- */
- public double getSignedClosestDifferenceToPoint(Hep3Vector point) {
- double tanLambda = sinLambda/cosLambda;
- double px = radius*cosPhi;
- double py = radius*sinPhi;
- double pz = radius*tanLambda;
- double rho = -cosLambda/radius;
- Hep3Vector pCrossB = new BasicHep3Vector(py, -px, 0);
- Hep3Vector xDiff = VecOp.sub(origin, point);
- double pMag = sqrt(px*px + py*py + pz*pz);
-
- double numerator = (-2*VecOp.dot(xDiff, pCrossB) + pMag*rho*(xDiff.magnitudeSquared()-xDiff.z()*xDiff.z()))/radius;
- double denominator = 1 + sqrt(1-2*rho*pMag*VecOp.dot(xDiff, pCrossB)/radius/radius + pMag*pMag*rho*rho*(xDiff.magnitudeSquared()-xDiff.z()*xDiff.z())/radius/radius);
- return numerator/denominator;
- }
+// /**
+// * Calculates the <em>signed</em> distance in mm between the Helix and an arbitrary point in Space
+// * @param point the point in space to calculate the distance to
+// * @return the distance in mm between the point and the helix at the point of closest approach
+// */
+// public double getSignedClosestDifferenceToPoint(Hep3Vector point) {
+// double tanLambda = sinLambda/cosLambda;
+// double px = radius*cosPhi;
+// double py = radius*sinPhi;
+// double pz = radius*tanLambda;
+// double rho = -cosLambda/radius;
+// Hep3Vector pCrossB = new BasicHep3Vector(py, -px, 0);
+// Hep3Vector xDiff = VecOp.sub(origin, point);
+// double pMag = sqrt(px*px + py*py + pz*pz);
+//
+// // translate along Z because algorithm can handle only numbers in the
+// // first quadrant
+// double zPos = xDiff.z();
+// // these are two mutually exclusive cases and two while loops
+// // may not be the best way to express this
+// while (zPos > radius*tanLambda*Math.PI/4) {
+// zPos -= radius*tanLambda*Math.PI/4;
+// }
+// while (zPos < -radius*tanLambda*Math.PI/4) {
+// zPos += radius*tanLambda*Math.PI/4;
+// }
+// xDiff = new BasicHep3Vector(xDiff.x(), xDiff.y(), zPos);
+//
+// double numerator = (-2*VecOp.dot(xDiff, pCrossB) + pMag*rho*(xDiff.magnitudeSquared()-xDiff.z()*xDiff.z()))/radius;
+// double denominator = 1 + sqrt(1-2*rho*pMag*VecOp.dot(xDiff, pCrossB)/radius/radius + pMag*pMag*rho*rho*(xDiff.magnitudeSquared()-xDiff.z()*xDiff.z())/radius/radius);
+// return numerator/denominator;
+// }
// Returns the "momentum" at the length s from the starting point.
// This uses the definition in http://www.phys.ufl.edu/~avery/fitting/transport.pdf
lcsim/src/org/lcsim/util/swim
diff -u -r1.6 -r1.7
--- HelixSwimmer.java 23 Aug 2005 05:41:07 -0000 1.6
+++ HelixSwimmer.java 26 Aug 2005 21:02:26 -0000 1.7
@@ -12,7 +12,7 @@
* A simple helix smimmer for use in org.lcsim. Uses standard lcsim units
* Tesla, mm, GeV. This swimmer works for charged and neutral tracks.
* @author tonyj
- * @version $Id: HelixSwimmer.java,v 1.6 2005/08/23 05:41:07 jstrube Exp $
+ * @version $Id: HelixSwimmer.java,v 1.7 2005/08/26 21:02:26 jstrube Exp $
*/
public class HelixSwimmer
{
@@ -90,6 +90,6 @@
}
public double getDistanceToPoint(Hep3Vector point) {
- return trajectory.getDistanceFromHelixToPoint(point);
+ return trajectory.getDistanceToPoint(point);
}
}
\ No newline at end of file
lcsim/src/org/lcsim/util/swim
diff -u -r1.3 -r1.4
--- Line.java 23 Aug 2005 05:41:07 -0000 1.3
+++ Line.java 26 Aug 2005 21:02:26 -0000 1.4
@@ -6,7 +6,7 @@
/**
* A straight line
* @author tonyj
- * @version $Id: Line.java,v 1.3 2005/08/23 05:41:07 jstrube Exp $
+ * @version $Id: Line.java,v 1.4 2005/08/26 21:02:26 jstrube Exp $
*/
public class Line implements Trajectory
{
@@ -46,7 +46,7 @@
}
// FIXME implementation
- public double getDistanceFromHelixToPoint(Hep3Vector point) {
+ public double getDistanceToPoint(Hep3Vector point) {
return Double.NaN;
}
}
lcsim/src/org/lcsim/util/swim
diff -u -r1.4 -r1.5
--- Trajectory.java 23 Aug 2005 05:41:07 -0000 1.4
+++ Trajectory.java 26 Aug 2005 21:02:26 -0000 1.5
@@ -5,7 +5,7 @@
/**
* A particle trajectory (either a Helix or a Line)
* @author tonyj
- * @version $Id: Trajectory.java,v 1.4 2005/08/23 05:41:07 jstrube Exp $
+ * @version $Id: Trajectory.java,v 1.5 2005/08/26 21:02:26 jstrube Exp $
*/
public interface Trajectory
{
@@ -33,5 +33,5 @@
* @return the distance of closest approach in mm
*/
// FIXME this should be a point rather than a vector
- public double getDistanceFromHelixToPoint(Hep3Vector point);
+ public double getDistanceToPoint(Hep3Vector point);
}
lcsim/test/org/lcsim/util/swim
diff -u -r1.5 -r1.6
--- HelixTest.java 26 Aug 2005 08:23:26 -0000 1.5
+++ HelixTest.java 26 Aug 2005 21:02:27 -0000 1.6
@@ -1,5 +1,7 @@
package org.lcsim.util.swim;
+import static java.lang.Math.cos;
+import static java.lang.Math.sin;
import hep.aida.ICloud2D;
import junit.framework.*;
import hep.physics.vec.BasicHep3Vector;
@@ -10,7 +12,7 @@
/**
*
* @author tonyj
- * @version $Id: HelixTest.java,v 1.5 2005/08/26 08:23:26 jstrube Exp $
+ * @version $Id: HelixTest.java,v 1.6 2005/08/26 21:02:27 jstrube Exp $
*/
public class HelixTest extends TestCase
{
@@ -41,15 +43,8 @@
assertTrue(Double.isNaN(circle.getDistanceToInfiniteCylinder(2)));
assertTrue(Double.isNaN(circle.getDistanceToInfiniteCylinder(0)));
assertEquals(0,circle.getDistanceToInfiniteCylinder(1),1e-14);
-
-// System.err.println(circle.getDistanceFromHelixToPoint(new BasicHep3Vector(1, 1, 1)));
-// System.err.println(circle.getDistanceFromHelixToPoint(new BasicHep3Vector(-.5, 0,0)));
-// System.err.println(circle.getDistanceFromHelixToPoint(new BasicHep3Vector(.5, 0,0)));
-// System.err.println(circle.getDistanceFromHelixToPoint(new BasicHep3Vector(0,.5,0)));
-// System.err.println(circle.getDistanceFromHelixToPoint(new BasicHep3Vector(0,-.5,0)));
-// System.err.println(circle.getDistanceFromHelixToPoint(new BasicHep3Vector(0,.5,1)));
-// System.err.println(circle.getDistanceFromHelixToPoint(new BasicHep3Vector(0,-.5,1)));
}
+
public void testCircle2()
{
Hep3Vector origin = new BasicHep3Vector(0,0,0);
@@ -67,6 +62,7 @@
assertEquals(radius*Math.PI,circle.getDistanceToInfiniteCylinder(2*radius),1e-14);
assertEquals(0,circle.getDistanceToInfiniteCylinder(0),1e-14);
}
+
public void testHelix()
{
Hep3Vector origin = new BasicHep3Vector(0,0,0);
@@ -136,70 +132,28 @@
}
public void testGetDistanceFromHelixToPoint() {
- Hep3Vector origin = new BasicHep3Vector(0,0,0);
- double radius = 1;
- double phi = Math.PI/2;
- double lambda = Math.PI/4;
- Helix helix = new Helix(origin,radius,phi,lambda);
-
- Hep3Vector pointOnHelix = helix.getPointAtDistance(2.1);
- System.out.println(helix.getDistanceFromHelixToPoint(pointOnHelix));
- System.out.println(helix.getDistanceFromHelixToPoint(helix.getPointAtDistance(21.1)));
- }
-
- public void testGetDistanceToPoint() {
- Hep3Vector origin = new BasicHep3Vector(0,0,0);
- double radius = 1;
- double phi = Math.PI/2;
+ Hep3Vector origin = new BasicHep3Vector(3, 6, 9);
+ double radius = 1.7;
+ double phi = Math.PI/4;
double lambda = Math.PI/4;
Helix helix = new Helix(origin, radius, phi, lambda);
-// assertEquals(helix.getDistanceFromHelixToPoint(origin), 0, 1e-14);
- for (int i=0; i<100; ++i) {
-// System.err.println("Trying to find" + i*.27);
- Hep3Vector pointOnHelix = helix.getPointAtDistance(i*.27);
- Hep3Vector pointToFind = new BasicHep3Vector(pointOnHelix.x()+1, pointOnHelix.y(), pointOnHelix.z());
-// System.err.println(helix.getDistanceFromHelixToPoint(pointToFind));
- }
-// System.err.println("Trying to find 2.66");
-// Hep3Vector pointOnHelix = helix.getPointAtDistance(2.66);
-// helix.getDistanceFromHelixToPoint(pointOnHelix);
-// System.err.println("Trying to find 20.66");
-// pointOnHelix = helix.getPointAtDistance(20.66);
-// helix.getDistanceFromHelixToPoint(pointOnHelix);
-// System.err.println("Trying to find 3*Pi/2");
- Hep3Vector pointOnHelix = helix.getPointAtDistance(radius*3*Math.PI/2);
- helix.getDistanceFromHelixToPoint(pointOnHelix);
-// System.err.println("Trying to find 1.3");
- pointOnHelix = helix.getPointAtDistance(1.3);
-// System.err.println(pointOnHelix);
-// assertEquals(helix.getDistanceFromHelixToPoint(pointOnHelix), 0, 1e-14);
- Hep3Vector newPoint = new BasicHep3Vector(0, 1.0, 0);
-// assertEquals(helix.getDistanceFromHelixToPoint(newPoint), 1.0, 1e-14);
- Hep3Vector centerOfCircle = new BasicHep3Vector(1, 0, 0);
-// assertEquals(helix.getDistanceFromHelixToPoint(centerOfCircle), radius, 1e-14);
- }
-
- public void testGetSignedClosestDifferenceToPoint() throws IOException {
- Hep3Vector origin = new BasicHep3Vector(0,0,0);
- double radius = 1;
- double phi = Math.PI/2;
- double lambda = Math.PI/4;
- Helix helix = new Helix(origin, radius, phi, lambda);
// test "random" points on helix
for (int i=0; i<500; ++i) {
- // System.err.println("Trying to find" + i*.37);
Hep3Vector pointOnHelix = helix.getPointAtDistance(i*.27);
- assertEquals(helix.getSignedClosestDifferenceToPoint(pointOnHelix), 0, 1e-10);
+// System.err.printf("%.3f, found: %s\n", i*.27, helix.getDistanceToPoint(pointOnHelix));
+ assertEquals(helix.getDistanceToPoint(pointOnHelix), i*.27, 1e-10);
+ }
+
+ // test that the distance to points has the right period along z
+ double period = 2*Math.PI*radius/Math.cos(lambda);
+ double offset = helix.getDistanceToPoint(new BasicHep3Vector(origin.x()+0.333, origin.y(), origin.z()));
+ for (int i=0; i<5; ++i) {
+ Hep3Vector newPoint = new BasicHep3Vector(origin.x()+0.333, origin.y(), origin.z()+i*radius*Math.tan(lambda)*Math.PI*2);
+ assertEquals(helix.getDistanceToPoint(newPoint), offset+i*period, 1e-10);
}
- Hep3Vector center = new BasicHep3Vector(radius, 0, 0);
- assertEquals(helix.getSignedClosestDifferenceToPoint(center), radius, 1e-10);
- Hep3Vector halfway = new BasicHep3Vector(radius/2, 0, 0);
- assertEquals(helix.getSignedClosestDifferenceToPoint(halfway), radius/2, 1e-10);
- Hep3Vector aboveCenter = new BasicHep3Vector(1, 2, 0);
- assertEquals(helix.getSignedClosestDifferenceToPoint(aboveCenter), -1, 1e-10);
}
-
+
private void assertEquals(Hep3Vector v1, Hep3Vector v2)
{
assertEquals(v1.x(),v2.x(), 1e-14);
CVSspam 0.2.8