java/branches/ecal-recon_HPSJAVA-93/src/main/java/org/hps/recon/ecal
--- java/branches/ecal-recon_HPSJAVA-93/src/main/java/org/hps/recon/ecal/EcalClusterIC.java 2014-06-19 20:56:41 UTC (rev 727)
+++ java/branches/ecal-recon_HPSJAVA-93/src/main/java/org/hps/recon/ecal/EcalClusterIC.java 2014-06-19 20:59:30 UTC (rev 728)
@@ -19,6 +19,7 @@
import org.lcsim.lcio.LCIOConstants;
import org.lcsim.util.Driver;
+
/**
* This Driver creates clusters from the CalorimeterHits of an
* {@link org.lcsim.geometry.subdetectur.HPSEcal3} detector.
@@ -38,15 +39,17 @@
*/
public class EcalClusterIC extends Driver {
// File writer to output cluster results.
- FileWriter writeHits;
+ FileWriter writeHits = null;
// LCIO collection name for calorimeter hits.
String ecalCollectionName;
// Name of the calorimeter detector object.
String ecalName = "Ecal";
// LCIO cluster collection name to which to write.
String clusterCollectionName = "EcalClusters";
+ // Collection name for rejected hits
+ String rejectedHitName = "RejectedHits";
// File path to which to write event display output.
- String outfile = "cluster-hit-IC.txt";
+ String outfile = null;
// Map of crystals to their neighbors.
NeighborMap neighborMap = null;
// Minimum energy threshold for hits; lower energy hits will be
@@ -54,10 +57,10 @@
double hitEnergyThreshold = 0.0075;
// Minimum energy threshold for seed hits; if seed hit is below
// cluster is excluded from output. Units in GeV.
- double seedEnergyThreshold = 0.2;
+ double seedEnergyThreshold = 0.1;
// Minimum energy threshold for cluster hits; if total cluster
// energy is below, the cluster is excluded. Units in GeV.
- double clusterEnergyThreshold = 0.4;
+ double clusterEnergyThreshold = 0.3;
// A Comparator that sorts CalorimeterHit objects from highest to
// lowest energy.
private static final EnergyComparator ENERGY_COMP = new EnergyComparator();
@@ -85,7 +88,19 @@
public void setEcalName(String ecalName) {
this.ecalName = ecalName;
}
+
+ /**
+ * Output file name for event display output.
+ * @param outfile
+ */
+ public void setOutfile(String outfile) {
+ this.outfile = outfile;
+ }
+ public void setRejectedHitName(String rejectedHitName){
+ this.rejectedHitName = rejectedHitName;
+ }
+
/**
* Minimum energy for a hit to be used in a cluster. Default of 0.0075 GeV
*
@@ -96,7 +111,7 @@
}
/**
- * Minimum energy for a seed hit. Default of 0.2 GeV
+ * Minimum energy for a seed hit. Default of 0.1 GeV
*
* @param seedEnergyThreshold
*/
@@ -105,7 +120,7 @@
}
/**
- * Minimum energy for a cluster. Default of 0.4 GeV
+ * Minimum energy for a cluster. Default of 0.3 GeV
*
* @param clusterEnergyThreshold
*/
@@ -152,12 +167,14 @@
throw new RuntimeException("The parameter ecalName was not set!");
}
- // Create a file writer and clear the output file, if it exists.
- try {
- writeHits = new FileWriter(outfile);
- writeHits.write("");
+ if (outfile!=null) {
+ // Create a file writer and clear the output file, if it exists.
+ try {
+ writeHits = new FileWriter(outfile);
+ writeHits.write("");
+ } catch (IOException e) {
+ }
}
- catch(IOException e) { }
}
public void detectorChanged(Detector detector) {
@@ -176,25 +193,29 @@
// Make sure the current event contains calorimeter hits.
if (event.hasCollection(CalorimeterHit.class, ecalCollectionName)) {
// Get the list of raw calorimeter hits.
- List<CalorimeterHit> hitList = event.get(CalorimeterHit.class, ecalCollectionName);
// Generate clusters from the calorimeter hits.
- List<HPSEcalCluster> clusterList = null;
- try { clusterList = createClusters(hitList); }
+ //List<HPSEcalClusterIC> clusterList = null;
+ try { createClusters(event); }
catch(IOException e) { }
-
- // If clusters were successfully created, put them in the event.
- if(clusterList != null) {
- int flag = 1 << LCIOConstants.CLBIT_HITS;
- event.put(clusterCollectionName, clusterList, HPSEcalCluster.class, flag);
- }
}
}
- public List<HPSEcalCluster> createClusters(List<CalorimeterHit> hitList) throws IOException {
+ public void createClusters(EventHeader event) throws IOException {
+
+ // Create a list to store the event hits in.
+ List<CalorimeterHit> hitList = new ArrayList<CalorimeterHit>();
+ List<CalorimeterHit> baseList = event.get(CalorimeterHit.class, ecalCollectionName);
+ for(CalorimeterHit r : baseList) {
+ hitList.add(r);
+ }
+
// Create a list to store the newly created clusters in.
- ArrayList<HPSEcalCluster> clusterList = new ArrayList<HPSEcalCluster>();
+ ArrayList<HPSEcalClusterIC> clusterList = new ArrayList<HPSEcalClusterIC>();
+ // Create a list to store the rejected hits in.
+ ArrayList<CalorimeterHit> rejectedHitList = new ArrayList<CalorimeterHit>();
+
// Sort the list of hits by energy.
Collections.sort(hitList, ENERGY_COMP);
@@ -204,7 +225,8 @@
for(int index = hitList.size() - 1; index >= 0; index--) {
// If the hit is below threshold or outside of time window, kill it.
if((hitList.get(index).getCorrectedEnergy() < hitEnergyThreshold)||
- (timeCut && (hitList.get(index).getTime() < minTime || hitList.get(index).getTime() > minTime + timeWindow))) {
+ (timeCut && (hitList.get(index).getTime() < minTime || hitList.get(index).getTime() > (minTime + timeWindow)))) {
+ rejectedHitList.add(hitList.get(index));
hitList.remove(index);
}
@@ -244,7 +266,7 @@
for (Long neighbor : neighbors) {
// Get the neighboring hit.
CalorimeterHit neighborHit = hitMap.get(neighbor);
-
+
// If it exists, add it to the list.
if(neighborHit != null) { neighborHits.add(neighborHit); }
}
@@ -257,10 +279,11 @@
// neighboring hits.
seedHitLoop:
for(CalorimeterHit neighbor : neighborHits) {
- if(hit.getCorrectedEnergy() <= neighbor.getCorrectedEnergy()) {
+ if(!equalEnergies(hit, neighbor)) {
isSeed = false;
- break seedHitLoop;
+ break seedHitLoop;
}
+
}
// If this hit is a seed hit, just map it to itself.
@@ -289,6 +312,7 @@
// common seeds.
commonHitList.add(neighborHit);
commonHitList.add(hitSeedMap.get(hit));
+
// Put the common seed list back into the set.
commonHits.put(hit, commonHitList);
@@ -299,7 +323,7 @@
// associate it with the neighboring seed and note
// that it has been clustered.
else {
- hitSeedMap.put(hit, neighborHit);
+ hitSeedMap.put(hit, neighborHit);
surrSeedSet.add(hit);
}
}
@@ -338,12 +362,27 @@
// If the neighboring hit is of lower energy than the
// current secondary hit, then associate the neighboring
// hit with the current secondary hit's seed.
- if (secondaryNeighborHit.getCorrectedEnergy() < secondaryHit.getCorrectedEnergy()) {
+
+ // if (secondaryNeighborHit.getCorrectedEnergy() < secondaryHit.getCorrectedEnergy()) {
+ if(!equalEnergies(secondaryNeighborHit, secondaryHit)) {
hitSeedMap.put(secondaryNeighborHit, hitSeedMap.get(secondaryHit));
}
+ else {continue;}
}
} // End component hits loop.
+
+ // This is a check to ensure ALL hits are either components or seeds.
+ for (CalorimeterHit check : hitList){
+ if(!hitSeedMap.containsKey(check)){
+ System.out.println("Something is not clustered or component!");
+ System.out.println("not clustered:"+"\t"+check.getIdentifierFieldValue("ix")+"\t"+
+ check.getIdentifierFieldValue("iy")+"\t"+check.getCorrectedEnergy());
+ }
+ }
+
+
+
// Performs second pass calculations for common hits.
commonHitsLoop:
for (CalorimeterHit clusteredHit : hitSeedMap.keySet()) {
@@ -376,9 +415,12 @@
// Check to make sure that the clustered neighbor hit
// is not already associated with the current clustered
// hit's seed.
- if (hitSeedMap.get(clusteredNeighborHit) != clusteredHitSeed) {
- if (clusteredHit.getCorrectedEnergy() < clusteredNeighborHit.getCorrectedEnergy()) {
- // Check and see if a list of common seeds
+
+ if (hitSeedMap.get(clusteredNeighborHit) != clusteredHitSeed){
+
+ //if (clusteredHit.getCorrectedEnergy() < clusteredNeighborHit.getCorrectedEnergy()) {
+ if(!equalEnergies(clusteredHit, clusteredNeighborHit)){
+ // Check and see if a list of common seeds
// for this hit already exists or not.
List<CalorimeterHit> commonHitList = commonHits.get(clusteredHit);
@@ -388,14 +430,17 @@
// Add the neighbors to the seeds to set of
// common seeds.
commonHitList.add(clusteredHitSeed);
- commonHitList.add(hitSeedMap.get(clusteredNeighborHit));
+ commonHitList.add(hitSeedMap.get(clusteredNeighborHit));
// Put the common seed list back into the set.
commonHits.put(clusteredHit, commonHitList);
}
}
+
+
}
} // End common hits loop.
+
// Remove any common hits from the clustered hits collection.
for(CalorimeterHit commonHit : commonHits.keySet()) {
@@ -403,7 +448,6 @@
}
-
/*
* All hits are sorted from above. The next part of the code is for calculating energies. Still
* needs implementation into new cluster collection so as to preserve shared hit energy
@@ -427,33 +471,32 @@
eEnergy += entry.getKey().getRawEnergy();
seedEnergy.put(eSeed, eEnergy);
}
-
+
//Distribute common hit energies with clusters
Map<CalorimeterHit, Double> seedEnergyTot = seedEnergy;
+
for (Map.Entry<CalorimeterHit, List<CalorimeterHit>> entry1 : commonHits.entrySet()) {
- CalorimeterHit commonCell = entry1.getKey();
- List<CalorimeterHit> commSeedList = entry1.getValue();
- CalorimeterHit seedA = commSeedList.get(0);
- CalorimeterHit seedB = commSeedList.get(1);
- double eFractionA = seedEnergy.get(seedA)/(seedEnergy.get(seedA)+seedEnergy.get(seedB));
+ CalorimeterHit commonCell = entry1.getKey();
+ CalorimeterHit seedA = entry1.getValue().get(0);
+ CalorimeterHit seedB = entry1.getValue().get(1);
+ double eFractionA = seedEnergy.get(seedA)/(seedEnergy.get(seedA)+seedEnergy.get(seedB));
double eFractionB = seedEnergy.get(seedB)/(seedEnergy.get(seedA)+seedEnergy.get(seedB));
- double currEnergyA = seedEnergyTot.get(seedA);
- double currEnergyB = seedEnergyTot.get(seedB);
- currEnergyA += eFractionA * commonCell.getCorrectedEnergy();
- currEnergyB += eFractionB * commonCell.getCorrectedEnergy();
+ double currEnergyA = seedEnergyTot.get(seedA);
+ double currEnergyB = seedEnergyTot.get(seedB);
+ currEnergyA += eFractionA * commonCell.getCorrectedEnergy();
+ currEnergyB += eFractionB * commonCell.getCorrectedEnergy();
- seedEnergyTot.put(seedA, currEnergyA);
- seedEnergyTot.put(seedB, currEnergyB);
+ seedEnergyTot.put(seedA, currEnergyA);
+ seedEnergyTot.put(seedB, currEnergyB);
}
+
+
-
-
-
/*
* Prints the results in event display format. Not analyzed
* for efficiency, as this will ultimately not be a part of
* the driver and should be handled by the event display output
- * driver instead.
+ * driver instead. Contains output loops to collection.
*/
// Only write to the output file is something actually exists.
if (hitMap.size() != 0) {
@@ -473,58 +516,86 @@
for (Map.Entry<CalorimeterHit, CalorimeterHit> entry2 : hitSeedMap.entrySet()) {
- if ((entry2.getKey() == entry2.getValue())&&(entry2.getKey().getCorrectedEnergy()>=seedEnergyThreshold)
- &&(seedEnergyTot.get(entry2.getKey())>=clusterEnergyThreshold)) {//seed and passes all thresholds
- int six = entry2.getKey().getIdentifierFieldValue("ix");
- int siy = entry2.getKey().getIdentifierFieldValue("iy");
- double energy = seedEnergyTot.get(entry2.getKey());
-// writeHits.append(String.format("Cluster\t%d\t%d\t%f%n", six, siy, energy));
+ if (entry2.getKey() == entry2.getValue()){
+ if((entry2.getKey().getCorrectedEnergy()<seedEnergyThreshold)
+ ||(seedEnergyTot.get(entry2.getKey())<clusterEnergyThreshold))
+ {
+ rejectedHitList.add(entry2.getKey());
+ }
- cluster.setParameters(entry2.getKey());
- cluster.addHit(entry2.getKey());
+ else{
+
+ int six = entry2.getKey().getIdentifierFieldValue("ix");
+ int siy = entry2.getKey().getIdentifierFieldValue("iy");
+ double energy = seedEnergyTot.get(entry2.getKey());
+// writeHits.append(String.format("Cluster\t%d\t%d\t%f%n", six, siy, energy));
+
+ HPSEcalClusterIC cluster = new HPSEcalClusterIC(entry2.getKey());
+ cluster.addHit(entry2.getKey());
- for (Map.Entry<CalorimeterHit, CalorimeterHit> entry3 : hitSeedMap.entrySet()) {
- if (entry3.getValue() == entry2.getValue()) {
- int ix = entry3.getKey().getIdentifierFieldValue("ix");
- int iy = entry3.getKey().getIdentifierFieldValue("iy");
-// writeHits.append(String.format("CompHit\t%d\t%d%n", ix, iy));
+ for (Map.Entry<CalorimeterHit, CalorimeterHit> entry3 : hitSeedMap.entrySet()) {
+ if (entry3.getValue() == entry2.getValue()) {
+ if(rejectedHitList.contains(entry2.getValue())){
+ rejectedHitList.add(entry3.getKey());
+ }
+ else{
+ int ix = entry3.getKey().getIdentifierFieldValue("ix");
+ int iy = entry3.getKey().getIdentifierFieldValue("iy");
+// writeHits.append(String.format("CompHit\t%d\t%d%n", ix, iy));
- cluster.addHit(entry3.getKey());
- }
- }
+ cluster.addHit(entry3.getKey());
+ }
+ }
+ }
+
for (Map.Entry<CalorimeterHit, List<CalorimeterHit>> entry4 : commonHits.entrySet()) {
if (entry4.getValue().contains(entry2.getKey())) {
int ix = entry4.getKey().getIdentifierFieldValue("ix");
int iy = entry4.getKey().getIdentifierFieldValue("iy");
// writeHits.append(String.format("SharHit\t%d\t%d%n", ix, iy));
- cluster.addHit(entry4.getKey());
+ // Added in shared hits for energy distribution between clusters, changed by HS 02JUN14
+// cluster.addHit(entry4.getKey());
+ cluster.addSharedHit(entry4.getKey());
}
}
-
- clusterList.add(cluster);
+ for(CalorimeterHit q : rejectedHitList)
+ {// This does not output in correct event display format, just for de-bugging
+// writeHits.append("Rejected"+q.getIdentifierFieldValue("ix")+"\t"+q.getIdentifierFieldValue("iy")+"\n");
+ }
+
+
+ cluster.setEnergy(energy);
+ clusterList.add(cluster);
+ }// End checking thresholds and write out.
}
- }
+
- // Write the event termination header.
+ } //End cluster loop
+ // Write the event termination header.
// writeHits.append("EndEvent\n");
- } //End event display output loop.
+// System.out.println("Number of clusters: "+clusterList.size());
+
+
+ } //End event display out loop.
+ int flag = 1 << LCIOConstants.CLBIT_HITS;
+ event.put(clusterCollectionName, clusterList, HPSEcalClusterIC.class, flag);
+ event.put(rejectedHitName, rejectedHitList, CalorimeterHit.class, flag);
-
-
-
- // Return the resulting cluster list.
- return clusterList;
}
public void endOfData() {
- // Close the event display output writer.
- try { writeHits.close(); }
- catch (IOException e) { }
+ if (writeHits != null) {
+ // Close the event display output writer.
+ try {
+ writeHits.close();
+ } catch (IOException e) {
+ }
+ }
}
- private static class EnergyComparator implements Comparator<CalorimeterHit> {
+ /* private static class EnergyComparator implements Comparator<CalorimeterHit> {
public int compare(CalorimeterHit o1, CalorimeterHit o2) {
// If the energies are equivalent, the same, the two hits
// are considered equivalent.
@@ -536,5 +607,58 @@
// Lower energy hits are ranked lower than higher energy hits.
else { return 1; }
}
+ }*/
+
+ // Also accounts for pathological case of cluster hits that are EXACTLY the same.
+ private static class EnergyComparator implements Comparator<CalorimeterHit> {
+ public int compare(CalorimeterHit o1, CalorimeterHit o2) {
+ // If the energies are equivalent, the same, the two hits
+ // are considered equivalent.
+ if(o1.getCorrectedEnergy() == o2.getCorrectedEnergy()) {
+ if(Math.abs(o1.getIdentifierFieldValue("iy")) < Math.abs(o2.getIdentifierFieldValue("iy"))){
+ return -1;
+ }
+ else if((Math.abs(o1.getIdentifierFieldValue("iy")) == Math.abs(o2.getIdentifierFieldValue("iy")))
+ && (o1.getIdentifierFieldValue("ix") < o2.getIdentifierFieldValue("ix"))){
+ return -1; }
+ else if (Math.abs(o1.getIdentifierFieldValue("iy")) > Math.abs(o2.getIdentifierFieldValue("iy"))){
+ return 1;
+ }
+ else{return 1;}
+ }
+ // Higher energy hits are ranked higher than lower energy hits.
+ else if(o1.getCorrectedEnergy() > o2.getCorrectedEnergy()) { return -1; }
+
+ // Lower energy hits are ranked lower than higher energy hits.
+ else { return 1; }
+ }
}
-}
+
+
+
+
+ // Handles pathological case where multiple neighboring crystals have EXACTLY the same energy.
+ private boolean equalEnergies(CalorimeterHit hit, CalorimeterHit neighbor){
+ boolean isSeed = true;
+
+ int hix = hit.getIdentifierFieldValue("ix");
+ int hiy = Math.abs(hit.getIdentifierFieldValue("iy"));
+ int nix = neighbor.getIdentifierFieldValue("ix");
+ int niy = Math.abs(neighbor.getIdentifierFieldValue("iy"));
+ double hE = hit.getCorrectedEnergy();
+ double nE = neighbor.getCorrectedEnergy();
+ if(hE < nE) {
+ isSeed = false;
+ }
+ else if((hE == nE) && (hiy > niy)) {
+ isSeed = false;
+ }
+ else if((hE == nE) && (hiy == niy) && (hix > nix)) {
+ isSeed = false;
+ }
+ return isSeed;
+ }
+
+
+
+}
\ No newline at end of file