Commit in lcsim-cal-calib/src/org/lcsim/cal/calib on MAIN | |||
PandoraPfoSamplingFractionAnalysislDriver.java | +21 | -5 | 1.2 -> 1.3 |
SamplingFractionAnalysisPolyCalDriver.java | +23 | -7 | 1.3 -> 1.4 |
+44 | -12 |
Fixed logic to identify the initial particle to also work on samples created by GEANT particle gun Corrected expected energy calculation
diff -u -r1.2 -r1.3 --- PandoraPfoSamplingFractionAnalysislDriver.java 4 Nov 2010 00:20:12 -0000 1.2 +++ PandoraPfoSamplingFractionAnalysislDriver.java 10 Feb 2012 15:26:26 -0000 1.3 @@ -1,7 +1,7 @@
/* * PandoraPfoSamplingFractionAnalysislDriver *
- * $Id: PandoraPfoSamplingFractionAnalysislDriver.java,v 1.2 2010/11/04 00:20:12 ngraf Exp $
+ * $Id: PandoraPfoSamplingFractionAnalysislDriver.java,v 1.3 2012/02/10 15:26:26 grefe Exp $
*/ package org.lcsim.cal.calib;
@@ -99,7 +99,19 @@
// organize the histogram tree by species and energy List<MCParticle> mcparts = event.getMCParticles();
- MCParticle mcpart = mcparts.get(mcparts.size() - 1);
+ //MCParticle mcpart = mcparts.get(mcparts.size() - 1); // this only works if particle was generated by stdhep, does not work with GEANT gps + MCParticle mcpart = null; + // Look for the most energetic final state particle + for (MCParticle myMCP : mcparts) { + if (myMCP.getGeneratorStatus() == MCParticle.FINAL_STATE) { + if (mcpart == null) { + mcpart = myMCP; + } else if (mcpart.getEnergy() < myMCP.getEnergy()) { + mcpart = myMCP; + } + } + } +
String particleType = mcpart.getType().getName(); double mcEnergy = mcpart.getEnergy(); long mcIntegerEnergy = Math.round(mcEnergy);
@@ -163,7 +175,11 @@
double clusterEnergy = c.getEnergy(); double mcMass = mcpart.getType().getMass(); // subtract the mass to get kinetic energy...
- double expectedEnergy = sqrt(mcEnergy * mcEnergy - mcMass * mcMass);
+ double expectedEnergy = mcEnergy; + // In case of protons and neutrons we expect only the kinetic energy in the shower + if (mcpart.getPDGID() == 2212 || mcpart.getPDGID() == 2112) { + expectedEnergy = mcEnergy - mcMass; + }
// System.out.println(mcpart.getType().getName()+" "+expectedEnergy); aida.cloud1D("measured - predicted energy").fill(clusterEnergy - expectedEnergy);
@@ -189,7 +205,7 @@
// System.out.println(detectorName); int caltype = 0; // FIXME Hard-coded name.
- if (detectorName.startsWith("Ecal"))
+ if (detectorName.toUpperCase().startsWith("ECAL"))
{ for (int i = 1; i < _ecalLayering.length + 1; ++i) {
@@ -198,7 +214,7 @@
} } // FIXME Hard-coded name.
- if (detectorName.startsWith("Hcal"))
+ if (detectorName.toUpperCase().startsWith("HCAL"))
{ caltype = 3; nHCalHits +=1;
diff -u -r1.3 -r1.4 --- SamplingFractionAnalysisPolyCalDriver.java 3 Nov 2010 15:19:51 -0000 1.3 +++ SamplingFractionAnalysisPolyCalDriver.java 10 Feb 2012 15:26:26 -0000 1.4 @@ -3,7 +3,7 @@
* * Created on May 19, 2008, 11:54 AM *
- * $Id: SamplingFractionAnalysisPolyCalDriver.java,v 1.3 2010/11/03 15:19:51 ngraf Exp $
+ * $Id: SamplingFractionAnalysisPolyCalDriver.java,v 1.4 2012/02/10 15:26:26 grefe Exp $
*/ package org.lcsim.cal.calib;
@@ -150,7 +150,19 @@
// organize the histogram tree by species and energy List<MCParticle> mcparts = event.getMCParticles();
- MCParticle mcpart = mcparts.get(mcparts.size() - 1);
+ //MCParticle mcpart = mcparts.get(mcparts.size() - 1); // this only works if particle was generated by stdhep, does not work with GEANT gps + MCParticle mcpart = null; + // Look for the most energetic final state particle + for (MCParticle myMCP : mcparts) { + if (myMCP.getGeneratorStatus() == MCParticle.FINAL_STATE) { + if (mcpart == null) { + mcpart = myMCP; + } else if (mcpart.getEnergy() < myMCP.getEnergy()) { + mcpart = myMCP; + } + } + } +
String particleType = mcpart.getType().getName(); double mcEnergy = mcpart.getEnergy(); long mcIntegerEnergy = Math.round(mcEnergy);
@@ -244,7 +256,11 @@
double clusterEnergy = c.getEnergy(); double mcMass = mcpart.getType().getMass(); // subtract the mass to get kinetic energy...
- double expectedEnergy = sqrt(mcEnergy * mcEnergy - mcMass * mcMass);
+ double expectedEnergy = mcEnergy; + // In case of protons and neutrons we expect only the kinetic energy in the shower + if (mcpart.getPDGID() == 2212 || mcpart.getPDGID() == 2112) { + expectedEnergy = mcEnergy - mcMass; + }
// System.out.println(mcpart.getType().getName()+" "+expectedEnergy); aida.cloud1D("measured - predicted energy").fill(clusterEnergy - expectedEnergy);
@@ -265,7 +281,7 @@
String detectorName = decoder.getSubdetector().getName(); int type = 0; // FIXME Hard-coded name.
- if (detectorName.startsWith("Ecal"))
+ if (detectorName.toUpperCase().startsWith("ECAL"))
{ if (layer >= firstEmStartLayer && layer < secondEmStartLayer) {
@@ -276,7 +292,7 @@
} } // FIXME Hard-coded name.
- if (detectorName.startsWith("Hcal"))
+ if (detectorName.toUpperCase().startsWith("HCAL"))
{ type = 2; }
@@ -287,10 +303,10 @@
// expectedEnergy = a*E1 + b*E2 +c*E3 for (int j = 0; j < 3; ++j) {
- _vec[j] += expectedEnergy * vals[j];
+ _vec[j] += vals[j];
for (int k = 0; k < 3; ++k) {
- _acc[j][k] += vals[j] * vals[k];
+ _acc[j][k] += (vals[j] * vals[k]) / expectedEnergy ;
} } } // end of single cluster cut
Use REPLY-ALL to reply to list
To unsubscribe from the LCD-CVS list, click the following link:
https://listserv.slac.stanford.edu/cgi-bin/wa?SUBED1=LCD-CVS&A=1