lcsim/src/org/lcsim/contrib/CosminDeaconu
diff -u -r1.1 -r1.2
--- ConnectTheDotsConverter.java 8 Dec 2007 00:23:31 -0000 1.1
+++ ConnectTheDotsConverter.java 18 Jan 2008 07:03:56 -0000 1.2
@@ -22,6 +22,7 @@
import java.util.Collections;
import java.util.Comparator;
import java.util.HashMap;
+import java.util.HashSet;
import java.util.List;
import org.lcsim.contrib.tracking.TrackerHitCheater;
import org.lcsim.event.EventHeader;
@@ -46,14 +47,10 @@
*
* USAGE:
*
- * This converter will operate on SimTrackerHitCollections called "AllSimHits"
- * (or whatever you change the variable COLLECTION_NAME to). You can use the
- * PutAllHitsInOneCollection class in my contrib folder to make this collection
- * for you.
+ * This converter will automatically operate on all the SimHits in the event
+ * (this is not really how these converters are supposed to work, but it's
+ * easier this way...)
- * If for some reason you want lines in between the SimTrackerHits instead of
- * helices, set FORCE_LINE_DRAW to true.
- *
* You may increase the resolution of the helices by modifying D_ALPHA...
*
*
@@ -84,8 +81,7 @@
* Note that the converter might fail if SimTrackerHits don't have momenta
* associated with them. I believe that this feature is enabled by default,
* but if you've changed THBIT_MOMENTUM in LCIOConstants, you might have
- * problems. In that case, you'll probably have to set FORCE_LINE_DRAW to true
- * and you'll only get lines connecting the SimTrackerHits...
+ * problems. .
*
* ISSUES:
*
@@ -106,7 +102,6 @@
private static final double D_ALPHA=5; //5 mm. Make sure to change MAX_SEX_PER_D appropriately too if you change this
private static final double MIN_ENERGY = 0; //allows cut on energy
- private static final String COLLECTION_NAME = "AllSimHits";
private static final String IDENTIFIER_STRING = "ConnectTheDotsMCParticles";
private static final int MAX_SEGS_PER_D = 10;//maximum segments per Euclidean distance before just drawing a line...
// this is sort of an upper limit on complexity and also a safeguard
@@ -118,13 +113,10 @@
private static final Color MUON_COLOR = Color.CYAN; //color for mu+/mu-
private static final Color PROTON_COLOR = Color.LIGHT_GRAY; //color for p, p_bar
private static final Color KAON_COLOR = Color.YELLOW; //color for K+,K-
+ private static final Color UNCHARGED_COLOR = Color.orange; // color for uncharged stuff...
private static final Color OTHER_COLOR = Color.MAGENTA; //color for anything else
private static final int MAX_POINTS = 3000; //this is the max number of points per instance in Wired4
-
- //always draw straight connecting lines, in case you don't have momentum info
- //in SimTrackerHits... or if you really like straight lines
- private static final boolean FORCE_LINE_DRAW = false;
private HepRepType type;
private HepRepType electronType;
@@ -132,20 +124,30 @@
private HepRepType muonType;
private HepRepType protonType;
private HepRepType kaonType;
+ private HepRepType unchargedType;
private HepRepType genericType;
+
+ private HashSet<Integer> event_hashcodes = new HashSet<Integer>();
+
/** Creates a new instance of ConnectTheDotsConverter */
public boolean canHandle(Class k) {
- return SimTrackerHit.class.isAssignableFrom(k);
+ return MCParticle.class.isAssignableFrom(k);
}
- public void convert(EventHeader event, List collection, HepRepFactory factory, HepRepTypeTree typeTree, HepRepInstanceTree instanceTree) {
+
+ //this is kind of cheating...it doesn't take the input collection (since it requires a collection with all of the SimTrackerHits in it)... it ignores it and takes all the SimHits from the event...
+ public void convert(EventHeader event, List ignoreme, HepRepFactory factory, HepRepTypeTree typeTree, HepRepInstanceTree instanceTree) {
- LCMetaData meta = event.getMetaData(collection);
+
+ if (event_hashcodes.contains(event.hashCode())) return; //dont' process the same event more than once in case there's more than one MCParticle collection...
+ event_hashcodes .add(event.hashCode());
+ List<SimTrackerHit> collection = new ArrayList<SimTrackerHit>();
- if (meta.getName()!=COLLECTION_NAME) return; //only process events with the right name
-
+ for (List<SimTrackerHit> l : event.get(SimTrackerHit.class)) collection.addAll(l);
+ LCMetaData meta = event.getMetaData(collection);
+
// set up the instance types
type = factory.createHepRepType(typeTree, IDENTIFIER_STRING);
@@ -157,7 +159,7 @@
type.addAttDef("numhits","Number of Hits","physics","");
type.addAttDef("pt","Initial Particle Transverse Energy","physics","GeV");
type.addAttDef("time","Particle Production Time","physics","nanoseconds");
- type.addAttDef("isCurler","Whether or not particle is a curler, by some definition","physics","");
+ // type.addAttDef("isCurler","Whether or not particle is a curler, by some definition","physics",""); //This doesn't work with the new detector model due to overlapping sensors '
type.addAttDef("particleStatus","the Simulator Status of the particle","physics","");
type.addAttDef("particleOriginRadius","The radius in mm from the ip to the origin of the MCParticle","physics","");
type.addAttDef("particleOriginZ","Z coordinate of particle origin","physics","");
@@ -178,6 +180,9 @@
kaonType = factory.createHepRepType(type,"Kaon");
kaonType.addAttValue("color",KAON_COLOR);
+ unchargedType = factory.createHepRepType(type, "Uncharged");
+ unchargedType.addAttValue("color",UNCHARGED_COLOR);
+
genericType = factory.createHepRepType(type,"Other Particle");
genericType.addAttValue("color",OTHER_COLOR);
@@ -202,7 +207,7 @@
//get Helix for swimming
HelixSwimmer helix = new HelixSwimmer(event.getDetector().getFieldMap().getField(new double[]{0,0,0})[2]);
- // loop through each particle
+ // loop through each particle (this is for hits that have simhits only)
for (MCParticle p : hits_by_particle.keySet()) {
if (p.getEnergy()<MIN_ENERGY) continue;
@@ -217,130 +222,174 @@
//loop through simhits
for (int i = 0; i< hitlist.size(); i++) {
SimTrackerHit h = hitlist.get(i);
-
- //if flag set, connect the hits with lines
- if (FORCE_LINE_DRAW) {
- //draw the first point if first iteration
- if (i==0) instance = makePoint(context,instance,p.getOrigin().x(),p.getOrigin().y(),p.getOrigin().z());
- instance = makePoint(context,instance,h.getPoint()[0],h.getPoint()[1],h.getPoint()[2]);
-
- //if last point, try to draw line to final point
- if (i==hitlist.size()-1) {
- try {
- double[] end = p.getEndPoint().v();
- instance = makePoint(context,instance,end[0],end[1],end[2]);
- }
- catch(RuntimeException re) {
- //for now, this does nothing... presumably we'd want to draw a line to the tracker boundary or something
+ double zStop; //swim the helix till we get here in z...
+ double charge = p.getCharge();
+
+ //first hit
+ if (i==0) {//draw from origin to first point if it's the first point
+ helix.setTrack(p.getMomentum(),new SpacePoint(p.getOrigin()),(int) charge);
+ zStop = hitlist.get(0).getPoint()[2];
+
+ int segs =1+(int) (helix.getDistanceToZ(zStop)/D_ALPHA);
+ double lDist = VecOp.sub(new BasicHep3Vector(hitlist.get(0).getPoint()),p.getOrigin()).magnitude();
+
+ //if segs exceeds MAX_SEGS...we just draw a straight line...
+ if (segs<MAX_SEGS_PER_D*lDist && segs>0) {
+ for (int j = 0; j<=segs; j++) {
+ SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
+ instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
}
}
+ else {
+ instance = makePoint(context,instance,p.getOriginX(),p.getOriginY(),p.getOriginZ());
+ instance = makePoint(context,instance,h.getPoint()[0],h.getPoint()[1],h.getPoint()[2]);
+ }
}
- //otherwise, connect hits with helices
- else {
-
- double zStop; //swim the helix till we get here in z...
- double charge = p.getCharge();
+ //case where we're not dealing with the last hit
+ if (i<hitlist.size()-1) {
- //case where we're not dealing with the last hit
- if (i<hitlist.size()-1) {
-
- if (i==0) {//draw from origin to first point if it's the first point
- helix.setTrack(p.getMomentum(),new SpacePoint(p.getOrigin()),(int) charge);
- zStop = hitlist.get(0).getPoint()[2];
-
- int segs =1+(int) (helix.getDistanceToZ(zStop)/D_ALPHA);
- double lDist = VecOp.sub(new BasicHep3Vector(hitlist.get(0).getPoint()),p.getOrigin()).magnitude();
-
- //if segs exceeds MAX_SEGS...we just draw a straight line...
- if (segs<MAX_SEGS_PER_D*lDist && segs>0) {
- for (int j = 0; j<=segs; j++) {
- SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
- instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
- }
- }
- else {
- instance = makePoint(context,instance,p.getOriginX(),p.getOriginY(),p.getOriginZ());
- instance = makePoint(context,instance,h.getPoint()[0],h.getPoint()[1],h.getPoint()[2]);
- }
+ zStop = hitlist.get(i+1).getPoint()[2];
+
+ Hep3Vector momentum = new BasicHep3Vector(h.getMomentum());
+ SpacePoint start = new SpacePoint(new BasicHep3Vector(h.getPoint()));
+
+ helix.setTrack(momentum,start,(int)charge);
+
+ int segs = 1+(int) (helix.getDistanceToZ(zStop)/D_ALPHA);
+ double lDist = VecOp.sub(new BasicHep3Vector(hitlist.get(i+1).getPoint()),start).magnitude();
+
+ if (segs<MAX_SEGS_PER_D*lDist) {
+ for (int j = 0; j<=segs; j++) {
+ SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
+ instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
}
-
- zStop = hitlist.get(i+1).getPoint()[2];
-
- Hep3Vector momentum = new BasicHep3Vector(h.getMomentum());
- SpacePoint start = new SpacePoint(new BasicHep3Vector(h.getPoint()));
-
- helix.setTrack(momentum,start,(int)charge);
-
- int segs = 1+(int) (helix.getDistanceToZ(zStop)/D_ALPHA);
- double lDist = VecOp.sub(new BasicHep3Vector(hitlist.get(i+1).getPoint()),start).magnitude();
-
- if (segs<MAX_SEGS_PER_D*lDist) {
+ }
+ else {
+ instance = makePoint(context,instance,h.getPoint()[0],h.getPoint()[1],h.getPoint()[2]);
+ instance = makePoint(context,instance,hitlist.get(i+1).getPoint()[0],hitlist.get(i+1).getPoint()[1],hitlist.get(i+1).getPoint()[2]);
+ }
+ }
+
+ //case where we ARE dealing with the last hit
+ else {
+ Hep3Vector momentum = new BasicHep3Vector(h.getMomentum());
+ SpacePoint start = new SpacePoint(new BasicHep3Vector(h.getPoint()));
+ helix.setTrack(momentum,start,(int)charge);
+
+ double r = event.getDetector().getConstants().get("tracking_region_radius").getValue();
+ double z = event.getDetector().getConstants().get("tracking_region_zmax").getValue();
+
+ try{
+ zStop = p.getEndPoint().z(); //this might thrown an exception if there is no endpoint defined
+
+ int segs = (int) Math.abs(helix.getDistanceToZ(zStop)/D_ALPHA);
+ double lDist = VecOp.sub(start,p.getEndPoint()).magnitude();
+
+ if (segs<MAX_SEGS_PER_D*lDist && segs>0) {
for (int j = 0; j<=segs; j++) {
SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
+
+ if (FORCE_STOP_AT_TRACKING_REGION_BOUNDARY) {if (pt.rxy()>r || Math.abs(pt.z()) > z) break;}
}
}
+
else {
instance = makePoint(context,instance,h.getPoint()[0],h.getPoint()[1],h.getPoint()[2]);
- instance = makePoint(context,instance,hitlist.get(i+1).getPoint()[0],hitlist.get(i+1).getPoint()[1],hitlist.get(i+1).getPoint()[2]);
+ instance = makePoint(context,instance,p.getEndPoint().x(),p.getEndPoint().y(),p.getEndPoint().z());
}
}
-
- //case where we ARE dealing with the last hit
- else {
- Hep3Vector momentum = new BasicHep3Vector(h.getMomentum());
- SpacePoint start = new SpacePoint(new BasicHep3Vector(h.getPoint()));
- helix.setTrack(momentum,start,(int)charge);
-
- double r = event.getDetector().getConstants().get("tracking_region_radius").getValue();
- double z = event.getDetector().getConstants().get("tracking_region_zmax").getValue();
-
- try{
- zStop = p.getEndPoint().z(); //this might thrown an exception if there is no endpoint defined
-
- int segs = (int) Math.abs(helix.getDistanceToZ(zStop)/D_ALPHA);
- double lDist = VecOp.sub(start,p.getEndPoint()).magnitude();
-
- if (segs<MAX_SEGS_PER_D*lDist && segs>0) {
- for (int j = 0; j<=segs; j++) {
- SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
- instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
-
- if (FORCE_STOP_AT_TRACKING_REGION_BOUNDARY) {if (pt.rxy()>r || Math.abs(pt.z()) > z) break;}
- }
- }
-
- else {
- instance = makePoint(context,instance,h.getPoint()[0],h.getPoint()[1],h.getPoint()[2]);
- instance = makePoint(context,instance,p.getEndPoint().x(),p.getEndPoint().y(),p.getEndPoint().z());
- }
- }
-
- //if no endpoint...
- catch (RuntimeException re) {
- int j = 0;
- int segs = 0;
- while (true) {
- SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
- double lDist = VecOp.sub(new BasicHep3Vector(h.getPoint()),new BasicHep3Vector(pt.v())).magnitudeSquared();
-
- if (segs>MAX_SEGS_PER_D*lDist) break;
- if (pt.rxy()>r || Math.abs(pt.z()) > z) break;
-
- instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
- j++;
- segs ++;
- }
+
+ //if no endpoint...
+ catch (RuntimeException re) {
+ int j = 0;
+ int segs = 0;
+ while (true) {
+ SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
+ double lDist = VecOp.sub(new BasicHep3Vector(h.getPoint()),new BasicHep3Vector(pt.v())).magnitudeSquared();
+
+ if (segs>MAX_SEGS_PER_D*lDist) break;
+ if (pt.rxy()>r || Math.abs(pt.z()) > z) break;
+
+ instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
+ j++;
+ segs ++;
}
- }//else last hit
- }//else draws helices
+ }
+ }//else last hit
}// loop through simhits
}//loop through particles
- }//end convert()
+
+ //get MCParticles that have no simhits...
+
+ List<MCParticle> allMCs = event.getMCParticles();
+ for (MCParticle p : allMCs) {
+
+ if (hits_by_particle.keySet().contains(p)) continue; //ignore particles that have hits
+
+ InstanceContext context = new InstanceContext(p,new ArrayList<SimTrackerHit>(),factory,instanceTree); //convenience class
+ HepRepInstance instance = getInstance(context); //create the instance
+ double charge = p.getCharge();
+ makePoint(context, instance,p.getOriginX(), p.getOriginY(), p.getOriginZ());
+
+ if (charge == 0) {
+ try {
+ makePoint(context, instance, p.getEndPoint().x(), p.getEndPoint().y(), p.getEndPoint().z());
+ continue;
+ }
+ catch (RuntimeException re) {}
+ }
+
+
+ helix.setTrack(p.getMomentum(),new SpacePoint(p.getOrigin()), (int) charge);
+ Hep3Vector start = p.getOrigin();
+
+ double r = event.getDetector().getConstants().get("tracking_region_radius").getValue();
+ double z = event.getDetector().getConstants().get("tracking_region_zmax").getValue();
+
+ try{
+ double zStop = p.getEndPoint().z(); //this might thrown an exception if there is no endpoint defined
+
+ int segs = (int) Math.abs(helix.getDistanceToZ(zStop)/D_ALPHA);
+ double lDist = VecOp.sub(start,p.getEndPoint()).magnitude();
+
+ if (segs<MAX_SEGS_PER_D*lDist && segs>0) {
+ for (int j = 0; j<=segs; j++) {
+ SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
+ instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
+
+ if (FORCE_STOP_AT_TRACKING_REGION_BOUNDARY) {if (pt.rxy()>r || Math.abs(pt.z()) > z) break;}
+ }
+ }
+
+ else instance = makePoint(context,instance,p.getEndPoint().x(),p.getEndPoint().y(),p.getEndPoint().z());
+ }
+
+ //if no endpoint...
+ catch (RuntimeException re) {
+ int j = 0;
+ int segs = 0;
+ while (true) {
+ SpacePoint pt = helix.getPointAtLength(D_ALPHA*j);
+ double lDist = VecOp.sub(start,new BasicHep3Vector(pt.v())).magnitudeSquared();
+
+ if (segs>MAX_SEGS_PER_D*lDist) break;
+ if (pt.rxy()>r || Math.abs(pt.z()) > z) break;
+
+ instance = makePoint(context,instance,pt.x(),pt.y(),pt.z());
+ j++;
+ segs ++;
+ }
+ }
+
+ }
+
+ }//end convert()
+ //THIS NO LONGER WORKS WITH THE NEW DETECTORS...
//int so you can cut on it...
int isCurler(ArrayList<SimTrackerHit> hitlist) {
@@ -373,7 +422,7 @@
return buff.toString();
}
- //wrapper class around createHepRepPoint that will create a new instance if too many points. Returns the instance for the next points to build on.
+ //wrapper function around createHepRepPoint that will create a new instance if too many points. Returns the instance for the next points to build on.
private HepRepInstance makePoint(InstanceContext context, HepRepInstance instance, double x, double y, double z) {
if (instance.getPoints().size()<MAX_POINTS) {
@@ -401,6 +450,7 @@
else if (ptype.equals("pi-")||ptype.equals("pi+")) instance = c.factory.createHepRepInstance(c.instanceTree,pionType);
else if (ptype.equals("p")||ptype.equals("p_bar")) instance = c.factory.createHepRepInstance(c.instanceTree,protonType);
else if (ptype.equals("K+")||ptype.equals("K-")) instance = c.factory.createHepRepInstance(c.instanceTree,kaonType);
+ else if (c.p.getCharge()==0) instance = c.factory.createHepRepInstance(c.instanceTree,unchargedType);
else instance = c.factory.createHepRepInstance(c.instanceTree,genericType);
double x = c.p.getMomentum().x();
@@ -414,7 +464,7 @@
instance.addAttValue("numhits",c.hitlist.size());
instance.addAttValue("time",c.p.getProductionTime());
instance.addAttValue("hashcode",c.p.hashCode());
- instance.addAttValue("isCurler",isCurler(c.hitlist));
+ // instance.addAttValue("isCurler",isCurler(c.hitlist)); //this doesn't work with the new detector due to overlapping sensors'
instance.addAttValue("particleOriginRadius",Math.sqrt(c.p.getOriginX()*c.p.getOriginX()+c.p.getOriginY()+c.p.getOriginY()));
instance.addAttValue("particleOriginZ",c.p.getOriginZ());
instance.addAttValue("particleStatus",getSimulatorString(c.p.getSimulatorStatus()));
@@ -457,6 +507,10 @@
return (l.number()==layerNumber && l.name().equals(detectorName));
}
+ public int hashcode() {
+ return layerNumber*10000000 + (int) detectorName.hashCode()/100;
+ }
+
public int number() {
return layerNumber;
}
@@ -467,9 +521,6 @@
}
-
-
-
class TimeCompare<SimTrackerHit> implements Comparator<SimTrackerHit>
{
public int compare(Object a, Object b)