Commit in hps-java/src/main/java/org/lcsim/hps/users/phansson on MAIN | |||
TwoTrackAnlysis.java | +191 | -71 | 1.2 -> 1.3 |
More branches to root tree
diff -u -r1.2 -r1.3 --- TwoTrackAnlysis.java 8 Jan 2013 00:37:47 -0000 1.2 +++ TwoTrackAnlysis.java 10 Jan 2013 18:35:53 -0000 1.3 @@ -9,6 +9,7 @@
import java.io.FileWriter; import java.io.IOException; import java.io.PrintWriter;
+import java.util.ArrayList;
import java.util.HashMap; import java.util.List; import java.util.Map;
@@ -19,9 +20,12 @@
import org.lcsim.event.Track; import org.lcsim.event.TrackerHit; import org.lcsim.event.util.ParticleTypeClassifier;
+import org.lcsim.fit.helicaltrack.HelicalTrackCross;
import org.lcsim.fit.helicaltrack.HelicalTrackHit;
+import org.lcsim.fit.helicaltrack.HelicalTrackStrip;
import org.lcsim.geometry.Detector; import org.lcsim.hps.analysis.ecal.HPSMCParticlePlotsDriver;
+import org.lcsim.hps.evio.TriggerData;
import org.lcsim.hps.recon.ecal.HPSEcalCluster; import org.lcsim.hps.recon.tracking.SvtTrackExtrapolator; import org.lcsim.hps.recon.vertexing.SimpleVertexer;
@@ -48,6 +52,7 @@
private boolean hideFrame = false; private String outputPlotFileName; private String ecalClusterCollectionName = "EcalClusters";
+ private String stereoHitCollectionName = "RotatedHelicalTrackHits";
private boolean _debug; private TwoTrackVertexer vertexer = new TwoTrackVertexer(); private TwoParticleVertexer particleVertexer = new TwoParticleVertexer();
@@ -109,6 +114,7 @@
@Override public void process(EventHeader event) {
+ if(this._debug) System.out.println(this.getClass().getSimpleName() + ": processing event "+totalEvents);
totalEvents++;
@@ -121,7 +127,7 @@
}
- if(tracklist.size()!=2) {
+ if(tracklist.size()<2) {
return; }
@@ -226,6 +232,7 @@
this.fillTextTuple(electron, positron, trk1, trk2, vtxPosMC, vtxPos, clusters, event);
+ if(this._debug) System.out.println(this.getClass().getSimpleName() + ": # two track events so far = "+totalTwoTrackEvents);
}
@@ -238,7 +245,7 @@
System.out.println(this.getClass().getSimpleName() + ": Total Number of Events = "+this.totalEvents); System.out.println(this.getClass().getSimpleName() + ": Total Number of Two-Track events Processed = "+this.totalTwoTrackEvents); System.out.println(this.getClass().getSimpleName() + ": Total Number of MCEvents = "+this.totalMCEvents);
- System.out.println(this.getClass().getSimpleName() + ": Total Number of Two-Track events Processed = "+this.totalTwoTrackMCEvents);
+ System.out.println(this.getClass().getSimpleName() + ": Total Number of Two-Track MC events Processed = "+this.totalTwoTrackMCEvents);
if (!"".equals(outputPlotFileName)) { try {
@@ -261,69 +268,7 @@
}
- private void makePlots() { - _vtxpos_x = aida.histogram1D("Vertex position X", 100, -1000, 1000); - _vtxpos_y = aida.histogram1D("Vertex position Y", 100, -200, 200); - _vtxpos_z = aida.histogram1D("Vertex position Z", 100, -200, 200); - _vtxpos_xy = aida.histogram2D("Vertex position Y vs X", 50, -200, 200, 50, -200, 200); - _partvtxpos_x = aida.histogram1D("Particle Vertex position X", 100, -1000, 1000); - _partvtxpos_y = aida.histogram1D("Particle Vertex position Y", 100, -200, 200); - _partvtxpos_z = aida.histogram1D("Particle Vertex position Z", 100, -200, 200); - _partvtxpos_xy = aida.histogram2D("Particle Vertex position Y vs X", 50, -200, 200, 50, -200, 200); - _nhits = aida.histogram1D("# hits", 3, 4, 7); - _nhits_vtxpospos = aida.histogram1D("# hits for pos vtx", 3, 4, 7); - _nhits_vtxposneg = aida.histogram1D("# hits for neg vtx", 3, 4, 7); - _layershit = aida.histogram1D("Layers hit", 5, 1, 6); - _layershit_vtxpospos = aida.histogram1D("Layers hit pos vtx", 5, 1, 6); - _layershit_vtxposneg = aida.histogram1D("Layers hit neg vtx", 5, 1, 6); - - - _plotterTrackVertex = aida.analysisFactory().createPlotterFactory().create(); - _plotterTrackVertex.createRegions(2,2); - _plotterTrackVertex.region(0).plot(_vtxpos_x); - _plotterTrackVertex.region(1).plot(_vtxpos_y); - _plotterTrackVertex.region(2).plot(_vtxpos_z); - _plotterTrackVertex.region(3).plot(_vtxpos_xy); - _plotterTrackVertex.region(3).style().setParameter("hist2DStyle", "colorMap"); - _plotterTrackVertex.region(3).style().dataStyle().fillStyle().setParameter("colorMapScheme", "rainbow"); - _plotterParticleVertex = aida.analysisFactory().createPlotterFactory().create(); - _plotterParticleVertex.createRegions(2,2); - _plotterParticleVertex.region(0).plot(_partvtxpos_x); - _plotterParticleVertex.region(1).plot(_partvtxpos_y); - _plotterParticleVertex.region(2).plot(_partvtxpos_z); - _plotterParticleVertex.region(3).plot(_partvtxpos_xy); - _plotterParticleVertex.region(3).style().setParameter("hist2DStyle", "colorMap"); - _plotterParticleVertex.region(3).style().dataStyle().fillStyle().setParameter("colorMapScheme", "rainbow"); - _plotterNhits = aida.analysisFactory().createPlotterFactory().create(); - _plotterNhits.createRegions(1,2); - _plotterNhits.region(0).plot(_nhits); - _plotterNhits.region(0).plot(this._nhits_vtxposneg,"mode=overlay"); - _plotterNhits.region(0).plot(this._nhits_vtxpospos,"mode=overlay"); - _plotterNhits.region(0).style().dataStyle().fillStyle().setVisible(false); - _plotterNhits.region(1).plot(this._layershit); - _plotterNhits.region(1).plot(this._layershit_vtxposneg,"mode=overlay"); - _plotterNhits.region(1).plot(this._layershit_vtxpospos,"mode=overlay"); - _plotterNhits.region(1).style().dataStyle().fillStyle().setVisible(false); - - _plotterParticleVertex.setTitle("MC particle Vertex"); - _plotterTrackVertex.setTitle("Two Track Vertex"); - _plotterNhits.setTitle("Hits on track"); - - for(int i=0;i<3;++i) { - ((PlotterRegion) _plotterParticleVertex.region(i)).getPlot().setAllowUserInteraction(true); - ((PlotterRegion) _plotterParticleVertex.region(i)).getPlot().setAllowPopupMenus(true); - ((PlotterRegion) _plotterTrackVertex.region(i)).getPlot().setAllowUserInteraction(true); - ((PlotterRegion) _plotterTrackVertex.region(i)).getPlot().setAllowPopupMenus(true); - - } - - if(!this.hideFrame) { - this._plotterParticleVertex.show(); - this._plotterTrackVertex.show(); - this._plotterNhits.show(); - } - } -
+
private boolean isFileEmpty(String fileName) { File f = new File(fileName);
@@ -337,13 +282,21 @@
br_line+="e_px/F:e_py/F:e_pz/F:p_px/F:p_py/F:p_pz/F:"; br_line+="trk1_q/I:trk1_chi2/F:trk1_px/F:trk1_py/F:trk1_pz/F:trk1_nhits/I:"; for(int iLayer=1;iLayer<=5;++iLayer) br_line+="trk1_hit"+iLayer+"_x/F:"+"trk1_hit"+iLayer+"_y/F:"+"trk1_hit"+iLayer+"_z/F:";
+ for(int iLayer=1;iLayer<=10;++iLayer) br_line+="trk1_strip"+iLayer+"_u/F:"; + for(int iLayer=1;iLayer<=10;++iLayer) br_line+="trk1_strip"+iLayer+"_n/F:";
br_line+="trk2_q/I:trk2_chi2/F:trk2_px/F:trk2_py/F:trk2_pz/F:trk2_nhits/I:"; for(int iLayer=1;iLayer<=5;++iLayer) br_line+="trk2_hit"+iLayer+"_x/F:"+"trk2_hit"+iLayer+"_y/F:"+"trk2_hit"+iLayer+"_z/F:";
+ for(int iLayer=1;iLayer<=10;++iLayer) br_line+="trk2_strip"+iLayer+"_u/F:"; + for(int iLayer=1;iLayer<=10;++iLayer) br_line+="trk2_strip"+iLayer+"_n/F:";
br_line+="vtx_truth_x/F:vtx_truth_y/F:vtx_truth_z/F:"; br_line+="vtx_x/F:vtx_y/F:vtx_z/F:"; br_line+="trk1_conv_x/F:trk1_conv_y/F:trk1_conv_z/F:"; br_line+="trk2_conv_x/F:trk2_conv_y/F:trk2_conv_z/F:";
- br_line+="ncl_top/I:ncl_bot/I";
+ br_line+="cl1_E/F:cl1_ix/I:cl1_iy/I:"; + br_line+="cl2_E/F:cl2_ix/I:cl2_iy/I:"; + br_line+="cl3_E/F:cl3_ix/I:cl3_iy/I:"; + br_line+="ncl_top/I:ncl_bot/I:"; + br_line+="trig_top/I:trig_bot/I";
printWriter.println(br_line); doPrintBranchInfoLine = false; }
@@ -360,7 +313,7 @@
for(TrackerHit hit : hitsOnTrack1) { HelicalTrackHit hth = (HelicalTrackHit) hit; hits1.put(hth.Layer(), hth);
- System.out.println("hit on layer " + hth.Layer());
+ //System.out.println("hit on layer " + hth.Layer());
} // stupid but I want to keep one line per event so default in case there is not hits in all layers for(int iLayer=0;iLayer<5;++iLayer) {
@@ -368,14 +321,44 @@
if (hitOnLayer != null) printWriter.format("%5.5f %5.5f %5.5f ", hitOnLayer.getPosition()[0],hitOnLayer.getPosition()[1],hitOnLayer.getPosition()[2]); else printWriter.format("%5.5f %5.5f %5.5f ", -9999999.9, -9999999.9, -9999999.9); }
+ List<HelicalTrackHit> stereoHits = new ArrayList<HelicalTrackHit>(); + if(event.hasCollection(HelicalTrackHit.class, stereoHitCollectionName)) { + stereoHits = event.get(HelicalTrackHit.class, stereoHitCollectionName); + } + HashMap<Integer,List<HelicalTrackStrip>> striphits1 = this.getStripHitsMap(hitsOnTrack1); + for(int iLayer=1;iLayer<=10;++iLayer) { + HelicalTrackStrip strip=null; + if(striphits1.containsKey(iLayer)) strip = striphits1.get(iLayer).get(0); + if(strip!=null) printWriter.format("%5.5f ", strip.umeas()); + else printWriter.format("%5.5f ", -99999999.9); + } + HashMap<Integer,List<HelicalTrackStrip>> allstriphits1 = this.getAllStripHitsMap(stereoHits); + for(int iLayer=1;iLayer<=10;++iLayer) { + if(striphits1.containsKey(iLayer)) printWriter.format("%5d ", allstriphits1.get(iLayer).size()); + else printWriter.format("%5d ", -99999999); + } +
List<TrackerHit> hitsOnTrack2 = trk2.getTrackerHits(); printWriter.format("%5d %5.5f %5.5f %5.5f %5.5f %5d ",trk2.getCharge(),trk2.getChi2(), trk2.getTrackStates().get(0).getMomentum()[0],trk2.getTrackStates().get(0).getMomentum()[1],trk2.getTrackStates().get(0).getMomentum()[2],hitsOnTrack2.size());
+ HashMap<Integer,List<HelicalTrackStrip>> striphits2 = this.getStripHitsMap(hitsOnTrack2);
HashMap<Integer,HelicalTrackHit> hits2 = this.getHitMap(hitsOnTrack2); for(int iLayer=0;iLayer<5;++iLayer) { HelicalTrackHit hitOnLayer = hits2.get(iLayer*2+1);// = this.getHitOnLayer(iLayer, hitsOnTrack); if (hitOnLayer != null) printWriter.format("%5.5f %5.5f %5.5f ", hitOnLayer.getPosition()[0],hitOnLayer.getPosition()[1],hitOnLayer.getPosition()[2]); else printWriter.format("%5.5f %5.5f %5.5f ", -9999999.9, -9999999.9, -9999999.9); }
+ for(int iLayer=1;iLayer<=10;++iLayer) { + HelicalTrackStrip strip=null; + if(striphits2.containsKey(iLayer)) strip = striphits2.get(iLayer).get(0); + if(strip!=null) printWriter.format("%5.5f ", strip.umeas()); + else printWriter.format("%5.5f ", -99999999.9); + } + HashMap<Integer,List<HelicalTrackStrip>> allstriphits2 = this.getAllStripHitsMap(stereoHits); + for(int iLayer=1;iLayer<=10;++iLayer) { + if(striphits2.containsKey(iLayer)) printWriter.format("%5d ", striphits2.get(iLayer).size()); + else printWriter.format("%5d ", -99999999); + } +
//Particle vtx if(vtxPosParticle!=null) printWriter.format("%5.5f %5.5f %5.5f ", vtxPosParticle.x(),vtxPosParticle.y(),vtxPosParticle.z() ); else printWriter.format("%5.5f %5.5f %5.5f ", -9999999., -9999999., -9999999. );
@@ -389,15 +372,47 @@
posAtConverter = this.vertexer.extrapolator().extrapolateTrack(SvtTrackExtrapolator.HARP_POSITION); printWriter.format("%5.5f %5.5f %5.5f ", posAtConverter.z(),posAtConverter.x(),posAtConverter.y()); //note rotation from JLab->tracking int ncl_t=0; int ncl_b=0;
- for(HPSEcalCluster cl : clusters) { - if(cl.getSeedHit().getIdentifierFieldValue("iy") > 0) ncl_t++; - else ncl_b++;
+ for(int i=0;i<3;++i) { + if(clusters.size()<=i) { + printWriter.format("%5.5f %5d %5d ",-999999.9,-999999,-999999); + } else { + //for(HPSEcalCluster cl : clusters) { + int iy = clusters.get(i).getSeedHit().getIdentifierFieldValue("iy"); + int ix = clusters.get(i).getSeedHit().getIdentifierFieldValue("ix"); + double E = clusters.get(i).getEnergy(); + printWriter.format("%5.5f %5d %5d ",E,ix,iy); + if( iy > 0) ncl_t++; + else ncl_b++; + }
}
- printWriter.format("%5d %5d",ncl_t,ncl_b);
+ printWriter.format("%5d %5d ",ncl_t,ncl_b);
+ printWriter.format("%5d %5d ",ncl_t,ncl_b);
+
+ TriggerData trigger = getTriggerInfo(event);
+ if(trigger==null) printWriter.format("%5d %5d",0,0);
+ else printWriter.format("%5d %5d",trigger.getTopTrig(),trigger.getBotTrig());
printWriter.println(); }
+ + private TriggerData getTriggerInfo(EventHeader event) { + String triggerDecisionCollectionName = "TriggerBank"; + List<TriggerData> t; + if(!event.hasCollection(TriggerData.class, triggerDecisionCollectionName)) { + if(_debug) System.out.println( "Event has NO trigger bank"); + return null; + } else { + List<TriggerData> triggerDataList = event.get(TriggerData.class, "TriggerBank"); + if(triggerDataList.isEmpty()) { + if(_debug) System.out.println( "Event has trigger bank exists but is empty"); + return null; + } else { + return triggerDataList.get(0); + } + } + } +
private HashMap<Integer,HelicalTrackHit> getHitMap(List<TrackerHit> hits) { HashMap<Integer,HelicalTrackHit> map = new HashMap<Integer,HelicalTrackHit>();
@@ -408,6 +423,49 @@
return map; }
+ private HashMap<Integer,HelicalTrackStrip> getStripHitMap(List<TrackerHit> hits) { + HashMap<Integer,HelicalTrackStrip> map = new HashMap<Integer,HelicalTrackStrip>(); + for(TrackerHit hit : hits) { + HelicalTrackHit hth = (HelicalTrackHit) hit; + HelicalTrackCross htc = (HelicalTrackCross) hth; + HelicalTrackStrip s1 = htc.getStrips().get(0); + HelicalTrackStrip s2 = htc.getStrips().get(1); + map.put(s1.layer(), s1); + map.put(s2.layer(), s2); + } + return map; + } + + private HashMap<Integer,List<HelicalTrackStrip>> getStripHitsMap(List<TrackerHit> hits) { + HashMap<Integer,List<HelicalTrackStrip>> map = new HashMap<Integer,List<HelicalTrackStrip>>(); + for(TrackerHit hit : hits) { + HelicalTrackHit hth = (HelicalTrackHit) hit; + HelicalTrackCross htc = (HelicalTrackCross) hth; + HelicalTrackStrip s1 = htc.getStrips().get(0); + HelicalTrackStrip s2 = htc.getStrips().get(1); + if(!map.containsKey(s1.layer())) map.put(s1.layer(), new ArrayList<HelicalTrackStrip>()); + if(!map.containsKey(s2.layer())) map.put(s2.layer(), new ArrayList<HelicalTrackStrip>()); + map.get(s1.layer()).add(s1); + map.get(s2.layer()).add(s2); + } + return map; + } + + private HashMap<Integer,List<HelicalTrackStrip>> getAllStripHitsMap(List<HelicalTrackHit> stereoHits) { + HashMap<Integer,List<HelicalTrackStrip>> map = new HashMap<Integer,List<HelicalTrackStrip>>(); + for(HelicalTrackHit hth : stereoHits) { + HelicalTrackCross htc = (HelicalTrackCross) hth; + HelicalTrackStrip s1 = htc.getStrips().get(0); + HelicalTrackStrip s2 = htc.getStrips().get(1); + if(!map.containsKey(s1.layer())) map.put(s1.layer(), new ArrayList<HelicalTrackStrip>()); + if(!map.containsKey(s2.layer())) map.put(s2.layer(), new ArrayList<HelicalTrackStrip>()); + map.get(s1.layer()).add(s1); + map.get(s2.layer()).add(s2); + + } + return map; + } +
private HelicalTrackHit getHitOnLayer(int layer, List<TrackerHit> hits) { HelicalTrackHit hitOnLayer = null; for(TrackerHit hit : hits) {
@@ -419,6 +477,68 @@
return hitOnLayer; }
+ private void makePlots() { + _vtxpos_x = aida.histogram1D("Vertex position X", 100, -1000, 1000); + _vtxpos_y = aida.histogram1D("Vertex position Y", 100, -200, 200); + _vtxpos_z = aida.histogram1D("Vertex position Z", 100, -200, 200); + _vtxpos_xy = aida.histogram2D("Vertex position Y vs X", 50, -200, 200, 50, -200, 200); + _partvtxpos_x = aida.histogram1D("Particle Vertex position X", 100, -1000, 1000); + _partvtxpos_y = aida.histogram1D("Particle Vertex position Y", 100, -200, 200); + _partvtxpos_z = aida.histogram1D("Particle Vertex position Z", 100, -200, 200); + _partvtxpos_xy = aida.histogram2D("Particle Vertex position Y vs X", 50, -200, 200, 50, -200, 200); + _nhits = aida.histogram1D("# hits", 3, 4, 7); + _nhits_vtxpospos = aida.histogram1D("# hits for pos vtx", 3, 4, 7); + _nhits_vtxposneg = aida.histogram1D("# hits for neg vtx", 3, 4, 7); + _layershit = aida.histogram1D("Layers hit", 5, 1, 6); + _layershit_vtxpospos = aida.histogram1D("Layers hit pos vtx", 5, 1, 6); + _layershit_vtxposneg = aida.histogram1D("Layers hit neg vtx", 5, 1, 6); + + + _plotterTrackVertex = aida.analysisFactory().createPlotterFactory().create(); + _plotterTrackVertex.createRegions(2,2); + _plotterTrackVertex.region(0).plot(_vtxpos_x); + _plotterTrackVertex.region(1).plot(_vtxpos_y); + _plotterTrackVertex.region(2).plot(_vtxpos_z); + _plotterTrackVertex.region(3).plot(_vtxpos_xy); + _plotterTrackVertex.region(3).style().setParameter("hist2DStyle", "colorMap"); + _plotterTrackVertex.region(3).style().dataStyle().fillStyle().setParameter("colorMapScheme", "rainbow"); + _plotterParticleVertex = aida.analysisFactory().createPlotterFactory().create(); + _plotterParticleVertex.createRegions(2,2); + _plotterParticleVertex.region(0).plot(_partvtxpos_x); + _plotterParticleVertex.region(1).plot(_partvtxpos_y); + _plotterParticleVertex.region(2).plot(_partvtxpos_z); + _plotterParticleVertex.region(3).plot(_partvtxpos_xy); + _plotterParticleVertex.region(3).style().setParameter("hist2DStyle", "colorMap"); + _plotterParticleVertex.region(3).style().dataStyle().fillStyle().setParameter("colorMapScheme", "rainbow"); + _plotterNhits = aida.analysisFactory().createPlotterFactory().create(); + _plotterNhits.createRegions(1,2); + _plotterNhits.region(0).plot(_nhits); + _plotterNhits.region(0).plot(this._nhits_vtxposneg,"mode=overlay"); + _plotterNhits.region(0).plot(this._nhits_vtxpospos,"mode=overlay"); + _plotterNhits.region(0).style().dataStyle().fillStyle().setVisible(false); + _plotterNhits.region(1).plot(this._layershit); + _plotterNhits.region(1).plot(this._layershit_vtxposneg,"mode=overlay"); + _plotterNhits.region(1).plot(this._layershit_vtxpospos,"mode=overlay"); + _plotterNhits.region(1).style().dataStyle().fillStyle().setVisible(false); + + _plotterParticleVertex.setTitle("MC particle Vertex"); + _plotterTrackVertex.setTitle("Two Track Vertex"); + _plotterNhits.setTitle("Hits on track"); + + for(int i=0;i<3;++i) { + ((PlotterRegion) _plotterParticleVertex.region(i)).getPlot().setAllowUserInteraction(true); + ((PlotterRegion) _plotterParticleVertex.region(i)).getPlot().setAllowPopupMenus(true); + ((PlotterRegion) _plotterTrackVertex.region(i)).getPlot().setAllowUserInteraction(true); + ((PlotterRegion) _plotterTrackVertex.region(i)).getPlot().setAllowPopupMenus(true); + + } + + if(!this.hideFrame) { + this._plotterParticleVertex.show(); + this._plotterTrackVertex.show(); + this._plotterNhits.show(); + } + }
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