diff --git a/tracking/src/main/java/org/hps/recon/tracking/DataTrackerHitDriver.java b/tracking/src/main/java/org/hps/recon/tracking/DataTrackerHitDriver.java index 222e6f1e72..2a418637f9 100644 --- a/tracking/src/main/java/org/hps/recon/tracking/DataTrackerHitDriver.java +++ b/tracking/src/main/java/org/hps/recon/tracking/DataTrackerHitDriver.java @@ -49,6 +49,8 @@ public class DataTrackerHitDriver extends Driver { private double threeClusterErr = clusterErrorMultiplier / 3.0; private double fourClusterErr = clusterErrorMultiplier / 2.0; private double fiveClusterErr = clusterErrorMultiplier / 1.0; + // weight the hits in a cluster by charge? (if not, all hits have equal weight) + private boolean useWeights = true; // Various data lists required by digitization. private List processPaths = new ArrayList(); private List processDEs = new ArrayList(); @@ -127,6 +129,10 @@ public void setFiveClusterErr(double fiveClusterErr) { this.fiveClusterErr = fiveClusterErr; } + public void setUseWeights(boolean useWeights){ + this.useWeights = useWeights; + } + /** * Creates a new instance of TrackerHitDriver. */ @@ -180,11 +186,18 @@ public void detectorChanged(Detector detector) { stripClusterer.setCentralStripAveragingThreshold(clusterCentralStripAveragingThreshold); // Set the cluster errors. - stripClusterer.SetOneClusterErr(oneClusterErr); - stripClusterer.SetTwoClusterErr(twoClusterErr); - stripClusterer.SetThreeClusterErr(threeClusterErr); - stripClusterer.SetFourClusterErr(fourClusterErr); - stripClusterer.SetFiveClusterErr(fiveClusterErr); + + DefaultSiliconResolutionModel model = new DefaultSiliconResolutionModel(); + + model.setOneClusterErr(oneClusterErr); + model.setTwoClusterErr(twoClusterErr); + model.setThreeClusterErr(threeClusterErr); + model.setFourClusterErr(fourClusterErr); + model.setFiveClusterErr(fiveClusterErr); + model.setUseWeights(useWeights); + + stripClusterer.setResolutionModel(model); + // Set the detector to process. processPaths.add(subdetectorName); diff --git a/tracking/src/main/java/org/hps/recon/tracking/DefaultSiliconResolutionModel.java b/tracking/src/main/java/org/hps/recon/tracking/DefaultSiliconResolutionModel.java new file mode 100644 index 0000000000..59df2cf8df --- /dev/null +++ b/tracking/src/main/java/org/hps/recon/tracking/DefaultSiliconResolutionModel.java @@ -0,0 +1,108 @@ +package org.hps.recon.tracking; + +import java.util.List; +import java.util.Map; + +import org.lcsim.detector.tracker.silicon.SiSensor; +import org.lcsim.detector.tracker.silicon.SiSensorElectrodes; +import org.lcsim.detector.tracker.silicon.SiStrips; + +import hep.physics.vec.BasicHep3Vector; +import hep.physics.vec.Hep3Vector; +import hep.physics.vec.VecOp; + +public class DefaultSiliconResolutionModel implements SiliconResolutionModel{ + + @Override + public double getMeasuredResolution(List cluster, SiSensorElectrodes electrodes) + + { + double measured_resolution; + + double sense_pitch = ((SiSensor) electrodes.getDetectorElement()).getSenseElectrodes(electrodes.getChargeCarrier()).getPitch(0); + + // double readout_pitch = electrodes.getPitch(0); + // double noise = + // _readout_chip.getChannel(strip_number).computeNoise(electrodes.getCapacitance(strip_number)); + // double signal_expected = (0.000280/DopedSilicon.ENERGY_EHPAIR) * + // ((SiSensor)electrodes.getDetectorElement()).getThickness(); // ~280 KeV/mm for thick Si + // sensors + if (cluster.size() == 1) { + measured_resolution = sense_pitch * _oneClusterErr; + } else if (cluster.size() == 2) { + measured_resolution = sense_pitch * _twoClusterErr; + } else if (cluster.size() == 3) { + measured_resolution = sense_pitch * _threeClusterErr; + } else if (cluster.size() == 4) { + measured_resolution = sense_pitch * _fourClusterErr; + } else { + measured_resolution = sense_pitch * _fiveClusterErr; + } + + return measured_resolution; + } + + public double getUnmeasuredResolution(List cluster, SiSensorElectrodes electrodes, Map strip_map) { + // Get length of longest strip in hit + double hit_length = 0; + for (FittedRawTrackerHit hit : cluster) { + hit_length = Math.max(hit_length, ((SiStrips) electrodes).getStripLength(strip_map.get(hit))); + } + return hit_length / Math.sqrt(12); + } + + + //perhaps the best values for these are .19, .12 and .2? + private double _oneClusterErr = 1 / Math.sqrt(12); + private double _twoClusterErr = 1 / 5.; + private double _threeClusterErr = 1 / 3.; + private double _fourClusterErr = 1 / 2.; + private double _fiveClusterErr = 1; + + private boolean _useWeights = true; + + public void setOneClusterErr(double err) { + _oneClusterErr = err; + } + + public void setTwoClusterErr(double err) { + _twoClusterErr = err; + } + + public void setThreeClusterErr(double err) { + _threeClusterErr = err; + } + + public void setFourClusterErr(double err) { + _fourClusterErr = err; + } + + public void setFiveClusterErr(double err) { + _fiveClusterErr = err; + } + + public void setUseWeights(boolean useWeights){ + _useWeights = useWeights; + } + + @Override + public Hep3Vector weightedAveragePosition(List signals, List positions) { + double total_weight = 0; + Hep3Vector position = new BasicHep3Vector(0, 0, 0); + + for (int istrip = 0; istrip < signals.size(); istrip++) { + double signal = signals.get(istrip); + + double weight = _useWeights ? signal : 1; + total_weight += weight; + position = VecOp.add(position, VecOp.mult(weight, positions.get(istrip))); + /*if (_debug) { + System.out.println(this.getClass().getSimpleName() + "strip " + istrip + ": signal " + signal + " position " + positions.get(istrip) + " -> total_position " + position.toString() + " ( total charge " + total_charge + ")"); + }*/ + + } + + return VecOp.mult(1 / total_weight, position); + } + +} diff --git a/tracking/src/main/java/org/hps/recon/tracking/SiliconResolutionModel.java b/tracking/src/main/java/org/hps/recon/tracking/SiliconResolutionModel.java new file mode 100644 index 0000000000..1d1b252e5c --- /dev/null +++ b/tracking/src/main/java/org/hps/recon/tracking/SiliconResolutionModel.java @@ -0,0 +1,14 @@ +package org.hps.recon.tracking; + +import java.util.List; +import java.util.Map; + +import org.lcsim.detector.tracker.silicon.SiSensorElectrodes; + +import hep.physics.vec.Hep3Vector; + +public interface SiliconResolutionModel { + public double getMeasuredResolution(List cluster, SiSensorElectrodes electrodes); + public double getUnmeasuredResolution(List cluster, SiSensorElectrodes electrodes, Map strip_map); + public Hep3Vector weightedAveragePosition(List signals, List positions); +} diff --git a/tracking/src/main/java/org/hps/recon/tracking/StripMaker.java b/tracking/src/main/java/org/hps/recon/tracking/StripMaker.java index 3860eb1f4d..078d4d7f7c 100644 --- a/tracking/src/main/java/org/hps/recon/tracking/StripMaker.java +++ b/tracking/src/main/java/org/hps/recon/tracking/StripMaker.java @@ -1,9 +1,7 @@ package org.hps.recon.tracking; 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.HashMap; @@ -46,14 +44,15 @@ public class StripMaker { SiTrackerIdentifierHelper _sid_helper; // Temporary map connecting hits to strip numbers for sake of speed (reset once per sensor) Map _strip_map = new HashMap(); - double _oneClusterErr = 1 / Math.sqrt(12); - double _twoClusterErr = 1 / 5; - double _threeClusterErr = 1 / 3; - double _fourClusterErr = 1 / 2; - double _fiveClusterErr = 1; + boolean _debug = false; + private SiliconResolutionModel _res_model = new DefaultSiliconResolutionModel(); + public void setResolutionModel(SiliconResolutionModel model){ + _res_model = model; + } + public StripMaker(ClusteringAlgorithm algo) { _clustering = algo; } @@ -154,25 +153,7 @@ public List makeHits(SiSensor sensor, SiSensorElectrodes electrode return hits; } - public void SetOneClusterErr(double err) { - _oneClusterErr = err; - } - - public void SetTwoClusterErr(double err) { - _twoClusterErr = err; - } - - public void SetThreeClusterErr(double err) { - _threeClusterErr = err; - } - - public void SetFourClusterErr(double err) { - _fourClusterErr = err; - } - - public void SetFiveClusterErr(double err) { - _fiveClusterErr = err; - } + public void setCentralStripAveragingThreshold(int max_noaverage_nstrips) { _max_noaverage_nstrips = max_noaverage_nstrips; @@ -242,20 +223,9 @@ private Hep3Vector getPosition(List cluster, SiSensorElectr System.out.println(this.getClass().getSimpleName() + " Calculate charge weighted mean for " + signals.size() + " signals"); } - double total_charge = 0; - Hep3Vector position = new BasicHep3Vector(0, 0, 0); - - for (int istrip = 0; istrip < signals.size(); istrip++) { - double signal = signals.get(istrip); - - total_charge += signal; - position = VecOp.add(position, VecOp.mult(signal, positions.get(istrip))); - if (_debug) { - System.out.println(this.getClass().getSimpleName() + "strip " + istrip + ": signal " + signal + " position " + positions.get(istrip) + " -> total_position " + position.toString() + " ( total charge " + total_charge + ")"); - } - - } - position = VecOp.mult(1 / total_charge, position); + + Hep3Vector position = _res_model.weightedAveragePosition(signals, positions); + if (_debug) { System.out.println(this.getClass().getSimpleName() + " charge weighted position " + position.toString() + " (before trans)"); } @@ -303,8 +273,8 @@ private double getTime(List cluster) { private SymmetricMatrix getCovariance(List cluster, SiSensorElectrodes electrodes) { SymmetricMatrix covariance = new SymmetricMatrix(3); - covariance.setElement(0, 0, Math.pow(getMeasuredResolution(cluster, electrodes), 2)); - covariance.setElement(1, 1, Math.pow(getUnmeasuredResolution(cluster, electrodes), 2)); + covariance.setElement(0, 0, Math.pow(_res_model.getMeasuredResolution(cluster, electrodes), 2)); + covariance.setElement(1, 1, Math.pow(_res_model.getUnmeasuredResolution(cluster, electrodes, _strip_map), 2)); covariance.setElement(2, 2, 0.0); SymmetricMatrix covariance_global = electrodes.getLocalToGlobal().transformed(covariance); @@ -329,57 +299,7 @@ private SymmetricMatrix getCovariance(List cluster, SiSenso // return new SymmetricMatrix((Matrix)covariance_global); } - private double getMeasuredResolution(List cluster, SiSensorElectrodes electrodes) // should - // replace - // this - // by - // a - // ResolutionModel - // class - // that - // gives - // expected - // resolution. - // This - // could - // be - // a - // big - // job. - { - double measured_resolution; - - double sense_pitch = ((SiSensor) electrodes.getDetectorElement()).getSenseElectrodes(electrodes.getChargeCarrier()).getPitch(0); - - // double readout_pitch = electrodes.getPitch(0); - // double noise = - // _readout_chip.getChannel(strip_number).computeNoise(electrodes.getCapacitance(strip_number)); - // double signal_expected = (0.000280/DopedSilicon.ENERGY_EHPAIR) * - // ((SiSensor)electrodes.getDetectorElement()).getThickness(); // ~280 KeV/mm for thick Si - // sensors - if (cluster.size() == 1) { - measured_resolution = sense_pitch * _oneClusterErr; - } else if (cluster.size() == 2) { - measured_resolution = sense_pitch * _twoClusterErr; - } else if (cluster.size() == 3) { - measured_resolution = sense_pitch * _threeClusterErr; - } else if (cluster.size() == 4) { - measured_resolution = sense_pitch * _fourClusterErr; - } else { - measured_resolution = sense_pitch * _fiveClusterErr; - } - - return measured_resolution; - } - - private double getUnmeasuredResolution(List cluster, SiSensorElectrodes electrodes) { - // Get length of longest strip in hit - double hit_length = 0; - for (FittedRawTrackerHit hit : cluster) { - hit_length = Math.max(hit_length, ((SiStrips) electrodes).getStripLength(_strip_map.get(hit))); - } - return hit_length / Math.sqrt(12); - } + private double getEnergy(List cluster) { double total_charge = 0.0;