lcsim-contrib/src/main/java/org/lcsim/contrib/SteveMagill
diff -N PFARPdijetm.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ PFARPdijetm.java 14 Feb 2012 19:33:14 -0000 1.1
@@ -0,0 +1,105 @@
+package org.lcsim.contrib.SteveMagill;
+
+import java.util.*;
+import org.lcsim.event.EventHeader;
+import org.lcsim.util.Driver;
+import hep.physics.vec.Hep3Vector;
+import org.lcsim.event.ReconstructedParticle;
+import org.lcsim.util.aida.AIDA;
+
+/*
+ * PFASelect includes fiducial cut to select region for PFA analysis
+ */
+
+public class PFARPdijetm extends Driver
+{
+ private String _scinjets;
+ private String _jdlinkmap;
+ private String _stjetmap;
+ private AIDA aida = AIDA.defaultInstance();
+
+ public PFARPdijetm()
+ {
+
+ }
+
+ protected void process(EventHeader event)
+ {
+ super.process(event);
+
+ double sEtot = 0;
+ double sjpx = 0;
+ double sjpy = 0;
+ double sjpz = 0;
+ double ujm2 = 0.;
+ // get E M of qqbar
+ List<Double> qqm = event.get(Double.class, "DJmass");
+ List<Double> qqE = event.get(Double.class, "DJEnergy");
+ // get corrected jets
+ List<ReconstructedParticle> scjets = event.get(ReconstructedParticle.class,_scinjets);
+ // get scintillator jet-track jet link map
+ Map<ReconstructedParticle, ReconstructedParticle> jtmap = (Map<ReconstructedParticle, ReconstructedParticle>) event.get(_stjetmap);
+ // get jet-deltam map
+ Map<ReconstructedParticle, Double> jdmap = (Map<ReconstructedParticle, Double>) event.get(_jdlinkmap);
+ if (scjets.size() == 2)
+ {
+ for (ReconstructedParticle scjet : scjets)
+ {
+ double jdm = 0.;
+ sEtot += scjet.getEnergy();
+ double jm = scjet.getMass();
+ aida.cloud1D("PFA Jet Mass before dM Corr").fill(jm);
+ if (jtmap.containsKey(scjet))
+ {
+ ReconstructedParticle jtl = jtmap.get(scjet);
+ if (jdmap.containsKey(jtl)) jdm = jdmap.get(jtl);
+ }
+ double jetm = jm-jdm;
+ if (jetm<0.) jetm = 0.;
+ aida.cloud1D("dM for PFA Jet").fill(jdm);
+ aida.cloud1D("PFA Jet Mass after dM Corr").fill(jetm);
+ double corrpmag = Math.sqrt(scjet.getEnergy()*scjet.getEnergy()-jetm*jetm);
+ double pmag = scjet.getMomentum().magnitude();
+ double cp = corrpmag/pmag;
+ Hep3Vector jp = scjet.getMomentum();
+ sjpx += jp.x()*cp;
+ sjpy += jp.y()*cp;
+ sjpz += jp.z()*cp;
+// ujm2 += scjet.getMass()*scjet.getMass();
+ ujm2 += jetm*jetm;
+ }
+ double sptot = Math.sqrt(sjpx*sjpx+sjpy*sjpy+sjpz*sjpz);
+ double sjmass = Math.sqrt(sEtot*sEtot-sptot*sptot+ujm2);
+ aida.cloud1D("PFA DiJet Mass dM").fill(sjmass);
+ aida.cloud1D("PFA DiJet ETot dM").fill(sEtot);
+ double jm = 0.;
+ for (Double qm : qqm)
+ {
+ jm = qm.doubleValue();
+ }
+ double jE = 0.;
+ for (Double qE : qqE)
+ {
+ jE = qE.doubleValue();
+ }
+ aida.cloud1D("Diff PFA DiJet Mass dM").fill(sjmass-jm);
+ aida.cloud1D("Diff PFA DiJet Energy dM").fill(sEtot-jE);
+ }
+ }
+
+ public void setInputScinJets(String scinjets)
+ {
+ _scinjets = scinjets;
+ }
+
+ public void setSTJetMap(String stjetmap)
+ {
+ _stjetmap = stjetmap;
+ }
+
+ public void setTJetDMMap(String jdlinkmap)
+ {
+ _jdlinkmap = jdlinkmap;
+ }
+
+}