diff -N FindSamplingFractions.java
--- /dev/null	1 Jan 1970 00:00:00 -0000
+++ FindSamplingFractions.java	9 Aug 2005 21:54:52 -0000	1.1
@@ -0,0 +1,189 @@

+import java.util.*;
+import hep.aida.*;
+import org.lcsim.event.MCParticle;
+import org.lcsim.event.Cluster;
+import org.lcsim.event.EventHeader;
+import org.lcsim.recon.cluster.cheat.*;
+import org.lcsim.recon.cluster.util.*;
+import org.lcsim.util.Driver;
+import org.lcsim.util.aida.AIDA;
+import org.lcsim.event.CalorimeterHit;
+import org.lcsim.event.SimCalorimeterHit;
+import org.lcsim.geometry.IDDecoder;
+import hep.physics.vec.*;
+
+/**
+ * Reconstruction: Clusters
+ */
+public class FindSamplingFractions extends Driver
+{
+   private AIDA aida = AIDA.defaultInstance();
+    IDDecoder _decoder;
+   String EMname = "EcalBarrHitsClusters";
+   String HADname = "HcalBarrHitsClusters";
+   String EMhitsname = "EcalBarrHits";
+   String HADhitsname = "HcalBarrHits";
+   CheatClusterer _clusterer;
+      IAnalysisFactory     af ;
+      ITree                tree ;
+      IDataPointSetFactory dpsf ;
+      IFunctionFactory     funcF ;
+      IFitFactory          fitF ;
+      IFitter              fitter;
+
+       IHistogramFactory hf;
+       ICloud2D c2g ;
+       ICloud2D c2n ;
+       ICloud2D c2nh ;
+      IDataPointSet dataPointSetg;
+      IDataPointSet dataPointSetn;
+      IDataPointSet dataPointSetnh;
+      int npg;
+      int npn;
+      int npnh;
+   
+   public FindSamplingFractions()
+   {
+      af     = IAnalysisFactory.create();
+      tree   = af.createTreeFactory().create();
+      dpsf   = af.createDataPointSetFactory(tree);
+      funcF  = af.createFunctionFactory(tree);
+      fitF   = af.createFitFactory();
+      fitter = fitF.createFitter("Chi2","jminuit");
+
+       hf = af.createHistogramFactory(tree);
+       c2g = hf.createCloud2D("gamma Edep vs Egen");
+       c2n = hf.createCloud2D("neutron Edep vs Egen");
+       c2nh = hf.createCloud2D("neutron #hits vs Egen");
+
+      // Create a two dimensional IDataPointSet.
+      dataPointSetg = dpsf.create("gamma dps","gamma Edep vs Egen",2);
+      dataPointSetn = dpsf.create("neutron dps","neutron Edep vs Egen",2);
+      dataPointSetnh = dpsf.create("neutron dps2","neutron #hits vs Egen",2);
+      npg = 0;
+      npn = 0;
+      npnh = 0;
+      
+        _clusterer = new CheatClusterer();
+   }
+   protected void process(EventHeader event)
+   {
+        List<List<SimCalorimeterHit>> collections = event.get(SimCalorimeterHit.class);
+        for (List<SimCalorimeterHit> collection : collections)
+        {
+            String name = event.getMetaData(collection).getName();
+            if( (name.compareTo(EMhitsname) == 0) |  (name.compareTo(HADhitsname) == 0))
+            {
+               Map<MCParticle,CheatCluster> result = _clusterer.findClusters(collection);
+               if (result.size() > 0) event.put(name+"Clusters",new ArrayList(result.values()));
+            }
+        }
+      List<List<Cluster>> clusterlists = event.get(Cluster.class);
+      for (List<Cluster> clcollection : clusterlists)
+      {
+         String name = event.getMetaData(clcollection).getName();
+         if(name.compareTo(EMname) == 0)
+         {
+            for(Cluster cl:clcollection)
+            {
+               CheatCluster chcl = (CheatCluster) cl;
+               MCParticle p = chcl.getMCParticle();
+               if(p.getPDGID() == 22)
+               {
+                  if(p.getEnergy() > 1.0)
+                  {
+                     Hep3Vector pp = p.getMomentum();
+                     double pperp = Math.sqrt(pp.x()*pp.x() + pp.y()*pp.y());
+                     if(pperp > Math.abs(pp.z()))
+                     {
+                        Hep3Vector vtx = p.getOrigin();
+                        double rsq = vtx.x()*vtx.x() + vtx.y()*vtx.y() + vtx.z()*vtx.z();
+                        if(rsq < 1.)
+                        {
+                           dataPointSetg.addPoint();
+                           IDataPoint dp = dataPointSetg.point(npg);
+                           dp.coordinate(0).setValue(p.getEnergy());
+                           dp.coordinate(1).setValue(chcl.getRawEnergy());
+                           double err = .2*Math.sqrt(p.getEnergy())*.012;
+                           dp.coordinate(1).setErrorPlus(err);
+                           dp.coordinate(1).setErrorMinus(err);
+                        
+                           c2g.fill(p.getEnergy(),chcl.getRawEnergy());
+                           npg++;
+                        }
+                     }
+                  }
+               }
+            }
+         }
+
+         if(name.compareTo(HADname) == 0)
+         {
+            for(Cluster cl:clcollection)
+            {
+               CheatCluster chcl = (CheatCluster) cl;
+               MCParticle p = chcl.getMCParticle();
+               if(p.getPDGID() == 2112)
+               {
+                  if((p.getEnergy()-p.getMass() > 1.0)&(p.getEnergy()-p.getMass() < 10.0))
+                  {
+                     Hep3Vector ep = p.getEndPoint();
+                     double eR = Math.sqrt(ep.x()*ep.x() + ep.y()*ep.y());
+                     if( (eR > 1390.) & (eR < 1810.) )
+                     {
+                        Hep3Vector pp = p.getMomentum();
+                        double pperp = Math.sqrt(pp.x()*pp.x() + pp.y()*pp.y());
+                        if(pperp > Math.abs(pp.z()))
+                        {
+                           Hep3Vector vtx = p.getOrigin();
+                           double rsq = vtx.x()*vtx.x() + vtx.y()*vtx.y() + vtx.z()*vtx.z();
+                           if(rsq < 1.)
+                           {
+                              double cutE = 0.;
+                              int cutN = 0;
+                              List<CalorimeterHit> chl = chcl.getCalorimeterHits();
+                              for(CalorimeterHit hit:chl)
+                              {
+                                 if(hit.getRawEnergy() > .0003)
+                                 {
+                                    cutN++;
+                                    cutE += hit.getRawEnergy();
+                                 }
+                              }
+                              dataPointSetn.addPoint();
+                              IDataPoint dp = dataPointSetn.point(npn);
+                              dp.coordinate(0).setValue(p.getEnergy() - p.getMass());
+                              dp.coordinate(1).setValue(cutE);
+                              double err = (Math.sqrt(p.getEnergy() - p.getMass())*.035)*.50;
+                              dp.coordinate(1).setErrorPlus(err);
+                              dp.coordinate(1).setErrorMinus(err);
+
+                              c2n.fill(p.getEnergy() - p.getMass(),cutE);
+                              npn++;
+/*
+                              List<CalorimeterHit> hits = chcl.getCalorimeterHits();
+                              int nhits = hits.size();
+                              double nhd = (double) nhits;
+*/
+                              double nhd = (double) cutN;
+                              dataPointSetnh.addPoint();
+                              dp = dataPointSetnh.point(npnh);
+                              dp.coordinate(0).setValue(p.getEnergy() - p.getMass());
+                              dp.coordinate(1).setValue(nhd);
+//                              err = (Math.sqrt(p.getEnergy() - p.getMass())*8.3 - 3.4)*.78;
+                              err = (Math.sqrt(p.getEnergy() - p.getMass())*12.)*.50;
+                              dp.coordinate(1).setErrorPlus(err);
+                              dp.coordinate(1).setErrorMinus(err);
+                              c2nh.fill(p.getEnergy() - p.getMass(),nhd);
+                              npnh++;
+                           }
+                        }
+                     }
+                  }
+               }
+            }
+         }
+
+      }
+   }
+}