Skip to content

Commit 9312b2b

Browse files
author
Matthew Solt
committed
tuple updates
1 parent 2898e31 commit 9312b2b

3 files changed

Lines changed: 271 additions & 11 deletions

File tree

analysis/src/main/java/org/hps/analysis/MC/MCFullDetectorTruth.java

Lines changed: 56 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -44,6 +44,62 @@ public static Map<MCParticle, List<SimTrackerHit>> BuildTrackerHitMap(List<SimTr
4444
return trackerHitMap;
4545
}
4646

47+
public static Map<MCParticle, Map<String, List<SimTrackerHit>>> BuildComboTrackerHitMap(List<SimTrackerHit> trackerHitsAct, List<SimTrackerHit> trackerHitsIn){
48+
Map<MCParticle,Map<String, List<SimTrackerHit>>> trackerNewHitMap = new HashMap<MCParticle,Map<String, List<SimTrackerHit>>>();
49+
for (SimTrackerHit hit : trackerHitsAct) {
50+
MCParticle p = hit.getMCParticle();
51+
if (p == null) {
52+
throw new RuntimeException("Tracker hit points to null MCParticle!");
53+
}
54+
if (trackerNewHitMap.get(p) == null) {
55+
trackerNewHitMap.put(p, new HashMap<String, List<SimTrackerHit>>());
56+
trackerNewHitMap.get(p).put("active", new ArrayList<SimTrackerHit>());
57+
}
58+
trackerNewHitMap.get(p).get("active").add(hit);
59+
}
60+
int i = 0;
61+
for (SimTrackerHit hit : trackerHitsIn) {
62+
if (i == trackerHitsIn.size()/2) break;
63+
MCParticle p = hit.getMCParticle();
64+
if (p == null) {
65+
throw new RuntimeException("Tracker hit points to null MCParticle!");
66+
}
67+
if (trackerNewHitMap.get(p) == null) {
68+
trackerNewHitMap.put(p, new HashMap<String, List<SimTrackerHit>>());
69+
trackerNewHitMap.get(p).put("inactive", new ArrayList<SimTrackerHit>());
70+
}
71+
if(trackerNewHitMap.get(p).get("inactive") == null){
72+
trackerNewHitMap.get(p).put("inactive", new ArrayList<SimTrackerHit>());
73+
}
74+
trackerNewHitMap.get(p).get("inactive").add(hit);
75+
i++;
76+
}
77+
return trackerNewHitMap;
78+
}
79+
80+
/*public static Map<MCParticle,List<SimTrackerHit>> BuildComboTrackerHitMap(Map<MCParticle, List<SimTrackerHit>> ActiveHitMap, Map<MCParticle, List<SimTrackerHit>> InActiveHitMap){
81+
Map<MCParticle, List<SimTrackerHit>> trackerHitMap = new HashMap<MCParticle, List<SimTrackerHit>>();
82+
for (Entry<MCParticle, List<SimTrackerHit>> entry_act : ActiveHitMap.entrySet()) {
83+
MCParticle p_act = entry_act.getKey();
84+
List<SimTrackerHit> hits_act = entry_act.getValue();
85+
for (Entry<MCParticle, List<SimTrackerHit>> entry_in : InActiveHitMap.entrySet()) {
86+
MCParticle p_in = entry_in.getKey();
87+
if(p_act != p_in) continue;
88+
List<SimTrackerHit> hits_in = entry_in.getValue();
89+
if (trackerHitMap.get(p_act) == null) {
90+
trackerHitMap.put(p_act, new ArrayList<SimTrackerHit>());
91+
}
92+
for (SimTrackerHit hit : hits_act) {
93+
trackerHitMap.get(p_act).add(hit);
94+
}
95+
for (SimTrackerHit hit : hits_in) {
96+
trackerHitMap.get(p_act).add(hit);
97+
}
98+
}
99+
}
100+
return trackerHitMap;
101+
}*/
102+
47103
public static Map<MCParticle, List<SimCalorimeterHit>> BuildCalHitMap(List<SimCalorimeterHit> calHits){
48104
// map particle to a list of its sim cal hits
49105
Map<MCParticle, List<SimCalorimeterHit>> calHitMap = new HashMap<MCParticle, List<SimCalorimeterHit>>();

analysis/src/main/java/org/hps/analysis/tuple/TridentTupleDriver.java

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -35,6 +35,7 @@ public class TridentTupleDriver extends TupleDriver {
3535
private final double tupleMaxSumCut = 1.3;
3636
private boolean getMC = false;
3737
private boolean getFullTruth = false;
38+
private int tupleevent = 0;
3839

3940
@Override
4041
protected void setupVariables() {
@@ -64,6 +65,8 @@ public void setGetFullTruth(boolean getFullTruth) {
6465
@Override
6566
public void process(EventHeader event) {
6667
/* make sure everything is there */
68+
tupleevent++;
69+
//tupleMap.put("tupleevent/I",(double) tupleevent);
6770
if (!event.hasCollection(ReconstructedParticle.class, unconstrainedV0CandidatesColName)) {
6871
return;
6972
}
@@ -115,6 +118,7 @@ public void process(EventHeader event) {
115118
tupleMap.put("minPositiveIso/D", minPositiveIso);
116119
tupleMap.put("minNegativeIso/D", minNegativeIso);
117120
tupleMap.put("minIso/D", minIso);
121+
tupleMap.put("tupleevent/I",(double) tupleevent);
118122

119123
if (getMC) {
120124
fillMCTridentVariables(event);

analysis/src/main/java/org/hps/analysis/tuple/TupleDriver.java

Lines changed: 211 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -95,7 +95,7 @@ public abstract class TupleDriver extends Driver {
9595
protected double ebeam = Double.NaN;
9696
private int nLay = 6;
9797
private int nEcalHit = 3;
98-
private int tupleevent = 0;
98+
//private int tupleevent = 0;
9999
FieldMap bFieldMap = null;
100100

101101
public void setNLay(int nLay) {
@@ -340,14 +340,14 @@ protected void addParticleVariables(String prefix) {
340340
protected void fillTruthEventVariables(EventHeader event) {
341341
tupleMap.put("run/I", (double) event.getRunNumber());
342342
tupleMap.put("event/I", (double) event.getEventNumber());
343-
tupleMap.put("tupleevent/I", (double) tupleevent);
344-
tupleevent++;
343+
//tupleMap.put("tupleevent/I", (double) tupleevent);
344+
//tupleevent++;
345345
}
346346
protected void fillEventVariables(EventHeader event, TIData triggerData) {
347347
tupleMap.put("run/I", (double) event.getRunNumber());
348348
tupleMap.put("event/I", (double) event.getEventNumber());
349-
tupleMap.put("tupleevent/I", (double) tupleevent);
350-
tupleevent++;
349+
//tupleMap.put("tupleevent/I", (double) tupleevent);
350+
//tupleevent++;
351351
List<ReconstructedParticle> fspList = event.get(ReconstructedParticle.class, finalStateParticlesColName);
352352
int npos = 0;
353353
int ntrk = 0;
@@ -1489,7 +1489,7 @@ protected void addFullTruthVertexVariables() {
14891489

14901490
protected void addMCParticleVariables(String prefix) {
14911491
String[] newVars = new String[] {"StartX/D", "StartY/D", "StartZ/D", "EndX/D", "EndY/D", "EndZ/D", "PX/D",
1492-
"PY/D", "PZ/D", "P/D", "M/D", "E/D","pdgid/I","parentID/I","HasTruthMatch/I"};
1492+
"PY/D", "PZ/D", "P/D", "M/D", "E/D","pdgid/I","parentID/I","HasTruthMatch/I","NTruthHits/I","NBadTruthHits/I","Purity/D"};
14931493

14941494
for (int i = 0; i < newVars.length; i++) {
14951495
newVars[i] = prefix + newVars[i];
@@ -1897,6 +1897,7 @@ protected void fillFullVertexTruth(EventHeader event, Track eleTrack, Track posT
18971897
Map<MCParticle, List<SimTrackerHit>> trackerHitMap = MCFullDetectorTruth.BuildTrackerHitMap(trackerHits);
18981898
Map<MCParticle, List<SimTrackerHit>> trackerInHitMap = MCFullDetectorTruth.BuildTrackerHitMap(trackerHits_Inactive);
18991899
Map<MCParticle, List<SimCalorimeterHit>> calHitMap = MCFullDetectorTruth.BuildCalHitMap(calHits);
1900+
Map<MCParticle, Map<String, List<SimTrackerHit>>> comboTrackerHitMap = MCFullDetectorTruth.BuildComboTrackerHitMap(trackerHits,trackerHits_Inactive);
19001901

19011902
TrackTruthMatching eleTruth = new TrackTruthMatching(eleTrack, hittomc);
19021903
TrackTruthMatching posTruth = new TrackTruthMatching(posTrack, hittomc);
@@ -1917,8 +1918,106 @@ protected void fillFullVertexTruth(EventHeader event, Track eleTrack, Track posT
19171918
if(ele != null){
19181919
String MCprefix = "ele";
19191920
tupleMap.put(MCprefix+"HasTruthMatch/I", (double) 1);
1921+
tupleMap.put(MCprefix+"NTruthHits/I", (double) eleTruth.getNHits());
1922+
tupleMap.put(MCprefix+"NBadTruthHits/I", (double) eleTruth.getNBadHits());
1923+
tupleMap.put(MCprefix+"Purity/D", eleTruth.getPurity());
19201924
fillMCParticleVariables("ele", ele);
1921-
for (Entry<MCParticle, List<SimTrackerHit>> entry : trackerInHitMap.entrySet()) {
1925+
for (Entry<MCParticle, Map<String, List<SimTrackerHit>>> entry : comboTrackerHitMap.entrySet()) {
1926+
MCParticle p = entry.getKey();
1927+
if (ele != p) continue;
1928+
List<SimTrackerHit> hits_act = entry.getValue().get("active");
1929+
List<SimTrackerHit> hits_in = entry.getValue().get("inactive");
1930+
1931+
Hep3Vector startPosition = p.getOrigin();
1932+
Hep3Vector startMomentum = p.getMomentum();
1933+
1934+
int i = 0;
1935+
int k = 0;
1936+
String layerprefix_1 = "";
1937+
// loop over particle's hits
1938+
boolean inactiveprev = false;
1939+
for (SimTrackerHit hit : hits_act) {
1940+
boolean inactive = false;
1941+
SimTrackerHit hit_act = hit;
1942+
do{
1943+
k++;
1944+
inactive = false;
1945+
trackerDecoder.setID(hit.getCellID64());
1946+
if(hits_in != null){
1947+
int j = 0;
1948+
for(SimTrackerHit hit_in:hits_in){
1949+
j++;
1950+
if(i >= j) continue;
1951+
if(hit_in.getLayer() >= hit_act.getLayer()) continue;
1952+
hit = hit_in;
1953+
inactive = true;
1954+
i = j;
1955+
break;
1956+
}
1957+
}
1958+
1959+
if(!inactive) hit = hit_act;
1960+
1961+
String layerprefix = MCFullDetectorTruth.trackHitLayer(hit);
1962+
//String layerprefix_1 = MCFullDetectorTruth.trackHitLayer_1(hit);
1963+
String prefix = "";
1964+
String prefix_1 = "";
1965+
if(!inactive){
1966+
prefix = MCprefix + layerprefix;
1967+
}
1968+
else{
1969+
prefix = MCprefix + layerprefix + "In";
1970+
}
1971+
if(!inactiveprev){
1972+
prefix_1 = MCprefix + layerprefix_1;
1973+
}
1974+
else{
1975+
prefix_1 = MCprefix + layerprefix_1 +"In";
1976+
}
1977+
1978+
Hep3Vector endPosition = hit.getPositionVec();
1979+
Hep3Vector endMomentum = new BasicHep3Vector(hit.getMomentum()[0], hit.getMomentum()[1], hit.getMomentum()[2]);
1980+
1981+
inactiveprev = inactive;
1982+
layerprefix_1 = layerprefix;
1983+
1984+
if(k == 1){
1985+
startPosition = endPosition;
1986+
startMomentum = endMomentum;
1987+
continue;
1988+
}
1989+
1990+
double q = p.getCharge();
1991+
Hep3Vector extrapPos = MCFullDetectorTruth.extrapolateTrackPosition(endPosition,endMomentum,startPosition,5,bFieldMap,q);
1992+
Hep3Vector extrapP = MCFullDetectorTruth.extrapolateTrackMomentum(endPosition,endMomentum,startPosition,5,bFieldMap,q);
1993+
1994+
if(extrapP == null) continue;
1995+
1996+
Hep3Vector startProt = VecOp.mult(beamAxisRotation, startMomentum);
1997+
Hep3Vector extrapProt = VecOp.mult(beamAxisRotation, VecOp.neg(extrapP));
1998+
Hep3Vector endPositionrot = VecOp.mult(beamAxisRotation, endPosition);
1999+
2000+
double thetaX = MCFullDetectorTruth.deltaThetaX(extrapProt,startProt);
2001+
double thetaY = MCFullDetectorTruth.deltaThetaY(extrapProt,startProt);
2002+
2003+
tupleMap.put(prefix+"svthitX/D", endPositionrot.x());
2004+
tupleMap.put(prefix+"svthitY/D", endPositionrot.y());
2005+
tupleMap.put(prefix+"svthitZ/D", endPositionrot.z());
2006+
tupleMap.put(prefix_1+"svthitPx/D", startProt.x());
2007+
tupleMap.put(prefix_1+"svthitPy/D", startProt.y());
2008+
tupleMap.put(prefix_1+"svthitPz/D", startProt.z());
2009+
tupleMap.put(prefix_1+"thetaX/D", thetaX);
2010+
tupleMap.put(prefix_1+"thetaY/D", thetaY);
2011+
tupleMap.put(prefix_1+"residualX/D", startPosition.x()-extrapPos.x());
2012+
tupleMap.put(prefix_1+"residualY/D", startPosition.y()-extrapPos.y());
2013+
2014+
startPosition = endPosition;
2015+
startMomentum = endMomentum;
2016+
} while(inactive);
2017+
}
2018+
break;
2019+
}
2020+
/*for (Entry<MCParticle, List<SimTrackerHit>> entry : trackerInHitMap.entrySet()) {
19222021
MCParticle p = entry.getKey();
19232022
if (ele != p) continue;
19242023
@@ -1929,7 +2028,6 @@ protected void fillFullVertexTruth(EventHeader event, Track eleTrack, Track posT
19292028
19302029
// loop over particle's hits
19312030
for (SimTrackerHit hit : hits) {
1932-
19332031
trackerDecoder.setID(hit.getCellID64());
19342032
19352033
String layerprefix = MCFullDetectorTruth.trackHitLayer(hit);
@@ -2030,7 +2128,7 @@ protected void fillFullVertexTruth(EventHeader event, Track eleTrack, Track posT
20302128
startMomentum = endMomentum;
20312129
}
20322130
break;
2033-
}
2131+
}*/
20342132
for (Entry<MCParticle, List<SimCalorimeterHit>> entry : calHitMap.entrySet()) {
20352133
MCParticle p = entry.getKey();
20362134
if (ele != p) continue;
@@ -2060,8 +2158,110 @@ protected void fillFullVertexTruth(EventHeader event, Track eleTrack, Track posT
20602158
if(pos != null){
20612159
String MCprefix = "pos";
20622160
tupleMap.put(MCprefix+"HasTruthMatch/I", (double) 1);
2161+
tupleMap.put(MCprefix+"HasTruthMatch/I", (double) 1);
2162+
tupleMap.put(MCprefix+"NTruthHits/I", (double) posTruth.getNHits());
2163+
tupleMap.put(MCprefix+"NBadTruthHits/I", (double) posTruth.getNBadHits());
2164+
tupleMap.put(MCprefix+"Purity/D", posTruth.getPurity());
20632165
fillMCParticleVariables("pos", pos);
2064-
for (Entry<MCParticle, List<SimTrackerHit>> entry : trackerInHitMap.entrySet()) {
2166+
2167+
for (Entry<MCParticle, Map<String, List<SimTrackerHit>>> entry : comboTrackerHitMap.entrySet()) {
2168+
MCParticle p = entry.getKey();
2169+
if (pos != p) continue;
2170+
List<SimTrackerHit> hits_act = entry.getValue().get("active");
2171+
List<SimTrackerHit> hits_in = entry.getValue().get("inactive");
2172+
2173+
Hep3Vector startPosition = p.getOrigin();
2174+
Hep3Vector startMomentum = p.getMomentum();
2175+
2176+
int i = 0;
2177+
int k = 0;
2178+
String layerprefix_1 = "";
2179+
// loop over particle's hits
2180+
boolean inactiveprev = false;
2181+
for (SimTrackerHit hit : hits_act) {
2182+
boolean inactive = false;
2183+
SimTrackerHit hit_act = hit;
2184+
do{
2185+
k++;
2186+
inactive = false;
2187+
trackerDecoder.setID(hit.getCellID64());
2188+
if(hits_in != null){
2189+
int j = 0;
2190+
for(SimTrackerHit hit_in:hits_in){
2191+
j++;
2192+
if(i >= j) continue;
2193+
if(hit_in.getLayer() >= hit_act.getLayer()) continue;
2194+
hit = hit_in;
2195+
inactive = true;
2196+
i = j;
2197+
break;
2198+
}
2199+
}
2200+
2201+
if(!inactive) hit = hit_act;
2202+
2203+
String layerprefix = MCFullDetectorTruth.trackHitLayer(hit);
2204+
//String layerprefix_1 = MCFullDetectorTruth.trackHitLayer_1(hit);
2205+
String prefix = "";
2206+
String prefix_1 = "";
2207+
if(!inactive){
2208+
prefix = MCprefix + layerprefix;
2209+
}
2210+
else{
2211+
prefix = MCprefix + layerprefix + "In";
2212+
}
2213+
if(!inactiveprev){
2214+
prefix_1 = MCprefix + layerprefix_1;
2215+
}
2216+
else{
2217+
prefix_1 = MCprefix + layerprefix_1 +"In";
2218+
}
2219+
2220+
Hep3Vector endPosition = hit.getPositionVec();
2221+
Hep3Vector endMomentum = new BasicHep3Vector(hit.getMomentum()[0], hit.getMomentum()[1], hit.getMomentum()[2]);
2222+
2223+
inactiveprev = inactive;
2224+
layerprefix_1 = layerprefix;
2225+
2226+
if(k == 1){
2227+
startPosition = endPosition;
2228+
startMomentum = endMomentum;
2229+
continue;
2230+
}
2231+
2232+
double q = p.getCharge();
2233+
Hep3Vector extrapPos = MCFullDetectorTruth.extrapolateTrackPosition(endPosition,endMomentum,startPosition,5,bFieldMap,q);
2234+
Hep3Vector extrapP = MCFullDetectorTruth.extrapolateTrackMomentum(endPosition,endMomentum,startPosition,5,bFieldMap,q);
2235+
2236+
if(extrapP == null) continue;
2237+
2238+
Hep3Vector startProt = VecOp.mult(beamAxisRotation, startMomentum);
2239+
Hep3Vector extrapProt = VecOp.mult(beamAxisRotation, VecOp.neg(extrapP));
2240+
Hep3Vector endPositionrot = VecOp.mult(beamAxisRotation, endPosition);
2241+
2242+
double thetaX = MCFullDetectorTruth.deltaThetaX(extrapProt,startProt);
2243+
double thetaY = MCFullDetectorTruth.deltaThetaY(extrapProt,startProt);
2244+
2245+
tupleMap.put(prefix+"svthitX/D", endPositionrot.x());
2246+
tupleMap.put(prefix+"svthitY/D", endPositionrot.y());
2247+
tupleMap.put(prefix+"svthitZ/D", endPositionrot.z());
2248+
tupleMap.put(prefix_1+"svthitPx/D", startProt.x());
2249+
tupleMap.put(prefix_1+"svthitPy/D", startProt.y());
2250+
tupleMap.put(prefix_1+"svthitPz/D", startProt.z());
2251+
tupleMap.put(prefix_1+"thetaX/D", thetaX);
2252+
tupleMap.put(prefix_1+"thetaY/D", thetaY);
2253+
tupleMap.put(prefix_1+"residualX/D", startPosition.x()-extrapPos.x());
2254+
tupleMap.put(prefix_1+"residualY/D", startPosition.y()-extrapPos.y());
2255+
2256+
startPosition = endPosition;
2257+
startMomentum = endMomentum;
2258+
} while(inactive);
2259+
}
2260+
break;
2261+
2262+
2263+
2264+
/*for (Entry<MCParticle, List<SimTrackerHit>> entry : trackerInHitMap.entrySet()) {
20652265
MCParticle p = entry.getKey();
20662266
if (pos != p) continue;
20672267
List<SimTrackerHit> hits = entry.getValue();
@@ -2171,7 +2371,7 @@ protected void fillFullVertexTruth(EventHeader event, Track eleTrack, Track posT
21712371
startPosition = endPosition;
21722372
startMomentum = endMomentum;
21732373
}
2174-
break;
2374+
break;*/
21752375
}
21762376
for (Entry<MCParticle, List<SimCalorimeterHit>> entry : calHitMap.entrySet()) {
21772377
MCParticle p = entry.getKey();

0 commit comments

Comments
 (0)