6 added + 4 modified, total 10 files
lcsim/src/org/lcsim/contrib/seedtracker
diff -N HitIdentifier.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ HitIdentifier.java 11 Dec 2007 20:58:28 -0000 1.1
@@ -0,0 +1,116 @@
+/*
+ * HitIdentifier.java
+ *
+ * Created on December 5, 2007, 3:04 PM
+ *
+ */
+
+package org.lcsim.contrib.seedtracker;
+
+import java.util.ArrayList;
+import java.util.List;
+
+import org.lcsim.detector.IDetectorElement;
+import org.lcsim.detector.DetectorIdentifierHelper;
+import org.lcsim.detector.identifier.IIdentifier;
+import org.lcsim.detector.identifier.IIdentifierHelper;
+import org.lcsim.detector.identifier.IExpandedIdentifier;
+import org.lcsim.detector.identifier.IIdentifierDictionary;
+import org.lcsim.detector.tracker.silicon.SiTrackerIdentifierHelper;
+import org.lcsim.event.TrackerHit;
+import org.lcsim.contrib.onoprien.tracking.hit.TrackerCluster;
+import org.lcsim.contrib.onoprien.tracking.hitmaking.OldTrackerHit;
+import org.lcsim.fit.helicaltrack.HelicalTrackHit;
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
+
+/**
+ *
+ * @author R. Partridge
+ * @version 1.0
+ */
+public class HitIdentifier {
+
+ /** Creates a new instance of HitIdentifier */
+ public HitIdentifier() {
+ }
+
+ public String Identifier(HelicalTrackHit hit) {
+ List<IDetectorElement> delist = getBaseDetElemList(hit);
+ IDetectorElement de = delist.get(0);
+ String detname = getName(de);
+ BarrelEndcapFlag beflag = getBarrelEndcapFlag(de);
+ int layer = getLayer(de, beflag);
+ if (delist.size() == 2) layer = layer / 2;
+ return detname+layer+beflag;
+ }
+
+ public String Identifier(String detname, int layer, BarrelEndcapFlag beflag) {
+ String identifier = detname + layer + beflag;
+ return identifier;
+ }
+
+ private List<IDetectorElement> getBaseDetElemList(HelicalTrackHit hit) {
+ List<IDetectorElement> delist = new ArrayList<IDetectorElement>();
+ TrackerHit thit = hit.getTrackerHit();
+ if (thit instanceof OldTrackerHit) {
+ OldTrackerHit othit = (OldTrackerHit) thit;
+ List<TrackerCluster> cluster_list = othit.getClusters();
+ for (TrackerCluster cluster : cluster_list) {
+ // Get detector element for this cluster
+ delist.add(cluster.getSensor().getDetectorElement());
+ }
+ }
+ return delist;
+ }
+
+ private String getName(IDetectorElement de) {
+ // Find the first level down from the top of the de tree
+ while (de.getParent().getParent() != null) de = de.getParent();
+ // Find the name of this detector
+ String detname = de.getName();
+ System.out.println(" System Name: "+detname);
+ return detname;
+ }
+
+ private int getLayer(IDetectorElement de, BarrelEndcapFlag beflag) {
+ int layer = 666;
+ // Are we in the barrel?
+ if (beflag == BarrelEndcapFlag.BARREL) {
+ // Find the detector element two levels below the top
+ while (de.getParent().getParent().getParent() != null) de = de.getParent();
+ } else {
+ // Find the detector element three levels below the top
+ while (de.getParent().getParent().getParent().getParent() != null) de = de.getParent();
+ }
+ // Get the DetectorIdentifierHelper
+ IIdentifierHelper hlp = de.getIdentifierHelper();
+ if (hlp instanceof DetectorIdentifierHelper) {
+ DetectorIdentifierHelper dehlp = (DetectorIdentifierHelper) hlp;
+ // Get the identifier
+ IIdentifier id = de.getIdentifier();
+ // Get the layer number
+ layer = dehlp.getLayerValue(id);
+ System.out.println("Layer: "+layer);
+ }
+ return layer;
+ }
+
+ private BarrelEndcapFlag getBarrelEndcapFlag(IDetectorElement de) {
+ BarrelEndcapFlag beflag = BarrelEndcapFlag.UNKNOWN;
+ // Find the second level down from the top of the de tree
+ while (de.getParent().getParent().getParent() != null) de = de.getParent();
+ // Get the DetectorIdentifierHelper
+ IIdentifierHelper hlp = de.getIdentifierHelper();
+ if (hlp instanceof DetectorIdentifierHelper) {
+ DetectorIdentifierHelper dehlp = (DetectorIdentifierHelper) hlp;
+ // Get the identifier
+ IIdentifier id = de.getIdentifier();
+ // Get the BarrelEndcapFlag
+ if (dehlp.isBarrel(id)) beflag = BarrelEndcapFlag.BARREL;
+ else if (dehlp.isEndcapPositive(id)) beflag = BarrelEndcapFlag.ENDCAP_NORTH;
+ else if (dehlp.isEndcapNegative(id)) beflag = BarrelEndcapFlag.ENDCAP_SOUTH;
+ System.out.println("BEFlag: "+beflag);
+ }
+ return beflag;
+ }
+}
\ No newline at end of file
lcsim/src/org/lcsim/contrib/seedtracker
diff -N MaterialCylinder.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ MaterialCylinder.java 11 Dec 2007 20:58:28 -0000 1.1
@@ -0,0 +1,59 @@
+/*
+ * MaterialCylinder.java
+ *
+ * Created on November 8, 2007, 2:10 PM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.contrib.seedtracker;
+
+import org.lcsim.detector.IPhysicalVolume;
+
+/**
+ *
+ * @author partridge
+ */
+public class MaterialCylinder {
+ private IPhysicalVolume _pv;
+ private double _radius;
+ private double _zmin;
+ private double _zmax;
+ private double _t_RL;
+
+ /** Creates a new instance of MaterialCylinder */
+ public MaterialCylinder(IPhysicalVolume pv, double radius, double zmin, double zmax, double t_RL) {
+ _pv = pv;
+ _radius = radius;
+ _zmin = zmin;
+ _zmax = zmax;
+ _t_RL = t_RL;
+ }
+
+ public double radius() {
+ return _radius;
+ }
+
+ public double zmin() {
+ return _zmin;
+ }
+
+ public double zmax() {
+ return _zmax;
+ }
+
+ public double ThicknessInRL() {
+ return _t_RL;
+ }
+
+ public String toString()
+ {
+ StringBuffer sb = new StringBuffer("Tracker MaterialCylinder for "+_pv.getName()+"\n");
+ sb.append("radius = "+_radius+"\n");
+ sb.append("Z min = "+_zmin+"\n");
+ sb.append("Z max = "+_zmax+"\n");
+ sb.append("Thickness (in RL) = "+_t_RL+"\n");
+ return sb.toString();
+ }
+}
\ No newline at end of file
lcsim/src/org/lcsim/contrib/seedtracker
diff -N MaterialDisk.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ MaterialDisk.java 11 Dec 2007 20:58:28 -0000 1.1
@@ -0,0 +1,60 @@
+/*
+ * MaterialDisk.java
+ *
+ * Created on November 8, 2007, 2:19 PM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.contrib.seedtracker;
+
+import org.lcsim.detector.IPhysicalVolume;
+
+/**
+ *
+ * @author partridge
+ */
+public class MaterialDisk {
+ private IPhysicalVolume _pv;
+ private double _rmin;
+ private double _rmax;
+ private double _z;
+ private double _t_RL;
+
+ /** Creates a new instance of MaterialCylinder */
+ public MaterialDisk(IPhysicalVolume pv, double rmin, double rmax, double z, double t_RL) {
+ _pv = pv;
+ _rmin = rmin;
+ _rmax = rmax;
+ _z = z;
+ _t_RL = t_RL;
+ }
+
+ public double rmin() {
+ return _rmin;
+ }
+
+ public double rmax() {
+ return _rmax;
+ }
+
+ public double z() {
+ return _z;
+ }
+
+ public double ThicknessInRL() {
+ return _t_RL;
+ }
+
+ public String toString()
+ {
+ StringBuffer sb = new StringBuffer("Tracker MaterialDisk "+_pv.getName()+"\n");
+ sb.append("Inner radius = "+_rmin+"\n");
+ sb.append("Outer radius = "+_rmax+"\n");
+ sb.append("Z = "+_z+"\n");
+ sb.append("Thickness (in RL) = "+_t_RL+"\n");
+ return sb.toString();
+ }
+
+}
lcsim/src/org/lcsim/contrib/seedtracker
diff -N MaterialLayer.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ MaterialLayer.java 11 Dec 2007 20:58:28 -0000 1.1
@@ -0,0 +1,72 @@
+/*
+ * MaterialLayer.java
+ *
+ * Created on August 12, 2007, 8:44 PM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.contrib.seedtracker;
+
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
+
+/**
+ *
+ * @author partridge
+ */
+public class MaterialLayer {
+ private double _rmin;
+ private double _rmax;
+ private double _zmin;
+ private double _zmax;
+ private double _X0;
+ private double _thickness;
+ private BarrelEndcapFlag _beflag;
+
+ /** Creates a new instance of MaterialLayer */
+ public MaterialLayer(double rmin, double rmax, double zmin, double zmax, double X0) {
+ _rmin = rmin;
+ _rmax = rmax;
+ _zmin = zmin;
+ _zmax = zmax;
+ _X0 = X0;
+ if (_rmax-_rmin < _zmax-_zmin) {
+ _beflag = BarrelEndcapFlag.BARREL;
+ _thickness = _rmax - _rmin;
+ }
+ else {
+ if (_zmax > 0) _beflag = BarrelEndcapFlag.ENDCAP_NORTH;
+ else _beflag = BarrelEndcapFlag.ENDCAP_SOUTH;
+ _thickness = _zmax - _zmin;
+ }
+ }
+
+ public double rmin() {
+ return _rmin;
+ }
+
+ public double rmax() {
+ return _rmax;
+ }
+
+ public double zmin() {
+ return _zmin;
+ }
+
+ public double zmax() {
+ return _zmax;
+ }
+
+ public double X0() {
+ return _X0;
+ }
+
+ public double thickness() {
+ return _thickness;
+ }
+
+ public BarrelEndcapFlag beflag() {
+ return _beflag;
+ }
+}
lcsim/src/org/lcsim/contrib/seedtracker
diff -N MaterialManager.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ MaterialManager.java 11 Dec 2007 20:58:28 -0000 1.1
@@ -0,0 +1,301 @@
+/*
+ * MaterialManager.java
+ *
+ * Created on August 13, 2007, 11:27 PM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.contrib.seedtracker;
+
+import java.util.ArrayList;
+import java.util.LinkedList;
+import java.util.List;
+import java.util.Map;
+import java.util.Set;
+
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.TrackerHit;
+import org.lcsim.geometry.Detector;
+import org.lcsim.detector.solids.ISolid;
+import org.lcsim.detector.solids.Tube;
+import org.lcsim.detector.ILogicalVolume;
+import org.lcsim.detector.material.IMaterial;
+import org.lcsim.detector.IPhysicalVolume;
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
+import org.lcsim.geometry.compact.Subdetector;
+
+/**
+ *
+ * @author partridge
+ */
+public class MaterialManager {
+ List<MaterialCylinder> _matcyl = new ArrayList<MaterialCylinder>();
+ List<MaterialDisk> _matdsk = new ArrayList<MaterialDisk>();
+ List<MaterialLayer> _matlyrs;
+
+ /** Creates a new instance of MaterialManager */
+ public MaterialManager(EventHeader event) {
+ // Build the model of tracker material
+ // Each volume defined in the compact.xml file is modelled as either
+ // a thin cylinder or disk with a thickness in radiation lengths
+ // that gives the correct total amount of material.
+ //
+ // First find the logical volume associated with the tracker
+ Detector det = event.getDetector();
+ ILogicalVolume ltrkr = det.getTrackingVolume().getLogicalVolume();
+ // Loop over the volumes defined at the compact.xml level
+ for (IPhysicalVolume pvtree : ltrkr.getDaughters()) {
+ // Flatten the geometry tree to get all daughters with material
+ List<IPhysicalVolume> pvflat = Flatten(pvtree);
+ // Calculate the total volume of material, skip this object if 0
+ double vtot = TotalVolume(pvflat);
+ if (vtot > 0.) {
+ // Calculate the average radiation length for this volume
+ double X0 = RadiationLength(pvflat);
+ // Determine if this volume should be modeled as a barrel or disk
+ if (isCylinder(pvtree, pvflat)) {
+ // Calculate the weighted radius of the elements
+ double rwgt = WeightedRadius(pvflat);
+ double zmin = zmin(pvtree, pvflat);
+ double zmax = zmax(pvtree, pvflat);
+ double zlen = zmax - zmin;
+ double thickness = vtot / (2. * Math.PI * rwgt * zlen * X0);
+ _matcyl.add(new MaterialCylinder(pvtree, rwgt, zmin, zmax, thickness));
+ } else {
+ double zwgt = WeightedZ(pvtree, pvflat);
+ double rmin = rmin(pvflat);
+ double rmax = rmax(pvflat);
+ double thickness = vtot / (Math.PI * (rmax*rmax - rmin*rmin) * X0);
+ _matdsk.add(new MaterialDisk(pvtree, rmin, rmax, zwgt, thickness));
+ }
+ }
+ }
+// for (MaterialCylinder cyl : _matcyl) {
+// System.out.println(cyl.toString());
+// }
+// for (MaterialDisk disk : _matdsk) {
+// System.out.println(disk.toString());
+// }
+ }
+
+ public double MSError(double p, double d0, double z0, double r, double z) {
+ double slope = (r-d0) / (z-z0);
+ double mserrorsq = 0.;
+ for (MaterialLayer lyr : _matlyrs) {
+ if (inside(d0,z0,lyr) != inside(r,z,lyr)) {
+ if (lyr.beflag() == BarrelEndcapFlag.BARREL) {
+ double zint = z0 + (lyr.rmax()-d0) * slope;
+ if ((zint > lyr.zmin()) && (zint < lyr.zmax())) {
+ double pathsq = Math.pow(r-lyr.rmax(),2)+Math.pow(z-zint,2);
+ double radlen = lyr.thickness() * Math.sqrt(1+Math.pow(slope,2)) / lyr.X0();
+ mserrorsq += msanglesq(p,radlen) * pathsq;
+ }
+ } else if (lyr.beflag() == BarrelEndcapFlag.ENDCAP_NORTH) {
+ double rint = d0 + (lyr.zmax()-z0) / slope;
+ if ((rint > lyr.rmin()) && (rint < lyr.rmax())) {
+ double pathsq = Math.pow(r-rint,2)+Math.pow(z-lyr.zmax(),2);
+ double radlen = lyr.thickness() * Math.sqrt(1+Math.pow(slope,-2)) / lyr.X0();
+ mserrorsq += msanglesq(p,radlen) * pathsq;
+ }
+ } else if (lyr.beflag() == BarrelEndcapFlag.ENDCAP_SOUTH) {
+ double rint = d0 + (lyr.zmin()-z0) / slope;
+ if ((rint > lyr.rmin()) && (rint < lyr.rmax())) {
+ double pathsq = Math.pow(r-rint,2)+Math.pow(z-lyr.zmin(),2);
+ double radlen = lyr.thickness() * Math.sqrt(1+Math.pow(slope,-2)) / lyr.X0();
+ mserrorsq += msanglesq(p,radlen) * pathsq;
+ }
+ }
+ }
+ }
+ return Math.sqrt(mserrorsq);
+ }
+
+ public double MSError(SeedCandidate seed, double bfield, TrackerHit hit) {
+ double[] pars = seed.getHelix().parameters();
+ double p = 0.3 * bfield * Math.abs(pars[2]) / Math.sqrt(1+Math.pow(pars[4],-2));
+ double[] pos = hit.getPosition();
+ double r = Math.sqrt(Math.pow(pos[0],2)+Math.pow(pos[1],2));
+ return MSError(p,pars[0],pars[3],r,pos[2]);
+ }
+
+ private List<IPhysicalVolume> Flatten(IPhysicalVolume xmlvol) {
+ // Flatten the geometry tree to find all the physical volumes for
+ // this tracker element that have no daughters
+ LinkedList<IPhysicalVolume> pvtree = new LinkedList<IPhysicalVolume>();
+ pvtree.add(xmlvol);
+ List<IPhysicalVolume> pvflat = new ArrayList<IPhysicalVolume>();
+ while (pvtree.size()>0) {
+ IPhysicalVolume pv = pvtree.poll();
+ if (pv.getLogicalVolume().getNumberOfDaughters()==0) {
+ pvflat.add(pv);
+ } else {
+ pvtree.addAll(pv.getLogicalVolume().getDaughters());
+ }
+ }
+ return pvflat;
+ }
+
+ private double TotalVolume(List<IPhysicalVolume> pvlist) {
+ double vtot = 0.;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // Only consider tubes until we can trust getCubicVolume
+ if (solid instanceof Tube) {
+ vtot += solid.getCubicVolume();
+ }
+ }
+ return vtot;
+ }
+
+ private double RadiationLength(List<IPhysicalVolume> pvlist) {
+ double X0wgt = 0.;
+ double invX0 = 0.;
+ double vtot = 0.;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // Only consider tubes until we can trust getCubicVolume
+ if (solid instanceof Tube) {
+ double vol = solid.getCubicVolume();
+ vtot += vol;
+ IMaterial mat = pv.getLogicalVolume().getMaterial();
+ double X0 = 10. * mat.getRadiationLength() / mat.getDensity();
+ if (X0 > 0.) invX0 += vol / X0;
+ }
+ }
+ if (invX0 > 0.) X0wgt = vtot / invX0;
+ return X0wgt;
+ }
+
+ private double rmin(List<IPhysicalVolume> pvlist) {
+ double rmin = 1.e10;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // For now, only consider elements that are tubes
+ if (solid instanceof Tube) {
+ Tube tube = (Tube) solid;
+ rmin = Math.min(rmin, tube.getInnerRadius());
+ }
+ }
+ return rmin;
+ }
+
+ private double rmax(List<IPhysicalVolume> pvlist) {
+ double rmax = 0.;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // For now, only consider elements that are tubes
+ if (solid instanceof Tube) {
+ Tube tube = (Tube) solid;
+ rmax = Math.max(rmax, tube.getOuterRadius());
+ }
+ }
+ return rmax;
+ }
+
+ private double zmin(IPhysicalVolume pvtree, List<IPhysicalVolume> pvlist) {
+ double zmin = 1.e10;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // For now, only consider elements that are tubes
+ if (solid instanceof Tube) {
+ Tube tube = (Tube) solid;
+ double z0 = pvtree.getTranslation().z();
+ zmin = Math.min(zmin, z0-tube.getZHalfLength());
+ }
+ }
+ return zmin;
+ }
+
+ private double zmax(IPhysicalVolume pvtree, List<IPhysicalVolume> pvlist) {
+ double zmax = -1.e10;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // For now, only consider elements that are tubes
+ if (solid instanceof Tube) {
+ Tube tube = (Tube) solid;
+ double z0 = pvtree.getTranslation().z();
+ zmax = Math.max(zmax, z0+tube.getZHalfLength());
+ }
+ }
+ return zmax;
+ }
+
+ private boolean isCylinder(IPhysicalVolume pvtree, List<IPhysicalVolume> pvlist) {
+ double rmin = rmin(pvlist);
+ double rmax = rmax(pvlist);
+ double zmin = zmin(pvtree, pvlist);
+ double zmax = zmax(pvtree, pvlist);
+ return (rmax - rmin) * Math.abs(zmax + zmin) < (zmax - zmin) * (rmax + rmin);
+ }
+
+ private double WeightedRadius(List<IPhysicalVolume> pvlist) {
+ double rwgt = 0.0;
+ double totwgt = 0.0;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // For now, only consider elements that are tubes
+ if (solid instanceof Tube) {
+ Tube tube = (Tube) solid;
+ double rmin = tube.getInnerRadius();
+ double rmax = tube.getOuterRadius();
+ double vol = solid.getCubicVolume();
+ IMaterial mat = pv.getLogicalVolume().getMaterial();
+ double X0 = 10.0 * mat.getRadiationLength() / mat.getDensity();
+ double wgt = vol / X0;
+ rwgt += 0.5 * (rmin + rmax) * wgt;
+ totwgt += wgt;
+ }
+ }
+ if (totwgt > 0.0) rwgt /= totwgt;
+ return rwgt;
+ }
+
+ private double WeightedZ(IPhysicalVolume pvtree, List<IPhysicalVolume> pvlist) {
+ double zwgt = 0.0;
+ double totwgt = 0.0;
+ for (IPhysicalVolume pv : pvlist) {
+ ISolid solid = pv.getLogicalVolume().getSolid();
+ // For now, only consider elements that are tubes
+ if (solid instanceof Tube) {
+ Tube tube = (Tube) solid;
+ double z0 = pvtree.getTranslation().z();
+ double vol = solid.getCubicVolume();
+ IMaterial mat = pv.getLogicalVolume().getMaterial();
+ double X0 = 10.0 * mat.getRadiationLength() / mat.getDensity();
+ double wgt = vol / X0;
+ zwgt += z0 * wgt;
+ totwgt += wgt;
+ }
+ }
+ if (totwgt > 0.0) zwgt /= totwgt;
+ return zwgt;
+ }
+
+
+ private boolean inside(double r, double z, MaterialLayer lyr) {
+ if (lyr.beflag() == BarrelEndcapFlag.BARREL) {
+ return r < lyr.rmin();
+ } else if (lyr.beflag() == BarrelEndcapFlag.ENDCAP_NORTH) {
+ return z < lyr.zmin();
+ } else {
+ return z > lyr.zmax();
+ }
+ }
+
+ private boolean outside(double r, double z, MaterialLayer lyr) {
+ if (lyr.beflag() == BarrelEndcapFlag.BARREL) {
+ return r > lyr.rmax();
+ } else if (lyr.beflag() == BarrelEndcapFlag.ENDCAP_NORTH) {
+ return z > lyr.zmax();
+ } else {
+ return z < lyr.zmin();
+ }
+ }
+
+ private double msanglesq(double p, double radlength) {
+ // Theta = (0.0136 GeV)*(z/p)*(sqrt(radLength))*(1+0.038*log(radLength))
+ return radlength * Math.pow(0.0136*(1.0 + 0.038 * Math.log(radlength))/p,2);
+ }
+}
\ No newline at end of file
lcsim/src/org/lcsim/contrib/seedtracker
diff -N ScatterLayer.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ ScatterLayer.java 11 Dec 2007 20:58:28 -0000 1.1
@@ -0,0 +1,60 @@
+/*
+ * ScatterLayer.java
+ *
+ * Created on August 30, 2007, 1:41 PM
+ *
+ * To change this template, choose Tools | Template Manager
+ * and open the template in the editor.
+ */
+
+package org.lcsim.contrib.seedtracker;
+
+import hep.physics.vec.Hep3Vector;
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
+
+/**
+ *
+ * @author partridge
+ */
+public class ScatterLayer {
+ private MaterialLayer _layer;
+ private double _s;
+ private Hep3Vector _pos;
+ private double _radlen;
+ private double _angle;
+
+ /** Creates a new instance of ScatterLayer */
+ public ScatterLayer(MaterialLayer layer, double p, double slope, double s, Hep3Vector pos) {
+ _layer = layer;
+ _s = s;
+ _pos = pos;
+ if (_layer.beflag() == BarrelEndcapFlag.BARREL) {
+ _radlen = _layer.thickness() * Math.sqrt(1+Math.pow(slope,2)) / _layer.X0();
+ } else if (_layer.beflag() == BarrelEndcapFlag.ENDCAP_NORTH || _layer.beflag() == BarrelEndcapFlag.ENDCAP_NORTH) {
+ _radlen = _layer.thickness() * Math.sqrt(1+Math.pow(slope,-2)) / _layer.X0();
+ }
+ // Theta = (0.0136 GeV)*(z/p)*(sqrt(radLength))*(1+0.038*log(radLength))
+ _angle = 0.0136 * (1.0 / p) * Math.sqrt(_radlen) * (1.0 + 0.038 * Math.log(_radlen));
+ }
+
+ public MaterialLayer layer() {
+ return _layer;
+ }
+
+ public double s() {
+ return _s;
+ }
+
+ public Hep3Vector position() {
+ return _pos;
+ }
+
+ public double radlen() {
+ return _radlen;
+ }
+
+ public double angle() {
+ return _angle;
+ }
+
+}
\ No newline at end of file
lcsim/src/org/lcsim/contrib/seedtracker
diff -u -r1.2 -r1.3
--- HitManager.java 11 Aug 2007 00:17:32 -0000 1.2
+++ HitManager.java 11 Dec 2007 20:58:28 -0000 1.3
@@ -17,6 +17,9 @@
import org.lcsim.event.EventHeader;
import org.lcsim.event.SimTrackerHit;
import org.lcsim.event.TrackerHit;
+import org.lcsim.fit.helicaltrack.HelicalTrackHit;
+import org.lcsim.fit.helicaltrack.HelicalTrack2DHit;
+import org.lcsim.fit.helicaltrack.HelicalTrack3DHit;
import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
/**
@@ -25,11 +28,12 @@
* @version 1.0
*/
public class HitManager {
- HashMap<String, List<TrackerHit>> _hitlist;
- HashMap<String, Double> _rmin;
- HashMap<String, Double> _rmax;
- HashMap<String, Double> _zmin;
- HashMap<String, Double> _zmax;
+ private HitIdentifier ID = new HitIdentifier();
+ private HashMap<String, List<HelicalTrackHit>> _hitlist;
+ private HashMap<String, Double> _rmin;
+ private HashMap<String, Double> _rmax;
+ private HashMap<String, Double> _zmin;
+ private HashMap<String, Double> _zmax;
/** Creates a new instance of HitManager */
public HitManager() {
@@ -41,30 +45,40 @@
* @param event EventHeader for the event to be organized
*/
public void OrganizeHits(EventHeader event) {
- _hitlist = new HashMap<String, List<TrackerHit>>();
+ _hitlist = new HashMap<String, List<HelicalTrackHit>>();
_rmin = new HashMap<String, Double>();
_rmax = new HashMap<String, Double>();
_zmin = new HashMap<String, Double>();
_zmax = new HashMap<String, Double>();
- List<List<TrackerHit>> hitcollections = event.get(TrackerHit.class);
- for (List<TrackerHit> hitcol : hitcollections) {
- for (TrackerHit hit : hitcol) {
- String identifier = LyrIdentifier(hit);
- if (!_hitlist.containsKey(identifier)) {
- _hitlist.put(identifier, new ArrayList<TrackerHit>());
- _rmin.put(identifier,9999999.);
- _rmax.put(identifier,0.);
- _zmin.put(identifier,9999999.);
- _zmax.put(identifier,-9999999.);
- }
- _hitlist.get(identifier).add(hit);
- double r = Math.sqrt(Math.pow(hit.getPosition()[0],2)+Math.pow(hit.getPosition()[1],2));
- double z = hit.getPosition()[2];
- if (r < _rmin.get(identifier)) _rmin.put(identifier,r);
- if (r > _rmax.get(identifier)) _rmax.put(identifier,r);
- if (z < _zmin.get(identifier)) _zmin.put(identifier,z);
- if (z > _zmax.get(identifier)) _zmax.put(identifier,z);
+ List<HelicalTrackHit> hitcol = (List<HelicalTrackHit>) event.get("HelicalTrackHits");
+// for (List<TrackerHit> hitcol : hitcollections) {
+ for (HelicalTrackHit hit : hitcol) {
+ String identifier = ID.Identifier(hit);
+ System.out.println("Hit Found and identified: "+identifier);
+ if (!_hitlist.containsKey(identifier)) {
+ _hitlist.put(identifier, new ArrayList<HelicalTrackHit>());
+ _rmin.put(identifier,9999999.);
+ _rmax.put(identifier,0.);
+ _zmin.put(identifier,9999999.);
+ _zmax.put(identifier,-9999999.);
}
+ List<HelicalTrackHit> hits = _hitlist.get(identifier);
+ hits.add(hit);
+// BarrelEndcapFlag beflag = getBarrelEndcapFlag(hit);
+// if (isPixel(hit)) {
+// hits.add(hit);
+// } else {
+// double zmin = getZMin(hit);
+// double zmax = getZMax(hit);
+// hits.add(new HelicalTrack2DHit(hit, zmin, zmax));
+// }
+// _hitlist.get(identifier).add(new HelicalTrack3DHit(hit, getBarrelEndcapFlag(hit)));
+ double r = Math.sqrt(Math.pow(hit.x(),2)+Math.pow(hit.y(),2));
+ double z = hit.z();
+ if (r < _rmin.get(identifier)) _rmin.put(identifier,r);
+ if (r > _rmax.get(identifier)) _rmax.put(identifier,r);
+ if (z < _zmin.get(identifier)) _zmin.put(identifier,z);
+ if (z > _zmax.get(identifier)) _zmax.put(identifier,z);
}
}
@@ -75,8 +89,9 @@
* @param beflag Barrel-endcap flag
* @return List of matching TrackerHits
*/
- public List<TrackerHit> getTrackerHits(String detname, int layer, BarrelEndcapFlag beflag) {
- return new ArrayList<TrackerHit>(_hitlist.get(Identifier(detname, layer, beflag)));
+ public List<HelicalTrackHit> getTrackerHits(String detname, int layer, BarrelEndcapFlag beflag) {
+// return new ArrayList<HelicalTrackHit>(_hitlist.get(Identifier(detname, layer, beflag)));
+ return _hitlist.get(ID.Identifier(detname, layer, beflag));
}
/**
@@ -84,28 +99,13 @@
* @param seedlayer Seedlayer to look at
* @return List of TrackerHits
*/
- public List<TrackerHit> getTrackerHits(SeedLayer seedlayer) {
- List<TrackerHit> hits = new ArrayList<TrackerHit>();
- String identifier = Identifier(seedlayer.getDetName(),seedlayer.getLayer(),seedlayer.getBarrelEndcapFlag());
- if (_hitlist.containsKey(identifier)) hits.addAll(_hitlist.get(identifier));
- return hits;
+ public List<HelicalTrackHit> getTrackerHits(SeedLayer seedlayer) {
+// List<HelicalTrackHit> hits = new ArrayList<HelicalTrackHit>();
+// String identifier = Identifier(seedlayer.getDetName(),seedlayer.getLayer(),seedlayer.getBarrelEndcapFlag());
+// if (_hitlist.containsKey(identifier)) hits.addAll(_hitlist.get(identifier));
+ return getTrackerHits(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
}
-
- private String LyrIdentifier(TrackerHit hit) {
- if (hit instanceof BaseTrackerHitMC) {
- BaseTrackerHitMC newhit = (BaseTrackerHitMC) hit;
- SimTrackerHit simhit = newhit.getSimHits().get(0);
- String detname = simhit.getSubdetector().getName();
- int layer = simhit.getLayer();
- BarrelEndcapFlag beflag = simhit.getIDDecoder().getBarrelEndcapFlag();
- return Identifier(detname, layer, beflag);
- } else return new String("Unknown");
- }
-
- private String Identifier(String detname, int layer, BarrelEndcapFlag beflag) {
- return new String(detname+layer+beflag.toString());
- }
-
+
/**
* Return the smallest radius hit for a particular SeedLayer
* @param seedlayer SeedLayer to consider
@@ -113,7 +113,7 @@
*/
public double getRMin(SeedLayer seedlayer) {
double rmin = 9999999.;
- String identifier = Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
+ String identifier = ID.Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
if (_rmin.containsKey(identifier)) rmin = Math.min(rmin, _rmin.get(identifier));
return rmin;
}
@@ -125,7 +125,7 @@
*/
public double getRMax(SeedLayer seedlayer) {
double rmax = 0.;
- String identifier = Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
+ String identifier = ID.Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
if (_rmax.containsKey(identifier)) rmax = Math.max(rmax, _rmax.get(identifier));
return rmax;
}
@@ -137,7 +137,7 @@
*/
public double getZMin(SeedLayer seedlayer) {
double zmin = 9999999.;
- String identifier = Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
+ String identifier = ID.Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
if (_zmin.containsKey(identifier)) zmin = Math.min(zmin, _zmin.get(identifier));
return zmin;
}
@@ -149,9 +149,52 @@
*/
public double getZMax(SeedLayer seedlayer) {
double zmax = -9999999.;
- String identifier = Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
+ String identifier = ID.Identifier(seedlayer.getDetName(), seedlayer.getLayer(), seedlayer.getBarrelEndcapFlag());
if (_zmax.containsKey(identifier)) zmax = Math.max(zmax, _zmax.get(identifier));
return zmax;
}
+ private String getName(TrackerHit hit) {
+ if (hit instanceof BaseTrackerHitMC) {
+ SimTrackerHit simhit = ((BaseTrackerHitMC) hit).getSimHits().get(0);
+ return simhit.getSubdetector().getName();
+ } else return new String("Unknown");
+ }
+
+ private int getLayer(TrackerHit hit) {
+ if (hit instanceof BaseTrackerHitMC) {
+ SimTrackerHit simhit = ((BaseTrackerHitMC) hit).getSimHits().get(0);
+ return simhit.getLayer();
+ } else return 999;
+ }
+
+ private BarrelEndcapFlag getBarrelEndcapFlag(TrackerHit hit) {
+ if (hit instanceof BaseTrackerHitMC) {
+ SimTrackerHit simhit = ((BaseTrackerHitMC) hit).getSimHits().get(0);
+ return simhit.getIDDecoder().getBarrelEndcapFlag();
+ } else return BarrelEndcapFlag.UNKNOWN;
+ }
+
+// private boolean isPixel(TrackerHit hit) {
+// if (hit instanceof BaseTrackerHitMC) {
+// return getName(hit).contains("vtx");
+// } else return true;
+// }
+
+// private double getZMin(TrackerHit hit) {
+// return -9999.;
+// }
+
+ private double getZMax(TrackerHit hit) {
+ return 9999.;
+ }
+
+ private String LyrIdentifier(TrackerHit hit) {
+ return ID.Identifier(getName(hit), getLayer(hit), getBarrelEndcapFlag(hit));
+ }
+
+// private String Identifier(String detname, int layer, BarrelEndcapFlag beflag) {
+// return new String(detname+layer+beflag.toString());
+// }
+
}
lcsim/src/org/lcsim/contrib/seedtracker
diff -u -r1.1 -r1.2
--- SeedCandidate.java 10 Aug 2007 15:50:26 -0000 1.1
+++ SeedCandidate.java 11 Dec 2007 20:58:28 -0000 1.2
@@ -14,6 +14,7 @@
import org.lcsim.event.TrackerHit;
import org.lcsim.fit.helicaltrack.HelicalTrackFit;
+import org.lcsim.fit.helicaltrack.HelicalTrackHit;
/**
* Candidate seed containing a list of hits that make up the SeedCandidate and the associated helix parameters
@@ -22,15 +23,15 @@
*/
public class SeedCandidate {
- private List<TrackerHit> _trackerhits = new ArrayList<TrackerHit>();
+ private List<HelicalTrackHit> _hits = new ArrayList<HelicalTrackHit>();
private HelicalTrackFit _helix;
/**
* Create a new SeedCandidate from a list of hits
* @param trackerhits List of tracker hits for the SeedCandidate
*/
- public SeedCandidate(List<TrackerHit> trackerhits) {
- _trackerhits.addAll(trackerhits);
+ public SeedCandidate(List<HelicalTrackHit> hits) {
+ _hits.addAll(hits);
}
/**
@@ -38,8 +39,8 @@
* @param trackerhits List of TrackerHits for the SeedCandidate
* @param helix HelicalTrackFit associated with the SeedCandidate
*/
- public SeedCandidate(List<TrackerHit> trackerhits, HelicalTrackFit helix) {
- this(trackerhits);
+ public SeedCandidate(List<HelicalTrackHit> hits, HelicalTrackFit helix) {
+ this(hits);
_helix = helix;
}
@@ -48,16 +49,16 @@
* @param seed Existing SeedCandidate to be cloned
*/
public SeedCandidate(SeedCandidate seed) {
- this(seed.getTrackerHits(), seed.getHelix());
+ this(seed.getHits(), seed.getHelix());
}
/**
* Assign a list of TrackerHits to the SeedCandidate
* @param trackerhits List of TrackerHits for the SeedCandidate
*/
- public void setTrackerHits(List<TrackerHit> trackerhits) {
- _trackerhits.clear();
- _trackerhits.addAll(trackerhits);
+ public void setTrackerHits(List<HelicalTrackHit> hits) {
+ _hits.clear();
+ _hits.addAll(hits);
return;
}
@@ -74,17 +75,25 @@
* Add a hit to the SeedCandidate
* @param hit TrackerHit to be added to the SeedCandidate
*/
- public void addHit(TrackerHit hit) {
- _trackerhits.add(hit);
+ public void addHit(HelicalTrackHit hit) {
+ _hits.add(hit);
return;
}
+ public List<HelicalTrackHit> getHits() {
+ return _hits;
+ }
+
/**
* Return the list of TrackerHits for the SeedCandidate
* @return List of TrackerHits for the SeedCandidate
*/
public List<TrackerHit> getTrackerHits() {
- return _trackerhits;
+ List<TrackerHit> tracker_hits = new ArrayList<TrackerHit>();
+ for (HelicalTrackHit hel_hit : _hits) {
+ tracker_hits.add(hel_hit.getTrackerHit());
+ }
+ return tracker_hits;
}
/**
lcsim/src/org/lcsim/contrib/seedtracker
diff -u -r1.5 -r1.6
--- SeedTracker.java 11 Aug 2007 01:27:15 -0000 1.5
+++ SeedTracker.java 11 Dec 2007 20:58:28 -0000 1.6
@@ -7,6 +7,8 @@
package org.lcsim.contrib.seedtracker;
+import hep.physics.vec.BasicHep3Vector;
+
import java.util.ArrayList;
import java.util.List;
import java.util.ListIterator;
@@ -21,6 +23,10 @@
//import org.lcsim.event.base.BaseTrackerHitMC;
import org.lcsim.fit.helicaltrack.HelicalTrackFit;
import org.lcsim.fit.helicaltrack.HelicalTrackFitter;
+import org.lcsim.fit.helicaltrack.HelicalTrackHit;
+import org.lcsim.fit.helicaltrack.HelicalTrack2DHit;
+import org.lcsim.fit.helicaltrack.HelicalTrack3DHit;
+import org.lcsim.spacegeom.SpacePoint;
import org.lcsim.util.Driver;
import hep.physics.vec.BasicHep3Vector;
@@ -55,6 +61,7 @@
*/
protected void process(EventHeader event) {
+ MaterialManager mm = new MaterialManager(event);
// Find the magnetic field
// _BField = event.getDetector().getFieldMap().getField(_IP).z();
_BField = event.getDetector().getFieldMap().getField(_OldIP)[2];
@@ -64,6 +71,7 @@
List<SeedCandidate> trackseeds = new ArrayList<SeedCandidate>();
// Loop over all strategies
for (SeedStrategy strategy : _strategylist) {
+ System.out.println(" Strategy "+strategy.getName());
// Find the 3 hits seeds for this strategy
List<SeedCandidate> seedlist = FindSeeds(strategy);
// Loop over the seeds
@@ -92,56 +100,45 @@
// Get the SeedLayers for this strategy
List<SeedLayer> seedlayerlist = strategy.getLayers(SeedType.Seed);
if (seedlayerlist.size() != 3) return candidates;
- // Get the hits for the seed layers
- List<List<TrackerHit>> seedhits = new ArrayList<List<TrackerHit>>();
- for (SeedLayer seedlyr : seedlayerlist) {
- seedhits.add(_HM.getTrackerHits(seedlyr));
- }
- // Order the seed layers from fewest hits to most hits
- for (int i=0; i<2; i++) {
- for (int j=i+1; j<3; j++) {
- if (seedhits.get(i).size() > seedhits.get(j).size()) {
- List<TrackerHit> temp = seedhits.get(i);
- seedhits.set(i,seedhits.get(j));
- seedhits.set(j,temp);
- }
- }
- }
- double RMin = strategy.getMinPT() / (0.0003 * _BField);
- double dMax = strategy.getMaxDCA();
-
- // Todo: add code to limit seeds to those that satisfy cuts
-
- for (TrackerHit s1 : seedhits.get(0)) {
- double[] pos1 = s1.getPosition();
- double r1 = Math.sqrt(pos1[0]*pos1[0] + pos1[1]*pos1[1]);
- double phi1 = Math.atan2(pos1[1],pos1[0]);
- double dphi1 = Math.asin((2*RMin*dMax - dMax*dMax-r1*r1)/(2*r1*(RMin-dMax)));
- for (TrackerHit s2 : seedhits.get(1)) {
- double[] pos2 = s2.getPosition();
- double r2 = Math.sqrt(pos2[0]*pos2[0] + pos2[1]*pos2[1]);
- double phi2 = Math.atan2(pos2[1],pos2[0]);
- double dphi2 = Math.asin((2*RMin*dMax - dMax*dMax-r2*r2)/(2*r2*(RMin-dMax)));
- if (Math.abs(phi2-phi1) > Math.abs(dphi2-dphi1)) continue;
- for (TrackerHit s3 : seedhits.get(2)) {
- double[] pos3 = s3.getPosition();
- double r3 = Math.sqrt(pos3[0]*pos3[0] + pos3[1]*pos3[1]);
- double phi3 = Math.atan2(pos3[1],pos3[0]);
- double dphi3 = Math.asin((2*RMin*dMax - dMax*dMax-r3*r3)/(2*r3*(RMin-dMax)));
- if (Math.abs(phi3-phi2) > Math.abs(dphi3-dphi2)) continue;
- List<TrackerHit> hitlist = new ArrayList<TrackerHit>();
- hitlist.add(s1);
- hitlist.add(s2);
- hitlist.add(s3);
- boolean status = _hf.fit(hitlist);
+// double RMin = strategy.getMinPT() / (0.0003 * _BField);
+// double dMax = strategy.getMaxDCA();
+
+ // Todo: add code to limit seeds to those that satisfy cuts
+ System.out.println("layer 0 hits: "+_HM.getTrackerHits(seedlayerlist.get(0)).size());
+ System.out.println("layer 1 hits: "+_HM.getTrackerHits(seedlayerlist.get(1)).size());
+ System.out.println("layer 2 hits: "+_HM.getTrackerHits(seedlayerlist.get(2)).size());
+
+ for (HelicalTrackHit s1 : _HM.getTrackerHits(seedlayerlist.get(0))) {
+// double[] pos1 = s1.getPosition();
+// double r1 = Math.sqrt(pos1[0]*pos1[0] + pos1[1]*pos1[1]);
+// double phi1 = Math.atan2(pos1[1],pos1[0]);
+// double dphi1 = Math.asin((2*RMin*dMax - dMax*dMax-r1*r1)/(2*r1*(RMin-dMax)));
+ for (HelicalTrackHit s2 : _HM.getTrackerHits(seedlayerlist.get(1))) {
+// double[] pos2 = s2.getPosition();
+// double r2 = Math.sqrt(pos2[0]*pos2[0] + pos2[1]*pos2[1]);
+// double phi2 = Math.atan2(pos2[1],pos2[0]);
+// double dphi2 = Math.asin((2*RMin*dMax - dMax*dMax-r2*r2)/(2*r2*(RMin-dMax)));
+// if (Math.abs(phi2-phi1) > Math.abs(dphi2-dphi1)) continue;
+ for (HelicalTrackHit s3 : _HM.getTrackerHits(seedlayerlist.get(2))) {
+// double[] pos3 = s3.getPosition();
+// double r3 = Math.sqrt(pos3[0]*pos3[0] + pos3[1]*pos3[1]);
+// double phi3 = Math.atan2(pos3[1],pos3[0]);
+// double dphi3 = Math.asin((2*RMin*dMax - dMax*dMax-r3*r3)/(2*r3*(RMin-dMax)));
+// if (Math.abs(phi3-phi2) > Math.abs(dphi3-dphi2)) continue;
+ List<HelicalTrackHit> helixhits = new ArrayList<HelicalTrackHit>();
+ helixhits.add(s1);
+ helixhits.add(s2);
+ helixhits.add(s3);
+ boolean status = _hf.fit(helixhits);
+ if (!status) System.out.println(" failed fit");
if (!status) continue;
HelicalTrackFit helix = _hf.getFit();
if (!HelixCuts(strategy, helix)) continue;
- for (TrackerHit hit : hitlist) {
- System.out.println(" Hit x "+hit.getPosition()[0]+" y "+hit.getPosition()[1]+" z "+hit.getPosition()[2]);
+ for (HelicalTrackHit hit : helixhits) {
+ System.out.println(" Hit x "+hit.x()+" y "+hit.y()+" z "+hit.z());
}
- candidates.add(new SeedCandidate(hitlist, helix));
+ candidates.add(new SeedCandidate(helixhits, helix));
}
}
}
lcsim/src/org/lcsim/contrib/seedtracker
diff -u -r1.3 -r1.4
--- SeedTrackerDriver.java 11 Aug 2007 00:17:32 -0000 1.3
+++ SeedTrackerDriver.java 11 Dec 2007 20:58:28 -0000 1.4
@@ -7,9 +7,16 @@
package org.lcsim.contrib.seedtracker;
-import org.lcsim.contrib.seedtracker.SeedTracker;
-import org.lcsim.contrib.seedtracker.SmearMCHits;
+import java.util.ArrayList;
+import java.util.List;
+
+import org.lcsim.contrib.onoprien.tracking.ExampleDriver2;
+import org.lcsim.fit.helicaltrack.HelicalTrackHitDriver.HitType;
+//import org.lcsim.contrib.seedtracker.SeedTracker;
+//import org.lcsim.contrib.seedtracker.HelicalTrackHitDriver;
+//import org.lcsim.contrib.seedtracker.SmearMCHits;
import org.lcsim.event.EventHeader;
+import org.lcsim.fit.helicaltrack.HelicalTrackHitDriver;
import org.lcsim.util.Driver;
/**
@@ -23,7 +30,11 @@
/** Creates a new instance of SeedTrackerDriver */
public SeedTrackerDriver()
{
- add(new SmearMCHits());
+// add(new SmearMCHits());
+ add(new ExampleDriver2());
+ HelicalTrackHitDriver hitdriver = new HelicalTrackHitDriver();
+ hitdriver.addCollection("StandardTrackerHits",HitType.VirtualSegmentation);
+ add(hitdriver);
add(new SeedTracker());
}
CVSspam 0.2.8