3 added + 6 modified, total 9 files
lcsim/src/org/lcsim/fit/helicaltrack
diff -N HelicalTrackCross.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ HelicalTrackCross.java 10 Jun 2008 19:35:09 -0000 1.1
@@ -0,0 +1,194 @@
+/*
+ * HelicalTrackCross.java
+ *
+ * Created on May 1, 2008, 11:02 AM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.fit.helicaltrack;
+
+
+import hep.physics.matrix.BasicMatrix;
+import hep.physics.matrix.Matrix;
+import hep.physics.matrix.MutableMatrix;
+import hep.physics.matrix.MatrixOp;
+import hep.physics.matrix.SymmetricMatrix;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+import javax.crypto.spec.DHGenParameterSpec;
+
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
+import org.lcsim.event.TrackerHit;
+
+/**
+ *
+ * @author richp
+ */
+public class HelicalTrackCross extends HelicalTrackHit {
+ private HelicalTrackStrip _strip1;
+ private HelicalTrackStrip _strip2;
+ private double _eps = 1e-6;
+ private double _separation;
+ private double _salpha;
+ private Hep3Vector _p1;
+ private Hep3Vector _p2;
+
+ /** Creates a new instance of HelicalTrackCross */
+ public HelicalTrackCross(TrackerHit trackerhit, HelicalTrackStrip strip1, HelicalTrackStrip strip2) {
+ super(trackerhit);
+ _strip1 = strip1;
+ _strip2 = strip2;
+
+ // Check if the sensors are parallel to each other
+ if (VecOp.cross(_strip1.w(),_strip2.w()).magnitude() > _eps) {
+ throw new RuntimeException("Trying to construct a stereo hit from non-parallel sensor planes!");
+ }
+
+ // Check that the normals point in the same direction
+ if (VecOp.dot(_strip1.w(),_strip2.w()) < 0.) {
+ throw new RuntimeException("Trying to construct a stereo hit using an in-consistent coordinate system!");
+ }
+
+ // Calculate the seperation between the sensor planes
+ _separation = VecOp.dot(_strip1.w(), VecOp.sub(_strip2.origin(), _strip1.origin()));
+
+ // Calculate v1hat . u2hat, which is equivalent to sin(alpha) where alpha is the stereo angle
+ _salpha = VecOp.dot(_strip1.v(),_strip2.u());
+ if (Math.abs(_salpha) < _eps) {
+ throw new RuntimeException("Trying to construct a stereo hit using parallel strips!");
+ }
+
+ // Calculate the hit locations for v = 0: p = Origin + u * u_hat
+ _p1 = VecOp.add(_strip1.origin(), VecOp.mult(_strip1.umeas(), _strip1.u()));
+ _p2 = VecOp.add(_strip2.origin(), VecOp.mult(_strip2.umeas(), _strip2.u()));
+
+ // Save the nominal position and covariance matrix for a track from the origin
+ Hep3Vector pos = PositionFromOrigin();
+ super.setPosition(pos);
+ super.setCovMatrix(CovarianceFromOrigin());
+
+ // Determine the BarrelEndcapFlag
+ if (Math.abs(_strip1.w().z()) < Math.abs(pos.z()) / pos.magnitude()) super.setBarrelEndcapFlag(BarrelEndcapFlag.BARREL);
+ else if (pos.z() > 0.) super.setBarrelEndcapFlag(BarrelEndcapFlag.ENDCAP_NORTH);
+ else super.setBarrelEndcapFlag(BarrelEndcapFlag.ENDCAP_SOUTH);
+ }
+
+ public void setTrackDirection(TrackDirection trkdir, HelicalTrackFit helix) {
+ super.setCorrectedPosition(PositionOnHelix(trkdir));
+ super.setCorrectedCovMatrix(CovarianceOnHelix(trkdir, helix));
+ }
+
+ public void resetTrackDirection() {
+ super.setCorrectedPosition(PositionFromOrigin());
+ super.setCorrectedCovMatrix(CovarianceFromOrigin());
+ }
+
+ private Hep3Vector PositionFromOrigin() {
+ // Assume the particle is coming from the origin, so x2 = gamma * x1
+ // gamma = Origin2 . w_hat / Origin1 . w_hat
+ double gamma = VecOp.dot(_strip2.origin(), _strip2.w()) / VecOp.dot(_strip1.origin(), _strip1.w());
+ // dp = (p2 - gamma * p1)
+ Hep3Vector dp = VecOp.sub(_p2, VecOp.mult(gamma, _p1));
+ // Unmeasured coordinate v1: v1 = (dp . u2_hat) / (gamma * sin(alpha))
+ double v1 = VecOp.dot(dp, _strip2.u()) / (gamma * _salpha);
+ // Position of hit on strip 1: r1 = p1 + v1 * v1_hat
+ Hep3Vector r1 = VecOp.add(_p1, VecOp.mult(v1, _strip1.v()));
+ // Take position to be the midpoint of r1 and r2: r = 0.5 * (1 + gamma) * r1
+ return VecOp.mult(0.5 * (1 + gamma), r1);
+ }
+
+ private Hep3Vector PositionOnHelix(TrackDirection trkdir) {
+ // Get the track direction unit vector
+ Hep3Vector dir = trkdir.Direction();
+ // Gamma is the distance the particle travels between the two sensors: gamma = separation / (what . dir)
+ double gamma = _separation / VecOp.dot(_strip1.w(), dir);
+ // dp = p2 - (p1 + gamma * dir)
+ Hep3Vector dp = VecOp.sub(_p2, VecOp.add(_p1, VecOp.mult(gamma, dir)));
+ // Unmeasured coordinate v1: v1 = (dp . u2hat) / (v1_hat . u2_hat)
+ double v1 = VecOp.dot(dp, _strip2.u()) / VecOp.dot(_strip1.v(), _strip2.u());
+ // Position of hit on strip 1: r1 = p1 + v1 * v1_hat
+ Hep3Vector r1 = VecOp.add(_p1, VecOp.mult(v1, _strip1.v()));
+ // Take position to be the midpoint of x1 and x2: r = r1 + 0.5 * gamma * dir
+ return VecOp.add(r1, VecOp.mult(0.5 * gamma, dir));
+ }
+
+ private SymmetricMatrix CovarianceFromOrigin() {
+ // Assume the particle is coming from the origin, so x2 = gamma * x1
+ // gamma = Origin2 . w_hat / Origin1 . w_hat
+ double gamma = VecOp.dot(_strip2.origin(), _strip2.w()) / VecOp.dot(_strip1.origin(), _strip1.w());
+ // Calculate the scale factor: factor = (1 + gamma) / 2 * sin(alpha)
+ double factor = (1 + gamma) / (2. * _salpha);
+ // Calculate the covariance matrix from resolution: cov = factor * (v2 * v2^T * du1^2 + v1 * v1^T * du2^2)
+ BasicMatrix cov = Vec2Matrix(VecOp.mult(_strip1.du() * factor, _strip2.v()));
+ MatrixAdd(Vec2Matrix(VecOp.mult(_strip2.du() * factor, _strip1.v())), cov);
+ // Add contributions from the unknown direction equal to the seperation along v1 and v2
+ factor = 0.5 * _separation / (Math.sqrt(12) * _salpha);
+ MatrixAdd(Vec2Matrix(VecOp.mult(factor, _strip1.v())), cov);
+ MatrixAdd(Vec2Matrix(VecOp.mult(factor, _strip2.v())), cov);
+ return new SymmetricMatrix (cov);
+ }
+
+ private SymmetricMatrix CovarianceOnHelix(TrackDirection trkdir, HelicalTrackFit helix){
+ Hep3Vector dir = trkdir.Direction();
+ Matrix dirderiv = trkdir.Covariance();
+ MutableMatrix hcov = helix.covariance();
+ // Calculate phat x v1
+ Hep3Vector pcrossv1 = VecOp.cross(dir, _strip1.v());
+ // Calculate phat x v2
+ Hep3Vector pcrossv2 = VecOp.cross(dir, _strip2.v());
+ // Calculate the scale factor: _separation / 2 * sin(alpha) * (phat . w)^2
+ double factor = _separation / (2 * _salpha * Math.pow(VecOp.dot(dir, _strip1.w()),2));
+ // Create the matrix of position derivatives: d_i,j = dr_i / dphat_j
+ // The matrix d is given by factor * (v1 * (phat x v2)^T + v2 * (phat x v1)^T where ^T means transpose
+ MutableMatrix d = Vec2Matrix(_strip1.v(), VecOp.mult(factor, pcrossv2));
+ MatrixAdd(Vec2Matrix(_strip2.v(), VecOp.mult(factor, pcrossv1)), d);
+ Matrix dh = MatrixOp.mult(d, dirderiv);
+ // Construct the transpose of dh
+ Matrix dht = MatrixTranspose(dh);
+ // Calculate the covariance contributions from the direction uncertainty: cov = d * dirderiv * d^T
+ MutableMatrix cov = (MutableMatrix) MatrixOp.mult(dh, MatrixOp.mult(hcov, dht));
+ // Add the contributions from measurement errors: cov += (v2 * v2^T * du1^2 + v1 * v1^T * du2^2) / sin(alpha)^2
+ MatrixAdd(Vec2Matrix(VecOp.mult(_strip1.du()/_salpha, _strip2.v())), cov);
+ MatrixAdd(Vec2Matrix(VecOp.mult(_strip2.du()/_salpha, _strip1.v())), cov);
+ // Convert to a symmetric matrix
+ return new SymmetricMatrix(cov);
+ }
+
+ private BasicMatrix Vec2Matrix(Hep3Vector vec1, Hep3Vector vec2) {
+ BasicMatrix mat = new BasicMatrix(3, 3);
+ for (int i = 0; i < 3; i++) {
+ for (int j = 0; j < 3; j++) {
+ mat.setElement(i, j, vec1.v()[i] * vec2.v()[j]);
+ }
+ }
+ return mat;
+ }
+
+ private BasicMatrix Vec2Matrix(Hep3Vector vec) {
+ return Vec2Matrix(vec, vec);
+ }
+
+ private void MatrixAdd(Matrix m1, MutableMatrix m2) {
+ // Add the contents of matrix m1 to matrix m2
+ if (m1.getNColumns() != m2.getNColumns()) throw new RuntimeException("HelicalTrackCross matrix dimension mismatch!");
+ if (m1.getNRows() != m2.getNRows()) throw new RuntimeException("HelicalTrackCross matrix dimension mismatch");
+ for (int i = 0; i < m1.getNRows(); i++) {
+ for (int j = 0; j < m1.getNColumns(); j++) {
+ m2.setElement(i, j, m1.e(i, j) + m2.e(i, j));
+ }
+ }
+ return;
+ }
+
+ private Matrix MatrixTranspose(Matrix m) {
+ MutableMatrix mt = new BasicMatrix(m.getNColumns(), m.getNRows());
+ for (int i = 0; i < m.getNRows(); i++) {
+ for (int j = 0; j < m.getNColumns(); j++) {
+ mt.setElement(j, i, m.e(i, j));
+ }
+ }
+ return mt;
+ }
+}
\ No newline at end of file
lcsim/src/org/lcsim/fit/helicaltrack
diff -N HelicalTrackStrip.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ HelicalTrackStrip.java 10 Jun 2008 19:35:09 -0000 1.1
@@ -0,0 +1,72 @@
+/*
+ * HelicalTrackHit1DHit.java
+ *
+ * Created on May 1, 2008, 10:06 AM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.fit.helicaltrack;
+
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+
+/**
+ *
+ * @author richp
+ */
+public class HelicalTrackStrip {
+ private Hep3Vector _origin;
+ private Hep3Vector _u;
+ private Hep3Vector _v;
+ private Hep3Vector _w;
+ private double _umeas;
+ private double _du;
+ private double _vmin;
+ private double _vmax;
+
+ /** Creates a new instance of HelicalTrackStripHit */
+ public HelicalTrackStrip(Hep3Vector origin, Hep3Vector u, Hep3Vector v, double umeas, double du, double vmin, double vmax) {
+ _origin = origin;
+ _u = VecOp.unit(u);
+ _v = VecOp.unit(v);
+ _w = VecOp.cross(u,v);
+ _umeas = umeas;
+ _du = du;
+ _vmin = vmin;
+ _vmax = vmax;
+ }
+
+ public Hep3Vector origin() {
+ return _origin;
+ }
+
+ public Hep3Vector u() {
+ return _u;
+ }
+
+ public Hep3Vector v() {
+ return _v;
+ }
+
+ public Hep3Vector w() {
+ return _w;
+ }
+
+ public double umeas() {
+ return _umeas;
+ }
+
+ public double du() {
+ return _du;
+ }
+
+ public double vmin() {
+ return _vmin;
+ }
+
+ public double vmax() {
+ return _vmax;
+ }
+}
lcsim/src/org/lcsim/fit/helicaltrack
diff -N TrackDirection.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ TrackDirection.java 10 Jun 2008 19:35:09 -0000 1.1
@@ -0,0 +1,36 @@
+/*
+ * TrackDirection.java
+ *
+ * Created on May 8, 2008, 10:28 AM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.fit.helicaltrack;
+
+import hep.physics.matrix.Matrix;
+import hep.physics.vec.Hep3Vector;
+
+/**
+ *
+ * @author Richard Partridge
+ */
+public class TrackDirection {
+ private Hep3Vector _dir;
+ private Matrix _cov;
+
+ /** Creates a new instance of TrackDirection */
+ public TrackDirection(Hep3Vector dir, Matrix cov) {
+ _dir = dir;
+ _cov = cov;
+ }
+
+ public Hep3Vector Direction() {
+ return _dir;
+ }
+
+ public Matrix Covariance() {
+ return _cov;
+ }
+}
lcsim/src/org/lcsim/fit/helicaltrack
diff -u -r1.1 -r1.2
--- HelicalTrack2DHit.java 11 Dec 2007 20:35:04 -0000 1.1
+++ HelicalTrack2DHit.java 10 Jun 2008 19:35:09 -0000 1.2
@@ -7,6 +7,9 @@
package org.lcsim.fit.helicaltrack;
+import hep.physics.matrix.SymmetricMatrix;
+import hep.physics.vec.Hep3Vector;
+
import org.lcsim.event.TrackerHit;
/**
@@ -31,35 +34,9 @@
_zmin = zmin;
_zmax = zmax;
}
-
- /**
- * Creates a new instance of HelicalTrack2DHit from a TrackerHit.
- * The z coordinate of the hit is assumed to be at the center of
- * the strip, and the error in z is assumed to be half the strip
- * length in z.
- * @param hit TrackerHit to be used
- */
- public HelicalTrack2DHit(TrackerHit hit) {
- super(hit);
- double z = this.z();
- double[] cov = hit.getCovMatrix();
- _zmin = z - Math.sqrt(cov[5]);
- _zmax = z + Math.sqrt(cov[5]);
- }
-
- /**
- * Create a new instance of HelicalTrack2DHit from coordinates. This
- * constructor is provided to support the legacy method of performing
- * HelicalTrackFits based on arrays of coordinates and errors.
- * @param x x coordinate of the hit
- * @param y y coordinate of the hit
- * @param z z coordinate of the hit
- * @param drphi uncertainty in the r-phi coordinate
- * @param zmin minimum z of the hit
- * @param zmax maximum z of the hit
- */
- public HelicalTrack2DHit(double x, double y, double z, double drphi, double zmin, double zmax) {
- super(x, y, z, drphi);
+
+ public HelicalTrack2DHit(TrackerHit hit, Hep3Vector pos, SymmetricMatrix cov, double zmin, double zmax) {
+ super(hit, pos, cov);
_zmin = zmin;
_zmax = zmax;
}
lcsim/src/org/lcsim/fit/helicaltrack
diff -u -r1.3 -r1.4
--- HelicalTrack3DHit.java 14 Mar 2008 17:49:10 -0000 1.3
+++ HelicalTrack3DHit.java 10 Jun 2008 19:35:09 -0000 1.4
@@ -7,6 +7,9 @@
package org.lcsim.fit.helicaltrack;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.matrix.SymmetricMatrix;
+
import org.lcsim.event.TrackerHit;
import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
@@ -16,7 +19,6 @@
* @version 1.0
*/
public class HelicalTrack3DHit extends HelicalTrackHit {
- private BarrelEndcapFlag _beflag;
private double _dz;
/**
@@ -26,25 +28,15 @@
*/
public HelicalTrack3DHit(TrackerHit hit, BarrelEndcapFlag beflag) {
super(hit);
- _beflag = beflag;
-// _dz = Math.sqrt(hit.getCovMatrix()[5]);
- _dz = 0.005;
+ super.setBarrelEndcapFlag(beflag);
+ _dz = Math.sqrt(hit.getCovMatrix()[5]);
}
- /**
- * Creates a HelicalTrack3DHit from scratch.
- * @param x x coordinate of the hit
- * @param y y coordinate of the hit
- * @param z z coordinate of the hit
- * @param drphi error in the bend (r-phi) coordinate
- * @param dz error in the non-bend coordinate
- */
- public HelicalTrack3DHit(double x, double y, double z, double drphi, double dz) {
- super(x, y, z, drphi);
- _beflag = BarrelEndcapFlag.BARREL;
- _dz = dz;
+ public HelicalTrack3DHit(TrackerHit hit, Hep3Vector pos, SymmetricMatrix cov) {
+ super(hit, pos, cov);
+ _dz = Math.sqrt(cov.e(2,2));
}
-
+
/**
* Return the uncertainty in the z coordinate
* @return uncertainty in the z coordinate
@@ -53,18 +45,5 @@
return _dz;
}
- /**
- * Returns true if the hit is in a barrel sensor, false for
- * endcap (disk-like) sensors.
- * @return true if the hit is a barrel sensor
- */
- public boolean isBarrel() {
- return _beflag == BarrelEndcapFlag.BARREL;
- }
-
- public BarrelEndcapFlag getBarrelEndcapFlag(){
- return _beflag;
- }
-
}
\ No newline at end of file
lcsim/src/org/lcsim/fit/helicaltrack
diff -u -r1.20 -r1.21
--- HelicalTrackFitter.java 23 May 2008 16:10:09 -0000 1.20
+++ HelicalTrackFitter.java 10 Jun 2008 19:35:09 -0000 1.21
@@ -4,11 +4,12 @@
*
* Created on March 25, 2006, 6:11 PM
*
- * $Id: HelicalTrackFitter.java,v 1.20 2008/05/23 16:10:09 partridge Exp $
+ * $Id: HelicalTrackFitter.java,v 1.21 2008/06/10 19:35:09 partridge Exp $
*/
import hep.physics.matrix.SymmetricMatrix;
import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
import java.util.ArrayList;
import java.util.Collections;
@@ -22,6 +23,7 @@
import org.lcsim.fit.line.SlopeInterceptLineFitter;
import org.lcsim.fit.zsegment.ZSegmentFit;
import org.lcsim.fit.zsegment.ZSegmentFitter;
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
/**
* Fit a helix to a set of space points. First, a circle is fit to the x-y coordinates.
@@ -49,6 +51,7 @@
private SlopeInterceptLineFit _lfit;
private ZSegmentFit _zfit;
private HelicalTrackFit _fit;
+ private int nsig = 3;
public enum FitStatus {Success, CircleFitFailed, LineFitFailed, ZSegmentFitFailed};
@@ -72,24 +75,24 @@
public FitStatus fit(double[] x, double[] y, double[] z, double[] drphi, double[] dz, int np) {
List<HelicalTrackHit> hitcol = new ArrayList<HelicalTrackHit>();
for(int i=0; i<np; i++) {
+ Hep3Vector pos = new BasicHep3Vector(x[i], y[i], z[i]);
if (dz[i] > 0.) {
- hitcol.add(new HelicalTrack3DHit(x[i], y[i], z[i], drphi[i], dz[i]));
+ hitcol.add(new HelicalTrack3DHit(null, pos, MakeCov(pos, drphi[i], 0., dz[i])));
} else {
double zmin = z[i] - Math.abs(dz[i]);
double zmax = z[i] + Math.abs(dz[i]);
- hitcol.add(new HelicalTrack2DHit(x[i], y[i], z[i], drphi[i], zmin, zmax));
+ hitcol.add(new HelicalTrack2DHit(null, pos, MakeCov(pos, drphi[i], 0., 0.), zmin, zmax));
}
}
- Map<HelicalTrackHit, MultipleScatter> msmap = new HashMap<HelicalTrackHit, MultipleScatter>();
- return fit(hitcol, msmap);
+ return fit(hitcol);
}
public FitStatus fit(List<HelicalTrackHit> hitcol) {
Map<HelicalTrackHit, MultipleScatter> msmap = new HashMap<HelicalTrackHit, MultipleScatter>();
- return fit(hitcol, msmap);
+ return fit(hitcol, msmap, null);
}
- public FitStatus fit(List<HelicalTrackHit> hitcol, Map<HelicalTrackHit, MultipleScatter> msmap) {
+ public FitStatus fit(List<HelicalTrackHit> hitcol, Map<HelicalTrackHit, MultipleScatter> msmap, HelicalTrackFit oldhelix) {
// Initialize the various fitter outputs
_cfit = null;
@@ -110,16 +113,21 @@
if (hit instanceof HelicalTrack3DHit) pixel_hits.add(hit);
// Strip hits
if (hit instanceof HelicalTrack2DHit) strip_hits.add(hit);
+ // Cross hits
+ if (hit instanceof HelicalTrackCross) pixel_hits.add(hit);
}
+ // Check to make sure we have at least 3 circle hits before proceeding
+ int nc = circle_hits.size();
+ if (nc < 3) return FitStatus.CircleFitFailed;
+
// Sort the hits to be monotonic in z so that we can follow a curling track
// It is assumed that the first hit on a track is closer to the origin than the last hit
// It is also assumed that the track won't curl through an angle > 180 degrees between
// neighboring points. This might be a problem for curlers with small dip angle.
- int nhits = circle_hits.size();
Collections.sort(circle_hits);
double zfirst = circle_hits.get(0).z();
- double zlast = circle_hits.get(nhits-1).z();
+ double zlast = circle_hits.get(nc-1).z();
if (Math.abs(zfirst) > Math.abs(zlast)) Collections.reverse(circle_hits);
// Create the arrays that will hold the fit output
@@ -129,8 +137,6 @@
SymmetricMatrix cov = new SymmetricMatrix(5);
// Setup for doing the circle fit
- int nc = circle_hits.size();
- if (nc < 3) return FitStatus.CircleFitFailed;
double[] x = new double[nc];
double[] y = new double[nc];
double[] wrphi = new double[nc];
@@ -155,16 +161,19 @@
chisq[0] = _cfit.chisq();
ndof[0] = nc - 3;
- par[0] = -_cfit.dca();
+ // The circle fitter has the opposite sign convention for d0 than has been adopted by org.lcsim (L3 Note 1666)
+ // Flip the sign of d0 and the corresponding covariance matrix elements below
+ par[0] = -1. * _cfit.dca(); // fix d0 sign convention
par[1] = _cfit.phi();
par[2] = _cfit.curvature();
- cov.setElement(0, 0, _cfit.cov()[0]);
- cov.setElement(0, 1, -_cfit.cov()[1]);
- cov.setElement(1, 1, _cfit.cov()[2]);
- cov.setElement(0, 2, -_cfit.cov()[3]);
- cov.setElement(1, 2, _cfit.cov()[4]);
- cov.setElement(2, 2, _cfit.cov()[5]);
-
+ // Covariance matrix from circle fitter is based on parameters in the order (curv, phi0, d0)
+ // This is opposite the org.lcsim convention of (d0, phi0, curv), so take care in filling covariance matrix!!
+ cov.setElement(2, 2, _cfit.cov()[0]); // cov(curv,curv)
+ cov.setElement(1, 2, _cfit.cov()[1]); // cov(phi0,curv)
+ cov.setElement(1, 1, _cfit.cov()[2]); // cov(phi0,phi0)
+ cov.setElement(0, 2, -1. * _cfit.cov()[3]); // cov(d0,curv), fix d0 sign convention
+ cov.setElement(0, 1, -1. * _cfit.cov()[4]); // cov(d0,phi0), fix d0 sign convention
+ cov.setElement(0, 0, _cfit.cov()[5]); // cov(d0,d0)
// Calculate the arc lengths from the DCA to each hit
Map<HelicalTrackHit, Double> smap = getPathLengths(_cfit, hitcol);
@@ -179,11 +188,23 @@
// Store the coordinate and errors for the line fit
for(int i=0; i<npix; i++) {
- HelicalTrack3DHit hit = (HelicalTrack3DHit) pixel_hits.get(i);
+ HelicalTrackHit hit = pixel_hits.get(i);
z[i] = hit.z();
+ // Get the multiple scattering error
double dz_ms = 0.;
if (msmap.containsKey(hit)) dz_ms = msmap.get(hit).dz();
- dz[i] = Math.sqrt(Math.pow(hit.dz(),2)+Math.pow(dz_ms,2));
+ // Get the resolution error
+ double dz_res = hit.getCovMatrix()[5];
+ // If the hit is in the endcap, inflate the error on r by 1/slope to approximate error on z
+ if (hit.getBarrelEndcapFlag() != BarrelEndcapFlag.BARREL) {
+ // Default slope assumes track from origin
+ double slope = Math.abs(hit.z() / hit.r());
+ // If we know the track direction, use it to get the local slope
+ if (oldhelix != null) slope = oldhelix.slope();
+ // Inflate the claimed z resolution to account for the fact we are measuring r, not z
+ dz_res = hit.dr() * slope;
+ }
+ dz[i] = Math.sqrt(dz_res*dz_res + dz_ms*dz_ms);
s[i] = smap.get(hit);
}
@@ -209,13 +230,24 @@
// If we have one pixel hit, turn it into a pseudo strip hit
if (npix == 1) {
- HelicalTrack3DHit hit = (HelicalTrack3DHit) pixel_hits.get(0);
+ HelicalTrackHit hit = pixel_hits.get(0);
double dz_ms = 0.;
if (msmap.containsKey(hit)) dz_ms = msmap.get(hit).dz();
- double dz = Math.sqrt(Math.pow(hit.dz(),2)+Math.pow(dz_ms,2));
- double zmin = hit.z() - dz;
- double zmax = hit.z() + dz;
- HelicalTrack2DHit fakestriphit = new HelicalTrack2DHit(hit.x(), hit.y(), hit.z(), hit.drphi(), zmin, zmax);
+ double dz_cov = hit.getCorrectedCovMatrix().diagonal(2);
+ double dz = Math.sqrt(dz_cov + dz_ms*dz_ms);
+ if (hit.getBarrelEndcapFlag() != BarrelEndcapFlag.BARREL) {
+ // Default slope assumes track from origin
+ double slope = Math.abs(hit.z() / hit.r());
+ // If we know the track direction, use it to get the local slope
+ if (oldhelix != null) slope = oldhelix.slope();
+ double dzres = hit.dr() * slope;
+ dz = Math.sqrt(dzres*dzres + dz_ms*dz_ms);
+ }
+ double zmin = hit.z() - nsig * dz;
+ double zmax = hit.z() + nsig * dz;
+ Hep3Vector fakepos = new BasicHep3Vector(hit.getPosition());
+ SymmetricMatrix fakecov = new SymmetricMatrix(3, hit.getCovMatrix(), true);
+ HelicalTrack2DHit fakestriphit = new HelicalTrack2DHit(null, fakepos, fakecov, zmin, zmax);
strip_hits.add(fakestriphit);
smap.put(fakestriphit, smap.get(hit));
}
@@ -323,4 +355,16 @@
return new CircleFit(oldfit.xref(), oldfit.yref(), curv, phi0, dca, oldfit.chisq(), cov);
}
+
+ private SymmetricMatrix MakeCov(Hep3Vector pos, double drphi, double dr, double dz) {
+ double x = pos.x();
+ double y = pos.y();
+ double r2 = x * x + y * y;
+ SymmetricMatrix cov = new SymmetricMatrix(3);
+ cov.setElement(0, 0, (y*y * drphi*drphi + x*x * dr*dr) / r2);
+ cov.setElement(0, 1, x * y * (dr*dr - drphi*drphi) / r2);
+ cov.setElement(1, 1, (x*x * drphi*drphi + y*y * dr*dr) / r2);
+ cov.setElement(2, 2, dz*dz);
+ return cov;
+ }
}
\ No newline at end of file
lcsim/src/org/lcsim/fit/helicaltrack
diff -u -r1.6 -r1.7
--- HelicalTrackHit.java 14 Mar 2008 18:21:50 -0000 1.6
+++ HelicalTrackHit.java 10 Jun 2008 19:35:09 -0000 1.7
@@ -8,16 +8,21 @@
package org.lcsim.fit.helicaltrack;
-import java.lang.Comparable;
-
import java.util.ArrayList;
+import java.lang.Comparable;
import java.util.List;
+
+import hep.physics.matrix.SymmetricMatrix;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
+
import org.lcsim.contrib.onoprien.tracking.hit.DigiTrackerHit;
import org.lcsim.contrib.onoprien.tracking.hit.TrackerCluster;
import org.lcsim.contrib.onoprien.tracking.hitmaking.OldTrackerHit;
import org.lcsim.event.MCParticle;
import org.lcsim.event.TrackerHit;
import org.lcsim.event.base.BaseTrackerHitMC;
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
/**
* Encapsulate the hit information needed by HelicalTrackFitter
@@ -27,9 +32,11 @@
public class HelicalTrackHit implements Comparable, TrackerHit {
private TrackerHit _trackerhit;
- private double _x;
- private double _y;
- private double _z;
+ private Hep3Vector _pos;
+ private Hep3Vector _poscor;
+ private SymmetricMatrix _cov;
+ private SymmetricMatrix _covcor;
+ private BarrelEndcapFlag _beflag;
private double _r;
private double _phi;
private double _drphi;
@@ -40,41 +47,18 @@
* @param hit Corresponding TrackerHit
*/
public HelicalTrackHit(TrackerHit hit) {
- _trackerhit = hit;
- double[] pos = hit.getPosition();
- _x = pos[0];
- _y = pos[1];
- _z = pos[2];
- setPolarCoordinates();
- double[] cov = hit.getCovMatrix();
-// _drphi = Math.sqrt(_y * _y * cov[0] + _x * _x * cov[2] - 2. * _x * _y * cov[1]) / _r;
-// _dr = Math.sqrt(_x * _x * cov[0] + _y * _y * cov[2] + 2. * _x * _y * cov[1]) / _r;
- if (_r > 100) {
- _drphi = 0.007;
- } else {
- _drphi = 0.005;
- }
- _dr = 0.0;
+ this(hit, new BasicHep3Vector(hit.getPosition()), new SymmetricMatrix(3, hit.getCovMatrix(), true));
}
- /**
- * Create a new instance of HelicalTrackHit from coordinates. This
- * constructor is provided to support the legacy method of performing
- * HelicalTrackFits based on arrays of coordinates and errors.
- * @param x x coordinate of the hit
- * @param y y coordinate of the hit
- * @param z z coordinate of the hit
- * @param drphi Uncertainty in the r-phi coordinate
- */
- public HelicalTrackHit(double x, double y, double z, double drphi) {
- _trackerhit = null;
- _x = x;
- _y = y;
- _z = z;
- _r = Math.sqrt(_x * _x + _y * _y);
- _phi = Math.atan2(_y, _x);
- _drphi = drphi;
- _dr = 0.;
+ public HelicalTrackHit(TrackerHit hit, Hep3Vector pos, SymmetricMatrix cov) {
+ _trackerhit = hit;
+ _pos = pos;
+ _poscor = pos;
+ _cov = cov;
+ _covcor = cov;
+ // Default to a barrel hit
+ _beflag = BarrelEndcapFlag.BARREL;
+ setPolarVariables();
}
/**
@@ -90,7 +74,7 @@
* @return x coordinate
*/
public double x() {
- return _x;
+ return _poscor.x();
}
/**
@@ -98,7 +82,7 @@
* @return y coordinate
*/
public double y() {
- return _y;
+ return _poscor.y();
}
/**
@@ -106,7 +90,7 @@
* @return z coordinate
*/
public double z() {
- return _z;
+ return _poscor.z();
}
/**
@@ -141,6 +125,10 @@
return _dr;
}
+ public BarrelEndcapFlag getBarrelEndcapFlag() {
+ return _beflag;
+ }
+
/**
* Implement comparable interface to allow hits to be sorted
* by their z coordinate
@@ -149,16 +137,16 @@
*/
public int compareTo(Object hit2) {
double zhit = ((HelicalTrackHit) hit2).z();
- if (_z < zhit) return -1;
- if (_z == zhit) return 0;
+ if (_poscor.z() < zhit) return -1;
+ if (_poscor.z() == zhit) return 0;
return 1;
}
- /**
- * Returns a list of MCParticles belonging to this hit.
- * null is returned if no list can be found.
- *
- * @return A list of MCParticles, or null if none can be retrieved.
+ /**
+ * Returns a list of MCParticles belonging to this hit.
+ * null is returned if no list can be found.
+ *
+ * @return A list of MCParticles, or null if none can be retrieved.
*/
public List<MCParticle> getMCParticles(){
@@ -166,66 +154,109 @@
if (_trackerhit instanceof BaseTrackerHitMC)
return ((BaseTrackerHitMC)_trackerhit).mcParticles();
- //Any of the SiStripSip BaseTrackerHits return a set of MCParticles, which must be converted to a list
+ //Any of the SiStripSip BaseTrackerHits return a set of MCParticles, which must be converted to a list
else if (_trackerhit instanceof org.lcsim.contrib.SiStripSim.BaseTrackerHit) {
- List<MCParticle> out = new ArrayList<MCParticle>();
- out.addAll(((org.lcsim.contrib.SiStripSim.BaseTrackerHit)_trackerhit).getMCParticles());
- return out;
+ List<MCParticle> out = new ArrayList<MCParticle>();
+ out.addAll(((org.lcsim.contrib.SiStripSim.BaseTrackerHit)_trackerhit).getMCParticles());
+ return out;
}
- //Getting this out of OldTrackerHit is slightly more convoluted...
+ //Getting this out of OldTrackerHit is slightly more convoluted...
else if (_trackerhit instanceof OldTrackerHit){
- OldTrackerHit h = (OldTrackerHit) _trackerhit;
+ OldTrackerHit h = (OldTrackerHit) _trackerhit;
List<MCParticle> out = new ArrayList<MCParticle>();
for (TrackerCluster cluster : h.getClusters()) {
for (DigiTrackerHit dhit : cluster.getDigiHits()){
for (DigiTrackerHit dhit2 : dhit.getElementalHits()){
- MCParticle p = dhit2.getMCParticle();
- if(p != null && !out.contains(p)) out.add(p);
+ MCParticle p = dhit2.getMCParticle();
+ if(p != null && !out.contains(p)) out.add(p);
}
}
}
- return out;
+ return out;
}
- return null;
+ return null;
}
public double[] getPosition() {
- return _trackerhit.getPosition();
+ return _pos.v();
}
-
+
+ public Hep3Vector getCorrectedPosition() {
+ return _poscor;
+ }
+
public double[] getCovMatrix() {
- return _trackerhit.getCovMatrix();
+ return _cov.asPackedArray(true);
}
-
+
+ public SymmetricMatrix getCorrectedCovMatrix() {
+ return _covcor;
+ }
+
public double getdEdx() {
- return _trackerhit.getdEdx();
+ if (_trackerhit != null) return _trackerhit.getdEdx();
+ else return 0.;
}
-
+
public double getTime() {
- return _trackerhit.getTime();
+ if (_trackerhit != null) return _trackerhit.getTime();
+ else return 0.;
}
-
+
public int getType() {
- return _trackerhit.getType();
+ if (_trackerhit != null) return _trackerhit.getType();
+ else return 0;
}
-
+
public List getRawHits() {
- return _trackerhit.getRawHits();
+ if (_trackerhit != null) return _trackerhit.getRawHits();
+ else return null;
+ }
+
+ protected void setBarrelEndcapFlag(BarrelEndcapFlag beflag) {
+ _beflag = beflag;
+ }
+
+ protected void setPosition(Hep3Vector pos) {
+ _pos = pos;
+ if (_poscor == null) setCorrectedPosition(pos);
+ return;
}
+ protected void setCorrectedPosition(Hep3Vector poscor) {
+ _poscor = poscor;
+ setPolarVariables();
+ return;
+ }
+
+ protected void setCovMatrix(SymmetricMatrix cov) {
+ _cov = cov;
+ if (_covcor == null) setCorrectedCovMatrix(cov);
+ return;
+ }
+
+ protected void setCorrectedCovMatrix(SymmetricMatrix covcor) {
+ _covcor = covcor;
+ setPolarVariables();
+ return;
+ }
/**
* Calculate the polar coordinates _r, _phi from the cartesian
* coordinates _x, _y and cache them in this object since we
* expect them to be used repeatedly by the track finding code.
*/
- private void setPolarCoordinates() {
- _r = Math.sqrt(_x * _x + _y * _y);
- _phi = Math.atan2(_y, _x);
+ private void setPolarVariables() {
+ double x = _poscor.x();
+ double y = _poscor.y();
+ _r = Math.sqrt(x*x + y*y);
+ _phi = Math.atan2(y, x);
if (_phi < 0.) _phi += 2. * Math.PI;
+ _drphi = Math.sqrt(y*y * _covcor.e(0,0) + x*x * _covcor.e(1,1)
+ - 2. * x * y * _covcor.e(0,1)) / _r;
+ _dr = Math.sqrt(x*x * _covcor.e(0,0) + y*y * _covcor.e(1,1)
+ + 2. * x * y * _covcor.e(0,1)) / _r;
return;
}
-
-
}
\ No newline at end of file
lcsim/src/org/lcsim/fit/helicaltrack
diff -u -r1.7 -r1.8
--- HelicalTrackHitDriver.java 14 Mar 2008 17:49:10 -0000 1.7
+++ HelicalTrackHitDriver.java 10 Jun 2008 19:35:09 -0000 1.8
@@ -7,8 +7,10 @@
package org.lcsim.fit.helicaltrack;
+import hep.physics.matrix.SymmetricMatrix;
import hep.physics.vec.BasicHep3Vector;
import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
import java.util.ArrayList;
import java.util.List;
@@ -20,6 +22,7 @@
import org.lcsim.contrib.onoprien.tracking.geom.Sensor;
import org.lcsim.contrib.onoprien.tracking.geom.SensorType;
import org.lcsim.contrib.onoprien.tracking.hit.TrackerCluster;
+import org.lcsim.contrib.onoprien.tracking.hitmaking.hitmakers.TrackerHitMakerBasic;
import org.lcsim.contrib.onoprien.tracking.hitmaking.OldTrackerHit;
import org.lcsim.detector.DetectorIdentifierHelper;
import org.lcsim.detector.identifier.IIdentifier;
@@ -27,6 +30,7 @@
import org.lcsim.event.TrackerHit;
import org.lcsim.fit.helicaltrack.HelicalTrackHit;
import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
+import org.lcsim.spacegeom.SpacePointVector;
import org.lcsim.util.Driver;
/**
@@ -34,14 +38,14 @@
* collections. The resulting HelicalTrackHits encapsulate the information
* needed to perform a helical track hit for either a segmented strip
* detector, a pixel detector, or cross hits from a stereo detector.
- *
+ *
* At this time, this driver supports the virtual segmentation hits
* produced by the packages in contrib.onoprien.tracking. Additional
* coding needs to be done to support fully segmented detectors.
- *
+ *
* The list of hit collections to be converted must be specified
* before the process method is executed.
- *
+ *
* Currently,
* @author Richard Partridge
* @version 1.0
@@ -57,16 +61,20 @@
*/
VirtualSegmentation,
/**
- *
+ *
* Digitized (SiTrackerHit)
*/
- Digitized
-
+ Digitized
+
}
+ private TrackerHitMakerBasic _hitmaker = new TrackerHitMakerBasic();
private List<String> _vscol = new ArrayList<String>();
private List<String> _smcol = new ArrayList<String>();
- private List<String> _digcol = new ArrayList<String>();
+ private List<String> _digcol = new ArrayList<String>();
private String _outname = "HelicalTrackHits";
+ private Hep3Vector _uloc = new BasicHep3Vector(1., 0., 0.);
+ private Hep3Vector _vloc = new BasicHep3Vector(0., 1., 0.);
+ private Hep3Vector _orgloc = new BasicHep3Vector(0., 0., 0.);
/** Creates a new instance of HelicalTrackHitDriver */
public HelicalTrackHitDriver() {
@@ -78,7 +86,7 @@
*/
public void process(EventHeader event) {
super.process(event);
-
+
// Initialize the list of HelicalTrackHits and vector with local z direction
List<HelicalTrackHit> helhits = new ArrayList<HelicalTrackHit>();
Hep3Vector lz = new BasicHep3Vector(0., 0., 1.);
@@ -121,10 +129,25 @@
// If there are more than one (i.e., two) clusters for this hit, we
// have a cross formed from two stereo layers
- if (clist.size() > 1) {
-
- // Stereo hit found - make a HelicalTrack3DHit
- helhits.add((HelicalTrackHit) new HelicalTrack3DHit(hit, beflag));
+ if (vshit.isStereo()) {
+ if (clist.size() != 2) throw new RuntimeException("More than 2 clusters on a stereo hit");
+ List<HelicalTrackStrip> strips = new ArrayList<HelicalTrackStrip>();
+ for (TrackerCluster cl : clist) {
+ if (cl.getTrackerHits().size() == 0) cl.addTrackerHit(_hitmaker.make(cl));
+ Sensor s = cl.getSensor();
+ Hep3Vector org = s.localToGlobal(_orgloc);
+ Hep3Vector uhat = VecOp.sub(s.localToGlobal(_uloc), org);
+ Hep3Vector vhat = VecOp.sub(s.localToGlobal(_vloc), org);
+ org.lcsim.contrib.onoprien.tracking.hit.TrackerHit dhit = cl.getTrackerHits().get(0);
+ double umeas = dhit.getLocalPosition().x();
+ double vmax = dhit.getSegment().magnitude() / 2.;
+ double vmin = -vmax;
+ double du = Math.sqrt(dhit.getLocalCovMatrix().diagonal(0));
+
+ Hep3Vector loc = dhit.getLocalPosition();
+ strips.add(new HelicalTrackStrip(org, uhat, vhat, umeas, du, vmin, vmax));
+ }
+ helhits.add((HelicalTrackHit) new HelicalTrackCross(vshit, strips.get(0), strips.get(1)));
} else {
// Check the number of readout dimensions (1 = strip, 2 = pixel)
@@ -163,13 +186,13 @@
//these are barrel strip hits
if (hit instanceof SiTrackerHitStrip1D) {
- //our current scheme only supports 2d hits in the barrel...
- //if there are 2d hits in the endcap, we'll ignore them
+ //our current scheme only supports 2d hits in the barrel...
+ //if there are 2d hits in the endcap, we'll ignore them
if(!id.isBarrel(detid)){
- continue;
+ continue;
}
- SiTrackerHitStrip1D h = (SiTrackerHitStrip1D) hit;
+ SiTrackerHitStrip1D h = (SiTrackerHitStrip1D) hit;
double z1 = h.getHitSegment().getPoints().get(0).z();
double z2 = h.getHitSegment().getPoints().get(1).z();
@@ -179,27 +202,24 @@
helhits.add((HelicalTrackHit) new HelicalTrack2DHit(h,zmin,zmax));
}
- //these are stereod strips and pixels
+ //these are stereod strips and pixels
else if (hit instanceof SiTrackerHitStrip2D || hit instanceof SiTrackerHitPixel) {
- BarrelEndcapFlag flag;
+ BarrelEndcapFlag flag;
//decide the flag
if(id.isBarrel(detid)) {
- flag = BarrelEndcapFlag.BARREL;
- }
- else if (id.isEndcapNegative(detid)) {
- flag = BarrelEndcapFlag.ENDCAP_SOUTH;
- }
- else if (id.isEndcapPositive(detid)) {
- flag = BarrelEndcapFlag.ENDCAP_NORTH;
- }
- else {
+ flag = BarrelEndcapFlag.BARREL;
+ } else if (id.isEndcapNegative(detid)) {
+ flag = BarrelEndcapFlag.ENDCAP_SOUTH;
+ } else if (id.isEndcapPositive(detid)) {
+ flag = BarrelEndcapFlag.ENDCAP_NORTH;
+ } else {
flag = BarrelEndcapFlag.UNKNOWN; //perhaps this should throw an error instead?
}
helhits.add((HelicalTrackHit) new HelicalTrack3DHit(hit,flag));
- }
+ }
}
}
@@ -217,8 +237,7 @@
_vscol.add(name);
} else if (type == HitType.Smeared) {
_smcol.add(name);
- }
- else if (type == HitType.Digitized)
+ } else if (type == HitType.Digitized)
_digcol.add(name);
return;
}
lcsim/src/org/lcsim/fit/helicaltrack
diff -u -r1.4 -r1.5
--- HelixUtils.java 7 Feb 2008 18:20:51 -0000 1.4
+++ HelixUtils.java 10 Jun 2008 19:35:09 -0000 1.5
@@ -9,6 +9,8 @@
import hep.physics.vec.BasicHep3Vector;
import hep.physics.vec.Hep3Vector;
+import hep.physics.matrix.BasicMatrix;
+import hep.physics.matrix.Matrix;
import java.util.ArrayList;
import java.util.List;
@@ -35,6 +37,8 @@
private double _yc;
private double _x0;
private double _y0;
+ private double _sphi0;
+ private double _cphi0;
private double _sth;
private double _cth;
private double _minslope = 1.e-6;
@@ -110,12 +114,46 @@
}
public Hep3Vector Direction(HelicalTrackFit helix, double s) {
+
+ // Cache the helix parameters
+ HelixCache(helix);
double phi = _phi0 - s / _RC;
double ux = Math.cos(phi) * _sth;
double uy = Math.sin(phi) * _sth;
double uz = _cth;
return new BasicHep3Vector(ux, uy, uz);
}
+
+ public Matrix DirectionDerivates(HelicalTrackFit helix, double s) {
+
+ // Cache the helix parameters
+ HelixCache(helix);
+
+ // Get the unit vector phat that is a "momentum" unit vector giving the track direction
+ Hep3Vector phat = Direction(helix, s);
+
+ // Create the matrix that will hold the derivatives
+ BasicMatrix deriv = new BasicMatrix(3,5);
+
+ // Calculate the non-zero derivatives of the direction with respect to the helix parameters
+ deriv.setElement(0, helix.curvatureIndex, (phat.x() - _cphi0 * _sth) / _omega); // dphat_x / domega
+ deriv.setElement(1, helix.curvatureIndex, (phat.y() - _sphi0 * _sth) / _omega); // dphat_y / domega
+ deriv.setElement(0, helix.dcaIndex, -_omega * _cphi0 * _sth); // dphat_x / ddca
+ deriv.setElement(1, helix.dcaIndex, -_omega * _sphi0 * _sth); // dphat_y / ddca
+ deriv.setElement(0, helix.phi0Index, -(1 - _dca * _omega) * _sphi0 * _sth); // dphat_x / dphi0
+ deriv.setElement(1, helix.phi0Index, (1 - _dca * _omega) * _cphi0 * _sth); // dphat_y / dphi0
+ deriv.setElement(0, helix.slopeIndex, -phat.x() * _sth * _cth); // dphat_x / dslope
+ deriv.setElement(1, helix.slopeIndex, -phat.y() * _sth * _cth); // dphat_y / dslope
+ deriv.setElement(2, helix.slopeIndex, Math.pow(_sth, 3)); // dphat_z / dslope
+
+ return deriv;
+ }
+
+ public TrackDirection getTrackDirection(HelicalTrackFit helix, double s) {
+ Hep3Vector dir = Direction(helix, s);
+ Matrix deriv = DirectionDerivates(helix, s);
+ return new TrackDirection(dir, deriv);
+ }
public Hep3Vector PointOnHelix(HelicalTrackFit helix, double s) {
@@ -161,10 +199,12 @@
_RC = 1.0 / _cfit.curvature();
_xr = _cfit.xref();
_yr = _cfit.yref();
- _xc = _xr + (_RC - _dca) * Math.sin(_phi0);
- _yc = _yr - (_RC - _dca) * Math.cos(_phi0);
- _x0 = _xr - _dca * Math.sin(_phi0);
- _y0 = _yr + _dca * Math.cos(_phi0);
+ _sphi0 = Math.sin(_phi0);
+ _cphi0 = Math.cos(_phi0);
+ _xc = _xr + (_RC - _dca) * _sphi0;
+ _yc = _yr - (_RC - _dca) * _cphi0;
+ _x0 = _xr - _dca * _sphi0;
+ _y0 = _yr + _dca * _cphi0;
}
return;
}
@@ -176,15 +216,18 @@
double[] param = _helix.parameters();
_dca = _helix.dca();
_phi0 = _helix.phi0();
- _RC = 1.0 / _helix.curvature();
+ _omega = _helix.curvature();
+ _RC = 1.0 / _omega;
_z0 = _helix.z0();
_slope = _helix.slope();
_xr = 0.;
_yr = 0.;
- _xc = _xr + (_RC - _dca) * Math.sin(_phi0);
- _yc = _yr - (_RC - _dca) * Math.cos(_phi0);
- _x0 = _xr - _dca * Math.sin(_phi0);
- _y0 = _yr + _dca * Math.cos(_phi0);
+ _sphi0 = Math.sin(_phi0);
+ _cphi0 = Math.cos(_phi0);
+ _xc = _xr + (_RC - _dca) * _sphi0;
+ _yc = _yr - (_RC - _dca) * _cphi0;
+ _x0 = _xr - _dca * _sphi0;
+ _y0 = _yr + _dca * _cphi0;
_sth = _helix.sth();
_cth = _helix.cth();
}
CVSspam 0.2.8