lcsim/src/org/lcsim/contrib/HansWenzel/Tracking
diff -N DigiTrackerHitsAccessDriver.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ DigiTrackerHitsAccessDriver.java 18 Oct 2007 22:09:06 -0000 1.1
@@ -0,0 +1,190 @@
+import java.util.*;
+
+import hep.aida.*;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.SimTrackerHit;
+import org.lcsim.recon.cat.util.NoSuchParameterException;
+import org.lcsim.units.clhep.SystemOfUnits;
+import org.lcsim.util.Driver;
+import org.lcsim.util.aida.AIDA;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+
+
+
+import org.lcsim.contrib.onoprien.tracking.geom.Sensor;
+import org.lcsim.contrib.onoprien.tracking.geom.SensorType;
+import org.lcsim.contrib.onoprien.tracking.hit.DigiTrackerHit;
+import org.lcsim.contrib.onoprien.tracking.hit.TrackerCluster;
+import org.lcsim.contrib.onoprien.tracking.hit.TrackerHit;
+import org.lcsim.contrib.onoprien.tracking.mctruth.MCTruth;
+import org.lcsim.contrib.onoprien.tracking.mctruth.MCTruthDriver;
+import org.lcsim.contrib.onoprien.tracking.mctruth.SimGroup;
+
+/**
+ * Driver for accessing the DigiTrackerHits in the Event
+ *
+ * @author Hans Wenzel
+ * @version $Id: DigiTrackerHitsAccessDriver.java,v 1.1 2007/10/18 22:09:06 wenzel Exp $
+ */
+public class DigiTrackerHitsAccessDriver extends Driver {
+
+// -- Constructor : -----------------------------------------------------------
+
+ public DigiTrackerHitsAccessDriver() {
+
+ _aida = AIDA.defaultInstance();
+ _hFactory = _aida.histogramFactory();
+
+ }
+
+// -- Processing event : ------------------------------------------------------
+
+ public void process(EventHeader event) {
+
+ System.out.println(" ");
+ System.out.println("Event "+event.getEventNumber());
+ super.process(event);
+ _iGroup = 0;
+
+ MCTruth mcTruth = (MCTruth) event.get("MCTruth");
+
+ Collection<SimGroup> simGroups = mcTruth.getSimGroups();
+ List<SimGroup> missedGroups = mcTruth.getMissedSimGroups();
+
+ System.out.println("Testing digitization, found "+simGroups.size()+" converted and "
+ +missedGroups.size()+" missed SimTrackerHit groups");
+
+ /* for (SimGroup group : simGroups) {
+ List<SimTrackerHit> simList = group.getSimTrackerHits();
+ List<DigiTrackerHit> digiList = group.getDigiTrackerHits();
+ plotGroup(simList, digiList);
+
+ }
+ */
+ }
+
+// -----------------------------------------------------------------------------
+
+ public void plotGroup(List<SimTrackerHit> simList, List<DigiTrackerHit> digiList) {
+
+ Sensor sensor = digiList.get(0).getSensor();
+ SensorType sensorType = sensor.getType();
+ boolean isStrip = (sensorType.getHitDimension() == 1);
+ String title = "Sensor" + sensor.getID() + "_" + _iGroup++ + "_";
+ String subDetName = simList.get(0).getSubdetector().getName();
+
+ // Calculate bin values and create histograms
+
+ double[] lowU = new double[digiList.size()];
+ double[] lowV = new double[digiList.size()];
+ double highU = - Double.MAX_VALUE;
+ double highV = - Double.MAX_VALUE;
+ int i=0;
+ for (DigiTrackerHit dHit : digiList) {
+ int channel = dHit.getChannel();
+ Hep3Vector channelSize = sensorType.getChannelDimensions(channel);
+ Hep3Vector locPos = sensorType.getChannelPosition(channel);
+ lowU[i] = locPos.x() - channelSize.x()/2.;
+ double uHigh = locPos.x() + channelSize.x()/2.;
+ if (uHigh > highU) highU = uHigh;
+ lowV[i] = locPos.y() - channelSize.y()/2.;
+ double vHigh = locPos.y() + channelSize.y()/2.;
+ if (vHigh > highV) highV = vHigh;
+ i++;
+ }
+ Arrays.sort(lowU);
+ double lastU = - Double.MAX_VALUE;
+ ArrayList<Double> binsU = new ArrayList<Double>();
+ for (i=0; i<lowU.length; i++) {
+ if (lowU[i] > lastU) {
+ lastU = lowU[i];
+ binsU.add(lastU);
+ }
+ }
+ binsU.add(highU);
+ int nU = binsU.size();
+ lowU = new double[nU];
+ for (i=0; i<nU; i++) {
+ lowU[i] = binsU.get(i);
+ }
+ IHistogram digiHist;
+ IHistogram simHist;
+ if (isStrip) {
+ digiHist = _hFactory.createHistogram1D(title+"digi", title+"digi", lowU);
+ simHist = _hFactory.createHistogram1D(title+"sim", title+"sim", lowU);
+ } else {
+ Arrays.sort(lowV);
+ double lastV = - Double.MAX_VALUE;
+ ArrayList<Double> binsV = new ArrayList<Double>();
+ for (i=0; i<lowV.length; i++) {
+ if (lowV[i] > lastV) {
+ lastV = lowV[i];
+ binsV.add(lastV);
+ }
+ }
+ binsV.add(highV);
+ int nV = binsV.size();
+ lowV = new double[nV];
+ for (i=0; i<nV; i++) {
+ lowV[i] = binsV.get(i);
+ }
+ digiHist = _hFactory.createHistogram2D(title+"digi", title+"digi", lowU, lowV);
+ simHist = _hFactory.createHistogram2D(title+"sim", title+"sim", lowU, lowV);
+ }
+ _eventHistList.add(digiHist);
+ _eventHistList.add(simHist);
+
+ // Fill histograms
+
+ Hep3Vector meanSimPos = new BasicHep3Vector();
+ double totCharge = 0.;
+ for (SimTrackerHit simHit : simList) {
+ double[] gPos = simHit.getPoint();
+ Hep3Vector globPos = new BasicHep3Vector(gPos);
+ Hep3Vector locPos = sensor.globalToLocal(globPos);
+ double charge = (simHit.getPathLength() > 0.) ? simHit.getPathLength()*simHit.getdEdx() : simHit.getdEdx();
+ meanSimPos = VecOp.add(meanSimPos, VecOp.mult(charge, sensor.localToGlobal(locPos)));
+ totCharge += charge;
+ if (isStrip) {
+ ((IHistogram1D)simHist).fill(locPos.x(), charge);
+ } else {
+ ((IHistogram2D)simHist).fill(locPos.x(), locPos.y(), charge);
+ }
+ }
+ meanSimPos = VecOp.mult(1./totCharge, meanSimPos);
+
+ Hep3Vector meanDigiPos = new BasicHep3Vector();
+ totCharge = 0.;
+ for (DigiTrackerHit digiHit : digiList) {
+ Hep3Vector locPos = sensorType.getChannelPosition(digiHit.getChannel());
+ Hep3Vector globPos = sensor.localToGlobal(locPos);
+ double charge = digiHit.getSignal();
+ meanDigiPos = VecOp.add(meanDigiPos, VecOp.mult(charge, globPos));
+ totCharge += charge;
+ if (isStrip) {
+ ((IHistogram1D)digiHist).fill(locPos.x(), charge);
+ } else {
+ ((IHistogram2D)digiHist).fill(locPos.x(), locPos.y(), charge);
+ }
+ }
+ meanDigiPos = VecOp.mult(1./totCharge, meanDigiPos);
+
+ double diffSimDigi = VecOp.sub(meanSimPos, meanDigiPos).magnitude();
+ _aida.cloud1D("Digi to Sim distance - " + subDetName).fill(diffSimDigi);
+ double stripLength = sensorType.getChannelDimensions(digiList.get(0).getChannel()).y();
+ if ((diffSimDigi > stripLength/2.) && (sensorType.getHitDimension() == 1)) {
+ System.out.println("Sensor "+digiList.get(0).getSensor().getID()+" "+_iGroup+" Dif: "+diffSimDigi+" Strip half length: "+(stripLength/2.));
+ }
+
+ }
+// -- Private parts : ---------------------------------------------------------
+
+ protected AIDA _aida;
+ protected IHistogramFactory _hFactory;
+ protected LinkedList<IBaseHistogram> _eventHistList;
+
+ protected int _iGroup;
+
+}