 |
 |
| Commit in hps-java/src/main/java/org/lcsim/hps/users/omoreno on MAIN | |||
| SvtHitEfficiency.java | +169 | -26 | 1.3 -> 1.4 |
Exclude tracks passing through bad channels from hit efficiency calculation
diff -u -r1.3 -r1.4 --- SvtHitEfficiency.java 8 Jan 2013 08:08:41 -0000 1.3 +++ SvtHitEfficiency.java 17 Jan 2013 16:41:52 -0000 1.4 @@ -9,9 +9,12 @@
import hep.aida.IPlotterStyle; import hep.physics.vec.BasicHep3Vector; import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
import java.util.ArrayList;
+import java.util.HashMap;
import java.util.List;
+import java.util.Map;
//--- lcsim ---// import org.lcsim.util.Driver;
@@ -21,18 +24,22 @@
import org.lcsim.detector.solids.Point3D; import org.lcsim.detector.solids.Polygon3D; import org.lcsim.detector.tracker.silicon.SiSensor;
+import org.lcsim.detector.tracker.silicon.SiStrips;
import org.lcsim.event.EventHeader; import org.lcsim.event.Track; import org.lcsim.event.TrackerHit; import org.lcsim.fit.helicaltrack.HelicalTrackHit; import org.lcsim.geometry.Detector;
+import org.lcsim.detector.tracker.silicon.ChargeCarrier;
-import org.lcsim.hps.monitoring.AIDAFrame; -import org.lcsim.hps.recon.ecal.HPSEcalCluster;
+import org.lcsim.hps.recon.tracking.HPSSVTCalibrationConstants;
//--- hps-java ---// import org.lcsim.hps.recon.tracking.SvtUtils; import org.lcsim.hps.recon.tracking.TrackUtils; import org.lcsim.hps.recon.tracking.SvtTrackExtrapolator;
+import org.lcsim.hps.recon.tracking.TrackerHitUtils; +import org.lcsim.hps.monitoring.AIDAFrame; +import org.lcsim.hps.recon.ecal.HPSEcalCluster;
public class SvtHitEfficiency extends Driver {
@@ -41,8 +48,10 @@
private List<IHistogram1D> histos1D = new ArrayList<IHistogram1D>();
private List<IHistogram2D> histos2D = new ArrayList<IHistogram2D>();
private List<IPlotter> plotters = new ArrayList<IPlotter>();
+ private Map<SiSensor, Map<Integer, Hep3Vector>> stripPositions = new HashMap<SiSensor, Map<Integer, Hep3Vector>>();
TrackUtils trkUtil = new TrackUtils();
SvtTrackExtrapolator extrapolator = new SvtTrackExtrapolator();
+ TrackerHitUtils trackerHitUtils = new TrackerHitUtils();
boolean debug = false;
boolean ecalClusterTrackMatch = false;
@@ -51,6 +60,7 @@
boolean enableMomentumPlots = true;
boolean enableChiSquaredPlots = true;
boolean enableTrackPositionPlots = true;
+ boolean maskBadChannels = false;
int plotterIndex = 0;
@@ -82,6 +92,20 @@
public SvtHitEfficiency(){
}
+ /**
+ * Enable/Disable debug
+ */
+ public void setDebug(boolean debug){
+ this.debug = debug;
+ }
+
+ /**
+ * Enable/Disable masking of bad channels
+ */
+ public void setMaskBadChannels(boolean maskBadChannels){
+ this.maskBadChannels = maskBadChannels;
+ }
+
public void detectorChanged(Detector detector){
// setup AIDA
@@ -99,6 +123,28 @@
IHistogram2D histo2D = null;
IHistogram1D histo1D = null;
ICloud2D cloud2D = null;
+
+ // Create a Map from sensor to bad channels and from bad channels to
+ // strip position
+ for(ChargeCarrier carrier : ChargeCarrier.values()){
+ for(SiSensor sensor : SvtUtils.getInstance().getSensors()){
+ if(sensor.hasElectrodesOnSide(carrier)){
+ stripPositions.put(sensor, new HashMap<Integer, Hep3Vector>());
+ SiStrips strips = (SiStrips) sensor.getReadoutElectrodes(carrier);
+ ITransform3D parentToLocal = sensor.getReadoutElectrodes(carrier).getParentToLocal();
+ ITransform3D localToGlobal = sensor.getReadoutElectrodes(carrier).getLocalToGlobal();
+ for(int physicalChannel = 0; physicalChannel < 640; physicalChannel++){
+ Hep3Vector localStripPosition = strips.getCellPosition(physicalChannel);
+ Hep3Vector stripPosition = parentToLocal.transformed(localStripPosition);
+ Hep3Vector globalStripPosition = localToGlobal.transformed(stripPosition);
+ //this.printDebug(SvtUtils.getInstance().getDescription(sensor) + " : Channel " + physicalChannel + " : Local Strip Position " + localStripPosition.toString());
+ //this.printDebug(SvtUtils.getInstance().getDescription(sensor) + " : Channel " + physicalChannel + " : Strip Position " + stripPosition.toString());
+ //this.printDebug(SvtUtils.getInstance().getDescription(sensor) + " : Channel " + physicalChannel + " : Strip Position " + globalStripPosition.toString());
+ stripPositions.get(sensor).put(physicalChannel, stripPosition);
+ }
+ }
+ }
+ }
//--- Momentum Plots ---//
//----------------------//
@@ -137,9 +183,11 @@
//--- Track Position Plots ---//
//----------------------------//
if(enableTrackPositionPlots){
+ int layerNumber = 1; + SiSensor sensor = null;
IPlotterStyle style = aida.analysisFactory().createPlotterFactory().createPlotterStyle();
for(int index = 1; index < 6; index++){
- plotters.add(PlotUtils.setupPlotter("Track Position - Layer " + index, 1, 2));
+ plotters.add(PlotUtils.setupPlotter("Track Position - Layer " + index, 2, 3));
title = "Track Position - Layer " + index + " - Tracks Within Acceptance";
cloud2D = aida.cloud2D(title);
PlotUtils.setup2DRegion(plotters.get(plotterIndex), title, 0, "x [mm]", "y [mm]", cloud2D, style);
@@ -149,14 +197,53 @@
title = "Track Position - Layer " + index + " - Difference";
cloud2D = aida.cloud2D(title);
PlotUtils.setup2DRegion(plotters.get(plotterIndex), title, 1, "x [mm]", "y [mm]", cloud2D, style);
- frames.get(1).addPlotter(plotters.get(plotterIndex));
+ sensor = SvtUtils.getInstance().getBottomSensor(layerNumber, 0); + title = SvtUtils.getInstance().getDescription(sensor) + " - Occupancy"; + histo1D = aida.histogram1D(title, 640, 0, 639); + histos1D.add(histo1D); + PlotUtils.setup1DRegion(plotters.get(plotterIndex), title, 2, "Channel #", histo1D); + sensor = SvtUtils.getInstance().getTopSensor(layerNumber, 0); + title = SvtUtils.getInstance().getDescription(sensor) + " - Occupancy"; + histo1D = aida.histogram1D(title, 640, 0, 639); + histos1D.add(histo1D); + PlotUtils.setup1DRegion(plotters.get(plotterIndex), title, 4, "Channel #", histo1D); + layerNumber++; + sensor = SvtUtils.getInstance().getBottomSensor(layerNumber, 0); + title = SvtUtils.getInstance().getDescription(sensor) + " - Occupancy"; + histo1D = aida.histogram1D(title, 640, 0, 639); + histos1D.add(histo1D); + PlotUtils.setup1DRegion(plotters.get(plotterIndex), title, 3, "Channel #", histo1D); + sensor = SvtUtils.getInstance().getTopSensor(layerNumber, 0); + title = SvtUtils.getInstance().getDescription(sensor) + " - Occupancy"; + histo1D = aida.histogram1D(title, 640, 0, 639); + histos1D.add(histo1D); + PlotUtils.setup1DRegion(plotters.get(plotterIndex), title, 5, "Channel #", histo1D); + layerNumber++; + //frames.get(1).addPlotter(plotters.get(plotterIndex));
plotterIndex++;
}
}
- for (AIDAFrame frame : frames) {
- frame.pack();
- frame.setVisible(true);
+ for (IPlotter plotter : plotters) {
+ plotter.show();
+ }
+ }
+
+ /**
+ * .
+ */
+ private Hep3Vector getStripPosition(SiSensor sensor, int physicalChannel){
+ return stripPositions.get(sensor).get(physicalChannel);
+ }
+
+ /**
+ * Print a debug message if they are enabled.
+ *
+ * @param debugMessage : message to be printed
+ */
+ private void printDebug(String debugMessage){
+ if(debug){
+ System.out.println(this.getClass().getSimpleName() + ": " + debugMessage);
}
}
@@ -168,23 +255,18 @@
// Get the list of tracks from the event
List<Track> tracks = event.get(Track.class, trackCollectionName);
- if(tracks.size() >= 2 ){
- if(debug)
- System.out.println(this.getClass().getSimpleName() + ": Two track event found! Skipping ...");
- return;
- }
+ if(tracks.size() >= 2 ) return;
// Get the list of Ecal clusters from the event
List<HPSEcalCluster> ecalClusters = event.get(HPSEcalCluster.class, ecalClustersCollectionName);
- -
+
for(Track track : tracks){
trkUtil.setTrack(track);
ecalClusterTrackMatch = false;
// Check if there is an Ecal cluster in the same detector volume as the track
- for(HPSEcalCluster ecalCluster : ecalClusters){
+ /*for(HPSEcalCluster ecalCluster : ecalClusters){
if(ecalCluster.getPosition()[1] > 0 && trkUtil.getZ0() > 0){
ecalClusterTrackMatch = true;
break;
@@ -193,12 +275,13 @@
ecalClusterTrackMatch = true;
break;
}
- }
+ }*/
+ /*
if(!ecalClusterTrackMatch){
if(debug) System.out.println(this.getClass().getSimpleName() + ": No matching Ecal cluster found");
continue;
- }
+ }*/
// Check that the track is associated with four hits only. This should
// be the case since the strategy is only requiring four hits to fit
@@ -254,22 +337,33 @@
aida.histogram1D("Track Momentum - Tracks With All Layers Hit").fill(momentum);
if(enableChiSquaredPlots)
aida.histogram1D("Chi Squared - Tracks With All Layers Hit").fill(track.getChi2());
- if(enableTrackPositionPlots){
- int layerNumber = (layer - 1)/2 + 1;
- String title = "Track Position - Layer " + layerNumber + " - Tracks With All Layers Hit";
- //aida.histogram2D(title).fill(trackPos.y(), trackPos.z());
- aida.cloud2D(title).fill(frontTrackPos.y(), frontTrackPos.z());
- }
+
return;
}
}
+ int layerNumber = (layer - 1)/2 + 1;
if(enableTrackPositionPlots){
- int layerNumber = (layer - 1)/2 + 1;
String title = "Track Position - Layer " + layerNumber + " - Difference";
//aida.histogram2D(title).fill(trackPos.y(), trackPos.z());
aida.cloud2D(title).fill(frontTrackPos.y(), frontTrackPos.z());
+
+ title = "Track Position - Layer " + layerNumber + " - Tracks With All Layers Hit";
+ //aida.histogram2D(title).fill(trackPos.y(), trackPos.z());
+ aida.cloud2D(title).fill(frontTrackPos.y(), frontTrackPos.z());
+ }
+
+ List<SiSensor> sensors = new ArrayList<SiSensor>();
+ if(trkUtil.getZ0() > 0){
+ sensors.add(SvtUtils.getInstance().getTopSensor(layer, 0));
+ sensors.add(SvtUtils.getInstance().getTopSensor(layer+1, 0));
+ } else {
+ sensors.add(SvtUtils.getInstance().getBottomSensor(layer, 0));
+ sensors.add(SvtUtils.getInstance().getBottomSensor(layer+1, 0));
}
+ aida.histogram1D(SvtUtils.getInstance().getDescription(sensors.get(0)) + " - Occupancy").fill(this.findIntersectingChannel(frontTrackPos, sensors.get(0))); + aida.histogram1D(SvtUtils.getInstance().getDescription(sensors.get(1)) + " - Occupancy").fill(this.findIntersectingChannel(rearTrackPos, sensors.get(1))); +
if(debug)
System.out.println(this.getClass().getSimpleName() + ": Stereo hit was not found.");
}
@@ -323,15 +417,65 @@
return true;
}
+
return false;
}
+
+ /**
+ *
+ */
+ public int findIntersectingChannel(Hep3Vector trackPosition, SiSensor sensor){
+
+ //--- Check that the track doesn't pass through a region of bad channels ---//
+ //--------------------------------------------------------------------------//
+
+ //Rotate the track position to the JLab coordinate system
+ this.printDebug("Track position in tracking frame: " + trackPosition.toString());
+ Hep3Vector trackPositionDet = VecOp.mult(VecOp.inverse(this.trackerHitUtils.detToTrackRotationMatrix()), trackPosition);
+ this.printDebug("Track position in JLab frame " + trackPositionDet.toString());
+ // Rotate the track to the sensor coordinates
+ ITransform3D globalToLocal = sensor.getReadoutElectrodes(ChargeCarrier.HOLE).getGlobalToLocal();
+ globalToLocal.transform(trackPositionDet);
+ this.printDebug("Track position in sensor electrodes frame " + trackPositionDet.toString());
+
+ // Find the closest channel to the track position
+ double deltaY = Double.MAX_VALUE;
+ int intersectingChannel = 0;
+ for(int physicalChannel= 0; physicalChannel < 639; physicalChannel++){
+ /*this.printDebug(SvtUtils.getInstance().getDescription(sensor) + " : Channel " + physicalChannel
+ + " : Strip Position " + stripPositions.get(sensor).get(physicalChannel));
+ this.printDebug(SvtUtils.getInstance().getDescription(sensor) + ": Channel " + physicalChannel
+ + " : Delta Y: " + Math.abs(trackPositionDet.x() - stripPositions.get(sensor).get(physicalChannel).x()));*/
+ if(Math.abs(trackPositionDet.x() - stripPositions.get(sensor).get(physicalChannel).x()) < deltaY ){
+ deltaY = Math.abs(trackPositionDet.x() - stripPositions.get(sensor).get(physicalChannel).x());
+ intersectingChannel = physicalChannel;
+ }
+ }
+
+ this.printDebug(SvtUtils.getInstance().getDescription(sensor) + ": Track intersects physical channel " + intersectingChannel);
+
+ return intersectingChannel;
+ }
/**
*
*/
public boolean sensorContainsTrack(Hep3Vector trackPosition, SiSensor sensor){
+
+
+ if(maskBadChannels){
+ int intersectingChannel = this.findIntersectingChannel(trackPosition, sensor);
+ if(intersectingChannel == 0 || intersectingChannel == 638) return false;
+
+ if(HPSSVTCalibrationConstants.isBadChannel(sensor, intersectingChannel)
+ || HPSSVTCalibrationConstants.isBadChannel(sensor, intersectingChannel+1)
+ || HPSSVTCalibrationConstants.isBadChannel(sensor, intersectingChannel-1)){
+ this.printDebug("Track intersects a bad channel!");
+ return false;
+ }
+ }
- ITransform3D localToGlobal = sensor.getGeometry().getLocalToGlobal();
+ ITransform3D localToGlobal = sensor.getGeometry().getLocalToGlobal();
Hep3Vector sensorPos = sensor.getGeometry().getPosition();
Box sensorSolid = (Box) sensor.getGeometry().getLogicalVolume().getSolid();
@@ -390,7 +534,6 @@
double area4 = this.findTriangleArea(vertices.get(3).x(), vertices.get(3).y(), vertices.get(0).x(), vertices.get(0).y(), trackPosition.y(), trackPosition.z());
if((area1 > 0 && area2 > 0 && area3 > 0 && area4 > 0) || (area1 < 0 && area2 < 0 && area3 < 0 && area4 < 0)) return true;
-
return false;
}
Use REPLY-ALL to reply to list
To unsubscribe from the LCD-CVS list, click the following link:
https://listserv.slac.stanford.edu/cgi-bin/wa?SUBED1=LCD-CVS&A=1