8 modified files
lcsim/src/org/lcsim/mc/fast/cluster/ronan
diff -u -r1.2 -r1.3
--- ClusterResolutionTables.java 15 Jul 2005 22:09:27 -0000 1.2
+++ ClusterResolutionTables.java 9 Aug 2005 18:34:45 -0000 1.3
@@ -8,7 +8,6 @@
private double EMConstantTerm;
private double EMPositionError;
private double EMResolution;
- private double EMmin;
private double EMOnset;
private double EMSharpness;
@@ -16,7 +15,6 @@
private double HADConstantTerm;
private double HADPositionError;
private double HADResolution;
- private double HADmin;
private double HADOnset;
private double HADSharpness;
@@ -27,7 +25,6 @@
ClusterResolutionTables(ConditionsSet set)
{
- EMmin = set.getDouble("EMmin");
EMOnset = set.getDouble("EMOnset");
EMSharpness = set.getDouble("EMSharpness");
PolarEMInner = set.getDouble("PolarEMInner");
@@ -38,7 +35,6 @@
EMPositionError = set.getDouble("EMPositionError");
EMAlignmentError = set.getDouble("EMAlignmentError");
- HADmin = set.getDouble("HADmin");
HADOnset = set.getDouble("HADOnset");
HADSharpness = set.getDouble("HADSharpness");
PolarHADInner = set.getDouble("PolarHADInner");
@@ -70,11 +66,6 @@
return EMResolution;
}
- public double getEMmin()
- {
- return EMmin;
- }
-
public double getEMOnset()
{
return EMOnset;
@@ -105,11 +96,6 @@
return HADResolution;
}
- public double getHADmin()
- {
- return HADmin;
- }
-
public double getHADOnset()
{
return HADOnset;
lcsim/src/org/lcsim/mc/fast/cluster/ronan
diff -u -r1.6 -r1.7
--- MCFastRonan.java 4 Aug 2005 02:23:23 -0000 1.6
+++ MCFastRonan.java 9 Aug 2005 18:34:45 -0000 1.7
@@ -9,10 +9,14 @@
import org.lcsim.conditions.ConditionsEvent;
import org.lcsim.conditions.ConditionsListener;
import org.lcsim.conditions.ConditionsSet;
+import org.lcsim.event.MCParticle;
import org.lcsim.event.Cluster;
import org.lcsim.event.EventHeader;
import org.lcsim.util.Driver;
+
+
+
/**
* Fast Monte Carlo cluster simulator
* @author M.Ronan Oct 2000 - Added "refined" cluster simulation
@@ -34,19 +38,20 @@
clusterParm = new ClusterResolutionTables(conditions);
}
List cl = new ArrayList();
- List p1 = new ArrayList();
+
boolean hist = getHistogramLevel() > 0;
- for (Particle p : event.getMCParticles())
- {
+ List<MCParticle> particles = event.getMCParticles();
+ for (MCParticle p : particles)
+ {
+
// filter for FINALSTATE
if (p.getGeneratorStatus() != p.FINAL_STATE)
{
continue;
}
- ParticleType ptype = p.getType();
- int PDGID = ptype.getPDGID();
+ int PDGID = p.getPDGID();
double charge = p.getCharge();
// filter neutrinos
@@ -57,10 +62,10 @@
}
double E = p.getEnergy();
- double pt2 = (p.getPX() * p.getPX()) + (p.getPY() * p.getPY());
+ double pt2 = p.getMomentum().magnitudeSquared()-p.getPZ()*p.getPZ();
double pt = Math.sqrt(pt2);
- double ptot = Math.sqrt(pt2 + (p.getPZ() * p.getPZ()));
- double cosTheta = p.getPZ() / ptot;
+ double ptot = p.getMomentum().magnitude();
+ double cosTheta = p.getPZ() / p.getMomentum().magnitude();
Random rand = getRandom();
@@ -68,7 +73,6 @@
if (PDGID == PhotonID)
{
// within acceptance
- double t = clusterParm.getEMOnset();
double thing = (1 - 1 / ( 1 + Math.exp( (E-clusterParm.getEMOnset())*clusterParm.getEMSharpness() ) ));
if (rand.nextDouble() > thing)
{
@@ -91,11 +95,11 @@
double thing = (1 - 1 / ( 1 + Math.exp( (E-clusterParm.getHADOnset())*clusterParm.getHADSharpness() ) ));
if (rand.nextDouble() > thing)
{
- continue;
+ //continue;
}
if (Math.abs(cosTheta) > clusterParm.getPolarHADOuter())
{
- continue;
+ //continue;
}
cl.add(new ReconHADCluster(clusterParm, rand, p, hist));
lcsim/src/org/lcsim/mc/fast/cluster/ronan
diff -u -r1.3 -r1.4
--- ReconCluster.java 4 Aug 2005 02:23:23 -0000 1.3
+++ ReconCluster.java 9 Aug 2005 18:34:45 -0000 1.4
@@ -13,7 +13,6 @@
protected double a = 0;
protected double b = 0;
protected double c = 0;
- protected double cosTh;
protected double d = 0;
protected double energy;
protected double phi;
@@ -49,10 +48,11 @@
protected void smearEnergy(Random rand, double E, boolean hist)
{
double sigma = ((a / Math.sqrt(E)) + b) * E;
- energy = E + (sigma * rand.nextGaussian());
- if (energy < mcp.getMass())
+
+ energy = 0;
+ while (energy <= mcp.getMass())
{
- energy = mcp.getMass() + 0.05;
+ energy = E + (sigma * rand.nextGaussian());
}
}
@@ -63,17 +63,14 @@
double Py = mcp.getPY();
double Pz = mcp.getPZ();
- double Pt2 = (Px * Px) + (Py * Py);
- double Pt = Math.sqrt(Pt2);
- double P = Math.sqrt(Pt2 + (Pz * Pz));
+ double P = Math.sqrt((Px * Px) + (Py * Py) + (Pz * Pz));
double Phi = Math.atan2(Py, Px);
if (Phi < 0)
{
Phi += (2 * Math.PI);
}
- double CosTh = Pz / P;
- double Theta = Math.acos(CosTh);
+ double Theta = Math.acos(Pz / P);
// Simulate position smearing on a sphere of radius 2 meters
radius = 2000.0;
@@ -82,19 +79,23 @@
double y = (radius * Py) / P;
double z = (radius * Pz) / P;
- double drPhi = transDist * CosTh;
- double transPhi = Math.PI * rand.nextDouble();
- x = x + (drPhi * Math.sin(transPhi));
- y = y + (drPhi * Math.cos(transPhi));
- z = z + (transDist * Math.sin(Theta));
+ //these vectors vt and vs (orthonorm.) span a plane perpendicular to the momentum vector,
+ //so smearing with a transdist will involve a lin. comb. of these
+ double[] vt = {-Math.cos(Theta)*Math.cos(Phi),-Math.cos(Theta)*Math.sin(Phi),Math.sin(Theta)};
+ double[] vs = {Math.sin(Phi),-Math.cos(Phi),0};
+
+ //restricted to [0,PI] since transdist can be negative
+ double alpha = rand.nextDouble() * Math.PI;
+ x = x + transDist * ( Math.cos(alpha) * vt[0] + Math.sin(alpha) * vs[0] );
+ y = y + transDist * ( Math.cos(alpha) * vt[1] + Math.sin(alpha) * vs[1] );
+ z = z + transDist * ( Math.cos(alpha) * vt[2] + Math.sin(alpha) * vs[2] );
phi = Math.atan2(y, x);
if (phi < 0)
{
phi += (2 * Math.PI);
}
- cosTh = z / radius;
- theta = Math.acos(cosTh);
+ theta = Math.acos(z / radius);
}
public Particle getMCParticle()
lcsim/src/org/lcsim/mc/fast/reconstructedparticle
diff -u -r1.7 -r1.8
--- MCFastReconstructedParticleDriver.java 18 Jul 2005 19:18:26 -0000 1.7
+++ MCFastReconstructedParticleDriver.java 9 Aug 2005 18:34:45 -0000 1.8
@@ -14,6 +14,7 @@
import java.util.ArrayList;
import java.util.List;
import java.util.Random;
+import org.lcsim.event.MCParticle;
import org.lcsim.event.Cluster;
import org.lcsim.event.EventHeader;
import org.lcsim.event.ReconstructedParticle;
@@ -55,6 +56,8 @@
protected void process(EventHeader event)
{
+ boolean hist = getHistogramLevel() > 0;
+
if (IDEff == null)
{
ConditionsSet conditions = getConditionsManager().getConditions("IDEfficiency");
@@ -80,23 +83,37 @@
// electrons and muons are special
if((abs(pdgid)== 11) && (rand.nextDouble() < IDEff.getElectronEff()))
{
- type = new Electron((int)Math.signum(pdgid));
+ type = new Electron(rt.getCharge());
}
else if((abs(pdgid)== 13) && (rand.nextDouble() < IDEff.getMuonEff()))
{
- type = new Muon((int)Math.signum(pdgid));
+ type = new Muon(rt.getCharge());
}
else
{
type = new Pion(rt.getCharge());
}
- aida.histogram1D("track-particle", 150, -10, 10).fill( (t.getPX()*t.getPX() + t.getPY()*t.getPY() + t.getPZ()*t.getPZ() + type.mass()*type.mass()) - p.getEnergy());
+ if ((p.getEnergy() > 10) && (hist))
+ {
+ if (Math.abs(t.getTrackParameter(4)) < 1)
+ {
+ aida.histogram1D("track-particle", 150, -0.0003, 0.0003).fill( ( Math.sqrt(t.getPX()*t.getPX() + t.getPY()*t.getPY() + t.getPZ()*t.getPZ()) - Math.sqrt(p.getPX()*p.getPX() + p.getPY()*p.getPY() + p.getPZ()*p.getPZ()) ) / ( p.getPX()*p.getPX() + p.getPY()*p.getPY() + p.getPZ()*p.getPZ()) );
+ }
+ if ((Math.abs(t.getTrackParameter(4)) < 2.16) && (Math.abs(t.getTrackParameter(4)) > 1.96) && (p.getMomentum().magnitude() > 80) && (p.getMomentum().magnitude() < 120))
+ {
+ aida.histogram1D("track-particle-cut", 150, -0.01, 0.01).fill( ( Math.sqrt(t.getPX()*t.getPX() + t.getPY()*t.getPY() + t.getPZ()*t.getPZ()) - Math.sqrt(p.getPX()*p.getPX() + p.getPY()*p.getPY() + p.getPZ()*p.getPZ()) ) / (Math.sqrt( p.getPX()*p.getPX() + p.getPY()*p.getPY() + p.getPZ()*p.getPZ())) );
+ }
+ }
// assume pion for remaining charged tracks
MCFastReconstructedParticle rp = new MCFastReconstructedParticle(t, type, p);
rpList.add(rp);
- aida.histogram1D("recon-particle", 150, -10, 10).fill(rp.getEnergy()-p.getEnergy());
+
+ if (hist)
+ {
+ aida.histogram1D("recon-particle", 150, -5, 5).fill((rp.getEnergy()-p.getEnergy())/(Math.sqrt(p.getEnergy())));
+ }
}
}
@@ -104,15 +121,18 @@
List<Cluster> clusters = event.getClusters();
for(Cluster c : clusters)
{
- ParticleType type = null;
Particle p = null;
+ ParticleType type = null;
// photons for EM
if( c instanceof ReconEMCluster)
{
ReconEMCluster emc = (ReconEMCluster) c;
p = emc.getMCParticle();
type = new Photon();
- aida.histogram1D("cluster-particle", 150, -10, 10).fill(emc.getEnergy()-emc.getMCParticle().getEnergy());
+ if (hist)
+ {
+ aida.histogram1D("photonCLS-particle", 150, -3, 3).fill((emc.getEnergy()-emc.getMCParticle().getEnergy())/(Math.sqrt(emc.getMCParticle().getEnergy())));
+ }
}
// assume a KZeroLong here for had cluster
else if(c instanceof ReconHADCluster)
@@ -120,7 +140,10 @@
ReconHADCluster emc = (ReconHADCluster) c;
p = emc.getMCParticle();
int pdgid = p.getPDGID();
- aida.histogram1D("cluster-particle", 150, -10, 10).fill(emc.getEnergy()-emc.getMCParticle().getEnergy());
+ if (hist)
+ {
+ aida.histogram1D("hadronCLS-particle", 150, -3, 3).fill((emc.getEnergy()-emc.getMCParticle().getEnergy())/(Math.sqrt(emc.getMCParticle().getEnergy())));
+ }
if ((abs(pdgid)==2112) && (rand.nextDouble() < IDEff.getNeutronEff()))
{
@@ -133,7 +156,11 @@
}
MCFastReconstructedParticle rp = new MCFastReconstructedParticle(c, type, p);
rpList.add(rp);
- aida.histogram1D("recon-particle", 150, -10, 10).fill(rp.getEnergy()-p.getEnergy());
+ if (hist)
+ {
+ aida.histogram1D("recon-particle", 150, -10, 10).fill(rp.getEnergy()-p.getEnergy());
+ }
+
}
// add the reconstructedparticles to the event
event.put(event.RECONSTRUCTEDPARTICLES, rpList, ReconstructedParticle.class, 0);
lcsim/src/org/lcsim/mc/fast/tracking
diff -u -r1.5 -r1.6
--- MCFastTracking.java 18 Jul 2005 19:18:26 -0000 1.5
+++ MCFastTracking.java 9 Aug 2005 18:34:45 -0000 1.6
@@ -20,7 +20,9 @@
public class MCFastTracking extends Driver implements ConditionsListener
{
private TrackResolutionTables parm;
+ private SimpleTables SmTbl;
private boolean beamSpotConstraint;
+ private boolean simple;
private final static double[] IP = { 0, 0, 0 };
public MCFastTracking()
@@ -32,6 +34,12 @@
{
this.beamSpotConstraint = beamSpotConstraint;
}
+
+ public MCFastTracking(boolean beamSpotConstraint, boolean simple)
+ {
+ this.beamSpotConstraint = beamSpotConstraint;
+ this.simple = simple;
+ }
public void setBeamSpotConstraint(boolean beamSpotConstraint)
{
@@ -68,6 +76,13 @@
parm = setTrackResolutionTables(conditions, beamSpotConstraint);
}
+ if (SmTbl == null)
+ {
+ ConditionsSet conditions = getConditionsManager().getConditions("SimpleTrack");
+ conditions.addConditionsListener(this);
+ SmTbl = new SimpleTables(conditions);
+ }
+
double bField = event.getDetector().getFieldMap().getField(IP)[2];
boolean hist = getHistogramLevel() > 0;
@@ -103,7 +118,7 @@
continue;
}
- trackList.add(new ReconTrack(bField, parm, getRandom(), p, hist));
+ trackList.add(new ReconTrack(bField, parm, SmTbl, getRandom(), p, hist, simple));
}
event.put(EventHeader.TRACKS, trackList, Track.class, 0);
}
@@ -111,6 +126,8 @@
public void conditionsChanged(ConditionsEvent event)
{
ConditionsSet conditions = getConditionsManager().getConditions("TrackParameters");
+ ConditionsSet simpleconditions = getConditionsManager().getConditions("SimpleTrack");
parm = setTrackResolutionTables(conditions, beamSpotConstraint);
+ SmTbl = new SimpleTables(simpleconditions);
}
}
lcsim/src/org/lcsim/mc/fast/tracking
diff -u -r1.5 -r1.6
--- ReconTrack.java 2 Jul 2005 06:40:45 -0000 1.5
+++ ReconTrack.java 9 Aug 2005 18:34:45 -0000 1.6
@@ -9,6 +9,7 @@
import java.util.Collections;
import java.util.List;
import java.util.Random;
+import org.lcsim.mc.fast.tracking.SimpleTables;
/**
@@ -17,7 +18,7 @@
* are provided. <br>
*
* @author Tony Johnson, Wolfgang Walkowiak
- * @version $Id: ReconTrack.java,v 1.5 2005/07/02 06:40:45 ngraf Exp $
+ * @version $Id: ReconTrack.java,v 1.6 2005/08/09 18:34:45 ngraf Exp $
*/
public class ReconTrack implements Track
{
@@ -40,20 +41,21 @@
transient private Particle mc;
private int m_tcharge;
- ReconTrack(double bField, TrackResolutionTables parm, Random rand, Particle mc, boolean hist)
+ ReconTrack(double bField, TrackResolutionTables parm, SimpleTables SmTbl, Random rand, Particle mc, boolean hist, boolean simple)
{
this.mc = mc;
// get original momentum from MCParticle
// convert to helical parameters
m_nosmear = new DocaTrackParameters(bField, mc);
+
+ double pt = m_nosmear.getPt();
if (hist)
{
double r = Math.abs(m_nosmear.getD0());
- double pt = m_nosmear.getPt();
-
+
AIDA aida = AIDA.defaultInstance();
aida.cloud1D("ptsqr").fill(pt*pt);
aida.cloud1D("pt").fill(pt);
@@ -72,8 +74,15 @@
m_nosmear.fillErrorMatrix(getErrorMatrixFromTable(table, abscth, ptot));
// smear tracks according to error matrix
- m_smear = (DocaTrackParameters) SmearTrack.smearTrack(bField, m_nosmear, rand);
-
+ if (simple == true)
+ {
+ m_smear = (DocaTrackParameters) SmearTrackSimple.smearTrackSimple(bField, m_nosmear, rand, SmTbl, pt, hist);
+ //double[] slice = {0, 0, 1, 0, 0};
+ //m_smear = (DocaTrackParameters) SmearTrackSimpleII.SmearTrackSimpleII(bField, m_nosmear, rand, SmTbl, slice, hist);
+ } else {
+ m_smear = (DocaTrackParameters) SmearTrack.smearTrack(bField, m_nosmear, rand);
+ }
+
if (hist)
{
AIDA aida = AIDA.defaultInstance();
lcsim/src/org/lcsim/mc/fast/tracking
diff -u -r1.1 -r1.2
--- SmearTrack.java 1 Feb 2005 19:42:48 -0000 1.1
+++ SmearTrack.java 9 Aug 2005 18:34:45 -0000 1.2
@@ -18,6 +18,7 @@
*/
static DocaTrackParameters smearTrack(double bField, TrackParameters noSmear, Random rand)
{
+
// get error matrix and do a sanity check
Matrix M = Matrix.constructWithCopy(noSmear.getErrorMatrix());
if (M.det() <= 0.)
lcsim/src/org/lcsim/mc/fast
diff -u -r1.2 -r1.3
--- MCFast.java 1 Jul 2005 23:51:57 -0000 1.2
+++ MCFast.java 9 Aug 2005 18:34:46 -0000 1.3
@@ -12,10 +12,15 @@
public class MCFast extends Driver
{
/** Creates a new instance of MCFast */
- public MCFast()
+ public MCFast(boolean beamSpotConstraint, boolean simple)
{
- add(new MCFastTracking());
+ add(new MCFastTracking(beamSpotConstraint, simple));
add(new MCFastRonan());
add(new MCFastReconstructedParticleDriver());
}
+
+ public MCFast()
+ {
+ this(false,false);
+ }
}
CVSspam 0.2.8