Print

Print


Commit in lcsim on MAIN
src/org/lcsim/contrib/JanStrube/tracking/NewTrack.java+2-21.5 -> 1.6
                                        /TransitionalTrack.java+2-21.2 -> 1.3
                                        /HelixSwimmer.java+2-21.5 -> 1.6
                                        /NewFastMCTrackFactory.java+2-21.7 -> 1.8
                                        /NewMCFastTrackDriver.java+6-41.1 -> 1.2
                                        /Constants.java-121.1 removed
src/org/lcsim/contrib/JanStrube/vtxFitter/Fitter.java+53-491.6 -> 1.7
                                         /vertexing.lyx+76-91.6 -> 1.7
                                         /vertexing.pdf[binary]1.3 -> 1.4
                                         /Vertex.java+4-41.2 -> 1.3
test/org/lcsim/contrib/JanStrube/tracking/FitterTest.java+4-31.1 -> 1.2
+151-89
1 removed + 10 modified, total 11 files
Removing my own value for the field conversion after establishing the the org.lcsim value has been correct since Sowmass.
This magically makes the fitter appear to converge to sensible results. More tests are still needed.

lcsim/src/org/lcsim/contrib/JanStrube/tracking
NewTrack.java 1.5 -> 1.6
diff -u -r1.5 -r1.6
--- NewTrack.java	7 Aug 2006 21:30:27 -0000	1.5
+++ NewTrack.java	15 Aug 2006 23:49:55 -0000	1.6
@@ -1,5 +1,5 @@
 /**
- * @version $Id: NewTrack.java,v 1.5 2006/08/07 21:30:27 jstrube Exp $
+ * @version $Id: NewTrack.java,v 1.6 2006/08/15 23:49:55 jstrube Exp $
  */
 package org.lcsim.contrib.JanStrube.tracking;
 
@@ -14,7 +14,7 @@
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.z0;
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.tanLambda;
 
-import static org.lcsim.contrib.JanStrube.tracking.Constants.fieldConversion;
+import static org.lcsim.constants.Constants.fieldConversion;
 
 import static java.lang.Math.sin;
 import static java.lang.Math.cos;

lcsim/src/org/lcsim/contrib/JanStrube/tracking
TransitionalTrack.java 1.2 -> 1.3
diff -u -r1.2 -r1.3
--- TransitionalTrack.java	18 Jul 2006 13:20:41 -0000	1.2
+++ TransitionalTrack.java	15 Aug 2006 23:49:55 -0000	1.3
@@ -1,5 +1,5 @@
 /**
- * @version $Id: TransitionalTrack.java,v 1.2 2006/07/18 13:20:41 jstrube Exp $
+ * @version $Id: TransitionalTrack.java,v 1.3 2006/08/15 23:49:55 jstrube Exp $
  */
 package org.lcsim.contrib.JanStrube.tracking;
 
