diff -N ClusterFinding.java
--- /dev/null	1 Jan 1970 00:00:00 -0000
+++ ClusterFinding.java	22 Aug 2011 20:00:07 -0000	1.1
@@ -0,0 +1,216 @@

+package org.lcsim.mcd.analysis;
+
+import java.util.Collections;
+import java.util.Comparator;
+import java.util.HashMap;
+import java.util.List;
+import java.util.Map;
+import org.lcsim.event.CalorimeterHit;
+import org.lcsim.event.Cluster;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.SimCalorimeterHit;
+import org.lcsim.recon.cluster.fixedcone.FixedConeClusterer;
+import org.lcsim.recon.cluster.fixedcone.FixedConeClusterer.FixedConeDistanceMetric;
+import org.lcsim.util.Driver;
+import org.lcsim.util.aida.AIDA;
+
+/**
+ * An example that shows how to find clusters and make some simple plots
+ * of the results.
+ * 
+ * @author Norman Graf
+ * @version $Id: ClusterFinding.java,v 1.1 2011/08/22 20:00:07 ngraf Exp $
+ * 
+ */
+public class ClusterFinding extends Driver
+{
+
+    private FixedConeClusterer _clusterer;
+    private String[] _collectionNames;
+
+    public class ClusterSortBySize implements Comparator<Cluster>
+    {
+
+        @Override
+        public int compare(Cluster cl1, Cluster cl2)
+        {
+            int diff = cl2.getSize() - cl1.getSize();
+            return diff;
+        }
+    }
+
+    public class ClusterSortByEnergy implements Comparator<Cluster>
+    {
+
+        @Override
+        public int compare(Cluster cl1, Cluster cl2)
+        {
+            double diff = cl2.getEnergy() - cl1.getEnergy();
+            if (diff > 0.0)
+            {
+                return 1;
+            } else if (diff < 0.0)
+            {
+                return -1;
+            } else
+            {
+                return 0;
+            }
+        }
+    }
+//
+// Use the "convenient" method of generating AIDA plots 
+//
+    private AIDA aida = AIDA.defaultInstance();
+
+    public ClusterFinding()
+    {
+
+        double radius = .4;
+        double seed = 0.;
+        double minE = 0.;
+        FixedConeDistanceMetric dm = FixedConeDistanceMetric.DOTPRODUCT;
+        _clusterer = new FixedConeClusterer(radius, seed, minE, dm);
+    }
+
+    public void setCollectionNames(String[] collectionNames)
+    {
+        _collectionNames = new String[collectionNames.length];
+        System.arraycopy(collectionNames, 0, _collectionNames, 0, collectionNames.length);
+    }
+//
+// Process an event
+//
+
+    @Override
+    protected void process(EventHeader event)
+    {
+        // which calorimeter hits should we cluster?
+        Map<Long, CalorimeterHit> hitmap = new HashMap<Long, CalorimeterHit>();
+        // which collections should we process?
+        if (_collectionNames != null)
+        {
+            // only process the collections we requested
+            for (int i = 0; i < _collectionNames.length; ++i)
+            {
+                List<SimCalorimeterHit> hits = event.getSimCalorimeterHits(_collectionNames[i]);
+                for (SimCalorimeterHit h : hits)
+                {
+                    // apply any time or energy cuts here...
+                    hitmap.put(h.getCellID(), h);
+                }
+            }
+        } else
+        {
+            // process all calorimetercollections...
+            List<List<SimCalorimeterHit>> hitList = event.get(SimCalorimeterHit.class);
+            for (List<SimCalorimeterHit> l : hitList)
+            {
+                for (SimCalorimeterHit h : l)
+                {
+                    // apply any time or energy cuts here...
+                    hitmap.put(h.getCellID(), h);
+                }
+            }
+        }
+        System.out.println("hitmap contains " + hitmap.size() + " hits");
+//
+// Make clusters
+//
+        List<Cluster> clusters = _clusterer.createClusters(hitmap);
+        System.out.println("found " + clusters.size() + " clusters.");
+        {
+//
+// Get the ClusterList name
+//
+//            String name = event.getMetaData(clusters).getName() + "/";
+            String name = "fixedCone";
+
+            if (clusters.size() > 0)
+                event.put(name + "_clusters", clusters);
+//
+// Histogram the number of clusters in the List
+//
+
+            aida.cloud1D(name + "clusters").fill(clusters.size());
+            Collections.sort(clusters, new ClusterSortByEnergy());
+
+            // Loop over all the clusters in a List
+            //
+            int count = 0;
+            for (Cluster cluster : clusters)
+            {
+                count++;
+                double[] pos1 =
+                {
+                    0.0, 0.0, 0.0
+                };
+                double[] pos2 =
+                {
+                    0.0, 0.0, 0.0
+                };
+                double[] pos3 =
+                {
+                    0.0, 0.0, 0.0
+                };
+                double e1 = 0.0;
+                double e2 = 0.0;
+//                System.out.println("Cluster size:  " + cluster.getSize() + "  Energy:  " + cluster.getEnergy());
+                if (count == 1)
+                {
+                    e1 = cluster.getEnergy();
+                    aida.cloud1D(name + "energy first cluster").fill(cluster.getEnergy());
+                    aida.cloud1D(name + "size first cluster").fill(cluster.getSize());
+                    pos1 = cluster.getPosition();
+                }
+                if (count == 2)
+                {
+                    e2 = cluster.getEnergy();
+                    aida.cloud1D(name + "energy second cluster").fill(cluster.getEnergy());
+                    aida.cloud1D(name + "size second cluster").fill(cluster.getSize());
+                    pos2 = cluster.getPosition();
+                    pos3[0] = pos2[0] - pos1[0];
+                    pos3[1] = pos2[1] - pos1[1];
+                    pos3[2] = pos2[2] - pos1[2];
+                    double delta = Math.sqrt(pos3[0] * pos3[0] + pos3[1] * pos3[1] + pos3[2] * pos3[2]);
+                    aida.cloud1D(name + "distance").fill(delta);
+                    aida.cloud1D(name + "deltaE").fill(e1 - e2);
+                    aida.cloud2D(name + "disvsdeltaE").fill(delta, e1 - e2);
+                }
+                if (count == 3)
+                {
+                    aida.cloud1D(name + "energy third cluster").fill(cluster.getEnergy());
+                }
+                if (count > 3)
+                {
+                    aida.cloud1D(name + "energy fourth and more cluster").fill(cluster.getEnergy());
+                }
+
+//
+// Histogram the "corrected" energy
+//
+                aida.cloud1D(name + "energy").fill(cluster.getEnergy());
+                aida.cloud1D(name + "size").fill(cluster.getSize());
+                aida.cloud2D(name + "sizevsenergy").fill(cluster.getSize(), cluster.getEnergy());
+//
+// Histogram the position as R vs Z
+//
+                double[] pos = cluster.getPosition();
+                double R = Math.sqrt(pos[0] * pos[0] + pos[1] * pos[1]);
+                aida.cloud2D(name + "Position:R vs Z").fill(pos[2], R);
+//
+// Histogram the computed direction
+//
+                aida.cloud1D(name + "Direction: theta").fill(cluster.getITheta());
+                aida.cloud1D(name + "Direction: phi").fill(cluster.getIPhi());
+//
+// Histogram the difference in direction and position theta,phi 
+//
+                double posphi = Math.atan2(pos[1], pos[0]);
+                aida.cloud1D(name + "delta phi").fill(posphi - cluster.getIPhi());
+                double postheta = Math.PI / 2. - Math.atan2(pos[2], R);
+                aida.cloud1D(name + "delta theta").fill(postheta - cluster.getITheta());
+            }
+        }
+    }
+}