diff -u -r1.31 -r1.32
--- ReclusterDTreeDriver.java	16 Jul 2008 21:04:31 -0000	1.31
+++ ReclusterDTreeDriver.java	17 Jul 2008 01:03:32 -0000	1.32
@@ -34,7 +34,7 @@

   * in this package, which uses the implementation in
   * org.lcsim.recon.cluster.directedtree developed by NIU).
   *

-  * @version $Id: ReclusterDTreeDriver.java,v 1.31 2008/07/16 21:04:31 mcharles Exp $

+  * @version $Id: ReclusterDTreeDriver.java,v 1.32 2008/07/17 01:03:32 mcharles Exp $

   * @author Mat Charles <[log in to unmask]>
   */
 

@@ -91,7 +91,7 @@

     }
     
     public ReclusterDTreeDriver(String dTreeClusterList, String trackList, String mcList) {

-	System.out.println("ReclusterDTreeDriver version 0.2");

+	System.out.println("ReclusterDTreeDriver unstable version");

 	initTrackMatch();
 	initCalibration();
 	initPlots();

@@ -1218,7 +1218,13 @@

 	// Calculate kinematic properties
 	double clusterEnergy = energy(clus, allSharedClusters, m_neutralCalib); // CHECK: Is this OK?
 	double clusterMomentumMagSq = clusterEnergy*clusterEnergy - mass_K0*mass_K0;

-	if (clusterMomentumMagSq<0.0) { clusterMomentumMagSq = 0.0; }

+	boolean treatAsZeroMass = false;
+	if (clusterMomentumMagSq<0.0) { 
+	    // Pathological case of E<M.
+	    // Revert to a zero mass particle here.
+	    clusterMomentumMagSq = 0.0; 
+	    treatAsZeroMass = true;
+	}

 	Hep3Vector pos = new BasicHep3Vector(clus.getPosition());
 	Hep3Vector threeMomentum = VecOp.mult(Math.sqrt(clusterMomentumMagSq), VecOp.unit(pos)); // assume it comes from the IP
 	HepLorentzVector fourMomentum = new BasicHepLorentzVector(clusterEnergy, threeMomentum);

@@ -1228,7 +1234,11 @@

 	// Set the PID and mass
 	part.addParticleID(pid_K0);
 	part.setParticleIdUsed(pid_K0);

-	part.setMass(mass_K0);

+	if (treatAsZeroMass) {
+	    part.setMass(0.0);
+	} else {
+	    part.setMass(mass_K0);
+	}

 	// All done
 	return part;
     }

@@ -2762,6 +2772,52 @@

 	    }
 	}
 

+	// There may be some "neutral hadron cluster cores" that are too small to be
+	// real -- their measured energy is below their mass. They will need special treatment.
+	Set<Cluster> lowEnergyNeutralShowers = new HashSet<Cluster>();
+	for (Cluster clus : neutralClusterCores) {
+	    double clusterEnergy = energy(clus, allSharedClusters, m_neutralCalib);
+	    if (clusterEnergy <= mass_K0) {
+		lowEnergyNeutralShowers.add(clus);
+	    }
+	}
+	// If they are entirely in the ECAL then assume they are photons.
+	for (Cluster clus : lowEnergyNeutralShowers) {
+	    List<CalorimeterHit> hits = clus.getCalorimeterHits();
+	    ListFilter filter = new ListFilter(new HitInECALDecision());
+	    List<CalorimeterHit> hitsInEcal = filter.filterList(hits);
+	    boolean allHitsInEcal = (hitsInEcal.size() == hits.size());
+	    if (allHitsInEcal) {
+		extraClustersToTreatAsPhotons.add(clus);
+	    } else {
+		// Not entirely in ECAL. For now, just add it to whatever non-tiny neutral
+		// cluster is nearest. (Not ideal, but not crazy.)
+		Cluster nearestCluster = null;
+		double nearestDistance = 0;
+		for (Cluster otherClus : neutralClusterCores) {
+		    if (lowEnergyNeutralShowers.contains(otherClus)) {
+			// Don't attach to another teeny shower
+			continue;
+		    } else {
+			double dist = proximity(clus, otherClus);
+			if (nearestCluster==null || dist<nearestDistance) {
+			    nearestCluster = otherClus;
+			    nearestDistance = dist;
+			}
+		    }
+		}
+		if (nearestCluster != null) {
+		    // Migrate to the other cluster
+		    ((BasicCluster)(nearestCluster)).addCluster(clus);
+		    neutralClusterCores.remove(clus);
+		} else {
+		    // Couldn't attach it to anything
+		    double clusterEnergy = energy(clus, allSharedClusters, m_neutralCalib);
+		    System.out.println("WARNING: Unable to attach teeny neutral hadron cluster with E="+clusterEnergy+" to a larger shower. It will have unphysical energy.");
+		}
+	    }
+	}
+

 	// Photons (...)
 	List<Cluster> photonsToUse = new Vector<Cluster>(modifiedPhotonClusters);
 	photonsToUse.addAll(extraClustersToTreatAsPhotons);

@@ -2794,7 +2850,7 @@

 	}
 
 	// Make reconstructed particles from the neutral hadron showers:

-

+	List<Cluster> neutralClusterCoresUsed = new Vector<Cluster>();

 	for (Cluster clus : neutralClusterCores) {
 	    if (extraClustersToTreatAsPhotons.contains(clus)) {
 		// Skip this (already used as a photon)

@@ -2802,6 +2858,7 @@

 	    }
 	    // write out neutral hadron
 	    BaseReconstructedParticle part = makeNeutralHadron(clus, allSharedClusters);

+	    neutralClusterCoresUsed.add(clus);

 	    // Add to the output list
 	    outputParticleList.add(part);
 	    outputNeutralHadronParticleList.add(part);

@@ -2865,8 +2922,8 @@

 	// Start by recording all used cluster cores...
 	List<Cluster> allUsedCoresList = new Vector<Cluster>();
 	Set<Cluster> allUsedCores = new HashSet<Cluster>();

-	allUsedCores    .addAll(neutralClusterCores);
-	allUsedCoresList.addAll(neutralClusterCores);

+	allUsedCores    .addAll(neutralClusterCoresUsed);
+	allUsedCoresList.addAll(neutralClusterCoresUsed);

 	if (allUsedCores.size() != allUsedCoresList.size()) { throw new AssertionError("Mis-counting of clusters: "+allUsedCores.size()+" vs "+allUsedCoresList.size()); }
 	allUsedCores    .addAll(photonsUsed);
 	allUsedCoresList.addAll(photonsUsed);