SlicDiagnostics/src/org/lcsim/slic/diagnostics
diff -u -r1.12 -r1.13
--- MCParticlePlots.java 1 Sep 2005 17:44:07 -0000 1.12
+++ MCParticlePlots.java 7 Sep 2005 05:47:03 -0000 1.13
@@ -5,20 +5,18 @@
*/
package org.lcsim.slic.diagnostics;
-
-import java.util.List;
import java.util.ArrayList;
import java.util.Map;
import java.util.HashMap;
import static java.lang.Math.sqrt;
import hep.aida.ICloud1D;
-import hep.aida.ICloud2D;
import hep.aida.IHistogram1D;
import hep.aida.IHistogram2D;
import hep.physics.vec.VecOp;
import hep.physics.particle.properties.ParticlePropertyManager;
import hep.physics.particle.properties.ParticlePropertyProvider;
import hep.physics.vec.Hep3Vector;
+import java.util.List;
import org.lcsim.event.MCParticle;
import org.lcsim.event.EventHeader;
@@ -28,7 +26,9 @@
*/
class MCParticlePlots extends AbstractPlots
{
- static ParticlePropertyProvider mcpLkp
+ Map<Integer, MCParticleInfo> mcpMap = new HashMap<Integer, MCParticleInfo>();
+
+ static ParticlePropertyProvider mcpLkp
= ParticlePropertyManager.getParticlePropertyProvider();
/* double values for sim status plot */
@@ -45,6 +45,7 @@
ICloud1D mcpTotalEnergyFinalState;
ICloud1D mcpCosTheta;
ICloud1D mcpPhi;
+ ICloud1D mcpBackscatterNotPEP;
IHistogram1D mcpSimStatus;
IHistogram1D mcpGenStatus;
IHistogram2D mcpSimVsGenStatus;
@@ -85,6 +86,7 @@
mcpVtxToEPDistGTMeV = aida.histogram1D(makePlotName("Vertex to Endpoint Distance for kE > 1 MeV"), 200, 0., 24000.0);
mcpKE = aida.histogram1D(makePlotName("Kinetic Energy"), 200, 0.0, 10.0);
mcpNoEndPoint = aida.cloud1D(makePlotName("No End Point"));
+ mcpBackscatterNotPEP = aida.cloud1D(makePlotName("Backscattered Particles without vertexIsNotEndPointOfParent flag"));
}
protected void fill(EventHeader event)
@@ -93,12 +95,24 @@
fill(event.get(MCParticle.class).get(0));
}
+ private MCParticleInfo getMCParticleInfo(int pdgid)
+ {
+ MCParticleInfo mcpInfo = mcpMap.get(pdgid);
+ if ( mcpInfo == null )
+ {
+ mcpMap.put(pdgid, new MCParticleInfo(pdgid));
+ mcpInfo = mcpMap.get(pdgid);
+ }
+ return mcpInfo;
+ }
+
private void fill(List<MCParticle> particles)
{
if ( particles.size() > 0 )
{
particleCounts.clear();
int FSCount = 0;
+ int backscatterNotPEPCount = 0;
int noEndPointCount = 0;
totE = 0;
for ( MCParticle particle : particles)
@@ -106,6 +120,15 @@
/* generator status */
mcpGenStatus.fill(particle.getGeneratorStatus());
+ /* backscatter w/o vertexIsNotEndpointOfParent flag */
+ if ( particle.getSimulatorStatus().isBackscatter() )
+ {
+ if ( !particle.getSimulatorStatus().vertexIsNotEndpointOfParent() )
+ {
+ backscatterNotPEPCount += 1;
+ }
+ }
+
/* make a list of doubles for sim statuses */
List<Double> statusList = makeSimulatorStatusList(particle.getSimulatorStatus());
@@ -166,12 +189,17 @@
/* fill plots for each final state particle type */
if ( particle.getGeneratorStatus() == MCParticle.FINAL_STATE )
{
+ MCParticleInfo info = getMCParticleInfo(particle.getPDGID());
+ info.addEventEnergy(particle.getEnergy());
+ info.incrNumHits();
+
FSCount += 1;
totE += particle.getEnergy();
mcpEnergyFinalState.fill(particle.getEnergy());
- String particleName = mcpLkp.get(particle.getPDGID()).getName();
+// String particleName = mcpLkp.get(particle.getPDGID()).getName();
+ String particleName = info.getName();
try
{
@@ -214,12 +242,14 @@
aida.cloud2D(makePlotName(particleName + ": Endpoint R vs Z")).fill(
HitUtil.getCylindricalRadius(particle.getEndPoint().v()),
particle.getEndPoint().z());
-
+
/* E vs E in EP daughters */
double dauE = 0.;
for ( MCParticle dau : particle.getDaughters() )
{
- if ( !dau.getSimulatorStatus().vertexIsNotEndpointOfParent() ) {
+ if ( !dau.getSimulatorStatus().vertexIsNotEndpointOfParent() &&
+ !dau.getSimulatorStatus().isBackscatter() )
+ {
dauE += dau.getEnergy();
}
}
@@ -228,30 +258,40 @@
aida.cloud2D(makePlotName(particleName + ": Final State Particle Energy vs Total Energy of Endpoint Daughters")).fill(
particle.getEnergy(), dauE);
}
- }
-
- /* increment particle count */
- if ( particleCounts.containsKey(particleName) )
- {
- int cnt = particleCounts.get(particleName);
- cnt += 1;
- particleCounts.remove(particleName);
- particleCounts.put(particleName, cnt);
- }
- else
- {
- particleCounts.put(particleName, 1);
}
}
tree.cd("..");
}
- }
+ }
- /* count per event for individual particle types */
- for ( String k : particleCounts.keySet() )
+ for ( MCParticleInfo info : mcpMap.values() )
{
+ String k = info.getName();
tree.cd(k);
- aida.cloud1D(makePlotName(k + ": Particle Count by Event")).fill(particleCounts.get(k));
+
+ int numHits = info.getEventNumHits();
+ double energy = info.getEventEnergy();
+
+ /* percentage of FS energy */
+ double percFSE = energy / totE;
+ if ( percFSE != 0 )
+ {
+ aida.cloud1D(makePlotName(k + ": Percentage of Total Final State Energy in Event")).fill(percFSE);
+ }
+
+ /* percentage of FS particles */
+ double percFSParticles = ((double)numHits) / ((double)FSCount);
+
+ if ( percFSParticles != 0 )
+ {
+ aida.cloud1D(makePlotName(k + ": Percentage of Final State Particle Count in Event")).fill(percFSParticles);
+ }
+
+ /* particle counts */
+ aida.cloud1D(makePlotName(k + ": Particle Count by Event")).fill(numHits);
+
+ info.resetTotals();
+
tree.cd("..");
}
@@ -266,6 +306,9 @@
/* # MCPs w/o endpoints */
mcpNoEndPoint.fill(noEndPointCount);
+
+ /* backscatter w/o vertexIsNotEndpointOfParent flag set */
+ mcpBackscatterNotPEP.fill(backscatterNotPEPCount);
}
}