lcsim/src/org/lcsim/recon/postrecon/leptonID/electron
diff -N AddBremPhotonsToElectrons.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ AddBremPhotonsToElectrons.java 27 May 2009 15:52:17 -0000 1.1
@@ -0,0 +1,127 @@
+package org.lcsim.recon.postrecon.leptonID.electron;
+import org.lcsim.util.Driver;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.MCParticle;
+import org.lcsim.event.Track;
+import org.lcsim.event.Cluster;
+import org.lcsim.event.ReconstructedParticle;
+import org.lcsim.event.base.BaseReconstructedParticle;
+import hep.physics.vec.*;
+import java.util.*;
+import java.io.IOException;
+public class AddBremPhotonsToElectrons extends Driver
+{
+ String rpl = "ReconstructedParticles";
+ String outrpn = "eehReconstructedParticles";
+//
+// Mass squared cut for associating photon with electron
+//
+ double msqcut = 0.49;
+//
+// # sigma cut on electron cluster match to ignore photon
+// correction
+//
+ double nsigematch = 4.;
+//
+// # sigma cut on combined cluster energy to ignore photon
+// correction
+//
+ double nsigtotignore = -999.;
+ public AddBremPhotonsToElectrons()
+ {
+ }
+ public void setMassSquaredCut(double x){msqcut = x;}
+ public void setOutputRPName(String x){outrpn = x;}
+ protected void process(EventHeader event)
+ {
+ List<ReconstructedParticle> el = new ArrayList<ReconstructedParticle>();
+ List<ReconstructedParticle> phl = new ArrayList<ReconstructedParticle>();
+ List<ReconstructedParticle> nrp = new ArrayList<ReconstructedParticle>(event.get(ReconstructedParticle.class,rpl));
+ List<ReconstructedParticle> rp = event.get(ReconstructedParticle.class,rpl);
+//
+// Find the reconstructed electrons and photons
+//
+ for(ReconstructedParticle p:rp)
+ {
+ if(p.getCharge() == 0.)
+ {
+ if(p.getMass() == 0.)phl.add(p);
+ }
+ else
+ {
+ if(p.getMass() < .001)
+ {
+ el.add(p);
+ }
+ }
+ }
+//
+// For each electron, loop over photons and use inv mass to
+// identify brem photons
+//
+ for(ReconstructedParticle pp:el)
+ {
+ Hep3Vector ep = pp.getMomentum();
+ double eE = pp.getEnergy();
+ List<ReconstructedParticle> addphl = new ArrayList<ReconstructedParticle>();
+ double addE = 0.;
+ for(ReconstructedParticle p:phl)
+ {
+ Hep3Vector php = p.getMomentum();
+ double phE = p.getEnergy();
+ double masssq = (eE+phE)*(eE+phE) - (ep.x()+php.x())*(ep.x()+php.x()) -
+ (ep.y()+php.y())*(ep.y()+php.y()) - (ep.z()+php.z())*(ep.z()+php.z());
+ if(masssq < msqcut)
+ {
+ addphl.add(p);
+ addE += phE;
+ }
+ }
+ double addf = 1. + addE/eE;
+ if(addE > 0.)
+ {
+//
+// Check for electron ID without an e/p match. In that case, assume combined brem and electron
+// shower, and just use the Track measurement.
+//
+ double clE = 0.;
+ for(Cluster c:pp.getClusters())
+ {
+ clE += c.getEnergy();
+ }
+ if(clE/eE < (1. - nsigematch*.18/Math.sqrt(eE)))addf = 1.;
+//
+// Also check that adding the brem photons exceed e/p match. If not,
+// assume the photons are electron fragments.
+//
+ if( (clE+addE)/eE < (1. + nsigtotignore*.18/Math.sqrt(eE)) )addf = 1.;
+//
+// Now modify the output ReconstructedParticle list
+//
+ BaseReconstructedParticle newe = new BaseReconstructedParticle(addf*eE,addf*ep.x(),addf*ep.y(),addf*ep.z());
+ newe.setMass(pp.getMass());
+ newe.setCharge(pp.getCharge());
+ newe.setReferencePoint(pp.getReferencePoint());
+ newe.setParticleIdUsed(pp.getParticleIDUsed());
+ newe.addTrack(pp.getTracks().get(0));
+ nrp.remove(pp);
+ nrp.add(newe);
+ for(Cluster c:pp.getClusters())
+ {
+ newe.addCluster(c);
+ }
+ newe.addParticle(pp);
+ for(ReconstructedParticle ph:addphl)
+ {
+ nrp.remove(ph);
+ newe.addParticle(ph);
+ for(Cluster c:ph.getClusters())
+ {
+ newe.addCluster(c);
+ }
+ }
+ }
+ }
+ event.put(outrpn,nrp,ReconstructedParticle.class,0);
+ }
+}
lcsim/src/org/lcsim/recon/postrecon/leptonID/electron
diff -N IsolatedHighPElectronIdentifier.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ IsolatedHighPElectronIdentifier.java 27 May 2009 15:52:17 -0000 1.1
@@ -0,0 +1,218 @@
+package org.lcsim.recon.postrecon.leptonID.electron;
+import org.lcsim.util.Driver;
+import hep.physics.particle.properties.ParticlePropertyManager;
+import hep.physics.particle.properties.ParticlePropertyProvider;
+import org.lcsim.event.ParticleID;
+import org.lcsim.event.base.BaseParticleID;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.MCParticle;
+import org.lcsim.event.Track;
+import org.lcsim.event.Cluster;
+import org.lcsim.event.ReconstructedParticle;
+import org.lcsim.event.base.BaseReconstructedParticle;
+import hep.physics.vec.*;
+import java.util.*;
+import java.io.IOException;
+public class IsolatedHighPElectronIdentifier extends Driver
+{
+ String rpl = "ReconstructedParticles";
+ double conect;
+ double maxE;
+ double minE;
+ final ParticlePropertyProvider dPPP;
+ ParticleID eppid;
+ ParticleID empid;
+ double emass;
+ Hep3Vector origin;
+ public IsolatedHighPElectronIdentifier()
+ {
+ conect = 0.90;
+ maxE = 125.;
+ minE = 10.;
+ dPPP = ParticlePropertyManager.getParticlePropertyProvider();
+ eppid = new BaseParticleID(dPPP.get(-11));
+ empid = new BaseParticleID(dPPP.get(11));
+ emass = dPPP.get(11).getMass();
+ origin = new BasicHep3Vector(0.,0.,0.);
+ }
+ public void setConeCosTheta(double x){conect = x;}
+ public void setMaxE(double x){maxE = x;}
+ public void setMinE(double x){minE = x;}
+ protected void process(EventHeader event)
+ {
+//
+// Match Tracks with Recon Particles
+//
+ List<Track> allTracks = event.get(Track.class,"Tracks");
+ List<Track> unusedTracks = new ArrayList<Track>(allTracks);
+ List<ReconstructedParticle> rp = event.get(ReconstructedParticle.class,rpl);
+ Map<Track,ReconstructedParticle> ttorp = new HashMap<Track,ReconstructedParticle>();
+ Map<Cluster,ReconstructedParticle> ctorp = new HashMap<Cluster,ReconstructedParticle>();
+ for(ReconstructedParticle p:rp)
+ {
+ if(p.getCharge() != 0.)
+ {
+ Track t = p.getTracks().get(0);
+ ttorp.put(t,p);
+ unusedTracks.remove(t);
+ }
+ else
+ {
+ if(p.getClusters().size() > 0)ctorp.put(p.getClusters().get(0),p);
+ }
+ }
+//
+// Map non-muon tracks > minE GeV to clusters within a 25 degree cone
+//
+ List<Track> xtrae = new ArrayList<Track>();
+ List<Cluster> cl = event.get(Cluster.class,"Clusters");
+ Map<Track,List<Cluster>> ttocll = new HashMap<Track,List<Cluster>>();
+ for(Track t:allTracks)
+ {
+ boolean eid = false;
+ Hep3Vector tp = new BasicHep3Vector(t.getMomentum());
+ if(tp.magnitude() < minE)continue;
+ if(tp.magnitude() > maxE)continue;
+ if(ttorp.containsKey(t))
+ {
+ double mass = ttorp.get(t).getMass();
+ if( Math.abs(mass - .105) < .005)continue;
+ if( Math.abs(mass - .0005) < .0001)eid = true;
+ }
+ List<Cluster> tcl = new ArrayList<Cluster>();
+ double HcalE = 0.;
+ for(Cluster c:cl)
+ {
+ Hep3Vector cp = new BasicHep3Vector(c.getPosition());
+ double ct = (tp.x()*cp.x()+tp.y()*cp.y()+tp.z()*cp.z())/tp.magnitude()/cp.magnitude();
+ double hec = c.getSubdetectorEnergies()[3] + c.getSubdetectorEnergies()[7];
+ if(ct > conect)
+ {
+ tcl.add(c);
+ HcalE += hec;
+ }
+ }
+ double f = HcalE/tp.magnitude();
+ if(!eid)
+ {
+ if(f < .04)xtrae.add(t);
+ }
+ ttocll.put(t,tcl);
+ }
+ if(xtrae.size() < 1)return;
+//
+// Have Tracks IDed as electrons. If we already have a ReconstructedParticle
+// for this Track, just replace it with clone with electron ID. Otherwise, we
+// need to make a new ReconstructedParticle and do some association
+// to determine which neutrals to remove.
+//
+ List<ReconstructedParticle> addlist = new ArrayList<ReconstructedParticle>();
+ List<ReconstructedParticle> removelist = new ArrayList<ReconstructedParticle>();
+ for(Track t:xtrae)
+ {
+ if(ttorp.containsKey(t))
+ {
+ ReconstructedParticle erp = makeElectronRP(ttorp.get(t));
+ addlist.add(erp);
+ removelist.add(ttorp.get(t));
+ }
+ else
+ {
+ BaseReconstructedParticle erp = makeElectronBRP(t);
+ if(ttocll.get(t).size() < 1)
+ {
+ addlist.add(erp);
+ }
+ else
+ {
+ List<Cluster> ocl = new ArrayList<Cluster>();
+ List<Cluster> mcl = new ArrayList<Cluster>(ttocll.get(t));
+ Hep3Vector eP = erp.getMomentum();
+ double eE = erp.getEnergy();
+ double cE = 0.;
+ for(int i=0;i<ttocll.get(t).size();i++)
+ {
+ double maxct = -1.;
+ Cluster c = null;
+ for(int j=0;j<mcl.size();j++)
+ {
+ Cluster tc = mcl.get(j);
+ Hep3Vector cp = new BasicHep3Vector(tc.getPosition());
+ double ct = (eP.x()*cp.x() + eP.y()*cp.y() + eP.z()*cp.z())/eP.magnitude()/cp.magnitude();
+ if(ct > maxct)
+ {
+ maxct = ct;
+ c = tc;
+ }
+ }
+ ocl.add(i,c);
+ mcl.remove(c);
+ }
+ List<Cluster> acl = new ArrayList<Cluster>();
+ for(Cluster c:ocl)
+ {
+ if( (cE + c.getEnergy()) < (eE + 3.*.2*Math.sqrt(eE)) )
+ {
+ cE += c.getEnergy();
+ acl.add(c);
+ }
+ }
+ addlist.add(erp);
+ double Erem = 0.;
+ int nrem = 0;
+ for(Cluster c:acl)
+ {
+ if(ctorp.containsKey(c))
+ {
+ erp.addCluster(c);
+ removelist.add(ctorp.get(c));
+ nrem++;
+ Erem += ctorp.get(c).getEnergy();
+ }
+ }
+ }
+ }
+ }
+ for(ReconstructedParticle p:removelist)
+ {
+ rp.remove(p);
+ }
+ for(ReconstructedParticle p:addlist)
+ {
+ rp.add(p);
+ }
+ }
+ public ReconstructedParticle makeElectronRP(ReconstructedParticle pi)
+ {
+ Hep3Vector m = pi.getMomentum();
+ double E = Math.sqrt(emass*emass + m.magnitudeSquared());
+ BaseReconstructedParticle p = new BaseReconstructedParticle(E,m);
+ p.setMass(emass);
+ p.setCharge(pi.getCharge());
+ p.setReferencePoint(pi.getReferencePoint());
+ if(pi.getCharge() > 0.)p.setParticleIdUsed(eppid);
+ else p.setParticleIdUsed(empid);
+ for(Cluster c:pi.getClusters())
+ {
+ p.addCluster(c);
+ }
+ for(Track c:pi.getTracks())
+ {
+ p.addTrack(c);
+ }
+ return p;
+ }
+ public BaseReconstructedParticle makeElectronBRP(Track t)
+ {
+ Hep3Vector m = new BasicHep3Vector(t.getMomentum());
+ double E = Math.sqrt(emass*emass + m.magnitudeSquared());
+ BaseReconstructedParticle p = new BaseReconstructedParticle(E,m);
+ p.setMass(emass);
+ p.setCharge((double) (t.getCharge()));
+ p.setReferencePoint(origin);
+ if(t.getCharge() > 0.)p.setParticleIdUsed(eppid);
+ else p.setParticleIdUsed(empid);
+ p.addTrack(t);
+ return p;
+ }
+}
