diff -u -r1.4 -r1.5
--- MiscUtilities.java	18 Oct 2005 22:32:57 -0000	1.4
+++ MiscUtilities.java	16 Dec 2005 21:09:28 -0000	1.5
@@ -16,9 +16,16 @@

 import org.lcsim.geometry.CylindricalSubdetector;
 //import org.lcsim.geometry.Subdetector;
 import org.lcsim.util.swim.Line;

+import org.lcsim.event.MCParticle;
+import org.lcsim.event.SimCalorimeterHit;
+import org.lcsim.event.EventHeader;

-class MiscUtilities

+import util.HitCountAssociator;
+
+public class MiscUtilities

 {

+    //MiscUtilities() {}
+

     /**
      * Get the minimum distance between hits in the two clusters.
      * If one or both of the clusters is empty, will return NaN.

@@ -45,26 +52,34 @@

 
     static protected double distance(Cluster clus, CalorimeterHit hit)
     {

-	// Loop over hits...

+	Hep3Vector point = new BasicHep3Vector(hit.getPosition());
+	return distance(clus, point);
+    }
+
+    static protected double distance(CalorimeterHit hit1, CalorimeterHit hit2)
+    {
+	Hep3Vector point2 = new BasicHep3Vector(hit2.getPosition());
+	return distance(hit1, point2);
+    }
+
+    static protected double distance(Cluster clus, Hep3Vector point)
+    {

 	boolean firstCheck = true;
 	double minDistance = Double.NaN; // Will stay NaN if clus is empty
 	List<CalorimeterHit> hits = clus.getCalorimeterHits();
 	for (CalorimeterHit hitInCluster : hits) {

-	    double dist = distance(hit, hitInCluster);

+	    double dist = distance(hitInCluster, point);

 	    if (firstCheck || dist<minDistance) {
 		minDistance = dist;
 		firstCheck = false;
 	    }
 	}

-

 	return minDistance;
     }
 

-    static protected double distance(CalorimeterHit hit1, CalorimeterHit hit2)
-    {
-	Hep3Vector vect1 = new BasicHep3Vector(hit1.getPosition());
-	Hep3Vector vect2 = new BasicHep3Vector(hit2.getPosition());
-	Hep3Vector displacement = VecOp.sub(vect1, vect2);

+    static protected double distance(CalorimeterHit hit, Hep3Vector point) {
+	Hep3Vector hitPoint = new BasicHep3Vector(hit.getPosition());
+	Hep3Vector displacement = VecOp.sub(hitPoint, point);

 	return displacement.magnitude();
     }
 

@@ -188,4 +203,213 @@

 	}
 	return matchingLayer;
     }

+
+    // Blah... belongs in Line.
+    static protected double findDOCAToPoint(Hep3Vector linePoint, Hep3Vector lineDirection, Hep3Vector point) 
+    {
+        // The first line is kind of ugly.
+        Hep3Vector displacement = new BasicHep3Vector(point.x() - linePoint.x(), point.y() - linePoint.y(), point.z() - linePoint.z());
+        double dotProduct = VecOp.dot(displacement, lineDirection);
+        double doca = Math.sqrt(displacement.magnitudeSquared() - dotProduct*dotProduct);
+        return doca;
+    }
+
+    /*
+      static protected boolean isHadronic(Cluster clus, EventHeader event)
+      {
+      throw new AssertionError("Not implemented, apparently");
+      //HitCountAssociator assoc = new HitCountAssociator(event);
+      //List<MCParticle> assocParticles = assoc.associateClusterToMCParticles(clus);
+      //MCParticle dominant = assocParticles.iterator().next();
+      //return ( Math.abs(dominant.getCharge()) < 0.001);
+      }
+
+      static protected double getCorrectedEnergy_sidaug05(Cluster clus, EventHeader event)
+      {
+      // Need to figure out the raw EM energy, the number of EM and HAD hits
+      double EMEraw = 0.0;
+      int nEMhits = 0;
+      int nHADhits = 0;
+      int nRejectedHitsHAD = 0;
+      for (CalorimeterHit hit : clus.getCalorimeterHits()) {
+      org.lcsim.geometry.Subdetector subdet = hit.getSubdetector();
+      if ( ! subdet.isCalorimeter() ) { throw new AssertionError("Cluster hit outside calorimeter"); }
+      String name = subdet.getName();
+      if (name.compareTo("EMBarrel") == 0 || name.compareTo("EMEndcap") == 0) {
+      // ECAL
+      nEMhits++;
+      EMEraw += hit.getRawEnergy();
+      } else if (name.compareTo("HADBarrel") == 0 || name.compareTo("HADEndcap") == 0) {
+      // HCAL
+      if (hit.getTime() < 100) {
+      // Within first 100 ns => OK
+      nHADhits++;
+      } else {
+      nRejectedHitsHAD++;
+      }
+      } else {
+      throw new AssertionError("DEBUG: Found non-calorimeterhit in calorimeter '"+subdet+"' with name '"+subdet.getName()+"'");
+      }
+      }
+
+      // Also need the "sangle", i.e. sin(theta)
+      double sangle = Double.NaN;
+      Hep3Vector[] posAndDir = MiscUtilities.getPositionAndDirection(clus);
+      Hep3Vector pos = null;
+      if (posAndDir != null) {
+      pos = posAndDir[0];
+      } else {
+      // Too few hits -- calculate energy-weighted position by hand.
+      double eRawSum = 0.0;
+      double[] energyWeightedPositionSum = new double[3];
+      for (int i=0; i<3; i++) { energyWeightedPositionSum[i] = 0.0; }
+      for (CalorimeterHit hit : clus.getCalorimeterHits()) {
+      double eRaw = hit.getRawEnergy();
+      eRawSum += eRaw;
+      double[] hitPos = hit.getPosition();
+      for (int i=0; i<3; i++) { energyWeightedPositionSum[i] += hitPos[i] * eRaw; }
+      }
+      for (int i=0; i<3; i++) { energyWeightedPositionSum[i] /= eRawSum ; }
+      pos = new BasicHep3Vector(energyWeightedPositionSum);
+      }
+      double mag = pos.magnitude();
+      double r = Math.sqrt(pos.x()*pos.x() + pos.y()*pos.y());
+      sangle = r/mag;
+      
+      // Now, is this an E/M shower or a hadronic shower?
+      boolean isHadronic = MiscUtilities.isHadronic(clus, event);
+      
+      if (isHadronic) {
+      // Hadronic -- use Ron's routine
+      double correctedEnergy = MiscUtilities.getCorrectedEnergy_sidaug05(EMEraw, nEMhits, nHADhits, sangle);
+      return correctedEnergy;
+      } else {
+      // Electromagnetic -- use default calibration
+      double defaultCalibrationEnergySum = 0.0;
+      for (CalorimeterHit hit : clus.getCalorimeterHits()) {
+      defaultCalibrationEnergySum += hit.getCorrectedEnergy();
+      }
+      return defaultCalibrationEnergySum;
+      }
+      }
+      
+      static protected double getCorrectedEnergy_sidaug05(double EMEraw, int nEMhits, int nHADhits, double sangle)
+      {
+      final double HADsf = 8.81;
+      final double HADoffset = 4.;
+      final double EMlin = .033897;
+      final double EMquad = -.00000781;
+      final double EMelin = 81.418;
+      final double EMequad = -8.65;
+      final int EMncut = 200;
+      final double AngleCoef = 1.;
+      final double offset = 0.;
+      
+      double Emeas = 0.;
+      double Emhad = 0.;
+      double Emem = 0.;
+      if(nEMhits < EMncut) {
+      Emem = nEMhits*EMlin + nEMhits*nEMhits*EMquad;
+      } else {
+      Emem = EMEraw*EMelin + EMEraw*EMEraw*EMequad;
+      }
+      double angcorr = 1. + AngleCoef*(Math.sqrt(sangle) -1.);
+      Emhad = (nHADhits+HADoffset)/HADsf/angcorr;
+      Emeas = Emem + Emhad;
+      //System.out.println("DEBUG: Emem="+Emem+", Emhad="+Emhad);
+      return (Emeas + offset);
+      }
+    */
+    
+    static protected double correctedEnergyInClusterECAL(MCParticle part, Cluster clus)
+    {
+	double energyECAL = 0.0;
+	for (CalorimeterHit hit : clus.getCalorimeterHits()) {
+	    boolean useThisHit = true;
+	    if (part != null) {
+		useThisHit = hitMatch(part, (SimCalorimeterHit)(hit));
+	    }
+	    if (useThisHit) {
+		org.lcsim.geometry.Subdetector subdet = hit.getSubdetector();
+		if ( ! subdet.isCalorimeter() ) { throw new AssertionError("Cluster hit outside calorimeter"); }
+		String name = subdet.getName();
+		if (name.compareTo("EMBarrel") == 0 || name.compareTo("EMEndcap") == 0) {
+		    // ECAL
+		    SimCalorimeterHit simHit = (SimCalorimeterHit) hit;
+		    double rawEnergy = hitMatchEnergy(part, simHit);
+		    double rawEnergyTotal = hit.getRawEnergy();
+		    double frac = (rawEnergy / rawEnergyTotal);
+		    energyECAL += frac * hit.getCorrectedEnergy();
+		}
+	    }
+	}
+	return energyECAL;
+    }
+
+    static protected int countHitsInClusterHCAL(MCParticle part, Cluster clus)
+    {
+	int countHitsHCAL = 0;
+	for (CalorimeterHit hit : clus.getCalorimeterHits()) {
+	    boolean useThisHit = true;
+	    if (part != null) {
+		useThisHit = hitMatch(part, (SimCalorimeterHit)(hit));
+	    }
+	    if (useThisHit) {
+		org.lcsim.geometry.Subdetector subdet = hit.getSubdetector();
+		if ( ! subdet.isCalorimeter() ) { throw new AssertionError("Cluster hit outside calorimeter"); }
+		String name = subdet.getName();
+		if (name.compareTo("HADBarrel") == 0 || name.compareTo("HADEndcap") == 0) {
+		    // HCAL
+		    if (hit.getTime() < 100) {
+			// Within first 100 ns => OK
+			countHitsHCAL++;
+		    }
+		}
+	    }
+	}
+	return countHitsHCAL;
+    }
+
+    static protected boolean hitMatch(MCParticle part, SimCalorimeterHit hit) {
+	int nContributingParticles = hit.getMCParticleCount();
+	for (int i=0; i<nContributingParticles; i++) {
+	    MCParticle hitPart = hit.getMCParticle(i);
+	    if (part == hitPart) {
+		return true;
+	    }
+	}
+	return false;
+    }
+
+    static protected double hitMatchEnergy(MCParticle part, SimCalorimeterHit hit) {
+	int nContributingParticles = hit.getMCParticleCount();
+	for (int i=0; i<nContributingParticles; i++) {
+	    MCParticle hitPart = hit.getMCParticle(i);
+	    if (part == hitPart) {
+		return hit.getContributedEnergy(i);
+	    }
+	}
+	return 0.0;
+    }
+
+    static protected int countHitsInCluster(MCParticle part, Cluster clus)
+    {
+	int countHits = 0;
+	for (CalorimeterHit hit : clus.getCalorimeterHits()) {
+	    if (hit instanceof SimCalorimeterHit) {
+		SimCalorimeterHit simHit = (SimCalorimeterHit) (hit);
+		int numContributingParticles = simHit.getMCParticleCount();
+		for (int i=0; i<numContributingParticles; i++) {
+		    MCParticle contributingParticle = simHit.getMCParticle(i);
+		    if (part == contributingParticle) {
+			// Match
+			countHits++;
+		    }
+		}
+	    } else {
+		throw new AssertionError("Non-simulated calorimeter hit");
+	    }
+	}
+	return countHits;
+    }

 }