Print

Print


Commit in hps-java/src/main/java/org/lcsim/hps/users/phansson on MAIN
TwoTrackAnlysis.java+191-711.2 -> 1.3
More branches to root tree

hps-java/src/main/java/org/lcsim/hps/users/phansson
TwoTrackAnlysis.java 1.2 -> 1.3
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();
+        } 
+    }
     
     
     
CVSspam 0.2.12


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