Author: [log in to unmask]
Date: Mon Jun 6 14:46:23 2016
New Revision: 4395
Log:
dqm beam energy can be set in the lcsim file
Modified:
java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/DataQualityMonitor.java
java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/EcalMonitoring.java
java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/FinalStateMonitoring.java
java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MollerMonitoring.java
java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MuonCandidateMonitoring.java
java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TridentMonitoring.java
java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/V0Monitoring.java
Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/DataQualityMonitor.java
=============================================================================
--- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/DataQualityMonitor.java (original)
+++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/DataQualityMonitor.java Mon Jun 6 14:46:23 2016
@@ -7,10 +7,13 @@
import java.util.Map;
import java.util.logging.Level;
import java.util.logging.Logger;
+
+import org.hps.conditions.beam.BeamEnergy.BeamEnergyCollection;
import org.hps.record.triggerbank.AbstractIntData;
import org.hps.record.triggerbank.TIData;
import org.lcsim.event.EventHeader;
import org.lcsim.event.GenericObject;
+import org.lcsim.geometry.Detector;
import org.lcsim.util.Driver;
import org.lcsim.util.aida.AIDA;
@@ -25,6 +28,14 @@
private static final Logger LOGGER = Logger.getLogger(DataQualityMonitor.class.getPackage().getName());
+ protected Double beamEnergy;
+ public void setBeamEnergy(double e){
+ this.beamEnergy = e;
+ }
+ public double getBeamEnergy(){
+ return this.beamEnergy;
+ }
+
protected AIDA aida = AIDA.defaultInstance();
protected DQMDatabaseManager manager;
protected String recoVersion = "v0.0";
@@ -37,6 +48,18 @@
protected boolean outputPlots = false;
protected String outputPlotDir = "DQMOutputPlots/";
+ @Override
+ protected void detectorChanged(Detector detector){
+ BeamEnergyCollection beamEnergyCollection =
+ this.getConditionsManager().getCachedConditions(BeamEnergyCollection.class, "beam_energies").getCachedData();
+ if(beamEnergyCollection != null && beamEnergyCollection.size() != 0)
+ beamEnergy = beamEnergyCollection.get(0).getBeamEnergy();
+ else
+ LOGGER.log(Level.WARNING, "warning: beam energy not found");
+
+ }
+
+
String triggerType = "all";//allowed types are "" (blank) or "all", singles0, singles1, pairs0,pairs1
public boolean isGBL = false;
Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/EcalMonitoring.java
=============================================================================
--- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/EcalMonitoring.java (original)
+++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/EcalMonitoring.java Mon Jun 6 14:46:23 2016
@@ -8,7 +8,6 @@
import java.util.logging.Logger;
import org.apache.commons.math.stat.StatUtils;
-import org.hps.conditions.beam.BeamEnergy.BeamEnergyCollection;
import org.hps.recon.ecal.cluster.ClusterUtilities;
import org.lcsim.event.CalorimeterHit;
import org.lcsim.event.Cluster;
@@ -102,9 +101,7 @@
@Override
protected void detectorChanged(Detector detector) {
- BeamEnergyCollection beamEnergyCollection =
- this.getConditionsManager().getCachedConditions(BeamEnergyCollection.class, "beam_energies").getCachedData();
- double beamEnergy = beamEnergyCollection.get(0).getBeamEnergy();
+ /*tab*/super.detectorChanged(detector);
//this.getConditionsManager().getCachedConditions(org.hps.conditions.EcalChannelCollection.class, "ecal_channels").
LOGGER.info("EcalMonitoring::detectorChanged Setting up the plotter");
aida.tree().cd("/");
@@ -126,7 +123,7 @@
twoclusterTotEnergy = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Two Cluster Energy Sum", 100, 0, beamEnergy*maxFactor);
twoclusterEnergyAsymmetry = aida.histogram1D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Two Cluster Energy Asymmetry", 100, 0, beamEnergy*maxFactor);
energyVsT = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Energy vs time", 400, 0.0, 200.0, 100, -0.1, beamEnergy*maxFactor);
- xVsY = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " X vs Y (NHits >1)", 200, -200.0, 200.0, 85, -85.0, 85.0);
+ xVsY = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " X vs Y (NHits >1)", 320, -270.0, 370.0, 90, -90.0, 90.0);
energyVsX = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Energy vs X", 50, 0, maxFactor*beamEnergy, 50, .0, 200.0);
energyVsY = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + " Energy vs Y", 50, 0, maxFactor*beamEnergy, 50, 20.0, 85.0);
pairsE1vsE2 = aida.histogram2D(plotClustersDir + triggerType + "/"+clusterCollectionName + "Pair E1 vs E2", 50, 0, beamEnergy*maxFactor, 50, 0, beamEnergy*maxFactor);
Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/FinalStateMonitoring.java
=============================================================================
--- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/FinalStateMonitoring.java (original)
+++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/FinalStateMonitoring.java Mon Jun 6 14:46:23 2016
@@ -17,7 +17,6 @@
import java.util.logging.Level;
import java.util.logging.Logger;
-import org.hps.conditions.beam.BeamEnergy.BeamEnergyCollection;
import org.hps.recon.ecal.cluster.ClusterUtilities;
import org.hps.recon.tracking.TrackType;
import org.hps.recon.tracking.TrackUtils;
@@ -104,14 +103,7 @@
@Override
protected void detectorChanged(Detector detector) {
- double beamEnergy = 6.6; //maximum possible beam energy is used as a default if the beam energy is unknown
- BeamEnergyCollection beamEnergyCollection =
- this.getConditionsManager().getCachedConditions(BeamEnergyCollection.class, "beam_energies").getCachedData();
- if(beamEnergyCollection != null && beamEnergyCollection.size() != 0)
- beamEnergy = beamEnergyCollection.get(0).getBeamEnergy();
- else
- LOGGER.log(Level.WARNING, "warning: beam energy not found. Assuming it is 6.6 GeV");
-
+ super.detectorChanged(detector);
double maxFactor = 1.5;
double feeMomentumCut = 0.75; //this number, multiplied by the beam energy, is the actual cut
Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MollerMonitoring.java
=============================================================================
--- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MollerMonitoring.java (original)
+++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MollerMonitoring.java Mon Jun 6 14:46:23 2016
@@ -9,7 +9,6 @@
import java.util.EnumSet;
import java.util.List;
import java.util.logging.Logger;
-import org.hps.conditions.beam.BeamEnergy.BeamEnergyCollection;
import org.hps.recon.particle.HpsReconParticleDriver;
import org.hps.recon.particle.ReconParticleDriver;
import org.hps.recon.tracking.TrackType;
@@ -224,16 +223,14 @@
this.beamPos[2] = beamPosY;
}
- double ebeam;
@Override
protected void detectorChanged(Detector detector) {
- LOGGER.info("MollerMonitoring::detectorChanged Setting up the plotter");
+ super.detectorChanged(detector);
+
+ /*tab*/LOGGER.info("MollerMonitoring::detectorChanged Setting up the plotter");
beamAxisRotation.setActiveEuler(Math.PI / 2, -0.0305, -Math.PI / 2);
- BeamEnergyCollection beamEnergyCollection
- = this.getConditionsManager().getCachedConditions(BeamEnergyCollection.class, "beam_energies").getCachedData();
- ebeam = beamEnergyCollection.get(0).getBeamEnergy();
aida.tree().cd("/");
String trkType = "SeedTrack/";
if (isGBL) {
@@ -258,16 +255,16 @@
// triTrackTimeDiff = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Track time difference", 100, -10, 10);
// triTrackTime2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Track time vs. track time", 100, -10, 10, 100, -10, 10);
- triTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Bottom vs. top momentum", 100, 0, v0PzMax * ebeam, 100, 0, v0PzMax * ebeam);
- triDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Moller: Bottom - top momentum", 100, -ebeam, ebeam);
- triSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Moller: Bottom + top momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- triPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Py(t) vs Py(b)", 50, -2 * v0PyMax * ebeam, 2 * v0PyMax * ebeam, 50, -2 * v0PyMax * ebeam, 2 * v0PyMax * ebeam);
- triPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Px(t) vs Px(b)", 50, -2 * v0PxMax * ebeam, 2 * v0PxMax * ebeam, 50, -2 * v0PxMax * ebeam, 2 * v0PxMax * ebeam);
-
- triMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Vertex mass vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, 0, 0.1 * ebeam);
- triZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Vertex Z vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, -v0VzMax, v0VzMax);
- triMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Moller: Vertex mass", 100, 0, 0.1 * ebeam);
- triZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Vertex Z vs. mass", 100, 0, 0.1 * ebeam, 100, -v0VzMax, v0VzMax);
+ triTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Bottom vs. top momentum", 100, 0, v0PzMax * beamEnergy, 100, 0, v0PzMax * beamEnergy);
+ triDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Moller: Bottom - top momentum", 100, -beamEnergy, beamEnergy);
+ triSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Moller: Bottom + top momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ triPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Py(t) vs Py(b)", 50, -2 * v0PyMax * beamEnergy, 2 * v0PyMax * beamEnergy, 50, -2 * v0PyMax * beamEnergy, 2 * v0PyMax * beamEnergy);
+ triPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Px(t) vs Px(b)", 50, -2 * v0PxMax * beamEnergy, 2 * v0PxMax * beamEnergy, 50, -2 * v0PxMax * beamEnergy, 2 * v0PxMax * beamEnergy);
+
+ triMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Vertex mass vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, 0, 0.1 * beamEnergy);
+ triZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Vertex Z vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, -v0VzMax, v0VzMax);
+ triMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Moller: Vertex mass", 100, 0, 0.1 * beamEnergy);
+ triZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Moller: Vertex Z vs. mass", 100, 0, 0.1 * beamEnergy, 100, -v0VzMax, v0VzMax);
// triX = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex X", 100, -v0VxMax, v0VxMax);
// triY = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Y", 100, -v0VyMax, v0VyMax);
// triZ = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Z", 100, -v0VzMax, v0VzMax);
@@ -282,16 +279,16 @@
// vertTrackTimeDiff = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Track time difference", 100, -10, 10);
// vertTrackTime2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Track time vs. track time", 100, -10, 10, 100, -10, 10);
- vertTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Bottom vs. top momentum", 100, 0, v0PzMax * ebeam, 100, 0, v0PzMax * ebeam);
- vertDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Bottom - top momentum", 100, -ebeam, ebeam);
- vertSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Bottom + top momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- vertPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Py(t) vs Py(b)", 50, -2 * v0PyMax * ebeam, 2 * v0PyMax * ebeam, 50, -2 * v0PyMax * ebeam, 2 * v0PyMax * ebeam);
- vertPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Px(t) vs Px(b)", 50, -2 * v0PxMax * ebeam, 2 * v0PxMax * ebeam, 50, -2 * v0PxMax * ebeam, 2 * v0PxMax * ebeam);
-
- vertMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, 0, 0.1 * ebeam);
- vertZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, -v0VzMax, v0VzMax);
- vertMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass", 100, 0, 0.1 * ebeam);
- vertZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. mass", 100, 0, 0.1 * ebeam, 100, -v0VzMax, v0VzMax);
+ vertTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Bottom vs. top momentum", 100, 0, v0PzMax * beamEnergy, 100, 0, v0PzMax * beamEnergy);
+ vertDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Bottom - top momentum", 100, -beamEnergy, beamEnergy);
+ vertSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Bottom + top momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ vertPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Py(t) vs Py(b)", 50, -2 * v0PyMax * beamEnergy, 2 * v0PyMax * beamEnergy, 50, -2 * v0PyMax * beamEnergy, 2 * v0PyMax * beamEnergy);
+ vertPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Px(t) vs Px(b)", 50, -2 * v0PxMax * beamEnergy, 2 * v0PxMax * beamEnergy, 50, -2 * v0PxMax * beamEnergy, 2 * v0PxMax * beamEnergy);
+
+ vertMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, 0, 0.1 * beamEnergy);
+ vertZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, -v0VzMax, v0VzMax);
+ vertMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass", 100, 0, 0.1 * beamEnergy);
+ vertZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. mass", 100, 0, 0.1 * beamEnergy, 100, -v0VzMax, v0VzMax);
vertX = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex X", 100, -v0VxMax, v0VxMax);
vertY = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Y", 100, -v0VyMax, v0VyMax);
// vertZ = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z", 100, -v0VzMax, v0VzMax);
@@ -427,7 +424,7 @@
bits.add(Cut.VTX_QUALITY);
}
- boolean vertexMomentumCut = v0MomRot.z() < v0PzMaxCut * ebeam && v0MomRot.z() > v0PzMinCut * ebeam && Math.abs(v0MomRot.x()) < v0PxCut * ebeam && Math.abs(v0MomRot.y()) < v0PyCut * ebeam;
+ boolean vertexMomentumCut = v0MomRot.z() < v0PzMaxCut * beamEnergy && v0MomRot.z() > v0PzMinCut * beamEnergy && Math.abs(v0MomRot.x()) < v0PxCut * beamEnergy && Math.abs(v0MomRot.y()) < v0PyCut * beamEnergy;
boolean vertexPositionCut = Math.abs(v0Vtx.x()) < v0UnconVxCut && Math.abs(v0Vtx.y()) < v0UnconVyCut && Math.abs(v0Vtx.z()) < v0UnconVzCut && Math.abs(bscVtx.x()) < v0BsconVxCut && Math.abs(bscVtx.y()) < v0BsconVyCut;
if (vertexMomentumCut && vertexPositionCut) {
bits.add(Cut.VERTEX_CUTS);
@@ -439,10 +436,10 @@
}
// boolean topBottomCut = top.getMomentum().y() * bot.getMomentum().y() < 0;
-// boolean pMinCut = top.getMomentum().magnitude() > minTrkPCut * ebeam && bot.getMomentum().magnitude() > minTrkPCut * ebeam;
- boolean pMaxCut = top.getMomentum().magnitude() < maxTrkPCut * ebeam && bot.getMomentum().magnitude() < maxTrkPCut * ebeam;
-// boolean pSumCut = top.getMomentum().magnitude() + bot.getMomentum().magnitude() < maxPSumCut * ebeam && top.getMomentum().magnitude() + bot.getMomentum().magnitude() < minPSumCut * ebeam;
-// boolean pMaxCut = top.getMomentum().magnitude() < maxTrkPCut * ebeam && bot.getMomentum().magnitude() < maxTrkPCut * ebeam;
+// boolean pMinCut = top.getMomentum().magnitude() > minTrkPCut * beamEnergy && bot.getMomentum().magnitude() > minTrkPCut * beamEnergy;
+ boolean pMaxCut = top.getMomentum().magnitude() < maxTrkPCut * beamEnergy && bot.getMomentum().magnitude() < maxTrkPCut * beamEnergy;
+// boolean pSumCut = top.getMomentum().magnitude() + bot.getMomentum().magnitude() < maxPSumCut * beamEnergy && top.getMomentum().magnitude() + bot.getMomentum().magnitude() < minPSumCut * beamEnergy;
+// boolean pMaxCut = top.getMomentum().magnitude() < maxTrkPCut * beamEnergy && bot.getMomentum().magnitude() < maxTrkPCut * beamEnergy;
if (pMaxCut) {
bits.add(Cut.TRACK_CUTS);
}
Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MuonCandidateMonitoring.java
=============================================================================
--- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MuonCandidateMonitoring.java (original)
+++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/MuonCandidateMonitoring.java Mon Jun 6 14:46:23 2016
@@ -7,7 +7,7 @@
import java.util.List;
import java.util.logging.Logger;
-import org.hps.conditions.beam.BeamEnergy.BeamEnergyCollection;
+// import org.hps.UnusedImportCheckstyleViolation
import org.hps.recon.tracking.TrackType;
import org.hps.recon.tracking.TrackUtils;
import org.lcsim.event.Cluster;
@@ -27,24 +27,15 @@
private String beamConV0CandidatesColName = "BeamspotConstrainedV0Candidates";
private String targetV0ConCandidatesColName = "TargetConstrainedV0Candidates";
private String clusterCollectionName = "EcalClustersCorr";
-
+
private int nRecoEvents = 0;
-
+
private String plotDir = "MuonCandidates/";
-
+
@Override
public void process(EventHeader event) {
/* make sure everything is there */
- if (!event.hasCollection(ReconstructedParticle.class, finalStateParticlesColName))
- return;
- if (!event.hasCollection(ReconstructedParticle.class, unconstrainedV0CandidatesColName))
- return;
- if (!event.hasCollection(ReconstructedParticle.class, beamConV0CandidatesColName))
- return;
- if (!event.hasCollection(ReconstructedParticle.class, targetV0ConCandidatesColName))
- return;
- if (!event.hasCollection(Cluster.class, clusterCollectionName))
- return;
+
//check to see if this event is from the correct trigger (or "all");
if (!matchTrigger(event))
@@ -54,10 +45,10 @@
RelationalTable hitToStrips = TrackUtils.getHitToStripsTable(event);
RelationalTable hitToRotated = TrackUtils.getHitToRotatedTable(event);
-
- List<ReconstructedParticle> unonstrainedV0List = event.get(ReconstructedParticle.class, unconstrainedV0CandidatesColName);
+
+ List<ReconstructedParticle> v0List = event.get(ReconstructedParticle.class, unconstrainedV0CandidatesColName);
List<Cluster> clusters = event.get(Cluster.class, clusterCollectionName);
- for (ReconstructedParticle uncV0 : unonstrainedV0List) {
+ for (ReconstructedParticle uncV0 : v0List) {
if (isGBL != TrackType.isGBL(uncV0.getType()))
continue;
Vertex uncVert = uncV0.getStartVertex();
@@ -88,89 +79,132 @@
muMinus = trks.get(0);
}
- if(uncVert.getAssociatedParticle().getEnergy() > .080) //conservative guess for twice minimum ionizing
- continue;
-
+
+
+
+
Hep3Vector pPlus = muPlus.getMomentum();
Hep3Vector pMinus = muMinus.getMomentum();
-
-
- double pp =pPlus.magnitude();
+
+ //cut out pairs that have both particles on the same half of SVT.
+ if(pPlus.y()*pMinus.y()>0){
+ continue;
+ }
+
+
+ double pp = pPlus.magnitude();
double pm = pMinus.magnitude();
-
+
double psum = pp + pm;
double pdiff = pp - pm;
if(psum > v0PSumMaxCut || psum < v0PSumMinCut)
continue;
if(pdiff > maxDiffCut || pdiff < -maxDiffCut)
continue;
- if(pp < minP|| pm < minP)
- continue;
-
-
+ if(pp < minP || pm < minP)
+ continue;
+ if(pp > maxP || pm > maxP)
+ continue;
+
+
Track trackPlus = muPlus.getTracks().get(0);
Track trackMinus = muMinus.getTracks().get(0);
-
+
double tPlus = TrackUtils.getTrackTime(trackPlus, hitToStrips, hitToRotated);
double tMinus = TrackUtils.getTrackTime(trackMinus, hitToStrips, hitToRotated);
-
-
-
-
- Hep3Vector xyPlus = TrackUtils.extrapolateTrack(trackPlus, ecalZPosition);
- Hep3Vector xyMinus = TrackUtils.extrapolateTrack(trackMinus, ecalZPosition);
-
-
-
-
+
+
+
+ {
+ double ecalZPosition = 1395;
+ Hep3Vector xyPlus = TrackUtils.extrapolateTrack(trackPlus, ecalZPosition);
+ Hep3Vector xyMinus = TrackUtils.extrapolateTrack(trackMinus, ecalZPosition);
+
+ xyAtEcalPlus.fill(xyPlus.x(), xyPlus.y());
+ xyAtEcalMinus.fill(xyMinus.x(), xyMinus.y());
+ if(!fid_ECal(xyPlus.x(), xyPlus.y()))
+ continue;
+ if(!fid_ECal(xyMinus.x(), xyMinus.y()))
+ continue;
+ }
+
+
//distance between the extrapolated track and the nearest cluster
double iso1 = 1000, iso2 = 1000;
-
- Cluster nc = getNearestCluster(clusters, xyPlus);
+
+ Cluster nc = getNearestCluster(clusters, trackPlus);
if(nc != null){
+ Hep3Vector xyPlus = TrackUtils.extrapolateTrack(trackPlus, nc.getPosition()[2]);
iso1 = Math.hypot(xyPlus.x()-nc.getPosition()[0], xyPlus.y()-nc.getPosition()[1]);
dxyNearestClusterPlus.fill(xyPlus.x()-nc.getPosition()[0], xyPlus.y()-nc.getPosition()[1]);
}
- nc = getNearestCluster(clusters, xyMinus);
+ nc = getNearestCluster(clusters, trackMinus);
if(nc != null){
+ Hep3Vector xyMinus = TrackUtils.extrapolateTrack(trackMinus, nc.getPosition()[2]);
iso2 = Math.hypot(xyMinus.x()-nc.getPosition()[0], xyMinus.y()-nc.getPosition()[1]);
dxyNearestClusterMinus.fill(xyMinus.x()-nc.getPosition()[0], xyMinus.y()-nc.getPosition()[1]);
}
-
+
if(iso1 < 20 || iso2 < 20)
continue;
-
+
tPlusVsTMinus.fill(tPlus, tMinus);
timeDiff.fill(tPlus- tMinus);
-
-
- xyAtEcalPlus.fill(xyPlus.x(), xyPlus.y());
- xyAtEcalMinus.fill(xyMinus.x(), xyMinus.y());
-
-
+
+
+
+
double mmu2 = .1057*.1057;
-
+
double mass = Math.sqrt(2*mmu2 + 2*Math.sqrt(pPlus.magnitudeSquared()+mmu2)*Math.sqrt(pMinus.magnitudeSquared()+mmu2)
-2*(pPlus.x()*pMinus.x()+pPlus.y()*pMinus.y()+pPlus.z()*pMinus.z()));
this.mass.fill(mass);
-
-
-
+
+ sumPxPy.fill(pPlus.x()+pMinus.x(), pPlus.y()+pMinus.y());
+
this.pPlus.fill(pp);
this.pMinus.fill(pm);
this.pPlusVsPMinus.fill(pp, pm);
this.pTot.fill(pm+pp);
-
+
for(Cluster c : clusters){
clusterXY.fill(c.getPosition()[0], c.getPosition()[1]);
}
}
}
- private Cluster getNearestCluster(List<Cluster> clusters, Hep3Vector xy) {
+ public static boolean fid_ECal(double x, double y)
+ {
+ y = Math.abs(y);
+
+ boolean in_fid = false;
+ double x_edge_low = -262.74;
+ double x_edge_high = 347.7;
+ double y_edge_low = 33.54;
+ double y_edge_high = 75.18;
+
+ double x_gap_low = -106.66;
+ double x_gap_high = 42.17;
+ double y_gap_high = 47.18;
+
+ y = Math.abs(y);
+
+ if( x > x_edge_low && x < x_edge_high && y > y_edge_low && y < y_edge_high )
+ {
+ if( !(x > x_gap_low && x < x_gap_high && y > y_edge_low && y < y_gap_high) )
+ {
+ in_fid = true;
+ }
+ }
+
+ return in_fid;
+ }
+
+ private Cluster getNearestCluster(List<Cluster> clusters, Track t) {
double rbest = 1000;
Cluster best = null;
for(Cluster c : clusters){
+ Hep3Vector xy = TrackUtils.extrapolateTrack(t, c.getPosition()[2]);
double r = Math.hypot(xy.x()-c.getPosition()[0], xy.y()-c.getPosition()[1]);
if(r<rbest){
rbest = r;
@@ -180,30 +214,28 @@
return best;
}
- private double ecalZPosition = 1300;
- private double v0PSumMinCut, v0PSumMaxCut, v0MaxPCut, beambeamCut, maxDiffCut;
+ private double v0PSumMinCut, v0PSumMaxCut, maxDiffCut;
private double minP;
-
+ private double maxP;
+
@Override
protected void detectorChanged(Detector detector) {
- BeamEnergyCollection beamEnergyCollection =
- this.getConditionsManager().getCachedConditions(BeamEnergyCollection.class, "beam_energies").getCachedData();
- double beamEnergy = beamEnergyCollection.get(0).getBeamEnergy();
- //feeMomentumCut = 0.75*beamEnergy; //GeV
+ super.detectorChanged(detector);
+ /*tab*///feeMomentumCut = 0.75*beamEnergy; //GeV
v0PSumMinCut = 0.2 * beamEnergy;
v0PSumMaxCut = 1.25 * beamEnergy;
minP = .1 *beamEnergy;
+ maxP = .85*beamEnergy;
maxDiffCut = .8*beamEnergy;
-
+
LOGGER.info("Setting up the plotter");
aida.tree().cd("/");
-
- String xtra = "Extras";
+
String trkType = "SeedTrack/";
if (isGBL)
trkType = "GBLTrack/";
@@ -215,37 +247,40 @@
mass = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Invariant Mass (GeV)", 100, 0, maxMass);
pPlus = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "P mu+", 100, 0, 1.5*beamEnergy);
pMinus = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "P mu-", 100, 0, 1.5*beamEnergy);
-
+
timeDiff = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Time Difference", 100, -10, 10);
tPlusVsTMinus = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "t+ versus t-", 100, -10, 10, 100, -10, 10);
-
+
pPlusVsPMinus = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "P mu+ vs P mu-", 100, 0, 1.5*beamEnergy, 100, 0, 1.5*beamEnergy);
-
-
-
-
+
+
+
+
pTot = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "P tot", 100, 0, 2*beamEnergy);
-
+
xyAtEcalMinus = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "XY at Ecal (mu-)", 200, -200.0, 200.0, 85, -85.0, 85.0);
xyAtEcalPlus = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "XY at Ecal (mu+)", 200, -200.0, 200.0, 85, -85.0, 85.0);
-
+
clusterXY = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Cluster XY", 200, -200.0, 200.0, 85, -85.0, 85.0);
-
-
+
+
dxyNearestClusterPlus = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "dXY from nearest cluster (mu+)", 200, -30.0, 30.0, 85, -30.0, 30.0);
dxyNearestClusterMinus = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "dXY from nearest cluster (mu-)", 200, -30.0, 30.0, 85, -30.0, 30.0);
- }
-
-
-
+ sumPxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "pair px vs py", 100, -.1, .1, 100, -.1, .1);
+ }
+
+
+
private IHistogram1D mass, pPlus, pMinus, pTot, timeDiff;
-
+
private IHistogram2D xyAtEcalPlus, xyAtEcalMinus, pPlusVsPMinus, tPlusVsTMinus;
-
+
private IHistogram2D dxyNearestClusterPlus, dxyNearestClusterMinus;
-
+
private IHistogram2D clusterXY;
-
+
+ private IHistogram2D sumPxPy;
+
}
Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TridentMonitoring.java
=============================================================================
--- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TridentMonitoring.java (original)
+++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/TridentMonitoring.java Mon Jun 6 14:46:23 2016
@@ -14,7 +14,7 @@
import java.util.List;
import java.util.Map.Entry;
import java.util.logging.Logger;
-import org.hps.conditions.beam.BeamEnergy.BeamEnergyCollection;
+// import org.hps.UnusedImportCheckstyleViolation
import org.hps.recon.ecal.cluster.ClusterUtilities;
import org.hps.recon.particle.HpsReconParticleDriver;
import org.hps.recon.particle.ReconParticleDriver;
@@ -327,22 +327,15 @@
this.beamPos[2] = beamPosY;
}
- public void setEbeam(double ebeam) {
- this.ebeam = ebeam;
- }
-
- private double ebeam = Double.NaN;
+
@Override
protected void detectorChanged(Detector detector) {
- LOGGER.info("TridendMonitoring::detectorChanged Setting up the plotter");
+ super.detectorChanged(detector);
+ /*tab*/LOGGER.info("TridendMonitoring::detectorChanged Setting up the plotter");
beamAxisRotation.setActiveEuler(Math.PI / 2, -0.0305, -Math.PI / 2);
- BeamEnergyCollection beamEnergyCollection
- = this.getConditionsManager().getCachedConditions(BeamEnergyCollection.class, "beam_energies").getCachedData();
- if (Double.isNaN(ebeam)) {
- ebeam = beamEnergyCollection.get(0).getBeamEnergy();
- }
+
aida.tree().cd("/");
String trkType = "SeedTrack/";
if (isGBL) {
@@ -367,16 +360,16 @@
// triTrackTimeDiff = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Track time difference", 100, -10, 10);
// triTrackTime2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Track time vs. track time", 100, -10, 10, 100, -10, 10);
- triTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Positron vs. electron momentum", 100, 0, v0PzMax * ebeam, 100, 0, v0PzMax * ebeam);
- triDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Positron - electron momentum", 100, -ebeam, ebeam);
- triSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Positron + electron momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- triPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Py(e) vs Py(p)", 50, -v0PyMax * ebeam, v0PyMax * ebeam, 50, -v0PyMax * ebeam, v0PyMax * ebeam);
- triPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Px(e) vs Px(p)", 50, -v0PxMax * ebeam, v0PxMax * ebeam, 50, -v0PxMax * ebeam, v0PxMax * ebeam);
-
- triMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Vertex mass vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, 0, 0.1 * ebeam);
- triZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Z vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, -v0VzMax, v0VzMax);
- triMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex mass", 100, 0, 0.1 * ebeam);
- triZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Z vs. mass", 100, 0, 0.1 * ebeam, 100, -v0VzMax, v0VzMax);
+ triTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Positron vs. electron momentum", 100, 0, v0PzMax * beamEnergy, 100, 0, v0PzMax * beamEnergy);
+ triDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Positron - electron momentum", 100, -beamEnergy, beamEnergy);
+ triSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Positron + electron momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ triPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Py(e) vs Py(p)", 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy, 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
+ triPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Px(e) vs Px(p)", 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy, 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy);
+
+ triMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Vertex mass vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, 0, 0.1 * beamEnergy);
+ triZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Z vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, -v0VzMax, v0VzMax);
+ triMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex mass", 100, 0, 0.1 * beamEnergy);
+ triZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Z vs. mass", 100, 0, 0.1 * beamEnergy, 100, -v0VzMax, v0VzMax);
// triX = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex X", 100, -v0VxMax, v0VxMax);
// triY = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Y", 100, -v0VyMax, v0VyMax);
// triZ = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Trident: Vertex Z", 100, -v0VzMax, v0VzMax);
@@ -392,40 +385,40 @@
triRadTrackTimeDiff = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Track time difference", 100, -10, 10);
triRadTrackTime2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Track time vs. track time", 100, -10, 10, 100, -10, 10);
- triRadTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Positron vs. electron momentum", 100, 0, v0PzMax * ebeam, 100, 0, v0PzMax * ebeam);
- triRadDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Positron - electron momentum", 100, -ebeam, ebeam);
- triRadSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Positron + electron momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- triRadPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Py(e) vs Py(p)", 50, -v0PyMax * ebeam, v0PyMax * ebeam, 50, -v0PyMax * ebeam, v0PyMax * ebeam);
- triRadPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Px(e) vs Px(p)", 50, -v0PxMax * ebeam, v0PxMax * ebeam, 50, -v0PxMax * ebeam, v0PxMax * ebeam);
+ triRadTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Positron vs. electron momentum", 100, 0, v0PzMax * beamEnergy, 100, 0, v0PzMax * beamEnergy);
+ triRadDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Positron - electron momentum", 100, -beamEnergy, beamEnergy);
+ triRadSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Positron + electron momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ triRadPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Py(e) vs Py(p)", 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy, 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
+ triRadPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Px(e) vs Px(p)", 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy, 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy);
// triRadMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex mass vs. vertex momentum", 100, v0PzMin, v0PzMax, 100, 0, 0.1);
// triRadZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Z vs. vertex momentum", 100, v0PzMin, v0PzMax, 100, -v0VzMax, v0VzMax);
- triRadMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex mass", 100, 0, 0.1 * ebeam);
- triRadZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Z vs. mass", 100, 0, 0.1 * ebeam, 100, -v0VzMax, v0VzMax);
+ triRadMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex mass", 100, 0, 0.1 * beamEnergy);
+ triRadZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Z vs. mass", 100, 0, 0.1 * beamEnergy, 100, -v0VzMax, v0VzMax);
// triRadX = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex X", 100, -v0VxMax, v0VxMax);
// triRadY = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Y", 100, -v0VyMax, v0VyMax);
// triRadZ = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Z", 100, -v0VzMax, v0VzMax);
// triRadXY = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Y vs. X", 100, -v0VxMax, v0VxMax, 100, -v0VyMax, v0VyMax);
// triRadZY = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Z vs. Y", 100, -v0VyMax, v0VyMax, 100, -v0VzMax, v0VzMax);
- triRadPx = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Px", 100, -v0PxMax * ebeam, v0PxMax * ebeam);
- triRadPy = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Py", 100, -v0PyMax * ebeam, v0PyMax * ebeam);
- triRadPz = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Pz", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- triRadPxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Py vs. Px", 100, -v0PxMax * ebeam, v0PxMax * ebeam, 100, -v0PyMax * ebeam, v0PyMax * ebeam);
+ triRadPx = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Px", 100, -v0PxMax * beamEnergy, v0PxMax * beamEnergy);
+ triRadPy = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Py", 100, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
+ triRadPz = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Pz", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ triRadPxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Py vs. Px", 100, -v0PxMax * beamEnergy, v0PxMax * beamEnergy, 100, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
triRadU = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Px over Ptot", 100, -0.1, 0.1);
triRadV = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative trident: Vertex Py over Ptot", 100, -0.1, 0.1);
// vertTrackTimeDiff = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Track time difference", 100, -10, 10);
// vertTrackTime2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Track time vs. track time", 100, -10, 10, 100, -10, 10);
- vertTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Positron vs. electron momentum", 100, 0, v0PzMax * ebeam, 100, 0, v0PzMax * ebeam);
- vertDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Positron - electron momentum", 100, -ebeam, ebeam);
- vertSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Positron + electron momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- vertPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Py(e) vs Py(p)", 50, -v0PyMax * ebeam, v0PyMax * ebeam, 50, -v0PyMax * ebeam, v0PyMax * ebeam);
- vertPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Px(e) vs Px(p)", 50, -v0PxMax * ebeam, v0PxMax * ebeam, 50, -v0PxMax * ebeam, v0PxMax * ebeam);
-
- vertMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, 0, 0.1 * ebeam);
- vertZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, -v0VzMax, v0VzMax);
- vertMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass", 100, 0, 0.1 * ebeam);
- vertZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. mass", 100, 0, 0.1 * ebeam, 100, -v0VzMax, v0VzMax);
+ vertTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Positron vs. electron momentum", 100, 0, v0PzMax * beamEnergy, 100, 0, v0PzMax * beamEnergy);
+ vertDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Positron - electron momentum", 100, -beamEnergy, beamEnergy);
+ vertSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Positron + electron momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ vertPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Py(e) vs Py(p)", 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy, 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
+ vertPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Px(e) vs Px(p)", 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy, 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy);
+
+ vertMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, 0, 0.1 * beamEnergy);
+ vertZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, -v0VzMax, v0VzMax);
+ vertMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex mass", 100, 0, 0.1 * beamEnergy);
+ vertZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z vs. mass", 100, 0, 0.1 * beamEnergy, 100, -v0VzMax, v0VzMax);
// vertX = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex X", 100, -v0VxMax, v0VxMax);
vertY = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Y", 100, -v0VyMax, v0VyMax);
// vertZ = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Vertex: Vertex Z", 100, -v0VzMax, v0VzMax);
@@ -441,25 +434,25 @@
vertRadTrackTimeDiff = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Track time difference", 100, -10, 10);
vertRadTrackTime2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Track time vs. track time", 100, -10, 10, 100, -10, 10);
- vertRadTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Positron vs. electron momentum", 100, 0, v0PzMax * ebeam, 100, 0, v0PzMax * ebeam);
- vertRadDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Positron - electron momentum", 100, -ebeam, ebeam);
- vertRadSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Positron + electron momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- vertRadPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Py(e) vs Py(p)", 50, -v0PyMax * ebeam, v0PyMax * ebeam, 50, -v0PyMax * ebeam, v0PyMax * ebeam);
- vertRadPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Px(e) vs Px(p)", 50, -v0PxMax * ebeam, v0PxMax * ebeam, 50, -v0PxMax * ebeam, v0PxMax * ebeam);
-
- vertRadMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex mass vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, 0, 0.1 * ebeam);
- vertRadZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Z vs. vertex momentum", 100, v0PzMin * ebeam, v0PzMax * ebeam, 100, -v0VzMax, v0VzMax);
- vertRadMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex mass", 100, 0, 0.1 * ebeam);
- vertRadZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Z vs. mass", 100, 0, 0.1 * ebeam, 100, -v0VzMax, v0VzMax);
+ vertRadTrackMomentum2D = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Positron vs. electron momentum", 100, 0, v0PzMax * beamEnergy, 100, 0, v0PzMax * beamEnergy);
+ vertRadDeltaP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Positron - electron momentum", 100, -beamEnergy, beamEnergy);
+ vertRadSumP = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Positron + electron momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ vertRadPyEleVsPyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Py(e) vs Py(p)", 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy, 50, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
+ vertRadPxEleVsPxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Px(e) vs Px(p)", 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy, 50, -v0PxMax * beamEnergy, v0PxMax * beamEnergy);
+
+ vertRadMassMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex mass vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, 0, 0.1 * beamEnergy);
+ vertRadZVsMomentum = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Z vs. vertex momentum", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy, 100, -v0VzMax, v0VzMax);
+ vertRadMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex mass", 100, 0, 0.1 * beamEnergy);
+ vertRadZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Z vs. mass", 100, 0, 0.1 * beamEnergy, 100, -v0VzMax, v0VzMax);
vertRadX = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex X", 100, -v0VxMax, v0VxMax);
vertRadY = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Y", 100, -v0VyMax, v0VyMax);
vertRadZ = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Z", 100, -v0VzMax, v0VzMax);
vertRadXY = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Y vs. X", 100, -v0VxMax, v0VxMax, 100, -v0VyMax, v0VyMax);
vertRadZY = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Z vs. Y", 100, -v0VyMax, v0VyMax, 100, -v0VzMax, v0VzMax);
- vertRadPx = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Px", 100, -v0PxMax * ebeam, v0PxMax * ebeam);
- vertRadPy = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Py", 100, -v0PyMax * ebeam, v0PyMax * ebeam);
- vertRadPz = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Pz", 100, v0PzMin * ebeam, v0PzMax * ebeam);
- vertRadPxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Py vs. Px", 100, -v0PxMax * ebeam, v0PxMax * ebeam, 100, -v0PyMax * ebeam, v0PyMax * ebeam);
+ vertRadPx = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Px", 100, -v0PxMax * beamEnergy, v0PxMax * beamEnergy);
+ vertRadPy = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Py", 100, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
+ vertRadPz = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Pz", 100, v0PzMin * beamEnergy, v0PzMax * beamEnergy);
+ vertRadPxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Py vs. Px", 100, -v0PxMax * beamEnergy, v0PxMax * beamEnergy, 100, -v0PyMax * beamEnergy, v0PyMax * beamEnergy);
vertRadU = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Px over Ptot", 100, -0.1, 0.1);
vertRadV = aida.histogram1D(plotDir + trkType + triggerType + "/" + "Radiative vertex: Vertex Py over Ptot", 100, -0.1, 0.1);
@@ -495,9 +488,9 @@
cutVertexZ[cut.ordinal()][i] = aida.histogram1D(String.format("%s%s%s/failed cut %d: %s/%s: Vertex Z position (mm)", plotDir, trkType, triggerType, cut.ordinal(), cut.name, i == VERTEX ? "vertex" : "trident"),
100, -v0VzMax, v0VzMax);
cutVertexMass[cut.ordinal()][i] = aida.histogram1D(String.format("%s%s%s/failed cut %d: %s/%s: Vertex mass (GeV)", plotDir, trkType, triggerType, cut.ordinal(), cut.name, i == VERTEX ? "vertex" : "trident"),
- 100, 0, 0.1 * ebeam);
+ 100, 0, 0.1 * beamEnergy);
cutVertexZVsMass[cut.ordinal()][i] = aida.histogram2D(String.format("%s%s%s/failed cut %d: %s/%s: Vertex Z vs. mass", plotDir, trkType, triggerType, cut.ordinal(), cut.name, i == VERTEX ? "vertex" : "trident"),
- 100, 0, 0.1 * ebeam, 100, -v0VzMax, v0VzMax);
+ 100, 0, 0.1 * beamEnergy, 100, -v0VzMax, v0VzMax);
}
}
}
@@ -695,7 +688,7 @@
bits.add(Cut.VTX_QUALITY);
}
- boolean vertexMomentumCut = v0MomRot.z() < v0PzMaxCut * ebeam && v0MomRot.z() > v0PzMinCut * ebeam && Math.abs(v0MomRot.x()) < v0PxCut * ebeam && Math.abs(v0MomRot.y()) < v0PyCut * ebeam;
+ boolean vertexMomentumCut = v0MomRot.z() < v0PzMaxCut * beamEnergy && v0MomRot.z() > v0PzMinCut * beamEnergy && Math.abs(v0MomRot.x()) < v0PxCut * beamEnergy && Math.abs(v0MomRot.y()) < v0PyCut * beamEnergy;
boolean vertexPositionCut = Math.abs(v0Vtx.x()) < v0UnconVxCut && Math.abs(v0Vtx.y()) < v0UnconVyCut && Math.abs(v0Vtx.z()) < v0UnconVzCut && Math.abs(bscVtx.x()) < v0BsconVxCut && Math.abs(bscVtx.y()) < v0BsconVyCut;
if (vertexMomentumCut && vertexPositionCut) {
bits.add(Cut.VERTEX_CUTS);
@@ -707,8 +700,8 @@
}
boolean topBottomCut = electron.getMomentum().y() * positron.getMomentum().y() < 0;
- boolean pMinCut = electron.getMomentum().magnitude() > minPCut * ebeam && positron.getMomentum().magnitude() > minPCut * ebeam;
- boolean pMaxCut = electron.getMomentum().magnitude() < beamPCut * ebeam && positron.getMomentum().magnitude() < beamPCut * ebeam;
+ boolean pMinCut = electron.getMomentum().magnitude() > minPCut * beamEnergy && positron.getMomentum().magnitude() > minPCut * beamEnergy;
+ boolean pMaxCut = electron.getMomentum().magnitude() < beamPCut * beamEnergy && positron.getMomentum().magnitude() < beamPCut * beamEnergy;
if (topBottomCut && pMaxCut && pMinCut) {
bits.add(Cut.TRACK_CUTS);
}
@@ -751,7 +744,7 @@
EnumSet<Cut> allButThisCut = EnumSet.allOf(Cut.class);
allButThisCut.remove(cut);
if (bits.containsAll(allButThisCut)) {
- if (uncV0.getMass() > plotsMinMass * ebeam && uncV0.getMass() < plotsMaxMass * ebeam && uncV0.getMomentum().magnitude() > radCut * ebeam) {
+ if (uncV0.getMass() > plotsMinMass * beamEnergy && uncV0.getMass() < plotsMaxMass * beamEnergy && uncV0.getMomentum().magnitude() > radCut * beamEnergy) {
switch (cut) {
case ISOLATION:
l1Iso.fill(minIso);
@@ -784,7 +777,7 @@
}
}
if (!bits.contains(cut)) {
- if (uncV0.getMass() > plotsMinMass * ebeam && uncV0.getMass() < plotsMaxMass * ebeam) {
+ if (uncV0.getMass() > plotsMinMass * beamEnergy && uncV0.getMass() < plotsMaxMass * beamEnergy) {
cutVertexZ[cut.ordinal()][VERTEX].fill(v0Vtx.z());
}
cutVertexMass[cut.ordinal()][VERTEX].fill(uncV0.getMass());
@@ -796,7 +789,7 @@
allTriCutsButThisCut.remove(cut);
if (bits.containsAll(allTriCutsButThisCut)) {
if (!bits.contains(cut)) {
- if (uncV0.getMass() > plotsMinMass * ebeam && uncV0.getMass() < plotsMaxMass * ebeam) {
+ if (uncV0.getMass() > plotsMinMass * beamEnergy && uncV0.getMass() < plotsMaxMass * beamEnergy) {
cutVertexZ[cut.ordinal()][TRIDENT].fill(v0Vtx.z());
}
cutVertexMass[cut.ordinal()][TRIDENT].fill(uncV0.getMass());
@@ -857,7 +850,7 @@
// triV.fill(pBestV0Rot.y() / pBestV0Rot.magnitude());
// triXY.fill(v0Vtx.x(), v0Vtx.y());
// triZY.fill(v0Vtx.y(), v0Vtx.z());
- if (bestCandidate.getMomentum().magnitude() > radCut * ebeam) {
+ if (bestCandidate.getMomentum().magnitude() > radCut * beamEnergy) {
triRadTrackTime2D.fill(tEle, tPos);
triRadTrackTimeDiff.fill(tEle - tPos);
// triRadZVsMomentum.fill(bestCandidate.getMomentum().magnitude(), v0Vtx.z());
@@ -926,7 +919,7 @@
// vertV.fill(pBestV0Rot.y() / pBestV0Rot.magnitude());
vertXY.fill(v0Vtx.x(), v0Vtx.y());
vertZY.fill(v0Vtx.y(), v0Vtx.z());
- if (bestCandidate.getMomentum().magnitude() > radCut * ebeam) {
+ if (bestCandidate.getMomentum().magnitude() > radCut * beamEnergy) {
BilliorVertexer vtxFitter = new BilliorVertexer(TrackUtils.getBField(event.getDetector()).y());
vtxFitter.setBeamSize(beamSize);
Modified: java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/V0Monitoring.java
=============================================================================
--- java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/V0Monitoring.java (original)
+++ java/trunk/analysis/src/main/java/org/hps/analysis/dataquality/V0Monitoring.java Mon Jun 6 14:46:23 2016
@@ -8,8 +8,8 @@
import hep.aida.IHistogram2D;
import hep.aida.IPlotter;
import hep.aida.IPlotterStyle;
+import hep.physics.vec.Hep3Vector;
import hep.physics.vec.BasicHep3Matrix;
-import hep.physics.vec.Hep3Vector;
import hep.physics.vec.VecOp;
import java.io.IOException;
@@ -18,8 +18,7 @@
import java.util.Map.Entry;
import java.util.logging.Level;
import java.util.logging.Logger;
-
-import org.hps.conditions.beam.BeamEnergy.BeamEnergyCollection;
+// import org.hps.UnusedImportCheckstyleViolation
import org.hps.recon.tracking.TrackType;
import org.hps.recon.tracking.TrackUtils;
import org.hps.recon.vertexing.BilliorTrack;
@@ -34,134 +33,529 @@
import org.lcsim.recon.tracking.seedtracker.SeedTrack;
/**
- * DQM driver V0 particles (i.e. e+e- pars) plots things like number of vertex position an mass
+ * DQM driver V0 particles (i.e. e+e- pars) plots things like number of vertex
+ * position an mass
*
* @author mgraham on May 14, 2014
+ *
*/
public class V0Monitoring extends DataQualityMonitor {
private static Logger LOGGER = Logger.getLogger(V0Monitoring.class.getPackage().getName());
- private static boolean hasSharedStrips(final ReconstructedParticle fs1, final ReconstructedParticle fs2,
- final RelationalTable hittostrip, final RelationalTable hittorotated) {
- return TrackUtils.hasSharedStrips(fs1.getTracks().get(0), fs2.getTracks().get(0), hittostrip, hittorotated);
- }
-
- private final BasicHep3Matrix beamAxisRotation = new BasicHep3Matrix();
- private final String beamConV0CandidatesColName = "BeamspotConstrainedV0Candidates";
- private IHistogram1D bsconChi2;
- private IHistogram2D bsconChi2VsTrkChi2;
+ private String finalStateParticlesColName = "FinalStateParticles";
+ private String unconstrainedV0CandidatesColName = "UnconstrainedV0Candidates";
+ private String beamConV0CandidatesColName = "BeamspotConstrainedV0Candidates";
+ private String targetV0ConCandidatesColName = "TargetConstrainedV0Candidates";
+ private String[] fpQuantNames = {"nV0_per_Event", "avg_BSCon_mass", "avg_BSCon_Vx", "avg_BSCon_Vy", "avg_BSCon_Vz", "sig_BSCon_Vx", "sig_BSCon_Vy", "sig_BSCon_Vz", "avg_BSCon_Chi2"};
+ //some counters
+ private int nRecoEvents = 0;
+ private int nTotV0 = 0;
+ private int nTot2Ele = 0;
+ //some summers
+ private double sumMass = 0.0;
+ private double sumVx = 0.0;
+ private double sumVy = 0.0;
+ private double sumVz = 0.0;
+ private double sumChi2 = 0.0;
+
+ /* V0 Quantities */
+ /* Mass, vertex, chi^2 of fit */
+ /* unconstrained */
+ private IHistogram1D unconMass;
+ private IHistogram1D unconVx;
+ private IHistogram1D unconVy;
+ private IHistogram1D unconVz;
+ private IHistogram1D unconChi2;
+ private IHistogram2D unconVzVsChi2;
+ private IHistogram2D unconChi2VsTrkChi2;
+ /* beamspot constrained */
+
+ private IHistogram1D nV0;
+
+ private IHistogram1D v0Time;
+ private IHistogram1D v0Dt;
+ private IHistogram2D trigTimeV0Time;
+ private IHistogram1D trigTime;
+
private IHistogram1D bsconMass;
private IHistogram1D bsconVx;
private IHistogram1D bsconVy;
private IHistogram1D bsconVz;
+ private IHistogram1D bsconChi2;
private IHistogram2D bsconVzVsChi2;
- private final String finalStateParticlesColName = "FinalStateParticles";
- private final String[] fpQuantNames = {"nV0_per_Event", "avg_BSCon_mass", "avg_BSCon_Vx", "avg_BSCon_Vy",
- "avg_BSCon_Vz", "sig_BSCon_Vx", "sig_BSCon_Vy", "sig_BSCon_Vz", "avg_BSCon_Chi2"};
-
- private final double maxFactor = 1.25;
- private IHistogram1D mollerHiP, mollerLoP, mollerEitherP, mollerPsum;
- private IHistogram1D mollerMass;
- private IHistogram1D mollerMassVtxCut;
- private IHistogram1D mollerUx;
- private IHistogram1D mollerUy;
- private IHistogram1D mollerVx;
-
- private IHistogram1D mollerVy;
-
- private IHistogram1D mollerVz;
- private IHistogram1D mollerVzVtxCut;
- private IHistogram2D mollerXVsVtxY;
- private IHistogram2D mollerXVsVtxZ;
-
- private IHistogram2D mollerYVsVtxZ;
- // some counters
- private int nRecoEvents = 0;
- private int nTotV0 = 0;
- private IHistogram1D numChargeHisto;
- private IHistogram1D nV0;
- private IHistogram1D pEle;
- private IHistogram2D pEleVspEle;
- private IHistogram2D pEleVspEleBeamBeam;
- private IHistogram2D pEleVspEleMoller;
- private IHistogram2D pEleVspEleNoBeam;
- private IHistogram2D pEleVspPos;
- private IHistogram2D pEleVspPosWithCut;
- private IHistogram2D pEleVsthetaBeamBeam;
- private IHistogram2D pEleVsthetaMoller;
-
- private IHistogram2D phiEleVsphiEle;
- private final String plotDir = "V0Monitoring/";
- private IHistogram1D pPos;
-
- private IHistogram2D pxEleVspxEle;
- private IHistogram2D pxEleVspxEleNoBeam;
- private IHistogram2D pxEleVspxPos;
-
- private IHistogram2D pyEleVspyEle;
- private IHistogram2D pyEleVspyEleNoBeam;
- private IHistogram2D pyEleVspyPos;
- private IHistogram1D sumChargeHisto;
- private double sumChi2 = 0.0;
- // some summers
- private double sumMass = 0.0;
- private double sumVx = 0.0;
- private double sumVy = 0.0;
- private double sumVz = 0.0;
-
- private IHistogram1D tarconChi2;
- private IHistogram2D tarconChi2VsTrkChi2;
-
+ private IHistogram2D bsconChi2VsTrkChi2;
/* target constrained */
private IHistogram1D tarconMass;
private IHistogram1D tarconVx;
private IHistogram1D tarconVy;
private IHistogram1D tarconVz;
+ private IHistogram1D tarconChi2;
private IHistogram2D tarconVzVsChi2;
- private final String targetV0ConCandidatesColName = "TargetConstrainedV0Candidates";
+ private IHistogram2D tarconChi2VsTrkChi2;
+
+ private IHistogram2D pEleVspPos;
+ private IHistogram1D pEle;
+ private IHistogram1D pPos;
+
+ private IHistogram2D pEleVspPosWithCut;
+ private IHistogram2D pyEleVspyPos;
+ private IHistogram2D pxEleVspxPos;
+
+ private IHistogram2D VtxZVsMass;
+ private IHistogram2D VtxYVsVtxZ;
+ private IHistogram2D VtxXVsVtxZ;
+ private IHistogram2D VtxXVsVtxY;
+ private IHistogram2D VtxXVsVtxPx;
+ private IHistogram2D VtxYVsVtxPy;
+ private IHistogram2D VtxZVsVtxPx;
+ private IHistogram2D VtxZVsVtxPy;
+ private IHistogram2D VtxZVsVtxPz;
+
+ private IHistogram2D VtxZVsL1Iso;
+ private IHistogram2D VtxZVsTrkChi2;
+
+ private IHistogram2D pEleVspEle;
+ private IHistogram2D phiEleVsphiEle;
+ private IHistogram2D pyEleVspyEle;
+ private IHistogram2D pxEleVspxEle;
+ private IHistogram2D pEleVspEleNoBeam;
+ private IHistogram2D pyEleVspyEleNoBeam;
+ private IHistogram2D pxEleVspxEleNoBeam;
+ private IHistogram2D pEleVspEleMoller;
+ private IHistogram2D pEleVsthetaMoller;
+ private IHistogram2D thetaEleVsthetaMoller;
+ private IHistogram2D pEleVspEleBeamBeam;
+ private IHistogram2D pEleVsthetaBeamBeam;
private IHistogram2D thetaEleVsthetaBeamBeam;
- private IHistogram2D thetaEleVsthetaMoller;
- private final double thetaMax = 0.06;
- private final double thetaMin = 0.015;
- private IHistogram1D trigTime;
- private IHistogram2D trigTimeV0Time;
- private IHistogram1D unconChi2;
-
- private IHistogram2D unconChi2VsTrkChi2;
- /* beamspot constrained */
- /* V0 Quantities */
- /* Mass, vertex, chi^2 of fit */
- /* unconstrained */
- private IHistogram1D unconMass;
- private final String unconstrainedV0CandidatesColName = "UnconstrainedV0Candidates";
- private IHistogram1D unconVx;
- private IHistogram1D unconVy;
- private IHistogram1D unconVz;
- private IHistogram2D unconVzVsChi2;
- private IHistogram1D v0Dt;
- private double v0ESumMinCut, v0MaxPCut, v0ESumMaxCut, molPSumMin, molPSumMax, beambeamCut;
-
- private IHistogram1D v0Time;
- private IHistogram2D VtxXVsVtxPx;
-
- private IHistogram2D VtxXVsVtxY;
- private IHistogram2D VtxXVsVtxZ;
-
- private IHistogram2D VtxYVsVtxPy;
-
- private IHistogram2D VtxYVsVtxZ;
-
- private IHistogram2D VtxZVsL1Iso;
-
- private IHistogram2D VtxZVsMass;
- private IHistogram2D VtxZVsTrkChi2;
-
- private IHistogram2D VtxZVsVtxPx;
-
- private IHistogram2D VtxZVsVtxPy;
-
- private IHistogram2D VtxZVsVtxPz;
+
+ private IHistogram1D mollerMass;
+ private IHistogram1D mollerMassVtxCut;
+ private IHistogram1D mollerVx;
+ private IHistogram1D mollerVy;
+ private IHistogram1D mollerVz;
+ private IHistogram1D mollerVzVtxCut;
+ private IHistogram2D mollerXVsVtxZ;
+ private IHistogram2D mollerYVsVtxZ;
+ private IHistogram2D mollerXVsVtxY;
+
+ private IHistogram1D mollerUx;
+ private IHistogram1D mollerUy;
+
+
+
+ private IHistogram1D sumChargeHisto;
+ private IHistogram1D numChargeHisto;
+
+ private final String plotDir = "V0Monitoring/";
+
+ private final BasicHep3Matrix beamAxisRotation = new BasicHep3Matrix();
+
+ private double maxFactor = 1.25;
+
+
+ private double thetaMax = 0.06;
+ private double thetaMin = 0.015;
+
+
+
+ private double feeMomentumCut, v0ESumMinCut, v0MaxPCut, v0ESumMaxCut,
+ molPSumMin, molPSumMax, beambeamCut;
+
+
+
+ @Override
+ protected void detectorChanged(Detector detector) {
+
+ super.detectorChanged(detector);
+
+ feeMomentumCut = 0.75*beamEnergy; //GeV
+
+ v0ESumMinCut = 0.8 * beamEnergy;
+ v0ESumMaxCut = 1.25 * beamEnergy;
+
+ v0MaxPCut = 1.05*beamEnergy;//GeV
+ molPSumMin = 0.80*beamEnergy;
+ molPSumMax = 1.25*beamEnergy;
+ beambeamCut = 0.80*beamEnergy;
+
+
+
+ beamAxisRotation.setActiveEuler(Math.PI / 2, -0.0305, -Math.PI / 2);
+
+ LOGGER.info("Setting up the plotter");
+ aida.tree().cd("/");
+ String xtra = "Extras";
+ String trkType = "SeedTrack/";
+ if (isGBL)
+ trkType = "GBLTrack/";
+
+ double maxMass = .2*beamEnergy;
+ double maxMassMoller = .1*Math.sqrt(beamEnergy);
+ /* V0 Quantities */
+ /* Mass, vertex, chi^2 of fit */
+ /* unconstrained */
+ unconMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Invariant Mass (GeV)", 100, 0, maxMass);
+ unconVx = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Vx (mm)", 50, -10, 10);
+ unconVy = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Vy (mm)", 50, -10, 10);
+ unconVz = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Vz (mm)", 50, -50, 50);
+ unconChi2 = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Chi2", 25, 0, 25);
+ unconVzVsChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Vz vs. Chi2", 25, 0, 25, 50, -50, 50);
+ unconChi2VsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/" + "Chi2 vs. total track chi2", 50, 0, 50, 50, 0, 25);
+ /* beamspot constrained */
+ bsconMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Mass (GeV)", 100, 0, maxMass);
+ bsconVx = aida.histogram1D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Vx (mm)", 50, -10, 10);
+ bsconVy = aida.histogram1D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Vy (mm)", 50, -10, 10);
+ bsconVz = aida.histogram1D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Vz (mm)", 50, -50, 50);
+ bsconChi2 = aida.histogram1D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Chi2", 25, 0, 25);
+ bsconVzVsChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Vz vs. Chi2", 25, 0, 25, 50, -50, 50);
+ bsconChi2VsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Chi2 vs. total track chi2", 50, 0, 50, 50, 0, 25);
+ /* target constrained */
+ tarconMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/" + "Mass (GeV)", 100, 0, maxMass);
+ tarconVx = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/" + "Vx (mm)", 50, -1, 1);
+ tarconVy = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/" + "Vy (mm)", 50, -1, 1);
+ tarconVz = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/" + "Vz (mm)", 50, -10, 10);
+ tarconChi2 = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/" + "Chi2", 25, 0, 25);
+ tarconVzVsChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/" + "Vz vs. Chi2", 25, 0, 25, 50, -50, 50);
+ tarconChi2VsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/" + "Chi2 vs. total track chi2", 50, 0, 50, 50, 0, 25);
+
+ nV0 = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Number of V0 per event", 10, 0, 10);
+ v0Time = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "V0 mean time", 100, -25, 25);
+ v0Dt = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "V0 time difference", 100, -25, 25);
+ trigTimeV0Time = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Trigger phase vs. V0 mean time", 100, -25, 25, 6, 0, 24);
+ trigTime = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Trigger phase", 6, 0, 24);
+
+ pEleVspPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "P(e) vs P(p)", 50, 0, beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
+
+
+ pEle = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "P(e)", 50, 0, beamEnergy * maxFactor);
+ pPos = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "P(p)", 50, 0, beamEnergy * maxFactor);
+
+ pEleVspPosWithCut = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "P(e) vs P(p): Radiative", 50, 0, beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
+ pyEleVspyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Py(e) vs Py(p)", 50, -0.04*beamEnergy, 0.04*beamEnergy, 50, -0.04*beamEnergy, 0.04*beamEnergy);
+ pxEleVspxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Px(e) vs Px(p)", 50, -0.04*beamEnergy, 0.04*beamEnergy, 50, -0.04*beamEnergy, 0.04*beamEnergy);
+ VtxZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Mass", 50, 0, maxMass, 50, -50, 80);
+ VtxXVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vx vs Vz", 100, -10, 10, 100, -50, 80);
+ VtxYVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vy vs Vz", 100, -5, 5, 100, -50, 80);
+ VtxXVsVtxY = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vx vs Vy", 100, -10, 10, 100, -5, 5);
+ VtxXVsVtxPx = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vx vs Px", 100, -0.1, 0.1, 100, -10, 10);
+ VtxYVsVtxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vy vs Py", 100, -0.1, 0.1, 100, -5, 5);
+ VtxZVsVtxPx = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Px", 100, -0.1, 0.1, 100, -50, 80);
+ VtxZVsVtxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Py", 100, -0.1, 0.1, 100, -50, 80);
+ VtxZVsVtxPz = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Pz", 100, 0.0, beamEnergy * maxFactor, 100, -50, 80);
+ VtxZVsL1Iso = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs L1 Isolation", 100, 0.0, 5.0, 50, -50, 80);
+ VtxZVsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Track Chi2", 50, 0, 50, 50, -50, 80);
+ phiEleVsphiEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/phi(e) vs phi(e)", 50, -Math.PI, Math.PI, 50, -Math.PI, Math.PI);
+ pyEleVspyEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Py(e) vs Py(e)", 50, -0.04*beamEnergy, 0.04*beamEnergy, 50, -0.04*beamEnergy, 0.04*beamEnergy);
+ pxEleVspxEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Px(e) vs Px(e)", 50, -0.02*beamEnergy, 0.06*beamEnergy, 50, -0.02*beamEnergy, 0.06*beamEnergy);
+
+ // electron vs electron momentum with different cuts
+ // 1) no cut
+ // 2) cut out FEE
+ // 3) cut out FEE and also cut on momentum sum
+ // 4) cut out everything except FEE coincidentals
+ pEleVspEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(e) vs P(e)", 50, 0, beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
+ pEleVspEleNoBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(e) vs P(e) NoBeam", 50, 0, beambeamCut, 50, 0, beambeamCut);
+ pEleVspEleMoller = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(e) vs P(e) Moller", 50, 0, beambeamCut, 50, 0, beambeamCut);
+ pEleVspEleBeamBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(e) vs P(e) BeamBeam", 50, beambeamCut, beamEnergy * maxFactor, 50, beambeamCut, beamEnergy * maxFactor);
+
+ pyEleVspyEleNoBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Py(e) vs Py(e) NoBeam", 50, -0.04*beamEnergy, 0.04*beamEnergy, 50, -0.04*beamEnergy, 0.04*beamEnergy);
+ pxEleVspxEleNoBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Px(e) vs Px(e) NoBeam", 50, -0.02*beamEnergy, 0.06*beamEnergy, 50, -0.02*beamEnergy, 0.06*beamEnergy);
+ sumChargeHisto = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Total Charge of Event", 5, -2, 3);
+ numChargeHisto = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Number of Charged Particles", 6, 0, 6);
+
+ pEleVsthetaMoller = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(e) vs Theta Moller", 50, 0, beambeamCut, 50, thetaMin, thetaMax);
+ thetaEleVsthetaMoller = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Theta vs Theta Moller", 50, thetaMin, thetaMax, 50, thetaMin, thetaMax);
+ pEleVsthetaBeamBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(e) vs Theta BeamBeam", 50, beambeamCut, beamEnergy * maxFactor, 50, thetaMin, thetaMax);
+ thetaEleVsthetaBeamBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Theta vs Theta BeamBeam", 50, thetaMin, thetaMax, 50, thetaMin, thetaMax);
+
+ mollerMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Mass (GeV)", 100, 0, maxMassMoller);
+ mollerMassVtxCut = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Mass (GeV): VtxCut", 100, 0, maxMassMoller);
+ mollerVx = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vx (mm)", 50, -10, 10);
+ mollerVy = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vy (mm)", 50, -2, 2);
+ mollerVz = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vz (mm)", 50, -50, 50);
+ mollerVzVtxCut = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vz (mm): VtxCut", 50, -50, 50);
+ mollerXVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vx vs Vz", 100, -5, 5, 100, -50, 50);
+ mollerYVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vy vs Vz", 100, -2, 2, 100, -50, 50);
+ mollerXVsVtxY = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vx vs Vy", 100, -5, 5, 100, -2, 2);
+
+ mollerUx = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Pair Momentum Direction Ux", 100, .015, .045);
+ mollerUy = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Pair Momentum Direction Uy", 100, -.01, .01);
+
+ mollerHiP = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(high)", 100, 0, beamEnergy*maxFactor);
+ mollerLoP = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(low)", 100, 0, beamEnergy*maxFactor);
+
+ mollerEitherP = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(either)", 100, 0, beamEnergy*maxFactor);
+ mollerPsum = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Psum", 100, 0, beamEnergy*maxFactor);
+
+ }
+
+ private IHistogram1D mollerHiP, mollerLoP, mollerEitherP, mollerPsum;
+
+ @Override
+ public void process(EventHeader event) {
+ /* make sure everything is there */
+ if (!event.hasCollection(ReconstructedParticle.class, finalStateParticlesColName))
+ throw new IllegalArgumentException("Collection missing: " + finalStateParticlesColName);
+ if (!event.hasCollection(ReconstructedParticle.class, unconstrainedV0CandidatesColName))
+ throw new IllegalArgumentException("Collection missing: " + unconstrainedV0CandidatesColName);
+ if (!event.hasCollection(ReconstructedParticle.class, beamConV0CandidatesColName))
+ throw new IllegalArgumentException("Collection missing: " + beamConV0CandidatesColName);
+ if (!event.hasCollection(ReconstructedParticle.class, targetV0ConCandidatesColName))
+ throw new IllegalArgumentException("Collection missing: " + targetV0ConCandidatesColName);
+
+ //check to see if this event is from the correct trigger (or "all");
+ if (!matchTrigger(event))
+ return;
+
+ nRecoEvents++;
+
+ RelationalTable hitToStrips = TrackUtils.getHitToStripsTable(event);
+ RelationalTable hitToRotated = TrackUtils.getHitToRotatedTable(event);
+
+ List<ReconstructedParticle> unonstrainedV0List = event.get(ReconstructedParticle.class, unconstrainedV0CandidatesColName);
+ for (ReconstructedParticle uncV0 : unonstrainedV0List) {
+ if (isGBL != TrackType.isGBL(uncV0.getType()))
+ continue;
+ Vertex uncVert = uncV0.getStartVertex();
+ Hep3Vector pVtxRot = VecOp.mult(beamAxisRotation, uncV0.getMomentum());
+ Hep3Vector vtxPosRot = VecOp.mult(beamAxisRotation, uncVert.getPosition());
+ double theta = Math.acos(pVtxRot.z() / pVtxRot.magnitude());
+ double phi = Math.atan2(pVtxRot.y(), pVtxRot.x());
+ unconVx.fill(vtxPosRot.x());
+ unconVy.fill(vtxPosRot.y());
+ unconVz.fill(vtxPosRot.z());
+ unconMass.fill(uncV0.getMass());
+ unconChi2.fill(uncVert.getChi2());
+ unconVzVsChi2.fill(uncVert.getChi2(), vtxPosRot.z());
+ unconChi2VsTrkChi2.fill(Math.max(uncV0.getParticles().get(0).getTracks().get(0).getChi2(), uncV0.getParticles().get(1).getTracks().get(0).getChi2()), uncVert.getChi2());
+
+ VtxZVsMass.fill(uncV0.getMass(), vtxPosRot.z());
+ VtxXVsVtxZ.fill(vtxPosRot.x(), vtxPosRot.z());
+ VtxYVsVtxZ.fill(vtxPosRot.y(), vtxPosRot.z());
+ VtxXVsVtxY.fill(vtxPosRot.x(), vtxPosRot.y());
+ VtxXVsVtxPx.fill(pVtxRot.x(), vtxPosRot.x());
+ VtxYVsVtxPy.fill(pVtxRot.y(), vtxPosRot.y());
+ VtxZVsVtxPx.fill(pVtxRot.x(), vtxPosRot.z());
+ VtxZVsVtxPy.fill(pVtxRot.y(), vtxPosRot.z());
+ VtxZVsVtxPz.fill(pVtxRot.z(), vtxPosRot.z());
+
+ //this always has 2 tracks.
+ List<ReconstructedParticle> trks = uncV0.getParticles();
+ // Track ele = trks.get(0).getTracks().get(0);
+ // Track pos = trks.get(1).getTracks().get(0);
+ // //if track #0 has charge>0 it's the electron! This seems mixed up, but remember the track
+ // //charge is assigned assuming a positive B-field, while ours is negative
+ // if (trks.get(0).getCharge() > 0) {
+ // pos = trks.get(0).getTracks().get(0);
+ // ele = trks.get(1).getTracks().get(0);
+ // }
+ // aida.histogram2D(plotDir + trkType + triggerType + "/" + "P(e) vs P(p)").fill(getMomentum(ele), getMomentum(pos));
+ // aida.histogram2D(plotDir + trkType + triggerType + "/" + "Px(e) vs Px(p)").fill(ele.getTrackStates().get(0).getMomentum()[1], pos.getTrackStates().get(0).getMomentum()[1]);
+ // aida.histogram2D(plotDir + trkType + triggerType + "/" + "Py(e) vs Py(p)").fill(ele.getTrackStates().get(0).getMomentum()[2], pos.getTrackStates().get(0).getMomentum()[2]);
+ ReconstructedParticle ele = trks.get(0);
+ ReconstructedParticle pos = trks.get(1);
+ //ReconParticles have the charge correct.
+ if (trks.get(0).getCharge() > 0) {
+ pos = trks.get(0);
+ ele = trks.get(1);
+ }
+ if (ele.getCharge() < 0 && pos.getCharge() > 0) {
+ VtxZVsTrkChi2.fill(Math.max(uncV0.getParticles().get(0).getTracks().get(0).getChi2(), uncV0.getParticles().get(1).getTracks().get(0).getChi2()), uncVert.getPosition().z());
+
+ Double[] eleIso = TrackUtils.getIsolations(ele.getTracks().get(0), hitToStrips, hitToRotated);
+ Double[] posIso = TrackUtils.getIsolations(pos.getTracks().get(0), hitToStrips, hitToRotated);
+ if (eleIso[0] != null && posIso[0] != null) {
+ double eleL1Iso = Math.min(Math.abs(eleIso[0]), Math.abs(eleIso[1]));
+ double posL1Iso = Math.min(Math.abs(posIso[0]), Math.abs(posIso[1]));
+ double minL1Iso = Math.min(eleL1Iso, posL1Iso);
+ VtxZVsL1Iso.fill(minL1Iso, uncVert.getPosition().z());
+ }
+
+ double pe = ele.getMomentum().magnitude();
+ double pp = pos.getMomentum().magnitude();
+ Hep3Vector pEleRot = VecOp.mult(beamAxisRotation, ele.getMomentum());
+ Hep3Vector pPosRot = VecOp.mult(beamAxisRotation, pos.getMomentum());
+
+ pEleVspPos.fill(pe, pp);
+ pEle.fill(pe);
+ pPos.fill(pp);
+
+
+ pxEleVspxPos.fill(pEleRot.x(), pPosRot.x());
+ pyEleVspyPos.fill(pEleRot.y(), pPosRot.y());
+ if (pe < v0MaxPCut && pp < v0MaxPCut && (pe + pp) > v0ESumMinCut && (pe + pp) < v0ESumMaxCut)//enrich radiative-like events
+
+ pEleVspPosWithCut.fill(pe, pp);
+ }
+
+ double eleT = TrackUtils.getTrackTime(ele.getTracks().get(0), hitToStrips, hitToRotated);
+ double posT = TrackUtils.getTrackTime(pos.getTracks().get(0), hitToStrips, hitToRotated);
+ double meanT = (eleT + posT) / 2.0;
+ v0Time.fill(meanT);
+ v0Dt.fill(eleT - posT);
+ trigTimeV0Time.fill(meanT, event.getTimeStamp() % 24);
+ trigTime.fill(event.getTimeStamp() % 24);
+ }
+
+ List<ReconstructedParticle> beamConstrainedV0List = event.get(ReconstructedParticle.class, beamConV0CandidatesColName);
+ //nV0.fill(beamConstrainedV0List.size());
+ int v0Count = 0;
+ for (ReconstructedParticle bsV0 : beamConstrainedV0List) {
+
+ if (isGBL != TrackType.isGBL(bsV0.getType()))
+ continue;
+ v0Count++;
+ nTotV0++;
+ Vertex bsVert = bsV0.getStartVertex();
+ Hep3Vector vtxPosRot = VecOp.mult(beamAxisRotation, bsVert.getPosition());
+ bsconVx.fill(vtxPosRot.x());
+ bsconVy.fill(vtxPosRot.y());
+ bsconVz.fill(vtxPosRot.z());
+ bsconMass.fill(bsV0.getMass());
+ bsconChi2.fill(bsVert.getChi2());
+ bsconVzVsChi2.fill(bsVert.getChi2(), vtxPosRot.z());
+ bsconChi2VsTrkChi2.fill(Math.max(bsV0.getParticles().get(0).getTracks().get(0).getChi2(), bsV0.getParticles().get(1).getTracks().get(0).getChi2()), bsVert.getChi2());
+ sumMass += bsV0.getMass();
+ sumVx += vtxPosRot.x();
+ sumVy += vtxPosRot.y();
+ sumVz += vtxPosRot.z();
+ sumChi2 += bsVert.getChi2();
+ }
+ nV0.fill(v0Count);
+ List<ReconstructedParticle> targetConstrainedV0List = event.get(ReconstructedParticle.class, targetV0ConCandidatesColName);
+ for (ReconstructedParticle tarV0 : targetConstrainedV0List) {
+
+ if (isGBL != TrackType.isGBL(tarV0.getType()))
+ continue;
+
+ Vertex tarVert = tarV0.getStartVertex();
+ Hep3Vector vtxPosRot = VecOp.mult(beamAxisRotation, tarVert.getPosition());
+ tarconVx.fill(vtxPosRot.x());
+ tarconVy.fill(vtxPosRot.y());
+ tarconVz.fill(vtxPosRot.z());
+ tarconMass.fill(tarV0.getMass());
+ tarconChi2.fill(tarVert.getChi2());
+ tarconVzVsChi2.fill(tarVert.getChi2(), vtxPosRot.z());
+ tarconChi2VsTrkChi2.fill(Math.max(tarV0.getParticles().get(0).getTracks().get(0).getChi2(), tarV0.getParticles().get(1).getTracks().get(0).getChi2()), tarVert.getChi2());
+ }
+ List<ReconstructedParticle> finalStateParticles = event.get(ReconstructedParticle.class, finalStateParticlesColName);
+ if (debug)
+ LOGGER.info("This events has " + finalStateParticles.size() + " final state particles");
+
+ ReconstructedParticle ele1 = null;
+ ReconstructedParticle ele2 = null;
+ int sumCharge = 0;
+ int numChargedParticles = 0;
+ for (ReconstructedParticle fsPart : finalStateParticles) {
+ if (isGBL != TrackType.isGBL(fsPart.getType()))
+ continue;
+ if (debug)
+ LOGGER.info("PDGID = " + fsPart.getParticleIDUsed() + "; charge = " + fsPart.getCharge() + "; pz = " + fsPart.getMomentum().x());
+ double charge = fsPart.getCharge();
+ sumCharge += charge;
+ if (charge != 0) {
+ numChargedParticles++;
+ if (charge < 1)
+ if (ele1 == null)
+ ele1 = fsPart;
+ else if (!hasSharedStrips(ele1, fsPart, hitToStrips, hitToRotated))
+ ele2 = fsPart;
+ }
+ }
+ sumChargeHisto.fill(sumCharge);
+ numChargeHisto.fill(numChargedParticles);
+
+ if (ele1 != null && ele2 != null) {
+ Hep3Vector p1 = VecOp.mult(beamAxisRotation, ele1.getMomentum());
+ Hep3Vector p2 = VecOp.mult(beamAxisRotation, ele2.getMomentum());
+ // Hep3Vector beamAxis = new BasicHep3Vector(Math.sin(0.0305), 0, Math.cos(0.0305));
+ // LOGGER.info(p1);
+ // LOGGER.info(VecOp.mult(rot, p1));
+
+ double theta1 = Math.acos(p1.z() / p1.magnitude());
+ double theta2 = Math.acos(p2.z() / p2.magnitude());
+ double phi1 = Math.atan2(p1.y(), p1.x());
+ double phi2 = Math.atan2(p2.y(), p2.x());
+ phiEleVsphiEle.fill(phi1, phi2);
+ pEleVspEle.fill(ele1.getMomentum().magnitude(), ele2.getMomentum().magnitude());
+ pyEleVspyEle.fill(ele1.getMomentum().y(), ele2.getMomentum().y());
+ pxEleVspxEle.fill(ele1.getMomentum().x(), ele2.getMomentum().x());
+ //remove beam electrons
+ if (ele1.getMomentum().magnitude() < beambeamCut && ele2.getMomentum().magnitude() < beambeamCut) {
+ pEleVspEleNoBeam.fill(ele1.getMomentum().magnitude(), ele2.getMomentum().magnitude());
+ pyEleVspyEleNoBeam.fill(ele1.getMomentum().y(), ele2.getMomentum().y());
+ pxEleVspxEleNoBeam.fill(ele1.getMomentum().x(), ele2.getMomentum().x());
+ }
+ //look at beam-beam events
+ if (ele1.getMomentum().magnitude() > beambeamCut && ele2.getMomentum().magnitude() > beambeamCut) {
+ pEleVspEleBeamBeam.fill(ele1.getMomentum().magnitude(), ele2.getMomentum().magnitude());
+ pEleVsthetaBeamBeam.fill(p1.magnitude(), theta1);
+ pEleVsthetaBeamBeam.fill(p2.magnitude(), theta2);
+ thetaEleVsthetaBeamBeam.fill(theta1, theta2);
+ }
+
+ //look at "Moller" events (if that's what they really are
+ if (ele1.getMomentum().magnitude() + ele2.getMomentum().magnitude() > molPSumMin
+ && ele1.getMomentum().magnitude() + ele2.getMomentum().magnitude() < molPSumMax
+ && (p1.magnitude() < beambeamCut && p2.magnitude() < beambeamCut)) {
+
+ Track ele1trk = ele1.getTracks().get(0);
+ Track ele2trk = ele2.getTracks().get(0);
+ SeedTrack stEle1 = TrackUtils.makeSeedTrackFromBaseTrack(ele1trk);
+ SeedTrack stEle2 = TrackUtils.makeSeedTrackFromBaseTrack(ele2trk);
+ BilliorTrack btEle1 = new BilliorTrack(stEle1.getSeedCandidate().getHelix());
+ BilliorTrack btEle2 = new BilliorTrack(stEle2.getSeedCandidate().getHelix());
+ BilliorVertex bv = fitVertex(btEle1, btEle2, TrackUtils.getBField(event.getDetector()).magnitude());
+ // LOGGER.info("ee vertex: "+bv.toString());
+ double invMass = bv.getParameters().get("invMass");
+ mollerMass.fill(invMass);
+ mollerVx.fill(bv.getPosition().x());
+ mollerVy.fill(bv.getPosition().y());
+ mollerVz.fill(bv.getPosition().z());
+ mollerXVsVtxZ.fill(bv.getPosition().x(), bv.getPosition().z());
+ mollerYVsVtxZ.fill(bv.getPosition().y(), bv.getPosition().z());
+ mollerXVsVtxY.fill(bv.getPosition().x(), bv.getPosition().y());
+
+ double ux = (ele1.getMomentum().x()+ele2.getMomentum().x())/(ele1.getMomentum().z()+ele2.getMomentum().z());
+ double uy = (ele1.getMomentum().y()+ele2.getMomentum().y())/(ele1.getMomentum().z()+ele2.getMomentum().z());
+ mollerUx.fill(ux);
+ mollerUy.fill(uy);
+
+ //higher and lower energy electrons in moller pair
+ double pt1 = ele1.getMomentum().magnitude();
+ double pt2 = ele2.getMomentum().magnitude();
+ double ph = (pt1>pt2) ? pt1 : pt2;
+ double pl = (pt1>pt2) ? pt2 : pt1;
+
+ mollerHiP.fill(ph);
+ mollerLoP.fill(pl);
+
+ mollerEitherP.fill(ph);
+ mollerEitherP.fill(pl);
+ mollerPsum.fill(pt1+pt2);
+
+
+ if (Math.abs(bv.getPosition().x()) < 2
+ && Math.abs(bv.getPosition().y()) < 0.5) {
+ mollerMassVtxCut.fill(invMass);
+ mollerVzVtxCut.fill(bv.getPosition().z());
+ }
+ pEleVspEleMoller.fill(p1.magnitude(), p2.magnitude());
+ pEleVsthetaMoller.fill(p1.magnitude(), theta1);
+ pEleVsthetaMoller.fill(p2.magnitude(), theta2);
+ thetaEleVsthetaMoller.fill(theta1, theta2);
+ }
+ }
+ }
+
+ @Override
+ public void printDQMData() {
+ LOGGER.info("V0Monitoring::printDQMData");
+ for (Entry<String, Double> entry : monitoredQuantityMap.entrySet())
+ LOGGER.info(entry.getKey() + " = " + entry.getValue());
+ LOGGER.info("*******************************");
+ }
/**
* Calculate the averages here and fill the map
@@ -169,29 +563,27 @@
@Override
public void calculateEndOfRunQuantities() {
- final IAnalysisFactory analysisFactory = IAnalysisFactory.create();
- final IFitFactory fitFactory = analysisFactory.createFitFactory();
- final IFitter fitter = fitFactory.createFitter("chi2");
- final double[] init = {50.0, 0.0, 0.2, 1.0, 0.0};
- final IFitResult resVx = this.fitVertexPosition(bsconVx, fitter, init, "range=\"(-0.5,0.5)\"");
- final double[] init2 = {50.0, 0.0, 0.04, 1.0, 0.0};
- final IFitResult resVy = this.fitVertexPosition(bsconVy, fitter, init2, "range=\"(-0.2,0.2)\"");
- final double[] init3 = {50.0, 0.0, 3.0, 1.0, 0.0};
- final IFitResult resVz = this.fitVertexPosition(bsconVz, fitter, init3, "range=\"(-6,6)\"");
+ IAnalysisFactory analysisFactory = IAnalysisFactory.create();
+ IFitFactory fitFactory = analysisFactory.createFitFactory();
+ IFitter fitter = fitFactory.createFitter("chi2");
+ double[] init = {50.0, 0.0, 0.2, 1.0, 0.0};
+ IFitResult resVx = fitVertexPosition(bsconVx, fitter, init, "range=\"(-0.5,0.5)\"");
+ double[] init2 = {50.0, 0.0, 0.04, 1.0, 0.0};
+ IFitResult resVy = fitVertexPosition(bsconVy, fitter, init2, "range=\"(-0.2,0.2)\"");
+ double[] init3 = {50.0, 0.0, 3.0, 1.0, 0.0};
+ IFitResult resVz = fitVertexPosition(bsconVz, fitter, init3, "range=\"(-6,6)\"");
if (resVx != null && resVy != null & resVz != null) {
- final double[] parsVx = resVx.fittedParameters();
- final double[] parsVy = resVy.fittedParameters();
- final double[] parsVz = resVz.fittedParameters();
-
- for (int i = 0; i < 5; i++) {
- LOGGER.info("Vertex Fit Parameters: " + resVx.fittedParameterNames()[i] + " = " + parsVx[i] + "; "
- + parsVy[i] + "; " + parsVz[i]);
- }
-
- final IPlotter plotter = analysisFactory.createPlotterFactory().create("Vertex Position");
+ double[] parsVx = resVx.fittedParameters();
+ double[] parsVy = resVy.fittedParameters();
+ double[] parsVz = resVz.fittedParameters();
+
+ for (int i = 0; i < 5; i++)
+ LOGGER.info("Vertex Fit Parameters: " + resVx.fittedParameterNames()[i] + " = " + parsVx[i] + "; " + parsVy[i] + "; " + parsVz[i]);
+
+ IPlotter plotter = analysisFactory.createPlotterFactory().create("Vertex Position");
plotter.createRegions(1, 3);
- final IPlotterStyle pstyle = plotter.style();
+ IPlotterStyle pstyle = plotter.style();
pstyle.legendBoxStyle().setVisible(false);
pstyle.dataStyle().fillStyle().setColor("green");
pstyle.dataStyle().lineStyle().setColor("black");
@@ -201,563 +593,73 @@
plotter.region(1).plot(resVy.fittedFunction());
plotter.region(2).plot(bsconVz);
plotter.region(2).plot(resVz.fittedFunction());
- if (outputPlots) {
+ if (outputPlots)
try {
plotter.writeToFile(outputPlotDir + "vertex.png");
- } catch (final IOException ex) {
+ } catch (IOException ex) {
Logger.getLogger(V0Monitoring.class.getName()).log(Level.SEVERE, null, ex);
}
- }
-
- // monitoredQuantityMap.put(fpQuantNames[2], sumVx / nTotV0);
- // monitoredQuantityMap.put(fpQuantNames[3], sumVy / nTotV0);
- // monitoredQuantityMap.put(fpQuantNames[4], sumVz / nTotV0);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[2], parsVx[1]);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[3], parsVy[1]);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[4], parsVz[1]);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[5], parsVx[2]);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[6], parsVy[2]);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[7], parsVz[2]);
+
+ // monitoredQuantityMap.put(fpQuantNames[2], sumVx / nTotV0);
+ // monitoredQuantityMap.put(fpQuantNames[3], sumVy / nTotV0);
+ // monitoredQuantityMap.put(fpQuantNames[4], sumVz / nTotV0);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[2], parsVx[1]);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[3], parsVy[1]);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[4], parsVz[1]);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[5], parsVx[2]);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[6], parsVy[2]);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[7], parsVz[2]);
}
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[0],
- (double) nTotV0 / nRecoEvents);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[1], sumMass
- / nTotV0);
- monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[8], sumChi2
- / nTotV0);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[0], (double) nTotV0 / nRecoEvents);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[1], sumMass / nTotV0);
+ monitoredQuantityMap.put(beamConV0CandidatesColName + " " + triggerType + " " + fpQuantNames[8], sumChi2 / nTotV0);
}
@Override
- protected void detectorChanged(final Detector detector) {
-
- final BeamEnergyCollection beamEnergyCollection = this.getConditionsManager()
- .getCachedConditions(BeamEnergyCollection.class, "beam_energies").getCachedData();
- final double beamEnergy = beamEnergyCollection.get(0).getBeamEnergy();
- v0ESumMinCut = 0.8 * beamEnergy;
- v0ESumMaxCut = 1.25 * beamEnergy;
-
- v0MaxPCut = 1.05 * beamEnergy;// GeV
- molPSumMin = 0.80 * beamEnergy;
- molPSumMax = 1.25 * beamEnergy;
- beambeamCut = 0.80 * beamEnergy;
-
- beamAxisRotation.setActiveEuler(Math.PI / 2, -0.0305, -Math.PI / 2);
-
- // LOGGER.info("Setting up the plotter");
- aida.tree().cd("/");
- final String xtra = "Extras";
- String trkType = "SeedTrack/";
- if (isGBL) {
- trkType = "GBLTrack/";
- }
-
- final double maxMass = .2 * beamEnergy;
- final double maxMassMoller = .1 * Math.sqrt(beamEnergy);
- /* V0 Quantities */
- /* Mass, vertex, chi^2 of fit */
- /* unconstrained */
- unconMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/"
- + "Invariant Mass (GeV)", 100, 0, maxMass);
- unconVx = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/"
- + "Vx (mm)", 50, -10, 10);
- unconVy = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/"
- + "Vy (mm)", 50, -10, 10);
- unconVz = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/"
- + "Vz (mm)", 50, -50, 50);
- unconChi2 = aida.histogram1D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/"
- + "Chi2", 25, 0, 25);
- unconVzVsChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName + "/"
- + "Vz vs. Chi2", 25, 0, 25, 50, -50, 50);
- unconChi2VsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + unconstrainedV0CandidatesColName
- + "/" + "Chi2 vs. total track chi2", 50, 0, 50, 50, 0, 25);
- /* beamspot constrained */
- bsconMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/"
- + "Mass (GeV)", 100, 0, maxMass);
- bsconVx = aida.histogram1D(
- plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Vx (mm)", 50, -10, 10);
- bsconVy = aida.histogram1D(
- plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Vy (mm)", 50, -10, 10);
- bsconVz = aida.histogram1D(
- plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Vz (mm)", 50, -50, 50);
- bsconChi2 = aida.histogram1D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/" + "Chi2",
- 25, 0, 25);
- bsconVzVsChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/"
- + "Vz vs. Chi2", 25, 0, 25, 50, -50, 50);
- bsconChi2VsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + beamConV0CandidatesColName + "/"
- + "Chi2 vs. total track chi2", 50, 0, 50, 50, 0, 25);
- /* target constrained */
- tarconMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/"
- + "Mass (GeV)", 100, 0, maxMass);
- tarconVx = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/"
- + "Vx (mm)", 50, -1, 1);
- tarconVy = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/"
- + "Vy (mm)", 50, -1, 1);
- tarconVz = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/"
- + "Vz (mm)", 50, -10, 10);
- tarconChi2 = aida.histogram1D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/"
- + "Chi2", 25, 0, 25);
- tarconVzVsChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName + "/"
- + "Vz vs. Chi2", 25, 0, 25, 50, -50, 50);
- tarconChi2VsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + targetV0ConCandidatesColName
- + "/" + "Chi2 vs. total track chi2", 50, 0, 50, 50, 0, 25);
-
- nV0 = aida
- .histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Number of V0 per event", 10, 0, 10);
- v0Time = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "V0 mean time", 100, -25, 25);
- v0Dt = aida
- .histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "V0 time difference", 100, -25, 25);
- trigTimeV0Time = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "Trigger phase vs. V0 mean time", 100, -25, 25, 6, 0, 24);
- trigTime = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Trigger phase", 6, 0, 24);
-
- pEleVspPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "P(e) vs P(p)", 50, 0,
- beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
-
- pEle = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "P(e)", 50, 0, beamEnergy
- * maxFactor);
- pPos = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "P(p)", 50, 0, beamEnergy
- * maxFactor);
-
- pEleVspPosWithCut = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "P(e) vs P(p): Radiative", 50, 0, beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
- pyEleVspyPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Py(e) vs Py(p)", 50,
- -0.04 * beamEnergy, 0.04 * beamEnergy, 50, -0.04 * beamEnergy, 0.04 * beamEnergy);
- pxEleVspxPos = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Px(e) vs Px(p)", 50,
- -0.04 * beamEnergy, 0.04 * beamEnergy, 50, -0.04 * beamEnergy, 0.04 * beamEnergy);
- VtxZVsMass = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Mass", 50, 0,
- maxMass, 50, -50, 80);
- VtxXVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vx vs Vz", 100, -10, 10,
- 100, -50, 80);
- VtxYVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vy vs Vz", 100, -5, 5, 100,
- -50, 80);
- VtxXVsVtxY = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vx vs Vy", 100, -10, 10,
- 100, -5, 5);
- VtxXVsVtxPx = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vx vs Px", 100, -0.1, 0.1,
- 100, -10, 10);
- VtxYVsVtxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vy vs Py", 100, -0.1, 0.1,
- 100, -5, 5);
- VtxZVsVtxPx = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Px", 100, -0.1, 0.1,
- 100, -50, 80);
- VtxZVsVtxPy = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Py", 100, -0.1, 0.1,
- 100, -50, 80);
- VtxZVsVtxPz = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Pz", 100, 0.0,
- beamEnergy * maxFactor, 100, -50, 80);
- VtxZVsL1Iso = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs L1 Isolation", 100,
- 0.0, 5.0, 50, -50, 80);
- VtxZVsTrkChi2 = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "Vz vs Track Chi2", 50,
- 0, 50, 50, -50, 80);
- phiEleVsphiEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/phi(e) vs phi(e)", 50, -Math.PI, Math.PI, 50, -Math.PI, Math.PI);
- pyEleVspyEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Py(e) vs Py(e)", 50, -0.04 * beamEnergy, 0.04 * beamEnergy, 50, -0.04 * beamEnergy,
- 0.04 * beamEnergy);
- pxEleVspxEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Px(e) vs Px(e)", 50, -0.02 * beamEnergy, 0.06 * beamEnergy, 50, -0.02 * beamEnergy,
- 0.06 * beamEnergy);
-
- // electron vs electron momentum with different cuts
- // 1) no cut
- // 2) cut out FEE
- // 3) cut out FEE and also cut on momentum sum
- // 4) cut out everything except FEE coincidentals
- pEleVspEle = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(e) vs P(e)",
- 50, 0, beamEnergy * maxFactor, 50, 0, beamEnergy * maxFactor);
- pEleVspEleNoBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/P(e) vs P(e) NoBeam", 50, 0, beambeamCut, 50, 0, beambeamCut);
- pEleVspEleMoller = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/P(e) vs P(e) Moller", 50, 0, beambeamCut, 50, 0, beambeamCut);
- pEleVspEleBeamBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/P(e) vs P(e) BeamBeam", 50, beambeamCut, beamEnergy * maxFactor, 50, beambeamCut,
- beamEnergy * maxFactor);
-
- pyEleVspyEleNoBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Py(e) vs Py(e) NoBeam", 50, -0.04 * beamEnergy, 0.04 * beamEnergy, 50,
- -0.04 * beamEnergy, 0.04 * beamEnergy);
- pxEleVspxEleNoBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Px(e) vs Px(e) NoBeam", 50, -0.02 * beamEnergy, 0.06 * beamEnergy, 50,
- -0.02 * beamEnergy, 0.06 * beamEnergy);
- sumChargeHisto = aida.histogram1D(
- plotDir + trkType + triggerType + "/" + xtra + "/" + "Total Charge of Event", 5, -2, 3);
- numChargeHisto = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "Number of Charged Particles", 6, 0, 6);
-
- pEleVsthetaMoller = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/P(e) vs Theta Moller", 50, 0, beambeamCut, 50, thetaMin, thetaMax);
- thetaEleVsthetaMoller = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Theta vs Theta Moller", 50, thetaMin, thetaMax, 50, thetaMin, thetaMax);
- pEleVsthetaBeamBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/P(e) vs Theta BeamBeam", 50, beambeamCut, beamEnergy * maxFactor, 50, thetaMin, thetaMax);
- thetaEleVsthetaBeamBeam = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Theta vs Theta BeamBeam", 50, thetaMin, thetaMax, 50, thetaMin, thetaMax);
-
- mollerMass = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Mass (GeV)", 100, 0, maxMassMoller);
- mollerMassVtxCut = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Mass (GeV): VtxCut", 100, 0, maxMassMoller);
- mollerVx = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vx (mm)",
- 50, -10, 10);
- mollerVy = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vy (mm)",
- 50, -2, 2);
- mollerVz = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Moller Vz (mm)",
- 50, -50, 50);
- mollerVzVtxCut = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Vz (mm): VtxCut", 50, -50, 50);
- mollerXVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Vx vs Vz", 100, -5, 5, 100, -50, 50);
- mollerYVsVtxZ = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Vy vs Vz", 100, -2, 2, 100, -50, 50);
- mollerXVsVtxY = aida.histogram2D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Vx vs Vy", 100, -5, 5, 100, -2, 2);
-
- mollerUx = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Pair Momentum Direction Ux", 100, .015, .045);
- mollerUy = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/"
- + "2 Electron/Moller Pair Momentum Direction Uy", 100, -.01, .01);
-
- mollerHiP = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(high)", 100, 0,
- beamEnergy * maxFactor);
- mollerLoP = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(low)", 100, 0,
- beamEnergy * maxFactor);
-
- mollerEitherP = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/P(either)",
- 100, 0, beamEnergy * maxFactor);
- mollerPsum = aida.histogram1D(plotDir + trkType + triggerType + "/" + xtra + "/" + "2 Electron/Psum", 100, 0,
- beamEnergy * maxFactor);
-
- }
-
- private BilliorVertex fitVertex(final BilliorTrack electron, final BilliorTrack positron, final double bField) {
- // Create a vertex fitter from the magnetic field.
- final double[] beamSize = {0.001, 0.2, 0.02};
- final BilliorVertexer vtxFitter = new BilliorVertexer(bField);
- // TODO: The beam size should come from the conditions database.
- vtxFitter.setBeamSize(beamSize);
-
- // Perform the vertexing based on the specified constraint.
- vtxFitter.doBeamSpotConstraint(false);
-
- // Add the electron and positron tracks to a track list for
- // the vertex fitter.
- final List<BilliorTrack> billiorTracks = new ArrayList<BilliorTrack>();
-
- billiorTracks.add(electron);
-
- billiorTracks.add(positron);
-
- // Find and return a vertex based on the tracks.
- return vtxFitter.fitVertex(billiorTracks);
- }
-
- IFitResult fitVertexPosition(final IHistogram1D h1d, final IFitter fitter, final double[] init, final String range) {
+ public void printDQMStrings() {
+ for (int i = 0; i < 9; i++)//TODO: do this in a smarter way...loop over the map
+ LOGGER.info("ALTER TABLE dqm ADD " + fpQuantNames[i] + " double;");
+ }
+
+ private IFitResult fitVertexPosition(IHistogram1D h1d, IFitter fitter, double[] init, String range
+ ) {
IFitResult ifr = null;
try {
ifr = fitter.fit(h1d, "g+p1", init, range);
- } catch (final RuntimeException ex) {
+ } catch (RuntimeException ex) {
LOGGER.info(this.getClass().getSimpleName() + ": caught exception in fitGaussian");
}
return ifr;
}
- @Override
- public void printDQMData() {
- LOGGER.info("V0Monitoring::printDQMData");
- for (final Entry<String, Double> entry : monitoredQuantityMap.entrySet()) {
- LOGGER.info(entry.getKey() + " = " + entry.getValue());
- }
- LOGGER.info("*******************************");
- }
-
- @Override
- public void printDQMStrings() {
- for (int i = 0; i < 9; i++) {
- LOGGER.info("ALTER TABLE dqm ADD " + fpQuantNames[i] + " double;");
- }
- }
-
- @Override
- public void process(final EventHeader event) {
- /* make sure everything is there */
- if (!event.hasCollection(ReconstructedParticle.class, finalStateParticlesColName)) {
- return;
- }
- if (!event.hasCollection(ReconstructedParticle.class, unconstrainedV0CandidatesColName)) {
- return;
- }
- if (!event.hasCollection(ReconstructedParticle.class, beamConV0CandidatesColName)) {
- return;
- }
- if (!event.hasCollection(ReconstructedParticle.class, targetV0ConCandidatesColName)) {
- return;
- }
-
- // check to see if this event is from the correct trigger (or "all");
- if (!this.matchTrigger(event)) {
- return;
- }
-
- nRecoEvents++;
-
- final RelationalTable<?, ?> hitToStrips = TrackUtils.getHitToStripsTable(event);
- final RelationalTable<?, ?> hitToRotated = TrackUtils.getHitToRotatedTable(event);
-
- final List<ReconstructedParticle> unonstrainedV0List = event.get(ReconstructedParticle.class,
- unconstrainedV0CandidatesColName);
- for (final ReconstructedParticle uncV0 : unonstrainedV0List) {
- if (isGBL != TrackType.isGBL(uncV0.getType())) {
- continue;
- }
- final Vertex uncVert = uncV0.getStartVertex();
- final Hep3Vector pVtxRot = VecOp.mult(beamAxisRotation, uncV0.getMomentum());
- final Hep3Vector vtxPosRot = VecOp.mult(beamAxisRotation, uncVert.getPosition());
- final double theta = Math.acos(pVtxRot.z() / pVtxRot.magnitude());
- final double phi = Math.atan2(pVtxRot.y(), pVtxRot.x());
- unconVx.fill(vtxPosRot.x());
- unconVy.fill(vtxPosRot.y());
- unconVz.fill(vtxPosRot.z());
- unconMass.fill(uncV0.getMass());
- unconChi2.fill(uncVert.getChi2());
- unconVzVsChi2.fill(uncVert.getChi2(), vtxPosRot.z());
- unconChi2VsTrkChi2.fill(
- Math.max(uncV0.getParticles().get(0).getTracks().get(0).getChi2(), uncV0.getParticles().get(1)
- .getTracks().get(0).getChi2()), uncVert.getChi2());
-
- VtxZVsMass.fill(uncV0.getMass(), vtxPosRot.z());
- VtxXVsVtxZ.fill(vtxPosRot.x(), vtxPosRot.z());
- VtxYVsVtxZ.fill(vtxPosRot.y(), vtxPosRot.z());
- VtxXVsVtxY.fill(vtxPosRot.x(), vtxPosRot.y());
- VtxXVsVtxPx.fill(pVtxRot.x(), vtxPosRot.x());
- VtxYVsVtxPy.fill(pVtxRot.y(), vtxPosRot.y());
- VtxZVsVtxPx.fill(pVtxRot.x(), vtxPosRot.z());
- VtxZVsVtxPy.fill(pVtxRot.y(), vtxPosRot.z());
- VtxZVsVtxPz.fill(pVtxRot.z(), vtxPosRot.z());
-
- // this always has 2 tracks.
- final List<ReconstructedParticle> trks = uncV0.getParticles();
- // Track ele = trks.get(0).getTracks().get(0);
- // Track pos = trks.get(1).getTracks().get(0);
- // //if track #0 has charge>0 it's the electron! This seems mixed up, but remember the track
- // //charge is assigned assuming a positive B-field, while ours is negative
- // if (trks.get(0).getCharge() > 0) {
- // pos = trks.get(0).getTracks().get(0);
- // ele = trks.get(1).getTracks().get(0);
- // }
- // aida.histogram2D(plotDir + trkType + triggerType + "/" + "P(e) vs P(p)").fill(getMomentum(ele),
- // getMomentum(pos));
- // aida.histogram2D(plotDir + trkType + triggerType + "/" +
- // "Px(e) vs Px(p)").fill(ele.getTrackStates().get(0).getMomentum()[1],
- // pos.getTrackStates().get(0).getMomentum()[1]);
- // aida.histogram2D(plotDir + trkType + triggerType + "/" +
- // "Py(e) vs Py(p)").fill(ele.getTrackStates().get(0).getMomentum()[2],
- // pos.getTrackStates().get(0).getMomentum()[2]);
- ReconstructedParticle ele = trks.get(0);
- ReconstructedParticle pos = trks.get(1);
- // ReconParticles have the charge correct.
- if (trks.get(0).getCharge() > 0) {
- pos = trks.get(0);
- ele = trks.get(1);
- }
- if (ele.getCharge() < 0 && pos.getCharge() > 0) {
- VtxZVsTrkChi2.fill(
- Math.max(uncV0.getParticles().get(0).getTracks().get(0).getChi2(), uncV0.getParticles().get(1)
- .getTracks().get(0).getChi2()), uncVert.getPosition().z());
-
- final Double[] eleIso = TrackUtils.getIsolations(ele.getTracks().get(0), hitToStrips, hitToRotated);
- final Double[] posIso = TrackUtils.getIsolations(pos.getTracks().get(0), hitToStrips, hitToRotated);
- if (eleIso[0] != null && posIso[0] != null) {
- final double eleL1Iso = Math.min(Math.abs(eleIso[0]), Math.abs(eleIso[1]));
- final double posL1Iso = Math.min(Math.abs(posIso[0]), Math.abs(posIso[1]));
- final double minL1Iso = Math.min(eleL1Iso, posL1Iso);
- VtxZVsL1Iso.fill(minL1Iso, uncVert.getPosition().z());
- }
-
- final double pe = ele.getMomentum().magnitude();
- final double pp = pos.getMomentum().magnitude();
- final Hep3Vector pEleRot = VecOp.mult(beamAxisRotation, ele.getMomentum());
- final Hep3Vector pPosRot = VecOp.mult(beamAxisRotation, pos.getMomentum());
-
- pEleVspPos.fill(pe, pp);
- pEle.fill(pe);
- pPos.fill(pp);
-
- pxEleVspxPos.fill(pEleRot.x(), pPosRot.x());
- pyEleVspyPos.fill(pEleRot.y(), pPosRot.y());
- if (pe < v0MaxPCut && pp < v0MaxPCut && pe + pp > v0ESumMinCut && pe + pp < v0ESumMaxCut) {
- pEleVspPosWithCut.fill(pe, pp);
- }
- }
-
- final double eleT = TrackUtils.getTrackTime(ele.getTracks().get(0), hitToStrips, hitToRotated);
- final double posT = TrackUtils.getTrackTime(pos.getTracks().get(0), hitToStrips, hitToRotated);
- final double meanT = (eleT + posT) / 2.0;
- v0Time.fill(meanT);
- v0Dt.fill(eleT - posT);
- trigTimeV0Time.fill(meanT, event.getTimeStamp() % 24);
- trigTime.fill(event.getTimeStamp() % 24);
- }
-
- final List<ReconstructedParticle> beamConstrainedV0List = event.get(ReconstructedParticle.class,
- beamConV0CandidatesColName);
- nV0.fill(beamConstrainedV0List.size());
- for (final ReconstructedParticle bsV0 : beamConstrainedV0List) {
-
- if (isGBL != TrackType.isGBL(bsV0.getType())) {
- continue;
- }
- nTotV0++;
- final Vertex bsVert = bsV0.getStartVertex();
- final Hep3Vector vtxPosRot = VecOp.mult(beamAxisRotation, bsVert.getPosition());
- bsconVx.fill(vtxPosRot.x());
- bsconVy.fill(vtxPosRot.y());
- bsconVz.fill(vtxPosRot.z());
- bsconMass.fill(bsV0.getMass());
- bsconChi2.fill(bsVert.getChi2());
- bsconVzVsChi2.fill(bsVert.getChi2(), vtxPosRot.z());
- bsconChi2VsTrkChi2.fill(
- Math.max(bsV0.getParticles().get(0).getTracks().get(0).getChi2(), bsV0.getParticles().get(1)
- .getTracks().get(0).getChi2()), bsVert.getChi2());
- sumMass += bsV0.getMass();
- sumVx += vtxPosRot.x();
- sumVy += vtxPosRot.y();
- sumVz += vtxPosRot.z();
- sumChi2 += bsVert.getChi2();
- }
-
- final List<ReconstructedParticle> targetConstrainedV0List = event.get(ReconstructedParticle.class,
- targetV0ConCandidatesColName);
- for (final ReconstructedParticle tarV0 : targetConstrainedV0List) {
-
- if (isGBL != TrackType.isGBL(tarV0.getType())) {
- continue;
- }
-
- final Vertex tarVert = tarV0.getStartVertex();
- final Hep3Vector vtxPosRot = VecOp.mult(beamAxisRotation, tarVert.getPosition());
- tarconVx.fill(vtxPosRot.x());
- tarconVy.fill(vtxPosRot.y());
- tarconVz.fill(vtxPosRot.z());
- tarconMass.fill(tarV0.getMass());
- tarconChi2.fill(tarVert.getChi2());
- tarconVzVsChi2.fill(tarVert.getChi2(), vtxPosRot.z());
- tarconChi2VsTrkChi2.fill(
- Math.max(tarV0.getParticles().get(0).getTracks().get(0).getChi2(), tarV0.getParticles().get(1)
- .getTracks().get(0).getChi2()), tarVert.getChi2());
- }
- final List<ReconstructedParticle> finalStateParticles = event.get(ReconstructedParticle.class,
- finalStateParticlesColName);
- if (debug) {
- LOGGER.info("This events has " + finalStateParticles.size() + " final state particles");
- }
-
- ReconstructedParticle ele1 = null;
- ReconstructedParticle ele2 = null;
- int sumCharge = 0;
- int numChargedParticles = 0;
- for (final ReconstructedParticle fsPart : finalStateParticles) {
- if (isGBL != TrackType.isGBL(fsPart.getType())) {
- continue;
- }
- if (debug) {
- LOGGER.info("PDGID = " + fsPart.getParticleIDUsed() + "; charge = " + fsPart.getCharge() + "; pz = "
- + fsPart.getMomentum().x());
- }
- final double charge = fsPart.getCharge();
- sumCharge += charge;
- if (charge != 0) {
- numChargedParticles++;
- if (charge < 1) {
- if (ele1 == null) {
- ele1 = fsPart;
- } else if (!hasSharedStrips(ele1, fsPart, hitToStrips, hitToRotated)) {
- ele2 = fsPart;
- }
- }
- }
- }
- sumChargeHisto.fill(sumCharge);
- numChargeHisto.fill(numChargedParticles);
-
- if (ele1 != null && ele2 != null) {
- final Hep3Vector p1 = VecOp.mult(beamAxisRotation, ele1.getMomentum());
- final Hep3Vector p2 = VecOp.mult(beamAxisRotation, ele2.getMomentum());
- // Hep3Vector beamAxis = new BasicHep3Vector(Math.sin(0.0305), 0, Math.cos(0.0305));
- // LOGGER.info(p1);
- // LOGGER.info(VecOp.mult(rot, p1));
-
- final double theta1 = Math.acos(p1.z() / p1.magnitude());
- final double theta2 = Math.acos(p2.z() / p2.magnitude());
- final double phi1 = Math.atan2(p1.y(), p1.x());
- final double phi2 = Math.atan2(p2.y(), p2.x());
- phiEleVsphiEle.fill(phi1, phi2);
- pEleVspEle.fill(ele1.getMomentum().magnitude(), ele2.getMomentum().magnitude());
- pyEleVspyEle.fill(ele1.getMomentum().y(), ele2.getMomentum().y());
- pxEleVspxEle.fill(ele1.getMomentum().x(), ele2.getMomentum().x());
- // remove beam electrons
- if (ele1.getMomentum().magnitude() < beambeamCut && ele2.getMomentum().magnitude() < beambeamCut) {
- pEleVspEleNoBeam.fill(ele1.getMomentum().magnitude(), ele2.getMomentum().magnitude());
- pyEleVspyEleNoBeam.fill(ele1.getMomentum().y(), ele2.getMomentum().y());
- pxEleVspxEleNoBeam.fill(ele1.getMomentum().x(), ele2.getMomentum().x());
- }
- // look at beam-beam events
- if (ele1.getMomentum().magnitude() > beambeamCut && ele2.getMomentum().magnitude() > beambeamCut) {
- pEleVspEleBeamBeam.fill(ele1.getMomentum().magnitude(), ele2.getMomentum().magnitude());
- pEleVsthetaBeamBeam.fill(p1.magnitude(), theta1);
- pEleVsthetaBeamBeam.fill(p2.magnitude(), theta2);
- thetaEleVsthetaBeamBeam.fill(theta1, theta2);
- }
-
- // look at "Moller" events (if that's what they really are
- if (ele1.getMomentum().magnitude() + ele2.getMomentum().magnitude() > molPSumMin
- && ele1.getMomentum().magnitude() + ele2.getMomentum().magnitude() < molPSumMax
- && p1.magnitude() < beambeamCut && p2.magnitude() < beambeamCut) {
-
- final Track ele1trk = ele1.getTracks().get(0);
- final Track ele2trk = ele2.getTracks().get(0);
- final SeedTrack stEle1 = TrackUtils.makeSeedTrackFromBaseTrack(ele1trk);
- final SeedTrack stEle2 = TrackUtils.makeSeedTrackFromBaseTrack(ele2trk);
- final BilliorTrack btEle1 = new BilliorTrack(stEle1.getSeedCandidate().getHelix());
- final BilliorTrack btEle2 = new BilliorTrack(stEle2.getSeedCandidate().getHelix());
- final BilliorVertex bv = this.fitVertex(btEle1, btEle2, TrackUtils.getBField(event.getDetector())
- .magnitude());
- // LOGGER.info("ee vertex: "+bv.toString());
- final double invMass = bv.getParameters().get("invMass");
- mollerMass.fill(invMass);
- mollerVx.fill(bv.getPosition().x());
- mollerVy.fill(bv.getPosition().y());
- mollerVz.fill(bv.getPosition().z());
- mollerXVsVtxZ.fill(bv.getPosition().x(), bv.getPosition().z());
- mollerYVsVtxZ.fill(bv.getPosition().y(), bv.getPosition().z());
- mollerXVsVtxY.fill(bv.getPosition().x(), bv.getPosition().y());
-
- final double ux = (ele1.getMomentum().x() + ele2.getMomentum().x())
- / (ele1.getMomentum().z() + ele2.getMomentum().z());
- final double uy = (ele1.getMomentum().y() + ele2.getMomentum().y())
- / (ele1.getMomentum().z() + ele2.getMomentum().z());
- mollerUx.fill(ux);
- mollerUy.fill(uy);
-
- // higher and lower energy electrons in moller pair
- final double pt1 = ele1.getMomentum().magnitude();
- final double pt2 = ele2.getMomentum().magnitude();
- final double ph = pt1 > pt2 ? pt1 : pt2;
- final double pl = pt1 > pt2 ? pt2 : pt1;
-
- mollerHiP.fill(ph);
- mollerLoP.fill(pl);
-
- mollerEitherP.fill(ph);
- mollerEitherP.fill(pl);
- mollerPsum.fill(pt1 + pt2);
-
- if (Math.abs(bv.getPosition().x()) < 2 && Math.abs(bv.getPosition().y()) < 0.5) {
- mollerMassVtxCut.fill(invMass);
- mollerVzVtxCut.fill(bv.getPosition().z());
- }
- pEleVspEleMoller.fill(p1.magnitude(), p2.magnitude());
- pEleVsthetaMoller.fill(p1.magnitude(), theta1);
- pEleVsthetaMoller.fill(p2.magnitude(), theta2);
- thetaEleVsthetaMoller.fill(theta1, theta2);
- }
- }
+ private BilliorVertex fitVertex(BilliorTrack electron, BilliorTrack positron, double bField) {
+ // Create a vertex fitter from the magnetic field.
+ double[] beamSize = {0.001, 0.2, 0.02};
+ BilliorVertexer vtxFitter = new BilliorVertexer(bField);
+ // TODO: The beam size should come from the conditions database.
+ vtxFitter.setBeamSize(beamSize);
+
+ // Perform the vertexing based on the specified constraint.
+ vtxFitter.doBeamSpotConstraint(false);
+
+ // Add the electron and positron tracks to a track list for
+ // the vertex fitter.
+ List<BilliorTrack> billiorTracks = new ArrayList<BilliorTrack>();
+
+ billiorTracks.add(electron);
+
+ billiorTracks.add(positron);
+
+ // Find and return a vertex based on the tracks.
+ return vtxFitter.fitVertex(billiorTracks);
+ }
+
+ private static boolean hasSharedStrips(ReconstructedParticle vertex, RelationalTable hittostrip, RelationalTable hittorotated) {
+ return hasSharedStrips(vertex.getParticles().get(0), vertex.getParticles().get(1), hittostrip, hittorotated);
+ }
+
+ private static boolean hasSharedStrips(ReconstructedParticle fs1, ReconstructedParticle fs2, RelationalTable hittostrip, RelationalTable hittorotated) {
+ return TrackUtils.hasSharedStrips(fs1.getTracks().get(0), fs2.getTracks().get(0), hittostrip, hittorotated);
}
}
