Commit in lcsim-contrib/src/main/java/org/lcsim/contrib/homer/vvHratio on MAIN | |||
MyDBDvvH.java | +39 | -5 | 1.1 -> 1.2 |
MyDBDvvHAnalysis.java | +133 | -51 | 1.1 -> 1.2 |
PFOJetFindingDriver.java | +5 | -3 | 1.1 -> 1.2 |
+177 | -59 |
progress
diff -u -r1.1 -r1.2 --- MyDBDvvH.java 11 Sep 2012 23:41:08 -0000 1.1 +++ MyDBDvvH.java 27 Sep 2012 09:00:41 -0000 1.2 @@ -1,22 +1,56 @@
import java.util.List;
+//import org.lcsim.event.*;
import org.lcsim.event.EventHeader; import org.lcsim.event.MCParticle; import org.lcsim.event.ReconstructedParticle; import org.lcsim.util.Driver; import org.lcsim.util.aida.AIDA;
+import org.lcsim.recon.postrecon.leptonID.electron.*; + +import java.util.*; +import java.io.IOException; + +
public class MyDBDvvH extends Driver {
+ + String rplnew = "eehReconstructedParticles"; + String rpl = "ReconstructedParticles"; + +
private AIDA aida = AIDA.defaultInstance(); public MyDBDvvH(){
+// +// Overwrite "ReconstructedParticles" with extra electron ID +// + add(new IsolatedHighPElectronIdentifier()); +// +// Make "eehReconstructedParticles" by adding brem photons to electrons +// +// add(new AddBremPhotonsToElectrons()); +
// add(new MyJetFindingDriver());
- add(new PFOJetFindingDriver()); - add(new MyDBDvvHAnalysis());
+ add(new PFOJetFindingDriver()); + add(new MyDBDvvHAnalysis()); + +
};
+ protected String mcParticlesSkimmedName = "MCParticlesSkimmed"; +
public void process(EventHeader event) {
- super.process(event); // this takes care that the child Drivers are loaded and processed.
+ List<MCParticle> particles = event.get(MCParticle.class,mcParticlesSkimmedName); + int doit=0; + // a nu as the first particle is sufficient to identify this as nunuH + // DON'T FORGET TO REMOVE FOR BACKGROUNDS! + for (MCParticle particle : particles) { + if (Math.abs(particle.getPDGID())==12) doit=1; + break; + } + System.out.println("doit = " + doit); + if (doit==1) + super.process(event); // this takes care that the child Drivers are loaded and processed.
// List<List<MCParticle>> jetList = (List<List<MCParticle>>) event.get("jetList");
- List<List<ReconstructedParticle>> jetList = (List<List<ReconstructedParticle>>) event.get("jetList"); - aida.cloud1D("nJets").fill(jetList.size());
+
} }
\ No newline at end of file
diff -u -r1.1 -r1.2 --- MyDBDvvHAnalysis.java 11 Sep 2012 23:41:08 -0000 1.1 +++ MyDBDvvHAnalysis.java 27 Sep 2012 09:00:41 -0000 1.2 @@ -1,5 +1,6 @@
import org.lcsim.util.aida.AIDA; import hep.physics.vec.VecOp;
+import hep.physics.vec.Hep3Vector;
import java.util.List; import org.lcsim.util.Driver; import org.lcsim.event.*;
@@ -13,13 +14,14 @@
import org.lcsim.spacegeom.SpaceVector;
+
public class MyDBDvvHAnalysis extends Driver {
+
// collections protected String inputPfoCollection = "PandoraPFOCollection"; protected String mcParticlesSkimmedName = "MCParticlesSkimmed";
-
// First get the reconstructed particles from the event List<ReconstructedParticle> PFOs = new ArrayList<ReconstructedParticle>();
@@ -35,38 +37,103 @@
// print(HLEVEL_DEFAULT, e.getMessage(), true); // }
+ // check for event with reconstruction errors + int skipEvent = 0; + for (ReconstructedParticle pfo : PFOs) { + boolean passPfoSelection = true; + + int type = pfo.getType(); +// double cosTheta = momentum.cosTheta(); + double energy = pfo.getEnergy(); + // System.out.println("energy = "+energy + "type = " + type); + if (energy>1000.) { + skipEvent = 1; + continue; + } + } + if (skipEvent == 1) return; // skip this bad event + +
List<Track> tracklist = event.getTracks();
- int nTrks = tracklist.size();
+ // int nTrks = tracklist.size();
+
+ List<List<ReconstructedParticle>> jetList = (List<List<ReconstructedParticle>>) event.get("jetList");
+ int nJets = jetList.size();
+ aida.cloud1D("Trk #all Tracks").fill(tracklist.size()); + aida.cloud1D("#Jets from PFO in ALL Events").fill(jetList.size());
- aida.cloud1D("Trk nTracks").fill(nTrks);
- -// List<MCParticle> particles = event.get(MCParticle.class,event.MC_PARTICLES);
+ // List<MCParticle> particles = event.get(MCParticle.class,event.MC_PARTICLES);
List<MCParticle> particles = event.get(MCParticle.class,mcParticlesSkimmedName);
- aida.cloud1D("Number of MC particles").fill(particles.size()); - for (MCParticle particle : particles) { - for (MCParticle mcp2 : particle.getParents()) { - int parentid = mcp2.getPDGID(); - if (Math.abs(parentid)==25) { - int id = particle.getPDGID(); - aida.cloud1D("All MC energy").fill(particle.getEnergy()); - aida.cloud1D("All MC cosTheta").fill(VecOp.cosTheta(particle.getMomentum())); - aida.cloud1D("All MC phi").fill(VecOp.phi(particle.getMomentum())); - aida.cloud1D("All MC H daughter ID").fill(id); - System.out.println("Parent type:"+mcp2.getType().getName()); - System.out.println("Daughter type:"+particle.getType().getName()); -// System.out.println("Starting point: "+mcp2.getOrigin().toString()); -// System.out.println("Momentum: "+mcp2.getMomentum()); -// System.out.println("Endpoint: "+mcp2.getEndPoint()); -// System.out.println("# Daughters: "+mcp2.getDaughters().size()); - System.out.println("Daughter id:"+id); - } - } - } -
+ aida.cloud1D("Number of MC particles").fill(particles.size()); + double MCM2vis = 0.0; + double EsumMC = 0.0; + double PsumMCX = 0.0; + double PsumMCY = 0.0; + double PsumMCZ = 0.0; + + int nid=0; + for (MCParticle particle : particles) { + // System.out.println("Sim stat: " + nid + " " + (particle.getSimulatorStatus()).getValue()); + if ((particle.getSimulatorStatus()).getValue()==67108864) { + // || + // (particle.getSimulatorStatus()).getValue()==16777216 ) { + double[] momentum = (particle.getMomentum()).v(); + double PX = momentum[0]; + double PY = momentum[1]; + double PZ = momentum[2]; + double PT = Math.sqrt(PX*PX + PY*PY); + double P = Math.sqrt(PT*PT+PZ*PZ); + aida.cloud1D("All MC particle cosTheta").fill(VecOp.cosTheta(particle.getMomentum())); + // System.out.println("Sim stat: " + nid + " " + (particle.getSimulatorStatus()).getValue() + " " + PX + " " + PY + " " + PZ + " " + P); + if (Math.abs(VecOp.cosTheta(particle.getMomentum()))<0.95) { + EsumMC += P; + PsumMCX += PX; + PsumMCY += PY; + PsumMCZ += PZ; + // System.out.println("Sums: " + nid + " " + (particle.getSimulatorStatus()).getValue() + " " + PsumMCX + " " + PsumMCY + " " + PsumMCZ + " " + EsumMC); + } + } + //don't want + //33554432 left + //1107296256 left + //134217728 tracker + //want + //67108864 cal + //67108864 cal + //16777216 stopped + + + for (MCParticle mcp2 : particle.getParents()) { + int parentid = mcp2.getPDGID(); + if (Math.abs(parentid)==25) { + int id = particle.getPDGID(); + aida.cloud1D("All MC H daughters energy").fill(particle.getEnergy()); + aida.cloud1D("All MC H daughters cosTheta").fill(VecOp.cosTheta(particle.getMomentum())); + aida.cloud1D("All MC H daughters phi").fill(VecOp.phi(particle.getMomentum())); + aida.cloud1D("All MC H daughter ID").fill(id); + System.out.println("Parent type:"+mcp2.getType().getName()); + System.out.println("Daughter type:"+particle.getType().getName()); + // System.out.println("Starting point: "+mcp2.getOrigin().toString()); + // System.out.println("Momentum: "+mcp2.getMomentum()); + // System.out.println("Endpoint: "+mcp2.getEndPoint()); + // System.out.println("# Daughters: "+mcp2.getDaughters().size()); + System.out.println("Daughter id:"+id); + } + } + nid++; + } + double EvisMC = EsumMC; + double PTvisMC = Math.sqrt(PsumMCX*PsumMCX+PsumMCY*PsumMCY); + double M2visMC = EsumMC*EsumMC - PTvisMC*PTvisMC - PsumMCZ*PsumMCZ; + aida.cloud1D("Evis from MC").fill(EvisMC); + aida.cloud1D("PTvis from MC").fill(PTvisMC); + aida.cloud1D("M2vis from MC").fill(M2visMC); + if (M2visMC > 0.) aida.cloud1D("Mvis from MC").fill(Math.sqrt(M2visMC));
// --------------------------------------------------
+ int nTrks = 0;
double EsumTrk = 0.0; double PsumTrkX = 0.0; double PsumTrkY = 0.0;
@@ -75,38 +142,42 @@
// aida.cloud1D("Trk cosTheta").fill(track.getMomentum()); TrackState trkstate = track.getTrackStates().get(0); double[] momentum = trkstate.getMomentum();
- System.out.println("mom = " + momentum[0] + "," + momentum[1] + "," + momentum[2]);
+ // System.out.println("mom = " + momentum[0] + "," + momentum[1] + "," + momentum[2]);
double PX = momentum[0]; double PY = momentum[1]; double PZ = momentum[2]; double PT = Math.sqrt(PX*PX + PY*PY); double P = Math.sqrt(PT*PT+PZ*PZ);
- EsumTrk += P; - PsumTrkX += PX; - PsumTrkY += PY; - PsumTrkZ += PY;
// System.out.println("PT = "+PT);
-// System.out.println("P = "+P);
+// System.out.println("P = "+P+" EsumTrk = "+EsumTrk);
double q = trkstate.getOmega()/Math.abs(trkstate.getOmega()); // System.out.println("q = "+q); double tlam = trkstate.getTanLambda(); // System.out.println("tlam = "+tlam); aida.cloud1D("Trk cosTheta").fill(PZ/P);
-
+ if (Math.abs(PZ/P)<0.95) { + if (P>0.5) nTrks++; + EsumTrk += P; + PsumTrkX += PX; + PsumTrkY += PY; + PsumTrkZ += PZ; + }
} double EvisTrk = EsumTrk; double PTvisTrk = Math.sqrt(PsumTrkX*PsumTrkX+PsumTrkY*PsumTrkY); double M2visTrk = EsumTrk*EsumTrk - PTvisTrk*PTvisTrk - PsumTrkZ*PsumTrkZ;
+ aida.cloud1D("#good Tracks").fill(nTrks);
aida.cloud1D("Evis from Tracks").fill(EvisTrk); aida.cloud1D("PTvis from Tracks").fill(PTvisTrk);
- aida.cloud1D("M2vis from Tracks").fill(M2visTrk); -
+ // aida.cloud1D("M2vis from Tracks").fill(M2visTrk); + if (M2visTrk>0.) aida.cloud1D("Mvis from Tracks").fill(Math.sqrt(M2visTrk));
- if ((PTvisTrk>20.0 && PTvisTrk<500.0) && - (EvisTrk>100.0 && EvisTrk<400.0) && - (nTrks>=6 && nTrks<=19)) {
+ if ((PTvisTrk>20.0 && PTvisTrk<250.0) && + (EvisTrk>100.0 && EvisTrk<500.0) && + (nTrks>=6 && nTrks<=45)) {
aida.cloud1D("Select M2vis from Tracks").fill(M2visTrk);
+ aida.cloud1D("#Jets from PFO in Track Selected Events").fill(jetList.size());
} // ----------------------------------------------------------------------------- double EsumPFO = 0.0;
@@ -120,39 +191,50 @@
int type = pfo.getType(); SpaceVector momentum = new SpaceVector(pfo.getMomentum()); double pT = momentum.rxy();
- double p = momentum.rxyz();
+ double P = momentum.rxyz();
double PX = momentum.x(); double PY = momentum.y(); double PZ = momentum.z(); double cosTheta = momentum.cosTheta(); double energy = pfo.getEnergy();
- EsumPFO += energy; - PsumPFOX += PX; - PsumPFOY += PY; - PsumPFOZ += PZ; - eTotalInput += energy; -// List<Cluster> pfoclusters = pfo.getClusters(); -// List<Track> pfotracks = pfo.getTracks();
+ // System.out.println("energy = "+energy + "type = " + type); + if (Math.abs(PZ/P)<0.95 && energy<1000.) { + EsumPFO += energy; + PsumPFOX += PX; + PsumPFOY += PY; + PsumPFOZ += PZ; + eTotalInput += energy; + // List<Cluster> pfoclusters = pfo.getClusters(); + // List<Track> pfotracks = pfo.getTracks(); + }
} double EvisPFO = EsumPFO; double PTvisPFO = Math.sqrt(PsumPFOX*PsumPFOX+PsumPFOY*PsumPFOY); double M2visPFO = EsumPFO*EsumPFO - PTvisPFO*PTvisPFO - PsumPFOZ*PsumPFOZ; aida.cloud1D("Evis from PFO").fill(EvisPFO); aida.cloud1D("PTvis from PFO").fill(PTvisPFO);
- aida.cloud1D("M2vis from PFO").fill(M2visPFO);
+ if (M2visPFO>0.) aida.cloud1D("Mvis from PFO").fill(Math.sqrt(M2visPFO));
// ------------------MC parent info for PFO selected events ----------------------------------------
- if ((PTvisPFO>20.0 && PTvisPFO<500.0) && - (EvisPFO>100.0 && EvisPFO<400.0) && - (nTrks>=6 && nTrks<=19)) { - aida.cloud1D("Select M2vis from PFO").fill(M2visPFO);
+ if ((PTvisPFO>20.0 && PTvisPFO<250.0) && + (EvisPFO>100.0 && EvisPFO<500.0) && + (nTrks>=6 && nTrks<=45)) { + aida.cloud1D("Select M2vis from PFO").fill(M2visPFO); + aida.cloud1D("#Jets from PFO in PFO Selected Events").fill(nJets); + if (nJets<4) aida.cloud1D("Select M2vis from PFO for #jets<4").fill(M2visPFO); + if (nJets<4 && M2visPFO>0.) { + aida.cloud1D("Select Mvisible from PFO for #jets<4").fill(Math.sqrt(M2visPFO)); +// aida.cloud1D("Select Mvisible from PFO for #jets<4 scaled using Evis").fill(Math.sqrt(M2visPFO)*500.0/EvisPFO); + }
for (MCParticle particle : particles) { for (MCParticle mcp2 : particle.getParents()) { int parentid = mcp2.getPDGID(); if (Math.abs(parentid)==25) { int id = particle.getPDGID(); aida.cloud1D("PFO Selected MC H daughter ID").fill(id);
+ if (nJets<4) aida.cloud1D("PFO Selected MC H daughter ID for #jets<4").fill(id); + if (nJets<4 && M2visPFO>0. && Math.abs(id)==5) aida.cloud1D("Select Mvisible from PFO for #jets<4 for H-bb").fill(Math.sqrt(M2visPFO));
} } }
diff -u -r1.1 -r1.2 --- PFOJetFindingDriver.java 11 Sep 2012 23:41:08 -0000 1.1 +++ PFOJetFindingDriver.java 27 Sep 2012 09:00:41 -0000 1.2 @@ -20,7 +20,7 @@
public class PFOJetFindingDriver extends Driver {
- private double yCut = 0.001;
+ private double yCut = 0.02;
protected String inputPfoCollection = "PandoraPFOCollection";
@@ -35,9 +35,11 @@
Map<HepLorentzVector, ReconstructedParticle> map = new HashMap<HepLorentzVector, ReconstructedParticle>(); for (ReconstructedParticle pfo : PFOs) {
- if (pfo.getCharge() != 0) {
+// if (pfo.getCharge() != 0) {
+ if (Math.abs(VecOp.cosTheta(pfo.getMomentum()))<0.95) {
map.put(pfo.asFourVector(), pfo);
- }
+ }
+// }
} JetFinder jetFind = new JadeEJetFinder(yCut); List<List<ReconstructedParticle>> jetList = new ArrayList<List<ReconstructedParticle>>();
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