4 added files
lcsim/test/org/lcsim/recon/tracking/trfzp
diff -N ZPlaneEloss_t.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ ZPlaneEloss_t.java 11 Sep 2007 20:23:50 -0000 1.1
@@ -0,0 +1,103 @@
+/*
+ * ZPlaneEloss_t.java
+ *
+ * Created on September 10, 2007, 11:05 AM
+ *
+ * $Id: ZPlaneEloss_t.java,v 1.1 2007/09/11 20:23:50 ngraf Exp $
+ */
+
+package org.lcsim.recon.tracking.trfzp;
+
+import junit.framework.TestCase;
+import org.lcsim.recon.tracking.trfbase.ETrack;
+import org.lcsim.recon.tracking.trfbase.PropDir;
+import org.lcsim.recon.tracking.trfbase.Surface;
+import org.lcsim.recon.tracking.trfbase.TrackError;
+import org.lcsim.recon.tracking.trfbase.TrackVector;
+import org.lcsim.recon.tracking.trfeloss.DeDx;
+import org.lcsim.recon.tracking.trfeloss.DeDxFixed;
+
+
+import static java.lang.Math.abs;
+import static java.lang.Math.sqrt;
+/**
+ *
+ * @author Norman Graf
+ */
+public class ZPlaneEloss_t extends TestCase
+{
+ private boolean debug;
+ /** Creates a new instance of ZPlaneEloss_t */
+ public void testZPlaneEloss()
+ {
+ String component = "ZPlaneEloss";
+ String ok_prefix = component + " (I): ";
+ String error_prefix = component + " test (E): ";
+
+ if(debug) System.out.println( ok_prefix
+ + "---------- Testing component " + component
+ + ". ----------" );
+
+ //********************************************************************
+ TrackVector trv = new TrackVector();
+ trv.set(0,1.0);
+ trv.set(1, 1.0);
+ trv.set(2, 0.1);
+ trv.set(3, 0.2);
+ trv.set(4, 0.0);
+
+ double density = 1.0; // g/cm^3
+ double thickness = 1.0; // cm
+
+ DeDx dedx = new DeDxFixed(density);
+ ZPlaneEloss passIt = new ZPlaneEloss(thickness, dedx);
+
+ TrackError initialError = new TrackError();
+
+ Surface srf = new SurfZPlane( 20.0 );
+
+ ETrack tmpTrack = new ETrack( srf, trv, initialError );
+
+ tmpTrack.setError( initialError );
+
+ TrackVector initialVector = tmpTrack.vector();
+ passIt.interact_dir( tmpTrack, PropDir.FORWARD );
+
+ TrackError finalError = tmpTrack.error();
+ TrackVector finalVector = tmpTrack.vector();
+
+ double particleMass = 0.13957;
+ double ptmax = 10000.; // GeV
+
+ double pinv = abs(trv.get(4));
+
+ // make sure pinv is greater than a threshold (1/ptmax)
+ // in this case assume q = 1, otherwise q = q/pt/abs(q/pt)
+
+ int sign = 1;
+ if(pinv < 1./ptmax)
+ pinv = 1./ptmax;
+ else
+ sign = (int) (trv.get(4)/abs(trv.get(4)));
+
+ double initialEnergy = sqrt(1./pinv/pinv + particleMass*particleMass);
+ double finalEnergy = initialEnergy;
+
+ double d = thickness * sqrt(1. + trv.get(2)*trv.get(2) + trv.get(3)*trv.get(3));
+
+ dedx.loseEnergy(finalEnergy, d);
+ double sigmaEnergy = dedx.sigmaEnergy(initialEnergy, d);
+
+ if(finalEnergy<particleMass) finalEnergy=initialEnergy;
+
+ // now evaluate the final q/p and error(4,4)
+ double finalQoverP = sign/sqrt(finalEnergy*finalEnergy - particleMass*particleMass);
+ double finalEr = sigmaEnergy*trv.get(4)*trv.get(4);
+ finalEr *= finalEr;
+
+ assertTrue(abs(finalVector.get(4)-finalQoverP)<0.00001);
+ assertTrue(abs(finalError.get(4,4)-finalEr)<0.00001);
+
+ }
+
+}
lcsim/src/org/lcsim/recon/tracking/trfxyp
diff -N XYPlaneEloss.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ XYPlaneEloss.java 11 Sep 2007 20:23:50 -0000 1.1
@@ -0,0 +1,167 @@
+/*
+ * XYPlaneEloss.java
+ *
+ * Created on September 10, 2007, 11:17 AM
+ *
+ * $Id: XYPlaneEloss.java,v 1.1 2007/09/11 20:23:50 ngraf Exp $
+ */
+
+package org.lcsim.recon.tracking.trfxyp;
+
+import org.lcsim.recon.tracking.trfbase.ETrack;
+import org.lcsim.recon.tracking.trfbase.Interactor;
+import org.lcsim.recon.tracking.trfbase.PropDir;
+import org.lcsim.recon.tracking.trfbase.Propagator;
+import org.lcsim.recon.tracking.trfbase.Surface;
+import org.lcsim.recon.tracking.trfbase.TrackError;
+import org.lcsim.recon.tracking.trfbase.TrackVector;
+import org.lcsim.recon.tracking.trfeloss.DeDx;
+import org.lcsim.recon.tracking.trfutil.Assert;
+
+import static org.lcsim.recon.tracking.trfxyp.SurfXYPlane.IV;
+import static org.lcsim.recon.tracking.trfxyp.SurfXYPlane.IZ;
+import static org.lcsim.recon.tracking.trfxyp.SurfXYPlane.IDVDU;
+import static org.lcsim.recon.tracking.trfxyp.SurfXYPlane.IDZDU;
+import static org.lcsim.recon.tracking.trfxyp.SurfXYPlane.IQP;
+
+import static java.lang.Math.abs;
+import static java.lang.Math.sqrt;
+
+/**
+ *
+ * @author Norman Graf
+ */
+public class XYPlaneEloss extends Interactor
+{
+ //attributes
+
+ private double _thickness;
+ private DeDx _dedx;
+
+ /** Creates a new instance of XYPlaneEloss */
+ public XYPlaneEloss(double thickness, DeDx dedx)
+ {
+ _thickness = thickness;
+ _dedx = dedx;
+ }
+
+
+ public void interact(ETrack tre)
+ {
+ }
+
+
+ public void interact_dir(ETrack theTrack, PropDir direction )
+ {
+ // This can only be used with xy-planes... check that we have one..
+
+ Surface srf = theTrack.surface();
+ Assert.assertTrue( srf instanceof SurfXYPlane);
+
+ // Reduced direction must be forward or backward.
+
+ // Reduced direction must be forward or backward.
+ Propagator.reduceDirection(direction);
+ Assert.assertTrue(direction == PropDir.FORWARD || direction == PropDir.BACKWARD);
+
+
+ TrackError cleanError = theTrack.error();
+
+
+ TrackError newError = new TrackError(cleanError);
+
+ TrackVector theVec = theTrack.vector();
+ TrackVector newVector = new TrackVector(theVec);
+
+ double pionMass = 0.13957; // GeV
+ double ptmax = 10000.; // GeV
+
+
+ double pinv = abs(theVec.get(IQP));
+ double dvdu = theVec.get(IDVDU);
+ double dzdu = theVec.get(IDZDU);
+ double lfac = sqrt(1. + dvdu*dvdu + dzdu*dzdu);
+
+ // make sure pinv is greater than a threshold (1/ptmax)
+ // in this case assume q = 1, otherwise q = q/p/abs(q/p)
+
+ double sign = 1.;
+ if ( pinv < 1./ptmax )
+ pinv = 1./ptmax;
+ else
+ sign = theVec.get(IQP) > 0. ? 1. : -1.;
+
+ // Evaluate the initial energy assuming the particle is a pion.
+
+ double trackEnergy = sqrt(1./(pinv*pinv) + pionMass*pionMass);
+
+ double trueLength = lfac * _thickness;
+
+ // assume the energy loss distribution to be Gaussian
+ double stdEnergy = _dedx.sigmaEnergy(trackEnergy, trueLength);
+
+ if(direction == PropDir.BACKWARD) trueLength = -trueLength;
+ _dedx.loseEnergy(trackEnergy, trueLength);
+
+ double newMomentum = trackEnergy>pionMass ?
+ sqrt(trackEnergy*trackEnergy - pionMass*pionMass) :
+ 1./pinv;
+
+ // Only vec(IQP) changes due to E loss.
+ newVector.set(IQP, sign/newMomentum);
+
+ // Also only error(IQP,IQP) changes.
+ double stdVec = theVec.get(IQP)*theVec.get(IQP)*stdEnergy;
+ stdVec *= stdVec;
+ newError.set(IQP, IQP, newError.get(IQP, IQP) + stdVec);
+
+ // Store modified track parameters.
+ theTrack.setVectorAndKeepDirection(newVector);
+ theTrack.setError( newError );
+
+ }
+
+
+ /**
+ * Return the thickness of material in the z plane.
+ *
+ * @return The thickness of the energy loss material.
+ */
+ public double thickness()
+ {
+ return _thickness;
+ }
+
+ /**
+ *Return the energy loss model used in this Interactor.
+ *
+ * @return The DeDx class representing energy loss.
+ */
+ public DeDx dEdX()
+ {
+ return _dedx; //cng shallow copy!
+ }
+
+ /**
+ *Make a clone of this object.
+ *
+ * @return A Clone of this instance.
+ */
+ public Interactor newCopy()
+ {
+ return new XYPlaneEloss(_thickness, _dedx);
+ }
+
+ /**
+ *output stream
+ *
+ * @return A String representation of this instance.
+ */
+ public String toString()
+ {
+ return "XYPlaneEloss with thickness "+_thickness+" and energy loss "+_dedx;
+ }
+
+
+
+}
lcsim/test/org/lcsim/recon/tracking/trfxyp
diff -N XYPlaneEloss_t.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ XYPlaneEloss_t.java 11 Sep 2007 20:23:51 -0000 1.1
@@ -0,0 +1,105 @@
+/*
+ * XYPlaneEloss_t.java
+ *
+ * Created on September 10, 2007, 11:33 AM
+ *
+ * $Id: XYPlaneEloss_t.java,v 1.1 2007/09/11 20:23:51 ngraf Exp $
+ */
+
+package org.lcsim.recon.tracking.trfxyp;
+
+import junit.framework.TestCase;
+import org.lcsim.recon.tracking.trfbase.ETrack;
+import org.lcsim.recon.tracking.trfbase.PropDir;
+import org.lcsim.recon.tracking.trfbase.Surface;
+import org.lcsim.recon.tracking.trfbase.TrackError;
+import org.lcsim.recon.tracking.trfbase.TrackVector;
+import org.lcsim.recon.tracking.trfeloss.DeDx;
+import org.lcsim.recon.tracking.trfeloss.DeDxFixed;
+import org.lcsim.recon.tracking.trfzp.ZPlaneEloss;
+
+
+import static java.lang.Math.abs;
+import static java.lang.Math.sqrt;
+
+/**
+ *
+ * @author Norman Graf
+ */
+public class XYPlaneEloss_t extends TestCase
+{
+ private boolean debug;
+
+ /** Creates a new instance of XYPlaneEloss */
+ public void testXYPlaneEloss()
+ {
+ String component = "XYPlaneEloss";
+ String ok_prefix = component + " (I): ";
+ String error_prefix = component + " test (E): ";
+
+ if(debug) System.out.println( ok_prefix
+ + "---------- Testing component " + component
+ + ". ----------" );
+
+ //********************************************************************
+ TrackVector trv = new TrackVector();
+ trv.set(0,1.0);
+ trv.set(1, 1.0);
+ trv.set(2, 0.1);
+ trv.set(3, 0.2);
+ trv.set(4, 0.0);
+
+ double density = 1.0; // g/cm^3
+ double thickness = 1.0; // cm
+
+ DeDx dedx = new DeDxFixed(density);
+ XYPlaneEloss passIt = new XYPlaneEloss(thickness, dedx);
+
+ TrackError initialError = new TrackError();
+
+ Surface srf = new SurfXYPlane( 20.0, 1. );
+
+ ETrack tmpTrack = new ETrack( srf, trv, initialError );
+
+ tmpTrack.setError( initialError );
+
+ TrackVector initialVector = tmpTrack.vector();
+ passIt.interact_dir( tmpTrack, PropDir.FORWARD );
+
+ TrackError finalError = tmpTrack.error();
+ TrackVector finalVector = tmpTrack.vector();
+
+ double particleMass = 0.13957;
+ double ptmax = 10000.; // GeV
+
+ double pinv = abs(trv.get(4));
+
+ // make sure pinv is greater than a threshold (1/ptmax)
+ // in this case assume q = 1, otherwise q = q/pt/abs(q/pt)
+
+ int sign = 1;
+ if(pinv < 1./ptmax)
+ pinv = 1./ptmax;
+ else
+ sign = (int) (trv.get(4)/abs(trv.get(4)));
+
+ double initialEnergy = sqrt(1./pinv/pinv + particleMass*particleMass);
+ double finalEnergy = initialEnergy;
+
+ double d = thickness * sqrt(1. + trv.get(2)*trv.get(2) + trv.get(3)*trv.get(3));
+
+ dedx.loseEnergy(finalEnergy, d);
+ double sigmaEnergy = dedx.sigmaEnergy(initialEnergy, d);
+
+ if(finalEnergy<particleMass) finalEnergy=initialEnergy;
+
+ // now evaluate the final q/p and error(4,4)
+ double finalQoverP = sign/sqrt(finalEnergy*finalEnergy - particleMass*particleMass);
+ double finalEr = sigmaEnergy*trv.get(4)*trv.get(4);
+ finalEr *= finalEr;
+
+ assertTrue(abs(finalVector.get(4)-finalQoverP)<0.00001);
+ assertTrue(abs(finalError.get(4,4)-finalEr)<0.00001);
+
+ }
+}
\ No newline at end of file
lcsim/src/org/lcsim/recon/tracking/trfzp
diff -N ZPlaneEloss.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ ZPlaneEloss.java 11 Sep 2007 20:23:51 -0000 1.1
@@ -0,0 +1,167 @@
+/*
+ * ZPlaneEloss.java
+ *
+ * Created on September 10, 2007, 10:40 AM
+ *
+ * $Id: ZPlaneEloss.java,v 1.1 2007/09/11 20:23:51 ngraf Exp $
+ */
+
+package org.lcsim.recon.tracking.trfzp;
+
+import org.lcsim.recon.tracking.trfbase.ETrack;
+import org.lcsim.recon.tracking.trfbase.Interactor;
+import org.lcsim.recon.tracking.trfbase.PropDir;
+import org.lcsim.recon.tracking.trfbase.Propagator;
+import org.lcsim.recon.tracking.trfbase.Surface;
+import org.lcsim.recon.tracking.trfbase.TrackError;
+import org.lcsim.recon.tracking.trfbase.TrackVector;
+import org.lcsim.recon.tracking.trfeloss.DeDx;
+import org.lcsim.recon.tracking.trfutil.Assert;
+
+import static java.lang.Math.abs;
+import static java.lang.Math.sqrt;
+
+import static org.lcsim.recon.tracking.trfzp.SurfZPlane.IX;
+import static org.lcsim.recon.tracking.trfzp.SurfZPlane.IY;
+import static org.lcsim.recon.tracking.trfzp.SurfZPlane.IDXDZ;
+import static org.lcsim.recon.tracking.trfzp.SurfZPlane.IDYDZ;
+import static org.lcsim.recon.tracking.trfzp.SurfZPlane.IQP;
+
+/**
+ * Class for modifying the covariance matrix of a track which
+ * has been propagated to an InteractingLayer containing
+ * a LayerZPlane.
+ * @author Norman Graf
+ */
+public class ZPlaneEloss extends Interactor
+{
+
+ //attributes
+
+ private double _thickness;
+ private DeDx _dedx;
+
+ /** Creates a new instance of ZPlaneEloss */
+ // constructor. The Interactor is constructed with the
+ // appropriate thickness and De/dx class.
+ public ZPlaneEloss(double thickness, DeDx dedx)
+ {
+ _thickness = thickness;
+ _dedx = dedx;
+ }
+
+ public void interact(ETrack tre)
+ {
+ }
+
+ public void interact_dir( ETrack theTrack, PropDir direction )
+ {
+
+ // This can only be used with z-planes... check that we have one..
+ Surface srf = theTrack.surface();
+ Assert.assertTrue( srf instanceof SurfZPlane);
+
+
+ // Reduced direction must be forward or backward.
+ Propagator.reduceDirection(direction);
+ Assert.assertTrue(direction == PropDir.FORWARD || direction == PropDir.BACKWARD);
+
+
+ TrackError cleanError = theTrack.error();
+
+
+ TrackError newError = new TrackError(cleanError);
+
+ TrackVector theVec = theTrack.vector();
+ TrackVector newVector = new TrackVector(theVec);
+
+ double pionMass = 0.13957; // GeV
+ double ptmax = 10000.; // GeV
+
+ double pinv = abs(theVec.get(IQP));
+ double dxdz = theVec.get(IDXDZ);
+ double dydz = theVec.get(IDYDZ);
+ double lfac = sqrt(1. + dxdz*dxdz + dydz*dydz);
+
+ // make sure pinv is greater than a threshold (1/ptmax)
+ // in this case assume q = 1, otherwise q = q/p/abs(q/p)
+
+ double sign = 1.;
+ if ( pinv < 1./ptmax )
+ pinv = 1./ptmax;
+ else
+ sign = theVec.get(IQP) > 0. ? 1. : -1.;
+
+ // Evaluate the initial energy assuming the particle is a pion.
+
+ double trackEnergy = sqrt(1./(pinv*pinv) + pionMass*pionMass);
+
+ double trueLength = lfac * _thickness;
+
+ // assume the energy loss distribution to be Gaussian
+ double stdEnergy = _dedx.sigmaEnergy(trackEnergy, trueLength);
+
+ if(direction == PropDir.BACKWARD) trueLength = -trueLength;
+ _dedx.loseEnergy(trackEnergy, trueLength);
+
+ double newMomentum = trackEnergy>pionMass ?
+ sqrt(trackEnergy*trackEnergy - pionMass*pionMass) :
+ 1./pinv;
+
+ // Only vec(IQP) changes due to E loss.
+ newVector.set(IQP, sign/newMomentum);
+
+ // Also only error(IQP,IQP) changes.
+ double stdVec = theVec.get(IQP)*theVec.get(IQP)*stdEnergy;
+ stdVec *= stdVec;
+ newError.set(IQP, IQP, newError.get(IQP, IQP)+stdVec);
+
+ // Store modified track parameters.
+ theTrack.setVectorAndKeepDirection(newVector);
+ theTrack.setError( newError );
+
+ }
+
+
+ /**
+ * Return the thickness of material in the z plane.
+ *
+ * @return The thickness of the energy loss material.
+ */
+ public double thickness()
+ {
+ return _thickness;
+ }
+
+ /**
+ *Return the energy loss model used in this Interactor.
+ *
+ * @return The DeDx class representing energy loss.
+ */
+ public DeDx dEdX()
+ {
+ return _dedx; //cng shallow copy!
+ }
+
+ /**
+ *Make a clone of this object.
+ *
+ * @return A Clone of this instance.
+ */
+ public Interactor newCopy()
+ {
+ return new ZPlaneEloss(_thickness, _dedx);
+ }
+
+ /**
+ *output stream
+ *
+ * @return A String representation of this instance.
+ */
+ public String toString()
+ {
+ return "ZPlaneEloss with thickness "+_thickness+" and energy loss "+_dedx;
+ }
+
+
+}
CVSspam 0.2.8