Author: [log in to unmask]
Date: Fri Aug 12 14:27:22 2016
New Revision: 4464
Log:
savepoint
Added:
java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/ECal_GBLRefitterDriver.java
Modified:
java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/ECAL_HelicalTrackHitDriver.java
java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/MakeGblTracks.java
Modified: java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/ECAL_HelicalTrackHitDriver.java
=============================================================================
--- java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/ECAL_HelicalTrackHitDriver.java (original)
+++ java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/ECAL_HelicalTrackHitDriver.java Fri Aug 12 14:27:22 2016
@@ -41,7 +41,13 @@
import org.lcsim.recon.tracking.digitization.sisim.TrackerHitType;
// For adding ECal hits to HelicalTrackHits
-import org.hps.recon.tracking.ecal.HelicalTrack3DHit;
+//import org.hps.recon.tracking.ecal.HelicalTrack3DHit;
+import org.lcsim.fit.helicaltrack.HelicalTrackHit;
+import org.lcsim.fit.helicaltrack.HelicalTrack3DHit;
+import org.lcsim.geometry.subdetector.HPSEcal3;
+import org.lcsim.geometry.subdetector.HPSEcal3.NeighborMap;
+import org.lcsim.event.CalorimeterHit;
+import org.lcsim.event.Cluster;
/**
* Driver used to create stereo hits from clusters.
@@ -63,6 +69,7 @@
private double maxDt = -99; // max time difference between the two hits in a cross
private double clusterAmplitudeCut = -99; // cluster amplitude cut
private String _subdetectorName = "Tracker";
+ // private String _subdetectorName2 = "";
private final Map<String, String> _stereomap = new HashMap<String, String>();
private List<SvtStereoLayer> stereoLayers = null;
private final List<String> _colnames = new ArrayList<String>();
@@ -71,6 +78,10 @@
private final String _axialname = "AxialTrackHits";
private final String _axialmcrelname = "AxialTrackHitsMCRelations";
private boolean rejectGhostHits = false;
+
+ // ECal cluster collection in event to create helicalTrackHits from
+ private final String _ecalclustername = "EcalClustersCorr";
+
public enum LayerGeometryType {
@@ -219,8 +230,8 @@
/////// For making HelicalTrackHits out of ECal hits
//*****
- List<HelicalTrack3DHit> ecalhits = new ArrayList<>();
- List<LCRelation> ecalmcrelations = new ArrayList<LCRelation>();
+ List<HelicalTrack3DHit> helicalecalhits = new ArrayList<>();
+ List<LCRelation> helicalecalmcrelations = new ArrayList<LCRelation>();
//*****
////////////////
@@ -288,9 +299,9 @@
} else {
// If not a 1D strip hit, make a pixel hit
// This should be removed as it is never used.
- HelicalTrackHit hit3d = this.makeDigi3DHit(hit);
- helhits.add(hit3d);
- hitrelations.add(new MyLCRelation(hit3d, hit));
+ // HelicalTrackHit hit3d = this.makeDigi3DHit(hit);
+ // helhits.add(hit3d);
+ // hitrelations.add(new MyLCRelation(hit3d, hit));
}
}
@@ -460,23 +471,66 @@
System.out.printf("%s: added %d stereo hits from %s collection \n", this.getClass().getSimpleName(), helicalTrackCrosses.size(), _colname);
}
- /*
+ /* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ *
* Look for ECal hits in the event to add to HelicalTrackHits
* These are HelicalTrack3DHits rather than HelicalTrack2DHits
*
**/
- // if (!event.hasCollection(CalorimeterHit.class, _colname)) {
- // if (_debug) {
- // System.out.println("Event: " + event.getRunNumber() + " does not contain the collection " + _colname);
- // }
- // continue;
- // }
-
-
-
-
-
-
+// if (!event.hasCollection(Cluster.class, _colname)) {
+// if (_debug) {
+// System.out.println("Event: " + event.getRunNumber() + " does not contain the collection " + _colname);
+// }
+// continue;
+// }
+
+ // List<Cluster> clusters = event.get(Cluster.class, _colname);
+
+// List<Cluster> clusters = event.get(Cluster.class, _ecalclustername);
+// if (_debug) {
+// System.out.printf("%s: found %d Clusters\n", this.getClass().getSimpleName(), clusters.size());
+// }
+
+// for (Cluster EcalCluster : clusters) {
+
+// HelicalTrackHit ecal_hit3d = this.makeEcal3DHit(EcalCluster);
+// helhits.add(ecal_hit3d);
+// hitrelations.add(new MyLCRelation(ecal_hit3d, EcalCluster));
+// }
+
+
+
+ /*
+ List<CalorimeterHit> ecalhitlist;
+ for (Cluster cluster : clusters) {
+
+ if(cluster.getCalorimeterHits().size() != isEmpty()){
+
+
+
+
+ }
+ }
+ List<CalorimeterHit> ecalhitlist = event.get(Cluster.class, _colname).getCalorimeterHits();
+ // List<Cluster> ecalhitlist = event.get(Cluster.class, _colname);
+ if (_debug) {
+ System.out.printf("%s: found %d Clusters\n", this.getClass().getSimpleName(), ecalhitlist.size());
+ }
+
+ for (Cluster ecal_hit : ecalhitlist) {
+
+ HelicalTrackHit ecal_hit3d = this.makeEcal3DHit(ecal_hit);
+ helhits.add(ecal_hit3d);
+ hitrelations.add(new MyLCRelation(ecal_hit3d, ecal_hit));
+ }
+ */
+
+ /*
+ *
+ * ECAL HITS !!!!!!!!!!!!!!!
+ *
+ *
+ **/
@@ -570,6 +624,30 @@
}
}
+
+ /* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ *
+ * Make HelicalTrack3DHits from ECal Clusters
+ *
+ * !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ */
+
+// private HelicalTrackHit makeEcal3DHit(Cluster c) {
+//private void makeEcal3DHit(Cluster c) {
+
+ // double z1 = h.getHitSegment().getEndPoint().x();
+ // double z2 = h.getHitSegment().getStartPoint().x();//x is the non-bend direction in the JLAB frame
+ // double zmin = Math.min(z1, z2); TODO: DETERMINE WHAT THIS NEEDS TO BE SET TO, USING THE ACTUAL GEOMETRY CLASSES
+ // double zmax = Math.max(z1, z2);
+// IDetectorElement de = c.getCalorimeterHits().get(0).getDetectorElement();
+
+//HelicalTrack3DHit ecal_hit = new HelicalTrack3DHit();
+// HelicalTrackHit ecal_hit = new HelicalTrackHit(c.getCalorimeterHits().get(0).getPositionVec(),
+// Matrix(c.getPositionError()), c.getEnergy(), c.getCalorimeterHits().get(0).getTime(),
+// null, _ID.getName(de), _ID.getLayer(de), _ID.getBarrelEndcapFlag(de));
+
+// return ecal_hit;
+// }
/*
* Make HelicalTrack2DHits from SiTrackerHitStrip1D...note that these HelicalTrack2DHits
Added: java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/ECal_GBLRefitterDriver.java
=============================================================================
--- java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/ECal_GBLRefitterDriver.java (added)
+++ java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/ECal_GBLRefitterDriver.java Fri Aug 12 14:27:22 2016
@@ -0,0 +1,184 @@
+package org.hps.recon.tracking.gbl;
+
+import java.util.ArrayList;
+import java.util.HashMap;
+import java.util.HashSet;
+import java.util.List;
+import java.util.Map;
+import java.util.Set;
+import org.apache.commons.math3.util.Pair;
+import org.hps.recon.tracking.MaterialSupervisor;
+import org.hps.recon.tracking.MultipleScattering;
+import org.hps.recon.tracking.TrackUtils;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.LCRelation;
+import org.lcsim.event.RelationalTable;
+import org.lcsim.event.Track;
+import org.lcsim.event.TrackerHit;
+import org.lcsim.event.base.BaseLCRelation;
+import org.lcsim.geometry.Detector;
+import org.lcsim.lcio.LCIOConstants;
+import org.lcsim.util.Driver;
+
+// 4 ECal
+import org.lcsim.event.Cluster;
+
+
+/**
+ * A Driver which refits tracks using GBL w/ ECal measurements (taken from ECalClustersCorr). Does not require GBL collections to
+ * be present in the event.
+ */
+public class ECal_GBLRefitterDriver extends Driver {
+
+ private String inputCollectionName = "MatchedTracks";
+ private String ECalCollectionName = "EcalClustersCorr";
+ private String outputCollectionName = "GBLTracks";
+ private String trackRelationCollectionName = "MatchedToGBLTrackRelations";
+
+ private double bfield;
+ private final MultipleScattering _scattering = new MultipleScattering(new MaterialSupervisor());
+ private boolean mergeTracks = false;
+
+ public void setInputCollectionName(String inputCollectionName) {
+ this.inputCollectionName = inputCollectionName;
+ }
+
+ public void setOutputCollectionName(String outputCollectionName) {
+ this.outputCollectionName = outputCollectionName;
+ }
+
+ /**
+ * Merge tracks with overlapping hit content. Right now nothing actually
+ * happens to the merged tracks; this is just for testing.
+ *
+ * @param mergeTracks default to false
+ */
+ public void setMergeTracks(boolean mergeTracks) {
+ this.mergeTracks = mergeTracks;
+ }
+
+ @Override
+ protected void detectorChanged(Detector detector) {
+ bfield = Math.abs(TrackUtils.getBField(detector).magnitude());
+ _scattering.getMaterialManager().buildModel(detector);
+ _scattering.setBField(bfield); // only absolute of B is needed as it's used for momentum calculation only
+ }
+
+ @Override
+ protected void process(EventHeader event) {
+ if (!event.hasCollection(Track.class, inputCollectionName)) {
+ return;
+ }
+
+ List<Track> tracks = event.get(Track.class, inputCollectionName);
+
+ RelationalTable hitToStrips = TrackUtils.getHitToStripsTable(event);
+ RelationalTable hitToRotated = TrackUtils.getHitToRotatedTable(event);
+
+ List<Track> refittedTracks = new ArrayList<Track>();
+ List<LCRelation> trackRelations = new ArrayList<LCRelation>();
+
+ List<GBLKinkData> kinkDataCollection = new ArrayList<GBLKinkData>();
+ List<LCRelation> kinkDataRelations = new ArrayList<LCRelation>();
+
+ Map<Track, Track> inputToRefitted = new HashMap<Track, Track>();
+
+
+
+
+ for (Track track : tracks) {
+
+ // If the event has ECal has ECal clusters, use them 4 track fitting!
+ if (event.hasCollection(Cluster.class, ECalCollectionName)) {
+
+ List<Cluster> clusters = event.get(Cluster.class, ECalCollectionName); // Get Clusters from the event to use in GBL
+
+ // alternate constructor to incorporate clusters into track refit
+ Pair<Track, GBLKinkData> newTrack = MakeGblTracks.refitTrack(TrackUtils.getHTF(track), TrackUtils.getStripHits(track, hitToStrips, hitToRotated), track.getTrackerHits(), clusters, 5, track.getType(), _scattering, bfield);
+
+ refittedTracks.add(newTrack.getFirst());
+ trackRelations.add(new BaseLCRelation(track, newTrack.getFirst()));
+ inputToRefitted.put(track, newTrack.getFirst());
+
+ kinkDataCollection.add(newTrack.getSecond());
+ kinkDataRelations.add(new BaseLCRelation(newTrack.getSecond(), newTrack.getFirst()));
+
+
+ }else{ /*..and if no clusters are in the event, just do the standard GBL track fitting~*/
+
+ List<Cluster> clusters;
+ Pair<Track, GBLKinkData> newTrack = MakeGblTracks.refitTrack(TrackUtils.getHTF(track), TrackUtils.getStripHits(track, hitToStrips, hitToRotated), track.getTrackerHits(), 5, track.getType(), _scattering, bfield);
+
+ refittedTracks.add(newTrack.getFirst());
+ trackRelations.add(new BaseLCRelation(track, newTrack.getFirst()));
+ inputToRefitted.put(track, newTrack.getFirst());
+
+ kinkDataCollection.add(newTrack.getSecond());
+ kinkDataRelations.add(new BaseLCRelation(newTrack.getSecond(), newTrack.getFirst()));
+
+ }
+// newTrack.getFirst().
+ }
+
+ if (mergeTracks) {
+ List<Track> mergedTracks = new ArrayList<Track>();
+
+ for (Track track : refittedTracks) {
+ List<TrackerHit> trackHth = track.getTrackerHits();
+ otherTrackLoop:
+ for (Track otherTrack : refittedTracks) {
+ if (track == otherTrack) {
+ continue;
+ }
+
+ Set<TrackerHit> allHth = new HashSet<TrackerHit>(otherTrack.getTrackerHits());
+ allHth.addAll(trackHth);
+// if (allHth.size() == trackHth.size()) {
+// continue;
+// }
+
+ boolean[] hasHit = new boolean[6];
+
+ for (TrackerHit hit : allHth) {
+ int layer = (TrackUtils.getLayer(hit) - 1) / 2;
+ if (hasHit[layer]) {
+ continue otherTrackLoop;
+ }
+ hasHit[layer] = true;
+ }
+ for (Track mergedTrack : mergedTracks) {
+ if (mergedTrack.getTrackerHits().containsAll(allHth)) {
+ continue otherTrackLoop;
+ }
+ }
+
+ Pair<Track, GBLKinkData> mergedTrack = MakeGblTracks.refitTrack(TrackUtils.getHTF(track), TrackUtils.getStripHits(track, hitToStrips, hitToRotated), allHth, 5, track.getType(), _scattering, bfield);
+ mergedTracks.add(mergedTrack.getFirst());
+// System.out.format("%f %f %f\n", fit.get_chi2(), inputToRefitted.get(track).getChi2(), inputToRefitted.get(otherTrack).getChi2());
+// mergedTrackToTrackList.put(mergedTrack, new ArrayList<Track>());
+ }
+ }
+
+ for (Track mergedTrack : mergedTracks) {
+ List<Track> subTracks = new ArrayList<Track>();
+ Set<TrackerHit> trackHth = new HashSet<TrackerHit>(mergedTrack.getTrackerHits());
+ for (Track track : refittedTracks) {
+ if (trackHth.containsAll(track.getTrackerHits())) {
+ subTracks.add(track);
+ }
+ }
+ System.out.format("%f:\t", mergedTrack.getChi2());
+ for (Track subTrack : subTracks) {
+ System.out.format("%f (%d)\t", subTrack.getChi2(), subTrack.getTrackerHits().size());
+ }
+ System.out.println();
+ }
+ }
+ // Put the tracks back into the event and exit
+ int flag = 1 << LCIOConstants.TRBIT_HITS;
+ event.put(outputCollectionName, refittedTracks, Track.class, flag);
+ event.put(trackRelationCollectionName, trackRelations, LCRelation.class, 0);
+ event.put(GBLKinkData.DATA_COLLECTION, kinkDataCollection, GBLKinkData.class, 0);
+ event.put(GBLKinkData.DATA_RELATION_COLLECTION, kinkDataRelations, LCRelation.class, 0);
+ }
+}
Modified: java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/MakeGblTracks.java
=============================================================================
--- java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/MakeGblTracks.java (original)
+++ java/branches/brad-dev/tracking/src/main/java/org/hps/recon/tracking/gbl/MakeGblTracks.java Fri Aug 12 14:27:22 2016
@@ -38,6 +38,10 @@
import org.lcsim.fit.helicaltrack.HelixUtils;
import org.lcsim.recon.tracking.digitization.sisim.SiTrackerHitStrip1D;
import org.lcsim.recon.tracking.digitization.sisim.TrackerHitType;
+
+// 4 ECal
+import org.lcsim.event.Cluster;
+
/**
* Utilities that create track objects from fitted GBL trajectories.
@@ -292,6 +296,7 @@
* (not covariance)
* @param stripHits Strip hits to be used for the GBL fit. Does not need to
* be in sorted order.
+ * @param clusters ECal clusters to be used in GBL fit
* @param hth Stereo hits for the track's hit list (these are not used in
* the GBL fit). Does not need to be in sorted order.
* @param nIterations Number of times to iterate the GBL fit.
@@ -299,6 +304,7 @@
* @param bfield B-field
* @return The refitted track.
*/
+// Old refitTrack (still used if no clusters are present)
public static Pair<Track, GBLKinkData> refitTrack(HelicalTrackFit helix, Collection<TrackerHit> stripHits, Collection<TrackerHit> hth, int nIterations, int trackType, MultipleScattering scattering, double bfield) {
List<TrackerHit> allHthList = TrackUtils.sortHits(hth);
List<TrackerHit> sortedStripHits = TrackUtils.sortHits(stripHits);
@@ -312,6 +318,21 @@
return mergedTrack;
}
+// Overloaded to fit using ECal clusters via 'doGBLFit'
+ public static Pair<Track, GBLKinkData> refitTrack(HelicalTrackFit helix, Collection<TrackerHit> stripHits, Collection<TrackerHit> hth, List<Cluster> clusters, int nIterations, int trackType, MultipleScattering scattering, double bfield) {
+ List<TrackerHit> allHthList = TrackUtils.sortHits(hth);
+ List<TrackerHit> sortedStripHits = TrackUtils.sortHits(stripHits);
+ FittedGblTrajectory fit = doGBLFit(helix, sortedStripHits, scattering, bfield, 0);
+ for (int i = 0; i < nIterations; i++) {
+ Pair<Track, GBLKinkData> newTrack = makeCorrectedTrack(fit, helix, allHthList, trackType, bfield);
+ helix = TrackUtils.getHTF(newTrack.getFirst());
+ fit = doGBLFit(helix, sortedStripHits, clusters, scattering, bfield, 0);
+ }
+ Pair<Track, GBLKinkData> mergedTrack = makeCorrectedTrack(fit, helix, allHthList, trackType, bfield);
+ return mergedTrack;
+ }
+
+// Old GBL fit
public static FittedGblTrajectory doGBLFit(HelicalTrackFit htf, List<TrackerHit> stripHits, MultipleScattering _scattering, double bfield, int debug) {
List<GBLStripClusterData> stripData = makeStripData(htf, stripHits, _scattering, bfield, debug);
double bfac = Constants.fieldConversion * bfield;
@@ -320,6 +341,16 @@
return fit;
}
+// Overloaded GBL fit using ECal clusters via 'makeStripData'
+ public static FittedGblTrajectory doGBLFit(HelicalTrackFit htf, List<TrackerHit> stripHits, List<Cluster> clusters, MultipleScattering _scattering, double bfield, int debug) {
+ List<GBLStripClusterData> stripData = makeStripData(htf, stripHits, _scattering, bfield, debug);
+ double bfac = Constants.fieldConversion * bfield;
+
+ FittedGblTrajectory fit = HpsGblRefitter.fit(stripData, bfac, debug > 0);
+ return fit;
+ }
+
+// Old makeStripData
public static List<GBLStripClusterData> makeStripData(HelicalTrackFit htf, List<TrackerHit> stripHits, MultipleScattering _scattering, double _B, int _debug) {
List<GBLStripClusterData> stripClusterDataList = new ArrayList<GBLStripClusterData>();
@@ -435,6 +466,185 @@
}
return stripClusterDataList;
}
+
+// Overloaded makeStripData to make "GBLStripClusterData" out of ECal clusters
+ public static List<GBLStripClusterData> makeStripData(HelicalTrackFit htf, List<TrackerHit> stripHits, List<Cluster> clusters, MultipleScattering _scattering, double _B, int _debug) {
+ List<GBLStripClusterData> stripClusterDataList = new ArrayList<GBLStripClusterData>();
+
+ // Find scatter points along the path
+ MultipleScattering.ScatterPoints scatters = _scattering.FindHPSScatterPoints(htf);
+
+ if (_debug > 0) {
+ System.out.printf("perPar covariance matrix\n%s\n", htf.covariance().toString());
+ }
+
+ for (TrackerHit stripHit : stripHits) {
+ HelicalTrackStripGbl strip;
+ if (stripHit instanceof SiTrackerHitStrip1D) {
+ strip = new HelicalTrackStripGbl(makeDigiStrip((SiTrackerHitStrip1D) stripHit), true);
+ } else {
+ SiTrackerHitStrip1D newHit = new SiTrackerHitStrip1D(stripHit);
+ strip = new HelicalTrackStripGbl(makeDigiStrip(newHit), true);
+ }
+
+ // find Millepede layer definition from DetectorElement
+ HpsSiSensor sensor = (HpsSiSensor) ((RawTrackerHit) stripHit.getRawHits().get(0)).getDetectorElement();
+
+ int millepedeId = sensor.getMillepedeId();
+
+ if (_debug > 0) {
+ System.out.printf("layer %d millepede %d (DE=\"%s\", origin %s) \n", strip.layer(), millepedeId, sensor.getName(), strip.origin().toString());
+ }
+
+ //Center of the sensor
+ Hep3Vector origin = strip.origin();
+
+ //Find intercept point with sensor in tracking frame
+ Hep3Vector trkpos = TrackUtils.getHelixPlaneIntercept(htf, strip, Math.abs(_B));
+ if (trkpos == null) {
+ if (_debug > 0) {
+ System.out.println("Can't find track intercept; use sensor origin");
+ }
+ trkpos = strip.origin();
+ }
+ if (_debug > 0) {
+ System.out.printf("trkpos at intercept [%.10f %.10f %.10f]\n", trkpos.x(), trkpos.y(), trkpos.z());
+ }
+
+ //GBLDATA
+ GBLStripClusterData stripData = new GBLStripClusterData(millepedeId);
+ //Add to output list
+ stripClusterDataList.add(stripData);
+
+ //path length to intercept
+ double s = HelixUtils.PathToXPlane(htf, trkpos.x(), 0, 0).get(0);
+ double s3D = s / Math.cos(Math.atan(htf.slope()));
+
+ //GBLDATA
+ stripData.setPath(s);
+ stripData.setPath3D(s3D);
+
+ //GBLDATA
+ stripData.setU(strip.u());
+ stripData.setV(strip.v());
+ stripData.setW(strip.w());
+
+ //Print track direction at intercept
+ Hep3Vector tDir = HelixUtils.Direction(htf, s);
+ double phi = htf.phi0() - s / htf.R();
+ double lambda = Math.atan(htf.slope());
+
+ //GBLDATA
+ stripData.setTrackDir(tDir);
+ stripData.setTrackPhi(phi);
+ stripData.setTrackLambda(lambda);
+
+ //Print residual in measurement system
+ // start by find the distance vector between the center and the track position
+ Hep3Vector vdiffTrk = VecOp.sub(trkpos, origin);
+
+ // then find the rotation from tracking to measurement frame
+ Hep3Matrix trkToStripRot = getTrackToStripRotation(sensor);
+
+ // then rotate that vector into the measurement frame to get the predicted measurement position
+ Hep3Vector trkpos_meas = VecOp.mult(trkToStripRot, vdiffTrk);
+
+ //GBLDATA
+ stripData.setMeas(strip.umeas());
+ stripData.setTrackPos(trkpos_meas);
+ stripData.setMeasErr(strip.du());
+
+ if (_debug > 1) {
+ System.out.printf("rotation matrix to meas frame\n%s\n", VecOp.toString(trkToStripRot));
+ System.out.printf("tPosGlobal %s origin %s\n", trkpos.toString(), origin.toString());
+ System.out.printf("tDiff %s\n", vdiffTrk.toString());
+ System.out.printf("tPosMeas %s\n", trkpos_meas.toString());
+ }
+
+ if (_debug > 0) {
+ System.out.printf("layer %d millePedeId %d uRes %.10f\n", strip.layer(), millepedeId, stripData.getMeas() - stripData.getTrackPos().x());
+ }
+
+ // find scattering angle
+ MultipleScattering.ScatterPoint scatter = scatters.getScatterPoint(((RawTrackerHit) strip.getStrip().rawhits().get(0)).getDetectorElement());
+ double scatAngle;
+
+ if (scatter != null) {
+ scatAngle = scatter.getScatterAngle().Angle();
+ } else {
+ if (_debug > 0) {
+ System.out.printf("WARNING cannot find scatter for detector %s with strip cluster at %s\n", ((RawTrackerHit) strip.getStrip().rawhits().get(0)).getDetectorElement().getName(), strip.origin().toString());
+ }
+ scatAngle = GblUtils.estimateScatter(sensor, htf, _scattering, _B);
+ }
+
+ //GBLDATA
+ stripData.setScatterAngle(scatAngle);
+ }
+
+ // *** Now loop over clusters and mock up GBLStripClusterData for them!~<3 ***
+ for (Cluster cluster : clusters) {
+
+ // This sets a BS Millipede ID (0) for the ECal...
+ GBLStripClusterData ECal_stripData = new GBLStripClusterData(666);
+ //Add ECal strips to the output list
+ stripClusterDataList.add(ECal_stripData);
+
+
+ // TRANSLATE ALL THIS GBLStripClusterData STUFF INTO ECAL, USING CLUSTERS!!!!
+
+ //path length to intercept
+ // double s = HelixUtils.PathToXPlane(htf, trkpos.x(), 0, 0).get(0);
+ // double s3D = s / Math.cos(Math.atan(htf.slope()));
+
+ //GBLDATA
+ // ECal_stripData.setPath(s);
+ // ECal_stripData.setPath3D(s3D);
+
+ double ecal_x = cluster.getPosition()[0];
+ double ecal_y = cluster.getPosition()[1];
+ double ecal_z = cluster.getPosition()[2];
+
+ // //GBLDATA
+ // ECal_stripData.setU(strip.u());
+ // ECal_stripData.setV(strip.v());
+ // ECal_stripData.setW(strip.w());
+
+ //Print track direction at intercept
+ // Hep3Vector tDir = HelixUtils.Direction(htf, s);
+ // double phi = htf.phi0() - s / htf.R();
+ // double lambda = Math.atan(htf.slope());
+
+ // //GBLDATA
+ // ECal_stripData.setTrackDir(tDir);
+ // ECal_stripData.setTrackPhi(phi);
+ // ECal_stripData.setTrackLambda(lambda);
+
+ //Print residual in measurement system
+ //Start by find the distance vector between the center and the track position
+ // Hep3Vector vdiffTrk = VecOp.sub(trkpos, origin);
+
+ // then find the rotation from tracking to measurement frame
+ // Hep3Matrix trkToStripRot = getTrackToStripRotation(sensor);
+
+ // then rotate that vector into the measurement frame to get the predicted measurement position
+ // Hep3Vector trkpos_meas = VecOp.mult(trkToStripRot, vdiffTrk);
+
+ //GBLDATA
+ // ECal_stripData.setMeas(strip.umeas());
+ // ECal_stripData.setTrackPos(trkpos_meas);
+ // ECal_stripData.setMeasErr(strip.du());
+
+ // No Scat for ECal
+ ECal_stripData.setScatterAngle(0);
+
+ }
+
+
+
+ return stripClusterDataList;
+ }
+
private static Hep3Matrix getTrackToStripRotation(SiSensor sensor) {
// This function transforms the hit to the sensor coordinates
|
