hps-java/src/main/java/org/lcsim/hps/recon/tracking
diff -u -r1.14 -r1.15
--- TrackUtils.java 23 Jan 2013 17:47:29 -0000 1.14
+++ TrackUtils.java 9 Feb 2013 00:34:16 -0000 1.15
@@ -3,6 +3,7 @@
//--- hep ---//
import hep.physics.matrix.SymmetricMatrix;
import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Matrix;
import hep.physics.vec.Hep3Vector;
import hep.physics.vec.VecOp;
import java.util.ArrayList;
@@ -14,16 +15,19 @@
import org.lcsim.detector.solids.Point3D;
import org.lcsim.detector.solids.Polygon3D;
import org.lcsim.detector.tracker.silicon.SiSensor;
+import org.lcsim.event.RawTrackerHit;
//--- org.lcsim ---//
import org.lcsim.event.Track;
import org.lcsim.fit.helicaltrack.*;
import org.lcsim.hps.users.phansson.WTrack;
+import org.lcsim.recon.tracking.seedtracker.SeedCandidate;
+import org.lcsim.recon.tracking.seedtracker.SeedTrack;
/**
*
* @author Omar Moreno <[log in to unmask]>
- * @version $Id: TrackUtils.java,v 1.14 2013/01/23 17:47:29 phansson Exp $
+ * @version $Id: TrackUtils.java,v 1.15 2013/02/09 00:34:16 phansson Exp $
* TODO: Switch to JLab coordinates
*/
@@ -264,7 +268,7 @@
/*
* Use code in WTrack to find the iterative solution to the interception
*/
- boolean debug = true;
+ boolean debug = false;
WTrack wtrack = new WTrack(helfit,bfield,true); // B-field sign is flipped so flip!
if(debug) System.out.printf("getHelixPlaneIntercept:find intercept between plane defined by point on plane %s, unit vec %s, bfield %.3f, and WTrack \n%s \n",point_on_plane.toString(),unit_vec_normal_to_plane.toString(), bfield,wtrack.toString());
Hep3Vector intercept_point = wtrack.getHelixAndPlaneIntercept(wtrack,point_on_plane, unit_vec_normal_to_plane, new BasicHep3Vector(0,0,1));
@@ -449,8 +453,91 @@
//System.out.printf("calculateTrackHitResidual: resz=%f eresz=%f dz_res=%f msdz=%f \n",resz,dz_res,Math.sqrt(dz_res2),msdz);
- return residuals;
- }
+ return residuals;
+ }
+
+
+ public static Map<String,Double> calculateLocalTrackHitResiduals(Track track, HelicalTrackHit hth, HelicalTrackStrip strip, double bFieldInZ) {
+
+ SeedTrack st = (SeedTrack) track;
+ SeedCandidate seed = st.getSeedCandidate();
+ HelicalTrackFit _trk = seed.getHelix();
+ Map<HelicalTrackHit,MultipleScatter> msmap = seed.getMSMap();
+ double msdrdphi = msmap.get(hth).drphi();
+ double msdz = msmap.get(hth).dz();
+ return calculateLocalTrackHitResiduals(_trk, strip, msdrdphi, msdz, bFieldInZ);
+ }
+
+ public static Map<String,Double> calculateLocalTrackHitResiduals(HelicalTrackFit _trk, HelicalTrackStrip strip, double msdrdphi, double msdz, double bFieldInZ) {
+
+
+ boolean debug = false;
+ boolean includeMS = true;
+
+
+ Hep3Vector u = strip.u();
+ Hep3Vector v = strip.v();
+ Hep3Vector w = strip.w();
+ Hep3Vector eta = VecOp.cross(u,v);
+ Hep3Vector corigin = strip.origin();
+ String side = SvtUtils.getInstance().isTopLayer((SiSensor)((RawTrackerHit)strip.rawhits().get(0)).getDetectorElement()) ? "top" : "bottom";
+
+ double bfield = bFieldInZ; //Math.abs(this._bfield.z());
+
+ //Find interception with plane that the strips belongs to
+ Hep3Vector trkpos = TrackUtils.getHelixPlaneIntercept(_trk, strip, bfield);
+
+ if(debug) {
+ System.out.printf("calculateLocalTrackHitResiduals: found interception point at %s \n",trkpos.toString());
+ }
+
+
+ if(Double.isNaN(trkpos.x()) || Double.isNaN(trkpos.y()) || Double.isNaN(trkpos.z())) {
+ System.out.printf("calculateLocalTrackHitResiduals: failed to get interception point (%s) \n",trkpos.toString());
+ System.out.printf("calculateLocalTrackHitResiduals: track params\n%s\n",_trk.toString());
+ System.out.printf("calculateLocalTrackHitResiduals: track pT=%.3f chi2=[%.3f][%.3f] \n",_trk.pT(bfield),_trk.chisq()[0],_trk.chisq()[1]);
+ trkpos = TrackUtils.getHelixPlaneIntercept(_trk, strip, bfield);
+ System.exit(1);
+ }
+
+ double xint = trkpos.x();
+ double phi0 = _trk.phi0();
+ double R = _trk.R();
+ double s = HelixUtils.PathToXPlane(_trk, xint, 0, 0).get(0);
+ double phi = -s/R + phi0;
+
+
+ Hep3Vector mserr = new BasicHep3Vector(msdrdphi * Math.sin(phi), msdrdphi * Math.sin(phi), msdz);
+ double msuError = VecOp.dot(mserr, u);
+
+ Hep3Vector vdiffTrk = VecOp.sub(trkpos, corigin);
+ TrackerHitUtils thu = new TrackerHitUtils(debug);
+ Hep3Matrix trkToStrip = thu.getTrackToStripRotation(strip);
+ Hep3Vector vdiff = VecOp.mult(trkToStrip, vdiffTrk);
+
+
+ double umc = vdiff.x();
+ double vmc = vdiff.y();
+ double wmc = vdiff.z();
+ double umeas = strip.umeas();
+ double uError = strip.du();
+ double vmeas = 0;
+ double vError = (strip.vmax() - strip.vmin()) / Math.sqrt(12);
+ double wmeas = 0;
+ double wError = 10.0/Math.sqrt(12); //0.001;
+
+ Map<String,Double> res = new HashMap<String,Double>();
+ res.put("ures", umeas-umc);
+ res.put("ureserr", includeMS ? Math.sqrt(uError * uError + msuError * msuError) : uError);
+ res.put("vres", vmeas-vmc);
+ res.put("vreserr", vError);
+ res.put("wres", wmeas-wmc);
+ res.put("wreserr", wError);
+
+ res.put("vdiffTrky",vdiffTrk.y());
+
+ return res;
+ }
}