Print

Print


Commit in lcsim-contrib/src/main/java/org/lcsim/contrib/jeremym/hps on MAIN
HPSEcalClusterer.java+206added 1.1
first implemention of basic HPS Ecal clustering algorithm

lcsim-contrib/src/main/java/org/lcsim/contrib/jeremym/hps
HPSEcalClusterer.java added at 1.1
diff -N HPSEcalClusterer.java
--- /dev/null	1 Jan 1970 00:00:00 -0000
+++ HPSEcalClusterer.java	28 Apr 2011 00:34:30 -0000	1.1
@@ -0,0 +1,206 @@
+package org.lcsim.contrib.jeremym.hps;
+
+import java.util.ArrayList;
+import java.util.HashMap;
+import java.util.List;
+import java.util.Map;
+import java.util.Set;
+
+import org.lcsim.event.CalorimeterHit;
+import org.lcsim.event.Cluster;
+import org.lcsim.event.EventHeader;
+import org.lcsim.geometry.Detector;
+import org.lcsim.geometry.IDDecoder;
+import org.lcsim.geometry.subdetector.HPSEcal;
+import org.lcsim.geometry.subdetector.HPSEcal.XYSide;
+import org.lcsim.geometry.util.IDEncoder;
+import org.lcsim.recon.cluster.util.BasicCluster;
+import org.lcsim.util.Driver;
+import org.lcsim.util.lcio.LCIOConstants;
+
+/**
+ * Creates clusters from CalorimeterHits in the HPSEcal detector.
+ * 
+ * Clustering algorithm is from pages 83 and 84 of the HPS Proposal. 
+ * 
+ * @author Jeremy McCormick <[log in to unmask]>
+ * @version $Id: HPSEcalClusterer.java,v 1.1 2011/04/28 00:34:30 jeremy Exp $
+ */
+public class HPSEcalClusterer extends Driver 
+{
+    HPSEcal ecal;
+    
+    String ecalCollectionName;
+    String ecalName;
+    String clusterCollectionName = "EcalClusters";
+        
+    // Minimum E for cluster seed.
+    double seedEMin = .05;
+    
+    // Minimum E to add hit to cluster.
+    double addEMin = .03;
+    
+    // Odd or even number of crystals in X.
+    boolean oddX;
+    
+    public HPSEcalClusterer()
+    {}
+    
+    public void setClusterCollectionName(String clusterCollectionName)
+    {
+        this.clusterCollectionName = clusterCollectionName;
+    }
+    
+    public void setSeedEMin(double seedEMin)
+    {
+        this.seedEMin = seedEMin;
+    }
+    
+    public void setAddEMin(double addEMin)
+    {
+        this.addEMin = addEMin;
+    }
+    
+    public void setEcalCollectionName(String ecalCollectionName)
+    {
+        this.ecalCollectionName = ecalCollectionName;
+    }
+    
+    public void setEcalName(String ecalName)
+    {
+        this.ecalName = ecalName;
+    }
+    
+    public void startOfData()
+    {
+        if (ecalCollectionName == null)
+            throw new RuntimeException("The parameter ecalCollectionName was not set!");
+        
+        if (ecalCollectionName == null)
+            throw new RuntimeException("The parameter ecalName was not set!");
+        
+        System.out.println(this.getClass().getCanonicalName());
+        System.out.println(" seedEMin="+seedEMin);
+        System.out.println(" addEMin="+addEMin);    
+        System.out.println();
+    }   
+    
+    public void detectorChanged(Detector detector)
+    {
+        ecal = (HPSEcal)detector.getSubdetector(ecalName);
+                        
+        System.out.println(ecal.getName());
+        System.out.println(" nx="+ecal.nx());
+        System.out.println(" ny="+ecal.ny());
+        System.out.println(" beamgap="+ecal.beamGap());
+        System.out.println(" dface="+ecal.distanceToFace());
+    }    
+    
+    public void process(EventHeader event)
+    {
+        //System.out.println("HPSEcalClusterer.process");
+                        
+        List<CalorimeterHit> hits = event.get(CalorimeterHit.class, ecalCollectionName);
+        if (hits == null)
+            throw new RuntimeException("Event is missing ECal hits collection!");
+            
+        IDDecoder dec = ecal.getIDDecoder();
+               
+        // Hit map.
+        Map<Long,CalorimeterHit> hitMap = new HashMap<Long,CalorimeterHit>();
+        
+        // Make map of x, y, and side to hit.
+        for (CalorimeterHit hit : hits)
+        {
+            dec.setID(hit.getCellID());
+            hitMap.put(hit.getCellID(), hit);            
+        }
+        
+        // Setup an encoder from the decoder.
+        IDEncoder enc = new IDEncoder(dec.getIDDescription());
+        
+        // Cluster list to be added to event.
+        List<Cluster> clusters = new ArrayList<Cluster>();
+        
+        // Loop over ECal hits.
+        for (CalorimeterHit hit : hits)
+        {
+            // Cut on min seed E.
+            if (hit.getRawEnergy() < seedEMin)
+                continue;
+            
+            // Set decoder's ID.
+            dec.setID(hit.getCellID());
+            
+            // Get ID values.
+            int ix = dec.getValue("ix");
+            int iy = dec.getValue("iy");
+            int side = dec.getValue("side");
+            
+            // Get neighboring crystals to hit.
+            Set<XYSide> neighbors = ecal.getNeighbors(ix, iy, side);
+            
+            //System.out.println("neighbors of " + ix + ", " + iy + ", " + side);
+           
+            int[] buffer = new int[dec.getFieldCount()];                        
+            dec.getValues(buffer);
+            
+            // Loop over neighbors.
+            List<CalorimeterHit> neighborHits = new ArrayList<CalorimeterHit>();
+            for (XYSide neighbor : neighbors)
+            {             
+                // debug print 
+                //System.out.println(" " + neighbor.x() + ", " + neighbor.y() + ", " + neighbor.side());
+                
+                buffer[dec.getFieldIndex("ix")] = neighbor.x();
+                buffer[dec.getFieldIndex("iy")] = neighbor.y();
+                buffer[dec.getFieldIndex("side")] = neighbor.side();                
+                
+                // Create ID for this crystal.
+                // FIXME This is redundant and inefficient.
+                long neighborId = enc.setValues(buffer);
+              
+                // Find neighbor hit if it exists.
+                CalorimeterHit neighborHit = hitMap.get(neighborId);
+                
+                // Got a hit?
+                if (neighborHit != null)
+                {                    
+                    // debug print
+                    /*
+                    dec.setID(neighborHit.getCellID());
+                    System.out.println("   found neighbor hit w/ x, y, side ( " 
+                            + dec.getValue("ix") + ", " 
+                            + dec.getValue("iy") + ", " 
+                            + dec.getValue("side") + ")");
+                            */
+                         
+                    // Neighbor is hotter so skip this hit.
+                    if (neighborHit.getRawEnergy() > hit.getRawEnergy())
+                    {
+                        //System.out.println("skipping hit because neighbor has more E");
+                        continue;
+                    } 
+                    
+                    // Add to cluster if above min E.
+                    if (neighborHit.getRawEnergy() >= addEMin)                    
+                        neighborHits.add(neighborHit);
+                }                                  
+            }
+            
+            // Make cluster.
+            BasicCluster cluster = new BasicCluster();
+            cluster.addHit(hit);
+            for (CalorimeterHit clusHit : neighborHits)
+            {
+                cluster.addHit(clusHit);
+            }            
+            clusters.add(cluster);
+            //System.out.println("cluster nhits = " + cluster.getCalorimeterHits().size());
+            //System.out.println("cluster E = " + cluster.getEnergy());
+        }
+        
+        int flag = 1 << LCIOConstants.CLBIT_HITS;        
+        event.put(clusterCollectionName, clusters, Cluster.class, flag);        
+    }         
+}
\ No newline at end of file
CVSspam 0.2.8