@@ -17,7 +17,7 @@
  */
 public class TransitionalTrack implements Track {
 
-    org.lcsim.contrib.JanStrube.tracking.Track _track;
+    public org.lcsim.contrib.JanStrube.tracking.Track _track;
     TransitionalTrack(org.lcsim.contrib.JanStrube.tracking.Track t) {
         _track = t;        
     }

lcsim/src/org/lcsim/contrib/JanStrube/tracking
HelixSwimmer.java 1.5 -> 1.6
diff -u -r1.5 -r1.6
--- HelixSwimmer.java	7 Aug 2006 21:30:27 -0000	1.5
+++ HelixSwimmer.java	15 Aug 2006 23:49:56 -0000	1.6
@@ -9,7 +9,7 @@
 import static java.lang.Math.cos;
 import static java.lang.Math.sqrt;
 import static java.lang.Math.signum;
-import static org.lcsim.contrib.JanStrube.tracking.Constants.fieldConversion;
+import static org.lcsim.constants.Constants.fieldConversion;
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.d0;
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.phi0;
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.omega;
@@ -23,7 +23,7 @@
  * For more info on swimming see <a href="doc-files/transport.pdf">this paper</a> by Paul Avery.
  * 
  * @author jstrube
- * @version $Id: HelixSwimmer.java,v 1.5 2006/08/07 21:30:27 jstrube Exp $
+ * @version $Id: HelixSwimmer.java,v 1.6 2006/08/15 23:49:56 jstrube Exp $
  */
 public class HelixSwimmer {
     public class SpatialParameters {

lcsim/src/org/lcsim/contrib/JanStrube/tracking
NewFastMCTrackFactory.java 1.7 -> 1.8
diff -u -r1.7 -r1.8
--- NewFastMCTrackFactory.java	7 Aug 2006 21:30:27 -0000	1.7
+++ NewFastMCTrackFactory.java	15 Aug 2006 23:49:56 -0000	1.8
@@ -1,5 +1,5 @@
 /**
- * @version $Id: NewFastMCTrackFactory.java,v 1.7 2006/08/07 21:30:27 jstrube Exp $
+ * @version $Id: NewFastMCTrackFactory.java,v 1.8 2006/08/15 23:49:56 jstrube Exp $
  */
 package org.lcsim.contrib.JanStrube.tracking;
 
@@ -9,7 +9,7 @@
 import static java.lang.Math.cos;
 import static java.lang.Math.sin;
 import static java.lang.Math.sqrt;
-import static org.lcsim.contrib.JanStrube.tracking.Constants.fieldConversion;
+import static org.lcsim.constants.Constants.fieldConversion;
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.d0;
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.omega;
 import static org.lcsim.contrib.JanStrube.tracking.Track.ParameterName.phi0;

lcsim/src/org/lcsim/contrib/JanStrube/tracking
NewMCFastTrackDriver.java 1.1 -> 1.2
diff -u -r1.1 -r1.2
--- NewMCFastTrackDriver.java	18 Jul 2006 13:20:41 -0000	1.1
+++ NewMCFastTrackDriver.java	15 Aug 2006 23:49:57 -0000	1.2
@@ -1,5 +1,5 @@
 /**
- * @version $Id: NewMCFastTrackDriver.java,v 1.1 2006/07/18 13:20:41 jstrube Exp $
+ * @version $Id: NewMCFastTrackDriver.java,v 1.2 2006/08/15 23:49:57 jstrube Exp $
  */
 package org.lcsim.contrib.JanStrube.tracking;
 
@@ -59,6 +59,7 @@
         boolean hist = getHistogramLevel() > 0;
 
         List<org.lcsim.event.Track> trackList = new ArrayList<org.lcsim.event.Track>();
+	List<Track> newTrackList = new ArrayList<Track>();
         for (Particle p : event.getMCParticles()) {
             // filter for FINAL_STATE
             if (p.getGeneratorStatus() != Particle.FINAL_STATE) {
@@ -87,10 +88,11 @@
             }
 
             Track t = factory.getTrack(p.getMomentum(), p.getOrigin(), (int)p.getCharge());
-            org.lcsim.event.Track passpartout = new TransitionalTrack(t);
-            trackList.add(passpartout);
+            // org.lcsim.event.Track passpartout = new TransitionalTrack(t);
+            // trackList.add(passpartout);
+	    newTrackList.add(t);
         }
-        event.put(EventHeader.TRACKS, trackList, Track.class, 0);
+        event.put(EventHeader.TRACKS, newTrackList, Track.class, 0);
     }
 
     public void conditionsChanged(ConditionsEvent event) {

lcsim/src/org/lcsim/contrib/JanStrube/tracking
Constants.java removed after 1.1
diff -N Constants.java
--- Constants.java	15 Jul 2006 09:12:45 -0000	1.1
+++ /dev/null	1 Jan 1970 00:00:00 -0000
@@ -1,12 +0,0 @@
-/**
- * @version $Id: Constants.java,v 1.1 2006/07/15 09:12:45 jstrube Exp $
- */
-package org.lcsim.contrib.JanStrube.tracking;
-
-/**
- * @author jstrube
- *
- */
-public final class Constants {
-    public static final double fieldConversion = 3335.641;
-}

lcsim/src/org/lcsim/contrib/JanStrube/vtxFitter
Fitter.java 1.6 -> 1.7
diff -u -r1.6 -r1.7
--- Fitter.java	7 Aug 2006 21:30:30 -0000	1.6
+++ Fitter.java	15 Aug 2006 23:50:01 -0000	1.7
@@ -1,11 +1,11 @@
 package org.lcsim.contrib.JanStrube.vtxFitter;
 /**
- * @version $Id: Fitter.java,v 1.6 2006/08/07 21:30:30 jstrube Exp $
+ * @version $Id: Fitter.java,v 1.7 2006/08/15 23:50:01 jstrube Exp $
  */
 
 import static java.lang.Math.atan2;
 import static java.lang.Math.sqrt;
-import static org.lcsim.contrib.JanStrube.tracking.Constants.fieldConversion;
+import static org.lcsim.constants.Constants.fieldConversion;
 import static org.lcsim.recon.vertexing.zvtop4.VectorArithmetic.subtract;
 import static org.lcsim.recon.vertexing.zvtop4.VectorArithmetic.add;
 import static org.lcsim.recon.vertexing.zvtop4.VectorArithmetic.dot;
@@ -27,15 +27,15 @@
 	double deltaChi2 = 1e-4;
 	
     double B_z;
-	// position
-	Matrix x_prev;
-	// covariance for x
-	Matrix C_prev;	
+    // position
+    Matrix x_prev;
+    // covariance for x
+    Matrix C_prev;	
     
     HelixSwimmer _swimmer;
 	
-	// chi2
-	double chi2_prev;
+    // chi2
+    double chi2_prev;
 		
 	/**
 	 * Definitions:
@@ -53,37 +53,40 @@
 	 * @param p
 	 */
 	public Fitter(Vertex p, double bz) {
-		this(p
+	    this(p
             , bz
             , new Matrix(new double[][] {{10000, 0, 0}, {0, 10000, 0}, {0, 0, 10000}})
         );
 	}
 
 	public Fitter(Vertex p, double bz, Matrix c0) {
-		unfittedParticle = p;
-        B_z = bz;
-		C_prev = c0;
-        x_prev = new Matrix(p.origin().getCartesianArray(), 3);
-        _swimmer = new HelixSwimmer(bz);
+	    unfittedParticle = p;
+            B_z = bz;
+	    C_prev = c0;
+	    x_prev = new Matrix(p.origin().getCartesianArray(), 3);
+	    _swimmer = new HelixSwimmer(bz);
 	}
 	
 	
 	public Vertex fit() {
-        particle = new Vertex();
-        particle.setOrigin(unfittedParticle.origin());
-        
-		for (Track t : unfittedParticle.getTracks()) {
-			filter(t);
-		}
-		for (Track t : particle.getTracks()) {
-			smoothe(t);
-		}
-        return particle;
+	    particle = new Vertex();
+	    particle.setOrigin(unfittedParticle.origin());
+	    for (int i=0; i<10; ++i) {
+	    	particle._tracks.clear();
+	    	chi2_prev = 0;
+	    	for (Track t : unfittedParticle.getTracks()) {
+	    		System.out.printf("filtering: chi2=%f\n", filter(t));
+	    	}
+	    	for (Track t : particle.getTracks()) {
+	    		System.out.printf("smoothing: %f\n", smoothe(t));
+	    	}
+	    }
+	    return particle;
 	}
     
     // The problem is that double[] is a column vector...
     // Turn everything into a Matrix
-	private double filter(Track t) {
+    private double filter(Track t) {
         // First, get the Matrices at the point of linearization x0, p0
         // a good choice for x0 is the previous best estimate and
         // p0 is the momentum of the track at the POCA to x0
@@ -115,6 +118,7 @@
 
         // cov(x)
         Matrix C_k = C_prev.inverse().plus(A_k.transpose().times(G_kB.times(A_k))).inverse();
+        System.out.printf("C_k: %g\t%g\t%g\n", C_k.get(0, 0), C_k.get(1, 1), C_k.get(2, 2));
         Matrix x_k = C_k.times(C_prev.inverse().times(x_prev).plus(
                  A_k.transpose().times(G_kB).times(m_k.minus(c_k0))));
         Matrix q_k = W_k.times(B_k.transpose().times(G_k.times(m_k.minus(c_k0.plus(A_k.times(x_k))))));
@@ -129,43 +133,43 @@
         double chi2_kf = r_k.transpose().times(G_k.times(r_k)).get(0, 0)
             + x_k.minus(x_prev).transpose().times(C_prev.inverse().times(x_k.minus(x_prev))).get(0, 0);
 		
-        
         double chi2 = chi2_prev + chi2_kf;
-		
         x_prev = x_k;
         C_prev = C_k;
         particle.addTrack(t, chi2_kf);
         particle._origin = new CartesianPoint(x_k.getRowPackedCopy());
+        System.out.println(particle._origin);
         chi2_prev = chi2;
         return chi2;
 	}
 	
-	private void smoothe(Track t) {
-        _swimmer.setTrack(t);
-        // at the moment, particle.origin is the best estimate so far
-        double alpha = _swimmer.getTrackLengthToPoint(particle.origin());
-        SpacePoint p = _swimmer.getMomentumAtLength(alpha);
-        // x is the POCA of the current track to the best estimate
-        Matrix q0 = new Matrix(p.getCartesianArray(), 3);
-
-        SpacePoint x = _swimmer.getPointAtLength(alpha);
-        
-        Matrix A_k = getSpatialDerivativeMatrix(particle.origin(), p, t.getCharge(), B_z);
-        Matrix B_k = getMomentumDerivativeMatrix(particle.origin(), p, t.getCharge(), B_z);
-        
-        Matrix G_k = t.getErrorMatrix().inverse();
-        Matrix W_k = B_k.transpose().times(G_k.times(B_k)).uminus();
-        // the measurement of the track is taken at the POCA to the best estimate
-        Matrix m_k = new Matrix(SpaceMomentum2Parameters(
+	private double smoothe(Track t) {
+	    _swimmer.setTrack(t);
+	    // at the moment, particle.origin is the best estimate so far
+	    double alpha = _swimmer.getTrackLengthToPoint(particle.origin());
+	    SpacePoint p = _swimmer.getMomentumAtLength(alpha);
+	    // x is the POCA of the current track to the best estimate
+	    Matrix q0 = new Matrix(p.getCartesianArray(), 3);
+            SpacePoint x = _swimmer.getPointAtLength(alpha);
+        
+	    Matrix A_k = getSpatialDerivativeMatrix(particle.origin(), p, t.getCharge(), B_z);
+	    Matrix B_k = getMomentumDerivativeMatrix(particle.origin(), p, t.getCharge(), B_z);
+        
+	    Matrix G_k = t.getErrorMatrix().inverse();
+	    Matrix W_k = B_k.transpose().times(G_k.times(B_k)).uminus();
+	    // the measurement of the track is taken at the POCA to the best estimate
+	    Matrix m_k = new Matrix(SpaceMomentum2Parameters(
                 x, p, particle.origin(), t.getCharge(), B_z).getValues()
                 , 5);
-        Matrix c_k0 =  m_k.minus(A_k.times(x_prev).plus(B_k.times(q0)));
+	    Matrix c_k0 =  m_k.minus(A_k.times(x_prev).plus(B_k.times(q0)));
 
-        Matrix q_kN = W_k.times(B_k.transpose().times(G_k)).times(m_k.minus(c_k0.plus(A_k.times(x_prev))));
-		Matrix D_kN = W_k.plus(W_k.times(B_k.transpose().times(G_k.times(A_k.times(C_prev.times(A_k.transpose().times(G_k.times(B_k.times(W_k)))))))));
-		Matrix E_kN = C_prev.times(A_k.transpose().times(G_k.times(B_k.times(W_k)))).uminus();
+	    Matrix q_kN = W_k.times(B_k.transpose().times(G_k)).times(m_k.minus(c_k0.plus(A_k.times(x_prev))));
+	    // Matrix D_kN = W_k.plus(W_k.times(B_k.transpose().times(G_k.times(A_k.times(C_prev.times(A_k.transpose().times(G_k.times(B_k.times(W_k)))))))));
+	    // Matrix E_kN = C_prev.times(A_k.transpose().times(G_k.times(B_k.times(W_k)))).uminus();
+	    return 0;
 	}
-	
+
+
 //	// returns y = A.x
 //	static double[] multiply(Matrix A, double[] x) {
 //		if (A.getColumnDimension() != x.length)

lcsim/src/org/lcsim/contrib/JanStrube/vtxFitter
vertexing.lyx 1.6 -> 1.7
diff -u -r1.6 -r1.7
--- vertexing.lyx	7 Aug 2006 21:30:30 -0000	1.6
+++ vertexing.lyx	15 Aug 2006 23:50:01 -0000	1.7
@@ -229,7 +229,11 @@
 \begin_inset Formula $\vec{P^{0}}$
 \end_inset
 
- is the point of closest approach to the reference point.
+ is the point of closest approach to the reference point 
+\emph on
+in the x-y plane
+\emph default
+.
  Since 
 \begin_inset Formula $|\vec{n}|\equiv1$
 \end_inset
@@ -278,17 +282,46 @@
 \end_layout
 
 \begin_layout Standard
-The unit conversion factor 
-\begin_inset Formula $a$
+The formula that relates the transverse momentum to field, charge and radius
+ is
+\begin_inset Formula \[
+p_{t}=aqB_{z}r,\]
+
 \end_inset
 
- by definition:
-\begin_inset Formula \[
-a=\frac{p_{t}}{qB_{z}r}=\text{const.}\]
 
+\end_layout
+
+\begin_layout Standard
+where 
+\begin_inset Formula $a$
+\end_inset
+
+ is a factor that takes care of the unit conversion.
+ Now in order to compute 
+\begin_inset Formula $a$
 \end_inset
 
+ one has to be careful: It is a unit conversion factor, 
+\emph on
+not
+\emph default
+ a physical quantity.
+ In our computations, we enter unitless numbers rather than physical quantities
+ into the formula above.
+ Therefore the constant 
+\begin_inset Formula $a$
+\end_inset
 
+ 
+\emph on
+collects
+\emph default
+ rather than 
+\emph on
+cancels
+\emph default
+ the units.
 \end_layout
 
 \begin_layout Standard
@@ -324,7 +357,7 @@
 
 \begin_layout Standard
 Now in LCIO, 
-\begin_inset Formula $[p]=\text{GeV}/c$
+\begin_inset Formula $[p_{t}]=\text{GeV}/c$
 \end_inset
 
 , 
@@ -345,9 +378,43 @@
 \end_inset
 
 .
- With this, we have
+\end_layout
+
+\begin_layout Standard
+In compuations we use the following:
+\begin_inset Formula \[
+\tilde{p}_{t}=\tilde{a}\tilde{q}\tilde{B_{z}}\tilde{r}\]
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+using 
+\begin_inset Formula $\tilde{p}_{t}=p_{t}/[p_{t}]$
+\end_inset
+
+, 
+\begin_inset Formula $\tilde{q}=q/[q]$
+\end_inset
+
+, 
+\begin_inset Formula $\tilde{B}_{z}=B_{z}/[B_{z}]$
+\end_inset
+
+ and 
+\begin_inset Formula $\tilde{r}=r/[r]$
+\end_inset
+
+.
+ Now it is easy to compute 
+\begin_inset Formula $\tilde{a}$
+\end_inset
+
+:
 \begin_inset Formula \[
-a=\frac{p_{t}}{qB_{z}r}=\frac{10^{9}\text{kg}\mathrm{m}^{2}}{c\mathrm{A}\mathrm{s}^{3}}\frac{\mathrm{A}\mathrm{s}^{2}}{\text{kg}}\frac{1000}{\mathrm{m}}=\frac{10^{12}\mathrm{m}/\mathrm{s}}{c}\approx\frac{10^{12}\mathrm{m}/\mathrm{s}}{3\times10^{8}\mathrm{m}/\mathrm{s}}=3335.641\]
+\tilde{a}=\frac{[q][B_{z}][r]}{[p_{t}]}=\frac{c\mathrm{A}\mathrm{s}^{3}}{10^{9}\text{kg}\mathrm{m}^{2}}\frac{\text{kg}}{\mathrm{A}\mathrm{s}^{2}}\frac{\mathrm{m}}{1000}=\frac{c}{10^{12}\mathrm{m}/\mathrm{s}}\approx2.9979\times10^{-4}\]
 
 \end_inset
 

lcsim/src/org/lcsim/contrib/JanStrube/vtxFitter
Vertex.java 1.2 -> 1.3
diff -u -r1.2 -r1.3
--- Vertex.java	7 Aug 2006 21:30:30 -0000	1.2
+++ Vertex.java	15 Aug 2006 23:50:02 -0000	1.3
@@ -16,12 +16,12 @@
 
 public class Vertex {
 	SpacePoint _origin;
-    private SpacePoint _originError;
-    private Map<Track, Double> _tracks;    
+    SpacePoint _originError;
+    Map<Track, Double> _tracks;    
     
     double _chi2;
-    private double _ndf;
-    private Matrix _covMatrix;
+    double _ndf;
+    Matrix _covMatrix;
     
 	public Vertex() {
 		_origin = new SpacePoint();

lcsim/test/org/lcsim/contrib/JanStrube/tracking
FitterTest.java 1.1 -> 1.2
diff -u -r1.1 -r1.2
--- FitterTest.java	7 Aug 2006 21:30:31 -0000	1.1
+++ FitterTest.java	15 Aug 2006 23:50:05 -0000	1.2
@@ -46,9 +46,10 @@
         SpacePoint mom1 = new CartesianPoint(1, 2, 3);
         SpacePoint mom2 = new CartesianPoint(2, -2, 1.5);
         SpacePoint mom3 = new CartesianPoint(-1, 2, 2);
-        particle.addTrack(fac.getTrack(mom1, pos, pos, -1, new Random(), false), 0);
-        particle.addTrack(fac.getTrack(mom2, pos, pos, 1, new Random(), false), 0);
-        particle.addTrack(fac.getTrack(mom3, pos, pos, -1, new Random(), false), 0);
+        SpacePoint ref = new CartesianPoint(3, 6, 9);
+        particle.addTrack(fac.getTrack(mom1, ref, ref, -1, new Random(), false), 0);
+        particle.addTrack(fac.getTrack(mom2, ref, ref, 1, new Random(), false), 0);
+        particle.addTrack(fac.getTrack(mom3, ref, ref, -1, new Random(), false), 0);
         particle.setOrigin(pos);
         Fitter f = new Fitter(particle, 5.0);
         f.fit();
CVSspam 0.2.8